EI、Scopus 收录
中文核心期刊

移动载荷作用下浮冰的非线性动力学响应

孟洋涵, 王展

孟洋涵, 王展. 移动载荷作用下浮冰的非线性动力学响应. 力学学报, 2022, 54(4): 862-871. DOI: 10.6052/0459-1879-22-040
引用本文: 孟洋涵, 王展. 移动载荷作用下浮冰的非线性动力学响应. 力学学报, 2022, 54(4): 862-871. DOI: 10.6052/0459-1879-22-040
Meng Yanghan, Wang Zhan. Nonlinear dynamic response of a floating ice sheet to a moving load. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(4): 862-871. DOI: 10.6052/0459-1879-22-040
Citation: Meng Yanghan, Wang Zhan. Nonlinear dynamic response of a floating ice sheet to a moving load. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(4): 862-871. DOI: 10.6052/0459-1879-22-040
孟洋涵, 王展. 移动载荷作用下浮冰的非线性动力学响应. 力学学报, 2022, 54(4): 862-871. CSTR: 32045.14.0459-1879-22-040
引用本文: 孟洋涵, 王展. 移动载荷作用下浮冰的非线性动力学响应. 力学学报, 2022, 54(4): 862-871. CSTR: 32045.14.0459-1879-22-040
Meng Yanghan, Wang Zhan. Nonlinear dynamic response of a floating ice sheet to a moving load. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(4): 862-871. CSTR: 32045.14.0459-1879-22-040
Citation: Meng Yanghan, Wang Zhan. Nonlinear dynamic response of a floating ice sheet to a moving load. Chinese Journal of Theoretical and Applied Mechanics, 2022, 54(4): 862-871. CSTR: 32045.14.0459-1879-22-040

移动载荷作用下浮冰的非线性动力学响应

基金项目: 中国科学院B类先导项目资助(XDB22040203)
详细信息
    作者简介:

    王展, 研究员, 主要研究方向: 水动力学. Email: zwang@imech.ac.cn

  • 中图分类号: O352

