Matlab task, please finish this as required with all the step with matlab

1, finish the question with all the steps, the file I already provided below.

2, when you finish the MATLAB code, please explain why you did this step with %% in MATLAB so that I can understand why you did that step.

3, important! Please finish on your own.4

4, pay attention to the ‘Hint 1 to 3’ in the file, this is the way you need to use.

5, you must use a quadratic spline polynomial, and you are not allowed to use built-in Matlab commands for interpolation or polynomial generation.

6, do a bonus if you can.

7, here is an example I provide to you. ”topic_3 page 51-67, problem 6 uses the similar method.

CHEE 3602, Topic 3: Interpolating polynomials
Fall 2023
Lab 6: Spline interpolation
Background
As discussed in the introductory lecture to interpolation, the T-xy diagram represents an
example of a function where it is relatively easy to calculate the fraction of a compound
(such as benzene) in the liquid or vapour phase given a temperature, but hard to do the
reverse (calculating the dew or bubble point of a mixture given only the liquid or vapour
fraction). One way to deal with this challenge is to calculate a set of liquid/vapour fraction
values over a range of temperatures, and then use an interpolating function to generate a
polynomial that expresses liquid/vapour fraction as a function of temperature.
In this lab, we will focus on just the liquid fraction of a binary benzene/toluene mixture.
Recall the process for calculating a liquid fraction from temperature:
1. Given an input pressure value 𝑃𝑡𝑜𝑡 , calculate the boiling point temperature for pure
components
– use Antoine’s law log10 (𝑃∗ ) = 𝐴 − 𝐵/(𝑇 + 𝐶)
2. With the overall temperature range established, choose a set of temperature values
within this range and calculate the vapour pressure of the two components for each
temperature
– once again, use Antoine’s law log10 (𝑃∗ ) = 𝐴 − 𝐵/(𝑇 + 𝐶)
– but reverse it to calculate 𝑃∗ from 𝑇
3. Then, calculate liquid fraction benzene (𝑥𝐵 ) based on these vapour pressures
∗
∗
– use Raoult’s and Dalton’s law, i.e., 𝑃𝑡𝑜𝑡 = 𝑃𝐵
𝑥𝐵 + 𝑃𝑇
𝑥𝑇
∗
∗
– recall that the mixture is binary, so 𝑥𝑇 = 1 − 𝑥𝐵 and 𝑃𝑡𝑜𝑡 = 𝑃𝐵
𝑥𝐵 + 𝑃𝑇
(1 − 𝑥𝐵 )
– and remember that 𝑃𝑡𝑜𝑡 is given as an input
4. The result is a set of temperature (𝑇) and liquid fraction benzene (𝑥𝐵 ) values, which
can be used to set up an interpolating function with 𝑇 as the 𝑦 values and 𝑥𝐵 as the x
values.
This lab’s task will be to estimate the bubble point temperature of a benzene/toluene mixture using a quadratic interpolating spline. Given a set of (𝑥𝑖 , 𝑓𝑖 ) points, recall the quadratic
spline fit process:
Page 1 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Fall 2023
1. Given 𝑛 points, there should be 𝑛 − 1 total spline polynomials
2. 𝑎𝑖 is set to the value of 𝑓𝑖
3. 𝑏𝑖 is solved from the following system of linear equations:
2( 𝑓2 − 𝑓1 )
⎡1
⎢0
⎢
⎢0
⎢
⎢⋮
⎢0
⎢
⎣0
4. 𝑐𝑖 is solved from
1
1
0
⋮
0
0
0
1
1
⋮
0
0
… 0
… 0
… 0
⋱ ⋮
… 1
… 0
⎡
⎤
ℎ1
0⎤⎡ 𝑏1 ⎤ ⎢ 2( 𝑓3 − 𝑓2 ) ⎥
⎢
⎥
ℎ2
⎥
0 ⎥⎢ 𝑏 2 ⎥ ⎢
⎥⎢
⎥ ⎢ 2( 𝑓4 − 𝑓3 ) ⎥
0 ⎥⎢ 𝑏 3 ⎥ ⎢
⎥
ℎ3
⎥⎢
⎥=⎢
⎥
⋮ ⎥⎢ ⋮ ⎥ ⎢
⎥
⋮
⎢𝑏
⎥ ⎢ 2( 𝑓 − 𝑓 ) ⎥
1⎥
⎥⎢ 𝑛−2 ⎥ ⎢ 𝑛−1 𝑛−2 ⎥
⎥
ℎ𝑛−2
1⎦⎣𝑏𝑛−1 ⎦ ⎢
⎢ 𝑓𝑛− 𝑓𝑛−1 ⎥
⎣
⎦
ℎ𝑛−1
𝑏𝑖+1 −𝑏𝑖
for [𝑖 = 1, 𝑛 − 2]
2ℎ𝑖
5. 𝑐𝑛−1 is set to 0
Once all polynomial terms are calculated, 𝑓(𝑥) is found by first, identifying which polynomial (𝑖) 𝑥 belongs to, and then plugging 𝑥 into the polynomial equation of 𝑓(𝑥) =
𝑎𝑖 + 𝑏𝑖 (𝑥 − 𝑥𝑖 ) + 𝑐𝑖 (𝑥 − 𝑥𝑖 )2 .
Task
Write a function to estimate the bubble point temperature of a benzene/toluene mixture
given total pressure, liquid fraction benzene, and the number of points that should be
used to build a quadratic interpolating polynomial. To make sure everyone uses the same
values, the spline function data must set the second derivative of the last polynomial to
zero (using pure toluene, i.e., no benzene, as the first point). Use the following data for
Antoine’s equation (with the temperature in Celsius and pressure in mmHg):
Compound
A
B
C
Benzene
6.89272
1203.531
219.888
Toluene
6.95805
1346.773
219.693
You must use a quadratic spline polynomial and you are not allowed to use built-in Matlab
commands for interpolation or polynomial generation.
Page 2 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Fall 2023
Core requirements
Function name
bt_txy
Input
Total pressure (in mmHg), liquid fraction benzene (𝑥𝐵 ), and the
total number of points used to perform the spline interpolation
(𝑛)
Output
A single temperature value corresponding to the bubble point of
the mixture
Errors
If pressure input is less than 0, fraction benzene is not a valid fraction, or the total number of points is less than 3
Example
bt_txy(760, 0.5, 4) should give 91.9104
There are several hints for this lab:
Hint 1
You will need to solve a system of linear equations to get all the 𝑏𝑖 values.
You do not need to program your own Gauss elimination to do this. Just
use Matlab’s linsolve command. Given a system of equations of the
form 𝐴𝑥 = 𝑏, you would use x = linsolve(A, b).
Hint 2
The system of equations that appears in the spline equation only has nonzero values on 2 diagonal lines. You can use Matlab’s diag command to
help construct this matrix. Check the diag help page for more information
Hint 3
It is possible to meet core requirements without any looping whatsoever.
You can use loops if you want, but remember Matlab’s array math.
Bonus requirements
Expand the original bt_txy function to also display a figure plotting the original data as
points and the interpolating spline as a line.
Page 3 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Function name
Fall 2023
bt_txy
Input
Total pressure (in mmHg), liquid fraction benzene (𝑥𝐵 ), and the
total number of points used to perform the spline interpolation
(𝑛)
Output
A single temperature value corresponding to the bubble point of
the mixture
Side-effect
Generate a figure plotting the T-xy data as points with the interpolating spline as a line as a pop-up
Errors
Example
If pressure input is less than 0, fraction benzene is not a valid fraction, or the total number of points is less than 3
bt_txy(760, 0.5, 4) should give 91.9104 as well as pop up a
plot similar to the example below
Submission checklist
Use the following checklist to help ensure you are meeting the core lab requirements:
You are submitting a single file called bt_txy.m which defines a function called bt_txy
Your function returns a value:
bt_txy(760, .5, 4)/2 should not cause problems
Your function should work correctly for any value of 𝑛 ≥ 3 and liquid fraction 0 ≥ 𝑥𝐵 ≥
1
Page 4 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Fall 2023
Lab 6: Spline interpolation
Background
As discussed in the introductory lecture to interpolation, the T-xy diagram represents an
example of a function where it is relatively easy to calculate the fraction of a compound
(such as benzene) in the liquid or vapour phase given a temperature, but hard to do the
reverse (calculating the dew or bubble point of a mixture given only the liquid or vapour
fraction). One way to deal with this challenge is to calculate a set of liquid/vapour fraction
values over a range of temperatures, and then use an interpolating function to generate a
polynomial that expresses liquid/vapour fraction as a function of temperature.
In this lab, we will focus on just the liquid fraction of a binary benzene/toluene mixture.
Recall the process for calculating a liquid fraction from temperature:
1. Given an input pressure value 𝑃𝑡𝑜𝑡 , calculate the boiling point temperature for pure
components
– use Antoine’s law log10 (𝑃∗ ) = 𝐴 − 𝐵/(𝑇 + 𝐶)
2. With the overall temperature range established, choose a set of temperature values
within this range and calculate the vapour pressure of the two components for each
temperature
– once again, use Antoine’s law log10 (𝑃∗ ) = 𝐴 − 𝐵/(𝑇 + 𝐶)
– but reverse it to calculate 𝑃∗ from 𝑇
3. Then, calculate liquid fraction benzene (𝑥𝐵 ) based on these vapour pressures
∗
∗
– use Raoult’s and Dalton’s law, i.e., 𝑃𝑡𝑜𝑡 = 𝑃𝐵
𝑥𝐵 + 𝑃𝑇
𝑥𝑇
∗
∗
– recall that the mixture is binary, so 𝑥𝑇 = 1 − 𝑥𝐵 and 𝑃𝑡𝑜𝑡 = 𝑃𝐵
𝑥𝐵 + 𝑃𝑇
(1 − 𝑥𝐵 )
– and remember that 𝑃𝑡𝑜𝑡 is given as an input
4. The result is a set of temperature (𝑇) and liquid fraction benzene (𝑥𝐵 ) values, which
can be used to set up an interpolating function with 𝑇 as the 𝑦 values and 𝑥𝐵 as the x
values.
This lab’s task will be to estimate the bubble point temperature of a benzene/toluene mixture using a quadratic interpolating spline. Given a set of (𝑥𝑖 , 𝑓𝑖 ) points, recall the quadratic
spline fit process:
Page 1 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Fall 2023
1. Given 𝑛 points, there should be 𝑛 − 1 total spline polynomials
2. 𝑎𝑖 is set to the value of 𝑓𝑖
3. 𝑏𝑖 is solved from the following system of linear equations:
2( 𝑓2 − 𝑓1 )
⎡1
⎢0
⎢
⎢0
⎢
⎢⋮
⎢0
⎢
⎣0
4. 𝑐𝑖 is solved from
1
1
0
⋮
0
0
0
1
1
⋮
0
0
… 0
… 0
… 0
⋱ ⋮
… 1
… 0
⎡
⎤
ℎ1
0⎤⎡ 𝑏1 ⎤ ⎢ 2( 𝑓3 − 𝑓2 ) ⎥
⎢
⎥
ℎ2
⎥
0 ⎥⎢ 𝑏 2 ⎥ ⎢
⎥⎢
⎥ ⎢ 2( 𝑓4 − 𝑓3 ) ⎥
0 ⎥⎢ 𝑏 3 ⎥ ⎢
⎥
ℎ3
⎥⎢
⎥=⎢
⎥
⋮ ⎥⎢ ⋮ ⎥ ⎢
⎥
⋮
⎢𝑏
⎥ ⎢ 2( 𝑓 − 𝑓 ) ⎥
1⎥
⎥⎢ 𝑛−2 ⎥ ⎢ 𝑛−1 𝑛−2 ⎥
⎥
ℎ𝑛−2
1⎦⎣𝑏𝑛−1 ⎦ ⎢
⎢ 𝑓𝑛− 𝑓𝑛−1 ⎥
⎣
⎦
ℎ𝑛−1
𝑏𝑖+1 −𝑏𝑖
for [𝑖 = 1, 𝑛 − 2]
2ℎ𝑖
5. 𝑐𝑛−1 is set to 0
Once all polynomial terms are calculated, 𝑓(𝑥) is found by first, identifying which polynomial (𝑖) 𝑥 belongs to, and then plugging 𝑥 into the polynomial equation of 𝑓(𝑥) =
𝑎𝑖 + 𝑏𝑖 (𝑥 − 𝑥𝑖 ) + 𝑐𝑖 (𝑥 − 𝑥𝑖 )2 .
Task
Write a function to estimate the bubble point temperature of a benzene/toluene mixture
given total pressure, liquid fraction benzene, and the number of points that should be
used to build a quadratic interpolating polynomial. To make sure everyone uses the same
values, the spline function data must set the second derivative of the last polynomial to
zero (using pure toluene, i.e., no benzene, as the first point). Use the following data for
Antoine’s equation (with the temperature in Celsius and pressure in mmHg):
Compound
A
B
C
Benzene
6.89272
1203.531
219.888
Toluene
6.95805
1346.773
219.693
You must use a quadratic spline polynomial and you are not allowed to use built-in Matlab
commands for interpolation or polynomial generation.
Page 2 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Fall 2023
Core requirements
Function name
bt_txy
Input
Total pressure (in mmHg), liquid fraction benzene (𝑥𝐵 ), and the
total number of points used to perform the spline interpolation
(𝑛)
Output
A single temperature value corresponding to the bubble point of
the mixture
Errors
If pressure input is less than 0, fraction benzene is not a valid fraction, or the total number of points is less than 3
Example
bt_txy(760, 0.5, 4) should give 91.9104
There are several hints for this lab:
Hint 1
You will need to solve a system of linear equations to get all the 𝑏𝑖 values.
You do not need to program your own Gauss elimination to do this. Just
use Matlab’s linsolve command. Given a system of equations of the
form 𝐴𝑥 = 𝑏, you would use x = linsolve(A, b).
Hint 2
The system of equations that appears in the spline equation only has nonzero values on 2 diagonal lines. You can use Matlab’s diag command to
help construct this matrix. Check the diag help page for more information
Hint 3
It is possible to meet core requirements without any looping whatsoever.
You can use loops if you want, but remember Matlab’s array math.
Bonus requirements
Expand the original bt_txy function to also display a figure plotting the original data as
points and the interpolating spline as a line.
Page 3 of 4
CHEE 3602, Topic 3: Interpolating polynomials
Function name
Fall 2023
bt_txy
Input
Total pressure (in mmHg), liquid fraction benzene (𝑥𝐵 ), and the
total number of points used to perform the spline interpolation
(𝑛)
Output
A single temperature value corresponding to the bubble point of
the mixture
Side-effect
Generate a figure plotting the T-xy data as points with the interpolating spline as a line as a pop-up
Errors
Example
If pressure input is less than 0, fraction benzene is not a valid fraction, or the total number of points is less than 3
bt_txy(760, 0.5, 4) should give 91.9104 as well as pop up a
plot similar to the example below
Submission checklist
Use the following checklist to help ensure you are meeting the core lab requirements:
You are submitting a single file called bt_txy.m which defines a function called bt_txy
Your function returns a value:
bt_txy(760, .5, 4)/2 should not cause problems
Your function should work correctly for any value of 𝑛 ≥ 3 and liquid fraction 0 ≥ 𝑥𝐵 ≥
1
Page 4 of 4
CHEE 3602 – Topic 3:
Interpolating polynomials
Stanislav Sokolenko
Fall 2023
Motivation
Interpolating polynomials
Motivation
Why interpolate?
Many important mathematical operations are only defined for functions (equations)
e.g. the concept of derivative/integral only makes sense for functions — what is “the rate of change” of a set of data points?
Collected empirical data must be converted into functional form so
that it can be processed/analyzed
e.g. calculating the rate of reaction given a set of concentrations
requires fitting a function to the concentration data
3/118
Interpolating polynomials
Motivation
Numerical methods
Interpolation is the process of fitting individual data points to
a function — enabling analysis that would not be possible using discrete points.
Many common numerical methods are derived using interpolation, even if there is no apparent interpolation taking place.
4/118
Interpolating polynomials
Motivation
Numerical methods
This topic is divided into two specific tasks:
Building interpolating polynomials
Applying interpolating polynomials to develop methods for differentiation/integration
We will start by going through a practical example of interpolation and
justifying the use of polynomials for interpolation.
5/118
Interpolating polynomials
Motivation
6/118
Defining interpolation
Interpolation refers to the approximation of a function
using a series of tabulated values.
𝑓 (𝑥)
Interpolating polynomials
Motivation
Reading tables
In effect, interpolation is used to estimate values between table rows.
T (°C)
P
*
(Hg)
90
525.76
92
566.99
94
610.91
96
657.62
98
707.27
100
760.00
Given the table above, what is the vapour pressure of water at 96.5 °C?
7/118
Interpolating polynomials
Motivation
8/118
Reading tables
Given that there is no entry for 96.5 °C, there are a couple of options
T (°C)
P
*
(Hg)
90
525.76
92
566.99
94
610.91
96
657.62
98
707.27
100
760.00
The simplest is to assume that
96.5 ≈ 96
Therefore,
P
(96.5) ≈ 657.62
*
…
Interpolating polynomials
Motivation
Reading tables
This is a very rough approximation. How can it be improved?
T (°C)
P
*
(Hg)
90
525.76
92
566.99
94
610.91
96
657.62
98
707.27
100
760.00
96.5 is between 96 and 98
(96.5 − 96)/(98 − 96) = 0.25
So the pressure should be higher
P
(96.5) ≈ 657.62 + (707.27 − 657.62)(0.25)
*
= 670.03
9/118
Interpolating polynomials
Motivation
Mathematically
Both approximations can be expressed as interpolation:
1.
The 657.62 estimate used a 0 order polynomial
2.
The 670.03 estimate used a 1st order polynomial
The 1st order polynomial captures more information around
the surrounding table entries than the 0 order.
the only option?
But is this is
10/118
Interpolating polynomials
Motivation
11/118
Interpolating functions
The goal of interpolation is to define a function that can fill in the
gaps between data table entries (or data points)
If the function that generated the original data is unknown, it makes
sense to choose something convenient
A good interpolating function should be easy to
generate (formulate)
evaluate
differentiate
integrate
Interpolating polynomials
Motivation
Polynomials
Polynomials satisfy every single one of the requirements!
12/118
Interpolating polynomials
Motivation
Polynomials
An nth-order polynomial is defined as
𝑃𝑛 (𝑥) = 𝑎0 + 𝑎1 𝑥 + 𝑎2 𝑥2 + … + 𝑎𝑛 𝑥𝑛
An nth-order polynomial
requires
𝑃𝑛 (𝑥) has 𝑛 + 1 coefficients and therefore
𝑛 + 1 points to estimate.
13/118
Interpolating polynomials
Motivation
Back to the table
Consider estimating vapour pressure at 96.5 one more time,
T (°C)
P
*
(Hg)
90
525.76
92
566.99
94
610.91
96
657.62
98
707.27
100
760.00
Three values can be used to fit a quadratic.
14/118
Exumple
Tefx
)
x *> ( x X 3 )
L
R 凶)
—
(x *) ( x ×7) fr
–
–
* *^) < L - x - , ⽕器 * - f — × ) 3 - t , 还 wn, + xx 0 0 ) )xc . y ( Lx t ⼀ - " . 3 )1 - c ) 9 ((1 0) t - JM L JC . } ) - 132 01 399 68 ^+ 284 , : - . = 13 66 x . ^ - . + 44 18x + 110 62 . . - c 世≈ qi 0 57 91 399 .68 . + t ) ω L ) ( 0 . 0 . 39 ) 391 . 96 ) ×+ 110 . 62 - Interpolating polynomials Motivation Formal justification A nice feature of using polynomials is that they can be used to approximate a subset of any continuous function. More formally, the Weierstrass approximation theorem states 𝑓 (𝑥) is a continuous function in the closed interval 𝑎 ≤ 𝑥 ≤ 𝑏, then for every 𝜀 > 0, there exists a polynomial
𝑃𝑛 (𝑥), where the value of 𝑛 depends on the value of 𝜀 such that
for all x in the closed interval 𝑎 ≤ 𝑥 ≤ 𝑏, |𝑃𝑛 (𝑥) − 𝑓 (𝑥)| < 𝜀.” that: “If 15/118 Interpolating polynomials Motivation Polynomial uniqueness It is only possible to define one polynomial of order through fined. 𝑛 that passes 𝑛 + 1 points, regardless of how the polynomial is de- This means that while some methods for polynomial construction may be more convenient than others, the final output is often identical. 16/118 Interpolating polynomials Motivation Chemical engineering example Temperature Recall the T-xy diagram: Vapour fraction Liquid fraction Fraction benzene Although there is an algorithm for getting straightforward way of getting 𝑇 = 𝑓 (𝑥, 𝑦) 𝑥, 𝑦 = 𝑓 (𝑇), there is no 17/118 Interpolating polynomials Motivation Interpolation One practical option is to generate a set of 𝑥, 𝑦 points using 𝑇 , and then interpolate to find the reverse relation: 𝑇 = 𝑓 (𝑥, 𝑦). different values of 18/118 Interpolating polynomials Motivation Interpolation Temperature Linear interpolation can get good results with a large number of points: Fraction benzene 19/118 Interpolating polynomials Motivation Interpolation Temperature But increasing interpolation order can yield dramatic improvements: Fraction benzene 20/118 Interpolating polynomials Motivation Nomenclature These slides reserve the term interpolation for fitting methods where the line of fit must pass through every single point. In contrast to regression, where the line of fit does not have to pass through every single data point. 21/118 Interpolating polynomials Motivation 22/118 Noise ... Note, passing through every single point is not always appropriate Good interpolation Bad interpolation Never interpolate noisy data! Good regression Interpolating polynomials Motivation Topic outline 1. Interpolating polynomials 2. Differentiation/integration 23/118 Interpolating polynomials Motivation 24/118 Learning outcomes By the end of this topic, you should be able to: 1. Apply Lagrange, divided difference, and Newton polynomial methods to generate interpolating functions 2. Estimate interpolation error 3. Understand the derivation of spline interpolation 4. Derive and apply finite difference and Newton-Cotes methods to differentiate or integrate data 5. Understand and apply adaptive integration methods Interpolating polynomials Interpolating polynomials Interpolating polynomials 26/118 Lagrange Lagrange polynomial The Lagrange interpolating polynomial is the most intuitive of all the interpolating polynomials. Given three points 𝑃2 (𝑥) = (𝑥1 , 𝑓1 ), (𝑥2 , 𝑓2 ), (𝑥3 , 𝑓3 ): (𝑥 − 𝑥2 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥1 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) 𝑓1 + 𝑓2 + 𝑓 (𝑥1 − 𝑥2 )(𝑥1 − 𝑥3 ) (𝑥2 − 𝑥1 )(𝑥2 − 𝑥3 ) (𝑥3 − 𝑥1 )(𝑥3 − 𝑥2 ) 3 ... It looks complicated, but this equation is straightforward to derive Interpolating polynomials Interpolating polynomials 27/118 Lagrange The derivation First, 𝑃2 (𝑥) must pass through all of the 𝑓𝑖 values, so it makes sense that they would be in the final equation: 𝑃2 (𝑥) = 𝑓1 + So what can be multiplied by when 𝑥 = 𝑥1 ? 𝑓2 + 𝑓3 𝑓2 to make sure that the 𝑓2 term disappears Interpolating polynomials Interpolating polynomials 28/118 Lagrange The derivation At 𝑥 = 𝑥1 , (𝑥 − 𝑥1 ) = 0, so the 𝑓2 and 𝑓3 terms disappear: 𝑃2 (𝑥) = But what about the 𝑓1 + (𝑥 − 𝑥1 ) 𝑓1 and 𝑓3 terms for 𝑥2 ? 𝑓2 + (𝑥 − 𝑥1 ) 𝑓3 Interpolating polynomials Interpolating polynomials 29/118 Lagrange The derivation Continuing: 𝑃2 (𝑥) = (𝑥 − 𝑥2 ) 𝑓1 + (𝑥 − 𝑥1 ) 𝑓2 + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) 𝑓3 Interpolating polynomials Interpolating polynomials 30/118 Lagrange The derivation And the numerators are done: 𝑃2 (𝑥) = Plug in (𝑥 − 𝑥2 )(𝑥 − 𝑥3 ) 𝑓1 + (𝑥 − 𝑥1 )(𝑥 − 𝑥3 ) 𝑓2 + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) 𝑥 = 𝑥1 , what should the denominator be to get 𝑓1 ? 𝑓3 Interpolating polynomials Interpolating polynomials Lagrange The derivation Divide by (𝑥1 − 𝑥2 )(𝑥1 − 𝑥3 ) to make sure 𝑃2 (𝑥1 ) = 𝑓1 : 𝑃2 (𝑥) = (𝑥 − 𝑥2 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥1 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) 𝑓1 + 𝑓2 + 𝑓3 (𝑥1 − 𝑥2 )(𝑥1 − 𝑥3 ) And follow the same approach for the other terms ... 31/118 Interpolating polynomials Interpolating polynomials Lagrange The derivation Overall result is therefore: 𝑃2 (𝑥) = (𝑥 − 𝑥2 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥1 )(𝑥 − 𝑥3 ) (𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) 𝑓1 + 𝑓2 + 𝑓 (𝑥1 − 𝑥2 )(𝑥1 − 𝑥3 ) (𝑥2 − 𝑥1 )(𝑥2 − 𝑥3 ) (𝑥3 − 𝑥1 )(𝑥3 − 𝑥2 ) 3 As illustrated before ... 32/118 Interpolating polynomials Interpolating polynomials Lagrange General form The overall equation can be concisely expressed as: 𝑥 − 𝑥𝑗 ⎞ ⎛ ⎟ 𝑓𝑖 𝑃𝑛 (𝑥) = ∑ ⎜∏ 𝑥 − 𝑥 𝑖 𝑗 𝑖=1 ⎝ 𝑗≠𝑖 ⎠ 𝑛+1 33/118 Interpolating polynomials Interpolating polynomials 34/118 Lagrange Polynomial evaluation Once a polynomial is constructed using a set of points, the denominator and the 𝑓𝑖 values do not change: 𝑃2 (𝑥) = (𝑥 − 𝑥2 )(𝑥 − 𝑥3 )𝑎1 + (𝑥 − 𝑥1 )(𝑥 − 𝑥3 )𝑎2 + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 )𝑎3 So it is more efficient to calculate the constant values of for a set of points and then store them for later use. ... But adding an extra point requires a lot of extra math 𝑎1 , 𝑎2 , 𝑎3 once Interpolating polynomials Interpolating polynomials Lagrange Example Problem 1 Interpolate the following T-xy temperature data as a function of fraction benzene using a Lagrange polynomial: i Temperature (°C) Fraction benzene 1 80.10 1.00 2 95.36 0.39 3 110.62 0.00 35/118 Exumple Tefx ) x *> ( x X 3 )
L
R 凶)
—
(x *) ( x ×7) fr
–
–
* *^) < L - x - , ⽕器 * - f — × ) 3 - t , 还 wn, + xx 0 0 ) )xc . y ( Lx t ⼀ - " . 3 )1 - c ) 9 ((1 0) t - JM L JC . } ) - 132 01 399 68 ^+ 284 , : - . = 13 66 x . ^ - . + 44 18x + 110 62 . . - c 世≈ qi 0 57 91 399 .68 . + t ) ω L ) ( 0 . 0 . 39 ) 391 . 96 ) ×+ 110 . 62 - Interpolating polynomials Interpolating polynomials Divided difference Divided difference polynomial The ideal method for polynomial construction should enable flexible calculation of different order polynomials using the same set of points. One way to do this is with divided difference polynomials. So, what is this divided difference? 36/118 Interpolating polynomials Interpolating polynomials Divided difference Divided difference A divided difference of order 1, the two points [ 𝑓1 , 𝑓2 ], is effectively the slope between (𝑥1 , 𝑓1 ), (𝑥2 , 𝑓2 ): [ 𝑓1 , 𝑓 2 ] = And a divided difference of order 2, two slopes between three points: [ 𝑓1 , 𝑓 2 , 𝑓 3 ] = 𝑓2 − 𝑓 1 𝑥2 − 𝑥 1 [ 𝑓1 , 𝑓2 , 𝑓3 ], is like a “slope” of the [ 𝑓2 , 𝑓 3 ] − [ 𝑓 1 , 𝑓 2 ] = 𝑥3 − 𝑥 1 𝑓3 − 𝑓 2 𝑓2 − 𝑓1 − 𝑥3 −𝑥2 𝑥2 −𝑥1 𝑥3 − 𝑥 1 37/118 Interpolating polynomials Interpolating polynomials 38/118 Divided difference Divided difference polynomial t ] 上 fi , : a* F, 涎 tniIr A polynomial of order 𝑛 can then be built using: 鸒 ! 𝑃𝑛 (𝑥) = [ 𝑓𝑖 ] + (𝑥 − 𝑥𝑖 )[ 𝑓𝑖 , 𝑓𝑖+1 ] + (𝑥 − 𝑥𝑖 )(𝑥 − 𝑥𝑖+1 )[ 𝑓𝑖 , 𝑓𝑖+1 , 𝑓𝑖+2 ] + ... Given three points would be: (𝑥1 , 𝑓1 ), (𝑥2 , 𝑓2 ), (𝑥3 , 𝑓3 ), a second order polynomial 𝑃2 (𝑥) = [ 𝑓1 ] + (𝑥 − 𝑥1 )[ 𝑓1 , 𝑓2 ] + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 )[ 𝑓1 , 𝑓2 , 𝑓3 ] Interpolating polynomials Interpolating polynomials Divided difference Divided difference polynomial 𝑃2 (𝑥) = [ 𝑓1 ] + (𝑥 − 𝑥1 )[ 𝑓1 , 𝑓2 ] + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 )[ 𝑓1 , 𝑓2 , 𝑓3 ] This may not look very convenient, but recall that: [ 𝑓1 , 𝑓 2 , 𝑓 3 ] = [ 𝑓2 , 𝑓 3 ] − [ 𝑓 1 , 𝑓 2 ] 𝑥3 − 𝑥 1 So polynomial order can be increased with one new term. 39/118 Interpolating polynomials Interpolating polynomials Divided difference Table of differences The divided difference polynomial can be built up using a table of pre-computed divided differences, with the total number of terms selected based on the required precision. This is a major improvement over the Lagrange polynomial, which must be almost entirely rebuilt if it is necessary to change the number of points used in the interpolation. 40/118 Interpolating polynomials Interpolating polynomials Divided difference Shorthand notation Since the square bracket notation [ 𝑓1 , 𝑓2 , 𝑓3 , ...] can quickly grow very long, it is common to use a shorthand notation: [ 𝑓1 ] = 𝑓1(0) = 𝑓1 [ 𝑓1 , 𝑓2 ] = 𝑓1(1) = 𝑓2(0) − 𝑓1(0) 𝑓2 − 𝑓 1 = 𝑥2 − 𝑥 1 𝑥2 − 𝑥 1 3 2 (1) (1) − 𝑥2−𝑥1 𝑓 − 𝑓 𝑥 −𝑥 2 1 1 [ 𝑓1 , 𝑓2 , 𝑓3 ] = 𝑓1(2) = 2 = 3 2 𝑥3 − 𝑥 1 𝑥3 − 𝑥 1 𝑓 −𝑓 𝑓 −𝑓 41/118 Interpolating polynomials Interpolating polynomials Divided difference General form So the overall equation can be concisely expressed as: 𝑛 𝑖 ⎜∏(𝑥 − 𝑥𝑗 )⎞ ⎟ 𝑓1(𝑖) 𝑃𝑛 (𝑥) = ∑ ⎛ 𝑖=0 ⎝ 𝑗=1 ⎠ 42/118 Interpolating polynomials Interpolating polynomials Divided difference Verification Problem 2 Take the definition of a divided difference polynomial: 𝑃𝑛 (𝑥) = [ 𝑓1 ] + (𝑥 − 𝑥1 )[ 𝑓1 , 𝑓2 ] + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 )[ 𝑓1 , 𝑓2 , 𝑓3 ] + ... and show that for 𝑥 = 𝑥2 , 𝑃𝑛 (𝑥2 ) = 𝑓2 . 43/118 Interpolating polynomials Interpolating polynomials Divided difference Example Problem 3 Interpolate the following T-xy temperature data as a function of fraction benzene using a divided difference polynomial: i Temperature (°C) Fraction benzene 1 80.10 1.00 2 95.36 0.39 3 110.62 0.00 44/118 " fi extid Px) : ) : t f," ) , ) t , 8. 1 : , ( - xp ( xx by alaulaticndiviliddiffereue start f x ×2 - 95 36 , : 87 (c - . - 25 . 0 z x, 0. 39 - 1 ti fi " fie ) i or ≈ 10. 10 1 0 2 0 } f . 14 ( ! 1 , y95 .36 - 34 - 13 … ) 011 .62 t ②] I ! .3 -25 " - 。 f … . … " , ⼀ , × 3 - ×, t 25 39 13 - . : o = if we ( - 02 , 1 4 11 , Ist wave P a) a 2 . ⼀ " 80 , 10 + order xL) L - ashimtc , 25 c 2 ) . udorder Prx) : 87 , 0 t L) ( - 2 i )+ L )( x- . n) 0 3 al 4 . 1) l Ervor atuate : )1 x 14 . ) 1 0 ( PRx ⼀⼼ 11 ) coupave Expoud to appractins . 4 . llx e ito 39 . tLi 8ol on 'xel 14 0 = 14 11x -44 ' . Tese of i 3) . . . 3Y , 6)x 4 ( 10 62 . RCC 39) . : 95 36 . , 2 Interpolating polynomials Interpolating polynomials Newton Newton polynomial Recall that one of the trade-offs considered when comparing numerical methods is robustness/generality vs. speed. A simplifying assumption can often be used to generate a simplified and potentially better version of an existing algorithm. One such simplification that is often made with interpolation is to choose equidistant points, making it possible to get rid of (𝑥𝑖+1 − 𝑥𝑖 ) terms in the divided difference formulation and replace them with a constant step size of ℎ. all 45/118 Interpolating polynomials Interpolating polynomials Newton Finite differences Since there is no need to keep track of all the Δ(0) 𝑓1 = 𝑓1 𝑥 terms, remove them: Δ(1) 𝑓1 = Δ(0) 𝑓2 − Δ(0) 𝑓1 = 𝑓2 − 𝑓1 Δ(2) 𝑓1 = Δ(1) 𝑓2 − Δ(1) 𝑓1 = ( 𝑓3 − 𝑓2 ) − ( 𝑓2 − 𝑓1 ) = 𝑓 3 − 2 𝑓2 + 𝑓 1 46/118 Interpolating polynomials Interpolating polynomials 47/118 Newton Newton polynomial (𝑥 − 𝑥𝑖 ) terms are replaced by the number of “steps” 𝑥 is from the first point 𝑥1 , where each step ℎ = 𝑥2 − 𝑥1 : All the 𝑠= 𝑥 − 𝑥1 ℎ Allowing the Newton polynomial to be defined as: 𝑃𝑛 (𝑠) = Δ(0) 𝑓1 + 𝑠Δ(1) 𝑓1 + 𝑠(𝑠 − 1) (2) 𝑠(𝑠 − 1)(𝑠 − 2) (3) Δ 𝑓1 + Δ 𝑓1 + ... 2! 3! Which is a considerable simplification over the divided difference approach that can be used when the data is on an equally spaced grid. Interpolating polynomials Interpolating polynomials Newton The kernel of many applications The Newton polynomial (finite difference) formulation forms the basis of topic 4. Start getting used to it now! 48/118 Interpolating polynomials Interpolating polynomials Newton Verification Problem 4 Substitute 𝑠 = (𝑥−𝑥1 )/ℎ into the divided difference equation: 𝑃2 (𝑥) = 𝑓1(0) + (𝑥 − 𝑥1 ) 𝑓1(1) + (𝑥 − 𝑥1 )(𝑥 − 𝑥2 ) 𝑓1(2) and show that it is equivalent to the Newton polynomial when 𝑥2 − 𝑥 1 = 𝑥 3 − 𝑥 2 = ℎ . 49/118 Interpolating polynomials Interpolating polynomials Newton Example Problem 5 Using the vapour pressure example data from the beginning of this topic, generate a 2nd order Newton polynomial from 𝑥1 = 94 to calculate 𝑃∗ (96.5). 50/118 Interpolating polynomials Interpolating polynomials Newton Different differences Newton polynomials were introduced with forward differences: Δ 𝑓𝑖 = 𝑓𝑖+1 − 𝑓𝑖 But backward differences are not much different: ∇ 𝑓𝑖 = 𝑓𝑖 − 𝑓𝑖−1 𝑃𝑛 (𝑠) = ∇(0) 𝑓𝑖 + 𝑠∇(1) 𝑓𝑖 + 𝑠= 𝑥 − 𝑥𝑖 ℎ 𝑠(𝑠 + 1) (2) 𝑠(𝑠 + 1)(𝑠 + 2) (3) ∇ 𝑓𝑖 + ∇ 𝑓𝑖 + ... 2! 3! 51/118 Interpolating polynomials Interpolating polynomials Spline interpolation Global vs. local fit 𝑛 points requires a polynomial of order 𝑛 − 1. However, working with high order polynomials (𝑛 > 20) can lead
Fitting a
to large round-off errors under some circumstances.
Instead of using a single high order polynomial to perform a
global fit, one solution is to combine many low order polynomial to perform a local fit.
52/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
Basic setup
Consider a small set of points to be fit:
(𝑥2 , 𝑓2 )
(𝑥1 , 𝑓1 )
(𝑥3 , 𝑓3 )
(𝑥4 , 𝑓4 )
53/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
1st order
Unique 1st order splines can be generated by forcing each spline to
touch neighbouring points.
(𝑥2 , 𝑓2 )
(𝑥1 , 𝑓1 )
(𝑥3 , 𝑓3 )
(𝑥4 , 𝑓4 )
54/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
1st order
𝑛 points can be fit by 𝑛 − 1 polynomials, so:
Parameters:
2 parameters per
𝑛 − 1 polynomial = 2(𝑛 − 1)
Available information:
Each of the
𝑛 − 1 polynomials touches 2 points
2(𝑛 − 1) − 2(𝑛 − 1) = 0
(sufficient to calculate all parameters)
55/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
2nd order
However, these conditions are not sufficient for 2nd order splines.
(𝑥2 , 𝑓2 )
(𝑥1 , 𝑓1 )
(𝑥3 , 𝑓3 )
(𝑥4 , 𝑓4 )
Data points alone do not specify a unique set of polynomials.
56/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
2nd order
𝑛 points can be fit by 𝑛 − 1 polynomials, so:
Parameters:
3 parameters per
𝑛 − 1 polynomial = 3(𝑛 − 1)
Available information:
Each of the
𝑛 − 1 polynomials touches 2 points
3(𝑛 − 1) − 2(𝑛 − 1) = 𝑛 − 1
(not enough to calculate parameters)
57/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
2nd order
Intuitively, it should be clear that this is a bad fit:
(𝑥2 , 𝑓2 )
(𝑥3 , 𝑓3 )
(𝑥1 , 𝑓1 )
But how can a good fit be described mathematically?
(𝑥4 , 𝑓4 )
58/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
2nd order
In general, a “good” fit is smooth
(𝑥2 , 𝑓2 )
(𝑥3 , 𝑓3 )
(𝑥1 , 𝑓1 )
And smoothness is defined by continuous derivatives
(𝑥4 , 𝑓4 )
…
59/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
2nd order
𝑛 points can be fit by 𝑛 − 1 polynomials, so:
Parameters:
3 parameters per
𝑛 − 1 polynomial = 3(𝑛 − 1)
Available information:
Each of the
𝑛 − 1 polynomials touches 2 points
Neighbouring polynomials have same 1st derivatives
at
𝑛 − 2 points
3(𝑛 − 1) − 2(𝑛 − 1) − (𝑛 − 2) = 1
(still not enough to calculate parameters)
60/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
2nd order
𝑛 points can be fit by 𝑛 − 1 polynomials, so:
Parameters:
3 parameters per
𝑛 − 1 polynomial = 3(𝑛 − 1)
Available information:
Each of the
𝑛 − 1 polynomials touches 2 points
Neighbouring polynomials have same 1st derivatives
at
𝑛 − 2 points
Force 1st or 2nd derivative to zero at 1 point
3(𝑛 − 1) − 2(𝑛 − 1) − (𝑛 − 2) − 1 = 0
(sufficient to calculate all parameters)
61/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
Higher orders
Increasing polynomial orders require more conditions
Typically, continuity conditions are extended to all higher orders
𝑛 − 2 points and have the
same 1st, 2nd, 3rd, etc., derivatives at 𝑛 − 2 points
neighbouring polynomials connect at
Extra conditions are provided by setting the derivatives of edge
polynomials to zero
62/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
In practice
Cubic polynomials are usually smooth enough for most applications.
Spline definition may not be applicable for special cases like
discontinuities. It may be necessary to implement custom spline
equation with special definitions or make use of advanced packages.
63/118
Interpolating polynomials
Interpolating polynomials
Spline interpolation
Some math
Problem 6
Derive a set of equations that can be used to calculate spline
parameters for a quadratic spline across
𝑛 points.
64/118
Example 6
.
sit
1
1
≈
i
∅
↓
xo
six]
ien
tie
is
tfinfisfi ,
di
:
need
we
Xi+ 1
Con
–
X0
≈
e
fod
eo
‘
bilx xi” t
bi ci
,
a,
( x xu)
C,
–
–
lne
fuall
iegnae ,
li
liunty candition
on
the
left
silxi)
xis
fl aitbi xix tcicxi
a
i
,
=
②
si (
)
fi):
ai
fi)
:
FEce
③
xit
Ist
–
,
=fi)
tbi ( xie
–
fls ebihiecihi
derivaeve
Si ‘(xlt 1
]
=Sit
xi
)
)
xi tCi ( xier
eo
1 ( Xit )
be
–
‘
equal
to
heighbauring
playoamials
,
bf t
2ci
hie
xity
lxour -xo
)
,
t
2 cieitxun –
( bit
靠 Jhi
fiur fitbihi t
=
fien fi
2
:
bo tbie
hi
↑
T
known
bitbit
,
unbow
2
=
)
州
hi
f)
2
bitr
bitlt
:
hitll 了
⼀
1
1
1
co
Δ
(
10
0
…
… …
–
last
our
”
s n-
[
| ]
…
…
…
carditichiwillsetund dorvaeve
uelaie
pilynnut
1 ( xn
)
2
: car
=
0
fuituy tbue hurs tCut
hnl ,
s bas
器
:
t
π
recap
!
ifwe
Cr
have
we
.
n
–
polynennele
l
.
I
l
I
[
?
=
pones
a
[
a =
b
hoave
ai
fiso
=
,
bi
is
fom
salved
asyslem
of
lnen
,
equalin
.
壑
i 鬭
囖
1
1
0
1
0
0
⼊
7
O
G
O
(
0
0
–
!
00
0
…
O
…
…
…
Cr
Eul ia
–
assumed
as
side
A
from biel bofer iilonmz
solue
2s
–
.
do
How
.
,
we
this
me
,
c
Λ
0
⇌
⼀
「
O
Ψ
⼄
f(
n
.
6
) ?
fiud si
xi = U
we
me
8
,
whochpoly namial
V= 2
cican 6
)
a
thi luo
su – ar
6
–
)
t
–
Interpolating polynomials
Interpolating polynomials
Method comparison
Algorithm selection
Algorithm selection depends on a large number of different factors, the
following is just a suggestion:
If you are dealing with more than 20 points: splines
If your data is on an evenly spaced grid: Newton polynomial
Otherwise: divided difference polynomial
The Lagrange interpolating polynomial was introduced because it is
easy to understand.
Its only advantage over the divided difference
method is that it is easier to program.
65/118
Interpolating polynomials
Interpolating polynomials
Method comparison
Caution
All interpolating techniques define a function that passes through
every single observed data point. You should not do this if your
data is noisy!
66/118
Interpolating polynomials
Interpolating polynomials
Method comparison
Further reading
Chapter 4 in the book covers everything in these slides and
more. Sections 4.1-4.2 serve as a good introduction to interpolation, sections 4.4-4.6 cover the three major global interpolation methods, and section 4.9 describes cubic splines.
67/118
Derivatives and integrals
Interpolating polynomials
Derivatives and integrals
Background
Derivatives and integrals
Beyond the direct application of interpolating polynomials to
“read between the table rows”, interpolating polynomials are
also used to develop methods for techniques such as differentiation and integration.
69/118
Interpolating polynomials
Derivatives and integrals
Background
Mathematical interpretation
Derivatives and integrals are often conceptualized geometrically
a derivative is the slope of the tangent line at a given point
an integral is the area under a section of a curve
This conceptualization can be used to determine how to calculate a
derivative or integral
but it does not explain why either is useful
70/118
Interpolating polynomials
Derivatives and integrals
Background
Physical interpretation
A derivative corresponds to a rate of change, while an integral
is an accumulation (typically over space or time).
These values may be useful on their own or as a means of calculating others.
71/118
Interpolating polynomials
Derivatives and integrals
72/118
Background
Rate of change
The flow rate (
𝑄) out of the pipe
is directly related to the rate of
change of fluid height in the ves-
H
Q
sel (
𝐻 ):
𝑄=𝐴
vessel
d𝐻
d𝑡
Interpolating polynomials
Derivatives and integrals
Background
Rate of change
Note that the above expression is not theoretical (we do not
necessarily know the rate of change in fluid height in advance).
However, fluid height is considerably easier to measure than
flow rate.
73/118
Interpolating polynomials
Derivatives and integrals
74/118
Background
Accumulation
Q
The overall fluid height in the ves-
𝐻 ) is directly related to the
flow rate (𝑄) out of the pipe over
sel (
time:
H
𝐻(𝑡) = 𝐻∣𝑡=𝑡 +
0
𝐴
1
vessel
𝑡
∫ 𝑄(𝑡)dt
𝑡0
Interpolating polynomials
Derivatives and integrals
Background
Accumulation
If we know the flow rate of fluid going into the pipe, we can
predict the height of the fluid in the vessel in advance.
75/118
Interpolating polynomials
Derivatives and integrals
Background
Another example
Derivatives and integrals are also used for monitoring reactions
Differentiation to calculate reaction rates from concentration data
𝑟=
d𝐶
d𝑡
Integration to calculate concentration from reaction rate
𝑡
𝐶(𝑡) = 𝐶0 + ∫ 𝑟dt
𝑡0
76/118
Interpolating polynomials
Derivatives and integrals
77/118
Background
Numerical vs. analytical
It is possible to solve simple derivatives and integrals analytically
typically using pattern recognition (to identify useful variable
substitutions or techniques like integration by parts)
Numerical methods are often simpler and allow the solution of
problems that cannot be solved analytically
Numerical methods can also be used to work with data (a series of
points) where the true function is unknown
Interpolating polynomials
Derivatives and integrals
Background
Simplifying equations
The general differentiation/integration technique is always the
same — interpolate input data into a function and find the derivative/integral of the function as normal.
However, this process can be derived into a simple algebraic
equation for common problems.
78/118
Interpolating polynomials
Derivatives and integrals
79/118
Differentiating points
Differentiation
If the input data points are evenly spaced, a derivative equation can be
calculated by fitting a Newton polynomial to points around
(𝑥𝑖 , 𝑓𝑖 )
(𝑥𝑖 , 𝑓𝑖 ).
Interpolating polynomials
Derivatives and integrals
Differentiating points
1st order
Recall the 1st order Newton polynomial:
𝑃1 (𝑠) = 𝑓𝑖 + 𝑠Δ(1) 𝑓𝑖
𝑠=
Where
𝑥 − 𝑥𝑖
ℎ
ℎ is the step size between 𝑥𝑖 , 𝑥𝑖+1 .
Calculating
d𝑃1 (𝑠)
requires using the chain rule:
d𝑥
d𝑃1 (𝑠) d𝑃1 (𝑠) d𝑠
=
d𝑥
d𝑠 d𝑥
80/118
Interpolating polynomials
Derivatives and integrals
Differentiating points
1st order
The two terms:
d𝑃1 (𝑠)
d
= ( 𝑓𝑖 + 𝑠Δ(1) 𝑓𝑖 )
d𝑠
d𝑠
= Δ(1) 𝑓𝑖 = 𝑓𝑖+1 − 𝑓𝑖
d𝑠
1
=
d𝑥 ℎ
Together:
d𝑃1 (𝑠)
𝑓𝑖+1 − 𝑓𝑖
=
d𝑥
ℎ
81/118
Interpolating polynomials
Derivatives and integrals
Differentiating points
Visually
Comparing the estimate and true value:
(𝑥𝑖 , 𝑓𝑖 )
The estimate is good if the function is smooth and the step size is small.
82/118
Interpolating polynomials
Derivatives and integrals
Differentiating points
Higher orders
( 𝑓𝑖+1 − 𝑓𝑖 )/ℎ is a first-order finite difference. The same process
can be used to define higher order finite difference equations.
83/118
Interpolating polynomials
Derivatives and integrals
Differentiating points
Higher orders
Problem 7
Derive a finite difference equation for the derivative at a point
𝑥𝑖 using a 2nd order Newton Polynomial equation.
84/118
Example
?
,
Phc
)
–
+
fc )
X
–
e
scfi
C
”
0
!
fc
)
Xi
H答装
:
ft
as
u
d
+
< e " , ÷ fc fit fi tI ( fiux - - 2 fie tfi ) , ) : cfit dd fo "ft : : 00) befere 2 “ : 6 fitn - h , fi o fier fie - fier ) fo , : - : fitv afie 1 tfi . - devvatve porhe 2 3 pozut derivatwa ←→←→ G ; ⾔ 。 … ⼀ ⼀ … ⼀ ⼀ !∴ … ②⼀ ⼀ … ⼀ 「 「 Δ ? Sco Three differae ae sio → our - fi fie afiei efo fier : at fon equations fit fi 4 + ti sal nfo - etie ae l - whrch is sel vs the most perfend accawate ( powes ! at clear to S= 13 ) … gul 5 : 2 Interpolating polynomials Derivatives and integrals 85/118 Integrating points Integration If the input data points are evenly spaced, an integral equation can be calculated by fitting a Newton polynomial to points around (𝑥𝑖 , 𝑓𝑖 ) (𝑥𝑖 , 𝑓𝑖 ). Interpolating polynomials Derivatives and integrals Integrating points 1st order Recall the same 1st order Newton polynomial: 𝑃1 (𝑠) = 𝑓𝑖 + 𝑠Δ(1) 𝑓𝑖 𝑠= Where 𝑥 − 𝑥𝑖 ℎ ℎ is the step size between 𝑥𝑖 , 𝑥𝑖+1 . Use the chain rule to go from 𝑥 to 𝑠: 𝑥𝑖+1 𝐼 = ∫ 𝑃1 (𝑠)d𝑥 𝑥𝑖 86/118 Interpolating polynomials Derivatives and integrals Integrating points 1st order First, convert 𝑥 to 𝑠 d𝑠 1 = d𝑥 ℎ d𝑥 = ℎ d𝑠 And do not forget about the limits of integration: 1 𝐼 = ∫ 𝑃1 (𝑠)ℎ d𝑠 0 With a 1st order estimate, use only two points 𝑥𝑖 , 𝑥𝑖+1 , or 𝑠 = 0, 1. 87/118 Interpolating polynomials Derivatives and integrals Integrating points 1st order 1 𝐼 = ∫ 𝑃1 (𝑠)ℎ d𝑠 0 1 = ℎ ∫ 𝑓𝑖 + 𝑠Δ(1) 𝑓𝑖 d𝑠 0 ∣ 𝑠2 (1) = ℎ ( 𝑓𝑖 𝑠 + Δ 𝑓𝑖 ) ∣∣ 2 ∣ 1 = ℎ ( 𝑓𝑖 + Δ(1) 𝑓𝑖 ) 2 = ℎ ( 𝑓𝑖 + = 𝑓𝑖+1 − 𝑓𝑖 ) 2 ℎ( 𝑓𝑖 + 𝑓𝑖+1 ) 2 1 𝑠=0 88/118 Interpolating polynomials Derivatives and integrals Integrating points Visually Comparing the estimate and true value: (𝑥𝑖 , 𝑓𝑖 ) The estimate is good if the function is smooth and the step size is small. 89/118 Interpolating polynomials Derivatives and integrals Integrating points Global estimate A single integral is typically calculated over all the points: 90/118 Interpolating polynomials Derivatives and integrals Integrating points Global estimate To integrate over 𝑛 points, sum over all 𝑛 − 1 sub-intervals: 𝑛−1 𝐼 = ∑ 𝐼𝑖 𝑖=1 But the formula allows some simplification: ℎ( 𝑓𝑖 + 𝑓𝑖+1 ) 2 𝑖=1 𝑛−1 𝐼=∑ 𝑛−1 ℎ = ( 𝑓1 + ∑ 2 𝑓𝑖 + 𝑓 𝑛 ) 2 𝑖=2 91/118 Interpolating polynomials Derivatives and integrals Integrating points Higher orders ℎ( 𝑓𝑖+1 + 𝑓𝑖 )/2 is a first-order Newton-Cotes equation (other- wise known as the trapezoid rule). The same process can be used to define higher order Newton-Cotes equations (such as the Simpson’s rule or Simpson’s 3/8 rule). 92/118 Interpolating polynomials Derivatives and integrals 93/118 Integrating points Higher orders Problem 8 Derive a Newton-Cotes equation for a global integral over points using a 2nd order Newton Polynomial equation. 𝑛 Interpolating polynomials Derivatives and integrals Integrating points Example Problem 9 Given the following set of data, use a 2nd order Newton polynomial to estimate the derivative at integral between 𝑥 = 0, 1 𝑥 = 0.5 and the global x 0.00 0.25 0.50 0.75 1.00 y 1.00 1.28 1.65 2.12 2.72 Compare the estimates to the true values given that 𝑦 = 𝑒𝑥 . 94/118 Interpolating polynomials Derivatives and integrals 95/118 Integrating functions Integrating functions When trying to find a derivative or integral from a set of points, the underlying function is unknown However, we may also be interested in finding the integral of a known function more common for integrals as derivatives are relatively easy to find analytically Recall the problem from topic 1 that could not be solved analytically: 𝑏 ∫ 𝑒−𝑥 d𝑥 𝑎 2 Interpolating polynomials Derivatives and integrals Integrating functions Integrating functions Numerically integrating functions is similar to points. Suggest a method to find the previously mentioned integral: 𝑏 ∫ 𝑒−𝑥 d𝑥 𝑎 2 96/118 Interpolating polynomials Derivatives and integrals 97/118 Integrating functions Beyond points If the function is known, it is easy to calculate function values at a series of points and directly apply the Newton-Cotes equations step size or order can be varied to estimate accuracy However, knowing the function enables a better estimate than possible with a fixed set of points Gauss-Legendre quadrature can be used to optimize which points to choose (outside the scope of the course) Richardson extrapolation can be used to generate a more accurate estimate from two estimates with different step sizes Interpolating polynomials Derivatives and integrals Integrating functions Beyond points By leveraging how truncation error works, so-called “adaptive” methods are able to get better estimates with less work. 98/118 Interpolating polynomials Derivatives and integrals Integrating functions Taylor series and truncation The Taylor series offers a convenient introduction to truncation error. Recall the Taylor series formulation: 𝑓 (𝑥) = 𝑓 (𝑥0 ) + 𝑓 ′ (𝑥0 )(𝑥 − 𝑥0 ) + 1 𝑓 ″ (𝑥0 )(𝑥 − 𝑥0 )2 + ... 2! Where the ellipses can be represented by a remainder term: 𝑓 (𝑥) = 𝑓 (𝑥0 ) + 𝑓 ′ (𝑥0 )(𝑥 − 𝑥0 ) + 1 𝑓 ″ (𝑥0 )(𝑥 − 𝑥0 )2 + 𝑅3 2! 99/118 Interpolating polynomials Derivatives and integrals 100/118 Integrating functions Taylor series and truncation It can be shown that the remainder can be defined as: 𝑓 (𝑥) = 𝑓 (𝑥0 ) + 𝑓 ′ (𝑥0 )(𝑥 − 𝑥0 ) + where 1 1 𝑓 ″ (𝑥0 )(𝑥 − 𝑥0 )2 + 𝑓 ‴ (𝜉 )(𝑥 − 𝑥0 )3 2! 3! 𝜉 is some value between 𝑥0 , 𝑥. Essentially, the maximum possible error is: 𝑅3 ≤ max ( 1 𝑓 ‴ (𝜉 )(𝑥 − 𝑥0 )3 ) 3! 𝜉 ∈ (𝑥0 , 𝑥) Interpolating polynomials Derivatives and integrals 101/118 Integrating functions Big O notation ℎ) from interpolation. It can be applied here as well to rewrite the Taylor series as an estimate of 𝑓 (𝑥 + ℎ) based on the estimate at 𝑓 (𝑥), a step ℎ away: Recall the idea of a step size ( Since ℎ2 ℎ3 𝑓 (𝑥 + ℎ) = 𝑓 (𝑥) + ℎ 𝑓 ′ (𝑥) + 𝑓 ″ (𝑥) + 𝑓 ‴ (𝜉 ) 2! 3! 𝑓 ‴ (𝜉 )/3! is a scalar term that modifies ℎ3 , it is simplified as: ℎ2 𝑓 (𝑥 + ℎ) = 𝑓 (𝑥) + ℎ 𝑓 ′ (𝑥) + 𝑓 ″ (𝑥) + 𝑂(ℎ3 ) 2! which just means that if the step size is decreased 2-fold, for example, the error will decrease 8-fold (or 23 ). Interpolating polynomials Derivatives and integrals Integrating functions Why do we care? Although finding some arbitrary relation for the error may seem pointless, understanding how error changes as a function of ℎ enables more accurate estimates. The following slides will demonstrate how to get this better estimate using a simple Taylor series expansion. 102/118 Interpolating polynomials Derivatives and integrals Integrating functions Richardson extrapolation Take a 1st order Taylor series expansion: ℎ2 ℎ3 𝑓 (𝑥 + ℎ) = 𝑓 (𝑥) + ℎ 𝑓 ′ (𝑥) + 𝑓 ″ (𝑥) + 𝑓 ‴ (𝑥) + ... 2! 3! And use it to estimate the derivative 𝑓 ′ (𝑥): 𝑓 (𝑥 + ℎ) − 𝑓 (𝑥) ℎ ℎ2 𝑓 ′ (𝑥) = − 𝑓 ″ (𝑥) − 𝑓 ‴ (𝑥) − ... 2! 3! ℎ 103/118 Interpolating polynomials Derivatives and integrals Integrating functions Richardson extrapolation The focus here is on ℎ so rename the 𝑓 ″ (𝑥) and 𝑓 ‴ (𝑥) terms: 𝑓 ′ (𝑥) = 𝑓 (𝑥 + ℎ) − 𝑓 (𝑥) − ℎ𝑟2 − ℎ2 𝑟3 − ... ℎ In fact, rename the other terms as well: 𝑓 ′ (𝑥) = 𝐴(ℎ) − ℎ𝑟2 − ℎ2 𝑟3 − ... where the sole focus is on ℎ. 104/118 Interpolating polynomials Derivatives and integrals Integrating functions Richardson extrapolation Starting with: 𝑓 ′ (𝑥) = 𝐴(ℎ) − ℎ𝑟2 − ℎ2 𝑟3 − ... What happens if the step size is halved? ℎ ℎ ℎ2 𝑓 ′ (𝑥) = 𝐴 ( ) − 𝑟2 − 𝑟3 − ... 2 2 4 Remember that all the other parts of the equation stay the same. 105/118 Interpolating polynomials Derivatives and integrals Integrating functions Richardson extrapolation Compare the two estimates: ℎ ℎ ℎ2 𝑓 ′ (𝑥) = 𝐴 ( ) − 𝑟2 − 𝑟3 − ... 2 2 4 𝑓 ′ (𝑥) = 𝐴(ℎ) − ℎ𝑟2 − ℎ2 𝑟3 − ... Now multiply the first equation by 2 and subtract the second: 𝑓 ′ (𝑥) = 2𝐴(ℎ/2) − 𝐴(ℎ) − (ℎ2 /2 + ℎ2 )𝑟3 − ... 106/118 Interpolating polynomials Derivatives and integrals Integrating functions Richardson extrapolation So what happened? There is now another estimate for 𝑓 ′ (𝑥): 𝑓 ′ (𝑥) = 2𝐴(ℎ/2) − 𝐴(ℎ) − (ℎ2 /2 + ℎ2 )𝑟3 − ... But with reduced error: 𝑓 ′ (𝑥) = 2𝐴(ℎ/2) − 𝐴(ℎ) + 𝑂(ℎ2 ) The error term went from 𝑂(ℎ) to 𝑂(ℎ2 )! 107/118 Interpolating polynomials Derivatives and integrals Integrating functions More generally Knowing the error structure makes it possible to combine multiple poor estimates to increase accuracy! 108/118 Interpolating polynomials Derivatives and integrals 109/118 Integrating functions Romberg integration Richardson extrapolation can be used equally well for calculating both derivatives and integrals The most common application is known as Romberg integration 1. a function interval is integrated using the trapezoid rule 2. the interval is then split in two and a new estimate is calculated 3. the two estimates are combined to get a higher order estimate 4. further subintervals are split until the approximation converges Interpolating polynomials Derivatives and integrals 110/118 Integrating functions Romberg integration The trapezoid rule error is 𝑂(ℎ2 ) and based on the structure of that er- 𝑂(ℎ2 ) estimates with different step sizes can be used to gen4 erate an 𝑂(ℎ ) estimate: ror, two 𝐴4 = 4𝐴2 (ℎ/2) − 𝐴2 (ℎ) 3 More generally, two estimates with order an improved estimate of order 𝑘 + 2: 𝑘 can be combined to generate 2𝑘 𝐴𝑘 (ℎ/2) − 𝐴𝑘 (ℎ) 𝐴𝑘+2 = 2𝑘 − 1 Interpolating polynomials Derivatives and integrals Integrating functions Example Problem 10 Perform 3 iterations of Romberg integration to find the inte- 𝑒𝑥 from 0 to 1. Compare the estimates to the true value 1 0 of 𝑒 − 𝑒 = 1.7182. gral of 111/118 Interpolating polynomials Derivatives and integrals 112/118 Method comparison Algorithm selection Algorithm selection depends on a large number of different factors, the following is just a suggestion: Derivatives from data on an even grid: finite-difference equations Integrals from data on an even grid: Newton-Cotes equations Calculating the integral of a known function: Romberg integration Note that these options do not cover all possible cases. The general process is to fit a function to your data and take the derivative/integral of that function. Interpolating polynomials Derivatives and integrals Method comparison Further reading The book separates differentiation and integration into Chapters 5 and 6, whereas the slides try to highlight the overarching similarities. Section 5.1 and 6.1 offer separate introduc- tions for the two topics. Finite difference and Newton-Coates equations are covered in 5.3 and 6.3. For adaptive integration, see sections 6.4 and 6.5. 113/118 Final thoughts Interpolating polynomials Final thoughts Interpolation is the core of numerical methods Analytical methods are generally continuous whereas most numerical methods work on discrete data. Interpolation fundamentally serves to convert discrete data into functional form. You will come across many different numerical “rules”, but most of these are derived directly from interpolation. It is necessary to develop a good understanding of interpolation to understand the principles behind other tasks like calculating derivatives and integrals. 115/118 Interpolating polynomials Final thoughts 116/118 Interpolation algorithm selection Algorithm selection depends on a large number of different factors, the following is just a suggestion: If you are dealing with more than 20 points: splines If your data is on an evenly spaced grid: Newton polynomial Otherwise: divided difference polynomial The Lagrange interpolating polynomial was introduced because it is easy to understand. Its only advantage over the divided difference method is that it is easier to program. Interpolating polynomials Final thoughts 117/118 Differentiation/integration algorithm selection Algorithm selection depends on a large number of different factors, the following is just a suggestion: Derivatives from data on an even grid: finite-difference equations Integrals from data on an even grid: Newton-Cotes equations Calculating the integral of a known function: Romberg integration Note that these options do not cover all possible cases. The general process is to fit a function to your data and take the derivative/integral of that function. Interpolating polynomials Final thoughts Caution All interpolating techniques define a function that passes through every single observed data point. You should not do this if your data is noisy! 118/118