extraction.m

来自「线性方程组的数值解法」· M 代码 · 共 41 行

function extraction(A,B)
%平方根法  Cholesky分解，参数说明(系数矩阵A,常数项B)
error(nargchk(2,2,nargin));%判断输入矩阵大小
d=max(size(A));%获取系数矩阵的阶
p=max(size(B));
if d==p
    dimention=d;
else
    fprintf('输入有误，不能构成线性方程组');
end

for j=1:d
    for k=1:j-1
        A(j,j)=A(j,j)-(A(j,k))^2;
    end
    A(j,j)=sqrt(A(j,j));
    for i=j+1:d
        for k=1:j-1        
            A(i,j)=A(i,j)-A(i,k)*A(j,k);
        end
        A(i,j)=A(i,j)/A(j,j);
    end
end
      
L=tril(A)

for i=1:d
    for k=1:i-1
        B(i)=B(i)-A(i,k)*B(k);
    end
    B(i)=B(i)/A(i,i);
end    
    
for i=d:-1:1   
    for k=i+1:d
        B(i)=B(i)-A(k,i)*B(k);
    end
    B(i)=B(i)/A(i,i);
end

X=B