脈沖燃燒風洞測力系統傳力特性研究

李世超， 高宏力， 劉勃鍇， 付國強， 魯采江

(西南交通大學 機械工程學院，成都 610000)

飛行器研發過程中需進行風洞測力試驗，以了解其氣動特性[1-5]。測試精度不但影響飛行器的研制風險與成本，而且還制約其氣動特性的評估與優化。脈沖燃燒風洞是用于揭示高超聲速飛行器氣動特性的地面試驗設備[6-8]，作用于飛行器模型表面的氣動力與力矩由測力系統中的天平測量。由于飛行器模型所受的氣動力要經過測力系統內部結構才能傳遞至天平上，因此測力系統的傳力特性會直接影響氣動力的測試精度[9]。

目前關于風洞測力系統測試精度的研究主要集中在以下幾個方面：新型測試原理的探索；慣性力補償技術研究；振動控制技術研究；在測試原理的研究方面，Sanderson等[10]基于應力波的傳遞原理，設計了一種應力波天平；羅也凡等[11]提出了一種能進行慣性力補償的加速度計應變天平；在慣性力補償方面，Bernstein等[12]通過在飛行器模型的適當位置安裝加速度傳感器來獲取其在試驗過程中因振動引起的慣性力，在此基礎上，對天平輸出信號進行補償；在振動控制方面，Capone等[13]采用被動減振技術，通過在測力系統的適當位置加裝阻尼元件，以降低飛行器模型的振動幅度；陳衛東等[14]基于主動減振技術，將激振器安裝在飛行器模型空腔內合適位置，并基于各種控制算法，通過控制激振器的輸出來達到降低飛行器模型振動幅度的目的；以上文獻均是針對改善測力系統的測試精度進行研究，然而關于測力系統傳力特性的評估卻鮮見報道。

本文提出了一種基于有限元仿真來評估測力系統傳力特性的方法。首先，通過參數辨識得到了各零部件材料屬性真實值，以此為基礎，建立了測力系統有限元模型，并對整機進行模態試驗，驗證了該有限元模型的準確性；其次，基于飛行動力學相關理論得出了飛行器模型重心與天平之間的傳力特性可以反映測力系統整機傳力特性的結論；對測力系統有限元模型進行諧響應分析，并提取天平中應變計粘貼處的位移幅頻曲線，結合位移幅頻曲線與天平測力公式獲得激勵力在飛行器模型重心與天平浮動框間的傳遞幅頻特性曲線；最后，將作用于飛行器表面的氣動力等效為飛行器模型重心處的集中力(激勵力)，并基于傅里葉變換展開成一系列諧波力信號，通過幅頻特性曲線計算出各諧波力傳遞至天平浮動框后的大小，以此為基礎，結合天平最低測試分辨率篩選出未被天平識別的諧波力信號，并計算未被識別的諧波力信號總能量占激勵力信號總能量的百分比，以此來評價測力系統傳力特性的優劣。

1 概 述

本文以某脈沖燃燒風洞測力系統為研究對象，其結構如圖1所示，由高超聲速飛行器模型、測力天平和支架三部分組成，天平浮動框與飛行器模型間采用螺釘固連，天平固定框與支架之間同樣采用螺釘緊固的方式連接。其中測力天平為盒式六分量應變天平，其結構簡圖如圖2所示，由浮動框、固定框和彈性測量元件組成，彈性測量元件將浮動框和固定框連接成一整體。工作時風洞產生的氣流作用在飛行器模型上產生氣動力，使得天平浮動框與固定框發生相對位移，粘貼于測量元件上的箔式應變計產生形變，從而導致電阻發生變化，惠斯通電橋將電阻的變化轉變為電壓輸出，此即測力系統的測試原理，圖3展示了箔式應變計A1～A10在天平中的安裝位置。

1-飛行器模型；2-應變天平；3-支架圖1 測力系統結構簡圖Fig.1 Measuring system structure diagram

圖2 應變天平結構簡圖Fig.2 Strain balance structure diagram

圖3 天平箔式應變計粘貼位置Fig.3 Foil strain gauge paste position

