相容拉格朗日-歐拉法求解黏性流體中彈性圓柱殼的振動

郝亞娟， 郭茜茜， 陳佳慧

(燕山大學 理學院， 河北 秦皇島 066004)

近年來，在海洋資源開發利用等方面的需求驅動下，許多科研人員著眼于研究水下生物的運動機理，指導水下機器人[1]的研究，所以仿生學[2]的研究具有了廣闊的應用前景和巨大潛在價值，涉及到流體力學、機械、材料、控制、生物等學科。

根據Lighthill[3]對水中生物運動方式的分類，波狀擺動推進是其中的一種，該推進方式生物身體做橫向扭曲、往復擺動，以橫波的方式由前向后或逆向傳播。早在20世紀50年代，Taylor[4]采用“靜態流體理論”分析計算微生物運動的流體力，該方法忽略慣性力影響，適用于雷諾數比較低的情況。Lighthill[5]提出一種應用于變形體的“細長體理論”，該理論指出細長魚類獲得較高推進效率的條件。隨后Wu[6]首先提出了“二維波動板理論”，該理論將魚體看作一塊彈性薄板，分析了二維柔性波動板的游動過程。Chopra等[7]又將“二維波動板理論”進行了擴展，提出了可用于大擺幅推進系統的理論。Triantafyllou[8]通過觀測證實了二維和三維波動板理論的真實性。Tian等[9]用不可壓縮的納維-斯托克斯方程數值解來了解行波表面推進機制，用靈活的箔來模擬行波運動游泳體模型，分析向前推進速度和流場。

國內外學者對相關理論開展廣泛的研究，取得了一些成果，但仍存在不足之處，對水中生物的推進機理還需進一步深入研究。本文嘗試將水中生物看作彈性圓柱殼體，假設圓柱殼長度足夠長并遠遠大于半徑，將其近似為無限長。另外假設速度很小，雷諾數Re<0.1，屬于低雷諾數流動，采用斯托克斯方程近似，這樣運動學問題和動力學問題相互獨立，因此速度與流體的黏性無關[10]。采用相容拉格朗日-歐拉法[11]求解無限長圓柱殼位于黏性流體中，當殼體表面發生行波振動，黏性流體的運動速度以及殼體的推進速度。相容拉格朗日-歐拉法的優勢在于，殼體采用拉格朗日法描述，流體采用歐拉法描述，在接觸面結合這兩種方法，這樣可直接利用流體力學和固體力學中的基本方程。相容拉格朗日-歐拉法求解彈性體理想流體繞流問題顯示了方法的有效性[12-13]，本文嘗試采用該方法研究黏性流體問題。

1 基本方程

對于不可滲透表面，應用黏性流體分子對表面的黏附條件和壓力矢量平衡條件，可以將表面接觸條件簡化。“相容”黏附條件簡化為

(1)

式中：u是彈性體固定點的位移矢量;V為與變形后彈性體重合的流體空間點處的速度矢量;t為時間。

在圓柱坐標系(r,θ,z)中，變形后空間點處沿單位矢量kl(l=r,θ,z)方向的速度矢量分量近似值Vl(l=r,θ,z)用變形前空間點處速度矢量和位移矢量沿坐標軸方向的分量vl和ul(l=r,θ,z)的泰勒級數展開式近似表示

(2)

由于殼的位移分量太小，故泰勒級數展開式的第三項不予考慮。

由式(1)得接觸面的運動條件為

(3)

Vz、Vθ和Vr的值可由式(2)計算得到。

對于不可壓縮黏性流體的納維-斯托克斯方程，在圓柱坐標系中形式如下[14]

(4)

其中

式中：μ表示流體動力黏性系數；ρ表示黏性流體質量密度，p為單位質量流體壓力；Fr、Fθ和Fz表示單位質量的質量力在三個坐標軸的投影。

假設流體速度很小，雷諾數Re< 0.1，式(4)左邊的慣性力同式(4)右邊的第三項黏性力相比可以忽略，而且運動是定常的，質量力可以略去，則與θ無關的二維流體運動的納維-斯托克方程簡化為斯托克斯方程[15]

(5)

其中

相應的連續方程為

(6)

2 問題的解

圓柱殼體表面行波振動形式如下

(7)

式中：α=2π/λ，ω=cα，λ和c為波的長度及速度，Zn和Rn為振幅。

引入流函數ψ滿足

rvz=?ψ/?r，rvr=-?ψ/?z

(8)

式(8)滿足連續方程式(6)。

由式(2)、式(3)和式(8)，接觸面運動方程式(3)可以寫成如下形式

(9)

(10)

式(5)中消去壓力p得到如下方程

(11)

方程式(11)具有如下形式的解

BnrK1(nαr)]sinn(αz-ωt)

(12)

該解滿足條件

vz=V0(r→∞)

式中：V0是與時間無關、相對于不動坐標系無窮遠處流體的運動速度，由給定問題及系數An和Bn同時確定，K0(x)，K1(x)是麥克唐納函數[16]。也就是說，相對于無窮遠處不動的流體而言，殼的推進速度等于-V0。

將式(7)和式(12)代入式(9)和式(10)，并利用如下關系式[16]

(13)

可得

(14)

(15)

(16)

由式(14)和式(15)得

(17)

(18)

由式(16) 、式(17)和式(18)可得

(19)

(20)

其中

(21)

由式(16)、式(19)和式(20)得速度V0

(22)

流體的速度vz和vr分別為

(23)

(24)

當nαR為大值時，利用函數漸進值

由式(22)，令圓柱殼的半徑R→∞，可以得到薄板振動時的無窮遠處流體速度V0為

(25)

3 結果分析

經編程分析可知，N取值不同，速度的具體數值不同，但變化規律一致，振幅和頻率對于速度的影響是顯著的，隨著振幅和頻率的增加，速度增加，而波長的變化對速度的增加幾乎沒有影響，該結果與文獻[11]一致。以下取N=1計算。

圖1 隨R增大圓柱殼與薄板值比較

對于彈性薄板，取N=1，表1給出了Z1和R1不同取值時V0/c的數值。從表1可以知道，純橫向振動(Z1=0，R1≠0)時無窮遠處流體的運動速度與波的方向相同，而純縱向振動(Z1≠0，R1=0)時二者方向相反，即純橫向振動時薄板的位移和波的方向相反，純縱向振動時二者方向相同。若兩種振動均有，由式(25)可知，當R1的數值位于-2.4Z1與0.4Z1之間時，無窮遠處流體的運動速度與波的方向相反，否則方向相同。

表1 對于薄板Z1和R1不同取值時的V0/c數值

表2 r=1.5R處和的最大值

4 結 論

(1) 接觸面條件中變形后的速度分量用變形前分量的泰勒級數展開式近似表示，簡單明了，適用性強；

表3 不同位置處的最大值

(2) 橫向振動對于黏性流體速度的影響是主要因素，縱向振動對結果影響很??；

(3) 圓柱殼表面做純縱向振動時，與時間無關的無窮遠處流體的運動速度與圓柱殼半徑無關；

(4) 彈性圓柱殼半徑較大(薄板)做純橫向振動時殼的位移方向與波的方向相反，純縱向振動時二者方向相同，兩種形式的振動均有，則方向可能相同，也可能相反。