Load motion on an ice cover in the presence of a liquid layer with velocity shear

Texto integral

В настоящее время ведутся активные работы по освоению Арктики. Задача о поведении ледяного покрова под действием движущейся нагрузки изучается, в целях исследования, с одной стороны, способов разрушения ледяного покрова с помощью судов на воздушной подушке, а с другой – возможности использования ледяного покрова в качестве переправ, плавающих платформ различного назначения. Для безопасности работ на льду необходимо изучать закономерности его поведения при механическом воздействии.

В настоящее время задача о поведении ледяного покрова под действием прямолинейно движущейся нагрузки изучена достаточно хорошо при ее стационарном и нестационарном движении. Имеется значительное число работ о влиянии различных физических факторов на характер распространения изгибно-гравитационных волн в плавающем ледяном покрове. Однако влияние подледного течения со сдвигом скорости изучено недостаточно.

Данная статья является продолжением работы [1], в которой исследованы критические скорости и прогиб ледяного покрова с учетом равномерного сжатия в ближнем и дальнем поле при наличии течения с постоянным сдвигом скорости по всей толщине жидкости. Однако в реальных условиях скорость течения жидкости может меняться произвольным образом и по величине, и по направлению. Поэтому необходимы дальнейшие исследования с учетом таких факторов. Обзор предыдущих источников по указанной теме можно найти в работе [1].

Достаточно хорошо изучено влияние течения с вертикальным сдвигом скорости общего вида на поверхностные волны. Аналитические результаты по волнам при наличии течения со сдвигом скорости получены только для нескольких типов распределения скорости течения по глубине: линейного [2], экспоненциального [3], степенного со степенью 1/7 [4], в виде суммы гиперболических функций [5]. Различные приближенные модели использовались для описания взаимодействия волн и течений с применением разложений по малому параметру [6–10].

Другой приближенный метод решения задач при сдвиговом течении с произвольным профилем по глубине жидкости основан на приближении его ломаной линией и делении области жидкости на слои [11–16]. В каждом слое завихренность постоянна, и решение строится аналитически.

Сходимость такого приближения к точному решению при N → ∞ (N – число слоев) доказана в работе [15]. Метод кусочно-линейной аппроксимации прост и эффективен, но он приводит к системе уравнений порядка N + 2 Дисперсионное соотношение в этом случае приближается полиномом порядка N + 1, которое имеет дополнительные корни. Смысл этих дополнительных собственных значений и векторов выяснен в работах [14, 15]. При этом непрерывное распределение скорости течения заменяется кусочно-линейным распределением с наличием N – 1 тангенциальных разрывов. Устойчивость такого течения к возмущениям необходимо дополнительно исследовать.

Далее приведено решение задачи о поведении плавающего ледяного покрова под действием прямолинейно, с постоянной скоростью движущейся нагрузки при наличии верхнего слоя жидкости с линейным сдвигом скорости.

1. ПОСТАНОВКА ЗАДАЧИ

Рассматривается бесконечный ледяной покров на поверхности идеальной несжимаемой жидкости конечной глубины H. Предполагается, что имеется два слоя жидкости, из которых нижний слой толщиной H₂ покоится, а верхний толщиной H₁ = H - H₂ имеет постоянный сдвиг скорости g. Ледяной покров моделируется упругой бесконечной пластиной постоянной толщины, плывущей по поверхности жидкости.

Если в неподвижной системе координат ось Ox ′ направлена вдоль течения жидкости, то в подвижной декартовой системе координат, связанной с пластиной, жидкость течет в противоположном направлении.

Введем декартову систему координат Oxyz, связанную с пластиной, с центром O на верхней границе жидкости и осью Oz, направленной вертикально вверх. Вектор скорости жидкости можно представить в виде

$v(x,y,z,t)=U\left(z\right)+u,v,w,U\left(z\right)=\left\{\begin{array}{cc}\gamma z,& -H_1<z<\mathrm{0,}\\ -\gamma H_1,& -H<z<-H_1.\end{array}\right.$

Основное течение жидкости U(z) не потенциально, имеет постоянную завихренность g в верхнем слое; u, v, w – возмущенные компоненты скорости, малые по сравнению с основным течением. Предполагается, что заданное внешнее давление p₀ равномерно распределено по области прямоугольной формы шириной 2b и длиной 2a. Область давления движется с постоянной скоростью V под углом y к оси Ox. Давление внешней нагрузки равно p₀ = gM/(4ab), где M — масса движущегося тела, g — ускорение свободного падения. Данная нагрузка моделирует судно на воздушной подушке.

Рассмотрим движение жидкости. Скорость ее течения удовлетворяет уравнениям Эйлера и неразрывности:

$\frac{dv}{dt}=-\nabla P/\rho -g{i}_z,\nabla \cdot v=\mathrm{0,}P=-\rho gz+p.$ (1.1)

Здесь d/dt — полная производная по времени, P — гидродинамическое давление, которое складывается из гидростатической части и возмущенного гидродинамического давления p, r – плотность жидкости, i_z – единичный орт оси Oz. Прогиб ледяного покрова h(x, y, t) описывается уравнением изгиба тонких упругих пластин Кирхгофа–Лява

$D{\Delta }^2\eta +Q\Delta \eta +{\rho }_{0}h{d}^2\eta /d{t}^2=-\rho g\eta +p(x,y\mathrm{,0,}t)-p_0(x,y,t),$ (1.2)

$\Delta ={\partial }^2\text{}/\partial x^2+{\partial }^2\text{}/\partial y^2,D=\frac{E{h}^3}{12(1-{\nu }^2)},$

где D — цилиндрическая жесткость пластины, Q — сжимающие усилия в пластине, ρ₀, h — плотность и толщина льда, E — модуль Юнга, v — коэффициент Пуассона.

Ставятся граничные условия: условие непротекания на дне

$w(x,y,-H,t)=\mathrm{0,}$ (1.3)

а на верхней границе жидкости – кинематическое условие

$w(x,y\mathrm{,0,}t)=d\eta /dt$ (1.4)

и динамическое условие (1.2). На границе между слоями ставятся условия непрерывности вертикальной скорости и давления.

2. ДИСПЕРСИОННЫЕ СООТНОШЕНИЯ

Введем безразмерные переменные, параметры и функции по формулам

$$\begin{gathered} (x^{\text{'}},y^{\text{'}},z^{\text{'}},H^{\text{'}},{H^{\text{'}}}_1,{H^{\text{'}}}_2,a^{\text{'}},b^{\text{'}})=(x,y,z,H,H_1,H_2,a,b)/L,t^{\text{'}}=t\sqrt{g/L},L=a, \\ F=\frac{V}{\sqrt{gL}},S=\gamma \sqrt{L/g},\beta =\frac{D}{\rho g{L}^4},q=\frac{Q}{\rho g{L}^2},\chi =\frac{{\rho }_{0}h}{\rho L}, \\ (U^{\text{'}},u^{\text{'}},v^{\text{'}},w^{\text{'}})=(U,u,v,w)/\sqrt{gL},{\eta }^{\text{'}}=\eta /\text{}L,(p^{\text{'}},{p^{\text{'}}}_0)=(p,p_0)/\left(\rho gL\right). \end{gathered}$$

