INSTRUCTION FOR RUNNING THE MATLAB CODE euler1.m


First download the codes and save them in files 
(with the same names) in the same directory, 
say in a directory called 
      NUMERICS

Then open MATLAB in the directory NUMERICS 
and type the commands (where >> is the MATLAB prompt) 
 
>> diary
>> format long

Recall that, after you type a command in MATLAB, you press "return" 
to execute it.  If you type a semicolon (;) after a command, 
the computer will execute it, but will not display the output.  


    EULER'S METHOD

Suppose that you want to solve the initial value problem (IVP) 

	dy
	-- = y - t^2 + 1 ,        0 <= t <= 2
	dt 

	y(0) = 0.5 .

To solve this IVP using Euler's method, 
you will use the function 
    euler1(f,a,b,ya,N)
which is saved in the file  euler1.m ; 
note that in MATLAB the name of the function 
and the name of the file must be the same 
(except, of course the standard MATLAB extension .m ).  
Here  f  is the name of the function  f(t,y) 
giving the right-hand side of the differential equation 
(in this case,  f(t,y) = y - t^2 + 1 ), 
a  and  b  are the initial and the final "times", 
ya  is the initial value  y(a)  of the function y, 
and  N  is the number of intervals in which 
you divide the interval  [a,b].  

The program  euler1.m  does the following (look at the file).  
First it computes the stepsize  h 
- note that at the end of the line there is a semicolon (;) 
which tells MATLAB not to print on the screen 
the result of executing this line.  
Then the code defines the vector-rows  T  and  Y , 
each of length  N+1 , and initializes them with zeros.  
If you type 

>> x = zeros(1,3)

and press "return", MATLAB will type 

x =

     0     0     0


The command that follows in the code  euler1.m , 

T=a:h:b;

initializes the mesh points as follows:  
    T(1)   = a
    T(2)   = a + h
    T(3)   = a + 2h
    ...
    T(N)   = a + (N-1)h
    T(N+1) = a + Nh = b

The next line in  euler1.m  defines the initial value 
of the unknown function:

Y(1)=ya;

The loop below computes the values of the approximations 
to the solution at the mesh points:  

for j=1:N
   Y(j+1)=Y(j)+h*feval(f,T(j),Y(j));
end

The command 
    feval(f,T(j),Y(j));
returns the value of the function  f(t,y)  
for  t = T(j)  and  y = Y(j) .

Finally, the line

E=[T' Y'];

defines the value E that the program returns.  

To run  euler1.m , you have to save the right-hand side 
of the differential equation, namely, 
     y - t^2 + 1 
in a file.  This function is saved in the file  rhs.m  
- again, the name of the file is the same as the value 
of the function it defines.  

To solve the IVP above by using Euler's method 
with  N = 10 , you have to type

>> euler1(@rhs,0,2,0.5,10)

and press "return" - you will see the values 
of the approximations to the true solutions.  

The symbol  @  in front of the function's name 
is called the "handle" to the function rhs 
- this way you make MATLAB call the function 
rhs  within the routine  euler1.m  every time 
it sees  f  in this code.  

Alternatively, you can type 

>> A = euler1(@rhs,0,2,0.5,10);

- this will assign the new variable A 
the output of the program (this is what the 
variable  E  on the last line of the program is for).  
Then type 

>> A

and see what happens.