cholesky楚列斯基分解代码

functionx=cholesky(A,b)%求正定对称矩阵A的楚列斯基分解%b为方程组的右端项%x为方程组的Ax=b的解[m,n]=size(A);ifm~=ndisp('ArgumentmatrixAmustbesquare!');return;elseifm~=length(b)disp('ThedimentionsofAandbdonotagree!');return;endL=zeros(n);x=zeros(n,1);y=zeros(n,1);%分别对L，x，y赋初始值fork=1:n;if(A(k,k)-L(k,1:k-1)*L(k,1:k-1)')=0disp('error');return;elseL(k,k)=sqrt(A(k,k)-L(k,1:k-1)*L(k,1:k-1)');%求下三角矩阵L对角线的元素L(k+1:n,k)=(A(k+1:n,k)-L(k+1:n,1:k-1)*L(k,1:k-1)')/L(k,k);%求L对角线下每列的元素endy(k)=(b(k)-L(k,1:k-1)*y(1:k-1))/L(k,k);%求下三角方程组Ly=b的解y；endx(n)=y(n)/L(n,n);fork=n:-1:1x(k)=(y(k)-L(k+1:n,k)'*x(k+1:n))/L(k,k);%求上三角方程组L'x=y的解xenda=10*ones(1,100);c=ones(1,99);d=ones(1,99);A=diag(a)+diag(c,1)+diag(d,-1);b=12*ones(100,1);b(1)=11;b(100)=11;x=cholesky(A,b)