[Fortran] 纯文本查看 复制代码
program main
implicit none
real(8)::a3(18,9),A3inv(9,18)
integer::log
a3=1.
call BGINV(18,9,A3,A3inv, Log, 1.e-15)
end
!!-----------------------------------------
!! general inverse
!!-----------------------------------------
SUBROUTINE BGINV(M,N,A,AA,L,EPS)
double precision :: A(M,N),U(M,M),V(N,N),AA(N,M),EPS
double precision,ALLOCATABLE :: S(:), E(:), WORK(:)
KA=MAX(M,N)+1
ALLOCATE(S(KA))
ALLOCATE(E(KA))
ALLOCATE(WORK(KA))
CALL BMUAV(A,M,N,U,V,L,EPS,KA,S,E,WORK)
IF (L.EQ.0) THEN
K=1
10 IF (A(K,K).NE.0.0) THEN
K=K+1
IF (K.LE.MIN(M,N)) GOTO 10
END IF
K=K-1
IF (K.NE.0) THEN
DO 40 I=1,N
DO 40 J=1,M
AA(I,J)=0.0
DO 30 II=1,K
30 AA(I,J)=AA(I,J)+V(II,I)*U(J,II)/A(II,II)
40 CONTINUE
END IF
END IF
DEALLOCATE(S)
DEALLOCATE(E)
DEALLOCATE(WORK)
RETURN
END
!!-----------------------------------------
!! singular decomposite
!!-----------------------------------------
SUBROUTINE BMUAV(A,M,N,U,V,L,EPS,KA,S,E,WORK)
!略去
end