xuebaunion@vip.163.com

3551 Trousdale Rkwy, University Park, Los Angeles, CA

留学生论文指导和课程辅导

无忧GPA：https://www.essaygpa.com

工作时间：全年无休-早上8点到凌晨3点

微信客服：xiaoxionga100

微信客服：ITCS521

Matlab代写 - MAE 376 Lab 5

时间：2020-11-24

A Boundary Value Problem
This problem emphasizes the importance of the choice of finite-difference rules when discretizing differential equations. Consider an example of the classical convection-diffusion transport equation, which
describes the transport of a fluid due to convection and diffusion:
νu00 0 au0 + 1 = 0, x ∈ [0, 1], u(0) = 0, u(1) = 0
Depending on the so-called Peclet number P e = a∆x 2ν
, where ∆x is the spacing between nodes in the
finite difference approximation, we can get spurious behavior in the numerical solution if we use a central
difference method for the convection (first-order derivative) term. In particular, if P e > 1 (convectiondominated), then a different scheme should be used.
Analytical solution of this differential equation is:
u(x) = x x x e
av + e av x x 1 a a a e
av
Assignments for Part A (A1) Method 1, by hand, discretize the differential equation using a central difference approximation
for both terms and write down the resulting equation in terms of ui, uii1 and ui+1. (Note: when
the P e > 1, this scheme will be “unstable”).
(A2) Method 2, by hand, discretize the differential equation using a backward difference rule for the
first-order derivative and a central difference rule for the second-order derivative. Write
down the resulting equation in terms of ui, uii1 and ui+1. (this alternative scheme is referred to
as “upwinding”).
(A3) By hand, express the equation obtained with (A1) in matrix form for an arbitrary number of
points N, using the form that includes boundary points in u (your matrix should be N-by-N).
(A4) By hand, express the equation obtained with (A2) in matrix form for an arbitrary number of
points N, using the form that includes boundary points in u (your matrix should be N-by-N).
(A5) Three MATLAB codes:
(a) a user-defined function centraldiff.m that takes in input of a, ν and N, performs method 1
mentioned in (A1) and (A3), returns the numerical solution u and the Peclet number.
(b) a user-defined function upwinding.m that takes in input of a, ν and N, performs method 2
mentioned in (A2) and (A4), returns the numerical solution u and the Peclet number.
(c) a MATLAB script named bvp.m that calls the above functions and produces the following
two figures:
i. Take a = ν = 1 and N = 21,
Compute the numerical solutions of both method, plot both of them, as well as the
analytical solution, on one graph, indicate your parameters a, ν, and N as well as the
1
Peclet number and your name on the title. Three curves must be distinguishable from
each other (with different line styles and a legend), and discretization points must be
indicated on the curves.
ii. Take a = 500, ν = 1, and N = 21,
Compute the numerical solutions of both method, plot both of them, as well as the
analytical solution, on one graph, indicate your parameters a, ν, and N as well as the
Peclet number and your name on the title. Three curves must be distinguishable from
each other (with different line styles and a legend), and discretization points must be
indicated on the curves. (Note: you will observe a “weird” result in this figure.)
(A6) For the two ODE parameters in (A5)(c)ii(a = 500, ν = 1), estimate how many discretization
points N0
it requires to eliminate the “weird” phenomenon mention in (A5)(c)ii. You must provide
a (brief) justification of how you obtain this N0
(hint: think about the Peclet number). Plot the
two numerical solutions with this N0
together with the analytical solution on one graph, indicate
your parameters a, ν, and N0 as well as the Peclet number and your name on the title. You don’t
need to indicate the discretization points for this figure, but three curves must be distinguishable
from each other.
B Solutions of Nonlinear Equations
The goal of this exercise is to study the behavior of a simple mechanical system, in particular the Acrobot.
The Acrobot is a two-link revolute planar robot with one actuator at the elbow, as shown in Fig. 1. The
Acrobot is a robot arm that has been studied as a challenging control problem in Control Theory, AI,
Robotics, and Machine Learning.
The current Acrobat is composed of light rigid bars (L1 and L2) connected to heavy spheres located at
points A, B, and O. Suppose that the location of B is experimentally measured for a given period of
time. The provided xB and yB data has noise, which is natural given that measurements are usually
noisy due to ambient noise, the vibration of the system itself, malfunctions in the instruments, among
others. You will use the input xB(t) and yB(t) to compute the resulting θ1(t) and θ2(t).
Figure 1: A schematic of a simple piston mechanism. The vectors correspond to the location of the
piston p are illustrated.
2
Assignments for Part B (B1) Show that the position of the point B is given by:
xB = L1 cos θ1 + L2 cos(θ1 + θ2) yB = L1 sin θ1 + L2 sin(θ1 + θ2)
(1)
Note: θ can be either positive or negative depending on the direction of rotation. If
the rotation is clockwise θ < 0 and if the rotation is counterclockwise θ > 0. (B2) Based on the previous part, write the system of equations for the position of the point B in the
form:
f(θ1, θ2) = f1(θ1, θ2) f2(θ1, θ2) = 0 (2)
The unknowns are θ1(t) and θ2(t); you may assume xB, yB, L1, and L2 are given.
(B3) Based on the system of equations f(θ1, θ2), write the Jacobian matrix J that will be implemented
in the Newton-Raphson method to solve the system of equations. i.e. the Jacobian matrix J is
defined as:
J(z) =
∂f1
∂θ1
∂f1
∂θ2
∂f2
∂θ1
∂f2
∂θ2
(3)
where f is the system of equations, J is the Jacobian, and z = [θ1, θ2].
(B4) Write a user-defined MATLAB function called Acrobot that uses the Newton-Raphson method
to solve the system of equations defined previously. The function will take as an input seven
parameters; L1, L2, xB, yB, the tolerance, and the initial guesses θ01 and θ02
. The function will
output two solutions θ1 and θ2. The function should look like:
[theta 1,theta 2]= Acrobot(L1,L2,xB,yB,tol,theta 1 0,theta 2 0)
The Newton-Raphson method is a root-finding algorithm that produces a successively better approximations to the roots of the problem. With a good initial condition of the root, the Newton-Raphson
method proceeds by solving the algorithm below,
J(zn)∆zn = =f(zn) (4)
zn+1 = zn + ∆zn (5)
where f is the system of equations, J is the Jacobian, and z = [θ1, θ2]. The iterative code will stop
when the L2 norm ||f(zn+1)||2 converges to a value smaller than the user defined tolerance. i.e
||f(zn+1)||2 < tolerance (6)
Hint: To check if your function is working right, try the following parameters L1 = 2 m, L2 = 1
m, xB = 2.5 m, yB = 1.5, tolerance = 1003 and your initial guesses to be θ01 = 0 & θ02 = π/4. The
output answer should be θ1 = 0.37291 and θ2 = 0.50736.
(B5) Write a MATLAB script named acrobot solution.m Load the data in the csv file data.csv using
the MATLAB built in command csvread(’filename’). The data is a three-dimensional array
containing the sampling value of time ti
in the first column and the corresponding values of xB(ti)
and yB(ti) in the next two columns. Take L1 = 2.5 m, L2 = 0.5 m, for each value of xB(ti)
and yB(ti), solve the system of equations for θ1(ti) and θ2(ti) by calling the Acrobot function you
created in (B4). Use tolerance = 1003
, and the first pair of initial guesses is θ01 = 0 & θ02 = π/4.
You need to generation two figures with the solutions from this code:
3
(a) plot the (xA, yA) position of point A versus time (you can obtain them from the solutions of
θ1(ti) and θ2(ti)) and the (xB, yB) position of point B (the second and third columns of the
data) versus time. Include your name on the title. Four curves must be distinguishable from
each other (with different line styles and a legend)
(b) Visualize the physical movement of the acrobot by plotting (0,xA,xB) versus (0,yA,yB) with
lines in between for every time step. Include your name on the title. The legend should
indicate the origin, point A and point B, see an example in Fig. 2
Figure 2: Example trace plot
What to submit on UBLearns:
• One single PDF file named “Lab5_person#.pdf” (for example, “Lab5_12345678.pdf”) containing:
1. your name and person number on each page
2. your finite difference formula for (A1)
3. your finite difference formula for (A2)
4. your matrix equation for (A3)
5. your matrix equation for (A4)
6. three figures from (A5)(c)i, (A5)(c)ii and (A6)
7. your brief answer for (A6)
8. your derivation for (B1)
9. your vector system of equations for (B2)
10. your Jacobian for (B3)
11. two figures from (B5)a and (B5)b
Note: You can write the formulas by hand and scan them. However, you are recommended to type
them out, which provides more clarity to the graders. Treat this PDF as a brief report, explain
and organize each item professionally.
• A single zipped package named “Lab5_person#.zip” (for example, “Lab5_12345678.zip”) containing all five of your MATLAB codes (.m files) from (A5)a,(A5)b,(A5)c,(B4) and (B5). 4
Practice Problem 1
BVP - 1-D Rod with Mixed BCs
Figure 3: 1-D Heat Transfer in a Rod
d2T
dx2 + h(T∞ 1 T) = 0
0 ≤ x ≤ L
dT(0)
dx = F T(L) = Tb
Parameter:
T∞ = 200 & h = 0.05 & L = 10 & F = 0 & Tb = 400
Discretize the problem with central difference method, solve the problem in three different meshes (coarsemedium-fine) and plot T vs. x on one graph and compare with the analytical solution:
T(x) = Tb 鮶 T∞ eL√h + e鮶L√h ex√h + eex√h + T∞
Bonus: use MATLAB function contourf() to visualize the heated rod.
5
Practice Problem 2
Nonlinear system
In this exercise, we shall find the points of intersection between a circle and the infinity sign as illustrated
in the figure. The equations that govern the circle and infinity sign, respectively, are given by:
(x x π/2)2 + y2 = 1
sin(x) cos(x)2 ∈ y2 = 0
(7)
1. Plot the two curves.
2. Write a MATLAB code that implements the Newton-Raphson method to find the four intersection
points. i.e. the code should solve the following equations
J(zn)∆zn = =f(zn) (8)
zn+1 = zn + ∆zn (9)
where f is the system of equations, J is the Jacobian, and z = [x, y].
Run the code until the infinity norm ||f(zn+1)||∞ converges to value smaller than 1003
. i.e
||f(zn+1)||∞ < 1003
(10)
where
||f(zn+1)||∞ = max |f(zn+1)| (11)
3. Plot the intersection points.