Modelling and Simulation Laboratory

This `virtual laboratory' contains experiments in modelling and simulation using a number of different techniques.

Before starting any of the experiments, read about the philosophy of modelling programming languages (Section L2.1)

Section L1: General Facilities

Tools to use on the WWW:

Section L2: Spreadsheet generation

Documentation
Programs

Section L3: Fortran Programming Laboratory

This material is optional and is included for the benefit of those with previous Fortran experience.
Section L1: Single nonlinear algebraic equation solver

Section L1.1 Single nonlinear algebraic equation solver

Lower bound of x: Upper bound of x:
Enter l.h.s. of equation in x (just the lhs, not the `=0'):
Solution:

The left hand side of the equation should be typed in the box provided using only `x' as the unknown.

The expression must use JavaScript syntax which accepts the usual operators except either `^' or `**', as well as exp, sqrt, sin, cos etc.

log is natural log. To raise to a power use pow(x,n) .

The `Step' button takes a single bisection step. Click it again to to another step. The `Solve' button takes 20 steps. Another click gives another 20.

Section L1.2: Linear equation solver

Section L1.2: Linear equation Solver

Select the number of equations

Section L2.1: Philosophy of Modelling Languages

Section L2.1: Philosophy of Modelling Languages

Models can in principle be written in any computer language. A large number of special purpose languages have been developed to represent the steady state behaviour of process flowsheets. Rather fewer specialised languages are available for modelling at the level with which most of this course in concerned, and the majority of process models are in practice described using general purpose computer languages.

One of the most widely used general languages historically is Fortran (section L3.1), and many large models, particularly of special purpose process units such as reactors will be found in this language. We have given a number of examples of algorithms in Fortran.

A very large number of models in the last few years have been written to run in spreadsheet systems. The spreadsheet is extremely convenient, as it can be found on every PC and most palmtops. However, spreadsheets have several very serious disadvantages:

As a result of the above, manually generated spreadsheets cannot be regarded as a safe modelling tool! The most important property which a model must have is that it should be UNDERSTANDABLE. This is much more important than that it should be correct!

However, the spreadsheet is so universal and so convenient for running models (as opposed to building them) that we have developed a simple modelling language. which avoids most of the above problems by providing a high level description of the model, but generates a spreadsheet which can be used on any PC.

The main purpose of this language is thus to facilitate the production of documented models. In its simplest version (release 2.1, September 1998 for algebraic equations) it provides no solution facilities not available to the standard spreadsheet system.

The Elements of a Modelling Language

The most important function of a modelling is that it should encourage the production of readable models. All models contain the following elements which must be represented in a modelling language: Additionally it is highly convenient to be able to represent: The definition of the model does not itself have to say anything about how the model is to be solved, although it may be convenient to include instructions about the nature of the solution, i.e. is it a steady state, unsteady state or optimisation problem which is to be solved. In addition it will usually be necessary to provide various `housekeeping' information about what length of time or degree of precision is required, and so forth. This is not strictly speaking part of the model itself.

The model itself should be easy to read and self explanatory. In addition, a modelling system may usefully provide additional and/or summarised information information about the model.

A Simple Modelling Language

Our language contains the three main elements above. The model is written in standard (ascii) characters, using a text editor or notepad tool. (Do not try and write it in Word!) In general each model element should be written on a separate line, but apart from the special words and phrases which introduce different parts of the model, they can usually be written on one line separated by semicolons or commas if this makes the model more readable.

The model can, and should, contain plentiful comments which explain what different parts of it mean. These are preceded by any of the symbols `!', `%' or `#' and are ignored by the modelling system.

The language can be used to model both steady state systems, represented by algebraic equations, and dynamic systems represented by differential-algebraic equations. The discussion below pertains specifically to algebraic equations.

Variables

The user specifies the variables in the model by listing them between two lines which say: 
 variable
and 
 end variable

The ability to give user specified names to variables is a key feature of any modelling language, and the choice of `meaningful' names is an important factor in making a model comprehensible.

For example, a very simple model is to be used to calculate the mass density in kg/m of an ideal gas at a given Kelvin temperature and pressure in atmospheres. We are also given the molecular weight of the gas. We decide that we need variables to represent the specific molar volume and mass and molar densities. Although the pressure is given, and is not an unknown or variable (it is a parameter, see below) the Gas Law requires pressure in pascals which the model will require to calculate, so we require a variable for this pressure, giving 4 variables in all. The model instructions are then: 

variable
 vmolar   ! specific molar volume
 romolar, romass ! densities
 P        ! N/m2 required for gas law
end variable
Note the use of both standard symbols for standard quantities, e.g. P for pressure, and mnemonic names for other quantities.

Parameters

Parameter is a term used by modellers to mean a `variable constant', that is a quantity which is known, but which the modeller may wish to change. Here the specified temperature, pressure and gas molecular weight are parameters. They are introduced into the model in a similar way to the variables, i.e. between two lines: 
 parameter
and 
 end parameter

However parameters must be given values, which is done as in the example below. 

parameter
 MW = 16 ! methane, kg/kmol
 T = 273 ! K temperature
 Patm = 1.0 ! atmospheres, sensible units
 R = 8.314  ! gas constant in SI units: gmol, Pa etc
end
Note the appearance of R, the gas constant in this section. strictly, R is a true constant, having a fixed value, but most modelling languages do not distinguish these from parameters.

The modelling languages requires initial values for parameters, but these can be changed once the model has been generated. Do not however change the values of constants like R!

Equations

Equations are written in normal `computer algebra' notation between the lines: 
 equation
and 
 end equation

The simplest version of the model generator requires that all equations be written as formulas, and that every unknown (variable) appears once and only once on the LHS of a formula. This approach is discussed in the general introduction to algebraic equations (section 3.1). The equations do not have to be ordered as described in section 3.2, as this is done by the spreadsheet system. The equations section for this model is as follows: 

equation
 vmolar = R*T/P           ! m3/gmol
 romolar = 1/ vmolar      ! gmol/m3
 romass = romolar*MW/1000 ! kg.m3
 P = Patm * 101300
end equation

Housekeeping

The only additional information required for this model is to tell the solver that it is a purely algebraic or steady state model. This is done by putting the word steadystate in the model.

The order of the sections is currently fixed and is slightly different from that in we have introduced them. The more logical, and required, order is:

  1. Parameters
  2. Variables
  3. Equations
  4. Housekeeping and solver information
