任意拉格朗日-歐拉描述的薄平板大幅扭轉振動氣動特性研究

祝志文 林君福 唐意 王欽華

摘要：采用任意拉格朗日一歐拉（Arbitrary Lagrangian-Eulerian，ALE）描述法，推導了二維流變區域不可壓流體流動的控制方程。基于剛性斷面邊界特點，簡化得到了薄平板強迫振動繞流場控制方程，采用有限差分法開展了薄平板大幅強迫扭轉振動的CFD模擬。研究發現，薄平板小幅扭轉振動氣動力呈線性特性，但大幅扭轉振動氣動力非線性顯著，且氣動力非線性隨來流折算風速的增大而變得顯著。從一個周期內氣流做功來看，雖然大幅扭轉振動時薄平板轉動軸到后緣是耗散氣流能量的氣動穩定區，但當折算風速較高后，上游側的半個薄平板表面變為氣動不穩定區，此時氣流對整個薄平板所作總功由負變為正，薄平板將從氣流中吸收能量。研究認為，當折算風速較高時，大幅扭轉振動的薄平板將進入氣動不穩定狀態。

關鍵詞：渦激振動;薄平板;大幅振動;氣動不穩定;任意Lagrange-Eulerian描述

中圖分類號：TU311.3 文獻標志碼：A 文章編號：1004-4523（2020）01-0012-12

DOI：10.16385/j.cnki.issn.1004-4523.2020.01.002

引言

處于大氣邊界層中的大跨度橋梁受自然風的作用。在某些條件下，橋梁或構件可能出現大幅度的振動，如拉索、吊桿和主梁的大幅渦激振動，風雨振和尾流馳振、鄰近或進人顫振狀態的主梁大幅顫振等。這些大幅振動都是構件表面作為流體邊界的流動域發生顯著變形的流體一固體耦合作用問題;對鈍體外形的橋梁構件，或構件因大幅振動產生較大的相對攻角效應，其氣動力系統可能會因這種大幅振動呈現顯著的非線性。

在連續介質力學中，流體微元是微觀上遠大于分子運動平均自由程而宏觀上遠小于所研究的流動域的一團流體，并在宏觀上將流體微元看作一個質點。流體力學描述流體微元集合的運動狀態，并通過三種方法描述流體的運動，其中一種是拉格朗日法（Lagrangian），這種方法跟蹤流場每一質點的運動，并獲得流場流動參數隨時問的變化，其特點是計算網格與流體質點固定并同步運動，不存在質點與網格問的相對運動。與歐拉描述法相比，這種方法不僅可跟蹤流體質點的運動軌跡，而且能準確描述流體的移動界面。但當流動出現大變形時，將直接導致計算網格的畸變，計算失敗。歐拉法（Euleri-an）通過將感興趣的流動區域離散為位置固定而數量有限的網格點，進而描述流動參數的空間分布隨時問的變化。因其計算網格不隨物體形位變化而改變，所以能容易地處理流動區域的大變形，但描述流變界面的運動需要引入非常復雜的映射關系，可能導致較大的計算誤差。任意拉格朗日一歐拉法就是將拉格朗日法和歐拉法組合起來，取長補短。其特點是描述流體運動的空問點位置既不固定，也不隨流體運動，而是以某種方式隨時問變化，這種質點的變化方式能反映流動域的實際變化，因而可準確地描述流體的運動界面，并維持單元的合理形狀，因而可解決只用歐拉法或只用拉格朗日法解決不了的復雜運動邊界的大變形流動問題。

橋梁構件大幅振動中，構件周圍空氣所占據的實際空問一般是隨時問變化的，并具有運動的界面。可見，離散點空問分布固定的歐拉方法顯然不適合描述此類問題，因流動區域扭曲而引起計算網格的畸形也會導致計算失敗，因而合理描述應該采用任意拉格朗日一歐拉法。文獻[11-12]通過CFD模擬獲得的振動剛性模型上的氣動力和強迫振動位移時程，建立氣動力離散時問氣動模型，開展了顫振導數識別，但模型豎彎和扭轉振動的幅值很小。文獻開展了薄平板豎彎和扭轉大幅振動CFD模擬和氣動力小波分析，發現在氣動力時程中出現了比薄平板振動頻率高的二階高頻成分，也即大幅振動的薄平板氣動力為非線性，但因CFD基于Fluent開展，而Fluent的動網格即使在均勻網格上也僅有一階精度，因此研究結果的準確性需要驗證。文獻[14]基于Volterra理論開展了薄平板非線性氣動力系統識別和CFD模擬，通過豎彎振動幅值為1%和5%薄平板寬度的薄平板強迫振動CFD分析，發現氣動力的非線性效應并不顯著，原因可能是薄平板的振動幅度偏小。下面就橋梁氣彈特性分析常采用的薄平板為研究對象，針對二維不可壓流動特征，首先推導能適應任意流變區域的流體運動控制方程，再根據剛性斷面振動的特點，給出合適的空問點運動方式和相應的流體運動控制方程。并以薄平板為例研究其大幅振動的氣動特征。

