在文献^[7]中,Faltinsen从物理分析的角度给出了静态情况下悬垂电缆线上各点的纵坐标$x$与 $y$, $T$的解析关系.在文献^[7]中假设条件下,本节的定理3.1利用数学的方法严格推导出了 与文献^[7]等价的$x$,$y$,$T$与弧长参数$s$的解析关系.并将定理3.1的结果延伸到更一般的情形 (定理 3.2),提出了相应的解的存在条件.同时,通过找出静态 解的解析式,我们可以利用(4.2)--(4.8)式计算相同条件下的数值解,通过数值解与解析式的比较, 一定程度上可以体现数值格式(4.2)--(4.8)式的高精度性.

现在我们来考虑(2.14)式的一种特殊情况下的解的存在性.

注 3.1 很明显(3.1)式是(2.14)式在 $$\frac{\partial x}{\partial t}\equiv 0,~~~\frac{\partial y}{\partial t}\equiv 0 ,$$ $\forall s \in [x_*,s_*],t\in [0,+\infty)$下的简化模型.

首先,考虑电缆线全局$(0\leq s \leq b(t))$光滑的情况,此时相应的边界条件如下

我们得到如下结论.

定理3.1在边界条件(3.2)下,方程组 (3.1 )存在唯一全局光滑解

证由(3.1)式的第一个方程关于$s$从$x_*$到$s$积分以及 (3.2)式,可得

由(3.1)式的第二个方程关于$s$从$x_*$到$s$积分以及(3.2)式,可得

由(3.1)式的第三个方程以及(3.4)--(3.5)式,可得

将(3.6)式代入(3.4)式,可得

将(3.6)式代入(3.5)式,可得

在施工过程中,最理想的状态是让电缆全局光滑,因为只有全局光滑的电缆线才能使海底接触点 $(a(t),0)$处的张力为$0$.然而,由(3.2)式可以看出,电缆线全局光滑的 条件是很苛刻的,而且电缆线施工过程中在船速太大的情况下, 电缆绳大多数不会保持全局光滑姿态,而是在$(a(t),0)$处产生尖点. 因此,研究电缆线局部光滑情况下 解的存在条件有着很重要的意义.

接下来,我们考虑电缆线局部光滑$(a(t)\leq s \leq b(t))$的情况.此时相应的边界条件如下

引理3.1假设$K(x_0-x_*)\geq\ln 2$,以及 $\frac{1}{2}{\rm e}^{K(x_0-x_*)}-\frac{1}{2}{\rm e}^{-K(x_0-x_*)}<\sqrt{K^2h^2+2Kh}$,则函数 $$\varphi(u)=\frac{1}{2}(\sqrt{1-u^2}+1){\rm e}^{\frac{K(x_0-x_*)}{u}} -\frac{u^2}{2(\sqrt{1-u^2}+1){\rm e}^{\frac{K(x_0-x_*)}{u}}}-\sqrt{(Kh+1)^2-u^2}$$ 在$(0,1)$上存在唯一的零点.

证 令$m=K(x_0-x_*)\geq \ln 2$,当$u\in (0,1)$时,通过简单计算可得 $$\begin{eqnarray*}\label{E:1.1} \varphi'(u)&=&-\frac{1}{2}\bigg[{\rm e}^{\frac{m}{u}}+\frac{u^2}{(\sqrt{1-u^2}+1)^2 {\rm e}^{\frac{m}{u}}}\bigg]\bigg[\frac{u}{\sqrt{1-u^2}}+\frac{(\sqrt{1-u^2}+1)m}{u^2}\bigg] \\ &&-\frac{u}{(\sqrt{1-u^2}+1){\rm e}^{\frac{m}{u}}} +\frac{u}{\sqrt{(Kh+1)^2-u^2}} \nonumber \\ &\leq & -\frac{1}{2}{\rm e}^{\frac{m}{u}}\frac{u}{\sqrt{1-u^2}}+\frac{u}{\sqrt{(Kh+1)^2-u^2}} \nonumber \\ &\leq &-\frac{1}{2}({\rm e}^{\frac{m}{u}}-2)\frac{u}{\sqrt{1-u^2}} \\ & <& -\frac{1}{2}({\rm e}^{m}-2)\frac{u}{\sqrt{1-u^2}} \\& <& -\frac{1}{2}({\rm e}^{ln2}-2)\frac{u}{\sqrt{1-u^2}} \\ & =& 0 . \end{eqnarray*}$$

