调用二分法子程序解方程出错，求助

Karmi

[Fortran] 纯文本查看 复制代码

  Program test
  Implicit none
  Integer :: i
  Real :: Maact1,Maact2
  Real,dimension(3) ::fact
  fact(1)=26.
  fact(2)=30.
  fact(3)=32.
  Do i=1,3                                      !输出的结果Maact1和Maact2结果一致，都为1.985840
    Call equation(fact(i),Maact1)               !但是把9-12行注释掉的话，
    Write(*,*) Maact1                           !输出的Maact2结果为2.1938，是正确的
  End Do
  Call equation(31.57,Maact2)
  Write(*,*) Maact2
  End Program test

  Subroutine equation(f,Ma)              !已知函数值，根据普朗特-迈耶函数求马赫数,1<=Ma<=3
  Implicit none
  Real,intent(in) :: f                   !输入普朗特-迈耶函数值
  Real,intent(out) :: Ma                 !输出马赫数
  Real :: eps=0.001                      !计算精度
  Real :: Ma1=1.0                        !区间下限
  Real :: Ma2=3.0                        !区间上限
  Real :: f_Ma1,f_Ma2,f_Ma               !对应函数值
  Real :: gamma=1.4
  Real,External :: P_Ma                  !普朗特迈耶函数
  !二分法解方程
  Ma=(Ma1+Ma2)/2
  f_Ma=P_Ma(Ma)-f
  f_Ma1=P_Ma(Ma1)-f
  f_Ma2=P_Ma(Ma2)-f
  Do while((abs(f_ma)>eps).and.((Ma2-Ma1)>eps))
    If(f_Ma*f_Ma1<0.) then
      Ma2=Ma
      Ma=(Ma1+Ma2)/2
      f_Ma2=f_Ma
    Else
      Ma1=Ma
      Ma=(Ma1+Ma2)/2
      f_Ma1=f_Ma
    End If
    f_Ma=P_Ma(Ma)-f
  End Do
  
  !f_Ma1=P_Ma(Ma1)-f
  !f_Ma2=P_Ma(Ma2)-f
  !Do
  !  Ma=(Ma1+Ma2)/2
  !  f_Ma=P_Ma(Ma)-f
  !  If ((abs(f_Ma)<=eps).or.abs(Ma1-Ma2)<=eps) exit
  !  If (f_Ma1*f_Ma<0) then
  !    Ma2=Ma
  !    f_Ma2=f_Ma
  !  Else
  !    Ma1=Ma
  !    f_Ma1=f_Ma
  !  End if
  !End do
  End subroutine equation

  Real Function P_Ma(Ma)
  Implicit None
  Real,Intent(in) :: Ma
  Real :: gamma=1.4
  P_Ma=sqrt((gamma+1)/(gamma-1))*atand(sqrt((gamma-1)/(gamma+1)*(Ma**2-1)))-atand(sqrt(Ma**2-1))
  End FUnction P_Ma

当我调用子程序的时候，输入如果是一个具体的数字的话，就能得到正确结果，但如果输入一个赋值的变量名的话就结果错误。
二分法也换了种写法，但还是同样的问题，是我的调用有什么问题吗？

Karmi

li913 发表于 2020-5-19 12:10
在定义变量时进行赋值，变量具有save属性，会保存前一次的计算结果。所以，把22,23行改一下就是。 ...

原来是这样，非常感谢，那么我把Ma1，Ma2作为参数输入来处理也是可以的了。

li913

在定义变量时进行赋值，变量具有save属性，会保存前一次的计算结果。所以，把22,23行改一下就是。