Fortran function return value polymorphism

Problem background

Suppose I need to develop a multi-precision math library and define a polymorphic function eye.
Is there any magic in Fortran to make x=eye(2,2) return the value of the eye function according to the type of x?

Change to subroutine?

subroutines have different advantages from functions. Functions can be written in formulas.

Overload equal sign(assignment(=))?

Fortran Class (*) in Function Result

module hello_bob
  interface assignment(=)
    module procedure int_equal_func_class
  end interface

contains

  subroutine int_equal_func_class(a,b)
    integer, intent(out) :: a(:)
    class(*), intent(in) :: b(:)

    select type (b)
      type is (integer)
        a = b
    end select
  end subroutine int_equal_func_class

   function func(a)
     class(*), intent(in) :: a(:)
     class(*), allocatable :: func(:)

     ! No intrinsic assignment supported, also see note about bounds
     allocate(func, source=a)
   end function func

end module hello_bob

program bob
  use hello_bob
  integer i(4)

  i=func([1,2,3,4])
  print*, i

end program bob

Is there any downside to customizing such overloads?
For example, print *, eye(2,2) will be invalid?

Or does Fortran have other magic?

(I’m just interested, and I don’t really need a solution. Thank you! :grinning:)

1 Like

There is a proposal open for this:

IMO, having this would allow programmers to create and use derived-types (multiprecision integers, higher precision reals, reals supporting auto-differentiation) that are fully compatible with all of the intrinsic functions. Currently, if you a have program that relies on explicit conversions like real(val,sp) (assume val is a double-precision), you cannot replace val with a custom numeric type.

The best workaround I can think of currently is to pass a dummy variable to the procedure.

function eye(n,i) result(res)
  integer, intent(in) :: n
  integer(kind=...) :: i   ! just here to get right output kind
  integer(kind=kind(i)) :: res(n,n)
end function
3 Likes

You can do this with deferred evaluation

module m
  implicit none
  
  integer, parameter :: dp = selected_real_kind(15,300)
  
  type :: Eye
    integer :: m
    integer :: n
  end type
  
  interface assignment(=)
    module procedure assign_integers_Eye
    module procedure assign_reals_Eye
  end interface
contains

subroutine assign_integers_Eye(output,input)
  integer, allocatable, intent(out) :: output(:,:)
  type(Eye),            intent(in)  :: input
  
  integer :: i
  
  allocate(output(input%m,input%n))
  output = 0
  do i=1,min(input%m,input%n)
    output(i,i) = 1
  enddo
end subroutine

subroutine assign_reals_Eye(output,input)
  real(dp), allocatable, intent(out) :: output(:,:)
  type(Eye),             intent(in)  :: input
  
  integer :: i
  
  allocate(output(input%m,input%n))
  output = 0
  do i=1,min(input%m,input%n)
    output(i,i) = 1
  enddo
end subroutine
end module

program p
  use m
  implicit none
  
  integer,  allocatable :: int_array(:,:)
  real(dp), allocatable :: real_array(:,:)
  
  integer :: i
  
  int_array = Eye(2,2)
  
  write(*,*)
  do i=1,2
    write(*,*) int_array(i,:)
  enddo
  
  real_array = Eye(2,2)
  write(*,*)
  do i=1,2
    write(*,*) real_array(i,:)
  enddo
end program

The general form for this kind of deferred evaluation is that a function foo(a,b,c) returns a type foo which simply contains a, b and c, and then the assignment bar = foo is overloaded to do what you’d expect the function foo(a,b,c) to do.

I wouldn’t recommend overloading integer foo(:) = class(*) bar(:), as this will conflict with all integer foo(:) assignments, including foo(:) = [1,2,3].

6 Likes

Your example is not runnable using the existing Fortran syntax, is it?

integer :: sp=4, dp = 8
print *, kind(sp), kind(dp)
   !! kind(4)=4, kind(8)=4

Apologies. Click below for a full example of the workaround I had in mind:

Click here
module eye_mod

  implicit none

  integer, parameter :: sp = kind(1.0e0)
  integer, parameter :: dp = kind(1.0d0)

  interface eye
    module procedure eye_default
    module procedure eye_sp
    module procedure eye_dp
  end interface

contains

  function eye_default(n) result(res)
    integer, intent(in) :: n
    real(sp) :: res(n,n)
    integer :: k
    res = 0
    do k = 1, n
      res(k,k) = 1
    end do
  end function

  function eye_sp(n,i) result(res)
    integer, intent(in) :: n
    real(sp), intent(in) :: i
    real(kind(i)) :: res(n,n)
    integer :: k
    res = 0
    do k = 1, n
      res(k,k) = 1
    end do
  end function

  function eye_dp(n,i) result(res)
    integer, intent(in) :: n
    real(dp), intent(in) :: i
    real(kind(i)) :: res(n,n)
    integer :: k
    res = 0
    do k = 1, n
      res(k,k) = 1
    end do
  end function

end module

The main program would look something like:

program test_eye

  use eye_mod, only: eye, wp => dp

  implicit none

  ! Dummy variable used to select precision
  real(wp), parameter :: dum_wp = 1.0_wp

  integer, parameter :: n = 3
  real(wp) :: I(n,n)
  integer :: k

  I = eye(n,dum_wp)

  do k = 1, n
    print *, I(k,:)
  end do

end program
1 Like

Thank you for providing a magic, eye has changed from a function to a type structure, but it achieves its function wonderfully!

Haha, If I have functions zeros(dim1), zeros(dim1,dim2), which have different dimensions, can I still use your method to construct such a polymorphic return value function?

Thank you very much! Clean and tidy!

Haha, yeah, “magic” is the right word. See the C++ Armadillo library for this kind of nonsense writ large (complete with incomprehensible template error messages yelling about “meat glue” for some reason).

You can make zeros work; you can have two types, zeros1d which is constructed by zeros(dim1) and zeros2d which is constructed by zeros(dim1,dim2). Then you just need to overload integer(:) = zeros1d and integer(:,:) = zeros2d.

You can get pretty far with this approach, although I can’t promise that anyone reading your code in the future (probably including future you) will thank you for it!

Can I achieve the effect of the official function, similar to real(x,8) instead of real(x, 8.0_wp)? The former may be more intuitive and concise. Is there such a method?

There is no method currently.

The workaround I presented relies on using an actual argument of given kind in place of its integer kind specifier.

:warning: Please avoid using integer literals like 4 and 8 as kind specifier at all costs. :warning: When you write real(x,8) it is not obvious that 8 is the kind specifier . One improvement is to be explicit with real(x,kind=8). This can still result in non-portability, as the integer kind specifiers are not universal across compilers!

When you write my_func(x,8.0_wp) it is again misleading, as the numerical value has no purpose other than to communicate the kind (e.g. you could pick -555.0_wp with the same end effect). That is why I defined a dummy parameter dum_wp in the example above.

1 Like