Numerical methods
Annie Fournier-Gagnoud
Phelma 2016-2017
OUTLINE
● Introduction
● Programming with Matlab
Introduction
Algorithm
Recommandations
Matlab software
● Ordinary Differential Equations ODE
● Partial Differential Equations PDE
A. Fournier-Gagnoud Numerical Methods 2016-2017 2
Introduction
To develop a program few steps are useful:
Choice of language,
Design the algorithm,
Construct the program,
Test the program step by step (sequence by
sequence).
A. Fournier-Gagnoud Numerical Methods 2016-2017 3
OUTLINE
● Introduction
● Matlab software
Introduction
Algorithm
Recommendations
Matlab
● Ordinary Differential Equations ODE
● Partial Differential Equations PDE
A. Fournier-Gagnoud Numerical Methods 2016-2017 4
ALGORITHM
Algorithm = sequence of finite instructions, used
for calculation and data processing
Algorithms are essential in computational
process
Algorithm must be rigorously defined
Algorithm is a precise list of precise steps
Order of instructions must to be listed explicitly
A. Fournier-Gagnoud Numerical Methods 2016-2017 5
ALGORITHM DESIGN
Algorithm design is a specific method to create a
mathematical process in solving problems.
In a first step we must analyze the program to
solve a given problem.
We have to identify the input data, the output
data, to choice the variables.
A. Fournier-Gagnoud Numerical Methods 2016-2017 6
….
Initialization of input data …
ALGORITHM DESIGN
Calculations
…..
….. Calculation 1
.….
Calculation 2
…
Saving results ….
…
Graphic representation
A. Fournier-Gagnoud Numerical Methods 2016-2017 7
ALGORITHM DESIGN
Main algorithm Function / subroutine
Main steps
of the calculation. Input parameters
Main program Output parameters
calls functions Successive Steps
or subroutines
A. Fournier-Gagnoud Numerical Methods 2016-2017 8
ALGORITHM DESIGN
Iterative algorithms: repetitive constructs like loops to solve
the given problems
An algorithm may be viewed as controlled logical deduction:
logical instruction
Serial Algorithms: computers execute one instruction of an
algorithm at a time
Parallel algorithms: take advantage of computer architectures
where several processors can work on a problem at the same
time
A. Fournier-Gagnoud Numerical Methods 2016-2017 9
OUTLINE
● Introduction
● Programming with Matlab
Introduction
Algorithm
Recommendations
Matlab software
● Ordinary Differential Equations ODE
● Partial Differential Equations PDE
A. Fournier-Gagnoud Numerical Methods 2016-2017 10
RECOMMENDATIONS
Planning the program
a series of smaller, independent tasks
each task as a separate function
General Coding Practices
variable names: easier to understand
< 80th column
Importance of comments %
Coding in Steps
Making Modifications in Steps
Testing the Final Program
A. Fournier-Gagnoud Numerical Methods 2016-2017 11
OUTLINE
● Introduction
● Programming with Matlab
Introduction
Algorithm
Recommendations
Matlab software
● Ordinary Differential Equations ODE
● Partial Differential Equations PDE
A. Fournier-Gagnoud Numerical Methods 2016-2017 12
MATLAB
MATrix LABoratory
Commercial software developed by the
MathWorks Inc.
Numerical calculations on vectors, matrices,
real or complex
Freeware : freemat, octave, scilab
Unix, Linux, Windows implementation
A. Fournier-Gagnoud Numerical Methods 2016-2017 13
MATLAB
Interactive software system for numerical
computations and graphics
Matlab is designed for Matrix computations:
Solving systems of linear equations
Computing eigenvalues - eigenvectors
Factoring matrices
A. Fournier-Gagnoud Numerical Methods 2016-2017 14
Beginning with Matlab
Start : double clicking the MATLAB shortcut
Desktop MATLAB appears, Including
Some menus
Help is very useful Getting Started
Three windows
Command window : Enter Matlab statements at the prompt
Showing you what you type and the result produced
At this level = Matlab is a pocket calculator
Files in the directory
A. Fournier-Gagnoud Numerical Methods 2016-2017 15
Beginning with Matlab
Command window
Variable Window
Click on help/MATLAB help
To reach on line help
Commands history
A. Fournier-Gagnoud Numerical Methods 2016-2017 16
MATLAB Help
Getting started
Introduction
Matrices and Array
Graphics
Programming
…….
A. Fournier-Gagnoud Numerical Methods 2016-2017 17
Matlab : variables
Variables are not declared prior
Easy – convenient
Careful users
Good use of Malab is a difficult task.
Basic data type : n-dimensional array of
double precision numbers
A. Fournier-Gagnoud Numerical Methods 2016-2017 18
Interactive system : vectors
>>a=[1 4 8 ]
a=
1 4 8
>>b=[3 ; 2 ]
b=
3
2
>>c=pi
A. Fournier-Gagnoud Numerical Methods 2016-2017 19
Matrix definition
>>A=[1 1 1; 2 4 8; 3 9 27 ]
A=
1. 1. 1. Matrix definition
2. 4. 8.
Each componant is separated by space
3. 9. 27.
Each line by ;
>>A' '
1. 2. 3.
Transposed Matrix
1. 4. 9. ' is the tranposition operator
1. 8. 27.
A. Fournier-Gagnoud Numerical Methods 2016-2017 20
Calculation: vectors-matrices
>>a=[1 4 8 ]
>>b=[3 ; 2 ;5]
>>a*c
??? Undefined function or variable ‘c’
>>a*b
ans=54
>>A=b*a
A=
3 12 24 3 3
2 8 16 a 1 4 8 b 2 A b a 2 1 4 8
5 20 40
5 5
A. Fournier-Gagnoud Numerical Methods 2016-2017 21
Calculation:
Matrix – Vector Product
>>b=[1 0 0 ]; No echo if ;
>>A*b
??? Error using ==> * |1 0 0|
Inner matrix dimensions must agree. | 1. 1. 1. |
>>A*b' | 2. 4. 8. |
1.
2.
| 3. 9. 27. | = ?
3.
Arithmetic operators + - * /
A. Fournier-Gagnoud Numerical Methods 2016-2017 22
Linear system: solver
--> b= [1 0 0]’;
transposed
--> x = A \ b
Resolution of A*x = b
x =
3.0 Verification that it runs !
- 2.5
0.5
A. Fournier-Gagnoud Numerical Methods 2016-2017 23
Matrix and Vectors building
[] [ ]
1 1 2 3
Defining a vector V: V = −2 a Matrix A:A= 2 −3 5
3 4 −3 1
• Display the 2nd component of V
• Calculate norm of V (3 methods)
• Build the vector W such as W=A . V and G=V ^ W
•Display the 2nd line, 3rd column element of A matrix
•Display first column of A
•Display 2nd line of A
•Calculate determinant of A
•Solve A X = V (2 methods)
A. Fournier-Gagnoud Numerical Methods 2016-2017 24
Matrix and Vectors building
>> A=[1 2 3;2 -3 5;4 -3 1]
>> V=[1 -2 3]’ or V=[1;-2;3]
>> A(2,3)
>> V(2) ans =
ans = 5
-2
>> A(:,1)
>> norm(V) ans =
>> sqrt(V’*V) 1
>> sqrt(dot(V,V)) 2
ans = 4
3.7417
>> A(2,:)
ans =
2 -3 5
A. Fournier-Gagnoud Numerical Methods 2016-2017 25
Matrix and Vectors building
>> W=A*V
W=
6
23
13
>> G=cross(V, W)
-30
-6
6
A. Fournier-Gagnoud Numerical Methods 2016-2017 26
Matrix and Vectors building
Build matrix T
Calculate with Matlab the dimension of T matrix
List of Matlab functions : zeros, ones, diag, eyes, size
2 -1 0 0 .... 0
-1 2 -1 0 0
0 -1 2 -1 0
T= 0 0
0 -1
0 ... .. -1 2
A. Fournier-Gagnoud Numerical Methods 2016-2017 27
First program « Hello World ! »
Program edition : File New M-file hello.m
Saving : File Save , write hello.m
Matlab code interpretation with editor
debug menu
Or by typing the name of the program file,
hello.m in the Matlab command window
A. Fournier-Gagnoud Numerical Methods 2016-2017 28
Matlab variables
All variables in Matlab are tables (matrices),
Matlab is case sensitive
variable names : I non equal to i
Auto-declaration of variables when initialized
Example : Matlab comments
t=0; % numeric variable
s=‘Hello World’; % table of characters variable
Type Overloading of a variable is possible
A. Fournier-Gagnoud Numerical Methods 2016-2017 29
Matlab variables
A. Fournier-Gagnoud Numerical Methods 2016-2017 30
Error ! Ok !
Look help informations for
zeros, ones et diag
A. Fournier-Gagnoud Numerical Methods 2016-2017 31
Variable: influence zone in program
Initialization of the variable N
Influence : from the
Initialization until end of file
No local declaration inside
! a bloc !!!
Danger : ‘:’ and , in a ‘for Loop’
have not the same effects !
for i=1, 10
i
end;
A. Fournier-Gagnoud Numerical Methods 2016-2017 32
Binary operators
binary operators
A<B
A>B
Operators : A && B
A <= B
A || B
A >= B
Arithmetic operators A == B
Logical operators A ~= B
Relational operators
A. Fournier-Gagnoud Numerical Methods 2016-2017 33
Control structures
compare.m
Condition :
if
Enter
true i & j
i<j ?
false Block of instructions
A. Fournier-Gagnoud Numerical Methods 2016-2017 34
When different instructions
Conditional connection have to be executed according
tochoice2.m
a variable value
choice1.m In replacement of
‘if, elseif, else, end’
control structure
no break instruction in a
block case as in C++!
A. Fournier-Gagnoud Numerical Methods 2016-2017 35
Loops
loop1.m
loop2.m
A. Fournier-Gagnoud Numerical Methods 2016-2017 36
Loops
loop3.m sum1.m
What is calculated ?
Break : loop breaking
A. Fournier-Gagnoud Numerical Methods 2016-2017 37
Mathematical functions
|x| abs(x)
arctg x atan(x)
cos x cos(x)
exp x exp(x)
ln x log(x)
log x log10(x)
x sqrt(x)
y
C&
x x^y C++
A. Fournier-Gagnoud Numerical Methods 2016-2017 38
Mathematical functions
: continuation following line
num2str : numeric conversion string
i ≡ sqrt(-1.)
A. Fournier-Gagnoud Numerical Methods 2016-2017 39
Writing a function
% Solving a linear system
% comments
function [x]=resolution(A, b) Entries are separated by a comma
x = A\b;
Outputs separated by a comma
• Save in a file with “resolution.m” name for example
• In matlab, select “resolution.m” with the menu : file open
• Call the resolution function from the Matlab command window :
x = resolution(A, b)
A. Fournier-Gagnoud Numerical Methods 2016-2017 40
Writing a function
function linsolver linsolver is the principal function
A=[1 2 3;2 -3 5;4 -3 1]
calling the resolution function
b=[1 -2 3]'
[x, ier]=resolution(A, b)
if (ier)
error('Pb resolution');
end;
% solving a linear system Conditional instruction
function [x, res]=resolution(A, b) if (condition)
if (det(A)==0) bloc;
res=1; determinant else
x=0; bloc;
else end
res=0;
x=A\b; What happens if det(A)=0 ?
end
A. Fournier-Gagnoud Numerical Methods 2016-2017 41
Graphical functions
For x[-50,50], calculate f=x*sin(x) and draw f(x)
List of Matlab functions : plot, xlabel, ylabel, title, grid on
x = linspace(-50,50,100);
y = x*sin(x);
plot(x,y)
hold on
title('Graph of f=x*Sine')
xlabel(' x ') % x-axis label
ylabel('f(x)=x*sin(x)') % y-axis label
A. Fournier-Gagnoud Numerical Methods 2016-2017 42
Graphical functions
For [0,2] et [0, ] calculate e=sin()^4*cos(2*)
Present e(, ) with a surface and isovalues of e(, ) as contour lines
List of Matlab functions : meshgrid, surf, contour
% Isovalues
clc
close all figure
colormap(jet)
[phi,tet]=meshgrid(0:0.02:2*pi,0:0.01:pi); [c,h]=contour(phi,tet,e,v);
e=sin(tet).^4.*cos(2*phi); set(h,'ShowText','on','TextStep',...
get(h,'LevelStep')*2);
% Fixed Isovalues
v=-1.:0.2:1; figure
colormap(jet)
colormap([0 0 0]) h=surf(phi,tet,e);
[c,h]=contour(phi,tet,e,v); shading interp
set(h,'ShowText','on','TextStep',...
get(h,'LevelStep')*2);
A. Fournier-Gagnoud Numerical Methods 2016-2017 43
Simplified Man-Machine-Interface
Example n°1
choice_list ={'banana', 'pear', 'apple', 'cherry'}
[rg, ok] = listdlg('PromptString','Select a fruit ?',...
'SelectionMode','single',...
'ListString', choice_list,...
'InitialValue',2,...
'listSize',[400 80],...
'Name','dialogue');
choice_list{rg}
ok % =1 if OK, =0 if CANCEL
A. Fournier-Gagnoud Numerical Methods 2016-2017 44
Simplified Man-Machine-interface
Example n°2
nb_lignes=1;
val_defaut={'293'};
txt = inputdlg('temperature' ,‘window_title, ...
nb_lignes, val_defaut);
T=str2double(txt) % convert a string into a real
A. Fournier-Gagnoud Numerical Methods 2016-2017 45
Structures
Point in 2D coordinates
M.x=input('abscisse du point : ');
x M.y=input('ordonnee du point : ');
y disp(M);
A. Fournier-Gagnoud Numerical Methods 2016-2017 46
M.x=input('abscisse du point : ');
M.y=input('ordonnee du point : ');
Structures
disp(M);
Point in 2D coordinates
Automatic declaration and
Dynamic memory allocation
of a table composed of point
structures
A. Fournier-Gagnoud Numerical Methods 2016-2017 47
Reading a Point table in a file
Objective: filling a data structure ‘tab’ from the data file ‘tab.dat’
Number of points ‘ tab.dat’ t=
63
Coordinates_of_points NP: 63
1 -3.1415927e+000 -1.2246468e-016 x: [1x63 double]
2 -3.0415927e+000 -9.9833417e-002
y: [1x63 double]
3 -2.9415927e+000 -1.9866933e-001
4 -2.8415927e+000 -2.9552021e-001
5 -2.7415927e+000 -3.8941834e-001
6 -2.6415927e+000 -4.7942554e-001
7 -2.5415927e+000 -5.6464247e-001
8 -2.4415927e+000 -6.4421769e-001
9 -2.3415927e+000 -7.1735609e-001
…
63 3.0584073e+000 8.3089403e-002
A. Fournier-Gagnoud Numerical Methods 2016-2017 48
function [tab, err]=reading(nomfich)
clear tab;
fich = fopen(nomfich);
if ( fich == -1 )
errordlg([‘Invalid file name : ' nomfich]); % !!does not stop
err = 1;
else
% reading the comment
comment = fscanf(fich, '%s', [1])
% reading the number of points
tab.NP = fscanf(fich,'%d', [1]) % [1] only one scalar value
% reading the comment function tab_point1
comment = fscanf(fich, '%s', [1]) clc
nom='tab.dat'
% reading point coordinates [t, ier]=reading(nom);
for np=1:tab.NP if (ier==1)
n = fscanf(fich,'%d', [1]) error(‘reading error);
end;
x = fscanf(fich,'%g', [1])
y = fscanf(fich,'%g', [1])
plot(t.x, t.y)
tab.x(n)=x;
tab.y(n)=y; t=
end; NP: 63
err = 0;
end; x: [1x63 double]
A. Fournier-Gagnoud Numerical Methods 2016-2017 49
fclose(fich); y: [1x63 double]
Reading a Point table in a file
Objective :
Filling a data structure ‘curve’ reading the data file ‘tab.dat’
curve= curve.NP = 63
NP: 63
M: [1x63 struct] curve.M =
1x63 struct array with fields:
x
y
curve.M(i) donne accès au ième point M(i) :
for ex. curve.M(1) = curve.M(1).x =
x: -3.1416
y: -1.2246e-016 -3.1416
A. Fournier-Gagnoud Numerical Methods 2016-2017 50
function [curve, err]=reading(nomfich)
clear curve;
fich = fopen(nomfich);
if ( fich == -1 )
errordlg([‘incorrect file name : ' nomfich]); % !! Does not stop
err = 1;
else
% reading a comment
comment = fscanf(fich, '%s', [1])
% reading the number of points =
curve.NP = fscanf(fich,'%d', [1]) % [1] only one scalar value
% reading a comment
comment = fscanf(fich, '%s', [1])
function tab_point2
% reading points coordinates clc
for np=1:curve.NP nom='tab.dat'
n = fscanf(fich,'%d', [1]) [curve, ier]=reading(nom);
x = fscanf(fich,'%g', [1]) if (ier==1)
y = fscanf(fich,'%g', [1]) error(‘error reading’);
curve.M(n).x=x; end;
curve.M(n).y=y;
end; for i=1:curve.NP
err = 0; x(i)=curve.M(i).x;
end; y(i)=curve.M(i).y;
fclose(fich); end;
A. Fournier-Gagnoud plot(x, y)
Numerical Methods 2016-2017 51
MATLAB
Website:
http://www.mathworks.com/academia/student_center/t
utorials/starting_matlab.html
On this web site read the different tutorials
and stop after: Creating scripts with
MATLAB Editor/Debugger
Two important points:
Visualizingdata
Creating scripts with MATLAB Editor/Debugger
A. Fournier-Gagnoud Numerical Methods 2016-2017 52
MATLAB
Indocumentation MATLAB
Getting started
IntroductionDocumentationProgramming
fundamental Functions and scripts
In this part the 5 sequences are very
important:
Program Development
Working with M-files
M-files Scripts and Functions
Calling Functions
Function Arguments
A. Fournier-Gagnoud Numerical Methods 2016-2017 53
