Interface mismatch - higher order functions

Question

I'm trying to 'reproduce' higher order functions in fortran.

module rk4

contains



 pure function f(t,x) result (fx)
    real, dimension(1), intent(in) :: x
    real, intent(in) :: t
    real, dimension(1) ::  fx

    fx = x
  end function f

  function step(x,f,dt) result(xn)
    real, intent(in) :: dt
    real, intent(in),  dimension(:) :: x
    real, dimension(:), allocatable :: k1,k2,k3,k4,xn
    real, external :: f
    integer :: N

    N = size(x)

    allocate(k1(N))
    allocate(k2(N))
    allocate(k3(N))
    allocate(k4(N))

    k1 = f(t,x)
    k2 = f(t+0.5*dt,x+0.5*k1*dt)
    k3 = f(t+0.5*dt,x+0.5*k2*dt)
    k4 = f(t+dt,x+dt*k3)


    allocate(xn(N))
    xn = x + (dt/6.)*(k1 + 2*k2 + 2*k3 + k4)


    deallocate(k1)
    deallocate(k2)
    deallocate(k3)
    deallocate(k4)

  end function step

end module rk4

When the module is called in the following way

 real, dimension(1) :: x0 = 2

 x0 = step(x0,f,0.01)

I get the following error

$gfortran -c test_rk4.f95 
test_rk4.f95:7.15:

  x0 = step(x0,f,0.01)
               1
Error: Interface mismatch in dummy procedure 'f' at (1): Type/rank mismatch in function result

What might be causing this?

Answer 1

The error message is complaining about the function f being incompatible with the dummy argument f.

You declare it as

   real, external :: f

which means it should return a scalar whereas in reality the function f returns an array.

The same names do not really help the understanding here. I changed the name of the dummy argument in the following code to g.

The easiest way to solve this is

 pure function f(t,x) result (fx)
    real, dimension(1), intent(in) :: x
    real, intent(in) :: t
    real, dimension(1) ::  fx

    fx = x
  end function f

  function step(x,g,dt) result(xn)
    real, intent(in) :: dt
    real, intent(in),  dimension(:) :: x
    real, dimension(:), allocatable :: xn
    procedure(f) :: g

    !here call g, not f!!!

The procedure statement comes from Fortran 2003 and causes the dummy argument procedure g to have the same interface as the procedure f.

Otherwise you can use an interface block:

  function step(x,g,dt) result(xn)
    real, intent(in) :: dt
    real, intent(in),  dimension(:) :: x
    real, dimension(:), allocatable :: xn

    interface
      pure function g(t,x) result (fx)
        real, dimension(1), intent(in) :: x
        real, intent(in) :: t
        real, dimension(1) ::  fx
      end function g
    end interface

The external statement should be used in modern code only in some exceptional cases.