Investigation of Heat Pipes for Drilling Applications

1
A brief introduction to FORTRAN
2
Numerical issues

• What numerical method ?
• What programming language ?
• Discretization of the Equations
• The pressure correction equation and the
  staggered grid
• How to tackle the non-linearities
• Solution of the linear system
• The conjugate gradient (CG), biconjugate gradient
  (biCG), Multigrid techniques, Parallel computing
  Domain decomposition

3

• What numerical method ?
• Finite differences ?
• Finite volume ?
• Finite elements ?

• What programming language ?
• Fortran or C ?

4
FORTRAN or C ? Why FORTRAN?

• There is a vast body of existing FORTRAN code
  (much of which is publicly available and of high
  quality). Numerical codes are particularly
  difficult to port, scientific establishments
  usually do not have large otherwise idle
  programming staffs, etc. so massive recoding into
  any new language is typically resisted quite
  strongly.
• FORTRAN 77 tends to be easier for non-experts to
  learn than C, because its 'mental model of the
  computer' is much simpler. For example, in
  FORTRAN 77 the programmer can generally avoid
  learning about pointers and memory addresses,
  while these are essential in C. More generally,
  in FORTRAN 77 the difference between (C notation)
  x, x, and often even x is basically hidden,
  while in C it's exposed. Consequently, FORTRAN 77
  is a much simpler language for people who are not
  experts at computer internals. Because of this
  relative simplicity, for simple programming tasks
  which fall within its domain, (say writing a
  simple least-squares fitting routine), FORTRAN 77
  generally requires much less computer science
  knowledge of the programmer than C does, and is
  thus much easier to use.
• Fortran 90 changes the picture somewhat, the new
  language is very rich and complex, but you don't
  have to use or even know about all this
  complexity.

5
FORTRAN or C ? Why FORTRAN?
The difference between Fortran and C, is the
difference between a language designed for
numerical computations, and a language designed
for other purposes (system programming).

• FORTRAN ADVANTAGES
• Scientifically oriented
• Better optimized code
• A lot of existing code
• d) Easier to learn
• e) More efficient mathematics
• f) Easier to use and more robust
• g) Better diagnostics

6
The Finite volume method
This method uses the integral form of the
conservation equation. The solution domain is
divided into a finite number of control volumes.
The governing equations are integrated over each
control volume and requires the conservation of
the variables for each of them. This clear
relationship between the numerical algorithm and
the underlying physical conservation principles
is one of the main attraction of this method.
Applying Gauss divergence theorem
7
Discretization of the Equations

• The discretization for one control volume

• After assembling all control volumes, the
  algebraic Equation System

Where ?p is the variable value at the center of
the control volume and ?l is the value of the
variable at the neighbor CVs. Coefficients Ap and
Al involve geometric quantities, fluid properties
and for non-linear equations, the variables
themselves. Qp contains all the terms which do
not contain unknown variables
8
How to tackle the non-linearities
Consider a 1D Conduction problem under steady
state condition
?xPW
?xEP
w
e
W
E
P
?x
And the diffusive term is evaluated as
9
How to tackle the non-linearities
1 - Start with a guess or estimates values of ?
at all grid points 2 - From this guessed ?s,
calculate values of the coefficients in the
discretization equation. 3 - Solve the nominally
linear set of algebraic equationns. To get new
values of ?. 4 - As these ?s values as better
guess, return to step 2 and repeat the process
until further iteration cease to produce a
significant changes in the values of ?.
10
Navier-Stokes equations
Continuity equation
Momentum in the x direction
Momentum in the y direction
Momentum equations in indicial notation
11
The Momentum Equation and the Pressure Correction
Equation

• Are the nonlinearities more difficult in the
  momentum equation?
• The answer is NO, in particular, the convection
  coefficient ?u being a function of the dependent
  variable u of the momentum equation is no
  different from the conductivity being a function
  of the temperature T.
• Why the momentum equation requires a special
  treatment ?

