How to pass variables to several levels up without going through headings (arguments) in functions

In fortran90, there are two ways to pass variables that is not through the function/subroutine heading/arguments. Those are common blocks (old Fortran77 style) or modules (new fortran90 style). This is extremely useful in large programs, where an argument (e.g. a parameter, but typically tens of them) may be used 5 levels up, but not in the intermediate calls. This facilitates good structuring (module names explain the types of variables contained and passed) and cleaner and safer code (no one wants one hundred parameters pased through the heading, level after level, to be used only at the end).

Is there a reasonable way to do this in Julia?

IMO it’s not exactly clear what are you asking for/missing (perhaps because I’m not familiar with Fortran). Do you know about Julia structs and modules?

You cannot use “using” in a function

What I want, I am missing, is to avoid passing tens of variables through all subtoutine/function levels, only to use them in the uppermost level. Imagine that at the 5th level you are using 100 parameters, that would mean that you need to pass all the 100 parameters through all the levels.

The option of declaring a variable “global” could solve it, but they recommend not using them because code is very slow then.

I think that maybe you just want composite types, a basic Julia feature? You can declare those with struct and mutable struct:

https://docs.julialang.org/en/v1/manual/types/

This is an example of the use of modules in a subroutine of a Finite Element Program (my program) in Fortran (every “use” have tens of variables, for example for dimensioning arrays according to the problem):

subroutine brck4(x0,     x,     d,     dpred,                    &  
&                matype, datam, ien,   lm,    x0l,   xl,         &  
&                xlint,  iprop, dpredl,bl,    bnl,   dmat,       &  
&                db,     strs,  strn,  cela,  imate,             &  
&                signl,  ekl,   eknl,  erf,   sbnl,              &  
&                shl,    shg,   shg0,  w,     det,   mpmat,      &
&                r  ,    shgm,  grad,  work,  tf,    det0,       &
&                vel,    vpred, apred, emass, vl,    al          )
!
!**** subroutine to compute the residual force vector and consistent
!     stifness matrix for the BRCK element
!
!     comments with # are the changes from QUAD subroutine
!
use c_consts
use c_dimens
use c_groupd
use c_apoint
use c_antype
use c_gsolve
use c_hhtpar
use c_iounit
use c_stepng
implicit none
real*8  x0(nsd,numnp)           !* original undeformed coordinates  
real*8  x(nsd,numnp)            !* current coordinates              
real*8  d(ndof,numnp)           !* incremental displacements        
real*8  dpred(ndof,numnp)       !* predictor values of the displacements
!
integer matype(2,nmat)          !* type of material of the material 
real*8  datam(ndatam,nmat)      !* material data                    
!
integer ien(nen,numel)          !* conectivity array
integer lm(nedof,nen,numel)     !* global equation locating matrix
real*8  x0l(nesd,nen)           !* initial coordinates of the element
real*8  xl(nesd,nen)            !* current coordinates of the element
real*8  xlint(nesd,nint)        !* id integration points
integer iprop(numel)            !* property user numbers for the elements
real*8  dpredl(nedof,nen)       !* predictor/incremental displacements
real*8  bl(nrowsb,nee)          !* linear B matix for each integration pt    
real*8  bnl(nrowsbnl,nee)       !* nonlinear B matrix for each int pt        
real*8  dmat(nrowsb,nrowsb)     !* material constitutive matrix              
real*8  db(nrowsb,nee)          !* matrix dmat times matrix bl               
real*8  strs(nsig,nint)         !* stress tensor in vector form              
real*8  strn(nsig,nint)         !* strain tensor in vector form              
real*8  cela(nrowsb,nrowsb)     !* elastic material matrix                  
integer imate(numel)            !* material label for the element            
real*8  signl(nsignl,nsignl)    !* stresses in special matrix form for bnl   
real*8  ekl(nee,nee)            !* element linear stiffness matrix           
real*8  eknl(nee,nee)           !* element nonlinear stifness matrix         
real*8  erf(nee)                !* element residual force                    
real*8  sbnl(nrowsbnl,nee)      !* nonlinear B matrix times special stress matrix signl
real*8  shl(nrowsh,nen,nint)    !* shape fcns derivatives and actual shape fcn (local)
real*8  shg(nrowsh,nen,nint)    !* shape fcns derivatives and actual shape fcn (global)
real*8  shg0(nrowsh,nen,nint)   !* shape fcns derivatives and actual shape fcn (global)
real*8  w(nint)                 !* weights for the integration points
real*8  det(nint)               !* determinant of jacobian at integr points
integer mpmat(mx1mpmat,numel)   !* pointer to the array of mat. data pointers
real*8  r(nint)                 !* radii at integr. points for axisim. case 
real*8  shgm(nrowsh,nen)        !* mean (b-bar) shape functions *unused*
real*8  grad(3,3)               !* deformation gradient tensor
real*8  work(*)                 !* a work array
real*8  tf(nptsltf,2,nltf)      !* time funcions                    
real*8  det0(nint)              !* initial determinant of jacobian at integr points
real*8  vel(ndof,numnp)         !* global nodal velocities
real*8  vpred(ndof,numnp)       !* global nodal velocity predictors
real*8  apred(ndof,numnp)       !* global nodal accelerations
real*8  emass(nemass,nemass)    !* element mass matrix
real*8  vl(nedof,nen)           !* element local velocities
real*8  al(nedof,nen)           !* element local accelerations
!
real*8  wd                      !* weight times jacobian
real*8  dens, dampk, dampm      !* density and stiffness and mass damping coefficients
real*8  grav(nesd)              !* prescribed constant body accelerations (gravities)
real*8  fcela                   !* factor multiplying elastic matrix (algo damping)
logical lformk                  !* if true, form stiffness matrix
logical lformm                  !* if true, form mass matrix
logical lformf                  !* if true, form residuan force
logical lmountk                 !* if true, mount stiffness matrix (indep. of lformk)
integer i 
!
real*8 gradd(nesd,nen,nint)
real*8 gradtempp(nesd,nen,nint)
!
real*8 eklj(nee,nee), eval(nee), evec(nee,nee), work1(nee), work2(nee)
!
call messages(50,'performing iterative phase tasks for BRCKs        ',-1)
!
!** loop on the elements of the group
!
do ielem = 1, numel
    !
    !** get element material label, element material type and properties material no
    !
    mat    = imate(ielem)
    matyp  = matype(1,mat)
    nprop  = iprop(ielem)
    !

