氣動薄膜執行機構模型優化及參數辨識

金縱橫，李世偉，尚群立

(浙江工業大學 信息工程學院，杭州 310023)

0 引言

摩擦存在于任何現實的機械設備中，其非線性特征是影響控制回路性能的主要因素之一[1-4]。深入研究氣動薄膜執行機構模型對開發新型智能電氣閥門定位器控制算法，提高生產效益有重要實用價值。

在先前的研究中，針對氣動閥中薄膜執行機構的研究偏少且主要研究其摩擦行為[5-6]。文獻[7]對氣動執行機構和比例方向閥組成的氣動系統進行機理建模。文獻[8]采集氣動執行機構輸入和輸出做基于數據驅動的系統辨識。對于新型摩擦模型不斷被提出。文獻[9]提出一種應用于切削過程的新型摩擦模型。文獻[10]通過改進LuGre摩擦模型(LuGre-Mod)，并驗證該模型可描述純滑動域中的順時針和逆時針磁滯回線。但新型摩擦模型并不能顯著提高實驗數據的預測精度，同時也大幅增加模型復雜度，提高工業生產成本和帶來更大的計算工作量[11]。在實際應用中需要在摩擦預測精度和模型復雜度之間做進一步權衡。對于摩擦模型的選擇，在工程應用中，仍普遍采用較簡單的摩擦模型，例如：Choudhury[12]，Kano[13]，He[14]。文獻[15]比較了應用于氣動閥的共八種基于物理原理和基于實驗數據的摩擦模型，在實際氣動閥上按ANSI/ISA標準進行測試，并與理論仿真結果進行比較，結果顯示Karnopp、LuGre、Kano三種摩擦模型具有最佳預期效果。對于摩擦模型參數辨識方法，文獻[16]提出將擬牛頓法應用于氣動閥摩擦模型辨識，進而得到該模型的Stribeck參數。文獻[17]針對Karnopp模型參數(m，k，Fc，Fv，Fs)，設計多種辨識算法，但每種算法僅能同時辨識一到兩個參數。為了降低測量噪聲的影響，需要對采集到的數據做濾波處理。文獻[18]對氣室氣壓和閥桿位移用零相位二階線性濾波器進行濾波。

上述研究從建模方法，提出新的摩擦模型，比較摩擦模型實驗效果，對摩擦模型進行更精準的參數辨識等方面去研究氣動執行機構等機械設備[19]。本文對Karnopp摩擦模型動靜摩擦切換條件進行修改，并將一種基于多元線性回歸和最小二乘法的參數辨識方法應用于估算氣動薄膜執行機構中移動部件質量，彈簧剛度系數，粘滯摩擦系數，庫倫摩擦，彈簧預緊力。通過仿真和現場實驗驗證參數辨識方法的有效性。結果表明，整合修改后摩擦模型的氣動薄膜執行機構動力學模型能準確模擬其動態過程。

1 氣動薄膜執行機構模型

工業生產中常見的氣開式調節閥，如圖1所示。

圖1 氣動閥示意圖

膜片將氣室氣壓轉換為力，并與方向相反的彈簧彈力共同決定閥桿的位移，進而改變閥芯和閥座之間的開合程度，控制管路流量。對圖1移動部件進行受力分析，得出氣動薄膜執行機構的動力學方程[15]：

(1)

式中，m是氣動薄膜執行機構移動部件的質量；Fe(t)=S·P(t)，其中Fe(t)是氣室壓力，S是膜片受力面積，P(t)是氣室氣壓；Fk(t)=k·x(t)，其中Fk(t)是彈簧彈力，k是彈簧剛度系數，x(t)是閥桿位移；Ff(t)是摩擦力；Fl是管道內液體施加在閥桿上的力，在氣動閥空載狀態下可以忽略；Fi是彈簧預緊力，當氣動薄膜執行機構不受氣室壓力時，用于壓緊閥芯與閥座，關閉管路，該力不會隨閥桿運動改變。

本次研究在氣動薄膜執行機構空載狀態下進行，省略后的動力學方程：

