（电工理论与新技术专业论文）蒙特卡罗方法在电磁场分析中的应用研究.pdf

华北电力人学硕十学位论文摘要 摘要 蒙特卡罗法( m o n t ec a r l om e t h o d ，简称m c m ) 是概率论学科中重要分支之 一，该方法既可以解决随机性问题，也可以解决确定性问题：随着电子计算机的飞 速发展，其应用范围越来越广泛，对m c m 研究具有一定的理论和实际应用价值。 本文应用m c m 解决静态电磁边值问题为主要目的，重点解决了三维场和多媒质场 边值问题的m c m ，填补了国内在这方面的空白。笔者在m a t l a b 环境下，自行设 计编写了求解静态电磁场边值问的m c m 程序。通过典型算例验证了算法和程序的 正确性，并对离散型和连续型m c m 进行了比较。最后分析了影响m c m 的主要因 素，得出了非常有用的结论。 关键词：电磁场，数值分析，蒙特卡罗法，随机游动 a b s t r a c t m o n t ec a r l om e t h o d ( m c m ) i sa ni m p o r t a n tb r a n c ho fp r o b a b i l i t yt h e o r y i tc a n s o l v en o to n l yt h er a n d o mp r o b l e m s ，b u ta l s ot h ed e f i n i t ep r o b l e m s a l o n gw i t ht h er a p i d d e v e l o p m e n to fe l e c t r o n i cc o m p u t e r s ，t h ea p p l i c a t i o ns c o p eo fm c m i sm o r ea n dm o r e w i d e s p r e a d s oi th a ss o m et h e o r e t i c a la n dp r a c t i c a la p p l i c a t i o nv a l u et os t u d ym c m t h e m a i n p u r p o s eo f t h i s p a p e ri s t os o l v e b o u n d a r yv a l u ep r o b l e m so ft h es t a t i c e l e c t r o m a g n e t i cf i e l d ，f o c u s i n go nt h e t h r e e d i m e n s i o n a lf i e l da n dm u l t i m e d i a f i e l d ，w h i c hh a sf i l l e dd o m e s t i cb l a n ki nt h i sa r e a t h ep r o j e c to fm c ms o l v i n gs t a t i c e l e c t r o m a g n e t i c f i e l di s d e s i g n e d a n dt h e p r o g r a mi sd e v e l o p e d i nm a t l a b e n v i r o n m e n t i na d d i t i o n ，t h ea c c u r a c yo fm c ma n dt h ep r o g r a mi sc o n f i r m e dt h r o u g h t y p i c a le x a m p l e s ，a n dt h ef i x e dr a n d o mw a l ki sc o m p a r e dw i t ht h ef l o a t i n gr a n d o m w a l k f i n a l l y , t h em a i nf a c t o r s w h i c hm a ya f f e c tm c ma r ed i s c u s s e da n ds o m eu s e f u l c o n c lu s i o n sa r ed r a w n w uw e i g e ( t h e o r ya n da d v a n c e dt e c h n o l o g yo fe l e c t r i c a le n g i n e e r i n g ) d i r e c t e db yp r o f l i uj i a n x i n k e yw o r d s ：e l e c t r o m a g n e t i cf i e l d ，n u m e r i c a la n a l y s i s ，m o n t ec a r l om e t h o d ， r a n d o mw a l k 声明户明 本人郑重声明：此处所提交的硕士学位论文蒙特卡罗方法在电磁场分析中 的应用研究，是本人在华北电力大学攻读硕士学位期间，在导师指导下进行的 研究工作和取得的研究成果。据本人所知，除了文中特别加以标注和致谢之处外， 论文中不包含其他人已经发表或撰写过的研究成果，也不包含为获得华北电力大 学或其他教育机构的学位或证书而使用过的材料。与我一同工作的同志对本研究 所做的任何贡献均己在论文中作了明确的说明并表示了谢意。 学位论文作者签名：日期：丝之：! 三：! ， 关于学位论文使用授权的说明 本人完全了解华北电力大学有关保留、使用学位论文的规定，即：学校有 权保管、并向有关部门送交学位论文的原件与复印件；学校可以采用影印、缩 印或其它复制手段复制并保存学位论文； 学校可允许学位论文被查阅或借阅； 学校可以学术交流为目的，复制赠送和交换学位论文；同意学校可以用不同 方式在不同媒体上发表、传播学位论文的全部或部分内容。 ( 涉密的学位论文在解密后遵守此规定) 作者签名： 绰 导师签名： 华北电力人学硕士学位论文 第一章引言 电磁场学科随着计算机技术、计算方法、实验研究技术的迅速发展和工业需 求的同益增加，其研究开发及工业应用逐步深化。电磁场领域在全球范围内的学术 交流频繁，许多国内外学者、专家、工程师投身于这一巨大的国际交流合作中，而 互联网技术的飞速发展极大地推动了它在地球村落内的交流，使得在天涯海角之间 的对话可在瞬间实现，电磁场领域的每两年一度、交替进行的高级国际会议 c o m p u m a g ( c o m p u t a t i o n o f e l e c t r o m a g n e t i cf i e l d ) 、c e f c ( c o n f e r e n c e o n e l e c t r o m a g n e t i cf i e l dc o m p u t a t i o n ) 等提供了交流合作的场合。 1 1 现代电磁技术的发展状况 经典电磁学理论为电磁技术的发展奠定了坚实的理论基础。进入2 0 世纪，随 着各个学科的发展，电磁技术与各相关学科相互交流渗透，形成了许多研究领域， 如：光磁电效应、磁光效应、磁电效应、磁致伸缩效应、磁共振效应等。 麦克斯韦方程组完整而充分地反映了电磁场客观运动规律，成为求解电磁场问 题的基本依据。工程上求解电磁场问题的基本任务是：首先根据物理场域和媒质特 性建立数学模型，即利用麦克斯韦理论导出控制方程，结合定解条件及源函数构成 边值问题( 或初值问题) ，求解出符合实际的场分布( 有时是反问题，求源分布) 。 电磁场边值问题的解法有图解法、模拟法、解析法和数值分析法。近2 0 年来， 由于高速度、大容量电子计算机的广泛应用和电算技术的进步，数值方法得到了迅 速的发展，其主要有：有限差分法、有限元法、积分方程法、边界元法和蒙特卡罗 法等基本类型，而且仍在继续发展，方兴未艾。 有限差分法( f d m ，f i n i t ed i f f e r e n c em e t h o d s ) 是最早应用于电磁场数值分析领 域的数值方法之一。它以概念清晰，方法简单、直观等特点，在电磁场数值计算领 域得到广泛应用。以此方法为基础发展起来的时域有限差分法 ( f d t d ，f i n i t e d i f f e r e n c et i m e d o m a i nm e t h o d ) 正在蓬勃发展，主要应用其对时 域场、瞬态场求解，在传输线、天线、电磁兼容预测和生物电磁分析等诸多领域得 到广泛的应用。 有限元法( f e m ，f i n i t ee l e m e n tm e t h o d ) 是上世纪6 0 年代提出，它特别适合于处 理场域内介质种类较多、交界面形状复杂的情况，其形成的代数方程具有系数矩阵 对称正定、稀疏等特点，所以该方法具有求解容易、收敛性好、占用计算机内存较 少等特点，成为静态场、工频涡流场求解的主要工具【1 1 。 积分方程法( v i e m ，v o l u m ei n t e g r a le q u a t i o nm e t h o d ) 的基础是麦克斯韦方程 1 华北电力人学硕士学位论文 的积分形式和毕奥一萨伐尔定律，对于线性问题具有较高的精度，特别适合于开域 场问题，但对于非线性问题，求解很复杂。 边界元法( b i e m ，b o u n d a r yi n t e g r a le q u a t i o nm e t h o d ) 是近1 0 年来发展形成的 数值计算方法。它计算精度高，数据处理工作量小，适合于解无限域问题，但在处 理多介质问题和复杂非线性问题等方面，不如有限元方便。 1 2 蒙特卡罗法的研究概况 蒙特卡罗法( m o n t ec a r l om e t h o d ，简称m c m ) ，也称统计模拟方法，或称计算 机随机模拟方法。这一方法源于美国在第一次世界大战进行的研制原子弹的“曼哈 顿计划”。该计划的主持人之一、数学家冯诺伊曼( v o nn e u m a n n ) 用驰名世界的赌 城一一摩纳哥的m o n t ec a r l o 来命名这种方法，为它蒙上了一层神秘面纱。 用蒙特卡罗法解数学问题是从相反关系出发，一旦某些概率满足所考虑的数学 方程，便作若干次随机抽样试验，产生随机变量，取其结果的平均值作为数学方程 解的近似估值。1 9 0 8 年，统计学家g o s s e t 在推导t 分布的同时，也猜出了样本方差 s 2 分布。为检验方差分布的正确性，他以3 0 0 0 个接近正态分布的人体测量数据作为 总体，用有放回抽样的方法抽样的方法抽出5 0 0 0 个n = 4 的样本，得出了样本方差 的分布，同时也验证了他理想中的t 分布的正确性。这就是数理统计中用模拟抽样 求样本统计量分布的一个例子。然而，在当时随机试验却受到一定限制，这是因为 要使计算结果准确度足够高，需要进行的试验次数相当大。因此人们都认为随机试 验要用人工进行冗长的计算，这实际上是不可能的。 随着电子计算机的出现和迅速发展，人们才有意识地、广泛地、系统地应用随 机抽样试验来解决数学物理问题，而且把蒙特卡罗法当作计算数学的一个新的重要 分支。随着计算机使用范围的日益广泛，它向各个学科的渗透也越来越深入。其中 计算物理学、计算概率统计学、计算机科学与统计学的界面学等边缘学科，它们都 和蒙特卡罗法有着密切的关系。 1 9 5 5 年以后，我国开展了蒙特卡罗方法的研究工作，在核科学、真空技术、地 质科学、医学统计、随机服务系统、系统模拟和可靠性等方面都解决了大量的实际 问题，并取得了一批理论和应用的成果i2 1 。 蒙特卡罗法在电磁场领域的发展如下： 1 9 4 4 年，k a k u t a n i 首先提出电位理论和随机游动之间的关系，文献1 3 7 j 进一步 扩展了k a k u t a n i 工作； 1 9 5 6 年，m 。e m u l l e r 在求解l a p l a c e 方程边值问题时，提出了一种不用差分近 似的连续型蒙特卡罗法【8 】； 1 9 6 2 年，王建华提出一般椭圆型方程第一边值问题蒙特卡罗法。还以二维 2 华北电力人学硕士学位论文 l a p l a c e 方程为例，说明了椭圆型方程第二、三边值问题蒙特卡罗法【9 j ； 1 9 7 1 年，g m r o y e r 提出了用连续型蒙特卡罗法解决静电场多介质边值问题的 处理办法。指出了蒙特卡罗法在求解部分节点数值时，是一种有效的方法； 1 9 7 8 年，s a n c h e z q u e s a d a 提出了应用离散型m c m 解决多媒质静电场边值问 题的解决方案； 1 9 8 8 年，r e i n h a r ds c h l o t t 弥补g m r o y e r 用连续型m c m 解决静电场多介质边 值问题的处理办法时存在的问题； 1 9 9 0 年，m n o s a d i k u 建议将蒙特卡罗法引入电磁场数值计算教学中，并用单 一媒质二维静电场典型的例子说明m c m 的基本思想； 1 9 9 2 年，m n o s a d i k u 和d h u n t 提出了一种不使用随机数发生器的新型 m c m ，用于计算静电场d i r i c h l e t 边值问题，相比离散型m c m ，该算法能在较 少的时间内达到更高的精度要求，但该方法最大缺点是编程非常困难； 1 9 9 3 年，j a y a n tn j e r e 和y a n n i c kl 运用提前计算转移概率的方法计算多介质 d i r i c h l e t 问题，提高了传统连续型m c m 的计算效率： 1 9 9 5 年，翟景春等解决了l a p l a c e 方程边值问题的离散型m c m ，建立了l a p l a c e 边值问题的概率模型，分析了蒙特卡罗法的误差和收敛速度，最后给出了在计算机 上实现的方法； 1 9 9 5 年，唐象能等首次应用离散型蒙特卡罗法计算了长方体区域内的静电场边 值问题，并提供了精度估算公式； 2 0 0 0 年，裴鹿成研究了蒙特卡罗法与问题维数的课题，指出一般数值方法的误 差不仅受问题维数影响，而且随着问题维数的增大按指数增加。与此情况截然不同 的是，蒙特卡罗法的误差与问题维数无关【1 0 j ； 2 0 0 1 年，马文淦首次在计算物理书中采用离散型m c m 对二维单一媒质静电场 边值问题为例，说明了蒙特卡罗法的应用。指出了蒙特卡罗法在求解部分节点数值 时，是一种有效的方法： 2 0 0 3 年，王娴等针对直接模拟的蒙特卡罗法的特点，提出了3 种计算程序的并 行方案，同时对3 种方案各自的优缺点进行了论述【1 l 】； 2 0 0 3 年，雷桂媛在其博士论文中指出，蒙特卡罗法是2 0 世纪对科学和工程计 算的发展和实践最具影响力的十大算法之一，蒙特卡罗法在新的2 1 世纪会更有光 明的前途【1 2 j ； 2 0 0 4 年，刘春光在其硕士论文中研究了m c m 在解微分方程的边值问题中的应 用，抛物方程以热传导方程为例，对齐次和非齐次方程的蒙特卡罗法分别做了较 详细的叙述，对于具体问题进行了数值试验和误差分析，研究得到的结果表明， 对于给定边值问题，不必一层层求解，可单独求解任意一点，并且可以利用增加试 验次数的方法来提高计算精度对高维情形，也做了简要的分析。文中还以泊松方 3 华北电力人学硕士学位论文 程为例，介绍了椭圆方程的蒙特卡罗法和拟蒙特卡罗法，并做了数值试验： 2 0 0 6 年，吕英华首次在计算电磁学书中采用离散型m c m 对二维单一媒质静电 场边值问题，说明了蒙特卡罗法在电磁场数值计算中的应用【l 引。 综上所述，蒙特卡罗法在电磁场领域中的应用在国内发展相对滞后，特别是在 多媒质电磁场、三维电磁场领域发展非常缓慢。 1 3 本课题研究意义和目的 与有限元法等常用数值方法相比，蒙特卡罗法在求解确定性问题时虽显得有些 “笨 ，但仍有显著的优点：( 1 ) 蒙特卡罗法可以单独求解计算区域内任意一点， 其值无需通过联解其它结点值求出，这使工作量大大减少；( 2 ) 蒙特卡罗法在求解 过程中只涉及局部几个结点的计算，不需要计算刚度矩阵，所用计算容量较少；有 限元法在求解过程中需计算刚度矩阵，所用计算容量随计算结点总数的增加而显著 增加；( 3 ) 在计算精度相当的情况下，蒙特卡罗法计算全场的速度比有限元法要慢 得多，但对于局部区域求解问题，其计算速度并不慢；( 4 ) 蒙特卡罗法对解决二维 和三维的问题来说，没有根本差别，这样容易把二维场问题推广到三维场问题；( 5 ) 蒙特卡罗法程序结构简单，在电子计算机上实现m c m 计算时，程序结构清晰简单， 便于编制和调试；( 6 ) 蒙特卡罗法便于并行计算。 m c m 最大弱点是收敛速度慢，误差大的概率性质。克服这个缺点可以通过并 行计算来增加计算次数，也可以与其他数值方法结合。随着计算机运行速度的不断 提高，蒙特卡罗法计算速度慢的缺点会逐渐消失，其用途会越来越广。因此对蒙特 卡罗法进行研究具有一定的理论和实际应用价值。 1 4 本论文主要研究内容 本文的主要研究内容有以下几个方面： ( 1 ) 回顾了电磁场数值分析领域的发展过程和国内外现状，比较了各种方法的 优缺点。阐述了m c m 基本原理、用m c m 求解问题的一般步骤及其误差分析。 ( 2 ) 在二维静电场边值问题的离散型和连续型m c m 基础之上，解决了三维静 电场边值问题的离散型和连续型m c m ，其内容包括区域和方程的离散、概率模型 的建立、随机游动点位置的确定等。 ( 3 ) 在单一媒质静电场边值问题的离散型和连续型m c m 基础之上，解决了多 媒质静电场边值问题的m c m ，重点解决了在随机游动点处于不同情况下，下一步 游动点位置的确定方法。 ( 4 ) 在m a t l a b 环境下设计编写了应用离散型和连续型m c m 计算二维、三维 电磁场边值问题计算程序和程序框图。 4 华j 匕电力人学硕+ 学位论文 ( 5 ) 分别用离散型和连续型m c m 分析了典型算例，并将计算结果与解析解进 行比较，同时对离散型和连续型m c m 进行了比较。 ( 6 ) 研究了影响离散型和连续型m c m 因素。内容包括：游动次数对计算精度、 计算时间的影响、待求点位置对计算时间的影响、网格长度对离散型m c m 的影响、 实际边界和假定边界问的距离对连续型m c m 的影响。通过实例计算分析，得到一 些有用结论。 5 华北电力入学硕+ 学位论文 第二章蒙特卡罗法基本原理 概率论是研究随机现象数量规律的科学，蒙特卡罗法是其中重要的一个分支。 它的思想虽然古老，但随着电子计算机的飞速发展，其应用范围越来越广泛。该方 法既可以解决随机性问题，也可以解决确定性问题。作为一种数值计算方法，其计 算简单灵活，结果形象直观，容易把二维问题推广到三维问题。 2 1m c m 计算步骤 蒙特卡罗法可解决各种类型的问题，但总的来说，视其是否涉及随机过程的性 态和结果可分为两类：第一类是确定性的数学问题。计算多重积分、求逆矩阵、解 线性代数方程、解积分方程、计算微分算子的特征值以及本文所研究的内容等都属 于这一类。第二类是随机性问题。例如中子在介质中的扩散等问题就属于随机性问 题，这是因为中子在介质内部不仅受到某些确定性的影响，而且更多的是受到随机 性的影响。对于这类问题，通常采用直接模拟方法，即根据实际物理情况的概率法 则，用电子计算机进行抽样试验。原子核物理问题、运筹学中的库存问题、随机服 务系统中的排队问题、动物的生态竞争和传染病的蔓延等都属于这一类。 蒙特卡罗法求解问题的一般步骤【1 4 , 1 5 j 如下： ( 1 ) 针对实际问题建立一个简单且便于实现的概率统计模型，确定随机变量享， 使得随机变量的数学期望e ( 宇) 等于所要求问题量妒。即 驴= e 倍) ( 2 - 1 ) 对于本身就具有随机性质的问题，主要是正确描述和模拟这个概率过程；对于 本来不具有随机性质的确定性问题，就必须事先构造一个人为的概率过程，使得它 的某些参量正好是所求问题的解。 ( 2 ) 对该模型的随机变量建立抽样方法，然后在计算机上进行模拟试验，抽取 足够的随机数岛，岛知。 ( 3 ) 用岛，岛知算术平均值多作为所求量9 的估计值，即 1n 驴一驴= 亩善舅 ( 2 屯) 根据大数定理有 p ( 妒2t 舞孑) 2 1 6 ( 2 - 3 ) 华北电力人学硕+ 学位论文 因此，当n 充分大时 驴= e ( 芋) 一驴= 百1 善n 皇 ( 2 - 4 ) 式( 2 - 4 ) 成立的概率等于1 ，这表明可以用算术平均值乒作为所求量驴的估计值。 2 2m c m 误差 由于多具有算术平均值的形式，袅，岛靠满足同分布相互独立的条件，所以 随机变量桨服从数学期望为0 ，方差为l 的标准正态分布。根据中心极限定理， 0 lq n 当n 充分大时，对于任意a 0 ，有 p 矧 a 卜去扩一倪 5 ， 其中口为显著水平，1 一口为置信水平，仃为均方差。这说明真实值与蒙特卡罗 法计算值之差在范围 驴一缈l 等 ( 2 - 6 ) 的概率为1 一口。 口与a 的关系可通过式( 2 - 5 ) 求得，也可以通过数学表格求得。例如取置信水 平1 一a 一9 9 ，可以查得a = 3 ，可以解释为不等式i 妒一驴i 斋成立的概率是9 9 。 v v 从上面的分析看到，蒙特卡罗法的误差与仃和n 有关。为了减小误差，就应 选择最优的随机变量，使其方差最小，对同一个问题，往往会有多个可以选择的随 机变量，这时就应当择优而用之。在方差固定时，增加模拟次数可以减小误差，把 试验次数增加1 0 0 倍，精度可以提高1 0 倍，但是这样会增加计算机的计算时间【“】。 2 3 本章小结 本章重点介绍了m c m 基本原理和解决问题的一般步骤，同时介绍了在方差固 定的情况下，m c m 误差与游动次数的关系，增加游动次数可以减小计算误差，但 同时也增加了计算机的计算时间。 7 华北电力人学硕十学位论文 第三章单一媒质电磁场边值问题的蒙特卡罗法 本章首先介绍单一媒质二维静电场的离散型和连续型m c m ，其次解决单一媒 质三维静电场的离散型和连续型m c m ，最后对离散型和连续型两种算法进行比较。 3 1 二维电磁场边值问题的m c m 设单一媒质二维静电场l a p l a c e 方程第一类边值问题是： j _ 弋7 2 0 9 = 警+ 窘一。洫。 争。， l 妒= ， o ns b 其中d 为求解区域，s b 是区域d 外边界。 式( 3 一1 ) 边值问题可以用两种m c m 进行求解。一种是离散型m c m ，它需要 把计算区域进行离散，其基础是差分方程；一种是连续型m c m 法，它不需要把计 算区域进行离散，其理论基础是电位平均值理论。 3 1 1 二维场离散型m c m 3 1 1 1 区域剖分和方程离散 7 | - i ly 、 【 ， 、一， 图3 一l 离散型m c m 区域剖分示意图 如图3 1 所示，以间隔缸一命一h 且平行于坐标轴的网格来覆盖区域d ，将区 域分割成许多的边长为h 的正方形。在接近边界s b 处，取与边界距离最近的结点， 8 华北电力人学硕十学位论文 然后将边界上的点连接成一个多边形区域，以s ：表示多边形区域的边界。剖分后 区域内任一结点g ，) ，) 与其相邻四个网格结点关系如图3 - 2 所示。 o ，y4 - j 1 1 ) 【一h ，y ) ( z ，y )b + ，y 0 ，) ，一 ) 图3 - 2 二维场离散型m c m 相邻结点示意图 下面的任务是用差分方程代替微分方程。其一阶导数表示为： 二阶导数表示为： 同理可得到： a 妒9 g + j i l ，y ) 一伊g ，y l a xh 生。翌：堡垒! 兰：z ! = 翌：堡二垒! 兰：兰! o x 2h 。翌堡垒! z ! 二兰翌垒：羔! 翌堡二垒：z2 h 2 旦二里。翌f 兰：z 垒! 二三翌垒! 芝2 翌垒! 羔二垒! o y 2 h 2 把式( 3 - 3 ) 、式( 3 4 ) 代入式( 3 - 1 ) 整理得到： 其中 够( z ，y ) = p 工+ c p ( x + h ，y ) + p 。一妒b h ，y ) + p y + 够& ，y + 危) + p ，一驴g ，y h ) l p x + 4p x 一| py + 2p y 一2 五 9 ( 3 - 2 ) ( 3 - 3 ) ( 3 - 4 ) ( 3 - 5 a ) ( 3 5 b ) 华北电力人学硕+ 学位论文 式( 3 - 6 ) 的概率意义可解释【1 7 - 2 0 】为：随机游动的粒子在点b ，y ) 向相邻结点 b + 庇，y l b 一无，y ) ，g ，y + j 1 1 ) ，b ，y 一无) 移动的概率相等，都是1 4 。 3 1 1 2 概率模型建立 q j i j 一。 j【 一 m ， j 1 m o m 卜斗 图3 - 3 离散型m c 法随机游动示意图 为求解内部结点m o 处电位值驴( m 。) ，现构造概率模型如下：考虑一质点在区域 d 网格结点上作随机游动，如图3 3 所示。设质点从结点m o 出发，第一步以1 4 的 概率向与结点m o 相邻的四个网格结点游动，到达其中的一个结点后，例如到达结 点m 1 ，如果结点m 1 仍是内部结点，第二步质点又继续以1 4 概率向与结点m l 相邻 的四个结点游动，直到质点游到某一个边界点q 时才停止游动，至此就完成了一次 随机游动。如果质点本来就在边界结点上，则质点不游动。 按照下面的方式定义一个随机变量亭。设随机游动任意一条游动路线为： y p m o - m l m 2 呻呻m x - l 呻qe s ： 定义随机变量亭的取值为： 亭；g - ，) = r ( o ) 若质点出发点是边界结点q ，那么就停留在边界点q 处， 亨= g o ) 一i ( o ) 可以证明【2 2 1 ，这样定义的随机变量亭的数学期望 ( 善) 一q g ( m 。) 1 n ( 3 - 6 ) 定义随机变量亭的取值为： ( 3 - 7 ) ( 3 - 8 ) 华北电力人学硕十学位论文 同上所述，一次游动可以得到一个随机变量亭的取值， 以得到n 个随机变量曼，邑昴，用平均值 手。专善氧g - 1 , 2 - n ) 作为驴。) 的近似值。 3 1 1 3 随机游动位置确定 经过次随机游动就可 ( 3 - 9 ) 设质点在游动过程中所在结点坐标是g ，) ，) ，则下一步游动位置的确定方法如 下：首先产生随机数口一v ( o ，1 ) ，然后利用式( 3 一1 0 ) 进行确定。 b ，y ) 一g + j l ，y ) g ，y ) 一g h ，y ) g ，y ) _ g ，y + 乃) g ，y ) 呻g ，y h ) 3 1 2 二维场连续型m c m 矿0 as 0 2 5 矿0 2 5 口s0 5 矿0 5 as 0 7 5 f ，0 7 5 口 0 时，至此完成一次随机游动，取边界s 。上与m 。距离最近的点a 。的边界 值f ( e ，) 作为随机变量的取值。同上述原则，从点m 。出发进行n 次游动， 这样我 们得到n 个随机变量的取粤厂( q ) ，厂( q z ) ，f ( e ) ，其算术平均值专著厂( q f ) 为的 驴( m 。) 近似值。 3 1 2 3 随机游动位置确定 设质点随机游动过程中所处点m i 的坐标是g ，y ；) ，则下一步游动点 m m o m ，y m ) 确定方法如下：首先产生随机数o t u ( o ，1 ) ，再令妒一2 z a ，最后用式 ( 3 1 2 ) 进行计算。 k l 。_ + 姚( 3 1 2 ) i y f + l 。y j + t s i n o i 从m c m 的求解过程可以看出：m c m 可以单独求解区域内任意感兴趣点的电 位值，其值不需要通过联解其它结点值求出，只需要边界点上的电位值即可，这是 m c m 的优势之一。另外m c m 在求解过程中只涉及局部几个结点的计算，不需要 计算刚度矩阵，所用计算容较少，相比之下，有限元法在求解过程中需计算刚度矩 阵，所用计算容量随计算结点总数的增加而显著增加。这也是m c m 的优势之一。 3 2 三维电磁场边值问题的m c m f v 2 t om 壹0 x 2 + 擎+ 壹0 2 2 一。 椭 i 驴= ， o nf 其中q 为空间求解区域，r 是区域q 的边界，m 是区域q 内任二点。 1 2 ( 3 1 3 a ) ( 3 1 3 b ) 华北电力人学硕+ 学位论文 3 2 1 三维场离散型m c m 3 2 1 1 区域剖分和方程离散 首先以等步长a x a yta z h 的正立方体剖分区域q ，在接近边界r 处，取与 边界距离最近的结点，剖分后区域的边界为f i 。剖分后区域内任一结点仁，y ，z ) 与其 相邻六个网格结点关系如图3 5 所示。式( 3 1 3 a ) 的差分形式为： 妒g ，y ，z ) 一p 。+ 9 b + ，y ，z ) + p ，一伊b h ，y ，z ) + p ，+ 驴b ，y + ，z ) + p ，一驴g ，y 一 ，z ) ( 3 1 4 a ) + p ：+ 妒b ，y ，z + j 1 ) + p ；一妒b ，y ，z h ) 其中 1 p x + 一p 工一p y + 一py 一一p ：+ 一p ：一t = ( 3 1 4 b ) o 式( 3 1 4 ) 的概率意义可解释为：随机游动的粒子在点g ，) ，z ) 向相邻结点 g + j i l ，y ，z ) ，g 一 ，) ，z ) g ，y + 庇，z ) ，g ，y j l ，z ) ，g ，y ，z + ) ，g ，y ，z 一 ) 移动的概率相等，都 是1 6 。 b ，y ，z + 办) g ，y + 矗， 产 ， ， ， t l , y ，z ) & + g ，y ，z ) ， ， ， ， ， ，一厶，、 h , y ，z ) 图3 5 三维场离散型m c 法相邻结点示意图 3 2 1 2 概率模型建立 为求解内部结点m o 电位值伊( m 。) ，现构造概率模型【2 2 ，2 3 1 如下：考虑一质点在区 域q 中的网格结点上作随机游动。设质点从结点m o 出发，第一步以1 6 的概率向与 结点m o 相邻的六个网格结点游动，到达其中的一个结点后，如果该结点仍是内部 结点，与第一步相似，质点继续以1 6 概率向与之相邻的六个结点游动，直到质点 1 3 华北电力火学硕士学位论文 一一_ 一 游到某一个边界点q 时才停止游动。如果质点本来就在边界结点上，则质点不游动。 按照下面的方式定义一个随机变量亭。设随机游动任意一条游动路线为： ，p ：m o 一肘l m 2 _ 一m x 一1 呻q r 定义随机变量亭的取值为： 宇= g 。) ；( o ) ( 3 一1 5 ) 若质点出发点是边界结点q ，则停留在边界点q 处，定义随机变量芋的取值为： 亭= g 忆) 一，( q ) ( 3 一1 6 ) 可以证明，这样定义的随机变量亭的数学期 e ( 占) = 驴。) ( 3 一1 7 ) 同上所述，一次游动可以得到一个随机变量宇的取值， 以得到n 个随机变量袅，邑昴，用平均值 ；一号耋氧g = l 2 ) 作为妒似o ) 的近似值。 3 2 1 3 随机游动位置确定 经过次随机游动就可 ( 3 - 1 8 ) 设质点在游动过程中所在结点坐标是b ，y ，z ) ，则下一步游动位置的确定方法如 下：首先产生随机数a u ( o ，1 ) ，然后利用式( 3 1 9 ) 进行确定。 b ，y ，z ) 呻g + j l l ，y ，z ) b ，y ，z ) - - b 一无，y ，z ) b ，y ，z ) 呻x ，y + | l ，z ) b ，y ，z ) 呻g ，y 一庇，z ) ( x ，y ，z ) 呻g ，y ，z + 元) b ，_ y ，z ) 呻g ，y ，z h ) 矿o 口s 丢 矿吉 口s 三 矿j 1 口s 三 矿j 1 as 吾 f ，2 亏 ，f ( q ) ，其算术平均值万1 善n 厂( q ) 为 缈。) 的近似值。 3 2 2 3 随机游动位置确定 设质点随机游动过程中所处点m i 的坐标是g ；，y ，z 。) ，则下一步游动点 m m “小y m ，z i + ，) 确定方法如下：首先产生两个随机数口，p u ( o ，1 ) ，然后令 痧。2 z a ，0aa r c c o s ( 1 2 ) ，最后用式( 3 2 1 ) 进行计算。 f x f + 1 一t + ，：fc o s 破s i n0 i y “l y f + s i n 驴fs i n o i ( 3 2 1 ) 【z f + l z i + 5c o s o i 从二维场到三维场m c m 的计算模型和计算过程来可以看出：用m c m 计算二 维和三维电磁场边值问题没有根本的差别，只要在二维场的基础上稍作改动，即可 实现三维场的计算，这也是m c m 优势之一。 1 5 华北电力人学硕十学位论文 3 3 离散型和连续型m c m 比较 从离散型和连续型m c m 计算过程可以看出：离散型m c m 需要把计算区域进 行剖分，在游动过程中，游动步长是固定的( 等于剖分网格的长度h ) 且游动的方 向都是平行于坐标轴方向。连续型m c m 不需要对计算区域进行剖分，在游动过程 中游动步长和游动方向都是不固定的，其游动一步可以相当离散型算法中游动多 步，能够较快地到达区域外边界而终止游动。因此连续型m c m 具有更高的计算效 率。两种算法详细比较请参见第六章有关内容。 3 4 本章小结 本章首先介绍了二维电磁场边值问题的m c m 。其中离散型m c m 内容包括区 域剖分和方程的离散、概率模型的建立、随机游动位置确定等；连续型m c m 内容 包括概率模型的建立、随机游动位置确定等；从m c m 的求解过程可以看出：m c m 可以单独求解计算区域内任意感兴趣的电位值，其电位值不需要通过联解其它结点 电位值求出，只需要边界点上的电位值即可，这是m c m 的优势之一。另外m c m 在求解过程中只涉及局部几个结点的计算，不需要计算刚度矩阵，所用计算容较少， 相比之下，有限元法在求解过程中需计算刚度矩阵，所用计算容量随计算结点总数 的增加而显著增加，这也是m c m 的优势之一。 其次在二维电磁场m c m 基础之上，解决了三维电磁场离散型和连续型m c m 。 从二维场到三维场的计算模型和计算过程可以看出：用m c m 求解二维和三维场没 有根本的差别，只要在二维场的基础上稍作改动，即可实现三维场的计算，这就是 m c m 优势之一。 、 最后对离散型和连续型m c m 进行了比较。指出连续型m c m 相比离散型m c m 具有更高的计算效率。 1 6 华北电力人学硕士学位论文 第四章多媒质电磁场边值问题的蒙特卡罗法 国内文献只是应用m c m 计算单一媒质电磁场边值问题，实际上各种电磁装置中大 量遇到的是多种媒质共存的状况，为此本章提出了多媒质电磁场边值问题的m c m 。 多媒质电磁场边值问题m c m 是在单一媒质电磁场边值问题基础上实现的，故与单 一媒质电磁场m c m 有许多相同的地方，本章重点针对多媒质电磁场边值问题m c m 中 特别的地方加以分析。 4 1 多媒质场数学模型 设求解区域d 由介质g ，和介质：组成的二维静电场狄里克莱边值问题是： v 2 蚜- 挚+ i a 2 q 丁。l - o 在区域d 1 中 v 2 仍a 矿：q d 2 + 争| 0 褪蜘z 1 ) ，亟。岛亟 在介厦分解面& 上 dd，li 仍。乏 在区域外边界咒上 仍一r 2i 式中，d 。为介质占。组成的区域，d ：为介质s ：组成的区域，s 。为区域d 的外边界， 吼为介质i 和介质s ：的分界面。 4 2 多媒质场离散型m c m 卜啼 ，一r _ _ - - 、 m 1mo - 一 竹，】 r ， 1 m i 一m 3 1 1 f 11_ y m ；1r e 2 l i ， 矽 l 图4 - 1 多媒质场离散型m c m 随机游动示意图 1 7 华j 匕电力火学硕十学位论文 如图4 - 1 所示，首先以缸一y = h 正方形剖分区域d 。其中咒为区域d 的外边 界，s d 为介质，和介质s ：的分界面，离散后区域边界为。 为叙述方便，假设在随机游动过程中某随机点的位置是m ；，则质点从点m ，b ，y ) 向相邻结点转移概率的计算方法分以下四种情况： ( 1 ) 当随机游动点m i 位于y c o n s t a n t 介质分界面上【2 4 】( 例如图4 1 中的点 m 3 ，m 4 ) o ，y + | z ) s 1 ：x h ，y ) g ，y )b + ，y ) e 2 + ) ， 0 ，y - h ) l o 图4 2 多媒质场离散型m c 法相邻结点示意图 如图4 2 所示，对介质分界面玩一d 知应用f l , d s = o 可以得到点m t g ，y ) 与 相邻结点电位之间关系是： 其中 驴g ，y ) - - p ：+ q g ( x + h ，y ) + p ；_ q d ( x 一 ，y ) + p y + 9 g ，y + j 1 1 ) + p y 一驴g ，y h ) ( 4 - 2 a ) p ，+ 。p ；一2 ：1 ip ，+ 一翻p ，一- 翻( 4 - 2 b ) 公式( 4 2 ) 概率意义表明：当随机游动点位于y = c o n s t a n t 介质分界面时，质点 向相邻结点游动概率并不相等，沿x 轴正、负方向的概率都是1 4 ，沿y 轴正方向游 动概率是，2 ( e 。+ 占：) ，沿) ，轴负方向游动概率是：2 ( e 。+ e ：) 。 设随机游动点& ，y ) 位于y c o n st a n t 介质分界面上，则下一步游动位置的确定 方法如下：首先产生随机数口u ( o ，1 ) ，然后利用式( 4 - 3 ) 进行确定。 1 8 华北电力人学硕+ 学位论文 b ，y ) 一x + j l z ，y ) ( x ，y ) _ x h ，y ) & ，y ) _ x ，y + h ) b ，y ) 一b ，y h ) 矿0 o r 0 2 5 矿0 2 5s a 0 5 i fo 5sa 稿 矿鹅鲫订 ( 4 - 3 ) ( 2 ) 当随机游动点m ：位于x = c o n st a n t 介质分界面上。与情况( 1 ) 相类似的推 导可以得出：粒子沿y 轴正、负方向游动的概率是1 4 ，沿石轴正方向游动概率是 毛2 ( f l + g ：) ，沿x 轴负方向游动概率是：2 ( ，+ ：) 。 ( 3 ) 当随机游动点m ；不在介质的分界面s 。上( 如图4 1 中的点 m o ，m ，m ：，m ；) 。这与单一媒质电磁场中的情况相同，即质点向相邻结点游动的概 率相同，都是1 4 。 4 3 多媒质场连续型m c m 4 3 1 理论基础 c 2 矽2 图4 3电位平均值理论示意图 多媒质电磁场边值问题连续型m c m 的理论基础是电位平均值理论另种表达 方式。如图4 - 3 所示，当圆心位于电介质分界面上，并且圆内电介质分界面是直线 时，圆周上电位平均值等于圆心的电位值。即： 荆。焘( d f ) + 丧( 瓢叫 治4 ) 1 9 华北电力人学硕+ 学位论文 4 3 2 随机游动位置确定 如图4 4 所示， s ：是由毛所包围并且接近s 。的假想边界，s ：、s d ”是接近介 质分界面s d 两侧假想分界面界，6 是s 6 - s ：、s d - s ：及s d 一醚之间的垂直距离。 为叙述方便，假设在随机游动过程中某随机点的位置是m ，以点m ；为圆心，半径 r i 构造最大圆为c ；，下一个随机点位置m ；n 产生在圆c ；上。 图4 - 4 多媒质场连续型m c m 随机游动示意图 随机点m g m ，y ) 确定办法分以下四种情况【2 5 - 2 8 】： ( 1 ) 当点m ，处在介质分界面s d 时( 如图4 4 中的点m 2 ) ； 首先把圆e 分成半圆c ；，和半圆c f 2 ，其中g 。处在介质s ，区域，e ：处在介质占2 区域。然后按照式( 4 5 ) 确定到达e ，和e 2 的转移概率，最后在选定的半圆上 均匀分布地产生下一点m m 。 p ，；上p ，。l 一l 。一生一 ? 12 赢见“一隶2 纛 ( 4 - 5 ) 因此在编制程序过程中，当点m ，处在介质分界面s d 时，需要两个随机数，一 个随机数用来确定随机点肘m 产生在半圆c ；，上，还是半圆c ；：上；另一个随机数用 来确定在半圆e 。或e ：上的具体位置。具体实现过程如下：首先产生一个随机数 口u ( o ，1 ) ，然后根据情况产生另一个随机数，即按照式( 4 6 ) 确定下一个随机点 位置。 2 0 华北电力人学硕+ 学位论文 旷a j l 占l + 2 矿