一类交叉非线性反应扩散方程组的数值模拟

Numerical Simulation of a Class of Cross Nonlinear Reaction-Diffusion Equations

摘要: 该文研究一类交叉非线性反应扩散方程组，而交叉扩散项在不稳定时可以创造图灵模式。为了寻找一种简单有效的非线性交叉反应扩散系统的数值方法，本文提出了采用这种新的空间谱插值配点方法模拟了一些数值算例，其结果和理论上的吻合度较好，结果表明了该方法的有效性。

Abstract: In this paper, we study a class of cross nonlinear reaction-diffusion equations whose cross diffusion terms can create Turing patterns when they are unstable. In order to find a simple and effective numerical method for nonlinear cross-reaction-diffusion systems, this paper presents a new spatial spectral interpolation method to simulate some numerical examples. The results agree well with the theory, and the results show that the method is effective.

1. 引言

交叉反应扩散问题在生物学 [1] [2]、医学 [3] [4] [5]、化学 [6] 等中都有广泛的应用，对于求解这类问题的数值方法也有很多，有限元分析法 [7]、无网格方法 [8] 等。在本文中，我们主要研究在理论生物学领域出现的模型，其模型如下：

$\{\begin{array}{l}\frac{\partial u}{\partial t}=d_{11}\Delta u+d_{12}\Delta v+f\left(u,v\right)\\ \frac{\partial v}{\partial t}=d_{21}\Delta u+d_{22}\Delta v+g\left(u,v\right)\end{array}$ (1)

其中： $u\left(x,y,t\right)$ 和 $v\left(x,y,t\right)$ 是未知函数， $d_{11},d_{12},d_{21},d_{22}$ 是扩散系数， $\Delta ={\partial }^2/\partial x^2+{\partial }^2/\partial y^2$ 是拉普拉斯算子， $\left(x,y\right)\in \Omega$， $t>0$ 是在光滑的边界 $\partial \Omega$ 上的齐次Neumann边界条件，即 $\partial u/\partial n=\partial v/\partial n=0$， $f\left(u,v\right),g\left(u,v\right)$ 是已知的光滑函数。

本文考虑了交叉反应扩散系统的以下两组非线性情况：

1) $f\left(u,v\right)=u\left(\alpha -\beta u\right)-\frac{cuv}{mv+1},g\left(u,v\right)=\frac{suv}{mv+1}-rv$

2) $f\left(u,v\right)=u\left(\sigma \left(1-\frac{u}{\beta }\right)-\left(1+\alpha v\right)v\right),g\left(u,v\right)=v\left(\left(1+\alpha v\right)u-1\right)$

2. 空间谱插值配点方法的描述

在本文中，我们使用空间谱插值配点方法来解决模型(1)。

从参考文献 [9] [10] 中，我们可以得到离散序列 $u_1,u_2,\cdots ,u_N$ 的插值函数 $I_{N}u\left(x\right)$ 可以写成：

$u\left(x\right)~I_{N}u\left(x\right)={\displaystyle \underset{m=1}{\overset{N}{\sum }}{u}_m{S}_N}\left(x-x_m\right)$ (2)

其中 $S_N\left(x\right)=\frac{\sin \left(\pi x/h\right)}{\left(2\pi /h\right)\tan \left(x/2\right)}$， $I_N$ 对于任何函数，都是这样的插值算子， $u\left(x\right)$ 在区间 $\left[0,2\pi \right]$ 定义为 $u_j=u\left(jh\right),j=1,\cdots ,N,x_j-x_m=\left(j-m\right)h$，插值空间为 $span\left\{S_N\left(x-jh\right),j=1,2,\cdots ,N\right\}$。由此推导出n阶导数的简单表达式并不难 $I_{N}u\left(x\right)$ 在 $x_j=jh$ 处的n阶导数的简单表达式并不难，为：

$I_N{u}^{\left(n\right)}\left(x_j\right)={\displaystyle \underset{m=1}{\overset{N}{\sum }}{u}_m}{S}_N^{\left(n\right)}\left(x_j-x_m\right)$ (3)