1二維流體運動的任意拉格朗日-歐拉描述方程

基于式（18）-（22）可在二維計算域上實現對任意流變區域內流體的運動描述，也即為描述任意流變區域不可壓流動的基本方程，其中映射關系式（1）需根據所研究的流動區域的具體變形形式來確定。當物理域變形特征簡單時，對應的映射關系可簡化，此時，式（18）-（22）可能會得到明顯簡化，比如橋梁斷面的氣動彈性問題，因橋梁斷面的絕對剛性假設，使得映射關系明顯簡化，如下文所述。

2薄平板強迫振動繞流場控制方程

式（1）和（18）-（22）能描述任意變形物理域的流體流動，比如帶移動邊界的流體動力學問題。移動邊界可為內流場或外流場流體外部邊界的變形或移動，也可為外流場繞流物體在流體力作用下其形位的改變，這些均能導致與繞流物體接觸的流體邊界的移動和變形。在進行橋梁斷面一類繞流場模擬時，由于作用在橋梁斷面上的氣動壓力和摩阻力很小，因此橋梁斷面的變形很小，可近似將橋梁斷面看作剛體。這樣，橋梁斷面繞流場內邊界的運動可用剛性橋梁斷面的整體運動來描述，也即采用剛體運動的描述方法。在橋梁氣動彈性分析中，通常用橋梁斷面的豎彎和繞某一參考軸的轉動來描述橋梁斷面上任意點的運動。

對振動的剛性斷面，如果其計算域的遠場邊界不變，將能恒定地給出遠場邊界條件，但由于計算域內、外邊界之問的物理區域在不斷發生變化，如果進行數值模擬，則需在每一時問步形成一次網格，如果計算的時問步數量很大，由網格生成所帶來的附加計算量將很大，且難以保證每一步上重新生成的網格的高質量。相比較，如果在離剛性斷面足夠遠處劃定外邊界，并對計算域劃分高質量網格，數值模擬時將計算域網格和剛性斷面一起同步振動，就可只在計算開始時形成一次網格，而在每一時問步上給出相應的邊界條件，如此處理，不僅使得數值模擬的網格質量明顯提高，也能省去每一時問步上網格的重新生成。這樣，需求解的物理區域只隨時問發生位置改變，形狀并不改變。由于計算域事先給定，計算域和物理域的對應關系也能明確。

取絕對坐標系OXY平面內任意點為物理域之參考點，剛性斷面的振動方式為相對于參考軸的平動和繞參考點的轉動。如取計算域的坐標原點在物理域的映射點為參考點，并在此參考點上建立計算域的相對坐標系oxy，當相對坐標系隨剛性斷面作同步運動時，計算域Dcomp（x，y）和物理域Dphys（x，y，t）存在下述映射關系

通過選擇合適的時問和空問步離散格式，對式（36）-（38）和式（41）-（42）進行離散，在合理計算域網格基礎上，基于給定的初始和邊界條件下，就可開展剛性斷面強迫振動的數值模擬。

3薄平板扭轉振動繞流場控制方程的數值求解

3.1CFD數值實現

對圖6所示的長度為26的無厚度扭轉振動剛性薄平板（為便于查看，圖示有厚度），其繞流控制方程式本文采用非定常不可壓流計算的二階Projec-tion算法，即將上述控制方程分裂為依次求解動量方程得到中問速度場，由中問速度場求解壓力場的Poisson型方程，再由獲得的壓力場更新中問速度場，得到修正的速度場方程，如此反復求解，反復修正，并依據修正速度場在連續方程式（37）中產生的殘余值的大小來判斷當前時問步流場變量的收斂情況。當修正速度對連續方程的殘余值滿足收斂條件時，求解進入下一時問步。