因此,$\varphi(u)$在$(0,1)$上单调递减.

此外还有 $$\lim \limits_{u\rightarrow 0}\varphi(u)=+\infty, \lim\limits_{u\rightarrow 1}\varphi(u)=\frac{1}{2}{\rm e}^{K(x_0-x_*)} -\frac{1}{2}{\rm e}^{-K(x_0-x_*)}-\sqrt{K^2h^2+2Kh}< 0.$$

由介值定理,可得$\varphi(u)$在$(0,1)$上存在唯一的零点.

定理 在引理3.1的条件以及边界条件(3.11)下,方程组(3.1)存在唯一局部光滑解

证由(3.1)式第一个方程以关于$s$从$x_*$到$s$及(3.11)式,可得

由(3.1)式第二个方程关于$s$从$x_*$到$s$以及(3.11)式,可得

由(3.1)式第三个方程以及(3.11)式,可得

将(3.15)式代入(3.13)式,可得

将(3.15)式代入(3.14)式,可得

由(3.17)式可得 $$x_0=x(s_*)=x_*+\frac{\frac{{\rm d}x(x_*)}{{\rm d}s}}{K}\ln\frac{\frac{{\rm d}y(x_*)}{{\rm d}s}+K(s_*-x_*)+\sqrt{[\frac{{\rm d}y(x_*)}{{\rm d}s}+K(s_*-x_*)]^2+[\frac{{\rm d}x(x_*)}{{\rm d}s}]^2}}{\frac{{\rm d}y(x_*)}{{\rm d}s}+1},$$ 则有

由(3.19)式,令 $y(s_{*})=h$,则有

令 $u=\frac{{\rm d}x(x_*)}{{\rm d}s},$则根据(3.1)式的第一个方程可得$\sqrt{1-u^2}=\frac{{\rm d}y(x_*)}{{\rm d}s}$. 再由(3.20)--(3.21)式可得 $$\varphi(u)=\frac{1}{2}(\sqrt{1-u^2}+1){\rm e}^{\frac{K(x_0-x_*)}{u}}-\frac{u^2}{2(\sqrt{1-u^2}+1){\rm e}^{\frac{K(x_0-x_*)}{u}}}-\sqrt{(Kh+1)^2-u^2}=0.$$

由注2.1以及引理3.1可知$\frac{{\rm d}x(x_*)}{{\rm d}s}$是$\varphi(u)$的唯一零点, $\frac{{\rm d}y(x_*)}{{\rm d}s}=\sqrt{1-(\frac{{\rm d}x(x_*)}{{\rm d}s} )^2}$.

4 数值模拟与实际应用 跳转到：本文顶部　1 引言 　2 数学模型 　3 静态解的存在性 　4 数值模拟与实际应用 　参考文献跳转到：本文顶部　1 引言 　2 数学模型 　3 静态解的存在性 　4 数值模拟与实际应用 　参考文献

自由边界问题(2.14)--(2.21)相当复杂,以至于目前还无法直接表达它的解. 然而,在电缆通信以及石油探测方面的应用要求我们至少能找到它的近似解. 因此,本节采用数值计算的方法^[13]来近似求解自由边界问题(2.14)--(2.21).

首先,我们将弧长参数$s$ 以及时间参数$t$ 离散化如下

$$\begin{equation}\label{eqn 4.1} 0=s_0<s_1<s_2<\cdots<s_{I_0+N},0=t_0<t_1<\cdots <t_n<\cdots <t_N; \end{equation}$$