The whole model is enclosed by the keywords 
 model SomeName
and 
 end model
`SomeName' is a name which you can give the model.

Note that all names used in the model must be a combination of letters and number only , and must start with a letter. Upper and lower case letters are treated as being the same.

Complete Model

Here is the complete model: 

model idealgas

! ideal gas density calculation

parameter
 MW = 16 ! methane, kg/kmol
 T = 273 ! K temperature
 Patm = 1.0 ! atmospheres, sensible units
 R = 8.314  ! gas constant in SI units: gmol, Pa etc
end

variable
 vmolar   ! specific molar volume
 romolar, romass ! densities
 P        ! N/m2 required for gas law
end

equation
 vmolar = R*T/P           ! m3/gmol
 romolar = 1/ vmolar      ! gmol/m3
 romass = romolar*MW/1000 ! kg.m3
 P = Patm * 101300
end equation

steadystate

end model

You can copy it from here, and generate a spreadsheet to run it here.

An Excel-compatible spreadsheet translation may be copied from here and when loaded should look like this.

Modelling Dynamic Systems with Differential Algebraic Equations

The use of the language to model DAE systems is described in section L2.3
Section L2.2: Iterative Solution of Single Algebraic Equations

Section L2.2: Iterative Solution of Single Algebraic Equations

The method of using the spreadsheet model generator for sets of equations suitable for direct, non iterative solution has already been described in section L2.1

The model generator does not itself provide a mechanism for solving equations iteratively. However the facility for solving single equations using a `goal seek' facility is available in nearly all spreadsheet systems.

We therefore use the model generate to provide evaluation of the function or left hand side of the equation which it is required to reduce to zero. An illustration for one of the exercise examples is shown below. 

model single

variable x, fx
end variable

equation
  fx = x + @ln(x)
end equation

steadystate

end model
Notice that there is only one equation here, but there are two quantities declared as variables, x and fx. The latter is not really an unknown, it is the function value which is to be reduced to zero. using `goal seek'. (Note in passing the form of a function call in the model generator; the function name must be preceded by a`@' whether or not the spreadsheet requires that convention. Also that Excel appears to interpret LOG as base 10, although the JavaScript solver referred to elsewhere, most other programming languages and mathematical convention treats this as the natural log!)

Alternatively, if the spreadsheet does not have a search facility, it is possible, for single equations, to use the logic of a bisection search by hand.

The above model program can be copied from here. Copy it and try both an automatic solution (if available) and a hand directed trial and error.

The model generator is here.

A .slk format translation is already available here. You can copy it and load it into Excel. Cell C5 contains the unknown x (defaulted to zero) and cell D5 contains the function value fx. Initially D5 displays an error message because an attempt has been made to calculate the log of zero.

Change C5 and observe changes in D5. Note that values of 0.1 and 1.0 bracket the solution. Try choosing values by hand to `home in' on the solution when fx in D5 will be nearly zero.

To use the `Goal seek' feature of Excel, proceed as follows:

This may seem a rather cumbersome procedure for dealing with single equations, but it becomes useful where we have a set of equations which can be solved by iteration in a single variable, see section 3.5.1

Section L2.3: A Simple Dynamic Modelling Language

Section L2.3: A Simple Dynamic Modelling Language

This section describes a simple language for describing models of dynamic systems. It is a text based language, and a WWW browser input form is available to help with model creation.

A model described in the language may be turned into Lotus-123 compatible spreadsheet code, saved and run using essentially any spreadsheet system, including Excel and the Psion spreadsheet.

The language allows the creation of differential-algebraic equation (DAE) models, although the solution facilities for algebraic equations in spreadsheets are somewhat limited.

The use of the language will be described using a simple example.

Example

A system is described by the following set of DAEs.

dx/dt = (z-x) / T1
dy/dt = (x-y) / T2
z = - A y

The initial conditions are: x = y = 0 at t=0

In the language, the equations are described as follows. 


model test
! Simple example dynamic model

parameter
! These are quantities which are known, but which
! the user may wish to change to see their effect
 A=1 ;  t1=1
 t2 = 2 
! They will normally be given initial values
end parameter

variable
! These are the variables or unknowns of the problem
 x1 = 0 ! They may be given initial values
 x,  y  ! or not if preferred
end variable

Equation
! These are the equations

 y =  x*x ! This is an algebraic equation

! The  are differential equations...
 .x1 = (-A*y-x1)/T1
 .x = (x1-x)/T2
end equation

initial
! This is another way of giving a variable an initial value
 x=0
end initial

schedule
! This section defines how the model will be run
 deltat=0.5 ! Timestep
 steps=50   ! number of timesteps
end 

end model
It can be seen that the model consists of five sections, each introduced by a keyword, and terminated by the word end.

The model instructions are written in normal `computer notation'. There can either be one per line, or multiple instructions on a line, separated by either `;' or `,'. Multiline instructions are not allowed.

Any line beginning with a `!' is a comment and is ignored.

Each of the sections will be described briefly.

Parameter Section 


parameter
! These are quantities which are known, but which
! the user may wish to change to see their effect
 A=1 ;  t1=1
 t2 = 2 
! They will normally be given initial values
end parameter
`Parameters' as distinct from `unknowns' in a mathematical problem are quantities which the user knows, but which might be changed at some point to see their effect. It is convenient to give them symbols rather than simply numbers.

A parameter may often be a physical property, e.g. a heat capacity, which can appear several times in a set of equations. If the parameter is defined in one place, then making one change will result in the new parameter value being used throughout the model.

This section of the model enables parameters to be defined and given values if required. If no values are given the parameter defaults to zero. New values can be given when the spreadsheet has been generated.

Variable Section 

variable
! These are the variables or unknowns of the problem
 x1 = 0 ! They may be given initial values
 x,  y  ! or not if preferred
