王光輝，谷湘潛

（中國(guó)氣象科學(xué)研究院 災(zāi)害天氣國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100081）

粘塑性海冰流變模型的穩(wěn)定性分析

王光輝，谷湘潛

（中國(guó)氣象科學(xué)研究院 災(zāi)害天氣國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100081）

對(duì)Hibler 1977提出的具有粘塑性本構(gòu)關(guān)系的海冰模型進(jìn)行了修改和穩(wěn)定性分析。在已有的工作中，應(yīng)用這樣一個(gè)粘塑性本構(gòu)關(guān)系進(jìn)行數(shù)值積分時(shí)大多都必須引入人工擴(kuò)散項(xiàng)。這里我們用數(shù)值逼近技巧提出計(jì)算相鄰浮冰之間相互作用力的新方法，并且證明修改后的模型無(wú)論是海冰的收斂還是發(fā)散流動(dòng)都是穩(wěn)定的，模式積分無(wú)需引入人工擴(kuò)散項(xiàng)。以二維結(jié)構(gòu)為例給出了詳細(xì)分析。

海冰；流變模型；穩(wěn)定性；橢圓屈服曲線；插值

1 前 言

海冰覆蓋地球海洋的大約7%，并且在這些海洋和大氣之間形成了盡管是部分但有時(shí)非常有效的屏障，阻止海－氣之間的相互交換。極地水冷卻到冰點(diǎn)形成海冰，其大小從浮在水中的小冰粒到幾米厚的大冰塊。雖然通常情況下海冰是結(jié)實(shí)的， 但冰層從來(lái)不會(huì)是完整的，因?yàn)楹Ｑ笾械木植考訜岷捅牟煌较虻倪\(yùn)動(dòng)，冰破裂露出水面來(lái)。冰經(jīng)常變形，將均勻一大片冰塊分裂成許多形狀不規(guī)則的大浮冰。與海洋大氣分布狀態(tài)明顯相反的是海冰的分布有非常大的季節(jié)和年度變化，因而對(duì)氣候模擬研究增加了復(fù)雜性[1]。長(zhǎng)期以來(lái)，人們認(rèn)識(shí)到在海洋，大氣和冰之間存在許多聯(lián)系。一方面，無(wú)論是海洋的上層還是海冰的覆蓋都將受到大氣有意義的影響。另一方面，海冰的存在有大量的氣候因果關(guān)系，影響大氣和海洋的溫度與環(huán)流形式。海冰減少吸收在海洋表面的太陽(yáng)輻射量。它是一個(gè)強(qiáng)大的絕緣體，限制著海洋與大氣之間熱量、質(zhì)量，動(dòng)量及化學(xué)要素的交換。此外，在全球氣候模式里，對(duì)海冰動(dòng)力過(guò)程的理解為大氣提供了臨界邊界條件。在海冰模擬中有兩個(gè)重要的變量，即：冰蓋率和冰的移速。冰蓋率對(duì)海洋和大氣之間的熱量交換有較大的作用。冰的移速對(duì)海洋冰的質(zhì)量，海溫及鹽度的重新分布有重要的影響。在收斂流動(dòng)中，相鄰的浮冰相互作用并向?qū)Ψ绞┘咏佑|力。而在發(fā)散流場(chǎng)內(nèi)浮冰彼此分離，接觸力消失。在維持聚合流動(dòng)期間，相鄰浮冰間不斷增加的接觸力通過(guò)破裂，變形和垂直位移形成在上面漂浮，在下傾覆的山脊式粘合冰。為了確定彈性響應(yīng)，應(yīng)變將被不確定的追蹤，呈現(xiàn)出豐富的理論[2]和數(shù)值問(wèn)題[3]。Hibler于1977年[4]通過(guò)在合適的時(shí)空尺度上取應(yīng)力和應(yīng)變率的平均值，提出了一種新的粘塑性流變理論。緊隨著這個(gè)工作，Hibler 1979年提出了動(dòng)力熱力學(xué)海冰模型[5]，在這個(gè)模型里，厚冰和薄冰被區(qū)別對(duì)待，厚冰的區(qū)域部分定義冰的聚合。這個(gè)理論作為一種經(jīng)典的海冰動(dòng)力學(xué)模型被廣泛接受。數(shù)值算法[6]形成了大多數(shù)其它海冰動(dòng)力模型的基礎(chǔ)。雖然計(jì)算的代價(jià)高昂，但仍有越來(lái)越多的近代模型包含復(fù)雜的流變學(xué)理論[7,8]。

