fmgv1d.m

数值计算中的多重网格计算程序 不过还需要做一些调试工作 基本思路是正确的

% Matlab code fmgv1D.m
% For "Applied Numerical Linear Algebra",  Question 6.16
% Written by James Demmel, Apr 22, 1995
%                Modified, Jun  2, 1997
%
% Full Multigrid V-Cycle for Poisson's equation on a 1D grid
% with zero Dirichlet boundary conditions.
% Algorithm from Brigg's ``Multigrid Tutorial''
% Include zero boundary values in arrays for ease of programming.
% Assume dimension n = 2^k + 1
% Inputs:
%   x = initial guess (n by 1 matrix with zeros on boundary)
%   b = right hand side (n by 1 matrix)
%   jac1, jac2 = number of weight Jacobi steps to do before and
%                after recursive call to mgv
% Outputs:
%   z = improved solution (n by n matrix with zeros on boundary)
%
% function z=fmgv1D(x,b,jac1,jac2);
%
function z=fmgv1D(x,b,jac1,jac2);
[nn,m]=size(b);
%  Compute residual
   r=zeros(nn,1);
   tmp = b(2:nn-1) - ( 2*x(2:nn-1) - x(1:nn-2) - x(3:nn) );
   r(2:nn-1) = tmp;
% Get right hand side for coarsest grid by repeated restriction of
% original right hand side to coarser grids
rc=mgvrhs1D(nn,3,r);
% Solve 1 by 1 problem
z=zeros(3,1);
% z(2)=rc(2)/4;
z(2)=rc(2)/2;
%  Assume nn = 2^k+1
k = round(log(nn-1)/log(2)); 
% Loop from the next-to-coarsest to finest grids
for i=2:k,
   n=2^i+1;
   m=2^(i-1)+1;
%  Interpolate coarser solution to this level
   zstrt=zeros(n,1);
   zstrt(3:2:n-2)=z(2:m-1);
   zstrt(2:2:n-1)=.5*(z(1:m-1)+z(2:m));
% Get right hand side for current grid by repeated restriction of
% original right hand side to coarser grids
   rhs=mgvrhs1D(nn,n,r);
%  Do Multigrid V-cycle
   z=mgv1D(zstrt,rhs,jac1,jac2);
end
%  plot(z),title(int2str(n)), pause
% Add final correction to initial guess
z=x+z;
return