Consider the time-domain differential equation:
d t d x ( t ) = x ( t )
Let us discretize the function and time itself. Function x 1 is the explicit Euler approximation, and x 2 is the implicit Euler approximation. Let the time step Δ t = 0 . 0 1
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
In the above equations, the k subscript denotes the present processing step, and the 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
When t = 3 , what is x 2 − x 1 , rounded to the nearest integer multiple of 0 . 1 ?
Bonus: Plot x 1 , x 2 , and the analytical solution together on a graph vs. time. Would the accuracy improve if we averaged x 1 and x 2 ? What would the resulting integration method be?
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.
Thanks for the solution. The average is indeed the trapezoidal method
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?
Log in to reply
Explicit Euler:
x k = x k − 1 + x ˙ k − 1 Δ t
Implicit Euler:
x k = x k − 1 + x ˙ k Δ t
Average of the two, which is the trapezoidal rule:
x k = x k − 1 + 2 1 x ˙ k − 1 Δ t + 2 1 x ˙ k Δ t
Log in to reply
@Steven Chase – I see. I have not used this method in practice. This method also appears to be implicit.
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 )
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.
@Karan Chatrath – In exchange for the more difficult solution, implicit methods are more numerically stable.
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.
Log in to reply
@Karan Chatrath – One thing I want to explore next is the case where x ˙ = x 2 . In the implicit case, you have to solve a quadratic for x k . Which root do you choose, since the quadratic equation gives two solutions?
Log in to reply
@Steven Chase – In this case, The roots of the quadratic equation lead to the iteration scheme:
x k + 1 = 2 Δ t 1 ± 1 − 4 Δ t x k
If: x k + 1 = 2 Δ t 1 + 1 − 4 Δ t x k
Then the limit as Δ t → 0 does not lead to x k + 1 = x k . However, using the second root and evaluating the same limit results in x k + 1 = x k which is expected.
x k + 1 = 2 Δ t 1 − 1 − 4 Δ t x k
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.
Problem Loading...
Note Loading...
Set Loading...
For the explicit Euler method:
x 1 ( k + 1 ) = x 1 ( k ) + Δ t x 1 ( k ) ⟹ x 1 ( k + 1 ) = ( 1 + Δ t ) x 1 ( k ) At the initial step: x 1 ( 1 ) = 1 . This gives the general formula on iteration:
x 1 ( k + 1 ) = ( 1 + Δ t ) k ⟹ x 1 ( k ) = ( 1 + Δ t ) k − 1
For the implicit Euler method:
x 2 ( k + 1 ) = x 2 ( k ) + Δ t x 2 ( k + 1 ) ⟹ x 2 ( k + 1 ) = ( 1 − Δ t ) x 2 ( k ) At the initial step: x 2 ( 1 ) = 1 . This gives the general formula on iteration:
x 2 ( k + 1 ) = ( 1 − Δ t ) k 1 ⟹ x 2 ( k ) = ( 1 − Δ t ) k − 1 1
Now, at t = 3 , k = 3 0 1 , which leads to the answer:
x 2 ( 3 0 1 ) − x 1 ( 3 0 1 ) = ( 1 + 0 . 0 1 ) 3 0 0 − ( 1 − 0 . 0 1 ) 3 0 0 1 ≈ 0 . 6
Averaging the two methods gives:
x 3 ( k ) = 2 1 ( ( 1 + Δ t ) k − 1 + ( 1 − Δ t ) k − 1 1 )
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 are:
∣ x 1 ( 3 0 1 ) − e 3 ∣ ≈ 0 . 2 9 7 1 ∣ x 2 ( 3 0 1 ) − e 3 ∣ ≈ 0 . 3 0 5 6 ∣ x 3 ( 3 0 1 ) − e 3 ∣ ≈ 4 . 2 6 8 9 × 1 0 − 3