Fortran 程序无法正确更新变量的值

Question

我正试图在Fortran中实现一个二分法子程序来解决一个计算科学程序，Fortran正在做一些奇怪的事情。因此，程序的目标是找到某个参数e0的超越方程的解，该参数在for循环的每一步更新并传递给子程序。

问题是:

• e0未按要求更新。e0的取值范围为0.1799…0.8999……在一次迭代中，理论上程序需要9次迭代才能完成。为什么会发生这种情况?
• 打印语句没有打印出来(见下面的输出)。我们期待一个"out"。打印，一些"打印(当子程序Bisection调用子程序f1时打印)，然后输出2;打印(用新的e0值)，一些"in"打印等。但是我们看到了第一个"1"。打印，一些"只打印"输出"打印。这是否意味着子程序f1只在"out "之间被调用?还有"out 2"?(它应该在每次"输出"之间调用;印刷)

我使用Fortran77进行数值求解已经有几年了，从来没有遇到过这样的事情，但是自从我没有编写任何程序以来已经将近6个月了，所以也许我错过了一个重要的事情。

代码:

program roots
implicit none
double precision A,B,eps,e0,e1
integer i,niter
external f1
A = 1d-1
B = 9d-1
eps = 10d-6
open(9,file='dades.dat',status='old')

do i=1,9
e0 = i*(B-A)/10 + A
c      Here is the first print. 'Out' meaning outside the subroutine
print *, 'out', i, e0
call Bisection(A,B,eps,f1,niter,e1,e0)
write (9,'(2(f10.5))') e0,e1
end do

close(9)
end program roots

subroutine Bisection(A,B,eps,f,niter,xroot,e0)
implicit none
double precision A,B,eps,xroot,fuc,fua,e0
integer niter,i
niter = nint(log((B-A)/eps)/log(2.))+1

do i=1,niter
xroot = (A+B)/2
c      Here the subroutine which uses e0 is called twice
call f(xroot,fuc,e0)
call f(A,fua,e0)
if (fuc .eq. 0) return
if (fuc*fua .lt. 0.) then
B = xroot
else
A = xroot
end if
end do
return
end subroutine Bisection

subroutine f1(x,f,e0)
implicit none
double precision x,f,e0
c      Here is the second print. 'In' meaning inside the subroutine
print *, 'in', e0
f = e0+1
end subroutine f1

输出:

out           1  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
in  0.17999999999999999
out           2  0.89999511718750003
out           3  0.89999572753906254
out           4  0.89999633789062505
out           5  0.89999694824218746
out           6  0.89999755859374997
out           7  0.89999816894531248
out           8  0.89999877929687500
out           9  0.89999938964843751

Answer 1

经过对代码的深入研究，我终于发现了这个问题。问题是变量A和B都是输入和输出，因为它们是在子例程中修改的。因此，当第一次调用子例程时，它会改变A和B的值(与代码开头设置的值不同)。其余的都是由此引发的问题。
解决方案是在每次迭代开始时重置值，如下所示:

program roots
implicit none
double precision A,B,eps,e0,e1
integer i,niter
external f1
eps = 10d-6
open(9,file='dades.dat',status='old')

do i=1,9
A = 1d-1
B = 9d-1
e0 = i*(B-A)/10 + A
c      Here is the first print. 'Out' meaning outside the subroutine
print *, 'out', i, e0
call Bisection(A,B,eps,f1,niter,e1,e0)
write (9,'(2(f10.5))') e0,e1
end do

close(9)
end program roots

subroutine Bisection(A,B,eps,f,niter,xroot,e0)
implicit none
double precision A,B,eps,xroot,fuc,fua,e0
integer niter,i
niter = nint(log((B-A)/eps)/log(2.))+1

do i=1,niter
xroot = (A+B)/2
c      Here the subroutine which uses e0 is called twice
call f(xroot,fuc,e0)
call f(A,fua,e0)
if (fuc .eq. 0) return
if (fuc*fua .lt. 0.) then
B = xroot
else
A = xroot
end if
end do
return
end subroutine Bisection

subroutine f1(x,f,e0)
implicit none
double precision x,f,e0
c      Here is the second print. 'In' meaning inside the subroutine
print *, 'in', e0
f = e0+1
end subroutine f1