# Euler's code for multiple ODEs

44 vues (au cours des 30 derniers jours)
Nicholas Moung le 8 Fév 2021
Commenté : James Tursa le 8 Fév 2021
Hello every one!
May I ask you a favour! I've written a small program for Euler's method to solve multiple ODEs. But I'm encountering an error. Could you please give me a pice of advice. In the following I would like to calculate (x and y). I'm sure I made some mistakes.
Function file for ODEs is
function dydx= FncTG(x,y,Wx_r,Wy_r,Wz_r)
dydx(1)=(Wy_r*sin(y))+(Wz_r*cos(y));
dydx(2)=(Wx_r-(tan (x)*Wy_r*cos(y))-(Wz_r*sin(y)));
dydx= dydx(:);
end
------------------------------------------------------------------------
Euler's method code is
function E = fncEg (f,tspan,y0,x,y,Wx_r,Wy_r,Wz_r)
h = 1/32;
t = 1:h:48;
N = length(t); % A vector to store the time values .
y = zeros (1,N); % Initialize the Y vector .
f=@(x,y)(FncTG);
for i = 1:(N-1)
y (i+1)= y(i)+ h*f(x(i),y(i)); % Update approximation y at t+h
end
E=[t' y'];
end
--------------------------------------------------------------------------
The calculation code is
%Initial values of angles
Theda(1)=0.0025/57.2958;
Gamma(1)=0.0013/57.2958;
h=1/32;
t = 1:h:48;
tspan=[1 48];
y0=[Theda(1),Gamma(1)];
%Importing data from STAT
STAT=importdata('STAT.txt');
nx= STAT(:,10);
ny= STAT(:,11);
nz= STAT(:,12);
Wx_d= STAT(:,7);
Wy_d= STAT(:,8);
Wz_d= STAT(:,9);
%Calculating by Euler method
[x,y]=fncEg(@(x,y)f,tspan,y0,x,y,Wx_r,Wy_r,Wz_r)
##### 0 commentairesAfficher -2 commentaires plus anciensMasquer -2 commentaires plus anciens

Connectez-vous pour commenter.

### Réponses (1)

James Tursa le 8 Fév 2021
Modifié(e) : James Tursa le 8 Fév 2021
You have two states, so your result should have two values per step. So this code:
y = zeros (1,N); % Initialize the Y vector .
f=@(x,y)(FncTG);
for i = 1:(N-1)
y (i+1)= y(i)+ h*f(x(i),y(i));
y = zeros (2,N); % Initialize the Y vector . 2 rows per step, not 1
y(:,1) = y0; % initialize first state
%f=@(x,y)(FncTG); % delete this line
for i = 1:(N-1)
y(:,i+1)= y(:,i) + h*f(x(i),y(:,i)); % work with columns of y, not elements of y
And you need to use the proper function handle. So this
[x,y]=fncEg(@(x,y)f,tspan,y0,x,y,Wx_r,Wy_r,Wz_r)
should be this
[x,y]=fncEg(@(x,y)FncTG(x,y,Wx_r,Wy_r,Wz_r),tspan,y0,x,y,Wx_r,Wy_r,Wz_r)
You might consider having your function use a tspan and h that is passed in instead of hard coding this inside the function.
##### 12 commentairesAfficher 10 commentaires plus anciensMasquer 10 commentaires plus anciens
Nicholas Moung le 8 Fév 2021
I've used x instead of y(1) and y instead of y(2). If I write the function file for derivatives like above, then will I need to change everywhere, such as
y = zeros (2,N); % Initialize the Y vector . 2 rows per step, not 1
y(:,1) = y0; % initialize first state
%f=@(x,y)(FncTG); % delete this line
for i = 1:(N-1)
y(:,i+1)= y(:,i) + h*f(x(i),y(:,i)); % work with columns of y, not elements of y
and
[x,y]=fncEg(@(x,y)FncTG(x,y,Wx_r,Wy_r,Wz_r),tspan,y0,x,y,Wx_r,Wy_r,Wz_r)
James Tursa le 8 Fév 2021
So, you are going to have to make a choice. Either you work with two different state variables x and y, or you work with one two-element state variable y where y(1) and y(2) are the states. If you work with x and y, then yes you will need to rewrite your propagation equations including both the x and the y propagation code explicitly. If you use the two-element y state vector, then you need to rewrite your derivative equation. You can't have it both ways simultaneously. So pick one method and make your code consistent.

Connectez-vous pour commenter.

### Catégories

En savoir plus sur Numerical Integration and Differential Equations dans Help Center et File Exchange

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by