(4.1)

在节点 $(s_i,t_{n+1})(n=1,2,\cdots,N-1;i=1,2,\cdots,I_0+n)$处,将方程组(2.14)离散化^[8]如下

$$\begin{eqnarray}\label{eqn 4.2} F_i(\mathbf{U})& =&T_{i+\frac{1}{2}}^{n+1}\alpha_{i+1}(x_{i+1}^{n+1}-x_{i}^{n+1})-T_{i-\frac{1}{2}}^{n+1}\beta_{i+1}(x_{i}^{n+1}-x_{i-1}^{n+1}) \nonumber\\ & & -K_1\bigg(\frac{x_{i}^{n+1}-x_{i}^{n}}{\Delta t_{n+1}}\bigg)^2 \frac{y_{i+1}^{n+1}-y_{i}^{n+1}}{\Delta s_{i+1}}-K_2 \bigg(\frac{y_{i}^{n+1}-y_{i}^{n}}{\Delta t_{n+1}}\bigg)^2\frac{y_{i+1}^{n+1}-y_{i}^{n+1}}{\Delta s_{i+1}} \nonumber\\ && -\rho[\gamma_{n+1}(x_{i}^{n+1}-x_{i}^{n})-\theta_{n+1}(x_{i}^{n}-x_{i}^{n-1})] \nonumber\\ &=&0, \end{eqnarray}$$

(4.2)

$$\begin{eqnarray}\label{eqn 4.3} G_i(\mathbf{U})& =&T_{i+\frac{1}{2}}^{n+1}\alpha_{i+1}(y_{i+1}^{n+1}-y_{i}^{n+1})-T_{i-\frac{1}{2}}^{n+1}\beta_{i+1}(y_{i}^{n+1}-y_{i-1}^{n+1}) \nonumber\\ && +K_1\bigg(\frac{y_{i}^{n+1}-y_{i}^{n}}{\Delta t_{n+1}}\bigg)^2 \frac{x_{i+1}^{n+1}-x_{i}^{n+1}}{\Delta s_{i+1}}+K_2 \bigg(\frac{x_{i}^{n+1}-x_{i}^{n}}{\Delta t_{n+1}}\bigg)^2\frac{x_{i+1}^{n+1}-x_{i}^{n+1}}{\Delta s_{i+1}} \nonumber\\ & & -(\rho-\rho_0)g-\rho[\gamma_{n+1}(y_{i}^{n+1}-y_{i}^{n})-\theta_{n+1}(y_{i}^{n}-y_{i}^{n-1})] \nonumber\\ &=&0, \end{eqnarray}$$

(4.3)

$$\begin{equation}\label{eqn 4.4} \Psi_i({\mathbf U})=(x_{i+1}^{n+1}-x_{i}^{n+1})^2+(y_{i+1}^{n+1}-y_{i}^{n+1})^2-(\Delta s_{i+1})^2=0 , \end{equation}$$

(4.4)

其中, $$\alpha_{i}=\frac {2}{\Delta s_{i+1}(\Delta s_{i+1}+\Delta s_{i})},~~\beta_{i}=\frac {2}{\Delta s_{i}(\Delta s_{i+1}+\Delta s_{i})},$$ $$\gamma_{n+1}=\frac {2}{\Delta t_{n+1}(\Delta t_{n+1}+\Delta t_{n})},~~\theta_{n+1}=\frac {2}{\Delta t_{n}(\Delta t_{n+1}+\Delta t_{n})},$$ $${\mathbf U}=(x_1^{n+1},x_2^{n+1},\cdots,x_{I_0+n}^{n+1},y_1^{n+1},y_2^{n+1},\cdots,y_{I_0+n}^{n+1},T_{\frac{3}{2}}^{n+1},T_{\frac{5}{2}}^{n+1},\cdots,T_{I_0+n+{\frac{1}{2}}}^{n+1})^{T} .$$

接下来,我们在初始节点$(s_i,t_0)$处,将(2.15)--(2.18)式离散化如下

