我数学不太好 所以看不懂这个求逆的subroutine,求出的矩阵值里有NaN 有无限 不知道问题出在哪里
SUBROUTINE matinv(n,a,ainv)
! inverts an (n x n) matrix
REAL*8 :: a(n,n),b(n,2*n),ainv(n,n),temp
INTEGER :: i,j,j1,j2,k,kp1,km1,l
do i=1,n
do j=1,n
b(i,j)=a(i,j)
END do
END do
j1=n+1
j2=2*n
do i=1,n
do j=j1,j2
b(i,j)=0.
END do
END do
do i=1,n
j=i+n
b(i,j)=1.
END do
do k=1,n
kp1=k+1
if (k == n) go to 6
l=k
do i=kp1,n
if (abs(b(i,k)) > abs(b(l,k))) l=i
END do
if (l == k) go to 6
do j=k,j2
temp=b(k,j)
b(k,j)=b(l,j)
b(l,j)=temp
END do
6 do j=kp1,j2
b(k,j)=b(k,j)/b(k,k)
END do
if (k == 1) go to 9
km1=k-1
do i=1,km1
do j=kp1,j2
b(i,j)=b(i,j)-b(i,k)*b(k,j)
END do
END do
if (k == n) go to 11
9 do i=kp1,n
do j=kp1,j2
b(i,j)=b(i,j)-b(i,k)*b(k,j)
END do
END do
END do
11 do i=1,n
do j=1,n
k=j+n
ainv(i,j)=b(i,k)
END do
END do
return
end subroutine matinv
SUBROUTINE matinv(n,a,ainv)
! inverts an (n x n) matrix
REAL*8 :: a(n,n),b(n,2*n),ainv(n,n),temp
INTEGER :: i,j,j1,j2,k,kp1,km1,l
do i=1,n
do j=1,n
b(i,j)=a(i,j)
END do
END do
j1=n+1
j2=2*n
do i=1,n
do j=j1,j2
b(i,j)=0.
END do
END do
do i=1,n
j=i+n
b(i,j)=1.
END do
do k=1,n
kp1=k+1
if (k == n) go to 6
l=k
do i=kp1,n
if (abs(b(i,k)) > abs(b(l,k))) l=i
END do
if (l == k) go to 6
do j=k,j2
temp=b(k,j)
b(k,j)=b(l,j)
b(l,j)=temp
END do
6 do j=kp1,j2
b(k,j)=b(k,j)/b(k,k)
END do
if (k == 1) go to 9
km1=k-1
do i=1,km1
do j=kp1,j2
b(i,j)=b(i,j)-b(i,k)*b(k,j)
END do
END do
if (k == n) go to 11
9 do i=kp1,n
do j=kp1,j2
b(i,j)=b(i,j)-b(i,k)*b(k,j)
END do
END do
END do
11 do i=1,n
do j=1,n
k=j+n
ainv(i,j)=b(i,k)
END do
END do
return
end subroutine matinv