Is creating nested subroutines/functions considered good practice in Fortran?

Help
17

The concept is known as a callback. In other words you have a subroutine that takes another subroutine or function as argument.

The classic example would be a numerical solver and a parametrized function of some sort. We can take root-solving as an example:

  ! find the root f(x) = 0, in the interval [ax,bx]
  x = zeroin(ax,bx,f,tol)

The function f passed as a parameter to zeroin must match a particular interface:

interface
  real(dp) function f(x)
    import dp
    real(dp), intent(in) :: x
  end function
end interface

To make this more concrete, say you had to find the friction factor f prescribed by the Colebrook-White equation:

\frac{1}{\sqrt{f}} = -2 \log \left( \frac{\epsilon}{3.7 D_\text{h}} + \frac{2.51}{\mathit{Re} \sqrt f} \right)

Finding f can be recast as finding the root of the equation:

G(f) = \frac{1}{\sqrt{f}} +2 \log \left( \frac{\epsilon}{3.7 D_\text{h}} + \frac{2.51}{\mathit{Re} \sqrt f} \right) = 0

In Fortran you might do this as follows:

real(dp) :: f, Re, Dh, eps
real(dp), parameter :: tol = 1.0e-9_dp

Re = 10000
Dh = 0.2_dp    ! in meters
eps = 0.002_dp ! 2 mm

f = zeroin(0.001_dp, 0.01_dp,friction_equation,tol)
print *, f

contains

  real(dp) function friction_equation(f) result(y)
    import, only: Re, Dh, eps  ! if supported by your Fortran compiler
    real(dp), intent(in) :: f
    y = 1.0_dp/sqrt(f) + 2*log(eps/(3.7_dp*Dh) + 2.51_dp/(Re*sqrt(f))) 
  end function

This could be made part of general subroutine taking Re, D_h and \epsilon as inputs, making sure they are valid and finally returning the friction factor f as output.

3 Likes