NONLINEAR DYNAMIC RESPONSE OF A FLOATING ICE SHEET TO A MOVING LOAD

  • 摘要: 本文考虑非线性、惯性和阻尼的影响, 研究了任意深度二维理想流体顶部浮冰的振动. 对相关的拟微分算子进行展开并将非线性项保留至三阶后, 完全非线性问题被简化为仅与自由面上的变量相关的三阶截断模型. 为了验证简化模型的准确性, 重点关注了自由孤立波解. 在不考虑阻尼的情况下, 采用多重尺度方法推导了三阶非线性薛定谔方程(NLS), 利用该方程预测了任意水深下原始欧拉方程中自由波包型孤立波解的存在性及三阶截断模型的准确性. 相比于Dinvay等所提出的二阶模型, 三阶截断模型的优势在于其对应的三阶NLS具有准确的非线性项系数, 能够在最小相速度附近更好地模拟冰层的动力学响应. 进一步地对自由孤立波解进行数值计算, 数值结果表明三阶截断模型在分岔曲线和孤立波波形上均与完全欧拉方程吻合良好, 准确性高于二阶截断模型. 基于三阶截断模型, 探究了匀速局域化载荷作用下的浮冰非线性动力学响应并将时间依赖解与实验测量数据进行比较, 数值计算结果与实验记录吻合良好.
    Abstract: Vibrations of a floating ice cover on top of a two-dimensional ideal fluid of arbitrary depth are studied when the effects of nonlinearity, inertia, and damping are all considered. We reduce the fully nonlinear problem to a cubic-truncation system involving variables on the free surface by expanding the relevant pseudo-differential operators and retaining nonlinear terms up to the third order. To validate the accuracy of the reduced model, we focus on the free wavepacket solitary wave solutions. In the absence of damping, the normal form analysis is performed to derive the cubic nonlinear Schrödinger equation, which predicts the existence of free wavepacket solitary waves in the primitive equations and the accuracy of the cubic-truncation model. The main advantage of the cubic-truncation approximation over the quadratic-truncation model is that the resultant NLS equation has correct coefficient of the nonlinear term, which allows a better approximation of dynamic responses of the ice cover near the phase speed minimum. Solitary waves are then numerically computed, and it is shown that the cubic-truncation approximation agrees well with the full Euler equations for bifurcation curves and wave profiles, indicating that the reduced model is more accurate than the quadratic truncation model. The nonlinear dynamic response of a floating ice sheet to a fully localized constant-moving load is investigated based on the cubic-truncation model. The time-dependent solutions are compared with the data from the field measurements, and good agreement is achieved between the numerical results and experimental records.
  • 作为经典水波理论的扩展, 水弹性波考虑了自由面的弹性效应, 可用于描述弹性薄板与流体之间的相互作用. 对水弹性波的研究起源于对高纬度地区大型冰盖的利用, 这些冰盖常被作为季节性的交通通道供车辆行驶及小型飞机起飞降落. 当车辆或小型飞机在冰盖上行驶时, 移动的压力源会激发一系列在冰层和流体界面处传播的波, 称之为水弹性波, 其涉及的回复力包括重力和冰层变形对流体施加的弹性. 与重力毛细波不同, 水弹性波的波长可达到数十米甚至上百米, 因此在实验中更易观测. Takizawa[1]和Squire等[2]分别在Lake Saroma和McMurdo Sound对移动载荷作用下的冰层响应进行了测量, 得到了不同运动速度下的冰层挠度及应变.

    为了更好地理解运动载荷下的冰层响应问题, 许多研究者从理论上对水弹性波进行了细致的研究, 但早期的研究多基于线性理论. Kheisin[3]在研究点载荷及线载荷作用下的冰层位移时发现线载荷存在两个临界速度, 在临界速度下冰层挠度趋于无穷. Nevel[4]扩展了Kheisin的分析方法以处理分布在圆形区域上的均匀载荷, 发现点载荷也存在对应的临界速度. 进一步的分析表明这个临界速度就是自由水弹性波的最小相速度[5]. 接着文献[6]运用渐近傅里叶展开研究了匀速运动载荷激发的远场稳定波形, 发现临界速度恰巧为水弹性波的群速度, 由此给出了临界速度处冰层响应趋于无穷的物理解释. 与此同时, 他们还在不同的系统参数下得到了各类复杂的波形包括焦散和零波响应区. Babaei等[7]以及van der Sanden和Short[8]利用卫星观测了车辆在冰面上行驶时所激发的波形, 观测结果验证了文献[6]对远场波形的理论预测. Schulkes等[9]在文献[6]工作的基础上, 进一步考虑了冰层压应力, 均匀流以及流体分层对冰层响应的影响. 对于临界速度处所出现的奇性, Kheisin[10]意识到可以从瞬态冰层动力学响应入手, 分析结果表明, 当线载荷以最小相速度$ c_{\mathrm{min}} $运动时水弹性波振幅与$ t^{1/2} $成正比, $ t $为时间参数. Schulkes和Sneyd[11]则在研究匀速线载荷作用下的时间依赖响应时得到了第二个临界速度$ \sqrt{gH} $ (g为重力加速度, $ H $为水深), 当载荷以此速度运动时水弹性波振幅与$ t^{1/3} $成正比. 然而Nugroho等[12]的研究结果却发现在三维情况下, 当载荷运动速度为$ \sqrt{gH} $时仍然可以得到稳态有界波形, 并不会出现水弹性波振幅无限增长的情况, 因此$ \sqrt{gH} $并不能称之为临界速度. Miles和Sneyd[13]探究了加速载荷作用下的冰层响应, 研究结果显示将载荷从静止逐渐加速到两个临界速度可以避免冰层响应在临界速度处的奇异.

    在对冰层动力学响应的实验研究中, Wilson[14]最早发现了观测点处冰层最大垂直位移对应的时刻与载荷经过观测点的时刻之间的差异. 随后Beltaos[15], Takizawa[1]以及Squire等[2]在他们的实验中证实了这种滞后效应的存在. 同时Takizawa[1]在控制方程中引入耗散项对实验中出现的滞后效应进行了解释, 但其所考虑的仅是稳态冰层响应而无法得到冰层响应的时历曲线. Hosking等[16], Wang等[17]利用记忆函数在线性理论中引入黏弹性得到了临界速度$ V = c_{\text{min}} $处的有界冰层响应以及滞后行为.

    黏性项的引入能够合理地解释实验中所出现的滞后效应, 但利用线性理论对临界速度附近的大振幅冰层动力学响应进行探究仍然存在较大的局限, 因此非线性对于描述移动载荷作用下的冰层效应同样重要. Părău和Dias[18]在考虑线性理论中的两个共振(奇点)情况时首次引入非线性效应并得到了临界速度附近的有界冰层位移. Dinvay等[19]考虑非线性、黏性以及惯性效应, 利用Dirichilet-Neumman算子建立了能够描述各类运动载荷下冰层响应的完全色散弱非线性模型并通过数值计算结果与实验结果之间的对比验证了模型的有效性. 但是需要指出的是由于非线性项仅保留至二阶, 该模型对应的非线性薛定谔方程是不准确的. 在利用二阶模型计算临界速度附近的冰层动力学响应时所得到的数值结果可能会与实验结果存在较大差异. 除此之外, 他们所关注的实验均为水深较小的情况, 而未对Squire等[2]在深水中的实验进行讨论. 因此仍需要发展更为准确的非线性黏弹性理论来描述不同水深情况下临界速度$ c_{\text{min}} $附近的冰层大幅度挠曲.

    通过对相关的拟微分算子进行展开并将非线性项保留到三阶, 本文将完全非线性二维问题转化为仅与自由面上变量相关的一维系统, 即三阶截断模型. 为验证三阶截断模型的准确性及其相对于二阶模型[19]的优势, 本文从理论和数值上对波包型孤立波解进行了探究. 理论上不考虑黏性和外加载荷的作用, 采用多重尺度方法推导三阶非线性Schrödinger方程, 基于此方程预测不同水深下孤立波解的存在性以及三阶截断模型的准确性. 数值上仍然暂不考虑黏性和外加载荷, 计算完全欧拉方程、三阶截断模型以及二阶截断模型在两个典型水深下的孤立波解及分岔曲线. 数值结果表明, 三阶截断模型能够较好地与欧拉方程吻合, 其精度远高于二阶截断模型. 最后将黏性和外加载荷考虑在内, 基于发展的三阶截断模型对两个工况中移动载荷作用下的浮冰动力学响应进行数值模拟并将数值结果与实验记录进行对比.

    考虑密度为$ \rho $, 水深为$ h $的不可压缩、无黏流体, 其顶部自由面覆盖密度为$ \rho_i $, 厚度为$ d $的无限大冰层, 可通过弹性形变对流体施加回复力. 假设冰层在受到移动载荷作用前没有受到压缩或拉伸. 流体底部为刚性边界, 不可穿透. 建立二维笛卡尔坐标系, 使$ x $轴与未发生形变时的冰层重合, $ y $轴垂直向上, 与重力加速度$ g $的方向相反, 冰层位移为$ \eta(x,t) $. 流体运动无旋, 因此引入速度势函数$ \phi $, 其满足Laplace方程

    $$ {\phi _{xx}} + {\phi _{yy}} = 0,\;\;\; - h < y < \eta (x,t) $$

    自由面上的运动学及动力学条件为

    $$ \begin{array}{l} \eta_{t}-\phi_{y}+\phi_{x}\eta_{x} = 0\\ \phi_{t}+\dfrac{1}{2}\left(\phi_{x}^2+\phi_{y}^2\right)+g\eta+\dfrac{P}{\rho} = p(x,t) \end{array} $$

    其中, $ P $为冰层弹性变形引起的回复力, $ p(x,t) $表示施加的外力载荷. 采用Toland[20]提出的水弹性模型, 同时考虑冰层的惯性及黏性, 于是自由面上的法向应力平衡方程为

    $$ P = P_{a}+D\left(\frac{1}{2}\kappa^3+\kappa_{ss}\right)+\rho_{i}\eta_{tt}+b\eta_{t} $$

    其中$ P_a $为大气压, 恒为常数, $ D = \dfrac{Ed^3}{12(1-\nu^2)} $为冰层的弹性刚度, $ \kappa $为冰层曲率, $ s $为弧长参数, 黏性参数$ b>0 $. 本文将冰层视为各向同性材料, 一方面是由于各向异性材料所涉及的弹性参数过于复杂, 在建立数学模型描述法向应力突变时会存在较大的困难; 另一方面, Takizawa[1]在实验中测量了冰层在横向和纵向上的弹性模量, 两者之间差异不大. Squire等[2]也指出真实弹性模量与冰层响应并不直接相关, 更为重要的是实验中的有效弹性模量, 因此本文采用了简化的各向同性假设[6, 19, 21, 22]并使用了实验测量得到的有效弹性模量. 此外底部边界处满足不可穿透条件${\phi _y} = 0,\;y = - h{\rm{ }} $.

    选取$\left(\dfrac{D}{g\rho}\right)^{\frac{1}{4}},\quad \left(\dfrac{D}{g\rho^5}\right)^{\frac{1}{8}}, \quad\left(\dfrac{Dg^3}{\rho}\right)^{\frac{1}{8}},\quad \left(\dfrac{D^3g}{\rho^3}\right)^{\frac{1}{8}} $为单位长度、单位时间、单位速度及单位势函数对原始欧拉方程进行无量纲化. 于是可得无量纲化后的欧拉方程

    $$ \phi_{xx}+\phi_{yy} = 0,\qquad -h<y<\eta(x,t) $$ (1)
    $$ \eta_t = \phi_y-\eta_{x}\phi_{x}, \qquad y = \eta(x,t) $$ (2)
    $$ \begin{aligned} {\phi _t} + \frac{1}{2}\left( {\phi _x^2 + \phi _y^2} \right) + & \eta + \left( {\frac{1}{2}{\kappa ^3} + {\kappa _{ss}}} \right) + \alpha {\eta _{tt}} + \beta {\eta _t} = p,\\ & y = \eta (x,t) \end{aligned}$$ (3)
    $$ \phi_y = 0,\qquad y = -h $$ (4)

    其中,惯性系数$ \alpha = \dfrac{{\rho_i}d}{\rho}\left(\dfrac{{\rho}g}{D}\right)^{\frac{1}{4}} $, 黏性系数$ \beta = \dfrac{b}{\rho}\left(\dfrac{Dg^3}{\rho}\right)^{1/8} $.

    为了计算载荷作用下的冰层响应并将数值结果与已有的实验结果进行对比, 引入Dirichilet-Neumman (DtN)算子对原始欧拉方程进行近似处理. 利用该算子可避免求解自由边界上的Laplace方程, 将自由面处的边界条件转化为以正则变量表达的形式, 消去边界条件对$ y $坐标的依赖, 将原始的二维问题转换为一维系统. 同时可对惯性项所引入的二阶时间导数项进行处理, 便于进行后续的时间依赖计算. 令自由面处的速度势函数$ \xi(x,t) = \phi(x,\eta(x,t),t) $, 定义DtN算子

    $$ {\rm{G}}(\eta)\xi = \phi_{y}-\eta_{x}\phi_{x} $$

    Craig和Sulem[23]已经证明当$ \eta $的范数小于一个确定的值, $ {\rm{G}}(\eta) $是解析的并且可以展开为级数形式

    $$ {\rm{G}}(\eta)\xi = \sum\limits_{j = 0}^{\infty}{\rm{G}}_{j}(\eta)\xi $$

    其前三项可写为

    $$ \begin{array}{l} {\rm{G}}_0 = {\rm{D}}\tanh(h{\rm{D}})\\ {\rm{G}}_1 = -{\rm{G}}_0{\eta}{\rm{G}}_0-\partial_{x}{\eta}\partial_{x}\\ {\rm{G}}_2 = \dfrac{1}{2}\partial_{xx}{\eta^2}{\rm{G}}_0+\dfrac{1}{2}{\rm{G}}_0{\eta^2}\partial_{xx}+{\rm{G}}_0{\eta}{\rm{G}}_0{\eta}{\rm{G}}_0 \end{array} $$

    其中$ {\rm{D}} = (-\partial_{xx})^{1/2} $. 深水情况下$ {\rm{G}}_0 = {\rm{D}} $, 其他两项形式不变. 由此, 惯性项$ \eta_{tt} $可用DtN算子表示

    $$ \begin{split} \eta_{tt} =& ({\rm{G}}(\eta)\xi)_t = {\rm{G}}(\eta_t)\xi+{\rm{G}}(\eta)\xi_t =\\ & ({\rm{G}}_0+{\rm{G}}_1+{\rm{G}}_2)\xi_t+{\rm{G}}_1({\rm{G}}_0\xi)\xi+{\rm{G}}_1({\rm{G}}_1\xi)\xi+\\ & {\rm{G}}_0({\eta}{\rm{G}}_0\xi)\xi_{xx}+\partial_{xx}\eta({\rm{G}}_0\xi)^2+\\ &{\rm{G}}_0({\rm{G}}_0\xi){\rm{G}}_0{\eta}{\rm{G}}_0\xi+{\rm{G}}_0{\eta}{\rm{G}}_0({\rm{G}}_0\xi)^2 \end{split} $$

    这里保留至关于$ \xi $的三阶项. 将用DtN算子表示的各项代入自由面的运动学和动力学条件, 可得

    $$ \eta_{t}-{\rm{G}}_0{\xi} = ({\rm{G}}_1+{\rm{G}}_2)\xi $$ (5)
    $$ \begin{split} {\rm{F}}\xi_t = & -\frac{1}{2}\left[{\xi_x}^2-({\rm{G}}_0\xi)^2-2({\rm{G}}_0\xi)({\rm{G}}_1\xi)-2\eta_{x}\xi_{x}{\rm{G}}_0\xi\right]-\eta-\\ &\alpha[{\rm{G}}_1({\rm{G}}_0\xi)\xi+{\rm{G}}_1({\rm{G}}_1\xi)\xi+{\rm{G}}_0({\eta}{\rm{G}}_0\xi)\xi_{xx}+\\ &\partial_{xx}\eta({\rm{G}}_0\xi)^2+{\rm{G}}_0({\rm{G}}_0\xi){\rm{G}}_0{\eta}{\rm{G}}_0\xi+{\rm{G}}_0{\eta}{\rm{G}}_0({\rm{G}}_0\xi)^2]-\\ &\left(\frac{1}{2}\kappa^3+\kappa_{ss}\right)-\beta\left({{\rm{G}}_0+{\rm{G}}_1+{\rm{G}}_2}\right)\xi+p \\[-15pt] \end{split}$$ (6)

    算子F定义为

    $$ {\rm{F}} = 1+\alpha{{\rm{G}}_0}+\alpha({\rm{G}}_1+{\rm{G}}_2) $$

    $ 1+\alpha{{\rm{G}}_0} $记为K. 对算子F取逆算子可得

    $$ {\rm{F}}^{-1} = {\rm{K}}^{-1}-\alpha{{\rm{K}}^{-1}}({\rm{G}}_1+{\rm{G}}_2){\rm{K}}^{-1}+\alpha^2{\rm{K}}^{-1}{{\rm{G}}_1}{\rm{K}}^{-1}{{\rm{G}}_1}{\rm{K}}^{-1} $$

    则动力学边界条件可进一步写为

    $$\begin{split} & {\xi _t} + {{\rm{K}}^{ - 1}}(1 + {\partial _{xxxx}})\eta = \\ &\;\;\;\; - \frac{1}{2}{{\rm{K}}^{ - 1}} \left\{ {{{\xi _x}^2 - {{\rm{G}}_0}\xi \left[({{\rm{G}}_0}\xi ) - 2{{\rm{G}}_0}\eta {{\rm{G}}_0}\xi - 2\eta {\xi _{xx}}\right]}} \right\} -\\ &\;\;\;\; {{\rm{F}}^{ - 1}}\left( {\frac{1}{2}{\kappa ^3} + {\kappa _{ss}}} \right) + {{\rm{K}}^{ - 1}}{\eta _{xxxx}} + \\ &\;\;\;\; \frac{\alpha }{2}{{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}\left[ {{\xi _x}^2 - {{({{\rm{G}}_0}\xi )}^2}} \right]-\;\;\;\;\;\;\;\\ &\;\;\;\; \; \alpha {{\rm{K}}^{ - 1}}[{{\rm{G}}_1}({{\rm{G}}_0}\xi )\xi + {{\rm{G}}_1}({{\rm{G}}_1}\xi )\xi + {{\rm{G}}_0}(\eta {{\rm{G}}_0}\xi ){\xi _{xx}} + \\ &\;\;\;\; {\partial _{xx}}\eta {({{\rm{G}}_0}\xi )^2} + {{\rm{G}}_0}({{\rm{G}}_0}\xi ){{\rm{G}}_0}\eta {{\rm{G}}_0}\xi + {{\rm{G}}_0}\eta {{\rm{G}}_0}{({{\rm{G}}_0}\xi )^2}] - \\ &\;\;\;\; {\alpha ^2}{{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}\eta + {\alpha ^2}{{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}{{\rm{G}}_1}({{\rm{G}}_0}\xi )\xi + \\ &\;\;\;\; \alpha {{\rm{K}}^{ - 1}}({{\rm{G}}_1} + {{\rm{G}}_2}){{\rm{K}}^{ - 1}}\eta - \beta {{\rm{K}}^{ - 1}}({{\rm{G}}_0} + {{\rm{G}}_1} + {{\rm{G}}_2})\xi + \\ &\;\;\;\; \alpha \beta {{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}({{\rm{G}}_0} + {{\rm{G}}_1})\xi + \alpha \beta {{\rm{K}}^{ - 1}}{{\rm{G}}_2}{{\rm{K}}^{ - 1}}{{\rm{G}}_0}\xi - \\ &\;\;\;\; {\alpha ^2}\beta {{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}{{\rm{G}}_1}{{\rm{K}}^{ - 1}}{{\rm{G}}_0}\xi + {{\rm{F}}^{ - 1}}p \end{split}$$ (7)

    这里将式(5)和式(7)称为三阶截断模型.

    本文将通过拟谱法在周期域上对三阶截断模型进行数值求解, 所有的导数项均在傅里叶空间中计算, 而非线性项则在物理空间中进行计算. 首先对方程式(5)和式(7)进行傅里叶变换

    $$ \begin{array}{l} \widehat{\eta}_t-k\tanh(kh)\widehat{\xi} = \widehat{N_1}\\ \widehat{\xi}_t+\dfrac{1+k^4}{1+\alpha{k\tanh(kh)}}\widehat{\eta} = \widehat{N_2} \end{array} $$

    $ {\cal{R}} = \sqrt{\dfrac{k\tanh(kh)(1+k^4)}{1+\alpha{k\tanh(kh)}}} $, 则原三阶截断模型可写为

    $$ \begin{array}{l} \widehat{p} = \widehat{\eta}+\dfrac{{ \text{i}}k\tanh(kh)}{{\cal{R}}}\widehat{\xi}\\ \widehat{q} = \widehat{\eta}-\dfrac{{ \text{i}}k\tanh(kh)}{{\cal{R}}}\widehat{\xi} \end{array} $$

    进一步地

    $$ \begin{array}{l} \widehat{p}_t+{ \text{i}}{\cal{R}}p = \widehat{N_1}+\dfrac{{ \text{i}}k\tanh(kh)}{{\cal{R}}}\widehat{N_2}\\ \widehat{q}_t-{ \text{i}}{\cal{R}}q = \widehat{N_1}-\dfrac{{ \text{i}}k\tanh(kh)}{{\cal{R}}}\widehat{N_2} \end{array} $$

    由于$ \xi $$ \eta $均为实数, 关于$ p $, $ q $的两方程实际上是等价的, 因此对三阶截断模型的求解可以简化为对下式的求解

    $$ \widehat{p_t}+{ \text{i}}{\cal{R}}p = \widehat{N}(p) $$ (8)

    在数值求解方程式(8)后, 可以进一步得到$ \xi $$ \eta $

    $$ \begin{split} & \widehat{\eta} = \frac{1}{2}\left[\widehat{p}(k)+{\widehat{p}(-k)}^*\right]\\ &\widehat{\xi} = \frac{{\cal{R}}}{2{ \text{i}}k\tanh(kh)}\left[\widehat{p}(k)-{\widehat{p}(-k)}^*\right] \end{split} $$

    在计算行波解时令$ \widehat{p_t} = -{ \text{i}}ck\widehat{p} $, 得到关于傅里叶系数的非线性代数方程, 继而利用牛顿迭代法对该非线性系统进行求解. 对于时间依赖的数值计算, 本文结合积分因子采用经典的四阶龙格-库塔方法进行时间积分[24], 同时在计算中通过对傅里叶模态加倍来解决由于混叠效应带来的数值误差.

    本节从自由波包型孤立波解入手, 通过比较完全欧拉方程、三阶截断模型和二阶截断模型的孤立波解, 验证三阶截断模型的准确性, 说明三阶截断模型相比于二阶截断模型[19]更为合理. 在弱非线性理论中, 波包型孤立波存在的条件是: 群速度与相速度在非色散点相等且在此波数下对应波包方程为焦聚型. 在不考虑惯性项的情况下, 水弹性波所对应的三阶非线性薛定谔方程(NLS)的性质会在临界水深$ h_c = 233 $发生变化[22]. 当水深小于临界深度时, NLS是焦聚型, 波包型孤立波从振幅无限小的周期波分岔而来. 而当水深大于临界深度时, NLS变为焦散型, 此时只存在有限振幅的波包型孤立波解. 为了更好地理解自由波包型孤立波解的存在性并对三阶截断模型的准确性进行验证, 首先令三阶截断模型中的黏性系数$ \beta = 0 $, 在无外力作用下推导考虑冰层惯性的三阶NLS方程. 对弱非线性水弹性波, 假设$ \eta\sim{O(\varepsilon)} $, $ \phi\sim{O(\varepsilon)} $, 这里$ \varepsilon $是用来衡量波陡的小参数. 由此, 自由面上的势函数可在$ y = 0 $处展开为

    $$ \phi(x,\eta,t) = \phi(x,0,t)+\eta\phi_{y}(x,0,t)+\frac{\eta^2}{2}\phi_{yy}(x,0,t)+O(\varepsilon^4) $$

    因此自由面上的运动学和动力学条件在$ y = 0 $处展开为

    $$ \begin{split} &\eta_{t}-\phi_{y} = \eta\phi_{yy}-\eta_{x}\phi_{x}+ \frac{\eta^2}{2}\phi_{yyy}-\eta\eta_{x}\phi_{xy}\\ & \phi_{t}+\eta+\eta_{xxxx}+\alpha\eta_{tt} = -\eta\phi_{ty}-\frac{\eta^2}{2}\phi_{tyy}-\\ & \qquad \frac{1}{2}\left(\phi_{x}^2+\phi_{y}^2\right)-\eta\left(\phi_{x}\phi_{xy}+\phi_{y}\phi_{yy}\right)\nonumber+\\ & \qquad \frac{5}{2}\eta_{x}^2\eta_{xxxx}+\frac{5}{2}\eta_{xx}^3+10\eta_{x}\eta_{xx}\eta_{xx} \end{split} $$

    为了得到波包的控制方程, 引入慢变量$ X = \varepsilon{x} $, $ T = \varepsilon{t} $$ \tau = \varepsilon^2{t} $. 设解的形式为

    $$ \begin{split} \eta = & \varepsilon{A_{11}(X-c_{\mathrm{g}}T,\tau)}{\rm{e}}^{{ \text{i}}\theta}+\\ &\sum\limits_{n = 2}^{\infty}\varepsilon^n\sum\limits_{j = 0}^{n}A_{nj}(X-c_{\mathrm{g}}T,\tau){\rm{e}}^{j{ \text{i}}\theta}+c.c. \end{split} $$ (9)
    $$ \phi = \sum\limits_{n = 1}^{\infty}\varepsilon^n\sum\limits_{j = 0}^{n}\phi_{nj}(X-c_{\mathrm{g}}T,y,\tau)+c.c. $$ (10)

    这里$ \theta = kx-\omega{t} $, $ c.c. $代表复共轭, $ \text{i} $为虚数单位. 将解式(9)和式(10)代入运动学和动力学边界条件中, 在$ O(\varepsilon) $阶得到

    $$ \begin{array}{l} { \text{i}}\omega{A_{11}}+k\sinh(kh)\varphi_{11} = 0\\ (1+k^4-\alpha\omega^2)A_{11}-{ \text{i}}w\cosh(kh)\varphi_{11} = 0 \end{array} $$

    从上式可知$ \varphi_{11} = \phi_{11}/\cosh(kh) $. 同时要使上述方程存在非零解需要满足可解性条件

    $$ \left[{1+\alpha{k\tanh(kh)}}\right]\omega^2 = (k+k^5)\tanh(kh) $$ (11)

    上式为水弹性波的色散关系式.

    $ O(\varepsilon^2) $阶, 可以得到

    $$ A_{11T}+c_{\mathrm{g}}A_{11X} = 0 $$ (12)

    其中

    $$ c_{\mathrm{g}} = \dfrac{(1+5k^4-\alpha\omega^2)\tanh(kh)+(1+k^4-\alpha\omega^2)kh\mathrm{sech}^2(kh)}{2\omega+2\alpha{k\tanh(kh)}} $$

    为水弹性波的群速度. 在$ O(\varepsilon^3) $阶时, 由可解性条件经过繁复的计算最终可以得到波包的控制方程, 即三阶非线性薛定谔方程

    $$ { \text{i}}A_{11\tau}+\frac{\omega_{{kk}}}{2}{A_{11XX}}+\frac{\alpha_1+\alpha_2+\alpha_3+\alpha_4}{2\omega^2\left[\coth(kh)+\alpha{k}\right]}{|A_{11}|^2A_{11}} = 0 $$ (13)

    其中

    $$ \begin{split} & \alpha_1 = -k\omega\left[\frac{h\omega^2}{c_g\sinh^2(kh)}+2\omega \coth(kh)\right]\cdot\\ &\qquad\bigg[\frac{2\omega}{\tanh(kh)}+\frac{c_g\omega^2}{\sinh^2(kh)}\bigg]\Big/(c_g^2-h)\\ & \alpha_2 = -\frac{2k\omega^4}{c_g\tanh(kh)\sinh^2(kh)}\\ & \alpha_3 = -\frac{k\omega^3}{\sinh^2(kh)}\frac{D_1}{D_0}-2k^2\omega^2\sinh(2kh)\frac{D_2}{{ \text{i}}D_0}+\\ &\qquad4k^2\omega^2\coth(kh)\cosh(2kh)\frac{D_2}{{ \text{i}}D_0}\\ &\alpha_4 = -4k^2\omega^3\coth(kh)+5k^7\omega \end{split} $$

    $ D_0 $, $ D_1 $$ D_2 $的具体表达式见附录. 为了后文讨论方便, 分别将色散项和非线性项的系数记作$\lambda,$ $ \gamma $. 令惯性项系数为0, 得到的三阶NLS与Milewski和Wang[25]所推导的方程退化到二维情形下的形式是一致的. 考虑深水情况$ h\to\infty $, 则此时色散项和非线性项的系数变为

    $$ \begin{array}{l} \lambda = \dfrac{1}{2\omega(1+\alpha{k})}\left[20k^3-4\alpha\omega\omega_{{k}}-2(1+\alpha{k})\omega_{{k}}^2\right]\\ \gamma = \dfrac{1}{2\omega(1+\alpha{k})}\left(\omega^2k^2\dfrac{-24\alpha{k^5}+6\alpha{k}-26k^4+4}{14k^4+12\alpha{k^5}-3\alpha{k}-1}+5k^7\right) \end{array} $$

    需要指出的是, 二阶截断模型[19]保留了原始欧拉方程的色散关系, 对应三阶NLS方程的色散项系数与本文所推导的是一致的, 但其在对原始方程进行近似时仅将非线性项保留至二阶, 因此对应的NLS中非线性项的系数是不准确的.

    为了验证三阶截断模型的准确性, 本文重点关注了水深$ h = 6.8\;\text{m} $$ h = 250\;\text{m} $时的NLS方程性质及孤立波解. 在不考虑惯性项影响时, 水深$ h = 6.8\;\text{m} $对应焦聚型NLS, 水深$ h = 250\;\text{m} $对应焦散型NLS. 现在在方程中加入惯性项, 在水深较小的情形下($ h = 6.8\;\text{m}<h_\mathrm{c} $)参考Takizawa[1]实验中冰层的物理参数取惯性系数$ \alpha = 0.069 $, 考虑惯性后的NLS仍然为焦聚型但非线性项的系数有所增大, 孤立波由振幅无限小的周期波分岔而来, 理论上此时的三阶截断模型会具有较高的精度. 进一步地, 数值计算波包型孤立波解以证明这一点. 对于完全欧拉方程, 采用保形映射将原不规则物理域映射为规则的矩形区域再结合牛顿迭代求解(关于此算法的详细介绍可见文献[22]). 将由数值计算得到的二阶截断模型和三阶截断模型的分岔曲线与完全欧拉方程的分岔曲线进行对比(图1). 在较大振幅范围内, 三阶截断模型所对应的分岔曲线均与完全欧拉方程的分岔曲线吻合得很好, 充分说明了三阶截断模型在该水深情况下具有较高的精度, 是一个好的近似模型. 而二阶截断模型的分岔曲线在分岔点附近与完全欧拉方程的分岔曲线相比仍然具有较为明显的差异, 这表明二阶截断模型的精度较差. 同时对相同波速下二阶截断模型、三阶截断模型以及完全欧拉方程的孤立波解进行比较后发现, 三阶截断模型的解与完全欧拉方程的解基本吻合, 两者之间的相对振幅差异约为$ 10^{-2} $, 而二阶截断模型的解则与完全欧拉的解相差较大.

    对水深较大的情况($ h = 250\;\mathrm{m}>h_\mathrm{c} $), 在参照Squire等[2]实验中的冰层参数后取惯性系数$\alpha = $$ 0.069$. 此时由于惯性项的影响, NLS方程由焦散型变为焦聚型, 因此三阶截断模型在水深较大的情况下仍然会具有较高的准确性. 同样地, 数值计算了该水深下二阶截断模型、三阶截断模型及完全欧拉方程的孤立波解及分岔曲线(图2). 从图中可以看出, 当孤立波振幅较小时三阶截断模型与完全欧拉方程的分岔曲线基本吻合, 但在中等振幅时两者会出现一定的差异, 这与水深$ h = 6.8\;\text{m} $时的情况略有不同. 产生这种现象的原因在于三阶NLS方程非线性项的系数: 当$ h = 250\;\text{m} $时, 虽然由于惯性项的加入使得方程由焦散型变为焦聚型, 但此时的非线性项系数$ \gamma $是一个非常接近于0的正数, 焦聚的特点表现得不明显, 这导致三阶截断模型的准确性相比于水深$ h = 6.8\;\mathrm{m} $时略有下降. 但在此水深情况下, 二阶截断模型与完全欧拉方程之间的差距较水深$ h = 6.8\;\text{m} $时更大, 完全无法反映出完全欧拉方程在对应波速处的波高, 其精度远小于三阶截断模型, 因此在水深较大的情况下使用二阶截断模型是不合适的. 图2(b)进一步展示了波速$c = 1.287\;8$时, 二阶截断模型、三阶截断模型和完全欧拉方程的孤立波解波形比对比. 因此综合两种水深下三阶截断模型和二阶截断模型与欧拉方程之间的对比, 三阶截断模型显然具有较高的准确性, 是一个更为合理的简化模型.

    图  2  (a) 水深$ h = 250\;\mathrm{m} $, 惯性$ \alpha = 0.069 $时, 孤立波解分岔图. 实线: 完全欧拉方程的分岔曲线, 短划线: 三阶截断模型的分岔曲线, 点线: 二阶截断模型的分岔曲线. (b) 波速$ c = 1.287\;8 $时, 欧拉方程、三阶截断模型和二阶截断模型的孤立波解. 实线: 完全欧拉方程, 短划线: 三阶截断模型, 点线: 二阶截断模型
    Figure  2.  The bifurcation curves of wavepacket solitary waves for $ h = 250\;\text{m} $ and $ \alpha = 0.069 $. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model. (b) The typical profiles for $ c = 1.287\;8 $. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model
    图  1  (a) 水深$ h = 6.8\;\text{m} $, 惯性系数$ \alpha = 0.069 $时, 孤立波解分岔图. 实线: 完全欧拉方程的分岔曲线, 短划线: 三阶截断模型的分岔曲线, 点线: 二阶截断模型的分岔曲线. (b) 波速$c = 1.287\;8$时, 欧拉方程、三阶截断模型和二阶截断模型的孤立波解. 实线: 完全欧拉方程, 短划线: 三阶截断模型, 点线: 二阶截断模型
    Figure  1.  (a) Bifurcation curves of wavepacket solitary waves for $ h = 6.8\;\text{m} $ and $ \alpha = 0.069 $. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model. (b) The typical profiles for $c = 1.287\;8$. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model

    基于三阶截断模型式(5)和式(7), 接下来将重点对匀速载荷作用下的冰层响应进行计算, 并将数值结果分别与实验结果进行对比. 实验中各物理参数的单位及具体取值如表1所示. 首先关注Takizawa[1]的实验结果. 该实验是在日本北海道的Lake Saroma完成的, 根据实验中所记录的参数可确定惯性项系数$\alpha = 0.069\;2$, 而外加载荷$ p(x-Ut) $则通过数值拟合静止载荷下的冰层挠度来确定. 在实验中Takizawa观察并记录了载荷经过监测点的时间与监测点处冰层最大垂直位移出现时间的差异, 这种滞后效应一般考虑是由于黏性效应造成的. 在实验中黏性的来源是多样的, 包括冰层本身的黏性, 覆盖在冰层上的积雪的影响以及冰层与流体界面处的边界层效应等. Hosking等[16]及Wang等[17]引入双参数记忆函数建立黏弹性模型对运动载荷下的浮冰动力学响应进行研究. 基于$ V<c_{\text{min}} $时滞后时间基本保持不变的特点, Hosking等[16]通过对多组数据进行试验并与实验记录进行对比确定了最佳的记忆函数参数$ A_0 $$ a_0 $. Dinvay等[19]假设阻尼与垂直速度成正比, 利用耗散项$ -b\eta_{t} $拟合黏性效应, 并将黏性系数视为待定参数, 在数值计算中通过反复试验优化滞后时间来确定$ b $. 除此之外, 这类确定待定系数的方法也在毛细重力波问题中得到了应用. Cho等[26]等在考虑黏性耗散对稳定及瞬时波形的影响时也将耗散项系数作为一个待定参数并通过将数值与实验结果进行拟合对比来确定具体取值. 本文采用较为简化的线性耗散模型$ -\beta\eta_{t} $描述实验中的黏性效应, $ \beta $的确定与上述参数的确定方法类似[16, 19, 26]: 选取多组黏性系数并数值计算同一黏性系数下不同载荷速度对应的波形和滞后时间, 将数值结果与实验记录对比直至两者能够较好吻合从而确定最佳黏性系数$ \beta $. 经过反复实验确定黏性系数$ \beta = 0.4 $, 此时数值结果与Takizawa[1]的实验结果吻合得最好.

    表  1  浅水[1]及深水[2]实验中各物理参数的取值及单位
    Table  1.  Values and units of the physical parameters in the shallow-water[1] and deep-water[2] experiments
    Physical parameterSymbolLake SaromaMcMurdo sound
    Young's modulus$E/( { {\rm{N} } \cdot { {\rm{m} }^{ - 2} } } )$$5.1\times{10^8}$$4.2\times{10^9}$
    Poisson ratio$\nu$$0.3$$0.3$
    ice thickness$d/{\rm{m}}$$0.17$$1.6$
    water depth$h/{\rm{m}}$$6.8$$\infty$
    flexural rigidity$D/({ {\rm{N} } } \cdot { {\rm{m} } } )$$2.35\times{10^5}$$1.8\times{10^9}$
    water density$\rho_{w}/(\mathrm{kg \cdot m^{-3} })$$1024$$1026$
    ice density$\rho_{i}/(\mathrm{kg \cdot m^{-3} })$$917$$917$
    下载: 导出CSV 
    | 显示表格

    图3展示了不同运动速度的载荷作用下的冰层动力学响应, 总体而言数值结果较为准确反映了实验结果的主要特征. 在准静态情形下, 即$ U = 2.2\;\mathrm{m/s} $$ U = 4.2\;\mathrm{m/s} $, 数值结果与实验结果吻合得非常好. 当载荷移动速度增大至$ U = 5.5\;\mathrm{m/s} $, 此时速度接近于最小相速度, 载荷两侧开始有波形成的趋势. 载荷速度为$ 6.2\;\mathrm{m/s} $时, 载荷前后出现明显的波动, 位于载荷前端的波长较短而后侧的波长较长, 三阶截断模型的数值计算结果较好地捕捉到了这两类波. 这说明三阶截断模型能够准确地反映出以最小相速度为分界点, 不同速度下冰层响应的差异. 进一步增大载荷移动速度至$ U = 8.9\;\mathrm{m/s} $, 此时载荷速度大于浅水中的重力波速度$ \sqrt{gH} $, 载荷后侧的重力波消失, 出现“声影区”[6]. 图3(c)~图3(e)中的冰层最大挠度与实验结果基本一致, 但载荷两侧的波形与实验记录存在一些差异. 这可能是因为在数值计算中所采用的黏性项较为简单而实验中黏性来源相对复杂, 随着载荷运动速度的增大, $ -\beta\eta_{t} $与实际黏性项产生的效果有所差异. 除此之外, 本文在建模时假设冰层厚度均匀, 但在实验中冰层厚度不均. 这些因素都可能会对载荷作用下产生的水弹性波波形造成影响. 除了对比不同载荷速度下的冰层垂直位移, 实验中出现的滞后效应在数值模拟中也得到了关注. 实验中Takizawa[1]记录了不同速度下载荷经过测量点的时刻及该时刻下载荷在$ y $方向上的位置, 在图3中以红点表示. 在数值计算中选取任意一点为观测点后可得到该点处冰层垂直位移的时历曲线, 同时记录移动载荷通过该点的时间及对应的垂向位置并在图3中以蓝点表示, 从图中可以看出数值计算得到的滞后时间在各个速度下都与实验结果较为吻合.

    图  3  三阶截断模型的数值结果与文献[1]实验中所记录的冰层形变对比. 其中虚线代表Takizawa的实验结果, 红点代表实验中的运动载荷经过冰层形变测量仪的垂直位置. 实线表示三阶截断模型式(5)和式(7)的数值结果, 蓝点代表数值计算中载荷经过监测点的垂向位置
    Figure  3.  Comparisons of the numerical results of Eq. (3) the numerical results of cubic-truncation model and the experimental records of Ref. [1]. The dashed lines represent the experimental data, and the red dots indicates the position of the load as it passed the deflectometer. The solid lines shows the numerical results of Eqs. (5) and (7), and the blue dots indicates the z-position of the load as it passes the point where the time series are obtained

    在Squire等[2]的实验中, 冰层厚度$ d = 1.6\;\mathrm{m} $, 冰层弹性刚度$D = 1.8 \times 10^9\;\mathrm{N \cdot m}$, 移动载荷重量为$ 2100\;\mathrm{kg} $. 与Takizawa[1]不同的是, Squire等[2]在实验中将多个应变计安装在载荷运动轨道两侧, 测量的是移动载荷作用下的冰层应变. 因此需要对数值计算得到的冰层垂直位移$ \eta $进一步处理以得到冰层线应变

    $$ \sigma = \frac{d}{2}\kappa^{*} $$

    其中, $ d $为冰层厚度, $ \kappa^{*} $为无量纲前的冰层曲率[27].

    考虑要与Squire的实验进行对比, 在基于三阶截断模型进行数值计算时, 取${\rm{G}}_0 = (-\partial_{xx})^{1/2}$以对应于深水情形.图4展示了三阶截断模型的数值计算结果与Squire等[2]实验结果之间的对比, 载荷速度分别为$ 4.5\;\mathrm{m/s} $, $ 17.5\;\mathrm{m/s} $, $ 18.4\;\mathrm{m/s} $, $ 20.8\;\mathrm{m/s} $, 前两个载荷速度为亚临界, 后两个速度为超临界速度. 在亚临界速度下, 数值结果与实验结果吻合得比较好, 不仅能够拟合冰层应变的幅值, 而且能够反映出随着载荷速度的增大应变曲线宽度逐渐变窄的趋势. 当载荷以较为接近最小相速度的速度$ 17.5\;\mathrm{m/s} $运动时, 载荷两侧逐渐升高并形成凹槽. 在图4(c)和图4(d)中, 数值预测的最大应变值略大于实验结果, 尤其是应变曲线中载荷前端的波形. 造成差异的原因除了以上在浅水实验中所提及的两个因素外, 也可能是由于本文所采用的是二维模型, 但在实验中所测得的应变是三维情形下冰层的应变, 从而导致两者所得到的结果存在差别. 但在超临界速度下, 数值计算仍然较为有效地捕捉到了载荷前后所形成的弹性波和重力波. 在数值计算中同样观察到了最大应变值点与载荷通过点之间的时间差异, 但由于Squire等[2]只是指出了滞后效应的存在, 没有明确记录差异时间, 因此无法对滞后的时间进行对比.

    图  4  三阶截断模型的数值结果与文献[2]实验中所测量的冰层应变对比, 图中应变值量级为$ 10^{-6} $. 虚线代表Squire的实验结果, 实线表示三阶截断模型式(5)和式(7)的数值结果
    Figure  4.  The comparison of microstrain between the numerical results of cubic-truncation model and numerical records of Ref. [2]. The microstrain is of $ 10^{-6} $. Dashed lines: experimental records, solid lines: numerical approximation of the cubic-truncation model Eqs. (5) and (7)

    本文主要考虑非线性, 惯性及黏性效应的影响, 研究了移动载荷作用下的浮冰动力学响应. 通过对Dirichilet-Neumman算子进行展开将原始的完全非线性问题简化为仅与自由面上的变量相关的三阶截断模型. 接着重点关注了此二维水弹性问题中的自由孤立波解以验证三阶截断模型的准确性. 在不考虑黏性和外力载荷下, 由多重尺度方法推导了三阶NLS方程并基于此方程对波包型孤立波解的存在性和三阶截断模型的准确性进行了探究. 惯性项的引入使得三阶NLS方程在浅水和深水中均表现为焦聚型, 因此理论上三阶截断模型在任意水深下都能够较好地近似完全欧拉方程. 另一方面, 对二阶截断模型、三阶截断模型和完全欧拉方程的分岔曲线和孤立波解进行数值求解并对比. 数值结果表明任意水深下三阶截断模型在一定振幅范围内都能够与完全欧拉方程较好地吻合, 是一个好的近似模型, 而二阶截断模型则在任意水深下都与完全欧拉方程存在较大的差异, 精度较差. 进一步地, 利用所得到的三阶截断模型对匀速运动载荷作用下的冰层响应进行了计算, 并将数值结果与Takizawa[1](浅水情况)和Squire等[2](深水情况)实验分别进行了对比, 结果表明此三阶截断模型能够较好地拟合移动载荷作用下冰层的垂直位移及应变情况.

    附录$ \;$

    在推导考虑惯性项的三阶非线性薛定谔方程中, 将解式(9)和式(10)代入原始欧拉方程中, 得到各阶方程. 这里给出推导的大致思路以及最终得到的三阶非线性薛定谔方程的系数表达式. 首先解Laplace方程, 可以得到

    $$ (\phi_{nj})_{yy}-(jk)^2\phi_{nj} = Q_{nj},\qquad (\phi_{nj})_{y} = 0 \qquad \mathrm{at} \quad y = -h $$

    这里$ n $表示阶数, $ j $表示模态, $ Q_{nj} $表示更低阶的项. 于是可以得到

    $$ \phi_{nj} = \varphi_{nj}(X-c_\mathrm{g}T,\tau)+\phi_{nj}^{Q} $$

    其中$ \phi_{nj}^{Q} $表示由于非齐次项$ Q_{nj} $而带来的特解. 将此解代入自由面边界条件, 即可得到

    $$ \begin{array}{l} -j{ \text{i}}\omega{A_{nj}}-jk\sinh(jkh)\varphi_{nj} = P_{nj}\\ \left[1+(j\text{i}k)^4-\alpha(j{ \text{i}}w)^2 \right] A_{nj}-j{ \text{i}}\omega\cosh(jkh)\varphi_{nj} = R_{nj} \end{array} $$

    这里$ P_{nj} $$ R_{nj} $仍然表示低阶项. 在$ O(\varepsilon) $阶, 由这两个方程的可解性条件, 可以得到色散关系式.

    $$ \left[1+\alpha{k\tanh(kh)}\right]\omega^2 = (k+k^5)\tanh(kh) $$

    $ O(\varepsilon^2) $阶, 可以得到群速度$ c_\mathrm{g} $

    $$ c_{\mathrm{g}} = \dfrac{(1+5k^4-\alpha\omega^2)\tanh(kh)+(1+k^4-\alpha\omega^2)kh\mathrm{sech}^2(kh)}{2\omega+2\alpha\omega{k\tanh(kh)}} $$

    不难验证$ c_\mathrm{g} = \omega_{k} $. 同时也可以得到$ A_{22} $$ \varphi_{22} $的表达式

    $$ A_{22} = \frac{D_1}{D_0},\qquad \varphi_{22} = \frac{D_2}{D_0} $$

    这里

    $$ \begin{split} & D_0 = -8{ \text{i}}k\omega^2\cosh(2kh)+2 \text{i}{k}\sinh(2kh)(2k+32k^5-8\alpha{k\omega^2})\nonumber\\ & D_1 = -8{ \text{i}}k^2\omega^2\cosh(2kh)\coth(kh)+\\ &\qquad 2{\rm{i}}{k^2}{\omega ^2}\sinh (2kh)\left[ {2 - \frac{1}{{{{\sinh }^2}(kh)}}} \right]\nonumber\\ & D_2 = -2k\omega\coth(kh)(2k+32k^5-8k\omega^2)+\\ &\qquad 2k{\omega ^3}\left[ {2 - \frac{1}{{{{\sinh }^2}(kh)}}} \right] \end{split} $$

    $ O(\varepsilon^3) $阶, 同样利用可解性条件, 在复杂的计算之后, 给出非线性薛定谔的系数表达式

    $$ \begin{array}{l} \lambda = \dfrac{1}{{2\omega \left[ {1 + \alpha k\tanh (kh)} \right]}}\left\{ {(10{k^3} - 2\alpha \omega {\omega _k})\tanh (kh) + } \right.\\ \qquad (1 + 5{k^4} - \alpha {\omega ^2})h{\rm{sec}}{{\rm{h}}^2}(kh) - 2\alpha kh\omega {\omega _k}{\rm{sec}}{{\rm{h}}^2}(kh) - \\ \qquad \left. {\;(1 + {k^4} - \alpha {\omega ^2})k{h^2}{\rm{sec}}{{\rm{h}}^2}(kh)\tanh (kh) - {\omega _k}\left[ {1 + \alpha k\tanh (kh)} \right]} \right\} \end{array} $$
  • 图  2   (a) 水深$ h = 250\;\mathrm{m} $, 惯性$ \alpha = 0.069 $时, 孤立波解分岔图. 实线: 完全欧拉方程的分岔曲线, 短划线: 三阶截断模型的分岔曲线, 点线: 二阶截断模型的分岔曲线. (b) 波速$ c = 1.287\;8 $时, 欧拉方程、三阶截断模型和二阶截断模型的孤立波解. 实线: 完全欧拉方程, 短划线: 三阶截断模型, 点线: 二阶截断模型

    Figure  2.   The bifurcation curves of wavepacket solitary waves for $ h = 250\;\text{m} $ and $ \alpha = 0.069 $. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model. (b) The typical profiles for $ c = 1.287\;8 $. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model

    图  1   (a) 水深$ h = 6.8\;\text{m} $, 惯性系数$ \alpha = 0.069 $时, 孤立波解分岔图. 实线: 完全欧拉方程的分岔曲线, 短划线: 三阶截断模型的分岔曲线, 点线: 二阶截断模型的分岔曲线. (b) 波速$c = 1.287\;8$时, 欧拉方程、三阶截断模型和二阶截断模型的孤立波解. 实线: 完全欧拉方程, 短划线: 三阶截断模型, 点线: 二阶截断模型

    Figure  1.   (a) Bifurcation curves of wavepacket solitary waves for $ h = 6.8\;\text{m} $ and $ \alpha = 0.069 $. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model. (b) The typical profiles for $c = 1.287\;8$. Solid line: full Euler equations, dashed line: third-truncation model, dotted line: quadratic-truncation model

    图  3   三阶截断模型的数值结果与文献[1]实验中所记录的冰层形变对比. 其中虚线代表Takizawa的实验结果, 红点代表实验中的运动载荷经过冰层形变测量仪的垂直位置. 实线表示三阶截断模型式(5)和式(7)的数值结果, 蓝点代表数值计算中载荷经过监测点的垂向位置

    Figure  3.   Comparisons of the numerical results of Eq. (3) the numerical results of cubic-truncation model and the experimental records of Ref. [1]. The dashed lines represent the experimental data, and the red dots indicates the position of the load as it passed the deflectometer. The solid lines shows the numerical results of Eqs. (5) and (7), and the blue dots indicates the z-position of the load as it passes the point where the time series are obtained

    图  4   三阶截断模型的数值结果与文献[2]实验中所测量的冰层应变对比, 图中应变值量级为$ 10^{-6} $. 虚线代表Squire的实验结果, 实线表示三阶截断模型式(5)和式(7)的数值结果

    Figure  4.   The comparison of microstrain between the numerical results of cubic-truncation model and numerical records of Ref. [2]. The microstrain is of $ 10^{-6} $. Dashed lines: experimental records, solid lines: numerical approximation of the cubic-truncation model Eqs. (5) and (7)

    表  1   浅水[1]及深水[2]实验中各物理参数的取值及单位

    Table  1   Values and units of the physical parameters in the shallow-water[1] and deep-water[2] experiments

    Physical parameterSymbolLake SaromaMcMurdo sound
    Young's modulus$E/( { {\rm{N} } \cdot { {\rm{m} }^{ - 2} } } )$$5.1\times{10^8}$$4.2\times{10^9}$
    Poisson ratio$\nu$$0.3$$0.3$
    ice thickness$d/{\rm{m}}$$0.17$$1.6$
    water depth$h/{\rm{m}}$$6.8$$\infty$
    flexural rigidity$D/({ {\rm{N} } } \cdot { {\rm{m} } } )$$2.35\times{10^5}$$1.8\times{10^9}$
    water density$\rho_{w}/(\mathrm{kg \cdot m^{-3} })$$1024$$1026$
    ice density$\rho_{i}/(\mathrm{kg \cdot m^{-3} })$$917$$917$
    下载: 导出CSV
  • [1]

    Takizawa T. Response of a floating sea ice sheet to a steadily moving load. Journal of Geophysical Research, 1988, 93: 5100-5112 doi: 10.1029/JC093iC05p05100

    [2]

    Squire VA, Robinson WH, Langhorne PJ, et al. Vehicles and aircraft on floating ice. Nature, 1988, 33: 159-161

    [3]

    Kheisin DY. Moving load on an elastic plate which floats on the surface of an ideal fluid. Izv. AN SSSR, OTN, Mekh. i Mashinostroenie, 1963, 1: 178-180 (in Russian

    [4]

    Nevel DE. Moving loads on a floating ice sheet. US Army CRREL Report, 1970: 261

    [5]

    Kerr AD. The critical velocities of a load moving on a floating ice plate that is subjected to implane force. Cold Regions Science and Technology, 1983, 6: 267-274 doi: 10.1016/0165-232X(83)90047-2

    [6]

    Davys JW, Hosking RJ. Waves due to a steadily moving source on a floating ice plate. Journal of Fluid Mechanics, 1985, 158: 269-287 doi: 10.1017/S0022112085002646

    [7]

    Babaei H, van der Sanden JJ, Short NH, et al. Lake ice cover deflection induced by moving vehicles: comparing theoretical results with satellite observations//New Research and Developments in Road Safety Session of the 2016 Conference of the Transportation Association of Canada Toronto, ON, 2016

    [8]

    van der Sanden JJ, Short NH. Radar satellites measure ice cover displacements induced by moving vehicles. Cold Regions Science and Technology, 2017, 533: 56-62

    [9]

    Schulkes RMSM, Hosking RJ, Sneyd AD. Waves due to a steadily moving source on a floating ice plate. Part 2. Journal of Fluid Mechanics, 1987, 180: 297-318

    [10]

    Kheisin DY. Some non-stationary problems of dynamics of the ice cove//Studies in Ice Physics and Ice Engineering (Iakolev ed.), Israel Program for Scientific Translations, 1971

    [11]

    Schulkes RMSM, Sneyd AD. Time-dependent response of floating ice to a steadily moving load. Journal of Fluid Mechanics, 1988, 186: 25-46 doi: 10.1017/S0022112088000023

    [12]

    Nugroho WS, Wang K, Hosking RJ, et al. Time-dependent response of a floating flexible plate to an impulsively started steadily moving load. Journal of Fluid Mechanics, 1999, 381: 337-355 doi: 10.1017/S0022112098003875

    [13]

    Miles J, Sneyd AD. The response of a floating ice sheet to an accelerating line load. Journal of Fluid Mechanics, 2003, 497: 435-439 doi: 10.1017/S002211200300675X

    [14]

    Wilson JT. Coupling between moving loads and flexural waves in floating ice sheets. US Army SIPRE Report, 1955

    [15]

    Beltaos S. Field studies on the response of floating ice sheets to moving loads. Canadian Journal of Civil Engineering, 1981, 8: 1-8 doi: 10.1139/l81-001

    [16]

    Hosking RJ, Sneyd AD, Waugh DW. Viscoelastic response of a floating ice plate to a steadily moving load. Journal of Fluid Mechanics, 1988, 196: 409-430 doi: 10.1017/S0022112088002757

    [17]

    Wang K, Hosking RJ, Milinazzo F. Time-dependent response of a floating viscoelastic plate to an impulsively started moving load. Journal of Fluid Mechanics, 2004, 521: 295-317 doi: 10.1017/S002211200400179X

    [18]

    Părău EI, Dias F. Nonlinear effects in the response of a floating ice plate to a moving load. Journal of Fluid Mechanics, 2002, 460: 281-305 doi: 10.1017/S0022112002008236

    [19]

    Dinvay E, Kalisch H, Părău EI. Fully dispersive models for moving loads on ice sheets. Journal of Fluid Mechanics, 2019, 876: 122-149 doi: 10.1017/jfm.2019.530

    [20]

    Toland JF. Heavy hydroelastic travelling waves. Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, 2007, 463: 2371-2397 doi: 10.1098/rspa.2007.1883

    [21] 浦俊, 卢东强. 三层流体中斜入射波作用下半无限板的水弹性响应. 力学学报, 2019, 51: 1614-1629 (Pu Jun, Lu Dongqiang. Hydroelastic response of a semi-infinite plate due to obliquely incident waves in a three-layer fluid. Chinese Journal of Theoretical and Applied Mechanics, 2019, 51: 1614-1629 (in Chinese) doi: 10.6052/0459-1879-19-081
    [22]

    Gao T, Vanden-Broeck JM, Wang Z. Numerical computations of two-dimensional flexural-gravity solitary waves on water of arbitrary depth. IMA Journal of Applied Mathematics, 2018, 83: 436-450 doi: 10.1137/S1064827597321532

    [23]

    Craig W, Sulem C. Numerical simulation of gravity waves. Journal of Computational Physics, 1993, 108: 73-83 doi: 10.1093/imamat/hxy007

    [24]

    Milewski PA, Tabak EG. A pseudospectral procedure for the solution of nonlinear wave equations with examples from free-surface flows. SIAM Journal on Scientific Computing, 1999, 21: 1102-1114 doi: 10.1006/jcph.1993.1164

    [25]

    Milewski PA, Wang Z. Three dimensional flexural-gravity waves. Studies in Applied Mathematics, 2013, 131: 135-148 doi: 10.1111/sapm.12005

    [26]

    Cho Y, Diorio JD, Akylas TR, et al. Resonantly forced gravity– capillary lumps on deep water. Part 2. Theoretical model. Journal of Fluid Mechanics, 2011, 672: 283-306

    [27]

    Squire VA. The breakp of shore fast sea ice. Cold Regions Science and Technology, 1993, 3: 211-218

  • 期刊类型引用(1)

    1. 倪宝玉,熊航,曾令东,孙林华. 移动压力面诱导的固壁冰航道内弯曲重力波传播特性研究. 力学学报. 2024(10): 2838-2853 . 本站查看

    其他类型引用(0)

图(4)  /  表(1)
计量
  • 文章访问数:  830
  • HTML全文浏览量:  269
  • PDF下载量:  111
  • 被引次数: 1
出版历程
  • 收稿日期:  2022-01-20
  • 录用日期:  2022-04-06
  • 网络出版日期:  2022-04-07
  • 刊出日期:  2022-04-17

目录

/

返回文章
返回