Fortran 95

Very quickly after F90 was released, the Fortran standards group started working on the next evolution: Fortran 95 (F95) [6]. The standard was officially published in 1997 and was oriented toward resolving some outstanding issues in F90. However, it did add a new extension for HPF [7]. Also, features that were noted as "obsolete" in F90 were removed in F95.

The HPF extension of F90 focused on parallel computing [8]. HPF uses a data parallel model to spread computations across multiple processors, allowing the use of multicore processors in a node that would otherwise be idle. This model works well for both single-instruction, multiple-data (SIMD) and multiple-instruction, multiple-data (MIMD) techniques.

A few HPF additions were the FORALL statement, compiler directives for array data distribution, a procedure interface for non-HPC parallel procedures, and additional library routines, including data scattering and gathering.

The FORALL construct allows writing array functions that look more like the original mathematical formula and makes it a bit easier to change the code for different ranges without worrying that the arrays are conformable. A quick example is the basic portion of the four-point stencil for the 2D Poisson equation problem [9] (Figure 1).

Figure 1: Four-point stencil formula [10]

Assume the domain ranges from 1 to n for both i and j . Therefore the iterations are over i = 2 to n -1 and j = 2 to n -1. You can write the iteration over the domain using array notation as follows:

a(2:n-1,2:n-1) = 0.25 *  (a(1:n-2,2:n) + a(3:n,2:n) + a(2:n,1:n-2) + a(2:n,3:n))

Using forall, the same can be written as:

forall (i=2:n-1, j=2:n-1) a(i,j) = 0.25*(a(i-1,j) + a(i+1,j) + a(i,j-1) + a(i,j+1))

Using forall makes the code a little easier to read. Moreover, it's very easy to change the domain without having to modify the formula itself.

Also introduced as part of HPF were compiler directives that deal with how data is distributed for processing. They are entered into the code as comments (e.g., !HPF$). Anything immediately after !HPF$ is considered an HPF directive [11]. The nice thing about HPF directives is that if the compiler is not HPF-aware, it views the directives as comments and just moves on. This is the same approach taken by OpenACC [12] and OpenMP [13].

The HPF extensions never seemed to be popular with users, and, in fact, not many compilers offered them. While the compilers were being written, a Message Passing Interface (MPI) standard for passing data between processors, even if they weren't on the same node, was becoming popular and provided a way to utilize other processors in the same system or processors on other systems.

Additionally, OpenMP was released, providing a threading model for applications. As with HPF, OpenMP allowed coders to insert directives that told the compiler how to break up the computations into threads and exchange data between threads. These directives were code comments, so compilers that were not OpenMP-capable just ignored them.

Some observed that these two developments led to the decline of HPF, although only a few F95 compilers implemented the extensions; for example, GFortran (GNU Fortran) [14] is not a true F95-compliant compiler because it does not implement HPF, although it implements almost all of the other F95 features.

Fortran 95 added some additional features outside of HPF [15], including FORALL and nested WHERE constructs to aid vectorization, user-defined PURE and ELEMENTAL procedures, and initialization of pointers to NULL(), among others.

In F95, some features or statements were deprecated (i.e., highly recommended not to be used before the next release of the standard) [15], including banning REAL and DOUBLE PRECISION index variables in DO statements and the removal of PAUSE, ASSIGN, the controversial assigned GOTO statement, and the H descriptor.

Fortran 2003

Both F90 and F95 introduced features that allowed object-oriented programming (OOP). Fortran 2003, the next standard, added even more features to enhance the object-oriented capability.

Before Fortran 2003, you could define the basics of classes in F90 [16]. Private and public functions could be used as needed, including constructors and destructors. Several resources are available online about OOP in F90 [17] and F95 [18]. One article even compares F90 and C++ [19], so you can understand the similarities between the two.

