Solving system of ODEs using Euler's method

Hello everyone,

I need to model a trajectory of a flying object and this process is described by a system of two 2nd-order ODEs. I have already reduced it to a system of four 1st-order ODEs:

with z1(0)=0, z2(0)=Vcosα, z3(0)=0, z4(0)=Vsin(α) while k is 0.1, m is the mass of the object, g is 9.8, V is the initial velocity and α is the launch angle.

I have only practcied solving a single 1st-order of ODEs using Euler's method before so I don't really get the hang of this problem. Could anyone please help me to solve this problem in the simplest way possible?

Many thanks!

[James Tursa]

What have you done so far? What does your single variable Euler code look like? You just need to add the code for the other three variables. Pretty much a cut & paste three times of your single variable Euler stepping code (the stuff inside your for-loop) and then modify as appropriate.

[Nobita Nobi]

Hi James,

Here is how my single variable Euler code looks like:

function [X,Y] = euler
a=0;
b=1;
n=10;
y0 = 1;
f=@(x,y)((-2*x*y)/(x+1));
h=(b-a)/n;
y = y0;
Y = []; 
X = a:h:b;

for x = a:h:b
   y = y + h*f(x,y);
   Y = [Y y];
end
plot(X,Y,'--pm'); title('Eulers  method')

[Eshan]

poora likhke bhejna pls 5 baje ke phle

Answer 1

[Samuele Sandrini]

Hi, i saw your code and the ODE's system, you have already done most of the problem. In fact, now you can bring yourself back to what you already know using matrix notation:

where is a vector containing your auxiliary variables and the function is equal to:

so in your code the function f will be replaced by:

f=@(t,z) [z(2);-z(2); z(4); -9.81;]; %z is the vector of auxiliary variables containing z1=z(1), z2=z(2)...

You will need to make other minor changes to your code, like the initial condition that will be a column vector like (i use z0 in according to your system notation but in your code is y0):

%z0=[z1_0; z2_0; z3_0; z4_0];           

and also the assignment of variable y ("z") during the cycle will be change in something like z(:)=... (and not y=y+.., but the second part h*f can have the same logic) beacuse in each column you will have ,... . Keep attention in the dimentions of the vectors z, z0 and use the same notation otherwise z(:) = z(:) + h * f can become a matrix and in this case is wrong. Also the construnction of vector Z will change and i think that is quite:

Z=[Z;z];    %instead of Y=[Y y];

it is not the most efficient way but it will works (more efficient is preallocate all matrix Z=zeros()..).

Let me know if you need more information and detail.

[shireesha myadari]

Hi, I just strated to lean differntial equations in matlab. I wnat to apply euler method to solve this equation as part of my thesis work. can some one help me. here we know the values of Δβ and γ and length of the fibre(z=50,100,150,200). I am confusing how can i write As and Ap while soving equation 1. i have asuumed intial conditions like ap(0)=1

I am insreting code what I have with me but I wrote equation with only ap values.

clear; clc; close('all');

c=3*10^8;%m

lamda_0=1559*10^-9;%units in m

lamda_p=1560*10^-9;%m

lamda_s=1540*10^-9;%m

pi=3.14;

dlamda = 0.03*10^-3;%dD/dlamda= s/m^2.m

deltabeeta =((((-2*pi*c)./((lamda_0)^2))).*(dlamda)*(lamda_p-lamda_0)*((lamda_p-lamda_s)^2));

gamma = 11*10^-3; %gamma=(2*pi*n2)./(lamda*A)=11 w^-1.m^-1

h=50;

z=0:h:300;

ap=zeros(size(z));

ap(1)=1;

n = numel(ap)

for i=1:n-1

f = (1j*gamma)*(((ap(i)^2)*ap(i)+(2*ap(i)*exp(1i*deltabeeta*z(i)))));

ap(i+1) = ap(i) + h * f;

end

plot(z,ap)

grid 'on'