代写AMME2000/BMET2960/BMET9960

AMME2000/BMET2960/BMET9960 - Assignment 2, 2023
Due: 11:59 pm Friday 19th May (Week 12) 2023
Assignment Informa on
Assignment 2 focuses on your understanding of the analy cal and numerical solu on to the Wave Equa on and
Laplace Equa on, and also the analy cal solu on to the Heat Equa on without boundary condi ons.
• Your report must be your own work. Please see the FAQ on Academic Integrity for informa on on
what is appropriate and what is not appropriate for Assignments.
Any students who use third-party services such as Chegg to post or view assignment ques ons will
receive 0 for the assignment and be reported. These services provide us with detailed data on user
viewing/pos ng because they are hos ng our copyrighted material. For more informa on, see the
first lecture where academic integrity, including penal es, was discussed at length.
• Present your assignment as a concise report in PDF format;
• The report must be typed, including all mathema cal working; handwriten assignments will
receive a mark of 0 and marks will be deducted for screenshots of equa ons and/or figures;
• Your MATLAB code must be included as an appendix to the assignment report. 10% of the
assignment marks are allocated for the readability and quality of your MATLAB code;
• The report should not exceed 10 pages (this limit does NOT include the code in the appendix);
• Structure the report using the appropriate Sec on and Ques on numbers (e.g. 1.1, 2.2). Ensure all
relevant working and final answers are in that sec on;
• All figures and tables in your report must be numbered and must be referred to and discussed;
• Submit your report to Turni n by the due date. Late submissions incur a penalty of 5%/day.
We have prepared an exemplary assignment to guide you in wri ng your own.
QUESTION 1: WAVE EQUATION [40 marks]
Relevant parts of the course for this ques on:
• Week 6 Lecture 1
• Week 6 Lecture 2
• Week 7 Pre-work video: ‘The starting problem’
• Week 7 tutorial and Week 8 tutorial
A torsional wave can propagate along a stiff rod if the rod experiences a torque. The resulting wave can be
modelled using the wave equation:
𝜕𝜕2𝜃𝜃
𝜕𝜕𝑡𝑡2 = 𝑐𝑐2 𝜕𝜕2𝜃𝜃
𝜕𝜕𝑥𝑥2 (1.1)
where 𝜃𝜃 = 𝜃𝜃(𝑥𝑥,𝑡𝑡) is the torsion angle of the wave in radians at location 𝑥𝑥 (along the length of the rod) and time 𝑡𝑡,
and 𝑐𝑐 is given by:
𝑐𝑐 =  𝐺𝐺
𝜌𝜌 (1.2)
where 𝐺𝐺 is the shear modulus of the rod in units of Pa, and 𝜌𝜌 is the density of the rod in units of kg.m-3. If the
boundary conditions are 𝜃𝜃(0,𝑡𝑡) = 𝜃𝜃(𝐿𝐿,𝑡𝑡) = 0, the general solution to Eqn (1.1) is given by:
𝜃𝜃(𝑥𝑥,𝑡𝑡) =  [𝐴𝐴𝑛𝑛 cos(𝜆𝜆𝑛𝑛𝑡𝑡) + 𝐵𝐵𝑛𝑛 sin(𝜆𝜆𝑛𝑛𝑡𝑡)]
∞
𝑛𝑛=0
sin  
𝑛𝑛𝑛𝑛𝑛𝑛
𝐿𝐿   (1.3)
where 𝜆𝜆𝑛𝑛 = 𝑛𝑛𝑛𝑛𝑛𝑛
𝐿𝐿 are the eigenvalues, 𝐿𝐿 is the length of the rod, and 𝐴𝐴𝑛𝑛 and 𝐵𝐵𝑛𝑛 are coefficients which depend on the
initial conditions. Consider a thin metallic filament of length L being used as a sensor in an engineering application
(Figure 1.1):
The filament is prone to torsional waves, which in turn impact the sensor function. To fully understand the impact,
you want to characterise the torsion 𝜃𝜃 = 𝜃𝜃(𝑥𝑥,𝑡𝑡) in the sensor filament. You may assume the boundary conditions
are 𝜃𝜃(0,𝑡𝑡) = 𝜃𝜃(𝐿𝐿,𝑡𝑡) = 0 so that the solution is given by Eqn (1.3).
1.1. If the initial torsion angle along the length of the rod is 0 everywhere and the initial torsional velocity along
the rod is given by 𝑔𝑔(𝑥𝑥), shown in Figure 1.2, write down the displacement initial condition and the velocity
initial condition mathematically. [4 marks]
1.2. Starting from Eqn (1.3), show why for this particular problem the solution reduces to:
𝜃𝜃(𝑥𝑥,𝑡𝑡) =   𝐵𝐵𝑛𝑛 sin(𝜆𝜆𝑛𝑛𝑡𝑡)
∞
𝑛𝑛=1
sin  
𝑛𝑛𝑛𝑛𝑛𝑛
𝐿𝐿   (1.4)
[2 marks]
1.3. Given the initial velocity profile (Figure 1.2), solve for the 𝐵𝐵𝑛𝑛 coefficients in Eqn (1.4). [4 marks]
Note: You may use the fact that the Fourier sine series for a triangular function of the form:
is given by:
𝑓𝑓(𝑥𝑥) =  2𝐴𝐴𝐴𝐴2𝑠𝑠𝑠𝑠𝑠𝑠  
𝑛𝑛𝑛𝑛
𝐾𝐾  
(𝐾𝐾 − 1)𝑛𝑛2𝜋𝜋2
∞
𝑛𝑛=1
𝑠𝑠𝑠𝑠𝑠𝑠  
𝑛𝑛𝑛𝑛𝑛𝑛
𝐿𝐿   (1.5)
Therefore, you do NOT need to perform any integrations yourself for this question.
1.4. Implement your analytical solution for the torsional wave in MATLAB. For this implementation use the
following parameter values:
• 𝐿𝐿 = 0.1 m
• 𝐺𝐺 = 30 GPa
• 𝜌𝜌 = 9 × 103 kg. m−3
• 𝑔𝑔0 = 50 rad. s−1
• 𝑛𝑛𝑥𝑥 = 101 (number of points in spatial domain)
In your report, provide a plot of your solution at 4 different time points, all overlaid on the same plot. Make
sure you indicate the time corresponding to each curve and provide a brief discussion of your results. Include
your MATLAB code in the appendix of your report. [10 marks]
1.5. Write down the form of the CTCS scheme to solve this problem numerically. (Note: You do not need to
derive the CTCS scheme.) Briefly justify why the CTCS scheme is a suitable choice to solve this problem.
(Hint: Consider the properties of the scheme and the expected physics of the problem. You do not need to
include any derivations for this part.) [2 marks]
1.6. Briefly describe what is meant by the ‘starting problem’ in the CTCS scheme and explain how you will
handle it for this problem. You should outline the specific steps and how you plan to implement those steps
in MATLAB. [4 marks]
1.7. Now implement the CTCS numerical scheme in MATLAB to solve the problem numerically. For the
numerical implementation you should use the following parameters:
• ∆𝑡𝑡 = the maximum stable time step
• 𝑛𝑛𝑥𝑥 = 101 (number of points in spatial domain)
In your report you should include the following:
• Explain how you computed the maximum stable time step. (Note that you are NOT required to
perform a von Neumann stability analysis.) [2 marks]
• Provide a plot showing the analytical and numerical solutions overlaid at the same 4 time points used
in 1.4. Briefly discuss your results. [6 marks]
Note: Make sure your implementation compares the analytical and numerical solutions at the same (not different)
time at each step.
1.8. Plot the maximum absolute percentage error between the analytical and numerical solution as a function of
time over one complete cycle of the wave up and down the rod. On the same plot overlay the mean absolute
percentage error as a function of time. Use 51 evenly spaced time points for your horizontal axis. (Note: In
your computation of the mean and maximum absolute percentage error, make sure you avoid any time points
or spatial locations that will introduce a divide-by-0 error.) Briefly explain how you generated the plot and
discuss the results. [6 marks]
QUESTION 2: LAPLACE EQUATION [40 marks]
Relevant parts of the course for this ques on:
• Week 7 Lecture 1
• Week 7 Lecture 2
• Week 9 tutorial
An engineering control experiment involves exposing a small rectangular wafer to several lasers. The laser light
results in a ‘photon pressure’ distribution on the wafer which we can model using the Laplace equation:
𝜕𝜕2𝑆𝑆
𝜕𝜕𝑥𝑥2 +
𝜕𝜕2𝑆𝑆
𝜕𝜕𝑦𝑦2 = 0 (2.1)
where 𝑆𝑆 = 𝑆𝑆(𝑥𝑥, 𝑦𝑦) is the 2-dimensional pressure distribution on the wafer in units of Pa. Figure 2.1 shows the
wafer dimensions and the pressure distribution at the boundary corresponding to 𝑥𝑥 = 𝑎𝑎. The other boundaries have
no pressure applied.
2.1 For the scenario shown in Figure 2.1, write down the 4 boundary conditions for the pressure S
mathematically. [4 marks]
2.2 Use separation of variables to show that the steady-state pressure distribution for the scenario in Figure 2.1 is
given by:
𝑆𝑆(𝑥𝑥, 𝑦𝑦) =  𝐴𝐴𝑛𝑛
∗ sinh  
𝑛𝑛𝑛𝑛𝑛𝑛
𝑏𝑏  
∞
𝑛𝑛=1
sin  
𝑛𝑛𝑛𝑛𝑛𝑛
𝑏𝑏   (2.2)
where 𝐴𝐴𝑛𝑛
∗ are coefficients yet to be determined. [4 marks]
2.3 Apply the non-zero boundary condition to solve for the 𝐴𝐴𝑛𝑛
∗ coefficients and hence write down the complete
solution for 𝑆𝑆(𝑥𝑥, 𝑦𝑦). [4 marks]
Note: For this question you may use Eqn (1.5) to help you solve for 𝐴𝐴𝑛𝑛
∗ – i.e. you should NOT perform any
integrations for this part.
2.4 Suppose the final arrangement of the lasers results in the boundary conditions shown in Figure 2.2:
Briefly explain how you would obtain the steady-state solution for 𝑆𝑆(𝑥𝑥, 𝑦𝑦) for this new scenario, and write
down this solution. Note: You do NOT need to perform any lengthy derivations, you should simply explain
what mathematical principles allow you to adapt your solution from 2.3 to this new solution. [4 marks]
2.5 Implement your analytical solution from 2.4 in MATLAB. For your implementation use the following
parameters:
• 𝑎𝑎 = 𝑏𝑏 = 0.1 m
• 𝑛𝑛𝑥𝑥 = 𝑛𝑛𝑦𝑦 = 101 (number of points in each dimension of the spatial domain)
• No. of Fourier terms = 100
• 𝑆𝑆1 = 𝑆𝑆2 = 2 × 10−3 Pa
In your report, provide a clearly labelled contour plot of the steady-state solution. Note that you can use the
pcolor function in MATLAB to generate this plot. Be sure to provide a brief discussion of your plot. Include
your MATLAB code in the appendix of your report. [8 marks]
2.6 Now implement the Gauss-Jacobi numerical scheme from Week 7 Lecture 2 to solve the same problem
numerically. In your report:
• Briefly outline the steps you used to implement the Gauss-Jacobi method, including your stopping
criteria. [4 marks]
• Present 2 plots side-by-side – the left plot should be your numerical solution to the steady-state
pressure distribution; the right plot should be the signed difference image between the numerical and
analytical solutions. Discuss your results. [8 marks]
2.7 The scientists want to know how much they can shorten the wafer before the pressure at the dead-centre
becomes too great. Using your analytical solution, determine at what wafer width, 𝑎𝑎𝑡𝑡ℎ, the pressure at the
centre first exceeds 75% of 𝑆𝑆1. You should determine 𝑎𝑎𝑡𝑡ℎ using a tolerance of 1 µm. Repeat your estimate
using your numerical solution. In your report:
• Briefly outline the method you used to compute 𝑎𝑎𝑡𝑡ℎ. [2 marks]
• Present your results using the analytical and numerical solutions and briefly discuss any differences.
[2 marks]
QUESTION 3: Fourier Integral [10 marks]
Relevant parts of the course for this ques on:
• Week 8 Lecture 2
• Week 9 Lecture 1
• Week 10 tutorial
A breach in a factory pipeline has caused a large amount of the insecticide permethrin to be released into a
long narrow waterway (Figure 3.1). The concentration C(x, t) of permethrin at location x and time t along
the waterway can be modelled using Fick’s second law, as the heat equation:
𝜕𝜕𝜕𝜕
𝜕𝜕𝑡𝑡 = 𝐷𝐷
𝜕𝜕2𝐶𝐶
𝜕𝜕𝑥𝑥2 (3.1)
where D is the mass diffusivity. Because the waterway is very long, there are no boundary conditions.
3.1 State the boundary conditions mathematically and briefly explain why you cannot solve this
problem analytically using Fourier series. [2 marks]
3.2 C(x, t) can be obtained using the Fourier Integral form of the heat equation solution, as shown in Eqn
(3.2):
𝐶𝐶(𝑥𝑥,𝑡𝑡) =
1
𝜋𝜋      𝑓𝑓(𝑣𝑣) cos(𝑝𝑝𝑝𝑝 − 𝑝𝑝𝑝𝑝)𝑒𝑒−𝐷𝐷𝑝𝑝2𝑡𝑡
∞
−∞
𝑑𝑑𝑑𝑑  𝑑𝑑𝑑𝑑 (3.2)
∞
0
where p is a continuous variable and f(x) = f(v) is the initial condition. Using the known integral on
the Data Sheet:
 𝑒𝑒−𝑠𝑠2
