Brusselator模型的空间谱插值配点方法数值模拟

Numerical Simulation of the Brusselator Model with Spatial Spectral Interpolation Coordination Method

摘要: 该文研究一类非线性反应扩散方程组Brusselator模型，因此寻找一种简单有效的非线性反应扩散系统的数值方法是非常重要的。基于这一问题，本文提出了采用这种新的空间谱插值配点方法模拟了一些数值算例，其结果和理论上的吻合较好，结果表明了该方法的有效性。

Abstract: This paper studies the Brusselator model of a class of nonlinear reaction-diffusion equations, so it is very important to find a simple and effective numerical method for nonlinear reaction-diffusion systems. Based on this problem, some numerical examples are simulated by using this new spatial spectral interpolation collocation method.

1. 引言

Brusselator模型是一个著名的模型 [1] [2]，它通常是用来描述化学反应过程中化学元素变化的一类反应扩散方程组 [3]。其形式为：

$\{\begin{array}{l}\frac{\partial u}{\partial t}=\Delta \left(u_{11}u+u_{12}v\right)-\left(\beta +1\right)u+u^{2}v+\alpha ,\left(x,t\right)\in \Omega \times \left[0,T\right]\\ \frac{\partial v}{\partial t}=\Delta \left(u_{21}u+u_{22}v\right)+\beta u-u^{2}v,\left(x,t\right)\in \Omega \times \left[0,T\right]\\ \frac{\partial u}{\partial n}=\frac{\partial v}{\partial n}=0,\left(x,t\right)\in \partial \Omega \times \left[0,T\right]\\ u\left(x,0\right)={\psi }_1\left(x\right),v\left(x,0\right)={\psi }_2\left(x\right),x\in \Omega \end{array}$ (1)

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

许多方法被用来解决Brusselator模型。这些方法包括有限差分法 [4]，Galerkin有限元法 [5]，最佳均匀(BURA)有理近似 [6]。本文基于一种新的空间谱插值配点方法，研究了Brusselator系统的动力学行为，分析了该系统的稳定性 [7] [8]。数值模拟结果和理论吻合较好，仿真结果表明了该方法的有效性。

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

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

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

$u\left(x\right)\sim 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)\sim 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)\sim 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)\sim 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)\sim 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)\sim 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)\sim 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)\sim 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)\sim 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)可以被写成一下矩阵形式：

$\begin{array}{l}{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.\end{array}$ (7)

$\begin{array}{l}{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,\end{array}$

其中 $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}{u}_{11}D& u_{12}D\\ u_{21}D& u_{22}D\end{array}\right]\left[\begin{array}{c}u\\ v\end{array}\right]=\left[\begin{array}{c}{f}_1\left(u,v\right)\\ f_2\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],$

$\begin{array}{l}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)},\end{array}$

$\begin{array}{l}\left[f_1\left(u,v\right),f_2\left(u,v\right)\right]\\ =\left[f_1\left(u_{11},v_{11}\right),\cdots ,f_1\left(u_{NN},v_{NN}\right),f_2\left(u_{11},v_{11}\right),\cdots ,f_2\left(u_{NN},v_{NN}\right),f_3\left(u_{11},v_{11}\right),\cdots ,f_3\left(u_{NN},v_{NN}\right)\right]\end{array}$

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

3. 动力系统分析

在这一部分，我们给出了Brusselator模型的图灵分叉分析 [11]。如果不考虑扩散现象，即假设扩散系数为零时，则方程(1)就可以写为：

$\{\begin{array}{l}\frac{\text{d}u}{\text{d}t}=\alpha +u^{2}v-\left(\beta +1\right)u\\ \frac{\text{d}v}{\text{d}t}=\beta u-u^{2}v\end{array}$ (10)

从动力学角度我们来分析系统(10)的解的性质 [12]。首先(10)只有唯一的平衡解 $\left(\alpha ,\beta /\alpha \right)$。其次，讨论这一平衡解的稳定性。为此考虑(10)的Jacobian矩阵：

$DF\left(u,v\right)=\left[\begin{array}{cc}2uv-\left(\beta +1\right)& u^2\\ \beta -2uv& -u^2\end{array}\right],$

其在平衡解 $\left(\alpha ,\beta /\alpha \right)$ 的Jacobian矩阵为：

$DF\left(\alpha ,\beta /\alpha \right)=\left[\begin{array}{cc}\beta -1& {\alpha }^2\\ -\beta & -{\alpha }^2\end{array}\right]$ (11)

显见矩阵(11)的迹 $\tau =tracDF\left(\alpha ,\beta /\alpha \right)=-{\alpha }^2+\beta -1$，判别式 $\Delta =-{\alpha }^2\left(\beta -1\right)+{\alpha }^2\beta ={\alpha }^2$。由于 $\Delta ={\alpha }^2>0$，所以(11)的特征根同号，故：

1) 若 $\tau <0$，即 $\beta >\alpha +1$，平衡解 $\left(\alpha ,\beta /\alpha \right)$ 是稳定的；

2) 若 $\tau >0$，即 $\beta <\alpha +1$，平衡解 $\left(\alpha ,\beta /\alpha \right)$ 是发散的。