(2)

在先前的研究中，不同的作者對Ff(t)有不同的描述，從而提出不同的摩擦模型，但新型摩擦模型包含更多的參數，更復雜，且并不能顯著提高仿真精度。行業中仍普遍使用較簡單的摩擦模型。文獻[10]實際測試八個摩擦模型中，Karnopp、LuGre、Kano三種摩擦模型具有最佳預期效果。結合實際需求和實驗可提供的測量數據，本次研究選擇Stribeck曲線作為閥桿相對運動速度不為零時所受摩擦力的大小，如圖2所示。該曲線中摩擦力的數學表達式可以用式(3)表示：

圖2 Stribeck曲線

(3)

(4)

(5)

當速度小于DV時，摩擦力：

Ff(Fr)=min(Fs，|Fr|)sign(Fr)

(6)

Fr=S·P(t)-t·x(t)-Fi

(7)

DV是臨界速度，如果閥桿速度小于臨界速度，則認為移動部件處于靜止狀態。式(5)的第一行為移動部件滑動時閥桿所受摩擦力，式(5)的第二行為移動部件靜止和動靜切換臨界時閥桿所受摩擦力。臨界速度的目的在于解決零速度檢測問題，界定閥桿是處于靜止狀態還是滑動狀態。文獻[10]提及臨界速度取值標準是觀察實驗和仿真的擬合度，這需多次調整臨界速度值并仿真。

針對臨界速度的取值標準需要基于實驗和仿真的擬合度這一問題，本文對簡化后的Karnopp摩擦模型提出進一步的修改，并對修改后的摩擦模型做離散化處理。

Ff[n]=

(8)

對移動部件處于滑動階段的氣動薄膜執行機構動力學方程做離散化處理：

(9)

式中，

(10)

(11)

n是自然數，T是采樣周期。

需要注意的是，在仿真模型中，需指定移動部件初始速度和初始加速度均為零。

2 參數辨識方法

2.1 基于多元線性回歸的方法

定義參數向量：

α=[FimFcFvk]

(12)

當閥桿處于滑動階段時，氣動薄膜執行機構離散化動力學方程：

(13)

式(13)和參數向量α是線性的，所以定義回歸向量β[n]：

(14)

參數向量α可以通過最小二乘法估計得到：

(15)

為了后續現場實驗得出的參數辨識結果作對比，接下來將引入傳統參數辨識方法。

2.2 傳統方法

剛度系數k：氣動薄膜執行機構做開環階躍測試時，記錄每個階躍測試下閥桿靜止時的平均位移和相應的平均氣室氣壓。

(16)

式中，ΔP，Δx分別為兩個階躍測試對應氣室氣壓差和閥桿位移差。

最大靜摩擦Fs：氣動薄膜執行機構做開環周期和幅值不變的三角或正弦測試時，記錄每個周期中兩個閥桿位移方向倒轉處的位移和對應氣室氣壓。

(17)

式中，ΔP，Δx分別為測試中每個周期兩個閥桿位移方向倒轉處對應氣室氣壓差和閥桿位移差。

庫侖摩擦Fc：氣動薄膜執行機構做開環周期三角或梯形斜坡測試時，記錄閥桿正向運動和反向運動期間，同一閥桿位移下的氣室氣壓。

(18)

式中，ΔP為閥桿正向運動和反向運動期間，同一閥桿位移下的氣室氣壓差。

粘滯摩擦系數Fv和預緊力Fi：氣動薄膜執行機構做開環周期三角或梯形斜坡測試時，記錄同一周期閥桿處于正向和反向運動狀態下各多組閥桿位移和相應氣室氣壓，并計算該周期正向和反向斜坡速度即為閥桿速度。

(19)

移動部件質量m通常由閥門供應商提供。

為保證結果正確性，所有參數辨識方法中每個測試需做多次，并取平均值。

3 參數辨識實驗及分析

3.1 仿真模擬分析

由氣動薄膜執行機構動力學方程和摩擦模型方程，通過MATLAB/Simulink搭建仿真模型，該模型如圖3所示。

