基于软组织的一种四阶微分方程零解稳定性

Zero-Solution Stability of a Fourth-Order Differential Equation Based on Soft Tissue

摘要: 本文主要讨论由软组织模型得到的一种四阶非线性微分方程零解的稳定性，通过使用能量度量算法构造该方程的Lyapunov函数，进而给出方程零解的稳定性的充分条件。

Abstract: In this paper, we mainly discuss the stability of the zero solution of a fourth-order nonlinear differ-ential equation obtained from the soft tissue model. By using the energy metric algorithm, we con-struct the Lyapunov function of the equation, and then give the sufficient conditions for the stability of the zero solution.

1. 研究背景

软组织通常是指人体的皮肤、皮下组织、肌肉、肌腱、韧带等。软组织弹性模型的研究可以追溯到上世纪80年代开始的形变模型的研究。本文将软组织看作是质量–弹簧模型，从而得到一个微分方程组。

李雅普诺夫稳定性理论是俄国数学家A. M. Lyapunov在1892年所创立的。该理论主要用于分析系统稳定性。李雅普诺夫稳定性理论在许多领域中都有广泛的应用，如生化反应、病毒动力学、生态生物、传染病模型研究等。从19世纪末以来，李雅普诺夫稳定性理论一直被用于指导稳定性的研究和应用。1967年，D.布肖对表征在集合与映射水平上的系统建立了李雅普诺夫第二方法。此时，Lyapunov函数被推广为向量形式，称为向量Lyapunov函数。李雅普诺夫稳定性理论主要指李雅普诺夫直接法。李雅普诺夫直接法可用于任意阶数的系统，运用这一方法可以不必求出微分方程的解而直接判定其解的稳定性。运用李雅普诺夫第二方法的关键是构造出一个合适的Lyapunov函数 [1] [2] ，一直以来，许多学者提出了各种不同的构造方法，如能量函数法、梯度法、能量度量算法、类比法 [3] [4] 等。

王联、王慕秋用类比法研究了一类三阶非线性系统的Lyapunov函数的构造问题 [5] 。秦宏立等学者应用类比法讨论了一类非线性系统的Lyapunov函数的构造 [6] 。柳苗应用能量度量算法和类比法讨论了两种四阶非线性微分方程的Lyapunov函数的构造方法 [7] 。

软组织可以看作是质量–弹簧模型，即单弹簧摆模型，根据牛顿第二定律，其运动方程如下：

$ma=K\left(L-\Vert r\Vert \right)\frac{r}{\Vert r\Vert }-K_d\stackrel{\dot{}}{r}+mg$

根据上述运动方程，可以得到其数学模型：

令 $\stackrel{\dot{}}{x}=v$ ， $\stackrel{\dot{}}{y}=z$ ， $\frac{K}{m}=\beta$ ， $K_d=\gamma$ ，则上述模型可简化为：