2 測力系統有限元模型的建立

2.1 測力系統幾何模型簡化

進行測力系統有限元建模時，由于系統中包含諸多圓角，倒角和孔等微小結構，它們對整機動態特性影響小，然而模型在進行有限元網格劃分時這些結構會明顯增加網格數量，增大計算量，因此，在保證有限元模型精度的前提下對整機模型進行合理的簡化能有效降低計算成本，提高計算效率，測力系統簡化原則如下：

(1) 忽略結構中的小倒角和圓角；

(2) 刪除結構中的小定位孔及螺紋孔等細小特征。

簡化后的模型如圖4所示。

2.2 測力系統有限元建模

測力系統是由各零部件通過結合部連接而成的裝配體。其中結合部指的是零部件間相互接觸并能傳遞力的部分，研究表明，裝配體整體剛度與阻尼的60%和90%以上由結合部貢獻[15]，因此結合部動態特性對整機動態特性有重要影響，在測力系統動力學建模時須正確描述結合部動態特性。因此，測力系統整機有限元建模包括各零部件有限元建模和結合部有限元建模，有限元建模采用Ansys軟件實現。

圖4 測力系統簡化模型Fig.4 Force measurement system simplified model

2.2.1 單一零部件建模

測力系統中的主要零部件包括飛行器模型，天平與支架，采用20節點的solid186單元和10節點的solid187單元分別對飛行器模型，天平和支架進行網格劃分, 共得到88 455個單元，181 626個節點。在建模過程中需要賦予各零部件相應的材料屬性，要想建立準確的各零部件有限元模型，需賦予有限元模型正確的材料屬性，然而各零部件材料屬性的參考值與真實值間存在一定的差異，因此，為了提高各零部件有限元模型的精度，需對零部件的材料屬性進行辨識。本文基于優化思想分別對飛行器模型與天平的彈性模量，密度，泊松比真實值進行辨識。

2.2.1.1 零部件材料屬性辨識

基于優化算法的材料屬性辨識思想為：采用優化策略不斷調整有限元模型中的材料屬性值，最小化零部件有限元模型前N階理論模態與前N階試驗模態間的綜合誤差，直至誤差小于收斂值ε時，即可認為此時有限元模型中的材料屬性值即為零部件材料屬性的真實值。

基于該辨識思想，具體的參數辨識流程描述如下：

式中：N為模態階數；ε為優化截止條件。

(2) 在MATLAB中編寫能與Ansys進行相互調用的優化程序，該程序基于鮑威爾優化理論，能夠實現單一目標函數、多變量的優化。首先，通過MATLAB調用Ansys進行模態分析并提取前四階模態頻率，基于Ansys的APDL語言將計算獲得的前四階模態頻率反饋到MATLAB中，通過優化算法，MATLAB將優化后各零部件彈性模量、泊松比與密度值傳遞給Ansys，以此流程自動循環迭代直至達到優化終止條件，輸出最終材料屬性值。

2.2.1.2 零部件模態試驗

在零部件參數辨識過程中需要各零部件的前N階試驗模態，本文采用單點激勵，多點響應數據采集(Single Input Multiple Output, SIMO)方式分別對零部件進行錘擊模態試驗。飛行器模型與天平均采用自由懸掛的方式，近似模擬零件自由狀態，各零部件上力錘激勵點和加速度測點的位置如圖5和圖6所示，試驗設備包括北京東方所的數據采集系統，力錘與三向壓電式加速度傳感器。

圖5 飛行器模型測點與激勵點的位置分布Fig.5 Distribution of locations of measurement points and excitation points in aircraft models

圖6 天平測點與激勵點的位置分布Fig.6 Distribution of locations of measurement points and excitation points in strain balance

2.2.1.3 零部件有限元模型驗證

經過參數辨識后，辨識得到的各零部件材料屬性值列于表1。對比各零部件參數辨識后的計算模態頻率與試驗模態頻率如表2所示，從表2可以發現，參數辨識后各零部件的計算模態值與試驗模態值間的最大誤差分別為8.5%與9.6%，其余誤差均在5%內，由此表明所建立的零部件有限元模型能準確描述其動力學特性。

