36-二维泊松方程的有限元法
多重网格在二维泊松方程有限元分析中的应用
１ 多重网格方法
多重 网 格 法 最 初 来 源 于 差 分 方 程 的 迭 代 求 解 数值解与差分方程精确解之 间 的 误 差 与 傅 立 ． 叶分析中不同频率 分 量 有 关 ， 不同的频率分量具有
［ ４ １］
犕 狌 犾 狋 犻 狉 犻 犱犳 犻 狀 犻 狋 犲 犲 犾 犲 犿 犲 狀 狋 犪 狀 犪 犾 狊 犻 狊 犳 狅 狉２ 犇犿 狅 犱 犲 犾 犻 狀 犳狆 狅 犻 狊 狊 狅 狀犲 狌 犪 狋 犻 狅 狀 犵 狔 犵狅 狇
， ， ＷＡＮＧ Ｑ ｉ ｎ ｉ ｎ Ａ ＩＷｕ ｍ ｉ ｎ ｏ ｎ ｌ ｉ ａ ｎ Ｂ ＷＡＮＧ Ｈ ｇ ｐ ｇ ｇ ｇ ｇ
Ｇ Ｓ Ｇ Ｓ 狉 狌 ｆ －犃 ｆ ｆ ， ｆ ＝犫 Ｇ Ｓ 然后 ， 将残差狉 ） 矩阵 犚 转换到 Ｒ ｅ ｓ ｔ ｒ ｉ ｃ ｔ ｉ ｏ ｎ ｆ 通过限制 （ Ｇ Ｓ 即狉 而在粗网格上求解犃 粗网格上 ， 狉 狌 ｃ ＝犚 ｆ ． ｃ ｃ ＝
．
多重网格方法通常建立在有限差分方法的基础 本文将多重网格方法引入到 ． 有限元方法中 ， 并以二维的泊松方程为例 ， 对单重网
多重网格在二维泊松方程有限元分析中的应用
王青平 ， 白武明 ， 王洪亮
（ 中国科学院地质与地球物理研究所 ， 北京 １ ） ０ ０ ０ ２ ９
摘 要 本文简要介绍多重网格 （ 算法的基本原理及基本步骤 ， 然后将多重网 格 算 法 引 入 有 限 单 元 中 ， 对二维泊 ＭＧ） 松方程进行求解 ． 单元数尺度从 ８×８ 逐次增加至 １ 并与单重网格中高斯 赛德尔迭代法 （ 、 共轭梯度法 ０ ２ ４×１ ０ ２ ４， Ｇ Ｓ） （ 结 果 表 明 ＭＧ 在 计 算 速 度 和 迭 代 次 数 都 明 显 优 于 Ｇ Ｃ Ｇ）在程序运行时间以及迭代次数方面进行比较 ． Ｓ、 Ｃ Ｇ 方 法． 在１ 比Ｃ 而且与理论解的误差更小 ． ０ ２ ４×１ ０ ２ ４ 网格中 ， ＭＧ 不仅比 Ｇ Ｓ快 ５ ０ ０ 多倍 ， Ｇ 快６ ０ 多倍 ， 关键词 多重网格 ， 共轭梯度法 ， 高斯 赛德尔迭代法 ， 泊松方程 ， 有限元分析 ： ／ ． ｉ ｓ ｓ ｎ ． １ ０ ０ ４ ２ ９ ０ ３． ２ ０ １ ０． ０ ４． ０ ３ ９ 中图分类号 Ｐ Ｄ Ｏ Ｉ １ ０． ３ ９ ６ ９ ６ ３ １ 文献标识码 Ａ ｊ
matlaB程序的有限元法解泊松方程
基于matlaB 编程的有限元法一、待求问题：泛定方程：2=x ϕ-∇边界条件：以（0,-1），（0,1），（1,0）为顶点的三角形区域边界上=0ϕ二、编程思路及方法1、给节点和三角形单元编号，并设定节点坐标画出以（0,-1），（0,1），（1,0）为顶点的三角形区域figure1由于积分区域规则，故采用特殊剖分单元，将区域沿水平竖直方向分等份，此时所有单元都是等腰直角三角形，剖分单元个数由自己输入，但竖直方向份数（用Jmax 表示）必须是水平方向份数（Imax ）的两倍，所以用户只需输入水平方向的份数Imax 。
采用上述剖分方法，节点位置也比较规则。
然后利用循环从区域内部（非边界）的节点开始编号，格式为NN(i,j)=n1，i ，j 分别表示节点所在列数与行数，并将节点坐标存入相应矩阵X(n1)，Y(n1)。
由于区域上下两部分形状不同因此，分两个循环分别编号赋值，然后再对边界节点编号赋值。
然后再每个单元的节点进行局部编号，由于求解区域和剖分单元的特殊性，分别对内部节点对应左上角正方形的两个三角形单元，上左，左上，下斜边界节点要对应三个单元，上左，左上，左下，右顶点的左下、左上，右上边界的左上，分别编号以保证覆盖整个区域。
2、求解泊松方程首先一次获得每个单元节点的整体编号，然后根据其坐标求出每个三角形单元的面积。
利用有限元方法的原理，分别求出系数矩阵和右端项，并且由于边界，因此做积分时只需对场域单元积分而不必对边界单元积条件特殊，边界上=0分。
求的两个矩阵后很容易得到节点电位向量，即泊松方程的解。
3、画解函数的平面图和曲面图由节点单位向量得到，j行i列节点的电位，然后调用绘图函数imagesc(NNV)与surf(X1,Y1,NNV')分别得到解函数的平面图figure2和曲面图figure3。
4、将结果输出为文本文件输出节点编号，坐标，电位值三、计算结果1、积分区域：2、f=1，x 方向75份，y 方向150份时，解函数平面图和曲面图20406080100120140102030405060700.0050.010.0150.020.0250.0320.0050.010.0150.020.0250.03对比：当f=1时，界函数平面图20406080100120140102030405060700.010.020.030.040.050.060.073、输出文本文件由于节点多较大，列在本文最末四、结果简析由于三角形区域分布的是正电荷，因此必定电位最高点在区域中部，且沿x 轴对称，三角形边界电位最低等于零。
有限元法求解二维Poisson方程的MATLAB实现
(x,y) e 9 0 ,
其中
— ax ay
e i 2(/3), 为 i?2 中的
有界凸区域，区 域 / 3 = { ( * ，；K) U2 + y2 < l }.
1 二 维 P oisson方程的有限元法
l .i 有限元方法的基本原理和步骤 有限元法是基于变分原理和剖分技术的一种数
值计算方法，把微分方法的定解问题转化为求解一
摘 要 ：文 章 讨 论 了 圆 形 区 域 上 的 三 角 形 单 元 剖 分 、有 限 元 空 间 ，通 过 变 分 形 式 离 散 得 到 有 限 元 方 程 .用 M A T L A B 编程求得数值解，并进行了误差分析. 关 键 词 ：Poisson方 程 ；有限元方法；M A T L A B 编 程 ；三角形单元剖分
U e l f +2( n , R m).
定理[7]1 (有 限 元 近 似 解 的 炉 模 估 计 ）假设 满足引理的条件，则 对 V f/ E 妒+1(/3,i T ) ，存在与 A 无 关 的常数C , 使得
W u - u . w, ^ chk \ u \ M
定理[7]2 (有 限 元 近 似 解 的 i 2 模 估 计 ）假设
1
0
0
0
2
3
0
0
細 !1[8]:
4
560ຫໍສະໝຸດ 中 图 分 类 号 ：0241.8
文 献 标 识 码 ：A
文章编号：1009 - 4 9 7 0 ( 2 0 1 8 ) 0 5 - 0015 - 04
0 引言
热 学 、流 体 力 学 、电 磁 学 、声 学 等 学 科 中 的 相
关 过 程 ，都 可 以 用 椭 圆 型 方 程 来 描 述 .最 为 典 型 的
泊松方程的解法及应用
泊松方程的解法及应用泊松方程是关于无限大区域内的某个标量势函数的二阶偏微分方程。
它在物理学和工程学中广泛应用，例如在电场、热传导、流体力学和弹性力学等领域。
本文将介绍泊松方程的解法及其在实践中的应用。
一、泊松方程的定义与基本性质泊松方程是具有如下形式的偏微分方程：∇²u = -ρ其中u是标量势函数，ρ是源项，∇²是拉普拉斯算子。
这个方程可以通过库伦定律推导出电力学中的几乎所有问题，是许多物理学领域研究的基础。
泊松方程有一些基本性质。
首先，它是线性的，也就是说，如果两个不同的源项ρ₁和ρ₂产生的标量势函数分别是u₁和u₂，那么对于常数a和b，它们的线性组合a u₁ + b u₂是对应于线性组合aρ₁ + bρ₂的标量势函数。
其次，它是反演对称的，也就是说，如果标量势函数u满足泊松方程，那么-u也满足泊松方程。
二、泊松方程的解法在实际应用中，我们需要求解泊松方程，以便计算出场的分布。
泊松方程的解法通常可以分为两种：1. 分离变量法分离变量法是将u(x, y, z)表示为三个独立变量x, y, z的函数的积的形式，即u(x, y, z) = X(x) Y(y) Z(z)，然后将泊松方程代入并对每个独立变量进行求导，最终得到连个常微分方程和一个初值问题，可由此得到标量势函数u的解析解。
2. 数值解法当求解泊松方程的解析解十分困难或不可能时，可以通过求解离散化的差分方程来得到数值解。
一般使用有限差分法、有限元法或谱方法，这些方法分别将无限大区域内的标量势函数划分为有限数量的子域，并在子域内使用数值技巧求解差分方程。
三、泊松方程在工程学中的应用泊松方程在物理学和工程学中的应用广泛，下面将介绍其中两个重要的应用：电势分布和热传导问题。
1. 电势分布在电场问题中，泊松方程描述了电场中的电势分布。
假设我们有一个电荷分布ρ(x, y, z)，根据库伦定律，它产生了电场。
泊松方程可以帮助我们计算出哪些区域具有高电势、低电势以及电压梯度等性质。
二维泊松方程的差分格式有限差分法
有限差分法（Finite Differential Method）是基于差分原理的一种
数值计算法。其基本思想：将场域离散为许多小网格，应用差分原理，将
求解连续函数的泊松方程的问题转换为求解网格节点上 的差分方程组的
问题。
1. 二维泊松方程的差分格式
二维静电场边值问题：
2
x 2
2
y 2
F
(1)
f (s)
(2)
L
通常将场域分成足够小的正方形网格， 网格线之间的距离为h，节点0,1,2,3,4上
的电位分别用0 ,1,和2 ,表3 示。4
设函数 在x0处可微 , 则沿x方向在 x0处的泰勒公式展开为
x
n (K )
Kn )
0
1 4
(1
2
3
4)
若场域离散为矩形网格, 差分格式为:
1•
2
1 h12
(1
2)
1 h2 2
( 2
4
)
(
1 h12
1 h2 2
)20
F
2.边界条件的离散化处理 ⑴第一类边界条件 给边界离散节点直接赋已知电位值。
⑵对称边界条件 合理减小计算场域, 差分格式为
•
0
1 4
(21
2
4
h2F)
⑶第二类边界条件 边界线与网格线相重合的差分格式:
(3)
将 x 和x1 分x别3 代入式(3),得
1
0
h(
x
)0
1 2!
h
2
(
2
x 2
)0
1 3!
h
3
(
3
x3
二维泊松方程很基础详细的求解过程
Topic 2: Elliptic Partial Differential EquationsLecture 2-4: Poisson’s Equation: Multigrid MethodsWednesday, February 3, 2010Contents1 Multigrid Methods2 Multigrid method for Poisson’s equation in 2-D3 Simple V −cycle algorithm4 Restricting the Residual to a Coarser Lattice 2 35 71 MULTIGRID METHODS5 Prolongation of the Correction to the Finer Lattice6 Cell-centered and Vertex-centered Grids and Coarsenings7 Boundary points8 Restriction and Prolongation Operators9 Improvements and More Complicated Multigrid Algorithms8 8 11 11 151 Multigrid MethodsThe multigrid method provides algorithms which can be used to accelerate the rate of convergence of iterative methods, such as Jacobi or Gauss-Seidel, for solving elliptic partial differential equations.Iterative methods start with an approximate guess for the solution to the differential equation. In each iteration, the difference between the approximate solution and the exact solution is made smaller.One can analyze this difference or error into components of different wavelengths, for example by using Fourier analysis. In general the error will have components of many different wavelengths: there will beshort wavelength error components and long wavelength error components.Algorithms like Jacobi or Gauss-Seidel are local because the new value for the solution at any lattice site depends only on the value of the previous iterate at neighboring points. Such local algorithms are generally more efficient in reducing short wavelength error components.The basic idea behind multigrid methods is to reduce long wavelength error components by updating blocks of grid points. This strategy is similar to that employed by cluster algorithms in Monte Carlo simulations of the Ising model close to the phase transtion temperature where long range correlations are important. In fact, multigrid algorithms can also be combined with Monte Carlo simulations.2 Multigrid method for Poisson’s equation in 2-DWith a small change in notation, Poisson’s equation in 2-D can be written:∂ 2u ∂x2 +∂ 2u∂y2= −f (x, y) ,where the unknown solution u(x, y) is determined by the given source term f (x, y) in a closed region. Let’s consider a square domain 0 ≤ x, y ≤ 1 with homogeneous Dirichlet boundary conditions u = 0 on the perimeter of the square. The equation is discretized on a grid with L + 2 lattice points, i.e., L interior points and 2 boundary points, in the x and y directions. At any interior point, the exact solution obeysu i,j = 14u i+1,j + u i−1,j + u i,j+1 + u i,j−1 + h2f i,j .The algorithm uses a succession of lattices or grids. The number of different grids is called the number of multigrid levels . The number of interior lattice points in the x and y directions is then taken to be 2 , so that L = 2 + 2, and the lattice spacing h = 1/(L − 1). L is chosen in this manner so that the downward multigrid iteration can construct a sequence of coarser lattices with2−1→ 2−2→ . . . → 20 = 1interior points in the x and y directions.Suppose that u(x, y) is the approximate solution at any stage in the calculation, and u exact(x, y) is the exact solution which we are trying to find. The multigrid algorithm uses the following definitions:· The correctionv = u exact− uis the function which must be added to the approximate solution to give the exact solution. · The residual or defect is defined asr = 2 u + f .Notice that the correction and the residual are related by the equation2 v = 2 u exact+ f − 2 u + f = −r .This equation has exactly the same form as Poisson’s equation with v playing the role of unknown function and r playing the role of known source function!3 SIMPLE V −CYCLE ALGORITHM3 Simple V −cycle algorithmThe simplest multigrid algorithm is based on a two-grid improvement scheme. Consider two grids:· a fine grid with L = 2 + 2 points in each direction, and· a coarse grid with L = 2−1 + 2 points.We need to be able to move from one grid to another, i.e., given any function on the lattice, we need to able to· restrict the function from fine → coarse, and· prolongate or interpolate the function from coarse → fine.Given these definitions, the multigrid V −cycle can be defined recursively as follows:· If = 0 there is only one interior point, so solve exactly foru1,1 = (u0,1 + u2,1 + u1,0 + u1,2 + h2f1,1)/4 .· Otherwise, calculate the current L = 2 + 2.3 SIMPLE V −CYCLE ALGORITHM· Perform a few pre-smoothing iterations using a local algorithm such as Gauss-Seidel. The idea is to damp or reduce the short wavelength errors in the solution.· Estimate the correction v = u exact− u as follows:– Compute the residualr i,j = 1h2 [u i+1,j + u i−1,j + u i,j+1 + u i,j−1− 4u i,j] + f i,j .–Restrict the residual r → R to the coarser grid.– Set the coarser grid correction V = 0 and improve it recursively.–Prolongate the correction V → v onto the finer grid.· Correct u → u + v.· Perform a few post-smoothing Gauss-Seidel interations and return this improved u. How does this recursive algorithm scale with L? The pre-smoothing and post-smoothing Jacobi or Gauss- Seidel iterations are the most time consuming parts of the calculation. Recall that a single Jacobi or Gauss-Seidel iteration scales like O(L2). The operations must be carried out on the sequence of grids with2 → 2−1→ 2−2→ . . . → 20 = 1interior lattice points in each direction. The total number of operations is of orderL2n=0122n≤ L211 .1 − 44 RESTRICTING THE RESIDUAL TO A COARSER LATTICEThus the multigrid V −cycle scales like O(L 2), i.e., linearly with the number of lattice points N!4Restricting the Residual to a Coarser LatticeThe coarser lattice with spacing H = 2h is constructed as shown. A simple algorithm for restricting the residual to the coarser lattice is to set its value to the average of the values on the four surrounding lattice points (cell-centered coarsening):6 CELL-CENTERED AND VERTEX-CENTERED GRIDS AND COARSENINGSR I,J = 14[r i,j + r i+1,j + r i,j+1 + r i+1,j+1] , i = 2I − 1 , j = 2J − 1 .5 Prolongation of the Correction to the Finer LatticeHaving restricted the residual to the coarser lattice with spacing H = 2h, we need to solve the equation2 V = −R(x, y) ,with the initial guess V (x, y) = 0. This is done by two-grid iterationV = twoGrid(H, V, R) .The output must now be i nterpolated or prolongated to the finer lattice. The simplest procedure is to copy the value of V I,J on the coarse lattice to the 4 neighboring cell points on the finer lattice: v i,j = v i+1,j = v i,j+1 = v i+1,j+1 = V I,J , i = 2I − 1, j = 2J − 1 .6 Cell-centered and Vertex-centered Grids and CoarseningsIn the cell-centered prescription, the spatial domain is partitioned into discrete cells. Lattice points are defined at the center of each cell as shown in the figure:The coarsening operation is defined by doubling the size of a cell in each spatial dimension and placing a coarse lattice point at the center of the doubled cell.Note that the number of lattice points or cells in each dimension must be a power of 2 if the coarsening operation is to terminate with a single cell. In the figure, the finest lattice has 23 = 8 cells in each dimension, and 3 coarsening operations reduce the number of cells in each dimension23= 8 → 22= 4 → 21= 2 → 20 = 1 .Note also that with the cell-centered prescription, the spatial location of lattice sites changes with each coarsening: coarse lattice sites are spatially displaced from fine lattice sites.A vertex-centered prescription is defined by partitioning the spatial domain into discrete cells and locating the discrete lattice points at the vertices of the cel ls as shown in the figure:The coarsening operation is implemented simply by dropping every other lattice site in each spatial dimension. Note that the number of lattice points in each dimension must be one greater than a power of 2 if the coarsening operation is to reduce the number of cells to a single coarsest cell. In the example in the figure the finest lattice has 23 + 1 = 9 lattice sites in each dimension, and 2 coarsening operations reduce the number of vertices in each dimension23+ 1 = 9 → 22+ 1 = 5 → 21 + 1 = 3 .The vertex-centered prescription has the property that the spatial locations of the discretization points are not changed by the coarsening operation.8 RESTRICTION AND PROLONGATION OPERATORS7 Boundary pointsLet’s assume that the outermost perimeter points are taken to be the boundary points. The behavior of these boundary points is different in the two prescriptions:· Cell-centered Prescription: The boundary points move in space towards the center of the region at each coarsening. This implies that one has to be careful in defining the “boundary values” of the solution.· Vertex-centered Prescription: The boundary points do not move when the lattice is coarsened.This make i t easier in principle to define the boundary values.These two different behaviors of the boundary points make the vertex-centered prescription a little more convenient to use in multigrid applications. However, there is no reason why the cell-centered prescription should not work as well.8 Restriction and Prolongation OperatorsIn the multigrid method it is necessary to move functions from a fine grid to the next coarser grid (Restric- tion), and from a coarse grid to the next finer grid (Prolongation). Many prescriptions for restricting andprolongating functions have been studied. Let’s consider two of the simplest prescriptions appropriate for cell- and vertex-centered coarsening:· Cell-centered Coarsening: In this prescription, a coarse lattice point is naturally associated with 2d neighboring fine lattice points in d-dimensions.· Suppose that f (x) is a function on the fine lattice at spatial position x, and F (X ) is the corresponding function on the coarse lattice, then this diagram suggests a simple prescription for restriction and prolongation.–Restriction: Average the function values at the 4 neighboring fine lattice sites x i:F (X ) = 144i=1f (x i) .– Prolongation: Inject the value of the function at the coarse lattice site to the 4 neighboring fine lattice sites:f (x i) = F (X ) , i = 1 . . . 4· Vertex-centered Coarsening: Consider a coarse lattice point and the 9 neighboring fine lattice points shown in the figure:2 2 · In this prescription, a coarse lattice point can naturally associated (in 2-D) with · the corresponding fine lattice point, or· the four nearest neighbor fine lattice points, left, right, up, and down, or · with the four diagonally nearest fine lattice points, etc.· It is a little more complicated here to define transfer operators. The problem is that the fine lattice points are associated with more than one coarse lattice point, unlike the cell-centered case: – The single red fine lattice point in the center coincides with an unique coarse lattice point. – Each of the 4 black fine lattice points however is equidistant from two coarse lattic e points. – Each of the 4 red fine lattice points is equidistant from four coarse lattice points. · This sharing of lattice points suggests the following prescriptions:· Prolongation: use bilinear interpolation in which the value of F at a coarse grid point is copied to 9 neighboring fine -grid points with the following weights:1 1 14 2 1 1 1 1 4 1 1.4 24This matrix is called the stencil for the prolongation.8 8 9 IMPROVEMENTS AND MORE COMPLICATED MULTIGRID ALGORITHMS· Restriction: The restriction operator is taken to be the adjoint of the prolongation operator:1 1 116 1 1 16 8 1 4 1 8 161 1 16.This choice of restriction operator is called full weighting.9Improvements and More Complicated Multigrid AlgorithmsThe algorithm implemented above is the simplest multigrid scheme with a single V-cycle. Section 19.6 of Numerical Recipes discusses various ways of improving this algorithm:· One can repeat the two-grid iteration more than once. If it is repeated twice in each multigrid level one obtains a W-cycle type of algorithm.· The Full Multigrid Algorithm starts with the coarsest grid on which the equation can be solved exactly. It then proceeds to finer grids, performing one or more V -cycles at each level along the way. Numerical Recipes gives a program mglin(u,n,ncycle) which accepts the source function −f in the first argument and implements the full multigrid algorithm with = log 2(n − 1) levels, performing ncycle V-cycles at each level, and returning the solution in the array parameter u. Note that this program assumes that the number of lattice points in each dimension L is odd, which leads to vertex centered coarsening:REFERENCESREFERENCESReferences[Recipes-C19-5] W.H. Press, S.A. Teukolsky, W. Vetterling, and B.P. Flannery, “Numerical Recipes in C”,Chapter 19 §6: Multigrid Methods for Boundary Value Problems, /a/bookcpdf/c19-6.pdf.。
二维泊松方程的有限元法
忠
2．单元网格划分
在二维情况下，单元可以是三角形和四边形。
具体要求是，三角形顶点连着顶点，
三角形的三条边长尽量接近
或三个内角尽量接近。
图示三角形的三个顶点，
i, j, k 的顺序按逆时针。
18／4／25
华北电力大学电气与电子工程学院
4
工程电磁场
主讲人： 王泽 忠
18／4／25
华北电力大学电气与电子工程学院
18／4／25
华北电力大学电气与电子工程学院
7
工程电磁场
主讲人： 王泽 忠
将 1 、 2 、 3 作为未知数，
求解上述方程组，并令
aaij
x j yk xk yi
xk yj xi yk
a
k
xi y j
x jyi
18／4／25
华北电力大学电气与电子工程学院
8
工程电磁场
1 xi
11 2
设第k 个节点是第一类边界上的节点，
其电位已知k k 0。
18／4／25
华北电力大学电气与电子工程学院
45
工程电磁场
主讲人： 王泽
忠
在总体系数矩阵和右端向量中，做如下处理：
（1） Akk 1；
（2） Ri Ri Aikk 0 （ i 1,2, , n ）；
（3） Rk k 0 ；
（4） Akj 0 （ i 1,2, , n ）；
N j • d N j ( )d N jd
e1 e
es1 es
e1 e
ne
nes
ne
Nk • d Nk ( )d Nk d
e1 e
es1 es
e1 e
18／4／25
二维泊松方程的求解公式
二维泊松方程的求解公式一、二维泊松方程的基本形式。
在二维直角坐标系中，泊松方程的形式为：∇^2u = f(x,y)其中，∇^2=frac{∂^2}{∂ x^2}+frac{∂^2}{∂ y^2}是二维拉普拉斯算子，u = u(x,y)是待求的函数，f(x,y)是已知的源函数。
二、求解方法之分离变量法（在特定边界条件下）1. 假设解的形式。
- 设u(x,y)=X(x)Y(y)，将其代入泊松方程frac{∂^2u}{∂ x^2}+frac{∂^2u}{∂y^2}=f(x,y)，得到：- X''(x)Y(y)+X(x)Y''(y)=f(x,y)- 两边同时除以X(x)Y(y)，得到：- (X''(x))/(X(x))+(Y''(y))/(Y(y))=(f(x,y))/(X(x)Y(y))2. 令(X''(x))/(X(x)) = - λ，(Y''(y))/(Y(y))=λ - (f(x,y))/(X(x)Y(y))（这里λ为分离常数）- 对于X(x)，我们得到方程X''(x)+λ X(x) = 0，其解的形式取决于λ的值。
- 当λ>0，设λ = k^2，k>0，则X(x)=Acos(kx)+Bsin(kx)。
- 当λ = 0，X(x)=Ax + B。
- 当λ<0，设λ=-k^2，k>0，则X(x)=Ae^kx+Be^-kx。
- 对于Y(y)，同样根据λ的值求解相应的二阶常微分方程。
3. 确定系数。
- 根据给定的边界条件（例如，在矩形区域0≤ x≤ a，0≤ y≤ b上的边界条件），确定系数A、B等。
三、格林函数法求解二维泊松方程。
1. 格林函数的定义。
- 对于二维泊松方程∇^2u = f(x,y)在区域Ω内，格林函数G(x,y;x_0,y_0)满足：- ∇^2G(x,y;x_0,y_0)=δ(x - x_0)δ(y - y_0)在Ω内，- G(x,y;x_0,y_0)=0在∂Ω上（∂Ω为区域Ω的边界）。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
9
工程电磁场
主讲人： 王泽忠
2．单元网格划分 在二维情况下，以三角形单元为例 网格划分就是把求解区域划分成有限个三角形。 具体要求是，三角形顶点连着顶点， 三角形三条边长或三个内角大小尽量接近。 图 显示了网格的一部分。 图 表示一个三角形的三个顶点，
2019/10/3
华北电力大学电气与电子工程学院
相应的待定常数为
u1, u2 , , un , unn
2019/10/3
华北电力大学电气与电子工程学院
3
工程电磁场
主讲人： 王泽忠
以 n 表示基函数序列通项的序号， nn 表示总项数。 u 的近似解（试探函数）表示为
nn
u M n (x, y)un
n 1
在伽辽金加权余量法中，权函数序列：
5
工程电磁场
主讲人： 王泽忠
代入第二类边界条件，得
aM m ud bM md M m f d
（ m 1, 2, , nn ）
将近似函数（试探函数）代入，得
nn
aM m ( M nun )d bMmd Mm f d
n 1
Ae ，Re ， Reb
2019/10/3
华北电力大学电气与电子工程学院
14
工程电磁场
主讲人： 王泽忠
单元系数矩阵和单元右端项的元素为
Ae,i, j (aNi N j )de
（ m 1, 2, , nn ）
以下为了书写方便，将 Mm (x, y) 写为 M m 。
对上述方程组应用格林公式，得
u
aM m
ud
M ma n d M m f d
（ m 1, 2, , nn ）
2019/10/3
华北电力大学电气与电子工程学院
7
工程电磁场
主讲人： 王泽忠
根据运算规则，将先求和在积分变为线积分在求和，
进一步整理得方程组
nn
[(a M m Mnd)un ] M m f d bM md
n1
（ m 1, 2, , nn ）
观察等号左右各项，可知
这是关于 u1, u2 , , un , unn 的代数方程组。 将 u1,u2, ,un , unn 写成列向量，
2019/10/3
华北电力大学电气与电子工程学院
2
工程电磁场
主讲人： 王泽忠
a 方程系数（表示材料参数），
u0 是第一类边界上已知位函数， b 是第二类边界上
已知场矢量的法向投影（表示边界上的源）。 设基函数序列为
M1(x, y), M 2 (x, y), , M n (x, y), , M nn (x, y)
（ m 1, 2, , nn ）
根据运算规则，将梯度运算移到求和运算之内，得
2019/10/3
华北电力大学电气与电子工程学院
6
工程电磁场
主讲人： 王泽忠
nn
aM m ( M nun )d bM md Mm f d
n 1
（ m 1, 2, , nn ）
10
工程电磁场
主讲人： 王泽忠
n1, n2 , n3 按逆时针排列。
2019/10/3
华北电力大学电气与电子工程学院
11
工程电磁场
主讲人： 王泽忠
3．单元系数矩阵和单元右端向量
进行离散化处理，
将整个区域的积分化为单元上的积分之和
ne
Am,n
(aM m M n )de
e1 e
13
工程电磁场
主讲人： 王泽忠
ne
Am,n
(aNi N j )de
e1 e
ne
neb
Rm
( Ni f )de
(bNi )de
e1 e
e 1 e
一次单元积分可以算出若干系数矩阵元素的单元贡献
和若干右端项元素的单元贡献
将其写成单元系数矩阵和单元右端向量
ne
neb
Rm
(M m f )de
(bM m )de
e1 e
e 1 e
式中， e 是单元编号， ne 华北电力大学电气与电子工程学院
12
工程电磁场
主讲人： 王泽忠
neb 是第二类边界上线段单元的总数。
在单元上，为将基函数、权函数更换为单元形状函数，
先将基函数和权函数表示为双重下标
ne
Am,n
(aM me,i M ne,j )de
e1 e
ne
neb
Rm
(M me ,i f )de
(bM me,i )de
e1 e
e 1 e
基函数和权函数可以用单元形状函数代替。表示为
2019/10/3
华北电力大学电气与电子工程学院
M1(x, y), M 2 (x, y), , M m (x, y), , Mnn (x, y)
以 m 表示权函数序列通项的序号， nn 表示总项数。
2019/10/3
华北电力大学电气与电子工程学院
4
工程电磁场
主讲人： 王泽忠
代入伽辽金加权余量方程，得如下方程组
M m (x, y)(a2u)d M m (x, y) f d
工程电磁场
主讲人： 王泽忠
工程电磁场
王泽忠
2019/10/3
华北电力大学电气与电子工程学院
1
工程电磁场
主讲人： 王泽忠
7.3 二维泊松方程的有限元法
1．二维泊松方程的伽辽金离散化
设位函数 u 满足泊松方程 a2u f
边界条件
u
1
u0
u
a n
2
b
式中， f 是已知函数（表示场源分布），
2019/10/3
华北电力大学电气与电子工程学院
8
工程电磁场
主讲人： 王泽忠
方程组写成矩阵的形式
A{u} {R}
式中，系数矩阵和右端列向量元素表达式
Am,n a M m M nd
Rm M m f d bM md
2019/10/3
华北电力大学电气与电子工程学院
求和运算的下标是 n ， M m 的下标 m 与求和无关，
可以将 M m 拿到求和运算之内。整理后得
nn
[(aM m M n )un ]d bM md M m f d
n 1
（ m 1, 2, , nn ）
2019/10/3
华北电力大学电气与电子工程学院