Implicit vs. Explicit Euler Integration

Calculus Level pending

Consider the time-domain differential equation:

d d t x ( t ) = x ( t ) \frac{d}{dt} x(t) = x(t)

Let us discretize the function and time itself. Function x 1 x_1 is the explicit Euler approximation, and x 2 x_2 is the implicit Euler approximation. Let the time step Δ t = 0.01 \Delta t = 0.01

x 1 k = x 1 k 1 + x ˙ 1 k 1 Δ t = x 1 k 1 + x 1 k 1 Δ t x 2 k = x 2 k 1 + x ˙ 2 k Δ t = x 2 k 1 + x 2 k Δ t t k = t k 1 + Δ t x_{1 k} = x_{1 \, k-1} + \dot{x}_{1 \, k-1} \Delta t = x_{1 \, k-1} + x_{1 \, k-1} \Delta t \\ x_{2 k} = x_{2 \, k-1} + \dot{x}_{2 \, k} \Delta t = x_{2 \, k-1} + x_{2 \, k} \Delta t \\ t_k = t_{k - 1} + \Delta t

In the above equations, the k k subscript denotes the present processing step, and the k 1 k-1 subscript denotes the previous processing step. The simulation is initialized as follows:

t 0 = 0 x 1 0 = 1 x 2 0 = 1 t_0 = 0 \\ x_{1 \, 0} = 1 \\ x_{2 \, 0} = 1 \\

When t = 3 t = 3 , what is x 2 x 1 x_2 - x_1 , rounded to the nearest integer multiple of 0.1 0.1 ?

Bonus: Plot x 1 x_1 , x 2 x_2 , and the analytical solution together on a graph vs. time. Would the accuracy improve if we averaged x 1 x_1 and x 2 x_2 ? What would the resulting integration method be?


The answer is 0.6.

This section requires Javascript.
You are seeing this because something didn't load right. We suggest you, (a) try refreshing the page, (b) enabling javascript if it is disabled on your browser and, finally, (c) loading the non-javascript version of this page . We're sorry about the hassle.

1 solution

Karan Chatrath
Mar 15, 2021

For the explicit Euler method:

x 1 ( k + 1 ) = x 1 ( k ) + Δ t x 1 ( k ) x_1(k+1) = x_1(k) + \Delta t x_1(k) x 1 ( k + 1 ) = ( 1 + Δ t ) x 1 ( k ) \implies x_1(k+1) = (1 + \Delta t) x_1(k) At the initial step: x 1 ( 1 ) = 1 x_1(1) = 1 . This gives the general formula on iteration:

x 1 ( k + 1 ) = ( 1 + Δ t ) k x_1(k+1) = (1 + \Delta t)^k x 1 ( k ) = ( 1 + Δ t ) k 1 \implies x_1(k) = (1 + \Delta t)^{k-1}

For the implicit Euler method:

x 2 ( k + 1 ) = x 2 ( k ) + Δ t x 2 ( k + 1 ) x_2(k+1) = x_2(k) + \Delta t x_2(k+1) x 2 ( k + 1 ) = x 2 ( k ) ( 1 Δ t ) \implies x_2(k+1) = \frac{x_2(k)}{(1 - \Delta t)} At the initial step: x 2 ( 1 ) = 1 x_2(1) = 1 . This gives the general formula on iteration:

x 2 ( k + 1 ) = 1 ( 1 Δ t ) k x_2(k+1) = \frac{1}{(1 - \Delta t)^k} x 2 ( k ) = 1 ( 1 Δ t ) k 1 \implies x_2(k) = \frac{1}{(1 - \Delta t)^{k-1}}

Now, at t = 3 t=3 , k = 301 k=301 , which leads to the answer:

x 2 ( 301 ) x 1 ( 301 ) = ( 1 + 0.01 ) 300 1 ( 1 0.01 ) 300 0.6 x_2(301) - x_1(301) = (1 + 0.01)^{300} - \frac{1}{(1 - 0.01)^{300}} \approx 0.6

Averaging the two methods gives:

x 3 ( k ) = 1 2 ( ( 1 + Δ t ) k 1 + 1 ( 1 Δ t ) k 1 ) x_3(k) = \frac{1}{2}\left((1 + \Delta t)^{k-1} + \frac{1}{(1 - \Delta t)^{k-1}}\right)