$$\begin{equation}\label{eqn 4.5} x_{i}^{0}=0 ,y_{i}^{0}=s_{i},x_{i}^{1}-x_{i}^{0}=0,y_{i}^{1}-y_{i}^{0}=0, T_{i+\frac{1}{2}}^{0}=(\rho-\rho_0)gs_{i+\frac{1}{2}} , \ i=0,1,\cdots,I_{0}. \end{equation}$$

(4.5)

设 $r$表示已经沉到海底的节点总数,则在边界节点处将(2.19)--(2.21)式离散化如下

$$\begin{equation}\label{eqn 4.6} y_{r}^{n+1}=h ,T_{r+\frac{1}{2}}^{n+1}=f_{0} ,x_{I_0+n}^{n+1}=xup(t_{n+1}) , y_{I_0+n}^{n+1}=0 ,T_{I_0+n+\frac{1}{2}}^{n+1}=f(t_{n+1}) , \ n=0,1,\cdots,N-1. \end{equation}$$

(4.6)

注 4.1 对于任意时刻 $t$,怎样获取$r$的值?首先, 当$t=0$时,令 $r=0$,并设定任意小的常量$\varepsilon >0$.假设在$t=t_{n}$时刻有 $r=r_{*}$.则在$t=t_{n+1}$时刻,我们比较 $y_i^{n+1}(i=r+1,r+2,\cdots,I_0+n)$ 与 $\varepsilon$ 的大小,并令 $l=\overline{\overline{\{i_{*}|y_{i_*}^{n+1}<\varepsilon, i_{*}=r+1,r+2,\cdots,I_0+n\}}}.$ 在$t=t_{n+1}$时刻 $r$的值为$r=r_{*}+l$. 其中,对任意集合$A$,$\overline{\overline{A}}$ 表示$A$中元素的个数.

由(4.2)--(4.4)式可知,在$t_{n+1}$时刻,当 $i=I_0+n$时, $x_{I_0+n}^{n-1},y_{I_0+n}^{n-1},z_{I_0+n}^{n-1}$ 是未知的(因为$n-1$时刻,第$I_0+n$个节点还未沉放入海底). 因此,我们需要给这几个变量取定一个合适的近似值,根据几何关系有

$$\begin{equation}\label{eqn 4.7} x_{I_0+n}^{n-1}=x_{I_0+n-1}^{n-1}+\frac {\Delta s_{I_0+n}}{\Delta s_{I_0+n-1}}(x_{I_0+n-1}^{n-1}-x_{I_0+n-2}^{n-1}), \end{equation}$$

(4.7)

$$\begin{equation}\label{eqn 4.8} y_{I_0+n}^{n-1}=y_{I_0+n-1}^{n-1}+\frac {\Delta s_{I_0+n}}{\Delta s_{I_0+n-1}}(y_{I_0+n-1}^{n-1}-y_{I_0+n-2}^{n-1}). \end{equation}$$

(4.8)

这种近似处理以后,$x_{I_0+n}^{n-1},y_{I_0+n}^{n-1}$就可以参与计算了.

由(2.17)和(2.18)式可知,在每个时刻 $t$,方程组(2.14)式的弧长参数取值范围为 $a(t) \leq s \leq b(t)$,其中$a(t)$ 随着沉入海底的部分越来越多而逐渐增大, $b(t)$随着探测船行走过程中释放缆绳也逐渐增大,这就涉及到一个自由边界计算的问题^[2]. 自由边界是个很复杂的问题.在这篇论文中,上边界$b(t)$由(2.18)式定义,下 边界$a(t)$是未知的.怎样取$a(t)$的值呢?本文采用如下方法: 为了让已经坠入海底的节点退出计算 ,我们在(4.2)--(4.4)式中将 $x_{i}^{n+1},y_{i}^{n+1},T_{i+\frac{1}{2}}^{n+1}~(0 \leq i\leq r)$作 为边界值,并令 $x_{i}^{n+1}=x_{i}^{n},y_{i}^{n+1}=y_{i}^{n}=0,$ $T_{i+\frac{1}{2}}^{n+1}=f_0~(0 \leq i\leq r).$ 同时,为了让刚刚从船上 放入海水中的节点参与计算,我们每一个时刻在缆绳顶端增加一个节点,即 增加$x_{I_0+n}^{n+1},y_{I_0+n}^{n+1},T_{I_0+n+\frac{1}{2}}^{n+1}$ 三个变量.除此以外,我们还用前面注释4.1的方法来取$r$的值,故有,$a(t_{n+1})=x_{r}^{n+1}$.