It is not about a type of variable, it is about passing variables not through heading (which has to be from the inmediate level below), but other ways (skipping levels)

“using” could be the solution, but it is not allowed inside functions, only in the lowest (main) level

This code look very odd to me – how many arguments does this function have?
And in addition it needs further configuration options?!?

Probably good to rethink the design a bit. Here some ideas that could be used in Julia:

  1. Collect related parameters into a struct and just pass it as a single argument.

  2. Partially apply a function and use that instead, e.g.,

    function myloss(param1, param2, argvec)
       # compute some stuff here ...
      (param1 * argvec[1] + param2 * argvec[2])^2
    end
    config = (; param1 = 1.2, param2 = 3.4)
    lossforopt(argvec) = myloss(config.param1, config.param2, argvec)
    Optim.optimize(lossforopt, [1.0, 2.0])
    

    This is just one example, but first-class functions allow to structure code very different from Fortran and can lead to very modular designs if used well.

  3. Use dynamic binding, e.g., via ScopedValues (in Base in Julia 1.11):

    using ScopedValues
    
    function mystuff(arg1, arg2)
        # Clever computation
       param1[] * arg1 + param2[] * arg2
    end
    #  Global config
    param1 = ScopedValue(1)
    param2 = ScopedValue(2)
    
    # Run with special config
    with(param1 => 3, param2 => 4) do
        mystuff(5, 6)
    end 
    

    While dynamic binding can be convenient at times, it does have a performance hit and is used sparingly in Julia.

3 Likes

I think the point that @nsajko was trying to make is:
If you group your 100 parameters into a structure/composite type, then you can pass a single variable…

2 Likes

It is not odd, it is quite typical in fortran90 (indeed even “modern”…in the 1990s :slight_smile: ). Sorry for being so old… :slight_smile:

A good finite element code is a very large code, with a library of elements, solvers (linear and nonlinear) and (often very complex) material models. It is very common to have tens of levels and tens of matrix arguments. Those matrices have dimensions which depend on the specific problems. Usually, matrices are allocated in a huge array (the “blank common”) and passed through headings. But passing dimensions through commons or modules, as well as flags and others is very usual, because otherwise the headings would be enormous.

Serious fortran90 programmers are used to that, as well as allocating most of the arrays in the “blank common” and just passing pointers (often also in modules or commons). “Allocate” and “Dealocate” is to be avoided.

