范景威 周偉燦 馮也騁 官元紅
臺(tái)風(fēng)強(qiáng)度變化是臺(tái)風(fēng)研究的一個(gè)重要課題,是當(dāng)前臺(tái)風(fēng)觀測(cè)和預(yù)報(bào)中的重難點(diǎn)與前沿問(wèn)題.影響臺(tái)風(fēng)強(qiáng)度變化的因子可以歸結(jié)為以下三類:環(huán)境氣流、下墊面和臺(tái)風(fēng)內(nèi)部結(jié)構(gòu)[1-2].隨著監(jiān)測(cè)手段的發(fā)展以及數(shù)值模式理論和技術(shù)的不斷成熟,國(guó)內(nèi)外學(xué)者針對(duì)臺(tái)風(fēng)強(qiáng)度變化的研究取得了一定成果,已有研究?jī)?nèi)容涉及諸多方面,包括環(huán)境風(fēng)垂直切變、冷空氣、高低空急流、西風(fēng)槽、地形以及海面溫度和噴沫等[3-9].
變分法是研究泛函極值的數(shù)學(xué)方法,其主要思想是構(gòu)造適當(dāng)?shù)姆汉?對(duì)泛函作變分并取變分為零,研究此時(shí)滿足的函數(shù)性態(tài).變分法廣泛地應(yīng)用于力學(xué)、數(shù)學(xué)物理反問(wèn)題、氣象資料同化和大氣運(yùn)動(dòng)的穩(wěn)定性問(wèn)題中.Arnold[10]和Vallis[11]運(yùn)用變分原理,指出流體運(yùn)動(dòng)和地轉(zhuǎn)平衡分別是在位渦守恒及勢(shì)能守恒的約束下,使總能量達(dá)到極?。闃s生[12]利用Finlayson[13]提出的限制性變分原理,求得了與非線性渦度方程對(duì)應(yīng)的Lagrange函數(shù)的近似表達(dá)式.Barth等[14]從線性淺水方程和系統(tǒng)能量變分出發(fā),推導(dǎo)了從高程和深度平均速度中去除慣性重力波的方法,指出變分濾波器在提前降低重力波方面的優(yōu)勢(shì).黃思訓(xùn)等[15]通過(guò)二次變分方法研究臺(tái)風(fēng)流場(chǎng)結(jié)構(gòu),將實(shí)測(cè)風(fēng)場(chǎng)分解成無(wú)旋場(chǎng)和有旋場(chǎng),再將有旋流場(chǎng)分解成對(duì)稱渦旋和非對(duì)稱渦旋對(duì),兩次提取到的都是最大渦旋.Wang等[16]用變分法推導(dǎo)了赤道電離層Rayleigh-Taylor(RT)穩(wěn)定和不穩(wěn)定的充分條件,并通過(guò)系統(tǒng)特征值計(jì)算了RT穩(wěn)定和不穩(wěn)定區(qū)域.Badin等[17]將變分原理應(yīng)用于波動(dòng)動(dòng)力學(xué)中,得到地表水波方程以及非線性問(wèn)題的波頻散關(guān)系.近年來(lái),隨著變分同化技術(shù)的發(fā)展與成熟,它被廣泛地應(yīng)用于大氣和海洋等不同領(lǐng)域[18-21],其目的是充分利用氣象四維觀測(cè)資料,在動(dòng)力預(yù)報(bào)模式中通過(guò)變分得到大氣海洋最優(yōu)狀態(tài)估計(jì),以獲得更好的模式初始場(chǎng)及預(yù)報(bào)效果.
目前對(duì)臺(tái)風(fēng)強(qiáng)度的研究多采用統(tǒng)計(jì)、診斷分析、數(shù)值模擬和敏感性試驗(yàn)等方法,缺少理論層面的研究.此外,前人的研究多采用直接分解的方法將臺(tái)風(fēng)流場(chǎng)分解成渦旋流場(chǎng)和無(wú)旋流場(chǎng)[15],這種分解無(wú)法得到有旋場(chǎng)中的最大渦旋.本文通過(guò)能量泛函變分,研究臺(tái)風(fēng)強(qiáng)度變化率最大時(shí)臺(tái)風(fēng)能量滿足的關(guān)系,在四力平衡關(guān)系的基礎(chǔ)上提取到有旋流場(chǎng)中的最大渦旋,給出臺(tái)風(fēng)強(qiáng)度變化最快時(shí)渦度的解析解.該結(jié)果為預(yù)報(bào)方程提供了一定的理論指引,對(duì)研究臺(tái)風(fēng)發(fā)展過(guò)程中速度的變化趨勢(shì)和臺(tái)風(fēng)的層次結(jié)構(gòu)具有實(shí)用價(jià)值.
假設(shè)大氣是均勻不可壓縮的流體,由連續(xù)方程有
(1)
其中,u=u(x,y,z,t)為緯向風(fēng)速,v=v(x,y,z,t)為經(jīng)向風(fēng)速,w=w(x,y,z,t)為垂直風(fēng)速.
大氣運(yùn)動(dòng)遵循牛頓第二運(yùn)動(dòng)定律.對(duì)于慣性坐標(biāo)系,牛頓第二運(yùn)動(dòng)定律可表示為
(2)
(3)
其中,V3=(u,v,w)為三維速度矢量,Ω為地球自轉(zhuǎn)角速度矢量,R為氣塊所在的緯圈平面內(nèi)從地軸到該氣塊的距離矢量,其大小為R=rcosφ.作用于空氣微團(tuán)的氣壓梯度力、分子黏性力、湍流黏性力和重力的表達(dá)式分別為
(4)
(5)
(6)
(7)
其中,ρ為密度,p=p(x,y,z,t)為氣壓,3為三維微分算子,μ為分子黏性系數(shù),為三維拉普拉斯算子,Az為湍流交換系數(shù),G為引力常數(shù),M為地球的質(zhì)量,r為空氣微團(tuán)的位置矢量.
將式(3)—(7)代入式(2),經(jīng)整理得到旋轉(zhuǎn)坐標(biāo)系中矢量形式的動(dòng)量方程:
(8)
其中,F=F2+F3為分子黏性力與湍流黏性力之和,即摩擦力.用V3點(diǎn)乘式(8)兩端得到旋轉(zhuǎn)坐標(biāo)系中的動(dòng)能方程:
(9)
定義能量泛函I(u,v,w,p)∈C2(Ω),
(10)
其中,Ω是任意臺(tái)風(fēng)區(qū)域,用該泛函表示臺(tái)風(fēng)強(qiáng)度變化率.假設(shè)在?Ω上u,v,w,p以及u,v,w關(guān)于x,y,z的所有一階偏導(dǎo)數(shù)均為常數(shù),此時(shí)臺(tái)風(fēng)在區(qū)域邊界上可視為基本氣流[22].對(duì)I作變分,有
(11)
令δI=0,利用δu,δv,δw的任意性得到Euler-Lagrange方程為
(12)
(13)
(14)
矢量形式為
2F+fp-g-3K=0,
(15)
其中,fp為氣壓梯度力,g=(0,0,g).式(15)表明,當(dāng)臺(tái)風(fēng)強(qiáng)度變化率達(dá)到最大時(shí),摩擦力、氣壓梯度力、重力和動(dòng)能梯度滿足四力平衡.因此,通過(guò)這四個(gè)力確定的向量可以作為臺(tái)風(fēng)強(qiáng)度的預(yù)報(bào)因子,更準(zhǔn)確地確定臺(tái)風(fēng)系統(tǒng)強(qiáng)度變化率達(dá)到最大的時(shí)間點(diǎn),研究臺(tái)風(fēng)強(qiáng)度的變化規(guī)律.此外,將式(15)變形成
(16)
則臺(tái)風(fēng)強(qiáng)度變化最快時(shí)的摩擦力表達(dá)式由式(16)給出,它為計(jì)算較困難且常常被忽略的摩擦力提供了一種新的計(jì)算方法.
為了進(jìn)一步研究臺(tái)風(fēng)強(qiáng)度變化最快時(shí)渦度的形式和流場(chǎng)的變化,采用黃思訓(xùn)等[15]的方法對(duì)臺(tái)風(fēng)風(fēng)場(chǎng)進(jìn)行分解.首先將風(fēng)場(chǎng)u(x,y,z),v(x,y,z)分解成平面場(chǎng)u′(x,y),v′(x,y)和垂直場(chǎng)H(z),即
u(x,y,z)=u′(x,y)H(z),
(17)
v(x,y,z)=v′(x,y)H(z),
(18)
其中,風(fēng)場(chǎng)u,v對(duì)應(yīng)臺(tái)風(fēng)強(qiáng)度變化率最大的時(shí)刻,分別滿足方程(12)和(13).由于臺(tái)風(fēng)的渦旋運(yùn)動(dòng)受環(huán)境氣流引導(dǎo),因此臺(tái)風(fēng)流場(chǎng)可以分解成有旋流場(chǎng)和無(wú)旋流場(chǎng).利用變分方法從風(fēng)場(chǎng)u′,v′中提取無(wú)旋氣流u1(x,y),v1(x,y)[15],使得泛函
(v′-v1)2]dΩ=min!
(19)
且u′,v′,u1,v1滿足
(20)
(21)
(22)
這是一個(gè)條件變分問(wèn)題,引入Lagrange乘子λ1(x,y,z),λ2(x,y,z),λ3(x,y,z),有
(23)
于是
(24)
假設(shè)在?Ω上λ1,λ2,λ3以及λ2,λ3關(guān)于x,y,z的所有一階偏導(dǎo)數(shù)均為0,則有
(25)
利用δu,δv,δw,δu1,δv的任意性,可得它們對(duì)應(yīng)的Euler-Lagrange方程分別為
(26)
(27)
(28)
(29)
(30)
(31)
(32)
(33)
(34)
(35)
根據(jù)式(29)和(30),有旋流場(chǎng)(u2,v2)和無(wú)旋流場(chǎng)(u1,v1)的表達(dá)形式分別為
(36)
(37)
(38)
(39)
在有旋場(chǎng)中引入流函數(shù)ψ(x,y),滿足
(40)
渦度為
(41)
將式(30)對(duì)x作微商,式(29)對(duì)y作微商,相減并結(jié)合式(33)得
(42)
將式(35)對(duì)x作微商,式(34)對(duì)y作微商,相減并結(jié)合式(42)得
(43)
將式(33)、(41)和(42)聯(lián)立得
ψ=λ1.
(44)
因此,在平面風(fēng)場(chǎng)中,二維Lagrange乘子λ1(x,y)即為有旋流場(chǎng)的流函數(shù)ψ(x,y).故方程(43)可改寫為
(45)
由于ζ=ζ(x,y),H=H(z),對(duì)方程(45)采用分離變量法,并設(shè):
(46)
則:
(47)
其中,λ為常數(shù).
由式(31)得
(48)
將式(46)、(48)代入式(26)、(27)得
(49)
(50)
給定方程(46)的邊界條件為H(0)=H,H(z1)=H1,方程(46)的解為
(51)
(52)
根據(jù)黃思訓(xùn)等[15],臺(tái)風(fēng)流場(chǎng)一般在圓域中經(jīng)變分分解能提取到最大渦旋,故將方程(47)轉(zhuǎn)化成極坐標(biāo)(r,θ)下的方程
0 (53) 為滿足齊次邊界條件,作變換 (54) 則η(0,θ)=η(r0,θ)=0.記ζ(r0,θ)=ζ0,方程(53)轉(zhuǎn)化成如下形式: (55) 令η(r,θ)=R(r)Φ(θ),ζ(r,θ)=ρ(r)Φ(θ),則有 (56) (57) 對(duì)方程(57)作變量分離得 (58) 其中,μ為常數(shù).結(jié)合η(r,θ)關(guān)于θ的周期性,可得 (59) (60) 式(59)是帶有周期性條件的特征值問(wèn)題,其特征值和特征函數(shù)分別為 μn=n2,Φn(θ)={cosnθ,sinnθ}. (61) 其中,n≥0.它的解為 (62) 其中,an,bn為常數(shù). (63) 當(dāng)λ=0時(shí),方程(63)為非齊次歐拉方程,它在邊界條件下的解為 (64) 此時(shí),結(jié)合式(56)得渦度ζ的表達(dá)式為 (65) 其中,A′n,B′n,c0為常數(shù). 當(dāng)λ>0時(shí),采用常數(shù)變易法[24]得方程(63)的通解為 (66) (67) 由邊界條件得 (68) (69) 因此,臺(tái)風(fēng)強(qiáng)度變化率最大時(shí)有旋流場(chǎng)渦度的解析解為 (70) 渦度擬能的表達(dá)式為 (71) 當(dāng)λ=0時(shí),流函數(shù)ψ(r,θ)滿足: (72) 根據(jù)ζ(r,θ)的表達(dá)式,ψ(r,θ)具有如下形式: E(r)sin2nθ+F(r)cos2nθ+G(r)sinnθcosnθ. (73) 由λ1|?Ω=0可知ψ|?Ω=0,故方程(72)的邊界條件為 C(r0)=D(r0)=E(r0)=F(r0)=G(r0)=0. (74) 將式(73)代入方程(72),整理得 (75) (76) (77) (78) (79) 方程(75)、(76)和(79)滿足邊界條件的解分別為 (80) (81) G(r)=0, (82) (83) (84) 且根據(jù)式(74),邊界條件為H(r0)=J(r0)=0.解得 (85) J(r)=0. (86) 由式(85)和(86)得 (87) 因此,流函數(shù)ψ(r,θ)的表達(dá)式為 (88) 根據(jù)式(40),有旋場(chǎng)u2(r,θ),v2(r,θ)滿足: (89) (90) 因此,有旋流場(chǎng)(u2,v2)的表達(dá)式為 (91) (92) 由變分的性質(zhì)可知,分解得到的有旋氣流(u2,v2)是最大的渦旋. 當(dāng)λ>0時(shí),為方便求解流函數(shù),對(duì)渦度ζ作如下簡(jiǎn)化: (93) (94) 同理,ψ(r,θ)具有如下形式: (95) 邊界條件為 H(r0)=M(r0)=0. (96) 將式(95)代入方程(94),整理得 (97) (98) H(r)=anGnr3-an(Gnr03-n+Fnr02n+2)rn+anFnr3n+2, (99) M(r)=bnGnr3-bn(Gnr03-n+Fnr02n+2)rn+bnFnr3n+2. (100) 因此,流函數(shù)ψ(r,θ)的表達(dá)式為 (101) 根據(jù)式(89)、(90),有旋流場(chǎng)(u2,v2)的表達(dá)式為 (3n+2)Fnr3n+1]sinθsinnθ+ an[-3Gnr2+n(Gnr03-n+Fnr02n+2)rn-1- (3n+2)Fnr3n+1]sinθcosnθ+ an[nGnr2-n(Gnr03-n+Fnr02n+2)rn-1+ nFnr3n+1]cosθsinnθ+ bn[-nGnr2+n(Gnr03-n+Fnr02n+2)rn-1- nFnr3n+1]cosθcosnθ, (102) nFnr3n+1]sinθsinnθ+bn[-nGnr2+n(Gnr03-n+ Fnr02n+2)rn-1-nFnr3n+1]sinθcosnθ+ bn[3Gnr2-n(Gnr03-n+Fnr02n+2)rn-1+ (3n+2)Fnr3n+1]cosθsinnθ+an[3Gnr2- n(Gnr03-n+Fnr02n+2)rn-1+ (3n+2)Fnr3n+1]cosθcosnθ. (103) 當(dāng)λ=0時(shí),根據(jù)有旋流場(chǎng)的表達(dá)式(91)、(92),取n=1,r0=1,A1=1,B1=1,c0=16,ρ(1)=2,臺(tái)風(fēng)強(qiáng)度變化最快時(shí)的流線圖如圖1所示.可以看出,流線圖在經(jīng)向和緯向上呈均勻梯度分布,符合實(shí)際流場(chǎng)分布[15].因此,在有旋流場(chǎng)的表達(dá)式(91)、(92)中,取n=1可得有旋流場(chǎng)及其渦度的解析解. 圖1 λ=0,n=1時(shí)臺(tái)風(fēng)強(qiáng)度變化率最大時(shí)刻的流線圖Fig.1 Plot of streamline at the maximum change rate of TC intensity when λ=0,n=1 圖2 λ>0,n=0時(shí)臺(tái)風(fēng)強(qiáng)度變化率最大時(shí)刻的流線圖Fig.2 Plot of streamline at the maximum change rate of TC intensity when λ>0,n=0 因此,臺(tái)風(fēng)強(qiáng)度變化最快時(shí)渦度、渦度擬能和流場(chǎng)的一個(gè)解析解分別為 (104) (105) (106) (107) 有旋流場(chǎng)(u2,v2)更直觀清晰地反映臺(tái)風(fēng)每一層的流場(chǎng)結(jié)構(gòu),渦度和渦度擬能的變化有助于研究臺(tái)風(fēng)的運(yùn)動(dòng)機(jī)理和風(fēng)場(chǎng)中能量的變化.臺(tái)風(fēng)的整體結(jié)構(gòu)取決于垂直廓線H(z),而λ是決定H(z)的主要參數(shù),因此λ決定整個(gè)臺(tái)風(fēng)的層次結(jié)構(gòu). 本文從理論層面研究臺(tái)風(fēng)強(qiáng)度變化率最大時(shí)臺(tái)風(fēng)能量滿足的關(guān)系,給出臺(tái)風(fēng)強(qiáng)度變化最快時(shí)渦度的解析解.結(jié)果表明,當(dāng)臺(tái)風(fēng)強(qiáng)度變化率達(dá)到最大時(shí),摩擦力、氣壓梯度力、重力和動(dòng)能梯度滿足四力平衡.通過(guò)這四個(gè)力確定的向量可以作為臺(tái)風(fēng)強(qiáng)度的預(yù)報(bào)因子,更準(zhǔn)確地確定臺(tái)風(fēng)系統(tǒng)強(qiáng)度變化率達(dá)到最大的時(shí)間點(diǎn),研究臺(tái)風(fēng)強(qiáng)度的變化規(guī)律.在四力平衡關(guān)系的基礎(chǔ)上,進(jìn)一步通過(guò)風(fēng)場(chǎng)變分分解提取到有旋場(chǎng)中的最大渦旋,得到臺(tái)風(fēng)強(qiáng)度變化最快時(shí)渦度和流場(chǎng)的一個(gè)解析解,從理論上證明了臺(tái)風(fēng)發(fā)展最快時(shí)前人的研究[15]結(jié)果,對(duì)研究臺(tái)風(fēng)發(fā)展過(guò)程尤其是平衡過(guò)程中的運(yùn)動(dòng)機(jī)理和臺(tái)風(fēng)的空間結(jié)構(gòu)具有一定的指導(dǎo)意義,為臺(tái)風(fēng)路徑和強(qiáng)度預(yù)報(bào)提供了一定的理論指引.此外,本文的結(jié)論為計(jì)算較困難且常常被忽略的摩擦力提供了一種新的計(jì)算方法.2.3 有旋流場(chǎng)的解析解
2.4 數(shù)值試驗(yàn)
3 結(jié)論