表1 零部件參數辨識值Tab.1 Parameter identification value of parts

表2 參數辨識后零部件理論模態值與試驗模態值對比Tab.2 Comparison of theoretical modal values and experimental modal values of component parts after parameter identification

2.2.2 結合部建模

測力系統中的主要結合部包括飛行器模型-天平螺釘固定結合部與天平-支架螺釘固定結合部，由于結合部表現出既存儲能量又釋放能量的特性[16]，因此采用彈簧-阻尼模型來描述結合部動態特性，對系統動態特性影響較小的結合部在Ansys中采用黏接處理。在螺釘固定結合部上的每個螺釘位置處用一個法向彈簧-阻尼單元和兩個切向彈簧-阻尼單元模擬，螺釘固定結合部彈簧-阻尼模型如圖7所示，其中圓圈代表彈簧-阻尼單元的位置，它與螺釘的中心軸位置重合。

圖7 結合部彈簧阻尼模型Fig.7 Spring damping model of joint

其中彈簧-阻尼單元的剛度系數與阻尼系數通過計算獲得，具體流程如下：

步驟1運用吉村允效法[17]獲得結合面上單位面積的剛度系數與阻尼系數；

步驟2對結合面上單位面積剛度系數與阻尼系數積分，求得結合面的總體剛度系數和阻尼系數；

步驟3將總體剛度系數與阻尼系數均分至每個彈簧-阻尼單元從而獲得每個彈簧-阻尼單元的剛度系數和阻尼系數；

步驟2中結合面總剛度系數和總阻尼系數的計算公式為

(1)

(2)

(3)

(4)

式中:Pn為結合部面壓；kn(Pn)，kt(Pn)，cn(Pn)，ct(Pn)分別為結合面單位接觸面積的法向剛度、阻尼系數和切向剛度、阻尼系數；Kn，Cn，Kt，Ct為結合面的總體法向剛度、總體法向阻尼系數、總體切向剛度、總體切向阻尼系數。結合部面壓可由式(5)、式(6)計算。

(5)

(6)

式中：Pn為結合面面壓；d2為螺紋中徑；φ為螺紋升角；ρv為螺紋當量摩擦角；μ為螺母與被連接件支撐面間的摩擦因數；T為螺釘的預緊力矩；Dw為六角螺母直徑；d0為螺釘中徑。

以測力系統中的飛行器模型-天平螺釘固定結合部、天平-支架螺釘固定結合部為例，闡述其彈簧-阻尼單元的剛度系數和阻尼系數的具體計算過程：首先由式(5)式(6)計算出飛行器模型-天平螺釘固定結合部、天平-支架螺釘固定結合部的接觸面壓分別為8.95 MPa，2.4 MPa。基于結合部面壓和結合面表面粗糙度，通過查詢文獻[18]得到結合面單位接觸面積的剛度與阻尼系數，如表3所示，以上結合部中結合面的粗糙度均為3.2 μm。再根據式(1)～式(4)求出結合面總體剛度系數和阻尼系數。最后，將結合部整體剛度系數與阻尼系數均分至該結合部上的所有彈簧-阻尼單元，從而獲得單個彈簧阻尼單元的剛度系數和阻尼系數如表4所示。

表3 結合部單位面積接觸剛度與阻尼Tab.3 The contact stiffness and damping of unit area

表4 彈簧-阻尼單元剛度系數與阻尼系數Tab.4 The stiffness coefficient and damping coefficient in spring-damping element

2.2.3 整機有限元模型驗證

將結合部彈簧-阻尼模型與零部件有限元模型綜合成測力系統整機有限元模型，其中，支架底面采用完全固定約束，如圖8所示。為驗證模型的準確性，對測力系統進行整機模態試驗，對比試驗模態頻率與計算模態頻率如表5所示，其中計算模態值與試驗模態值間的最大誤差為10.3%，其余誤差均在10%內，證明測力系統整機有限元模型能準確描述其動態特性。

圖8 測力系統整機有限元模型Fig.8 Finite element model of the force measuring system

表5 模型驗證Tab.5 Model validation

