扩展Fisher-Kolmogorov方程的格子Boltzmann方法.pdf

第4 4卷 第4期 吉首大学学报(自然科学版)V o l.4 4 N o.4 2 0 2 3年7月J o u r n a l o f J i s h o u U n i v e r s i t y(N a t u r a l S c i e n c e s E d i t i o n)J u l.2 0 2 3 文章编号:1 0 0 7 2 9 8 5(2 0 2 3)0 4 0 0 1 9 1 2扩展F i s h e r-K o l m o g o r o v方程的格子B o l t z m a n n方法*冯颖欣,戴厚平,汪 辰,魏雪丹(吉首大学数学与统计学院,湖南 吉首 4 1 6 0 0 0)摘 要:基于格子B o l t z m a n n方法构造了一种新的求解扩展F i s h e r-K o l m o g o r o v(E F K)方程的数值格式.利用T a y l o r展开和C h a p m a n-E n s k o g多尺度展开技术,得到一维、二维E F K方程格子B o l t z m a n n模型的局部平衡态分布函数和修正函数,进一步恢复对应宏观方程.在此基础上,推导出模型迭代格式的稳定性条件,并利用模型求解E F K方程的初边值问题.关键词:扩展F i s h e r-K o l m o g o r o v方程;格子B o l t z m a n n方法;数值计算;稳定性中图分类号:O 2 4 2.1 文献标志码:AD O I:1 0.1 3 4 3 8/j.c n k i.j d z k.2 0 2 3.0 4.0 0 3F i s h e r-K o l m o g o r o v(F K)方程是描述生物、化学等学科中反应扩散现象的一类重要方程.1 9 8 7年,C o u l l e t等1为了进一步扩大方程的模拟优势,在F K方程中添加了一项稳定的四阶导数项,得到扩展F i s h e r-K o l m o g o r o v(E x t e n d e d F i s h e r-K o l m o g o r o v,E F K)方程.E F K方程在许多领域表现出巨大的模拟优势,如双稳态系统中发生的时空湍流和图像形成12、族群遗传3、液晶体中磁畴壁的扩散现象4及L i f-s h i t z点附近发生的二阶相变5等.但由于源项函数为0时,E F K方程的精确解暂时无法求得6,因此学者采用有限差分法78、正交三次样条配置法9、五次B样条微分求积方法1 0等对E F K方程进行了数值研究.而这些方法的差分格式较复杂,在求解方程初边值问题时,不利于灵活处理边界,也无法进行大规模的并行计算.近年来,格子B o l t z m a n n方法(L a t t i c e B o l t z m a n n M e t h o d,L BM)作为一种数值模拟方法,广泛应用于各种复杂流体系统1 11 3中.这种方法也被越来越多的学者用于非线性偏微分方程的数值求解,如变系数反应扩散方程1 4、N a v i e r-S t o k e s方程1 5和分数阶对流扩散方程1 6等.与传统的数值计算方法1 71 9相比,L BM更加简单、灵活,在处理复杂边界条件时效率更高,代码编译过程也较系统化.因此,笔者拟采用L BM求解一维、二维E F K方程,为其数值解的研究提供一种新途径.1 E F K方程的格子B o l t z m a n n模型1.1 格子B o l t z m a n n模型考虑以下形式的一维、二维E F K方程2 0:ut+uX X X X-uX X+g(u)=F(X,t).(1)其中:为常数;g(u)=u3-u;F(X,t)为源项函数.鉴于(1)式的四阶特性,笔者分别采用D 1 Q 5和D 2 Q 1 7离散速度模型建立一维、二维E F K方程的格子*收稿日期:2 0 2 3 0 3 1 9基金项目:湖南省自然科学基金资助项目(2 0 2 1 J J 3 0 5 4 8);湖南省教育厅科学研究重点项目(2 1 A 0 3 2 9);吉首大学校级科研项目(J D Y 2 2 0 0 6);湖南省教育厅科学研究一般项目(2 2 C 0 2 7 9)作者简介:冯颖欣(1 9 9 8),女,陕西铜川人,吉首大学数学与统计学院硕士研究生,主要从事微分方程数值解研究通信作者:戴厚平(1 9 7 9),男,湖南邵阳人,吉首大学数学与统计学院副教授,博士,主要从事微分方程数值解研究.B o l t z m a n n(L a t t i c e B o l t z m a n n,L B)模型.该模型的演化方程为fi(X+c eit,t+t)-fi(X,t)=-1(fi(X,t)-f(e q)i(X,t)+t Ri(X,t).(2)其中:fi(X,t)为粒子在(X,t)处的密度分布函数;f(e q)i(X,t)为粒子在(X,t)处的平衡态分布函数;t Ri(X,t)为修正项函数;为弛豫时间;c为网格速度,c ei=ci,这里ci为不同方向上的格子离散速度.利用T a y l o r公式将(2)式的左侧展开到t4,可得t(t+ci?)fi+t22(t+ci?)2fi+t36(t+ci?)3fi+t42 4(t+ci?)4fi+O(t5)=-1(fi-f(e q)i)+t Ri.(3)对(3)式进行如下形式的C h a p m a n-E n s k o g多尺度展开:t=t0+t1+2t2+3t3,fi=f(0)i+f(1)i+2f(2)i+3f(3)i+4f(4)i+O(5),Ri=3R(1)i,D=t0+ci?.(4)其中:为K n u d s e n(克努森)数,可表示为粒子平均自由程与特征长度的比值;t与同数量级.将(4)式代入(3)式,并比较的各阶系数,可得O(0):f(0)i=f(e q)i,(5)O(1):C1D f(0)i=-1f(1)i,(6)O(2):C2D2f(0)i+t1f(0)i=-1f(2)i,(7)O(3):C3D3f(0)i+2C2t1D f(0)i+t2f(0)i=-1f(3)i,(8)O(4):C4D4f(0)i+3C3t1D2f(0)i+2C2t2D f(0)i+t3f(0)i+C22t21f(0)i=-1f(4)i+R(1)i.(9)其中:C1=1;C2=12-;C3=2-+16;C4=-3+322-71 2+12 4.定义宏观量u(X,t)为点(X,t)处各方向的局部粒子分布函数之和,即u(X,t)=ifi(X,t).对(5)式关于i求和,并由两者质量守恒可知,局部平衡态分布函数应满足if(e q)i=if(0)i=u.结合(4)式,有if(k)i=0(k=1,2,3,4).对(6)(9)式关于i求和,可得iD f(0)i=0,(1 0)t1if(0)i+C2iD2f(0)i=0,(1 1)t2if(0)i+C3iD3f(0)i+2C2t1iD f(0)i=0,(1 2)t3if(0)i+C4iD4f(0)i+3C3t1iD2f(0)i+2C2t2iD f(0)i+C22t21if(0)i=iR(1)i.(1 3)1.2 恢复一维E F K方程基于D 1 Q 5模型,构建一维E F K方程的L B模型,其离散速度为(e0,e1,e2,e3,e4)=(0,1,-1,2,-2).02吉首大学学报(自然科学版)第4 4卷定义平衡态分布函数f(0)i(i=0,1,2,3,4)的k(k=0,1,2,3,4)阶矩为if(0)i=A1u,icif(0)i=0,ic2if(0)i=A2u,ic3if(0)i=0,ic4if(0)i=A3u,iRi=-g(u)+F(x,t).(1 4)将(1 4)式代入(1 0)(1 3)式,可得A1ut0=0,(1 5)A1ut1+A2C22ux2=0,(1 6)A1ut2=0,(1 7)A1ut3+A3C44ux4+A1C22ut21+3A2C33ut1x2=iR(1)i.(1 8)将(1 5)式+(1 6)式+2(1 7)式+3(1 8)式,可得A1ut+A2C22ux2+3A3C44ux4+g(u)+3A1C22ut21+33A2C33ut1x2=F(x,t).(1 9)令A1=1,A2C2=-1,3A3C4=,E1=-3A1C22ut21-33A2C33ut1x2,则(1 9)式可恢复为如下一维E F K方程:ut+4ux4-2ux2+g(u)=F(x,t)+E1,其中E1为误差项.1.3 恢复二维E F K方程基于D 2 Q 1 7模型,构建二维E F K方程的L B模型,其离散速度为(e0,e1,e2,e3,e4,e5,e6,e7,e8,e9,e1 0,e1 1,e1 2,e1 3,e1 4,e1 5,e1 6)=0 1 0-101-1-112 0-202-2-220 0 10-1 11-1-1 0 20-2 22-2-2 .定义平衡态分布函数f(0)i(i=0,1,1 6)的k(k=0,1,2,3,4)阶矩为if(0)i=A1u,icif(0)i=(0,0),ic2if(0)i=(A2u,0,A2u),ic3if(0)i=(0,0,0,0),ic4if(0)i=(A3u,0,0,0,A3u),iRi=-g(u)+F(x,y,t).(2 0)12第4期 冯颖欣,等:扩展F i s h e r-K o l m o g o r o v方程的格子B o l t z m a n n方法将(2 0)式代入(1 0)(1 3)式,可得A1ut0=0,(2 1)A1ut1+A2C22ux2+A2C22uy2=0,(2 2)A1ut2=0,(2 3)A1ut3+A3C44ux4+A3C44uy4+A1C22ut21+3A2C33ut1x2+3A2C33ut1y2=iR(1)i.(2 4)将(2 1)式+(2 2)式+2(2 3)式+3(2 4)式,可得A1ut+A2C2(2ux2+2uy2)+3A3C4(4ux4+4uy4)+g(u)+A1C232ut21+3A2C333ut1x2+3A2C333ut1y2=F(x,t).(2 5)令A1=1,A2C2=-1,3A3C4=,E2=-(A1C232ut21+3A2C33(3ut1x2+3ut1y2),则(2 5)式可恢复为如下二维E F K方程:ut+4ux4+4uy4-2ux2-2uy2+g(u)=F(x,y,t)+E2,其中E2为误差项.2 稳定性分析接下来证明一维、二维E F K方程的L B模型是稳定的.由(1 4)式,求得如下一维E F K方程的L B模型的平衡态分布函数:f(0)0=(1+54c2 C2+4c43C4)u,f(0)1=-(23c2 C2+6c43C4)u,f(0)2=-(23c2 C2+6c43C4)u,f(0)3=(12 4c2 C2+2 4c43C4)u,f(0)4=(12 4c2 C2+2 4c43C4)u,Ri=-g(u)5=15(u-u3).采用经典有限差分法重写L B模型的迭代格式.定义fni,j=fi(xj,tn),unj=u(xj,tn),tn=nt,xj=jx,于是迭代格式可改写为22吉首大学学报(自然科学版)第4 4卷fn+10,j=(1-1)fn0,j+1(1+54c2 C2+4c43C4)unj+t5(unj-(unj)3),fn+11,j+1=(1-1)fn1,j-1(23c2 C2+6c43C4)unj+t5(unj-(unj)3),fn+12,j-1=(1-1)fn2,j-1(23c2 C2+6c43C4)unj+t5(unj-(unj)3),fn+13,j+2=(1-1)fn3,j+1(12 4c2 C2+2 4c43C4)unj+t5(unj-(unj)3),fn+14,j-2=(1-1)fn4,j+1(12 4c2 C2+2 4c43C4)unj+t5(unj-(unj)3),(2 6)且迭代格式(2 6)符合如下初始条件:f00,j=(1+54c2 C2+4c43C4)u0j,f01,j=-(23c2 C2+6c43C4)u0j,f02,j=-(23c2 C2+6c43C4)u0j,f03,j=(12 4c2 C2+2 4c43C4)u0j,f04,j=(12 4c2 C2+2 4c43C4)u0j.在t=(n+1)t时,有un+1j=fn+10,j+fn+11,j+fn+12,j+fn+13,j+fn+14,j.假设初始值u0j=u0(xj)是有界且足够光滑的,接下来证明L B模型迭代格式在L1L上是稳定的.假设|u0(xj)|1,定义un+1j=fn+10,j+fn+11,j+fn+12,j+fn+13,j+fn+14,j=(1-1)(fn0,j+fn1,j-1+fn2,j+1+fn3,j-2+fn4,j+2)+1(1+54c2 C2+4c43C4)unj-1(23c2 C2+6c43C4)unj-1-1(23c2 C2+6c43C4)unj+1+1(12 4c2 C2+2 4c43C4)unj-2+1(12 4c2 C2+2 4c43C4)unj+2+t5(unj-(unj)3)+(unj-1-(unj-1)3)+(unj+1-(unj+1)3)+(unj-2-(unj-2)3)+(unj+2-(unj+2)3)H(fn0,j,fn1,j-1,fn2,j+1,fn3,j-2,fn4,j+2,unj,unj-1,unj+1,unj-2,unj+2)H(v1,v2,v3,v4,v5,v6,v7,v8,v9,v1 0).命题1 假设n0,1212+66,|vi|1(i=6,7,1 0),且42(-2 54+c2(-12)2+5(-12)c2)c2(2-+11 2)5,(2 7)则H(v1,v2,v3,v4,v5,v6,v7,v8,v9,v1 0)是单调的.证明 在命题条件下,可得32第4期 冯颖欣,等:扩展F i s h e r-K o l m o g o r o v方程的格子B o l t z m a n n方法Hv1=Hv2=Hv3=Hv4=Hv5=1-10,Hv6=1(1+54c2 C2+4c43C4)+t5(1-3(unj)2)1(1+54c2 C2+4c43C4)+t50,Hv7=-1(23c2 C2+6c43C4)+t5(1-3(unj-1)2)1(-23c2 C2-6c43C4)+t50,Hv8=-1(23c2 C2+6c43C4)+t5(1-3(unj+1)2)1(-23c2 C2-6c43C4)+t50,Hv9=1(12 4c2 C2+2 4c43C4)+t5(1-3(unj-2)2)1(12 4c2 C2+2 4c43C4)+t50,Hv1 0=1(12 4c2 C2+2 4c43C4)+t5(1-3(unj+2)2)1(12 4c2 C2+2 4c43C4)+t50.因此结论明显成立.证毕.命题2 若1212+66,|u0(xj)|1,n0,且单调条件(2 7)成立,则对于jZ,lZ,有m i nl u0l+T F(u0m i n)un+1jm a xl u0l+T F(u0m a x).证明 令unm i n=m i nl unl,unm a x=m a xl unl(lZ),由初始迭代格式可知,对于jk(k=0,1,2,3,4),当n=1时,有4k=0f1k,jk=H(f00,j0,f01,j1-1,f02,j2+1,f03,j3-2,f04,j4+2,u0j0,u0j1-1,u0j2+1,u0j2+1,u0j4+2)(1-1)u0m a x+1u0m a x+t F(u0m a x)=u0m a x+t F(u0m a x).同理,有4k=0f1k,jk=H(f00,j0,f01,j1-1,f02,j2+1,f03,j3-2,f04,j4+2,u0j0,u0j1-1,u0j2+1,u0j2+1,u0j4+2)u0m i n+t F(u0m i n).由数学归纳法,假定对于n,有u0m i n+nt F(u0m i n)4k=0fnk,jku0m a x+nt F(u0m a x),即u0m i n+nt F(u0m i n)unju0m a x+nt F(u0m a x),则对于n+1,有4k=0fn+1k,jk=H(fn0,j0,fn1,j1-1,fn2,j2+1,fn3,j3-2,fn4,j4+2,unj0,unj1-1,unj2+1,unj3-2,unj4+2)(1-1)(u0m a x+nt F(u0m a x)+1(u0m a x+nt F(u0m a x)+t F(u0m a x)=u0m a x+(n+1)t F(u0m a x)=u0m a x+(1+1n)nt F(u0m a x)n+u0m a x+nt F(u0m a x)=u0m a x+T F(u0m a x).同理,有4k=0fn+1k,jk=H(fn0,j0,fn1,j1-1,fn2,j2+1,fn3,j3-2,fn4,j4+2,unj0,unj1-1,unj2+1,unj3-2,unj4+2)(1-1)(u0m i n+nt F(u0m i n)+1(u0m i n+nt F(u0m i n)+t F(u0m i n)=u0m i n+(n+1)t F(u0m i n)=u0m i n+(1+1n)nt F(u0m i n)n+u0m i n+nt F(u0m i n)=u0m i n+T F(u0m i n).42吉首大学学报(自然科学版)第4 4卷证毕.以下给出稳定性的证明.假设u(x,t)是满足初始条件u(x,0)=u0(x),|u0(x)|1,并采用相同迭代格式(2 6)和单调条件(2 7)的一维E F K方程的另一个解.命题3 假设满足命题2的条件,则存在如下不等式:jm a xun+1j,un+1jjm a xu0j,u0j+(n+1)jm a xt F(u0j),t F(u0j),jm i nun+1j,un+1jjm i nu0j,u0j+(n+1)jm i nt F(u0j),t F(u0j).证明 由j(m a xf00,j,f00,j+m a xf01,j,f01,j+m a xf02,j,f02,j+m a xf03,j,f03,j+m a xf04,j,f04,j)=jm a xu0j,u0j,有jm a xu1j,u1jj(m a xf10,j,f10,j+m a xf11,j,f11,j+m a xf12,j,f12,j+m a xf13,j,f13,j+m a xf14,j,f14,j)(1-1)j(m a xf00,j,f00,j+m a xf01,j,f01,j+m a xf02,j,f02,j+m a xf03,j,f03,j+m a xf04,j,f04,j)+1jm a xu0j,u0j+jm a xt F(u0j),t F(u0j)=jm a xu0j,u0j+jm a xt F(u0j),t F(u0j).采用归纳法,可得jm a xun+1j,un+1j(1-1)(jm a xu0j,u0j+njm a xt F(u0j),t F(u0j)+1(jm a xu0j,u0j+njm a xt F(u0j),t F(u0j)+jm a xt F(u0j),t F(u0j)jm a xu0j,u0j+(n+1)jm a xt F(u0j),t F(u0j).同理,可得jm i nun+1j,un+1jjm i nu0j,u0j+(n+1)jm i nt F(u0j),t F(u0j).证毕.定理1 L B模型的迭代格式(2 6)在L1上是稳定的.证明 利用|un+1j-un+1j|=m a xun+1j,un+1j-m i nun+1j,un+1j及命题3,可得j|un+1j-un+1j|j|u0j-u0j|+(n+1)t|F(u0j)-F(u0j)|=j|u0j-u0j|+(1+1n)nt|F(u0j)-F(u0j)|n+j|u0j-u0j|+nt|F(u0j)-F(u0j)|=j|u0j-u0j|+T|F(u0j)-F(u0j)|,即L B模型的迭代格式(2 6)在L1上是稳定的.证毕.同理,假设120,jZ,得到如下二维E F K方程L B模型的稳定性条件:2(-2 52+c2(-12)2+5(-12)c2)2c2(2-+11 2)5.52第4期 冯颖欣,等:扩展F i s h e r-K o l m o g o r o v方程的格子B o l t z m a n n方法3 数值算例现分别给出一维、二维E F K方程的数值算例,用以验证本研究构建的L B模型对求解一维、二维E F K方程初边值问题的有效性.引入全局相对误差(G l o b a l R e l a t i v e E r r o r,用G R E表示)和最大绝对误差(M a x i m u m A b s o l u t e E r r o r,用MA E表示),来分析初边值问题精确解和数值解的误差.G R E和MA E的定义如下:G R E=i|u(Xi,t)-u*(Xi,t)|i|u*(Xi,t)|,MA E=m a xi|u(Xi,t)-u*(Xi,t)|.其中:u(Xi,t)和u*(Xi,t)分别为(Xi,t)处E F K方程的数值解和精确解,u(Xi,t)通过uk+1j迭代计算求得;k为迭代次数.此外,选取非平衡外推格式进行边界处理.算例1 考虑如下形式的一维E F K方程的初边值问题:ut+ux x x x-ux x+g(u)=F(x,t)x-8,8,t0,T,u(x,0)=s i n(4(x-4)u(-8,t)=u(8,t)=0.(2 8)其中:g(u)=u3-u;源项函数F(x,t)=12 5 6e-ts i n(4(x-4)(4+2 5 6 e-2ts i n(4(x-4)2+1 62-5 1 2),其精确解u(x,t)=e-ts i n(4(x-4).取定参数=0.0 1,x=0.1 6,t=0.0 0 1,=0.5 6.表1给出了不同时刻(T)下,u的全局相对误差和最大绝对误差.表1 不同时刻下u的全局相对误差和最大绝对误差T a b l e 1 G l o b a l R e l a t i v e E r r o r a n d M a x i m u m A b s o l u t e E r r o r o f u w i t h D i f f e r e n t T i m eTG R E/1 0-3MA E/1 0-30.13.6 6 11 5.3 6 90.52.8 0 31 0.4 9 50.62.6 3 69.5 0 11.01.9 3 46.3 7 41.51.5 0 13.8 6 22.03.4 6 82.3 3 02.58.4 0 71.3 9 4 由表1可知,G R E明显保持在1 0-3数量级,MA E也在T=0.5之后保持在1 0-3数量级,数值解与精确解基本一致,表明一维E F K方程的L B模型求解问题(2 8)是有效的.T=0.1,0.5,1时,问题(2 8)的数值解和精确解如图1所示.图1 不同时刻下问题(2 8)的数值解和精确解F i g.1 N u m e r i c a l S o l u t i o n a n d E x a c t S o l u t i o n o f P r o b l e m(2 8)w i t h D i f f e r e n t T i m e62吉首大学学报(自然科学版)第4 4卷由图1可见,数值解和精确解吻合度较高,进一步显示了L BM在一维E F K方程的初边值问题求解中的适用性.表2示出了T=1时,不同格子数(N)下的全局相对误差和最大绝对误差.表2 不同格子数下的全局相对误差和最大绝对误差T a b l e 2 G l o b a l R e l a t i v e E r r o r a n d M a x i m u m A b s o l u t e E r r o r w i t h D i f f e r e n t L a t t i c eNxG R E/1 0-3MA E/1 0-3C P U/s6 00.2 71.6 3 72.5 0 90.0 4 2 1 87 00.2 30.9 1 51.3 0 00.0 2 6 8 08 00.2 01.4 0 12.7 9 90.0 2 5 5 29 00.1 81.5 9 94.4 1 00.0 2 5 8 91 0 00.1 61.9 3 46.3 7 40.0 2 8 1 01 1 00.1 52.7 8 69.0 5 80.0 3 4 0 21 2 00.1 34.6 7 41.3 5 70.0 3 5 4 0 由表2可知,G R E和MA E在空间方向上的收敛速度极快.算例2 考虑如下形式的一维E F K方程的初边值问题:ut+ux x x x-ux x+g(u)=0 x-4,4,t0,T,u(x,0)=1 0-3e x p(-x2)u(-4,t)=u(4,t)=1.(2 9)由于问题(2 9)暂时无法求得精确解1 0,因此笔者采用一维E F K方程的L B模型对其进行数值求解,并与文献1 0进行图像对比.取定参数=0.0 0 0 1,x=0.1,t=0.0 0 1.图2示出了N=8 0时,不同时刻T下问题(2 9)的数值解.图2 不同时刻下问题(2 9)的数值解F i g.2 N u m e r i c a l S o l u t i o n o f P r o b l e m(2 9)w i t h D i f f e r e n t T i m e对比图2与文献1 0中的图3可知,两图形基本一致,即在时间增加的条件下数值解衰变并达到1的稳定状态.这说明,一维E F K方程的L B模型对求解问题(2 9)具有较强的适用性.算例3 考虑如下形式的二维E F K方程的初边值问题:ut+4ux4+4uy4-2ux2-2uy2+g(u)=0 (x,y)-4,4-4,4,t0,T,u(x,y,0)=s i n(4(x-4)s i n(4(y-4),u(-4,y,t)=u(4,y,t)=0,u(x,-4,t)=u(x,4,t)=0.(3 0)其中:g(u)=u3-u;源项函数72第4期 冯颖欣,等:扩展F i s h e r-K o l m o g o r o v方程的格子B o l t z m a n n方法F(x,y,t)=e-ts i n(4(x-4)s i n(4(y-4)(1 2 8 e-2ts i n(4(x-4)2s i n(4(y-4)2+4+1 6 2-2 5 6)1 2 8,其精确解u(x,y,t)=e-ts i n(4(x-4)s i n(4(y-4).取定参数=0.0 1,x=41 5,t=0.0 0 1,=0.5 8 7.表3给出了不同时刻(T)下,u的全局相对误差和最大绝对误差.表3 不同时刻下u的全局相对误差和最大绝对误差T a b l e 3 G l o b a l R e l a t i v e E r r o r a n d M a x i m u m A b s o l u t e E r r o r o f u w i t h D i f f e r e n t T i m eTG R E/1 0-3MA E/1 0-30.15.4 6 66.9 6 70.54.9 6 84.3 5 10.64.8 9 53.9 0 51.04.7 3 22.5 7 13.04.7 6 20.3 5 25.05.5 6 90.0 5 37.09.7 8 70.0 1 0 由表3可知,G R E和MA E都明显保持在1 0-3数量级,数值解与精确解基本一致,表明二维E F K方程的L B模型求解问题(3 0)的实用性较高.T=1,N=3 0时,图3示出了问题(3 0)的数值解和精确解,图4示出了问题(3 0)的数值解和精确解的绝对误差.图3 T=1,N=3 0时问题(3 0)的数值解和精确解F i g.3 N u m e r i c a l S o l u t i o n a n d E x a c t S o l u t i o n o f P r o b l e m(3 0)w i t h T=1,N=3 0图4 T=1,N=3 0时问题(3 0)的数值解和精确解的绝对误差F i g.4 A b s o l u t e E r r o r f o r N u m e r i c a l S o l u t i o n a n d E x a c t S o l u t i o n o f P r o b l e m(3 0)w i t h T=1,N=3 082吉首大学学报(自然科学版)第4 4卷由图3,4可见,数值解和精确解相似度较高,进一步显示了L BM在二维E F K方程的初边值问题求解中的适用性.表4给出了T=1时,不同格子数(N)下的全局相对误差和最大绝对误差.表4 不同格子数下的全局相对误差和最大绝对误差T a b l e 4 G l o b a l R e l a t i v e E r r o r a n d M a x i m u m A b s o l u t e E r r o r w i t h D i f f e r e n t L a t t i c eNG R E/1 0-3MA E/1 0-3C P U/s3 04.7 2 82.5 8 90.3 1 1 1 34 07.4 9 84.5 8 30.4 9 0 9 05 08.5 7 99.8 8 90.7 3 4 0 86 03 3.9 5 01 9.9 3 61.1 9 0 7 7 由表4可知,随着格子数的减小,G R E和MA E都逐渐减小,且在空间方向上的收敛速度较快.4 结语笔者采用L BM构建一维、二维E F K方程的L B模型,并用其求解了一维、二维E F K方程的初边值问题.L B模型可以灵活地处理非线性微分方程中的高阶项和非线性项,模型数值格式较简单.此外,笔者给出了L B模型迭代格式的稳定条件,并选择合适的参数进行了数值实验.实验结果表明,初边值问题的数值解与精确解一致,全局相对误差保持在1 0-3数量级左右,模型的可靠性高.接下来,笔者考虑将本研究结果推广到其他高阶非线性偏微分方程的数值求解中,进一步拓宽L BM的应用范围.参考文献:1 C OU L L E T P,E L P H I C K C,R E P AUX D.N a t u r e o f S p a t i a l C h a o sJ.P h y s i c a l R e v i e w L e t t e r s,1 9 8 7,5 8:4 3 1 4 3 4.2 D E E G T,VAN S AA R L OO S W I M.B i s t a b l e S y s t e m s w i t h P r o p a g a t i n g F r o n t s L e a d i n g t o P a t t e r n F o r m a t i o nJ.P h y s i-c a l R e v i e w L e t t e r s,1 9 8 8,6 0(2 5):2 6 4 1 2 6 4 4.3 A R ON S ON D G,WE I N B E R G E R H F.M u l t i d i m e n s i o n a l N o n l i n e a r D i f f u s i o n A r i s i n g i n P o p u l a t i o n G e n e t i c sJ.A d-v a n c e s i n M a t h e m a t i c s,1 9 7 8,3 0(1):3 3 7 6.4 Z HU G u o z h e n.E x p e r i m e n t s o n D i r e c t o r W a v e s i n N e m a t i c L i q u i d C r y s t a l sJ.P h y s i c a l R e v i e w L e t t e r s,1 9 8 2,4 9(1 8):1 3 3 2 1 3 3 5.5 HO R N R E I CH R M,L U B AN MA R S HA L L,S HT R I KMAN S.C r i t i c a l B e h a v i o r a t t h e O n s e t o f k-S p a c e I n s t a b i l i t y o n t h e L i n eJ.P h y s i c a l R e v i e w L e t t e r s,1 9 7 5,3 5(2 5):1 6 7 8 1 6 8 1.6 李 娟,高广花.二维扩展F i s h e r-K o l m o g o r o v方程的线性化紧差分格式的最大模误差分析J.西南师范大学学报(自然科学版),2 0 1 7,4 2(3):1 2 2 1.7 SWE I L AMN H,E L S AKOUT D M,MUT T A R D I M M.N u m e r i c a l S o l u t i o n f o r S t o c h a s t i c E x t e n d e d F i s h e r-K o l m o g o r o v E q u a t i o nJ.C h a o s,S o l i t o n s a n d F r a c t a l s,2 0 2 1,1 5 1:1 1 1 2 1 3.D O I:1 0.1 0 1 6/j.c h a o s.2 0 2 1.1 1 1 2 1 3.8 L I S h u g u a n g,X U D a,Z HAN G J i e,e t a l.A N e w T h r e e-L e v e l F o u r t h-O r d e r C o m p a c t F i n i t e D i f f e r e n c e S c h e m e f o r t h e E x-t e n d e d F i s h e r-K o l m o g o r o v E q u a t i o nJ.A p p l i e d N u m e r i c a l M a t h e m a t i c s,2 0 2 2,1 7 8:4 1 5 1.9 D ANUM J AYA P,P AN I AM I YA K.O r t h o g o n a l C u b i c S p l i n e C o l l o c a t i o n M e t h o d f o r t h e E x t e n d e d F i s h e r-K o l m o g o r o v E q u a t i o nJ.J o u r n a l o f C o m p u t a t i o n a l a n d A p p l i e d M a t h e m a t i c s,2 0 0 5,1 7 4(1):1 0 1 1 1 7.1 0 M I T T A L R C,D AH I YA S UM I TA.A S t u d y o f Q u i n t i c B-S p l i n e B a s e d D i f f e r e n t i a l Q u a d r a t u r e M e t h o d f o r a C l a s s o f S e m i-L i n e a r F i s h e r