numerical technics in engineering

courseENGT5140_2018_20191.pdf

Home >Engineering homework help >Mechanical Engineering homework help >numerical technics in engineering

COURSEWORK on

Numerical Techniques in Engineering (ENGT5140 )

This is a course work on Numerical Techniques in Engineering (ENGT5140), which is the part of the module assessment. This course work determines 50% of your grade for the module.

All coursework will be marked against the following criteria:

• Correctness of results (60%), including mathematical derivations and ac- curacy of simulation

• Elegance and clarity of presentation (20%). Your work should be presented as if writing for an audience of fellow students, and explain everything they would need to reproduce your work.

• Elegance and style of computer programming work (20%)

• Any course work must be conducted individually, cheating of any kind will not be tolerated. You may be asked for a presentation or an oral examination based on your course work, if necessary.

In your submitted coursework, it should be included:

• Details of your coursework report, such as any necessary materials such as plots, tables and diagrams, in particular, any discussions and analysis about your results.

• Programming codes in Matlab

• Please copy the question and then follow the answers

You are required to submit an electronic copy of the coursework report through the blackboard. Submitted report will be checked for plagiarism and provide a submission record. The deadline for the electronic report submission is on Thursday 10th January 2018, before 2:00pm . Failure to submit elec- tronic copy on time might result in your report not being marked according to the University policies.

A brief guideline on electronic submission through the Blackboard is given as follows:

• Login BB and into ENGT5140 shell

• Click Assessments

• Click Coursework Submission

• In Section ASSIGNMENT SUBMISSION please choose Browse My Computer to submit your work

Your report should be in PDF format except MATLAB codes. Any other format is not allowed.

Coursework Questions

There are 15 questions in full, however, only four questions are your coursework, which will be chosen (determined) by the computer (by a random choice).

There is no mark if you do the questions selected by yourself.

1. Write Matlab codes to find the roots of the functions

(a) Consider the function

f(x) = 14xex−2 − 12ex−2 − 7x3 + 20x2 − 26x + 12

on the interval [0, 3]

i. Plot the function on the interval

ii. Apply Newton’s method to find both roots

iii. Determine for which roots Newton’s method convergences linearly and for which the convergence is quadratic

(b) Given the function

f(x) = 3

√ 1 −

i. Find the root of f(x)

ii. Apply Newton’s method and plot the first 50 iterates. This is another way Newton’s method can fail, by producing a chaotic trajectory. Please give a good discussion for this problem.

2. Newton’s method in an application problem.

(a) The idea gas law for a gas at low temperature and pressure is

PV = nRT

where P is pressure (in atm), V is volume (in L), T is temperature (in K), n is the number of moles of the gas, and R = 0.0820578 is the molar gas constant. Van der Waal’s equation (

P + n2a

V 2

) (V −nb) = nRT

covers the non-ideal case where these assumptions do not hold. Use the ideal gas law to compute an initial guess, followed by Newton’s method applied to van der Waals’ equation to find the volume of one mole of oxygen at 320 deg.K. and a pressure of 15 atm. For oxygen, a = 1.36 L2 atm/mole2 and b = 0.003183 L/mole. State your initial guess and solution with three significant digits.

(b) Use the same code, to find the volume of 1 mole of benzene vapor at 700 deg.K. under a pressure of 20 atm. For benzene, a = 18.0L2.atm/mole2 and b = 0.1154 L/mole

3. Numerical integrations

(a) Apply the Simpson’s method to the integrals with 20 < n < 200

(a)

∫ 1 0 ex

2 dx, (b)

∫ π 2

ex − 1 sin x

dx (c)

∫ π 0 ecosxdx

You must choose at least three different ns to calculate the solutions and make your discussions for your numerical results.

(b) Approximate the integrals using n = 3 Gaussian quadrature and compare with the correct value and give the approximation error

(a)

∫ 1 −1

(x3 + 2x)dx, (b)

∫ 1 −1

x2 + 1 dx, (c)

∫ 3 −3 xe−x

2/2dx

4. Solve linear systems Write Matlab codes for both Jacobi and Gauss-Seidel methods, respectively, to solve the following equations. Make sure at the beginning you solve each equation for the variable that has the largest coefficient, why?

(a) Apply your codes to the system of equations

x1 + 9x2 − 2x3 = 36 2x1 −x2 + 8x3 = 121 6x1 + x2 + x3 = 107

with an initial x1 = 1, x2 = 1, x3 = 1

(b) Apply your codes to solve the system

A(t)x = b

starting from x = [0,0,0]T, where

A(t) =

 1 t tt 1 t t t 1

  , b =

 22

 

For t = 0.2, 0.5, 0.8, 0.9, 0.99, determine the number of steps (or iterations) to obtain the solution to 6 significant digits in the computations. Plot the number of steps as functions of t and make your comments.