In Fortran 2003, the ability to create type-bound procedures [20] was introduced, which makes calling functions within a class much easier. Fortran 2003 also introduced type extensions and inheritance (key aspects of OOP), as well as polymorphism and dynamic type allocation. This collection of additions really completed support for abstract data types, making OOP much easier.

To round out Fortran 2003 [21], derived type enhancements (parameterized derived types, improved control of accessibility, improved structure constructors, and finalizers) and procedure pointers were added.

Classic object-oriented languages such as C++ or Python allow you to define classes that include data and methods (functions) that operate on that data. You can instantiate (create an instance) of the class that has its own data and methods or several instances of the class that each have their own unique data and methods.

In Fortran, modules are roughly equivalent to classes. Modules can contain data and methods (functions and subroutines), but it really doesn't support the concept of separate instances of a module. To get around this, you create a module to contain the methods for the specific class with a derived type that contains the data. From this combination, you can create separate "instances" of the type that pass data to the methods in the modules.

I'll look at a quick F90 example [16] for creating and using classes by creating a module that has two public functions: One computes the area of a circle given the radius, and the other prints out the area of the circle. In the main program, I can create a derived type using the module and then use methods in the module using the derived type (Listing 8).

01  MODULE CLASS_CIRCLE
02     IMPLICIT NONE
03     PRIVATE
04     PUBLIC :: CIRCLE, CIRCLE_AREA, CIRCLE_PRINT
05
06     REAL :: PI = 3.1415926535897931d0 ! Class-wide private constant
07
08     TYPE CIRCLE
09         REAL :: RADIUS
10     END TYPE CIRCLE
11
12 CONTAINS
13
14     FUNCTION CIRCLE_AREA(THIS) RESULT(AREA)
15         TYPE(CIRCLE), INTENT(IN) :: THIS
16         REAL :: AREA
17         AREA = PI * THIS%RADIUS**2.0
18     END FUNCTION CIRCLE_AREA
19
20     SUBROUTINE CIRCLE_PRINT(THIS)
21         TYPE(CIRCLE), INTENT(IN) :: THIS
22         REAL :: AREA
23         AREA = CIRCLE_AREA(THIS)   ! Call the circle_area function
24         PRINT *, 'Circle: r = ', this%radius, ' area = ', AREA
25     END SUBROUTINE CIRCLE_PRINT
26 END MODULE CLASS_CIRCLE
27
28 PROGRAM CIRCLE_TEST
29     USE CLASS_CIRCLE
30     IMPLICIT NONE
31
32     TYPE(CIRCLE) :: C    ! Declare a variable of type Circle.
33     C = CIRCLE(1.5)      ! Use the implicit constructor, radius = 1.5
34     CALL CIRCLE_PRINT(C) ! Call a class subroutine
35
36 END PROGRAM CIRCLE_TEST

Notice (lines 32-24) that the derived type is used to create the variable C; then, I can "store" the radius in C by calling the "constructor" of the class (all it does in this case is store the value), and the program calls the "method" to print out the value by passing the derived type instance C.

In Fortran 2003, you can take advantage of some new features to make the same code look more object oriented. In the example in Listing 9, type-bound procedures can be used within the module.

