Your Federal Quarterly Tax Payments are due April 15th

# Tutorial introduction to Fortran 90 - PDF by msh36617

VIEWS: 0 PAGES: 6

• pg 1
```									    Tutorial introduction to Fortran 90

A simple program

The study of physics usually begins with mechanics, so it makes sense to start our study
of computational physics in the same place.

The intent of this tutorial is to allow you to get some familiarity with the programming
process — creating a code ﬁle, compiling and debugging a program — and not to teach
you to program. The program examples can be typed in as is, and along the way you will
get some idea of the language.

Of course, there are other good reasons. Mechanics has evolved methods of extreme
sophistication in the 350 years since Newton. One of the main impetuses has been the
study of orbital motion in the solar system. But these sophisticated methods are there to
provide a good starting point for perturbation treatments; often obscured in the beauty
of the formalism is the fact that there are few exact solutions. So numerical treatments
— starting already with Euler — have provided a powerful alternative.

Let’s start oﬀ with something very simple: dropping a ball from, say, the leaning tower
of Pisa. The equation of motion is
d2 y
m        = −mg
dt2
with y the height, t the time, m the mass of the ball and g the acceleration due to gravity.

Of course, this has the solution at time t, assuming an initial height y0 and velocity zero
at time t = 0,
1
y = y 0 − g t2
2

We can begin our programming in the simplest way, using the computer to evaluate this
mathematical expression. We will evaluate the height for a set of times 0.1, 0.2, . . . , 1.0
s. We will do this pretty much in the same way we would do the task with a scientiﬁc
calculator, so think of the program as a list of instructions for using such a tool.

Here is a simple program in Fortran to calculate the height as a function of time:

program drop0
2   ! drop0.f90 - rwf
! Simple Fortran program
4   ! straightforward evaluation of expressions

6   implicit none

1
8   real : :   g =9.8         !    accel due to gravity (m/s^2)
real : :   m =1.0         !    mass of ball ( kg)
10   real : :   y0 = 1 0 . 0   !    initial height 10 m
real : :   v0 = 0 . 0     !    initial vertical velocity component (m/s)
12   real : :   tmax = 1 . 0   !    drop time (s)
real : :   dt = 0 . 1     !    time step
14   integer    : : Nsteps , step
real : :   y
16
! determine number of time steps:
18   Nsteps=tmax / dt+1
! calculation loop
20   do step =1, Nsteps
t = step ∗ dt
22       y = y0 −0.5∗ g∗t ∗∗2
print ∗ , t , y
24   end do

26   end program drop0

(Note that the line numbers printed here are just for reference below, and are not part of
the program!)

• Comments, which are there to provide explanations of what is going on, start from
an exclamation mark and continue to the end of the line.
• Note that the program is “bracketed” by program ... and end program ... state-
ments. The name given in the statement is arbitrary, and does not have to be the
same as the ﬁle.
• implicit none is a statement that turns oﬀ some features of older versions of
Fortran.
• Note that all variables used in the program are declared in the ﬁrst section of the
program; some are deﬁned to have initial values.
• Lines 8-15 deﬁne various constants and variables. Numerical values could simply
be used instead in the appropriate places, but it is usually better to give names to
• This program uses a do loop to iterate over a range of integers which are used to
calculate the time and hence the height of the ball.
• In line 18, the number of steps to be taken is calculated.
• Line 20 introduces the do loop: The variable step is given values in succession from
1 to Nsteps and the indented portion below is executed. Note that the statement
ends with a statement end do. Do loops can be contained in other do loops, but
cannot overlap.

2
• Lines 21-23 constitute the body of the loop. Note that the indentation is not required,
but is useful for seeing the structure of the code.

• Note the exponentiation operator is **.

• The print statement will print the values of the variables in its argument list on the
terminal. The asterisk * is a placeholder for a format statement, and implies a free
format. You can specify the format of the data by using a proper format statement
instead of the asterisk; this allows you to specify number of signiﬁcant digits, etc.
Note the comma after the *.

Type this into a ﬁle. The comments are unimportant from the compiler’s point of view
and can be omitted — they are there for you to keep track of things. Save the ﬁle with a
.f90 extension, e.g. something like drop0.f90.

Now, before we can run the program, it must be compiled. We will use the program ifc
(Intel Fortran compiler) for this. The command for this is

ifc -o drop0 drop0.f90

This takes the ﬁle drop0.f90 and creates an executable program drop0 (that’s the -o
drop0

If the command prompt returns after this with no messages all went well. On the other
hand it is probable that you got some stream of error messages and warnings. You will
have to see what they say and try and ﬁx your code, maybe with some help.

OK, it compiles cleanly. Now run it:

./drop0

This line runs the program drop0 that exists in the current directory (that’s what the ./
indicates). (Technical note: the Unix PATH, which lists where executable programs exist,
usually excludes the current directory for security reasons, so you need to explicitly set
the path).

You program runs! You get a list of numbers: do they seem sensible?

Now make the problem a little less trivial by including the drag force, so that the equation
of motion is
d2 y             1
m 2 = −mg + CρAv 2
dt              2
Here, C is a constant; ρ the density of air and A the cross-sectional area of the body.

3
This actually has the solution,

1            1 + (γv0 )2
y = y0 +        ln
γ2g           cosh γgt

where γ =        CρAv 2 /2mg and y0 is the initial height.

This expression has a number of mathematical functions in it, including the square root.
There are suitable versions of these supplied by the system. Here is a simple program in
Fortran to evaluate this expression. You can enter this, but it’s quicker to copy your ﬁrst
program to a new ﬁle (cp drop0.f90 drop1.f90) and modify the new ﬁle. Use what
you have.

program drop1
2   ! drop1.f90 - rwf
! Simple Fortran program
4   ! straightforward evaluation of expressions

6   implicit none

8   real : :    g =9.8         !        accel due to gravity (m/s^2)
real : :    m =1.0         !        mass of ball ( kg)
10   real : :    y0 = 1 0 . 0   !        initial height 10 m
real : :    v0 = 0 . 0     !        initial vertical velocity component (m/s)
12   real : :    tmax = 2 . 0   !        drop time (s)
real : :    dt = 0 . 1     !        time step
14   real : :    C =0.5         !        drag coefficient
real : :    rho = 1 . 0    !        density of air ( kg/m^3)
16   real : :    R =0.1         !        radius of ball (m)
real : :    pi =3.1415926
18   integer     : : Nsteps , step
real : :    y , t , gamma
20
! determine number of time steps:
22   Nsteps=tmax / dt+1
! calculation loop
24   do step =1, Nsteps
t = step ∗ dt
26       gamma=sqrt ( C∗ rho ∗ pi∗R ∗ ∗ 2 / ( 2 ∗ m∗g ) )
y = y0 + ( 1 . 0 / ( g∗ gamma ∗ ∗ 2 ) ) ∗ log ( sqrt ( 1 . 0 + ( gamma ∗v0 ) ∗ ∗ 2 ) / ( cosh ( gamma ∗g∗t ) )
)
28       print ∗ , t , y
end do
30
end program drop1

• In many languages you will need to import the math functions from a library or
module. Fortran has them built into the language. (You can of course import other
things from a library or module).

4
• Note the way of calling a function — this is very similar to obvious mathematical
notation. A function typically returns a value.

• Note that trigonometric functions expect their argument in radians.

• The system does not supply a value for π.

The distance fallen in a given time is presumably smaller than in the original code, where
there was no air resistance. How much smaller? Calculate it! Let’s add a line to calculate
the vacuum value of free fall, say yvac = y0-0.5*g**2. Add this line to the do loop, and
alter your program to print it out. Note that you have a new variable name.

How about a graph? The simplest way to do this is to use a program like gnuplot. This
requires the data to be plotted to be in a ﬁle. The simplest way to do this is to cut and
paste from the screen to a text editor, and save the ﬁle, for example as drop1.out. The
start gnuplot and plot it:

gnuplot> plot ’drop1.out’

Hmm. This only plots one set of data. You have to tell the program which columns to
plot:

gnuplot> plot ’drop1.out’ using 1:2, ’drop1.out’ using 1:3

This can be abbreviated; you might like to use lines instead of points.

gnuplot> plot ’drop1.out’ u 1:2 with lines, ’drop1.out’ u 1:3 w l

Now let’s change the program so that it runs until one of the objects hits the ground at
y = 0. We can do this by replacing the do loop with a do-while loop:

do while ( yvac > 0 . 0 )
2       step = step + 1

While we’re at it change the distance fallen to 50 m.
Now, when we run the program, the data scrolls out of view. It’s time to write it to a ﬁle
from within the program. We can do so by ﬁrst opening a ﬁle

open(1, file=’drop2.out’)

5
Note that the ﬁle name is a character string enclosed in quotes. the “1” is a logical unit
used to refer to the ﬁle. To write to the ﬁle we can replace the print statements with:

write(1,*) t,y,yvac

Here the “1” stands for the ﬁle; the “*” is still a placeholder for a format. Note that there
is no comma before the t. Finally, we can close the ﬁle (after we have written everything
to it.

close(1)

Plot the results using gnuplot.

6

```
To top