微分方程有限元法的边值与初值问题相容性处理
2013-11-12李纯清龚兴厚张高文
李纯清, 龚兴厚, 张高文
(湖北工业大学材料科学与工程学院, 湖北 武汉 430068)
有限元法是一种离散的变分法或基于变分原理的离散方法.它吸收了古典Ritz-Galerkin方法和有限差分法的精髓,既以变分原理为基础,又对解空间进行离散近似,同时还利用分段插值多项式特点,从而克服了微分方程数值求解许多刺手问题1-2.
1 数值求解中的边值及初值问题处理原理与方法
以比较典型的一类微分方程及其各类边界条件为例[3]:
)+q(x)y=
f(x),x∈(a,b);
(1)
p(x)∈C1[a,b],q(x),f(x)∈C[a,b].
式(1)的Galerkin变分形式略推
(Ly,v)=(f,v),
即
(2)
对式(2)左边第一项进行分部积分:
v(a)ga-v(a)αay(a)-v(b)gb+v(b)αby(b)+
将分部积分结果代入式(2)并整理后得
)dx-
v(a)αay(a)+v(b)αby(b)=
(3)
将上面等式中两个积分项中y,v进行有限维化,也就是构造连续函数集合S的有限N+1维子集SN+1,使得y,v∈SN+1=span{φ0,φ1,φ2,…,φn},记为yN+1,vN+1(称为Galerkin连续变分问题离散近似化),并将积分区间[a,b]分段化a,x1,x2,…,xn-1,b,在每个结点上定义插值基函数φ0,φ1,…,φn(即有限元化),得到所谓总体刚度矩阵和总体载荷向量间的方程,取插值基函数φi为一次项式,其一般形式如
结合式(3),可将式(1)中的边值及初值条件的处理
(4)
2 实例验证
例 1: 初值问题
其准确解为y(x)=ex-e-x+1.
将求解区间分段化,取6个节点分别为0,0.2,0.4,0.6,0.8,1,离散为5个子区间[0,0.2]、[0.2,0.4],[0.4,0.6],[0.6,0.8],[0.8,1],在各个区间的节点上取一次插值基函数,用MatLab容易实现例1对应的Galerkin变分形式的线性方程组为
对例1给定初值条件进行处理
1×2+α0×1=g0, 1×dy5+α1×y5=g1
令α0=0;则g0=2
令α1=0;则g1=dy5
将α0=0,g0=2 ,α1=0,g1=dy5(待定)按式(4)处理方法可得到例1的初值处理后的线性方程组
对上述方程组利用MatLab的符号运算功能易求取其解
其中已知y0=y(0)=1,即-1.6193+0.8476×dy5=1,可得到dy5=3.0901,并代入上式得到例1在节点0,0.2,0.4,0.6,0.8,1的数值解(表1),与理论解十分接近.
表1 例1数值解
本例也可令α0=2,则g0=4; 令α1=1,则g1=dy5+y5.但需要建立方程组求出待定的dy5与y5,具体过程
(5)
建立方程组求解dy5,y5
解之得
代入式(5),也可求出在节点0,0.2,0.4,0.6,0.8,1的数值解,结果与前述方法基本相同.
例2:两点边值问题:
y'(1)=1,y'(2)=4
用MatLab可求取其理论解表达式,取x={1,1.25,1.5,1.75,2}处的理论值,则y(x)分别为7.6515,8.0707,8.7289,9.5469,10.4904.
同理,插值基函数取二次项式,用MatLab容易实现例2对应的Galerkin变分形式的线性方程组为
对例2给定边值条件进行处理如
1×1+α1×y0=g1
4×4+α2×y4=g2
令α1=0;则g1=1
令α2=0;则g2=16
将α1=0,g1=1 ,α2=0,g2=16按式(4)处理方法可得到例2的边值问题处理后的线性方程组
其解与理论解比较非常接近(表2).
表2 例2数值解
例3:两点边值问题:
y=-xy(1)=1,y(2)=4.
其理论解在x为1,1.25,1.5,1.75,2值时,y(x)依次对应为1.0000,1.9435,2.6945,3.3626,4.0000.
对例3给定边值条件进行处理
1×dy0+α1×1=g1
4×dy2+α2×4=g2
令α1=1,则g1=1+dy0
令α2=2,则g2=8+4×dy4.
将α1=1,g1=1+dy0;α2=2,g2=8+4×dy4(其中dy0,dy4待求),按式(4)处理方法可得到例3的边值问题处理后的线性方程组
同样,借助MatLab符号运算功能,求出上述方程的解为
(6)
由于y0=y(1)=1,y4=y(2)=4,可如下建立方程组,求出dy0,dy4.
再将dy0,dy4结果代入式(6),可计算出例3在节点x=1,1.25,1.5,1.75,2处的数值解(表3),与理论值也非常接近.
表3 例3数值解
3 讨论与总结
1)无论是微分方程的初始条件,还是边界条件,本处理方法都适用.尤其是非齐次边界条件,比打靶法优越.
2)本文提出的对微分方程的初值和边值问题的处理方法与技巧主要是针对它的Galerkin变分形式进行讨论的,但也同样适用于微分方程的Ritz变分形式,因为Ritz变分形式只是Galerkin变分形式的一种特殊情况.
3)处理的依据是微分方程对应的等价泛函变分形式.在寻求变分形式过程中,微分方程特定的初始条件和(或)边界条件的也映射为泛函变分形式另一种特定形式,第三类边界条件可以说是两者间的映射公式.
4)一般而言,并不是所有的边值与初值问题存在对应的变分形式,这就使得有限元法应用受到限制.其实,如果给定的边值或初值条件能与微分方程相兼容,认为有限元法就可应用,比如实际物理问题(温度场,应变场等)的微分方程描述.
5)将微分方程不进行变分方法,而直接应用有限元法进行数值求解,此问题有待探索.因为有限元思想是独立,当与Ritz-Galerkin变分近似解法相结合,才有用有限元法数值求解微分方程.
[参考文献]
[1] 施妙根. 科学和工程计算基础[M]. 北京:清华大学出版社, 2002(7):239-266.
[2] 老大中. 变分法基础[M]. 北京:国防工业出版社, 2004(9):285-323.
[3] 姚端正. 数学物理方法[M]. 北京:科学出版社, 2010(3):(220-235.