Штрихи далее опускаем. Среди всех безразмерных параметров задачи параметр c, который характеризует инерцию пластины, является малым и несущественным, инерция пластины мала по сравнению с инерцией жидкости. Для упрощения выражений можно положить его равным нулю.

Введем вертикальную координату в каждом слое:

$z_1=z,z_2=z+H_1,$

так что на верхней границе каждого слоя z _j = 0, а на нижней границе z _j = -H_j . Скорость течения в каждом слое имеет вид

$U_1\left(z_1\right)=S_1{z}_1,U_2=-S_1{H}_1.$

Рассмотрим волновые линейные возмущения вида

$$\begin{gathered} (\eta ,u,v,w,p)=(\overline{\eta },\overline{u}(z),\overline{v}(z),\overline{w}(z),\overline{p}(z\left)\right)e{}^{i(k\cdot r-\Omega t)}, \\ \Omega =\omega \sqrt{L/g},k=(k_x,k_y),r=(x,y), \end{gathered}$$

где w — частота в размерных переменных. Тогда система уравнений (1.1) линеаризуется и записывается в каждом слое в виде

$$\begin{gathered} i(-\Omega +k_x{U}_j(z_j\left)\right){\overline{u}}_j+S_j{\overline{w}}_j=-i{k}_x{\overline{p}}_j, \\ i(-\Omega +k_x{U}_j(z_j\left)\right){\overline{v}}_j=-i{k}_y{\overline{p}}_j, \\ i(-\Omega +k_x{U}_j(z_j\left)\right){\overline{w}}_j=-{{\overline{p}}^{\text{'}}}_j, \\ i{k}_x{\overline{u}}_j+i{k}_y{\overline{v}}_j+{{\overline{w}}^{\text{'}}}_j=\mathrm{0,} \end{gathered}$$ (2.1)

где j = 1, 2; S₂ = 0, штрих означает дифференцирование по z.

После преобразований в каждом слое получаем уравнение

