Writing procedures that operate on both scalars and arrays

When a procedure is defined to operate on scalar arguments, it’s relatively straightforward to make it work with array arguments as well. For example, recall our pure function sum from the previous subsection:

pure integer function sum(a, b)
  integer, intent(in) :: a, b
  sum = a + b
end function sum

Invoking this function as, say, sum(3, 5) will evaluate to 8. Is there a way to pass array arguments to this function such that it returns an array as a result? For example, if we called sum([1, 2, 3], [2, 3, 4]), we’d get [3, 5, 7] as a result. One approach would be to declare another function that receives arrays as arguments, and we’d invoke that function instead. Fortran offers a much more elegant approach to this. You can declare the procedure as elemental, which automatically allows the scalar dummy arguments to be treated as arrays, if the arguments passed in are arrays. The result of the procedure then takes the same shape as the input arrays.

Consider the following definition:

pure elemental integer function sum(a, b)     ❶
  integer, intent(in) :: a, b
  sum = a + b
end function sum

❶ The elemental attribute allows receiving array arguments in place of scalars.

With the function defined in this way, you can pass an array as an argument to either a or b, or both. If more than one argument passed is an array, then all the array arguments have to be of the same shape.

Specifically, you can call the function like this:

print *, sum(3, 5)                     ❶
print *, sum([1, 2], 3)                ❷
print *, sum(1, [2, 3, 4])             ❸
print *, sum([1, 2, 3], [2, 3, 4])     ❹
print *, sum([1, 2], [2, 3, 4])        ❺

❶ Both arguments are scalars; evaluates to 8.

❷ Only the first argument is an array; evaluates to [4, 5].

❸ Only the second argument is an array; evaluates to [3, 4, 5].

❹ Both arguments are arrays; evaluates to [3, 5, 7].

❺ Arrays are not of the same shape–this is illegal!

These elemental snippets demonstrate different ways you can invoke an elemental function. If you try to compile the program with the last line in there (sum([1, 2], [2, 3, 4])), the compiler will raise an error. For example, gfortran reports

sum_function_elemental.f90:9:23:
 
  print *, sum([1, 2], [2, 3, 4])
                       1
Error: Different shape for elemental procedure at (1) on dimension 1 (3 and 2)

This is an important restriction of elemental procedures to keep in mind. If you pass multiple arrays as arguments to an elemental procedure, they all have to be of conforming shape.

The cold front function that we worked on earlier in this chapter is the perfect candidate for an elemental function. Try doing the exercise in the sidebar to redefine that function with the elemental attribute, and call it by passing arrays to it.

Exercise 2: Writing an elemental function that operates on both scalars and arrays

In subsection 3.2.1 we wrote a function that calculates the cold front temperature given five real scalar arguments as input (see listing 3.3). In the main program of that same listing, we used a do loop to iterate over different values of time increment to invoke the cold front temperature function at different times.

Can you use the elemental feature to redefine that function, and call it from the main program with an array of times (instead of a do loop)? For example, you could invoke the function like this:

real :: dt(8)
 
dt = [6, 12, 18, 24, 30, 36, 42, 48]
cold_front_temperature(12., 24., 20., 960., dt)

As a result, the function should return an array of the same length as dt:

22.5000000       21.0000000       19.5000000       18.0000000
   16.5000000       15.0000000       13.5000000       12.0000000

You can find the solution to this exercise in the “Answer key” section near the end of this chapter.

When you use the elemental attribute to define a procedure, it’s automatically defined as pure, even if pure is not explicitly specified. It is, however, good practice to specify both attributes for clarity.

Impure elemental?

I mentioned that an elemental procedure is automatically promoted to a pure procedure, even if the pure attribute is omitted from the procedure definition. Is the pure attribute really necessary for elemental properties? Since the Fortran 2008 standard, you can define the procedure as impure elemental. This feature is specifically designed to allow elemental behavior for nonpure procedures. In practice, you’d want to use impure elemental whenever you have a function that operates on both scalars and arrays but needs functionality that’s not permitted in pure procedures. These include I/O to and from screen or external files, calling C functions, or exchanging data with other parallel processors.