然而，粘塑性模型同被采用的相應(yīng)橢圓屈服曲線在它的表達(dá)式里是否存在內(nèi)在的不穩(wěn)定性仍是一個(gè)有爭(zhēng)議的問(wèn)題。Gray和Killworth 1995年作了穩(wěn)定性分析，證明模型是不穩(wěn)定的，而且問(wèn)題根本就是不適定的[7]。以前的研究工作表明數(shù)值模式的穩(wěn)定性依賴(lài)于人工擴(kuò)散項(xiàng)的引進(jìn)，而且，模式積分的數(shù)值結(jié)果是（至少部分）人工數(shù)值擴(kuò)散項(xiàng)的函數(shù)。在粘塑性海冰模擬過(guò)程中，動(dòng)量平衡方程和連續(xù)方程構(gòu)成一組基本的方程組。不同于Gray和Killworth的處理方法，我們?cè)谶@里提出：用數(shù)值逼近的方法去計(jì)算相鄰浮冰之間的相互作用力。 在這個(gè)模型里，由于計(jì)算相鄰海冰相互作用力的應(yīng)力張量不連續(xù)，采用一個(gè)連續(xù)函數(shù)序列去逼近不連續(xù)的應(yīng)力張量函數(shù)。然后，通過(guò)這個(gè)序列，重新計(jì)算相鄰海冰間的相互作用力。由此得到的模型是穩(wěn)定的，無(wú)需為了穩(wěn)定而引入人工擴(kuò)散項(xiàng)。

文章的其余部分安排如下：第2→部分，敘述了粘塑性模型和相鄰浮冰之間作用力F的基本計(jì)算方法→。第3部分，主要對(duì)二維水平情況的相互作用力F進(jìn)行近似計(jì)算。第4部分，用擾動(dòng)法[9?11]對(duì)修改后的海冰模型的解的穩(wěn)定性進(jìn)行了分析，證明了在不引進(jìn)人工擴(kuò)散項(xiàng)的情況下，模型是穩(wěn)定的。最后給出本文的結(jié)論。

2 粘塑性海冰模型

以 Gray 和 Morland 1994年提出的相互作用的連續(xù)框架描述大片浮冰的質(zhì)量、動(dòng)量及能量平衡關(guān)系[8]為基礎(chǔ)。用浮冰和前導(dǎo)水組成的水平二維混合公式對(duì)Hilbler 1979年提出的經(jīng)典框架推導(dǎo)類(lèi)似的平衡定理[5]。

2.1 海冰的密集度與厚度的連續(xù)方程

設(shè)xoy是固定在旋轉(zhuǎn)地球的Cartesian 直角坐標(biāo)系，并且定義一個(gè)非慣性框架。冰的密集度和厚度滿(mǎn)足的連續(xù)方程是[7]

2.2 海冰的動(dòng)量平衡方程

考慮海冰在二維平面內(nèi)移動(dòng)。在海冰運(yùn)動(dòng)平面的Cartesian 直角坐標(biāo)系內(nèi)，動(dòng)量平衡方程是[7?8]：

對(duì)于一個(gè)擁有橢圓縱橫比為e的橢圓屈服表面，由文獻(xiàn)[7]，粘塑性模型將通過(guò)法向流規(guī)則將應(yīng)變張量D和應(yīng)力張量Ne聯(lián)系起來(lái)，即：

