流動線性穩定性方程推導那點事兒(一)

引言

申請開通這個專欄有一段時間了,一直沒有更新。從專欄的名字可以看出來,本意是想寫一點與流動穩定性以及轉捩相關的內容。流動穩定性與轉捩雖然是流體力學中的一個很大的研究方向,但是相對從事CFD國內流體力學大同行來說,了解並且對其有所認識的同志還是數量有限的。特別是最近十年的時間,國內的轉捩相關的研究需求也逐漸擴張,實際工作中也確實發現很多大同行是有相關的基礎知識的需求的,為此打算利用業餘時間在這個專欄裡面寫點自己所了解的流動穩定性以及轉捩。我個人是希望更多的大同行以及廣大的本科生研究生對流動穩定性以及轉捩問題產生興趣的,並且介入到這個研究領域,大家共同來提高我國的研究水平。同時呢,寫這個的時候我也會盡量把問題說的詳細一點點,包括會設置一些練習題,有興趣的同志可以自己做做試試,算是很好的入門過程。再有呢,我也想通過這個渠道告訴大家,流動穩定性雖然看起來門檻略高,但是實際上還好,而且因為其相對比較倚重數學推導,因此它的邏輯性以及嚴謹性都非常漂亮,實際接觸下來是非常有趣的,有種從基本假設一層一層的建模、推導、分析與理解的過程中的學習快感,是非常非常好玩的事情,而且真心花點時間來研究的話,其實入門本身並不是很難的。所以這也是變相的給我們課題組,即天津大學機械學院高速空氣動力學研究室做個廣告,有興趣的同學可以來我們室交流或者讀研究生。

當然,我自身的認識水平也不大高,難免有所紕漏或者不對的地方,也請讀者在閱讀中能夠指出,大家一起探討,共同提高。

第一次呢,我想先寫一點關於線性穩定性方程的事情。本來么,理論上來說,應該先寫一個引論的,即闡明流動穩定性與轉捩預測到底是怎麼回事,它有哪些研究內容等等相關的內容。但是呢,因為是專欄,類似於筆記性質的合集的東西,所以往往是我想到哪裡,就寫到哪裡。為此,從科普講義的角度來看存在一些順序錯亂,還請見諒。適時的時候我會出一個閱讀目錄給出一個我所認為的合理的閱讀順序的。 好吧,廢話不說了,我先開始講講線性穩定性方程的推導那些事兒。


1. 以單波方程為例

我們先以單波方程為例。已知一維單波方程如下:

 frac {partial u} {partial t}+cfrac {partial u} {partial x}=0 	ag{1}

它的解 u(x, t) 是一個單波,可以寫作為 u(x, t)=hat u {
m{e}}^{i(alpha x-omega t)} 的形式,其中 hat u 我們叫它擾動的形狀。將其代入(1)式,可以簡化成:

 (omega-calpha) hat u=0 	ag{2}

顯然,方程若有非平凡解(即 hat u 
eq 0 ),則需要滿足 omega=calpha 。根據單波方程的性質,波傳播的相速度為 c ,這時該關係即對應於擾動波的色散關係(即波數與頻率之間的數學關係)。通常情況下,頻率 omega 和波數 alpha 二者知其一。因此,需要根據色散關係求出未知量。在這個求解過程中,方程(2)對應於求解一個特徵值問題,即根據方程係數行列式等於零確定色散關係,根據特徵向量來確定具體的基本擾動形狀 hat u

以上即為一個最簡單的問題模型。

2. 基本問題背後的數學概念

在明確了單波方程這一最為簡單的一維問題的情況後,我們來簡要分析一下該問題的基本數學概念。

從泛函問題的角度來考慮,則描述小擾動的任一線性偏微分方程,如果它的解可以設成單波形式,即滿足下式:

u(t, x_1, x_2, ldots, x_n)=hat u{
m e}^{i(-omega t + alpha_1 x_1 + alpha_2 x_2+cdots +alpha_n x_n)}.	ag{3}

其中, t 對應時間坐標, x_i 對應空間坐標,通常為三維問題。將該式帶入到原偏微分方程中,那麼則構成如下特徵值問題:

L(omega, alpha_1, alpha_2, ldots, alpha_n)hat u =0 	ag{4}

顯然,對應於非平凡解,需要滿足 left | L 
ight |=0 。該式參數 (omega, alpha_1, alpha_2, ldots, alpha_n) 中通常有一個未知,其它均為給定。通過求解該特徵值問題,即可求得這個未知量。這樣也就得到了基本的求解擾動波(即微分方程的解)對應的色散關係,並且同時根據特徵函數確定擾動的形狀向量 hat u

需要指出,對於該問題來說,如果偏微分方程每個坐標都可以寫成波數的形式從而對方程進行降維處理,那麼特徵向量 hat u 為一個值或者為一個向量,否則 hat u 即為未能寫成波數形式的那個坐標的函數,因此我們常稱 hat u 為形函數(或特徵函數)。

例如,如果對於擾動在 x_2 方向不能寫成波數形式,那麼它只能寫成如下形式:

u(t, x_1, x_2, ldots, x_n)=hat u(x_2) hat u{
m e}^{I(-omega t + alpha_1 x_1 + alpha_3 x_3+cdots +alpha_n x_n)}.	ag{5}

此時特徵值問題對應於

 L(omega, alpha_1, alpha_3, ldots, alpha_n)hat u(x_2) =0 	ag{6}

顯然,這裡特徵值問題的解矢量對應於 hat u(x_2) 是一個特徵函數。

所以,從數學本質上講,所謂線性穩定性問題對應於一個特徵值問題的解。因此,其解集對應的就是 L 運算元所對應的子空間。此時特徵值問題中 L 運算元的特徵值所滿足的方程 |L|=0 對應於擾動的色散關係,而特徵值問題的特徵向量則對應於子空間的線性無關向量。也就是說,子空間內的任一元素都可以用用特徵值問題給出的這些特徵函數的線性組合來表示。通常情況下會將這些線性無關的特徵函數歸一化後作為基向量使用。

註:由於對應流動線性穩定性問題,通常情況下運算元 L 不是自伴隨的,因此這些線性無關的特徵函數雖然線性無關但互相之間並不正交。這也是進一步的對流動線性穩定性問題存在伴隨問題的數學依據。

運算元 L 不是自伴隨對應於如果 L 是矩陣的話,則 L 並非Hermitian矩陣(即共軛轉置是它本身的矩陣)

當寫成矩陣形式時,大多數時候未必能寫成經典特徵值問題 Ax=lambda x 的形式,這時依然可以寫成廣義特徵值 Ax=lambda Bx 的形式。

3. 流體力學線性穩定性方程的一般推導方法

以下介紹流體力學線性穩定性方程的一般推導方法

3.1 小擾動方程

以下先以一維Burgers方程為例進行推導。

設某問題存在以下控制方程:

 frac {partial u} {partial t}+ufrac {partial u} {partial x}=0	ag{7}

那麼,我們可以先認為瞬時狀態下的 u 可以寫成如下形式: u=U_0+u^prime 。其中 u 為瞬時量, U_0 為基本狀態下的量(對應流動穩定性問題則為待分析動力穩定性的基本流動),而 u^prime 為擾動量。顯然,對應於線性動力系統即線性穩定性問題時, u^prime 是小量。因此,我們也可以將瞬時狀態下做如下分解以利整理:

 u=U_0+epsilon u^prime 	ag{8}

這裡, epsilon 是小參數。我們將上式代入方程(7),那麼就有:

 frac {partial{(U_0+epsilon u^prime)}} {partial t}+ (U_0+epsilon u^prime) frac {partial{(U_0+epsilon u^prime)}} {partial x}=0 	ag{9}

對該方程整理,顯然有:

 frac {partial U_0 } {partial t}+U_0frac {partial U_0} {partial x}+ epsilonleft( frac {partial u^prime } {partial t}+ U_0frac {partial u^prime} {partial x}+ u^primefrac {partial U_0} {partial x}
ight)+ epsilon^2 u^primefrac{partial u^prime}{partial x}=0	ag{10}