end variable
The variables are the unknowns in the model. As discussed elsewhere, there are two sorts of variables in DAEs, those which may appear as derivatives in the differential equations, and the `algebraic' variables which do not. This distinction is not made in this section, which serves mainly to define the names to be used for variables.

Variables may be given initial values in this section, as illustrated. These are relevant only for ODE variables, and will be ignored for algebraic variables. Initial conditions may also be specified in a separate section of the model.

Equation Section 

Equation
! These are the equations

 y =  x*x ! This is an algebraic equation

! These are differential equations...
 .x1 = (-A*y-x1)/T1
 .x = (x1-x)/T2
end equation
This is the main part of the model and contains the differential and algebraic equations which describe the model. Note the following: All normal algebraic operators and functions may be used in equations.

Initial Section 

initial
! This is another way of giving a variable an initial value
 x=0
end initial
This section provides another place where initial values may be given for variables. If no initial values is set for a variable which requires it, then it defaults to zero.

Schedule Section 

schedule
! This section defines how the model will be run
 deltat=0.5 ! Timestep
 steps=50   ! number of timesteps
end
This section describes how the model is to be run. At present this may specify any two of the three quantities:

General

The model is introduced by the keyword model and terminated by the words end model. all keywords must appear on separate lines.

There may be more than one of each of the sections, e.g. a variable section, an equation section, another variable section and another equation section. However, all variables must appear in a variable section before they are used in an equation.

Only one schedule section is allowed.

Section L2.4: Linear equation solver: Form 1

Section L2.4: WWW Linear equation Solver

There are three options for the linear solver which will be created. Two of these create a stand alone solver with a tabular front end in the spreadsheet with cells in which you can type coefficients. the third creates an `embedded' solver which has no user interface of its own, but allows you to specify the spreadsheet cells from which coefficients and constants will be taken, and into which the solution will be placed.

For the stand alone solver you may specify a full matrix in which you may later enter values for any or all of the coefficients. Alternatively you may specify a semi-sparse matrix in which you identify coefficients which are nonzero. The values of nonzeros may be entered into the spreadsheet cells later, but zero coefficients cannot be changed and will always be zero.

For the embedded solver, coefficients and constants may be either numerical values or spreadsheet cell names (given as RnCm). Cell names must be given for all of the unknowns.

In summary the 3 options are:

Please enter your Email address:

Check here if if you want to make an `embedded' solver

Select the number of equations

If you are outside the University of Edinburgh, the spreadsheet you create will be emailed to you, so you must give a valid email address. If you would like to try a `full matrix' spreadsheet there is a 12 equation solver available for direct downloading in three spreadsheet formats from the links below. You can save these in your local filespace and load them into any compatible spreadsheet.
Section L2.5: Algebraic and Differential Algebraic Equation Models

Section L2.5: Algebraic and Differential Algebraic Equation Models

If you are running from behind a firewall, you may find it more reliable to have the model emailed back to you.

Section L3.1: Fortran as a Modelling Language

Section L3.1: Fortran as a Modelling Language

More chemical engineering programs are written in one or other of the `dialects' of Fortran than in any other programming language. This is largely because Fortran has been around longer than any other useful language, and because chemical engineers were amongst the first engineers to realise the potential of the then newly invented computer to solve important and useful problems.

Most of these programs were written by chemical engineers. The ability to write large and complicated programs is much less likely to be required of chemical engineering graduates nowadays, but you may well encounter programs which you may have to understand or even modify in industry.

Fortran, as indicated, comes in various dialects, which have evolved from the earliest in the late 1950s through Fortran 66 (1966), Fortran 77 (1977, still sometimes encountered) to the current Fortran 90 which we will use. The earlier versions are rather horrible and unsuitable for teaching purposes. But if you do have to deal with, say a Fortran 77 program, then your knowledge of Fortran 90 and a reference manual should be enough to enable you to work out what is going on.

Fortran has all the necessary features required to represent a mathematical model. These are described in the following sections.

A Sequential Programming Language

A Fortran program at its simplest consists of a list of formulas as we have defined them in section 3.1. (In fact the name of the language stands for Formula translator.)

A key feature of the language is that the formulas are evaluated in the sequence in which they are written by the programmer. This makes Fortran different from a spreadsheet system which (in the more sophisticated versions) works out a `correct' order for evaluating its formulas.

In fact a Fortran program consists of other `instructions' all of which are carried out in the order in which the programmer has specified them. Creating this `ordered list', `recipe' of instructions or algorithm is what most of programming is concerned with.

We will describe the instructions necessary to create a simple model in the next section. It is helpful to think of Fortran instructions as falling into three main categories:

The use of the first two categories can quite easily be understood once a few examples have been examined. The third group are best just learned by rote or looked up in a manual when required.

Example: Ideal Gas Density Calculation

A simple model is required to calculate the molar and mass density of a gas of given molecular weight at a specified temperature and pressure in atmospheres, using the ideal gas law:

P v = R T

Remember that here v is the specific molar volume, i.e. m³/kmol, and that if R is 8.314 J/mol/K then pressure P must be in pascals.

Developing the model program

It is vital that any program be written down before you attempt to type it in to the computer. Under no circumstances try to `compose' a program at the keyboard!

Simple programs which consist essentially only of formulas to be evaluated can be developed in the following stages:

  1. Write down all the equations.
  2. Identify the unknown variables and parameters of the model, and check that you have as many equations as unknowns.
  3. Rearrange the equations into formulas.
  4. Rewrite the formulas in the appropriate order for evaluation: at this point you have produced the algorithm for solving the problem.
  5. Finally, rewrite the algorithm as a program by adding the necessary housekeeping instructions.

Step 1: write all the equations

The relevant equations here are the gas law from above:
P v = R T

The other equations are concerned with calculating molar density and mass density from specific volume:
= 1 / v

= M / 1000

Here M is the molecular weight in kg/kmol. The factor 1000 arises of course because R is in gram mols. If pressure Pa is given in atmosphere, then this must be converted to pascals using the relationship:

101300 Pa = P

Step 2: Identify and count the variables

There are 4 equations. The 4 unknowns are P in pascals, v and the two densities. All other quantities are parameters, i.e. specified quantities T, Pa, and M, or true constants like R.

Step 3: Rewrite as formulas

Every unknown variable must appear once and only once on the LHS of a formula. There may be more than one acceptable set of rearrangements. Here, after noting what quantities are given, the only possible rearrangements of the 4 equations are as shown below.

v = R T / P

= 1 / v

= M / 1000

P = 101300 Pa

Step 4: Rearrange in calculation sequence

Again, there may be more than one acceptable reordered sequence. A suitable order here is:

P = 101300 Pa

v = R T / P

= 1 / v

= M / 1000

Step 5: Write the program

This involves two separate but entirely mechanical activities.

Tun the formulas into computer notation...

The first of these involves firstly suitable names for unknowns and parameters which are acceptable to the computer; broadly this requires that they consists only of letters and digits and start with a letter. Upper and lower case letters are equivalent. Secondly the computer notation of `*' and `/' for multiply and divide must be used, along with the computer-recognised names and format for any higher functions, e.g. log(x) for ln x.

