Sine equation in Euler Method ODE

2 views (last 30 days)
Chris Horne
Chris Horne on 31 Mar 2022
Commented: Chris Horne on 1 Apr 2022
I am trying to code the following equation
y = SQRT(1+x)*cos(x^2)
The exact solution is so much different from the Euler result. I am njot sure if it is axes or code typo. Any help is appreciated. Thanks
% 2D SYSTEM SIN-COS EQUATIONS USING FORWARD EULER METHOD
% % Initial conditions and setup
% x(t)= sqrt(t+1)*cos(t).^2
t(1)=0; % initial t
x(1)=1; y(1)=0; ; % initial x,y
dt=0.005; % time step
nn=10000; % number of time steps in radians
for i=1:nn-1 % loop counter
dx=(sqrt(i+1)*cos(i.^2))./2*(i+1) - 2*i*(sqrt(i+1)*sin(i.^2)); %x' equation
dy=(sqrt(i+1)*sin(i.^2))./2*(i+1) - 2*i*(sqrt(i+1)*cos(i.^2)); %y' equation
x(i+1)=x(i)+dt*dx; % Find new x
y(i+1)=y(i)+dt*dy; % Find new y
t(i+1)=t(i)+dt; % Find new t
end
% exact solution without Euler
tt = linspace(0,10000,length(x));%set the t,x vectors to same lengthtt=0:0.1:10000;
x_exact = sqrt(1+tt).*cos(tt.^2); %????? typo?
% Plot x vs. t with exact solution
figure
t = linspace(0,10000,length(x)); %set the t,x vectors to same length
plot(t,x,'b')
%axis([0 10000 -5000 5000])
hold on
title('2D System sqrt(1+t)cos equation vs exact soln')
xlabel('t')
ylabel('x')
plot(t,x_exact,'r')
  2 Comments
Chris Horne
Chris Horne on 31 Mar 2022
Yes, dx/dt = (sqrt(t+1)*cos(t^2)) / 2*(t+1) ) - 2*t*(sqrt(t+1)*sin(t^2))
with exact solution x = sqrt(1+t)*cos(t^2)
Thanks for your interest

Sign in to comment.

Accepted Answer

Torsten
Torsten on 1 Apr 2022
% 2D SYSTEM SIN-COS EQUATIONS USING FORWARD EULER METHOD
% solve u'(t) = F(t,u(t)) where u(t)= sqrt(t+1)*cos^2(t),sqrt(t+1)*sin^2(t) ;
% % Initial conditions and setup
% x(t)= sqrt(t+1)*cos(t).^2
% y(t)= sqrt(t+1)*sin(t).^2
t(1)=0; % initial t
dt=0.005; % time step
nn=1000; % number of time steps in radians
t = zeros(1,nn);
x = zeros(1,nn);
y = zeros(1,nn);
x(1) = 1;
y(1) = 0; % initial x,y
for i=1:nn-1 % loop counter
dx=sqrt(t(i)+1)*cos(t(i)^2)./(2*(t(i)+1)) - 2*t(i)*(sqrt(t(i)+1)*sin(t(i)^2));%x' equation
dy=sqrt(t(i)+1)*sin(t(i)^2)./(2*(t(i)+1)) + 2*t(i)*(sqrt(t(i)+1)*cos(t(i)^2)); %y' equation
x(i+1)=x(i)+dt*dx; % Find new x
y(i+1)=y(i)+dt*dy; % Find new y
t(i+1)=t(i)+dt; % Find new t
end
% x Exact value without Euler
tt = linspace(t(1),t(end),numel(t));%set the tt,x vectors to same length
x_exact = sqrt(1+tt).*cos(tt.^2);
% Plot x vs. t with exact solution
figure
%t = linspace(0,1000,length(x));%set the t,x vectors to same length
plot(t,x,'b') % plot approximation in blue
hold on
title('2D System sqrt(1+t)cos equation vs exact soln')
xlabel('t')
ylabel('x')
plot(tt,x_exact,'r') % plot exact in red
legend('x','exact')
hold off;
figure
% y Exact value without Euler
tt = linspace(t(1),t(end),numel(t));%set the tt,y vectors to same length
y_exact = sqrt(1+tt).*sin(tt.^2);
% Plot y vs. t with exact solution
plot(t,y,'b')
hold on
title('2D System sqrt(1+t)sin equation vs exact soln')
xlabel('t')
ylabel('y')
%t = linspace(0,1000,length(y_exact));%set the t,y vectors to same length
plot(tt,y_exact,'r')
hold on
legend('y','exact')
hold off
  3 Comments

