Some of the finer details of Fortran

Fortran is well described in a number of books, but not all aspects of it will become clear just from reading and studying them. This page is intended to highlight a few details that can be easily overlooked and to describe a few styles of programming that may be unfamiliar to most Fortran programmers.

Compiler extensions and limitations
Modules and generic source code
Functional programming
Using work arrays
Another look at generic source code

Compiler extensions and limitations

The Fortran standard specifies the Fortran language in detail but it also intentionally leaves a number of things up to the compiler writers. Furthermore, compiler writers can add options or features that are not covered by the standard or that relax some restrictions.

This section is intended to list a few of the aspects that can cause surprises.

First of all, not all compilers implement the complete standard. The features listed below are not supported by at least one compiler:

Then there are a number of extensions - features that are accepted by some compilers, but that would flag errors with others:

Modules and generic source code

The FLIBS collection of Fortran source files contains several files that are intended for generic purposes, most notably the files in the datastructures directory.

While "true" genericity may not be possible with the current standard (Fortran 90/95), we can achieve a lot with just a few devices that are available:

Both these techniques are used extensively in the examples/test programs. (I even wrote an article for the Fortran Forum journal using these techniques to show how design patterns, well-known from the Java and C++ world, can be applied in a Fortran programs.)

Here is an example (from "test_queue.f90"): 

module MYDATA_MODULE

type MYDATA
    character(len=20) :: string
end type MYDATA

end module

module MYDATA_QUEUES
    use MYDATA_MODULE, QUEUE_DATA => MYDATA

    include "queues.f90"
end module MYDATA_QUEUES
The code in the file "queues.f90" is completely independent of the type of data we want to store in the queues.

The main limitation of this style of genericity is that you can only store data of the same type in such queues. When you need queues for more than one data type, you can follow the method of example "two_lists.f90" in the tests directory.

Another limitation is that you can only store derived types in this way - basic types like "INTEGER" or "REAL" can not be renamed (and the syntax gets in the way: type(QUEUE_DATA) versus REAL for instance.)

Whether these limitations are prohibitive, depends entirely on the problem at hand: if you do need to store data of different types in a single "container", then you also need some way of identifying the type. And type safety is becoming an issue - will you always use the contained data in the right way (for instance pass it to routines that can handle such data)?

In principle it is possible to store data of different types in a single container like the queues above: via the TRANSER() function. But if that can be avoided, your program will probably be simpler to understand.

Importing from modules

PM

Combining interfaces

PM

Character strings with arbitrary lengths

PM

Returning a pointer?

PM

Using work arrays

Often you need a work array of a size related to the size of the problem that you are trying to solve, for instance the transpose of a matrix or a few arrays that hold the intermediate steps in a numerical integration method for systems of differential equations.

As long as the work arrays are "small", you can use automatic arrays, but the problem is that these can cause a stack overflow - resulting in an abrupt abortion of the program without a chance to catch it.

The alternative is to allocate an array (alocatable or pointer) of the appropriate size, but that has a performance penalty - memory allocation is time-consuming. You also need to take care that the memory is deallocated again (though in Fortran 95 this is automatically taken care of by the compiler under certain circumstances - allocatable arrays that do not have the target or the save attribute). The advantage, however, is that runtime errors can be caught, so that you can write a decent error message to the screen or log file.

Can these two approaches be combined? Yes, they can, as the following example shows: 

