Recall the parachute problem where the velocity satisfied v g k m v where g was 32 and k m was 1 5 We implemented a method Euler s method by putting v v t h v t h and we showed different that for different values of h In simulink there was an integrator In Matlab itself there is an integrator which you can program that picks an appropriate step size They are built on Runge Kutta methods More on this later Create an M file function velocityp fall time velocity velocityp 32 1 5 velocity save in fall This just evaluates the function at time and velocity Now in command window type time vel ode23 fall 0 10 0 plot time vel Means that we are starting at time 0 and ending at time 50 and we are starting at velocity 0 What if we also wanted the distance also The nice thing is that ode23 can take the second component as a vector and return a vector If velocity height we can use 2 components one for the velocity and one for the height We thus take a vector y velocity height a column vector and write function vp2 fall2 t y vp2 32 1 5 y 1 y 1 In command window time y2 ode23 fall2 0 10 0 1000 plot time y2 2 To plot time versus distance What is the problem Set 1000 to 180 and see what happens Again we need to stop We can set stopping criteria using options that are set using odeset This such a common occurrence that they have figured out something called an event Do an event subroutine function gstop isterminal direction g t y gstop y 2 isterminal 1 direction which says when the second variable becomes 0 from any direction stop Save it and then in the command window issue opts odeset events g time y ode23 fall2 0 10 0 180 opts In the football problem we need 3 parameters y y x The function is now function vp3 football t y vp3 32 100 sin pi 6 y 1 100 cos pi 6 save this in football We can use the same g function and issue time y3 ode23 football 0 10 50 0 0 0 opts And to plot it plot y3 3 y3 2 The simplest are of the form of R k 1 Take a half step get the new slope and then take a full step with that slope starting at the beginning midpoint Or 2 Take a step get a new slope take the average of the slopes at the beginning and end and use that average from the beginning trapezoidal rule Often what is done is use a rule and then do it again with step size halved If the 2 are close use the same step size else halve the step size Now let us do another problem predator prey model or chemical reaction model sometimes on exams for graduate students in chemistry Rabbits have an infinite supply of food and foxes prey on rabbits for food The more rabbits the more foxes which produces fewer rabbits which produces fewer foxes which produces more rabbits Ode r 2r a r f f f a r f If a is 0 they do not interact and rabbits multiply and foxes die off of starvation With a 0 an encounter results in a rabbit being eaten so rabbits decrease and foxes increase Set a 0 and initially have 300 rabbits and 150 foxes and go for 10 units of time Let us have the function take a vector called animals besides the time t The user function is function rf predator t animal a 0 01 rf 2 animal 1 a animal 1 animal 2 animal 2 a animal 1 animal 2 In the command window use time rabfox ode23 predator 0 10 300 150 plot rabfox 1 rabfox 2 Now as a function of time plot time rabfox 2 plot time rabfox 1 r We can do all these in simulink