${{{\overline{w}}^{\text{'}}}^{\text{'}}}_j-k^2{\overline{w}}_j=\mathrm{0,}{k}^2=k_x^2+k_y^2.$

Общее решение имеет вид

${\overline{w}}_j=A_{j}sh\left(k\right(z_j+H_j\left)\right)+B_{j}ch\left(k\right(z_j+H_j\left)\right).$

Для определения четырех неизвестных констант A_j, B_j имеем два уравнения на границе между слоями, а также условия на дне (1.3) и верхней границе (1.2). Из условия (1.3) и условия непрерывности вертикальной скорости следует

$B_2=\mathrm{0,}{B}_1=A_{2}shk{H}_2.$

Для давления получаем формулу

$\begin{array}{c}{k}^2{p}_j=\left[i\right(\Omega -k_x{U}_j\left(z_j\right)\left)kchk\right(z_j+H_j)+i{k}_x{S}_{j}shk(z_j+H_j\left)\right]A_j+\\ +\left[i\right(\Omega -k_x{U}_j\left(z_j\right)\left)kshk\right(z_j+H_j)+i{k}_x{S}_{j}chk(z_j+H_j\left)\right]B_j.\end{array}$ (2.2)

Из условия непрерывности давления на границе слоев находим ($z_1=-H_1, z_2=\mathrm{0,} S_2=\mathrm{0,} S_1=S$)

$i(\Omega +k_{x}S{H}_1)A_{1}k+i{k}_{x}S{B}_1=i(\Omega +k_{x}S{H}_1)k{A}_{2}chk{H}_2.$

Подставляем выражение для A₂ через B₁:

$i(\Omega +k_{x}S{H}_1)k{A}_1+i{k}_{x}S{B}_1-i(\Omega +k_{x}S{H}_1)kcthk{H}_2{B}_1=0.$ (2.3)

Уравнение (1.2) на верхней границе с помощью выражения для гидродинамического давления (2.2) и условия (1.4) при p₀ = 0 преобразуется к виду

$\begin{array}{c}\left[k^2\right(\beta k^4-q{k}^2+1)thk{H}_1-{\Omega }^{2}k-\Omega k_{x}Sthk{H}_1]A_1+\\ +\left[k^2\right(\beta k^4-q{k}^2+1)-{\Omega }^{2}kthk{H}_1-\Omega k_{x}S]B_1=0.\end{array}$ (2.4)

Уравнения (2.3), (2.4) образуют однородную систему линейных алгебраических уравнений. Чтобы она имела нетривиальное решение, определитель должен быть равен нулю. Этот определитель и является дисперсионным соотношением, определяющим зависимость частоты от волнового числа k и направления распространения. В данном случае это полином третьего порядка. Приравнивая его нулю и вводя цилиндрические координаты в плоскости волновых чисел $k_x=k\cos \alpha ,$ после преобразований получаем уравнение

$\begin{array}{c}{\Omega }^3+{\Omega }^{2}S\cos \alpha \left(k{H}_1+\frac{shk{H}_{1}chk{H}_2}{chkH}\right)-\\ -\Omega \left[k(\beta k^4-q{k}^2+1)thkH-S^2{\cos }^2\alpha \left(k{H}_{1}thkH-\frac{shk{H}_{1}shk{H}_2}{chkH}\right)\right]-\\ -k\left(\beta k^4-q{k}^2+1\right)S\cos \alpha \left(k{H}_{1}thkH-\frac{shk{H}_{1}shk{H}_2}{chkH}\right)=0.\end{array}$ (2.5)

Можно эту систему свести к обобщенной задаче на собственные значения для системы удвоенного порядка.

Проведенные расчеты (исходные данные приведены в разд. 4) показали, что в системе есть две частоты, соответствующие изгибно-гравитационным волнам в жидкости под упругой пластиной, а также дополнительная частота малой величины, соответствующая волнам на границе слоев. Различить частоты можно по их собственным векторам [14].

Собственные векторы, соответствующие изгибно-гравитационным волнам, имеют максимум на верхней поверхности, а собственный вектор дополнительной частоты имеет максимум на границе слоев. Пример таких векторов показан на рис. 1. Из рисунка видно, что векторы с номерами 1, 2 соответствуют изгибно-гравитационным волнам, а с номером 3 — волнам на границе слоев.

Рис. 1. Пример собственных векторов дисперсионного соотношения для системы уравнений (2.3), (2.4).

Исследовалась линейная устойчивость двуслойного течения к малым возмущениям в зависимости от толщины сдвигового слоя. Если кубическое уравнение (2.5) имеет только один действительный корень при некоторых параметрах задачи, то два других корня являются комплексно-сопряженными, одно из них дает растущее со временем решение.

3. РЕШЕНИЕ ЗАДАЧИ

Решение строится аналогично работе [14]. Поскольку мы ищем стационарное решение в системе координат, движущейся вместе с нагрузкой, в исходной системе координат, связанной с ледяной пластиной, возмущенные скорости движения жидкости имеют вид

$(u,v,w)=\left(u\right(r_1,z),v(r_1,z),w(r_1,z\left)\right),r_1=r-Ft,F=(F\cos \psi ,F\sin \psi ).$

Применяем преобразование Фурье по переменным x, y. Неизвестные функции u, v, w, h, p ищем в виде

$$\begin{gathered} (u,v,w)(r_1,z)=\frac{1}{4{\pi }^2}\underset{-\infty }{\overset{\infty }{\int }}\underset{-\infty }{\overset{\infty }{\int }}(\overline{u},\overline{v},\overline{w})(k,z)e{}^{ik\cdot r_1}d{k}_{x}d{k}_y,k=(k_x,k_y), \\ \eta \left(r_1\right)=\frac{1}{4{\pi }^2}\underset{-\infty }{\overset{\infty }{\int }}\underset{-\infty }{\overset{\infty }{\int }}\overline{\eta }\left(k\right)e{}^{ik\cdot r_1}d{k}_{x}d{k}_y, \\ p(r_1,z)=\frac{1}{4{\pi }^2}\underset{-\infty }{\overset{\infty }{\int }}\underset{-\infty }{\overset{\infty }{\int }}\overline{p}(k,z)e{}^{ik\cdot r_1}d{k}_{x}d{k}_y. \end{gathered}$$

Функцию  также представим в виде

$p_0\left(r_1\right)=\frac{1}{4{\pi }^2}\underset{-\infty }{\overset{\infty }{\int }}\underset{-\infty }{\overset{\infty }{\int }}{\overline{p}}_0\left(k\right)e{}^{ik\cdot r_1}d{k}_{x}d{k}_y.$

Решение строим таким же образом, как в разд. 2, только в уравнениях (2.1)–(2.3)  Уравнение движения пластины теперь содержит в правой части движущуюся область давления

$\begin{array}{c}\left[k^2(\beta k^4-q{k}^2+1)thk{H}_1-{(k\cdot F)}^{2}k-(k\cdot F)k_{x}Sthk{H}_1\right]A_1+\\ +\left[k^2(\beta k^4-q{k}^2+1)-{(k\cdot F)}^{2}kthk{H}_1-(k\cdot F)k_{x}S\right]B_1=i(k\cdot F)k^2{\overline{p}}_0/chk{H}_1.\end{array}$

После преобразований, вводя полярные координаты   получаем уравнение

$\begin{array}{c}\left[(\beta k^4-q{k}^2+1)\left(k-\frac{S\cos \alpha S_H}{S{H}_1\cos \alpha +F\cos (\alpha -\psi )}\right)-k^2{F}^2{\cos }^2(\alpha -\psi )cthkH+\right.\\ +\frac{SF\cos \alpha \cos (\alpha -\psi )}{S{H}_1\cos \alpha +F\cos (\alpha -\psi )}\left[kF\cos (\alpha -\psi )+S\cos \alpha thk{H}_1\right]C_H-\\ \left.\begin{array}{c}\\ \end{array}-kSF\cos \alpha \cos (\alpha -\psi )\right]B_1=\frac{i{k}^{2}F\cos (\alpha -\psi )shk{H}_2}{shkH}{\overline{p}}_0,\end{array}$

где

$S_H=\frac{shk{H}_{1}shk{H}_2}{shkH},C_H=\frac{chk{H}_{1}shk{H}_2}{shkH}.$

Перейдем в систему координат Ox₁ y₁z, движущуюся вместе с нагрузкой, повернутую на угол y относительно исходной системы, с осью Ox₁ направленной вдоль линии движения, осью Oy₁ — перпендикулярной к ней. Получаем

$f(k,\sigma )B_1=\frac{i{k}^{2}F\cos \sigma shk{H}_2}{shkH}{\overline{p}}_0,$ (3.1)

$\begin{array}{c}\sigma =\alpha -\psi ,\\ {\overline{p}}_0\left(k\right)=4p_0\frac{\sin \left(ka\cos \sigma \right)\sin \left(kb\sin \sigma \right)}{k^2\cos \sigma \sin \sigma },\end{array}$

$\begin{array}{c}f(k,\sigma )=\left(\beta k^4-q{k}^2+1\right)\left(k-\frac{S\cos (\sigma +\psi )S_H}{S{H}_1\cos (\sigma +\psi )+F\cos \sigma }\right)-k^2{F}^2{\cos }^2\text{}\sigma cthkH-\\ -kSF\cos \sigma \cos (\sigma +\psi )+\frac{SF\cos \sigma \cos (\sigma +\psi )}{S{H}_1\cos (\sigma +\psi )+F\cos \sigma }\left[kF\cos \sigma +S\cos (\sigma +\psi )thk{H}_1\right]C_H.\end{array}$

Из условий (1.4), (2.3) находим

$\overline{\eta }=\frac{A_{1}shk{H}_1+B_{1}chk{H}_1}{-i\left(kF\right)}=-\frac{{\overline{p}}_0\left[kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H)\right]}{f(k,\sigma )\left[F\cos \sigma +S{H}_1\cos (\sigma +\psi )\right]}.$