圖3 氣動薄膜執行機構Simulink模型

該仿真模型能夠模擬氣動薄膜執行機構的開環運行狀態，模型中各參數設定值：m=2 kg，Fv=5 000 N·s/m，Fc=200 N，Fs=200 N，k=200 000 N/m，Fi=1 800 N。為驗證摩擦模型修改前后的區別及參數辨識方法的可靠性，進行了如下仿真實驗。

仿真實驗一：

在開環無測量噪聲情況下，原摩擦模型與修改后的摩擦模型在同一正弦激勵下的閥桿位移過程的對比。

實驗目的：觀察原模型與修改后模型之間的區別，驗證優化結果。

采樣頻率：10 000 Hz

圖4中，實線表示模型修改后DV=0.000 5 m/s在的閥桿位移(m)隨時間(s)變化曲線，虛線代表原模型在DV=0.000 5 m/s的閥桿位移(m)隨時間(s)變化曲線。

圖4 摩擦模型修改前后閥桿位移對比

圖5 加入噪聲的仿真結果

實驗結果顯示，摩擦模型修改前，當閥桿運動臨近靜止時，閥桿位移呈階梯式運動。修改后的模型與原模型對比，修改后的模型閥桿位移更平滑。模型修改前后閥桿位移差異并不顯著，原因在于在氣動薄膜執行機構動力學模型中，相對于彈簧彈力，摩擦力占比并不顯著。仿真實驗二選擇修改后的模型進行實驗。

仿真實驗二：

在開環無測量噪聲和有測量噪聲情況下，修改后的氣動薄膜執行機構模型在同一激勵下的參數辨識結果對比。

實驗目的：模擬驗證在實際工況下，參數辨識方法的可行性。

采樣頻率：10 000 Hz

實驗具體數據如表1所示。

表1 仿真實驗二參數辨識結果

結果表明，在有噪聲的情況下，參數辨識方法得出的辨識值誤差在可接受范圍內，該方法可行。

3.2 現場實驗

本章對前文所述的參數辨識方法在實驗臺架上進行驗證，實驗臺架結構如圖6所示。

實驗裝置：氣動閥，I/P轉換器，閥位變送器，壓力變送器，其中氣動閥閥桿行程16 mm，填料為聚四氟乙烯，膜片有效受力面積0.038 m2，閥門移動部件質量1.6 kg，上述氣動閥參數均由閥門供應商提供。數據采集裝置為美國國家儀器有限公司的NI-DAQ-9332機箱。

為了獲取參數辨識所需數據，確保回歸向量β[n]中所有系數都能影響氣動薄膜執行機構模型動力學方程，通過上位機給定正弦激勵信號，同時為了保證數值推算的準確性，設定采樣頻率為1 000 Hz。實驗采集的氣室氣壓，閥桿位移，如圖7所示。由于I/P轉換器的動態特性，采集的氣室壓力與給定的正弦激勵信號在波形上有些許差異。

圖7 實驗結果

為了獲取氣動薄膜執行機構模型參數的最優值，下面將引入三個統計學參數作為評判標準。

式(20)為實際氣室壓力和估計氣室壓力的標準誤差。

(20)

式(21)為實際氣室壓力和估計氣室壓力的相關系數。

(21)

最后一個評判標準是參數向量α中每個參數對應的p-value。當p-value小于0.05時，對應參數符合大于95%置信區間標準。

3.3 實驗結果分析

對實驗平臺得出的實驗數據進行分析處理。

其中在不同Δυ下的各項參數辨識值變化曲線如圖8所示。

圖8 不同Δυ下各參數辨識值

估計外力與實際外力的標準誤差和相關系數如圖9所示。

圖9 估計外力與實際外力的標準誤差和相關系數

不同Δυ下的各項參數p-value變化曲線如圖10所示。

圖10 不同Δυ下各參數p-value

