Engineering at Alberta Courses

Ordinary Differential Equations: Problems

The current code given in the explicit Euler method assumes that $h$ is chosen such that $\frac{t_{max}-t_0}{h}$ is a natural number. Modify the code so that if this is not the case, the solution for the dependent variable at $t_{max}$ is available.
Consider the following IVP:
$\frac{\mathrm{d}x}{\mathrm{d}t}=xt+4x$

with the initial condition .
- Use the built-in DSolve function to find the exact solution and plot $x$ for $t\in[0,1]$ .
- Use the explicit Euler method with $h=0.1$ to find the values of $x_i$ . Plot the exact solution overlapping the obtained numerical solution.
- Use the implicit Euler method with $h=0.1$ to find the values of $x_i$ . Plot the exact solution overlapping the obtained numerical solution.
- Use Heun’s method with $h=0.1$ to find the values of $x_i$ . Plot the exact solution overlapping the obtained numerical solution.
- Use the midpoint method with $h=0.1$ to find the values of $x_i$ . Plot the exact solution overlapping the obtained numerical solution.
- Use the classical Runge-Kutta method with $h=0.1$ to find the values of $x_i$ . Plot the exact solution overlapping the obtained numerical solution.
Consider a mass-spring system with units, units with the initial conditions .
- Use the built-in DSolve function to find the exact solution and plot the position $x$ and velocity $x'$ for $t\in[0,2\pi]$ .
- Use the explicit Euler method with $h=0.2$ to find the values of $x_i$ and $x'_i$ . Plot the exact solutions overlapping the obtained numerical solutions.
- Use the classical Runge-Kutta method with $h=0.2$ to find the values of $x_i$ and $x'_i$ . Plot the exact solution overlapping the obtained numerical solution.
Consider the following IVP:
$y''(t)+9y=0$

with the initial conditions $y(0)=1$ and $y'(0)=0$ . Use the implicit Euler method to find a numerical solution to the given IVP for $t\in[0,5]$ with $h=0.1$ . Compare with the exact solution.
Consider the following IVP:
$y''(t)+2y'(t)+6y=\cos{(2t)}$

with the initial conditions $y(0)=1$ and $y'(0)=0$ . Use the classical 4th-order Runge-Kutta method to find a numerical solution to the given IVP for $t\in[0,10]$ with $h=0.1$ . Compare with the exact solution.
Consider the following nonlinear system of IVPs:
$\begin{split}\frac{\mathrm{d}x}{\mathrm{d}t}&=10(y-x)\\\frac{\mathrm{d}y}{\mathrm{d}t}&=x(28-z)-y\\\frac{\mathrm{d}z}{\mathrm{d}t}&=xy-\frac{8}{3}z\\\end{split}$

with the initial conditions $x(0)=-10$ , $y(0)=-10$ , and $z(0)=30$ . Use the classical 4th-order Runge-Kutta method to find a numerical solution to the given IVP for $t\in[0,20]$ with $h=0.01$ . Then, using list plots, plot the three curves showing $y$ versus $x$ , $z$ versus $x$ , and $y$ versus $z$ .
Note: The above system is an example of a Lorenz system of equations which can have chaotic solutions depending on the parameters used and the initial conditions.
Use the finite difference method to find a numerical solution to the following BVP:
$4y''(x)-2y-x=0$

with the boundary conditions $y(0)=0$ and $y(6)=0$ . Use $n=6$ intervals, and compare with the exact solution.
The deflection $y$ as a function of the position $x$ along the length of a simply supported Euler-Bernoulli beam uniformly loaded is given by:
$EI\frac{\mathrm{d}^2y}{\mathrm{d}x^2}=\frac{wLx}{2}-\frac{wx^2}{2}$

Where $E=200,000$ MPa is the Young’s modulus, $I=30,000 cm^4$ is the moment of inertia, $w=15 kN/m$ is the uniformly distributed load, and $L=3$ m is the length of the beam. The boundary conditions are given by $y(0)=y(L)=0$ . Use the finite difference method to find a numerical solution to the deflection $y$ using $n=5$ intervals. Compare with the exact solution.
A disgruntled golfer hits a small stone with all his strength. The stone’s acceleration, affected by the drag of the air and by the gravity force, can be modelled using the following pair of differential equations:
(1) $\frac{du}{dt}=-C_d\left(u\sqrt{u^2+v^2}\right)$

(2) $\frac{dv}{dt}=-C_d\left(v\sqrt{u^2+v^2}\right)-g$

where $C_d = 0.056$ is the drag coefficient, $g=9.81 \text{ m/s}^2$ in the gravitational acceleration.
When the stone (assumed to be spherical) leaves the club face, it has an initial speed of $U_0$ and an initial angle of $\theta_0$ with respect to the horizontal plane.
1. Calculate precisely the horizontal distance traveled before the rock first strikes the ground (L), for an initial speed of $U_0=265.7 \text{ km/h}$ and for an initial angle of $\theta_0=36^\circ$ . Use the Euler Explicit Method and select the time step to obtain a solution with a Relative Error Tolerance of $10^{-5}$ . Describe in detail in your Main Report the procedure used to find L.
2. Use the same function (or similar) to find accurately and efficiently, the angle corresponding to the maximum horizontal distance that the rock will travel ( $\theta_{max}$ ). Report it together with the corresponding distance ( $L_{max}$ ). Describe the procedure used to find $\theta_{max}$ .

