海啸条件下潜艇运动的数学建模与仿真

赵冬冬， 张　栋， 吴猛猛

(海军潜艇学院， 山东 青岛　266042)

海啸是一种具有强大破坏力的海浪，水下地震、火山爆发或水下海底塌陷和滑坡等大地活动都可能引发海啸。海啸按照受灾地点与海啸源点距离的远近，分为越洋海啸和局地海啸。

越洋海啸对中国大陆沿海地区产生的影响很小，但是局地海啸对我国的影响较大，虽然渤海和黄海爆发局地海啸的可能性很小，但东海和南海却具备海啸产生的条件，历史上这两个海域也产生过海啸[1],那么这就要求开展近场海啸的危险性分析和海啸预警等工作，本文依据近场海啸数值模拟理论和计算模式建立海啸数学仿真模型,研究海啸对潜艇运动的影响。

1　海啸数学仿真模型的建立

1.1　数学模拟的计算方法

由于海洋深度的不规格变化、海岛的不均匀分布、海啸本身的复杂性等不确定性问题，给海啸的模拟带来众多困难。利用传统的理论方法对海啸方程进行求解，计算相当繁琐。

随着计算机技术的发展，关于海啸数值模拟的方法也越来越普及。目前求解海啸数值模拟的方法大多采用有限差分法[2]，也有一些数值计算模型采用有限单元法和有限体积法。

海啸波浪不同于海面波浪，其波长远大于其发生地和传播地的海水深度，因此海啸波属于长波范畴，质点运动的垂直加速度相对于重力加速度来说非常小，可以忽略不计[3-4]。基于这些假设，可用非线性浅水方程表示海啸波的运动：

(1)

(2)

(3)

其中，g是重力加速度，x和y为平面坐标，h是海区深度，t是时间，u和v分别为x和y的分速度，η是自由表面高度，τx/ρD和τy/ρD分别是x和y的底部摩擦力。

底部摩擦力可表示成：

(4)

(5)

其中，D是实际深度(D=h+η)，n为曼宁粗糙系数。M、N分别为X和Y方向的流量关于u和v的表达式如下：

M=u(h+η)=uD

(6)

N=ν(h+η)=νD

(7)

将式(4)～(7)分别代入上式得：

(8)

(9)

(10)

因此可知地震海啸海区任意一点的水深之后，波浪传播到此点的速度为

(11)

海啸在X和Y方向的传播速度为

(12)

(13)

其中，θ为海啸传播方向与x轴的夹角。

1.2　地震初始位移场的确定

首先进行断层参数的确定，主要有两种情况：

1)对于已发生的地震海啸，由于数据较多，且经过大量研究工作者的研究，可以直接引用。

2)如果地震海啸未发生，可采用日本气象厅给出的经验公式[5],见图1。

logL=0.5P-1.9

(14)

logW=0.5P-2.2

(15)

logQ=0.5P-3.2

(16)

图1　中心差分断层参数示意图

其中，P为地震震级，L为断层长度，W为断层宽度，Q为滑动距离。

本文在研究日本地震算例时采用第二种方法来确定地震初始位移场。

1.3　算例分析

2011年3月11日，日本东部海域发生高强度地震，测得地震震级强度里氏9.0级，震源深度10km，距离海岸190km。

首先计算其初始垂直位移场。断层参数采用经验公式计算，即式(14)～式(16)，结果为

L=708km,W=354km,Q=20.0km。

需要说明的是本文实际采用的滑移量Q=14.3m是作了修正的，因为采用式(16)计算得到的结果偏大，理由是用此式计算智利海啸和苏门答腊海啸的滑移量，结果偏大[4]。由于缺乏数据，故而只选取智利海啸和苏门答腊海啸来修正此次日本海啸，采用取平均值的方法，修正过程见表1。

表1　海啸滑移量计算修正

本文计算采用的距离步长为0.1km，海区深度步长为1m，代入方程编程得出海啸传播过程仿真图，如图2所示。

图2　海啸产生的波浪与离岸距离关系示意图

通过计算可得，采用本文模型模拟3·11地震海啸的最大浪高为22.7m，在海区水深205m时达到最大值，因此海啸从海区水深205m传播至岸边时，海啸产生的浪高度都会达到最大高度22.7m。而在对日本海啸的实际观测中，测得海啸的最大高度24m，与本文模型计算的数值相差1.3m，分析原因有：一是本文采用的海底地形是缓变地形，即海区水深变化缓慢，不存在地形突变现象；二是海区水文条件复杂，而数值模拟相对简单，因此与实际情况存在一定的差异[6]。

2　潜艇在海啸波浪中的纵、横倾

