（计算数学专业论文）不规则区域上抛物问题的配置有限元法.pdf

山东大学硕士学位论文 不规则区域上抛物问题的配置有限元法 陈圆圆 ( 山东大学数学学院，济南，山东2 5 0 1 0 0 ) 中文摘要 配置法是近几十年发展起来的一种数值求解方法，它是以满足纯插值约 束条件的方式，寻求算子方程近似解的方法配置法不必计算数值积分，逼近 方程容易形成，计算简便而且收敛精度高，因此在数值求解椭圆型方程、抛物 型方程、双曲型方程中得到广泛应用，但是主要研究工作大都针对规则区域展 开 本文建立了解一类不规则区域上抛物方程初边值问题的配置有限元法 该方法结合了正交配置法的简单性和可变有限元的多样性在每个结点使用 了具有4 个自由度的双三次h e r m i t e 元 本文共分两章： 第一章简单介绍了f r i n d 和p i n d e r 在f 1 】中建立的一种解一类边值问题 的潜力巨大的数值方法该方法结合了正交配置法的简单性和可变有限元的 多样性，在每个结点使用了具有4 个自由度的双三次h e r m i t e 元其中使用 的亚参数变换允许最高精度配置点的准确定位以及不规则边界的唯一表示 本章分为三节：第一节是引言。简单介绍了配置法的发展概况；第二节详 细介绍了配置方程的推导该方法对一类不规则区域上的势能问题： f 貉+ 寄= 0 t = 缸co nr 1 【爱= 一q c 觎r 2 提出了一种配置有限元解法，并给出了与g a l e r k i n 有限元法计算效率的比较 第三节对这种方法的优势做了简单总结 山东大学硕士学位论文 第二章针对初边值问题展开讨论对抛物型方程应用第一章的方法。给出 了在扇环形区域上一类热传导方程初边值问题的相应的配置解法 本章分为四节，第一节是引言，简单介绍了抛物问题的配置法发展；第二 节给出了在扇环区域上的热传导问题t i 象( z ，耖，t ) 一a u ( x ，y ，t ) = 0 ，( z ，耖) q ，t ( o ，t 1 缸l r 。= 0 ，爱i r 。= g ， ( z ，y ) 锄，t ( o ，t l it i t ：o = u o ( x ，y ) ， ( z ，y ) q 并应用第一章介绍的方法得出配置方程；第三节把前面的结论推广到更一般 的情形，指出本文所介绍方法具有重要的应用价值；第四节指出了本文存在的 不足 关键词：不规则区域。抛物方程，双三次h e r m i t e 元，配置有限元法 山东大学硕士学位论文 ac o l l o c a t i o nf i n i t ee l e m e n t m e t h o df o rp a r a b o l i cp r o b l e m s i ni r r e g u l a rdo m a i n s y u a n y u a nc h e n ( s c h o o lo fm a t h e m a t i c s ，s h a n d o n gu n i v e r s i t y , j i n a n ，s h a n d o n g2 5 0 1 0 0 ，p r 。c h i n a ) a b s t r a c t t h ec o l l o c a t i o nm e t h o di san u m e r i c a lm e t h o dw h i c hh a sd e v e l o p e df o rs e v e r a l d e c a d e s i ti sam e t h o dw h i c hs e a r c hf o rt h ea p p r o x i m a t es o l u t i o n so ft h eo p e r a t o r f u n c t i o n sb ys a t i s f y i n gp u r ei n t e r p o l a t i o nc o n d i t i o n t h ec o l l o c a t i o nm e t h o dh a s m a n ya d v a n t a g e s ，f o re x a m p l e ，i tn e e d n tc a l c u l a t en u m e r i c a li n t e g r a l ；i tc a nf o r m t h ea p p p r o x i m a t ee q u a t i o ne a s i l y ；i t sc o m p u t ei ss i m p l ea n dc o n v e n i e n t ，a n di th a s h i g h - o r d e ra c c u r a c y h e n c e ，t h ec o l l o c a t i o ni sw i d e l yu s e df o rs o l i n ge l l i p t i ee q u a - t i o n s ，p r a r b o l i ee q u a t i o n sa n dh y p e r b o l i ce q u a t i o n s b u tm o s ts 乞u d yw a sr i s t r i c t e d i nr e g u l a rd o m a i n s i nt h i sp a p e r ，ac o l l o c a t i o nf i n i t ee l e m e n tm e t h o df o rs o l v i n gc e r t a i ni n i t i a l b o u n d a r yv a l u ep r o b l e m si ni r r e g u l a rd o m a i n s t h em e t h o dc o m b i n e st h es i m p l i c i t y o fo r t h o g o n a lc o l l o c a t i o nw i t ht h ev e r s a t i l i t yo fd e f o r m a b l ef i n i t ee l e m e n t s b i c u b i c h e r m i t ee l e m e n t sw i t hf o u rd e g r e e s - o f - f r e e d o mp e rn o d ea r eu s e d a n a l y s i so ft h i s p a p e rs h o wt h a tb yt a k i n ga d v a n t a g eo ft h eb o u n d a r yc o n d i t i o n s ，am i n i m u ms u m b e to fc o l l o c a t i o np o i n t sc a nb eu s e d t h i sa r t i c l ei sd i v i d e di n t ot w oc h a p t e r s c h a p t e r1i n t r o d u c e sap o t e n t i a l l yp o w e r f u ln u m e r i c a lm e t h o df o rs o l v i n go e r - t a i nb o u n d a r yv a l u ep r o b l e m sw h i c hi sd e v e l o p e di n 【1 】t h em e t h o dc o m b i n e st h e s i m p l i c i t yo fo r t h o g o n a lc o l l o c a t i o nw i t ht h ev e r s a t i l i t yo fd e f o r m a b l ef i n i t ee l e m e n t s b i c u b i ch e r m i t ee l e m e n t sw i t hf o u rd e g r e e s - o f - f r e e d o mp e rn o d ea r eu s e d as u b - i l l 山东大学硕士学位论文 p a r a m e t r i ct r a n s f o r m a t i o np e r m i t st h ep r e c i s ep o s i t i o n i n go ft h ec o l l o c a t i o np o i n t s f o rm a x i m u m a c c u r 警y 黟w e l la 8 & u n i q u er e p r e s e n t a t i o no fi r r e g u l a rb o u n d a r i e s t h i sc h a p t e ri sd e v i d e di n t ot h r e es e c t i o n s t h ef i r s ts e c t i o ni si n t r o d u c t i o n ， w h i c hi n t r o d u c e st h eg e n e r a ls i t u a t i o no ft h ec o l l o c a t i o nm e t h o ds i m p l y t h es e c o n d s e c t i o ng i v e st h ei n t r o d u c t i o no ft h ea r i t h m e t i c t h i sm e t h o ds e t su pac o l l o c a t i o n f i n i t em e t h o df o rp o t e n t i a lp r o b l e m si ni r r e g u l a rd o m a i n s ： f 器+ 务= o u5 o nf i i 赛= - q co nr 2 a n dg i v e st h ec o m p a r e 诵t ht h eg a l e r k i nf i n i t em e t h o di nt h ec o m p u t a t i o n a le f l i - c i e n c y t h et h i r ds e c t i o nc o n c l u d e st h em e t h o ds i m p l y c h a p t e r2a p p l i e st h em e t h o dt oap a r a b o l i ce q u a t i o n ，a n dg i v e st h ec o r r e s p o n d - i n gc o l l o c a t i o nm e t h o do fc e r t a i nh e a te x c h a n g i n gp r o b l e mi ni r r e g u l a rd o m a i n s t h i sc h a p t e ri sd e r i d e di n t of o u rs e c t i o n s t h ef i r s ts e c t i o ni st h ei n t r o d u c t i o n s w h i c hi n t r o d u c e st h ed e v e l o p m e n to ft h ec o l l o c a t i o nm e t h o df o rp a r a b o l i cp r o b l e m s s i m p l y , a n di ta l s og i v e st h em a i np r o b l e m so ft h i sa r t i c l e 够w e l l t h es e c o n ds e c t i o n g i v e sah e a te x c h a n g i n gp r o b l e mi ns e c t o r - c i r q u e - s h a p e dd o m a i n ： i 甓( z ，奢，t ) 一a u ( x ，y ，) = o ，( z ，s ，) f t ，t ( o ，t i 训r 。= 0 ，舞l r 。= 9 ， ( 霸y ) 舶，( o ，t i 【u k = o = 伽( z ，可) ， ( z ，耖) q u s i n gt h em e t h o di n t r o d u c e di nt h ec h a p t e rl ，w eg a i nt h ec o l l o c a t i o nm e t h o d t h et h i r ds e c t i o ne x t e n d st h em e t h o do ft h i sa r t i c l e ，a n di n d i c a t e st h a tt h em e t h o d h a sv e r yi m p o r t a n tv a l u e ；t h ef o u r t hs e c t i o ng i v e st h es h o r t a g eo ft h i sp a p e r k e y w o r d s ：i r r e g u l a rd o m a i n s ，t h ep a r a b o l i ce q u a t i o n s ，b i c u b i ch e r m i t ee l e - m e r i t s ，c o l l o c a t i o nf i n i t em e t h o d i v 原创性声明 本人郑重声明：所呈交的学位论文，是本人在导师的指导下，独 立进行研究所取得的成果。除文中已经注明引用的内容外，本论文不 包含任何其他个人或集体已经发表或撰写过的科研成果。对本文的研 究作出重要贡献的个人和集体，均已在文中以明确方式标明。本人完 全意识到本声明的法律责任由本人承担。 论文作者签名：盐园圃 日 关于学位论文使用授权的声明 本人完全了解山东大学有关保留、使用学位论文的规定，同 意学校保留或向国家有关部门或机构送交论文的复印件和电子 版，允许论文被查阅和借阅；本人授权山东大学可以将本学位论 文的全部或部分内容编入有关数据库进行检索，可以采用影印、 缩印或其他复制手段保存论文和汇编本学位论文。 ( 保密论文在解密后应遵守此规定) 论文作者签名：脚导师签名槛日 期：狸芝：竺兰 山东大学硕士学位论文 第一章不规则区域上势能问题的配置有限元法 本章主要回顾了f r i n d 和p i n d e r 在【1 】中建立的一种解一类边值问题的 潜力巨大的数值方法该方法结合了正交配置法的简单性和可变有限元的多 样性，在每个结点使用了具有4 个自由度的双三次h e r m i t e 元其中使用的 亚参数变换允许最高精度配置点的准确定位以及不规则边界的唯一表示 1 1 引言 数值分析的配置方法在近几十年有了很大发展它是以满足纯插值约束 条件的方式，寻求算子方程近似解的方法简单地说，就是对算子方程l u = ， 寻找形如 口( z ) = a j c j ( z ) j f f i l 的近似解，其中勺，j = l ，2 ，是待定系数，咖，j = 1 ，2 ，是线性无 关的基函数 我们选择n 个配置点戤，满足插值约束f i ( x ) = 札( z ) ，得到配置方程 即可求出系数n j 配置法具有无需计算数值积分，计算简便及收敛精度高等优点，因此在工 程技术和计算数学的许多领域得到广泛的应用 对椭圆，抛物、双曲等各种类型方程应用配置法( 【4 】，【5 1 ，【6 】) 都取得了很好 的结果我们已经知道配置点取l e g e n d r e 多项式的根时，即为正交配置法，可 得最高精度的配置格式有限元也可以引入到配置法中，形成有限元配置法。 但配置点的准确定位要求限制了区域形状的选择，因此大部分研究都针对矩 形区域展开( f 3 】，f 7 1 ) 本章指出，在曲h e r m i t e 四边形上，对椭圆方程可以得 到一种配置有限元格式，兼有有限元和配置法的优点 磊 = 0 戤 ，j = 以 咖 己 吩 m = 一 z ul 山东大学硕士学位论文 1 2 配置方程 考虑势能方程的边值问题： f 蠹+ 酽o - u = 。 ( 1 - 2 1 ) 。 让= 缸c o nr 1 i l 姿：一口c o nr 2 其中口c 是法向通量，r l + n 是区域a ( x ，y ) 的边界假设区域是可划分为如 图l ( a ) 所示的四边形元的形状不规则的曲h e r m i t e 四边形，图x ( a ) 中的四边 形元有成直角的隅角，而且可以按照随着元边界正交的曲线叽7 坐标系方便 地定义，也很容易映射到f ，t 7 坐标下的基本2 2 正方形( 图1 ( b ) ) 这个变换 除了使不规则区域的处理更容易，还允许正交配置点( 实际上就是图1 ( b ) 中 基本正方形的g a u s s 结点) 的准确定位 图1 ：等参四边形元；( a ) 曲边形，( b ) 基本正方形 我们选择每个结点处有4 个自由度的双三次h e r m i t e 元，并在每个元内 取4 个高斯结点定义在基本正方形( 图l ( b ) ) 上的检验函数取为 蛾川= 善旧时妨。簧+ 田1 哿+ 谚1 蒜】：( 1 2 j = l 一 山东大学硕士学位论文 其札，簧，筹，意是结点未髓咖3 0 0 钙j ，钙1 ，锡1 廷双三次h e m i t e 基 函数 利用变换t ( 1 2 3 a ) ( 1 2 3 b ) a 2 u 伊盯抛a 2 丁钆，拈升曲升、a 2 u曲曲a 2 乞 折a 7 - a 2 t i 一= = = 一一+ 一一+ i 一一一一l 一+ l 却必却抛。越却升。、必却。卸诺7 a 仃丹必却a 盯2 。箧却折2 ( 1 2 3 c ) 我们可以把检验函数方程( 1 2 2 ) 变换到仉7 坐标把上述关系应用到结点值 并注意 到，在元边界处有赛= 鬻= o ，我们可以用在整体叮，7 - 坐标系中定义的结点 未知量来表示检验函数，得到如下形式： 狮) = 妻咖螂等m 鲁十孵7 器】 ( 1 2 4 ) 其中变换基函数定义为 四= 妒 ( 1 2 5 a ) 蜉= 呓。簧十钙1 器 蝣= 谚1 面a 1 - j 丽a 2 r j 妒= 弼1 警鲁 兵甲的绪点导数另有讨论 这样，检验函数方程( 1 2 4 ) 就可以代入到微分方程( 1 2 2 ) 此时( 1 2 2 ) 是写在所有配置点上的，对于被划分为勺个元素的网格，就得到如下 n = 蛾勺个方程。 ；t c 第十第mc 箬十等，鲁+ c 萼a x + 等，鲁 + ( 磐+ 纂2o , r 2 ) 黑祧，2 ， ( 1 2 6 ) 沈2 a 3 ，7 曲升 v p r 吖 讥一升c蔷一升 丹一武打一却 + 十 让一盯 u 一盯孔一衍讥一曲 丝西丝却 i j j | 丝武丝却 、，、，、i， 的 钯 配 2 2 2 1 l l ，l，l，k 山东大学硬士学位论文 其中求和j 足对配置点所在元的四个结点进行，基函数导数在配置点i 上求 值把已知量移到右边。可得矩阵形式： c u = f ( 12 7 ) 其中c 是系数矩阵，u = ( u 。，等，鲁，嘉，让：，) t 包含除去边界值的所 有结点未知量，f 包含已知的边界信息 方程( 1 2 7 ) 就是我们要求的配置方程，可求出结点未知量 上述方法产生的配置矩阵是一个稀疏矩阵，其实际结构取决于配置点和 结点未知量的编号方式最紧凑的编号方式是对配置点按元连续编号，对未 知量按结点连续编号( 如图2 ( a ) 所示) ，这时产生的矩阵形式如图2 ( b ) 所示 在这种方式下，每一个元上与坐标变换相关的所有计算可以一次完成，但矩 阵的主对角线几乎全部由0 组成因此我们要对编号方式进行调整。采取的 方法是：结点未知量的编号方式不变，而把配置点都分配给最近的结点。并按 结点连续编号( 如图3 ( a ) 所示) 这样得到的配置矩阵如图3 ( b ) 所示。可以用 t h o m a s 型算法或高斯消去法等多种方法求解 剑世趔 堕到阻 劁 l l 1 2 2 32 4 3 53 6 4 7 4 8 91 02 l2 z3 33 44 54 6 司厩 圃团两两 碉衙 剑殴纠睦到逝 蚓睦 78 1 92 03 13 24 3 “ 56 1 71 82 93 04 14 2 3l 1 1 1习晒羽瞬碉雨 到巴倒也到 型幽 341 51 62 72 83 9 4 0 l2 1 31 4 2 52 6 3 73 8 司厅 司丽酮丽 习丽 4 图2 ：配置矩阵：( a ) 配置点的按元编号网格，( b ) 结果配置矩阵 山东大学硕士学位论文 剑虹型隆 型坠 纠盥 6 1 7 1 82 93 04 14 24 8 51 5 1 62 72 83 94 04 7 司两困厨 碉两 碉而 剑也 幽隆型逝 到鲢 4 1 31 42 5 2 63 7 3 84 6 3l l l z2 32 43 53 64 5 司匝习嘶2 ：3 5词丽 到怛劂也 到殴 则唑 2 9 1 02 l2 23 33 44 4 l 7 8 1 9 2 03 1 3 24 3 曰厅 司两 两匝硇匝 图3 ：配置矩阵：( a ) 配置点的按结点编号网格。( b ) 结果配置矩阵 1 3 结论 f r i n d 和p i n d e r 在1 1 】中还给出了配置有限元法与g a l e r k i n 有限元法 计算效率的比较。指出求解同样的势能问题时，本文介绍的方法的运算量比 g a l e r k i n 有限元法在公式表述部分少接近9 0 ，而在求解部分少5 0 同时 1 1 1 中给出的算例还显示了两f , - ) r 法所得解的精度相同，配置解甚至具有更高 精度 可见，配置有限元法作为边值问题解的方法有相当可观的前景由于可变 形元的使用，它实际上与传统的g a l e r k i n 有限元法一样有很大的灵活性，产 生的精度也至少相同。但计算量要少的多 因此【1 】中得出结论：对于要求c 1 连续解的势能问题，配置法具有很大 的优势；对( 严连续解，优势则相对较弱一般来说，当所使用远的自由度数 目较大时，配置有限元法要优于g a l e r k i n 方法 5 山东大学硕士学位论文 第二章不规则区域上热传导问题的配置有限元法 上一章指出，配置有限元法解椭圆边值问题的前景广阔本章则针对初边 值同题展开讨论对抛物型方程，我们应用第一章的方法，给出了在扇环形区 域上一类热传导方程初边值问题的相应的配置解法 2 1 引言 我们已经知道，配置法是以满足纯插值约束条件的方式。寻求算子方程近 似解的方法，具有不必计算数值积分，逼近方程容易形成，计算简便而且收敛 精度高等优点配置法对椭圆型、抛物型、双曲型等各种方程都可以适用，在 规则区域上可以和有限元法自由结合我们还知道，以l e g e n d r e 多项式的根 为配置点的配置法即正交配置法，可以得到最高阶精度但以上研究主要集中 在规则区域，尤其是矩形区域上 对于一般的抛物型方程， h 舢) 鲁- d ( 础嚣- 6 瓴扯，爱攀，o z 1 ，0 螂 1札(z，0)=，(z)，0z1 i 铭( o ，d = o o ( t ) ，u ( 1 ，t ) = g l ( t ) ，0 t t 一般有时闻连续和时间离散两种配置格式( 1 2 1 ) 所谓时间连续的配置方法。 就是寻找一个关于时间可微的映射u ：( 0 ，7 1 _ a 1c r ，占) ，满足： i ) u ( x ，0 ) 一，( z ) 是小量； i i ) ( c ( ( ，) 警_ o ( ( ，) 等“( 仉蚴驴0 歹= 1 ，m ；k = 1 ，r 一1 ；0 t t ； i i i ) u ( o ，t ) = g oc t ) ，u ( 1 ，t ) = 9 1 ( t ) ，0 t t 其中a l ( r ，6 ) = u c 1 ( 驯”只( 五) ，t = i ，m ) ，r 3 ，耳( e ) 表示i 上 次数不超过r 的多项式函数空间在e 上的限制i = 【0 ，1 1 ，h = k i 一1 ，反l ，t = 1 ，2 ，m ，其中 孔 丝。是区间i 的剖分结点 时间离散即全离散配置格式是指用有限差分在时间方向离散：令 个 u n = u “ ) = ，( z ，t n ) ，k = n a t ，t2 素； 6 山东大学硕士学位论文 u - + = 三( 严+ u n + - ) ，五u 竹二竿 则方程的c r a n k - n i c o l s o n 配置格式为：寻找映射u ：t o ，l ，t _ 人l ( r ，6 ) 满足： ) u o 一，是小量； i ) c ( 妙n + ) d t u - a ( u n + ) 矽譬一b ( u n + ，瑶+ ；) 代巧) ：0 ， j = 1 ，m ：七= 1 ，r l ；仡= 0 ，n 一1 ； i i ) u “( 0 ) = g o ( t 。) ，u “( 1 ) = g l ( “) ，r i = 0 ， 其中a l ( r ，6 ) 的定义同上 上一章我们介绍了在一类不规则区域上处理椭圆问题的一种结合了有限 元灵活性的配置方法，本章的主要目的是将这种方法推广到抛物问题上 考虑热传导方程初边值问题， f 象( z ，耖一u ( z ，秒，) = o ， ，可) q ，t ( o ，卅 ( 2 2 1 口) r 1 o ，关h _ g ， ( 砌) 砌艇( o ，刁( 2 姗) o 、 新0 图4 ：扇环形区域 7 山东大学硕士学位论文 其中q 为在极坐标中用n r b ，口0 p 表示的扇环。如图4 所 示魂为q 的边界，n + r 2 = 撇，其中r l 表示曲线r = b 加上0 = o 和0 = p 部分的边界，而r 2 表示r = o 部分的边界 给定极坐标下的两个剖分t 口= t o t 1 r 2 7 r a = b ， o l = o o 0 1 0 2 靠= p 则区域q 被划分为图5 所示的楔形网格，显然这种网格所划分的四边形元具 有成直角的隅角，符合第一章对区域和网格的要求 图5 ：楔形网格 首先我们对( 2 2 1 a ) 式的时间t 进行有限差分设u ( x ，y ，t ) 为方程的近似 解，记 扩= 矽“( z ，分) = u ( x ，y ，k ) ，“= n a t ，t2 亩； c ，n + 吾= 三( u n + u n + 1 ) ，也u 拜= u n 斗矿l _ g n 则方程( 2 2 i a ) 可写成如下离散形式： 磊沪一沪+ 考= 0 ( 2 2 ，2 ) 类似地，我们选择每个结点处有4 个自由度的双三次h e r m i t e 元定义在基 本正方形( 图1 ( b ) ) 上的检验函数仍取为 8 谶叩) = 妾f 妒十科。簧+ 鳄1 等+ 鸩1 客1 ( 2 2 3 ) 山东大学硕士学位论文 其中，簧，筹，意是结点未知量双三夹h e r m i t e 。基函数定义为 钙o = 去( f 十岛) 2 ( 乓岛一2 ) + 仍) 2 ( m 一2 ) = 一击白( f 十6 ) 2 ( 鸽一1 ) ( 露+ 仍) 2 ( 栅一 = 一去恁+ 岛) 2 岛一2 ) 仍( ，7 + ) 2 ( 印仍一 = 击岛( + 岛) 2 ( 鹃一1 ) 协+ 仍) 2 ( 叩仍一卜 其中岛，仍是结点j 的坐标值这些函数都是二阶可微的，因此适合代入方程 ( 2 2 1 ) 利用坐标变换( 1 2 3 ) 把检验函数方程( 2 2 3 ) 变换到o r ，7 坐标。 砬( 仉r ) = 其中变换基函数 其中的结点导数 斟4 o u ：十蝣鲁十蛎等增嘉1 后面会专门讨论 得 ( 2 2 5 ) 把( 2 2 5 ) 代入方程( 2 2 2 ) ，就得到配置方程对于划分成m 个楔形 元的扇环形区域，每个元内取4 个高斯结点，我们得到的配置方程是写在所 有m = 4 m 忌个配置点上的为了便于书写，我们记 屯= ( 鹋，蝣，孵，妒) ， 叼娟务磐，等，筹n 9 、l，、l，、l，、l， 硒 & 耐 2 2 2 2 仁 q 偿 生却铲一越 以巧一研弓一田 + 铲一垮a；c罄藩 o 0 i 0 秽 鳄 却 = = i l r = 旦却铲一武生却伊一武耽两堕西 山东大学硕士学位论文 则铲( 吼r ) = 呜叼，于是配置方程可以写成 ( 轰一学一学) 旷l 】( g ) = 耻差+ 学+ 学 ( 2 2 7 ) 其中 c ：c 等20 ，箬，箬，警， c 蚶= c 第，等，等，等， q ，i = i ，2 ，m 为全部配置点，求和是对每个配置点c 所在元的4 个结点 而作，基函数导数以及之在配置点处求值 把方程组( 2 2 7 ) 写成矩阵形式为 c u n + 1 = c u ，叼、 其中c 为系数矩阵- 扩2 【? j ， 将u 1 中已知边界值连同相应系数移到右边。记为f n ，并记余下的结点未 知量为【2 “，相应的系数矩阵变为a ，则有 a 叼“= g u n + p ( 2 2 8 ) 这就是我们要求的配置方程，可以求出结点未知量 由上面的分析可知，虽然元素的具体表示形式不同。我们对热传导问题应 用第一章中介绍的算法产生的配置矩阵与势能问题的配置矩阵具有相同的形 式这样我们就可以求出热传导问题( 2 2 1 ) 的配置解 在求解之前。我们还要解决两个问题t 一确定基函数o ，护，矿，吖的 二阶导数；二确定( 2 2 6 ) 中的结点导数 这两个问题1 中都有详细分析对本例的讨论，第一个问题可以得出和 中完全相同的结论，这里我们只简要给出结果 对于标准基函数妒= 多恁，7 ) ，二阶导数可展成 象= 篆簧+ 蠹筹+ 2 蹇塞器c 墓，2 雾+ c 瓦o n ) - i - 22 翥2 c 2 2 9 口， 一= 一一十一一- 1 rz i j 一- 1 i 一 iz zm ，z - 如2如2 必。如却一如沈陡锄如7 吠2 却 r 一7 山东大学硕士学位论文 玳= 喜( 如埘簧聊筹y x i ， z ( ，7 ) = ( 钙。巧+ 钙。鬻+ 钙1 丽j 3 ，( f ，7 ) ：声( 罗耽+ 旌。霉+ 旌t 挚) 其中出现的结点导数由旋转可得； 簧= 篑哪 一= 一f n s r ，二 必必叮 蔷一鲁咖岛却却一 簧= 簧咖易 一= ；s l nh ： 爱必一7 = o 儿1 - - - a | c o s o j 8 nu ，h 其中p 是从水平方向到盯轴的转角 对方程( 2 2 1 0 ) 求导，形成j a c o b i 矩阵，由其逆 厂1 = ( 2 2 9 b ) 此只需要 的基本正 ( 2 2 1 0 a ) ( 2 2 1 0 b ) ( 2 2 1 l 口) ( 2 2 1 l b ) ( 2 2 1l c ) ( 2 2 1l d ) ) 一1 一高( 一娄一霎) = ( 蓁萎) c 2 2 ，2 ， 可得一阶导数篡，塞，荔，舅的值二阶导数的计算在的附录中有详细介 绍，我们只给出结果： 萎；i)(一三)2兰4c 赛，2 八。衢 却一曲元 剐增蝴湘 黝一孵惭 删驹咖m蝴啉 瑶篡警硼一却剔挫把赎 棚枷一椭莩一 争脯硼虾 两芘西一卸如一必如一却 封 n l ，j、ilj，2 3 t 5 阳r d 口d 一 一 却一却缸一却以一越珏一鹫胁一暇c一假 a a h 溉秀塑却 一卅 萨锄矿鱼群塑妒 我丽欢妒 山东大学硕士学位论文 其中。从= i 竺：宠i ，p m 。c m ，n = ，2 ，3 ，4 ，5 ，是矩阵 p ：2 a z 戳 a z 却 伊z 延却 伊z 越2 0 2 z 却2 a s 必 却 却 伊 武却 0 2 耖 民2 0 2 y 砌2 0 0 o0 00 ，沈西o z 西、如o z 西o u 、必卸。却茂7 麓却嚣却 2 赛赛c 褰) 2c 赛，2。a 越、a 7 、a 2 骞岛( 骞) 2 ( 赛) 2。却却。却7、却7 下面我们根据【1 1 中的讨论来确定本例中的结点导数警，玺，亳，和 凌0 2 丽1 - j ，即确定第二个问题 为了简化运算、减少运算量，我们在算例中使用等距剖分- 令r = 与，p = 丝n 2 ，贝0n = a 十t 7 ，l = o ，l ，n l ；以= a + j a o ，j = o ，1 ，他 在【l l 中由于四边形元的不规则性，为保证梯度的连续性。对一阶导数的 逼近使用了线性插值由于本算例及其剖分的特殊性，我们直接使用曲边元和 基本正方形相应边长的比，即 ( 氟= 鱼2 ， c 飘= 争， 由于区域的特殊性以及选择了等距剖分。事实上我们有 ( 氟= 警， ( 2 2 讹) ( 务l 】；= 百h - ， ，( 2 2 ，1 4 b ) 这样对所有的葶，歹，葡o r 只需计算一次t 而砸0 0 也只需要计算死z 次。大大减少 1 2 山东大学硕士学位论文 对于交叉导数蠢，蠢，和【1 】中一样。我们沿着元边界把它们取为常 色( 飘- 州塾r 鱼l ! 垂竺! 二! 亟! ! 、炎铆”。 厶乇 丝( 孰一( 孰，盟l 。! 互! ! 二：互! 兰 、d 毛拶叼”。上乙 利用( 2 2 9 ) 式我们有 ( 翥k ? 警圳( 2 2 1 5 a ) ( 意k2 喾= 。( 2 2 1 5 b ) 至此，我们可以给出配置方程为： c u 。v ：2 + 蜘l = g u 竹+ f “ n = 0 1 - 其中g ，职，c ，泸以及f “定义如前 由于系数矩阵在各时间层完全相同，只需计算一次，所以对于抛物方程 配置有限元法只增加了时间方面的运算，计算量仍然比较少 我们把货：法简要概括如下： ( 1 ) 计算所有配置点已，i = 1 ，2 ，m 以及网格结点； ( 2 ) 根据( 2 2 1 4 ) 和( 2 2 1 5 ) 式计算相应的结点导数从而求出( 2 2 1 1 ) 中的 一o x j 堕塑和塑 必却必1 ”却 ( 3 ) 根据( 2 2 6 ) 在配置点上计算( 2 2 7 ) 中的罟； ( 4 ) 计算( 2 2 4 的导数簧，鬻，意 以及象，髻； ( 5 ) 根据( 2 2 t 2 ) 和( 2 2 t 3 ) 计算墓，塞，喾，舅以及嘉，象，雾，雾； ( 6 ) 把( 2 ) ，( 4 ) ，( 5 ) 步的结果代入( 2 2 。9 ) 计算的二阶导数； 1 3 山东大学硕士学位论文 ( 7 ) 根据( 2 2 7 ) 式确定系数矩阵c ，并确定相应的g ； ( 8 ) 计算u o = u o ，并置r , = 0 ； ( 9 ) 确定f n ； ( 1 0 ) 根据配置方程g u 2 + 1 = c u n + p ，求解 “； ( 1 1 ) 令n ：= n + l ，转( 8 ) 我们可以根据以上算法编写程序来证明算法的有效性由于作者编程基 础较弱，程序存在不足，运行出的结果不理想附录中给出了前面所提到的求 解前应解决的两个问题部分的程序 我们已经知道，按照图3 ( a ) 的编号方式，a 的结构仍为图3 ( b ) 所示。 不同的只是元素的具体形式，因此，1 1 l 中对该方法与g a l e r k i n 有限元法计算 有效性的比较仍然成立同时由【1 】中的分析可知，本文所介绍的配置有限元 法在空间上的精度可达d ( 酽) ，而对时间采用一阶差分。时间精度为o ( ) 因 此，结合了可变元的配置法在解抛物方程时也具有相当高的计算效率以及理 想的精度 2 3 区域的简单推广 上一节我们的算例区域为扇环，如果令o l = 0 ，口= 2 1 r ，原有的两条半径 部分的边界重合，即取口，b ，则区域变成圆环；如果再把a 取成一个极 小的正数，区域就接近整圆，这对实际问题，尤其是热传导问题无疑 力口t 具有重要意义这时，我们可以把边界条件改为u i r = 0 ，兰i r 2 = g ，( z ，耖) 0 f l ，t ( 0 ，r i ，其中r l 表示外圆，r 2 表示内圆，r i 十r 2 = a q 这时。由于区域不再足四边形，结点未知量的编号有所变化。当然所形成 的矩阵也稍有不同对划分为3x4 网格的圆环。我们把结点未知量的编号方 式用图6 ( a ) 表示，对于配置点的编号，仍然是分配给最近的结点并按结点连 续编号图中省略了配置点的编号而相应的配置矩阵结构用图6 ( b ) 表示， 这个矩阵容易求解因此，我们的推广是有意义的 1 4 山东大学硕士学位论文 嚣：；嚣 # # 4 1 n l t j t n l l ：l 一 嚣： 幸 _ 援攒# 4 t 瓣 4 4 t 材 材 i h t ； 蚺 t # 掇援 图6 ：配置矩阵t( a ) 结点未知量的编号网格，( b ) 圆环的配置矩阵 2 4 本文的不足 由于作者的水平有限，本文还存在很多不足之处其中最主要的不足在于 以下两个方面： 一、本文所介绍的方法摆脱了配置法在矩形区域上的限制，结合了可变元 的灵活性，是一大突破；但所选区域的局限性很大，不能实现广泛应用 二，本文没能给出完整有效的程序运行出理想的结果来证明算法的有效 性，在实践方面有待提高 1 5 山东大学硕士学位论文 附录 对2 2 的后半部分提出的两个问题，即：一确定基函数扩，妒，矿，矿r 的二阶导数；二确定( 2 2 6 ) 中的结点导数当网格数为2 0x2 0 时，这一部 分的程序整理如下( 由于个人能力的原因，其中部分语句使用有误，但从中可 体现出实现思路) ： 1 6 x l = z e r o s ( 1 ，2 1 ) ； f o ri = 1 ：2 1 ； x z ( i ) = 0 i + 0 - i ) 0 8 2 0 ； e n d y l - - - - z e r o s ( 1 ，2 1 ) ； f o ri = l ：2 l ： y 1 ( i ) - - - - p i ( i 一1 ) 4 0 ； e n d d x l = o 9 2 0 ； d y l = p i 4 0 ； m = ( 3 - s q r t ( 3 ) ) 6 d y l ； b l = ( 3 + | q r t ( 3 ) ) 6 d y l ； a 2 = ( 3 - s q r t ( 3 ) ) 6 d x l ； b 2 = ( 3 + s q r t ( 3 ) ) 6 + d x l ； y 2 ( p ) = z e r o s ( 1 ，1 6 0 0 ) ； f o rp = l ：4 0 ： y 2 ( p ) = y l ( o ) + a l ； e n d f o ri - - 0 ：1 8 ； f o rj = 4 1 ：2 ：11 9 ； p = j + i ： y 2 ( p ) - - - y l ( i + 1 ) + b l ； e n d f o rj = 4 2 ：2 ：1 2 0 ； p = j + 8 0 。i ； y 2 ( p ) - - - - y l ( i + 1 ) + a l ； e n d e n d f o ri = 1 5 6 1 ：1 6 0 0 ； y 2 ( p ) - - - - y 1 ( 2 0 ) + b l ； e n d x 2 ( p ) 司娜( 1 ，1 6 0 0 ) ； f o rp = 1 ：2 ：3 9 ； x 2 ( p ) = x l ( ( p - 1 ) 2 ) + a 2 ； e n d f o rp = 2 ：2 ：4 0 ； x 2 ( p ) - - x l ( p 2 i ) + b 2 ； e n d f o fi = l ：4 0 ： f o rj = 3 9 ：4 0 ：1 4 7 9 ； p 司+ 2 + i ； x 2 ( p ) = ) 【2 ( i ) ； e n d f o rj = 4 0 ：4 0 ：1 4 8 0 ； p _ - - - j + 2 i ： x 2 ( p ) = 垃( i ) ； e n d e n d f o rp = 1 5 6 1 ：1 6 0 0 ； x 2 ( p ) = x 2 ( p - 1 5 6 0 ) ； e n d x 1 = z e r m ( 2 1 ，2 1 ) ； f o ri = 1 ：2 1 ； f o rj = 1 ：2 1 ； x l ( i j ) = x l o ) * o