# [Numerical] Euler Backward method

#### Nick89

Hi,

I need to solve the following simple ODE with both the Euler Forward and Euler Backward numerical methods. I also need to answer for which values of T this can still be calculated:
$$y'(t) = \frac{-1}{2y(t)}$$
$$y(0) = 2$$
$$t \in (0, T]$$

Obviously the analytical solution is
$$y(t) = \sqrt{4 - t}$$
So it would seem T must be between 0 and 4 for the root to be real.

The Euler Forward method, if you are not familiar with it, follows from doing a taylor expansion and reads:
$$y_{n+1} = y_n + h y'_n$$
Here, h is the 'step size'.
This method is called the explicit method because y_n+1 can be calculated explicitly by simply evaluating the right hand side, which is known.
If an initial value is known, any next value can be calculated from this. Because y(0) is given (=2) I have been able to put this into a Matlab script and graph the numerical solution, which matched the analytical solution perfectly as far as I could tell.

Next, I also need to graph the solution using the Euler Backward method, which looks similar but is not quite the same:
$$y_n = y_{n+1} - h y'_{n+1}$$
This method is also called the implicit method because y_n+1 cannot be calculated explicitly by evaluating the right hand side. Instead it needs to be solved for y_n+1.

Unfortunately, I have no clue how to build 'solve for y_n+1' into Matlab, and I'm pretty sure we are not required to do this.
So I tried just entering the ODE to get:
$$y_n = y_{n+1} + \frac{h}{2 y_{n+1}}$$

Solving this for y_n+1 by hand, I obviously get a quadratic equation for y_n+1, for which there are two solutions:
$$y_{n+1} = \frac{ -y_n \pm \sqrt{y_n^2 - 2h}}{2}$$

Now finally my question... Of course I can use the equation above to evaluate y_n+1, but how do I know whether to use the + or the - ?
I really have no idea..? I could use the + all the time, or the - all the time, but I have no reason to believe that I should not use +, -, +, - etc or something like that. It does not sound logical to require it to use +,-,+,- etc but still, I simply don't know...

How can I determine whether to use the + or -?

Thanks!

Related Differential Equations News on Phys.org

#### matematikawan

Just wondering. What advantages does the Euler Backward method has over the Euler Forward?

#### D H

Staff Emeritus
Solving this for y_n+1 by hand, I obviously get a quadratic equation for y_n+1, for which there are two solutions:
$$y_{n+1} = \frac{ -y_n \pm \sqrt{y_n^2 - 2h}}{2}$$
Note well: You have a sign error. The two solutions are

$$y_{n+1} = \frac{y_n \pm \sqrt{y_n^2 - 2h}}{2}$$

Which to choose? First, it helps to rewrite the above as

$$y_{n+1} = y_n\left(\frac{1 \pm \sqrt{1 - 2h/y_n^2}}{2}\right)$$

Examine what happens as the step size becomes very small. One choice makes a lot of sense. The other choice is complete nonsense.

Just wondering. What advantages does the Euler Backward method has over the Euler Forward?
Note well: Euler techniques almost always yield very poor results.

That said, the advantage of using implicit integration techniques is stability (but typically at the cost of increased complexity and sometimes decreased accuracy). A simple example is exponential decay with a time constant $\tau$:

$$\frac{dx}{dt} = -{x(t)}/{\tau}$$

Here $\tau>0$. Forward Euler integration with step size $\Delta t$ leads to

$$x_{n+1} = x_n(1 - {\Delta t}/{\tau})$$

If you choose a time step that is larger than the time constant the above will diverge. Even if the time step is less than the time constant, the results can be quite bad. You choose a time step that is much less than the time constant.

In comparison, backwards Euler leads to

$$x_{n+1} = \frac{x_n}{1 + {\Delta t}/{\tau}}$$

which is always stable, even with a bad choice of a step size.

#### Nick89

Note well: You have a sign error. The two solutions are

$$y_{n+1} = \frac{y_n \pm \sqrt{y_n^2 - 2h}}{2}$$
Oops, I think you're right.

Which to choose? First, it helps to rewrite the above as

$$y_{n+1} = y_n\left(\frac{1 \pm \sqrt{1 - 2h/y_n^2}}{2}\right)$$

Examine what happens as the step size becomes very small. One choice makes a lot of sense. The other choice is complete nonsense.
When I use the + sign it converges to y_n, when I use the - sign it becomes zero. That makes sense, thanks.

I ended up using the Newton iteration method to solve the equation numerically (just because it is more general and can be applied to other functions as well), but it yields the same result (just a little slower).

Thanks for the help!

#### matematikawan

Note well: Euler techniques almost always yield very poor results.

That said, the advantage of using implicit integration techniques is stability (but typically at the cost of increased complexity and sometimes decreased accuracy). A simple example is exponential decay with a time constant $\tau$:

$$\frac{dx}{dt} = -{x(t)}/{\tau}$$

Here $\tau>0$. Forward Euler integration with step size $\Delta t$ leads to

$$x_{n+1} = x_n(1 - {\Delta t}/{\tau})$$

If you choose a time step that is larger than the time constant the above will diverge. Even if the time step is less than the time constant, the results can be quite bad. You choose a time step that is much less than the time constant.

In comparison, backwards Euler leads to

$$x_{n+1} = \frac{x_n}{1 + {\Delta t}/{\tau}}$$

which is always stable, even with a bad choice of a step size.
This interesting. I have only basic knowledge of Euler's method. Never been introduced to the term forward or backward method before.