考慮 epsilon 為任一小量。式(10)忽略掉二階小量,剩餘部分有如下形式:

frac {partial U_0 } {partial t}+U_0frac {partial U_0} {partial x}=0 label{u0}	ag{11}

 frac {partial u^prime } {partial t}+ U_0frac {partial u^prime} {partial x}+ frac {partial U_0} {partial x} u^prime=0label{u}	ag{12}

方程(11)描述了基本狀態量的數學關係。事實上對於基本狀態量 U_0 ,顯然它依然滿足於基本方程(7),兩個方程形式完全相同。而對於方程(12),它描述了一階擾動量的數學關係,即是線性擾動方程,又稱小擾動方程。

3.2 線性穩定性方程

得到小擾動方程之後,我們就可以考慮解的形式,即將解寫成(3)的形式。但是這裡有一個問題——小擾動解的形式實際上是與 U_0 的特徵有關係的。我們注意到(3)的形式中暗含著波動解的概念,因此就需要看在什麼坐標下解可以寫作為波動的形式。

首先,通常 U_0 是定常的,因此小擾動在時間方向上是波動的,這樣時間項就可以寫成e指數的形式,對應擾動波的頻率。而在空間上,就需要考慮 U_0 的變化了,即小擾動是否可以在空間上寫成波動函數的形式。

這裡我們引入一個假設,即平行性假設。我們認為小擾動在某一方向上存在變化尺度,即在該方向上的波長 lambda 。而基本流 U_0 在同樣的方向上也存在著一個變化尺度,稱之為 l 。如果在某一方向上 l>>lambda ,那麼就可以認為相對於擾動波的變化,基本流流場的變化可以忽略,而此時擾動在該方向上就可以寫作為波動函數的形式。顯然這是一種假設,因為除非基本流 U_0 在某一方向上絕對平行(即不變),否則波動函數的形式一定存在著誤差。我們將這種假設稱之為平行性假設。相反,如果基本流 U_0 在某一方向上是快變的,即基本流變化尺度與擾動波尺度相當,那麼就不能在該方向上寫成擾動波函數的形式,即只能採用方程(5)的形式。這時就要將該方向寫在特徵函數的坐標中,並不寫到e指數所表徵的相位上去。

還是以3.1節的Burgers方程為例,如果小擾動解 u^primex 方向上相對 U_0 是慢變的,或者

U_0x 方向上不變,那麼方程(12)解的形式則可以寫成:

 u^prime=hat u {
m{e}}^{i(alpha x-omega t)};

否則,如果 U_0x 方向上也是快變的,那麼小擾動解只能寫成:

 u^prime=hat u(x) {
m{e}}^{-iomega t}.

同樣的,與此類似,對於二維平板邊界層問題,由於其在流向方向 x 上是緩慢增長的,而在垂直壁面方向 y 上變化比較劇烈,顯然,此時邊界層內的擾動就要寫成:

u^prime=hat u(y) {
m{e}}^{alpha x+ eta z-iomega t}.

這樣,我們就得到了小擾動解的基本形式。接下來只需要將其帶入到小擾動方程中進行整理,就可以寫成特徵值問題的形式了。這樣也就得到了線性穩定性方程。

以上是用標量方程舉例。如果未知數為矢量,其推導過程和方法是完全一樣的。這時方程化簡為如下特徵值問題:

Lhatphi=0

其中, L=L(omega, alpha_1, alpha_2,ldots) 。對應於特徵值問題, |L|=0 給出色散關係,特徵函數給出

擾動波的形狀。

習題1:試推導不可壓縮邊界層方程

 
abla cdot mathbf{u}=0 \ 
hofrac{Dmathbf{u}}{Dt}=-
abla p+frac{1}{Re}Delta mathbf{u}

的線性穩定性方程。

(提示:未知數為 phi^prime=(	extbf{u}^prime, p^prime)^T 。)

習題2: 利用二維Euler方程推導可壓縮均勻來流中的小擾動形式。

(提示:設來流馬赫數為$M$,特徵值問題的特徵值為 u+a , u-a 以及 u )

推薦閱讀:

TAG:流體力學 | 計算流體力學CFD | 湍流 |