4 穩(wěn)定性分析

式中上面有橫杠的變量表示平衡態(tài)，上面有波浪線的量為擾動(dòng)變量。將上述表達(dá)式代入(15)－(18)并舍棄高階小量，得到其線性化方程組

5 結(jié) 論

本文通過(guò)插值逼近方法重新計(jì)算浮冰之間的相互作用力，修改了粘塑性流變模型，并用擾動(dòng)理論對(duì)其解的穩(wěn)定性進(jìn)行了分析。證明修改后的模型在某些條件下，其解無(wú)任是浮冰收斂時(shí)的聚集還是發(fā)散時(shí)的分離，其解都是穩(wěn)定的。克服了修改前模型的解存在不穩(wěn)定的現(xiàn)象。

致謝：感謝美國(guó)Vermont 大學(xué)Yu Jun教授對(duì)本文的幫助。

[1]Washington M W, Claire L P.An introduction to Three-Dimensional Climate Modeling [M].New York: Oxford University Press, 1986.

[2]Pritchard R S.An elastic-plastic constitutive law for sea ice [J].Journal of Applied Mechanics, 1975(42): 379-384.

[3]Colony R, Pritchard R S.Integration of elastic-plastic constructive laws [J].AIDJEX- Bull, 1975(30): 55-80.

[4]Hibler W D.A viscous sea ice law as a stochastic average of plasticity [J].Journal of Geophysical Research, 1977, 82(27): 3932-3938.

[5]Hibler W D.A dynamic Thermodynamic sea ice model [J].Journal of Physical Oceanography, 1979(9): 815-845.

[6]Hibler W D.Modeling pack ice as a viscous-plastic continuum: Preliminary results.Sea Ice Processes and Models [M].R.S.Pritchard ED: University of Washington Press, 1980.

[7]Gray J M, Killworth P D.Stability of the viscous plastic sea ice Rheology [J].Journal of Physical Oceanography, 1995(25): 971-978.

[8]Gray J M, Morland L.A two-dimensional model for the dynamics of sea ice [J].Phil.Trans.Roy.Soc, 1994, 34(7): 219-290.

[9]Yu J, Kevorkian J, Haberman R.Weak Nonlinear Long Waves I in Channel Flow with internal Dissipation [J].Studies in applied mathematics, 2000(105): 143-163.

[10]Yu J, Kevorkian J.Nonlinear evolution of small disturbances into roll wave in an inclined open channel [J].Journal of Fluid Mechanics, 1992(243): 575-594.

[11]Horst L.Stability Theory: An Introduction to the Stability of Dynamic Systems and Rigid Bodies [M].Stuttgart: John Wiley & Sons Ltd.and B.G.Teubner, 1987.

Stability analysis of the viscous plastic sea ice rheology model

WANG Guang-hui, GU Xiang-qian

(State Key Laboratory of Severe Weather, Chinese Academy of Meteorological Sciences, Beijing 10081, China)

Revision and stability analysis are carried out for a sea ice model with a viscous-plastic constitutive relation of Hibler (1977) in this paper.In the existing work, numerical model integration in many studies using this viscous-plastic constitute relation are only made possible by introducing an artificial diffusion.Here, we present a new way to compute the interaction force of adjacent sea ice floes by using a mathematical approximate method.The revised model is shown to be stable in all cases including both convergent and divergent flows.The need of artificial diffusion in model integration has now been removed.As an example, we present detailed analysis for a two-dimensional configuration.

sea ice; rheology model; stability; elliptical yield curve; interpolation

P731.15

A

1001-6932(2011)06-0625-08

2011-07-07；

2011-09-26

中國(guó)自然科學(xué)基金 (40875066, 41075080 )。

王光輝(1962－), 男, 研究員, 主要從事數(shù)值天氣預(yù)報(bào)研究。電子郵箱：earth54@126.com。