Explicit vs Implicit 1st order Euler method for Numerical Integration

In this project, I will discuss the necessity for an implicit numerical scheme and its advantages over an explicit one.

For this demonstration, I will use the first order Euler Schemes for Numerical Integration as it is the easiest to use and understand,

The first order Euler Numerical scheme is derived from the Taylors Series, in which only the first term is considered, and the rest are ignored, hence the error of both these schemes will be of second order

From the equation of the Taylor series as shown here,

we just take into consideration the first term, h is the step size, Rearranging the equation,

`x_(i+1) = x_(i) + hcdotf(x_(i))`

This is called the Explicit Euler method, where we use data available at (i)th point to calculate the unknown value at the (i+1)th point.

The other alternative for this method is called the Implicit Euler Method, here converse to the other method we solve the non-linear equation which arises by formulating the expression in the below-shown way, using numerical root finding methods.

`x_(i+1) = x_(i) + hcdotf(x_(i+1))`

The second method is substantially more expensive to compute than the explicit method, but is it useful?? I shall try to explain the usefulness and effectiveness of a implicit method over a conventional explicit method using an example.

Let us take us take the example equation,

`dx/dt=-1000*x - e^-t`

This equation is solvable by hand, if done so we get the answer as follows if we use the initial condition 0 at time 0 we get,

`x = e^(1000*t)*(1-e^(999*t))/999`

Now, we will write a program in Python to test the integrity of the solutions of the implicit, explicit methods against the exact analytical solution.

import math 
import matplotlib.pyplot as plt
import numpy as np


def f(x,t):
	
	eqn = -1000*x - np.exp(-t)

	return eqn

def anal(tspan,dt):
	
	t = np.arange(0,tspan,dt)
	x = np.zeros(len(t))
	
	for i in range(0,len(t)):

		#term1 = 	
		x[i] = np.exp(-1000*t[i])*(1-np.exp(999*t[i]))/999

	return [t,x]

def explicit_euler(intx,tspan,dt):
	
	t = np.arange(0,tspan,dt)
	x = np.zeros(len(t))
	x[0] = intx
	
	for i in range(1,len(t)):
	
		x[i] = x[i-1] + dt*f(x[i-1],t[i-1]) 

	return [t,x]

def implicit_euler(intx,tspan,dt):
	
	t = np.arange(0,tspan,dt)
	x = np.zeros(len(t))
	x[0] = intx
	
	for i in range(1,len(t)):
	
		x[i] = (x[i-1] - np.exp(-t[i])*dt)/(1+1000*dt)

	return [t,x]


t_start = 0
t_end = 0.1
n = 100

dt = (t_end - t_start)/(n-1)

ta,anl = anal(t_end,dt)
ti,imp = implicit_euler(0,t_end,dt)
te,exp = explicit_euler(0,t_end,dt)

plt.plot(ti,imp,'r')
plt.plot(te,exp,'b')
plt.plot(ta,anl,'g')
plt.legend(["Implicit Solution","Explicit Solution","Analytical Solution"])
plt.show()

Here we write the program to find the solution between the time interval, 0 to 0.1 seconds and we assume the solution is calculated at 100 points between the set time interval,

The solution when we have 100 steps in the time interval,

As you can see both the schemes (implicit and explicit) seem to have a very similar result to the Exact Analytical solution. Though some unregister can be seen in the explicit scheme.

To further investigate this issue, let us compute the same program but with 75 steps in the time interval, the results are shown below,

As you can see here the error of the explicit scheme has increased in terms of greater oscillations, this is called and instability error, and because of the lower number of points we have a higher error in the explicit scheme,however after some time the error does converge and matches the solution with the analytical solution, whereas the implicit scheme still remains very consistent.

Now, let us further push this and make the number of points between the interval as 50, Now the graphs obtained is as follows,

As you can see here as expected the explicit solution still deviates, but unlike the last two cases, it doesn't consistenly converge with the analytical solution after some time period, it just consistenly grows in scale to point where the solution set offered by the explicit scheme is completely wrong.

This is the issue with the explicit scheme, in certain cases, like the system portrayed by this equation, the solution rapidly changes in a set interval of the solution domain, in this case at the beginning of the solution.

When the solution to an ODE changes rapidly in a defined interval, this system or ODE is referred to as a stiff ODE/system,

Explicit schemes are pretty much useless in these stiff regions as they can completely mess up the solution if the step size isn't low enough, This is where the implicit scheme shines through,

Becasue the implicit shceme computes the slope at the next point, it doesnt encounter this error.