С помощью обратного преобразования Фурье из (3.1) получаем

$\eta (x_1,y_1)=-\frac{2p_0}{{\pi }^2}\underset{0}{\overset{\infty }{\int }}\underset{-\pi /2}{\overset{\pi /2}{\int }}\frac{\sin \left(ka\cos \sigma \right)\sin \left(kb\sin \sigma \right)}{k\cos \sigma \sin \sigma f_1(k,\sigma )}\cos \left(k\right(x_1\cos \sigma +y_1\sin \sigma \left)\right)g(k,\sigma )dkd\sigma ,$                       (3.2)

$$\begin{gathered} g(k,\sigma )=kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H), \\ \begin{array}{c}{f}_1(k,\sigma )=(\beta k^4-q{k}^2+1)\left[kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H)\right]-\\ -\left[k^2{F}^2{\cos }^2\sigma cthkH+kSF\cos \sigma \cos (\sigma +\psi )\right]\left[F\cos \sigma +S{H}_1\cos (\sigma +\psi )\right]+\\ +SF\cos \sigma \cos (\sigma +\psi )\left[kF\cos \sigma +S\cos (\sigma +\psi )thk{H}_1\right]C_H.\end{array} \end{gathered}$$

В частных случаях $\psi =\mathrm{0,}\pi$ формулы упрощаются: при $\psi =0$

$$\begin{gathered} g(k,\sigma )=kF+S(k{H}_1-S_H), \\ {f}_1(k,\sigma )=(\beta k^4-q{k}^2+1)\left[kF+S(k{H}_1-S_H)\right]-{\cos }^2\sigma \left[(k^2{F}^{2}cthkH+kSF)(F+S{H}_1)-SF(kF{C}_H+S{S}_H)\right]; \end{gathered}$$

при $\psi =\pi$

$g(k,\sigma )=kF-S(k{H}_1-S_H),$

$f_1(k,\sigma )=(\beta k^4-q{k}^2+1)\left[kF-S(k{H}_1-S_H)\right]-{\cos }^2\sigma \left[(k^2{F}^{2}cthkH-kSF)(F-S{H}_1)+SF(kF{C}_H-S{S}_H)\right].$

Известно, что существует критическая скорость движения нагрузки. При движении со скоростью, меньшей критической, функция $f_1(k,\sigma )$ не имеет корней, соответствующих изгибно-гравитационным модам. Волны в ледяном покрове не возникают, деформации льда локализованы вблизи области нагрузки. Для определения критической скорости исследуем корни функции f₁. Рассмотрим случай y = 0. Из условия $f_1(k,\sigma )=0$ получаем

${\cos }^2\sigma =\frac{(\beta k^4-q{k}^2+1)[kF+S(k{H}_1-S_H\left)\right]}{(k^2{F}^{2}cthkH+kSF)(F+S{H}_1)-SF(kF{C}_H+S{S}_H)}\le 1.$                                                         (3.3)

Критическое значение параметра $F_{*}$ достигается при ${\cos }^2\sigma =1.$ Получаем

$\begin{array}{c}G(k,F_{*})=\mathrm{0,}G(k,F)=F^3{k}^{2}cthkH+k{F}^{2}S(k{H}_{1}cthkH+1-C_H)+\\ +F\left[S^2(k{H}_1-S_H)-k(\beta k^4-q{k}^2+1)\right]-S(\beta k^4-q{k}^2+1)(k{H}_1-S_H).\end{array}$                                                 (3.4)

Поскольку $F_{*}$ является наименьшим значением, для которого выполнено уравнение (3.4), то $\partial F/\partial k=0.$ Дифференцируем уравнение (3.4) по k:

$\partial G/\partial k(k,F_{*})=0.$                                                                                                                                  (3.5)

Уравнения (3.4) и (3.5) с условием $\partial F\text{}/\partial k=0$ составляют систему для определения критических значений $F_{*},k_{*}.$ Система решалась методом Брауна [17].

При $\psi =\pi$ получаем аналогично:

${\cos }^2\sigma =\frac{(\beta k^4-q{k}^2+1)[kF-S(k{H}_1-S_H\left)\right]}{(k^2{F}^{2}cthkH-kSF)(F-S{H}_1)+SF(kF{C}_H-S{S}_H)}\le 1.$                                                            (3.6)