Another advantage is to have access to any array in the problem from any subroutine (only need to “use” the blank common and the pointer). This is crucial in many tasks.

Frankly, it is a real handicap not be able to do so. I have seen this question asked in other sites, but no one was able to answer. And this is very important for those programming very large codes.

Probably the struct is a workaround…

…yes, struct is probably the closest thing in readability (avoiding passing 100 parameters)…only that it has to be passed through all the levels.

Thank you everyone for the help!

In general, passing parameters explicitly makes it easier to reason about the code. If this is of no concern then on could of course also utilize global or module-level variables, like so:

julia> module Mod

       my_options = (; a=1, b=2)

       foo(c) = my_options.a + my_options.b + c

       end
Main.Mod

julia> Mod.foo(3)
6

julia> Mod.my_options = (; a=10, b=20)
(a = 10, b = 20)

julia> Mod.foo(3)
33

Naturally, this approach can also be combined with structs.

1 Like

It sounds to me like what you’re looking for is “splatting”:

julia> fn(a, b, c) = a + b + c
fn (generic function with 1 method)

julia> t = (1, 5, 7)
(1, 5, 7)

julia> fn(t...)
13

You can also splat a Vector, a Tuple is the more efficient choice if the values or number of values never changes, which is likely in your case.

You can use a NamedTuple to help keep track of which is which:

julia> t2 = (a=1, b=5, c=7)
(a = 1, b = 5, c = 7)

julia> fn(t2...)
13
1 Like

“Passing through commons or modules” is more typically known as global variables. Yes, you can have global variables in Julia too, and yes, you can put the globals into a module (to give them a separate namespace).

There are lots of software-engineering reasons why this is a bad idea, however (in any language, including Fortran — google “global variables bad” to find many discussions). As other people alluded to, it’s generally recommended (not just in Julia) to create and pass around data structures (that wrap all the data you need) instead of using globals (or passing zillions of individual parameters).

And yes, I know that lots of people use COMMON blocks and other global data in Fortran (and in many other languages). After all, Fortran 77 didn’t have user-defined struct types. (And a lot of computational-science code is written by people with little formal training in software engineering.) But just because it’s a widespread practice doesn’t make it a good idea.

But again, if you really want this in Julia, knock yourself out — for example, just make a module Foo and define some global variables x, y, z in it, then access them elsewhere via Foo.x etcetera.

Just make sure they are const or typed in order to avoid the performance pitfall of untyped globals.

Yes, pre-allocating and re-using large arrays is good (in any language). When you pass an Array variable in Julia (whether as an explicit argument or inside a data structure) it is effectively passed as a pointer (the array data is shared, not copied).

PS. What you are calling “passing through headings” is usually called “function arguments” (in all of computer science, not just Julia).

PPS. In Julia, passing array dimensions explicitly is rarely needed, because the dimensions are stored in the Array data structure. And if you want to pass a subarray, you can use a view.

8 Likes

Thank you! Very good comments. I will try Foo.x. It looks like a good option…

…wrt “function arguments” versus “passing through headings”, I used it because fortran has functions and subroutines (and subroutines change arguments); I just did not want to refer just to functions…but probably it was even less clear :slight_smile:

In both cases they are called “arguments”, even in Fortran. Calling it “passing through headings” will just confuse people.

Well, I guess we have quite different views on “good code”, i.e., I hate when I have to scroll while reading through a single function.
Julia and Fortran are quite different languages and many Fortran idioms/solutions can look very different in Julia:

  • In Julia it is easy to define custom data types and these are as efficient as build-in ones.
  • A function can return multiple values (as a tuple) at basically no cost, i.e., mutating arguments in order to pass multiple values to the caller is not required.
  • Functions are first class and can be passed as arguments, stored in vectors etc.
  • All functions are generic and support multiple dispatch – which is probably the killer feature of Julia. In particular, code can specialize on different arguments and thus, different methods, implementations can be cleanly separated:
    struct MeshLayout
        # Specific parameter for geometry go here ...
     end
     struct Material
         # Parameters concerning material go here
     end
     struct Algorithm
         # Parameters controlling the solver go here
     end
    
     function solveFEM(geometry::Mesh, stuff::Material, solver::Algorithm;
                                     preallocated = nothing)
         # Specific implementation for this combination of geometry, stuff and solver goes here
         # Also takes preallocated scratch space if available
     end
    
    (This is just an example – as I have no idea how an FEM-solver actually looks like – but should illustrate the idea)
