直接三角分解LU法

silence123

%直接三角分解法LU;
5 n" ^4 u& O8 y9 m3 A1 Q; tformat rational
disp('请输入矩阵的阶数');
N=input('N=');
- C0 M; c8 D7 a3 ?for i=1:N
3 O0 }! i4 @  X0 P' o    for j=1:N
" i9 P* `; a; A, @7 u# R" y4 h- W7 F        A(i,j)=input('A(i,j)=');
    end
end
1 E$ b; u* j) a, @% T7 ?% A/ D9 y; Y0 k( QL=eye(N);
$ }. D3 r3 M0 M5 W! D' q; o8 N' cfor i=1:N
    U(1,i)=A(1,i);
    L(i,1)=A(i,1)/U(1,1);
& M8 e8 }! Y! E3 m+ ^$ kend
for r=2:N
6 M% y* Y! t+ L& Z+ d; J    for i=r:N
: H3 ]  B( x6 Q4 n% \        U(r,i)=A(r,i);
4 o6 p  @0 Y0 s+ D6 i; 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);
; \0 T% [' ^6 J# q8 D8 n        end
, b: ]! W0 f4 m( r6 O  i6 \        L(i,r)=L(i,r)/U(r,r);
    end
end
; Z# }3 j, d% g5 zdisp('A=');
disp(A);
* G- F7 G0 D- U5 a% o1 Bdisp('L=');
5 B. U2 W. `( D% [5 M$ xdisp(L);
8 d) Z+ q( u6 R9 Fdisp('U=')
# x' C* b: Y' C$ [8 ?disp(U);

0.9清1.8清2.7清