User:EGM6341.s10.Team1.Kumanchik/HW7

=Problem 5: Determine the limit on the step size in the Hermite-Simpson method, Forward Euler, Backward Euler=

Statement
Hermite-Simpson algorithm was used to find $$\hat{h}=0.125$$ to give an error of order $$10^{-6}$$ from the true value where $$\hat{h}$$ is the step size between nodes in the numerical solution. The function to be integrated was the population logistic equation,


 * $$\dot{x}=rx\left ( 1-\frac{x}{x_{max}} \right )$$

The numerical integration using the Hermite-Simpson algorithm was compared to the true value for x,


 * $$x\left ( t \right )=\frac{x_0x_{max}e^{rt}}{x_{max}+x_0\left ( e^{rt}-1 \right )}$$

For this problem, the step size is increased at a rate of $$h_{new}=2^k\hat{h}$$ where $$k=1,2,3,...$$

Solution
The differential equation is first integrated using ODE45 (at MATLAB command).

Next the Hermite-Simpson algorithm is run using increasing values of step size. The max step size is limited by the size of the interval $$t\in [0,10]$$. The step size maxes out at 8. A step size of 16 is not allowed. Therefore this algorithm is stable at integrating this function.

Next the Forward Euler integration is used. The step size maxes at 2 at which point it wildly oscillates about $$x_{max}=10$$

Finally, backward Euler integration is used. The step size is allowed to go to 8 so it is stable like the Hermite-Simpson algorithm. Notice, though, that the answer is less accurate for the same number of steps compared to the Hermite-Simpson method

=Problem 6: Interceptor Simulation using Hermite-Simpson algorithm=

Statement
Use the Hermite-Simpson algorithm to determine the trajectory of an interceptor given the following inputs.


 * $$T(t)=\left\{\begin{matrix}

6000 N, & 0\leq t< 27, & 33\leq t\leq 40\\ 1000 N, & 27\leq t< 33 & \end{matrix}\right.$$


 * $$\alpha(t)=\left\{\begin{matrix}

0.03 rad & 0\leq t< 21\\ \frac{0.13-0.09}{27-21}(t-21)+0.09 rad & 21\leq t< 27\\ \frac{-0.2+0.13}{33-27}(t-27)-0.13 rad & 27\leq t< 33\\ \frac{-0.13+0.2}{40-33}(t-33)-0.2 rad & 33\leq t\leq 40 \end{matrix}\right.$$


 * $$\mathbf{f}=

\begin{bmatrix} \dot{\gamma}\\ \dot{v}\\ \dot{x}\\ \dot{y} \end{bmatrix}=\begin{bmatrix} \frac{T-D}{mv}sin(\alpha)+\frac{L}{mv}cos(\alpha)-\frac{g}{v}cos(\gamma)\\ \frac{T-D}{m}cos(\alpha)+\frac{L}{m}sin(\alpha)-gsin(\gamma)\\ vcos(\gamma)\\ vsin(\gamma) \end{bmatrix}$$


 * $$\mathbf{z}=\begin{bmatrix}

\gamma\\ v\\ x\\ y \end{bmatrix}$$


 * $$DL=\left\{\begin{matrix}

D(\alpha,v,y)=0.5C_d\rho v^2S_{ref}\\ L(\alpha,v,y)=0.5C_c\rho v^2S_{ref} \end{matrix}\right. $$


 * $$\left\{\begin{matrix}

C_d(\alpha)=A_1\alpha^2+A_2\alpha+A_3\\ C_l(\alpha)=B_1\alpha+B_2\\ \rho(y)=C_1y^2+C_2y+C_3 \end{matrix}\right.$$


 * $$\left\{\begin{matrix}

m=1005 kg\\ g=9.81 \frac{m}{s^2}\\ S_{ref}=0.3376 m^2\\ A_1=-1.9431\\ A_2=-0.1499\\ A_3=0.2359\\ B_1=21.9\\ B_2=0\\ C_1=3.312e^{-9} \frac{kg}{m^5}\\ C_2=-1.142e^{-4} \frac{kg}{m^4}\\ C_3=1.224 \frac{kg}{m^3} \end{matrix}\right.$$


 * $$t\in [0,40]$$ seconds
 * $$h=1$$ seconds

Solution
The control inputs are the thrust, T, and the angle of attack, $$\alpha$$. They are plotted below.

The Hermite-Simpson algorithm takes the solution to a differential equation and subdivides it into nodes. Each node is related to the next by,


 * $$z_{i+1}=z_{i}+\frac{h}{6}\left [ f_i+4f_{i+\frac{1}{2}}+f_{i+1} \right ]$$
 * $$z_{i+\frac{1}{2}}=\frac{1}{2}\left [ z_i+z_{i+1} \right ]+\frac{h}{8}\left [ f_i-f_{i+1} \right ]$$

However, this equation is dependent on the $$z_{i+1} \,$$ node which appears in the $$f_{i+1} \,$$ and $$f_{i+\frac{1}{2}}$$ terms. Therefore, the equation is solved for by the Newton Raphson method,


 * $$F=-z_{i+1}+z_{i}+\frac{h}{6}\left [ f_i+4f_{i+\frac{1}{2}}+f_{i+1} \right ]=0$$

Intelligently choosing values for $$z_{i+1} \,$$ allows $$F \,$$ to approach 0. The intelligent algorithm is,


 * $$z^{(k+1)}=z^{(k)}-\frac{dF(z^{(k)})}{dz}^{-1}F(z^{(k)})$$

where $$k $$ is the previous guess and $$k+1$$ is the next guess. This algorithm can search for $$z_{i+1}$$ by starting with $$z^{k}=z_i$$. Using $$F$$ as a metric, the algorithm can stop when $$F\approx0$$ or for this problem we choose $$F$$ of $$O(10^{-9})$$.

This solution can be extended to this interceptor problem which solves a system of differential equations. Placing the differential equations in matrix form, the solution remains identical. Keep in mind that the differential of a matrix is the Jacobian. The result of the interceptor problem is plotted below. The time step is 1 second between 0 and 40 seconds.



The source code of the interceptor algorithm

Thrust and angle functions (control variables)

=Problem 7: Integrate the logistic equation using an inconsistent trapezoidal-Simpson algorithm=

Statement
Integrate the logistic equation,
 * $$\dot{x}=rx\left ( 1-\frac{x}{x_{max}} \right )$$

Using the Hermite-Simpson algorithm and the inconsistent Trapezoidal-Simpson algorithm, compare their errors for the same step size, $$h=0.125$$. Both methods require the Newton-Raphson method so use a error for each node of $$O(10^{-9})$$. Compare the integrations to the true integral of the logistic equation,
 * $$x\left ( t \right )=\frac{x_0x_{max}e^{rt}}{x_{max}+x_0\left ( e^{rt}-1 \right )}$$

Solution
The Hermite-Simpson algorithm was adjusted to create the inconsistent trapezoidal-Simpson algorithm. Both solutions were subtracted from the true solution. The error in the Hermite-Simpson algorithm was designed from problem 5 as less than $$O(10^{-6})$$. Using the same parameters on the inconsistent trapezoidal-Simpson algorithm yields worse error by about 3 orders of magnitude. The errors from the true value are plotted below. The line on the bottom is the Hermite-Simpson error. It is much smaller than the inconsistent trapezoidal-Simpson error.