$\{\begin{array}{l}\stackrel{\dot{}}{x}=v\\ \stackrel{\dot{}}{y}=z\\ \stackrel{\dot{}}{v}=\beta \left(\frac{L}{\sqrt{x^2+y^2}}-1\right)x-\gamma y\\ \stackrel{\dot{}}{z}=\beta \left(\frac{L}{\sqrt{x^2+y^2}}-1\right)y-\gamma +g\end{array}$ (1.1)

运用能量度量算法的 [8] [9] 过程如下：

先将方程变换成等价的微分方程组，

${\stackrel{\dot{}}{x}}_i=F_i\left(x\right)$ (1.2)

其中 $x=col\left(x_1,x_2,\cdots ,x_n\right),i=1,\cdots ,n$ 。

从而有

$\frac{\text{d}{x}_i}{\text{d}{x}_j}=\frac{F_i\left(x\right)}{F_j\left(x\right)},i<j;i,j=1,\cdots ,n.$ (1.3)

于是对上述等式变换成下面形式，

$F_j\left(x\right)\text{d}{x}_i-F_i\left(x\right)\text{d}{x}_j=0,i<j;i,j=1,\cdots ,n.$ (1.4)

对(1.4)式做适当的变换得，

$w_1\left(x\right)\text{d}{x}_1+w_2\left(x\right)\text{d}{x}_2+\cdots +w_n\left(x\right)\text{d}{x}_n=0$

其中

$w_1\left(x\right)=w_1\left(x_1,0,\cdots ,0\right),w_2\left(x\right)=w_2\left(x_1,x_2,\cdots ,0\right),\cdots ,w_n\left(x\right)=w_n\left(x_1,x_2,\cdots ,x_n\right)$

假设

$V\left(x\right)={\displaystyle {\int }_0^{x_1}{w}_1\left({\tau }_1,0,\cdots ,0\right)\text{d}{\tau }_1}+{\displaystyle {\int }_0^{x_2}{w}_2\left(x_1,{\tau }_2,\cdots ,0\right)\text{d}{\tau }_2}+\cdots +{\displaystyle {\int }_0^{x_n}{w}_n\left(x_1,x_2,\cdots ,{\tau }_n\right)\text{d}{\tau }_n}$ (1.5)

求出 $V\left(x\right)$ 沿着方程(1.2)的全导数

$\frac{\text{d}V\left(x\right)}{\text{d}t}=\frac{\partial V\left(x\right)}{\partial x_1}{\stackrel{\dot{}}{x}}_1+\frac{\partial V\left(x\right)}{\partial x_2}{\stackrel{\dot{}}{x}}_2+\cdots +\frac{\partial V\left(x\right)}{\partial x_n}{\stackrel{\dot{}}{x}}_n$

根据常微分方程稳定性的相关理论以及 $V\left(x\right)$ 和 $\frac{\text{d}V\left(x\right)}{\text{d}t}$ 的相关性质，可以讨论出方程(1.2)的零解是否

具有稳定性 [10] [11] [12] [13] 。

定理1：设系统的状态方程为

$\stackrel{\dot{}}{x}=f(\; x\; )$

其中 $x_e=0$ ，为其平衡态.若存在一个有连续一阶偏导数的正定函数 $V\left(x\right)$ ，满足下述条件： $\frac{\text{d}V\left(x\right)}{\text{d}t}$ 为非正定(半负定)的，则该系统在原点处的平衡态是一致稳定的 [14] [15] 。

2. 主要结果

本文先通过使用能量度量算法构造出合适的李雅普诺夫函数 $V\left(x\right)$ ，再给出方程的零解的稳定性的充分条件，以及充分条件的证明过程。

考虑方程

$\{\begin{array}{l}m\ddot{x}=K\left(L-\sqrt{x^2+y^2}\right)\frac{x}{\sqrt{x^2+y^2}}-K_d\stackrel{\dot{}}{x}\\ m\ddot{y}=K\left(L-\sqrt{x^2+y^2}\right)\frac{y}{\sqrt{x^2+y^2}}-K_d\stackrel{\dot{}}{y}+mg\end{array}$ (2.1)

的等价微分方程组：

$\{\begin{array}{l}{\stackrel{\dot{}}{x}}_1=x_3\\ {\stackrel{\dot{}}{x}}_2=x_4\\ {\stackrel{\dot{}}{x}}_3=f\left(x_1,x_2\right)x_1-\gamma x_3\\ {\stackrel{\dot{}}{x}}_4=f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g\end{array}$ (2.2)

其中 $f\left(x_1,x_2\right)=\beta \left(\frac{L}{\sqrt{x_1^2+{\left(x_2+c\right)}^2}}-1\right)$ ， $c=\frac{L\beta +g}{\beta }$ 。

等价微分方程组按照(1.3)式两两相除，得到：

$\frac{\text{d}{x}_1}{\text{d}{x}_2}=\frac{x_3}{x_4}$

$\frac{\text{d}{x}_1}{\text{d}{x}_3}=\frac{x_3}{f\left(x_1,x_2\right)x_1-\gamma x_3}$

$\frac{\text{d}{x}_1}{\text{d}{x}_4}=\frac{x_3}{f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g}$

$\frac{\text{d}{x}_2}{\text{d}{x}_3}=\frac{x_4}{f\left(x_1,x_2\right)x_1-\gamma x_3}$

$\frac{\text{d}{x}_2}{\text{d}{x}_4}=\frac{x_4}{f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g}$

$\frac{\text{d}{x}_3}{\text{d}{x}_4}=\frac{f\left(x_1,x_2\right)x_1-\gamma x_3}{f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g}$

即

$x_3\text{d}{x}_2=x_4\text{d}{x}_1$ (2.3)

$x_3\text{d}{x}_3=\left(f\left(x_1,x_2\right)x_1-\gamma x_3\right)\text{d}{x}_1$ (2.4)

$x_3\text{d}{x}_4=\left[f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g\right]\text{d}{x}_1$ (2.5)

$x_4\text{d}{x}_3=\left(f\left(x_1,x_2\right)x_1-\gamma x_3\right)\text{d}{x}_2$ (2.6)

$x_4\text{d}{x}_4=\left[f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g\right]\text{d}{x}_2$ (2.7)

$\left(f\left(x_1,x_2\right)x_1-\gamma x_3\right)\text{d}{x}_4=\left[f\left(x_1,x_2\right)\left(x_2+c\right)-\gamma x_4+g\right]\text{d}{x}_3$ (2.8)

做如下适当变换：

$\left(2.3\right)\times cf\left(x_1,x_2\right)x_3{x}_4+\left(2.4\right)\times x_3{x}_4^2-\left(2.6\right)\times x_4{x}_3^2-\left(2.5\right)\times x_3{x}_4^2+\left(2.7\right)\times x_4{x}_3^2+\left(2.8\right)\times x_1{x}_4{x}_3^2$

可以得到，

$\begin{array}{l}\left(\gamma x_3^2{x}_4^2-\gamma x_3{x}_4^3+g{x}_3{x}_4^2\right)\text{d}{x}_1+\left(\gamma x_3^2{x}_4^2-\gamma x_3^2{x}_4-g{x}_4{x}_3^2\right)\text{d}{x}_2\\ +\left(\gamma x_1{x}_3^2{x}_4^2-g{x}_1{x}_4{x}_3^2-f\left(x_1,x_2\right)x_1^2{x}_3^2{x}_4-cf\left(x_1,x_2\right)x_1{x}_4{x}_3^2\right)\text{d}{x}_3\\ +\left(f\left(x_1,x_2\right)x_1^2{x}_3^2{x}_4-\gamma x_1{x}_4{x}_3^3\right)\text{d}{x}_4=0\end{array}$ (2.9)

由于

${\displaystyle {\int }_0^{x_1}{w}_1\left({\tau }_1,0,0,0\right)\text{d}{\tau }_1}=0$

${\displaystyle {\int }_0^{x_2}{w}_2\left(x_1,{\tau }_2,0,0\right)\text{d}{\tau }_2}=0$

${\displaystyle {\int }_0^{x_3}{w}_3\left(x_1,x_2,{\tau }_3,0\right)\text{d}{\tau }_3}=0$

$\begin{array}{c}{\displaystyle {\int }_0^{x_4}{w}_4\left(x_1,x_2,x_3,{\tau }_4\right)\text{d}{\tau }_4}={\displaystyle {\int }_0^{x_4}\left(f\left(x_1,x_2\right)x_1^2{x}_3^2{\tau }_4-\gamma x_1{\tau }_4{x}_3^3\right)\text{d}{\tau }_4}\\ =\frac{1}{2}\left(f\left(x_1,x_2\right)x_1^2{x}_3^2-\gamma x_1{x}_3^3\right)x_4^2\end{array}$

对(2.9)式进行积分，

$\begin{array}{c}V\left(x_1,x_2,x_3,x_4\right)={\displaystyle {\int }_0^{x_1}{w}_1\left({\tau }_1,0,0,0\right)\text{d}{\tau }_1}+{\displaystyle {\int }_0^{x_2}{w}_2\left(x_1,{\tau }_2,0,0\right)\text{d}{\tau }_2}\\ +{\displaystyle {\int }_0^{x_3}{w}_3\left(x_1,x_2,{\tau }_3,0\right)\text{d}{\tau }_3}+{\displaystyle {\int }_0^{x_4}{w}_4\left(x_1,x_2,x_3,{\tau }_4\right)\text{d}{\tau }_4}\\ =\frac{1}{2}\left(f\left(x_1,x_2\right)x_1^2{x}_3^2-\gamma x_1{x}_3^3\right)x_4^2\end{array}$

函数 $V\left(x_1,x_2,x_3,x_4\right)$ 沿着(2.2)式的全导数为：

$\begin{array}{c}\frac{\text{d}V\left(x\right)}{\text{d}t}=\frac{\partial V\left(x\right)}{\partial x_1}{\stackrel{\dot{}}{x}}_1+\frac{\partial V\left(x\right)}{\partial x_2}{\stackrel{\dot{}}{x}}_2+\frac{\partial V\left(x\right)}{\partial x_3}{\stackrel{\dot{}}{x}}_3+\frac{\partial V\left(x\right)}{\partial x_4}{\stackrel{\dot{}}{x}}_4\\ =\left[2f\left(x_1,x_2\right)x_1{x}_3^3+\left(\frac{\partial f\left(x_1,x_2\right)}{\partial x_1}{x}_3^2-\gamma f\left(x_1,x_2\right)\right)x_3{x}_1^2+f{\left(x_1,x_2\right)}^2{x}_3{x}_1^3\right]x_4^2\\ +\gamma x_1{x}_4\left(\gamma x_4-g-cf\left(x_1,x_2\right)-x_{1}f\left(x_1,x_2\right)\right)x_3^2\\ +\left[\left(\frac{1}{2}{x}_4\frac{\partial f\left(x_1,x_2\right)}{\partial x_2}-2\gamma f\left(x_1,x_2\right)\right)x_4^2{x}_3^2+\left(f\left(x_1,x_2\right)x_1+c+g\right)f\left(x_1,x_2\right)x_3^2{x}_4\right]x_1^2\end{array}$

定理2.1 若系统(2.2)满足以下条件：

(1) 函数 $f\left(x_1,x_2\right)$ 连续可微。

(2) $f\left(x_1,x_2\right)x_1>\gamma x_3$ 。

(3) $\gamma x_4-g>\left(c+x_1\right)f\left(x_1,x_2\right)$ 。

(4) $x_1<0,x_2\le 0,x_3>0,x_4>0$ 。

(5) $\gamma >\max \left\{\frac{x_3^2}{\min \left\{f\left(x_1,x_2\right)\right\}}\frac{\partial f\left(x_1,x_2\right)}{\partial x_1},\frac{x_4}{4\min \left\{f\left(x_1,x_2\right)\right\}}\frac{\partial f\left(x_1,x_2\right)}{\partial x_2}\right\}$ 。

则系统(2.2)的零解是渐进稳定的。

证明：根据实际情况，可知系统(2.2)满足 $\gamma >0$ ，并且

$\sqrt{x_1^2+{\left(x_2+c\right)}^2}>\frac{L\beta +g}{\beta }$

所以有

$\frac{L}{\sqrt{x_1^2+{\left(x_2+c\right)}^2}}<\frac{L\beta }{L\beta +g}<1$

因此

$f\left(x_1,x_2\right)<0.$

又由于 $f\left(x_1,x_2\right)x_1>\gamma x_3$ ，则

$f\left(x_1,x_2\right)x_1^2-\gamma x_1{x}_3>0$

$\frac{1}{2}\left(f\left(x_1,x_2\right)x_1^2-\gamma x_1{x}_3\right)x_4^2{x}_3^2>0$

因此

$V\left(x_1,x_2,x_3,x_4\right)>0$

即 $V\left(x_1,x_2,x_3,x_4\right)$ 为正定函数。

又由于 $\gamma x_4-g>\left(c+x_1\right)f\left(x_1,x_2\right)$ ，则

$\gamma x_4-g-cf\left(x_1,x_2\right)-x_{1}f\left(x_1,x_2\right)>0$

$\gamma x_1{x}_4\left(\gamma x_4-g-cf\left(x_1,x_2\right)-x_{1}f\left(x_1,x_2\right)\right)<0$

$\left(f\left(x_1,x_2\right)x_1+c+g\right)f\left(x_1,x_2\right)x_3^2{x}_4<0$

由于 $$ ，所以有

$\left(\frac{1}{2}{x}_4\frac{\partial f\left(x_1,x_2\right)}{\partial x_2}-2\gamma f\left(x_1,x_2\right)\right)x_4^2{x}_3^2<0$

$\left(\frac{\partial f\left(x_1,x_2\right)}{\partial x_2}{x}_3^2-\gamma f\left(x_1,x_2\right)\right)x_1^2{x}_3<0$

显然， $f\left(x_1,x_2\right)<-\frac{2x_3^2}{x_1^2}$ 成立

所以有

$2f\left(x_1,x_2\right)x_1{x}_3^3+f{\left(x_1,x_2\right)}^2{x}_1^3{x}_3<0$

因此

$\stackrel{\dot{}}{V}\left(x\right)<0$

根据定理1.1，系统(2.2)的零解是渐进稳定的。

进而可知模型(1.1)在弹性系数K，阻尼系数K_d存在的情况下，在受到外力作用时，运动是稳定的。

本文创新性地使用能量度量算法构造出了质量–弹簧模型的Lypunov函数，证明了模型在有弹性系数和阻尼系数存在时，其运动是渐进稳定的。后期将对模型进行改进，改为正六边形弹簧网络模型，使其更加接近现实中的软组织结构。

参考文献

NOTES

^*通讯作者。