技术领域
[0001] 本发明涉及水下声源激发的海水表面波预报方法技术领域，具体地涉及实际海洋波导中点声源激发的海水表面波预报方法。
相关背景技术
[0002] 由于水和空气的声阻抗差异极大，海水‑空气界面不是一个理想的压力释放界面。因此在声场作用下，海水表面会形成一层微扰的薄层，并产生沿海水表面传播的微波。这种存在于海水和空气交界面的微波，就是海水表面波。
[0003] 海面的波动还存在由于自然条件干扰所形成的自然表面波，但水下声源激发的海水表面波频率和声源的频率相同，而风力扰动等环境因素激励的自然水表面波，频率一般小于1Hz。因此，可以从频谱上将这种微波和自然水表面波相区分。这为研究人员观测海水表面波，并以此来探测水下目标提供了充分的依据。
[0004] 现有技术对于海水表面波以及探测水下目标的研究和应用具有如下缺陷：
[0005] 1.由于现有技术对于水下声源激励海面形成微幅波的理论研究目前尚处于起步阶段，从而尚无可工业化的对于海水表面波预报解决方案；
[0006] 2.由于现有海水表面波理论的边界是平整海面，但实际上海面是起伏不定的，且具有一定的时空特性；从而导致现有技术的理论基础实际上是不符合实际情形，无法实际应用；
[0007] 3.由于现有技术在研究中是以平面波作为激励波，以计算无限大半空间条件下海水表面波特性，从而并未考虑实际浅海波导条件(水深有限、海底弹性) 对海水表面波的影响，以至于进一步降低了现有技术实际应用可能性。
具体实施方式
[0084] 下面结合具体实施例，进一步阐明本发明，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。
[0085] 一种实际海洋波导中点声源激发的海水表面波预报方法，包含以下步骤：
[0086] S100.建立粗糙海面海浪谱模型。
[0087] 本具体实施例中，粗糙海面海浪谱模型采用高斯谱的功率谱，按式(1)表达：
[0088]
[0089] 其中：S(K)为关于海浪的空间波数的功率谱密度；σ为高斯粗糙面的髙度起伏均方根；l为表面相关长度；K为海浪的空间波数。
[0090] 需要说明的是，对粗糙海面的仿真，目前比较成熟的方法是通过海浪谱进行粗糙海面建模。海浪谱定义为海面高度起伏相关函数的Fourier变换，它给出了海面各谐波分量相对于波数的关系。通过海浪谱可以反映出海面中各谐波分量相对于空间波数和方位的分布，同时也是描述海面的最基本二阶统计量。由海浪谱也可以求出海面一、二阶统计参量，比如海面的波高相关函数和有效波高等参量。常见的海浪有PM谱、高斯谱等，由于高斯谱能直接反映粗糙海面的均方波高，因此这里采用常见的高斯谱。
[0091] S200.建立海面粗糙时无限大半空间下点源激发的海水表面波模型；具体包含以下步骤：
[0092] S210.建立海水表面波对应速度势函数的偏微分方程。
[0093] 如图1所示，对于海水表面波，同样可假设其波面方程为η(x,t)，其中粗糙海面为h(x,t)，亦即粗糙海面的高度函数，用竖直方向的位移来表征粗糙海面相对时间的波形。对入射声波做远场近似，即声源到声波入射点的距离 时 (其中λ为海水中声波的波长)，可将入射的球面波看做平面波。
[0094] 于是，可以将速度势函数和波面函数按照小参数ε进行幂级数展开，并将其偏导数在y＝h(x,t)处对y做泰勒级数展开，按式(2)表达：
[0095]
[0096] 然后，再由线性化自由表面运动学条件 得到海水表面波对应速度势函数的偏微分方程。
[0097] 本具体实施例中，海水表面波对应速度势函数的偏微分方程按式(3)表达：
[0098]
[0099] 其中：t为时间；x为横向距离；y为纵向波高；g为重力加速度；ρ为海水密度；T为海水表面波张力系数；h(x,t)为粗糙海面的高度函数，按式(4) 表达：
[0100]
[0101] 其中：L为粗糙海面的长度；N为长度L内的离散点；Kn为海浪的离散波数，且Kn＝2πn/L；n为计数器；ΔK为相邻谱域之间的谐波样本的海浪空间波数差；N(0,1)为正态分布；e为自然对数；i为虚数单位。
[0102] F(Kn)为粗糙海面的谐波谱函数，按式(5)表达：
[0103]
[0104] 为速度势函数，按式()表达：
[0105]
[0106] 需要说明的是，速度势函数的得到原理如下：
[0107] 首先，粗糙海面的水面微波模型和平整海面的模型相比，边界条件发生了改变，所以可将速度势函数写成式(6)的表现形式：
[0108]
[0109] 其中，A、B、A′之间的关系按式(7)表达：
[0110]
[0111] 然后，将速度势函数代入海面的边界条件中，可得式(8)：
[0112]
[0113] 进一步的，为了使波面方程的解满足初始位置条件，将波面方程写成式(9) 的表现形式：
[0114] η＝Aej(kx‑ωt) (9)
[0115] 其中：A为海水表面波幅值，按式(10)表达：
[0116]
[0117] 其中：A′为求解过程中的待定系数，按式(11)表达：
[0118]
[0119] 需要说明的是，A、A′是通过将式(8)带入式(9)得到的。
[0120] k为海水表面波波数；H(x)为海底随水平距离x变换的海底高度函数；ω为粗糙海面下的频散关系，按式(12)表达：
[0121]
[0122] 这样，就得到了速度势函数的表达式。
[0123] S220.建立海面粗糙时点声源激发的海水表面波模型。
[0124] 本具体实施例中，海面粗糙时点声源激发的海水表面波模型按式(13)表达：
[0125]
[0126] S230.建立海面粗糙时点声源激发的海水表面波位移函数。
[0127] 本具体实施例中，将速度势函数代入到线性化自由表面运动学条件中，从而得到海面粗糙时点声源激发的海水表面波位移函数按式(14)表达：
[0128]
[0129] 海面粗糙时点声源激发的海水表面波位移函数即为S200最终建立并输出的海面粗糙时无限大半空间下点源激发的海水表面波模型。
[0130] S300.建立实际海洋波导条件下点源激发海水表面波模型；具体包含以下步骤：
[0131] S310.解算出海面粗糙时点声源激发的海水表面波反射系数。
[0132] 首先如图2所示，当海面粗糙时，声波以θ的掠射角入射，反射声波pr通过平均反射系数求解。
[0133] 然后利用小斜率近似理论，得到声波在海水中的波数为k1的声波以掠射角为θ的角度入射到海面粗糙时的海面粗糙时点声源激发的海水表面波反射系数；其中：k1为声波在海水中的波数，按式(15)表达：
[0134]
[0135] 本具体实施例中，海面粗糙时点声源激发的海水表面波反射系数按式(16) 表达：
[0136]
[0137] 其中：θ为声波相对于粗糙海面的掠射角；R(sinθ)为海面反射系数；gs(K)为中间变量。
[0138] S320.解算出监测点入射声压。
[0139] 本具体实施例中，监测点入射声压按式(17)表达：
[0140]
[0141] 其中：Q为点声源振幅； 和 分别为经不同反射途径到达监测点的掠射角；γ为海底反射系数；Rf为海面反射系数；pi为监测点入射声压；rn1为直达声波的传播距离；rn2为海底反射声波的传播距离；k1为海水中波数；k为海水表面波波数。
[0142] 需要说明的是，如图3所示，此时监测点入射声压按式(18)表达：
[0143]
[0144] 需要进一步说明的是，由于声波是向各个方向发出的，从声源出发的声波有多条途径到达监测点。这里取n为声波从声源到监测点过程经过粗糙海面反射的次数，则将式(18)解算出来，即得式(17)。
[0145] S330.解算出弹性海底的反射系数。
[0146] 如图4所示，本具体实施例中，弹性海底的反射系数按式(19)表达：
[0147]
[0148] k1 cosθn2＝kp sinβ＝ks sinα
[0149] 其中：γ海底； 为纵波在竖直方向上的阻抗； 为横波在竖直方向上的阻抗； 为入射声波在竖直方向上的阻抗。
[0150] 其中：ρ2为弹性海底的密度；cp为弹性海底的纵波波速；cs为弹性海底的横波波速；为弹性海底纵波波数； 为弹性海底横波波数；θn2为声波入射海底的掠射
角；θn2为横波透射角；β为纵波透射角。
[0151] S340.利用镜像法，将所有的反射声压进行叠加。
[0152] 首先，根据镜像法原理，水中声源出发的声波到达海面监测点主要包括以下声路径：直达波、海底一次反射波、海面‑海底一次反射波，以及其他多次反射波。
[0153] 如图5所示。由监测点入射声压可以得到在该监测点的反射声压按式(20) 表达：
[0154]
[0155] 然后利用线性化欧拉公式可知：
[0156] 本具体实施例中，被叠加后的所有的总振速按式(21)表达：
[0157]
[0158] 其中：V为被叠加后的所有的总振速。
[0159] S350.计算监测点处质点振幅，然后根据监测点处质点振幅，获取表面波模型。
[0160] 本具体实施例中，监测点处质点振幅按式(22)表达：
[0161]
[0162] 需要说明的是，此处的监测点为界面监测点。
[0163] 本具体实施例中，表面波模型按式(23)表达：
[0164]
[0165] 其中：η(x,t)为表面波模型。
[0166] S360.在考虑表面波衰减的条件下，计算海面粗糙时点声源激发的表面波幅值，得到实际海洋波导条件下点源激发海水表面波模型；实际海洋波导条件下点源激发海水表面波模型即为本实际海洋波导中点声源激发的海水表面波预报方法的最终结果。
[0167] 考虑海水的粘性对表面波衰减影响的，增加相应的衰减系数项，可得海面粗糙时点声源激发的表面波的模型按式(24)表达：
[0168]
[0169] 由于点源激励的水面微幅波在水平面的各个方向都有相同的传播规律，因此可将模型推广为三维模型：
[0170] 本具体实施例中，实际海洋波导条件下点源激发海水表面波模型按式(25) 表达：
[0171]
[0172] 其中：α为衰减系数，按式(26)表达，由现有参考文件可查：
[0173]
[0174] 为了进一步说明本发明的原理，发明人在本具体实施例中还进一步对实际海洋波导条件下点源激发海水表面波模型进行了实验验证，具体如下：
[0175] 建立水下声源激发粗糙海面海水表面波的有限元模型，通过自由液面(海面) 的流体力学本构和声学计算域的耦合进行数值模拟，对理论进行验证。流体的自由液面会依据流体的本构发生变形，从而产生海浪。同时划分的网格会跟随液面变形，使声学计算域的水面边界变形，形成声学动边界。并以此计算点声源激发声场，从而得到海水表面波。海面为自由液面，自由液面上方为空气，声源为脉动点源。
[0176] 如图6所示，本具体实施例中，流体计算域的长度为16m，声源深度为1.5m，声源频率1000Hz，完美匹配层厚度为1m，海水边界为开放边界。
[0177] 经模拟，通过所划分的有限元网格如图7所示，初步验证结果如图8所示。
[0178] 为了进一步证明本发明的有效性，发明人在本具体实施例中还上述实际海洋波导条件下点源激发海水表面波模型实验进行了，具体如下：
[0179] 海水表面波沿着粗糙海面传播，但它们的传播规律有着很大区别。
[0180] 设表面波位移曲线频率为100Hz，声源深度为30m，海水深度为100m，风浪为四级风浪。粗糙海面长度为1000m，相关长度为10m，均方波高为四级风浪对应的1.2672。可以看到海水表面波逐渐衰减，而海浪谱则没有呈现有规律的衰减。这是由于粘滞力的存在，使得表面波在沿着海面传播的过程中其幅值将逐渐减小，幅值的最大值位于表面波入射中心处。这些区别表明海水表面波的传播独立于海浪，这为后续海水表面波的探测提供了理论基础。
[0181] 如图9所示，为给出了不同风浪下的表面波幅值随监测点位置的变化曲线。表面波的位移幅值随着距离的增加发生衰减，且随着风浪增加表面波的位移幅值也会逐渐减小。
[0182] 不同的海水深度会对经海底反射后的声波的大小产生影响。如图10所示，为声源频率为100Hz时，不同风浪下表面波幅值随海水深度的变化曲线。可以看到表面波幅值随着海水深度的增大整体上呈现衰减的趋势，但在一些局部海水深度，表面波幅值存在峰值。表面波幅值随着风浪的增大而减小，且在声源频率较高时，2～6级风浪下的表面波曲线发生重合。这是由于声源频率较大时，平均反射系数随均方波高的增大几乎不发生变化。表面波幅值会随海水深度的变化呈现起伏，是声波通过海底反射到达海面引起的。从图中可以看到，随着海水深度的增大，表面波幅值的起伏逐渐减小，并趋近于一定值。
[0183] 在上述的详细描述中，各种特征一起组合在单个的实施方案中，以简化本公开。不应该将这种公开方法解释为反映了这样的意图，即，所要求保护的主题的实施方案需要比清楚地在每个权利要求中所陈述的特征更多的特征。相反，如所附的权利要求书所反映的那样，本发明处于比所公开的单个实施方案的全部特征少的状态。因此，所附的权利要求书特此清楚地被并入详细描述中，其中每项权利要求独自作为本发明单独的优选实施方案。
[0184] 为使本领域内的任何技术人员能够实现或者使用本发明，上面对所公开实施例进行了描述。对于本领域技术人员来说；这些实施例的各种修改方式都是显而易见的，并且本文定义的一般原理也可以在不脱离本公开的精神和保护范围的基础上适用于其它实施例。因此，本公开并不限于本文给出的实施例，而是与本申请公开的原理和新颖性特征的最广范围相一致。
[0185] 上文的描述包括一个或多个实施例的举例。当然，为了描述上述实施例而描述部件或方法的所有可能的结合是不可能的，但是本领域普通技术人员应该认识到，各个实施例可以做进一步的组合和排列。因此，本文中描述的实施例旨在涵盖落入所附权利要求书的保护范围内的所有这样的改变、修改和变型。此外，就说明书或权利要求书中使用的术语“包含”，该词的涵盖方式类似于术语“包括”，就如同“包括，”在权利要求中用作衔接词所解释的那样。此外，使用在权利要求书的说明书中的任何一个术语“或”是要表示“非排它性的或者”。
[0186] 以上所述的具体实施方式，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施方式而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
[0187] 最后，应当指出，以上实施例仅是本发明较有代表性的例子。显然，本发明不限于上述实施例，还可以有许多变形。凡是依据本发明的技术实质对以上实施例所做的任何简单修改、等同变化及修饰，均应认为属于本发明的保护范围。
[0188] 在此，需要说明的是，上述技术方案的描述是示例性的，本说明书可以以不同形式来体现，并且不应被解释为限于本文阐述的技术方案。相反，提供这些说明将使得本发明公开将是彻底和完整的，并且将向本领域技术人员充分传达本说明书所公开的范围。此外，本发明的技术方案仅由权利要求的范围限定。
[0189] 用于描述本说明书和权利要求的各方面公开的形状、尺寸、比率、角度和数字仅仅是示例，因此，本说明书和权利要求的不限于所示出的细节。在以下描述中，当相关的已知功能或配置的详细描述被确定为不必要地模糊本说明书和权利要求的重点时，将省略详细描述。
[0190] 在使用本说明书中描述的“包括”、“具有”和“包含”的情况下，除非使用否则还可以具有另一部分或其他部分，所用的术语通常可以是单数但也可以表示复数形式。
[0191] 应该指出，尽管在本说明书可能出现并使用术语“第一”、“第二”、“顶部”、“底部”、“一侧”、“另一侧”、“一端”、“另一端”等来描述各种不同的组件，但是这些成分和部分不应受这些术语的限制。这些术语仅用于区分一个成分和部分和另一个成分和部分。例如，在不脱离本说明书的范围的情况下，第一部件可以被称为第二部件，并且类似地，第二部件可以被称为第一部件，顶部和底部的部件在一定情况下，也可以彼此对调或转换；一端和另一端的部件可以彼此性能相同或者不同。
[0192] 此外，在构成部件时，尽管没有其明确的描述，但可以理解必然包括一定的误差区域。
[0193] 在描述位置关系时，例如，当位置顺序被描述为“在...上”、“在...上方”、“在... 下方”和“下一个”时，除非使用“恰好”或“直接”这样的词汇或术语，此外则可以包括它们之间不接触或者接触的情形。如果提到第一元件位于第二元件“上”，则并不意味着在图中第一元件必须位于第二元件的上方。所述部件的上部和下部会根据观察的角度和定向的改变而改变。因此，在附图中或在实际构造中，如果涉及了第一元件位于第二元件“上”的情况可以包括第一元件位于第二元件“下方”的情况以及第一元件位于第二元件“上方”的情况。在描述时间关系时，除非使用“恰好”或“直接”，否则在描述“之后”、“后续”、“随后”和“之前”时，可以包括步骤之间并不连续的情况。本发明的各种实施方案的特征可以部分地或全部地彼此组合或者拼接，并且可以如本领域技术人员可以充分理解的以各种不同地构造来执行。本发明的实施方案可以彼此独立地执行，或者可以以相互依赖的关系一起执行。