下面讨论(10)的分岔。固定 $\alpha$，让 $\beta$ 变动。由 $\tau =0$ 得到 ${\beta }_H={\alpha }^2+1$。下面我们来证明 ${\beta }_H$ 是(10)的Hopf分岔点 [13] [14] [15]，设 $\lambda$ 为(11)的特征根，为此必须证明：(a) (11)在 ${\beta }_H$ 点只有一对纯虚根。(b) $\partial \mathrm{Re}\left(\lambda \right)/\partial \beta \ne 0$。事实上，(11)的特征根为 ${\lambda }_{1,2}=\left(\tau \pm \sqrt{{\tau }^2-4\Delta }\right)/2$，由 $\tau =0$ 及 $\Delta ={\alpha }^2$ 知 ${\lambda }_{1,2}=\pm \alpha i$。所以(a)成立。其次，由 ${\lambda }_{1,2}=\left(\tau \pm \sqrt{{\tau }^2-4\Delta }\right)/2$ 可得 $\mathrm{Re}\left(\lambda \right)=\left(-{\alpha }^2+\beta -1\right)/2$，所以 $\partial \mathrm{Re}\left(\lambda \right)/\partial B=1/2\ne 0$，故(b)也成立。综上所述， ${\beta }_H$ 为系统(10)的Hopf分岔点，且由于 $\tau <0$ 时系统的平衡不稳定，所以由 $\tau <0$ 变化到 $\tau =0$ 系统在Hopf分岔点处得到稳定性 [14]。

4. 数值实验

在这一节中，我们给出了一些数值例子，以便更好地解释使用不同的初始条件和参数的上述分析结果。

例1 考虑以下形式的Brusselator模型 [3]：