2.1　水面状态潜艇在海啸波浪中的横倾及纵倾

由于潜艇的航速相对于海啸的传播速度非常小，因此潜艇的速度可以忽略为零，在这里只考虑最严重的情况，即潜艇无航速受海啸波浪的影响，且海啸波浪从潜艇的正横方向传过来，此时产生的横倾最大。由于水表面不是静水面，而是海啸波浪的波面，其中α表示波倾角。因潜艇的宽度远小于海啸所产生的波长，因此可认为船宽范围内波面是直线。如图3所示。

图3　潜艇相对波面横倾

如图3所示，潜艇相对于静水面横倾θ角，相对于波面则横倾θ-α角度，此时船所受的扶正力矩可写成：

M(θ)=-Dh(θ-α)=-Dhθ+Dhα

(17)

再加上惯性力矩的阻尼力矩，按照牛顿第二定律可写成在规则波上的横倾微分方程为

(18)

或者

(19)

由于上面分析时相对理想化、简单化，因此考虑到实际情况引入一个修正因子∂eθ。一般可取

∂eθ=e-kd/2

(20)

式中,k为波数，d为潜艇的吃水，相当于取二分之一吃水处的次波面来计算波浪的扰动力矩。

式(19)是数学上典型的二阶常系数线性非齐次微分方程，它的稳定解为

θ=θαsin(ωt-εθ-α)

(21)

式中：θα为横倾幅值，εθ-α为横倾运动与波倾角的相位差。它们的表达式分别为：

(22)

(23)

上式还可改写为

(24)

海啸波的波数与海啸波的波长有关，因此需要求解海啸的波长。海啸的波长也是与深度密切相关的，海啸波是重力波，在深海位置的波长特别长，运动到近岸处，波长变短，频率变大。在本文中根据优于线性插值更佳的经验公式：

λ=const1+h2×const2

(25)

其中，const1为基准波长，在深海中取200km-400km，近海岸中取10km-50km，h为海区水深，单位为m，const2为调整系数，本文取1.4×10-5。

如图4所示，由于潜艇水上状态纵稳性高较高，因此潜艇遭受海啸时纵倾变化相对于横倾来说数值较小。本文在计算潜艇水上纵稳性高时，由于浮心和重心间距离较短，远远小于纵稳心与重心间的距离，可以忽略不计，因此本文中将纵稳心半径作为纵稳心高。

图4　潜艇浮心、重心及纵稳心高示意图

图4中为M点：纵稳心；C点为浮心；g点为重心；R点为纵稳心半径，即为纵稳心与浮心间的垂向距离。

2.2　水下状态潜艇在海啸波浪中的横倾及纵倾

与研究水面状态潜艇在海啸波浪中的横倾一样，当处于平衡状态的水下潜艇遭受海啸袭击时，使其产生横倾和纵倾。由于潜艇在水下状态时，无论倾角多大，其排水容积和形状不变，所以水下排水量和对应的水下浮心位置不变。另外认为潜艇在水下倾斜过程中艇上载荷没有增减和移动，故潜艇水下重量和重心也不变，此时，潜艇的有效水线面积为零，所以水线面积惯性矩为零，即稳定中心半径为零，也就是说潜艇在水下状态时，浮心、横稳定中心和纵稳定中心三点重合。可见，当不考虑液体载荷的自由液面影响时，潜艇水下的纵稳度与横稳度相等。

如图5所示，潜艇水下纵、横稳心高相等，则有

H↓=h↓=Zc↓-Zg↓

(26)

水下扶正力矩应有如下公式：

mθ↓=P↓(Zc↓-Zg↓)sinθ

(27)

MΨ↓=P↓(Zc↓-Zg↓)sinΨ

(28)

当纵倾与横倾相等时，其扶正力矩也相等。

图5　潜艇水下相对波面横倾

水下潜艇横倾1°，其横倾力矩如下式：

mθ↓=P↓H↓sin1°=0.0175P↓H↓

(29)

或者

(30)

水下潜艇纵倾1°，其纵倾力矩如下式：

MΨ↓=P↓H↓sin1°=0.0175P↓H↓

(31)

或者：

(32)

这时，只要能求出外力矩即可求得倾角。外力矩的求解由水上横倾计算公式求得。本章所建立的模型在海洋环境模拟装置中进行了实验验证，得出的实验结果与计算结果接近，说明所建立的模型较为准确[7]。

本文选取国内外某4型现役潜艇为例(参数如表2所示)，将潜艇性能参数代入模型求得潜艇在海啸波浪中所产生的横、纵倾(以日本3·11海啸为模型)，计算仿真结果如图6～图9所示。

