Friday, 10 March 2017

1-D Parabolic Equation PDE schemes

1-D Heat conduction equation (with no internal heat generation)

$\frac{\delta T}{\delta t} = \alpha \left ( \frac{\delta ^2 T}{\delta x^2} \right )$

Before we dwell any further in the discretization of the above PDE we need to know 3 important terms.

  1. Consistency
  2. Stability 
  3. Convergence
Here a finite difference scheme is said to be consistent when as $\Delta x,\Delta t \rightarrow 0$ , then the scheme  $\rightarrow$ to the original PDE.

A scheme is said to be stable when the ratio of the error propagation from one time step to another is less than $1$.

A scheme is said to converge when the above to conditions are met as well as the solution is close to the actual solution (analytic or numerical) of the PDE .
  
Note : Implicit schemes are mostly unconditionally stable ,while the explicit schemes are conditionally stable or unconditionally unstable
Explicit schemes
FTCS scheme

$\frac{T_{i} ^{n+1} -T_{i} ^{n}}{\Delta t} = \frac{T_{i+1} ^{n} -2T_{i} ^{n}+T_{i-1} ^{n}}{(\Delta x)^2}  + O(\Delta t,(\Delta x)^2)$

This scheme is conditionally stable , consistent ,convergent ,second order accurate in space,1st order accurate in time.
For more information on the stability criterion please check the tutorial 1 under the CFD tab.

CTCS scheme

$\frac{T_{i} ^{n+1} -T_{i} ^{n-1}}{2 \Delta t} = \frac{T_{i+1} ^{n} -2T_{i} ^{n}+T_{i-1} ^{n}}{(\Delta x)^2}  + O((\Delta t)^2,(\Delta x)^2)$

This scheme is unconditionally unstable. Hence this scheme is absolutely useless .Therefore the Dufort -Frankel scheme emerged.

Dufort Frankel scheme

$\frac{T_{i} ^{n+1} -T_{i} ^{n-1}}{2 \Delta t} = \frac{T_{i+1} ^{n} -2\left (\frac{T_{i} ^{n+1} +T_{i} ^{n-1}}{2}\right )+T_{i-1} ^{n}}{(\Delta x)^2}  + O((\Delta t)^2,(\Delta x)^2,(\frac{\Delta t}{\Delta x})^2)$

This scheme is unconditionally stable. Hence this scheme is brilliant . The error term can be easily derived , do this as an exercise.

Implicit schemes
BTCS scheme (laasonen scheme)

$\frac{T_{i} ^{n+1} -T_{i} ^{n}}{\Delta t} = \frac{T_{i+1} ^{n+1} -2T_{i} ^{n+1}+T_{i-1} ^{n+1}}{(\Delta x)^2}  + O(\Delta t,(\Delta x)^2)$

Unconditionally stable scheme . This scheme forms a tridiagonal matrix , which can be solved using matrix inversion ,or the Gauss elimination method or the TDMA (thomas) algorithm. For more information on the TDMA algorithm check tutorial 2 under the CFD tab .


Crank Nicolson scheme 

$\frac{T_{i} ^{n+1} -T_{i} ^{n}}{\Delta t} =  \frac{1}{2}\left (\frac{T_{i+1} ^{n+1} -2T_{i} ^{n+1}+T_{i-1} ^{n+1}}{(\Delta x)^2} +\frac{T_{i+1} ^{n} -2T_{i} ^{n}+T_{i-1} ^{n}}{(\Delta x)^2}    \right ) + O((\Delta t )^2,(\Delta x)^2)$

This is a really interesting scheme centered around the time step $n + 1/2$ . As can be seen the above scheme is second order accurate in time as well.
Unconditionally stable scheme . This scheme forms a tridiagonal matrix , which can be solved using matrix inversion ,or the Gauss elimination method or the TDMA (thomas) algorithm


Beta Formulation
Finally comes the $\beta $ formulation , this combines all the above mentioned schemes in a single shot .

$\frac{T_{i} ^{n+1} -T_{i} ^{n}}{\Delta t} =  \beta \left (\frac{T_{i+1} ^{n+1} -2T_{i} ^{n+1}+T_{i-1} ^{n+1}}{(\Delta x)^2}  \right )+(1-\beta) \left ( \frac{T_{i+1} ^{n} -2T_{i} ^{n}+T_{i-1} ^{n}}{(\Delta x)^2}   \right )  $

This scheme behaves like an FTCS scheme for $\beta < 1/2$
This scheme behaves like an BTCS scheme for $\beta > 1/2$
This scheme behaves like an Crank Nicolson scheme for $\beta = 1/2$ 


2D Parabolic PDEs