其中 $D_N^{\left(n\right)}={\left[S_N^{\left(n\right)}\left(x_j-x_m\right)\right]}_{i,j=1,2,\cdots ,N}$ 被称为第n个谱微分矩阵 [9]。

这里我们考虑有限空间域 $\Omega =\left[0,2\pi \right]\times \left[0,2\pi \right]$。我们定义 $N^2$ 个等间距的网格点在区域 $\Omega$ 上，因此

$\left(x_i,y_j\right)=\left(ih,jh\right),i,j=1,2,3,\cdots ,N.$

其中 $h=\frac{2\pi }{N}$ 对于给定的有限自然数 $N\in \mathbb{N}$。使用(2)则配点函数 $I_{N}u\left(x,y,t\right)$ 和 $I_{N}v\left(x,y,t\right)$ 关于函数 $u\left(x,y,t\right)$ 和 $v\left(x,y,t\right)$ 可以被写为：

$\begin{array}{l}u\left(x,y,t\right)~I_{N}u\left(x,y,t\right)={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N\left(x-x_i\right)S_N\left(y-y_i\right)u\left(x_i,y_j,t\right)}}\\ v\left(x,y,t\right)~I_{N}v\left(x,y,t\right)={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N\left(x-x_i\right)S_N\left(y-y_i\right)v\left(x_i,y_j,t\right)}}\end{array}$ (4)

其中 $u_{i,j}=u\left(x_i,y_j,t\right),v_{i,j}=v\left(x_i,y_j,t\right),i,j=1,2,\cdots ,N$。因此，下列关系在搭配点 $\left(x_p,y_q\right)$ 处成立：

$\begin{array}{l}u\left(x_p,y_q,t\right)~I_{N}u\left(x_p,y_q,t\right)={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N\left(x_p-x_i\right)S_N\left(y_q-y_i\right)u\left(x_i,y_j,t\right)}}\\ v\left(x_p,y_q,t\right)~I_{N}v\left(x_p,y_q,t\right)={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N\left(x_p-x_i\right)S_N\left(y_q-y_i\right)v\left(x_i,y_j,t\right)}}\end{array}$ (5)

$\begin{array}{l}{u}^{\left(2,0\right)}\left(x_p,y_q,t\right)~I_N{u}^{\left(2,0\right)}\left(x_p,y_q,t\right)=\frac{{\partial }^{2}u\left(x_p,y_q,t\right)}{\partial x^2}={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N^{\left(2\right)}}}\left(x_p-x_i\right)S_N\left(y_q-y_i\right)u\left(x_i,y_j,t\right)\\ u^{\left(0,2\right)}\left(x_p,y_q,t\right)~I_N{u}^{\left(0,2\right)}\left(x_p,y_q,t\right)=\frac{{\partial }^{2}u\left(x_p,y_q,t\right)}{\partial y^2}={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N}}\left(x_p-x_i\right)S_N^{\left(2\right)}\left(y_q-y_i\right)u\left(x_i,y_j,t\right)\\ v^{\left(2,0\right)}\left(x_p,y_q,t\right)~I_N{v}^{\left(2,0\right)}\left(x_p,y_q,t\right)=\frac{{\partial }^{2}v\left(x_p,y_q,t\right)}{\partial x^2}={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N^{\left(2\right)}}}\left(x_p-x_i\right)S_N\left(y_q-y_i\right)v\left(x_i,y_j,t\right)\\ v^{\left(0,2\right)}\left(x_p,y_q,t\right)~I_N{v}^{\left(0,2\right)}\left(x_p,y_q,t\right)=\frac{{\partial }^{2}v\left(x_p,y_q,t\right)}{\partial y^2}={\displaystyle \underset{i=1}{\overset{N}{\sum }}{\displaystyle \underset{j=1}{\overset{N}{\sum }}{S}_N}}\left(x_p-x_i\right)S_N^{\left(2\right)}\left(y_q-y_i\right)v\left(x_i,y_j,t\right)\end{array}$

注意的是

