Best way to create a routine that can handle different precisions of input/output

Say I’m writing a simple helper routine within a module, for example

subroutine add_reals(result, x, y)
    real, intent(out) :: result
    real, intent(in ) :: x,y
    result = x + y
end subroutine add_reals

However, within the code that will use it, there might be a case where I want:

real(kind = 8), intent(out) :: result
real(kind = 8), intent(in ) :: x,y

or

real, intent(out) :: result
real, intent(in ) :: x,y
! i.e. I want to rely on the default real kind

What is the best way to set this up? Is the only way to define a generic interface like

interface add_reals
   module procedure add_reals_default_kind
   module procedure add_reals_r8
end interface

and then effectively repeat the routine?

I’m hoping to avoid repeating the contents, for maintenance reasons, of the actual routine I’m writing (which is notably more complex), as they will be identical outside of the specification section. However, if its necessary, its necessary! Just want to check with others on here!

Thanks!

It seems to me, that due to strong type checking, the only way I could do this would be via some clever use of includes?

include (either the intrinsic include of Fortran or the Fortran preprocessor include) is currently the most concise and energy-efficient way to implement generic Fortran interfaces. I recommend putting your generic interfaces in a module and the implementations in a submodule of that module, which includes a third file containing the generic type/kind-agnostic algorithm. This hierarchy of modules/submodules has the extra benefit of 1) avoiding long compilation cascades if you change the algorithm later and 2) hiding the implementation if needed and exposing only the interface to the end user.

1 Like

Use of INCLUDE is a useful way of generating generics. For completeness there is an ongoing discussion of templating proposals for F2023 that should help address this issue; and preprocessors like prep(1) allow for a basic templating capability (also see m4(1), …) and also that you can cast or promote as GitHub - urbanjost/M_anything: Use polymorphism to allow promoting, casting and molding intrinsic types demostrates.

Casting or promoting can simplify routines with a lot of options that you want to take any type but should be used sparingly, and probably only for procedures that are used just a few times. See the discussions about what troubles the old C language caused by promoting essentially all float values
for some basics on why that is not without peril before making heavy use of it (ie something like “anyscalar_to_double” in M_anything).

INCLUDE and generics should be considered first, but there are alternatives (with caveats).

Admittedly I am biased, but I use prep and $BLOCK/$SET/$POST a lot for complex generics, but I basically only post the resulting standard expanded Fortran code, but virturally all my modules are actually built with prep(1).

So the alternative approach illustrated here shows how the example “sum two numbers” can be implemented as “sum up to ten intrinsic scalars” very generically. Of course this simple example
is not too exciting, but it still illustrates something that would take a huge number of procedures if done as a generic; where the number of arguments can vary and the arguments can be of many different types.

module m_anything
use, intrinsic :: iso_fortran_env, only : int8, int16, int32, int64       
use, intrinsic :: iso_fortran_env, only : real32, real64, real128         
implicit none
private
integer,parameter :: dp=kind(0.0d0)
public anyscalar_to_double   
public generic_sum
contains

impure elemental function anyscalar_to_double(valuein) result(d_out)
class(*),intent(in) :: valuein
real(kind=dp)       :: d_out
   select type(valuein)
   type is (integer(kind=int8));   d_out=dble(valuein)
   type is (integer(kind=int16));  d_out=dble(valuein)
   type is (integer(kind=int32));  d_out=dble(valuein)
   type is (integer(kind=int64));  d_out=dble(valuein)
   type is (real(kind=real32));    d_out=dble(valuein)
   type is (real(kind=real64));    d_out=dble(valuein)
   type is (real(kind=real128));   d_out=dble(valuein)
   type is (complex);              d_out=abs(valuein)
   type is (logical);              d_out=merge(0.0d0,1.0d0,valuein)
   type is (character(len=*));     read(valuein,*) d_out
   class default; stop '*anyscalar_to_double* <ERROR> unknown type'
   end select
end function anyscalar_to_double

