Sobolev方程的半有限元方法
2015-03-21辜继明瞿少成
辜继明, 俞 辉*, 瞿少成
(1.三峡大学 理学院, 湖北 宜昌 443000; 2.华中师范大学 信息技术系, 武汉 430079)
辜继明1, 俞 辉1*, 瞿少成2
(1.三峡大学 理学院, 湖北 宜昌 443000; 2.华中师范大学 信息技术系, 武汉 430079)
对Sobolev方程采用半有限元法进行数值模拟.通过将空间变量和时间变量分离,得到Sobolev方程的离散格式.首先对空间变量应用有限元方法进行离散化,得到常微分方程组的初值问题;再对时间变量应用有限差分法进行离散化,得到一系列线性方程组,求解可得到Sobolev方程的数值解.本文从理论上推导出了本文所讨论的Sobolev方程半有限元算法的矩阵算法格式,分析了其可行性.在最后给出了数值例子,从数值例子中进一步验证了半有限元方法的可行性.
半有限元法; Sobolev; 方程有限元法; 有限差分法
讨论以下Sobolev方程
其中,Ω⊂Rd(d=1,2,3)为方程的有界区域,∂Ω是方程的边界,u0(x),f(x,t)是已知函数,a是关于x,t的函数a(x,t)的连续函数,且a(x,t)>0,b(x,t)≥0.
Sobolev方程是偏微分方程,主要是研究不定常微分方程.对Sobolev方程的研究具有重大的物理工程意义.通过对Sobolev方程的研究,可以分析流体穿过裂缝岩石的理论,可得出土壤湿气迁移的过程,同样也可解决不同介质的热传导问题.人们对Sobolev方程很早就开始研究了,如标准混合有限元法应用于Sobolev方程[1],Galerkin-有限元逼近法[2-3],Sobolev方程的H1-Galerkin混合有限元方法[4-5].
本文从解抛物型方程的角度,提出半有限元的一种算法.由于Sobolev方程牵涉到时间和空间变量,是一个不定常的偏微分方程.半有限元方法的基本思想是将时间和空间变量进行分开处理.第1步、对空间变量应用有限元方法进行离散化,可以得到一个常微分方程组的初值问题,这一个过程称为半离散化过程;第2步、对第1步得到的初值问题中的时间变量应用有限差分法再进行离散化,得到一系列线性代数方程组,这一过程称为全离散化过程.
很多的物理或工程问题,都可以用微分方程来表示,求解微分方程的数值解最早是用有限差分法,有限差分法在求解微分方程时效率很高,尤其是在求解一、二维问题时,但是在高维复杂的问题时,因为跳跃误差存在,导致迭代出现误差累积,从而会有不收敛现象.有限元法是求解微分方程数值解最重要的方法之一,在高维复杂的问题求解时,可以得到较高的精度,但是所需计算的时间往往较大.
本文结合有限差分法和有限元法,发挥各自的优势,在空间变量上,应用有限元求解,得到精度较高的关于时间的微分方程,再用有限差分法,可以高效的求解动态的偏微分方程.这种半有限元法比单纯使用有限元方法求解时效率更高,比单纯使用有限差分法精度更高,具有深远的研究意义,尤其是推广到多维动态微分方程问题时.
本文的主要内容:在第1部分主要推导出了半有限元方法的空间变量离散格式,得到关于时间变量的一个微分方程组的初值问题,最后以矩阵的形式表达出来;第2部分主要对微分方程进行时间离散化,得到全离散格式,最后以矩阵的形式表达出来;第3部分通过数值例子说明方法的有效性.
1 半有限元法的半离散格式
考虑如下经典的一维Sobolev方程初边值问题:
引理1[6]设F是一个微分算子,则
Fu=0∀v.
从引理1出发,可以把微分方程转化为积分方程,具体做法是,将该微分方程两边乘以v并积分,得
(1)
对式(1)左端第2项应用分部积分法得
其中,
Bx(1,0,t)=a(1,t)ux′(1,t)v(1)-a(0,t)ux′(0,t)v(0),
Bt(1,0,t)=b(1,t)ut′(1,t)v(1)-b(0,t)ut′(0,t)v(0),
于是
(2)
比较式(1)和式(2)可得,原问题中要求u二阶导数存在,而在积分形式中,由于积分和分部积分技巧的运用,使得对解光滑性的要求降低,只要u以及其的一阶导数u′可积即可.
记
则式(2)可写成
考虑到原问题中的边界条件,即u,v应满足
u(0,t)=u(1,t)=v(0)=v(1)=0,
因此Bx(1,0,t)=0,Bt(1,0,t)=0,并且,由初始条件u(x,0)=u0(x)可得
(u(x,0),v)=(u0(x),v),
于是原问题转化为:
(3)
式(3)就是Sobolev方程初边值的积分形式.
引进Sobolev函数空间
半离散过程是对空间变量应用有限元方法进行处理.具体的做法如下:
对[0,1]进行任意划分:
0=x0 取多项式样条函数空间 n;uh(0)=uh(1)=0}, 则有 (4) 将式(4)写成矩阵形式得 (设为Aξ(t)) (设为b) 即 (5) 考虑边界条件,由(uh(x,0),vh)=(u0(x),vh)可得 (6) 将式(6)成矩阵形式,得 (设为D1ξ(0)) (设为g) 即 D1ξ(0)=g. (7) 由式(5)和式(7)可得如下式子 (8) 则式(8)即是将空间变量进行有限元离散化,得到的一常微分方程组的初值问题.可以得出D1,D2,A具有对称性和正定性,则式(8)具有唯一解,表示有限元的半离散化是可行的. 考虑一维经典初值问题: 差分法的基本思想是:以差商代替微商,最简单最经典是Euler法[7],本文构造的差分法是在Euler法的基础上,基于线性组合思想的线性多步法(LMS法). Euler算法思想如下: 对时间变量t等距剖分:0=t0 在[tm,tm+1]上,作 其中,u(tm)是微分方程在tm处的精确解;um是差分方程在tm的精确解. 在Euler法基础上,利用线性组合思想推广如下: 本文采用的LMS法是Adams二步外插法,此差分法是一种稳定的,并且二阶收敛的方法,具体参考文献[9].Adams二步外插法具体算法格式如下: (9) 将Adams二步外插法运用到有限元半离散形式上,得 或写成矩阵格式 即 通过差分法对时间变量进行离散化,得到Sobolev方程全离散格式如下: 容易验证此时方程的真解为 本文有限元用的是Lagrange二次插值,在计算数值积分时用到了局部高斯积分技术[8-10],有限差分法用的是Adams二步外插法,所有的程序是在Matlab2012编程实现的.本文空间步长取的是0.1,时间步长取的是0.25. 通过Matlab程序运行,得到真实值和半有限元法解,如图1所示. 图1 半有限元数值解和真实解的比较Fig.1 Solution compared with the real numerical solutions of half finite element method 图1中,横坐标是x,纵坐标是u.左上角的图表示在t=1时,通过半有限元得到的x,u的关系;右上角的图表示在t=1时,真实函数x,u的关系.左下角的图表示在t=0.75时,通过半有限元得到的x,u的关系;右下角的图表示在t=0.75时,真实函数x,u的关系.通过程序得出误差: 通过图1和误差可以得出半有限元方法在解决一维Sobolev方程是可行的,故推广到二维和三维中是非常有意义的.图2具体给出了t=0.00,0.25,0.50,0.75,1.00半有限元解和真实解的关系. 图2 半有限元解逼近真实解Fig.2 Half finite element solution to approximate the real solution Sobolev方程的研究具有重大的工程意义,本文用半有限元方法对Sobolev微分方程进行了数值解的研究.半有限元方法通过对空间变量进行离散化,得到了一组常微分方程;再对时间变量进行离散化,得到了矩阵差分格式,从而求解得到Sobolev方程的数据解.从数值结果可以得出半有限元方法的可行性,并且在一维计算的精度非常高,故将该方法推广到二维、三维甚至多维上具有重大意义. [1] Jiang Z W,Chen H Z. Error estimates for mixed finite element methods for sobolev equation[J]. Northeast Math J, 2001, 17(3): 301-314. [2] Ewing R E. Time-stepping Galerkin methods for nonlinear sobolev partial differential equations[J]. SIAM J Nmer Anal, 1978,15:1125-1150. [3] Nakao M T. Error estimate of a Galerkin method for some nonlinear Sobolev equations in one dimention[J]. Mumer Math,1985, 47:139-157. [4] Douglas J R J, Dupont T F,Wheeler M F.H1-Galerkin Methods for the Laplace and Heat Equations [M]. New York:Academic Press, 1975. [5] Pani A K. AnH1-Galerkin mixed finite element method for parabolic partial difference equations[J]. SIAM J Nmer Anal, 1998, 35: 721-727. [6] 林 群.微分方程数值解法基础教程[M].第2版.北京:科学出版社, 2003. [7] Batina J T. A Gridless Euler/Navier-Stokes Solution Algorithm for Complex Aircraft Applications[R]. AIAA Paper， 1993. [8] Li Jian. Investigations on two kinds of two-level stabilized finite element methods for the stationary Navier-Stokes equations[J]. Appl Math Comput, 2006, 182: 1470-1481. [9] Bochev P B, Dohrmann C R, Gunzburger M D. Stabilization of low-order mixed finite elements for the stokes equations[J]. SIAM J Numer Anal, 2006, 44: 82-101. [10] Li Jian, He Yinnian. A stabilized finite element method based on local Gauss integral technique for the stationary Navier-Stokes equations [J]. J Comp Appl Math, 2008, 214: 58-65. The half finite element method for Sobolev equation GU Jiming1, YU Hui1, QU Shaocheng2 (1.College of Science, China Three Gorges University, Yichang, Hubei 443000;2.Department of Information Technology, Central China Normal University, Wuhan 430079) In this paper， a half finite element method was proposed to solve the Sobolev equation. Discretization of Sobolev equation was presented by separating space and time variables. Firstly， applying finite element method to get discretization of the spatial variables， the problem was transformed into the initial value problems of ordinary differential equations. Secondly， applying finite-difference method to get discretization of the time variables， the problem was transformed into a series of linear equations. Thirdly， Numerical solution of the Sobolev equation was acquired by solving the series of linear equations. The method was proved to be feasible by theoretical analysis the Matrix algorithm format of the Sobolev equation was presented. Finally， numerical results were provided to illustrate the efficiency of our method． half finite element method; Sobolev equation; finite element method; finite-difference method 2014-10-18. 国家自然科学基金项目(61273183,61174216,61374028,61304162);湖北省自然科学基金项目(2011CDB187);湖北省高等学校优秀中青年科技创新团队项目(T201103). 1000-1190(2015)03-0334-05 O175 A *通讯联系人. E-mail: yuhui@ctgu.edu.cn.
2 半有限元法的全离散格式
3 数值例子
4 结论