$\begin{array}{l}u={\left[u_{11},u_{21},\cdots ,u_{N1},u_{12},u_{22},\cdots ,u_{N2},u_{1N},\cdots ,u_{NN}\right]}^{\text{T}},\\ v={\left[v_{11},v_{21},\cdots ,v_{N1},v_{12},v_{22},\cdots ,v_{N2},v_{1N},\cdots ,v_{NN}\right]}^{\text{T}},\end{array}$ (6)

因此，公式(5)可以被写成一下矩阵形式：

$u^{\left(2,0\right)}=D_N^{\left(2,0\right)}u,u^{\left(0,2\right)}=D_N^{\left(0,2\right)}u,v^{\left(2,0\right)}=D_N^{\left(2,0\right)}v,v^{\left(0,2\right)}=D_N^{\left(0,2\right)}v.$ (7)

$D_N^{\left(2,0\right)}u=D_N^{\left(2\right)}\otimes E_N,D_N^{\left(0,2\right)}=E_N\otimes D_N^{\left(2\right)},D_N^{\left(0,0\right)}=E_N\otimes E_N,$

其中 $E_N$ 是N阶单位矩阵， $\otimes$ 是矩阵的克罗内克积，相反地。二阶谱微分矩阵为：

$D_N^{\left(2\right)}=\left(\begin{array}{ccccc}-\frac{{\pi }^2}{3h^2}-\frac{1}{6}& \ddots & ⋮& \ddots & \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)\\ \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)& \ddots & -\frac{1}{2}{\csc }^2\left(\frac{2h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{2h}{2}\right)\\ \frac{1}{2}{\csc }^2\left(\frac{2h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{3h}{2}\right)\\ ⋮& \ddots & -\frac{{\pi }^2}{3h^2}-\frac{1}{6}& \ddots & ⋮\\ \frac{1}{2}{\csc }^2\left(\frac{3h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{2h}{2}\right)\\ \frac{1}{2}{\csc }^2\left(\frac{2h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)& \ddots & \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)\\ \frac{1}{2}{\csc }^2\left(\frac{1h}{2}\right)& \ddots & ⋮& \ddots & -\frac{{\pi }^2}{3h^2}-\frac{1}{6}\end{array}\right)$ (8)

结合等式(6)和等式(7)，Equation (1)可以写成以下系统：

$\frac{\partial }{\partial t}\left[\begin{array}{c}u\\ v\end{array}\right]=\left[\begin{array}{cc}{d}_{11}D& d_{12}D\\ d_{21}D& d_{22}D\end{array}\right]\left[\begin{array}{c}u\\ v\end{array}\right]=\left[\begin{array}{c}f\left(u,v\right)\\ g\left(u,v\right)\end{array}\right]$ (9)

这里

$\left[u,v\right]=\left[u_{11},\cdots ,u_{N1},u_{12},\cdots ,u_{N2},u_{1N},\cdots ,u_{NN},v_{11},\cdots ,v_{N1},v_{12},\cdots ,v_{N2},v_{1N},\cdots ,v_{NN}\right],$

$D=D_N^{\left(2,0\right)}+D_N^{\left(0,2\right)}=D_N^{\left(2\right)}\otimes E_N+E_N\otimes D_N^{\left(2\right)},$

$\left[f\left(u,v\right),g\left(u,v\right)\right]=\left[f\left(u_{11},v_{11}\right),\cdots ,f\left(u_{NN},v_{NN}\right),g\left(u_{11},v_{11}\right),\cdots ,g\left(u_{NN},v_{NN}\right)\right]$

利用Matlab中的ode45求解器求解系统(9)，得到系统(1)的数值解。

3. 系统分叉分析

在这一部分，我们给出了交叉反应扩散问题的图灵分叉分析。我们假设非扩散系统(1)的平衡点是 $E_{*}=\left(u_0,v_0\right)$，所以

$\{\begin{array}{l}f\left(u_0,v_0\right)=0\\ g\left(u_0,v_0\right)=0\end{array}$ (10)

现在我们对平衡点 $E_{*}$ 进行线性稳定性分析。我们在E处对系统的雅可比矩阵 $A_0$ 进行

$A_0=\left[\begin{array}{cc}\frac{\partial f}{\partial u}\left(u_0,v_0\right)& \frac{\partial f}{\partial v}\left(u_0,v_0\right)\\ \frac{\partial g}{\partial u}\left(u_0,v_0\right)& \frac{\partial g}{\partial v}\left(u_0,v_0\right)\end{array}\right]$ (11)

$A_0$ 的特征根由

${\lambda }_{1,2}\left(0\right)=\frac{1}{2}\left[tr\left(A_0\right)\pm \sqrt{{\left[tr\left(A_0\right)\right]}^2-4\det \left(A_0\right)}\right]$

其中 $tr\left(A_0\right)=\frac{\partial f}{\partial u}\left(u_0,v_0\right)+\frac{\partial g}{\partial v}\left(u_0,v_0\right)$ 和 $\det \left(A_0\right)=\left[\frac{\partial f}{\partial u}\left(u_0,v_0\right)\right]\times \left[\frac{\partial g}{\partial v}\left(u_0,v_0\right)\right]-\left[\frac{\partial g}{\partial u}\left(u_0,v_0\right)\right]\times \left[\frac{\partial f}{\partial v}\left(u_0,v_0\right)\right]$。

一般情况下，如果特征值的实部 ${\lambda }_{1,2}\left(0\right)$ 为负值，则非扩散系统(1)是稳定的。

接下来，我们导出了交叉扩散驱动的不稳定性条件，并证明了这些是在没有交叉扩散的情况下经典扩散驱动不稳定性条件的推广。

$A_k=\left[\begin{array}{cc}{b}_{11}& b_{12}\\ b_{21}& b_{22}\end{array}\right]-k^2\left[\begin{array}{cc}{d}_{11}& d_{12}\\ d_{21}& d_{22}\end{array}\right]$ (12)

其中 $b_{11}=\frac{\partial f}{\partial u}\left(u_0,v_0\right),b_{12}=\frac{\partial f}{\partial v}\left(u_0,v_0\right),b_{21}=\frac{\partial g}{\partial u}\left(u_0,v_0\right),b_{22}=\frac{\partial g}{\partial v}\left(u_0,v_0\right)$。

特征根由 $A_k$ 给出：

${\lambda }_{1,2}\left(k\right)=\frac{1}{2}\left[tr\left(A_k\right)\pm \sqrt{{\left[tr\left(A_k\right)\right]}^2-4\det \left(A_k\right)}\right]$

其中 $tr\left(A_k\right)=b_{11}+b_{22}-k^2\left(d_{11}+d_{22}\right)$ 和 $\det \left(A_k\right)=\left(d_{11}{d}_{22}-d_{12}{d}_{21}\right)k^4-\left(b_{11}{d}_{22}+b_{22}{d}_{11}-b_{12}{d}_{21}-b_{21}{d}_{12}\right)k^2+\left(b_{11}{b}_{22}-b_{12}{b}_{21}\right)$。

图灵分岔发生在没有非扩散时稳定的等分态，但在交叉扩散时变得不稳定。因此， $E_{*}$ 成为交叉扩散系统(1)不稳定点的唯一途径是特征值 ${\lambda }_{1,2}\left(k\right)$ 的实部是正的，那么交叉扩散系统(1)是不稳定的。

4. 数值实验

在这一部分中，我们给出了一些数值插图，以更好地解释上述分析结果使用不同的初始条件和参数。如果非线性项为 $f\left(u,v\right)=u\left(\alpha -\beta u\right)-\frac{cuv}{mv+1},g\left(u,v\right)=\frac{suv}{mv+1}-rv$ 时，则反应扩散体系(1)则变为例1 [10] 形式。其中 [10] $E_{*}=\left(u_0,v_0\right)$， $u_0=\frac{s\left(m\alpha -c\right)+\sqrt{s^2{\left(m\alpha -c\right)}^2+4cmrs\beta }}{2ms\beta }$，

$v_0=\frac{s{u}_0-r}{mr}$. $A_0=\left[\begin{array}{cc}-\beta u_0& -\frac{c{u}_0}{{\left(m{v}_0+1\right)}^2}\\ \frac{s{v}_0}{\left(m{v}_0+1\right)}& -\frac{ms{u}_0{v}_0}{{\left(m{v}_0+1\right)}^2}\end{array}\right]$.

例1 考虑以下形式的交叉反应扩散模型 [10]：

$\{\begin{array}{l}\frac{\partial u}{\partial t}-d_{11}\Delta u-d_{12}\Delta v=u\left(\alpha -\beta u\right)-\frac{cuv}{mv+1},\left(x,t\right)\in \Omega \times \left[0,\infty \right]\\ \frac{\partial v}{\partial t}-d_{21}\Delta u-d_{22}\Delta v=\frac{suv}{mv+1}-rv,\left(x,t\right)\in \Omega \times \left[0,\infty \right]\\ \frac{\partial u}{\partial n}=\frac{\partial v}{\partial n}=0,\left(x,t\right)\in \partial \Omega \times \left[0,\infty \right]\\ u\left(x,0\right)=u_0\left(x\right)\ge 0,v\left(x,0\right)=v_0\left(x\right)\ge 0,x\in \Omega \end{array}$ (13)

这里 $d_{11},d_{12},d_{21},d_{22}>0$， $\alpha ,\beta ,\gamma ,s,m,c$ 为常数， $t\in R_{+}$ 为时间变量， $x,y\in R$ 是空间变量， $u,v$ 是t与 $x,y$ 的函数。其中参数1 $\beta =1,c=1,\gamma =1,s=1,m=4,\alpha =0.8,d_{11}=1,d_{22}=10,d_{12}=1,d_{21}=10$ 和参数2 $\beta =1,c=1,\gamma =2,s=1,m=4,\alpha =1.2,d_{11}=1,d_{22}=1,d_{12}=0.4,d_{21}=2$，固定一个初始条件 $v=0.5ones\left(N\right)$，当u以及参数变化时模拟的结果见图1~11，具体初始条件见表1。

如果非线性项为 $f\left(u,v\right)=u\left(\sigma \left(1-\frac{u}{\beta }\right)-\left(1+\alpha v\right)v\right),g\left(u,v\right)=v\left(\left(1+\alpha v\right)-1\right)$ 时，则反应扩散体系(1)则变成例2形式。其中 [11] $u_0=\beta \left(1-\frac{v_0\left(1+\alpha v_0\right)}{\sigma }\right)$，而 $v_0$ 是下列三次多项式方程 $\beta {\alpha }^2{v}^3+2\beta \alpha v^2+\beta \left(1-\alpha \sigma \right)v+\sigma \left(1-\beta \right)=0$ 的根。 $A_0=\left[\begin{array}{cc}{a}_{11}& a_{12}\\ a_{21}& a_{22}\end{array}\right]$， $a_{11}=\sigma -\frac{2\sigma u_0}{\beta }-v_0-\alpha v_0^2$， $a_{12}=-u_0-2\alpha u_0{v}_0$， $a_{21}=v_0+\alpha v_0^2$， $a_{22}=u_0+2\alpha u_0{v}_0-1$。

例2 考虑以下形式的交叉反应扩散模型[11]：

$\{\begin{array}{l}\frac{\partial u}{\partial t}-d_{11}\Delta u-d_{12}\Delta v=u\left(\sigma \left(1-\frac{u}{\beta }\right)-\left(1+\alpha v\right)v\right),\left(x,t\right)\in \Omega \times \left[0,\infty \right]\\ \frac{\partial v}{\partial t}-d_{21}\Delta u-d_{22}\Delta v=v\left(\left(1+\alpha v\right)u-1\right),\left(x,t\right)\in \Omega \times \left[0,\infty \right]\\ \frac{\partial u}{\partial n}=\frac{\partial v}{\partial n}=0,\left(x,t\right)\in \partial \Omega \times \left[0,\infty \right]\\ u\left(x,0\right)=u_0\left(x\right)\ge 0,v\left(x,0\right)=v_0\left(x\right)\ge 0,x\in \Omega \end{array}$ (14)

这里 $d_{11},d_{12},d_{21},d_{22}>0$， $\alpha ,\beta ,\sigma$ 为常数， $t\in R_{+}$ 为时间变量， $x,y\in R$ 是空间变量， $u,v$ 是t与 $x,y$ 的函数。其中 $\delta =4,\alpha =0.3,\beta =1.5,d_{11}=0.5,d_{22}=10$，参数1为 $d_{12}=0.2,d_{21}=10.5$ ；参数2为 $d_{12}=1,d_{21}=0.2$，固定一个初始条件 $v=0.5ones\left(N\right)$，当v固定时改变参数和u得到的模拟结果见图12~19，固定 $u=\cos \left(x^2+y^2\right)$ 时得到的模拟结果见图20、图21，具体的初始条件见表2。

Table 1. The different initial conditions corresponding to the numerical solution and pattern in example 1 are shown in Figures 1~11

表1. 例1中数值解和斑图所对应的不同的初始条件在图1~11

Figure 1. Shows the numerical solution of $u=ones\left(N\right)+\text{sech}\left(x^2/0.02-9+2y/0.01\right)$ in parameter 1 and initial condition of example 1

图1. 在参数1和初始条件为 $u=ones\left(N\right)+\text{sech}\left(x^2/0.02-9+2y/0.01\right)$ 的数值解

Figure 2. Shows the pattern of $u=ones\left(N\right)+\text{sech}\left(x^2/0.02-9+2y/0.01\right)$ in parameter 1 and initial condition of example 1

图2. 例1在参数1和初始条件为 $u=ones\left(N\right)+\text{sech}\left(x^2/0.02-9+2y/0.01\right)$ 的斑图

Figure 3. Shows the pattern of $u=ones\left(N\right)+\text{sech}\left(x^2/0.02-9+2y/0.01\right)$ in parameter 2 and initial condition of example 1

图3. 例1在参数2和初始条件为 $u=ones\left(N\right)+\text{sech}\left(x^2/0.02-9+2y/0.01\right)$ 的斑图

Figure 4. Shows the numerical solution of $u=\sin \left(x^2/0.02-9+2y/0.01\right)$ in parameter 1 and initial condition of example 1

图4. 例1在参数1和初始条件为 $u=\sin \left(x^2/0.02-9+2y/0.01\right)$ 的数值解

Figure 5. Shows the pattern of $u=\sin \left(x^2/0.02-9+2y/0.01\right)$ in parameter 1 and initial condition of example 1

图5. 例1在参数1和初始条件为 $u=\sin \left(x^2/0.02-9+2y/0.01\right)$ 的斑图

Figure 6. Shows the numerical solution of $u=\sin \left(\exp \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left(\cos \left({\left(x+0.4\right)}^3+{\left(y-0.4\right)}^2\right)\right)$ in parameter 1 and initial condition of example 1

图6. 例1在参数1和初始条件为 $u=\sin \left(\exp \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left(\cos \left({\left(x+0.4\right)}^3+{\left(y-0.4\right)}^2\right)\right)$ 的数值解

Figure 7. Shows the pattern of $u=\sin \left(\exp \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left(\cos \left({\left(x+0.4\right)}^3+{\left(y-0.4\right)}^2\right)\right)$ in parameter 1 and initial condition of example 1

图7. 例1在参数1和初始条件为 $u=\sin \left(\exp \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left(\cos \left({\left(x+0.4\right)}^3+{\left(y-0.4\right)}^2\right)\right)$ 的斑图

Figure 8. Shows the numerical solution of $u=\text{sech}\left({\left(-x^3-y^5\right)}^3\right)$ in parameter 1 and initial condition of example 1

图8. 例1在参数1和初始条件为 $u=\text{sech}\left({\left(-x^3-y^5\right)}^3\right)$ 的数值解

Figure 9. Shows the pattern of $u=\text{sech}\left({\left(-x^3-y^5\right)}^3\right)$ in parameter 1 and initial condition of example 1

图9. 例1在参数1和初始条件为 $u=\text{sech}\left({\left(-x^3-y^5\right)}^3\right)$ 的斑图

Figure 10. Shows the numerical solution of $u=\cos \left(\sin \left(-x^3+y^5\right)\right)$ in parameter 1 and initial condition of example 1

图10. 例1在参数1和初始条件为 $u=\cos \left(\sin \left(-x^3+y^5\right)\right)$ 的数值解

Figure 11. Shows the pattern of $u=\cos \left(\sin \left(-x^3+y^5\right)\right)$ in parameter 1 and initial condition of example 1

图11. 图为例1在参数1和初始条件为 $u=\cos \left(\sin \left(-x^3+y^5\right)\right)$ 的斑图

Figure 12. Shows the numerical solution of $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)-\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ in parameter 1 and initial condition of example 2

图12. 为例2在参数1和初始条件为 $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)-\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ 的数值解

Figure 13. Shows the pattern of $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)-\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ in parameter 1 and initial condition of example 2

图13. 为例2在参数1和初始条件为 $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)-\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ 的斑图

Figure 14. Shows the numerical solution of $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)+\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ in parameter 2 and initial condition of example 2

图14. 图为例2在参数2和初始条件为 $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)+\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ 的数值解

Figure 15. Shows the pattern of $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)+\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ in parameter 2 and initial condition of example 2

图15. 为例2在参数2和初始条件为 $u=\sin \left(\sin \left({\left(x-0.4\right)}^3{\left(y+0.4\right)}^2\right)\right)+\cos \left(\exp \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)\right)$ 的斑图

Figure 16. Shows the numerical solution of $u=\text{sech}\left(\pi \left(-x^2+y^2\right)\right)$ in parameter 1 and initial condition of example 2

图16. 图为例2在参数1和初始条件为 $u=\text{sech}\left(\pi \left(-x^2+y^2\right)\right)$ 的数值解

Figure 17. Shows the pattern of $u=\text{sech}\left(\pi \left(-x^2+y^2\right)\right)$ in parameter 1 and initial condition of example 2

图17. 为例2在参数1和初始条件为 $u=\text{sech}\left(\pi \left(-x^2+y^2\right)\right)$ 的斑图

Figure 18. Shows the numerical solution with the initial condition of $u=\cos \left(\exp \left(\sin \left(-x^2-y^2\right)\right)\right)$ example 2

图18. 图为例2在参数1和初始条件为 $u=\cos \left(\exp \left(\sin \left(-x^2-y^2\right)\right)\right)$ 的数值解

Figure 19. Shows the pattern of $u=\cos \left(\exp \left(\sin \left(-x^2-y^2\right)\right)\right)$ in parameter 1 and initial condition of example 2

图19. 为例2在参数1和初始条件为 $u=\cos \left(\exp \left(\sin \left(-x^2-y^2\right)\right)\right)$ 的斑图

Figure 20. Shows the numerical solution of $u=\cos \left(x^2+y^2\right),v=\sin \left(x^2/0.02-9+2y/0.01\right)$ in parameter 2 and initial condition of example 2

图20. 图为例2在参数2和初始条件为 $u=\cos \left(x^2+y^2\right),v=\sin \left(x^2/0.02-9+2y/0.01\right)$ 的数值解

Figure 21. Shows the pattern of $u=\cos \left(x^2+y^2\right),v=\sin \left(x^2/0.02-9+2y/0.01\right)$ in parameter 2 and initial condition of example 2

图21. 为例2在参数2和初始条件为 $u=\cos \left(x^2+y^2\right),v=\sin \left(x^2/0.02-9+2y/0.01\right)$ 的斑图

Table 2. The different initial conditions corresponding to the numerical solution and pattern in example 2 are shown in Figures 12~21

表2. 例2中数值解和斑图所对应的不同的初始条件在图12~21

5. 结论

本文用谱插值配点方法求解一类非线性交叉反应扩散模型，数值结果表明了该方法与理论吻合较好。本文所有程序由matlab2017b算得。

致谢

感谢王玉兰老师的支持与帮助。

参考文献