how to solve a second order differential equation using Euler's method?

How to solve a second order differential equation (boundary value problem) using Euler's Method without using inbuilt matlab functions such as ode45?

[James Tursa]

What instructions were you given in class? What have you done so far? What specific problems are you having with your code?

[Remston Martis]

I have a second order differential equations with 2 boundary conditions.

I have to use Euler's method(the shooting method) to solve the equation.

I am able to code for a first order differential equation but not for a second order differential equation.

    function Eout = Eulers(F, yint,h,yfinal,x0)
% F is the desired differential equation
% yint and yfinal are the boundary value conditions
% h is the step size
% x0 is the initial guess

      x = x0;
      Eout = x;
      for y = yint : h : tfinal-h
         s = F(y,x);
         x = x + h*s;
         Eout = [Eout; x];
      end

The code works for a simple first order differential equation.

My differential equation is d^2y/dx^2 = 2y +8x(9-x)

I have managed to solve it mathematically on paper but am struggling with the code. I can provide you with the mathematical solution if you require it.

Answer 1

[Torsten]

If you define

F=@(y,x)[y(2);2*y(1)+8*x*(9-x)]

you can work with the code from above.

Best wishes

Torsten.

[Remston Martis]

When I do that I receive the error "Index exceeds array bounds"

[Torsten]

You've set x0 as a (2x1) vector, haven't you ?

[Jan]

@Remston Matris: Please post the complete error message, such that the readers do not have to guess, in which line the error occurs.

[Remston Martis]

Sorry for being excessively terse.

error message "Error in @(y,x)[y(2);2*y(1)+8*x*(9-x)]

Error in Eulers(line 14) s=F(y,x)

[Remston Martis]

I have not set x0 as a vector. It is an integer.

x0 is my initial "shooting" guess. I do not understand why it should be a vector.

[Torsten]

You have a second order differential equation.

This has to be written as a system of two first-order ODEs:

y1'=y2
y2'=y1+8*x*(9-x)

So you will need two initial values.

That's the reason why x0 (and x) are (2x1).

But I think, according to your notation y is the independent variable and x is the dependent variable to be solved for. Thus, F has to read

F=@(y,x)[x(2);2*x(1)+8*y*(9-y)]

[Remston Martis]

No, y is my dependent variable. I'll try your code and revert if it works. I think I understand where you're getting at. Thank you .

[Torsten]

y can't be your dependent variable because you update x as

x = x + h*s;

where s is obtained by

s = F(y,x);

[Jan]

@Remston Martis: Using y, x and t causes confusion here. You have e.g. "yfinal" and "tfinal". Use either "x and y" or "x and t". Then:

function dx = F(t, x)
dx = [x(2); x(1) + 8 * t * (9 - t)];
end

Now "yint and yfinal are the boundary value conditions" confuses me, but I assume you want something like:

x    = x0;
tvec = tint:h:tfinal-h;
Eout = zeros(2, size(tvec));
Eout(:, 1) = x;

for idx = 1:length(tvec)
   t  = tvec(idx);
   dx = F(t,x);
   x  = x + h * dx;
   Eout(:, idx) = x;
end

Answer 2

[Jan]

It is worth to be nitpicking:

 % x0 is the initial guess

No, x0 is the initial value of the trajectory when you consider the integration. To solve a boundary value problem, you need an additional layer around the integration: e.g. a single shooting or multiple shooting method. Then this x0 is the initial guess of the shooting method.

To solve your problem, convert the 2nd order equation to a system of two equations of order 1. Then y has 2 components: The initial position and velocity. Converting higher order equations to order 1 is the first step for almost all integrators.

[Remston Martis]

Hi, Thank you for the information. I have actually solved the entire question on paper but mathematically so I guess I'm on the same page with what you have stated.

However, I do not know how to apply the 2nd order equation conversion into my code. I am unable to translate it into the code.

Would greatly appreciate if you could just help with those couple lines of code that I could add into the above.

Thank you!