分布式离散中多重网格预处理求解微分方程

打开实时脚本

本例展示了如何使用预条件迭代求解器和分布式数组求解泊松方程。通过使用分布式阵列，可以使用机器集群的内存来扩展计算，而不仅仅是单个机器的内存。您可以在不更改代码的情况下进行扩展。

本例继续介绍的主题用分布式阵列迭代法求解线性方程组．在[1]的基础上，本例通过使用泊松方程来模拟房间内的热分布，其形式被称为齐次稳态热方程。稳态意味着热量不随时间变化，均匀意味着没有外部热源。

$- Δ u ＝ - （ \frac{\partial^{2} u}{\partial x^{2}} + \frac{\partial^{2} u}{\partial y^{2}} + \frac{\partial^{2} u}{\partial z^{2}} ）＝ 0$

在方程中， $u$ 表示每一点的温度 $（ x ， y ， z ）$ 房间的。要解这个方程，首先要用有限差分离散化方法用线性方程组近似它。然后，使用预条件共轭梯度(pcg)方法求解方程组。预处理对问题进行转换，以提高数值求解器的性能。通过使用分布式阵列，您可以利用机器集群的组合内存，并允许更精细的离散化。

学习如何:

建立离散三维网格和边界条件。
定义一个离散化和一个多重网格预处理。
应用预处理数值求解器求解三维体积的热方程。

离散空间维度

在这个例子中，一个立方体的边1给房间做模型。第一步是用三维网格离散化它。

本例中的预处理方法使用了几个具有不同粒度级别的网格。每一层都将网格粗化一倍2在每个维度。定义多重网格级别的数量。

multigridLevels = 2;

定义每个维度的点数，X，Y,Z，最好的网格。预处理方法要求网格中的点数能被整除2 ^ multigridLevels．在这种情况下，点的数量必须能被4整除，因为多重网格级别的数量为2．

numPoints。X = 32;numPoints。Y = 32;numPoints。Z = 32;

控件对三维网格的空间维度进行离散化meshgrid函数。按点数均匀划分各维，用linspace．请注意，要包括立方体的边界，必须添加两个附加点。

[X,Y,Z] = meshgrid(linspace(0,1,numPoints.X+2)，.．.linspace (0, 1, numPoints.Y + 2),.．.linspace (0, 1, numPoints.Z + 2));

定义边界条件

假设这个房间有一扇窗和一扇门。墙壁和天花板的恒温为0度，窗户的恒温为16度，门的恒温为15度。地板的温度是0.5度。目标是确定整个房间内部的温度分布。

使用关系操作符定义地板、窗户和门的坐标，并在这些边界元素上定义温度。边界是立方体的面，因此是其中之一X，Y,或Z必须0或1．将立方体的其余边界和内部设置为0．

地板= (0.0 < = X和X < = 1.0)和(0.0 < = Y & Y < = 1) & (Z = = 0.0);窗口= (X = = 1) & (0.2 < = Y & Y < = 0.8)和(0.4 < = & Z < = 0.6);门= (0.4 < = X & X < = 0.6) & (Y = = 1.0)和(0.0 < = & Z < = 0.6);u = 0(大小(X));U(楼层)= 0.5;U(窗口)= 16;U(门)= 15;

这些边界条件指定了解必须沿着域边界取的常数值。这种类型的边界条件被称为狄利克雷边界条件。

图形化边界条件片函数。使用位于立方体边界的切片来显示非零边界条件。