(c) Solve the sparse system within six correct decimal places for n = 100, n = 1000 and n = 10000 by using your codes which you made for solving (a) and (b). The correct solution is xe = [1, 1, · · · , 1]T . Report the number of step iterations and the error max |xn −xe| for all n, here xn is the numerical solutions and xe is the exact solution. The system is

 3 −1 −1 3 −1

. . . . . .

. . .

−1 3 −1 −1 3

 

  x1 x2 ... ... xn

  =

 

2 1 ... 1 2

 

5. A cylindrical pipe has a hot fluid flowing through it. Because the pressure is very high, the walls of the pipe are thick. For such a situation, the temperature distribution in the metal wall related to radial distance can be modelled by a second order differential equation

r d2θ

dr2 + dθ

dr = 0

where θ is temperature, and r is the radial distance from the centreline.

In the following questions, you must use both the shooting method and finite difference method to solve the above equation, and compare the solutions with different mesh steps, h, by each method. Plot solutions and make your judgment on solutions for each method.

(a) Solve for the temperature θ(r) within a pipe whose radius is 1 cm and whose outer radius is 2 cm if the fluid is at 500oC and the temperature of the outer circumference is 25oC.

(b) The pipe is insulated to reduce the heat loss. The insulation used has proper-

ties such that the gradient dθ

dr at the outer circumference is proportional to the

difference in temperature from the outer wall to the surroundings, which can be modelled by

dθ

∣∣∣∣ r=2

= 0.083(θ(2) − 20)

Solve the equation with this boundary condition.

6. A cable hung at its two ends as shown in Fig. 1 by its own weight will have a catenary shape described by the equation

y = Tx w

[ cosh

Tx − 1 ]

(1)

where w is the weight per unit length and Tx is the horizontal, x-component of the tension of the cable. Equation (1) is for the case when both w and Tx are constant

Figure 1: A diagram for a cable.

throughout the cable. In fact, Equation (1) is the solution of the differential equation

d2y

dx2 =

dx =

[ 1 +

( dy

)2]1/2 (2)

where s is a variable along the length of the cable.

(a) Applying the Runge-Kutta method to solve the equation (2), where w = 0.12 kN/m, xA = yA = 0, xB = 200 m, and yB = 50 m. Let the initial conditions be y = 0 at x = xA = 0, iterate Tx value until yB is within 99.9% of 50 m.

(b) How could the problem be solved if xA = 100 m and yA = 25 m by application of the Runge-Kutta method (w remains equal to 0.12 kN/m)?

7. The equation y′ = 1 + y2, y(0) = 0

has the analytical solution y(t) = tan t. The tangent function is infinite at t = π/2.

(a) Write your program by using both Euler and Runge-Kutta method to solve this initial value problem between t = 0 and t = 1.6, and compare the results of the program with the analytical solution

(b) What is the behaviour of the Runge-Kutta fourth order method when used be- tween t = 0 and t = 1.6

8. Given the following non-linear boundary value problem

y′′ − 18y2 = 0, y(1) = 1

3 , y(2) =

(a) Use the shooting method to approximate solution

(b) Use finite difference to approximate solution

and discuss your results with both methods

9. Here is a typical steady-state heat flow problem. Consider a thin steel plate to be a 10 × 20 (cm)2 rectangle. If one side of the 10 cm edge is held at 1000C and the other three edges are held at 00C, what are the steady-state temperature at interior points? We can state the problem mathematically in this way if we assume that heat flows only in the x and y directions: Find u(x,y) (temperature) such that

∂2u

∂x2 + ∂2u

∂y2 = 0 (3)

with boundary conditions

u(x, 0) = 0

u(x, 10) = 0

u(0,y) = 0

u(20,y) = 100

We replace the differential equation by a difference equation

h2 [ui+1,j + ui−1,j + ui,j+1 + ui,j−1 − 4ui,j] = 0 (4)

which relates the temperature at the point (xi,yj) to the temperature at four neigh- bouring points, each the distance h away from (xi,yj). An approximation of Equation (3) results when we select a set of such points (these are often called as nodes) and find the solution to the set of difference equations that result.

(a) If we choose h = 5 cm , find the temperature at interior points.

(b) Write a program to calculate the temperature distribution on interior points with h = 2.5, h = 0.25, h = 0.025 and h = 0.0025 cm. Discuss your solutions and examine the effect of grid size h.

∂2u

∂x2 + ∂2u

∂y2 = xy(x− 2)(y − 2)

on the region 0 ≤ x ≤ 2, 0 ≤ y ≤ 2

with boundary condition u = 0 on all boundaries except for y = 0, where u = 1.0. Write and run the program with different grid sizes h and discuss your numerical results.

10. If a anti level beam of length L, which bends due to a uniform load of w lb/ft, is also subject to an axial force P at its free end (see Figure 2). The equation of its elastic curve is

EI d2y

dx2 = Py −

wx2

2 (5)

For this equation, the origin O has been taken at the free end. I is the moment of

