SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
én Sevilla
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Description of the module
▪ Lectures
• Lectures on Mondays
• FIRST THREE WEEKS: 10.00h-12.00h
• REST OF THE SEMESTER: 10.00h-11.00h
• PC Labs on Tuesday 16.00h-18.00h
• Students 689000 to 954000: Engineering Central C109, PC Room
• Students 954000 to 980000: Engineering Central C104, PC Room
• 1 demonstrator per room (Photos will be made available in Blackboard so that you can recognise the demonstrator in your room)
• Office Hours on Mondays 12.00h-13.00h
• Room A128, Engineering Central
• Matlab support
• PC Labs. This is a unique opportunity to get help from the demonstrators
• Math and CAE cafes. http://www.swansea.ac.uk/engineering/courses/mathematics-resources/
2
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Content of the module
▪ Syllabus
• Introduction to Matlab
• Introduction to programming 1: functions and scripts
• Introduction to programming 2: control flow and loops
• Problem solving using numerical methods in Matlab
• Interpolation and approximation
• Root finding
• Numericalintegration
• Ordinary differential equations ▪ Learning outcomes
• At the end of the semester you must be able to use Matlab to solve problems that you will be given in Level 3 modules and the project
• To gain knowledge in using Matlab as a numerical calculator to solve engineering problems
• To gain an ability to design and test computer codes for solving engineering problems
3
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Material
▪ Course material (in Blackboard)
• Lecture notes (plus extra Matlab resources, videos, etc)
• PC Labs notes
• Model answers
▪ Books
• More than 40 books of Matlab in the library
• One of them is recommended in Blackboard
▪ Online Resources
• Go to Blackboard→College of Engineering Mathematics Resources • Videos
• Manuals • Notes
• Examples • etc.
4
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Assessment
▪ Group Coursework 1 (15%). Submission on week 25 • Show knowledge on basic Matlab programming concepts
▪ Group Coursework 2 (25%). Submission on week 29
• Show knowledge on how to use numerical methods within Matlab
▪ Class test (60%). Week 29 (before Easter break) • Class test will take place in the PC Lab
Assessment calendar in College Intranet
5
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Academic Misconduct
▪ Plagiarism
• Using, without acknowledgment, another person’s work and
submitting it for assessment as though it were one’s own work
▪ Collusion
• Two or more students or other persons working together without prior
authorisation
▪ Commissioning
• Paying for, or arranging for another to produce, a piece of work which is then submitted for assessment
6
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Academic Misconduct
▪ AM is on the rise and I take it extremely serious
• Last 3 years in EG-228 only
• ~50 students penalised including 5 expelled from the University
• If you don’t feel confident with Matlab
• Work hard from day 1
• Come to all PC labs
• Useallavailableresources
• Office hours
• Demonstrators • Math Cafés
• And think about your future
7
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Assessment
▪ Group Coursework
• Defining the groups
• The maximum number of group members is 3
• A list of empty groups is available in Blackboard
• Self-enrol before the end of teaching week 3
• At the end of week 3, self-enrolment will automatically close
• If you do not have a group by then, you will automatically be joined into a group and no changes will be allowed
• If you are in a group of 1 or 2 members, your group might be completed with other members
• Submission
• One submission per group. ONLY electronic submissions
• All students are responsible of the submission
• Further instructions will be provided in due course
8
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Assessment
▪ Group Coursework
• All coursework questions will be answered ONLY using a Blackboard
forum (not by email) so everybody can see Q&A
• You can use the demonstrators and the Math Cafes to
• Clarify questions related to the use of Matlab
• Solve questions about the list of exercises proposed in the PC lab notes
• Demonstrators do not know the marking scheme. Following their advice to complete the coursework does not guarantee full marks. You should not ask them to help you do the coursework in the first place
• You can use office hours to clarify questions regarding the
• Theoretical content of the module
• The use of Matlab
• The list of exercises proposed in the PC lab notes
• I will not answer if a solution is correct/incorrect before submitting (will you ask this in an exam?!)
• Each coursework is must reflect your own work
• If during the marking any doubts arise, the whole group can be called for an interview that will be used to decide the final mark
9
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
én Sevilla
PART 1 Introduction to Matlab Using Matlab as a calculator
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ The Matlab environment
Help
Current folder
Command window
Command prompt >>
EDU>
Workspace
Command History
11
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Matlab as a numerical calculator
• Basic arithmetic operations • + Addition
• – Subtraction
• * Multiplication • / Division
• Hierarchy of arithmetic operations
1. Content of brackets is evaluated starting from the innermost bracket
2. Exponentials, from left to right (The symbol for exponentiation is ^)
3. Multiplications and divisions, from left to right
4. Additions and subtractions, from left to right
12
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Definition of variables
• To create a variable assign a value
to a name using =
>> r = 2.131
• Use ; at the end to suppress the output
>> r = 2.131;
• The name of variables is case sensitive
>> r = 2.131; >> R = 90.121;
• Variables can be created using expressions involving other variables: >> L = 2*pi*r;
>> A = pi*r^2
• The variables appear in the workspace.
• Note that “pi” is a reserved word. There is no need to define it
13
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Matlab built-in functions
• A great number of functions is available
• Trigonometric:
• sin, cos, tan, … (angle in radians)
• sind, cosd, tand, … (angle in degrees)
• asin, acos, sinh, …
• Exponential:
• exp
• Logarithms
• log, log2, log10
• Squareroot • sqrt
• Rounding:
• ceil, floor, round,…
• Check the Matlab help for these functions
>> sin(pi/2) >> cosd(45) >> asin(-1)
>> exp(2.1) >> log10(100) >> sqrt(9.5) >> ceil(7.1)
14
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Use Maltab help!
• Extensive use of Matlab product help is strongly recommended, specially at the early stages of your learning!
15
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Vectors
• Row vector: using square brackets and separating the components by commas or spaces.
• Example: w = [0, 7.8 9, -1.1];
• Column vector: using square brackets and separating the components by semicolon
• Example: >> v = [-6.2; 5.1; 0.1];
16
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Vectors
• 1.
Ways to build large vectors without writing all the components: Function linspace:
>> v = linspace(0, 2, 11); Row vector with 11 equally
spaced points between 0 and 2 Syntax initial:step:final
>> w = [0: 0.4 : 1.2];
Row vector starting at 0 and,
increasing by 0.4, ending at 1.2 Function zeros:
>> x = zeros(1,4);
Vector of 1 row and 4 columns fill in with zeros
2.
3.
17
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Vectors
• 4.
5.
Ways to build large vectors without writing all the components: Function ones:
>> y = ones(3,1) ;
Vector of 3 rows and 1 column fill in with ones
Concatenate vectors of the same type (row or column)
Row vectors:
>> v = [1, 8, 4];
>> w = [7, 9, 1, -4]; >> z = [v, w];
Column vectors: >> v = [1; 4];
>> w = [7; 9; -4]; >> z = [v; w];
18
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Operations with vectors
• Addition and subtraction of vectors of the same type (row or column) and length
>> v = [8, 2, 9]; >> w = [-2, 0.1, 7]; >> v + w
• Scalar product of a row vector and a column vector >> v = [6, 3, 1];
>> w = [4; 2; 0]; >> v*w
• Functions applied to vectors >> v = [30 60 90];
>> w = sind(v)
This means w = [ sind(30), sind(60), sind(90)];
19
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Operations with vectors – Usual mistakes
• Addition and subtraction with vectors of different type (row or column) and/or length
>> v = [6, 3, 1];
>> w = [4; 2; 0];
>> x = [-1; 0; 10; 3]; >> v + w
>> w + x
• Scalar product of non matching vectors >> v = [6, 3 1];
>> w = [4 2 0];
>> v*w
>> v*transpose(w)
• Note the use of the function transpose
• When something is not working: read the error message!
20
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Basic 2D plots
Linear plot using two vectors
• Plot the function cos(x) for the angle x between [0,720] degrees
>> x = [0:45:720]; >> y = cosd(x); >> plot(x,y)
• Increase the number of sample points to get more “resolution”.
>> x = [0:5:720]; >> y = cosd(x); >> plot(x,y)
21
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Basic 2D plots
Improving the look of the plots by
changing line and colours >>plot(x,y,’LineWidth’,4)
4 is the thickness of the line
>>plot(x,y, ‘r–o’)
r indicates red colour,
— indicates a discontinuous line o indicates symbol=circle
22
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Basic 2D plots
Adding information
>> title(‘Function plots’) >> xlabel(‘x (angle)’) >> ylabel(‘y’)
>> grid on
Several plots in one figure: >> hold on
>> z = sind(x);
>> plot(x,z, ‘k-s’)
>> legend(‘cosinus’,’sinus’)
“hold on” Allows you to plot data in an existing figure without losing the previous plot. If not, only the last plot is shown
23
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ The name of the vector followed by the component/s to access
in between brackets
>> v = [10, -12, 35, 7.5, 9 , 12, 4]; >> v(3)
▪ The result is a scalar
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
24
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ The name of the vector followed by the component/s to access
in between brackets
>> v = [10, -12, 35, 7.5, 9 , 12, 4]; >> v(3)
▪ We can access several components using the semicolon >> v(1:5)
▪ The result is a row vector
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
25
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ The name of the vector followed by the component/s to access
in between brackets
>> v = [10, -12, 35, 7.5, 9 , 12, 4]; >> v(3)
▪ We can access several components using the semicolon >> v(1:5)
>> v(1:2:5)
▪ The result is a row vector
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
26
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ The name of the vector followed by the component/s to access
in between brackets
>> v = [10, -12, 35, 7.5, 9 , 12, 4]; >> v(3)
▪ We can access several components using the semicolon >> v(1:5)
>> v(5:-2:1)
▪ The result is a row vector
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
27
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ The name of the vector followed by the component/s to access
in between brackets
>> v = [10, -12, 35, 7.5, 9 , 12, 4]; >> v(3)
▪ We can access several components using the semicolon >> v(1:5)
>> v(1:2:5), v(5:-2:1)
▪ Even if the components don’t follow a pattern >> v([1,4,6])
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
28
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ But sometimes, we don’t know a priori where the information is • Obtain the maximum value in a vector
>> v = [10, -12, 35, 7.5, 9 , 12, 4]; >> a = max(v)
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
29
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from vectors
▪ But sometimes we don’t know a priori where the information is
• Obtain the maximum value in a vector
>> v = [10, -12, 35, 7.5, 9 , 12, 4];
>> a = max(v)
• But where is the maximum?
• >> [a, position] = max(v)
1
2
3
4
5
6
7
10
-12
35
7.5
9
12
4
30
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Matrices
• Matrices defined element by element
• The elements in between square brackets
• The comma (or space) separates elements in the same row
• The semicolon separates rows
>> A = [1, 2, 6; -1 0 4; 3 7 9];
• Matrices defined by multiplying
a column vector by a row vector
>> v = [1; 5; 2]; >> w = [4, 3]; >> B = v*w
31
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Matrices
• Matrices defined using special functions:
• Functionzeros:
>> A = zeros(3,4)
Matrix of dimension 3×4 filled with zeros
• Functionones:
>> B = ones(8,10);
Matrix of dimension 8×10 filled with ones
• Functioneye:
>> C = eye(3,3)
Identity matrix of dimension 3×3
• If you call zeros, ones and eye with one number only it produces a square matrix
>> D = ones(4)
Matrix of dimension 4×4 filled with ones
32
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Matrices
• Matrices defined concatenating existing vectors and/or matrices
• The dimensions must be consistent
>> v = [1 2]; >> w = [-1 0]; >> A = [v; w]
>> v = [7; 2]; >> w = [0; 5]; >> B = [v, w]
>> v = [1 2];
>> w = [-1 0];
>> x = [3, 7];
>> y = [6; 4; 9];
>> A = [[v; w; x], y]
33
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Matrices
• Matrices defined concatenating existing vectors and/or matrices
• The dimensions must be consistent
>> A = [1 0; 7 2];
>> B = [-1 3; 4 2; 9 6]; >> C = [A; B]
>> A = [1 0; 7 2];
>> B = [-1 3; 4 2; 9 6];
>> C = [[A; B], [transpose(B); eye(3)]]
34
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ The name of the matrix followed by the components to access
in between brackets and separated by a comma
▪ First index indicate row, second index indicate column >> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6];
>> A(2 , 3)
▪ The result is a scalar
1
2
3
4
1
10
21
13
1
2
0
35
4
-8
3
-1
7
9
2
4
-4
4
5
6
35
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ The name of the matrix followed by the components to access
in between brackets and separated by a comma
▪ We can access to several components using the colon >> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6];
>> A(2 , : )
▪ The result is a row vector
Colon means “ALL COLUMNS”
1
2
3
4
1
10
21
13
1
2
0
35
4
-8
3
-1
7
9
2
4
-4
4
5
6
36
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ The name of the matrix followed by the components to access
in between brackets and separated by a comma
▪ We can access to several components using the colon >> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6];
>> A(2 , 2:4)
▪ The result is a row vector
1
2
3
4
1
10
21
13
1
2
0
35
4
-8
3
-1
7
9
2
4
-4
4
5
6
37
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ The name of the matrix followed by the components to access
in between brackets and separated by a comma
▪ We can access to several components using the colon >> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6];
>> A([1,3] , 4)
▪ The result is a column vector
1
2
3
4
1
10
21
13
1
2
0
35
4
-8
3
-1
7
9
2
4
-4
4
5
6
38
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ The name of the matrix followed by the components to access
in between brackets and separated by a comma
▪ We can access to several components using the colon >> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6];
>> A([1,4] , 3:4)
▪ The result is a matrix
1
2
3
4
1
10
21
13
1
2
0
35
4
-8
3
-1
7
9
2
4
-4
4
5
6
39
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ But sometimes we don’t know a priori where the information is • Extract the maximum value in a matrix
>> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6]; >> z = max(A)
Result: z = [10 35 13 6];
▪ The result is a vector containing the maximum of each column operating along all rows
1
2
3
4
1
10
21
13
1
2
0
35
4
-8
3
-1
7
9
2
4
-4
4
5
6
40
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Accessing information from matrices
▪ But sometimes we don’t know a priori where the information is • Extract the maximum value in a matrix
>> A = [10, 21, 13, 1; 0, 35, 4, -8; -1 7 9 2; -4 4 5 6]; >> z = max(A)
Result: z = [10 35 13 6];
▪ To obtain the maximum of the entire matrix we need to use, again, the function max
>> max(z)
▪ Or, in a single instruction, >> max(max(A))
41
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Matlab
▪ Solving linear systems
▪ The slash and backslash operators
>> A = [1, 2, 6; -1 0 4; 3 7 9] >> b = [2; 1; 0];
>> x = A\b
• Solves the system Ax = b using the most appropriate method (equivalent to x= A-1*b without inverting A, which is very inefficient)
>> A = [1, 2, 6; -1 0 4; 3 7 9] >> B = [4, 7, 2; 8 1 1; 0 1 -1]; >> A\B
• Returns A-1*B • >>A/B
• Returns A*B-1 (by solving linear systems)
42
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
én Sevilla
PART 1 Introduction to Programming 1- Functions and scripts
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Outline
• The Matlab editor
• Scripts
• How to create a script?
• How to execute a script?
• Functions
• Inputs (arguments)
• Outputs (results)
• How to create a function?
• How to execute a function?
• Tips for designing good programs
• Presenting the results on the screen
44
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ The Matlab editor
For simple computations it is possible to work in the command
window but, in general we must use scripts and functions ▪ Script
• Matlab M-file containing a sequence of Matlab instructions ▪ Function
• Matlab M-file capable of taking input variables, operating with them and returning the result as outputs
• Syntax:
• Examples:
• >> x = sind(30)
x is the output, sind is the name of the function, 30 is the input
• >> v = linspace(0,1,10);
0, 1 and 10 are the three inputs and the vector v is the output
45
function [output1, output2,…]=functionName(input1, input2,…)
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ To create a script
• Create a new M-file by choosing “New script”
• The Matlab editor will open (blank page)
Write a set of Matlab instructions that will be executed in order
• Save the M-file by choosing “Save” in the menu and name the file
The extension of the file must be “.m”
• Use comments (lines starting with %) in order to describe what is the program doing
SAVE THE FILE IN YOUR WORKING DIRECTORY!
46
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
• • •
▪ Example of a script
Comments (lines starting with %), in green, are not Matlab instructions
Note that the semicolon has been forgotten in lines 15 and 18.
This must be avoided!
▪ The Matlab editor is interactive. It uses colours to highlight Warnings and Errors
▪ In fact, forgetting the semicolon at the end of a instruction is highlighted with a warning colour on the equal symbol of lines 15 and 18
47
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Script execution
• Once the file is saved, you can execute a script by
• Writing the name in the command window (with no extension)
• Clicking “Run” in the Matlab editor
48
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ To create a function
• Create a new M-file by choosing “New script”**
• The Matlab editor will open (blank page)
• Write the function header •
Starts always with the word “function” that will be highlighted in blue in the Matlab editor
• Write a set of Matlab instructions to be
executed in order.
The instructions are designed to obtain the outputs using the inputs
• Save the M-file by choosing “Save” and name as functionName.m
**Note that you can also create a Matlab function by choosing “New Function”. The function header will be created but you still have to name
the function and provide the inputs and outputs
function outputs = functionName(inputs)
THE NAME OF THE M-FILE MUST BE EQUAL TO THE NAME OF THE FUNCTION
49
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Function example – One input and one output
• A function that, given the radius of a sphere, computes its area • Name:areaSphere.m
• Input:Radius(R)
• Output: Area (A)
50
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Function execution – One input and one output
• Once the file is saved write in the command window
answer=functionName(input)
• Note that when a function is used, there
is no need to use the same name of the inputs and outputs in the command window
51
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Function example – Several inputs and outputs
• A function that, given the height and radius of a cylinder return its
area and volume
• Name: areaVolumeCylinder.m
• Inputs: Radius (R) and Height (H)
• Outputs: Area (A) and Volume (V)
• Note that:
• Inputs are separate by commas and
in between brackets (R,H)
• Outputs are separated by commas and
in between square brackets [A,V]
▪ Following the same syntax, functions with an arbitrary number of inputs/outputs can be defined
52
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Function execution – Several inputs and outputs
• Once the file is saved, you can execute (use) a function by
• Writing the name in the command window providing names for the inputs and outputs
• Note that the order is important!!
• 1.123 will be assigned to R (first input) and 0.9129 to H (second input)
• Analogously, the result A is assigned to area (first output) and the result V to the
variable volume (second output)
53
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Tips for designing good programs (1) – MARKING SCHEME
• Sketch the structure of your program before starting Matlab
• Decide how many functions and scripts you need
• Is never a good idea to repeat formulas. Define functions instead and use them whenever is necessary
• Decide the name of the files. Do not underestimate this part, it is important. The names should already indicate what the program is meant to do
• The sketch should indicate the workflow: which functions are called by other Matlab files
• Write each script or function in a separate Matlab file (M-file)
• It is possible to write several functions inside a single M-file but you should
avoid this practice (at least in this module)
• Use meaningful names for the variables
• Regularly introduce comments that help to understand what the program is doing
54
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Tips for designing good programs (2) – MARKING SCHEME
• Structure for scripts:
• Start scripts with a short description or title
• Clear memory and the screen to avoid interference with existing data
• Define the given data or user parameters
• Perform the computations (usually calling functions)
• Present the results on the screen. Avoid multiple outputs on the screen by suppressing or forgetting the semicolon. User disp and fprintf
• Structure for functions:
• Start functions with the header (syntax, inputs and outputs)
• Provide a description of the inputs and outputs
• Define any parameters required
• Perform the computations
▪ Your program must not only produce the correct answer, but also to be understood by anyone else in your absence
55
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Presenting the results on the screen
▪ To write a sentence on the screen use disp
• Example: >> disp(‘Hello world!’)
• The sentence between apostrophes
▪ To write a sentence that contains variable values use fprintf
• Example:
>> a=436;
>> fprintf(‘The current value of a is %d. It is an integer.’, a)
• The sentence between apostrophes
• To write a variable in the sentence use %d if the variable is integer
• After the apostrophe, introduce a comma and the name of the variable
• When the sentence is written, Matlab will change the symbol %d by the current value of the variable a
56
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Presenting the results on the screen
▪ To force a line break after the use fprintf use \n. Examples:
>> fprintf(‘The current value of a is %d. It is an integer.’, a)
>> fprintf(‘The current value of a is %d. It is an integer.\n’, a)
>> fprintf(‘The current value of a is %d.\nIt is an integer.\n’, a)
57
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Presenting the results on the screen
▪ To introduce several variables, introduce them in order • Example:
>> a=2;
>> b=15;
>> fprintf(‘The current value of a is %d.\nThe value of b is %d.\n’, a, b)
▪ Note that the order is important. The first symbol %d will be substituted by the first variable found after the apostrophe (a), The second symbol %d will be substituted by the second variable found after the apostrophe (a) and so on.
58
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 1
▪ Presenting the results on the screen
▪ If the variable is a real number use %f. • Example:
>> a=13.489;
>> fprintf(‘The current value of a is %f.\n’, a)
▪ By default %f will show 6 decimal places. The number of decimal places can be controlled by using %A.Bf being A the number of total digits and B the number of decimal digits
▪ Examples:
59
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
én Sevilla
PART 1 Introduction to Programming 2- Control flow and loops
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Outline
• Control flow
• Take decisions
• Relational and logical operators
• Loops
• Repeat instructions many times
• Fill in vectors
• Iterate until a condition is satisfied
61
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Control flow
▪ So far, we have designed simple programs where all the
instructions written in the M-file are executed in order
▪ What happens when, depending on the value of a variable, different instructions should be executed?
▪ Example:
• Using your mobile: If you write the correct PIN you can access
otherwise you receive a message and you cannot access ▪ Matlab implementation:
• Create a function that, given a pin number, compares it to the correct password and returns 1 if the PIN is correct and 0 otherwise
62
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Control flow
▪ How do we check if a condition is satisfied?
if condition
instructions to be executed
end
▪ The instructions inside the if statement are only executed if the condition is true
▪ == means “equal to”
Two equals are used to compare One equal is used to assign variables
63
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Control flow
▪ If we have several conditions and different instructions to be
executed in each case we use
How does it work?
• If condition1 is true, instructions1 are executed and no more conditions are checked
• If condition1 is false, condition2 is checked. If conditions2 is true, instructions2 are executed and no more conditions are checked
• If conditions1 and condition2 are false, instructions3 are executed
if (condition1) instructions1
elseif (condition2) instructions2
else
instructions3
end
64
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Control flow
▪ How do we check if a condition is satisfied?
▪ Relational and logical operators:
• == Equal to
• > Greater than
• < Less than
• >= Greater than or equal to
• <= Less than or equal to
• ~= NOT equal to
65
SWANSEA UNIVERSITY - College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Control flow ▪ Example:
• Create a function that, given a number returns one of the following three messages
• The number is less than 10.
• The number is greater than (or equal to) 10 but less than 20.
• The number is greater than (or equal to) 20.
66
SWANSEA UNIVERSITY - College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Control flow
▪ If we have to check two (or more) conditions at the same time
use the following symbols for AND and OR
• && AND operator
• || OR operator
▪ Example:
• Define the following function in Matlab
67
SWANSEA UNIVERSITY - College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Loops
▪ for loops allow us to repeat certain instructions reducing the
amount of writing significantly ▪ Syntax:
for counter=initial:step:final instruction
end
• Write a Matlab program that writes on the screen the odd numbers from 1 to 23
• Remember: 1:N is equivalent to “1:1:N” (initial:step:final)
68
SWANSEA UNIVERSITY - College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Loops
▪ Example 1: Write a function in Matlab that, given a positive
integer n computes the sum
▪ Note how the result is first initialised to be zero
▪ Then each term is added to the previous value, so we accumulate the final value of the sum
69
SWANSEA UNIVERSITY - College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Loops
▪ for loops are commonly used to fill in vectors (or matrices) and
to iterate until something happens
▪ Fill in a vector example: Plot the function below between [-3,0]
▪ Strategy
• Define a vector with enough points between [-3,0]
• x = [-3, -2.95, -2.90, -2.85, ....., -0.10, -0.05, 0]
• In Matlab: >> x = [-3:0.05:0]; or >> x = linspace(-3,0,61);
• Compute the value of G for each component of x and store it in another vector
• In Matlab: initialize a vector of zeros and fill in each component with a for loop from i=1 to 61
70
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Loops
▪ Fill in a vector example: Plot the function below between [-3,0]
71
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Loops
▪ Iterate until something happens
▪ Exercise: The pressure (p) on a can depends on the temperature (T) according to the law
p(T)=0.00105T2 + 0.0042T + 1.352
▪ Create a program that, starting with a temperature of 15oC, and increasing the temperature by 0.5oC, checks when the pressure is more than 3.2atm. The program should end when this value is detected
▪ Use break to stop the program when the desired condition is true
72
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Introduction to Programming 2
▪ Loops
▪ Model answer
▪ The command break stops the execution of the loop, so we are not actually iterating until T=Tfinal
73
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
PART 2 Interpolation and approximation
én Sevilla
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Objectives
▪ Interpolation and approximation
• Types of approximation
• Approximation criteria
▪ How to build an approximation?
▪ Criteria for deciding which is the best way to approximate a function or a set of data
75
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ After conducting experiments in a wind tunnel we have obtained the following data for the drag coefficient of a sphere as a function of the Reynolds number
▪ How can we obtain the drag for Re=30?
76
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ The temperature f of a certain gas in terms of the electric charge x is known to be given by
▪ If we need to evaluate repeatedly the derivatives or integrals of this function it is desirable to approximate the function by a simpler expression
▪ A polynomial will be
• Easy to differentiate
• Easy to integrate
77
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ We need to specify
▪ The function to be approximated
• Discrete: set of points (table, usually measurements)
• Continuous: analytical expression
78
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ We need to specify
▪ The function to be approximated (Discrete – Continuous)
▪ Type of approximant
• Polynomials: cheap to evaluate, easy to differentiate and to integrate
79
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ We need to specify
▪ The function to be approximated (Discrete – Continuous)
▪ Type of approximant (Polynomials)
• Usually the function is approximated with a single polynomial, but
sometimes it is better to split the problem in smaller problems:
Splines: They allow more local control, i.e. modifying the data in one region does not affect the whole description of the approximation, only locally
80
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ We need to specify
▪ The function to be approximated (Discrete – Continuous) ▪ Type of approximant (Polynomials – Polynomial Splines)
▪ Approximation criteria
• Interpolation: find the polynomial that passes through the data
• Approximation: find the polynomial that “best fits” the data
81
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
How is the function
to be approximated
defined? (Continuous – Discrete)
Which type of function do
we use as an
approximation? (Polynomial – Polynomial Splines)
In which sense is the
approximation defined? (Interpolation – Approximation)
▪ We need to specify
▪ The function to be approximated (Discrete – Continuous) ▪ Type of approximant (Polynomials – Polynomial Splines) ▪ Approximation criteria (Interpolation – Approximation)
82
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Polynomial interpolation
▪ Interpolation is useful
• When we want the approximation to pass exactly through the data.
• To approximate the value of the solution at intermediate points
▪ Interpolation is not useful
• When the data is obtained from experimental results. In this situation the results contain some errors and, therefore, it is not a good idea to ask for the approximation to satisfy exactly the conditions
• To extrapolate results (out of the domain)
▪ An alternative is to use a polynomial approximation
• The degree of the polynomial will be defined by us, not by the amount of data
• We say that the polynomial is the best fit of the data, but what does “best fit” means?
83
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Polynomial least-squares approximation
▪ Linear regression: approximation with a straight line
▪ Given a set of data we look for the straight line that better approximates the data in a least-squares sense
▪ The equation of the line is
with the coefficients obtained solving
di (xi,yi)
84
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Polynomial least-squares approximation ▪ Least-squares approximation is useful
• When the data is obtained from experimental results. In this situation the results contain some errors and, therefore, it is not a good idea to ask for the approximation to satisfy exactly the conditions
• To extrapolate results (out of the domain) because the least-squares approximation gives the trend of the solution
▪ Least-squares approximation not useful
• When we want the approximation to pass exactly through the data.
85
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Polynomial interpolation with splines ▪ With polynomial interpolation
• The degree of the polynomial is given by the amount of data. A big amount of data means a very high-order polynomial
• The approximation is global (a small change on the data produce global changes on the approximation)
▪ With least-squares approximation
• The approximation does not passes through the data
• The approximation is global (a small change on the data produce global changes on the approximation)
▪ Al alternative is to define a piecewise approximation (splines)
86
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Polynomial interpolation with splines
▪ The most simple and intuitive option: Linear C0 spline
• Linear because a polynomial of degree 1 is used on each subinterval
• C0 because the approximation is continuous between intervals but the derivative is not
87
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Polynomial interpolation with splines ▪ Sometimes we require smoothness
• Design of a roller coaster
▪ The most common splines are
• Cubic C1 spline. We impose
• Interpolation of the data
• Continuity of the approximation and the first derivative between subintervals (control oscillations)
• Cubic C2 spline. We impose
• Interpolation of the data
• Continuity of the approximation, the first and second derivatives between subintervals (control smoothness)
88
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Interpolation and approximation
▪ Summary
▪ With polynomial interpolation
• The degree of the polynomial is given by the amount of data. A big amount of data means a very high-order polynomial
• The approximation is global (a small change on the data produce global changes on the approximation)
▪ With least-squares approximation
• The approximation does not passes through the data
• The approximation is global (a small change on the data produce global changes on the approximation)
▪ Al alternative is to define a piecewise approximation (splines)
• We can do interpolation or approximation to build each piece
• The approximation is local (a small change on the data produce changes on just some pieces)
89
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
én Sevilla
PART 2 Root finding
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ Objectives
▪ What is root?
▪ How to approximate the roots of a nonlinear equation?
▪ Numerical methods: • Bisection
• Newton
▪ How to obtain the solution with the desired accuracy
91
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ An aluminium column is subjected to an eccentric axial load P
▪ The maximum stress is given by
▪ We want to obtain the maximum load P that can be supported
• the cross sectional area is A=25800mm2, the length is L=7100mm, the eccentricity of the load is e=85mm, the half depth of the column is c=170mm, the radius of gyration is r=142mm and the Young modulus is E=71·109 and the maximum stress is 120·106Pa
92
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ We can plot the function using Matlab to have an very crude approximation but, how can we obtain an accurate value?
93
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ Root finding consists on determining the values of x for which f(x)=0
▪ Numerical methods for root finding are iterative. This means that an initial estimate is required. After that, a succession is constructed in such a way that it converges to the root
▪ Note that the function f(x) may have several roots, so the initial estimate is crucial. A bad estimate can produce convergence to the “wrong” result
▪ In order to obtain a good initial approximation we can use
• Physical knowledge of the problem
• A plot of the function
94
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ First, we need to define the function for which the root has to be found
▪ For instance, in the column example, if we are looking for the P corresponding to the maximum stress,
this means that we want to obtain the root of the function
95
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ Plotting the function f(P) we observe that the root is in the interval
▪ The value at the extremes is
▪ How do we know that the root
is inside the interval?
• Continuous function
• Change of sign
▪ We can evaluate the function at the midpoint:
▪ What can we say about the root?
96
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ We know that the root has to be in the interval ▪ The value at the extremes is
▪ We can evaluate the function at the midpoint:
▪ What can we say about the root?
97
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ We know that the root has to be in the interval
▪ The value at the extremes is
▪ We can evaluate the function at the midpoint:
▪ What can we say about the root?
▪ Repeating the process we can determine the root with the
desired accuracy
▪ This process is called the BISECTION METHOD
98
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ When do we stop the computation?
• We can monitor the value of the function at the approximated root at
each iteration.
• Note also that if the difference between two approximations is not physically relevant, there is no need to continue iterating
99
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ When do we stop the computation? We introduce in the code two stopping criteria. They have to be used together
• Criteria on f: when the value of the function at the approximate root is small enough
• Criteria on x: when the difference between two approximations of the root is small enough
▪ What does small enough means?
• Depends on the problem. You need to use your physical knowledge
of the problem. For instance:
If we are computing a distance in a problem with dimensions in meters, we should know if a precision of millimetres is enough
100
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ Another approach using the geometric interpretation of derivatives
• Select an initial guess
• Evaluate the derivative of the function at this point
• Compute the tangent line and the intersection of the tangent and the
x axis
• This is the new approximation of the root
• Repeating the process we build a succession
101
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ How is obtained in terms of ?
▪ The equation of the tangent is
where slope (m) is given by the derivative and the ordinate in the origin (n) is obtained imposing that the line passes through
▪ The tangent is
▪ Then, we can obtain by imposing y=0 in the previous equation:
This is the NEWTON’S METHOD
102
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ Results for the aluminium column problem
▪ Don’t forget to use the two stopping criteria:
• Criteria on f: when the value of the function at the approximate root is
small enough
• Criteria on x: when the difference between two approximations of the
root is small enough
103
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ Discussion
▪ The Newton’s method is faster.
It needs less iterations to compute the root with a desired accuracy.
104
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Root finding
▪ In Matlab
▪ In the PC Lab you will learn to solve root finding problems using
a technique that uses the ideas presented here
▪ The key function is called fzero. It is not bisection or Newton but it is similar as
• Requires an initial approximation
• Iterates until the solution is found
105
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
én Sevilla
PART 2 Numerical integration
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Objectives
▪ Numerical integration
• Simple trapezoidal rule
• Simple Simpson’s rule
• Composite rules
• Adaptive rules
107
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Analytical integration is not always feasible
▪ Debye’s formula for the heat capacity of a solid is , where N is the number of particles in the solid, k is the Boltzmann constant, u is the ratio between the absolute and Debye’s temperature and
▪ In fact, we do not require the analytical expression of the integral, we just need a way to evaluate the integral for each value of u
108
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Sometimes we don’t even know the analytical expression
▪ The pressure exerted by the wind is measured at various heights on a vertical wall, as shown in the figure. Find the height of the pressure centre, defined as
▪ How can we integrate a function if we don’t know its analytical expression?
109
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ How can we integrate a function if we don’t know its analytical expression?
▪ We know how to interpolate/approximate a function, so we can approximate the integral of a function by the integral of the approximation of the function
▪ If we interpolate the function using a polynomial of degree one
110
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ How can we integrate a function if we don’t know its analytical expression?
▪ We know how to interpolate/approximate a function, so we can approximate the integral of a function by the integral of the approximation of the function
▪ If we interpolate the function using a polynomial of degree one
▪ TRAPEZOIDAL RULE
111
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ How can we integrate a function if we don’t know its analytical expression?
▪ We know how to interpolate/approximate a function, so we can approximate the integral of a function by the integral of the approximation of the function
▪ If we interpolate the function using a polynomial of degree two
112
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ How can we integrate a function if we don’t know its analytical expression?
▪ We know how to interpolate/approximate a function, so we can approximate the integral of a function by the integral of the approximation of the function
▪ If we interpolate the function using a polynomial of degree two
▪ SIMPSON RULE
113
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ In general is not enough to employ the simple trapezoidal or Simpson rules because the variation of the function makes impossible to properly capture it using a single line or parabola
Trapezoidal rule Simpson’s rule
114
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ A simple alternative is to split the interval [a,b] in several subintervals and apply the trapezoidal or the Simpson’s rule in each subinterval.
▪ Composite trapezoidal rule
115
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ A simple alternative is to split the interval [a,b] in several subintervals and to apply the trapezoidal or the Simpson’s rule in each subinterval
▪ Composite trapezoidal rule
116
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite trapezoidal rule with equally-spaced points
117
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite trapezoidal rule with equally-spaced points
118
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite trapezoidal rule with equally-spaced points
119
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite trapezoidal rule with equally-spaced points • Adding all contributions
120
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite Simpson’s rule with equally-spaced points
121
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite Simpson’s rule with equally-spaced points
122
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Composite rules
▪ Composite Simpson’s rule with equally-spaced points
123
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Numerical integration
▪ Comparison
▪ The efficiency is measured by the number of points used, so we are interested in the
quadrature that uses
the minimum number
of points to give
a desired accuracy
124
SWANSEA UNIVERSITY College of Engineering
Problem Solving in Engineering with MATLAB
EG-228
Dr. Rub
PART 2 Ordinary Differential Equations (ODEs)
én Sevilla
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Objectives
▪ Ordinary differential equations (ODEs)
• What is an ODE?
• Order of an ODE
• Initial conditions
• From a high-order ODE to a system of first-order ODEs
126
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ ODEs involve a function and its derivatives. They appear in a large number of engineering applications
▪ The vertical deflection y(x) of a cantilever beam due to its own weight is given by the ODE
where
127
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ The order of an ODE is the highest derivative that appears in the equation.
128
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Initial conditions
▪ Can you guess the solution of this ODE?
▪ The following function satisfies the ODE:
▪ But also any other function such as , where C is
a constant
▪ We need one condition in order to determine the constant, for instance the solution at one extreme of the solution interval
▪ The constant is then determined
INITIAL CONDITION
129
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Initial conditions
▪ Can you guess the solution of this ODE?
▪ The following function satisfies the ODE:
▪ But also any other function such as , where C is
a constant
▪ In fact, any function such as D are constants
, where C and
▪ We need now two conditions in order to determine the solution 130
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Initial Value Problems (IVP)
▪ To determine the solution of an ODE we need as many initial
conditions as the order of the ODE
▪ An ODE with initial condition/s is called an initial value problem
(IVP)
Example:
Pendulum with no friction:
131
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Numerical solution of ODEs
▪ Numerical methods for ODEs are only designed to solve first-
order equations or systems written in the form First-order ODE First-order ODE system
▪ If we have a higher order ODE, the first step is to write the high order ODE as a system of first-order ODEs
132
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Numerical solution of IVPs
▪ The velocity of a skydiver in free fall is given by the differential
equation
where the mass is m=80 Kg, the gravity acceleration is taken as
g=9.81 m/s2 and the drag coefficient is k=0.2028 Kg/m.
▪ First we need to write the equation in the general form
▪ In this case
133
This example is solved in the PC Lab notes
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ If we have a ODE of order n, we will rewrite it as a system of n
first-order ODEs
▪ Movement of a pendulum with no friction:
• Second-order ODE • System of two first-order ODEs
with initial conditions with initial condition
▪ We need to detail the expression of , and
134
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ The strategy consists on taking y as a vector containing the unknown function and its n-1 first derivatives
▪ Then we deduce the expression of f component by component • First component:
In terms of y
135
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ The strategy consists on taking y as a vector containing the unknown function and its n-1 first derivatives
▪ Then we deduce the expression of f component by component • First component:
136
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ The strategy consists on taking y as a vector containing the unknown function and its n-1 first derivatives
▪ Then we deduce the expression of f component by component • First component:
137
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ The strategy consists on taking y as a vector containing the unknown function and its n-1 first derivatives
▪ Then we deduce the expression of f component by component
• First component:
• Second component:
Using the ODE
138
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ The strategy consists on taking y as a vector containing the unknown function and its n-1 first derivatives
▪ Then we deduce the expression of f component by component
• First component:
• Second component:
In terms of y
139
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ The strategy consists on taking y as a vector containing the unknown function and its n-1 first derivatives
▪ Then we deduce the expression of f component by component
• First component:
• Second component:
140
SWANSEA UNIVERSITY – College of Engineering Problem Solving in Engineering with MATLAB · EG-228 Dr. Rubén Sevilla
Ordinary Differential Equations (ODEs)
▪ Higher order ODEs
▪ Movement of a pendulum w?ith no friction:
▪ So, we have
▪ Finally, do not forget the initial conditions
▪ Once the high-order equation is translated to a system, we can apply standard numerical methods for solving it
141