■LU Decompositionのプログラミング: 対称バンドマトリクスの場合■
次に対称バンドマトリクスに対し三角行列分解を応用する場合を紹介します。バンドマトリクスですから正方マトリクスを
下図の左に示すようにマトリクスを押し潰します。更に対称バンドマトリクスですから、左半分を捨てます。
=== | 対称であるから 左半分を削除 | ===> |
SUBROUTINE LLDECOMPVT ( A, NNODE, NBW, MXW, MXN ) IMPLICIT REAL*8(A-H,O-Z) DIMENSION A(MXN,MXW) C-------I = 1 DO J = 2, NBW A(1,J) = A(1,J)/A(1,1) END DO C------ I >= 2 DO I = 2 , NNODE SUM = 0.D0 DO M = 2, MIN0(I,NBW) II = I - M + 1 SUM = SUM + A(II,M)**2*A(II,1) END DO A(I,1) = A(I,1) - SUM DO J = 2, MIN0(NBW,NNODE-I+1) SUM = 0.D0 DO M = 2, MIN0( I,NBW, NBW-J+1 ) II = I - M + 1 L = J + M - 1 SUM = SUM + A(II,1)*A(II,M)*A(II,L) END DO A(I,J) = (A(I,J)-SUM)/A(I,1) END DO END DO RETURN ENDNBWは、ハーフバンド幅を示します。MXNとMXWは、[A]のメモリー宣言のための整数です。 また、連立方程式をとく手順は以下のようになります。
CALL LLDECOMPVT ( A, NNODE, NBW, MXW, MXN ) C DO I=2,NNODE DO J = 2, MIN0(I,NBW) B(I) = B(I) - A(I-J+1,J)*B(I-J+1) END DO END DO C DO I = 1 , NNODE B(I) = B(I) / A(I,1) END DO C DO I = NNODE-1 , 1 , -1 DO J = 2 , MIN0(NNODE-I+1,NBW) B(I) = B(I) - A(I,J)*B(I+J-1) END DO END DO
BACK | NEXT |
---|
Menu | LU Decompo | Stiff | 3D Solid | 3D Fluid | Eigen&Lanczos | Sound Eigen | Solid Eigen | Solid Axisym |