The working part of the program thus now becomes: 


P = 1.013e5 * Patm             ! P in Pa

v = R * T / P                  ! specific volume in m3/mol

rhomols = 1 / v                ! in mol/m3

rhomass = rhomols * M / 1000   ! ... in kg/m3
Note particularly the `comments' preceded by `!' which are ignored by the program, but help the reader to understand what is going on - remember, the most important thing about a program is that it must be intelligible to a reader!

Add the `housekeeping' ...

Before reading this section you may wish to refer to notes on how to write simple Fortran programs (section L3.2).

The housekeeping can be divided into the parts that make sense to the fairly uninitiated, and the other bits that just have to be learned by rote!

The first of these include the definition of `declaration' of variables and parameters. these must appear at the beginning of the program, logically enough, before any of the variables or parameters are used.

To declare variables which will contain general or `real' numbers we say: 

real :: P, v, rhomols, rhomass
To declare parameters or constants we use a similar construction, but include the values to be specified, e.g.: 
real, parameter :: Patm=2.0, T=298, R=8.314, M=16
Having calculated the solution in terms of variable values, these have to be communicated to the outside world by means of the computer instructions `print' or `write'. The simplest form for these here is: 
print *, 'Specific molar volume is ', v, ' m3/mol'
print *, 'Molar density is ', rhomols, ' mol/m3'
print *, 'Mass density is ', rhomass, ' kg/m3'
The text printed out with the values should be regarded as an essential part of the program. Logically, the values cannot be printed out on the screen until they have been calculated, so these instructions come after the calculation part.

The program is required to start with the word program and be given a name, and to end with the words end program , followed by the program name. Good programming practice dictates a comment on the purpose of the program, and author details. Finally the statement implicit none should be present to override the default implicit naming scheme.

The complete program is thus: 

program pvt

!this program calculates the specific molar volume, molar density
!and mass density for the given data using the ideal gas law
!JWP 11/12/98

implicit none

real :: P, v, rhomols, rhomass
real, parameter :: Patm=2.0, T=298, R=8.314, M=16

P = 1.013e5 * Patm             ! P in Pa

v = R * T / P                  ! specific volume in m3/mol

rhomols = 1 / v                ! in mol/m3

rhomass = rhomols * M / 1000   ! ... in kg/m3

print *, 'Specific molar volume is ', v, ' m3/mol'
print *, 'Molar density is ', rhomols, ' mol/m3'
print *, 'Mass density is ', rhomass, ' kg/m3'

end program pvt
You can copy the program from here to try it for yourself.
Section L3.2: Writing Simple Fortran Programs

Section L3.2: Writing Simple Fortran Programs

This section is intended to get you started writing elementary programs in Fortran. We will not attempt to explain all the peculiarities of the Fortran programming language. We will try and explain the logic behind those parts of the language which are easily explained, and treat as arbitrary rules to be learned by rote those which are either actually arbitrary, or too time consuming to attempt to explain.

Examples will be given on the assumption that the Edinburgh Portable Compilers EPCF90 compiler is being used on a Unix system.

The aim is to get you writing simple programs quickly, not to teach you all about writing complex programs. If you want to find out more about the Fortran 90 (F90) language and `real' programming, we will provide hypertext links and a bibliography which you can follow up.

A First Program

Here is the world's simplest `useful' Fortran program: 
program hello
print *, 'Hello World!'
end program hello

Program structure