Sign in to comment.

More Answers (1)

James Tursa
James Tursa on 31 Mar 2022
Edited: James Tursa on 31 Mar 2022
You've only got one scalar differential equation, so I don't understand why you think you need two variables x and y along with independent variable t to handle this. You only need the x and t that are present in your original equation. E.g., something like this outline for the Euler part based on your posted code:
% % Initial conditions and setup
nn = 10000; % number of time steps in radians
t = zeros(1,nn);
x = zeros(1,nn);
t(1) = 0; % initial t
x(1) = 1; % initial x
dt = 0.005; % time step
for i=1:nn-1 % loop counter
dx = sqrt(t(i)+1)*cos(t(i)^2)/(2*(t(i)+1)) - 2*t(i)*(sqrt(t(i)+1)*sin(t(i)^2)); %x' equation
x(i+1) = x(i) + dt*dx; % Find new x
t(i+1) = t(i) + dt; % Find new t
end
So, only one derivative equation for the x. And all of those i's get replaced with t(i)'s. Also I have added a pair of parentheses to ensure that the entire 2*(t(i)+1) is in the denominator where it belongs.
For the exact solution, again you need only the one scalar x equation. And you need to make sure you are using the same t scale. E.g.,
tt = linspace(t(1),t(end),numel(t));
Finally, it is not clear to me what your intended final time is. In your Euler code it seems to be 9999*dt, but in your exact solution code it seems to be 10000. You need to make sure these are the same in order for your plots to match up properly.
  1 Comment
Chris Horne
Chris Horne on 31 Mar 2022
James, Thanks for the tip on the t(i). Then I did not expect the plot results to be so much different. The 'cos' code is same as the 'sin' code, more or less. The exact solution for y is much different than x.
% 2D SYSTEM SIN-COS EQUATIONS USING FORWARD EULER METHOD
% solve u'(t) = F(t,u(t)) where u(t)= sqrt(t+1)*cos^2(t),sqrt(t+1)*sin^2(t) ;
% % Initial conditions and setup
% x(t)= sqrt(t+1)*cos(t).^2
% y(t)= sqrt(t+1)*sin(t).^2
t(1)=0; % initial t
x(1)=1; y(1)=0; % initial x,y
dt=0.005; % time step
nn=1000; % number of time steps in radians
t = zeros(1,nn);
x = zeros(1,nn);
for i=1:nn-1 % loop counter
dx=sqrt(t(i)+1)*cos(t(i)^2)./(2*(t(i)+1)) - 2*t(i)*(sqrt(t(i)+1)*sin(t(i)^2));%x' equation
dy=sqrt(t(i)+1)*sin(t(i)^2)./(2*(t(i)+1)) - 2*t(i)*(sqrt(t(i)+1)*cos(t(i)^2)); %y' equation
x(i+1)=x(i)+dt*dx; % Find new x
y(i+1)=y(i)+dt*dy; % Find new y
t(i+1)=t(i)+dt; % Find new t
end
% x Exact value without Euler
tt = linspace(t(1),t(end),numel(t));%set the tt,x vectors to same length
x_exact = sqrt(1+tt).*cos(tt.^2);
% Plot x vs. t with exact solution
figure
t = linspace(0,1000,length(x));%set the t,x vectors to same length
plot(t,x,'b') % plot approximation in blue
hold on
title('2D System sqrt(1+t)cos equation vs exact soln')
xlabel('t')
ylabel('x')
plot(t,x_exact,'r') % plot exact in red
legend('x','exact')
hold off;
figure
% y Exact value without Euler
tt = linspace(t(1),t(end),numel(t));%set the tt,y vectors to same length
y_exact = sqrt(1+tt).*sin(tt.^2);
% Plot y vs. t with exact solution
plot(t,y,'b')
hold on
title('2D System sqrt(1+t)sin equation vs exact soln')
xlabel('t')
ylabel('y')
t = linspace(0,1000,length(y_exact));%set the t,y vectors to same length
plot(t,y_exact,'r')
hold on
legend('y','exact')
hold off;

Sign in to comment.

Products


Release

R2021b

Community Treasure Hunt

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

Start Hunting!