直接三角分解LU法

silence123

%直接三角分解法LU;
9 n, d5 W" }) r% }" Z) c: w: Jformat rational
7 X7 l* ^! [: [, K( Ldisp('请输入矩阵的阶数');
: Z" [: K. G7 K( a3 Q  U# \1 N4 H  oN=input('N=');
for i=1:N
& h9 f' _0 [4 s2 v3 t    for j=1:N
        A(i,j)=input('A(i,j)=');
# `3 d1 c! w& x9 L: I+ d2 I    end
end
5 x" M. Y! W" r  @L=eye(N);
/ x0 `' P  p* j( K# L5 afor i=1:N
    U(1,i)=A(1,i);
! R& G0 m8 d# k# P    L(i,1)=A(i,1)/U(1,1);
* |7 e0 Q! Y2 D* O% ]+ I4 Rend
for r=2:N
/ u+ F! H8 {, I) w3 P2 f7 O    for i=r:N
$ \3 C/ I  f/ D  b7 Z        U(r,i)=A(r,i);
        for k=1:r-1
. @+ V4 C& ?) v4 o            U(r,i)=U(r,i)-L(r,k)*U(k,i);
        end
. t/ b: w& h+ D( t. [        L(i,r)=A(i,r);
# U. g( D1 B- I1 o4 }! y        for k=1:r-1
3 Q0 w$ w! E, B7 p            L(i,r)=L(i,r)-L(i,k)*U(k,r);
        end
        L(i,r)=L(i,r)/U(r,r);
    end
end
! \" S, T" h5 V2 idisp('A=');
. ~) n- U. d1 ]  X8 P5 O0 Cdisp(A);
disp('L=');
disp(L);
" Z/ p' Q+ \$ v! T5 B2 J& Ddisp('U=')
disp(U);

0.9清1.8清2.7清