12
The Momentum Equation and ...
WHERE IS THEN THE REAL DIFFICULTY? The real
difficulty lies in the unknown PRESSURE FIELD. We
need to develop an additional equation for the
pressure that lead to proper, converged solution.
That is, we need to construct a formula for the
pressure which is simple a numerical artifice
designed to expedite the convergence of the
velocity field that satisfy continuity Using this
idea, there are different algorithms to achieve
this, like SIMPLE Semi Implicit Method for
Pressure-Linked Equation SIMPLER Semi Implicit
Method for Pressure-Linked Equation
Revised SIMPLEC Semi Implicit Method for
Pressure-Linked Equation Consistent PISO
Pressure Implicit with Splitting of Operators
13
The under relaxation factors
14
Discretization of the U-momentum equation
Using a power law scheme for the convection terms
and a central difference scheme for the diffusion
terms, we obtain
The deferred correction scheme
Momentum equation for Vz , V? and T are obtained
in the same fashion
15
Discretization of the P-correction equation
16
Comments about the solution of the linear system
P is a preconditioning matrix
17
The conjugate gradient (CG)
We transform the problem to a minimization
problem. The CG is applicable only to symmetric
system.
18
The Biconjugate Gradient (BiCG)
The convection/diffusion equation generates
non-symmetric systems
Other variant of BiCG CGS (CG Squared) CGSTAB
(CG Stabilized) GMRES (Generalized Minimum
Residual)
19
Numerical Test non-linear heat conduction in 3D
20
The Grid
Heat input zone
Grid 70 x 70 x 16 78,400 nodes
21
Convergence history for the different numerical
schemes
22
Multigrid techniques
Multigrid is more a strategy to accelerate
convergence than a particular method
23

• For a two grid level
• On the fine grid, perform iterations with a
  method that gives a smooth error
• Compute the residual on the fine grid
• Restrict the residual to the coarse grid
• Perform iterations of the correction equation on
  the coarse grid
• Interpolate the correction to the fine grid
• Update the solution on the fine grid
• Repeat the entire procedure until the residual
  is reduced to the desired level

24
Parallel computing Domain decomposition
Structure of a global coefficient matrix when a
square subdomain is subdivided into 4 subdomains
25
Architecture of the computer program
Pre - Processing
Geometrical and physical parameters from an input
file. Grid generation.
Kernel algorithm CALL CALCT ? CALL
SOLVER CALL CALCU ? CALL SOLVER CALL
CALCV ? CALL SOLVER CALL CALCP
? CALL SOLVER CALL CALCW ?
CALL SOLVER
Processing
Post-Processing
Write results to an output file
26
Basic concepts in FORTRAN
FORTRAN (Formula Translation)
High-level language and low level language
(machine language) Source program (text file)
Compilation
Object program (binary .obj file) Executable
(.exe file) style Syntax error and execution
errors
PROGRAM SUMA ! This is my second fortran
program REAL A, B, C READ (,) A, B C A
B PRINT 5, C 5 FORMAT(1X, 'La suma es ',
F5.2)      END
PROGRAM MYNAME CHARACTER NAME20 !
This is my first fortran program PRINT, 'WHAT
IS YOUR NAME?' READ, NAME PRINT, 'HI, THERE
', NAME END
27
Fortran Evolution first compiler 1957. first
official standard 1972 FORTRAN 66'. updated
in 1980 to FORTRAN 77. updated further in 1991
to Fortran 90. upgraded in 1996 - remove
obsolescent features, correct mistakes and add
limited basket of new facilities. Fortran is
now an ISO/IEC and ANSI standard.
28
Main characteristics of Fortran 90 1) Free
Format in the source code. In Fortran 90, we can
use the format of Fortran 77 or the free format.
If the free format is used, the extension .f90
should be used in the name of the file. 2)
Dynamic assignment. It is possible to use a
dynamic assignment, then all the arrays "work"
independently of their size. 3) User-defined
data. Different types of user-defined data can
be defined, of the same form used in C with
struct. 3) Modules. The modules permit to do a
programming in an oriented objects style similar
to C . The modules can also be used to hide
global variables, then the common sentence fall
in disuse. 4) Functions recursive. Now as a
part of the language. 5) Operations with arrays
built internally. The sentences as A 0 and C
A B are now valid when A, B and C are arrays.
Also, there is a function for the multiplication
of matrices (matmul). 6) Definition of
operators An own meaning for operators can be
defined like and for the user-defined types
of data (objects). 7) Fortran 90 is very
different to the first versions of Fortran. But
guard compatibility with Fortran 77. Fortran 77
has been included as subset of Fortran 90.
29
Introduction The FORTRAN (FORMULATES TRANslator)
is a very powerful programming language. It is
the first scientific language, being developed at
the end of the 50s, but still extensively used
in scientific applications and in engineering.
The last versions of FORTRAN at present are the
versions FORTRAN 90 and 95, that replaced the
version FORTRAN 77 during the 90s. FORTRAN 90
includes new structures and characteristics in
its syntax, though still is compatible with the
previous version. Originally, any calculation
that implied math formulas evaluation had to be
done by means of complicated programs, that
translated those formulas to the language of
computers, in a very primitive way. FORTRAN was
the first high-level language, oriented to
facilitate the things to the user, and to allow
to write programs in an almost identical way as
formulas are written. In FORTRAN the
instructions that the user gives to the computer
are found in a text file (in code ASCII). The
name of the file follows the conventions of the
operating system of the computer. In WINDOWS a
program FORTRAN usually has the extension .FOR',
for example, 'CALCULATE.FOR' in Unix is
traditional the extension is .f', for example,
'calculate.f'. Contrary to BASIC, that is a
interpreted language, FORTRAN is a compiled
language. This means that, once a program is
written, the program should be translated in
block to machine language, which is the language
that the processor of the computer understands,
by means of a process called compilation in
contrast, in an interpreted language, the lines
of the program are translated according to the
flow goes passing, which turns the execution
usually slows.
30
Writing a program if FORTRAN 77
The instructions of the program are written as a
text file in lines. In FORTRAN 77, columns from
the 1 to the 5 are reserved to numbers to
indicates the flow of the program goes to that
line from another part of the program. The
column 6 is reserved to a symbol that indicates
the compiler that that line is a continuation of
the previous line. Therefore, the instructions
of the program begins from the column 7 to 80).
An example of simple program is
Programs must finish with the instruction END.
The instruction that indicates the beginning of
the program, PROGRAM, is not necessary, although
it is recommended to make clear the structure of
the program in the case that subroutines and
functions exist.
31
VARIABLE The first important concept is the
concept of a variable. A variable is a zone
reserved in the memory of the computer, to which
a name is assigned. In this zone we can keep a
number, a sequence of numbers or character
string. The variables in which numbers can be
stored are called numerical, and they can have as
names a combination of characters, though the
first one cannot be a number the number of
characters varies according to the compiler. In
the previous program-example, A, B and C are
variables. Other examples can be SUM, REMAINDER,
RESUL1, etc. The variables are defined at the
beginning of the program (it is convenient,
although not necessary in FORTRAN 77, to define
all the variables of a program at the beginning).
The numerical variables can be integer, real or
complex. We will refer here to the entire and to
the real. The integer variables contain integer
numbers (numbers without a fractional part), for
instance, -5, and they are defined a the
beginning of the program by means of the
instruction INTEGER. The real variables contain
real numbers, (3.14159, 1.457E5, which means
1,457.105(exponential notation)), and they are
defined at the beginning of the program by means
of the REAL instruction. Depending on the
precision we want to represent the numbers, the
real and integer variables can be of simple
precision (four bytes) or of double precision
(eight bytes). This is indicated in the
definition of the variable REAL4 SUMA
                         REAL8 RESTO If the
