Week 5 - Literature review: ODE Stability

Aim: Write a program to check explicit and implicit schemes, and check the stability of ODE.

Objective:

• To determine which scheme is stable.
• Plot for unstability condition.
• Define how can we derive stability for all problems.

Introduction:

When we are deriving finite difference approximation of a exact derivative we are introducing error in a system this error evolves with time according to governing equation of ordinary differential equation, and this error goes on increasing exponentially which can be seen as sharp peak in the plots.

There are two techniques to calculate the ODE :

• Implicit method
• Explicit method

Implicit method:In this method the values at the right handside of any function are unknown.This scheme is more stable.

Explicit method:In this method the values at the right side of the functions are known.scheme is less stable.

To define this schemes we will write a python program for both schemes for a ODE and plot a results.

The ODE is dy/dt = -1000y - exp(-t)

Python code for stability analysis:

import math
import matplotlib.pyplot as plt
import numpy as np 
#function for analytic solution
def analytic(t_span,dt):
	t = np.arange(0,t_span,dt)
	y = np.zeros(len(t))
	for i in range(1,len(t)):
		y[i] = (np.exp(-1000*t[i])*(1 - np.exp(999*t[i])))/999
	return(t,y)

#function for explicit method
def explicit(t_span,dt):
	t = np.arange(0,t_span,dt)
	y = np.zeros(len(t))
	for i in range(1,len(t)):
		y[i] = y[i-1]- dt*(1000*y[i-1] + np.exp(-t[i-1]))
	return(t,y)	

#function for implicit method
def implicit(t_span,dt):
	t = np.arange(0,t_span,dt)
	y = np.zeros(len(t))
	for i in range(1,len(t)):
		y[i] = (y[i-1] - np.exp(-t[i-1])*dt)/(1+1000*dt)
	return(t,y)	

#simulation time
t_span = 0.1
#time step
dt = 1.0e-3

#calling function to store results in variables
t, y_explicit= explicit(t_span,dt)
t, y_implicit= implicit(t_span,dt)
t, y_analytic= analytic(t_span,dt)

#plotting
plt.plot(t , y_implicit, label = 'implicit')
plt.plot(t , y_explicit, label = 'explicit')
plt.plot(t , y_analytic, label = 'analytic')
plt.legend(bbox_to_anchor =(0.85, 1.15), ncol = 3) 
plt.show()

PLOTS

Plot for different schemes at dt = 1e-3:

From above plot we can give a conclusion that the implicit scheme results have less deviation from analytic solution compared to explicit scheme.

Plot for different schemes at dt = 0.2e-2:

After seeing above plot it can be visualized that the explicit scheme results are totally blown up giving no solution where as the implicit scheme still provide a solution which is near to analytic solution.

When does a ODE become unstable?

--> At large values of time step the ode generally get unstable, And the unstability is greatly seen when we use explicit method.

How does stability affect engineering simulations?

--> It disturbs whole solution leading to erroneius results.

Can stability condition be derived for all types of problems?

--> The stability condition can be dervied by using von nieman stability analysis and eigen value analysis

In von nieman analysis we use courant fredrich lewis number(CFL) number which defines whether solution will be stable or not.

Generally CFL<=1.