Learn from Julia

[jacobwilliams]

Just wanted to put this comment here: https://fortran-lang.discourse.group/t/is-there-a-good-fortran-package-for-sde-dde-and-ode/1756/49

You could relatively easily rewrite the basic solvers in Fortran. That said, Fortran isn’t powerful enough to easily deal with things like sparsity of different structures, different input types (eg DoubleFloat, Float16), symbolic simplification, and all the fancy stuff like that. To rival the Julia ecosystem, you would have to impliment a good type system with generics and dispatch, a symbolic algebra library, and basically the rest of a more versatile programming language (Greenspun's tenth rule - Wikipedia).

I hope the end result of this generics activity will be to enable us to write libraries like DifferentialEquations.jl in Fortran. Dare we hope for such a thing?

[tclune]

@jacobwilliams I'll definitely take a look. As you say, there are some built in limitations of what Fortran is that will prevent going down the same path as Julia. I would hope that we'll be able to write an analogous package to DifferentialEquations.jl, but it will certainly be more verbose and probably will require more effort per extension when done.

There are days when I wish the codes that I support were all just simple F77 so that I could just wrap the old code in Julia and do all further development in Julia. Too much advanced Fortran in the code base for there to be a good solution there though. (Must admit though that my knowledge of Julia is rather thin.)

[wclodius2]

The type declarations of Fortran are not in a form that makes it straightforward to be as flexible as most generics. Ideally you want to be able declare something like:

    TYPE( ARRAY( ARRAY( REAL(DP), ALLOCATABLE ), ALLOCATABLE ) ) :: x(:)(:)

but I don't see how to get there without a lot of controversy, and redundancy in type declarations.

[tclune]

My (vague) hope is that in any one layer you won't need to be explicitly aware of the nesting. E.g., the template does not need to know that the type name parameter is itself an array, and should not be accessing the elements at that level. The F90 approach in this situation is to have a wrapper type:

type :: wrapper
      real(dp), allocatable :: inner(:)
end type

The template would then just have an array of wrapper objects.

So the (vague) idea is to greatly simplify the introduction of such a wrapper and as much as possible make the wrapper object act like the thing it wraps.

I agree that at best there will be a lot of controversy. I'm still not sure that the set of acceptable generic solutions is non-empty. :-)