表2　国内外某4型现役潜艇参数表

下面针对不同艇型计算其在海啸波中的倾斜角度。

图6　某Ⅰ型潜艇不同状态倾角与离岸距离关系图

图7　某Ⅱ型潜艇不同状态倾角与离岸距离关系图

图8　某Ⅲ型潜艇不同状态倾角与离岸距离关系图

图9　某Ⅳ型潜艇不同状态倾角与离岸距离关系图

从图中可以看出，相同条件下Ⅳ型潜艇在遭遇海啸时其稳性要比其他类型核潜艇好，主要是由于其不论是从艇体尺度还是排水量、稳心高等等，都是较大的。同时，Ⅳ型潜艇水下横倾要比水上横倾小，原因除了潜艇水下排水量大之外，还由于潜艇水下稳心高要比水上稳心高大，潜艇的扶正力矩较大。

3　结果分析

通过前面的仿真结果，可以得出以下结论：

1)同一艇型的潜艇，水上横倾和水下横倾的数值不同，但是水上横倾数值和水下横倾数值哪个大，与艇型有直接关系，除Ⅳ型潜艇水下横倾值比水上横倾值小之外其余艇型潜艇，水下横倾值均比水上横倾值大，因此，潜艇一旦遭遇海啸时应尽量在海面进行航行规避。

2)潜艇的横倾值、纵倾值与潜艇的基本性能参数有关，其角度大小取决于排水量、稳心高等，不同艇型的潜艇，随着排水量的增大以及稳心高的增大，其最大横倾值和纵倾值不断减小。

3)潜艇在水面状态主要是考虑横倾，但是在水下状态的失稳主要来自纵稳性，尤其是潜艇纵向较长，载荷移动力臂大，因此移动不大的载荷都有可能造成可观的纵倾力矩，因此潜艇在水下应将活动物品固定牢靠，防止额外力矩增大潜艇水下纵倾。

4)通过以上示意图可以看出，不同艇型在同一海啸中所产生的横倾也不同，与潜艇的性能参数有关。同时也可以看出当潜艇遭遇海啸时，距离海岸越远所产生的横倾越小，潜艇不易察觉，但是近岸所产生的横倾较大，且随着海区水深的急剧减小而骤然增大，因此潜艇在海上一旦遇到海啸，应采取远离近岸航行的措施应对。

5)通过以上示意图的潜艇水上纵倾曲线可以看出，遭受海啸时潜艇水上纵倾变化不大，最大值为9.2°，相对于横倾30°、40°以上，纵倾角度非常小，同时随着潜艇排水量的增加、纵稳心高的增大等，潜艇的水上纵倾角会变小。

6)潜艇在近岸水下航行时，当潜艇横倾瞬间增大，潜艇瞬间偏离主航向，难以恢复航向时，如若查明潜艇本身机械一切正常，即可认为潜艇水下遭遇巨大波浪，有可能遭遇海啸，此时应当尽量增速，改变潜艇的航向，使之顶浪航行，并且注水，以稳住潜艇，只要条件允许尽快向远海处航行。

4　结束语

研究地震海啸对潜艇的影响，需要我们精确模拟地震海啸从震源至岸边的整个传播过程，但是地震海啸受海底地形、海洋环境等的影响较大，本文建立的仿真数学模型经日本“311”地震海啸检验与实际监测数据基本吻合，但是仍然需要加强对地震海啸的模拟。

此外，研究地震海啸对潜艇的影响是为后面研究制定潜艇应对措施奠定基础的，因此研究潜艇应对地震海啸的措施是后续要进行的主要工作。

参考文献:

[1]王培涛，高义，于福江，等.基于数值模拟的渤海海域地震海啸危险性定量化研究[J].海洋学报,2014,36(1):56-64.

[2]姚远，蔡树群，王盛安，等.海啸波数值模拟的研究现状[J].海洋科学进展,2007,25(4):488-494.

[3]Meteorology and Hydrology Center,PLA Headquarters of General Staff . Focus on t sunamis[M].Beijing:Military Science PublishingHouse, 2005:5-28.

[4]任叶飞.基于数值模拟的我国地震海啸危险性分析研究[D].哈尔滨:中国地震局工程力学研究所,2007.

[5]付庆军.越洋和局地海啸数值模拟的理论研究及其应用[D].天津:天津大学,2006.

[6]潘文亮，王盛安，孙璐，等.2010年智利和2011年日本海啸在华南沿岸的实测海啸波形和特征[J].热带海洋学报,2014,33(6):17-23.

[7]刘常波，张栋，等.潜艇海上应对地震海啸有关问题研究[R].青岛:海军潜艇学院,2012.