# Solving a system of second order ODE backwards

25 views (last 30 days)
Commented: John D'Errico on 4 Mar 2023
Dear all,
I am trying to solve a simple second order ODE but I was hoping to solve it backwards. With that I mean, that I have info at t=1 and I want to get the value of the solution for t in (0,1]. I believe that myODE will solve it forward...meaning you need an initial condition for t=0. Could you help me out here?
function dv=myODE(h,y)
% v = v_tt + 1/t v_t+u
%where u is known and v_t=0 at t=1.
t = (0+h :h: 1-h);
u=1+cos(2*pi*t);
dy(1) = y(2);
dy(2)=y(1)-h*y(2)-u;
end
Also, when calling this function, I am not sure why it is giving me an error
[t,v]=ode45(@myode((t,v),h:T,[1;0]) (Here I was trying to put some intiial conditions because Im not sure how to do it backwards, but it still gives me errors forward in time.)

John D'Errico on 2 Mar 2023
Edited: John D'Errico on 2 Mar 2023
Easy, peasy. For example, solve the ODE
y' = sin(t)
subject to the initial condition, that at t==1 we have y(1)=1/2. Now solve it moving from 1 to 0.
odefun = @(t,y) sin(t);
y0 = 1/2;
tspan = [1 0]; % it will solve with a negative time step
[tout,yout] = ode45(odefun,tspan,y0);
plot(tout,yout,'r-',1,y0,'bo') So it started at t==1, and then used a negative time step, moving down to t==0.
Again, the initial condition is indeed y(1)==1/2. You can see the solution passes through that point.
##### 2 CommentsShow 1 older commentHide 1 older comment
John D'Errico on 4 Mar 2023
ODE45 is an adaptive solver, so your time steps are not truly the steps used by ODE45. It may need to go smaller steps, or it may interpolate as needed. But if you specify specific time steps, ODE45 will give them to you.
odefun = @(t,y) sin(t);
y0 = 1/2;
tspan = 1:-0.1:0;
[tout,yout] = ode45(odefun,tspan,y0);
[tout,yout]
ans = 11×2
1.0000 0.5000 0.9000 0.4187 0.8000 0.3436 0.7000 0.2755 0.6000 0.2150 0.5000 0.1627 0.4000 0.1192 0.3000 0.0850 0.2000 0.0602 0.1000 0.0453

Torsten on 2 Mar 2023
An example:
fun = @(t,y) y;
[t1,y1] = ode45(fun,0:0.01:1,1);
[t2,y2] = ode45(fun,1:-0.01:0,exp(1));
hold on
plot(t1,y1)
plot(t2,y2)
hold off
grid on I see! thanks!

William Rose on 2 Mar 2023
First, try to get your script to work in the normal forward direction.
The funciton needs to be adjusted. You have
function dv=myODE(h,y)
% v = v_tt + 1/t v_t+u
%where u is known and v_t=0 at t=1.
t = (0+h :h: 1-h);
u=1+cos(2*pi*t);
dy(1) = y(2);
dy(2)=y(1)-h*y(2)-u;
end
That will not work for several reasons.
1. The output variable name, dv, does not match dy(1) and dy(2) used in the function body.
2. dy must be a column vector, not a row vector.
3. The internal definition of t as a vector is not needed, and causes trouble.
4. You need to pass t and y and, possibly, h, where h is a constant of the model. Alternatively, you could define h inside the function, which is what I do below. I assume u is an external forcing function.
A better version of the function would be
function dy=myODE(t,y)
% v = v_tt + 1/t v_t+u
% where u is known and v_t=0 at t=1.
dy=[0;0]; % assure dy is a column vector
u=1+cos(2*pi*t); % external forcing
dy(1) = y(2);
dy(2) = y(1)-h*y(2)-u;
end
Your call to ode45() includes "myode" but its name is "myODE". The case matters.
tspan=[0,1]; % time span, going forward
y0=[1,0]; % initial conditions
[t,v]=ode45(@myODE,tspan,y0);
Try that. Once it is working in the forward direction, then we can think about running it backward.
One option for going backwards is to define t2=1-t. Then reverse the signs of the derivatives in myODE (since that is what will happen when you do the substitution of variables), then integrate t2 forward, from t2=0 to 1. Specify the initial condition (at time t2=0) as whatever the state is at t=1. When you get the answer [t2,y] from ode45(), compute t=1-t2, and plot t versus y.
##### 2 CommentsShow 1 older commentHide 1 older comment
William Rose on 2 Mar 2023

Thank you all! This is very helpful!