Transcript An Introduction to FORTRAN Lecturer: Dr Dan Moore

An Introduction to FORTRAN
Lecturer: Dr Dan Moore
(50 years of structureless programming!)
Lecture #3
Arithmetic IF vs. Logical IF (revisited)
FORTRAN II has the arithmetic IF statement:
IF(number)negative label, zero label, positive label
This leads to very unstructured programming:
IF(D)20,30,40
20 DR=sqrt(-D)
PRINT *,'Complex Roots'
GO TO 10
30 PRINT *,'Repeated Roots, R1=R2=',-B/(2*A)
GO TO 10
40 DR=sqrt(D)
PRINT *,'Roots = ‘,(-B+DR)/(2*A), (-B-DR)/(2*A)
GO TO 10
Note the many GO TO statements. It is not clear when one
code block ends and another code block begins.
You must jump out of each separate code section explicitly.
17/07/2015
Fortran 2014 Lecture #3
2
(better is the) Logical IF
• IF(logical expression)THEN
statement (do if logical expression is true)
statement
. . .
ELSE
statement (do if logical expressionis false)
statement
. . .
ENDIF
This is much more Structured Programming as you can see that the
THEN, ELSE and ENDIF statements that top and tail each separate
code block.
Advice: Use the arithmetic IF statement if and only if its 3 way switch
leads to more ‘transparent’ (easy to understand) code.
17/07/2015
Fortran 2014 Lecture #3
3
Indexed DO loops
• Formal Syntax for the Start:
DO [label] int_var=start,stop[,stride]
Formal syntax for the end of the loop: Classic FORTRAN
label CONTINUE or
label executable statement
Formal syntax for the end of the loop: Modern FORTRAN
END DO
•
The number of iterations, which is evaluated before execution of the loop begins,
is calculated as: MAX(INT((< stop>-< start >+< stride>)/< stride>),0)
If this is zero or negative then the loop is not executed.
If < stride > is absent it is assumed to be equal to 1.
•
You can exit a DO loop by simply jumping out of it: IF(err.lt.1.e-10)EXIT
will jump to the next statement after the END DO statement.
•
YOU MAY NOT CHANGE THE DO LOOP VARIABLE INSIDE THE DO LOOP!
17/07/2015
Fortran 2014 Lecture #3
4
Nested DO loops (loops inside loops)
Classic FORTRAN Nested Loop Structure (Matrix Multiply!):
DO 10 I = 1,L
DO 10 J = 1,M
C(I,J) = A(I,1)*B(1,J)
DO 10 K=2,N
10
C(I,J)=C(I,J)+A(I,K)*B(K,J)
(Several DO loops may terminate on the same labelled Fortran statement!
Do this to annoy the Computer Scientists!)
Modern FORTRAN recommended Nested Loop Structure:
DO I = 1,N
DO J = 1,N
C(I,J) = A(I,1)*B(1,J)
DO K=2,N
C(I,J)=C(I,J)+A(I,K)*B(K,J)
END DO
END DO
END DO
(Each DO loop has its own END DO statement. Can make debugging easier!)
17/07/2015
Fortran 2014 Lecture #3
5
Flow Control: Classic vs. Modern
•
Classic FORTRAN DO loops
DO 10 I=1,101,5
X(I)=Y(I)
DO 10 J=2,5
10 X(I)=X(I)+Y(J)
DO I=1,101,5
X(I)=Y(I)
DO J=2,5
X(I)=X(I)+Y(J)
END DO
END DO
Computed GO TO Statement
GO TO (200,300,400),NTRA
(Code for NTRA <1 or >4)
GO TO 500
200 (Code for NTRA =1) )
...
GO TO 500
300 (Code for NTRA =2 )
...
GO TO 500
400 (Code for NTRA =3 )
...
500 CONTINUE
(the indentation is for appearances only!)
Note that you must jump out of each code block!
Computer scientist blame
ALL OF THE EVILS OF PROGRAMMING
on the GO TO statement !!
All of the languages they develop ban the
GO TO
Statement.
17/07/2015
Modern Fortran – DO loops
Case Statement
SELECT CASE (NTRA)
CASE(l)
(Code for NTRA =l)
CASE(2)
(Code for NTRA=2)
CASE(3,4,5)
(Code for NTRA=3,4,5)
CASE(DEFAULT)
(Code for NTRA=anything else)
END SELECT
Note that you can execute only one code block!
(Although you can have multiple case values for
one code block.)
Fortran 2014 Lecture #3
6
Modern Flow Control:
• DO WHILE (logical expression)
. . .
END DO
or
DO 200 WHILE (logical expression)
. . .
200 CONTINUE
The (lable) CONTINUE statement can be used to terminate any
DO (lable) WHILE type loop block.
The code inside of the loop must change at least one of the variables in the
(logical expression) to cause the loop to terminate.
Otherwise, you have written an INFINITE LOOP!
(And your program will never stop running!)
17/07/2015
Fortran 2014 Lecture #3
7
Escaping from a Loop
• The EXIT command will cause the program to jump to
the first statement after the end of the current DO loop
(END DO) This is the ‘escape’ command.
• The GO TO label will ‘escape’ from the current
loop iff label is a statement outside of the current
loop. You can jump to anywhere in the current program
block, even to a label before the GO TO label
statement.
(The backwards jump, again anathema to
our Computer Scientist Colleagues!)
17/07/2015
Fortran 2014 Lecture #3
8
More Loop Control:
• The CYCLE command interrupts the execution of the
DO code block and begins a new cycle of the
DO loop with the next value of the DO loop control
(goes back to the DO statement!).
I=0
DO
I = I+1
IF(I >= 3 .AND. I <= 8)CYCLE
IF(I > 10) EXIT
PRINT *, ‘ I = ‘,I
END DO
Prints:
I = 1
I = 2
I = 9
I = 10
17/07/2015
Fortran 2014 Lecture #3
9
New Chapter: FORTRAN Arrays:
• An Array (or matrix) holds a collection of different values of
the same data type as a single named variable, accessible by
an address.
• The rules for array names are the same as for other
FORTRAN Variables, including any IMPLICIT and explicit
typing already declared in the program block.
• Arrays may be indexed with one, two, three up to seven
indices. These are referred to as the Dimensions of the array.
(My current record is 5 indices for a 3-D velocity field with
two time levels: VEL(3,NX,NY,NZ,2)
• All arrays must be declared before they are used!
17/07/2015
Fortran 2014 Lecture #3
10
Array Terminology
• Rank – Number of indices or dimensions
• Bounds – permitted lower and upper limits for
each index
• Extent – total number of possible index values
in each dimension
• Size – total number of elements
• Shape – rank and extents
• Conformable – two arrays are said the be
conformable when they are the same shape.
17/07/2015
Fortran 2014 Lecture #3
11
Array Declaration
Arrays are declared at the start of each program block and are of a fixed size!
(This is why FORTRAN Code is FAST! Static Memory Allocation)
Classic:
DIMENSION A(15),B(5,3),C(3,5)
REAL*8 D(-100:100),E(201)
REAL :: A(15),B(5,3),C(3,5)
REAL, DIMENSION(10) :: G,H,P
REAL*8, DIMENSION(-100:100) :: A
REAL*8, DIMENSION(:,:),ALLOCATABLE :: PSI
The last statement declares PSI(:,:) to be a variable size array!
The ultimate size of Array PSI(:,:) is to be determined by a subsequent
ALLOCATE statement such as:
Modern:
READ *,isz1,isz2
ALLOCATE(PSI(isz1,isz2),STAT=ierr)
IF(ierr.NE.0)THEN
PRINT *,’PSI allocation fails!’
STOP
ENDIF
We know PSI is a Two Dimension array but we don’t know what size when we write the code.
17/07/2015
Fortran 2014 Lecture #3
12
Array Declaration - II
• The array indices can begin (lower bound) and
end (upper bound) with any literal integer or integer
parameters, but the lower bound must be less than
or equal to the upper bound.
REAL*8 D(-100:100)
• If no lower bound is specified, it is assumed to be 1.
DIMENSION A(15)
DIMENSION A(1:15) are equivalent!
• Arrays can be zero-sized (but then they can’t be used!)
17/07/2015
Fortran 2014 Lecture #3
13
Accessing array elements
• Individual array elements are accessed by integer
index values:
DO I=1,5
B(I)= ATAN(DBLE(I)/2)
C(2*I-1)=COS(B(I))
END DO
• Whole arrays may be ‘assigned’ to a literal
or to a single variable value:
A=0
• A range of array elements may be assigned in a
similar manner: A(5:10)=1
17/07/2015
Fortran 2014 Lecture #3
14
Accessing array elements - II
• Referencing an array element with an address outside of
the declared range for that array has totally
unpredictable results!
IT IS A VERY BAD IDEA !!!
If your compiler is in ‘Debug’ mode,
Array Bound Checking is generally switched on and an
error message and a program stop generally results.
In ‘Release’ mode Array Bound Checking is generally
switched off and serious problems can arise.
17/07/2015
Fortran 2014 Lecture #3
15
Array Operations
• Whole array names may appear on the Right Hand side of an
assignment statement if their use is appropriate:
DIMENSION A(5,3),B(5,3),C(5,3)
…
A = B + C
is allowed.
(But so is A = B * C and
A = B * sin(C) !! What nonsense!!)
the compiler translates this last statement as
A(1,1) = B(1,1)*sin(C(1,1))
A(2,1) = B(2,1)*sin(C(2,1))
etc…
B * C is
NOT MATRIX MULTIPLICATION AS WE KNOW IT!
17/07/2015
Fortran 2014 Lecture #3
16
Array Operations - II
• Element based array operations are the safest:
• Array subsections can be used on both sides of =:
DIMENSION A(15)
…
A(:)=1.
! set whole array to 1.
A(3:6)=2.
! set A(3),A(4),..A(6) to 2.
A(1:7:2)=3. ! set A(1),A(3),..A(7) to 3.
A(1:11:2)=B(1:6) ! set A(1),A(3),…A(11) to B(1),B(2), …B(6)
<start>:<stop>[:<stride>] is the format and all three
must be integer expressions. If <start> is omitted, the
index begins at the lowest possible bound. If <stop> is
omitted the index ranges to the highest possible bound.
17/07/2015
Fortran 2014 Lecture #3
17
How Arrays are stored in memory
• FORTRAN doesn’t know or care !!!
(Formally, there is no guaranteed storage association.)
• However, for the purposes of I/O, the elements of an
array are read in or printed out with the first index
varying the fastest, then the next index and then the
next index:
DIMENSION A(5,3)
PRINT *,’A = ‘, A
will print out the array in the order:
A(1,1), A(2,1),A(3,1),A(4,1),A(5,1), A(1,2), A(2,2),A(3,2),
… A(5,3)
READ *, A
expects the data for A to be entered (typed in) in that order.
17/07/2015
Fortran 2014 Lecture #3
18
Array initialization (Constructors)
• Initialization of a 1 D array:
REAL, DIMENSION(4) :: heights
heights=(/5.10, 5.6, 4.0, 3.6/)
All array initializations must be 1-D
• Initialization of a 2 D array: (is very awkward!)
REAL, DIMENSION(10,10) :: T
T=RESHAPE([(i,i=1,100)],[10,10])
The initialization can be done at declaration:
REAL,DIMENSION(4)::heights=(/5.10,5.6,4.0,3.6/)
• Named Array Parameters can be created:
REAL, DIMENSION(3,3), PARAMETER :: Ident_3 =
RESHAPE([1,0,0,0,1,0,0,0,1],[3,3])
17/07/2015
Fortran 2014 Lecture #3
19
&
Arrays as Arguments in SUBROUTINEs and FUNCTIONs
• Array names are passed by NAME in SUBROUTINEs:
The sizes of the arrays can be dummy arguments.
•
REAL, DIMENSION(:,:) :: A(5,4),B(4,5),C(5,5),D(4,4)
DO 10 I=1,5
DO 10 J=1,4
A(I,J)=4*I+J
10
B(J,I)=5*J+I
CALL myMATMUL(A,B,C,5,4,5)
CALL myMATMUL(B,A,D,4,5,4)
10
SUBROUTINE myMATMUL(A,B,C,L,M,N)
DIMENSION A(L,M),B(M,N),C(L,N)
DO 10 I = 1,L
DO 10 J = 1,M
C(I,J) = A(I,1)*B(1,J)
DO 10 K=2,N
C(I,J)=C(I,J)+A(I,K)*B(K,J)
RETURN
END
17/07/2015
Fortran 2014 Lecture #3
20
Automatic Arrays
• In SUBROUTINE and FUNCTION program blocks, it is
possible to declare local arrays whose size depends on the
Dummy Arguments.
Such arrays are called Automatic Arrays.
They cannot have the SAVE attribute or be initialized:
SUBROUTINE MAGIC(A,M,N)
DIMENSION A(M,N)
DIMENSION ROWSUM(M),COLSUM(N)
DO I=1,M
ROWSUM(M)=A(M,1)
DO J=2,N
ROWSUM(M)=ROWSUM(M)+A(M,N)
END DO
END DO
. . .
ROWSUM and COLSUM are Automatic Arrays
17/07/2015
Fortran 2014 Lecture #3
21
Classic Fortran Bad Habits
• Memory was extraordinarily expensive and limited before
1995! And Old FORTRAN did not allow dynamic arrays.
• So some ‘old’ FORTRAN Programs did ‘Bad’ things such as:
i. Creating large 1-D arrays in the Main Program
block (e.g. DIMENSION A(1000) )
and then passing an array element as an argument
to a subroutine expecting a 2-D array as an argument!
CALL KE(A(1),A(257),16,16,Energy)
SUBROUTINE KE(Psi,Omega,NX,NY,E)
DIMENSION Psi(NX,NY),Omega(NX,NY)
This was Legal and commonplace 30 years ago!
Modern FORTRAN Compilers should complain about this
but generally can be persuaded to compile it as expected.
17/07/2015
Fortran 2014 Lecture #3
22
Classic Fortran Bad Habits - II
• There is no formal requirement in Classic FORTRAN that the
number of DIMENSIONs of arrays passed as SUBROUTINE or
FUNCTION arguments agree between program blocks!.
An array can be 1-D in the main program block,
5-D in a SUBROUTINE and then 3-D or 2-D in a FUNCTION
called by the SUBROUTINE.
• If a FUNCTION or a SUBROUTINE expects an array as a dummy
argument, the calling program can pass any array address as that
argument. The FUNCTION/SUBROUTINE will expect the ‘local’
array to start at that location and that the array supplied as an
argument will be declared big enough to keep any address
references with the array bounds:
EXECUTIVE SUMMARY! Turn on Array Bounds Checking!
17/07/2015
Fortran 2014 Lecture #3
23
Allocatable Arrays
• All arrays used in Classic FORTRAN had to be declared explicitly in the Main
Program Block or before they were used in a FUNCTION or SUBROUTINE
This is referred to as a STATIC Memory Model.
Modern Fortran introduced DYNAMIC arrays whose size can be determined when
the program is run. These are ALLOCATABLE arrays.
Their Dimension (number of indices) is declared before use, but not their size:
REAL*8, DIMENSION(:,:),ALLOCATABLE :: PSI
The ultimate size of array PSI is determined by an ALLOCATE statement:
READ *,isz1,isz2
ALLOCATE(PSI(isz1,isz2),STAT=ierr)
IF(ierr.NE.0)THEN
PRINT *,’PSI allocation fails!’
STOP
ENDIF
• The IF(ierr.NE.0)THEN loop is good coding practice to warn you if the
attempt to allocate memory for the array of the required size fails!
17/07/2015
Fortran 2014 Lecture #3
24
Allocatable Arrays - II
• Allocated arrays can be de-allocated to reuse the Computer
memory for other things:
IF(ALLOCATED(PSI))DEALLOCATE(PSI)
• It is an error to attempt to de-allocated an array that has not
been declared to be ALLOCATABLE and ALLOCATED !
• The built in function ALLOCATED(array)
returns a value .TRUE. if the array can be DEALLOCATEd
• If a program block containing an allocatable array which has
not been declared with a SAVE is exited without the array
being DEALLOCATEd then this memory becomes inaccessible
and you have a ‘Memory Leak’.
17/07/2015
Fortran 2014 Lecture #3
25
Character Declarations and access
• CHARACTER strings are declared in a similar manner
to 1-D arrays:
• CHARACTER(LEN=10)
:: StrNum,Postcode
CHARACTER
:: sex
CHARACTER (LEN=32)
:: address,street
CHARACTER (LEN=80),DIMENSION(60)
:: PAGE
• Individual characters and substrings in a character
string can be accessed with the usual subarray notation:
address=‘180 Queens Gate, SW7 2BZ’
StrNum=address(1:3)
Postcode=address(19:25)
• Two character strings are joined together to make a
single string with the // operator
17/07/2015
Fortran 2014 Lecture #3
26
Derived Types
• Modern FORTRAN allows compound entities or derived types to be defined:
TYPE POINT
REAL :: x,y,z
END TYPE POINT
• TYPE(POINT)
:: p1,p2,p3 p4
TYPE TRIANGLE
TYPE(POINT) :: p1,p2,p3
END TYPE TRIANGLE
Values can be assigned by component:
p1%x=0.0;p1%y=0.0;p1%z=0.0
p2%x=1.0;p2%y=0.0;p2%z=0.0
p3%x=0.0;p1%y=1.0;p1%z=0.0
Values can be assigned as an object:
p4=POINT(0.0,0.0,1.0)
t1=TRIANGLE(p1,p2,p3)
T3D1=TETRA(p1,p2,p3,p4)
• TYPE(TRIANGLE):: t1,t2,t3,t4
Derived types can be arguments to
FUNCTIONs and subroutines.
TYPE TETRA
TYPE(POINT) :: p1,p2,p3,p4
END TYPE TETRA
Derived types can be returned by
FUNCTIONs
• TYPE(TETRA) :: T3D1
17/07/2015
Fortran 2014 Lecture #3
27
Derived Types-II
•
Derived Type I/O :
PRINT *, p1
is the same as
PRINT *, p1%x, p1%y, p1%z
Derived Type definitions should be put in a MODULE program block.
Derived Types may be arguments to FUNCTIONs and may be the returned value of a FUNCTION.
These FUNCTIONs should be put into the MODULE as well.
REAL FUNCTION TETRAvolume(T)
USE TETRAdef
TYPE TETRA :: T
DIMENSION V12(3),V13(3),V14(3)
V12(1)=T%p2%x-T%p1%x
V12(2)=T%p2%y-T%p1%
…
TETRAvolume=ABS(DOT(V12,CROSS(V13,V14)))/6.
RETURN
END
REAL FUNCTION DOT(V1,V2)
DIMENSION V1(3),V2(3)
DOT=V1(1)*V2(1)+V1(2)*V2(2)+V1(3)*V2(3)
RETURN
END
17/07/2015
Fortran 2014 Lecture #3
28
Interactive Developer Environment
• Are you content with the current
a. Compose/Edit
b. Compile/Link
c. Run
cycle to create a working FORTRAN Program?
If(YES)THEN  Start Practical Work!
So you think this edit/compile/run cycle can be done
more easily? Other people have thought so too.
17/07/2015
Fortran 2014 Lecture #3
29
Fortran Compilers and IDEs
• What is an IDE ?
• An Interactive Development Environment is a program to let you write,
edit, compile, link and run a program in a specific computer language.
It was invented by Borland in the mid 1980's to support their PASCAL
compiler. Several IDEs for FORTRAN are currently in widespread use.
• The most common IDE in Windows is Microsoft Visual Studio, currently in
version 2013. Intel and SilverFrost FORTRAN can integrate with this IDE.
• The Netbeans and Eclipse IDEs are also available in LINUX
(and Apple MAC) as well as Windows. These work with GNU gfortran
• On a Mac, there is a program called Xcode. Many commercial FORTRAN
packages are sold with their own proprietary IDE. See the course Web site
for the details on how to download and install various free IDEs and the
Fortran Compilers to go with them.
17/07/2015
Fortran 2014 Lecture #1
30/137
Your choices on these ICT PCs:
• You can choose to use:
– (i) MS VS 2013 + Intel Fortran Compiler
or
– (ii) Netbeans + MinGW GNU gfortran
to write and run your FORTRAN. However, neither of
these are likely to be available within your own
department unless they have specifically requested
them from ICT.
17/07/2015
Fortran 2014 Lecture #1
31/137
Fortran on your own IT:
• Students may download the Silverfrost FORTRAN95 at no cost and install
it in Windows. The link is on the Course Website. It comes with its own
parochial IDE : Plato. The Silverfrost FORTRAN95 Compiler may be
integrated into MS Visual Studio 2013, available at no cost to students
from Microsoft. See the Course Web Page for details on how to get Visual
Studio at no cost.
• Anyone may download and install the GNU Compiler Suite including
gfortran. See the Course Web Page for the relevant link. You can then
choose to use:
(i) Netbeans as your IDE,
(ii) Eclipse as your IDE and or
(iii) run directly in the MinGW Shell in command line mode or in the
cygwin terminal shell.
• Again, the links (and hopefully, the installation and use instructions) are
available on the Course Web Page.
17/07/2015
Fortran 2014 Lecture #1
32/137