3. Explain why the angle for the maximum distance is not $\theta_{max}=45^\circ$ , You will find the plotting of your results useful to obtain the explanation.
Your company manufactures coil springs of high quality steel and you are asked to program a process controller based on estimated heat treatment times, because actual temperatures cannot be measured. First, the straight steel wires need to be heated until they become ductile, a process called “annealing”. Then they are wound around a cylinder into their coil shape and, finally, they are tempered to give them strength by immersing them in oil (“quenching”) for rapid cooling. You need to estimate the times of annealing, of quenching and total process time.
In summary, the process is as follows:
1. Annealing: Steel wires are heated starting from $T_{amb}$ by radiant heaters at $T_{rad}$ that fully surround the wires, while gaining/loosing some heat by free convection to the air at $T_{oven}$ , until $T = T_{anneal}$ is reached.
2. Forming: When $T = T_{anneal}$ is reached, the wires are removed from the oven and automatically wound around a calibrated cylinder. The forming takes 50 s, during which the wires are cooling by radiation to the surroundings and by free convection to the air at $T_{amb}$ .
3. Tempering: Finally, the springs are quenched in oil at $T_{oil}$ , an intense free convection process (no radiation), until $T = T_{cool}$ is reached.

The equation governing the wire temperature $T$ (assumed constant through the cross-section) is
$mc_p\frac{dT}{dt}=\dot{Q}_c+\dot{Q}_r$

where $m$ is the mass of the wire, $A$ is the surface area of the wire (lateral surface only, neglect ends), the convective heat gain/loss is:
$\dot{Q}_c=\begin{cases}-h_{c_{air}}A(T-T_{oven})\quad\text{during annealing}\\-h_{c_{air}}A(T-T_{amb})\quad\text{during forming}\\-h_{c_{coil}}A(T-T_{oil})\quad\text{during quenching}\end{cases}$

and the radiative heat gain/loss is:
$\dot{Q}_r=\begin{cases}-\epsilon\sigma A(T^4-T_{rad}^4)\quad\text{during annealing}\\-\epsilon\sigma A(T^4-T_{amb}^4)\quad\text{during forming}\end{cases}$

The values for the convection heat transfer coefficient $h_c$ are:
$h_{c_{air}}=5\rm W/(m^2K)$

$h_{c_{oil}}=250\rm W/(m^2K)$

The values for temperatures are:
$\begin{aligned}T_{amb}&=21^\circ \rm{C} = 294\rm{K}\\ T_{rad}&=850^\circ \rm{C} = 1123\rm{K}\\ T_{amb}&=250^\circ \rm{C} = 523\rm{K}\\ T_{amb}&=650^\circ \rm{C} = 923\rm{K}\\ T_{amb}&=40^\circ \rm{C} = 313\rm{K}\\ T_{amb}&=45^\circ \rm{C} = 318\rm{K}\\ \end{aligned}$

The dimensions and properties for the spring wires are: length $l= 0.05\rm m$ , diameter $d= 0.003 \rm m$ , density $\rho=8131\rm kg/m^3$ , $c_p=434 \rm J/(kg K)$ , emissivity $\epsilon = 0.32$ , and the Boltzmann constant is $\sigma = 5.67\times 10^{-8}\rm W/m^2K^4$ .

1. Calculate precisely the time required to complete the annealing ( $t_{anneal}$ ) and to finish the entire process ( $t_{total}$ ) for the conditions described above. Use the Euler Explicit Method and select the time step to obtain a solution with a Relative Error of $t_{total}$ below the Tolerance of $10^{-5}$ . Plot the temperature history for the final curve.
Important: Use $T_{anneal}$ and the precisely calculated time of the end of the annealing process ( $t_{anneal}$ ) as Initial Condition and initial time for the forming process. Use the temperature and time at the end of the forming process as Initial Condition and initial time for the tempering process.

2. To minimize the energy consumption of the process, you need to reduce the temperature of the radiant heaters as low as possible, while ensuring the process still completes in a reasonable time. Use the same function above (or similar) to find accurately and efficiently, the minimum temperature for the radiant heaters $T_{{rad}_{min}}$ that ensures a total time for completion of less than $250 \rm s$ ( $t_{total}\le 250$ ). Describe the procedure used to find $T_{{rad}_{min}}$ .

3. Explain (give physical reason) why does the process fail if $T_{rad}= 953 \rm K$ ?
(Hint: compare intermediate values)

Open Educational Resources

Ordinary Differential Equations: Problems

Leave a Reply Cancel reply