直接三角分解LU法

silence123

%直接三角分解法LU;
/ d4 P3 z: d0 u, x4 A( L4 @format rational
5 Z: r% e$ r+ A  N, t5 cdisp('请输入矩阵的阶数');
N=input('N=');
& @4 \* \0 y2 x7 _for i=1:N
# r/ J/ c0 X( w6 {8 v5 E: Q    for j=1:N
8 }) [9 Y) s& D% ^% ?        A(i,j)=input('A(i,j)=');
5 ?+ c* V/ t& _  B    end
; T" w* z5 [4 }2 _6 Y/ Iend
L=eye(N);
for i=1:N
, I" F: G3 s( d* @+ R    U(1,i)=A(1,i);
1 s# z1 j- d0 j0 d! G+ ]+ U6 o; m+ F    L(i,1)=A(i,1)/U(1,1);
5 g2 Z/ J( d. U! i9 ]1 gend
) x# F4 z/ y6 ^# afor r=2:N
- e8 ^! K( \; y+ h" A* o    for i=r:N
  J5 K; R4 t9 k) D) |3 {! o        U(r,i)=A(r,i);
        for k=1:r-1
            U(r,i)=U(r,i)-L(r,k)*U(k,i);
        end
        L(i,r)=A(i,r);
        for k=1:r-1
            L(i,r)=L(i,r)-L(i,k)*U(k,r);
+ Q  v6 ?! ?/ `, }$ `9 N/ b# G; F0 v8 E        end
. C  W0 |$ R0 i        L(i,r)=L(i,r)/U(r,r);
4 ]* t' F3 ?+ m    end
end
, j: k, Y1 j& a& H0 {5 R1 q7 |& {+ _disp('A=');
0 H0 I7 |  ]& R; c5 `! Xdisp(A);
( v) M+ C2 s6 r+ Bdisp('L=');
( s8 `! {3 l. j% Jdisp(L);
, C- u: k/ j8 S, ~$ l2 v1 e4 z9 ddisp('U=')
disp(U);

0.9清1.8清2.7清