该楼层疑似违规已被系统折叠 隐藏此楼查看此楼
function [x]=gaijin(A,b,n) %用改进平方根法求解Ax=b
L=zeros(n,n); %L为n*n矩阵
D=diag(n,0); %D为n*n的主对角矩阵
S=L*D;
for i=1:n %L的主对角元素均为1
L(i,i)=1;
end
for i=1:n
for j=1:n %验证A是否为对称正定矩阵
if (eig(A)<=0)||(A(i,j)~=A(j,i)) %A的特征值小于0或A非对称时,输出wrong
disp('wrong');break;
end
end
end
D(1,1)=A(1,1); %将A分解使得A=LDLT
for i=2:n
for j=1:i-1
S(i,j)=A(i,j)-sum(S(i,1:j-1)*L(j,1:j-1)');
L(i,1:i-1)=S(i,1:i-1)/D(1:i-1,1:i-1);
end
D(i,i)=A(i,i)-sum(S(i,1:i-1)*L(i,1:i-1)');
end
y=zeros(n,1); % x,y为n*1阶矩阵
x=zeros(n,1);
for i=1:n
y(i)=(b(i)-sum(L(i,1:i-1)*D(1:i-1,1:i-1)*y(1:i-1)))/D(i,i); %通过 LDy=b解得y的值
end
for i=n:-1:1
x(i)=y(i)-sum(L(i+1:n,i)'*x(i+1:n)); %通过LTx=y解得x的值
end
>> A=[2 -1 1,-1 -2 3,1 3 1]
b=[4 5 6]
n=3
gaijin(A,b,n)
错误使用 eig
对于标准特征值问题 EIG(A),A 必须为方阵。
出错 gaijin (line 15)
if (eig(A)<=0)||(A(i,j)~=A(j,i))
%A的特征值小于0或A非对称时,输出wrong