Question

求fortran高手进来！！！就帮忙解决一下问题！！感激涕零，重赏

!**********************************************************************************************
subroutine nihe(x,y,pset,edp)
real*8 pset,edp
integer n,k
common/group3/ n,k
real*8 x(n),y(n),A(k+1,k+1),B(k+1,1)
!A为系数矩阵，(A,B)为增广矩阵
integer i,j,m
real*8 pp
pp=0

do m=0,k
do j=0,k
do i=1,n
pp=pp+x(i)**(m+j)
!print*,pp
end do
A(m+1,j+1)=pp
pp=0
! print*,A(m+1,j+1)
end do
end do

pp=0
do m=0,k
do i=1,n
pp=pp+y(i)*x(i)**m
end do
B(m+1,1)=pp
pp=0
!print*,b(m+1,1)
end do
!print*,B
call solveeq(A,B,pset,edp)

end subroutine nihe
!**********************************************************************************
subroutine solveeq(A,B,pset,edp)
real*8 pset,edp
integer n,k
common/group3/ n,k
real*8 A(k+1,k+1),B(k+1,1),Ainv(k+1,k+1),X(k+1,1)
!A为系数矩阵，(A,B)为增广矩阵,Ainv为矩阵A的逆矩阵。X为未知数矩阵
integer i,j
real*8 pp

pp=0edp=0
call qiuni(A,Ainv,k+1)
!print*,Ainv
do i=1,k+1
do j=1,k+1
pp=pp+Ainv(i,j)*B(j,1)
end do
X(i,1)=pp
pp=0
end do

do i=1,k+1
edp=edp+X(i,1)*pset**(i-1)
end do
!print*,edp
end subroutine solveeq
!**********************************************************************************
subroutine qiuni(aa,b,n)
!这是一个求逆矩阵的子程序。 aa为原矩阵，b为存放aa的逆矩阵，n为矩阵aa的维数
integer:: n,i,j,k
real*8 aa(n,n),b(n,n),a(n,n)
a=aa
do i=1,n
b(i,i)=1
end do
do i=1,n
b(i,:)=b(i,:)/a(i,i)
a(i,i:n)=a(i,i:n)/a(i,i)
do j=i+1,n
do k=1,n
b(j,k)=b(j,k)-b(i,k)*a(j,i)
end do
a(j,i:n)=a(j,i:n)-a(i,i:n)*a(j,i)
end do
end do
do i=n,1,-1
do j=i-1,1,-1
do k=1,n
b(j,k)=b(j,k)-b(i,k)*a(j,i)
end do
end do
end do
end
!**********************************************************************************

这里只列出了几个子程序，主程序调用的是call nihe(x,y,pset,edp)
x,y里面的元素数值都很小10的负3次方吧，但均不为零。我是这输出了一下系数矩阵A的逆Ainv,但是Ainv里面的元素大部分都是NaN,导致接下来的工作没法做，我想知道这是什么原因，该怎么改？？？？感激涕零！！！

Answer 1

1.建议你给出调用nihe的代码，及其参数的值。否则其他人无法帮助你。
2.某些算法可能具有局限性。不能适合于某些病态的算例。或者容易引起误差。
3.学会下断点，单步 debug。这是程序员的基本功