由二階Projection算法得到的偏微分方程，空問離散采用基于交錯網格的有限差分法。在關于中問速度場的動量方程求解上，左端動量項采用二階外插的Adams-Bashforth格式，右端黏性項采用二階中心差分格式;壓力場方程左端項采用二階中心差分格式;時問導數采用二階中心格式，得到的差分方程均采用三對角矩陣算法求解，并通過Fortran程序編程實現。

因剛性薄平板厚度為零，計算域網格采用直角坐標系下的結構化網格，物面上下第一個網格高度為0.001b，并沿垂直物面方向網格尺度逐漸增大，但增長率不大于1.05，以捕捉邊界層的復雜流動，并適度控制網格規模和非定常計算量。薄平板沿流向均勻劃分32個網格，但從薄平板前緣往計算域入口，以及從薄平板后緣往計算域出口，采用不大于1.05的網格增長率逐漸增大網格。計算域總網格數為168×168，圖7是計算域整體和薄平板周圍的網格劃分，其中紅色方框內有模擬的薄平板。

設圖6所示的薄平板繞垂直OXY平面的薄平板轉動軸轉動，振動方式為正弦振動，即

式中fa為薄平板扭轉振動的頻率;a（t）和ao分別為扭轉振動位移角和位移角幅值（°），以薄平板前緣抬頭為正向，如圖6所示。

從控制方程式（34）-（36）和式（39），（40）的推導可知，網格系統和薄平板剛性固定，同步振動。因此在每一時問步，基于式（41）可得到網格和薄平板的同步振動速度;薄平板表面采用無滑移條件，因此薄平板物面流體的運動速度等于所在薄平板表面扭轉速度;對計算域的人口和上、下邊界，因來流平行于絕對參考系OXY的OX軸，需在每一時問步上分別投影到運動參考系oxy的兩個坐標軸方向，成為計算域對應邊界的速度邊界條件;計算域出口設置為紐曼邊界條件。對壓力邊界條件，在薄平板物面和計算域人口、上下邊界也采用紐曼邊界條件。式中B=26。考慮作用在振動薄平板上的氣動力對雷諾數不敏感，從降低不可壓流非定常迭代求解的計算量考慮，本文所有模擬的雷諾數均為1200。

3.2CFD方法驗證

薄平板顫振導數有Theodorsen解析解，通常用于檢驗氣彈分析CFD方法和程序的正確性。對單自由度扭轉強迫振動的薄平板，其非定常氣動力與相應顫振導數的表達式為：

分別開展不同折算風速的CFD模擬，得到作用在扭轉振動薄平板上的氣動力，利用式（46）和（47）基于最小二乘法識別薄平板的4個顫振導數，并將結果和Theodorsen解析解作了比較，如圖8所示。可見二者具有一致的趨勢線，高折算風速二者出現較小偏差，可能是Theodorsen解析解是基于無黏、無流動分離的勢流理論，但本文流動是真實的有黏空氣，且薄平板振動會導致流動分離，而黏性效應和流動分離可能隨折算風速的增大而增大。研究表明，本文CFD方法能較準確得到振動薄平板上的氣動力，數值方法是正確的。

4薄平板大幅扭轉強迫振動的能量特征

4.1振動能量的定義

在來流中繞扭轉軸作扭轉強迫振動的薄平板，定義薄平板表面的壓力系數為

為對比薄平板大幅振動的氣動力特征，下面給出了薄平板小幅扭轉振動（振幅0.2°）的氣動力時程。圖9（a）和（b）分別是折算風速Vr=4和24時，一個周期內的氣動力系數時程曲線。可見無論是升力還是扭矩系數時程，其曲線具有明顯的諧波特性，其波形頻率等于扭轉強迫振動頻率，可知當薄平板小幅扭轉振動時，其氣動力系統具有良好的線性特性;另外，從不同折算風速下的曲線來看，升力和扭矩問的相位差隨折算風速不同而有變化。

