直接三角分解LU法

silence123

%直接三角分解法LU;
! X0 o( o8 C' s5 k, qformat rational
8 @, {7 p8 p2 ?- |disp('请输入矩阵的阶数');
* x( R0 t) E4 ?7 BN=input('N=');
for i=1:N
    for j=1:N
        A(i,j)=input('A(i,j)=');
    end
end
L=eye(N);
/ \' Y& Q$ r5 j+ A) H. W, J% c! zfor i=1:N
9 D' p5 f$ |! b! w1 a. ?/ X% L    U(1,i)=A(1,i);
) v6 |: F) b9 O5 u& V2 ]    L(i,1)=A(i,1)/U(1,1);
6 \& u0 H+ f' s# Qend
0 ^- P. `; I- D* R( z9 D  R) i; Afor r=2:N
    for i=r:N
        U(r,i)=A(r,i);
        for k=1:r-1
' ?) J4 Z( D$ n% u            U(r,i)=U(r,i)-L(r,k)*U(k,i);
        end
        L(i,r)=A(i,r);
        for k=1:r-1
- a9 X4 o0 ?/ _: A! s            L(i,r)=L(i,r)-L(i,k)*U(k,r);
- j! _. }# [8 p1 a! o+ i% \  M        end
        L(i,r)=L(i,r)/U(r,r);
+ w' g5 K# X3 I1 F5 e2 R4 {    end
0 j7 y0 d4 a' L- rend
disp('A=');
disp(A);
7 K3 |- ]+ Y5 |+ P0 ^4 ?2 Sdisp('L=');
: c' }) n" D6 N+ k/ D; f; ydisp(L);
  s& }  b2 N0 m2 vdisp('U=')
8 E, Q+ H) s& v" I$ z% r' S8 Xdisp(U);

0.9清1.8清2.7清