kac.m

来自「这是在网上下的一个东东」· M 代码 · 共 70 行

%
% x=kac(A,b,tolx,maxiter)
%
% Inputs:
%      A        Constraint matrix.
%      b        right hand side.
%      tolx     Terminate when the relative diff between two iterations
%               is less than tolx.
%      maxiter  Stop after maxiter iterations.
%
% Outputs:
%      x        Solution.
%
function x=kac(A,b,tolx,maxiter)
%
% First, find the size of the matrix.
%
[m,n]=size(A);
%
% Make a copy of A' to speed up some accesses.
%
AP=A';
%
% Setup an initial solution of all zeros.
%
x=zeros(n,1);
iter=0;
%
%  Precompute the row norms squared.
%
n2=zeros(m,1);
for i=1:m
  n2(i)=norm(AP(:,i),2)^2;
end
%
%  The main loop performs iterations of Kaczmarz algorithm until 
%  maxiters is exceeded or successive iterates differ by less 
%  than tolx.  
%
while (iter <= maxiter)
%
% Update the iteration count.
%
  iter=iter+1;
%
%  Start the update cycle with the current solution.
%
  newx=x;
%
%  Perform a cycle of m updates.
%
  for i=1:m
    newx=newx-((newx'*AP(:,i)-b(i))/(n2(i)))*AP(:,i);
  end
%
%  Check for convergence to fixed solution.
%
  if (norm(newx-x)/(1+norm(x)) < tolx)
    x=newx;
    return;
  end
%
%  Update x for the next major iteration.
%
  x=newx;
end

disp('Max iterations exceeded.');
return;