圖10為小幅扭轉振動薄平板在一個完整周期內4個關鍵時刻的壓力等高線，對應從平衡位置出發后達到最大正攻角位置（圖10（a）），然后回到平衡位置（圖10（b）），再達到最大負攻角位置（圖10（c）），最后回到零攻角的平衡位置（圖10（d））。可見從一個狀態到另外一個狀態，壓力等高線呈現有規律地變化，薄平板前緣為高風壓區，后緣為低風壓區。等高線分布規律的改變反映了薄平板的扭轉振動方向。

圖11（a）和（b）分別為在一個扭轉振動周期內，氣流在振動薄平板上表面各點、以及整個薄平板所作的無量綱總功。可見薄平板前緣（L.E.）局部為穩定區外，其后到轉動軸（A.R.）為不穩定區，也即氣流對薄平板作正功，薄平板從氣流中吸收能量，因此為激勵區;扭轉軸到薄平板后緣（T.E.）氣流做負功，因而是耗散氣流能量的穩定區，因此薄平板上激勵區與穩定區共存。圖11（b）給出了Vr=4-40，一個扭轉周期內氣流對薄平板所作的總功。可見在所計算的折算風速范圍內，一個周期內氣流對扭轉振動的薄平板均做負功，也即振動薄平板消耗氣流的能量，因此小幅扭轉振動的薄平板是穩定的。

4.3大幅扭轉振動

強迫薄平板發生大幅度的扭轉振動，如扭轉振幅角為20°，圖12（a）-（d）分別為折算風速Vr=4，24，44和68時，一個周期內的氣動力系數時程曲線。與圖9相比，可見隨著折算風速的提高，升力和扭矩系數時程逐漸失去諧波特性，特別是高折算風速。由于扭轉振動為正弦強迫位移，可見當薄平板大幅扭轉振動時，其氣動力系統不再為線性，且隨著折算風速的提高，其非線性特性越來越顯著。

圖13為不同折算風速下大幅扭轉振動的薄平板，在一個振動周期內氣流在振動薄平板上表面各點所作的無量綱總功。與圖11對比可見，在計算的所有折算風速范圍上，轉動軸到薄平板后緣氣流做負功，因而是薄平板上耗散氣流能量的穩定區，且在低折算風速（Vr=4-16）下，薄平板表面的能量消耗特征與小幅度扭轉振動類似。但當折算風速高于20后，迎風側的半個薄平板上每個點在一個周期內氣流做的功均為正，即為不穩定區，氣流對薄平板作正功，薄平板從氣流中吸收能量，因此薄平板上的激勵區與穩定區共存，而一個周期內氣流對整個薄平板表面所有點的總功之和決定了一個周期內氣流對薄平板所作的總功，該總功的正負反映了扭轉振動薄平板的氣動穩定性。

圖14給出了一個扭轉振動周期內氣流對薄平板所作總功隨折算風速的變化，可見當折算風速小于16時，一個周期內氣流對薄平板所作總功為負，也即薄平板消耗氣流的能量;但當折算風速大于20后，氣流對薄平板所作的總功由負變為正，此時一個周期內扭轉振動薄平板從氣流中吸收能量，這表明大幅度扭轉振動的薄平板已進入氣動不穩定狀態。

5結論

本文推導了二維流變區域的任意拉格朗日一歐拉描述方程，基于剛性薄平板簡化得到了薄平板強迫振動繞流場控制方程，采用數值方法開展了薄平板大幅強迫扭轉振動的氣動特征研究，可得到下述結論：

1）基于任意拉格朗日一歐拉描述的二維流體運動方程，通過建立計算域和物理域的隨時問變化的映射關系，可在不隨時問變化的二維計算域上，采用歐拉描述法實現對二維物理域上任意時變區域的流體流動的描述。對薄平板扭轉強迫振動，相應的映射關系和控制方程明顯簡化。

2）采用有限差分法開展了薄平板扭轉振動的CFD模擬。從小幅扭轉與大幅扭轉振動氣動力時程的對比發現，大幅扭轉振動氣動力呈現非線性特征，折算風速越高，氣動力的非線性特征越顯著。

3）大幅扭轉振動時，薄平板轉動軸到后緣是能量耗散區。但當折算風速高于20后，迎風側的半個薄平板從氣流中吸收能量，為氣動不穩定區;扭轉振動一個周期內氣流對整個薄平板作的總功由負變為正，此時薄平板從氣流中吸收能量，表明大幅扭轉強迫振動的薄平板已進入氣動不穩的狀態。