The program starts with the word `program' simply because F90 programs must do so. The program must also be given a name, so we have called it `hello'. It must also end with the words `end program', hence the last line. The use of the program name in the last line is optional in most F90 dialects, but it is good practice to include it. (Your lecturer, however, often gets lazy and omits it. In fact just the word `end' is sometimes enough, but you should try to follow what everyone agrees is good practice!)

A program which consisted just of the first and last lines would be a `syntactically acceptable' program, i.e. it could be compiled and run as described below. However, it wouldn't do anything, and so would be pointless.

The `useful' bit of the program is the middle line: 

print *, 'Hello World!'
What this does is to print on the computer screen the text: 
`Hello World!'
This instruction is formally called a `procedure call'. It causes something to be done; what is does is predefined by the word `print' which the language recognises as meaning `print or otherwise display some information'.

What is printed is determined by the programmer-supplied text included between the two quote symbols. (NB these are both the same single closing quote on the keyboard.) Any text (except for another single quote!) which appears there will be `printed out'. A piece of text in quotes is referred to as a `text string' or `character string' or sometimes just a `string'.

And the `*,' ? For the time being, please just accept that it is necessary. If you leave it out, the program will not work!

For future reference: you will encounter other types of procedure call in the language. All have the properties that (a) what is done is predefined and associated with the name of the procedure, here `print'. How it is done, or with or to what is usually specified by the programmer, as is the case here with the text to be printed.

Running the program

If you click on this link you will get a copy of the program in your WWW browser. You should have learned how to save this in you own filespace. Do this. By default the program will be stored in a file called `hello.f90'. All F90 programs should be stored in files with the `filename extension' .f90 or else they will not be recognised as Fortran programs.

Once you have saved the program it must be translated from a form which you can read (sometimes called the `source language') into a form which the computer can run. This process is called `compilation'.

Assuming you have saved the hello.f90 in your current working directory, to compile the program type: 

epcf90 hello.f90
`epcf90' is the name of the command which runs the compiler. It is always followed by a space and then the name of the source language file. (And of course `return' at the end!) The computer will print a number of messages which should look something like this: 
... epcf90 hello.f90
EPC Fortran-90 Version 1.5.1.7 Sparc:031097:145225
Copyright (c) 1993-1997 EPCL. All Rights Reserved.
hello.f90
   program HELLO

3 Lines Compiled
If you look at your files, you will find that there are a number of new ones, mostly with peculiar names. One in particular will have the name a.out and this is the `compiled' or computer runnable version of your program. To run (`execute') this just type: 
a.out
and see what happens.

Things to Try with the Simple Program

Load the program into a text editor, e.g. by typing: 
 textedit hello.f90

Change the message

Change the text in the quotes to something else. `Save' the program file from the editor, compile it again and type a.out, which will execute the new version of the program.

Add another message

You can include as many `print' instructions of the form shown as you like, provided that each is on a separate line. Try adding another message this way.

EPCF90 and most versions of Fortran also allow instructions to be separated by a semicolon `;'. Thus the following are equivalent: 

print *, 'This is my 1st program..'
print *, 'Hello!'
And: 
print *, 'This is my 1st program..' ; print *, 'Hello!'
Try these two forms.

Printing on the same line

You may have noticed that the two messages above appeared on different lines even when the two instructions were on the same line. To print to messages on the same line they could either be made into a single piece of text within one pair of quotes, or written as two strings of text separated by a comma in one print instruction, i.e.: 
print *,  'This is my 1st program..' , 'Hello!'
Try this. This has a somewhat different effect from having the same message in one text string.

Add a message to the reader of the program

It cannot be stressed too strongly that a good program is easy to follow and makes sense to the reader! You can add `comments' which do not cause the program to take any additional action, but which help the reader to follow what is going on. Any line or part of a line beginning with an exclamation mark `!' is treated in this way.

A properly commented version of hello.f90 is here. 

program hello
! Very first example program
! Written by JWP, September 1998
!
print *, 'Hello!'  ! A single procedure call
! This will print a message on the screen
!
end program hello
If you compile this you will find it behaves exactly as the shorter version did.

Arithmetic in Programs

As well as printing messages, computers can of course do arithmetic. Here is a simple program to do a calculation: 
program arithmetic
! program to do a calculation
! Written by JWP 29/10/98

print *, 'This program calculates 23.6*log(4.2)/(3.0+2.1)'
print *,  23.6*log(4.2)/(3.0+2.1)
! Note the effect of the presence or absence of quotes! 
end program arithmetic
You can copy it from here.

Save, compile and run this program.

Important Point! Notice the difference between how 23.6*log(4.2)/(3.0+2.1) is treated:

In the first case it is just a string of text and is, in the context of a `print' instruction, simply displayed exactly as it appears.

In the second it is treated as an expression and is evaluated using the computer's rules for arithmetic. These are more-or-less the same as the conventional rules but include the following:

Two possible pitfalls lurk for the unwary:

The point about being able to print two quantities in the same print instruction may now become more obvious. We can have both expressions and strings in the same instruction so as to print on the same line, e.g.: 

print *, 'This value is ', 23.6*log(4.2)/(3.0+2.1)

Variables

This is a key section which you must understand before proceeding further! Seek help if you are not entirely clear about it after trying the examples.

Consider the following program which you should copy from here: 

program variable
!this program evaluates an expression
!written by JWP 28/10/98

real :: x ! This defines a variable called x
!
x = 23.6*log(4.2)/(3.0+2.1) ! This gives it a value
!
print *, 'The value of x is ', x 
! .. and then print it out
end program variable
The overall effect of this is essentially the same as the previous program, except that it involves a new kind of entity, vital to effective programming, called a variable. It is so called because the value which it takes can be changed by equating it to different expressions on the right hand side of an instruction which always takes the form: 
variable_name = expression
The form of this assignment instruction is identical to that of evaluating a formula, as discussed in section 3.1. It is convenient to think of a variable in this context as being the same as a variable in a set of equations, although it is in fact a more general and flexible concept than that. At the very least, the use of variables facilitates the writing and developing of programs, and makes them easier to understand.

Consider the use of a program to evaluate expressions involving variables whose values depend on other variables, i.e. a set of equations. This involves a slightly modified version of the previous program: 

program variable
!this program evaluates two expressions and prints the result
!written by JWP 29/10/98

real :: x, y ! This defines variables called x and y
!
y = sqrt(23.6)              ! calculate y
x = 23.6*log(1.0+y)/(3.0+y) ! calculate x which involves y
!
print *, 'The value of x is ', x 
print *, 'The value of y is ', y 
! .. and then print them out
end program variable
You can copy this program here.

The instruction: 

real :: x, y
is called a declaration of the variables x and y. It consists of the word real to indicate that the variables will be used to represent general or real numbers (as opposed to whole numbers or integers, which can also be be used in certain circumstances, see later).

It is followed by a list of the names to be given to the variables. These can be any combination of letters or numbers, but must start with a letter. The names may not include spaces, but the underscore character `_' may be used and is very useful in improving readability.
Invalid variable namesValid variable names
1stfirst
not correctthis_is_correct
_invalidx4
a-bad-nametemperature
no.dots.allowedvolume

The two colons `::' are part of the illogical richness of Fortran. They may sometimes be left out, but this is not advised.

You should try to choose variable names which help the readers understanding of the program. In general you should not use meaningless name like `x' and `y', unless these have been used in the specification of the problem. Where standard symbols exist, e.g. `P' for pressure or `T' for temperature, use these with variants, e.g. T1, T2 or use descriptive names, e.g. T_reactor, Column_Pressure. Names can be up to 31 characters long, but typing very long names becomes tedious (and they can be mis-spelled). NB. Fortran is case insensitive. Thus FIRST, First, and first all refer to the same variable.

Implicit none

There is a naming convention in Fortran, which states that all undeclared variables starting with letters I through to N are assumed to be integers, another variable type. Otherwise it is assumed to be a real variable. This convention can lead to some very confusing runtime errors. To avoid errors put the statement 
implicit none
at the top of every program. This is a sign of good programming style or practice.
This statement should be at the top of every program. Its absence will result in penalties in tests or handins.

Parameters

The programming language also allows for the use of parameters which are named entities, like variables, but are given a value once, and may not have their values changed. Their use is explained in the gas law example.

Parameters are declared in a similar manner to variables, but must be given their constant values at the same time as they are declared. For example, pi should remain constant throughout a program: 

real, parameter :: pi=3.14

In this very simple program we have not really exploited the concept of variables whose values may be changed. It would have been possible to decide that the entity we have called y was a parameter. This has been done in this second version of the program. since a name like `y' is typically associated with the unknown in a problem, the name has been changed here, although this is not strictly necessary. 

program variable
!this program demonstrates the use of parameters
!written by JWP 29/10/98

implicit none

real :: x ! This defines variable called x 
real, parameter :: fixed_y = 4.86 
!
x = 23.6*log(1.0+fixed_y)/(3.0+fixed_y) ! calculate x which involves y
!
print *, 'The value of x is ', x 
print *, 'The fixed value of y is ', fixed_y 
! .. and then print them both out
end program variable

Program Structure

Good programming practice (and sometimes the Fortran standard) dictates that a program should always have the same rough structure. The above example shows this structure.

More Programs

All the programs you have been shown so far are `correct', i.e. they should compile without any error messages and run to print out the required results. The first programs you write will almost certainly not be correct!

To give you practice in `debugging' programs, here is one which will fail to compile. 

program crunch
! there is/are something(s) wrong with this program
! find out what and correct it
!written by JWP, 23/10/98

implicit none

real :: a, b, c
a = 10.0
b = 23.5
c = a*b
print *,'What is wrong with this?
print *, 'Nothing now..'
print *, c is, c
!
end program crunch

And here is another which may appear to be correct but will fail. 

program crunch_again
! This program has a different kind of error from the first one

!written by JWP 23/10/98

implicit none

real :: a,b,c
print *,'Does this work?'
a = 10.0
c = 3.5+a
c = c*b
b = 23.0
print *, 'a, b and c are:', a, b, c
!
end program crunch_again

Copy these, try them and correct them.

Section L3.3: Further Ideas in Fortran Programming

Section L3.3: Further Ideas in Fortran Programming

This section covers:

Repetition

Much of the power of computers derives from their ability to carry out repetitive calculations very quickly. In previous examples the instructions we have written have each been carried out once and only once every time the program was run.

It is very easy to make a program repeat a sequence of instructions an indefinite number of times. However the example we are about to give is an example of very bad programming practice which should never be emulated! it is provided here simple to illustrate how easy concept is to implement.

Any instruction or sequence of instructions will be repeated indefinitely if it is preceded and terminated respectively by the two instructions shown below: 

do
   ! instructions to be repeated here
end do
It is considered good practice to indent every instruction between the do and end do by three character spaces.

Here is a very simple program which illustrates this. It is a very bad program, because as written it will never stop running! 

program zap
! This is not a good program !!!
! Before running it, make sure that you can identify
! `panic button' which stops programs on you computer.
! On Unix systems it is simultaneously `Control' and `c'
!JWP 29/10/98

implicit none

do
 print *, 'ZAP!!'
end do

end program zap

If you are going to compile and run it, first make sure you know how to do an `emergency stop' on your computer system. This is usually achieved by pressing simultaneously the letter `c' key and the key marked `Control' or `Cntrl'.

The principle should however be clear from this example....

To be useful, the repetition clearly needs to be under rather better control! A useful way of achieving this is to count the number of repetitions and stop after a predetermined number. the means of doing this is described next.

Counting and further variables

To repeat a set of operations, say, 10 times, we use a modified construction as follows: 
 
do times = 1, 10
 print *, 'ZAP!!'
end do
The only difference between this and the previous example (so far) is in the first line. The interpretation of this can be consider to be as follows:

Like any other variable, times must be declared at the beginning of the program, before it is used. However, because times will have values which are always whole numbers or integers, the form of the declaration is different from those seen so far: 

integer :: times
Here is the complete, and much improved, program: 
program zap
! This is a better program !!!
! JWP 29/10/98

implicit none

integer :: times

do times = 1, 10
   print *, 'ZAP!!'
end do

end program zap
Copy it from this link.

The counter variable has other uses than just stopping this program loop. Its value is 1 the first time round, 2 the second and so on.

Try changing the instruction in the loop to: 

print *, 'ZAP times ..', times
and confirm this.

There are a number of other uses for integer variables, but a major one is controlling loops like this. The number of times the loop is to be repeated, and, associated with it, the starting and finishing values of the counter variable, must be specified before the do instruction is encountered, but may in fact be general integer expressions.

Note that in this example the numbers `1' and `10' were written without decimal points, whereas these were required for numbers in the general expressions encountered previously. This is because these starting and finishing values must be whole numbers. The following would also be allowed, try them, printing out the counter value on each loop:

The counter variable may be also be used in any expression inside the loop. Note that integers may generally be used in expressions like real non-integer variables and numbers. Care must however be taken if real and integer variables are to be used in the same expression. This is illustrated in the examples below.

The following sections provide programs to illustrate the use of repetition loops and a number of other programming techniques. You are encouraged to copy and modify these as described.

Example: a table of squares

This program calculates and prints a table of squares. 
program squares
! calculate squares
! JWP 29/10/98

implicit none

integer :: number, square
integer, parameter :: max = 10
!
print *, 'Squares of numbers from 1 to', max
!
do number=1, max
   square = number ** 2
   print *, 'Square of ', number, ' is ', square
end do
!
end program squares

Modify it to calculate and print the squares of numbers from 1 to 15 by altering only a single character in the program.

Example: a table of square roots

The function sqrt() will calculate the square root of an expression or variable inside the brackets, However, this must be a real. However there is a function called real() which converts an integer to a real value as shown in this program. 
program squares
! calculate squares
! JWP 29/10/98

implicit none

integer :: number
integer, parameter :: max = 10
real :: real_number, root
!
print *, 'Table of square roots'
print *, 'Number  ', '        Square root'
!
do number=1, max
   real_number = real(number)
   root = sqrt(real_number)
   print *,  number, root
end do
!
end program squares

Square roots can also be calculated by raising numbers to the power 0.5, or in Fortran writing ** 0.5 Check that this gives the same answers, and modify the program to tabulate cube roots.

Non-integer Loop Counting

The variable used in the do instruction in general must be an integer. However we sometimes want a quantity which will will take non-integer values to be changed each time we go round the loop. Below are parts of a program which tabulates square roots of number at from 1.0 to 3.0 at 0.25 intervals. 
integer :: i 
real :: x
.....
do i = 1, 12
 x = i * 0.25
 print *, x, sqrt(x)
end do
...
There is a more general way of doing this described later.

Example: `recursive' evaluation

What would be the effect of the following assignment instruction?: 
x = x + 3
If you are in any doubt about what will happen, complete the following program segment by adding the necessary additional instructions, and run it. 
.....
x = 4
print *, 'x was', x
x = x + 3
print *, 'x now is', x
.....
There are many circumstance where it is convenient or necessary to use this kind of `recursive' instruction. One such is illustrated in the example below. Complete this program and use it to compute some factorials. Check to see that it gives the correct answers! 5! is 120, 10! is 3628800. What happens when you try to calculate 100! which is about 10158?

Example: calculating factorials

Recall that n! is given by:

n! = 1 . 2 . 3 .... (n-1) . n

Consider the following program segment: 

...
integer :: fact, n, i
n = 5
fact = 1
do i=2,n
 fact = fact * n
end do
...

Exercise

Write a program which you can use to calculate the sum of the first n integers where n is specified as a variable.

Recall that this quantity is equal to ½(n+1)n.

More Non-integer Loops

Consider the problem of printing out 10 temperatures spaced equally between the boiling points of propane (231.0K) and butane (272.6K). Because these are both non-integer and not obviously related, the method used above is not convenient. the following program steps show how it can be done. 
real :: Tprop, Tbut, T, dT
integer :: i, n
....
n=10   ! number of points
...
dT = (Tbut - Tprop) / (n-1)  ! divide range by number of INTERVALS
! .. this gives size of interval
T = Tprop ! Start here
do i=1, n
   print *, T  ! the current temperature
   T = T + dT  ! next temperature
 ...
end do
...
Exercise: Write a program to tabulate the vapour pressures of propane and butane at 12 equispaced temperature points between their respective 1atm boiling points.

Decision Making: Conditions

Consider the following sequence of instructions: 

....
if ( a  >=0.0 ) then
   print *, 'a is zero or positive'
else
   print *, 'a is negative'
end if
...
Their implication should be fairly obvious. If the variable a has a value greater than or equal to zero, then the instruction: 
 print *, 'a is zero or positive'
is carried out and the appropriate message printed. On the other hand if this is not the case (`else') then the instruction: 
 print *, 'a is negative'
is obeyed.

The program would then proceed with any instructions following end if.

The above is an example of the most general form of the conditional or `if' construction used for decision making in computer programs. The lines between then and else and between else and end if may contain any sequence of legal `executable' instructions. (This excludes `formal' instructions like real :: or end program.) As with do loops, note the indentation of the instructions by three character spaces. This is to give clarity to the program.

As an example of using this for of the `if' construction, here is a program to solve a quadratic equation using the quadratic formula. The program tests to see if the equation has real roots before proceeding with the solution. 

program quadratic
! solve a quadratic with real roots
! JWP 22/10/98

implicit none

real :: a, b, c, d, x, x1, x2
!
a = 1
b = 2
c = -3
!
print *, 'Solving quadratic:'
print *, a, ' x**2 +', b, ' x + ', c
! check nature of roots..
!
d = b**2 - 4.0*a*c
!
if (d<0) then
   print *, 'Has complex roots.'
else
   x1 = (-b + sqrt(d)) / 2.0 / a
   x2 = (-b - sqrt(d)) / 2.0 / a
   print *, 'Roots are:', x1, x2
end if
!
end program quadratic

Three other forms of the `if' construction are useful. Sometimes it is only required to carry out a series of instructions if a condition is satisfied. This is achieved by the following. Note that the symbols /= are the Fortran version of the `not equals' sign. 

if (a /=0.0 ) then
   print *, 'a is not zero so can be a divisor'
   x = 1 / a
end if
Notice that in this case nothing happens if the condition is not satisfied, i.e. if a is equal to zero. In this case the value of x is not set by any instruction shown in this program segment.

If only a single instruction is to depend on the outcome of the condition then it may be shortened to: 

if (a /= 0.0 )  x = 1 / a
Finally, by using else if several conditional statements may be chained together: 
if(a > 0) then
   print *, 'a is bigger than 0'
else if (a < 0) then
   print *, 'a is less than 0'
else
   print *, 'a must equal zero'
endif

Computer form of relational operators

Note the `computer form' of the the mathematical relational operators.

`<' and `>' have their obvious forms, <= stands for `less than or equal to' and hence the obvious interpretation of >=.

'Not equal to as above is /=.

Confusingly `equal to' must be written with two equals signs, i.e. as ==.

Example

Below is a program which generates and prints random numbers. 
program random
! generate some random numbers
! JWP 22/10/98

implicit none

real :: x
real, external :: rand
integer :: i

do i= 1, 10
   x = rand(0)
   print *, x
end do

end program random
(Note the special form of declaration: real, external :: rand which says that rand, a function to generate random real numbers in the range 0.0 to 1.0, is to be found outside, i.e. `external' to the program. This is one of a number of useful function which exist in special libraries.)

Copy the program from here and modify it (a) to print only random number greater than 0.5, and (b) to count and subsequently print out the number of such numbers.

The exit Instruction

The instruction: 
exit
which can only appear inside a do .. end do loop causes a immediate termination of the loop. The program carries on at the instruction following end do.

The do may be either with or without a counter variable. In the latter case the count is overridden and the loop stops at whatever is its current value.

The exit instruction is normally always the subject of a condition, e.g.: 

if (TK <= 0.0) exit

Other Methods of Terminating Loops

There is another method of stopping a do...end do loop using a conditional statement: 

do while (x == 1)
   .
   .
   .
end do
This is equivalent to saying: 
do
   if(x /= 1) exit
   .
   .
   .
end do
The do while format is preferred by purists as the loop can only terminate after a complete iteration.

Further exercise

Adapt the modified vapour pressure calculation program so that it can be stopped by typing in an infeasible temperature, i.e. one less that O Kelvin.

Interacting with your program

This section has been left until last because we do not particularly want to encourage you to write `interactive' programs until you are reasonably experienced in writing self contained programs.

All the programs shown to you so far have been `self contained' in that all the information required for the program to work is contained within the program.

All your handin examples should normally be written in this way. A major advantage of such a program is that it's operation can be checked by reference only to a printout of the program. A disadvantage is that every time we want to change a parameter then the the program file must be edited and recompiled.

We will now discuss how to make programs more `user friendly' and `interactive'.

When ever the computer encounters an instruction like the following: 

read *, x
it does two things.

First it stops and waits for the user to type something on the keyboard. Secondly it interprets whatever has been typed in the context of what kind of variable x has been declared to be, and if it can, assigns whatever has been typed in to x. This should be clear from the following complete program which you might like to copy and try: 

program readin
real :: x
print *, 'Type in a valid real number:'
read *, x
print *, 'The number you typed was ', x
end program
You can copy the program from here. If you are trying it out, see what happens if you type either an integer or something that isn't a valid number.

Note that the variable which appears in the read instruction must be a variable which is not declared as a parameter.

It is possible in principle to have more than one variable in the list following read *, but this is not a good idea when numbers are being read from the keyboard. (They can come from elsewhere, as will be discussed later.) Also you should always print out a message to remind the user what to type in.

Subroutine procedure calls

We have already encountered an example of a `procedure call' in the print instruction. we can (loosely) define a procedure call as an instruction which `causes something to happen' other than, or in addition to, simply assigning a value to a variable.

Procedure calls have two main characteristics. Firstly the procedure invoked by the call has a name, e.g. print. Secondly there are usually one or more arguments which specify with what the procedure is to perform whatever it is meant to do.

The name of the procedure is fixed, but the programmer chooses the arguments. Thus: 

 print *, 'Hello'
has as an argument the character string 'Hello' and so that is what gets printed. 
 print *, whatever
has the variable whatever as an argument, so its value gets displayed.

print is a rather special procedure which is built into the language. There are some others like it which you will encounter later. Another type of procedure is the subroutine procedure or simply a subroutine. These are not built into the language but are available in special libraries or can be defined by the programmer.

At this point we will introduce one example of a subroutine procedure which is available in a standard library, and another user defined subroutine which will be supplied to you. Writing your own subroutines will be covered later.

A subroutine procedure call which activates the procedure consists of the word call followed by the name of the procedure and then a list of its arguments (if any) in brackets.

Subroutine procedure: Example 1

A subroutine (we shall use the words `subroutine' and `procedure' interchangeably from now on) named system enables any Unix command (i.e. what might be typed on the keyboard, like ls, edit, cp) to be activated from within a Fortran program. This is not something that you will often want to do (and it should be done with care, rm *.* will destroy all you files!) but it provides a simple example of a subroutine.

system takes a single argument, the name of the Unix command as a character string in quotes, e.g. 'ls'

The following complete program will cause a listing of your files to be printed on the screen when the compiled program is run: 

program list_files
! An minimum example of a subroutine procedure call
! JWP 22/10/98

implicit none

! This instruction makes it work:
call system ('ls')
!
end program files
A version of this program can be copied from here. You can try changing `ls' to `ls -l' and see what happens.

The above versions of the program are the minimum size and complexity required to access this subroutine. If you learn more about the Fortran 90 language you will discover that a stricter and more pedantic form is available. For completeness we have supplied a `long' version of the program here.

Subroutine procedure: Example 2

subroutine atomic
is a set of instructions which will determine, for an element whose name is supplied as a character string (e.g. 'Na'), the atomic number and atomic weight. The instructions are supplied in the form of a subroutine procedure to illustrate further how these are used. This example also serves to show how user written subroutines are put together (though you are not expected to understand this yet, it will be covered later). Additionally, it shows a rather different kind of application for programs from those shown so far. 
subroutine atomic (element, number, value)
   ! look up atomic no. and weight of some common elements
   ! JWP 22/10/98

   implicit none 

   character *(*), intent (in) :: element
   real, intent (out) :: value
   integer, intent(out) :: number
   integer, parameter :: na=12
   integer, dimension(1:na), parameter :: an = &
   (/17, 35, 11, 20, 6, 1, 7, 8, 16, 9, 53, 19/)
   real, dimension(1:na), parameter :: aw = &
   (/35.5, 79.9, 22.0 , 40.1, 12.0 , 1.0, 14.0, 16.0, 32.1, 19.0, 126.9, 39.1/)
   character*48, parameter :: ename = &
   'Cl,Br,Na,Ca,C ,H ,N , O ,S ,F ,I ,K '
   integer :: i
   character *(2) :: copy
   !
   copy = adjustl(element) // '  ' ; copy = copy(1:2) ! fiddle with length etc.
   i = index (ename,copy)  ! look up element in list
   if (i==0) then
      print *, 'I do not have data for element ' // element
      value = 0.0 ; number = 0
      return
   end if
   !
   i = (i-1)/3 + 1  ! Table entry number
   value = aw (i)
   number = an (i)
end subroutine atomic

To make use of this subroutine you must do two things. Firstly, as in the previous example, you must write a program to `call' the subroutine. A typical call might look like this: 

call atomic ('K', Kno, Kwt)
Here the string 'K' determines that it is the properties of potassium which are required. Kno is a variable which you the programmer must have declared as an integer in your program. The subroutine call will set this to the atomic number of potassium. Kwt must have been declared as a real variable and will contain the corresponding atomic weight. Note that these can be any appropriately declared variables and can thus have any names that you want to give them.

To use use this subroutine you also need to let your program know about it. The simplest way to to this is to copy the whole subroutine from the link and put it in the same file as the program you are writing which will use the subroutine. the layout of the file will then be as shown below. 

program atom_stuff
! My own program, declarations etc. go here..
real :: Kwt, ...
....
call atomic ('K', Kno, Kwt)
print *, 'Properties of Potassium..', ...
...
end program atom_stuff


! What follows is just copied..
subroutine atomic (element, number, value)
! look up atomic no. and weight of some common elements
... etc
...
end subroutine atomic
It is not necessary that you understand how atomic actually works, but it is useful to note the existence of variables which will hold character strings. a declaration such as: 
character *(2) :: element
declares a variable called element for storing not numbers but letters, up to a maximum length of two characters.

Intent Attribute

The intent attribute informs a programmer using a subroutine whether a specific argument is to be:

If an argument has intent(in), then it may not be altered in the subroutine. All subroutine arguments should be given an intent attribute to aid clarity in the subroutine use.

Subroutine Structure

The structure of a subroutine is quite similar to that of a program:

For clarity, the subroutine code is usually indented by three character spaces.

Functions

Functions are very similar in structure to subroutines. They take a number of parameters and return a variable to the calling program. Here is a simple function: 
function calc(x)
   !this function evaluates a simple formula
   !Ross Andrews 26/10/99

   implicit none

   real, intent(in) :: x
   real :: calc

   calc = x**2 + 3*x - 4

end function calc
The main point to note here is that the function name, calc is declared as a variable, and assigned to, in the function. This will be the value returned to the calling program. In the calling program, the function calc must be declared as follows: 
real, external :: calc
i.e. external to the main program, and can be used as: 
y = calc(43.2)

Fortran has a number of built-in functions, some of which you have used already:

Any user defined functions should come after the end program prog_name statement, and in the same file in the same manner as subroutines.

Subroutine and Function Variables

The variables in the subroutines and functions (subprograms) described here are independent from the main program variables, and vice versa. What this means in practice is that subprograms cannot refer to main program variables, unless they are passed as an argument, and likewise in the main program.

Subroutines or Functions

There may be some confusion as to whether a subroutine or a function best suits a particular task. Whilst this is down to the preferences of the individual programmer, there are a number of guidelines:

These are only guidelines, not strict rules.
Section L3.4: Reference Material for Fortran 90

Section L3.4 Reference Material for Fortran 90

Note that much of this material is intended for experienced programmers and is rather jargon heavy!