Of course in practise people seldom use Euler forward method. It is only for illustration of numerical DE mthod. I do not know about Euler backward method.

So it is all about stability in using Euler backward method. But it look like that this method only work if we could algebraically solve yn=yn+1-hy'n+1 for yn+1. That a big problem. Not many DE can be solve with this method.

Another curious question. Do they also have such thing as Modified Euler Backward Method?

#### Nick89

But it look like that this method only work if we could algebraically solve yn=yn+1-hy'n+1 for yn+1. That a big problem. Not many DE can be solve with this method.
Actually, you can just enter the DE for y'n+1 which will probably be a function of yn+1. The 'problem' is now reduced simply to solving a (possibly non-linear) equation. Many equations like this can be solved analytically, and if they cannot there are always numerical methods.

#### matematikawan

Still need to be convinced on using Euler backward method. Please explain how to solve this IVP.

y' = sin y subject to let say y(0)=1.

yn = yn+1 - h sin(yn+1)

Then what next? It is not easy to express yn+1 in term of yn.

#### D H

Staff Emeritus
The typical approach to using implicit methods for non-linear ODEs is to linearize the derivative function, y'=f(x,y), where x is the independent variable. The method is no longer purely implicit and is better called semi-implicit. Implicit Euler integration uses

$$y_{n+1}=y_n+h f(x_{n+1},y_{n+1})$$

Suppose the derivative function does not depend on the dependent variable. Linearizing the derivative function about yn,

$$f(y_{n+1}) \approx f(y_n) + (y_{n+1}-y_n)\left.\frac{\partial f}{\partial y}\right|_{y_n}$$

With this, the semi-implicit Euler integration scheme becomes

$$y_{n+1} = y_n + h \left( f(y_n)+(y_{n+1}-y_n)\left.\frac{\partial f}{\partial y}\right|_{y_n} \right)$$

Solving for yn+1,

$$y_{n+1} = y_n + h\left(1-h\left.\frac{\partial f}{\partial y}\right|_{y_n}\right)^{-1} f(y_n)$$

In the problem at hand, y'=sin y, this becomes

$$y_{n+1} = y_n + h\left(1-h\cos y_n\right)^{-1} \sin y_n$$

#### matematikawan

Thanks DH for the explaination.

So there can be two type of Euler Backward Method.
(i) the implicit yn=yn+1 - h y'n+1 and
(ii) the semi-implicit $$y_{n+1}=y_n+h f(x_{n+1},y_{n+1})$$ .

Or is there others?

I hope they have the theorem which state that this method is always stable. I will try later to write a matlab script later to see what result this method gives. Has anyone seen the analytic solution for
$$\frac{dy}{dx} = \sin(y)$$
before?

#### matematikawan

I have written the matlab script. If only I know the analytic solution to compare.

% Solving IVP y'=sin(y) , y(0)=1
% Using Euler Backward Method

h=0.1; % step size
y(1)=1; % matlab index start from 1. y0=y(1)
xi=0; xf=10; % domain [xi , xf] num_pt=(xf-xi)/h+1 = 101
for n=2:101
y(n)=y(n-1)+h*sin(y(n-1))/(1-h*cos(y(n-1)));
end
x=xi:h:xf;
plot(x,y)

#### Attachments

• 8.5 KB Views: 501

#### matematikawan

OK got it.

$$\int \csc (y) dy = x + c$$
$$\log (\tan (\frac{y}{2})) = x + c$$

ref: http://integrals.wolfram.com/index.jsp?expr=Csc[x]&random=false

$$y = 2\tan^{-1} (\tan (0.5) e^x)$$ #### ametista82

I have the same problem
I have this force....
F= k*l-ln/ln*ppj/||p*pj|| where p and pj is the position of particle p and pj, l and ln are the length of a spring. how can use the implicit integration in this ode?
and If i have a costant force that depends only by an constant acceleration? thank you for your help

#### matematikawan

Not sure whether I could still remember this stuff. If you could please write down the related differential equation (if possible in latex) and initial condition it will be more clearer your problem.

#### ametista82

Sorry...I read a Baraff's paper about implicit integration..(sorry but I don't use Latex..)...so by this paper:

deltax = h(v0 +deltav) where deltax=x(t0+h)-x(t0)

(I -hdf/dv -h^2df/dx)deltav=h(f0+hdf/dxv0)

so the problem is to solve for deltav, so if my force is

F=fg + k* l-ln/ln ppj/|ppj| where fg is gravitational force, l is the natural length of the spring and ppj are the position of two particle betwen the string.

So how can use implicit integration? And another question if I have only a costant force as F=ma how can use an implicit integration?

#### klunkerz

Solving this for y_n+1 by hand, I obviously get a quadratic equation for y_n+1, for which there are two solutions:
$$y_{n+1} = \frac{ -y_n \pm \sqrt{y_n^2 - 2h}}{2}$$

Hi,
How to solve the above equation by hand?
Can anyone please explain it to me?

#### matematikawan

I though everybody has access to at least a calculator to compute square root of a number.

Referring back to DH post (the 3rd post)
$$y_{n+1} = y_n\left(\frac{1 + \sqrt{1 - 2h/y_n^2}}{2}\right)$$

Apply the binomial expansion
$$y_{n+1} = \frac{y_n}{2}\left(1 + 1-\frac{2h}{2y_n^2}+...\right)$$
$$\approx y_n\left(1 -\frac{h}{2y_n^2}\right)$$

which can be calculated using a 'slide rule'/table 