E7 – Prof. Alam Fall 2019, UC Berkeley
Homework Assignment 11
Interpolation, Root
Finding, Numerical
Differentiation and
Integration
This assignment will introduce you to interpolation, root finding and numerical differentia-
tion and integration in MATLAB.
The due time is 11:59am 4 Dec 2019. NO LATE HOMEWORK WILL BE AC-
CEPTED. Please upload the following files through the bCourses website:
• The function files myInterp.m, interpthermo.m, refineimage.m, damping.m MyDer.m
and mySimpson.m.
Directions to upload files can be found here.
1 Interpolating Data from an Experiment
You are doing an experiment that finds the relation between two quantities t and d. t is an
independent variable which you vary and you measure d for each such value of t in order to
obtain some readings. In this example, your independent variable t is stored in the variable
t and your dependent variable displacement is stored in the variable d. Now that you have
a set of t and d data points, you want to interpolate a curve to fit the data.
Write a function myInterp with the function declaration line
func t i on D = myInterp ( t , d ,T, opt ion )
The inputs are:
• t: an n-by-1 column vector containing the t data points.
• d: an n-by-1 column vector containing the d data points.
• T: an m-by-1 column vector containing the t instances for which an interpolation is
desired. In general, m ≥ n.
1
https:
guides.instructure.com/m/4212/l/54353-how-do-i-upload-a-file-to-my-assignment-submission
E7 – Prof. Alam Fall 2019, UC Berkeley
• option: a string, containing either ’nearest’, ’linear’, ’cubic’, ’regression’, o
’cubicpoly’.
The output is:
• D: an m-by-1 column vector containing the interpolated d data.
Use the MATLAB built-in function interp1 to do the interpolation for the options ’nearest’,
’linear’, and ’cubic’. Additionally, if the content of option is invalid, generate an e
o
’Option selection is inco
ect.’ using the e
or function.
For the options ’regression’ and ’cubicpoly’, perform regression of the data values on
a straight line and a cubic polynomial respectively. That means, you should solve for the
coefficients of the line by representing the problem in the Ax = b form, and evaluate T on
the equation representing the line. Use the backslash operator to solve for the coefficients.
You can check the following cases to test your code:
t = [ XXXXXXXXXX 0 ] ’ ;
d = [ XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX XXXXXXXXXX ] ’ ;
D1 = myInterp ( t , d , l i n s p a c e (0 , 60 , 100 ) ’ , ’ nearest ’ ) ;
D2 = myInterp ( t , d , l i n s p a c e (0 , 60 , 100 ) ’ , ’ l i n e a r ’ ) ;
D3 = myInterp ( t , d , l i n s p a c e (0 , 60 , 100 ) ’ , ’ cubic ’ ) ;
D3 = myInterp ( t , d , l i n s p a c e (0 , 60 , 100 ) ’ , ’ r e g r e s s i o n ’ ) ;
D3 = myInterp ( t , d , l i n s p a c e (0 , 60 , 100 ) ’ , ’ cubicpoly ’ ) ;
You can generate a graph of your data points (t versus d)) as well as the interpolated line
(T versus D) to visually check the output of your function.
2 Thermodynamic Property Tables
Consider a thermodynamic system that operates with superheated water vapour. In analysis,
it is common to refer to tables that contain properties of the operating fluid.
In this problem, you are going to work with a cell a
ay that stores data for internal energy
as a function of temperature. You are going to use linear interpolation to find internal energy
values for temperatures located in between those in the table.
(a) Create a function interpthermo which takes in two input arguments T1 and newT. T1
is an input cell a
ay of properties. newT is a double a
ay of new temperatures. You
function should call a subfunction interpthermo2 with the appropriate arguments and
2
E7 – Prof. Alam Fall 2019, UC Berkeley
create a new cell a
ay T2 which contains properties at new temperatures in addition
to those already in T1.
(b) Create a subfunction interpthermo2 which takes in two input arguments T1 and newT.
The input T1 is a cell a
ay of properties stored in ascending order of temperatures and
newT is a temperature which does not exist in the cell a
ay T1. The output argument
is T2, a cell a
ay with the new interpolated internal energy for temperature newT with
the same names for properties, also stored in ascending order of temperatures. You
should make use of the function interp1 in your computations.
You can use the following example to check your code:
T = { XXXXXXXXXX ; XXXXXXXXXX ; XXXXXXXXXX ; XXXXXXXXXX ; XXXXXXXXXX}
5 x 2 c e l l a
ay
T =
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
Tnew = interpthermo (T, [ XXXXXXXXXX ] )
ans =
9 x 2 c e l l a
ay
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
XXXXXXXXXX
3
E7 – Prof. Alam Fall 2019, UC Berkeley
3 Image Interpolation
In this problem you will use interpolation to refine a image. From bcourses, you can download
imdata.mat (class double), which is a 2-D a
ay data storing a image that co
esponds to a
sample image. You can show this color image with the following code
load imdata . mat
f i g u r e
imagesc ( data )
If you want to display the image in black, white, and gray, just add the following command.
colormap gray
Write a function with the following header:
f unc t i on re f ined imagedata = re f i n e image ( imagedata )
where input is the name of the image data file (string), and output is refined image data,
which is generated by halving the grid size of the image.
In the function, use built in function interp2 to create interpolated image data. Use the
”spline” method to do the 2-D interpolation.
You can use the commands below to check your code:
load imdata . mat
f i g u r e (1 )
imagesc ( data )
colormap gray
a x i s image
a x i s o f f
e f ined imagedata = re f i n e image ( data )
f i g u r e (2 )
imagesc ( re f ined imagedata ) ;
colormap gray
shading i n t e r p
a x i s image
a x i s o f f
4 Root Finding
The spring shock system in automobiles is an example of a real mechanical system which
involves the deflection of nonlinear springs that would demonstrate a damping effect ove
4
E7 – Prof. Alam Fall 2019, UC Berkeley
time (see behavior response in Figure 1 below). The resistance force of the spring F (N),
with respect to time t (s), is given by the following equation
F (t) = k1(cos (πt))(e
− t
k2π )
where k1 and k2 are defined as spring constants (N/s)
Figure 1: Damping effect with spring constants k1 = 1000 N/s and k2 = 0.6
We will use the above equation to describe the force on a car’s shock system following a
speed bump encountered at t = 0 seconds. We want to find out a time instance near t = 5.5
seconds when the resistance force F is equal to zero. In order to do that, we need to find the
oots of F (t).
You will find the roots using three different methods. In part (a), you will use the MATLAB
uilt-in function fzero. In part (b), you will use the Newton-Raphson method. In part
(c), you will use the bisection method.
Write a function damping with the function declaration line
func t i on [ S ] = damping ( k1 , k2 , i n i t i a l )
The inputs are:
• k1, k2: the spring constants defined above.
• initial: the time instance of interest, near which we want to find a root.
5
E7 – Prof. Alam Fall 2019, UC Berkeley
The output is:
• S: a 1-by-1 structure a
ay. Its fieldnames and content will be described in detail in the
following subproblems.
In your function, represent F (t) as a function handle stored in variable F. Also, create a
variable tol = 1e-6 which represents the value of the tolerance required for the bisection
method and the Newton-Raphson method. Note that the tolerance is defined such that if
the absolute value of the difference between the actual function value F(t) and the expected
function value (0, in this case) is ≤ tol, then the algorithm should be stopped.
(a) Use fzero to determine the root of F with initial as the guess location. Store the
oot to S in the fieldname fzero root.
(b) Within damping.m, write a subfunction myNewton with function declaration line
func t i on newton root = myNewton( f , df , x0 , t o l )
which will recursively calculate the root nearest to x0 given the function handle f,
the function handle of the derivative of f (which you need to calculate on your own)
df, the initial guess x0, and the tolerance tol.
Store the root found using Newton-Raphson method to S in the fieldname newton.
(c) Within damping.m, write a subfunction bisection with function declaration line
func t i on b i s e c t r o o t = b i s e c t i o n ( f , a , b , t o l )
which will recursively calculate the root located between a and b given the function
handle f, the endpoints of the range a and b, and the tolerance tol.
Use a range of ±0.5 (make sure to check the signs of the function evaluated at the
endpoints) and store the root found using bisection method to S in the fieldname
isect root.
The expected