1 Like

It is not odd, it is quite typical in fortran90 (indeed even “modern”…in the 1990s :slight_smile: )

As someone who’s been writing Fortran code for 15 years I’m going to have to push back on that. Common blocks are not typical in Fortran 90. They are typical in Fortran 77, and strongly discouraged in Fortran 90 (to the point where I don’t even mention the syntax in my Modern Fortran Reference Card). In Fortran 90, you would replace them with types or (preferably not, since not thread-save) module variables. Both of these are available in pretty much exactly the same way in Julia.

I would even say that Julia and modern Fortran are so similar that I’ve been able to port large swaths of Fortran code by copying it to a .jl file and doing some very basic search-and-replace to adjust the syntax.

A good finite element code is a very large code

More importantly, old code (Fortran 77). But if this code is still maintained, you’d want to move that over bit by bit to more modern constructs, and those would most definitely not use common blocks.

Another advantage is to have access to any array in the problem from any subroutine (only need to “use” the blank common and the pointer). This is crucial in many tasks.

That is not crucial. That is a remnant of unstructured programming from the 60s, before we had languages that support encapsulation. I understand that there is still code being used in production that uses that kind of approach, but it is something you actively want to move away from. It requires extreme discipline not to shoot yourself in the foot this way. And even then: the code will be difficult to understand, to say the least.

Frankly, it is a real handicap not be able to do so

It really isn’t. It just requires changing your mindset from Fortran 77. Trust me, you’ll be doing yourself a favor :wink:

Probably the struct is a workaround…

Yes, it is. Put those variables in a big struct and pass that around. Not just in Julia, but in Fortran as well!

4 Likes

You look young, 35 years of fortran here… :frowning: since the Digital VMS era… Used C++ and Matlab, never got along with Python, but liking Julia a lot…

But I must admit that I say COMMON but they really are modules of variables (but basically we think of them as commons). We still call the entry lines of a model as “cards” :slight_smile:

For example from my code (it is fortran90):

module blank_common
     !
     !**** This is the typical main blank common where dynamic memory is stored
     !     max_dim_bc  is a parameter that contains the maximum dimension of
     !     the blank common in words of 4 bytes.
     !     The definition of ia and ca avoids errors in some compilers
     !
     use maxdimm
     implicit none
     !
     real        a(maxdimbc)	        !* 1 real = 4 bytes, real*8 = 8 bytes
     integer     ia(maxdimbc)           !* 1 integer = 4 bytes
     character*4 ca(maxdimbc)           !* 4 characters = 4 bytes
     equivalence (a, ia, ca)            !* same stuff
     !common a                           !* passed in the blank common
end module
!
module c_antype
    !
    !** contains some option variables
    !
    implicit none
    integer iexec       !* if = 0, stop after printing dictionary, do not execute
    integer lagrange    !* formulation flag - changed by the element group
    integer iantype     !* type of analysis:
                        !*     1 = static implicit, 2 = dynamic implicit, 3 = modal
end module
module c_apoint
     !
     !** pointers to main arrays except eqn. system (array name without the first mp, see main subroutines)
     !
     implicit none
     integer         mpitseq, mprtseq,                                                 &
     &               mpnpar, mpnumelg, mpmpel, mpid, mpnelno, mpx0, mpx, mpd, mpdpred, &
     &               mplbc, mprforces, mpvel, mpapred, mpvpred, mpacel, mpdinc,        &
     &               mpsumrfor, mpsumsfor,                                             &
     &               mpf, mpfv, mpndfv, mptf, mpitf, mpietype, mpmatype, mpdatam,      &
     &               mprhs , mpihist, mphist, mphistt, mpixypl, mpidebug
end module
!
module c_bpoint
     !
     !** mpoint routine pointers for dynamic memory handling
     !
     implicit none
     integer,parameter:: n4bdic = 8     !* number of *4 bytes in dictionary per array
     integer mfirst     !* location for the next array
     integer mlast      !* position where starts array info storage
     integer mtot       !* total available space
     integer iprec      !* precision of the program
     integer mtrestart  !* total memory printed at restart files
     integer idcomm     !* if = 1, include comments on variables in memory for output
     integer nallocs    !* number of variables allocated in blank common array a
end module
!
!

Certainly, that is most probably the best solution nowadays. Thank you very much!!