xSlices = 1;ySlices = 1;zSlices = 0;f = slice(X,Y,Z,u,xSlices,ySlices,zSlices，“最近的”）;标题(“常数非零边界条件”),包含(“x”), ylabel (“y”), zlabel (“z”）;colorbar, colormap很酷的；阴影插值函数；集(f,“EdgeColor”，[0 0 0]);

微分方程的离散化与求解

本例利用有限差分逼近法将微分方程离散为线性系统，并使用多重网格预处理来提高迭代求解器的性能。对于本例，在支持函数中使用离散化和预处理金宝appdiscretizePoissonEquation而且multigridPreconditioner．对于其他问题，选择适合您的应用程序的离散化和预处理。

在这个例子中，discretizePoissonEquation采用7点模板有限差分法将泊松方程离散为多个不同粒度的网格。该函数创建了离散化的多重网格结构，具有预先计算的三角分解和在粗级和细级之间映射的操作符。预条件器使用该多重网格信息以有效的方式近似解。

在其他技术中，这个预处理应用平滑来最小化误差与一系列的逼近。定义平滑步骤的数量。使用更多的步骤可以使近似更准确，但也需要更多的计算。然后对微分方程进行离散化，建立预条件。

numberOfSmootherSteps = 1;[A,b,multigridData] =离散poissonequation (numPoints,multigridLevels,numberOfSmootherSteps,u);

级别0:问题的维度是32768，有223232个非零。第1级:问题的维数是4096，有27136个非零。第2级:问题的维数是512，有3200个非零。

pre冷气= setuppre冷气(multigridData);

利用预条件共轭梯度求解线性方程组。

Tol = 1e-12;maxit = number (A);pcg (A, b,托尔,麦克斯特,预调节器);

PCG在迭代45时收敛到相对残差为5.4e-13的解。

使用分布式阵列进行扩展

如果需要更多的计算资源(比如内存)，可以使用分布式数组进行扩展，而不需要更改代码。分布式数组将数据分布到多个工作线程上，它们可以利用机器集群的计算性能和内存。

启动一个并行工作池。默认情况下,parpool使用默认集群。在MATLAB中检查默认集群配置文件首页选项卡，在环境区,在平行>选择默认集群．

parpool;

启动并行池(parpool)使用'MyCluster'配置文件…连接到并行池(工人数:12)。

分布温度变量u方法在集群中的工作内存中调用分布式函数。

distU =分布式(u);

您可以使用与以前相同的代码;不需要更改，因为如果输入是分布式数组，则离散化和预处理函数将创建分布式数组。许多MATLAB函数都针对分布式数组进行了增强，因此您可以使用与内存中数组相同的方式使用它们。

请注意,discretizePoissonEquation返回包含分布式数据的结构。若要以分布式方式在结构中使用分布式数据，必须在对象中创建结构spmd块。类中使用它的任何函数也必须调用spmd块。

spmd[A,b,multigridData] = discretizePoissonEquation(numPoints,multigridLevels,numberOfSmootherSteps,distU);pre冷气= setuppre冷气(multigridData);结束

分析文件并将文件传输给工作人员…完成。实验1:0级:问题的维度是32768，有223232个非零。第1级:问题的维数是4096，有27136个非零。第2级:问题的维数是512，有3200个非零。

使用pcg在一个spmd以分布式的方式求解线性系统的块。

spmdx = pcg(A,b,tol,maxit, preconditioning);结束

实验1:pcg在迭代45时收敛到相对残差5.4e-13的解。

阴谋的结果

求解器的解是一个适合内存的向量。通过使用将工人的数据发送到客户端收集．将数据重新塑造成3-D数组并重新排列尺寸以产生最终解决方案。设置内部的部分u对这个解。外部部分，即边界，已经包含了边界条件的值。

x3D =重塑(gather(x)，numPoints.X,numPoints.Y,numPoints.Z);u (2: end-1, 2: end-1, 2: end-1) =排列(x3D (2, 1, 3));

方法可视化解决方案片函数。添加额外的切片来绘制立方体内部的温度。您可以使用旋转工具，或改变切片的位置，以探索解决方案。

xSlices = [.5,1];ySlices = [.5,1];zSlices = [0，.5];f = slice(X,Y,Z,u,xSlices,ySlices,zSlices，“最近的”）;标题(的热量分布),包含(“x”), ylabel (“y”), zlabel (“z”）;colorbar, colormap很酷的；阴影插值函数；

你可以尝试不同的值numPoints在本例中测试不同级别的离散化。使用较大的值会增加分辨率，但需要更多内存。此外，规模较大multigridLevels是，预调节器的内存效率就越高。然而，一个更大的multigridLevels暗示一个不太准确的预处理条件，因为粗化会降低每个级别的准确性。因此，求解器可能需要更多的迭代才能达到相同的精度水平。

定义预调节器

定义一个用于预处理共轭梯度方法的多网格预处理。这种类型的预处理使用几个具有不同粒度级别的离散网格来更有效地近似线性方程组的解。本例中的预处理方法基于[2]，并遵循以下主要阶段:

采用高斯-塞德尔近似方法进行预平滑。
在较粗的水平上计算剩余解。
递归的前提条件在一个较粗的级别，或直接解决如果在最粗的级别。
用更粗的网格解决方案更新解决方案。
Postsmooth使用高斯-塞德尔近似方法。

函数x = multigridPreconditioner(mgData,r,level)如果(level < mgData(level).MaxLevel) x = 0 (size(r)，“喜欢”, r);使用高斯-赛德尔预平滑为i = 1: mgData(水平)。NumberOfSmootherSteps x = mgData(level). matrices。LD \ (-mgData(level). matrices。U*x + r);x = mgData(level). matrices。DU \ (-mgData(level).矩阵。L*x + r);结束%在较粗的级别上计算残差Axf = mgData(level).Matrices.A*x;rc = r(mgData(level).Fine2CoarseOperator)- Axf(mgData(level).Fine2CoarseOperator);%递归调用，直到达到最粗糙的级别xc = multigridPreconditioner(mgData, rc, level+1);%更新溶液与延长粗网格溶液x(mgData(level).Fine2CoarseOperator) = x(mgData(level).Fine2CoarseOperator)+xc;% Postsmooth使用Gauss-Seidel为i = 1:mgData(level)。NumberOfSmootherSteps x = mgData(level). matrices。LD \ (-mgData(level). matrices。U*x + r);x = mgData(level). matrices。DU \ (-mgData(level).矩阵。L*x + r);结束其他的在最粗糙的水平上获得精确的溶液。x = mgData(level). matrices。A \ r;结束结束

创建一个函数，该函数接受多重网格数据并返回一个函数句柄，该句柄将预处理条件应用于输入数据。在本例中，此函数句柄是的前置条件输入pcg．必须创建此函数，因为不可能在其中定义匿名函数spmd块。

函数pre冷气= setuppre冷气(multigridData)如果~isempty(multigridData)预处理条件= @(x,varargin) multigrid预处理条件(multigridData,x,1);其他的预处理条件= [];结束结束

参考文献

[1]唐加拉，J. M. A. Heroux, P. Luszczek。“HPCG基准:高性能计算系统排名的新度量。”诺克斯维尔，田纳西州:田纳西大学，2015.

[2]埃尔曼，H. C.， D. J.西尔维斯特，A. J.瓦特。有限元和快速迭代求解:在不可压缩流体动力学中的应用。牛津，英国:牛津大学出版社，2005年，第2.5节。

另请参阅

pcg|分布式|spmd|收集