（应用数学专业论文）用径向基配点法求解地下水稳定流问题.pdf

摘要 内容摘要：无网格方法是求解微分方程定解问题的一类新的数值方法，它采用基 于点的近似，可以彻底或部分地消除网格，不仅可以保证计算的精度，而且可以 减少计算难度。径向基函数配点法实施过程简单，是一种纯无网格方法。本文首 先介绍了无网格法的基本原理，然后讲述了m q 配点法的基本原理，并将这种 方法应用于求解非均质多孔介质中的稳定流问题。全文共分五章。第一章主要总 结了无网格法产生的背景、研究现状和发展趋势。第二章介绍了无网格法的基本 原理，并详细介绍了径向基函数和配点法离散。第三章引入径向基函数配点法， 并用其解泊松方程。第四章将m q 配点法用于地下水数值模拟中，通过数值试验， 分析了用m q 配点法解决非均质多孔介质中的二维稳定流问题时节点分布与形 状参数对精度的影响，最后将m q 配点法应用于求解非均质多孔介质中的三维 稳定流问题。第五章对全文进行了总结与展望。 关键词：径向基函数，配点法，多孔介质，稳定流，地下水数值模拟 a b s t r a c t c o n t e n t ：m e s h l e s sm e t h o di san e w l y d e v e l o p i n gn u m e r i c a lm e t h o di ns o l v i n gp a r t i a l d i f f e r e n t i a ie q u a t i o n s b a s e do nt h ep o i n ta p p r o x i m a t i o n ，t h em e s h l e s sm e t h o dc a n e l i m i n a t em e s h e sp a r t i a l l yo rc o m p l e t e l y , t h u sr e d u c et h ec o m p u t a t i o n a lt i m e m e a n w h i l e , i tm a i n t a i n st h es o l u t i o na c c u r a c y r a d i a lb a s i sf u n c t i o nc o l l o c a t i o n m e t h o di sat r u l ym e s h l e s sm e t h o da n de a s y - t o p r o g r a m i nt h i sd i s s e r t a t i o n ，w ef i r s t l y i n t r o d u c et h eb a s i ct h e o r yo ft h em e s h l e s sm e t h o d a n dt h e n ，w en a r r a t et h e f u n d a m e n t a lt h e o r yo ft h em u l t i q u a d r i c s ( m q ) c o l l o c a t i o nm e t h o da n da p p l yt h i s m e t h o dt os o l v es t e a d yf l o wp r o b l e m se n c o u n t e r e di nh e t e r o g e n e o u sp o r o u sm e d i a t h eo u t l i n eo ft h ed i s s e r t a t i o ni sa sf o l l o w s i nc h a p t e ro n e ，w es u m m a r i z et h e b a c k g r o u n d ，c u r r e n tr e s e a r c ha n df u t u r et e n d e n c yo ft h em e s h l e s sm e t h o d i nc h a p t e r t w o ，w eo u t l i n et h em qc o l l o c a t i o nm e t h o da n du s et h i sm e t h o dt os o l v ep o i s s o n p r o b l e m s f o l l o w e db yc h a p t e r4 ，t h em qc o l l o c a t i o nm e t h o di sa p p l i e dt os i m u l a t e g r o u n d w a t e rp r o b l e m s a c c o r d i n gt oo u rn u m e r i c a le x a m i n a t i o ni ns o l v i n gt w o d i m e n s i o n a ls t e a d yf l o wp r o b l e m si nh e t e r o g e n e o u sp o r o u sm e d i a , w ea n a l y z et h e s e n s i t i v i t yo ft h es o l u t i o na c c u r a c yt ot h en o d e sd i s t r i b u t i o na n ds h a p ep a r a m e t e r s t h en u m e r i c a ls i m u l a t i o nf o rt h r e ed i m e n s i o n a l s t e a d y f l o wp r o b l e m si n h e t e r o g e n e o u sp o r o u sm e d i ai sa l s oi n v e s t i g a t e d f i n a l l y , c o n c l u s i o n sa n dp r o s p e c ta r e g i v e ni nc h a p t e rf i v e k e yw o r d s ：r a d i a lb a s i sf u n c t i o n , c o l l o c a t i o nm e t h o d ，p o r o u sm e d i a , s t e a d yf l o w , s i m u l a t i o ni ng r o u n d w a t e r 用径向基配点法求解地下水稳定流问题 学位论文独创性声明 本人承诺：所呈交的学位论文是本人在导师指导下所取得的研究成果。论文中除特 别加以标注和致谢的地方外，不包含他人和其他机构已经撰写或发表过的研究成果，其 他同志的研究成果对本人的启示和所提供的帮助，均已在论文中做了明确的声明并表示 谢意。 学位论文作者签名：圣坚! 至自 学位论文版权的使用授权书 本学位论文作者完全了解辽宁师范大学有关保留、使用学位论文的规定，及学校有 权保留并向国家有关部门或机构送交复印件或磁盘，允许论文被查阅和借阅。本文授权 辽宁师范大学，可以将学位论文的全部或部分内容编入有关数据库并进行检索，可以采 用影印、缩印或扫描等复制手段保存、汇编学位论文，并且本人电子文档的内容和纸质 论文的内容相一致。 保密的学位论文在解密后使用本授权书。 学位论文作者签名：兰坚! 熏鱼指导教师签名：2 訇兰幺毫 一 签名同期： 2 7 1 年明如日 用径向基配点法求解地下水稳定流问题 用径向基配点法求解地下水稳定流问题 1 引言 1 1 无网格法产生的背景 大量科学和工程问题的数学模型都可以归结为偏微分方程( p a r t i a ld i f f e r e n t i a l e q u a t i o n s ，p d e ) 的定解问题。这类问题往往比较复杂，而且计算量大，所以对偏微分 方程数值方法的研究已经成为当前计算数学的主流方向，它也是大规模科学和工程计算 中的核心内容。 有限元法( f e m ) 是数值方法领域最重要的方法之一。它使得数值模拟不仅是计算 数学的重要分支，而且成为工程科学领域的重要组成部分，并深深地影响了工程物理学 科的各个分支。尽管有限元法从理论基础到实际应用都十分成熟，但其内在固有的局限 性使它在一些领域难以应用。例如，在求解金属冲压成形、高速撞击等涉及特大变形的 问题时，有限元需要网格重构甚至会产生严重扭曲。为了达到满足工程要求的计算结果， 需要不断地重新划分网格，这大大地增加了网格剖分的工作量，从而增加了计算成本， 而且也严重地影响解的精度。所以要想彻底解决有限元法面临的这些问题，就应该避免 使用网格的思想。因此，无网格法被提出来并在近年来得到了迅速发展。 1 2 无网格法研究现状及发展趋势 从1 9 7 7 年至今，国内外学者已经提出了几十种不同的无网格法。1 9 7 7 年l u c y 提 出的光滑粒子流体动力学法( s p h ) 被普遍认为是最早的无网格法。归一化光滑函数算 法是由j o h n s o n 等提出的，该方法提高了s p h 法的精度。近年来，s p h 法被应用于一 些数值模拟领域，国内学者对s p h 法也做了很多研究【l 】。s p h 方法优点是格式简单、 计算效率高，但精度较差，在一些对精度要求不高的领域仍然被广泛使用。 1 9 9 4 年b e l y t s c 址。提出了e f g 法，掀起了无网格法研究的热潮。e f g 法是在d e m 方法的基础上发展而来，随后b e l y t s c h k o t 2 】等对d e m 进行了两点改进，提出了无单元 g a l e r k i n 法；c h u n g 等对e f g 方法进行了误差分析；b e l y t s c h k o 等在动态裂纹扩展的数 值模拟中成功的使用了e f g 方法，研究了e f g 法中的积分方案和近似函数的计算方法； k r y s l 等在板壳分析中使用了e f g 法；张雄等将e f g 法的思想应用于节理岩体的分析中。 研究表明，e f g 法精度和收敛速度都比有限元法高，虽然其具有较好的稳定性但e f g 法计算量大，并且需要背景网格进行积分，并不是真正的无网格法。 l i u 等在g a l e r k i n 法的基础上提出了重构核质点法( r k p m ) ，并结合了小波( w a v e l e t s ) 的思想，构造了多尺度重构核质点法( m u l t i s c a l er e p r o d u c i n gk e r n e lp a r t i c l em e t h o d ： 用径向基配点法求解地下水稳定流问题 m r k p m ) 。r k p m 方法精度高，收敛快，推动了无网格法的发展。但r k p m 法需要背 景网格积分，且效率较低。为了解决这个缺点，c h e i l 提出了一种稳定的点积分方案。 b a b u s k a 等提出了单位分解有限元法。a t l u r i 等提出了两种都是用移动最小二乘法建立 场函数的近似方法，分别是局部边界积分方程法( l o c a lb o u n d a r yi n t e g r a le q u a t i o n m e t h o d ，l b i e ) 和无网格局部p e 仃o v g a l e r k i n 法( m e s h l e s sl o c a lp e 廿o v - g a l e r k i nm e t h o d ， m l p g ) 。后者建立无网格形式积分时，不需要背景网格。l b i e 可以看成是m l p g 的一 种特殊情况，但需要进行奇异积分计算。张见明等提出了杂交边界点法，该方法计算精 度高，而且收敛性好。 与上述方法不同，径向基函数配点法不需要背景网格和数值积分，是一种纯无网格 方法。在径向基函数法中，m q ( m u l t i q u a d r i c ) 径向基函数法由于高精确度、形式简洁、 求解过程及程序简单等受到国内外学者的青睐。m q 径向基函数最初是由h a r d y 在 1 9 7 1 年提出的，1 9 8 2 年，f r a n k 分析了2 9 种用于散乱数据拟合的二维的插值方法， 通过比较发现m q 径向基函数综合性能是最好的。1 9 9 0 年，k a n s a 将径向基函数引入 配点法中，发现m q 径向基函数具有指数收敛性，且具有高的精度，实现过程非常简 单，计算效率也很高。吴宗敏等人证明了用径向基函数进行离散数据插值和求解偏微分 方程的优越性，并给出了误差估计。h o n 等人采用径向基函数求解了两相流问题。因此 有必要对m q 径向基函数进行更深入的研究。 无网格法的发展虽然比有限元和边界元晚，但它有自己的优越性和灵活性，它不依 赖于网格，有利于分析复杂的三维问题。目前无网格法仍处在初级阶段，还有很大的发 展空问，随着无网格法的发展，更多数值方法会被应用于更广泛的领域，有效的解决实 际工程问题。 1 3 地下水数值模拟的意义 近年来，随着工农业的发展，城市的扩张，人口的增长，全世界范围内出现了不同 程度的供水危机问题。2 0 0 6 年3 月联合国教科文组织发出通告，全球约有1 5 的人口无 法获得安全的饮用水。地下水是主要的供水水源，由于水资源短缺，长期过量开采地下 水，造成地下水位持续下降，地下水资源逐渐枯竭，同时，还引发了许多环境地质问题。 从保证可持续发展的资源条件来看，这已经成为制约经济和社会发展的主要因素之一。 因此有必要对地下水资源进行更精确的评价，以便能合理开发利用地下水资源。从目前 情况看，地下水数值模拟中存在着一些问题：如何加强参数的研究，提高地下水资源数 值模拟的精度仍是亟待解决的问题；三维问题的解决是值得重视的工作。 在水文地质领域内常用的数值模拟方法有：有限差分法、有限元法、边界元法和有 限体积法等。上述几种方法各有优缺点，但任一种方法都不可能很有效地解决所有数学 物理方程，并且他们共同的缺点是耗费的成本都很高。为消除这种现象，将m q 径向基 2 用径向基配点法求解地下水稳定流问题 函数法应用到地下水数值模拟中是非常有意义的尝试，因为这种方法的易用性和精度都 较传统有限元法有了较大提高，可以更有效的解决多种问题。 1 4 本文的主要工作 本文介绍了无网格法的基本原理，并详细介绍了径向基函数和配点法离散，然后讲 述了m q 配点法的基本原理，求解并重点分析了用m q 配点法解决非均质多孔介质中 的二维稳定流问题时节点分布与形状参数对精度的影响，最后将m q 配点法用于求解 非均质多孔介质中的三维稳定流问题。 3 用径向基配点法求解地下水稳定流问题 2 无网格法的基本原理 2 1 近似函数口1 无网格法的近似函数是通过一组离散点x ，( ，= l ，2 ，以) 来建立的，不依赖于网格。 待求函数“( x ) 可以近似表示为 “( x ) “( x ) = m ( x ) 蜥= n r ( x ) u ( 2 1 1 ) l = l 其中蜥= u ( x ，) 是函数u ( x ) 在节点x ，处的值，( x ) 为节点x ，的紧支函数，也叫形函数， 以为形函数的节点总数， n ( x ) = 【l ( x ) ，n a x ) ，圯( x ) 】7 ，u = 【“。，“：，】2 在多维问题中，近似函数式( 2 1 1 ) 可以改写为 吩( x ) “? ( x ) = m ( x ) = n ( x ) u ( 2 1 2 ) 1 = 1 其中= u t ( x ，) ， n ( x ) = 【l ( x ) ，2 ( x ) ，m ( x ) ，】，u = 【“n ，“：。，”。：，鸭l 一，u n 。】2 ，为 单位矩阵。不同的无网格近似函数具有不同的形函数，与有限元法不同，大多数无网格 近似函数不具有插值特性，因此“，一般不再是试函数“6 ( x ) 在节点x ，处的值，即 “4 ( x ，) u ，m ( x ，) 如。 2 1 1 径向基函数 径向基函数( r a d i a lb a s i sf u n c t i o n ，简记为r b f ) 又称距离基函数，是指以径向坐标 ，为中间变量的复合函数矽( ，) = 矽( 1 i x x 川) ，其中，= u x - x 川是结点x 与弓的欧几里德距 离。对d 维实空间中定义在q 上的任意连续可微函数“( x ) 可以近似为 u ( x ) - - u 6 ( x ) = 吩力( x ) = m 7 ( x ) a x qcr d ( 2 1 1 1 ) j = l 其中x = “，恐，吃) ，a = ( q ，a 2 ，a n ) ，n 为插值结点数，x ，为插值结点的坐标 向量， 吩为待定系数 一般说来，径向基函数可以分为紧支径向基函数( c o m p a c t l ys u p p o r t e dr a d i a lb a s i s f u n c t i o n ，c s r b f ) 和全局径向基函数( g l o b a l l ys u p p o r t e dr a d i a lb a s i sf u n c t i o n ，g s r b f ) 。 后者在插值计算中有较好的精度，但计算量要比前者大。常见的g s r b f 有 m u l t i q u a d r i c s ( m q ) ：国，( x ) = ( c 2 + d ? ) “2( 2 2 ) r e c i p r o a lm u l t i - q u a d r i c s ( r m q ) ：c o , ( x ) = ( c 2 + 彳) 。1 佗( 2 1 1 3 ) g a u s s i a n s ( e x p ) ：q ( x ) = e x p ( 一叫)( 2 1 1 4 ) t h i n - p l a t es p l i n e s ( t p s ) ：c o , ( x ) = 刃芦l o gd t( 2 1 1 5 ) 4 用径向基配点法求解地下水稳定流问题 正定紧支径向基函数有 c s r b f l ： c s r b f 2 ： c s r b f 3 ： q ( x ) = ( 1 一r ) 4 + ( 4 + 1 6 r + 1 2 r 2 + 3 r 3 ) q ( x ) = ( 1 - r ) 6 + ( 6 + 3 6 r + 8 2 r 2 + 7 2 r 3 + 3 0 r 4 + 5 r 5 ) 姒x ) = ；+ ，2 - 詈r 3 + 2 一n ， c s r b f 4 ： 州x ) = 1 l5 + 1 9 6r 2 - 学，3 + j 41 1 6 5r 5 + + 2 ，2 l n ， c s r b f 5 ： q ( x ) = ( 1 一r ) ：( 3 + 1 8 r + 3 5 r 2 ) c s r b f 6 ： 国，( x ) = ( 1 ，) ：( 1 + 18 r + 2 5 r 2 + 3 2 r 3 )( 2 1 1 6 ) 式中c 为大于零的常数，为整数。呸= 0 x 一写0 为计算点x 到节点x ，的距离，：d ，d 。， 厶是定义在节点x ，处的径向基函数的支撑域半径。( 1 - 厂) + 定义为 ”，) + = 一，。i a 砌o 绑 _ 蜓1 叫7 ) 在二维问题中，上述说明对，的定义是假设支撑域为圆形域的情况。关于，更一般 的定义为 ( 2 1 1 8 ) 如果九九，则对应的支撑域为椭圆。本文取九= d 订= 。 径向基函数是一类特殊的函数，它以空间距离，为自变量，具有形式简单、空间维 数无关、各向同性等优点，适合在数值计算中应用，也适用于大型科学与工程问题的求 解。在径向基函数法中，m q ( m u l t i q u a d r i c ) 径向基函数法由于高精确度、形式简洁、求 解过程及程序简单等更受人们关注。m q 径向基函数只与空间点的分布有关，而无需任 何网格剖分。本文采用的m q 函数可表示为 红( x ) = 豚 其中e 为参数。 ( 2 1 1 9 ) 2 2 加权残量法 加权残量法是一种用于求解常微分方程( o d e s ) 或偏微分方程( p d e s ) 龇a 解的通用 的方法，许多数值方法均基于加权残量法。本文所用的配点法也是基于加权残量法的。 许多工程及物理问题都可以用带一组边界条件的( o d e s ) 或( p d e s ) 来描述，例如 5 用径向基配点法求解地下水稳定流问题 怒三厂g 嚣 2 【g ( “) = 在边界r 上 其中f ，g 为微分算子，r 是域q 的边界， u 是待求函数。 在实际问题求解过程中，式( 2 2 1 ) 无法精确求解，只能求其数值解。假设”6 ( x ) 为 “( x ) 的一个近似解，函数u 近似的表示为 u ( x ) m u 6 ( x ) = 口，谚( x ) = 西r ( x ) a ( 2 2 2 ) 其中矽，( x ) 是第项基函数或试函数，a y 是所待求的第歹项基函数的系数， ，l 是基函数 的项数。一般情况下，以越大，近似解“6 ( x ) 越接近精确解u ( x ) 。 将近似函数“6 ( x ) 代入微分方程时，通常会出现 jf ( “| i ：) 一f o 在问譬域? 内 ( 2 2 3 ) 【g ( u “) 一g 0 在边界r 上 因此一般会产生残量r ，和r 6 ： r。s：=f(uh卜)-rg ( u 嚣 旺2 舢 【6 =6 ) 一g 在边界r 上 由于所选择的近似函数的不同，将导致残量函数随之变化。所以应选择适当的近似 函数以使其相应的残量尽可能的小。我们可令其在平均意义上的残量加权积分为零，即 i w , r ，d r 2 + i v r 6 d f = 0 ( 2 2 5 ) 其中i = l ，2 ，o ，n ；w ：和v 为我们选取的权函数。 将( 2 2 4 ) 代入( 2 2 5 ) 得 ，w ie f ( h ) 一厂 d q + ve g ( u h ) 一g l a r = 0 ( 2 2 6 ) 再将( 2 2 2 ) 代入得 f w , f ( 7 ( x ) a ) - f i d e + g ( 西7 ( x ) a - g a r = o ( 2 2 7 ) 将上式对i = 1 ，2 ，z 展开可得： 6 用径向基配点法求解地下水稳定流问题 w l f ( 西7 ( x ) a ) 一厂 d q + v l g ( m7 ( x ) a ) 一g d r = 0 q l w , e f ( 叭x ) a ) 卅船m g ( 叭x ) a ) 一g p ( 2 l2 8 )qr 二o ， f ( 西7 ( x ) a ) 一厂 d q + f v g ( 辔7 ( x ) a ) 一g d r = 0 由此我们可获得以个未知量q ( f = 1 ，2 ，聆) 的刀个方程。求解这些方程可获得q 及 其相对应的近似解，使得残量r 。和r 6 在平均意义上为零。 以上所描述的便为一般形式的加权残量法。加权残量法中权函数的选取将形成不同 的数值近似方法，常见的有：重构核质点法，局部p e t r o v - g a l e r k i n 法，配点法等。 2 3 配点法离散 使微分方程在问题域上的一组特定点上得到满足，这就是所谓的配点法。配点型无 网格法是纯无网格法，不需要借助于任何网格和计算积分。配点法的早期发展和应用可 参见s l a t e r ( 1 9 3 4 ) ，b a r t a ( 1 9 3 7 ) ，f r a z e r ( 1 9 3 7 ) ，l 雒c z o s ( 1 9 3 8 ) 等学者所做的工作【4 1 。基于 配点法的有光滑质点流体动力学方法( s p h ) ，有限点法( f p m ) ，无网格配点法( p c m ) ， h p 无网格云法，最小二乘配点无网格法等。 由于d i r a c 6 函数a ( x x ，) 具有如下性质： ia ( x x i ) = 0 x x i 竹砸1 ) 出= l c _ o q a d 【葺：c 取d i r a e 万函数a ( x x ，) 作为权函数，令 w ：= 烈x i 毛耋i 茎；筹，v = 万。x 三x ，耄i 三；荨 ( 2 3 2 ) 其中i = l ，2 ，刀 将( 2 3 2 ) 代入 即 w r ，d q + v f r 6 d f = 0 q r ( 2 3 3 ) 由加权残量法公式可导出配点法 ，3 ( x - x i ) r ，d q + 1 3 ( x - x i ) r 6 d f = s 3 ( x - x , ) f ( m7 ( x ) a ) 一厂 d q orq + f 3 ( x - x , ) g ( m7 ( x ) a ) 一g d f = 0 ( 2 3 4 ) r 7 用径向基配点法求解地下水稳定流问题 j 6 ( x x ) r ，d q + p ( x x i ) r 6 d f = r 。( x i ) + r 6 ( x i ) - o ( 2 3 5 ) qf 将上式应用于问题域中的玎个特定点，即配点法迫使在问题域中的n 个特定点 x 。0 = l ，2 ，h ) 处的残量为零。 8 用径向基配点法求解地下水稳定流问题 3 径向基函数配点法 3 1m q 配点法基本原理 径向基函数是距离的函数，而配点法是通过配点的形式确定系数的方法。径向基函 数配点法不需要背景网格和数值积分，是一种纯无网格方法。在径向基函数法中，m q ( m u l t i q u a d r i c ) 径向基函数法由于高精确度、形式简洁、求解过程及程序简单等受到国内 外学者的青睐。m q 径向基函数最初是由h a r d y 5 】在1 9 7 1 年提出的，1 9 8 2 年，f r a n k 6 】 分析了2 9 种用于散乱数据拟合的二维的插值方法，通过比较发现m q 径向基函数综合 性能是最好的。1 9 9 0 年，k a n s a 7 】将径向基函数引入配点法中，发现m q 径向基函数具 有指数收敛性，且具有高的精度，实现过程非常简单，计算效率也很高。r b f 在计算中， 不再需要网格剖分，同时又由于r b f 仅与距离相关，不依赖于空间维数，因而特别适用 于高维问题的处理。因此本文采用了m q 径向基函数。 由于径向基函数的连续可微性，将由其近似表示的场函数及各阶导数直接代入边值 问题的微分方程或边界条件即可建立求解待定系数向量a 的代数方程组。配点法的思路 是在域内和边界上选定个点( 配点) ，在这个点上满足微分方程或边界条件，建立 个代数方程，求解个待定系数。 边值问题的一般提法可以表示为( 2 2 1 ) ，将近似函数用( 2 1 1 9 ) 所表示的径向基 函数表示，代入( 2 2 1 ) ，则有 跗锄：兰州厄万) ：厂x q g 。甜。x ，：兰i = 1q g 。j i 再币，：g x r 3 1 2 ) 矩形域内的节点 边界上的节点 图1二维问题的节点分布图 二维问题的节点分布如图l 所示，在求解域布置结点，在区域q 内部布置m 个结 点x 。，x ：，x m ，在边界r 布置结点x ，x 。在个结点处使用方程( 3 1 2 ) ，得到 阶线性方程组 6 ( x t ) ) = 6 ( x 。) ) = 兰州厄习阿mx 。q ) = gx r k = 1 ，2 ，n ( 3 1 3 ) 可求得未知系数( 口l ，a ：，) ，再把( 口。，a ：，) 代入( 2 2 2 ) ，就得到了任意结 点x 处的近似解u 6 ( x ) 。 其中对于二维m q 径向基函数 对其求一阶偏导数 对其求二阶偏导数 缈，( x ) = 厄i 再瓦万 ，：= 9j ) - - - 孵习 型：! 兰： - = 一= = 融 驴，( x ) a 缈，( x ) y y j - ：一= 一=砂 缈，( x ) x x j a 2 缈，( x ) ( y y ，) 2 + c 2 = 一= ：一= c o x 2 妒砖) 3 塑+ 塑= ( 3 x 2 。 加2 一 ( 7 一y ，) 2 + c 2( x x ，) 2 + c 2 ( 彬+ c 2 ) 3( 彬+ c 2 ) 3 ( 3 1 4 ) ( 3 1 5 ) ( 3 1 6 ) ( 3 1 7 ) ( 3 1 8 ) ( 3 1 9 ) ( 3 1 1 0 ) ( 3 1 1 1 ) 对于三维m q 径向基函数 c j ( x ) = 肛i 再石瓦万f ( 3 1 1 2 ) 令_ = 对其求一阶偏导数 可得 伊小) = 仃万 l o ( 3 1 1 3 ) 竖厅 用径向基配点法求解地下水稳定流问题 d 缈f ( x )x x ，x x ， = 一= 二一= = = = = = = = = 三= 8 x c j ( x ) 心t + o a 缈，( x ) y 一少，少一y ， = _ 一= 二= = = = = = = ，_ 砂 缈小) 彬+ c 2 a c j ( x ) ：三三：丝 瑟 驴砖) 彬+ c 2 对其求二阶偏导数 0 2 缈一j ( x ) ：业羔生墼掣：业羔4 竺坐 ( 3 1 1 7 ) 一=-jo：一=：一 li ，j 缸2 ( o j ( x ) 3 ( 呼+ c 2 ) 3 a 2 缈一j ( x ) = 坚二丛生掣望= 堡二兰! ：型! ： ( 3 1 1 8 ) 一=：=j：一=：=r：一 t i lxj c 少 缈，【x j ) 。 ( 哆+ c 2 ) 3 塑掣：坠型丛掣丝生：坠型善垒丛生 ( 3 9 ) 二一= ：- = _ _ - 一= _ = 二= = = = = = = = = - = _ 一 i llvj o | z 缈，( x ) 。 ( 哆+ c 2 ) 3 可0 2 缈j ( x ) + 学+ 学= 丽2 r + 3 c 2 m 。， 从上面我们可以看出，m q 配点法求解稳定问题的过程时，计算量小且形式非常简 洁。从式( 3 1 4 ) 到式( 3 1 2 0 ) 在程序的运算中得以应用。 式( 3 1 3 ) 形成的线性方程组的矩阵形式及其解为 lf ( 矽1 a - - fi l g 坳) jl g j a = 嘲。1 厂l 厂g j m 3 ) l g ( ) j f ( 矽) = g ( 矽) = ，( 办( 五) ) f ( q k l ( x 2 ) ) ，( 唬( 一) ) ，( 政( 屯) ) ，( 九( 五) ) ，( 九( 屯) ) f ( q k l ( x u , ) ) ，( 唬( h ) ) f ( 九( h ) ) g ( 识( h + ) )g ( 唬( h + - ) ) g ( 九( x m + 。) ) g ( 么( + ：) ) g ( 欢( h + ：) ) g ( 九( h + z ) ) ：。： g ( 办( h ) )g ( 唬( h ) ) g ( 九( x d ) ( 3 1 2 4 ) ) ) ) 4 5 6 1 l ，i 1 l ，lil 1 j 1 j 2 j ( ( ( 用径向基配点法求解地下水稳定流问题 从上述m q 配点法的求解过程可以看出，其形式非常简洁。径向基函数只与空间距 离有关，与空间维数无关，且无需任何网格剖分，这使得径向基配点法的研究与应用具 有广阔的前景。 3 2 算例 考虑矩形区域上的二维泊松方程 j “( z ，y ) = 一2 ( x + y x 2 一y 2 ) ( x ，y ) e - o 5 ，0 5 】 ( 3 2 1 ) 【u ( x ，y ) = 0( z ，y ) a q 上式的解析解为 u ( x ，y ) = ( x x 2 ) ( y y 2 ) ( 3 2 2 ) 图2 泊松问题节点分布图 将矩形区域在x 方向和1 ，方向各均匀划分为2 0 段，因此所取的节点数目为2 1 2 1 。 图2 给出了节点分布的情况，图3 给出t m q 配点法数值解。用m q 配点法计算时，取 c = 1 0 。图4 给出了用m q 配点法计算时所有计算点处的误差结果。结果表明，使用m q 配点法的最大相对误差小于1 5 1 0 。 这里相对误差定义为 西( x f ) =l u ( x ，) - u 6 ( x u ( x f ) l i = l ，2 ，m ( 3 2 3 ) 其中“( x ，) 是精确解“( x ) 在计算点处的值，u h ( x ，) 是m q 近似解，x ，是区域内的计算点， m 是计算点数。 1 2 月径向基配点击求解地t 水穗定流月趣 05 0 5 图3p o s s i o n 问题数值解 误差 - 0505 图4 p o s s i o n 问题相对误差 。瞻 盛犁鞴 用径向基配点法求解地下水稳定流问题 4 数值算例 4 1 应用径向基函数法求解非均质多孔介质中的二维稳定流问题 自然界的多孔介质几乎都是非均质的，基于上面对径向基函数法的分析，下面通过 径向基函数法求解非均质多孔介质中的水流问题，来说明这种方法的优越性。所用算例 来自文献【8 】。 假定研究区为正方形，在正方形内部满足控制方程 丢( k 罢) + * 等 = 。 ， 其中k ( x ，y ) = x 2 。 未给定边界条件时，式( 4 1 1 ) 有多解，这里令解析解 h = x 2 3 y 2 ( 4 1 2 ) 在边界上满足第一类边界条件，可由式( 4 1 2 ) 确定。 4 1 1 参数连续变化的二维稳定流问题 假定研究区域对应的x ，y 区间均为【5 ，5 】。k 随x 区间的确定而确定。研究区内 任一点的水头值可用式( 4 1 2 ) 算出。为与传统有限元法进行对比，对计算区域按n n 等分构造网格，在此相同的网格节点上分别用径向基函数法和传统有限元法进行计算。 本文以下各例研究区均用等距的规则网格剖分，单元内子单元的剖分也用等距的规则网 格，如图5 所示。 一专一j _ 一毒一一毒一一专一一毒一毒_ 一丰一毒一一善一 v 节点 卡计算点 月往向基配 法隶解地t 水稳定流目题 c m q 径向基形状参数 x ，区域内的计算点 盯计算点数 精确解“( 砷中的最大值 ”( x ，) m q 近似解 西击 最大误差 西k均方根误差 西矗= 西k ( 2 ) 参数连续变化的二维稳定流问题解析解与数值解 图6 参数连续变化的二维稳定流问题解析解和数值解的图象对比 图6 为参数连续变化的二维稳定流问题解析解和数值解的图形对比。从中可以看出 应用m q 配点法拟合效果很好。 ( 3 ) 节点分布与形状参数对精度的影响 首先给出两种不同的误差定义：最大值误差和均方根误差 ：坚蛑掣 托， 坐u 竺 用径向基配点法求解地下水稳定流问题 西= u m a x i ( 4 1 4 ) 其中x ；是区域内的计算点，m 是计算点数，甜一是精确解甜( x ) 中的最大值，u h ( 墨) 是 m q 近似解。 取n = 1 0 ，1 5 ，2 0 ；c = 0 5 ，1 0 1 5 ，2 0 ，实验结果如表1 。 表1n = 1 0 ，1 5 ，2 0 ；c = 0 5 ，1 0 ，1 5 ，2 0 实验结果 将西缸看作是c 和n 的函数，即历赢= t i c ，n ) ，如图7 所示。历k 曲面，形如峡 谷，沿n 增大的方向，函数递减；沿c 增大的方向，函数先减后增。 1 6 m 径向基配点法衷杆地t 水稳定洼目题 08 08 枷04 媸 02 0 图7 节点分布与形状参数对精度的影响 表2 节点分布与形状参数对精度的影响 节点个数0 5 0 7091 11 3 1 51 71 92 1 l on 4 9 9 60 2 3 4 5n 1 4 2 30 0 2 6 90 0 6 7 8n 1 1 3 60 2 9 7 60 4 5 3 40 7 6 7 6 1 2 n 3 6 8 702 0 9 9 仉1 1 0 50 0 2 0 9 0 0 4 7 80 1 0 0 902 0 8 70 3 7 8 60 7 8 1 40 3 0 6 l01 6 8 9n 0 9 0 400 1 5 60 0 3 9 80 0 9 9 701 3 7 90 3 0 9 80 6 6 5 8 1 60 2 6 7 80 1 2 9 80 0 7 8 400 1 2 90 0 1 2 90 0 8 3 60 1 0 9 80 1 9 8 705 7 6 4 1 802 1 0 901 0 9 300 6 7 8n 0 0 9 70 0 1 0 70 0 7 9 5n 0 8 6 70 1 6 7 80 4 9 8 7 2 0 01 5 3 4n d 9 9 80 0 5 5 90 0 0 5 2 0 0 0 9 80 0 5 2 1n 0 5 6 90 1 3 4 90 3 4 6 5 2 2 0 1 5 4 80 1 0 0 5 仉0 4 5 90 5 4n 0 0 8 30 0 4 8 70 0 4 9 801 2 9 80 3 0 9 4 2 4 0 1 5 5 90 1 0 2 30 0 4 8 70 0 0 6 40 0 1 0 500 5 3 100 5 9 401 3 0 90 3 4 8 7 ( 4 ) 参数c 最佳取值的确定 为确定参数c 的最佳取值，我们取2 0 2 0 的节点分布，多次试验，得到图8 ，西k 随 c 的变化，呈凸性，c 的最佳取值近似为1 2 。在这里，图9 给出了不同参数c 对实验结 果的影响。 甩枉向基配点法求解地t 扣穗定渣目题 图8 参数c 的优选 图9不同参数c 对实验结果的影响( 这里n - 2 0 ) 孵_】 用径向基配点法求解地下水稳定流问题 ( 5 )饥。，随节点数的增加变化趋势 节点的分布可归结为节点均匀分布和节点随机分布。通过大量实验，合理选择形状 参数c 时，两种节点分布所得数值解均有较好精度，而节点均匀分布比节点随机分布得 出的数值解的精度要高出一个数量级。总体来说，在m q 径向基配点法求解偏微分方 程时，节点均匀分布时得到数值解的精度要比节点随机分布时得出的数值解的精度几乎 处处都高。因此，我们这里选用的是节点均匀分布。图1 0 给出了节点均匀分布时不同 节点数对误差的影响。 图1 0 节点数对误差的影响 理论上，误差随着节点数的增加而递减。多次实验中发现：当误差达到一定精度要 求后，节点数再增加，不会带来更好的精度，有时甚至会引起误差略微反弹。 ( 6 ) 结论 从参数连续变化的二维稳定流问题我们可以看出，m q 配点法在求解稳定问题时， 计算量很小、形式也非常简洁。同时，如何分布节点和如何选择径向基形状参数对于用 径向基配点法求解偏微分方程的精度均有直接影响。相对于节点随机分布而言，节点均 匀分布所得到的数值解精度更高一些，而且形状参数所取的范围更广一些。如何选择形 状参数直接影响了求解偏微分方程的精度，这点与用径向基函数去拟合函数所得到的规 律是一致的。 配点法存在不稳定问题的缺陷，即计算中需要对b 阵求逆，而b 阵有时会是奇异、 或接近奇异的矩阵。不稳定问题的解决方法：通常是在b 阵中加入扰动项，但实际计算 时，扰动项的大小不易给出。大的会带来较大误差，小了又没有效果。我们在实验发现， 出现不稳定问题时，把步长取为非整数，比如8 9 ，11 1 2 而不取l ，也可以解决这个问题， 但通常能够达到更好的精度，有时精度竞能提高一个数量级。这只是我们在实验中发现 1 9 用径向基配点法求解地下水稳定流问题 的一个现象，至于它的理论依据和工程价值，没有作深入的研究。 ( 7 ) 传统有限元法( l f e m ) 和径向基函数配点法( 1 沿f c m ) 误差对比 假定研究区域对应的x ，y 区间均为 5 ，2 0 】。k 的取值范围随x 区间确定而确定；研 究区间内的任意一点的水头值可用式( 4 1 2 ) 算出，并且还用了数值方法求出水头的近 似解。本文以下各个例予研究区域均用等距的规则网格剖分。三种方法均将研究区域划 分为7 2 个单元，4 9 个结点；上述不同数值法所求结果的精度见图1 l 。图1 1 是研究区域 内截面y = 1 2 5 m 处各种数值方法所得水头值的误差对比，横坐标为截面上各点的x 坐 标，纵坐标为数值解与解析解间误差的绝对值。径向基函数配点法不仅计算量小，且从 图1 1 中可以看出径向基函数配点法的计算结果更精确。 图1 1参数连续变化条件下二维稳定流在y = 1 2 5 m 处截面上数值解对比 4 1 2 参数渐变的二维稳定流问题 、 研究区为正方形，对应的x ，y 区间为l o k m x l o k m ，以( o