Modeling Basics

                                 Population Models

8/20             8/27            9/3            9/10               9/17             9/24                10/1                 10/8

8/20: Welcome and Introdunction to first order initial value problems. (ODE, IVP)

8/27: Solution techniques ; Use Maple to solve and plot.

Example:  y' = x^2 - 2 y/x, use Maple to:
  (1)  Find the general solution.
   (2)   Plot the directional fiels and analyze the solution.
   (3)   Plot the solutions s.t. y(1)=2,  y(0.5)=2;
   (4)   Plot  all the solution together.

9/3: Solution techniques , continued;  Matlab basics; Why numerical solution?

       See:  Computer Practice 9/3
 

9/10: Solution techniques; Approximate solution using Matlab.

                y_i+1 = y_i + h f ( x_i, y_i) ,           i = 0, 1, 2, ...

Ex 1:     y' = x + y,       y(0) = -0.5,   h = 0.1.
 
Euler's method:    y0 = -0.5
                             y1 = y0 + h f(x0,y0) = -0.5 + 0.1 ( 0 +(-0.5)) = -0.55,
                             y2 = y1+ h f(x1,y1) = -0.55 + 0.1 (0.1 + (-0.55)) = -0.595,
                             y3 = y2 + h f(x2,y2) = -0.595 + 0.1 (0.2+ (-0.595)) = -0.6345,
                                            .................................................................................................

Tedious and easy to make errors. But it is important to understand the method.
 

Use Matlab to solve it.

        function [x,y] = euler_simple(x0,y0,xfinal,n)

        h = (xfinal - x0)/n;
        x(1) = x0;
        y(1) = y0;

        for i=1:n,
           y(i+1) = y(i) + h*f(x(i),y(i));
           x(i+1) = x(i) + h;
        end

Step 2:  Create a file called f.m  whose contents are:

           function yp = f(x,y)
                yp = x + y;

Step 3: Solve the problem using Matlab by typing

        x0 =0; y0=-0.5; xfinal=2; n=40;
       [x,y] = euler_simple(x0,y0,xfinal,n);
       plot(x,y)                                                       % We also see the solution plot!

A better way is to put every thing above in the file called main.m.
 

 9/17:  Euler's method and Matlab implementation, Modified Euler's method.

         function [x,y] = euler_simple_f (x0,y0,xfinal,n, f )

        h = (xfinal - x0)/n;
        x(1) = x0;
        y(1) = y0;

        for i=1:n,
           y(i+1) = y(i) + h*feval(f, x(i),y(i));
           x(i+1) = x(i) + h;
        end
 
In calling routine (script file), can be put into a text file, for example, main.m

        x0 =0; y0=-0.5; xfinal=2; n=40;
       [x,y] = euler_simple_f(x0,y0,xfinal,n, 'f' );
       plot(x,y, 'o' )                                                       % We also see the solution plot!

Script file is a collection of Matlab commands., something like .main.m  drive.m etc.
 
Modified Euler's method:

       y_i+1 = y_i + h f ( x_I + h/2, y_i + h f ( x_i, y_i  )  ),           i = 0, 1, 2, ...

Matlab Code:

        function [x,y] = euler_simple_b(x0,y0,xfinal,n,f)
 
        h = (xfinal - x0)/n;
        x(1) = x0;
        y(1) = y0;
 
        for i=1:n,
           y(i+1) = y(i) + h*feval(f,x(i)+h/2,y(i) + h*feval(f,x(i),y(i) ) );
           x(i+1) = x(i) + h;
        end

Test problem,  main4.m

            x0=0; y0=1; xfinal=2*pi; n1=40;
            [x1,y1] = euler_simple_b(x0,y0,xfinal,n1,'f2');
 
            plot(x1,y1,'o')
            hold
            [x1,y1] = euler_simple_b(x0,y0,xfinal,160,'f2');
            plot(x1,y1)

9/24:  Lab #1, HA237-238, and HA269.