Here let us now consider a plate whose initial conditions are specified .
FTCS scheme
So Let us discretize the governing 2 -Dimensional Parabolic PDE using FTCS scheme .
$$\frac{T_{i,j}^{n+1}-T_{i,j}^{n}}{\Delta t} = \alpha \left (\frac{T_{i+1,j}^{n}-2T_{i,j}^{n}+T_{i-1,j}^{n}}{\Delta x ^2}  +\frac{T_{i,j+1}^{n}-2T_{i,j}^{n}+T_{i,j-1}^{n}}{\Delta y ^2} \right ) $$
Since this is an explicit scheme it becomes easy to calculate the solutions of $T$ as a function of time.(ie)

$$T_{i,j}^{n+1} = T_{i,j}^{n}+ \Delta t \alpha \left (\frac{T_{i+1,j}^{n}-2T_{i,j}^{n}+T_{i-1,j}^{n}}{\Delta x ^2}  +\frac{T_{i,j+1}^{n}-2T_{i,j}^{n}+T_{i,j-1}^{n}}{\Delta y ^2} \right )$$
For computational ease $\Delta x = \Delta y$ can be used. Now let us see the stability criterion of the scheme below.
Stability Condition
It can be derived to show that the stability for the 2 D heat conduction equation is similar to that of the 1-D heat equation ,hence 

$$ \alpha \Delta (t) \left (\frac{1}{(\Delta x )^2} + \frac{1}{(\Delta x )^2}  \right  )\leq \frac{1}{2}$$
Check Computational Fluid Dynamics by  Klaus A. Hoffmann (Author), Steve T. Chiang (Author) for reference.

This can be further reduced to 
$$ \alpha \Delta (t) \left  (\frac{1}{(\Delta x )^2} \right )\leq \frac{1}{4}$$ when $\Delta x = \Delta y$   

BTCS scheme
The same can be done in BTCS scheme to obtain 
$$\frac{T_{i,j}^{n+1}-T_{i,j}^{n}}{\Delta t} = \alpha \left (\frac{T_{i+1,j}^{n+1}-2T_{i,j}^{n+1}+T_{i-1,j}^{n+1}}{\Delta x ^2}  +\frac{T_{i,j+1}^{n+1}-2T_{i,j}^{n+1}+T_{i,j-1}^{n+1}}{\Delta y ^2} \right ) $$

The above equation as can be seen that it forms a pentadiagonal matrix which can be solved by matrix inversion which is computationally expensive ,hence this scheme wont be used most of the time. Hence emerges the ADI scheme

ADI scheme
The ADI (alternate direction implicit method) is the execution of the implicit schemes at one direction at time steps  $n,n +1/2$
$$\frac{T_{i,j}^{n+1/2}-T_{i,j}^{n}}{\Delta t} = \alpha \left (\frac{T_{i+1,j}^{n+1/2}-2T_{i,j}^{n+1/2}+T_{i-1,j}^{n+1/2}}{\Delta x ^2}  +\frac{T_{i,j+1}^{n}-2T_{i,j}^{n}+T_{i,j-1}^{n}}{\Delta y ^2} \right ) $$
$$\frac{T_{i,j}^{n+1}-T_{i,j}^{n+1/2}}{\Delta t} = \alpha \left (\frac{T_{i+1,j}^{n+1}-2T_{i,j}^{n+1}+T_{i-1,j}^{n+1}}{\Delta x ^2}  +\frac{T_{i,j+1}^{n+1/2}-2T_{i,j}^{n+1/2}+T_{i,j-1}^{n+1/2}}{\Delta y ^2} \right ) $$

The above equations as can be seen that it forms a tridiagonal matrix which can easily be solved by TDMA algorithm . The first equation is solved for finding the solutions for intermediate time step $n +1/2$ then this is used to find the solution at time step $n+1$ .

CFD - Finite Difference schemes (Intro)

First Order Derivatives


Forward Difference

$\frac{\delta u}{\delta x} =  \frac{u_{i+1} - u_{i}}{\Delta x} + O(\Delta x)$

Derivation :
  
$u_{i+1}=u_{i} + \left .\frac{\Delta x}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + \left .\frac{(\Delta x)^2}{2 !} \frac{\delta ^2u}{\delta x ^2} \right \vert _{i} + \left .\frac{(\Delta x)^3}{3 !} \frac{\delta ^3u}{\delta x ^3} \right \vert _{i} +O((\Delta x)^4) $
$\implies \frac{u_{i+1}-u_{i}}{\Delta x} = \left .\frac{1}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + \left .\frac{\Delta x}{2 !} \frac{\delta ^2u}{\delta x ^2} \right \vert _{i} + \left .\frac{(\Delta x)^2}{3 !} \frac{\delta ^3u}{\delta x ^3} \right \vert _{i} +O((\Delta x)^3) $
$\implies \frac{u_{i+1}-u_{i}}{\Delta x} = \left .\frac{1}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + O((\Delta x)) $


