与 quadpack 的多维集成导致对非递归过程的递归调用 'dqage'

Question

我正在尝试使用 quadpack double library. Since this website does not have a Latex engine integrated I attach a picture 和数学计算二维积分。

我还在下面附上我在 test_module 模块中实现此操作。

module test_module
  use f90_kind ! defines dp
  use physconst ! defines pi
  use global_variables ! defines the variables for quadpack
  use csv_file ! allows me to write to a csv file

  implicit none

  real(dp) :: a, b, c ! parameters I want to set

contains

  real(dp) function multi_val_func(theta)
    real(dp), intent(in) :: theta

    multi_val_func = exp(theta**2 + c*theta)

  end function multi_val_func

  real(dp) function theta_integral(phi)
    real(dp), intent(in) :: phi

    c = phi*(a+b)

    ! Quadpack variables
    epsabs = 1.0E-14_dp
    epsrel = 0.0E0_dp
    key    = 6
    ilow   = 0.0E0_dp
    ihigh  = 2.0E0_dp*pi

    call dqage(multi_val_func, ilow, ihigh, epsabs, epsrel, key, limit, &
               res, abserr, neval, ier, alist, blist, rlist, elist, iord, last)

    theta_integral = res

  end function theta_integral

  subroutine phi_integral(final_result)
    real(dp),intent(out) :: final_result

    a = 1.0E0_dp
    b = 1.0E0_dp

    ! Quadpack variables
    epsabs = 1.0E-14_dp
    epsrel = 0.0E0_dp
    key    = 6
    ilow   = 0.0E0_dp
    ihigh  = pi

    call dqage(theta_integral, ilow, ihigh, epsabs, epsrel, key, limit, &
               res, abserr, neval, ier, alist, blist, rlist, elist, iord, last)

    final_result = res

    print *, final_result

  end subroutine phi_integral

end module test_module

在主文件中，我定义了 result 变量，然后是 call phi_integral(result)。我使用以下 makefile

编译此代码

FFLAGS = -O0 -fcheck=all -ffree-line-length-none

debug:
    gfortran -c $(FFLAGS) $(srcdir)/f90_kind.f90
    gfortran -c $(FFLAGS) $(srcdir)/physconst.f90
    gfortran -c $(FFLAGS) $(srcdir)/global_variables.f90
    gfortran -c $(FFLAGS) $(extlib_qp)/quadpack_double.f90
    gfortran -c $(FFLAGS) $(extlib_csv)/csv_file.f90
    gfortran -c $(FFLAGS) $(srcdir)/multi_var_integration.f90
    gfortran -c $(FFLAGS) $(srcdir)/main.f90
    gfortran *.o -o debug -lblas -llapack
    rm -f *.o *.mod

代码可以编译，但是当我运行它时，我收到以下错误：

Fortran runtime error: Recursive call to nonrecursive procedure 'dqage'

Error termination. Backtrace:
#0  0x103a28d3d
#1  0x103a299f5
#2  0x103a29fd6
#3  0x1039f94b6
#4  0x1039dee59
#5  0x1039fad68
#6  0x1039f9948
#7  0x1039dec7d
#8  0x1039def61
#9  0x1039defa3

非常感谢任何处理此问题的建议。

Answer 1

并非 Fortran 中的每个子例程或函数都是可重入的，或者像 Fortran 所说的那样是“递归的”。在 Fortran 2018 之前，可以调用的子例程，如果已经在其自身内部（即使是间接调用），则必须用 recursive 属性标记。

即使这样标记子程序，如果它包含静态（保存的）数据或使用全局（公共块、模块）变量，它也可能无法按预期工作。

我无法立即看到 http://www.netlib.no/netlib/quadpack/dqage.f 中使用的任何此类静态数据，因此它可能在添加 recursive 后实际工作。

但是，正如 francescalus 指出的那样，您最好使用专用的 2d 正交库。

Answer 2

theta积分可以通过误差函数erf完成。参见 wolfram alpha。

其结果是一个变量 phi 的函数，可以通过 quadpack 对其进行评估。因此，你根本不需要多维积分例程。

与 quadpack 的多维集成导致对非递归过程的递归调用 'dqage'

Multi dimensional integration with quadpack results in recursive call to nonrecursive procedure 'dqage'

fortran