01 MODULE CLASS_CIRCLE
02     IMPLICIT NONE
03     PRIVATE
04     PUBLIC :: CIRCLE, CIRCLE_AREA, CIRCLE_PRINT
05
06     REAL :: PI = 3.1415926535897931d0 ! Class-wide private constant
07
08     TYPE, PUBLIC :: CIRCLE
09         REAL :: RADIUS
10         CONTAINS
11             PROCEDURE :: AREA => CIRCLE_AREA
12             PROCEDURE :: PRINT => CIRCLE_PRINT
13     END TYPE CIRCLE
14
15 CONTAINS
16
17     FUNCTION CIRCLE_AREA(THIS) RESULT(AREA)
18         CLASS(CIRCLE), INTENT(IN) :: THIS
19         REAL :: AREA
20         AREA = PI * THIS%RADIUS**2.0
21     END FUNCTION CIRCLE_AREA
22
23     SUBROUTINE CIRCLE_PRINT(THIS)
24         CLASS(CIRCLE), INTENT(IN) :: THIS
25         REAL :: AREA
26         AREA = THIS%AREA()    ! Call the type_bound function
27         PRINT *, 'Circle: r = ', this%radius, ' area = ', AREA
28     END SUBROUTINE CIRCLE_PRINT
29 END MODULE CLASS_CIRCLE
30
31 PROGRAM CIRCLE_TEST
32     USE CLASS_CIRCLE
33     IMPLICIT NONE
34
35     TYPE(CIRCLE) :: C ! Declare a variable of type Circle.
36     C = CIRCLE(1.5)   ! Use the implicit constructor, radius = 1.5
37     CALL C%PRINT      ! Call a type-bound subroutine
38
39 END PROGRAM CIRCLE_TEST

In the module definition, I type-bound the procedures CIRCLE_AREA and CIRCLE_PRINT to the module (lines 11-12). In the main program, once I instantiate the derived type C = CIRCLE(1.5) (line 36), I can call methods (functions) as I would a derived type (CALL C%PRINT) (line 37). This looks a bit more like OOP.

In addition to object-oriented characteristics, Fortran 2003 introduced new features to enhance its ability to interact with C and C++ libraries or code. If you have tried to integrate C and Fortran in the past, you understand the fundamental differences between the two. For one, multidimensional arrays in C/C++ are in the opposite order in Fortran. Additionally, all arguments in Fortran are passed by reference and not by value, which means C must pass Fortran arguments as pointers. Finally, by default, C arrays start with 0 (zero) and Fortran arrays start with 1. Modern Fortran allows you to define arrays starting with 0 (zero) – or any integer, including negative numbers.

Integrating Fortran and C has been notoriously difficult, although possible, and linkers will blissfully create a binary for you. Moreover, it has never been easy to create portable C/Fortran integration. Sometimes, you have to test the compilers to determine exactly how they integrate.

Fortran 2003 added a standardized way to generate procedures, derived type declaration, and global variables that are interoperable with C. Because of standardization, Fortran 2003-compliant compilers make portability of C and Fortran integration much better.

The basic model for Fortran/C integration comprises three basic parts: (1) Define Fortran "kinds" that map to C types. Remember that in F90, data types started defining the precision of variables via the kind definition. (2) Use an intrinsic module to define names: USE, INTRINSIC :: ISO_C_BINDING. (3) Use a BIND(C) attribute to specify C linkage. Discussing the details is beyond the scope of this article, but look around the web if C and Fortran integration is important.

Another set of big changes in Fortran 2003 has to do with I/O. Before Fortran 2003, all I/O was record based, which worked fine for data that had been produced by other Fortran code or for reading ASCII files. However, it doesn't work well when reading data from a file that was created by an instrument or a spreadsheet or a database. Fortran I/O had to become much more general, including the ability to do "stream" I/O. With Fortran 2003, you could now do both formatted and unformatted stream I/O.

Using stream file access has several benefits. The first is that the file will likely be smaller because there are no record markers. However, this depends on how you were performing I/O before you started using streams. Second, streams allow you to read data from files you normally would not be able to read, such as databases. Third, new I/O functions allow you to move to different "positions" within the file to perform I/O.

Unformatted stream writes have no record markers. Subsequent writes just append data to the file [22]. For example, consider two WRITE() statements that write first and then second . The first write adds 5 bytes to the file, and the second write adds 6 bytes to the file. Therefore, the next write should start at byte 12. To get the current position in the file, use INQUIRE(UNIT=<fn>, POS=mypos).

If the code did not use streams, then the writes would add a record marker after each write. Therefore, the position would be something greater than 12 bytes. If you compile the code with a compiler that supports streams, the resulting position will be 12.