3 天平標定

天平測力公式是傳力特性分析的基礎，在Ansys中采用虛擬靜態校準方法，獲得天平輸入與輸出之間的關系，天平測力公式為

(7)

式中：Fx，Fy，Fz，Mx，My，Mz分別為作用于天平浮動框沿x，y，z方向的軸向力與轉矩；εn,x，εn,y，εn,z分別為編號為n的箔式應變計在x，y，z方向上的位移；kn,n為天平靜校矩陣中的系數。

采用表6所示的校準方案，校準力的加載位置為浮動框上表面，如圖9所示，靜校過程中每次加載只施加一種類型的力，其余力均為零，提取每次加載后所有箔式應變計粘貼處x，y，z三個方向上的位移，將每次校準過程中的校準力向量與箔式應變計的位移向量帶入到式(7)中，可得標定矩陣中的一列，重復此過程，直至獲得標定矩陣中所有系數為止。

為驗證天平測力公式的準確性，在天平浮動框上表面同時施加表7所列的6種力，對比真實載荷與辨識載荷如表8所示，誤差幾乎為零，由此表明，天平測力公式精度完全滿足要求。

圖9 天平校準力加載位置Fig.9 The position of balance calibration loading

表6 天平校準方案Tab.6 The calibration program of balance

表7 天平測力公式精度驗證施力方案Tab.7 The force scheme of accuracy verification of balance force formula

表8 天平測力公式驗證Tab.8 Accuracy verification of balance force formula

4 測力系統傳力特性分析

在風洞測力試驗中，飛行器模型受到的氣動力需要經過測力系統內部結構才能傳遞至天平上，因此測力系統傳力特性的優劣性將直接影響其測試精度，為了保障測力數據的可信度，需要對測力系統的傳力特性進行分析與評價。

天平對力的測量有分辨率限制，當作用于天平上的力低于其分辨率時，天平將無信號輸出，此力便未能測得，本文中的應變天平分辨率為2 N，以此作為測力系統傳力特性的判斷依據，傳力特性分析流程大致如下，流程圖如圖10所示：

1) 選取能反映測力系統整機傳力特性的傳力路徑；

2) 基于有限元模型進行諧響應分析，提取激振力F(ω)作用于傳力路徑輸入點時天平中所有應變計粘貼處的位移幅頻特性曲線X1(ω)，X2(ω)，…,Xn(ω)，其中n表示應變計編號；

3) 將應變計粘貼處位移幅頻特性曲線X1(ω)，X2(ω)，…,Xn(ω)代入天平測力式(7)中，求出激振力F(ω)傳遞至天平浮動框后的力幅頻特性曲線F′(ω)；

4) 結合激振力幅頻特性曲線F(ω)及激振力傳遞至天平浮動框后的力幅頻特性曲線F′(ω)，獲得力在傳力路徑輸入點與天平浮動框間的傳遞特性，計算公式如(8)所示；

5) 將作用于飛行器表面的氣動力等效成集中力F并施加于傳力路徑的輸入點；

6) 基于傅里葉變換將力F展開成一系列諧波信號并代入到傳遞函數式(8)計算出氣動力F的各諧波力信號從傳力路徑的輸入點傳遞至天平浮動框后的大小，傅里葉變換的展開公式如式(9)所示；

7) 以天平分辨率為依據篩選出為被識別的諧波力信號；

8) 計算未被識別的諧波力信號的能量E′占輸入力信號總能量E的百分比ε，ε的大小可以反映測力系統的傳力特性，計算公式如式(10)所示，能量E與E′的計算公式如式(11)與式(12)所示

(8)

(9)

(10)

(11)

(12)

式中：E為原始信號的總能量；F(ω)為激勵信號傅里葉變換展開后頻率為ω的諧波力信號的幅值；Δω為各諧波力信號頻率間的差值；E′為未被天平識別的諧波力信號所占的能量；ωn表示第n個未識別的諧波力信號的頻率；N為未識別的諧波力信號的總個數。

圖10 傳力特性分析流程圖Fig.10 The flow chart of force characteristics analysis

4.1 傳遞路徑的選擇