precision is not indicated, the compiler
understands that is a real variables in simple
precision. Some compilers allow the use of real
variables of quadruple precision, REAL16. For
the integer variables, we have INTEGER4 FIRST
INTEGER8 SECOND
32
Variables in which characters can be stored are
called character variables its name follows the
same conventions of the numerical variables, but
they are defined by means of the instruction
CHARACTER, followed by n, being n the number of
characters that can store in the variable. Thus,
CHARACTER10 NAME Defines a variable called
NAME that can contain up to ten characters. The
characters stored in the character variables can
be any and they should be included among simple
quotation marks, 1', ' The result is', etc.
Arrays of variables can be defined, by means of
components of a vector or a matrix. These
variables have a name and an index. For example,
A(I) is a component of the vector A the variable
I is an integer number and is an index. Also we
can have matrices, B(I,J), with two indices, I
and J. The variable vectors and matrices should
also be declared at the beginning of the program,
with the instruction DIMENSION Thus, we can
write DIMENSION A (10) Which indicates that the
vector A has ten components, A(1), A(2),...,
A(10).
33
INPUT and OUTPUT The input and the output of
data are carried out in FORTRAN through the
instructions READ and WRITE. We can also use the
instruction PRINT to write to the screen. The
syntax of the instructions READ and WRITE
is READ (unidad, formato) lista de variables
WRITE (unidad, formato) lista de variables A
common situation is the one that the data are
introduced from the keyboard, and the data are
written to in the screen, with free format. In
these cases the syntax is READ (,) lista de
variables WRITE (,) lista de variables We
can also use PRINT , lista de variables
   PROGRAM SUMA    REAL A, B, C    READ (,)