$\begin{array}{c}G(k,F_{*})=\mathrm{0,}G(k,F)=F^3{k}^{2}cthkH-k{F}^{2}S(k{H}_{1}cthkH+1-C_H)+\\ +F\left[S^2(k{H}_1-S_H)-k(\beta k^4-q{k}^2+1)\right]+S(\beta k^4-q{k}^2+1)(k{H}_1-S_H).\end{array}$                                              (3.7)

Из уравнений (3.5), (3.7) с условием $\partial F\text{}/\text{}\partial k=0$ находим критические значения $F_{*},k_{*}.$ В общем случае произвольного угла y уравнение $f_1(k,\sigma )=0$ преобразуется к виду

$a_1\left(k\right)tg{}^3\sigma +a_2\left(k\right)tg{}^2\sigma +a_3\left(k\right)tg\sigma +a_4\left(k\right)=\mathrm{0,}$ (3.8)

$$\begin{gathered} a_1\left(k\right)=S\sin \psi (\beta k^4-q{k}^2+1)(k{H}_1-S_H), \\ {a}_2\left(k\right)=F{S}^2{\sin }^2\psi (k{H}_1-S_H)-(\beta k^4-q{k}^2+1)\left[kF+S\cos \psi (k{H}_1-S_H)\right], \\ {a}_3\left(k\right)=S\sin \psi \left[(\beta k^4-q{k}^2+1)(k{H}_1-S_H)-k^2{F}^2{H}_{1}cthkH-kF(F+2S{H}_1\cos \psi )+F(kF+2S\cos \psi thk{H}_1)C_H\right]. \\ \begin{array}{c}{a}_4\left(k\right)=-(\beta k^4-q{k}^2+1)\left[kF+S\cos \psi (k{H}_1-S_H)\right]+F{S}^2{\cos }^2\psi (k{H}_1-S_H)+\\ +k^2{F}^{2}cthkH(F+S{H}_1\cos \psi )+k{F}^{2}S\cos \psi (1-C_H).\end{array} \end{gathered}$$

Поскольку это кубическое уравнение с вещественными коэффициентами, то оно всегда имеет действительный корень. При F > F_* имеется интервал значений параметра k, где уравнение (3.8) имеет три действительных корня. На границах этого интервала комплексно-сопряженные корни становятся действительными и соответственно кратными. Запишем это уравнение в виде

$a_1\left(k\right)t^3+a_2\left(k\right)t^2+a_3\left(k\right)t+a_4\left(k\right)=\mathrm{0,}t=tg\sigma .$

Тогда производная по t тоже равна нулю на границах интервала

$3a_1\left(k\right)t^2+2a_2\left(k\right)t+a_3\left(k\right)=\mathrm{0,}$

$t_1=\frac{-a_2+\sqrt{a_2^2-3a_1{a}_3}}{3a_1},t_2=-\frac{a_2+\sqrt{a_2^2-3a_1{a}_3}}{3a_1}.$ (3.9)

Подставляем значения t₁, t₂ в уравнение (3.8) и получаем два трансцендентных уравнения. Для t₁ корней нет, а для t₂ при F > F_* находятся два корня. При F = F_* эти два корня совпадают. Из условия кратности корней производная по k тоже обращается в нуль:

${a^{\text{'}}}_1\left(k\right)t_2^3+{a^{\text{'}}}_2\left(k\right)t_2+{a^{\text{'}}}_3\left(k\right)t_2+{a^{\text{'}}}_4\left(k\right)=0.$ (3.10)

Таким образом, получили опять систему двух уравнений: (3.8) при tg s = t₂ и (3.10), которая решается методом Брауна. Начальные приближения для поиска F_*, k_* находятся с помощью численных расчетов.

Множество всех корней функции f₁(k,s) при F > F_* образует в плоскости волновых чисел две кривые: C_k¹ и $C_k^2.$ C_k¹ соответствует изгибно-гравитационным волнам, а $C_k^2$ — волнам на границе слоев. Анализ функции f₁(k, s) показывает, что при F > F_* имеется интервал (s₁, s₂) угла s, где функция f₁ имеет два действительных корня соответствующих изгибно-гравитационным волнам. А также есть небольшой интервал (s₃, s₄), где функция f₁ имеет дополнительный корень k₀(s), соответствующий волнам на границе слоев. При $\psi =\mathrm{0,}\pi$ дополнительной моды нет.

При y = 0,p значения s₁,s₂ находим из формул (3.3), (3.6) как минимальные и максимальные величины. При произвольном угле y и F > F_* сначала находим c помощью уравнений (3.8), (3.9) интервал значений $(k_1^{*},k_2^{*}),$ в котором существуют изгибно-гравитационные волны, т.е. уравнение (3.8) имеет три корня. На этом интервале находим значения параметра k, при которых $\partial \sigma /\partial k=0$ для верхней и нижней кривой. Значения угла s в этих точках определяют (s₁, s₂). На рис. 2 изображены для примера кривые $C_k^1,C_k^2$ для ледяного покрова толщиной h = 1 м, H₁ = 5 м, H = 50 м, Q = 0, y = p/3, S = 0.4; 0.8; 1.2 (остальные исходные данные приведены в разд. 4). Из рисунка видно, что при S = 1.2 кривая C_k¹ проходит через начало координат.

Рис. 2. Кривые волновых чисел $C_k^1,C_k^2$ для ледяного покрова толщиной 50 м, Q=0, $\psi =\pi /3;$  1–3 – S=0.4, 0.8, 1.2.

Поскольку дополнительная мода соответствует действительному корню уравнения (3.8), который существует при любых значениях параметра k, то σ₃ и s₄ определяются из асимптотики при k → 0 и k → ∞. При k → 0 имеем

$$\begin{gathered} a_1\left(k\right)\approx kS{H}_1^2\sin \psi /H,a_2\left(k\right)\approx k(-F-S\cos \psi H_1^2\text{}/\text{}H+F{S}^2{\sin }^2\psi H_1^2\text{}/\text{}H), \\ {a}_3\left(k\right)\approx kS\sin \psi H_1[H_1-2F(F+S{H}_1\cos \psi \left)\right]/\text{}H, \\ {a}_4\left(k\right)\approx k[-F-S\cos \psi H_1^2\text{}/\text{}H+F{(F+S{H}_1\cos \psi )}^2\text{}/\text{}H]. \end{gathered}$$