subroutine mkwork( size_problem )
    integer :: size_problem
    real, dimension(min(100,size_problem), target :: work1
    real, dimension(:), pointer                   :: work

    !
    ! Decide which of the two to use
    ! (really requires Fortran 95, because of the
    ! status of the pointer)
    !
    if ( size_problem <= 100 ) then
        work => work1
    else
        allocate( work(1:size_problem) )
    endif

    ... Do the job

    !
    ! Free the work array, when it was allocated
    !
    if ( .not. associated( work, work1 ) ) then
        deallocate( work )
    endif
end subroutine mkwork
Another way to solve the problem, is simply to use allocatable/pointer arrays with the save attribute: 
subroutine mkwork( size_problem )
    integer :: size_problem
    real, dimension(:), allocatable, save :: work

    if ( .not. allocated(work) ) then
        allocate( work(size_problem) )
    else
        if ( size_problem >= size(work) ) then
           deallocate( work )
           allocate( work(size_problem) )
        endif
    endif

    ... Do the job

    !
    ! We do not free the memory!
    !
end subroutine mkwork
That way you keep the cost of memory allocation down to a minimum: you simply allocate a new array when encountering a problem of a larger size. Any smaller problem can be solved with the previously allocated memory.

The drawback is that you have a small memory leak to worry about, but it is limited as there is only one chunk of memory per subroutine.

Measuring performance

PM

Distinguising "missing values"

PM

The perils of double precision

PM

Quantum computing

PM

Functional programming

One well-known feature that Fortran 90 introduced to the Fortran language is the ability to handle arrays as a whole, instead of merely via their elements: 
real, dimension(1:100) :: a, b

b = 2*a
This is in fact a simple form of functional programming. In languages like LISP, Haskell or ML the fundamental data type for a collection of data is the list. While not at all equivalent to an array in Fortran, these lists are often used via operations, such as map and filter that are quite easy to simulate in Fortran, or are simply present in another guise: Here are two practical examples of a functional programming style in Fortran:

Checking integer values

Here is a small function to reliably check if a string is a valid representation of a single integer value: 
logical function okay(string)
    character(len=*) :: string
    character(len=20)          :: form
    integer                    :: ierr
    integer                    :: ival

    write( form, '(a,i0,a)' ) '(i', len(string), ')'
    read( string, form, iostat = ierr ) ival

    okay = ierr == 0

end function okay
As discussed above, it is not enough to read the string with list-directed input (the very powerful "read(...,*)", you need to use formatted input.

If the string contains several integers, like "12, 34, ...", you can check the first value like this: 

logical function okay(string)
    character(len=*)            :: string
    character(len=len(string))  :: first
    character(len=20)           :: form
    integer                     :: ierr
    integer                     :: ival

    okay = .false.
    read(string, *, iostat = ierr ) first

    if ( ierr == 0 ) then
        write( form, '(a,i0,a)' ) '(i', len(string), ')'
        read( first, form, iostat = ierr ) ival
    endif

    okay = ierr == 0

end function okay

Another look at generic source code

One vexing problem of modern programming practice is that type safety and flexibility do not always go together. Here is one such dilemma:

Solve a system of ordinary differential equations: 

    x'  =  f(t,x,a),  t = 0,...,T

    x(t=0) = x0
where x is a vector of N dependent variables, f a function describing the derivative of x w.r.t. time and a is a set of parameters.

We would like to make a general library for solving this type of systems. The start is simple: Define a routine that takes a parameter f to represent the right-hand side (the name of a function subprogram) and the routine can be any convenient numerical method.

But then we run into problems: how to represent the set of parameters a? For one problem an array of reals might do, for another we would need a more complicated data structure.

The numerical method is independent of this parameter - just pass it around. But Fortran does not allow us to define a generic type. Mucking about with the transfer() function is one solution but you can not hide it from the user then.

Another problem, perhaps more mundane, is that we would like to print (output, save) intermediate results. But that leads to a second routine we would like to pass to the generic, mathematical algorithm.

Can we solve this problem in a more general, indeed generic way? Yes, we can: using reverse communication. But that leads to a complicated interface:

It gives us full flexibility in implementing the computation of the derivative and the printing of the intermediate result, we do not have to interfer with the library of numerical algorithms in any way.

But the user does have to program this loop. And we would like to hide the gory details.

This is actually possible and here is a sketch of how to do that:

-- Arjen Markus, dd. april 2008