Fortran 90 IntensiveWorkshopDr Stephen So([email protected])Griffith School of Engineering

Brief Outline Fortran versus MATLAB / CFortran program development and syntaxIntrinsic data types and operatorsArrays and array manipulationFile reading and writingSubprograms – functions and subroutinesModulesDerived data typesFunction and operating overloadinggfortran compiler optimisationCalling LAPACK and FFTPACKCalling Fortran 90 subprograms from MATLAB

Introduction Essentialelements of the Fortran 90 (90/95/2003) programming language will becovered Prior programming experience (e.g. MATLAB)is assumed Working in the Windows environment See

What do I need in order towrite Fortran programs? Fortran compilerconverts your source coded into an executableproblem Texteditor or integrated developmentenvironment (IDE) used to write your source code into a text file Somenumerical libraries (LAPACK, BLAS,FFTPACK, RKSUITE, etc.) contains subprograms written (and tested) byothers

Fortran compilers There are many commercial Fortran compilers onthe market (cost a lot of ) Intel Visual FortranPGI Visual FortranNAG FortranLahey FortranAbsoft FortranCompaq Visual Fortran (discontinued, replaced byIntel VF)We will use a free and open-source compiler GNU Fortran (gfortran)

IDEs Commercial Microsoft Visual Studio (used by Intel VF, PGIVF, Lahey, Compaq VF)Custom IDEs (used by Absoft and older versionof Lahey Fortran 95) Free IDEsand open-source IDEsCodeBlocksEclipse (Photran tran

Numerical Libraries Commercial Intel Math Kernel Library (MKL)International Mathematics and StatisticsLibrary (IMSL) Free maths librariesmaths librariesAMD Core Math Library (ACML) (no sourcecode)BLAS/ATLAS, LAPACK, FFTPACK, etc. fortran

Fortran versus MATLAB IBM Mathematical Formula Translation System(Fortran II, III, IV, 66, 77, 90, 95, 2003, 2008) Advantages: Very fast (compiled versus interpreted)Efficient with memory (can deallocate arrays)More control over data representation andformattingVery popular in scientific and engineeringcommunityLots of legacy Fortran 77 code available (BLAS,LAPACK, etc.)

Fortran versus MATLAB Disadvantages: Tedious and difficult to learnMore verboseVery limited matrix manipulationNo built-in routines or toolboxesNo built-in visualisation tools

Fortran versus C Cstrengths are in systems programming Fortran's strengths are in number crunchingand numerical data processing Fortran compilers are more sophisticatedthan C compilers for numerical optimisation Modern Fortran compilers have support forparallel processing and GPU processing(CUDA) as well

Compiling and linking Fortran90 source code is typed into a textfile ( *.f90, *.95) To run the Fortran program: Step 1: Compile source code object file(*.o)Step 2: Link in other object files (othersubprograms) and stub executable file(*.exe)

Compiling and linkingStatic func1.oWindows

Using gfortran (fromcommand-line) Ina command line, to compile and link asingle Fortran source file (hello.f90) gfortran hello.f90 -o hello.exe To gfortran hello.f90 -c To link multiple object filesgfortran hello.o func1.o func2.o -o hello.exe To perform compiling only (no linking)include static librariesgfortran hello.o liblapack.a libblas.a -o hello.exe

Fortran program structureprogram program namespecification statementsexecutable statements[contains][internal subprograms]end program program name Notethat [ ] means this section is optional

Simple Fortran programThis statement should always beprogram helloincludedimplicit none! This is a commentprint *, 'Hello, world'end program hello Savethis into a text file called hello.f90 Compile in command line gfortran hello.f90 -o hello.exe

Comments and continuation Anyline of code after ! is treated as acomment A statement can be broken into multiplelines by appending a & (in MATLAB it is .)print *, 'Length of the two sides are', &side 1, 'and', side 2, &'metres' Multiplestatements can be put onto a singleline by separating them with semi-colons (;)side 1 3; side 2 4; hypot 10

Printing to the screen Wecan see the print statement allows us todisplay text on the screen Syntax: print fmt, list fmtspecifies how the text is to be formatted(more about this later) If fmt is an asterisk (*), formatting isautomatic (list-directed I/O) list can consist of strings (text) or variables

Reading from keyboard Toread data from the keyboard, we use theread statement Syntax: read fmt, list Iffmt is an asterisk (*), the format of inputdata is automatic list consists of variables to store the data

Variables LikeMATLAB, variables are memorylocations used to store data Unlike MATLAB: we must declare every variable that we usein the program (unlike MATLAB)Variable names are not case-sensitive (e.g. pi,pI, Pi, PI are the same name) Syntax for declaring variablesdata type [, specifiers ::] variable name

Examples of variabledeclarations integer:: counter, i, j real :: mass real, parameter :: pi 3.1416 logical :: flag .false. complex :: impedance character(20) :: name

Intrinsic numeric types Integer A whole number that has no fractional part ordecimal point e.g. 3, -44, 120Largest: 2147483647 Real (4 bytes / 32 bits)(4 bytes / 32 bits)A number that can have a fractional part or decimalpoint e.g. 3.0, 2.34, -98.12, 3.3e8IEEE 754 single precision: 6 decimal digitsSmallest, Largest: 1.17549435E-38, 3.40282347E 38

Intrinsic numeric types Doubleprecision (8 bytes / 64 bits)A real number with twice the precision e.g. 3.3d8 IEEE 754 double precision: 15 decimal digits Smallest: 2.22507385850720138E-308 Largest: 1.79769313486231571E 308Complex Complexnumber consisting of a real part andimaginary part e.g. 3.2 – 4.98 i (3.2, 4.98)

Intrinsic non-numeric types Character ASCIIcharacter or text e.g. 'a', 't', ' ' character(N) to declare string of N chars Can use either “ or ' to delimit strings // for concatenation e.g. 'abc'//'def' len() function returns the number ofcharacters Logical Either.true. or .false.

Kind parameter Therange of an integer (max 2147483647)may not be enough The precision of double precision (15) maynot be enough Fortran 90 allows us to specify the kind ofinteger and the kind of real The kind number is usually (but not always)the number of bytes integer(kind 8) or integer(8) – 8 byte integerreal(8) – 8 byte real (double precision)

Kinds for different compilers Compilerspecific (refer to kinds.f90 andkindfind.f90) gfortran supports: integer(1), integer(2), integer(4), integer(8),integer(16) (64 bit only)real(4), real(8) (double precision),real(10) – precision of 18 (extended precision) LaheyFortran 95 and Intel Visual Fortransupports: integer(1), integer(2), integer(4), integer(8)real(4), real(8) (double precision),real(16) – precision of 33 (quad precision)

Intrinsic kind functions kind(a)will return the kind of variable a Use the following functions for portability(since the kind number may not be equal tothe number of bytes): selected int kind(p) returns the kind ofinteger with p significant digitsselected real kind(p, r) returns the kindnumber of a real with p digit precision anddecimal exponent range of -r and r Thesefunctions return -1 if not supported bythe compiler, generating a compile error

Specifying the kind of literalconstant Inorder to specify the kind for a literalconstant, we use an underscore followed bythe kind number 231 2 – integer(2) literal constant 23.13 8 – real(8) literal constant More convenient to use parametersinteger, parameter :: sp selected real kind(6)integer, parameter :: dp selected real kind(15) Thento specify a double precision literal, wetype: 3.1416 dp

Kind number exampleprogram kind exampleimplicit noneinteger, parameter :: dp selected real kind(15)integer, parameter :: ep selected real kind(18)real :: adouble precision :: b, creal(ep) :: da 1.234567890123456789012345678901234567879b 1.234567890123456789012345678901234567879c 1.234567890123456789012345678901234567879 dpd 1.234567890123456789012345678901234567879 epprint *, 'real(4) ', aprint *, 'real(4) ', bprint *, 'real(8) ', cprint *, 'real(10) ', dend program kind example

Intrinsic functions forcomplex variables real(z)– returns the real part of z aimag(z) – returns the imaginary part of z conjg(z) – return complex conjugate of z abs(z) – returns magnitude of z z cmplx(a, b) – form complex numberfrom two real variables a and b

Initialising variables Thereare three ways to initialise variableswith starting values: During the declaration:real :: mass1 50.3, mass2 100.0 In a data statement (version 1):real :: mass1, mass2,data mass1, mass2 /50.3, 100.0/ In a data statement (version 2):real :: mass1, mass2,data mass1 /50.3/, mass2 /100.0/

Repeating initialising data Ifa list of variables need to be initialised withthe same value, we can repeat them Instead of:integer :: a, b, c, ddata a, b, c, d /0, 0, 0, 0/ We can use:integer :: a, b, c, ddata a, b, c, d /4 * 0/ Works for initialising arrays too

Arithmetic operators Addition Subtraction Multiplication * Division / Exponentiation **312 3 ** 12 For integer powers, it is faster to multiplymanually e.g.

Relational operators Greaterthan .gt. or Greater than or equal .ge. or Less than .lt. or Less than or equal .le. or Equal .eq. or Not equal .ne. or /

Logical operators Not(unary) .not. And .and. Or .or. Logical equivalence .eqv. Logical non-equivalence .neqv.

Simple example programprogram age programimplicit noneinteger :: year, age, present year 2011character(10) :: my name ! string of 10 charsprint *, 'What is your name?'read *, my nameprint *, 'What year were you born?'read *, yearage present year – yearprint *, 'Hello', my name, 'you are of age', ageend program age program

Data type conversions Needto be careful when mixing differentdata types in expressions Fortran will sometimes automatically convertthe variables to preserve precision (but stillneed to check)integer :: a 3, c 7real :: b 4.3, dc a b ! truncation will occurd a b ! no truncationprint *, c, d

Explicit data type conversion Wecan explicitly force data type conversions,just to be certain int(x) converts to integer real(x) convert to real dble(x), dfloat(x) convert to doubleprecision cmplx(x) converts real/integer to complex(no imaginary part) cmplx(x, y) converts real/integer tocomplex (real x and imaginary y)

Exampleinteger :: a 3, c 7real :: b 4.3, dcomplex :: ec a bd real(a) be cmplx(b, c) ! not e (b, c)print *, c, d, e

Do loop Syntax:do [var start, end, step]statementsend do varmust be an integer Example:do i 1, 10, 2total total iend do

Do loop (no counter) Infinitedo (no var) can be stoppedwith exitdoif (condition) thenexitend ifend do

Implied Do loop Can be used for:Initialising arraysReading a list of values (useful for 2D arrays)Printing a list of values Syntax: (statement, var start, stop, step) print*, (i, i 0, 20, 2) read *, (table(i), i 1, 10) real :: values(20)values (/ (0.2 * i, i 1, 20) /)

Reading with implied do Comparethis code:do i 1, 10read *, value(i)end do With this code:read *, (value(i), i 1, 10) Is there a difference?

While loop Introducedin Fortran 90 (some Fortran 77compilers supported it) Syntax:do while (expression)statementsend do Example:do while (total .lt. 100)total total 1end do

If then else. Example:if (discrim .gt. 0) thenprint *, 'There are two roots'else if (discrim .eq. 0) thenprint *, 'There is a single root'elseprint *, 'The roots are complex'end if

Select statementselect case (option)case(1)print *, 'Gardening'case(2:10)print *, 'Clothes'case(11, 12, 15, 20:30)print *, 'Restaurant'case defaultprint *, 'Carpark'end select

Arrays Arraysare variables that consist of multipleelements of the same type Fortran supports up to 7 dimensions By default, array indices start at 1 Examples:integer :: table1(20), vector(5, 5)complex, dimension(3, 3) :: matrix1, matrix2table1(1) 3vector(2, 3) -10matrix1(1, 2) (-2.4, 9.12)

Custom array indices andinitialisation We can redefine the index rangeinteger :: number(0 : 10), cost(-10 : 20, 3 : 9) Wecan initialise (and set) 1D arrays usingarray literal constants (/ . /)real :: values(3) (/3.2, 4.7, -1.0/)integer :: table(2)table (/2, 4/) Notethat the number of elements (on theright) must match the array dimensions

Array initialisation using adata statement Wecan also initialise arrays using a datastatement (from Fortran 77)real :: readings(5), count(3), voltage(8)data readings /1.2, 0.9, -2.0, 56.0, -0.89/, &count /2.0, 6.0, 7.0/, voltage /8 * 0/ 2Darrays are internally stored as 1D arraysinteger :: matrix(2, 2)data matrix /2, 5, 3, 10/ Youcan also reshape() function (see later)

How multidimensional arraysare stored UnlikeC, Fortran stores multidimensionalarrays using column-order indices to the left change fastest Imagine a 2D array: integer :: a(2, 3)a(1,1)a(1,2)a(1,3)a(2,1)a(2,2)a(2,3) InFortran, this array is stored in this order:[a(1,1), a(2,1), a(1,2), a(2,2), a(1,3), a(2,3)]

Array slicing LikeMATLAB, we can select sub-sections orslices of an array Syntax: array name(start:end:stride)integer :: a(10), b(5), totaldata a /1, 2, 3, 4, 5, 6, 7, 8, 9, 10/b a(3:7) ! a(3),a(4),a(5),a(6),a(7)b a(:5)! a(1),a(2),a(3),a(4),a(5)b a(6:)! a(6),a(7),a(8),a(9),a(10)b a(1:10:2) ! a(1),a(3),a(5),a(7),a(9)b a((/ 2, 5, 7, 8, 9 /)) ! a(2),a(5),a(7),a(8),a(9)total sum(a(2:7)) ! Sum is an intrinsic function

Array operators Notethat the arrays must conform with eachother Addition and subtraction integer :: a(10), b(10), c(10)c a badds each corresponding element Multiplication and divisionc a*bmultiplies each corresponding element Exponentiation c a ** 4.4Raises every element of a to power of 4.4

Where.elsewhere.endwhere Allowsus to perform an operation only oncertain elements of an array (based on amasking condition)e.g. where (a 0.0) a 1.0 / a This performs a reciprocal only on theelements of array a that are positive

Where.elsewhere.endwhere Moregeneral form:where (logical expression)array operationselsewherearray operationsend where Forexample:where (a 0.0)a log(a) ! log of positive elementselsewherea 0.0 ! set negative elements to zeroend where

Intrinsic vector/matrixfunctions Transpose a transpose(b) Dot productc dot product(a, b) Matrix multiplication (matrices must conform)c matmul(a, b)

Matrix example Giventwo matrices A and B[ ]A 1 23 4B [ ]5 67 8 Wewant to write a Fortran program thatcomputes and prints the following:C ABAT

Array terminology Thenumber of dimensions is called the rank The number of elements along a singledimension is called the extent in thatdimension Sequence of extents is the shape of thearray Example: integer :: a(3, 4, 10:15) Rank is 3Extent of 1st dim 3, 2nd dim 4, 3rd dim 6Shape is (/ 3, 4, 6 /)

Inquiring and reshapingarrays s shape(a) returns the shape of array a n size(a) returns the total number ofelements in array a n size(a, i) returns the total number ofelements along dimension i of array a c reshape(b, s) function to change shapeof array b to shape of a (stored as s) e.g. c reshape(b, (/3, 4/))

Example array programprogram reshape1implicit noneinteger :: a(4), b(2, 2), c(2, 2), i, j, d(3, 5)data a /1, 2, 3, 4/, b /1, 2, 3, 4/print *, 'Shape of a ', shape(a)print *, 'Shape of b ', shape(b)print *, 'size of dim 2 ', size(d, 2)c reshape(a, (/2, 2/))print *, 'Matrix b'do i 1, 2print *, (b(i, j), j 1, 2)end doprint *, 'Matrix c'do i 1, 2print *, (c(i, j), j 1, 2)end doend program reshape1


Allocatable arrays Normal arrays need to have their size fixed andknown at compile-timeUse allocatable arrays if their size is not knownreal, allocatable :: array(:)In the code, to allocate this array to n elementsallocate(array(n))When not needed, we can free the memorydeallocate(array)You can allocate multidimensional arraysreal, allocatable :: array(:, :, :)

Formatting screen output Wehave used print *, var1, var2, var3 The asterisk * means use list-directed output(compiler automatically formats thevariables) For more fine formatting controls, there aretwo methods: Including a format specifier stringUsing the format statement

Formatting screen output Format specifier string: print '(1x, i5, 5x, f10.8)', index, densityFormat statement:print 10, index, density10 format(1x, i5, 5x, f10.8) These examples both use the followingformatting for each line: One space (1x), thenInteger number with width of 5 places (i5), thenFive spaces (5x), thenReal number with width of 10 places and 8 decimalplaces (f10.8)

Formatting specifiers Nx– N blank spaces iW – integer number with W places fW.D – real number with W total places(includes decimal point) and D decimalplaces eW.DeE – real number in exponential formatwith W places, D decimal places, Eexponential places (e.g. 0.xxxezz) enW.D – engineering notation (e3, e6, e9,.) esW.D – scientific notation