From this formula, I cannot really identify a method, but I am inclined to think it is the Trapezoidal method given the average involved in the formula. But I am unsure.

The average plot clearly appears to be overlapping the exact solution. For a little more perspective, the absolute errors between the exact solution and each method evaluated at t = 3 t=3 are:

x 1 ( 301 ) e 3 0.2971 \lvert x_1(301) - \mathrm{e}^{3} \rvert \approx 0.2971 x 2 ( 301 ) e 3 0.3056 \lvert x_2(301) - \mathrm{e}^{3} \rvert \approx 0.3056 x 3 ( 301 ) e 3 4.2689 × 1 0 3 \lvert x_3(301) - \mathrm{e}^{3} \rvert \approx 4.2689 \times 10^{-3}

Thanks for the solution. The average is indeed the trapezoidal method

Steven Chase - 2 months, 4 weeks ago

Log in to reply

Thanks. I am unable to derive the iteration formula for the trapezoidal rule from the steps I performed. Could you show it, if possible?

Karan Chatrath - 2 months, 4 weeks ago

Log in to reply

Explicit Euler:

x k = x k 1 + x ˙ k 1 Δ t x_k = x_{k-1} + \dot{x}_{k-1} \Delta t

Implicit Euler:

x k = x k 1 + x ˙ k Δ t x_k = x_{k-1} + \dot{x}_{k} \Delta t

Average of the two, which is the trapezoidal rule:

x k = x k 1 + 1 2 x ˙ k 1 Δ t + 1 2 x ˙ k Δ t x_k = x_{k-1} + \frac{1}{2} \dot{x}_{k-1} \Delta t + \frac{1}{2} \dot{x}_{k} \Delta t

Steven Chase - 2 months, 4 weeks ago

Log in to reply

@Steven Chase I see. I have not used this method in practice. This method also appears to be implicit.

Karan Chatrath - 2 months, 4 weeks ago

Log in to reply

@Karan Chatrath Usually what happens is you have the numerical integration method, and then you have a state-space of the form:

x ˙ = f ( x , u ) \dot{x} = f(x,u)

So you can replace the time derivative terms with some function of the state variable and the forcing signal. And then in the explicit case, the right side consists entirely of history terms, and is trivial to solve. In the implicit case, you have to solve a more complicated equation for the present value, since it is present on both sides. The trap rule is implicit in that sense.

Steven Chase - 2 months, 4 weeks ago

@Karan Chatrath In exchange for the more difficult solution, implicit methods are more numerically stable.

Steven Chase - 2 months, 4 weeks ago

Log in to reply

@Steven Chase Thanks for the comments. I recall some of your exercises involving the computation of the stability region areas of implicit and explicit Euler schemes.

Karan Chatrath - 2 months, 4 weeks ago

Log in to reply

@Karan Chatrath One thing I want to explore next is the case where x ˙ = x 2 \dot{x } = x^2 . In the implicit case, you have to solve a quadratic for x k x_k . Which root do you choose, since the quadratic equation gives two solutions?

Steven Chase - 2 months, 4 weeks ago

Log in to reply

@Steven Chase In this case, The roots of the quadratic equation lead to the iteration scheme:

x k + 1 = 1 ± 1 4 Δ t x k 2 Δ t x_{k+1} = \frac{1 \pm \sqrt{1 - 4 \Delta t x_k}}{2\Delta t}

If: x k + 1 = 1 + 1 4 Δ t x k 2 Δ t x_{k+1} = \frac{1 + \sqrt{1 - 4 \Delta t x_k}}{2\Delta t}

Then the limit as Δ t 0 \Delta t \to 0 does not lead to x k + 1 = x k x_{k+1}=x_{k} . However, using the second root and evaluating the same limit results in x k + 1 = x k x_{k+1}=x_{k} which is expected.

x k + 1 = 1 1 4 Δ t x k 2 Δ t x_{k+1} = \frac{1 - \sqrt{1 - 4 \Delta t x_k}}{2\Delta t}

Hence, the above root is the appropriate choice. This gets tricky for highly nonlinear equations and would perhaps require Newton-Raphson or some other numerical root searching scheme.

Karan Chatrath - 2 months, 4 weeks ago

0 pending reports

×

Problem Loading...

Note Loading...

Set Loading...