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 >>
原文地址:https://www.cnblogs.com/guliangt/p/12119299.html
时间: 2024-10-11 08:55:19