Is there a way to implicitly define a function midprogram in fortran'

Question

I am trying to solve some eigenvalue problems. I am using gfortran. I have written a module that contains the bisection_method subroutine. In my program, I defined the secular(m,x) function that takes a 2D array m and a float x and outputs the characteristic polynomial of m at x; i.e. det (m - x Id). One of the arguments of the bisection_method subroutine is a function f, which was intended as a function of one real input and real output. I would want to input, however, the 'partially applied' function secular(m,_) to this subroutine. Is there a way to do this, without explicitly defining this function in a module?

I can't define this functions explicitly, for I am to do this procedure for several matrices m. ALso, I cannot modify the body of the bisection_method because I use it for function of one real argument as well. Is there a way around this in Fortran?

Do you mean like a [_closure_](https://stackoverflow.com/q/9692697/3157076)? — francescalus, Nov 19 '20 at 15:06
isn't it possible to pass a wrapper internal procedure of the form 'f(x)' (that references 'm' via host association) to 'bisection_method()' as an argument? — roygvib, Nov 19 '20 at 17:27
@roygvib I solved it defining a wrapper function with parameter `m`, and added a setter subroutine. — Victor Gustavo May, Nov 19 '20 at 17:42

score 2 · Answer 1 · answered Nov 22 '20 at 12:38

As @francescalus pointed out, you want to have a closure for your problem.

Closures are partially supported in Fortran by using internal procedures, because the internal procedure has access to all variables in the surrounding scope.¹

Assuming that you want to find the eigenvalues of M with your code it could be structured like this.²

module bisection_and_linalg
    use iso_fortran_env, only: real64
    integer, parameter :: wp = real64
    implicit none(type, external)

    abstract interface
        real(wp) pure function real_function(x)
            real(wp), intent(in) :: x
        end function
    end interface

contains

    !> Find the root of f in the interval I
    real(wp) pure function bisect(f, I)
        procedure(real_function) :: f
        real(wp) :: I(2)
        ...
    end function

    !> Evaluate the characteristic polynomial of m at x
    real(wp) pure function secular(M, x)
        real(wp), intent(in) :: M(:, :), x

        ...
    end function

    !> Get eigenvalues
    real(wp) pure function eigenvalues(M)
        real(wp), intent(in) :: M(:, :)

        ...

        ! n.b. here you can use the bisection method to
        !   find your eigenvalues.
        bisect(f, ...)

    contains

        real(wp) pure function f(x)
            ! n.b. here you have your closure.
            !   M is captured from surrounding scope.
            f = secular(M, x)
        end function
    end function

end module

¹ The only caveat is, that internal procedures only exist as long the surrounding scope exists. So it is unfortunately not possible to write a generic function that takes a function and returns a function pointer to a partially applied version of it by using internal procedures. But this is not a concern for your problem.

² Perhaps it would be better to return an actual characteristic polynomial. Then you could derive it, to use e.g. Newton-Raphson instead of bisection.

Is there a way to implicitly define a function midprogram in fortran'

1 Answers1