Fortran
Fortran
Fortran
FORTRAN (FORmula TRANslation) is a third-generation ( 3GL ) programming language that was
designed for use by engineers, mathematicians, and other users and creators of scientific algorithms.
It has a very succinct and spartan syntax. Today, the C language has largely displaced FORTRAN.
History
One of the oldest programming languages, the FORTRAN was developed by a team of programmers at
IBM led by John Backus, and was first published in 1957. The name FORTRAN is an acronym for
FORmula TRANslation, because it was designed to allow easy translation of math formulas into code.
Often referred to as a scientific language, FORTRAN was the first high-level language, using the first
compiler ever developed. Prior to the development of FORTRAN computer programmers were
required to program in machine/assembly code, which was an extremely difficult and time consuming
task, not to mention the dreadful chore of debugging the code. The objective during it's design was to
create a programming language that would be: simple to learn, suitable for a wide variety of
applications, machine independent, and would allow complex mathematical expressions to be stated
similarly to regular algebraic notation. While still being almost as efficient in execution as assembly
language. Since FORTRAN was so much easier to code, programmers were able to write programs
500% faster than before, while execution efficiency was only reduced by 20%, this allowed them to
focus more on the problem solving aspects of a problem, and less on coding.
FORTRAN was so innovative not only because it was the first high-level language, but also because of
it's compiler, which is credited as giving rise to the branch of computer science now known
as compiler theory. Several years after it's release FORTRAN had developed many
different dialects, (due to special tweaking by programmers trying to make it better suit their personal
needs) making it very difficult to transfer programs from one machine to another.
Despite this standardization, a few years later, various new dialects began to surface again, requiring
the Standards Association review the language again. This version is known as FORTRAN '77. This
version was released in 1978 (it was called '77 because the Association began it's review in 1977),
with several new features. Some of the more notable properties were; new error handling methods,
and mechanisms for managing large-scale programs. The latest version; Fortran '90 (released in 1990,
using the new capitalization scheme) added even more new features, such as support for: recursion,
pointers, and for programmer-defined data types. {Fortran 90's future - Current research in
complier theory involves equipping compilers to generate object code, that is able to exploit the
capabilities of massively parallel computers. Thr Fortran 90 compilers are key targets of such
research}
Areas of Application
FORTRAN is useful for a wide variety of applications, some of the more outstanding ones are as
follows:
Number crunching - due to the more natural (like it's true algebraic form) way of expressing
complex mathematical functions and it's quick execution time, FORTRAN is easy and efficient
at processing mathematical equations.
Scientific, mathematical, statistical, and engineering type procedures - due to it's rapid
number-crunching ability FORTRAN is a good choice for these type of applications.
Basically FORTRAN is most useful for applications that are "computational-bound" rather than "I/O
bound".
Command Description
g95 -c h1.f90 h2.f90 h3.f90 Compiles multiple source files. If all goes well, object
files h1.o, h2.o and h3.o are created
g95 -o hello h1.f90 h2.f90 h3.f90 Compiles multiple source files and links them together
to an executable file named 'hello'
Each program contains one main program and may or may not contain other program units. The
syntax of the main program is as follows:
program program_name
implicit none
! type declaration statements
! executable statements
end program program_name
program addNumbers
! This simple program adds two numbers
implicit none
! Type declarations
real :: a, b, result
! Executable statements
a = 12.0
b = 15.0
result = a + b
print *, 'The total is ', result
end program addNumbers
When you compile and execute the above program, it produces the following result:
The total is 27.0000000
Please note that:
All
Fortran programs start with the keyword program and end with the keyword end
program, followed by the name of the program.
The implicit
none statement allows the compiler to check that all your variable types are
declared properly. You must always use implicit none at the start of every program.
Comments in Fortran are started with the exclamation mark (!), as all characters after this
(except in a character string) are ignored by the compiler.
Fortranallows both uppercase and lowercase letters. Fortran is case-insensitive, except for
string literals.
Basics
The basic character set of Fortran contains:
Identifier
An identifier is a name used to identify a variable, procedure, or any other user-defined item. A name
in Fortran must follow the following rules:
It must
be composed of alphanumeric characters (all the letters of the alphabet, and the digits
0 to 9) and underscores (_).
Keywords
Keywords are special words, reserved for the language. These reserved words cannot be used as
identifiers or names.
The following table, lists the Fortran keywords:
end interface end module end program end select end subroutine
Where While
Integer type
Real type
Complex type
Logical type
Character type
Integer Type
The integer types can hold only integer values. The following example extracts the largest value that
can be held in a usual four byte integer:
program testingInt
implicit none
integer :: largeval
print *, huge(largeval)
end program testingInt
When you compile and execute the above program it produces the following result:
2147483647
Note that the huge() function gives the largest number that can be held by the specific integer data
type. You can also specify the number of bytes using the kind specifier. The following example
demonstrates this:
program testingInt
implicit none
!two byte integer
integer(kind=2) :: shortval
!four byte integer
integer(kind=4) :: longval
!eight byte integer
integer(kind=8) :: verylongval
!sixteen byte integer
integer(kind=16) :: veryverylongval
!default integer
integer :: defval
print *, huge(shortval)
print *, huge(longval)
print *, huge(verylongval)
print *, huge(veryverylongval)
print *, huge(defval)
end program testingInt
When you compile and execute the above program, it produces the following result:
32767
2147483647
9223372036854775807
170141183460469231731687303715884105727
2147483647
Real Type
It stores the floating point numbers, such as 2.0, 3.1415, -100.876, etc.
Traditionally there are two different real types, the default real type and double precision type.
However, Fortran 90/95 provides more control over the precision of real and integer data types
through thekindspecifier, which we will study in the chapter on Numbers.
program division
implicit none
! Define real variables
real :: p, q, realRes
! Define integer variables
integer :: i, j, intRes
! Assigning values
p = 2.0
q = 3.0
i = 2
j = 3
! floating point division
realRes = p/q
intRes = i/j
print *, realRes
print *, intRes
end program division
When you compile and execute the above program it produces the following result:
0.666666687
0
Complex Type
This is used for storing complex numbers. A complex number has two parts, the real part and the
imaginary part. Two consecutive numeric storage units store these two parts.
For example, the complex number (3.0, -5.0) is equal to 3.0 – 5.0i
Logical Type
There are only two logical values: .true. and .false.
Character Type
The character type stores characters and strings. The length of the string can be specified by len
specifier. If no length is specified, it is 1.
For example,
character (len=40) :: name
name = “Zara Ali”
The expression, name(1:4) would give the substring “Zara”.
Implicit Typing
Older versions of Fortran allowed a feature called implicit typing, i.e., you do not have to declare the
variables before use. If a variable is not declared, then the first letter of its name will determine its
type.
Variable names starting with i, j, k, l, m, or n, are considered to be for integer variable and others are
real variables. However, you must declare all the variables as it is good programming practice. For
that you start your program with the statement:
implicit none
This statement turns off implicit typing.
Fortran – Variables
A variable is nothing but a name given to a storage area that our programs can manipulate. Each
variable should have a specific type, which determines the size and layout of the variable's memory;
the range of values that can be stored within that memory; and the set of operations that can be
applied to the variable.
The name of a variable can be composed of letters, digits, and the underscore character. A name in
Fortran must follow the following rules:
It cannot be longer than 31 characters.
It must
be composed of alphanumeric characters (all the letters of the alphabet, and the digits
0 to 9) and underscores (_).
Based on the basic types explained in previous chapter, following are the variable types:
Type Description
Variable Declaration
Variables are declared at the beginning of a program (or subprogram) in a type declaration
statement.
type-specifier :: variable_name
For example,
integer :: total
real :: average
complex :: cx
logical :: done
character(len=80) :: message ! a string of 80 characters
Later you can assign values to these variables, like,
total = 20000
average = 1666.67
done = .true.
message = “A big Hello from Softecks”
cx = (3.0, 5.0) ! cx = 3.0 + 5.0i
You can also use the intrinsic function cmplx, to assign values to a complex variable:
cx = cmplx (1.0/2.0, -7.0) ! cx = 0.5 – 7.0i
cx = cmplx (x, y) ! cx = x + yi
Example
The following example demonstrates variable declaration, assignment and display on screen:
program variableTesting
implicit none
! declaring variables
integer :: total
real :: average
complex :: cx
logical :: done
character(len=80) :: message ! a string of 80 characters
!assigning values
total = 20000
average = 1666.67
done = .true.
message = "A big Hello from Softecks"
cx = (3.0, 5.0) ! cx = 3.0 + 5.0i
Print *, total
Print *, average
Print *, cx
Print *, done
Print *, message
end program variableTesting
When the above code is compiled and executed, it produces the following result:
20000
1666.67004
(3.00000000, 5.00000000 )
T
A big Hello from Softecks
Fortran – Constants
The constants refer to the fixed values that the program cannot alter during its execution. These fixed
values are also called literals.
Constants can be of any of the basic data types like an integer constant, a floating constant, a
character constant, a complex constant, or a string literal. There are only two logical
constants : .true. and .false.
The constants are treated just like regular variables, except that their values cannot be modified after
their definition.
Literal constants
Named constants
A literal constant have a value, but no name.
Type Example
Named constants should be declared at the beginning of a program or procedure, just like a variable
type declaration, indicating its name and type. Named constants are declared with the parameter
attribute. For example,
real, parameter :: pi = 3.1415927
Example
The following program calculates the displacement due to vertical motion under gravity.
program gravitationalDisp
! this program calculates vertical motion under gravity
implicit none
! gravitational acceleration
real, parameter :: g = 9.81
! variable declaration
real :: s ! displacement
real :: t ! time
real :: u ! initial speed
! assigning values
t = 5.0
u = 50
! displacement
s = u * t - g * (t**2) / 2
! output
print *, "Time = ", t
print *, 'Displacement = ',s
end program gravitationalDisp
When the above code is compiled and executed, it produces the following result:
Time = 5.00000000
Displacement = 127.374992
Fortran – Operators
An operator is a symbol that tells the compiler to perform specific mathematical or logical
manipulations. Fortran provides the following types of operators:
Arithmetic Operators
Relational Operators
Logical Operators
Let us look at all these types of operators one by one.
Arithmetic Operators
Following table shows all the arithmetic operators supported by Fortran. Assume variable A holds 5
and variable B holds 3 then:
** Exponentiation Operator, raises one operand to the power A ** B will give 125
of the other.
Relational Operators
Following table shows all the relational operators supported by Fortran. Assume variable A holds 10
and variable B holds 20, then:
> .gt. Checks if the value of left operand is greater than (A > B) is
the value of right operand, if yes then condition not true.
becomes true.
< .lt. Checks if the value of left operand is less than the (A < B) is
value of right operand, if yes then condition true.
becomes true.
>= .ge. Checks if the value of left operand is greater than or (A >= B) is
equal to the value of right operand, if yes then not true.
condition becomes true.
<= .le. Checks if the value of left operand is less than or (A <= B) is
equal to the value of right operand, if yes then true.
condition becomes true.
For example, x = 7 + 3 * 2; here, x is assigned 13, not 20 because operator * has higher precedence
than +, so it first gets multiplied with 3*2 and then adds into 7.
Here, operators with the highest precedence appear at the top of the table, those with the lowest
appear at the bottom. Within an expression, higher precedence operators will be evaluated first.
Fortran – Decisions
Decision making structures require that the programmer specify one or more conditions to be
evaluated or tested by the program, along with a statement or statements to be executed, if the
condition is determined to be true, and optionally, other statements to be executed if the condition is
determined to be false.
Following is the general form of a typical decision making structure found in most of the
programming languages:
Statement Description
nested select case construct You can use one select case statement inside
another select case statement(s).
The optional else is placed at the end and it is executed when none of the above conditions hold true.
[name:]
if (logical expression 1) then
! block 1
else if (logical expression 2) then
! block 2
else if (logical expression 3) then
! block 3
else
! block 4
end if [name]
Example
program ifElseIfElseProg
implicit none
! local variable declaration
integer :: a = 100
! check the logical condition using if statement
if( a == 10 ) then
! if condition is true then print the following
print*, "Value of a is 10"
else if( a == 20 ) then
! if else if condition is true
print*, "Value of a is 20"
else if( a == 30 ) then
! if else if condition is true
print*, "Value of a is 30"
else
! if none of the conditions is true
print*, "None of the values is matching"
end if
print*, "exact value of a is ", a
end program ifElseIfElseProg
When the above code is compiled and executed, it produces the following result:
Fortran – Loops
There may be a situation, when you need to execute a block of code several number of times. In
general, statements are executed sequentially : The first statement in a function is executed first,
followed by the second, and so on.
Programming languages provide various control structures that allow for more complicated
execution paths.
A loop statement allows us to execute a statement or group of statements multiple times and
following is the general form of a loop statement in most of the programming languages:
Fortran provides the following types of loop constructs to handle looping requirements. Click the
following links to check their detail.
do while loop Repeats a statement or group of statements while a given condition is true. It
tests the condition before executing the loop body.
nested loops You can use one or more loop construct inside any other loop construct.
Fortran supports the following control statements. Click the following links to check their detail.
exit If the exit statement is executed, the loop is exited, and the execution
of the program continues at the first executable statement after the
end do statement.
stop If you wish execution of your program to stop, you can insert a stop
statement
Fortran – Numbers
Numbers in Fortran are represented by three intrinsic data types:
Integer type
Real type
Complex type
Integer Type
The integer types can hold only integer values. The following example extracts the largest value that
could be hold in a usual four byte integer:
program testingInt
implicit none
integer :: largeval
print *, huge(largeval)
end program testingInt
When you compile and execute the above program it produces the following result:
2147483647
Please note that the huge() function gives the largest number that can be held by the specific integer
data type. You can also specify the number of bytes using the kind specifier. The following example
demonstrates this:
program testingInt
implicit none
!two byte integer
integer(kind=2) :: shortval
!four byte integer
integer(kind=4) :: longval
!eight byte integer
integer(kind=8) :: verylongval
!sixteen byte integer
integer(kind=16) :: veryverylongval
!default integer
integer :: defval
print *, huge(shortval)
print *, huge(longval)
print *, huge(verylongval)
print *, huge(veryverylongval)
print *, huge(defval)
end program testingInt
When you compile and execute the above program it produces the following result:
32767
2147483647
9223372036854775807
170141183460469231731687303715884105727
2147483647
Real Type
It stores the floating point numbers, such as 2.0, 3.1415, -100.876, etc.
Traditionally there were two different real types : the default real type and double precision type.
However, Fortran 90/95 provides more control over the precision of real and integer data types
through the kind specifier, which we will study shortly.
program division
implicit none
! Define real variables
real :: p, q, realRes
! Define integer variables
integer :: i, j, intRes
! Assigning values
p = 2.0
q = 3.0
i = 2
j = 3
! floating point division
realRes = p/q
intRes = i/j
print *, realRes
print *, intRes
end program division
When you compile and execute the above program it produces the following result:
0.666666687
0
Complex Type
This is used for storing complex numbers. A complex number has two parts : the real part and the
imaginary part. Two consecutive numeric storage units store these two parts.
For example, the complex number (3.0, -5.0) is equal to 3.0 – 5.0i
The generic function cmplx() creates a complex number. It produces a result who’s real and
imaginary parts are single precision, irrespective of the type of the input arguments.
program createComplex
implicit none
integer :: i = 10
real :: x = 5.17
print *, cmplx(i, x)
end program createComplex
When you compile and execute the above program it produces the following result:
(10.0000000, 5.17000008)
The following program demonstrates complex number arithmetic:
program ComplexArithmatic
implicit none
complex, parameter :: i = (0, 1) ! sqrt(-1)
complex :: x, y, z
x = (7, 8);
y = (5, -7)
write(*,*) i * x * y
z = x + y
print *, "z = x + y = ", z
z = x - y
print *, "z = x - y = ", z
z = x * y
print *, "z = x * y = ", z
z = x / y
print *, "z = x / y = ", z
end program ComplexArithmatic
When you compile and execute the above program it produces the following result:
(9.00000000, 91.0000000)
z = x + y = (12.0000000, 1.00000000)
z = x - y = (2.00000000, 15.0000000)
z = x * y = (91.0000000, -9.00000000)
z = x / y = (-0.283783793, 1.20270276)
The following table displays the number of bits and range for integers:
64 9,223,372,036,854,774,807 (2**63)–1
32 2,147,483,647 (2**31)–1
The following table displays the number of bits, smallest and largest value, and the precision for real
numbers.
program rangePrecision
implicit none
real:: x, y, z
x = 1.5e+40
y = 3.73e+40
z = x * y
print *, z
end program rangePrecision
When you compile and execute the above program it produces the following result:
x = 1.5e+40
1
Error : Real cos its kind at (1)
main.f95:5.12:
y = 3.73e+40
1
Error : Real constant overflows its kind at (1)
Now let us use a smaller number:
program rangePrecision
implicit none
real:: x, y, z
x = 1.5e+20
y = 3.73e+20
z = x * y
print *, z
z = x/y
print *, z
end program rangePrecision
When you compile and execute the above program it produces the following result:
Infinity
0.402144760
Now let’s watch underflow:
program rangePrecision
implicit none
real:: x, y, z
x = 1.5e-30
y = 3.73e-60
z = x * y
print *, z
z = x/y
print *, z
end program rangePrecision
When you compile and execute the above program it produces the following result:
y = 3.73e-60
1
Warning : Real constant underflows its kind at (1)
Executing the program....
$demo
0.00000000E+00
Infinity
The Kind Specifier
In scientific programming, one often needs to know the range and precision of data of the hardware
platform on which the work is being done.
The intrinsic function kind() allows you to query the details of the hardware’s data representations
before running a program.
program kindCheck
implicit none
integer :: i
real :: r
complex :: cp
print *,' Integer ', kind(i)
print *,' Real ', kind(r)
print *,' Complex ', kind(cp)
end program kindCheck
When you compile and execute the above program it produces the following result:
Integer 4
Real 4
Complex 4
You can also check the kind of all data types:
program checkKind
implicit none
integer :: i
real :: r
character*1 :: c
logical :: lg
complex :: cp
print *,' Integer ', kind(i)
print *,' Real ', kind(r)
print *,' Complex ', kind(cp)
print *,' Character ', kind(c)
print *,' Logical ', kind(lg)
end program checkKind
When you compile and execute the above program it produces the following result:
Integer 4
Real 4
Complex 4
Character 1
Logical 4
Fortran – Characters
The Fortran language can treat characters as single character or contiguous strings.
Characters could be any symbol taken from the basic character set, i.e., from the letters, the decimal
digits, the underscore, and 21 special characters.
The intrinsic data type character stores characters and strings. The length of the string can be
specified by len specifier. If no length is specified, it is 1. You can refer individual characters within a
string referring by position; the left most character is at position 1.
Character Declaration
Declaring a character type data is same as other variables:
type-specifier :: variable_name
For example,
character :: reply, sex
you can assign a value like,
reply = ‘N’
sex = ‘F’
The following example demonstrates declaration and use of character data type:
program hello
implicit none
character(len=15) :: surname, firstname
character(len=6) :: title
character(len=25)::greetings
title = 'Mr. '
firstname = 'Rowan '
surname = 'Atkinson'
greetings = 'A big hello from Mr. Bean'
print *, 'Here is ', title, firstname, surname
print *, greetings
end program hello
When you compile and execute the above program it produces the following result:
Here is Mr. Rowan Atkinson
A big hello from Mr. Bean
Concatenation of Characters
The concatenation operator //, concatenates characters.
program hello
implicit none
character(len=15) :: surname, firstname
character(len=6) :: title
character(len=40):: name
character(len=25)::greetings
title = 'Mr. '
firstname = 'Rowan '
surname = 'Atkinson'
name = title//firstname//surname
greetings = 'A big hello from Mr. Bean'
print *, 'Here is ', name
print *, greetings
end program hello
When you compile and execute the above program it produces the following result:
Here is Mr.Rowan Atkinson
A big hello from Mr.Bean
Function Description
verify(string, chars) It scans the "string" from left to right (unless back=.true.) for the first
occurrence of any character not contained in "chars". It returns an
integer giving the position of that character, or zero if only the
characters in "chars" have been found
repeat(string,ncopy) It returns a string with length equal to "ncopy" times the length of
"string", and containing "ncopy" concatenated copies of "string"
Example 1
This example shows the use of the index function:
program testingChars
implicit none
character (80) :: text
integer :: i
text = 'The intrinsic data type character stores characters and strings.'
i=index(text,'character')
if (i /= 0) then
print *, ' The word character found at position ',i
print *, ' in text: ', text
end if
end program testingChars
When you compile and execute the above program it produces the following result:
Example 2
This example demonstrates the use of the trim function:
program hello
implicit none
character(len=15) :: surname, firstname
character(len=6) :: title
character(len=25)::greetings
title = 'Mr.'
firstname = 'Rowan'
surname = 'Atkinson'
print *, 'Here is', title, firstname, surname
print *, 'Here is', trim(title),' ',trim(firstname),' ', trim(surname)
end program hello
When you compile and execute the above program it produces the following result:
Here is Mr. Rowan Atkinson
Here is Mr. Rowan Atkinson
Example 3
This example demonstrates the use of achar function
program testingChars
implicit none
character:: ch
integer:: i
do i=65, 90
ch = achar(i)
print*, i, ' ', ch
end do
end program testingChars
When you compile and execute the above program it produces the following result:
65 A
66 B
67 C
68 D
69 E
70 F
71 G
72 H
73 I
74 J
75 K
76 L
77 M
78 N
79 O
80 P
81 Q
82 R
83 S
84 T
85 U
86 V
87 W
88 X
89 Y
90 Z
Function Description
lle(char, char) Compares whether the first character is lexically less than or equal to the
second
lge(char, char) Compares whether the first character is lexically greater than or equal to
the second
lgt(char, char) Compares whether the first character is lexically greater than the second
llt(char, char) Compares whether the first character is lexically less than the second
Example 4
The following function demonstrates the use:
program testingChars
implicit none
character:: a, b, c
a = 'A'
b = 'a'
c = 'B'
if(lgt(a,b)) then
print *, 'A is lexically greater than a'
else
print *, 'a is lexically greater than A'
end if
if(lgt(a,c)) then
print *, 'A is lexically greater than B'
else
print *, 'B is lexically greater than A'
end if
if(llt(a,b)) then
print *, 'A is lexically less than a'
end if
if(llt(a,c)) then
print *, 'A is lexically less than B'
end if
end program testingChars
When you compile and execute the above program it produces the following result:
Fortran – Strings
The Fortran language can treat characters as single character or contiguous strings.
A character string may be only one character in length, or it could even be of zero length. In Fortran,
character constants are given between a pair of double or single quotes.
The intrinsic data type character stores characters and strings. The length of the string can be
specified by len specifier. If no length is specified, it is 1. You can refer individual characters within a
string referring by position; the left most character is at position 1.
String Declaration
Declaring a string is same as other variables:
type-specifier :: variable_name
For example,
Character(len=20) :: firstname, surname
you can assign a value like,
character (len=40) :: name
name = “Zara Ali”
The following example demonstrates declaration and use of character data type:
program hello
implicit none
character(len=15) :: surname, firstname
character(len=6) :: title
character(len=25)::greetings
title = 'Mr.'
firstname = 'Rowan'
surname = 'Atkinson'
greetings = 'A big hello from Mr. Beans'
print *, 'Here is', title, firstname, surname
print *, greetings
end program hello
When you compile and execute the above program it produces the following result:
Here is Mr. Rowan Atkinson
A big hello from Mr. Bean
String Concatenation
The concatenation operator //, concatenates strings.
program hello
implicit none
character(len=15) :: surname, firstname
character(len=6) :: title
character(len=40):: name
character(len=25)::greetings
title = 'Mr.'
firstname = 'Rowan'
surname = 'Atkinson'
name = title//firstname//surname
greetings = 'A big hello from Mr. Beans'
print *, 'Here is', name
print *, greetings
end program hello
When you compile and execute the above program it produces the following result:
Here is Mr. Rowan Atkinson
A big hello from Mr. Bean
Extracting Substrings
In Fortran, you can extract a substring from a string by indexing the string, giving the start and the
end index of the substring in a pair of brackets. This is called extent specifier.
The following example shows how to extract the substring ‘world’ from the string ‘hello world’:
program subString
character(len=11)::hello
hello = "Hello World"
print*, hello(7:11)
end program subString
When you compile and execute the above program it produces the following result:
World
Example
The following example uses the date_and_time function to give the date and time string. We use
extent specifiers to extract the year, date, month, hour, minutes and second information separately.
program datetime
implicit none
character(len = 8) :: dateinfo ! ccyymmdd
character(len = 4) :: year, month*2, day*2
character(len = 10) :: timeinfo ! hhmmss.sss
character(len = 2) :: hour, minute, second*6
call date_and_time(dateinfo, timeinfo)
! let’s break dateinfo into year, month and day.
! dateinfo has a form of ccyymmdd, where cc = century, yy = year
! mm = month and dd = day
year = dateinfo(1:4)
month = dateinfo(5:6)
day = dateinfo(7:8)
print*, 'Date String:', dateinfo
print*, 'Year:', year
print *,'Month:', month
print *,'Day:', day
! let’s break timeinfo into hour, minute and second.
! timeinfo has a form of hhmmss.sss, where h = hour, m = minute
! and s = second
hour = timeinfo(1:2)
minute = timeinfo(3:4)
second = timeinfo(5:10)
print*, 'Time String:', timeinfo
print*, 'Hour:', hour
print*, 'Minute:', minute
print*, 'Second:', second
end program datetime
When you compile and execute the above program, it gives the detailed date and time information:
Date String: 20140803
Year: 2014
Month: 08
Day: 03
Time String: 075835.466
Hour: 07
Minute: 58
Second: 35.466
Trimming Strings
The trim function takes a string, and returns the input string after removing all trailing blanks.
Example
program trimString
implicit none
character (len=*), parameter :: fname="Susanne", sname="Rizwan"
character (len=20) :: fullname
fullname=fname//" "//sname !concatenating the strings
print*,fullname,", the beautiful dancer from the east!"
print*,trim(fullname),", the beautiful dancer from the east!"
end program trimString
When you compile and execute the above program it produces the following result:
The function adjustr takes a string and returns it by removing the trailing blanks and appending
them as leading blanks.
Example
program hello
implicit none
character(len=15) :: surname, firstname
character(len=6) :: title
character(len=40):: name
character(len=25):: greetings
title = 'Mr. '
firstname = 'Rowan'
surname = 'Atkinson'
greetings = 'A big hello from Mr. Beans'
name = adjustl(title)//adjustl(firstname)//adjustl(surname)
print *, 'Here is', name
print *, greetings
name = adjustr(title)//adjustr(firstname)//adjustr(surname)
print *, 'Here is', name
print *, greetings
name = trim(title)//trim(firstname)//trim(surname)
print *, 'Here is', name
print *, greetings
end program hello
When you compile and execute the above program it produces the following result:
Here is Mr. Rowan Atkinson
A big hello from Mr. Bean
Here is Mr. Rowan Atkinson
A big hello from Mr. Bean
Here is Mr.RowanAtkinson
A big hello from Mr. Bean
Example
program hello
implicit none
character(len=30) :: myString
character(len=10) :: testString
myString = 'This is a test'
testString = 'test'
if(index(myString, testString) == 0)then
print *, 'test is not found'
else
print *, 'test is found at index: ', index(myString, testString)
end if
end program hello
When you compile and execute the above program it produces the following result:
test is found at index: 11
Fortran – Arrays
Arrays can store a fixed-size sequential collection of elements of the same type. An array is used to
store a collection of data, but it is often more useful to think of an array as a collection of variables of
the same type.
All arrays consist of contiguous memory locations. The lowest address corresponds to the first
element and the highest address to the last element.
Numbers(1) Numbers(2) Numbers(3) Numbers(4)
Arrays can be one- dimensional (like vectors), two-dimensional (like matrices) and Fortran allows
you to create up to 7-dimensional arrays.
Declaring Arrays
Arrays are declared with the dimension attribute.
For example, to declare a one-dimensional array named number, of real numbers containing 5
elements, you write,
real, dimension(5) :: numbers
The individual elements of arrays are referenced by specifying their subscripts. The first element of
an array has a subscript of one. The array numbers contains five real variables –numbers(1),
numbers(2), numbers(3), numbers(4), and numbers(5).
integer, dimension (5,5) :: matrix
You can also declare an array with some explicit lower bound, for example:
real, dimension(2:6) :: numbers
integer, dimension (-3:2,0:4) :: matrix
Assigning Values
You can either assign values to individual members, like,
numbers(1) = 2.0
or, you can use a loop,
do i=1,5
numbers(i) = i * 2.0
end do
One dimensional array elements can be directly assigned values using a short hand symbol, called
array constructor, like,
numbers = (/1.5, 3.2,4.5,0.9,7.2 /)
please note that there are no spaces allowed between the brackets ‘( ‘and the back slash ‘/’
Example
The following example demonstrates the concepts discussed above.
program arrayProg
real :: numbers(5) !one dimensional integer array
integer :: matrix(3,3), i , j !two dimensional real array
!assigning some values to the array numbers
do i=1,5
numbers(i) = i * 2.0
end do
!display the values
do i = 1, 5
Print *, numbers(i)
end do
!assigning some values to the array matrix
do i=1,3
do j = 1, 3
matrix(i, j) = i+j
end do
end do
!display the values
do i=1,3
do j = 1, 3
Print *, matrix(i,j)
end do
end do
!short hand assignment
numbers = (/1.5, 3.2,4.5,0.9,7.2 /)
!display the values
do i = 1, 5
Print *, numbers(i)
end do
end program arrayProg
When the above code is compiled and executed, it produces the following result:
2.00000000
4.00000000
6.00000000
8.00000000
10.0000000
2
3
4
3
4
5
4
5
6
1.50000000
3.20000005
4.50000000
0.899999976
7.19999981
Term Meaning
Rank It is the number of dimensions an array has. For example, for the array named
matrix, rank is 2, and for the array named numbers, rank is 1.
Exten It is the number of elements along a dimension. For example, the array numbers has
t extent 5 and the array named matrix has extent 3 in both dimensions.
Shape The shape of an array is a one-dimensional integer array, containing the number of
elements (the extent) in each dimension. For example, for the array matrix, shape is
(3, 3) and the array numbers it is (5).
Size It is the number of elements an array contains. For the array matrix, it is 9, and for
the array numbers, it is 5.
program arrayToProcedure
implicit none
integer, dimension (5) :: myArray
integer :: i
call fillArray (myArray)
call printArray(myArray)
end program arrayToProcedure
subroutine fillArray (a)
implicit none
integer, dimension (5), intent (out) :: a
! local variables
integer :: i
do i = 1, 5
a(i) = i
end do
end subroutine fillArray
subroutine printArray(a)
integer, dimension (5) :: a
integer::i
do i = 1, 5
Print *, a(i)
end do
end subroutine printArray
When the above code is compiled and executed, it produces the following result:
1
2
3
4
5
In the above example, the subroutine fillArray and printArray can only be called with arrays with
dimension 5. However, to write subroutines that can be used for arrays of any size, you can rewrite it
using the following technique:
program arrayToProcedure
implicit none
integer, dimension (10) :: myArray
integer :: i
interface
subroutine fillArray (a)
integer, dimension(:), intent (out) :: a
integer :: i
end subroutine fillArray
subroutine printArray (a)
integer, dimension(:) :: a
integer :: i
end subroutine printArray
end interface
call fillArray (myArray)
call printArray(myArray)
end program arrayToProcedure
subroutine fillArray (a)
implicit none
integer,dimension (:), intent (out) :: a
! local variables
integer :: i, arraySize
arraySize = size(a)
do i = 1, arraySize
a(i) = i
end do
end subroutine fillArray
subroutine printArray(a)
implicit none
integer,dimension (:) :: a
integer::i, arraySize
arraySize = size(a)
do i = 1, arraySize
Print *, a(i)
end do
end subroutine printArray
Please note that the program is using the size function to get the size of the array.
When the above code is compiled and executed, it produces the following result:
1
2
3
4
5
6
7
8
9
10
Array Sections
So far we have referred to the whole array, Fortran provides an easy way to refer several elements,
or a section of an array, using a single statement.
To access an array section, you need to provide the lower and the upper bound of the section, as well
as a stride (increment), for all the dimensions. This notation is called a subscript triplet:
array ([lower]:[upper][:stride], ...)
When no lower and upper bounds are mentioned, it defaults to the extents you declared, and stride
value defaults to 1.
program arraySubsection
real, dimension(10) :: a, b
integer:: i, asize, bsize
a(1:7) = 5.0 ! a(1) to a(7) assigned 5.0
a(8:) = 0.0 ! rest are 0.0
b(2:10:2) = 3.9
b(1:9:2) = 2.5
!display
asize = size(a)
bsize = size(b)
do i = 1, asize
Print *, a(i)
end do
do i = 1, bsize
Print *, b(i)
end do
end program arraySubsection
When the above code is compiled and executed, it produces the following result:
5.00000000
5.00000000
5.00000000
5.00000000
5.00000000
5.00000000
5.00000000
0.00000000E+00
0.00000000E+00
0.00000000E+00
2.50000000
3.90000010
2.50000000
3.90000010
2.50000000
3.90000010
2.50000000
3.90000010
2.50000000
3.90000010
Reduction
Inquiry
Construction
Reshape
Manipulation
Location
For example,
real, dimension (:,:), allocatable :: darray
The rank of the array, i.e., the dimensions has to be mentioned however, to allocate memory to such
an array, you use the allocate function.
allocate ( darray(s1,s2) )
After the array is used, in the program, the memory created should be freed using
the deallocate function
deallocate (darray)
Example
The following example demonstrates the concepts discussed above.
program dynamic_array
implicit none
!rank is 2, but size not known
real, dimension (:,:), allocatable :: darray
integer :: s1, s2
integer :: i, j
print*, "Enter the size of the array:"
read*, s1, s2
! allocate memory
allocate ( darray(s1,s2) )
do i = 1, s1
do j = 1, s2
darray(i,j) = i*j
print*, "darray(",i,",",j,") = ", darray(i,j)
end do
end do
deallocate (darray)
end program dynamic_array
When the above code is compiled and executed, it produces the following result:
data variable / list / ...
Example
The following example demonstrates the concept:
program dataStatement
implicit none
integer :: a(5), b(3,3), c(10),i, j
data a /7,8,9,10,11/
data b(1,:) /1,1,1/
data b(2,:)/2,2,2/
data b(3,:)/3,3,3/
data (c(i),i=1,10,2) /4,5,6,7,8/
data (c(i),i=2,10,2)/5*2/
Print *, 'The A array:'
do j = 1, 5
print*, a(j)
end do
Print *, 'The B array:'
do i = lbound(b,1), ubound(b,1)
write(*,*) (b(i,j), j = lbound(b,2), ubound(b,2))
end do
Print *, 'The C array:'
do j = 1, 10
print*, c(j)
end do
end program dataStatement
When the above code is compiled and executed, it produces the following result:
The A array:
7
8
9
10
11
The B array:
1 1 1
2 2 2
3 3 3
The C array:
4
2
5
2
6
2
7
2
8
2
Example
The following example demonstrates the concept:
program whereStatement
implicit none
integer :: a(3,5), i , j
do i = 1,3
do j = 1, 5
a(i,j) = j-i
end do
end do
Print *, 'The A array:'
do i = lbound(a,1), ubound(a,1)
write(*,*) (a(i,j), j = lbound(a,2), ubound(a,2))
end do
where( a<0 )
a = 1
elsewhere
a = 5
end where
Print *, 'The A array:'
do i = lbound(a,1), ubound(a,1)
write(*,*) (a(i,j), j = lbound(a,2), ubound(a,2))
end do
end program whereStatement
When the above code is compiled and executed, it produces the following result:
The A array:
0 1 2 3 4
-1 0 1 2 3
-2 -1 0 1 2
The A array:
5 5 5 5 5
1 5 5 5 5
1 1 5 5 5
Derived data types are used to represent a record. E.g. you want to keep track of your books in a
library, you might want to track the following attributes about each book:
Title
Author
Subject
Book ID
Defining a Derived data type
To define a derived data type, the type and end type statements are used. . The type statement
defines a new data type, with more than one member for your program. The format of the type
statement is this:
type type_name
declarations
end type
Here is the way you would declare the Book structure:
type Books
character(len=50) :: title
character(len=50) :: author
character(len=150) :: subject
integer :: book_id
end type Books
type(Books) :: book1
The components of the structure can be accessed using the component selector character (%):
book1%title = "C Programming"
book1%author = "Nuha Ali"
book1%subject = "C Programming practice"
book1%book_id = 123456
Note that there are no spaces before and after the % symbol.
Example
The following program illustrates the above concepts:
program deriveDataType
!type declaration
type Books
character(len=50) :: title
character(len=50) :: author
character(len=150) :: subject
integer :: book_id
end type Books
!declaring type variables
type(Books) :: book1
type(Books) :: book2
!accessing the components of the structure
book1%title = "C Programming"
book1%author = "Nuha Ali"
book1%subject = "C Programming practice"
book1%book_id = 123456
book2%title = "Telecom Billing"
book2%author = "Zara Ali"
book2%subject = "Telecom Billing practice"
book2%book_id = 654321
!display book info
Print *, book1%title
Print *, book1%author
Print *, book1%subject
Print *, book1%book_id
Print *, book2%title
Print *, book2%author
Print *, book2%subject
Print *, book2%book_id
end program deriveDataType
When the above code is compiled and executed, it produces the following result:
C Programming
Nuha Ali
C Programming practice
123456
Telecom Billing
Zara Ali
Telecom Billing practice
654321
Array of Structures
You can also create arrays of a derived type:
type(Books), dimension(2) :: list
Individual elements of the array could be accessed as:
list(1)%title = "C Programming"
list(1)%author = "Nuha Ali"
list(1)%subject = "C Programming practice"
list(1)%book_id = 123456
The following program illustrates the concept:
program deriveDataType
!type declaration
type Books
character(len=50) :: title
character(len=50) :: author
character(len=150) :: subject
integer :: book_id
end type Books
!declaring array of books
type(Books), dimension(2) :: list
!accessing the components of the structure
list(1)%title = "C Programming"
list(1)%author = "Nuha Ali"
list(1)%subject = "C Programming practice"
list(1)%book_id = 123456
list(2)%title = "Telecom Billing"
list(2)%author = "Zara Ali"
list(2)%subject = "Telecom Billing practice"
list(2)%book_id = 654321
!display book info
Print *, list(1)%title
Print *, list(1)%author
Print *, list(1)%subject
Print *, list(1)%book_id
Print *, list(1)%title
Print *, list(2)%author
Print *, list(2)%subject
Print *, list(2)%book_id
end program deriveDataType
When the above code is compiled and executed, it produces the following result:
C Programming
Nuha Ali
C Programming practice
123456
C Programming
Zara Ali
Telecom Billing practice
654321
Fortran – Pointers
In most programming languages, a pointer variable stores the memory address of an object.
However, in Fortran, a pointer is a data object that has more functionalities than just storing the
memory address. It contains more information about a particular object, like type, rank, extents, and
memory address.
integer, pointer :: p1 ! pointer to integer
real, pointer, dimension (:) :: pra ! pointer to 1-dim real array
real, pointer, dimension (:,:) :: pra2 ! pointer to 2-dim real array
A pointer can point to:
1
5
You should empty the allocated storage space by the deallocate statement when it is no longer
required and avoid accumulation of unused and unusable memory space.
You associate a pointer variable with a target variable using the association operator (=>).
program pointerExample
implicit none
integer, pointer :: p1
integer, target :: t1
p1=>t1
p1 = 1
Print *, p1
Print *, t1
p1 = p1 + 4
Print *, p1
Print *, t1
t1 = 8
Print *, p1
Print *, t1
end program pointerExample
When the above code is compiled and executed, it produces the following result:
1
1
5
5
8
8
A pointer can be:
Undefined
Associated
Disassociated
In the above program, we have associated the pointer p1, with the target t1, using the => operator.
The function associated, tests a pointer’s association status.
Nullify does not empty the targets as there could be more than one pointer pointing to the same
target. However, emptying the pointer implies nullification also.
Example 1
The following example demonstrates the concepts:
program pointerExample
implicit none
integer, pointer :: p1
integer, target :: t1
integer, target :: t2
p1=>t1
p1 = 1
Print *, p1
Print *, t1
p1 = p1 + 4
Print *, p1
Print *, t1
t1 = 8
Print *, p1
Print *, t1
nullify(p1)
Print *, t1
p1=>t2
Print *, associated(p1)
Print*, associated(p1, t1)
Print*, associated(p1, t2)
!what is the value of p1 at present
Print *, p1
Print *, t2
p1 = 10
Print *, p1
Print *, t2
end program pointerExample
When the above code is compiled and executed, it produces the following result:
1
1
5
5
8
8
8
T
F
T
952754640
952754640
10
10
Please note that each time you run the code, the memory addresses will be different.
Example 2
program pointerExample
implicit none
integer, pointer :: a, b
integer, target :: t
integer :: n
t= 1
a=>t
t = 2
b => t
n = a + b
Print *, a, b, t, n
end program pointerExample
When the above code is compiled and executed, it produces the following result:
2 2 2 4
read(*,*) item1, item2, item3...
print *, item1, item2, item3
write(*,*) item1, item2, item3...
However the formatted I/O gives you more flexibility over data transfer.
read fmt, variable_list
print fmt, variable_list
write fmt, variable_list
Where,
Format specification defines the way in which formatted data is displayed. It consists of a string,
containing a list of edit descriptors in parentheses.
An edit descriptor specifies the exact format, for example, width, digits after decimal point etc., in
which characters and numbers are displayed.
For example:
Print "(f6.3)", pi
The following table describes the descriptors:
I This is used for integer output. This takes the form print "(3i5)", i, j, k
‘rIw.m’ where the meanings of r, w and m are
given in the table below. Integer values are right
justified in their fields. If the field width is not
large enough to accommodate an integer then the
field is filled with asterisks.
F This is used for real number output. This takes the print "(f12.3)",pi
form ‘rFw.d’ where the meanings of r, w and d are
given in the table below. Real values are right
justified in their fields. If the field width is not
large enough to accommodate the real number
then the field is filled with asterisks.
A This is used for character output. This takes the print "(a10)", str
form ‘rAw’ where the meanings of r and w are
given in the table below. Character types are right
justified in their fields. If the field width is not
large enough to accommodate the character string
then the field is filled with the first ‘w’ characters
of the string.
X This is used for space output. This takes the form print "(5x, a10)", str
‘nX’ where ‘n’ is the number of desired spaces.
/ Slash descriptor – used to insert blank lines. This print "(/,5x, a10)", str
takes the form ‘/’ and forces the next data output
to be on a new line.
Following symbols are used with the format descriptors:
Symbol Description
c Column number
d Number of digits to right of the decimal place for real input or output
w Field width – the number of characters to use for the input or output
Example 1
program printPi
pi = 3.141592653589793238
Print "(f6.3)", pi
Print "(f10.7)", pi
Print "(f20.15)", pi
Print "(e16.4)", pi/100
end program printPi
When the above code is compiled and executed, it produces the following result:
3.142
3.1415927
3.141592741012573
0.3142E-01
Example 2
program printName
implicit none
character (len=15) :: first_name
print *,' Enter your first name.'
print *,' Up to 20 characters, please'
read *,first_name
print "(1x,a)",first_name
end program printName
When the above code is compiled and executed, it produces the following result: (assume the user
enters the name Zara)
Example 3
program formattedPrint
implicit none
real :: c = 1.2786456e-9, d = 0.1234567e3
integer :: n = 300789, k = 45, i = 2
character (len=15) :: str="Tutorials Point"
print "(i6)", k
print "(i6.3)", k
print "(3i10)", n, k, i
print "(i10,i3,i5)", n, k, i
print "(a15)",str
print "(f12.3)", d
print "(e12.4)", c
print '(/,3x,"n = ",i6, 3x, "d = ",f7.4)', n, d
end program formattedPrint
When the above code is compiled and executed, it produces the following result:
45
045
300789 45 2
300789 45 2
Tutorials Point
123.457
0.1279E-08
n = 300789 d = *******
program productDetails
implicit none
character (len=15) :: name
integer :: id
real :: weight
name = 'Ardupilot'
id = 1
weight = 0.08
print *,' The product details are'
print 100
100 format (7x,'Name:', 7x, 'Id:', 1x, 'Weight:')
print 200, name, id, weight
200 format(1x, a, 2x, i3, 2x, f5.2)
end program productDetails
When the above code is compiled and executed, it produces the following result:
In the last chapter, you have seen how to read data from, and write data to the terminal. In this
chapter you will study file input and output functionalities provided by Fortran.
You can read and write to one or more files. The OPEN, WRITE, READ and CLOSE statements allow
you to achieve this.
open (unit = number, file = "name").
However, the open statement may have a general form:
open (list-of-specifiers)
The following table describes the most commonly used specifiers:
Specifier Description
[UNIT=] u The unit number u could be any number in the range 9-99 and it indicates
the file, you may choose any number but every open file in the program
must have a unique number
IOSTAT= ios It is the I/O status identifier and should be an integer variable. If the open
statement is successful then the ios value returned is zero else a non-zero
value.
ERR = err It is a label to which the control jumps in case of any error.
STATUS = sta It shows the prior status of the file. A character string and can have one of
the three values NEW, OLD or SCRATCH. A scratch file is created and
deleted when closed or the program ends.
ACCESS = acc It is the file access mode. Can have either of the two values, SEQUENTIAL
or DIRECT. The default is SEQUENTIAL.
FORM= frm It gives the formatting status of the file. Can have either of the two values
FORMATTED or UNFORMATTED. The default is UNFORMATTED
After the file has been opened, it is accessed by read and write statements. Once done, it should be
closed using the close statement.
close ([UNIT=]u[,IOSTAT=ios,ERR=err,STATUS=sta])
Please note that the parameters in brackets are optional.
Example
This example demonstrates opening a new file for writing some data into the file.
program outputdata
implicit none
real, dimension(100) :: x, y
real, dimension(100) :: p, q
integer :: i
! data
do i=1,100
x(i) = i * 0.1
y(i) = sin(x(i)) * (1-cos(x(i)/3.0))
end do
! output data into a file
open(1, file='data1.dat', status='new')
do i=1,100
write(1,*) x(i), y(i)
end do
close(1)
end program outputdata
When the above code is compiled and executed, it creates the file data1.dat and writes the x and y
array values into it. And then closes the file.
The read and write statements respectively are used for reading from and writing into a file
respectively.
read ([UNIT=]u, [FMT=]fmt, IOSTAT=ios, ERR=err, END=s)
write([UNIT=]u, [FMT=]fmt, IOSTAT=ios, ERR=err, END=s)
Most of the specifiers have already been discussed in the above table.
The END=s specifier is a statement label where the program jumps, when it reaches end-of-file.
Example
This example demonstrates reading from and writing into a file.
In this program we read from the file, we created in the last example, data1.dat, and display it on
screen.
program outputdata
implicit none
real, dimension(100) :: x, y
real, dimension(100) :: p, q
integer :: i
! data
do i=1,100
x(i) = i * 0.1
y(i) = sin(x(i)) * (1-cos(x(i)/3.0))
end do
! output data into a file
open(1, file='data1.dat', status='new')
do i=1,100
write(1,*) x(i), y(i)
end do
close(1)
! opening the file for reading
open (2, file='data1.dat', status='old')
do i=1,100
read(2,*) p(i), q(i)
end do
close(2)
do i=1,100
write(*,*) p(i), q(i)
end do
end program outputdata
When the above code is compiled and executed, it produces the following result:
0.100000001 5.54589933E-05
0.200000003 4.41325130E-04
0.300000012 1.47636665E-03
0.400000006 3.45637114E-03
0.500000000 6.64328877E-03
0.600000024 1.12552457E-02
0.699999988 1.74576249E-02
0.800000012 2.53552198E-02
0.900000036 3.49861123E-02
1.00000000 4.63171229E-02
1.10000002 5.92407547E-02
1.20000005 7.35742599E-02
1.30000007 8.90605897E-02
1.39999998 0.105371222
1.50000000 0.122110792
1.60000002 0.138823599
1.70000005 0.155002072
1.80000007 0.170096487
1.89999998 0.183526158
2.00000000 0.194692180
2.10000014 0.202990443
2.20000005 0.207826138
2.29999995 0.208628103
2.40000010 0.204863414
2.50000000 0.196052119
2.60000014 0.181780845
2.70000005 0.161716297
2.79999995 0.135617107
2.90000010 0.103344671
3.00000000 6.48725405E-02
3.10000014 2.02930309E-02
3.20000005 -3.01767997E-02
3.29999995 -8.61928314E-02
3.40000010 -0.147283033
3.50000000 -0.212848678
3.60000014 -0.282169819
3.70000005 -0.354410470
3.79999995 -0.428629100
3.90000010 -0.503789663
4.00000000 -0.578774154
4.09999990 -0.652400017
4.20000029 -0.723436713
4.30000019 -0.790623367
4.40000010 -0.852691114
4.50000000 -0.908382416
4.59999990 -0.956472993
4.70000029 -0.995793998
4.80000019 -1.02525222
4.90000010 -1.04385209
5.00000000 -1.05071592
5.09999990 -1.04510069
5.20000029 -1.02641726
5.30000019 -0.994243503
5.40000010 -0.948338211
5.50000000 -0.888650239
5.59999990 -0.815326691
5.70000029 -0.728716135
5.80000019 -0.629372001
5.90000010 -0.518047631
6.00000000 -0.395693362
6.09999990 -0.263447165
6.20000029 -0.122622721
6.30000019 2.53026206E-02
6.40000010 0.178709000
6.50000000 0.335851669
6.59999990 0.494883657
6.70000029 0.653881252
6.80000019 0.810866773
6.90000010 0.963840425
7.00000000 1.11080539
7.09999990 1.24979746
7.20000029 1.37891412
7.30000019 1.49633956
7.40000010 1.60037732
7.50000000 1.68947268
7.59999990 1.76223695
7.70000029 1.81747139
7.80000019 1.85418403
7.90000010 1.87160957
8.00000000 1.86922085
8.10000038 1.84674001
8.19999981 1.80414569
8.30000019 1.74167395
8.40000057 1.65982044
8.50000000 1.55933595
8.60000038 1.44121361
8.69999981 1.30668485
8.80000019 1.15719533
8.90000057 0.994394958
9.00000000 0.820112705
9.10000038 0.636327863
9.19999981 0.445154816
9.30000019 0.248800844
9.40000057 4.95488606E-02
9.50000000 -0.150278628
9.60000038 -0.348357052
9.69999981 -0.542378068
9.80000019 -0.730095863
9.90000057 -0.909344316
10.0000000 -1.07807255
Fortran – Procedures
A procedure is a group of statements that perform a well-defined task and can be invoked from your
program. Information (or data) is passed to the calling program, to the procedure as arguments.
Functions
Subroutines
Function
A function is a procedure that returns a single quantity. A function should not modify its arguments.
The returned quantity is known as function value, and it is denoted by the function name.
Syntax:
Syntax for a function is as follows:
function name(arg1, arg2, ....)
[declarations, including those for the arguments]
[executable statements]
end function [name]
The following example demonstrates a function named area_of_circle. It calculates the area of a circle
with radius r.
program calling_func
real :: a
a = area_of_circle(2.0)
Print *, "The area of a circle with radius 2.0 is"
Print *, a
end program calling_func
function area_of_circle (r)
! this function computes the area of a circle with radius r
implicit none
! function result
real :: area_of_circle
! dummy arguments
real :: r
! local variables
real :: pi
pi = 4 * atan (1.0)
area_of_circle = pi * r**2
end function area_of_circle
When you compile and execute the above program, it produces the following result:
The area of a circle with radius 2.0 is
12.5663710
Please note that:
You must specify implicit none in both the main program as well as the procedure.
function name(arg1, arg2, ....) result (return_var_name)
[declarations, including those for the arguments]
[executable statements]
end function [name]
Subroutine
A subroutine does not return a value, however it can modify its arguments.
Syntax
subroutine name(arg1, arg2, ....)
[declarations, including those for the arguments]
[executable statements]
end subroutine [name]
Calling a Subroutine
You need to invoke a subroutine using the call statement.
The following example demonstrates the definition and use of a subroutine swap, that changes the
values of its arguments.
program calling_func
implicit none
real :: a, b
a = 2.0
b = 3.0
Print *, "Before calling swap"
Print *, "a = ", a
Print *, "b = ", b
call swap(a, b)
Print *, "After calling swap"
Print *, "a = ", a
Print *, "b = ", b
end program calling_func
subroutine swap(x, y)
implicit none
real :: x, y, temp
temp = x
x = y
y = temp
end subroutine swap
When you compile and execute the above program, it produces the following result:
Before calling swap
a = 2.00000000
b = 3.00000000
After calling swap
a = 3.00000000
b = 2.00000000
program calling_func
implicit none
real :: x, y, z, disc
x= 1.0
y = 5.0
z = 2.0
call intent_example(x, y, z, disc)
Print *, "The value of the discriminant is"
Print *, disc
end program calling_func
subroutine intent_example (a, b, c, d)
implicit none
! dummy arguments
real, intent (in) :: a
real, intent (in) :: b
real, intent (in) :: c
real, intent (out) :: d
d = b*b - 4.0*a*c
end subroutine intent_example
When you compile and execute the above program, it produces the following result:
Recursive Procedures
Recursion occurs when a programming languages allows you to call a function inside the same
function. It is called recursive call of the function.
When a procedure calls itself, directly or indirectly, is called a recursive procedure. You should
declare this type of procedures by preceding the word recursive before its declaration.
Following is an example, which calculates factorial for a given number using a recursive procedure:
program calling_func
implicit none
integer :: i, f
i = 15
Print *, "The value of factorial 15 is"
f = myfactorial(15)
Print *, f
end program calling_func
recursive function myfactorial (n) result (fac)
! computes the factorial of n (n!)
implicit none
! function result
integer :: fac
! dummy arguments
integer, intent (in) :: n
select case (n)
case (0:1)
fac = 1
case default
fac = n * myfactorial (n-1)
end select
end function myfactorial
Internal Procedures
When a procedure is contained within a program, it is called the internal procedure of the program.
The syntax for containing an internal procedure is as follows:
program program_name
implicit none
! type declaration statements
! executable statements
. . .
contains
! internal procedures
. . .
end program program_name
The following example demonstrates the concept:
program mainprog
implicit none
real :: a, b
a = 2.0
b = 3.0
Print *, "Before calling swap"
Print *, "a = ", a
Print *, "b = ", b
call swap(a, b)
Print *, "After calling swap"
Print *, "a = ", a
Print *, "b = ", b
contains
subroutine swap(x, y)
real :: x, y, temp
temp = x
x = y
y = temp
end subroutine swap
end program mainprog
When you compile and execute the above program, it produces the following result:
Before calling swap
a = 2.00000000
b = 3.00000000
After calling swap
a = 3.00000000
b = 2.00000000
Fortran – Modules
A module is like a package where you can keep your functions and subroutines, in case you are
writing a very big program, or your functions or subroutines can be used in more than one program.
Modules provide you a way of splitting your programs between multiple files.
module name
[statement declarations]
[contains [subroutine and function definitions] ]
end module [name]
use name
Please note that
You
can add as many modules as needed, each will be in separate files and compiled
separately.
The variables declared in a module specification part, are global to the module.
The
variables declared in a module become global variables in any program or routine where
the module is used.
The
use statement can appear in the main program, or any other subroutine or module which
uses the routines or variables declared in a particular module.
Example
The following example demonstrates the concept:
module constants
implicit none
real, parameter :: pi = 3.1415926536
real, parameter :: e = 2.7182818285
contains
subroutine show_consts()
print*, "Pi = ", pi
print*, "e = ", e
end subroutine show_consts
end module constants
program module_example
use constants
implicit none
real :: x, ePowerx, area, radius
x = 2.0
radius = 7.0
ePowerx = e ** x
area = pi * radius**2
call show_consts()
print*, "e raised to the power of 2.0 = ", ePowerx
print*, "Area of a circle with radius 7.0 = ", area
end program module_example
When you compile and execute the above program, it produces the following result:
Pi = 3.14159274
e = 2.71828175
e raised to the power of 2.0 = 7.38905573
Area of a circle with radius 7.0 = 153.938049
However, you can control the accessibility of module code using the private and public attributes.
When you declare some variable or subroutine as private, it is not available outside the module.
Example
The following example illustrates the concept:
In the previous example, we had two module variables, e and pi. Let us make them private and
observe the output:
module constants
implicit none
real, parameter,private :: pi = 3.1415926536
real, parameter, private :: e = 2.7182818285
contains
subroutine show_consts()
print*, "Pi = ", pi
print*, "e = ", e
end subroutine show_consts
end module constants
program module_example
use constants
implicit none
real :: x, ePowerx, area, radius
x = 2.0
radius = 7.0
ePowerx = e ** x
area = pi * radius**2
call show_consts()
print*, "e raised to the power of 2.0 = ", ePowerx
print*, "Area of a circle with radius 7.0 = ", area
end program module_example
When you compile and execute the above program, it gives the following error message:
ePowerx = e ** x
1
Error: Symbol 'e' at (1) has no IMPLICIT type
main.f95:19.13:
area = pi * radius**2
1
Error: Symbol 'pi' at (1) has no IMPLICIT type
Since e and pi, both are declared private, the program module_example cannot access these variables
anymore.
module constants
implicit none
real, parameter,private :: pi = 3.1415926536
real, parameter, private :: e = 2.7182818285
contains
subroutine show_consts()
print*, "Pi = ", pi
print*, "e = ", e
end subroutine show_consts
function ePowerx(x)result(ePx)
implicit none
real::x
real::ePx
ePx = e ** x
end function ePowerx
function areaCircle(r)result(a)
implicit none
real::r
real::a
a = pi * r**2
end function areaCircle
end module constants
program module_example
use constants
implicit none
call show_consts()
Print*, "e raised to the power of 2.0 = ", ePowerx(2.0)
print*, "Area of a circle with radius 7.0 = ", areaCircle(7.0)
end program module_example
When you compile and execute the above program, it produces the following result:
Pi = 3.14159274
e = 2.71828175
e raised to the power of 2.0 = 7.38905573
Area of a circle with radius 7.0 = 153.938049
Numeric Functions
Mathematical Functions
Numeric Inquiry Functions
Floating-Point Manipulation Functions
Bit Manipulation Functions
Character Functions
Kind Functions
Logical Functions
Array Functions
We have discussed the array functions in the Arrays chapter. In the following section we provide
brief descriptions of all these functions from other categories.
ANINT (A [, KIND]) It returns a real value, the nearest integer or whole number.
CEILING (A [, KIND]) It returns the least integer greater than or equal to number A.
CMPLX (X [, Y, KIND]) It converts the real variables X and Y to a complex number X+iY; if Y
is absent, 0 is used.
FLOOR (A [, KIND]) It provides the greatest integer less than or equal to number A.
INT (A [, KIND]) It converts a number (real or integer) to integer, truncating the real
part towards zero.
MAX (A1, A2 [, A3,...]) It returns the maximum value from the arguments, all being of
same type.
MIN (A1, A2 [, A3,...]) It returns the minimum value from the arguments, all being of same
type.
Example
program numericFunctions
implicit none
! define constants
!define variables
real :: a, b
complex :: z
!values for a, b
a = 15.2345
b = -20.7689
write(*,*) 'abs(a): ',abs(a),' abs(b): ',abs(b)
write(*,*) 'aint(a): ',aint(a),' aint(b): ',aint(b)
write(*,*) 'ceiling(a): ',ceiling(a),' ceiling(b): ',ceiling(b)
write(*,*) 'floor(a): ',floor(a),' floor(b): ',floor(b)
z = cmplx(a, b)
write(*,*) 'z: ',z
end program numericFunctions
When you compile and execute the above program, it produces the following result:
abs(a): 15.2344999 abs(b): 20.7688999
aint(a): 15.0000000 aint(b): -20.0000000
ceiling(a): 16 ceiling(b): -20
floor(a): 15 floor(b): -21
z: (15.2344999, -20.7688999)
Mathematical Functions
Function Description
ACOS (X) It returns the inverse cosine in the range (0, π), in radians.
ASIN (X) It returns the inverse sine in the range (-π/2, π/2), in radians.
ATAN (X) It returns the inverse tangent in the range (-π/2, π/2), in radians.
ATAN2 (Y, X) It returns the inverse tangent in the range (-π, π), in radians.
Example
The following program computes the horizontal and vertical position x and y respectively of a
projectile after a time, t:
program projectileMotion
implicit none
! define constants
real, parameter :: g = 9.8
real, parameter :: pi = 3.1415927
!define variables
real :: a, t, u, x, y
!values for a, t, and u
a = 45.0
t = 20.0
u = 10.0
! convert angle to radians
a = a * pi / 180.0
x = u * cos(a) * t
y = u * sin(a) * t - 0.5 * g * t * t
write(*,*) 'x: ',x,' y: ',y
end program projectileMotion
When you compile and execute the above program, it produces the following result:
x: 141.421356 y: -1818.57861
Function Description
EPSILON (X) It returns the number that is almost negligible compared to one. In
other words, it returns the smallest value such that REAL( 1.0,
KIND(X)) + EPSILON(X) is not equal to REAL( 1.0, KIND(X)).
RRSPACING (X) It returns the reciprocal of the relative spacing of model numbers
near given number
SPACING (X) It returns the absolute spacing of model numbers near given
number
Character Functions
Function Description
ACHAR (I) It returns the Ith character in the ASCII collating sequence.
ADJUSTL (STRING) It adjusts string left by removing any leading blanks and
inserting trailing blanks
CHAR (I [, KIND]) It returns the Ith character in the machine specific collating
sequence
IACHAR (C) It returns the position of the character in the ASCII collating
sequence.
INDEX (STRING, SUBSTRING [, It returns the leftmost (rightmost if BACK is .TRUE.) starting
BACK]) position of SUBSTRING within STRING.
SCAN (STRING, SET [, BACK]) It returns the index of the leftmost (rightmost if BACK is
.TRUE.) character of STRING that belong to SET, or 0 if none
belong.
Kind Functions
Function Description
SELECTED_REAL_KIND ([P, R]) Real kind type parameter value, given precision and range
Logical Function
Function Description
LOGICAL (L [, KIND]) Convert between objects of type logical with different kind type
parameters
However, Fortran 90/95 provides more control over the precision of real and integer data types
through the kind specifie.
real, kind = 2 :: a, b, c
real, kind = 4 :: e, f, g
integer, kind = 2 :: i, j, k
integer, kind = 3 :: l, m, n
In the above declaration, the real variables e, f and g have more precision than the real variables a, b
and c. The integer variables l, m and n, can store larger values and have more digits for storage than
the integer variables i, j and k. Although this is machine dependent.
Example
program kindSpecifier
implicit none
real(kind = 4) :: a, b, c
real(kind = 8) :: e, f, g
integer(kind = 2) :: i, j, k
integer(kind = 4) :: l, m, n
integer :: kind_a, kind_i, kind_e, kind_l
kind_a = kind(a)
kind_i = kind(i)
kind_e = kind(e)
kind_l = kind(l)
print *,'default kind for real is', kind_a
print *,'default kind for int is', kind_i
print *,'extended kind for real is', kind_e
print *,'default kind for int is', kind_l
end program kindSpecifier
When you compile and execute the above program it produces the following result:
default kind for real is 4
default kind for int is 2
extended kind for real is 8
default kind for int is 4
For example, the bit_size(i) intrinsic function specifies the number of bits used for storage. For real
numbers, the precision(x) intrinsic function, returns the number of decimal digits of precision, while
the range(x) intrinsic function returns the decimal range of the exponent.
Example
program getSize
implicit none
real (kind = 4) :: a
real (kind = 8) :: b
integer (kind = 2) :: i
integer (kind = 4) :: j
print *,'precision of real(4) =', precision(a)
print *,'precision of real(8) =', precision(b)
print *,'range of real(4) =', range(a)
print *,'range of real(8) =', range(b)
print *,'maximum exponent of real(4) =' , maxexponent(a)
print *,'maximum exponent of real(8) =' , maxexponent(b)
print *,'minimum exponent of real(4) =' , minexponent(a)
print *,'minimum exponent of real(8) =' , minexponent(b)
print *,'bits in integer(2) =' , bit_size(i)
print *,'bits in integer(4) =' , bit_size(j)
end program getSize
When you compile and execute the above program it produces the following result:
precision of real(4) = 6
precision of real(8) = 15
range of real(4) = 37
range of real(8) = 307
maximum exponent of real(4) = 128
maximum exponent of real(8) = 1024
minimum exponent of real(4) = -125
minimum exponent of real(8) = -1021
bits in integer(2) = 16
bits in integer(4) = 32
selected_int_kind (r)
selected_real_kind ([p, r])
The selected_real_kind function returns an integer that is the kind type parameter value necessary
for a given decimal precision p and decimal exponent range r. The decimal precision is the number of
significant digits, and the decimal exponent range specifies the smallest and largest representable
number. The range is thus from 10-r to 10+r.
For example, selected_real_kind (p = 10, r = 99) returns the kind value needed for a precision of 10
decimal places, and a range of at least 10-99 to 10+99.
Example
program getKind
implicit none
integer:: i
i = selected_real_kind (p = 10, r = 99)
print *,'selected_real_kind (p = 10, r = 99)', i
end program getKind
When you compile and execute the above program it produces the following result:
selected_real_kind (p = 10, r = 99) 8
Readability
Proper logical structure
Self-explanatory notes and comments
For example, if you make a comment like the following, it will not be of much help:
! loop from 1 to 10
do i=1,10
However, if you are calculating binomial coefficient, and need this loop for nCr then a comment like
this will be helpful:
Self-checking codes to ensure there will be no numerical errors like division by zero, square root of a
negative real number or logarithm of a negative real number.
Including codes that ensure variables do not take illegal or out of range values, i.e., input validation.
Not putting checks where it would be unnecessary and slows down the execution. For example:
real :: x
x = sin(y) + 1.0
if (x >= 0.0) then
z = sqrt(x)
end if
Clearly written code using appropriate algorithms.
Splitting the long expressions using the continuation marker ‘&’.
Making meaningful variable names.
A debugger program steps through the code and allows you to examine the values in the variables
and other data objects during execution of the program.
It loads the source code and you are supposed to run the program within the debugger. Debuggers
debug a program by:
Setting breakpoints,
Stepping through the source code,
Setting watch points
Breakpoints specify where the program should stop, specifically after a critical line of code. Program
executions after the variables are checked at a breakpoint.
Watch points are the points where the values of some variables are needed to be checked,
particularly after a read or write operation.
Comman Purpose
d
next Executes only the next line of source code, without stepping into any function
call
step Execute the next line of source code by stepping into a function in case of a
function call.
Command Purpose
next Executes only the next line of source code, without stepping into any
function call.
step Execute the next line of source code by stepping into a function in case of a
function call.
Plato is a "programming environment". Within Plato, you can create and edit programs and get them
to run. Plato's editor is special – it understands the syntax of various programming languages. We tell
Plato which language we are using when we create our empty file and save it with a .f95 (FORTRAN
95) extension. Provided you have given your file the appropriate extension, Plato's editor will be able
to check the syntax of the program, highlighting the various keywords that it knows about using a
colour code to distinguish between the various elements of the language.
Always ensure that your program files have a .f95 extension
Running your first FORTRAN 95 Program
Exercise 1.1
Type in the following exactly as shown:
!My first program
program first
print *,'This is my first program'
end program first
Click on edit
Click Select all
Click copy
Open a new document in Word or Notepad and click paste.
Making Decisions
Assignment
When we start programming, the similarity between mathematical equations and FORTRAN
statements can be confusing.
Consider the following FORTRAN statements:
x = 2 Store the value 2 in memory location x
y = 3 Store the value 3 in memory location y
z = x + y Add the values stored in memory location
x and y and store the result
in memory location z
In mathematics, “x = 2” means that the variable x is equal to 2. In FORTRAN it means “store the value 2
in the memory location that we have given the name x”.
The significance of this is made clearer by the following equation in mathematics:
x + y =z
In mathematics, this means that the left hand side of the equation is equal to the right hand side.
In FORTRAN, this expression is meaningless: there is no memory location "x+y" and so it would lead
to a compiler error.
Rule – there can only ever be ONE variable name on the left hand side of an equals sign
Exercise 2.1
Write a program which reads in two numbers a and b. Get the program to swap the values around so
that the value that was in a is now in b, and print out the result. Hint you need to declare a third
variable for intermediate storage of the data.
Arithmetic
The arithmetic operators are
+,- plus and minus
*,/ multiply and divide
** exponentiation (raise to the power)
() brackets
The order of precedence in FORTRAN is identical to that of mathematics.
Unlike algebra, the operator must always be present xy is not the same as x*y
Where operations are of equal precedence they are evaluated left to right
Consecutive exponentiations are evaluated right to left
We can override the order of evaluation by use of brackets
Exercise 2.2
The following program is an example of the use of arithmetic.
program calculate
implicit none
! a simple calculator
real :: x,y,z,answer
x=1.5
y=2.5
z=3.5
answer=x+y/z
print *,'result is ',answer
end program calculate
Explore the use of arithmetic operators by modifying program calculate. Use it to calculate the values:
Intrinsic Functions
FORTRAN is especially useful for mathematical computation because of its rich library of inbuilt
functions (intrinsic functions). We shall mention a few briefly here:
Exercise 5.3
To illustrate these limits copy file magnitude.f95 and run the program. Take a while to get a feel for
what is going on. Try inputting various values for the variable maxpower (eg 400). Can you confirm
that the table above is correct?
One interesting construct is
print *,i,2.0_ikind**i
Here, we are telling the compiler that the real constant 2.0 is also using extended precision. Check
what happens if you select extended precision (option 3) and enter a value of maxpower of 400. See
what happens if you rewrite the line to be
print *,i,2.0**i
Run the program again and enter the same values. Can you explain what is going on?
Convergence – exiting loops on a condition
In the program extended.f95, we found the sum of
It is useful to determine at what point such sums converge to a steady value – otherwise we may make
arbitrary assumptions about the summation range. Clearly, a point will be reached, within the
precision range that can be handled on our computer, that the term
will be too small to contribute to the sum. At this point we should exit the loop otherwise the program
will do more computation than is required.
One way to do this is to compare the value of the variable sum with its previous value, and if the
difference between the two is very small, then exit the loop.
program whileloop
implicit none
integer, parameter :: ikind=selected_real_kind(p=15)
real (kind=ikind) :: sum,previoussum,x,smallnumber,error
integer :: i
sum = 0.0
previoussum = 0.0
smallnumber = 10.0_ikind**(-15.0)
do i=1,1000
x=i
sum = sum + 1.0 /(x**6)
error=abs(sum-previoussum)
if (error<smallnumber) then
print *,'sum ',sum,' number of loops ',i
exit
end if
previoussum = sum
end do
end program whileloop
IMPORTANT NOTE
In the real world, we have to make choices about the amount of precision we need to work to. It is
pointless asking for 15 digits of precision if, for example, we can only take a measurement to + or – 1%
accuracy!
It is not necessary to always use a loop counter in a do loop. If we don't actually specify a counter, the
program will loop forever. Constructs like this are OK:
smallnumber = .0000001_ikind
do
print *, 'enter a positive number '
read *, number
if (number <= smallnumber) exit
end do
The disadvantage is that, if you get the code wrong, you run the risk of the program looping forever –
generally it's safer to use a loop counter!
Arrays and Formatted I/O
Arrays
Let us imagine that we want to find the average of 10 numbers. One (crude) method is shown in the
next program.
program av
real :: x1,x2,x3,x4,x5,x6,x7,x8,x9,x10,average
read *, x1,x2,x3,x4,x5,x6,x7,x8,x9,x10
average = (x1 + x2 + x3 + x4 + x5 + x6 + x7 + x8 + x9 + x10)/10
print *, 'the average is ',average
print *, 'the numbers are:'
print *, x1
print *, x2
print *, x3
print *, x4
print *, x5
print *, x6
print *, x7
print *, x8
print *, x9
print *, x10
end program av
This approach is messy, involves a lot of typing and is prone to error. Imagine if we had to deal with
thousands of numbers!
The way around this is to use arrays. An array is a list that we can access through a subscript. To
indicate to FORTRAN that we are using an array, we just specify its size when we declare it.
real, dimension(100) ::x
.
.
x(1) = 3
x(66) = 4
This snippet of code allocates 100 memory locations to the array x. To access an individual location,
called an array element, we use a subscript – here we are assigning the number 4 to the 66th element
of array x and 3 to the 1st element.
Now let's return to program av at the start of this worksheet, we'll re-write it using an array.
program av2
implicit none
real ,dimension(10) :: x
real :: average,sum
integer :: i
print *, 'enter 10 numbers'
sum=0.0
do i=1,10
read *, x(i)
sum=sum+x(i)
end do
average=sum/10
print *, 'the average is ',average
print *, 'the numbers are'
print *,x
end program av2
Notice that if we type
print*, x
the program will print out the entire contents of the array.
The additional benefit of this program is that with very few changes, we could make it deal with any
number of items in our list. We can improve on this still further by making use the parameter data
type:
program av3
!just change the value of the parameter to change the size of the !array
implicit none
integer, parameter :: imax = 10
real,dimension(imax) :: x
real :: average,sum
integer :: i
print *, 'enter’ ,imax, ’ numbers'
sum=0.0
do i=1,imax
read *, x(i)
sum=sum+x(i)
end do
average=sum/imax
print *, 'the average is ',average
print *, 'the numbers are'
print *,x
end program av3
Note this is an example of good programming. The code is easily maintainable – all we have to do to
find an average of a list of numbers of any size is just to change the size of the parameter imax. We can
also allocate the size of the array at run time by dynamically allocating memory.
The following program demonstrates the use of arrays where we do not know the size of the array.
program alloc
implicit none
integer, allocatable,dimension(:):: vector
!note syntax - dimension(:)
integer :: elements,i
print *,'enter the number of elements in the vector'
read *,elements
allocate(vector(elements))
!allocates the correct amount of memory
print *,' your vector is of size ',elements,'. Now enter each element'
do i=1,elements
read *,vector(i)
end do
print *,'This is your vector'
do i=1,elements
print *,vector(i)
end do
deallocate(vector)
!tidies up the memory
end program alloc
The program is called alloc.f95 and can be copied from the web page. Note in particular the bolded
lines. The new way of declaring the array vector tells the compiler that it is allocatable – ie the size
will be determined at run time.
We shall look at this further in Section 7.
Exercise 5.1
Write a program that asks the user how many numbers they want to enter, call this value imax.
Allocate imax elements to two arrays, a and b. Read in imax numbers to a and do the same to b. Print
out the arrays a, b and print out the sum of a and b. Compare your attempt with sumalloc.f95.
Array magic
One of the benefits of arrays is that you can easily do operations on every element by using simple
arithmetic operators.
program ramagic
implicit none
real ,dimension(10) :: a,b,c,d
open(10,file='data.txt')
read(10,*) a
b=a*10
c=b-a
d=1
print *, 'a= ',a
print *, 'b= ',b
print *, 'c= ',c
print *, 'd= ',d
end program ramagic
Exercise 5.2
Copy program ramagic.f95 and file data.txt to your own computer. Run the program and examine the
output.
Exercise 5.3
Write a program that fills a 10 element array x with values between 0 and .9 in steps of .1. Print the
values of sin(x) and cos(x) using the properties of arrays to simplify your program. Compare your
answer with ramagic2.f95.
Multi-dimensional arrays
The arrays we have looked at so far have been one dimensional, that is a single list of numbers that
are accessed using a single subscript. In concept, 1 dimensional arrays work in a similar way to
vectors. We can also use two dimensional arrays which conceptually are equivalent to matrices.
So, for example,
integer, dimension(5,5) :: a
sets up a storage space with 25 integer locations.
The next program creates a 2 dimensional array with 2 rows and 3 columns. It fills all locations in
column 1 with 1, columns 2 with 2, column 3 with 3 and so on.
program twodra
implicit none
integer,dimension(2,3) :: a
integer ::row,col,count
count = 0
!creates an array with 3 cols and 2 rows
!sets col 1 to 1, col2 to 2 and so on
do row=1,2
count=0
do col =1,3
count=count+1
a(row,col)=count
end do
end do
do row=1,2
do col =1,3
print *,a(row,col)
end do
end do
end program twodra
FORTRAN actually allows the use of arrays of up to 7 dimensions, a feature which is rarely needed. To
specify a extended precision 3 dimensional array b with subscripts ranging from 1 to 10, 1 to 20 and 1
to 30 we would write:
real (kind=ikind),dimension(10,20,30) :: b
Exercise 5.4
Using a 4*4 array create an identity matrix, that is, a matrix of the form:
1 0 0 0
0 1 0 0
0 0 1 0
0 0 0 1
and output it. Wouldn't it be nice if we could actually output the matrix elements in rows and
columns? At the end of this section we shall see exactly how to do this.
Formatting your output
You may now be wondering if there is any way to have better control over what your output looks
like. So far we have been using the default output option – that's what the *'s are for in
the write and print statements:
write(10,*) x,y,z
print *, 'program finished'
Exercise 5.5
Copy format.f95, and run it
program format
implicit none
!demonstrates use of the format statement
integer, parameter :: ikind=selected_real_kind(p=15)
real , dimension(4) :: x
integer, dimension(4) :: nums
integer :: i
real(kind=ikind),dimension(4) :: computed
!fill up the arrays with something
do i = 1,4
nums(i) = i * 10
computed(i) = cos(0.1*i)
x(i) = computed(i)
end do
print *,'nums - integer'
write(*,1) nums
1 format(2i10)
print *, 'x - real'
write(*,2) x
2 format(f6.2)
print *, 'computed - double precision'
write(*,3) computed
3 format(f20.7)
end program format
You can see that the write and format statements come in pairs.
write(output device,label) variable(s)
label format(specification)
We are using in this example a * as the output device – in other words, the screen.
The format statement can actually go anywhere in the program, but by convention we usually place
them just after the associated write or all together at the end of the program. It's just a matter of taste.
The tricky part here is the specification. There are different specifications for integer, real, and
character variables.
Integer Specification
General form : nim
Right justified
m is the number of character spaces reserved for printing (including the sign if there is one)
If the actual width is less than m, blanks are printed
n is the number of integers to output per line. If omitted, one number is output per line.
Floating point Specification
General form : nfm.d
Right justified
m is the number of character spaces reserved for printing (including the sign if there is one),
and the decimal point.
If the actual width is less than m, blanks are printed
n is the number of real numbers to output per line. If omitted, one number is output per line.
d is the number of spaces reserved for the fractional part of the number – filled with 0's if
fewer spaces are needed. If the fractional part is too wide it is rounded.
If the total width for output (m) is too small, FORTRAN will just output *'s.
Rule m >= width of the integer part plus
d plus
1 (space for decimal point) plus
1 (space for sign – if negative)
Essentially, make m nice and wide and you won't have any trouble!
Exponential Specification
General form nEm.d
Alternative specification for outputting real
d is the number of decimal places
m is the total width of the field including the sign (if any), the character E and its sign, the
decimal point and the number of places of decimals. Again make m nice and wide to ensure the
field is properly printed out.
n is the number of exponential numbers to output per line. If omitted, one number is output
per line.
Example
real :: a,b
a = sqrt(5.0)
b = -sqrt(a)
write(*,10) a,b
10 format(2E14.5)
produces:
0.22361E+01 -0.14953E+01
Character Specification
General form nAm
n is the number of strings to print
m is the maximum number of characters to output
Example:
program chars
implicit none
character ::a*10,b*10
a='hello'
b='goodbye'
write(*,10) a,b
10 format(2a10)
end program chars
Exercise 5.6
Using the format specifications in format.f95 as a guide, produce a table of
x ex
where
,
for values of x in increments of 0.1. Write your output to a file called myoutput. Ensure that your
output lines up neatly in columns. An example program is neatoutput.f95 is available on the website.
Implied Do Loop to write arrays
So far, the method we have used for input and output of arrays is:
integer :: col,row
real :: ra(10,10)
!using do loop
do row = 1,10
do col = 1,10
read *, ra(row,col)
write(*,*) ra(row,col)
end do
end do
The trouble with this method is that the rows and columns are not preserved on output. An
alternative, and neater method is to use an implied do loop in the write statement.
real :: ra(10,10)
integer :: row,col
!use implied do
do row = 1,10
do col = 1,10
read *, ra(row,col)
end do
end do
do row=1,10
write(*,10) (ra(row,col),col=1,10)
end do
10 format(10f5.1)
Subroutines and Functions
Re-using code – the subroutine
Examine the following program
program output
implicit none
real,dimension(3) :: a,b,c
character :: answer*1
!initialise arrays
a = 1.5
b = 2.5
c = 3.5
write(*,1) 'a',a
print *, 'type y to continue or any other key to finish'
read *, answer
if (answer /= 'y') stop
write(*,1) 'b',b
print *, 'type y to continue or any other key to finish'
read *, answer
if (answer /= 'y') stop
write(*,1) 'c',c
print *, 'type y to continue or any other key to finish'
read *, answer
if (answer /= 'y') stop
write(*,1) 'a*b*c',a * b * c
1 format(a,3f8.3)
end program output
The program sets up some arrays and then outputs them. At three stages in the program (bolded), it
asks whether it should continue; it stops if the answer is not 'y'. Notice that the three bolded parts of
the code are identical.
Simple enough – but look at the amount of code! Most of it is the same – wouldn't it be nice to re-use
the code and cut down on the typing? The answer is to use subroutines.
program output1
implicit none
real,dimension(3) :: a,b,c
!initialise arrays
a = 1.5
b = 2.5
c = 3.5
write(*,1) 'a',a
call prompt()
write(*,1) 'b',b
call prompt()
write(*,1) 'c',c
call prompt()
write(*,1) 'a*b*c',a * b * c
1 format(a,3f8.3)
end program output1
!++++++++++++++++++++++++++++++++++++++++++++++
subroutine prompt()
!prompts for a keypress
implicit none
character answer*1
print *, 'type y to continue or any other key to finish'
read *, answer
if (answer /= 'y') stop
end subroutine prompt
Examine the code, each time we use type
call prompt()
the program jumps to the line
subroutine prompt()
then executes each line of the code it finds in the subroutine until it reaches the line
end subroutine prompt
and then returns to the main program and carries on where it left off.
The program is much easier to understand now. All the code for prompting is in one place. If we ever
need to change the code which prompts the user to continue, we will only ever need to change it once.
This makes the program more maintainable.
Arguments to subroutines
We have seen that subroutines are very useful where we need to execute the same bit of code
repeatedly.
The subroutine can be thought of as a separate program which we can call on whenever we wish to do
a specific task. It is independent of the main program – it knows nothing about the variables used in
the main program. Also, the main program knows nothing about the variables used in the subroutine.
This can be useful – we can write a subroutine using any variable names we wish and we know that
they will not interfere with anything we have already set up in the main program.
This immediately poses a problem – what if we want the subroutine to do calculations for us that we
can use in the main program? The following program uses arguments to do just that.
Example: a program that calculates the difference in volume between 2 spheres.
program vols
!Calculates difference in volume of 2 spheres
implicit none
real :: rad1,rad2,vol1,vol2
character :: response
do
print *, 'Please enter the two radii'
read *, rad1,rad2
call volume(rad1,vol1)
call volume(rad2,vol2)
write(*,10) 'The difference in volumes is, ',abs(vol1-vol2)
10 format(a,2f10.3)
print *, 'Any more? - hit Y for yes, otherwise hit any key'
read *, response
if (response /= 'Y' .and. response /= 'y') stop
end do
end program vols
!________________________________________________
subroutine volume(rad,vol)
implicit none
real :: rad,vol,pi
!calculates the volume of a sphere
pi=4.0*atan(1.0)
vol=4./3.*pi*rad*rad*rad
!It's a little quicker in processing to do r*r*r than r**3!
end subroutine volume
When the program reaches the lines
call volume(rad1,vol1)
It jumps to the line
subroutine volume(rad,vol)
The values, rad1 and vol1 are passed to the subroutine. The subroutine calculates a value for the
volume and when the line :
end subroutine volume
is reached, the value of the volume is returned to the main program
Points to notice – these are very important – please read carefully
You may have several subroutines in your program. Ideally, a subroutine should do a specific
task – reflected by its name.
All the variables in subroutines, apart from the ones passed as arguments, are 'hidden' from
the main program. That means that you can use the same names in your subroutine as in the
main program and the values stored in each will be unaffected – unless the variable is passed
as an argument to the subroutine.
It is very easy to forget to declare variables in subroutines. Always use implicit none in your
subroutines to guard against this.
All the variables in the subroutine must be declared.
The positioning of the arguments (in this case, rad and vol) is important. The subroutine has
no knowledge of what the variables are called in the main program. It is vital that the
arguments agree both in position and type. So, if an argument to the subroutine is real in the
main program, it must also be real in the subroutine.
If an argument to the subroutine is an array, it must also be declared as an array in the
subroutine.
Exercise 6.1
Write a program that calculates the difference in area between two triangles. Your program should
prompt the user for the information it needs to do the calculation. Use a subroutine to calculate the
actual area. Pass information to the subroutine using arguments.
User Defined Functions
We have already met FORTRAN intrinsic functions like abs, cos, sqrt. We can also define our own
functions – they work in a similar way to subroutines.
As an example, let's write a program (func.f95) that does some trigonometry. As you know, the trig
routines in FORTRAN use radians, not degrees - so it would be nice to write a function that does all the
conversion for us.
print *,'Enter a number'
read *, a
pi=4.0*atan(1.0)
print *,'the sine of ',a,' is ',sin(a*pi/180)
In this snippet, we are having to code the conversion from degrees to radians directly into the main
part of the program. That's OK for a 'one-off', but what if we needed to do the conversion several
times. Now look at this:
program func
!demonstrates use of user defined functions
implicit none
integer, parameter :: ikind=selected_real_kind(p=15)
real (kind=ikind):: deg,rads
print *, 'Enter an angle in degrees'
read *, deg
write(*,10) 'sin = ',sin(rads(deg))
write(*,10) 'tan = ',tan(rads(deg))
write(*,10) 'cos = ',cos(rads(deg))
10 format(a,f10.8)
end program func
!_____________________________________________
function rads(degrees)
implicit none
integer, parameter :: ikind=selected_real_kind(p=15)
! returns radians
real (kind=ikind) :: pi,degrees,rads
pi=4.0_ikind*atan(1.0_ikind)
rads=(degrees*pi/180.0_ikind)
end function rads
What we have done, in effect, is to create our own function rads, which is used in an identical way to
the intrinsic ones you have used already like sqrt, cos, and abs.
When the line
write(*,10) 'sin = ',sin(rads(deg))
is reached, the program jumps to
function rads(degrees)
the value, degrees, is passed to the function. The function does some computation, then finally returns
the calculated value to the main program with the line
rads=(degrees*pi/180.0_ikind)
Note carefully that it doesn't return the value in the argument list (as does a subroutine) but
actually assigns the value to its own name rads.
The function rads converts the value of the argument, degrees, to radians.
Notice that we must declare the data type of the function both in the main program, and in the
function itself as if it were a variable.
Functions return one value. This value, when calculated, is assigned to the name of the
function as if it were a variable –
rads=(degrees*pi/180.0_ikind)
Exercise 6.2
Write a program that includes a function called
real function average(n,list)
where n is integer and is the number of items in the list, and list is a real array.
Write suitable code for reading the numbers from a file (or keyboard), and output the average of the
numbers.
Exercise 6.3
Write a program that allows a user to enter the size of a square matrix. In the program write a
subroutine to compute a finite difference matrix. Ensure your output is neatly formatted in rows and
columns.
So, for a 10 by 10 matrix, we expect output to look like this
2 -1 0 0 0 0 0 0 0 0
-1 2 -1 0 0 0 0 0 0 0
0 -1 2 -1 0 0 0 0 0 0
0 0 -1 2 -1 0 0 0 0 0
0 0 0 -1 2 -1 0 0 0 0
0 0 0 0 -1 2 -1 0 0 0
0 0 0 0 0 -1 2 -1 0 0
0 0 0 0 0 0 -1 2 -1 0
0 0 0 0 0 0 0 -1 2 -1
0 0 0 0 0 0 0 0 -1 2
Advanced Topics
Array Functions
FORTRAN provides a number of intrinsic functions that are useful for working with arrays. Among
these are some which are specifically aimed at working with matrices and vectors.
Matrix
multiplication of
MATMUL Matrix/vector two matrices or a
matrix and a
vector.
Scalar (dot)
DOT_PRODUCT Vector product of two
vectors
Transpose of a
TRANSPOSE Matrix
matrix
Maximum value of
an array, or of all
the elements along
MAXVAL Any array
a specified
dimension of an
array.
Minimum value of
an array, or of all
the elements along
MINVAL Any array
a specified
dimension of an
array.
Where in 2005 you bought a single CPU with a high clock rate, you now get a machine with
multiple cores. Most machines you can get these days have a minimum of 2 cores, with quad-core
machines becoming more and more common. But, there is always a but, even though you now have
access to multiple times the processing power of 2005, this does not mean that your own code will be
able to use it. Unfortunately there is no simple compiler switch which makes your code parallel (like
the -m64 switch which makes your code 64-bit), you have to do this yourself (the free lunch is over).
Two commonly used frameworks for this task are OpenMPand MPI. The former mainly focuses on
shared memory configurations (laptops, desktops, single nodes in a cluster), while the latter focuses
of large distributed memory setups (multi-node clusters) and is thus well-suited for creating codes
that need to run on hundreds or even thousands of CPU’s. The two frameworks differ significantly in
their complexity, fortunately for us, OpenMP is both the easier one, and the one most suited for a
modern multi-core computer. The OpenMP framework consists of pragma’s (or directives) which can
be added in an existing code as comment lines, and tell a compiler knowledgeable of OpenMP how
to parallelize the code. (It is interesting to note that MPI and OpenMP are inteded for parallel
programming in either C, C++ or fortran … a hint that what the important programming languages
are.)
B. Machine properties
OpenMP contains several functions which allow you to query and set several environment variables
(check out these cheat-sheets for OpenMP v3.0 and v4.0).
o omp_get_num_procs() : returns the number of processors your code sees (in hyper-threaded
CPU’s this will be double of the actual number of processor cores).
o omp_get_thread_num() : returns the index of the specific thread you are in [0..I[
Simple OpenMP functions
1. subroutine OpenMPTest1()
2. use omp_lib;
3.
14.
15. end subroutine OpenMPTest1
Notice in the example code above that outside the parallel section indicated with the directives $OMP
PARALLEL and $OMP END PARALLEL, the program only sees a single thread, while inside the
parallel section 8 threads will run (independent of the number of cores available).
C. Simple parallelization
The OpenMP frameworks consists of an set of directives which can be used to manage the
parallelization of your code (cheat-sheets for OpenMP v3.0 and v4.0). I will not describe them in detail
as there exists several very well written and full tutorials on the subject, we’ll just have a look at a
quick and easy parallelization of a big for-loop. As said, OpenMP makes use of directives (or Pragma’s)
which are placed as comments inside the code. As such they will not interfere with your code when it
is compiled as a serial code (i.e. without the -fopenmp compiler flag). The directives are preceded by
what is called a sentinel ( $OMP ). In the above example code, we already saw a first
directive: PARALLEL. Only inside blocks delimited by this directive, can your code be parallel.
OpenMP second test
1. subroutine OMPTest2()
2. use omp_lib;
3.
5. doubleprecision, allocatable :: A(:,:,:)
6. doubleprecision :: RD(1:1000)
8.
9. call random_seed()
10. call random_number(RD(1:1000))
12. allocate(A(1:IDT,1:IDT,1:IDT))
13.
15. read(*,*) NT
16. call omp_set_num_threads(NT)
17. startT=omp_get_wtime()
19. stt=omp_get_wtime()
20.
21. !$OMP DO
22. do nrz=1,IDT
23. do nry=1,IDT
24. do nrx=1,IDT
25. A(nrx,nry,nrz)=RD(modulo(nrx+nry+nrz,1000)+1)
26. end do
27. end do
28. end do
32. TTime=(omp_get_wtime()-startT)/omp_get_wtick()
34.
35. deallocate(A)
36. end subroutine RunTest2
The program above fills up a large 3D array with random values taken from a predetermined list. The
user is asked to set the number of threads (lines 14-16), and the function omp_get_wtime() is used to
obtain the number of seconds since epoch, while the function omp_get_wtick() gives the number of
seconds between ticks. These functions can be used to get some timing data for each thread, but also
for the entire program. For each thread, the starting time is stored in the variable stt. To protect this
variable of being overwritten by each separate thread, this variable is declared as private to the
thread (line 18: PRIVATE(stt) ). As a result, each thread will have it’s own private copy of the stt
variable.
The DO directive on line 21, tells the compiler that the following loop needs to be parallelized. Putting
the !$OMP DOpragma around the outer do-loop has the advantage that it minimizes the overhead
produced by the parallelization (i.e. resources required to make local copies of variables, calculating
the distribution of the workload over the different threads at the start of the loop, and combining the
results at the end of the loop).
As you can see, parallelizing a loop is rather simple. It takes only 4 additional comment lines (!$OMP
PARALLEL , !$OMP DO, !$OMP END DO and !$OMP END PARALLEL) and some time figuring out
which variables should be private for each thread, i.e. which are the variables that get updated during
each cycle of a loop. Loop counters you can even ignore as these are by default considered private. In
addition, the number of threads is set on another line giving us 5 new lines of code in total. It is of
course possible to go much further, but this is the basis of what you generally need.
Unfortunately, the presented example is not that computationally demanding, so it will be hard to see
the full effect of the parallelization. Simply increasing the array size will not resolve this as you will
quickly run out of memory. Only with more complex operations in the loop will you clearly see the
parallelization. An example of a more complex piece of code is given below (it is part of the phonon-
subroutine in HIVE):
Parallel example: Phonon spectra
2. N = this%DimDynMat
3. LWORK = 2*N - 1
4. call omp_set_num_threads(this%nthreads)
5. chunk=(this%nkz)/(this%nthreads*5)
6. chunk=max(chunk,1)
8. allocate(DM(N,N))
10. !the write statement only needs to be done by a single thread, and the other threads do
not need to wait for it
16. do nrz=1,this%nkz
17. do nry=1,this%nky
18. do nrx=1,this%nkx
19. if (this%kpointListBZ1(nrx,nry,nrz)) then
21. WORK = 0.0_R_double
22. RWORK = 0.0_R_double
23. DM(1:this%DimDynMat,1:this%DimDynMat)=this%dynmatFIpart(1:this
%DimDynMat,1:this%DimDynMat)! make a local copy
24. do nri=1,this%poscar%nrions
25. do nrj=1,this%poscar%nrions
26. Rpart=cmplx(0.0_R_double,0.0_R_double)
27. do ns=this%vilst(1,nri,nrj),this%vilst(2,nri,nrj)
28. Rpart=Rpart + exp(i*(dot_product(this%rvlst(1:3,ns),this
%kpointList(:,nrx,nry,nrz))))
29. end do
30. Rpart=Rpart/this%mult(nri,nrj)
31. DM(((nri-1)*3)+1:((nri-1)*3)+3,((nrj-1)*3)+1:((nrj-1)*3)+3) = &
32. & DM(((nri-1)*3)+1:((nri-1)*3)+3,((nrj-1)*3)+1:((nrj-1)*3)+3)*Rpart
33. end do
34. end do
35. call MatrixHermitianize(DM,IOS=IO)
37. this%FullPhonFreqList(:,nrx,nry,nrz)=sign(sqrt(abs(W)),W)*fac
38. end if
39. end do
40. end do
41. end do
In the above code, a set of equations is solved using the LAPACK eigenvalue solver ZHEEV to obtain
the energies of the phonon-modes in each point of the Brillouin zone. As the calculation of the
eigenvalue spectrum for each point is independent of all other points, this is extremely well-suited for
parallelization, so we can add !$OMP PARALLELand !$OMP END PARALLEL on lines 7 and 47. Inside
this parallel section there are several variables which are recycled for every grid point, so we will
make them PRIVATE (cf. line 7, most of them are work-arrays for the ZHEEV subroutine).
Lines 12 and 44 both contain a write-statement. Without further action, each thread will perform this
write action, and we’ll end up with multiple copies of the same line (Although this will not break your
code it will look very sloppy to any user of the code). To circumvent this problem we make use of the !
$OMP SINGLE directive. This directive makes sure only 1 thread (the first to arrive) will perform the
write action. Unfortunately, the SINGLE block will create an implicit barrier at which all other
threads will wait. To prevent this from happening, the NOWAIT clause is added at the end of the
block. In this specific case, the NOWAIT clause will have only very limited impact due to the location of
the write-statements. But this need not always to be the case.
On line 15 the !$OMP DO pragma indicates a loop will follow that should be parallelized. Again we
choose for the outer loop, as to reduce the overhead due to the parallelization procedure. We also tell
the compiler how the work should be distributed using the SCHEDULE(TYPE,CHUNK) clause. There
are three types of scheduling:
1. STATIC: which is best suited for homogeneous workloads. The loop is split in equal pieces (size
given by the optional parameter CHUNK, else equal pieces with size=total size/#threads)
2. DYNAMIC: which is better suited if the workload is not homogeneous.(in this case the central if-
clause on line 19 complicates things). CHUNK can again be used to define the sizes of the
workload blocks.
From this real life example, it is again clear that OpenMP parallelization in fortran can be very simple.
D. Speedup?
On my loyal sidekick (with hyper-threaded quad-core core i7) I was able to get following speedups for
the phonon-code (the run was limited to performing only a phonon-DOS calculation):
2. There needs to be sufficient blocks of work for all threads (indicated by nkz in the plot)
In case of the DM nored calculations, the parallelized loop clearly takes the biggest part of the
calculation-time, while for the DM red calculation, also the section generating the q-point grid takes a
large fraction of the calculation time limiting the effect of parallelization. An improvement here would
be to also parallelize the subroutine generating the grid, but that will be for future work. For now, the
expensive DM nored calculations show an acceptable speedup.
To x64 or not to x64: Installing a 64-bit fortran compiler
Current day computers generally have 64-bit processors, and most even have 64-bit operating
systems. On such systems, 32-bit programs will run fine, but 64-bit programs can make more efficient
use of the underlying system. When we installed a fortran compiler and the code::blocks IDE, the
default fortran compiler generated 32-bit programs. This generally is not an issue, unless you need a
large amount of memory, for example to store a temporary array with 4003 double precision
coordinates (as I did for a project I’m currently working on). You may first start to look for ways
of increasing the stack-size of your program, but you will soon discover that the problem is more
profound: a 32-bit program cannot access address spacing beyond 4Gb. (In practice, generally you will
not even reach 4Gb before running into problems.) This is because the memory address of your data is
stored as a 32-bit value (232 = 4 294 967 296 = 4Gb) so the only way out of this predicament is a
“larger address” aka 64-bit. So you need to install a new compiler capable of providing 64-bit
programs.
A. Installing minGW64 for code::blocks
Just unzip the folders in the location of your choice, et voilà, the compiler has been installed.
Now select Toolchain Executables, enter the path to the minGW64 installation and set the compilers
and linker to x86_64-w64-mingw32-gfortran.exe. (Don’t use gfortran.exe as this might lead to
confusion, both for yourself and the IDE.) In addition, set the ar-tool to x86_64-w64-mingw32-gcc-
ar.exe, otherwise the compiler will not be able to link libraries made with this 64-bit compiler.
Furthermore, under the tab compiler settings in the tab Other compiler options add the -m64 flag.
Finally, you may also want to update possible libraries in the linker settings tab (think for example of
libgomp.dll.a if you use openmp for parallelisation).
At this point you can compile and run your fortran code using this compiler. The resulting 64-bit
programs will now be able to access more than 2-4Gb of memory, and might turn out to run faster if
you use much double precision arithmetic.
x86_64-w64-mingw32-gfortran.exe : use the explicit name of the 64-bit compiler, not the
generic “gfortran”. This will prevent the make-program of accidentally using the wrong version
when multiple versions are installed. When they are all included in the windows PATH
variable, the make-program will just use the first one it encounters. (Use gfortran -v in the
Command Prompt to find out which version is your current default.) Not using the explicit
name may result in the following error message: “ gfortran: error: CreateProcess: No such file
or directory “ when running makelibs.bat.
2. private
3. logical :: init
4. logical :: running
5. logical :: done
6. character(len=255) :: message
7. character(len=30) :: progressString
8. character(len=20) :: bar
9. real :: progress
10. contains
11. private
12. procedure,pass(this),public :: initialize
13. procedure,pass(this),public :: reset
14. procedure,pass(this),public :: run
15. procedure,pass(this),private:: printbar
16. procedure,pass(this),private:: updateBar
All properties of the class are private (data hiding), and only 3 procedures are available to the
user: initialize, runand reset. The procedures, printbar and updatebar are private, because we
intend the class to be smart enough to decide if a new print and/or update is required. The reset
procedure is intended to reset all properties of the class. Although one might consider to make this
procedure private as well, it may be useful to allow the user to reset a progress bar in mid progress.
(The same goes for the initialize procedure.)
Run procedure of the TProgressBar class
1. subroutine run(this,pct,Ix,msg)
2. class(TProgressBar) :: this
3. real::pct
5. character(len=*),intent(in),optional :: msg
6.
9. this%running=.true.
10. this%progress=pct
11. call this%updateBar(Ix)
12. call this%printbar()
14. this%done=.true.
16. end if
17. end if
18.
In practice, the run procedure is the heart of the class, and the only procedure needed in most
applications. It takes 3 parameters: The progress (pct), the number of digits to print of pct (Ix),and
the <string> message (msg). The later two parameters are even optional, since msg may already have
been provided if the initialize procedure was called by the user. If the class was not yet initialized it
will be done at the start of the procedure. And while the progress bar has not yet reached 100%
(within 1 millionth of a %) updates and prints of the bar are performed. Using a set of Boolean
properties (init, running, done), the class keeps track of its status. The update and print procedures
just do this: update the progress bar data and print the progress bar. To print the progress bar time
and time again on the same line, we need to make use of the carriage return character (character 13 of
the ASCII table):
write(*,trim(fm), advance='NO') achar(13), trim(this%message),trim(adjustl(this
%progressString)),'%','[',trim(adjustl(this%bar))
The advance=’NO‘ option prevents the write statement to move to the next line. This can sometimes
have the unwanted side-effect that the write statement above does not appear on the screen. To force
this, we can use the fortran 2003 statement flush(OUTPUT_UNIT), where “OUTPUT_UNIT” is a
constant defined in the intrinsic fortran 2003 module iso_fortran_env. For older versions of fortran,
several compilers provided a (non standard) flush subroutine that could be called to perform the
same action. As such, we now have our class ready to be used. The only thing left to do is to turn it into
a dll or shared library.
Windows
The windows approach is very easy. We let Code::Blocks do all the hard work.
If there are modules included in the dll, the *.mod files need to be present for the compiler to access
upon compilation of the user program. (Which results in a limitation with regard to distribution
of the dll.)
Add the library to the linker settings of the program (project>build options>linker settings), and
then add the .dll file.
Upon running the program you only need the program executable and the dll.
static library
The entire setup is the same as for the shared library. This time, however, choose the “Fortran
Library” option instead of Fortran dll. As the static library is included in the executable, there is no
need to ship it with the executable, as is the case for the dll.
Unix
For the unix approach we will be working on the command line, using the intel compiler, since this
compiler is often installed at HPC infrastructures.
Special thanks also to the people of stack-exchange for clarifying some of the issues with the modules
2. private
4. contains
5. private
A Fortran class consist out of two parts: the first part lists all the properties, while the second
(optional) part, following the containskeyword, lists the methods of the object. Note that I start each
part with the private keyword, to indicate that by default everything should be hidden for the outside
world. This is part of the OOP concept of encapsulation and data hiding. In data hiding one allows the
outside world access to a property of an object through accessor and mutator methods (set and get).
This way the programmer can make sure modification of a given property is always performed in the
intended fashion. Furthermore, the use of accessors allows for properties to be internally stored in a
different fashion, or not at all, and calculated on the fly when requested. As such, the opinion variable
oi is outside the agent-object only accessable through the method getOpinion, and can only be
modified through the method setOpinion. The method updateOpinion will perform the opinion
update. Each of these methods has the public keyword, which means that they are accessible from
outside the TAgentClass-object.
There are two more important keywords at play in the above example: procedure and pass. As you
may have noted, our description of the methods does not contain any information on which
parameters they take, if they are functions or subroutines. The only information given is their name,
and that they are procedures, which just means that it are either functions or subroutines. The
advantage of this approach is that you can easily change your function/subroutine calling sequence
without the need for updating your class definition. The pass(this) keyword is a way of telling Fortran
that this procedure will take at least one variable (which doesn’t need to be provided during the
procedure call) called ‘this’(you may name it whatever you want) which is your actual TAgentClass-
object. This is a major difference with languages like C++, where the this-pointer (or self in case
of pascal) is implicitly present.
The setup of our module containing this class will then look like this:
General setup of the agent class module.
As is clear in the code above, the “procedures” are just functions and subroutines known from
standard Fortran 90/95, only the use of the class keyword is new. The class indication is similar to the
type indication, but extends on it. The concept of inheritance allows one to derive one class from
another. For example the classes TCar and TPlane could both be child classes of the class TVehicle,
they contain all methods and properties of the TVehicle class plus some specific to the TCar or TPlane
class. If a function now has a parameter of the TVehicle class, then it means that also objects of the
types TCar and TPlane are acceptable input. We will come back to this in a future tutorial when we
discuss inheritance.
A simple agent class with several properties and methods.
Returning to the paper of Sobkowicz,we find that there are other variables involved in the simulation:
Tolerance threshold d: This parameter defines if opinions update to a new value or remain
the same. In this paper, this parameter is specific for each agent, and a function of it’s opinion:
Exponent L : Another global variable which defines how d varies through the range.
Constant parameter µ : which controls the speed of convergence during the opinion updates,
and is a real value between 0 and 0.5
From these, it is clear that tolerance d is a property of our TAgentClass, and should as such be added.
All other variables are the same for all agents in the simulation, but will change for different
simulations (with the possible exception of µ). So these variable can either be defined as global
variables, which conflicts with our goal of OOP, or as properties of our TAgentClass. The later also
allows us to modify our TAgentClass, allowing for all agents to have different personal values for these
parameters. Another added bonus is the fact that these variables do not need to be given to functions
that update the opinion. Now that we have this “large” set of properties we need to think of their
initialization as well. Furthermore, since the tolerance d depends on the agent’s opinion, this also
means that it needs to be updated each time the opinion of the agent changes. Our new improved
TAgentClass looks as follows:
Extended TAgentClass
2. private
9. contains
10. private
With our TAgentClass defined (and all procedures filled out), we can now setup the actual simulation
and run it. This is done in the RunTutorial1 subroutine of the previous session. The first thing to do is
to set up our population of agents. Generally one uses an array for such a purpose.
Schematic representation of an agent based opinion dynamics program.
Following the paper, we discover that the layout of an actual simulation is a set of nested loops.
There exist two levels of time-steps. The smallest one represent a single communication between
two agents. We also learn that all communications are sequential and random, which means that on
average each agent has two communications, one as a listener, and one as a talker. Because of this
sequential nature, there is no need to have a temporary array with the previous opinion status of each
agent, since it is assumed that on a second communication event the agent will be acting with his
updated opinion. This is quite different from a mean field type approach where all communication
happens at the same instant.
The second, larger level of time-step defines a full simulation step. In this all agents have had their two
communications.
These time-steps define the two inner loops of our simulation. In addition, there is a third outer loop.
Because of the random nature of the communications, the same simulation is run multiple times, and
the evolution of the opinions over time (i.e. the large time-steps) is averaged over these runs. When
all this is finished, we just need to write down the results, and clean up the memory (e.g. deallocate
the allocatable arrays with all our agents).
Where will we be calling our TAgentClass objects?
With this in hand, we now have a fully working agents program, which can run a simulation of 20K
agents over 500 time steps in about 13 seconds on a 1GHz intel atom CPU (or 40-45 minutes for a 200
run full simulation, as presented in the paper). With 20 million opinion updates per run, this is already
entering the realm of number-crunching or scientific computing.
2. Provide a “Project title” and a location to store the files of the new project.
3. Choose your compiler (GNU Fortran Compiler, make sure this is the one selected) and Finish.
4. Now Code::Blocks will present you a new project, with already one file present: main.f90
5. First we are going to rename this main file, e.g. to main.f95, to indicate that we will be using
Fortran 95 code in this file. By right-clicking on the file name (indicated with the red arrow) you
see a list of possible op w this option while the file is open in the right-hand pane. Close the file
there, and now you can rename main.f95.
6. Although we will be following the OOP paradigm for the main internal workings of our program, a
procedural setup will be used to set up the global lines of the program. I split the code in three
levels: The top-level being the main program loop which is limited to the code below, the middle
level which contains what you may normally consider the program layout (which you can easily
merge with the top level) and the bottom level which contains the OO innards of our program. All
subroutines/functions belonging to the middle level will be placed in one external module
(Tutorial1) which is linked to the main program through the use-statement at line 2 and the
single subroutine call (RunTutorial1).
1. program AgentTutorials
2. use Tutorial1;
3. implicit none
4.
6. call RunTutorial1()
Tutorial1.f95
1. module Tutorial1
2. implicit none
3. private
4.
5. public :: RunTutorial1
6.
7. contains
8. !+++++++++++++++++++++++++++++++++++++
10. !!
12. !<------------------------------------
13. subroutine RunTutorial1()
17.
b. private : The private statement on line 3 makes sure that everything inside the module remains
hidden for the outside world, unless you explicitly make it visible. In addition
to data/information-hiding this also helps you keep your namespace in check.
c. public :: Runtutorial1 : Only the function Runtutorial1 is available outside this module.(Because
we want to use it in program AgentTutorials.)
d. contains : Below this keyword all functions and subroutines of the module will be placed.
e. “!” : The exclamation marks indicate a comment block. In this case it is formatted to be parsed
by doxygen.
f. subroutine XXX / end subroutine XXX : The subroutine we will be implementing/using as the
main loop of our agents program. Note that in recent incarnations of fortran the name of the
subroutine/function/module/… is repeated at their end statement. For small functions this may
look ridiculous, but for larger subroutines it improves readability.
At this point or program doesn’t do that much, we have only set up the framework in which we will be
working. Compiling and running the program (F9 or Build > Build and run) should show you:
In the next session of this tutorial we will use OO-Fortran 2003 to create the agent-class and create a
first small simulation
Start to Fortran
If you are used to programming in C/C++, Java or Pascal, you probably do this using an Integrated
Development Environment (IDE’s) such as Dev-Cpp/Pascal, Netbeans, Eclipse, … There are dozens of
free IDE’s for each of these languages. When starting to use fortran, you are in for a bit of a surprise.
There are some commercial IDE’s that can handle fortran (MS Visual Studio, or the Lahey IDE). Free
fortran IDEs are rather scarce and quite often are the result of the extension of a C++ focused IDE.
This, however, does not make them less useful. Code::Blocks is such an IDE. It supports
several programming and scripting languages including C and fortran, making it also suited for mixed
languages development. In addition, this IDE has been developed for Windows, Linux and Mac Os X,
making it highly portable. Furthermore, installing this IDE combined with for example the gcc
compiler can be done quickly and without much hassle, as is explained in this excellent tutorial. In
5 steps everything is installed and you are up and running:
1. Get a gfortran compiler at https://gcc.gnu.org/wiki/GFortran
Go for binaries and get the installer if you are using Windows. This will provide you with the
latest build. Be careful if you are doing this while upgrading from gfortran 4.8 to 4.9 or 4.10. The
latter two are known to have a recently fixed compiler-bug related to the automatic finalization of
objects. A solution to this problem is given in this post.
UPDATE 03/02/2017: As the gcc page has changed significantly since this post was written, I
suggest to follow the procedure described here for the installation of a 64bit version of the
compiler.
2. Get the Code::Blocks IDE
at http://www.codeblocks.org/ or http://cbfortran.sourceforge.net/ (preferred)
Since version 13.12 the Fortranproject plugin is included in the Code::Blocks installation.
3. Setup gfortran
Run the installer obtained at step 1…i.e. keep clicking OK until all is finished.
4. Setup Code::Blocks for fortran
ii. Run Code::Blocks and set your freshly installed GNU fortran compiler as default.
iii. Associate file types with Code::Blocks. If you are not using other IDE’s this may be an
interesting idea
iv. Go to settings, select “Compiler and Debugger”, click on “Toolchain executables” and set the
correct paths.
iv. Make sure the compiler is set to “GNU Fortran Compiler”, and click Finish.
v. A new project is now being created, containing a main file named “main.f90”
Of course any real project will contain many files, and when you start to create fortran 2003/2008
code you will want to use “.f2003” or “.f03” instead of “.f90” . The Code::Blocks IDE is well suited for
the former tasks, and we will return to these later. Playing with this IDE is the only way to learn about
all its options. Two really nice plugins are “Format Fortran Indent” and “Code statistics”. The first one
can be used to auto-indent your Fortran code, making it easier to find those nasty missing “end”
statements.
Fortran: Tales of the Living Dead?
Fortran, just like COBOL (not to be confused with cobold), is a programming language which is most
widely known for its predicted demise. It has been around for over half a century, making it the oldest
high level programming language. Due to this age, it is perceived by many to be a language that should
have gone the way of the dinosaur a long time ago, however, it seems to persist, despite futile
attempts to retire it. One of the main concerns is that the language is not up to date, and there are so
many nice high level languages which allow you to do things much easier and “as fast”. Before we look
into that, it is important to know that there are roughly two types of programming languages:
1. Compiled languages (e.g. Fortran, C/C++, Pascal)
The former languages result in binary code that is executed directly on the machine of the user. These
programs are generally speaking fast and efficient. Their drawback, however, is that they are not very
transferable (different hardware, e.g. 32 vs. 64 bit, tend to be a problem). In contrast, interpreted
languages are ‘compiled/interpreted’ at run-time by additional software installed on the machine of
the user (e.g. JVM: Java virtual machine), making these scripts rather slow and inefficient since they
are reinterpreted on each run. Their advantage, however, is their transferability and ease of use. Note
that Java is a bit of a borderline case; it is not entirely a full programming language like C and Fortran,
since it requires a JVM to run, however, it is also not a pure scripting language like python or bash
where the language is mainly used to glue other programs together.
Fortran has been around since the age of the dinosaurs. (via Onionesque Reality)
Now let us return to Fortran. As seen above, it is a compiled language, making it pretty fast. It was
designed for scientific number-crunching purposes (FORTRAN comes from FORmula TRANslation)
and as such it is used in many numerical high performance libraries (e.g. (sca) LAPACK and BLAS).
The fact that it appeared in 1957 does not mean nothing has happened since. Over the years the
language evolved from a procedural language, in FORTRAN II, to one supporting also
modern Object Oriented Programming (OOP) techniques, in Fortran 2003. It is true that new
techniques were introduced later than it was the case in other languages (e.g. the OOP concept), and
many existing scientific codes contain quite some “old school” FORTRAN 77. This gives the impression
that the language is rather limited compared to a modern language like C++.
So why is it still in use? Due to its age, many (numerical) libraries exist written in Fortran, and due to
its performance in a HPC environment. This last point is also a cause of much debate between Fortran
and C adepts. Which one is faster? This depends on many things. Over the years compilers have been
developed for both languages aiming at speed. In addition to the compiler, also the programming
skills of the scientist writing the code are of importance. As a result, comparative tests end up showing
very little difference in performance [link1, link2]. In the end, for scientific programming, I think the
most important aspect to consider is the fact that most scientist are not that good at programming as
they would like/think (author included), and as such, the difference between C(++) and Fortran
speeds for new projects will mainly be due to this lack of skills.
However, if you have no previous programming experience, I think Fortran may be easier and safer to
learn (you can not play with pointers as is possible with C(++) and Pascal, which is a good thing, and
you are required to define your variables, another good coding practice (Okay, you can use implicit
typing in theory, but more or less everybody will suggest against this, since it is bad coding practice)).
It is also easier to write more or less bug free code than is possible in C(++) (remember defining a
global constant PI and ending up with the integer value 3 instead of 3.1415…). Also its long standing
procedural setup keeps things a bit more simple, without the need to dive into the nitty gritty details
of OOP, where you should know that you are handling pointers (This may be news for people used to
Java, and explain some odd unexpected behavior) and getting to grips with concepts like inheritance
and polymorphism, which, to my opinion, are rather complex in C++.
In addition, Fortran allows you to grow, while retaining ‘old’ code. You can start out with simple
procedural design (Fortran 95) and move toward Object Oriented Programming (Fortran 2003)
easily. My own Fortran code is a mixture of Fortran 95 and Fortran 2003. (Note for those who think
code written using OOP is much slower than procedural programming: you should set the relevant
compiler flags, like –ipo )
In conclusion, we end up with a programming language which is fast, if not the fastest, and contains
most modern features (like OOP). Unlike some more recent languages, it has a more limited user base
since it is not that extensively used for commercial purposes, leading to a slower development of the
compilers (though these are moving along nicely, and probably will still be moving along nicely when
most of the new languages have already been forgotten again). Tracking the popularity of
programming languages is a nice pastime, which will generally show you C/C++ being one of the most
popular languages, while languages like Pascal and Fortran dangle somewhere around 20th-
40th position, and remain there over long timescales.
The fact that Fortran is considered rather obscure by proponents of newer scripting languages
like Python can lead to slightly funny comment like:“Why don’t you use Python instead of such an old
and deprecated language, it is so such easier to use and with the NumPy and SciPy library you can also
do number-crunching.”. First of all, Python is a scripting language (which in my mind unfortunately
puts it at about the same level as HTML and CSS ), but more interestingly, those libraries are merely
wrappers around C-wrappers around FORTRAN 77 libraries like MINPACK. So people suggesting to
use Python over Fortran 95/2003 code, are actually suggesting using FORTRAN 77 over more recent
Fortran standards. Such comments just put a smile on my face. With all this in mind, I hope to show
in this blog that modern Fortran can tackle all challenges of modern scientific programming.