Consider the following ODE$-$IVP:

$\frac{dy}{dt}=y(-2-5y-y^2)+0\mathrm{.}5cos(t),y(0)=3,$tin$[0,5]$

that represents a CSTR with a fluctuating feed stream and a second$-$order reaction. This is very closely related to the problem introduced in

Section $4\mathrm{.}4$ in Dorfman and Daoutidis.

Part $($a$)$

Write a MATLAB program that solves this problem numerically using RK$4.$ This function should take as input the initial condition $(t_0,y_0),$

stepsize $h,$ and the final time $t_f,$ and automatically generate a labeled plot of $y(t).$ Test this function for $h=0\mathrm{.}01$ and report your results.

$\%$ Implement your main RK$4$ function at the bottom, call your function and

$\%$ plot the results using code here

$\%$ Common constants:

t$0=$ i$\%$ initial time to start integrating from

y$0=$ i$\%$ initial condition

tf $=$ ; final time to stop integrating

h$=$ ;$\%$ given stepsize

$\%$ call the RK$4$ function to solve the problem numerically

$\%$ plot the results and comment on them

Part $($b$)$

Write a function that applies linear stability analysis to determine the maximum stepsize required for stability at the current time $(t,y(t)).$ Using

this function, determine the maximum stepsize $h^{**}$ so that the system is stable at the first step.

$\%$ Implement your main stability function at the bottom, call your function

$\%$ here to calculate h$*$

Part $($c$)$

Using your code RK$4$ code from Part $($a$),$ write a program that generates two plots: $(1)$ a plot of the maximum stepsize versus $t,$ and $(2)$ a plot

of the numerical solution of $y$ versus $t$ for $h=0\mathrm{.}02,h=h^{**},$ and $h=h^{**}+0\mathrm{.}035($on a single plot$).$ For full credit, don't forget to comment on

your findings.

$\%$ call the RK$4$ function for the three different h values and plot the

$\%$ results on the same plot.

$\%$ call the stability function at each integration step to calculate the

$\%$ maximum stepsize at each timestep over the horizon and plot the results

Comment on results:

Common Functions

function $[$t$,$y$]=$ RK$4($t$0,$tf$,$y$0,$h$)$

$\%$ implement the RK$4$ algorithm and output a vector of t and y for plotting

end

function hs $=$ stability$()$

$\%$ calculates the maximum stepsize h$*$ based on the linear stability

$\%$ criterion

end

function f $=$ fun$($t$,$y$)$

$\%$ evaluates the right$-$hand side function for the RK$4$ algorithm

end

Problem $3$ Consider the following ODE$-$IVP:

$\frac{dy}{dt}=y(-2-5y-y^2)+0\mathrm{.}5cos(t),y(0)=3,$tin$[0,5]$

that represents a CSTR with a fluctuating feed stream and a second$-$order reaction.

a$)$ Write a MATLAB program that solves this problem numerically using RK$4.$ This function

should take as input the initial condition $(t_0,y_0),$ stepsize $h,$ and the final time $t_f,$ and

automatically generate a labeled plot of $y(t).$ Test this function for $h=0\mathrm{.}01$ and report your

results.

b$)$ Write a function that applies linear stability analysis to determine the maximum stepsize

required for stability at the current time $(t,y(t)).$ Using this function, determine the maximum

stepsize $h^{**}$ so that the system is stable at the first step.

c$)$ Using your code from Part $($a$),$ write a program that generates two plots: $(1)$ a plot of the

maximum stepsize versus $t,$ and $(2)$ a plot of the numerical solution $y$ versus $t$ for $h=0\mathrm{.}02,$

$h=h^{**},$ and $h=h^{**}+0\mathrm{.}035($on a single plot$).$ For full credit, don't forget to comment on

your findings.