$\{\begin{array}{l}\frac{\partial u}{\partial t}=\Delta u_{1}u-\left(\beta +1\right)u+u^{2}v+\alpha ,\left(x,t\right)\in \Omega \times \left[0,T\right]\\ \frac{\partial v}{\partial t}=\Delta u_{2}v+\beta u-u^{2}v,(x,t)\in \Omega \times \left[0,T\right]\\ \frac{\partial u}{\partial n}=\frac{\partial v}{\partial n}=0,\left(x,t\right)\in \partial \Omega \times \left[0,T\right]\\ u\left(x,0\right)={\psi }_1\left(x\right),v\left(x,0\right)={\psi }_2\left(x\right),x\in \Omega \end{array}$ (12)

这里 $\alpha ,\beta >0$， $u_1,u_2$ 为常数， $t\in R_{+}$ 为时间变量， $x,y\in R$ 是空间变量， $u,v$ 是t与 $x,y$ 的函数。其中 $u_1=0.0025,u_2=0.0032,\alpha =2,\beta =1,v=ones\left(N\right)$，固定一个初始条件，当u变化时模拟的结果见图1~10，具体初始条件见表1。

例2：考虑以下形式的Brusselator模型 [3]：

$\{\begin{array}{l}\frac{\partial u}{\partial t}=\Delta \left(u_{11}u+u_{12}v\right)-\left(\beta +1\right)u+u^{2}v+\alpha ,\left(x,t\right)\in \Omega \times \left[0,T\right]\\ \frac{\partial v}{\partial t}=\Delta \left(u_{21}u+u_{22}v\right)+\beta u-u^{2}v,\left(x,t\right)\in \Omega \times \left[0,T\right]\\ \frac{\partial u}{\partial n}=\frac{\partial v}{\partial n}=0,\left(x,t\right)\in \partial \Omega \times \left[0,T\right]\\ u\left(x,0\right)={\psi }_1\left(x\right),v\left(x,0\right)={\psi }_2(x),x\in \Omega \end{array}$ (13)

这里 $\alpha ,\beta >0$， $u_{11},u_{12},u_{21},u_{22}$ 为常数， $t\in R_{+}$ 为时间变量， $x,y\in R$ 是空间变量， $u,v$ 是t与 $x,y$ 的函数。其中 $u_{11}=0.4,u_{12}=22.2665,u_{21}=0.02,u_{22}=2,\alpha =6,\beta =1$ ；本文方法得到的数值结果见图11~22，具体的初始条件见表2。

Figure 1. Shows the numerical solution with the initial condition of $u=\sin \left(1000000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ of example 1

图1. 例1初始条件为 $u=\sin \left(1000000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ 的数值解

Figure 2. Shows the pattern with the initial condition of $u=\sin \left(1000000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ example 1

图2. 例1初始条件为 $u=\sin \left(1000000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ 的斑图

Figure 3. Shows the numerical solution with the initial condition of $u=\sin \left(1000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ of example 1

图3. 例1初始条件是 $u=\sin \left(1000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ 的数值解

Figure 4. Shows the pattern with the initial condition of $u=\sin \left(1000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ example 1

图4. 例1初始条件是 $u=\sin \left(1000\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ 的斑图

Figure 5. Shows the numerical solution with the initial condition of $u=\text{sech}\left(100\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ of example 1

图5. 例1初始条件是 $u=\text{sech}\left(100\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ 的数值解

Figure 6. Shows the pattern with the initial condition of $u=\text{sech}\left(100\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ example 1

图6. 例1初始条件是 $u=\text{sech}\left(100\left({\left(x-1/3\right)}^2-{\left(y-1/2\right)}^2\right)\right)$ 的的斑图

Figure 7. Shows the numerical solution with the initial condition of $u=\text{sech}\left(x^2/0.02-2\pi y^2/0.01\right)$ of example 1

图7. 例1初始条件是 $u=\text{sech}\left(x^2/0.02-2\pi y^2/0.01\right)$ 的数值解

Figure 8. Shows the pattern with the initial condition of $u=\text{sech}\left(x^2/0.02-2\pi y^2/0.01\right)$ example 1

图8. 例1初始条件是 $u=\text{sech}\left(x^2/0.02-2\pi y^2/0.01\right)$ 的斑图

Figure 9. Shows the numerical solution with the initial condition of $u=12/25+\sin \left(\pi x^2-510\right)+\cos \left(\pi y-220\right)$ of example 1

图9. 例1初始条件是 $u=12/25+\sin \left(\pi x^2-510\right)+\cos \left(\pi y-220\right)$ 的数值解

Figure 10. Shows the pattern with the initial condition of $u=12/25+\sin \left(\pi x^2-510\right)+\cos \left(\pi y-220\right)$ example 1

图10. 例1初始条件是 $u=12/25+\sin \left(\pi x^2-510\right)+\cos \left(\pi y-220\right)$ 的斑图

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

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

Figure 11. Shows the numerical solution with the initial condition of $u=\cos \left(-\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)+{\text{e}}^{-\sin \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)}$ of example 2

图11. 图为例2初始条件是 $u=\cos \left(-\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)+{\text{e}}^{-\sin \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)}$ 的数值解

Figure 12. Shows the pattern with the initial condition of $u=\cos \left(-\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)+{\text{e}}^{-\sin \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)}$ example 2

图12. 例2初始条件是 $u=\cos \left(-\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)+{\text{e}}^{-\sin \left({\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2\right)}$ 的斑图

Figure 13. Shows the numerical solution with the initial condition of $u=\sin \left(\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left({\text{e}}^{{\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2}\right)$ of example 2

图13. 例2初始条件是 $u=\sin \left(\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left({\text{e}}^{{\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2}\right)$ 的数值解

Figure 14. Shows the pattern with the initial condition of $u=\sin \left(\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left({\text{e}}^{{\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2}\right)$ example 2

图14. 图为例2初始条件是 $u=\sin \left(\pi \left({\left(x-0.4\right)}^2+{\left(y+0.4\right)}^2\right)\right)-\sin \left({\text{e}}^{{\left(x+0.4\right)}^2+{\left(y-0.4\right)}^2}\right)$ 的斑图

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

图15. 例2初始条件是 $u=\text{sech}\left(\pi \left(-x^2+y^2\right)\right)$ 的数值解

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

图16. 图为例2初始条件是 $u=\text{sech}\left(\pi \left(-x^2+y^2\right)\right)$ 的斑图

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

图17. 例2初始条件是 $u=\sin \left(\pi \left(-x^2-y^2\right)\right)$ 的数值解

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

图18. 图为例2初始条件是 $u=\sin \left(\pi \left(-x^2-y^2\right)\right)$ 的斑图

Figure 19. Shows the numerical solution with the initial condition of $u=12/25-\sin \left(10x^2-510\right)+\sin \left(\pi y-220\right)$ of example 2

图19. 例2初始条件是 $u=12/25-\sin \left(10x^2-510\right)+\sin \left(\pi y-220\right)$ 的数值解

Figure 20. Shows the pattern with the initial condition of $u=12/25-\sin \left(10x^2-510\right)+\sin \left(\pi y-220\right)$ example 2

图20. 图为例2初始条件是 $u=12/25-\sin \left(10x^2-510\right)+\sin \left(\pi y-220\right)$ 的斑图

Figure 21. Shows the numerical solution with the initial condition of $u=12/25-\sin \left(\pi x^2-510\right)+\sin \left(\cos \left(x^2+y^2-22\right)\right)$ of example 2

图21. 例2初始条件是 $u=12/25-\sin \left(\pi x^2-510\right)+\sin \left(\cos \left(x^2+y^2-22\right)\right)$ 的数值解

Figure 22. Shows the pattern with the initial condition of $u=12/25-\sin \left(\pi x^2-510\right)+\sin \left(\cos \left(x^2+y^2-22\right)\right)$ example 2

图22. 图为例2初始条件是 $u=12/25-\sin \left(\pi x^2-510\right)+\sin \left(\cos \left(x^2+y^2-22\right)\right)$ 的斑图

Table 2. The different initial conditions corresponding to the numerical solution and pattern in example 2 are shown in Figures 11-22

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

5. 结论

本文用谱插值配点方法求解Brusselator模型，数值结果表明了该方法与理论吻合较好。本文所有程序由matlab2017b算得。

致谢

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

参考文献