Решаем уравнение (3.8) с этими коэффициентами, деленными на k, и получаем значение tg s₄. Если это уравнение имеет три действительных корня, то это означает, что кривая C_k¹ проходит через начало координат, и два других корня дают углы, под которыми кривая волновых чисел C_k¹ подходит к началу координат.

При k → ∞ для tg s₃ находим

$$\begin{gathered} a_1\left(k\right)\approx \beta k^{5}S{H}_1\sin \psi ,a_2\left(k\right)\approx -\beta k^5(F+S{H}_1\cos \psi ), \\ {a}_3\left(k\right)\approx \beta k^{5}S{H}_1\sin \psi ,a_4\left(k\right)\approx -\beta k^5(F+S{H}_1\cos \psi ). \end{gathered}$$

Получаем

$tg{\sigma }_3=\frac{F+S{H}_1\cos \psi }{S{H}_1\sin \psi }.$

При y → 0, p s₃ → p/2.

При F > F_* интеграл (3.2) является интегралом в смысле главного значения. Чтобы правильно выполнить условие излучения на бесконечности и найти стационарное решение, применялся метод, изложенный в работе [18]. В области нагрузки вводится нестационарное давление с малым параметром e: $p(x,y,t)=p_{0}e{}^{\epsilon t},\epsilon >0.$ Тогда полюсы сдвигаются с действительной оси в верхнюю или нижнюю полуплоскость. Этот сдвиг и определяет правило обхода. Стационарное решение определяется как предел решения, полученного при e → 0. Получено, что при вычислении интеграла (3.2) корни k₀, k₁ обходятся сверху, а корень k₂ – снизу.

Вклад дополнительной моды в прогиб ледяного покрова мал. Как видно из рис. 2, интервал (s₃, s₄) мал, и собственный вектор, соответствующий дополнительной моде, дает значение на верхней границе, близкое к нулю (см. рис. 1).

Как и в работе [1], картину изгибно-гравитационных волн в дальнем поле можно построить с помощью асимптотических методов. При этом дополнительную моду можно не учитывать. Гребни волн в дальнем поле на основе кривых C_k¹ строятся аналогичным образом, как в работе [1], и имеют такой же характер.

Волновое сопротивление F₁ и боковая сила F₂, действующие на движущееся тело, определяются по формуле

$(F_1,F_2)=-p_0\underset{-a}{\overset{a}{\int }}\underset{-b}{\overset{b}{\int }}\left(\frac{\partial \eta }{\partial x_1},\frac{\partial \eta }{\partial y_1}\right)d{x}_{1}d{y}_1.$

Подставляем выражение для прогиба (3.2) и интегрируем:

$$\begin{gathered} F_1=i\frac{4p_0^2}{{\pi }^2}\underset{0}{\overset{2\pi }{\int }}\underset{0}{\overset{\infty }{\int }}\frac{{\displaystyle {\sin }^2\left(ka\cos \sigma \right){\sin }^2\left(kb\sin \sigma \right)}}{k^2{\sin }^2\sigma \cos \sigma f_1(k,\sigma )}\left[kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H)\right]dkd\sigma , \\ {F}_1=i\frac{4p_0^2}{{\pi }^2}\underset{0}{\overset{2\pi }{\int }}\underset{0}{\overset{\infty }{\int }}\frac{{\displaystyle {\sin }^2\left(ka\cos \sigma \right){\sin }^2\left(kb\sin \sigma \right)}}{k^2{\sin }^2\sigma \cos \sigma f_1(k,\sigma )}\left[kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H)\right]dkd\sigma , \end{gathered}$$

В силу симметрии вклад в эти интегралы дают только вычеты в нулях знаменателя. Получаем

