www.gusucode.com > pde 案例源码 matlab代码程序 > pde/DomainDecompositionProblemExample.m

    %% Domain Decomposition Problem
% This example shows how to perform one-level domain decomposition for
% complicated geometries, where you can decompose this geometry into the
% union of more subdomains of simpler structure. Such structures are often
% introduced by the PDE app.
%
% Assume now that $\Omega$ is the disjoint union of some subdomains
% $\Omega_1, \Omega_2, \ldots, \Omega_n$. Then you could renumber the nodes
% of a mesh on $\Omega$ such that the indices of the nodes of each
% subdomain are grouped together, while all the indices of nodes common to
% two or more subdomains come last. Since |K| has nonzero entries only at
% the lines and columns that are indices of neighboring nodes, the
% stiffness matrix is partitioned as follows:
%
% $$ K = 
% \left(
% \begin{array}{ccccc}
% K_1 & 0 & \cdots & 0 & B_1^T \\
% 0 & K_2 & \cdots & 0 & B_2^T \\
% \vdots & \vdots & \ddots & \vdots & \vdots \\
% 0 & 0 & \cdots & 0 & B_n^T \\
% B_1 & B_2 & \cdots & B_n & C \\
% \end{array}
% \right) $$
%
% while the right side is
%
% $$ F =
% \left(
% \begin{array}{c}
% f_1 \\
% f_2 \\
% \vdots \\
% f_n \\
% f_c \\
% \end{array}
% \right) $$
%
% The Partial Differential Equation Toolbox(TM) function |assempde| can
% assemble the matrices $K_j$, $B_j$, and $C$ separately. You have full
% control over the storage and further processing of these matrices.
%
% Furthermore, the structure of the linear system
%
% $$Ku = F$$
%
% is simplified by decomposing $K$ into the partitioned matrix.
%
% Now consider the geometry of the L-shaped membrane. You can plot the
% geometry of the membrane by typing
%
%  pdegplot('lshapeg')
%%
% Notice the borders between the subdomains. There are three subdomains.
% Thus the matrix formulas with $n = 3$ can be used. Now generate a mesh
% for the geometry:
[p,e,t] = initmesh('lshapeg'); 
[p,e,t] = refinemesh('lshapeg',p,e,t); 
[p,e,t] = refinemesh('lshapeg',p,e,t);
%%
% So for this case, with $n = 3$, you have
%
% $$ \left(
% \begin{array}{cccc}
% K_1 & 0 & 0 & B_1^T \\
% 0 & K_2 & 0 & B_2^T \\
% 0 & 0 & K_3 & B_3^T \\
% B_1 & B_2 & B_3 & C \\
% \end{array}
% \right)
% \left(
% \begin{array}{c}
% u_1 \\
% u_2 \\
% u_3 \\
% u_c \\
% \end{array}
% \right) = 
% \left(
% \begin{array}{c}
% f_1 \\
% f_2 \\
% f_3 \\
% f_c \\
% \end{array}
% \right) $$
%
% and the solution is given by block elimination:
%
% $$ \left( C - B_1 K_1^{-1}B_1^T - B_2 K_2^{-1} B_2^T - B_3 K_3^{-1} B_3^T
% \right) = f_c - B_1 K_1^{-1} f_1 - B_2 K_2^{-1} f_2 - B_3 K_3^{-1} f_3 $$
%
% $$ u_1 = K_1^{-1} \left( f_1 - B_1^T u_c \right) $$
%
% $$ \cdots $$
%%
% In the following MATLAB(TM) solution, a more efficient algorithm using
% Cholesky factorization is used:
time = [];
np = size(p,2);
% Find common points
c = pdesdp(p,e,t);

nc = length(c);
C = zeros(nc,nc);
FC = zeros(nc,1);

[i1,c1] = pdesdp(p,e,t,1);
ic1 = pdesubix(c,c1);
[K,F] = assempde('lshapeb',p,e,t,1,0,1,time,1);
K1 = K(i1,i1);
d = symamd(K1);
i1 = i1(d); 
K1 = chol(K1(d,d));
B1 = K(c1,i1);
a1 = B1/K1;
C(ic1,ic1) = C(ic1,ic1)+K(c1,c1)-a1*a1';
f1 = F(i1);
e1 = K1'\f1;
FC(ic1) = FC(ic1)+F(c1)-a1*e1;

[i2,c2] = pdesdp(p,e,t,2);
ic2 = pdesubix(c,c2); 
[K,F] = assempde('lshapeb',p,e,t,1,0,1,time,2);
K2 = K(i2,i2);
d = symamd(K2);
i2 = i2(d); 
K2 = chol(K2(d,d));
B2 = K(c2,i2);
a2 = B2/K2;
C(ic2,ic2) = C(ic2,ic2)+K(c2,c2)-a2*a2';
f2 = F(i2);
e2 = K2'\f2;
FC(ic2) = FC(ic2)+F(c2)-a2*e2;

[i3,c3] = pdesdp(p,e,t,3);
ic3 = pdesubix(c,c3);
[K,F] = assempde('lshapeb',p,e,t,1,0,1,time,3);
K3 = K(i3,i3);
d = symamd(K3);
i3 = i3(d);
K3 = chol(K3(d,d));
B3 = K(c3,i3);
a3 = B3/K3;
C(ic3,ic3) = C(ic3,ic3)+K(c3,c3)-a3*a3'; 
f3 = F(i3);
e3 = K3'\f3;
FC(ic3) = FC(ic3)+F(c3)-a3*e3;

% Solve 
u = zeros(np,1);
u(c) = C\ FC;
u(i1) = K1\(e1-a1'*u(c1));
u(i2) = K2\(e2-a2'*u(c2));
u(i3) = K3\(e3-a3'*u(c3));
%%
% The problem can also be solved by typing

% Compare with solution not using subdomains
[K,F] = assempde('lshapeb',p,e,t,1,0,1);
u1 = K\F;
norm(u-u1,'inf')
pdesurf(p,t,u)
%%
% For the entire example, see
% <http://www.mathworks.com/help/pde/examples/poisson-s-equation-using-domain-decomposition.html Poisson's Equation Using Domain Decomposition>.