數值分析實驗之平方根法解線性方程組（MATLAB代碼）

平方根法解如下方程組

  

 

• 法一： 

 

clear
clc
A=input('輸入對稱正定矩陣A=')
B=input('輸入自由項B=')
n=length(A(:,1));
for k=1:n
   if (det(A(1:k,1:k))<=0)
       input('矩陣不是正定矩陣,請重新運行程序')
   end
end
%分解A=L*L'
for i=1:n
   t=0;
   for s=1:i-1
       t=t+L(i,s)^2;
   end
   L(i,i)=sqrt(A(i,i)-t);
   for k=i+1:n
       u=0;
       for s=1:i-1
           u=u+L(i,s)*L(k,s);
       end
       L(k,i)=(A(k,i)-u)/L(i,i);
   end
end
%分解AX=B為Ly=B  Lx=y

%求y
for i=1:n
   r=0;
   for k=1:i-1
       r=r+L(i,k)*y(k);
   end
   y(i)=(B(i)-r)/L(i,i);
end
%求x
for i=n:-1:1
   q=0;
   for k=i+1:n
       q=q+L(k,i)*x(k);
   end
   x(i)=(y(i)-q)/L(i,i);
end
x' %若要輸出行向量，此處改為x即可

　　　　運行結果：

       

 

• 法二：

function pngfanggen(n)
A=[4 2 -4 0 2 4 0 0;
   2 2 -1 -2 1 3 2 0;
  -4 -1 14 1 -8 -3 5 6;
   0 -2 1 6 -1 -4 -3 3;
   2 1 -8 -1 22 4 -10 -3;
   4 3 -3 -4 4 11 1 -4;
   0 2 5 -3 -10 1 14 2;
   0 0 6 3 -3 -4 2 19];
B = [0 -6 20 23 9 -22 -15 45];
L=zeros(8,1);
for j=1:n
    L(j,j)=sqrt(A(j,j)-sum(L(j,:)*L(j,:)'));
    for i=j+1:n
        L(i,j)=(A(i,j)-sum(L(i,:)*L(j,:)'))/L(j,j);
    end
end


Y=zeros(1,n);
for i=1:n
    Y(i)=(B(i)-sum(L(i,:)*Y'))/L(i,i);
end

X(n)=Y(n)/L(n, n);
for i=n-1:-1:1
    X(i)=(Y(i)-sum(L(:,i)'*X'))/L(i,i);
end

disp('平方根法求得L=');
disp(L);
disp('該線性方程組的解為X=');disp(X');
end

　　運行結果：