线性方程组的迭代解法——共轭梯度法

Posted 2020-11-26 guliangt

　　1.代码

%%共轭梯度法（用于求解正定对称方程组）
%%线性方程组M*X = b，M是方阵，X0是初始解向量,epsilon是控制精度
function CGM = Conjugate_gradient_method(M,b,X0,epsilon)
m = size(M);up = 1000;e = floor(abs(log(epsilon)));
X(:,1) = X0;
r(:,1) = b-M*X0;p(:,1) = r(:,1);
for k = 1:up
    alpha = Inner_product(r(:,k),r(:,k))/Inner_product(p(:,k),M*p(:,k));
    X(:,k+1) = X(:,k)+alpha*p(:,k);
    r(:,k+1) = r(:,k)-alpha*M*p(:,k);
    beta(:,k) = Inner_product(r(:,k+1),r(:,k+1))/Inner_product(r(:,k),r(:,k));
    p(:,k+1) = r(:,k+1)+beta(:,k)*p(:,k);
    X_delta(:,k) = X(:,k+1)-X(:,k);
    if sqrt(Inner_product(X_delta(:,k),M*X_delta(:,k))) < epsilon
        break;
    end
end
disp(‘迭代次数为：‘);
k-1
CGM = vpa(X(:,k),e);
    %%内积
    function IP = Inner_product(M1,M2)
        MAX = max(size(M1));
        sum = 0;
        for i = 1:MAX
            sum = sum+M1(i)*M2(i);
        end
        IP = sum;
    end
end

　　2.例子

迭代次数为：
ans =
     2
S =
  -2.12121212
 -0.454545455
   1.21212121
   2.87878788
ans =
   -2.1212
   -0.4545
    1.2121
    2.8788
>>