（计算数学专业论文）用间断galerkin有限元求解线性和非线性微分方程.pdf

用间断g a l e r k i n 有限元求解线性和非线性微分方程 i o 1中文摘要， 本文使用间断g a l e r k i n 有限元( d g ) 来求解二维奇异摄动问题， 并结合插值系数有限元方法求解一维、二维半线性微分方程多解问 题。 在求解一维奇异摄动问题时，间断有限元在节点处具有比连续有 限元更强的超收敛性，能够很好模拟解的剧烈变化【2 9 】。受其启发， 本文我们将用d g 方法来求解二维奇异摄动问题。在一致网格和两 种s h i s h i k i n 网格的剖分下，我们的数值例子表明d g 方法不仅没有 产生任何振荡，而且具有超收敛性质。 对于半线性椭圆型方程的多解问题，当m o r s ei n d e x 很大时，这些 解的图形含有很多陡峭的边界峰和内部峰，而d g 方法正是处理边 界峰和内部峰最有效的方法之一。我们对一维和二维问题的数值结 果证实了该方法的有效性。 关键词：间断g a l e r k i n 有限元，奇异摄动问题，半线性微分方程， 多解，插值系数有限元 i i 湖南师范大学2 0 0 6 届硕士学位论文 0 2a b s t r a c t i nt h i sp a p e r ，d i s c o n t i n u o u sg a l e r k i nm e t h o d ( d g ) i sd i s c u s s e dt os o l v e t h et w o - d i m e n s i o n a ls i n g u l a r l yp e r t u r b e dp r o b l e ma n ds e m i l i n e a rd i f f e r e n t i a l e q u a t i o n s f o rt h el a t t e r ，t h es o - c a l l e di n t e r p o l a t e dc o e f f i c i e n tf i n i t ee l e m e n t m e t h o di sc o m b i n e dt od e a lw i t ht h en o n l i n e a rt e r m a sw ec a l c u l a t e dt h ea p p r o x i m a t es o l u t i o n so fs i n g u l a r l yp e r t u r b e dc o n v e c t i o n d i f f u s i o ne q u a t i o n si nt h eo n e d i m e n s i o n a ls e t t i n g ，w ef i n dt h a tt h es u p e r c o n - v e r g e n c eo fd gm e t h o di ss t r o n g e rt h a nt h a to ft r a d i t i o n a lf i n i t ee l e m e n t m e t h o da tn o d e s m o r e o v e r ，t h ed gm e t h o da l s oc a ns i m u l a t et h ea c u t ev a r i a t i o no fs o h i t i o n 2 9 i n s p i r e db yt h o s er e s u l t s ，w eu s et h ed gm e t h o dt os o l v e t w o d i m e n s i o n a ls i n g u l a r l yp e r t u r b e dp r o b l e m u n d e ru n i f o r mm e s ha n dt w o k i n d so fs h i s h i k i nm e s h ，o u rn u m e r i c a lr e s u l t so b t a i n e db yd gm e t h o dv e r i f y t h es u p e r c o n v e r g e n c ea n dn oo s c i l l a t i o ni so b s e r v e d w i t hr e g a r dt ot h em u l t i p l es o l u t i o n so fs e m i l i n e a re l l i p t i ce q u a t i o n s ，d g m e t h o di se f f i c i e n tt od e a lw i t ht h es h a r pp e a ko fs o l u t i o n s ，w h i c hw i l lo c c u r w h e nt h em o r s ei n d e xi sh i g h ，i nt h ed o m a i no ro nt h eb o u n d a r y o u rn u m e r i c a l r e s u l t sa r ea l s op r e s e n t e dt os h o wt h ee f f i c i e n c yo ft h em e t h o d k e yw o r d s ：d i s c o n t i n u o u sg a l e r k i nm e t h o d ，s i n g u l a r l yp e r t u r b e dp r o b - l e m ，s e m i l i n e a rd i f f e r e n t i a le q u a t i o n s ，m u l t i p l es o l u t i o n s ，i n t e r p o l a t e dc o e f f i c i e n tf i n i t ee l e m e n tm e t h o d 用间断g a l e r k i n 有限元求解线性和非线性微分方程 l 第一章绪论 1 1引言 自然科学中许多问题的数学模型都是非线性微分方程，这些方 程往往有多个甚至无穷多个解，它们的结构相当复杂。这些解具有 何种性质和结构，如何在数值上求解它们，尤其是用数值方法求非 稳定解，是具有挑战性的重要课题。对这个问题的深入研究，对非线 性科学是有重要意义的。 在现代非线性科学中多解现象广为出现。细杆受压是材料力学 中最简单而重要的问题。早在2 0 0 多年前e u l e r 提出了细杆纵向受压 的稳定性问题。设杆在z 轴的( 0 ，f ) 上，一端固定，另一端加纵向压 力p 。此力不大时，杆只有微小压缩变形，不会弯曲。当p 大到一定 程度时，会发生戏剧性变化。设杆有微小横向位移u 0 ，而成为挠曲 轴，其弧长为5 ，挠曲轴的切线与z 轴的交角是日，挠曲轴的曲率为 肛= 五d o 。在压力p 作用下产生弯矩p u ，细杆弯曲变形与其曲率p = 忑d o 成正比，则细杆的弯曲方程为 上式对弧长s 微分一次，坐d s = 甭d 2 0 ，并利用u ，( s ) = s i n 0 得到一个半线 性方程( 限于考虑两端简支情形) 目( s ) + s i n 曰= 0 ，口，( 0 ) = p 懈) = 0 ， 准确求解此问题较复杂，要用到第1 类椭圆积分，其结果是：当p p l ( 临界压力) 时，它有两个非零解，其幅度的大小是确定的，与p 有 关。这就解决了压杆的弯曲问题，因此e u l e r 提出的细杆轴向受压的 弯曲不稳定问题是历史上第1 个多解的典型例子，它也是因参数p 的变化使解的个数发生改变的一种分叉现象。 但是 8 j 中第9 章对此问题的数值计算表明，当p 继续增大时， 此非线性方程还有多波的变号解，它们可能很不稳定，但毕竟是存 湖南师范大学2 0 0 6 届硕士学位论文 在的，在实际生活中，n 1 的多波弯曲现象似乎难以理解，但它在很 特定的条件下存在，可能更不稳定，在2 0 世纪对薄壁圆筒所做的轴 向加压的破坏试验中，就观察到薄壁圆筒产生的多波皱折型变形， 很像是形如s i n n x ( n 1 ) 的解。 这种分叉和多解现象不仅在各种复杂的力学结构设计中广泛出 现，而且也存在于天体物理、量子物理等领域。如天体物理中出现3 维空间r 3 中的半线性椭圆问题 乱+ k ( x ) f ( u ) = 0 ( 在r 3 内) ，弘( ) = 0 早在1 8 6 9 年h l a n e 计算太阳表面的温度和质量密度时导出了此方 程，这里k = 1 ，f ( u ) = l u r l 让，u 是引力位势，p 是多方指数，称为l a n e - e m d e n 方程。它是研究星球内部结构在自由磁场作用下的磁流体力学 方程。为考虑星球团的动力学，a s e d d i n g t o n ( 1 9 5 1 ) 和t m a t u k u m a ( 1 9 3 0 ) 先后分别改进为 k ( x ) f ( u ) = e 2 u ( i + h 2 ) ，k ( x ) f ( u ) = l u l p - 1 u ( 1 + h 2 ) 对l a n e - e m d e n 方程，c h a n d r a s e k h a r ( 1 9 6 7 ) 发现了3 种可解情形：p = 0 ，1 ，5 。 数十年后，经过一些学者的一系列数学研究工作发现，e d d i n g t o n 方程 没有任何整体正解，显然是非线性项按指数增长坏了事，而对m a t u k u m a 方程在某些限制下证明了多解的猜想。 在量子物理中最基本的方程是著名的s c h r s d i n g e r 方程 i w = x w + j w i p 一1 w n ( z ) 叫， 这里w = u + v i 是电场强度( 复函数) ，它描述例如非线性光学、晶 体、凝聚态物质等现象，这里用非线性项模拟了多质点间的相互作 用。仿照s i l b e r b e r g 寻找谐波形特解u = u ( r ) e 印( 一i a ( + i c o ) ，这里可 认为u ( r ) 20 是非负函数，它满足如下非线性具有奇异系数的特征 常微分方程 u r r + _ d - - 1 珥一2 a u + 2 3 = o ， u ( o ) = 1 ，u ( o o ) = o 这种立方非线性方程正是产生多解的重要来源之一。 用间断g a l e r k i n 有限元求解线性和非线性微分方程 3 我们讨论半线性椭圆边值问题 a u + f ( x ，札) = 0i nq ，“= 0o na q ( 1 1 ) 的多解及其数值计算，这里= 磋+ d ：是l a p l a c e 算子，御是区域 g t 的边界。对f ( x ，t ) 的正则性和增长性假设如下： ( a 1 ) f ( x ，t ) 在豆r 上局部l i p s c h i t z 连续； ( a 2 ) 当n 3 时，存在常数c ，和g 使得i f ( x ，) i c 1 + c 2 1 t l p ， 这里0 p 2 ，m 0 使得对于i t i m 有 0 0 ，那么 问题( 1 1 ) 存在无数多对解。当奇非线性项有轻微摄动时，由m o r s e 理论知( 1 1 ) 有无穷多个解。为了研究非线性方程的解( 包括多解) ， 在过去3 0 多年间，现代变分学( 临界点理论) 利用各种极大极小定理 湖南师范大学2 0 0 6 届硕士学位论文 或爬山定理，证明了许多重要的结果 1 8 。受这些极大极小定理的启 发，近十年来出现了3 种相关的计算方法：山路算法( m p a ) 2 4 ，高环 绕算法( h l a ) 2 6 ，最大最小算法( m n a ) 2 5 。对于奇非线性情形，由 于高m o r s e 指标临界点的多重性、退化性和不稳定性，m p a ，h l a 和 m n a 的应用遇到了固有的困难，在某些限制下它们都只计算出若干 个多解。 陈传淼和谢资清在【9 ，3 0 中提出了一种全新的多解计算理论和 方法一搜索延拓法( s e m ) ，它与爬山定理及各种极小极大方法无 关，是用少数几个特征基的线性组合搜索所有解的初值，然后利用 一种谨慎而有效的延拓法和有限元法迭代完成精密化计算。它在三 层子空间& 。c 品c8 n 中进行，首先在& 。中搜索所有的初值，然后 用延拓法( 甚至多启动延拓法) 进一步逼近某个所需的解u 2 & ，最 后在最大的子空间踟中以前一步中得到的解作为初值用n e w t o n 法 求解目标方程组。搜索延拓法理论上能够计算出任意高m o r s ei n d e x 的解，但事实上，随着m o r s ei n d e x 的提高，解的峰变得越来越多，并 且峰尖和峰谷也越来越细，当m o r s ei n d e x 很大时，这些解的图形含 有很多陡峭的边界峰和内部峰，数值计算的精度不是特别高，只有 当网格充分密时，才能得到没有振荡的数值解。 而近三十年来，间断有限元方法在双曲型方程的数值处理上运 用得相当广泛。间断有限元方法是1 9 7 3 年由r e e d 和h i l l 2 1 1 首先提 出，并应用于求解中子输运方程，但这种方法长期以来一直没有得 到很好的研究和应用。直到2 0 世纪8 0 年代后期和9 0 年代，c o e k b u r n 和s h u 等 1 ，2 ，3 ，4 ，5 ，1 2 ，1 5 结合r u n g e - k u t t a 方法，将间断有限元方 法推广到非线性一维守恒律方程和方程组，高维守恒律方程和方程 组，并给出了部分收敛性理论证明后6 ，1 1 ，1 4 1 ，这一方法才引起人们 的注意，并逐渐开始应用于流体力学计算领域【7 ，1 3 ，1 6 】。 间断有限元方法能够应用于流体力学计算领域，并受到人们的 关注，主要在于它保持了通常有限元方法的优点：能够处理复杂的 区域边界和具有复杂的边界条件问题，并获得与区域内部一致的计 算精度；易于网格加密和高精度处理边界条件，实现自适应计算；可 以得到任意阶精度的格式，同时又具有很好的局部紧致性，可以更 用间断g a l e r k i n 有限元求解线性和非线性微分方程 5 - 好的模拟解的剧烈改变。对于带小参数的奇异摄动问题，传统的数 值计算方法，如连续有限元、有限差分法在求解时遇到了很大的困 难。当小参数趋于0 时，解在某些局部区域变化很剧烈，形成内部层 和边界层，数值上很难模拟，用一般的连续有限元计算，数值上振荡 非常剧烈，结果很不准确。 求解这类问题的理想的数值方法，应该无闭锁，无物理振荡，并 且关于奇异参数一致最优阶收敛甚至超收敛。而间断有限元方法可 以很好的模拟解的剧烈变化，具有适合解这类问题的优点。【2 0 】中对 常微分方程初值问题的间断有限元解得到了2 p + 1 阶的超收敛估计。 【2 9 中也得到了一维奇异摄动问题的间断有限元解的2 p + 1 阶的超收 敛估计，并观察到，即使在均匀网格下，当奇异参数非常小时，也 没有任何振荡现象，这对于传统的数值方法是不可想象的。鉴于d g 方法在求解一维区域奇异摄动问题时的高效性，在第二章中我们试 图用该方法来求解二维正方形区域在c a r t e s i a n 网格下的奇异摄动问 题，据我们所知，这是d g 方法首次用于求解二维区域上的奇异摄动 问题，数值结果表明，在均匀网格下，d g 方法所产生的数值解没有 任何振荡，并且在s h i s h i k i n 网格下可以观察到p + 2 阶一致超收敛估 计。 对于非线性椭圆型方程，当解的m o r s e 指标很高时，解含有多个 越来越陡峭的边界峰和内部峰，而d g 方法正是求解含边界峰和内 部峰的解的一个有效的方法，所以本文中我们将用d g 方法处理弱 形式中乱所对应的项，同时结合插值系数有限元方法 2 8 来简化非 线性项，( u ) 的处理，以求解一维和二维区域上半线性椭圆型方程的 多解。 1 2预备知识：插值系数有限元 有限元法完全可用于解非线性问题，且具有解线性问题相同的 好精度。但数值求解非线性方程组时由于n e w t o n 法的迭代过程中要 多次对不同的值计算切矩阵，其工作量巨大。 湖南师范大学2 0 0 6 届硕士学位论文 讨论一种特定的半线性椭圆边值问题 一d i ( a i j ( x ) d j u ) + f ( u ) = 9 ( z ) i nq ，u = 0o na q ， ( 1 2 ) q 是具有l i p s c h i t z 边界勰的二维有界区域。它相应的弱形式为 q ( u ，u ) = a ( u ，u ) 十( 厂( 札) ，v ) = 0 ，v v s o ， 其中双线型a ( u ，u ) = 厶a i j ( z ) n 乱功u 有界且岛强制，s o = u h 1 ( q ) ，u = 0o na q - 。 将区域q 作剖分j ，q 被剖分成有限个小单元丁。设玩= x j f f 表 示所有内节点的集合。记瓯cs o 为n 次有限元子空间， 屿( z ) ) r 为瓯 中的一组基，则( 1 2 ) 有限元解u 可以表示为u h ( x ) = 丝。g j ( z ) s ，u s = u h ( z f ) ，且满足 a ( u h ，v ) + ( f ( u h ) ，v ) = 0 ，v v s h ( 1 3 ) 取v = n i ，i = 1 ，2 ，m ，代入u ，( 1 3 ) 可转化成一非线性方程组 mm a ( n j ，n , ) u j ( ，( n j ( x ) u j ) ，g i ) = o ，i = 1 2 m ( 1 4 ) 5 = 15 = 1 我们一般采用牛顿法解此非线性方程组。众所周知，牛顿法中最让 人关心的是它的j a c o b i 矩阵。( 1 4 ) 的j a c o b i 矩阵为 彳 j = a ( n j ，m ) 一，7 ( 肌巩) 码，g d m ( 1 5 ) k = l 可以看出，( 1 5 ) 中的第二项是个很大的矩阵，它在牛顿法的每一步 迭代过程中都将重复计算，因此会耗费大量的时间和代价。 现在我们介绍解决( 1 2 ) 的插值系数有限元 2 8 。设厶，( 仳h ) = 丝，( z ) 厂( u s ) 为，( u ) 的分片插值多项式。用厶厂( u ) 代替，( 札 ) 代 入到( 1 3 ) 中，插值系数有限元解札 仍记为u h ( x ) = 丝。屿( z ) ，这 时可得到一个新的有限元方程 a ( u h ，v ) + ( i h f ( u h ) ，v ) = 0 ，v v 瓯 用间断g a l e r k i n 有限元求解线性和非线性微分方程 7 将u 九， h f ( u ) 的表达式代入上式，最终将得到一个非线性代数方程组 m ( b 屿+ m i j f ( u j ) ) = 0 ，i = l 2 ，m ， ( 1 6 ) j = l 其中刚度矩阵= l ( 屿，f ) ，质量矩阵m 玎= ( g j ，m ) 均为常数矩 阵，在牛顿迭代过程中仅需计算一次。而( 1 6 ) 的j a c o b i 矩阵为 j = + m i j f ( u j - ) ) m m 当忌i j 和m 巧给定后，j a c o b i 矩阵只需将m 玎与厂n ，) 简单相乘即可， 因此解( 1 6 ) 与解( 1 4 ) 相比计算量大大减少。 很长时期里，人们不知道插值系数有限元法的逼近精度到底如 何。首先z l a m a l 3 1 】讨论热传导问题及线性元，假定有限元的最大模 m 觚i u ( z ) i 关于h 一致有界，他证明了最佳误差估计一札h 1 | = o ( h z ) 。 以后l a r s s o n t h o m e e ( 1 9 8 8 ) 不作此假定证明了i i u u _ i i l = o ( ) 。不 久c h e n l a r s s o n z h a n g 1 0 利用分片均匀三角剖分上椭圆有限元 梯度的超收敛估计，证明了几乎最佳阶i i 让一t t h i i = o ( h 。i n h ) 。文【2 8 】 中谢资清、陈传淼证明了，对拟一致剖分网格，m 次插值系数有限元 法仍有最佳阶误差il u u h li = o ( h 仇+ ) 。文【2 3 中熊之光、陈传淼进一 步证明了，上述问题有限元的经典结果对插值系数有限元法仍然保 持。 湖南师范大学2 0 0 6 届硕士学位论文 第二章用间断有限元方法求解二维区域上奇异摄动问题 2 1二维奇异摄动问题 我们讨论如下二维对流扩散方程的d i r i c h l e t 问题 j e a u + 卢v u - t - c u = ，i n q = ( o ，1 ) ( o ，1 f 2 1 ) 【乱= u 0 o na q j 其中万= 洒，仍) ( o l ，a ) ( 0 ，o ) ，并且卢。，仍均为常数，是一个小参 数。 设q = v u ，那么( 2 1 ) 可以改写为 | ，g = 乳 i nq 一v q + p v u + c u = f i nq ( 2 2 ) 【u = u 。 。n 独 将区域q 在z 方向上作剖分：0 = x i l z z ； z 。+ = 1 ， 步长h = x i + 一x i _ ，在y 方向上作剖分：0 = 可 可； 可g 可卅 = 1 ，步长b = + 一y j 一 记单元= k 一；，+ 一；，+ 】。 记q = u k = l , , m x a k k 。定义u + 和u 一分别为u 在单元内部、外部的 极限值。设e 为两个相邻单元的公共边，则对于此边上的任意一点 z ，记矿和礼一为该点的正、负法向向量。引入节点处的跃度为【u 】- u + n + + u 一几一， q 】- q + 矿+ q 一n 一，平均为【u ) = u + + 广u - ， q ) = 里2 。定 义有限元空间坻，y u 为 m n = q l 2 ( q ) 2 ：q l k s ( 后) 2 ，v k t ) ， v r = 札l 2 ( q ) ：乱i 七s ( 七) ，y k t ) ， 其中s ( 忌) 为单元k 上每个变量次数不超过p 的多项式空间q p ( 南) 。 分别用试探函数r l 2 ( q ) 2 ，v l 2 ( q ) ，乘以方程( 2 2 ) 中的前两个 方程，并在每个单元上分部积分，得到( 2 2 ) 的弱形式 | q r d x = 一| ，润r d x + | 埘n k d s 、 用间断g a l e r k i n 有限元求解线性和非线性微分方程 9 fq v v e x - j ( 让厅v v d x + f k c u v d x - e f o o 七 j n k d s + z 。u u 矽n k d s = f k f u d z ， 其中仡t 是边界上的外法向向量。 基于以上弱形式，我们引入相应的d g 方法。设有限元解q m n ，似满足上述弱形式，则对于任意的r m n ，u ，有 j ( q g r d x = - zu n v r d x + z 七咖舢， ( 2 3 ) v v d x - - 卜j n 如+ kc u n v d x - - e 8 n k d s + l o o k 雹n v 石n k d s = f kf 诎， ( 2 4 ) 式中的 ， 为数值通量，即让，q 在单元边界上的近似值。 本文中我们采用的通量表示如下 轴= 协- c l - m - c l z 【口( 2 5 ) 【 = 札) - i - c 1 2 u 、7 如果e 是q 的外边界，则通量就是 j = q + 一c l l ( l l u o ) n + 1 一u n = u o 本章的数值计算中我们取1 2 1 = 1 ，c z = 划2 ， = ( 1 ，1 ) ，且采用双线 性有限元基函数，故单元k 上的有限元解可以表示为 44 札= 乱磬妒，q n = 毋妒乎) ， n = l n = i 其中让磬n = 1 ，4 ) 为单元上各插值点的值，妒磬( n = 1 ，4 ) 为 q 1 ( 尼) 中的一组基。 将札、q 的表达式代入( 2 3 ) ，可以局部的解出q 的系数，再 将q g 的表达式代入( 2 4 ) ，最后将得到一个关于乱的系数趾譬( 咒= 1 ，2 ，3 ，4 ，k = 1 ，m 礼) 的线性代数方程组 g u = p 1 0 湖南师范大学2 0 0 6 届硕士学位论文 其中 g = c ( 1 ) b ( 2 ) a ( 3 ) d ( 1 ) c ( 2 ) b ( 3 ) e ( 1 ) d ( 2 ) c ( 3 ) e ( 2 ) d ( 3 )e ( 3 ) - a ( m n 一2 ) b ( ”m 一2 ) c ( ”m 一2 ) d ( m r 一2 ) e ( ，n 礼一2 ) a ( m n 一1 ) b ( m n 一1 ) c ( m n 一1 ) d ( , n n 一1 ) a ( ”t n )b ( m n )c ( m , o a ( m ，b ( 舢，c ( m ，d ( m ，e ( ) 是与对流矩阵尬和质量矩阵趔有关的 4 4 矩阵，且 尬= 凫) ，) 4 4 ，硝钟= 川m ，妒渤。4 ，i ，j = 1 ，4 ， f ：( f ( 1 ) t ，f ( 2 ) 丁，f ( m n ) 丁】，f ( ：【( 厂，v ! y ) 1 1 x 4 ，礼= 1 ，4 u ： u 咒u 孵，钆t 一 2 2数值例子 我们考虑用上节介绍的双线性间断有限元方法来解决二维的奇 异摄动问题。考虑方程 k 寰舡乳棚叫篇，1 ) ( 0 1 ) ( 2 6 ) 其中选取p ( z ，可) = ( 1 ，1 ) ，c ( x ，) = 0 和f ( x ，可) = + y ) ( 1 一e 一字e 一字) + ( z 一夕) ( e 一孚一e 一孚) ，其真解为“( z ，可) = x y ( i e 一下1 - - x ) ( 1 一e 一宰) ( 1 ) 均匀网格下求解 用间断g a l e r k i n 有限元求解线性和非线性微分方程 1 1 图2 。l ：均匀网格下的真解n =图2 2 ：均匀网格下的有限元解 1 6 e = 1 0 1 u ，n = 1 6 ，= 1 0 1 e = 1 = 0 1e = 0 0 1 n 0 u e l l 。n o r d e r i | 札一训n o r d e r l i u 一也i i 。 o r d e r 44 9 0 e - 0 4一 1 0 4 e - 0 1 6 3 1 争0 1 84 3 7 e - 0 53 4 82 1 0 e 一0 2 2 3 17 7 8 e - 0 10 3 0 1 6 4 3 5 e 0 63 3 22 1 6 e 一0 33 2 8 4 7 9 e - o l0 7 0 3 24 7 2 e - 0 7 3 2 02 4 8 e - 0 43 1 21 7 5 e 一0 11 4 6 6 4 5 4 5 e - 0 83 1 22 9 6 e - 0 53 0 73 1 0 e 一0 22 4 9 1 2 8 6 5 2 e 一0 93 0 63 5 9 e - 0 63 0 43 2 5 e 0 33 2 5 e = 1= 0 1= 0 0 1 n 岛 o r d e r a 口 o r d e r 乏q o r d e r 4 1 0 9 e 一0 1 8 8 3 e 一0 1 9 5 0 e 一0 0 8 2 1 4 e - 0 22 3 63 2 1 e 一0 11 4 61 6 3 e + 0 1一o 7 8 1 64 8 4 e 一0 32 1 57 1 7 e 0 2 2 1 61 4 4 e + 0 10 1 7 3 2 1 1 7 争0 32 0 51 9 6 争0 21 8 72 5 0 e - 0 02 5 3 6 4 2 8 8 e - 0 42 0 25 2 8 争0 31 8 92 5 7 e 一0 13 2 8 1 2 8 7 1 6 e - 0 52 0 01 3 8 e 一0 31 9 47 0 7 e - 0 21 8 6 1 2 湖南师范大学2 0 0 6 届硕士学位论文 图2 3 ：均匀网格下真解的z 方向 导数，n = 3 2 ，e = 1 0 1 图2 5 ：均匀网格下真解的y 方向 导数，n = 3 2 ，= 1 0 1 图2 4 ：均匀网格下有限元解的x 方向导数，n = 3 2 ，= 1 0 - 1 图2 6 ：均匀网格下有限元解的y 方向导数，n = 3 2 ，= 1 0 - 1 用间断g a l e r k i n 有限元求解线性和非线性微分方程 1 3 图2 7 ：均匀网格下的有限元解误 差，n = 1 6 ，= 1 0 2 图2 8 ：均匀网格下的有限元解误 差，n = i 6 ，= 1 0 3 表2 - 1 和表2 - 2 分别列出了间断有限元解和其导数在一致均匀网 格下的最大模误差，其中分别取= 1 ，= 0 1 和= o 0 1 ，x , y 方 向上的剖分数均为n 。计算结果表明有限元解在x , y 方向上的导数误 差值相同。从表中我们得到当用一次元去逼近时，有限元解可达到3 阶收敛，导数相应的可达到2 阶收敛。在我们的数值实验中，当e 非 常小时，导数在边界层内变化非常快，所以表中的导数误差是相对 误差，即瞥。 图2 1 和图2 2 分别是当= 1 6 ，= 1 0 - 1 时的真解和有限元解， 从图中可看到，间断有限元解很好地逼近了真解。 图2 3 和图2 5 分别是真解在x 方向和y 方向的导数图，与之对 应，图2 4 和图2 6 分别是有限元解在z 方向和y 方向的导数分布图。 同有限元解一样，导数即使在边界层附近也没有产生振荡现象。在 n = 1 6 剖分的均匀网格下，我们作出了当e = 1 0 - 2 和e = 1 0 - 3 时的 有限元解的误差分布图，见图2 7 和图2 8 。从图中我们可以明显地 看出误差主要出现在边界层内，特别是在两个方向的边界层的交汇 处，误差最大处形成了一个尖形，但误差也没有任何振荡现象，这是 相对连续元的间断有限元的最大优势之一。 湖南师范大学2 0 0 6 届硕士学位论文 ( 2 ) 第一种s h i s h i k i n 网格下求解，7 - = ( 2 p + 1 ) el n ( n + 1 ) 下面我们转到局部加密网格，即s h i s h i k i n 网格上求解。所谓s h i s h i k i n 网格，即选取一个0 丁 1 2 ，使得丁与边界层的宽度几乎一致。 在我们的数值实例中，边界层宽度= o ( ) 在z 方向和y 方向同时将 ( 0 ，1 一下) 和( 1 一下，1 ) 两个区间等分为g 2 个小区间，因此所得到单元 数为。个 = 1 0 2e = 1 0 4e = 1 0 6 n i i u a l l 。n o r d e r | i u 一刮n o r d e r 一训。n o r d e r 41 0 1 e - 0 1 9 9 2 e - 0 2 9 9 1 e 一0 2 8 3 6 4 e 一0 21 4 83 4 9 e 一0 21 5 13 4 8 e - 0 2 1 5 1 1 68 7 3 e - 0 32 0 68 1 l e 0 32 1 08 1 0 e 一0 32 1 0 3 21 9 5 e 0 32 1 61 8 6 e 一0 32 1 21 8 6 e 0 3 2 1 2 6 44 1 2 e 一0 42 2 44 0 0 e 一0 42 2 23 9 9 e 0 42 2 2 1 2 88 3 4 e 一0 52 3 08 1 l e 一0 52 3 08 1 l e 0 52 3 0 e ：1 0 2：1 0 4e ：1 0 6 n 色口 o r d e r a 口 o r d e r 芭q o r d e r 4 8 0 4 e - 0 1 6 6 7 e 一0 1 6 6 6 e - 0 1 82 9 0 争0 11 4 72 6 4 e 一0 11 3 42 。6 3 e 一0 11 3 4 1 61 2 2 e 一0 l1 2 51 1 5 e - 0 l1 2 01 1 5 e - 0 11 2 0 3 2 5 2 2 e 一0 21 2 24 9 9 e 0 21 2 04 9 8 争0 21 2 0 6 42 1 0 e 一0 21 3 12 0 2 e 一0 21 3 12 0 2 e 一0 21 3 1 1 2 87 7 9 e 一0 31 4 37 4 9 e - 0 31 4 37 4 9 e 一0 31 4 3 首先我们选取7 - = ( 2 p + 1 ) l n ( n + 1 ) 。表2 - 3 和表2 4 分别列出了 有限元解及其导数在= 1 0 ，= 1 0 一和= 1 0 6 时的最大模误差和 收敛阶。从表2 - 3 和表2 4 可以得到，解及其导数都有一致的超收敛 误差估计 蚓j 。c ( 掣) 卿 蚓i 。c ( 警) 州 一星盥鲤鱼兰竖：兰! 呈查墼垂查蟹垡堡塑韭垡堡丝坌查堡 ：望： _ _ - _ - - - - _ _ - _ l _ _ - - - - _ _ _ - - _ _ _ _ - _ _ - _ _ _ - - _ _ _ _ - - - _ _ - - _ _ _ _ - - - _ _ - _ - _ _ _ - - - _ _ 一一一一一 图2 。9 ：s h i s h i k i n 网格下的有限元 解误差，r = ( 2 p + l 弦l n ( + 1 ) ，n 一 1 6 = l o 一2 图2 i i ：s h i s h h l d n 网格下有限元 解z 方向的导数误差，丁= ( 印+ 1 ) s l n ( n + 1 ) ，n = 1 6 ，e = 1 0 2 图2 1 0 ：s h i s h i k i n 网格下的有限元 解误差，r = ( 印+ 1 ) l n ( + 1 ) ，n = 1 6 ，= 1 0 6 图2 1 2 ：s h i s h i l d n 网格下有限元 解。方向的导数误差，f = ( 2 p + 1 ) l n ( n + 1 ) ，n = 1 6 ，= 1 0 6 1 6 湖南师范大学2 0 0 6 届硕士学位论文 图2 1 3 ：s h i s h i k i n 网格下有限元 解y 方向的导数误差，丁= ( 2 p + 1 ) e l n ( n + 1 ) ，n = 1 6 ，e = 1 0 2 图2 1 4 ：s h i s h i l 【i n 网格下有限元 解y 方向的导数误差，丁= ( 2 p + 1 ) e l n ( n + 1 ) ，n = 1 6 ，e = 1 0 6 这里c 不依赖于奇异系数。这与用间断有限元解1 维奇异摄动问题 的结果一致。 图2 9 至图2 1 4 都是取7 _ = ( 2 p + 1 ) l n ( n + 1 ) 的s h i s h i k i n 网格情 形，其中图2 9 和图2 1 0 分别是取= 1 0 ，e = 1 0 川的有限元解的误 差图。从图中可以看出同均匀网格的情形一样，误差主要分布在边 界层，当越小时，边界层也越窄，但误差没有增大，有一致的收敛 性。 图2 1 0 和图2 1 1 是解在z 方向导数的误差情形。误差主要出现 在与y 轴平行的边界层内，类比对y 方向导数误差主要出现在与z 轴方向平行的边界层内，具体见图2 1 3 和图2 1 4 。因方程( 2 6 ) 的真解 在x ，y 方向具有对称性，所以在z 和y 方向导数误差也是对称的。 ( 3 ) 第二种s h i s h i k i n 网格下求解，7 = 一( 2 p + 1 ) 1 n 下面选取另一种类型s h i s h i k i n 网格即7 = 一( 2 p + 1 ) ei n 的情形。 同第一种s h i s h i k i n 网格一样，表2 5 和表2 6 分别列出了有限元解和 其导数在= 1 0 ，= 1 0 3 和= 1 0 一时的最大模误差和收敛阶，其 两个方向上的导数误差相同。当丁= 一( 2 p + 1 ) i n 时也有一致的超收 用间断g a l e r k i n 有限元求解线性和非线性微分方程 1 7 且似乎有如下的误差估计 l i a 。c n 一( p + 2 ) 驯。c n 一( p + 1 ) 图2 1 5 至图2 2 0 都是取7 - = 一( 2 p + 1 ) e l n c 的s h i s h i k i n 网格情 形，与第一种s h i s h i k i n 网格情形相对应，图2 1 5 和图2 1 6 分别是取 e = 1 0 ，e = 1 0 - 4 的有限元解误差图。图2 1 7 和图2 1 8 是分别取 e = 1 0 ，e = 1 0 - 3 时的z 方向的导数误差图，y 方向的导数误差图见 图2 1 9 和图2 2 0 。与取丁= ( 2 p + 1 ) e l n ( n + 1 ) 的s h i s h i k i n 网格相比， 第二种s h i s h i k i n 网格即取丁= 一( 2 p + 1 ) el ne 时所得到的误差没有第 一种稳定，但具有相同的收敛阶。 = 1 0 2e = 1 0 3 = 1 0 4 n i i u e l l 。n o r d e r i i u e l l 。n o r d e r i i u e l l 。n o r d e r 45 2 9 e 0 1 7 2 l e - 0 1 8 4 8 e - 0 1 8 2 1 0 e - 0 11 3 33 8 5 e - 0 l0 9 15 2 5 e 一0 10 6 9 1 64 1 7 e - 0 22 3 3 1 1 6 e 一0 11 7 32 0 8 e - 0 l1 3 4 3 24 3 3 争0 33 2 71 6 2 e - 0 22 8 34 0 0 e - 0 22 3 8 6 4 5 5 6 e - 0 42 9 61 8 0 e 0 33 1 74 1 0 争0 33 2 9 1 2 87 0 9 争0 52 9 72 2 9 e - 0 42 9 7 5 3 4 e - 0 42 9 4 ：1 0 2 e = 1 0 3= 1 0 4 n 岛 o r d e r 岛 o r d e r 岛 o r d e r 4 1 5 2 e + 0 1 1 5 7 e + 0 2 1 5 1 e + 0 3 84 3 8 e + 0 01 8 0 2 6 9 e + 0 12 5 41 6 0 e + 0 23 2 4 1 6 3 2 2 e 一0 13 7 68 6 2 e - 0 14 9 62 7 2 e 一0 05 8 8 3 2 8 4 2 e - 0 21 9 41 6 5 e - 0 12 3 82 9 2 e 一0 13 2 2 6 42 5 0 融0 21 7 5 4 9 0 e - 0 21 7 58 0 l e 0 21 8 6 1 2 8 7 0 5 e 一0 31 8 31 4 4 e - 0 21 7 72 4 0 e - 0 2 1 7 4 湖南师范大学2 0 0 6 届硕士学位论文 图2 1 5 ：s h i s h i k i n 网格下的有限 元解误差，丁= 一( 2 p + 1 ) v l n ，n = 1 6 e = 1 0 2 图2 1 7 ：s h i s h i m n 网格下的有限 元解z 方向导数误差，r = 一( 2 p + 1 ) z l n e ，n = 1 6 ，= 1 0 2 图2 1 6 ：s h i s h n i n 网格下的有限 元解误差，丁= 一( 2 p + 1 ) i n ，n =