本文中的力傳遞路徑指的是氣動力在測力系統中任意兩點之間傳遞的路徑，然而氣動力是以面力的形式作用在飛行器模型上的，因此氣動力在測力系統中的傳遞路徑有無數多條，需從中找出一條傳遞路徑能夠描述整機的傳力特性。由于測力系統所測得力指的是飛行器模型重心處的集中力，這是因為做飛行器飛行動力學分析時需要將其看成是剛體，并將面力形式的氣動力等效成作用在飛行器重心處的集中力，所以，風洞測力試驗的目的是給出飛行器模型在各種飛行條件下其重心處所受的集中力，基于此在設計時給予了飛行器模型較高的剛度，盡量讓其傳力特性接近于剛體的傳力特性，所以，飛行器模型上任意一點至天平的傳力特性可以等效成飛行器模型重心處到天平的傳力特性，因此選用飛行器模型重心處與天平間的傳力特性來描述測力系統整機傳力特性。通過計算，飛行器模型重心在笛卡爾坐標系中的位置如表9所示。

表9 飛行器模型重心位置坐標Tab.9 The center of gravity location about aircraft model

4.2 測力系統諧響應分析

為了分析飛行器模型重心處與天平浮動框表面間的傳力特性，首先需要獲得激振力作用在飛行器模型重心處時，天平各應變計處的位移幅頻響應曲線，因此對其進行諧響應分析。由于飛行器模型的重心并不在飛行器上，因此將重心沿x,y,z方向投影至飛行器模型上，如圖11中的a，b，c三點。分別在此三點處施加幅值為1 000 N、激振頻率為0～3 000 Hz的簡諧力，其中a點施加x方向的激振力，b點施加y方向的激振力，c點施加z方向上的激振力，由于篇幅的限制這里只展示天平中編號為a1的箔式應變計粘貼處在x，y，z方向上的位移幅頻特性曲線，如圖12所示。

圖11 激振力在飛行器模型上的投影Fig.11 Exciting force projection on aircraft model

(a) x方向的位移幅頻特性曲線

(b) y方向的位移幅頻特性曲線

(c) z方向的位移幅頻特性曲線圖12 應變片a1在x，y，z方向的位移幅頻特性曲線Fig.12 The displacement amplitude-frequency characteristic curve of a1 in x, y, z direction

4.3 測力系統軸向傳力特性

通過諧響應分析得到了激振力在a點沿x方向激振時各應變計粘貼處位移幅頻特性曲線，分別將該曲線中各頻率成分的幅值除以與之頻率相對應的激振力幅值，得到單位激振力作用下的位移幅頻特性曲線，并將其帶入到天平測力式(7)中進行計算，得到作用在a點的單位激振力傳遞至天平浮動框表面后的力幅頻特性曲線，并根據式(8)計算獲得激勵力輸入a點與天平浮動框間的傳力特性，如圖13所示。

圖14為某飛行器模型在進行脈沖風洞測力試驗過程中風洞測試段總壓和天平軸向輸出信號圖，虛線表示風洞試驗段總壓變化規律，實線表示天平軸向輸出電壓信號。從圖中可以看出，總壓上升階段為風洞啟動過程，總壓穩定階段為有效試驗階段，持續約500 ms。因此，在進行傳力特性分析時可以將作用于飛行器模型表面的氣動力等效成作用于模型a,b,c處分別沿x，y，z方向的階躍力，如圖15所示。基于傅里葉變換將階躍力展開成諧波力信號，將模型a點-天平浮動框力傳遞特性曲線(如圖13所示)中各頻率的幅值和與之相對應的諧波力幅值相乘，就可以得到作用在a點處的各諧波力信號傳遞至天平后的幅值，再以天平測試分辨率為評判依據篩選出輸入力中未被天平所感應到的諧波力成分，即找出幅值低于2 N的諧波力，分析結果表明，輸入力中未被感受到的諧波力共有99個，其頻率值離散分布在51～3 000 Hz，由式(10)計算出力在測力系統軸向方向上的傳遞誤差ε為1.681 9×10-5，當a點激勵力的幅值不斷改變時，測力系統軸向方向上的傳遞誤差ε如表10所示。由此表明，大部分的氣動力在測力系統軸向上可以很好的傳遞。