inertia; here I = bh3

12 . At the point x = L,

dx = 0 and at x = 0, y = 0.

O X

Figure 2: Diagram for a cantilever beam.

(a) Solve the boundary value problem by any numerical method for a wooden beam, 2 in × 4 in × 10 ft, where E = 12 × 105 lb/in2. Find y versus x when the beam has the 4 in dimension vertical with w = 25 lb/ft and a tension force of P = 500 lb.

(b) And also solve for the deflections if the beam is turned so that the 4 in dimension is horizontal.

[1 + (y′)2]3/2 d2y

dx2 = Py −

wx2

2 (6)

If the defection of the beam is small, the difference is negligible, but in some cases this is not true. Furthermore, if there is considerable bending of the beam, the horizontal distance from the origin to the wall is less than L, the original length of the beam. Solve Equation (6) with the same conditions, and determine by how much the deflections differ from those previous calculations.

11. Consider a laterally insulated metal bar of length L = 1, which satisfies the heat equation

∂u

∂t = ∂2u

∂x2

Suppose that the ends of the bar are kept at temperature u = 0oC and the temperature in the bar at t = 0, u(0,x) = sin πx.

(a) Using an explicit method with different meshes h to solve the equation and find the numerical solutions, then plot the solutions for t = 0, t = 1, t = 2, t = 3, t = 4 and t = 5 in the same figure.

(b) Using an implicit method with different meshes h to solve the equation and find the numerical solutions, then plot the solutions for t = 0, t = 1, t = 2, t = 3, t = 4 and t = 5 in the same figure.

(c) Using Crank-Nicolson method with different meshes h to solve the equation and find the numerical solutions, then plot the solutions for t = 0, t = 1, t = 2, t = 3, t = 4 and t = 5 in the same figure.

(d) Making your discussion for these three methods with different size of h

The exact solution for this problem is

u(t,x) = sin(πx)e−π 2t

which can be used to compare with your numerical solutions.

12. Condon and Odishaw (1967) discuss Duffing equation for the flux φ in a transformer. This nonlinear differential equation is

d2φ

dt2 + ω20φ + bφ

3 = ω

N E cos(ωt), φ(0) = φ′(0) = 0

In this equation , E cos ωt is the sinusoidal source voltage and N is the number of turns in the primary winding, where ω0 and b are parameters of the transformer design. Make a plot of φ versus t (and compare to the source voltage) if E = 165, ω = 120π, N = 600, ω20 = 83 and b = 0.14. For approximate calculations, the nonlinear term bφ3 is sometimes neglected. Evaluate your result to determine whether this makes a significant error in the results. (Condon E. and Odishaw H., (1967) Handbook of Physics, New York, McGraw-Hill)

13. An 80kg paratrooper is dropped from an airplane at a height of 600m. After 5 second the chute opens. The paratrooper’s height as a function of time, y(t), is given by

d2y

dt2 = −g +

α(t)

m y(0) = 600 (m)

y′(0) = 0 (m/s)

where g = 9.81m/s2 is the acceleration due to gravity and m = 80kg is the para- trooper’s mass. The air resistance α(t) is proportional to the square of the velocity, with different proportionality constants before and after the chute opens,

α(t) =

{ K1y

′(t)2, t < 5s K2y

′(t)2, t ≥ 5s

(a) Find the analytical solution for free fall, K1 = 0 and K2 = 0. At what height does the chute open? How long does it take to reach the ground? What is the impact velocity? Plot the height versus time and label the plot appropriately.

(b) Consider the case K1 = 1 15

. and K2 = 4 15

. At what height does the chute open? How long does it take to reach the ground? What is the impact velocity? Make a plot of the height versus time and label the plot appropriately.

14. Apply the shooting method to the boundary value problem

y′1 = 4 −y2 t3

y′2 = −e y1

y1(1) = 0, y2(2) = 0, 1 ≤ t ≤ 2

Note that if the initial condition y2(1) were given, this would be an initial value problem. You may apply the shooting method to determine the unknown y2(1), using matlab routine ODE45 to solve the initial problem. Please give the details about using the shooting method to solve this problem.

The exact solutions are y1(t) = ln t and y2(t) = 2 − t2/2. You can use it to compare with your numerical solutions.

15. The shooting method for ODE:

(a) Consider the second order differential equation

y′′ + 2y′ − 3y = 0, y(0) = e2, y(2) = 1

First, find the exact solution for this equation, and then write a matlab code to use the shooting method to find the approximate solutions. In this question, you should give the detail about choosing the initial guesses. Please read lecture note p68-69. You should compare your numerical results with the exact solutions.

(b) Then use your code to solve the following nonlinear second order differential equation

y′′ + 0.1(y′)2 + 0.6y = sin t, y(0) = 1, y(2) = 2

Clearly, you could not find the exact solution, however, you should plot your numerical solutions and make a good discussion with different guesses.