A, B    C A B    WRITE (,) 'La suma es
', C    END
34
Intrinsic functions
Numbers and numerical variables can be combined
utilizing mathematical operators as (sum), -
(minus), (multiplication) and / (division) in
the case of division of integer, the division is
the integer division, ex., 3/2 1. Another
mathematical operation of interest is the
exponentiation, which is represented by means of
ex., 3 2 (32). Some functions are
included in FORTRAN. They are intrinsic
functions, whose argument (or arguments) are
included among parenthesis. Intrinsic functions
are the square root, SQRT, trigonometric sine
SIN, cosine, COS, tangent, TAN, arcsine, ASIN,
arccosine ACOS, arctangent ATAN, the exponential
function, EXP and some other. The syntax of
these functions is, ex., SQRT(3.0), SIN(2.3),
EXP (-1.0). If the arguments are in double
precision one must place the letter D in front of
the name of the function. Thus, DSQRT(3.0D0),
DSIN(2.3D0), DEXP (-1D0).
35
Loops
A loop is a structure that permits to perform
repetitions of a group of instructions. The
commands to control a loop are DO-END DO.
Associated to the loop there is a variable,
normally an integer that works as an index. The
general syntax of the loop is DO index initial
value, final value, increment... END DO In each
iteration the index changes from the initial
value to the final value, increasing by the value
of the increment. If the increment doesnt
appear it is considered implicitly that its value
is one. The following example computes a series
of N terms     PROGRAM SUM SUM 0.0
     READ (,) N       DO I 1, N
            SUM SUM 1.0/I2       END DO
     PRINT , SUM      END The way in which
the structure DO-END DO functions is the
following one. In the first iteration the
control variable, in this case I, assumes the
initial value, 1. Instructions inside the loop
are executed, and when the the instruction END DO
is reached the flow is transfer to the beginning
of the loop. The index I takes now the value 2,
and thus the loop is repeated until the value of
I is N, in which the iteration is performed one
more time and the flow continues with the
instructions that come after the END DO.
36
Conditional structures
A condition is a situation in which two variables
are compared and, depending on the result of that
comparison, an instruction or group of
instructions is executed. A conditional structure
has the following structure IF
(condición) THEN        ...        ELSE       
...        END IF The conditional operators
are .EQ. (EQual to) .LT. (Less Than)
.GT. (Greater Than) .LE. (Less to or Equal
than) .GE. (Greater to or Equal than)       
 
PROGRAM  CONDITIONAL      REAL SUM
     INTEGER M, N      READ(,) N,M
          IF (N .LE. M) THEN            
PRINT,N IS LESS OR EQUAL THAN M
ELSE            PRINT,N IS LARGER THAN
M            END IF        END
37
Subroutines and functions
Many times is necessary in a program to perform
an group of instructions in a repetitive way.
Subroutines and functions avoid the writing of
the same piece of code several times. Besides,
they help to make the program more clear and
intelligible. The subroutines and extrinsic
functions are pieces of program that are defined
separately of the main program, and their
variables need to be defined in the subroutine
since are local variables, not shared with the
main program or with other subroutines.
Subroutines are defined with the instruction
SUBROUTINE. A subroutine should end with the
instruction END. The syntax is SUBROUTINE
name (list of arguments)       ...       END
A subroutine is called with the instruction
CALL. This is an example PROGRAM 
POLAR_CARTES       REAL X, Y, RHO, THETA      
READ (,) RHO, THETA       CALL POLARES (X, Y,
RHO, THETA)       PRINT , X, Y       END
      SUBROUTINE POLARES (X, Y, RHO, THETA)
      REAL X, Y, RHO, THETA       X RHO COS
(THETA)       Y RHO SIN (THETA)       END
38
FUNCTIONS
The syntax is very similar to the syntax of a
subroutines but, differently than in a
subroutine, in which there can be any number of
arguments, in a function there is only one a
variable as an output, which is name of the
function itself. Functions are defined by the
instruction FUNCTION, and end also with END. The
general structure is FUNCTION name (list of
arguments)... END As an example PROGRAM
CALCULATE_FACT INTEGER N , FACTORIAL,
FACT READ (,) N FACT FACTORIAL
(N) PRINT , FACT END
FUNCTION FACTORIAL (N) INTEGER I, N,
FACTORIAL FACTORIAL 1 DO I 2, N
FACTORIAL FACTORIAL I
END DO END Note that the name of the
function should be declared as if it were another
variable, in the main program that call the
function as well as in the function itself.
39
Modules

• A module differs from a subprogram in two
  important ways
• a module may contain more than one
  subprogram (called module subprograms)
• a module may contain declaration and
  specification statements which are accessible to
  all program units which use the module.
• Modules are also compiled separately.
• The general form of a module is
• MODULE name
• declaration statements
• CONTAINS
• module subprograms
• END MODULE name

MODULE MyUtils REAL, PARAMETER Pi
3.1415927 CONTAINS SUBROUTINE SWOP( X, Y )
REAL Temp, X, Y Temp X X Y Y
Temp END SUBROUTINE SWOP END MODULE MyUtils
USE MyUtils IMPLICIT NONE REAL A, B READ, A B
Pi CALL SWOP( A, B ) PRINT, A, B END