Backward Difference

$\frac{\delta u}{\delta x} =  \frac{u_{i} - u_{i-1}}{\Delta x} + O(\Delta x)$

Derivation :

  $u_{i-1}=u_{i} - \left .\frac{\Delta x}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + \left .\frac{(\Delta x)^2}{2 !} \frac{\delta ^2u}{\delta x ^2} \right \vert _{i} - \left .\frac{(\Delta x)^3}{3 !} \frac{\delta ^3u}{\delta x ^3} \right \vert _{i} +O((\Delta x)^4) $
$\implies \frac{u_{i}-u_{i-1}}{\Delta x} = \left .\frac{1}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + O((\Delta x)) $


Central Difference

$\frac{\delta u}{\delta x} =  \frac{u_{i+1} - u_{i-1}}{ 2\Delta x } + O((\Delta x)^2)$

Derivation :

$-u_{i+1}+u_{i} + \left .\frac{\Delta x}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + \left .\frac{(\Delta x)^2}{2 !} \frac{\delta ^2u}{\delta x ^2} \right \vert _{i} + \left .\frac{(\Delta x)^3}{3 !} \frac{\delta ^3u}{\delta x ^3} \right \vert _{i} +O((\Delta x)^4)  =0      -(1) $ 

  $-u_{i-1}+u_{i} - \left .\frac{\Delta x}{1 !} \frac{\delta u}{\delta x} \right \vert _{i} + \left .\frac{(\Delta x)^2}{2 !} \frac{\delta ^2u}{\delta x ^2} \right \vert _{i} - \left .\frac{(\Delta x)^3}{3 !} \frac{\delta ^3u}{\delta x ^3} \right \vert _{i} +O((\Delta x)^4) =0-(2)$ 

$(\alpha eqn_1 +\beta eqn_2) = ..... + 0*\left .\frac{(\Delta x)^2}{2 !} \frac{\delta ^2u}{\delta x ^2} \right \vert _{i}  +O((\Delta x)^3)  =0 -(3)$
Solving for $\alpha,\beta $ substituting in equation 3 , rearranging we get the above result. 


Similarly using the above method, we can obtain higher order accuracies for  forward, backward difference schemes as well.

Second Order Derivatives 

Forward Difference

$\frac{\delta^2 u}{\delta x^2} =  \frac{u_{i+2}-2u_{i+1} + u_{i}}{ \Delta x ^2} + O((\Delta x))$

Central Difference

$\frac{\delta^2 u}{\delta x^2} =  \frac{u_{i}-2u_{i-1} + u_{i-2}}{ \Delta x ^2} + O((\Delta x))$

Central Difference

$\frac{\delta^2 u}{\delta x^2} =  \frac{u_{i+1}-2u_{i} + u_{i-1}}{ \Delta x ^2} + O((\Delta x)^2)$

Note : Derive the above equations for your own practice.

Other methods of derivation

There are many methods using which the above schemes can be derived . One is that of the polynomial method ,  the other one is discussed for a simple case below :

 $$\frac{\delta ^n f}{\delta x ^n} = \frac{\delta ^{n-1} }{\delta x ^{n-1}} \left (\frac{\delta f}{\delta x }\right )    $$
Hence for a second order derivative (using forward difference) we get

$$\frac{\delta ^2 f}{\delta x ^2} = \frac{\delta  }{\delta x  } \left (\frac{f_{i+1} - f_{i}}{\Delta x} \right )   $$
$$\frac{\delta ^2 f}{\delta x ^2} = \left (\frac{\delta f_{i+1}  }{\delta x}-\frac{\delta f_{i}  }{\delta x}\right )\frac{1}{\Delta x}   $$
$$\frac{\delta ^2 f}{\delta x ^2} = \left (\frac{ f_{i+1}  -f_{i+1}  }{\Delta x}-\frac{ f_{i+1} - f_{i}  }{\Delta x}\right )\frac{1}{\Delta x}   $$

$$\frac{\delta ^2 f}{\delta x ^2} = \left (\frac{ f_{i+1}  - 2f_{i+1} +f_{i} }{\Delta x ^2} \right )$$

Polynomial method

Let there be a polynomial $f(x)   = ax^2+bx+c$ .
$\frac{\delta f}{\delta x} = 2ax+b    -(1)$
$\frac{\delta ^2 f}{\delta x ^2} = 2a  -(2)$

Let at $x =0 ,f(x) = f_{i},x = \Delta x ,f(x) = f_{i+1} ,x =2 \Delta x ,f(x) = f_{i+2} $.
Now solve for $a , b,c$ in terms of $f_{i},f_{i+1},f_{i+2},\Delta x$ 

Substituting $a , b$ back in equation 1, 2 will give you the appropriate derivatives.

Try the above calculations for practice.