差分格式的离散化已经完成,下面我们采用迭代算法来求解非线性方程组(4.2)--(4.4)式,具体过程如下.

在$t=t_{n+1}$时刻,取如下初始迭代值: $x_{i}^{(0)}=x_{i}^{n},$ $y_{i}^{(0)}=y_{i}^{n},$ $T_{i+\frac{1}{2}}^{(0)} =T_{i+\frac{1}{2}}^{n}$,$i=r+1,\cdots,$ $I_0+n.$

(1) $x$迭代计算

$$\begin{eqnarray}\label{eqn 4.9} &&T_{i+\frac{1}{2}}^{(s)}\alpha_{i+1}x_{i+1}^{(s+1)}-(T_{i+\frac{1}{2}}^{(s)}\alpha_{i+1}+T_{i-\frac{1}{2}}^{(s)}\beta_{i+1}+\rho\gamma_{n+1})x_{i}^{(s+1)}+T_{i-\frac{1}{2}}^{(s)}\beta_{i+1}x_{i-1}^{(s+1)} \nonumber\\ &=&K_1\bigg(\frac{x_{i}^{(s)}-x_{i}^{n}}{\Delta t_{n+1}}\bigg)^2\frac{y_{i+1}^{(s)} -y_{i}^{(s)}}{\Delta s_{i+1}}+K_2 \bigg(\frac{y_{i}^{(s)}-y_{i}^{n}}{\Delta t_{n+1}}\bigg)^2\frac{y_{i+1}^{(s)}-y_{i}^{(s)}}{\Delta s_{i+1}} \nonumber\\ &&-\rho\gamma_{n+1}x_{i}^{n}-\rho\theta_{n+1}(x_{i}^{n}-x_{i}^{n-1}). \end{eqnarray}$$

(4.9)

$x_{i}^{(s)},y_{i}^{(s)},T_{i+\frac{1}{2}}^{(s)}(i=r+1,\cdots,I_0+n)$ 的值已经在前面的迭代中算得,那么(4.9)式就是关于 $x_{i}^{(s+1)}$ $(i=r+1,\cdots,I_0+n)$ 的线性方程组.求解这样一个代数方程组, 可以得到$x_{i}^{(s+1)}$ $(i=r+1,\cdots,I_0+n)$的值.

(2) $y$迭代计算

$$\begin{equation}\label{eqn 4.10} y_{i+1}^{(s+1)}-y_{i}^{(s+1)}=\sqrt{|(\Delta s_{i+1})^2-(x_{i+1}^{(s+1)}-x_{i}^{(s+1)})^2|}. \end{equation}$$

(4.10)

$x_{i}^{(s+1)}$的值已由 (4.9)式获得,则(4.10) 式就成为关于$y_{i}^{(s+1)}$ $(i=r+1,\cdots,I_0+n)$的线性方程组.这样就可以计算出$y_{i}^{(s+1)}$ $(i=r+1,\cdots,I_0+n)$ .

(3) $T$迭代计算

