mb_45方法.m

求解大规模对称矩阵的CG算法

clear all;
format long;  %修改默认数据类型
tic;   %开始计时
n=1000;   %矩阵规模
e=0.00000001;  %精度
A=hilb(n);   %生成系数矩阵
K=cond(A);  %条件数
M=zeros(n);
sum=0;
for i=1:1:45  %一共45个子矩阵
   for j=1:1:i
       for k=1:1:i
           line=j+sum;
           low=k+sum;            
           if(line<=n)&&(low<=n)
               M(line,low)=A(line,low);
           end
       end
   end
   sum=sum+i;
end
x=ones(n,1);   %生成精确解
b=A*x;   %生成精确解对应的b
X=zeros(n,1);  %初始
ek1=sqrt((X-x)'*A*(X-x)); %开始是不会知道精确解的！
sulv=0;  % 评价指标
r=b-A*X;
ni=inv(M);
z=ni*r;
p=z;
k=1;%迭代次数
fan=norm(r);
while (k<=100)&&(fan>e)
zr=r'*z;
afa=zr/(p'*(A*p));
X=X+afa*p;
ek2=sqrt((X-x)'*A*(X-x));
sulv=sulv+ek2/ek1;
ek1=ek2;
r=r-afa*(A*p);
z=ni*r;
fan=norm(r);
beta=r'*z/zr;
p=z+beta*p;
k=k+1;
end
en=norm(X-x);
sulv=sulv/(k-1);
B=A*X;
toc;