圖13 激勵力輸入a點與天平浮動框間的傳力特性Fig.13The force characteristics between point a and balance float box

圖14 脈沖風洞測力試驗總壓圖Fig.14 The total pressure diagram of pulse wind tunnel force test

圖15 作用于a，b，c三點的氣動力Fig.15 Air force acting on a, b, c three points

4.4 測力系統縱向傳力特性

采用如上方法，對測力系統縱向進行傳力特性分析，b點處作用激振力時，得到單位激振力傳遞至天平浮動框表面后的力幅頻特性曲線如圖16所示，分析結果表明，輸入力中未被感受到的諧波力共有304個，其頻率值離散分布在51～3 000 Hz，由式(10)計算出力在測力系統軸向方向上的傳遞誤差ε為5.920 1×10-5，當b點激勵力的幅值不斷改變時，測力系統縱向方向上的傳遞誤差ε如表10所示。由此表明，氣動力在測力系統縱向的傳遞損失較小，大部分的氣動力在測力系統軸向上可以很好的傳遞。

4.5 測力系統徑向傳力特性

同理，對測力系統徑向的傳力特性進行分析，得到單位激振力作用在c點時激振與天平浮動框間的力傳遞幅頻特性曲線如圖17所示，分析結果表明，輸入力中未被感受到的諧波力共有105個，其頻率值離散分布在51～3 000 Hz，由式(10)計算出力在測力系統軸向方向上的傳遞誤差ε為6.253 5×10-5，當c點激勵力的幅值不斷改變時，測力系統徑向力的傳遞誤差ε如表10所示。說明氣動力在測力系統縱向的傳遞損失很小，大部分的氣動力在測力系統軸向上可以很好的傳遞。

圖16 激勵力輸入b點與天平浮動框間的傳力特性Fig.16 The force characteristics between point b and balance float box

圖17 激勵力輸入c點與天平浮動框間的傳力特性Fig.17 The force characteristics between point c and balance float box

表10 測力系統傳力特性 Tab.10 The force characteristics of force measurement system

5 結 論

本文圍繞風洞測力系統傳力特性評估開展研究，建立了測力系統整機有限元模型；得到了飛行器模型重心與天平間的傳遞特性可以描述整機傳力特性的結論；獲得了力在飛行器模型重心與天平浮動框間傳遞的幅頻特性曲線；揭示了測力系統軸向、縱向與徑向傳力特性；通過對仿真結果分析，得出以下結論：

(1) 對測力系統進行了整機模態試驗，提取了系統前4階模態并與有限元仿真值進行了對比，其中第四階計算模態與試驗模態的誤差達到10.3%，其余均在10%以內，由此驗證了該模型的可靠性，從而保證了基于有限元模型的傳力特性分析的可信度；

(2) 研究發現，基于傅里葉變化將作用飛行器模型重心處的階躍力展開成一系列諧波力，并根據力在飛行器模型重心-天平浮動框間的傳遞特性計算各諧波力傳遞至天平浮動框后的大小，以天平的最小測試分辨率篩選出未被天平識別到的諧波力成分，計算未被辨識信號能量占原始信號總能量的百分比，以此來評估測力系統傳力特性優劣的方法是可行的；

(3) 測力系統傳力特性分析表明，1 000 N大小的氣動力沿軸向、縱向與徑向傳遞時，未被天平識別的諧波力信號能量占輸入信號總能量的百分比分別為1.681 9×10-5，5.920 1×10-5，6.235 3×10-5，并且隨著軸向力幅值的不斷變大，傳遞誤差均減小，氣動力幅值由1 000 N逐漸增大至6 000 N時，軸向、縱向與徑向的傳遞誤差變化范圍分別為：1.681 9×10-5～0，5.920 1×10-5～0，6.235 3×10-5～2.747 8×10-5，其中軸向傳力誤差最小，由此表明氣動力在測力系統的軸向、縱向與徑向能夠很好的傳遞，測力系統具有良好的測試精度。