$$\begin{eqnarray}\label{eqn 4.11} &&\alpha_{i+1}(y_{i+1}^{(s+1)}-y_{i}^{(s+1)})T_{i+\frac{1}{2}}^{(s+1)}-\beta_{i+1}(y_{i}^{(s+1)}-y_{i-1}^{(s+1)})T_{i-\frac{1}{2}}^{(s+1)} \nonumber\\ &=&-K_1\bigg(\frac{y_{i}^{(s+1)}-y_{i}^{n}}{\Delta t_{n+1}}\bigg)^2\frac{x_{i+1}^{(s+1)} -x_{i}^{(s+1)}}{\Delta s_{i+1}}-K_2 \bigg(\frac{x_{i}^{(s+1)}-x_{i}^{n}}{\Delta t_{n+1}}\bigg)^2\frac{x_{i+1}^{(s+1)}-x_{i}^{(s+1)}}{\Delta s_{i+1}} \nonumber\\ &&+(\rho-\rho_0)g+\rho[\gamma_{n+1}(y_{i}^{(s+1)}-y_{i}^{n})-\theta_{n+1}(y_{i}^{n}-y_{i}^{n-1})]. \end{eqnarray}$$

(4.11)

$x_{i}^{(s+1)},y_{i}^{(s+1)}$ 已由(4.9)--(4.10)式获得,则 (4.11)式是关于 $T_{i+\frac{1}{2}}^{(s+1)}$ $(i=r+1,\cdots,I_0+n)$ 的线性方程组. 因此,可以计算出 $T_{i+\frac{1}{2}}^{(s+1)}$,$i=r+1,\cdots,I_0+n$.

现在 $x_{i}^{(s+1)},y_{i}^{(s+1)},T_{i+\frac{1}{2}}^{(s+1)}$的值已获得, 重复 (4.9)--(4.11)式 依此取$s=0,1,2,\cdots ,M$, 经过有限次迭代,得到$t=t_{n+1}$时刻的数值解如下 $$x_{i}^{n+1}=x_{i}^{(M)},y_{i}^{n+1}=y_{i}^{(M)}, T_{i+\frac{1}{2}}^{n+1}=T_{i+\frac{1}{2}}^{(M)},~ i=r+1,r+2,\cdots,I_0+n.$$

由定理3.1和定理3.2,我们知道,在解的存在性条件下,(3.3)和 (3.12)式是 (3.1)式 加上某些特定边界条件的解.现在为了检验我们提出的算法的有效性, 我们用此格式来数值模拟静态方程 (3.1)的解,并与其解析解(3.3),(3.12)进行对比, 通过下面的两组数值实验可以看出格式(4.2)--(4.8) 有很高的精度.

表 1

表 1 模拟参数 参数 变量(试验1) 变量(试验2) 参数 变量(应用) 参数 变量(应用) $x_0$ -- 5m $C_d$ 0.008 $g$ 9.8m/s$^2$ $h$ 10m 10m $C_L$ 0.008 $f_0$ 10N $g$ 9.8m/s$^2$ $9.8 m/s^2$ $d$ 0.025m $\sigma$ 1.8m/s $\rho_c$ 1500kg/m$^3$ 1500kg/m$^3$ $\rho_c$ 2000kg/m$^3$ $a$ 0.05m/s$^2$ $\rho_w$ 1000kg/m$^3$ 1000kg/m$^3$ $\rho_w$ 1000kg/m$^3$ $h$ 20m

表 1 模拟参数

试验1 考虑方程组(3.1)在边界条件(3.2)下的解.各物理参数取值如表 1所示.模拟结果见图 4.

图 4

图 4 方程组(3.1)在边界条件(3.2)下的数值解与解析解的对比

试验2 考虑方程组(3.1)在边界条件(3.11)下的解. 各物理参数取值如表 1所示.模拟结果见图 5.

图 5

图 5 方程组(3.1)在边界条件(3.11)下的数值解与解析解的对比

图 6

图 6 二维电缆短时间实际应用

二维电缆短时间实际应用模拟. 由图 4--图 5可以看出,我们提出的数值算法有很高的精度. 下面我们利用该数值算法来动态模拟二维电缆铺设过程中电缆线的运动轨迹 (定解问题(2.14)--(2.21)). 各参数取值如表 1 所示.模拟结果如图 6. 其中,彩色线表示电缆线.