$\begin{array}{c}{F}_1=\frac{8p_0^2}{\pi }{\left.\left\{\underset{{\sigma }_1}{\overset{{\sigma }_2}{\int }}\left(\frac{{\sin }^2\left(ka\cos \sigma \right){\sin }^2\left(kb\sin \sigma \right)}{k^2{\sin }^2\sigma \cos \sigma {f^{\text{'}}}_1(k,\sigma )}\left[kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H)\right]\right)\right.\right|}_{k_2}^{k_1}\text{}d\sigma +\\ +\underset{{\sigma }_3}{\overset{{\sigma }_4}{\int }}\left.\frac{{\sin }^2\left(k_{0}a\cos \sigma \right){\sin }^2\left(k_{0}b\sin \sigma \right)}{k_0^2{\sin }^2\sigma \cos \sigma {f^{\text{'}}}_1(k_0,\sigma )}\left[k_{0}F\cos \sigma +S\cos (\sigma +\psi )(k_0{H}_1-S_H)\right]d\sigma \right\}.\end{array}$ (3.11)

$$\begin{gathered} F_2=\frac{8p_0^2}{\pi }{\left\{\underset{{\sigma }_1}{\overset{{\sigma }_2}{\int }}\left.\left(\frac{{\sin }^2\left(ka\cos \sigma \right){\sin }^2\left(kb\sin \sigma \right)}{k^2\sin \sigma {\cos }^2\sigma {f^{\text{'}}}_1(k,\sigma )}\left[kF\cos \sigma +S\cos (\sigma +\psi )(k{H}_1-S_H)\right]\right)\right|\right.}_{k_2}^{k_1}\text{}d\sigma + \\ +\underset{{\sigma }_3}{\overset{{\sigma }_4}{\int }}\left.\frac{{\sin }^2\left(k_{0}a\cos \sigma \right){\sin }^2\left(k_{0}b\sin \sigma \right)}{k_0^2\sin \sigma {\cos }^2\sigma {f^{\text{'}}}_1(k_0,\sigma )}\left[k_{0}F\cos \sigma +S\cos (\sigma +\psi )(k_0{H}_1-S_H)\right]d\sigma \right\}. \end{gathered}$$ (3.12)

Безразмерные коэффициенты волновых сил находятся по формуле

$(A_1,A_2)=\frac{g\rho }{2b{p}_0^2}(-F_1,F_2).$

3. ЧИСЛЕННЫЕ РЕЗУЛЬТАТЫ

Численные расчеты проводили для ледяного покрова при следующих входных параметрах задачи: E = 5 ГПа, r = 1025 кг/м³, ρ₀ = 922.5 кг/м³, v = 1/3, a = 20 м, b = 10 м, толщина льда h = 1 м, глубина жидкости H = 50 м. Сжимающие усилия в пластине заданы в виде $Q=\overline{Q}\sqrt{\rho gD},\overline{Q}=0$ и 1. Глубина сдвигового слоя H₁, скорость и направление движения нагрузки менялись.

Проведены расчеты на исследование линейной устойчивости возникающих изгибно-гравитационных волн к малым возмущениям при значениях параметра S в пределах от 0 до 1 и различных значенииях угла a, коэффициента сжатия $\overline{Q}$, толщины сдвигового слоя H₁ и волнового числа k с помощью уравнения (2.5).

Оказалось, что при $\overline{Q}=0$ решение устойчиво при всех значениях указанных параметров. С увеличением коэффициента cжатия Q появляются зоны неустойчивости при больших значениях глубины сдвигового слоя. Так, при $\overline{Q}=1$ неустойчивость обнаружена при S=0.2, H₁≥45 м; S=0.3, H₁≥35 м;  S=0.4, H₁≥29 м;  S=0.5, H₁≥24 м и т. д. С ростом параметров $\overline{Q}$ и S область неустойчивости растет. При $\overline{Q}=1.5$ решение неустойчиво c любой толщиной сдвигового слоя.

Проведенные расчеты показали, что наибольшее отличие критических скоростей наблюдается при углах y = 0 и y = p. При y = p/2 критические скорости при различных значениях параметра сдвига S близки к значению при S = 0 (в отсутствие течения). На рис. 3 приведены зависимости критических скоростей изгибно-гравитационных волн в ледяном покрове от глубины сдвигового слоя H₁ при y = 0, p; $\overline{Q}=\mathrm{0,1}$ и различных значениях параметра сдвига S.

Рис. 3. Зависимость критических скоростей от толщины сдвигового слоя при $\overline{Q}=0$ (а, б) и $\overline{Q}=1$ (в, г), Ψ = 0 (а, в) и ψ = π (б, г): кривые 1–5 соответствуют значениям S=0.1, 0.2, 0.3, 0.4, 0.5.

При $\overline{Q}=1$ кривые представлены в пределах устойчивости. Минимальные критические скорости достигаются при y = 0. Максимальные критические скорости наблюдаются при y = p. При увеличении коэффициента сжатия льда критические скорости значительно уменьшаются. Критические скорости уменьшаются с ростом толщины сдвигового слоя при y = 0 и растут при y = p. Видно, что сдвиг скорости, направление движения, коэффициент сжатия льда и толщина сдвигового слоя оказывают существенное влияние на критические скорости изгибно-гравитационных волн.

На рис. 4 показаны зависимости прогиба ледяного покрова вдоль линии движения нагрузки  от координаты x₁ при V = 20 м/с, в отсутствие сжатия $\overline{Q}=0$ (а, б) и $\overline{Q}=1$ (в, г), при y = 0 (а, в) и y = p (б, г) для случаев однослойной (H₁ = H) и двуслойной жидкости с H₁ = 5 м, S = 0, 0.3, 0.5.

Рис. 4. Прогиб ледяного покрова при $\overline{Q}$ = 0 (а, б) и $\overline{Q}$ = 1 (в, г), ψ = 0 (а, в) и ψ = π (б, г): 1 – S = 0; 2, 3 – S = 0.3, 0.5 для двух слоев H₁ = 5 м, H₂ = 45 м; 4, 5 – S = 0.3, 0.5 для одного слоя H₁ = H = 50 м.

Как видно из рисунка, волны перед нагрузкой для всех случаев практически одинаковой амплитуды. Только на рис. 4б для однослойной жидкости при S=0.5 (кривая 5) амплитуды волн и перед нагрузкой, и за нагрузкой значительно выше остальных. Это объясняется тем, что для этого случая скорость движения нагрузки близка к критической скорости V_* = 19.3 м/с.

При y = 0 волны перед нагрузкой укорачиваются, а за нагрузкой удлиняются в обоих случаях. При y = p волны перед нагрузкой удлиняются, а за нагрузкой укорачиваются. При y = 0 с увеличением параметра S амплитуда волн за нагрузкой становится меньше в обоих случаях: двуслойной и однослойной жидкости. При y = p с ростом параметра S все происходит наоборот – амплитуда волн за нагрузкой увеличивается в обоих случаях. Сжатие льда приводит к увеличению прогиба ледяного покрова. Увеличение толщины сдвигового слоя дает такой же результат, как увеличение параметра S – при y = 0 амплитуды прогиба за нагрузкой уменьшаются, а при y = p растут. При движении нагрузки под произвольным углом к направлению течения симметрия прогиба нарушается.

На рис. 5 изображены зависимости безразмерных коэффициентов волновых сил A₁, A₂ от скорости движения нагрузки при –Q = 0, y = 0, p/2, p для однослойной (H₁ = H) и двуслойной жидкости с H₁ = 5 м, S = 0, 0.3, 0.5.

Рис. 5. Безразмерные коэффициенты волновых сил A₁, A₂ в зависимости от скорости движения нагрузки при $\overline{Q}$ = 0, ψ = 0 (а), ψ = π/2 (б) и ψ = π (в): 1 – S = 0; 2, 3 – S = 0.3, 0.5 для двух слоев: H₁ = 5 м, H₂ = 45 м; 4, 5 – S = 0.3, 0.5 для одного слоя: H = 50 м.

Как видно из рисунка, графики сил подобны, основное отличие связано с различием критических скоростей. Поскольку при y = p/2 критические скорости мало меняются с изменением параметра S, все кривые для коэффициента продольной волновой силы A₁ практически сливаются. Коэффициент боковой силы A₂ максимален при y = p/2 и растет c увеличением параметра S.

На рис. 6 представлены аналогичные графики для случая –Q = 1. Видно, что коэффициенты волновых сил в этом случае в 1,5 раза выше, чем при –Q = 0. Продольная волновая сила уменьшается с ростом параметра S при y = 0 и увеличивается при y = p в области сверхкритических скоростей. При y = 0 увеличение толщины сдвигового слоя приводит к уменьшению пика продольной силы, а при y = p – к увеличению.

Рис. 6. Безразмерные коэффициенты волновых сил A₁, A₂ в зависимости от скорости движения нагрузки при $\overline{Q}$ = 1, y = 0 (а), ψ = p/2 (б) и ψ = π (в): 1 – S = 0; 2, 3 – S = 0.3, 0.5 для двух слоев: H₁ = 5 м, H₂ = 45 м; 4, 5 – S = 0.3, 0.5 для одного слоя: H = 50 м.

Исследован вклад дополнительных корней в волновые силы, который характеризуется отношениями

$R_1=F_1^0\text{}/\text{}{F}_1,R_2=F_2^0\text{}/\text{}{F}_2,$

где $F_1^0,F_2^0$, определяются интегралами по интервалу $({\sigma }_3,{\sigma }_4)$ в формулах (3.11), (3.12). На рис. 7 показаны зависимости отношений R₁, R₂ от скорости движения нагрузки V при y = p/2, –Q = 0 (а) и –Q = 1 (б), S = 0.3, 0.5. Как видно из графиков, вклад дополнительных корней в волновые силы мал, вклад их в боковую силу значительно выше, чем в продольную. Сжатие ледяного покрова приводит к увеличению интервала $({\sigma }_3,{\sigma }_4),$ и вклад дополнительных корней в волновые силы увеличивается.

Рис. 7. Вклад дополнительных корней в волновые силы в зависимости от скорости движения нагрузки при H₁ = 5 м, ψ = π/2, $\overline{Q}$ = 0 (а), $\overline{Q}$ = 1 (б): 1–2 – S = 0.3, 0.5.

ЗАКЛЮЧЕНИЕ

Исследовано поведение ледяного покрова при равномерном прямолинейном движении нагрузки при наличии в верхнем слое сдвигового потока с учетом однородного сжатия льда. Поле скоростей жидкости не потенциально. Показано, что сдвиг скорости, толщина сдвигового слоя, направление движения нагрузки и коэффициент сжатия льда оказывают существенное влияние на критические скорости изгибно-гравитационных волн и прогиб ледяного покрова. При увеличении коэффициента сжатия льда критические скорости значительно уменьшаются. Сжатие льда приводит к неустойчивости возмущений при большой толщине сдвигового слоя.

При увеличении коэффициента сжатия амплитуды прогиба льда значительно увеличиваются за нагрузкой и уменьшаются перед нагрузкой. Увеличение параметра сдвига S оказывает качественно такое же влияние, как увеличение скорости движения нагрузки. При движении нагрузки под произвольным углом к направлению течения симметрия прогиба нарушается. Боковая сила максимальна при движении нагрузки перпендикулярно течению. При увеличении коэффициента сжатия льда волновые силы значительно растут.

Sobre autores

L. Tkacheva

Autor responsável pela correspondência
Email: tkacheva@hydro.nsc.ru
Rússia, Новосибирск

Bibliografia

Arquivos suplementares

Arquivos suplementares

Ação

1. JATS XML

Baixar

2. Fig. 1. An example of eigenvectors of the dispersion relation for the system of equations (2.3), (2.4).

Baixar (62KB)

Metadados

3. Fig. 2. Wave number curves for an ice cover of 50 m thickness, 1–3 –

Baixar (105KB)

Metadados

4. Fig. 3 Dependence of critical velocities on the thickness of the shear layer at (a, b) and (c, d), y = 0 (a, c) and y = p (b, d): curves 1–5 correspond to the values

Baixar (263KB)

Metadados

5. Fig. 4. Ice cover deflection at –Q = 0 (a, b) and –Q = 1 (c, d), y = 0 (a, c) and y = p (b, d): 1 – S = 0; 2, 3 – S = 0.3, 0.5 for two layers H1 = 5 m, H2 = 45 m; 4, 5 – S = 0.3, 0.5 for one layer H1 = H = 50 m.

Baixar (498KB)

Metadados

6. Fig. 5. Dimensionless wave force coefficients A1, A2 depending on the load movement speed at –Q = 0, y = 0 (a), y = p/2 (b) and y = p (c): 1 – S = 0; 2, 3 – S = 0.3, 0.5 for two layers: H1 = 5 m, H2 = 45 m; 4, 5 – S = 0.3, 0.5 for one layer: H = 50 m.

Baixar (178KB)

Metadados

7. Fig. 6. Dimensionless wave force coefficients A1, A2 depending on the load movement speed at = 1, y = 0 (a), ψ = p/2 (b) and ψ = π (c): 1 – S = 0; 2, 3 – S = 0.3, 0.5 for two layers: H1 = 5 m, H2 = 45 m; 4, 5 – S = 0.3, 0.5 for one layer: H = 50 m.

Baixar (175KB)

Metadados

8. Fig. 7. Contribution of additional roots to wave forces depending on the speed of load movement at H1 = 5 m, ψ = π/2,  = 0 (a),  = 1 (b): 1–2 – S = 0.3, 0.5.

Baixar (141KB)

Metadados