直接三角分解LU法

silence123

%直接三角分解法LU;
format rational
+ ~% s4 y0 l" y2 h9 vdisp('请输入矩阵的阶数');
N=input('N=');
for i=1:N
    for j=1:N
        A(i,j)=input('A(i,j)=');
    end
) r! ^+ M+ R$ ~7 v: k: Y1 b; Mend
: i; t1 b0 n) d- gL=eye(N);
  g. G- e' x- t& w; _# A2 v7 M9 y; Mfor i=1:N
, _0 y  ?, t( T) X: b" ^    U(1,i)=A(1,i);
/ s% w  u3 c" `4 t  A( _    L(i,1)=A(i,1)/U(1,1);
end
7 ^; r( N0 V) a9 V. d$ f) hfor r=2:N
    for i=r:N
) u* k) N: P8 ?$ O9 p. O6 F& z7 q        U(r,i)=A(r,i);
- f5 B5 A/ Q6 Z% T        for k=1:r-1
3 Y, m2 K' Z9 p/ m# d& l            U(r,i)=U(r,i)-L(r,k)*U(k,i);
9 B7 V# p! T" E% L7 l        end
4 G3 O& n2 H' b" u5 r6 l) }        L(i,r)=A(i,r);
        for k=1:r-1
            L(i,r)=L(i,r)-L(i,k)*U(k,r);
9 _5 e! t7 P& L; F% b* t6 B# q3 L# B        end
        L(i,r)=L(i,r)/U(r,r);
) O1 v1 n9 W; @2 \; N% z3 `* m* O- Q    end
end
disp('A=');
# x, x, y" ?% @- J2 G- Edisp(A);
2 M) S$ Y; V% {7 ]. i. tdisp('L=');
disp(L);
; A$ @, J1 {4 P2 B* ?1 T1 g. [; `disp('U=')
disp(U);

0.9清1.8清2.7清