Using streams, you can specify a position within a file for a READ() or WRITE() function. For example, READ(11, POS=12) string reads a string at position 12 (i.e., 12 bytes into a file). The start of a file is position 0.

A couple of quick pieces of advice when using stream output: If you use it to create an output file, older Fortran code will not be able to read it (no record markers). Only code built with compilers that implement streams will be able to read the file. Conversely, if you are reading a file written by an older version of Fortran, you can't read it using streams – you have to read it as classic Fortran with record markers.

Other features added in Fortran 2003 [21] include (1) data manipulation enhancements (allocatable components, deferred type parameters, VOLATILE attribute, explicit type specification in array constructors and allocate statements, pointer enhancements, extended initialization expressions, and enhanced intrinsic procedures). (2) I/O enhancements: asynchronous transfer, stream access (previously discussed), user-specified transfer operations for derived types, user-specified control of rounding during format conversions, named constants for preconnected units, the FLUSH statement, regularization of keywords, and access to error messages. (3) Procedure pointers. (4) Support for IEEE floating-point arithmetic and floating-point exception handling. (5) Support for international usage: access to ISO 10646 4-byte characters and choice of decimal or comma in numeric-formatted I/O. (6) Enhanced integration with the host operating system: access to command-line arguments, environment variables, and processor error messages.

The Fortran developers were not done, however. Fortran 2008 incorporates some new features to improve performance, and Fortran 2015, which, although still under development, shows some serious promise.

Fortran 2008

Fortran 2003 is to Fortran 2008 much as F90 is to F95. The revision added some corrections and clarifications to Fortran 2003 while introducing new capabilities. Probably the biggest added feature was the idea of concurrent computing [23] via Coarray Fortran [24]. Concurrent computing executes several computations during overlapping periods of time, allowing users to run sections of code at the same time. In contrast, sequential computations occur one after the other and have to wait for the previous computations to finish.

As a point of a clarification, concurrent computing is related to, but different from, parallel programming, which is so prevalent in HPC. In parallel computing, various processes run at the same time. On the other hand, concurrent computing has processes that overlap in terms of duration, but their execution doesn't have to happen at the same instant. A good way to explain the differences between serial, sequential, concurrent, and parallel is shown in Table 1 [23], which assumes you have two tasks, T1 and T2.

"Simultaneous" means executed in the same physical instant, and "sequential" is an antonym of "concurrent" and "parallel." In this case, sequential typically means something is performed or used serially in a particular order.

Coarray Fortran (CAF) is a set of extensions for Fortran 95/2003 that was developed outside of the Fortran standard so that experimentation could take place quickly. The extensions were adopted in Fortran 2008, with the syntax varying a little bit relative to the original CAF definition. CAF is an example of a Partitioned Global Address Space (PGAS) [25], which assumes a global address space that is logically partitioned. Each process has a portion of the address space, which can have an affinity for a particular process.

CAF uses a parallel execution model, so code performance can be improved. A draft of the Fortran 2008 standard [26] stated that "A Fortran program containing coarrays is interpreted as if it were replicated a fixed number of times and all copies were executed asynchronously. Each copy has its own set of data objects and is called an image. The array syntax of Fortran is extended with additional trailing subscripts in square brackets to give a clear and straightforward representation of access to data on other images."

Data references without square brackets are local data (local to the image). If the square brackets are included, then the data might need to be communicated between images. CAF uses a single-program, multiple data (SPMD) model.

When a coarray-enabled application starts, a copy of the application is executed on each processor. However, the images (each copy of the application) are executed asynchronously. The images are distinguished from one another by an index between 1 and n , the number of images. Notice it starts with 1 and not 0, which is perhaps influenced by the Fortran roots.

