（计算数学专业论文）基于多水平残量空间方法.pdf

摘要 偏微分方程的数值解法作为科学工程计算的核心问题，而求解大规模的线 性方程组义是偏微分方程数值解法研究的一个核心课题之一，在当前的求解大 规模的线性方程组的方法中，多重网格技术占据着举足轻重的地位。对于特殊 问题构造特殊算法一直是计算数学家和学者的研究热点，本篇硕士论文研究和 介绍了基于多水平残量空间的求解大规模代数方程组的一种新提出的方法。 本文共分三章第一章首先介绍了多重网格技术产牛的背景，然后通过一 个简单的模型问题介绍二重网格技术的主要思想。简单的说明了多重网格技术 里涉及到的插值算子、限制算子以及粗网格上的离散算子构造问题。在对二重 网格技术的介绍基础之后，给出了一般形式的几何多重网格方法，例如v 循 环，w 循环以及完全多重网格( f m g ) 等几种形式的算法描述。 第二章主要介绍了基于多水平残量空间的方法。首先介绍了多水平残量空 间这一概念，然后对我们的残量空间利用g r a m s c h m i d t 止交化方法，得到一 组相互共轭的基。然后我们把p e t r o v g a l e r k i n 方法( 投影方法) 相结合构造出 基于多水平残量空间的p e t r o v g a l e r k i n 方法，同时我们把这种方法从代数方程 组的系数矩阵为对称正定矩阵形式推广到系数矩阵可逆但不要求对称正定的情 况，而且还对残量空间进行了扩充，得到扩充形式的算法。本章最后一部分介 绍了多重网格共轭梯度法，而且给出了另外不同的部分正交化算法。 第三章中我们通过几个数值算例介绍了基于多水平残量空间方法在以及其 扩充形式的算法求解偏微分方程数值解中的应用，我们主要把这些方法运用在 p o s s i o n 方程以及对流扩散方程上面并与其他的一些算法进行比较，通过这些算 例来说明我们介绍的多水平残量空间方法在求解偏微分方程数值解中的应用， 通过比较会发现我们的扩充形式多水平残量空间算法以及相应的改进技术对某 些特殊问题还是比较有效的。 关键词：几何多重网格方法，多水平残量空间，扩充多水平残量空间算 法，p e t r o v - g a l e r k i n 条件 a b s t r a c t n u m e r i c a lm e t h o do fp d ei st h ec o r eo fs c i e n t i f i cc o m p u t i n g ，a tt h es a m e t i m eh o wt os o l v el a r g e - s c a l el i n e a rs y s t e mi so n eo ft h em o s ti m p o r t a n tt h i n g s i nt h er e s e a r c ho fn u m e r i c a lm e t h o df o rp d e t h e r ea r em a n ye f f i c i e n tm e t h o d s t os o l v el a r g e - s c a l el i n e a rs y s t e m ，a m o n gt h e s em e t h o d s ，m u l t i g r i dm e t h o di so n e o ft h em o s te f f e c t i v ea n ds c a l a b l em e t h o d s t oc o n s t r u c ta ne f f i c i e n tm e t h o df o r s p e c i a lp r o b l e m si sp o p u l a ra m o n ge x p e r t i s em a t h e m a t i c i a no fs c i e n t i f i cc o i n - p u t i n g i nt h i sp a p e r w ei n t r o d u c ean e wt e c h n o l o g yo fs o l v i n gl a r g e - s c a l el i n e a r s y s t e m s t h ep a p e ri s d i 、i d e di n t ot h r e ep a r t s 、v ei n t r o d u c et h eb a c k g r o u n do f m u l t i - g r i dm e t h o di nt h ef i r s tc h a p t e r ，a n dt h e nu s eas i m p l em o d e lp r o b l e mt o i n t r o d u c et h eb a s i ci d e ao ft w o - g r i dm e t h o d j 胎a l s oi n t r o d u c eh o wt oc o n s t r u c t a n dc h o o s et h ei n t e r p o l a t eo p e r a t o r ，r e s t r i c t i o no p e r a t o ra n dd i s c r e t eo p e r a t o r o nc o a r s e rg r i d a f t e rc o m p l e t i n gt h e s et h i n g s ，w ei n t r o d u c em o r eg e n e r a lm u l t i - g r i dm e t h o da n dt h e i rd e s c r i p t i o n ，s u c ha sv - c y c l e ，w c y c l ea n df u l lm u l t i g r i d m e t h o d i nc h a p t e rt w o ，w ef i r s ti n t r o d u c et h eb a s i cc o n c e p to fm u l t i l e v e lr e s i d u a l s p a c ea n dt h eu s es c h m i d to r t h o g o n a l i z a t i o nm e t h o dt oo b t a i nag r o u po fo r t h o g o n a lb a s e t h e nw ec o m b i n et h em u l t i l e v e lr e s i d u a ls p a c ea n dp r o j e c t i o nm e t h o d t oc o n s t r u c tam e t h o dc a l l e dm u l t i g r i d p e t r o v e - g a l e r k i nm e t h o d t h e nw ee x - t e n dt h i sm e t h o df r o mt h ec a s ew h i c hi st h a tt h em a t r i xo fs y s t e mi ss p dt om o r e g e n e r a lc a s et h a ti st h em a t r i xh a si n v e r s ee v e ni fi ti sn o ts p d t h ew ee x p a n d t h eb a s eo ft h es p a c ea n do b t a i nt h ee x p a n d e da l g o r i t h mo ft h i sm e t h o d w ba l s o i n t r o d u c et h eb a s i ci d e aa n dc o n c e p to ft h em u l t i g r i dc o n ju g a t em e t h o da n dg i v e a n o t h e ra l g o r i t h mo fp a r t i a lo r t h o g o n a l i z a t i o na tl a s t i nc h a p t e rt h r e e ，w ea p p l yt h em e t h o d si n t r o d u c e di nt h el a s tp a r tt os e v e r a l p d es u c ha sp o i s s o ne q u a t i o n c o n v e c t i o n - d i f f u s i o ne q u a t i o na n de q u a t i o nw i t h l a r g ej u m p i n gc o e f f i c i e n t s t h r o u g ht h e s en u m e r i c a le x a m p l ea n dc o m p a r i s o nw i t h o t h e rm e t h o d sw ew i l lf i n dt h a tt h em e t h o d sb a s e do nm u l t i l e v e lr e s i d u a ls p a c e a b s t r a c ti u i sm o r ee f f i c i e n tf o rs o m es p e c i a lp r o b l e m s k e y w o r d s ：g e o m e t r i cm u l t i g r i d ，e x p a n d e dm u l t i g r i d p e t r o v g a l e r k i nm e t h o d ， m u l t i l e v e lr e s i d u a ls p a c e ，p e t r o v e - g a l e r k i nc o n d i t i o n 第一章几何多重网格综述 1 1引言 在日常生活中，当我们打开电视看天气预报节目的时候，看到预报员对着 卫星云图给我们进行预报天气的时候；当我们看到神六升空的时候，了解到嫦 娥号卫星在太空轨道飞行的时候；这些事件或者现象中涉及到求解复杂的偏微 分方程( 组) ，例如空气动力学、流体力学、三维立体仿真、热传导等等，求解 这些复杂的方程( 组) 就是我们通常所说的科学工程计算。毫无疑问，椭圆和双 曲偏微分方稗模型在工程计算巾占据着核心地位，很多情况下，这些问题的计 算超过了最强大计算机的计算性能或者需要的时间太大以致_ 丁在仪器设备的设 计中不能包含高等的数学模型，使得优化设计更加困难，多重网格方法就是在 算法效率上有晕要的优点，它能够在求解n 个未知量的代数方程组时，只需要 o ( n ) 的运算和存储，这个方法不仪仪对某类问题，而且对很多类问题都适合。 多重网格技术真正的出现是在1 9 6 4 年由f e d o r e n k o 给出的在一个正方形 区域上的求解p o s s i o n 方程的五点有限差分方法时提到的，这个方法在1 9 6 6 年 由b a c h v a l o v 在q = ( 0 ，1 ) ( 0 ，1 ) 求解一般性椭圆方程： l u ：= 一( n 。z u n ) p + ( 6 q 札) 。+ c a = s 利用中心差分离散时加以运用。但是这段时间内，关于多重网格技术在理论上 和实践上都没有得到重大的发展，没有引起足够的重视。第一个将多重网格方 法运用到实践中的先行者是以色列的数学家b r a n d t ，他在1 9 7 3 年应用多重网 格方法求解实际问题，在1 9 7 7 年他又写了一篇论文，这篇文章中他清楚的给出 了多重网格技术的主要原则和实际应用中的效果，正因为这篇文章，使得了多 重网格方法得到了广泛的关注和快速的发展。 h a c k b u s c h 在1 9 7 6 年独立的发现了多重网格方法，而且给出了多重网格 方法的坚实的数学理论，也给出了可靠的算法，他对多重网格技术的完善主 要是在1 9 7 8 、1 9 8 0 、1 9 8 1 年这三年期间完善的。f r e d e r i c k s o n 给出的一个关于 p o s s i o n 方程求解的方法使得w e s s e l i n g 在1 9 7 7 年给出了一个关于涡流形式的 n a v i e r s t o c k s 方程的有效求解方法。 第一章几何多重网格综述2 在多重网格技术的发展过程中，形成了两个分支：几何多重网格和代数多 重网格，几何多重网格主要是在粗网格粗化或者矫正的过程中依赖于已经给出 的网格或者剖分，形式比较固定；或者说设计几何多重网格方法的时候，限制 算子与插值算子是同定的，更多的专注在迭代格式的选择和设计。关于多重网 格更加详细的介绍和说明见【1 ，2 ，3 】，这些专著中对这些问题有专门的讨论和说 明。代数多重网格方法则纯粹是从代数观点出发，完全依赖于矩阵的本身信息， 在设计代数多重网格方法的过程中，我们通常选择固定而又简单的的迭代格式 例如g a u s s - s e i d e l 迭代或者j a c o b i 迭代，而主要专注于插值算子的构造，插值 算子的构造是代数多重网格方法的核心，关于此方面的内容可以见【7 ，3 ，8 ，1 4 1 。 下面我们将主要介绍几何多重网格的一些基本内容。 1 2几何多重网格方法 我们现在通过一个一维模型问题来介绍二重网格的主要思想。我们的模型 问题定义如下： u - 。l 。u ，：= 牡。：j 7 三丢：z ( 。，1 ) i c 1 1 ， 对于模型问题，我们采用中心差分的格式进行离散，采用的网格剖分为： 巧= j h ，j = 1 ，2 ，m ，x o = 0 ，x m + i = 1 ，h = 上m + i 离散后的方程组为： fl 呦：= 一丝生l 二挚= 乃，办：= ，( ) ，j - - - - 1 , 2 , , m ( 1 2 ) l u 02u m + 12 0 写成代数方程组的形式就是： 4 乱= b 第一章几何多重网格综述 3 其中a 为j 对角矩阵，札= 【u 1 ，u 2 ，u m + 1 】t ，b = b l ，5 2 ，b m + l 】r 1 a = 商 凡 一12o 一121 o 0 00 00 00 121 012 对于这个线性方程组，我们采用阻尼j a c o b i 迭代法求解，。司以得到如f 的迭代 格式： j ，吗【舛1 】= 譬办+ ( 啦。+ 掣。) 、母+ u = u 吗【p + 1 】+ ( 1 一u ) 够1 与之对应的迭代矩阵为： g j ，= u 称为松弛因子，当u = 1 时就是通常情况卜的j a c o b i 迭代法。我们记u 为方 程组的准确解，则迭代误差e u + l = u u 1 】满足如下的齐次线性方程组： e p + 1 】：g ，e mour 容易求得g ，的m 个特征值以及相应的特征向量，如下： 入七( g ，) = 1 2 _ vs i n 2 ( 警) ，1 k m v k = 【s i n ( k h 丌) ，s i n ( m k h t r ) r ，1 后m 显然这m 个线性无关的特征向量构成了m 维线性空间的一组基，所以我们可 以将初始误差e 【o 】表示为它们的线性组合，即 m e 【o 】= f 口七 j _ - _ 一 七= 1 由误差迭代方程，得 e u + 1 1 = 口+ 1 ) v k u 0 0 一 生2 1 2 u 一2 u 0 o 一 0 l u 型2 一 型2 n u u 丝2 一 一 一 一 1 u一 竺2 0 o l 第一章几何多重网格综述 4 这说明误差迭代矩阵的特征向量可以来表示误差分量，而对应的特征值表示误 差分量的增长系数。从特征值的表达式可以知道，仪当0 u 1 时，才能保证 p ( a j ，) 1 ，这说明只有当0 u 1 时，阻尼j a c o b i 迭代才收敛。 我们把误差分量饥，的分成两组：低频部分1 k m 2 和高频部分 m 2 k m 。由于k 表示半波数，可以知道k 越小则半波的个数就少，从 而就比较平滑；反之k 越大，则半波个数就大，半波就多从而在几何上来看就 比较振荡。因而低频和高频部分也分别成为光滑误差分量和振荡误差分量。下 图给出了m = 2 4 时几个不同频率的误差分量的波形，从几何上来看低频部分 ( k = 1 ，2 ，3 ，4 ) 是相对平滑的，而对于高频部分( k = 7 ，9 ，1 0 ，1 1 ) 是相对振荡的。 图1 1 ：误差分量v k ，及其衰减凶子九( g j ，) 现在我们来看看误差分量对应的衰减因子九( g ，) 随着频葺夏的变化曲线图 ( 从上至下依次对应于u = ，石1 ，；，1 ) 从图中的曲线可以知道，对于低频分量衰 减因子都不会很小，特别地当k = 1 时，误差分量a 1 ( g ，) = 1 2 u s i n 2 ( h 7 r 2 ) 1 一华相对最大，而对于高频分量，我们则可以通过对u 进行选择而加以控制 衰减囚子。我们选择u 使得入m 2 ( g j ，) = 一a m ( g j ，) ，可以得到= ；，此时有 i , 入m 2 ( g j ) i = i a m ( g j ，) i = 吾1 ，l 入k l 吾1 ，m z k m ，也就是说当u = ；时， 一次阻尼j a c o b i 迭代就可以使得高频误差分量衰减到初始误差分量的盲1 。这些 结果表明低频误差分量是阻尼j a c o b i 迭代收敛慢的原因，而高频误差分量则在 迭代中衰减很快。换句话说，随着迭代次数的增加，误差的光滑性也进一步增 强。对于这种现象，其他的迭代方法也是存在的。 假设我们的网格剖分为q m ，q l ，对应的网格步长分别为h 饥，h l ，h o 圈园匿区霸弱 第一章几何多重网格综述 5 而且网格步长满足关系k 一1 = 2 h m ( 在利用多重网格方法来求解偏微分方程数 值解中通常的网格步长也是这样取的，因为便于程序设计而且即使是满足其他 的关系也不会得到更好的结果) 。我们来看看相邻两层网格之间的误差分量的关 系，此时我们分别记u 乏2 f 以及 恐来表示细网格偶数节点编号误差分量和粗网 格误差分量。当1 k 警时有 u 乏巧= s i n ( 等) = s i n ( 镌) = 恐 这表明细网格上的低频误差分量变成了粗网格上的误差分量，有可能是光滑的 也可能是振荡的。但是从几何上来看，细网格上的低频分量在粗网格上来看比 在细网格上要振荡。从这个关系，就可以得到启示，可以将细网格上的光滑误 差分量限制到粗网格上，然后在粗网格上进行迭代处理就可以消除部分误差分 量，基于这种思想我们可以一层一层做下去直到消除所有的误差分量为止，最 后还需要将求解问题一层一层的返同，直到最初的细网恪上。这就是几何多重 网格的基本思想。 现在我们来介绍多晕网格技术中最基本的二重网格方法，双罩网格技术是 一切多重网格方法的基础或者说复杂的多重网格循环都是建立在双重网格方法 之上的。我们就模型问题来说明双重网格方法。假设存在对区间的两个嘲格剖 分q 1 ，q 2 ，对应的网格步长分别为h 1 ，h 2 而且满足h 2 = 2 h 1 ，那么对于模型问题 的双重网格方法描述如下： ( i ) 在细网格q 1 上用迭代法求解离散问题 j l h 。牡_ f l 。，j = 。，j1 j 慨 1 乱 1 ，o ：u h l ，胁= o 假设迭代了1 步，对应的近似解记为v h 。，j ，并计算相应的残量 r l j = 厶1 j l h l , j v h l j ，l j 尬 则迭代耽后的近似解逼近差分解的误差e 慰= u h 。歹一y h 。j 满足残量方程 l h l e h l = l l h l v h l ( i i ) 将细网格上的误差转移到粗网格q 2 上，也就是将残量“。限制到q 2 上，通常的做法之一是对残量，加权，例如 jr 2 j = j ( 7 _ l l ，2 j 一1 + 2 r l i l l ，2 j + r h l ，2 j + 1 ) ，l 歹m 2 ir h 2 ，02r h 2 ，m 2 + 120 第一章几何多重网格综述6 ( i i i ) 在粗网格上求解离散问题 j l h 2 e h 2 = r h 2 le k o 2 e 2 拖+ 150 ( i v ) 将粗网格上的误差插值或延拓到细网格上 ( v ) 对近似解做修正 l j 卜口i l l l j + e h l j ，l 歹a 矗 然后再以，作为初值在细网格上迭代沈次，从而义得到新的解，我们又回到 第( i ) 步，可以一直循环直到得到满足条件的解。从上面的二重网格算法来看， 我们有三个问题还没有给出明确的答案，它们分别是：插值算子、限制算子以 及粗网格上的离散算j rl ：。对于这三个问题，下面我们分别加以简单的说明。 对于一维问题，插值算子通常按照如一卜方式定义： ju 易+ - = ( 矽+ 噼。) 【呜= 嵋 ，0 j 一1 这也就是我们通常所说的分段线性插值。如果要求高精度的插值，那么我们也 可以选择精度更高的插值方式，例如二次插值以及三次插值。对于更加详细的 插值形式可以见【1 ，2 】。对于二维情况下的插值，通常按照如下方式定义： i 吃，巧= 黪 lu 袅+ 1 ，巧= 虱1 “v t 2 j h + 2 + h l ，j ) 壤，巧+ - = ( u 芍+ 呓+ 。) i 吃+ 1 ，巧+ l = j ( 嵋；+ 2 + h 1 j + 蝣+ 1 + 噼1 ，j + 1 ) 【，0 i ，歹5 詈一1 对于一维问题的限制算子，最为直接和明显的方式就是直接限制方式，因 为粗网格节点必是细网格上的节点，因而我们直接令粗网格节点上的取值就为 在细网格上的节点值：u = 咄；另外一个常用的限制算子就是完全加权算子， 定义如下： 矽= 去( 易一1 + 2 吃+ u 乞+ 1 ) ，1 j 等一l 埘 n 1 现在我们定义捅值算子如下： 露一1 ：舻_ 舻t ，i = z ，2 定义限制算子如下： 第二章基于多水平残量空间方法 1 0 鬈1 ：彤t _ 俨，e q = ( 难1 ) t ，i = “，2 定义对角占优对称正定矩阵： 舰：舻l _ r 啦，i = f ，1 在有限元离散中，我们通常取坛为质量矩阵，在实际应用中我们也经常取尬 为单位矩阵因为便于计算，而且有时候也能取得很好效果。现在我们进一步定 义下面的算子： g ：= 伍l f 一- 霹+ 1 ，歹 s 巧：= 毋1 ”i j j + 2 鬈，歹 1 ：砖- 1j - ad ；：i j 冗( + 呱一l = 呱+ i 联一l 如果要求砖上钟，v i ，j 虽然可以满足d 七l a d 扣1 ，但是在这种情况下前 露( 彤t ) ，i = 1 ，l 一般是彳i 会成