function generic_sum(a,b,c,d,e,f,g,h,i,j) result (dvalue)
class(*),intent(in),optional :: a,b,c,d,e,f,g,h,i,j
doubleprecision              :: dvalue
   dvalue=0.0d0
   if(present(a))dvalue=dvalue+anyscalar_to_double(a)
   if(present(b))dvalue=dvalue+anyscalar_to_double(b)
   if(present(c))dvalue=dvalue+anyscalar_to_double(c)
   if(present(d))dvalue=dvalue+anyscalar_to_double(d)
   if(present(e))dvalue=dvalue+anyscalar_to_double(e)
   if(present(f))dvalue=dvalue+anyscalar_to_double(f)
   if(present(g))dvalue=dvalue+anyscalar_to_double(g)
   if(present(h))dvalue=dvalue+anyscalar_to_double(h)
   if(present(i))dvalue=dvalue+anyscalar_to_double(i)
   if(present(j))dvalue=dvalue+anyscalar_to_double(j)
end function generic_sum

end module m_anything

program demo_anyscalar_to_double
use m_anything, only : gsum=>generic_sum
use, intrinsic :: iso_fortran_env, only : int8, int16, int32, int64
use, intrinsic :: iso_fortran_env, only : real32, real64, real128
implicit none
   ! call same function with different type and number of  scalar input types
   write(*,*)gsum(2,3.0,4.0d0,(3.0,4.0))
   write(*,*)gsum(1,2,3,4,5,6,7,8,9)
   write(*,*)gsum('123',456,.true.,.false.,1>=2,10==5*2)
   ! something more interesting
   if (gsum(1>=2,10==5*2,10>3,3+3==6) < 2) then
      write(*,*)'close enough, less than 2 conditions not met'
   else
      write(*,*)'two or more conditions not met'
   endif
end program demo_anyscalar_to_double
1 Like

Yes, it’s not exciting; nonetheless it’s a good illustration of something that you feel like doing because it can be done but which you should seriously resist doing!

1 Like

I agree. I recommend trying out the standard INCLUDE first and see if that can work for you.

1 Like

As I mentioned with just a few parameters I would start with INCLUDE. I have noted in the past that CLASS(*) has pitfalls, but REAL_VALUE1 == REAL_VALUE2 has more and it is used everywhere. The example is not too exciting because it just takes the simple OP-supplied subroutine a little farther. The technique itself is powerful, especially when the alternative would require hundreds of procedures in a generic. CLASS(*) is here to stay; but as I mentioned promotion and casting have their pitfalls and this is a relatively new one; but Fortran is full of them, as in I=1/3, X==Y, A=1234567890123456789.12345, DO i=1,huge(0); and on and on.

Since it is now more than ever easy to do – is there a good article online anyone knows of that discusses casting and promotion in particular? There is some mention of it in literature discussing the perils of floating point. A good description of why and when to use CLASS(*) and what it was most likely designed to do versus what it can do would be a great Fortran Wiki article, I think.

All that being said there are times where this usage pattern is great. Even this trivial routine creates a function that can convert strings to numeric values, lets you operate with logical values numerically, could easily be extended to do runtime checks for overflow, allows for conditionless expressions like A=GSUM(B<0)*B to avoid summing negative numbers, and so on. There will be more and more uses of CLASS(*), a good guide would be useful.

1 Like

Great thanks everyone - yeah, I’ve seen some pretty concise examples of setting something like this up with includes. I’ll try that route first. I thought of adding many optional input variables, but it can get ugly pretty fast.

It can be the briefest of Wiki articles:

Don’t use class(*)

Nothing further is needed.

:+1:

Re: “it can get ugly pretty fast,” agree.

Depending on the complexity in your actual routine(s), there may be some “tricks” using ASSOCIATE / BLOCK constructs that might help. If you notice things start to get too “ugly”, if you can post here some mockup of your routine(s) [presuming you may not be able to share actual], readers can try to give you pointers to mitigate the “ugliness”.

Or otherwise, you can turn toward the preprocessor route.

1 Like

I am not seeing this whole “simple and elegant” argument for using include “generic_implementation” to implement generic functions on many data types. Sure, maybe you can actually type this out pretty quick, because so little is typed… But maintenance of this type of structure would be terrible. You have a module that defines a bunch of interfaces, telling you what the input/outputs are and their types, then a separate file submodule which basically does nothing but have repeated 3 line

module procedure func_type
include "generic implementation"
end module procedure func_type

chunks like that many times over. The submodule file tells you basically nothing, other than the name of the implementation you now have to go open to see what the function is even doing. Except, problem here is you can’t see what any of the variables that have been genericized are… I don’t really see how an approach like this could scale beyond the most elementary math operations. Any logic beyond 100 lines would be dreadful.