In the application code, you define a coarray with a trailing []. The coarray then exists in each image, with the same name and having the same size. They can be scalar, array, static, or dynamic and of intrinsic or derived type. Synchronization statements can be used to maintain program correctness. Table 2 shows examples of the ways you can define a coarray. Notice that scalars can be coarrays, fixed-size arrays, allocatable arrays, and derived types, and you can have coarrays with different coranks. However, you cannot have a coarray that is a constant or a pointer.

What is a "corank"? Variables in Fortran code can have rank , bounds , extent , size , and shape , which are defined inside the parentheses of an array declaration, as specified from F90 onward. For coarrays, the corank , cobounds , and coextents are given by data in square brackets. The cosize of a coarray is always equal to the number of images specified when the application is executed. A coarray has a final coextent , a final upper cobound , and a coshape that depends on the number of images.

A simple example declaration is integer :: a(4)[*]. Here, each image ("process") is an integer array of 4. You can assemble a 2D (rank = 2) coarray data structure from these 1D (rank = 1) arrays. If using four images with the previous declaration, the coarray looks like the illustration in Figure 2. The arrays are stacked on top of each other because of the 1D coarray specification [*].

Figure 2: Two-dimensional array from 1D arrays: integer :: a(4)*.

If you specify the coarray in a different manner, but still use four images, a 2D coarray can be declared as integer :: a(4)[2,*]. Each image is again a 1D integer with length (extent) 4. However, the coarray specifies that the 2D array is assembled differently (as in Figure 3). Remember that each image has it's own array that is the same as the others but can be combined using coarrays.

Figure 3: Two-dimensional array from 1D arrays: integer :: a(4)2,*.

The use of coarrays can be thought of as opposite the way distributed arrays are used in MPI. With MPI applications, each rank or process has a local array; then, the process needs to be mapped from the local array to the global array so that local data can be mapped to the larger global array. The starting point is with global arrays and then to local arrays.

Coarrays are the opposite. You can take local arrays and combine them into a global array using a coarray. You can access the local array (local to the image) using the "usual" array notation. You can also access data on another image almost in the same way, but you have to use the coarray notation.

In another simple example, each image has a 2D array that is combined into a coarray to create a larger 2D array. Again, assume four images with the coarray declaration integer :: a(4,4)[2,*]. Each image has a 2D integer array a. With the coarray definition given by the square bracket notation, the "local" arrays are combined into a coarray (globally accessible) that looks like Figure 4.

Figure 4: Two-dimensional array from 2D arrays: integer :: a(4,4)2,*.

If the array is local to the image, you access the array as you normally would. For example, for image 3 to access element (2,2) from the array, the statement would be something like b = a(2,2). You can always use coarray notation if needed, but in this case, you know the data is local, so you can access it using local notation. If image 1, 2, or 4 wanted to access that element (global access), then the statement would have to be b = a(2,2)[1,2].

You still have to pay attention to what image holds what data, but writing the statements to access the data is fairly easy. The key is the coarray subscripts. The declaration integer, dimension(10,4) :: a is a simple local variable with rank 2 (i.e., two indices) and has lower bounds 1 and 1. The upper bounds are 10 and 4. The shape is [10,4].

The declaration integer :: a(10,4)[3,*] would convert the array to a coarray, which adds a corank of 2 (two indices). The lower cobounds are 1 and 1 . The upper cobounds are 3 and m , and the coshape is , where m = ceiling(num_images()/3) .

Using coarrays, you no longer have to use MPI_SEND() and MPI_RECV() or their non-blocking equivalents for sending data from one process to another. You just access the remote data using the coarray syntax and let the underlying coarray library do the work, making multiprocess coding easier.

CAF for Fortran 2008 can be implemented in different ways. A common way is to implement it using MPI, as in GFortran [27]. Except for a small number of functions, GFortran provides coarray functionality starting in Fortran 2008 v5.1 using the coarray capability in OpenCoarrays [28]. The compiler translates the coarray syntax into library calls that then use MPI functions underneath.