由第二節參數辨識方法，結合圖8～10，在Δυ介于[0.004，0.005 8]區間內，圖8中參數向量α中各參數辨識值變化趨于平穩，圖9中標準誤差達到最小范圍，相關系數接近于1，圖10中各參數p-value趨于0，說明該區間內模型中各參數是顯著的。在Δυ=5.72×10-3m/s時，標準誤差趨于最小值，相關系數達到最大值，各參數對應p-value均小于0.05。圖9所示，當Δυ<4×10-3m/s時，標準誤差和相關系數均顯著高于區間[0.004，0.005 8]，且Fc出現負數的情況，原因在于閥桿相對運動速度較小的時間段多集中在閥桿位移方向反轉的過度區間，該區間氣室氣壓抖動劇烈，如圖7氣室氣壓曲線所示，這極大影響辨識的精度。當Δυ>5.8×10-3m/s時，m，Fv，Fc，k，Fi辨識結果值均隨著Δυ的遞增波動加劇，原因在于隨著Δυ的遞增，符合參數辨識的樣本越少。綜上所述，當Δυ=5.72×10-3m/s時，各參數的辨識結果值最接近真實值。

其中基于多元線性回歸方法和傳統方法進行參數辨識的結果數據對比如表2所示。

表2 現場實驗參數辨識結果

表2中粘滯摩擦系數的辨識結果差異較大，原因有以下兩點：

1)在氣動薄膜執行機構的動力學方程中，相對于占比更高的彈簧彈力，粘滯摩擦力占比極小，對氣動薄膜執行機構運動的貢獻并不顯著。

2)式(19)在估算粘滯摩擦系數和預緊力時，還有彈簧剛度系數和庫倫摩擦參與計算，兩個變量的估算本身就存在誤差。

4 氣動薄膜執行機構模型驗證

為了驗證氣動薄膜執行機構模型的可靠性，需要比對仿真和現場實驗結果。仿真模型輸入信號均為現場采集后經濾波的氣室氣壓。

測試一：開環階躍輸入響應測試，驗證模型能否重現死區和評估氣動薄膜執行機構動態響應。圖11為該激勵信號作用下仿真與實驗的閥桿位移對照結果。

圖11 階躍輸入響應測試

測試二：輸入信號為隨機信號。選擇該信號的原因在于當閥桿處于滑動狀態時候，閥桿位移，速度和加速度都隨時間不斷變化，此時，氣動閥模型中所有參數都會影響移動部件的力平衡。圖12為該激勵信號下仿真與實驗的閥桿位置對照結果。

圖12 正弦掃頻測試

需要注意的是，由于I/P轉換器的動態特性，在閥桿位移變小，I/P轉換器對外排氣時，會出現圖11和圖12中排氣過量導致位移出現尖峰的情況。

圖11和圖12中仿真和實驗閥桿位移沒有完全重合的原因有以下兩點：

1)實際生產中，氣動薄膜執行機構閥桿不同部位與填料接觸部分是不均勻的，摩擦系數可變，而人為量化的摩擦系數均為定值。

2)I/P轉換器在給氣室充排氣瞬間，氣源供氣壓力不穩定且氣路和氣室存在一定程度漏氣。

3)實驗使用的設備存在不同程度誤差。

5 結束語

本研究對應用于氣動薄膜執行機構的Karnopp摩擦模型動靜摩擦切換條件進行修改。在傳統方法的基礎是，將多元線性回歸和最小二乘法對氣動薄膜執行機構參數進行辨識，通過仿真和現場實驗驗證模型修改前后區別和參數辨識方法的可行性，得出如下結論：

1)應用于氣動薄膜執行機構的傳統Karnopp摩擦模型存在的實際問題進行分析，并對其動靜摩擦切換條件進行修改，解決了DV取值問題。

2)通過帶測量噪聲的仿真和現場采集氣動薄膜執行機構輸入輸出的真實實驗數據，驗證基于多元線性回歸和最小二乘法的參數辨識方法是有效的。該方法相對于傳統方法，極大減少建模工作量。

3)修改后的模型和參數辨識方法并不局限于本研究，經適當修改可應用于其它領域，提高其實用性和普適性。