cos(2𝑏𝑏𝑏𝑏) .
∞
0
𝑑𝑑𝑑𝑑 = √𝜋𝜋
2 𝑒𝑒−𝑏𝑏2
(3.3)
show how Eqn (3.2) can be simplified to:
𝐶𝐶(𝑥𝑥,𝑡𝑡) = 1
2√𝐷𝐷𝐷𝐷𝐷𝐷   𝑓𝑓(𝑣𝑣)𝑒𝑒
−
(𝑥𝑥−𝑣𝑣)2
4𝐷𝐷𝐷𝐷  
𝑑𝑑𝑑𝑑
∞
−∞
(3.4)
[2 marks]
3.3 Assume the initial condition is given by:
𝐶𝐶(𝑥𝑥, 0) = 𝑓𝑓(𝑥𝑥) =   𝐶𝐶0 if − 10 < 𝑥𝑥 < 10
0 otherwise
as shown in Figure 3.2. By making a suitable variable substitution, solve the remaining integral in
Eqn (3.4) to show that the final solution for C(x,t) can be expressed as the following sum of error
functions:
𝐶𝐶(𝑥𝑥,𝑡𝑡) = 𝐶𝐶0
2  erf  
(10 + 𝑥𝑥)
2√𝐷𝐷𝐷𝐷   + erf  
(10 − 𝑥𝑥)
2√𝐷𝐷𝐷𝐷    (3.5)
[4 marks]
3.4 There is an endangered frog species 300 m from the initial site of deposition of the permethrin. The
frogs will not survive if they encounter a concentration of permethrin greater than 1% of C0. Assuming
that D=0.38 m2
/s and using the table of error function values in the Data Sheet to help you, determine
if the frogs will be able to survive for 24 hours before the habitat relocation team can rescue them.
Show all working and explain all reasoning. [2 marks]