Implicit Compact Pade (ICP) Scheme

Approximation of First Derivative

Given the values of a function on a set of nodes the finite difference approximation to the derivative of the function is expressed as a linear combination of the given function values. For simplicity consider a uniformly spaced mesh where the nodes are indexed by \(i\). The independent variable at the nodes is \(x_{i}=h(i-1)\) for \(1\le i\le N\) and the function values at the nodes \(f_{i}=f(x_{i})\) are given. The finite difference approximation \({f}'_{i}\) to the first derivative \((df/dx)(x_{i})\) at the node \(i\) depends on the function values at nodes near \(i\). For second- and fourth-order central differences the approximation \({f}'_{i}\) depends on the sets \((f_{i-1},f_{i+1})\) and \((f_{i-2},f_{i-1},f_{i+1},f_{i+2})\), respectively. In the spectral methods, however, the value off \({f}'_{i}\) depends on all the nodal values. The Pade or compact finite difference schemes mimic this global dependence. The schemes presented here are generalizations of the Pade scheme.

These generalizations are derived by writing approximations of the form:

\[\begin{split}\begin{align} &\beta{f}'_{i-2}+\alpha{f}'_{i-1}+{f}'_{i}+\alpha{f}'_{i+1}+\beta{f}'_{i+2}\\ = &c\cfrac{{f}_{i+3}-{f}_{i-3}}{6h}+b\cfrac{{f}_{i+2}-{f}_{i-2}}{4h}+a\cfrac{{f}_{i+1}-{f}_{i-1}}{2h} \end{align}\end{split}\]

Approximation of Second Derivative

The derivation of compact approximations for the second derivative proceeds exactly analogous to the first derivative. Once again the starting point is a relation of the form

\[\begin{split}\begin{align} &\beta{f}''_{i-2}+\alpha{f}''_{i-1}+{f}''_{i}+\alpha{f}''_{i+1}+\beta{f}''_{i+2}\\ = &c\cfrac{{f}_{i+3}-2{f}_{i}+{f}_{i-3}}{9h^{2}}+b\cfrac{{f}_{i+2}-2{f}_{i}+{f}_{i-2}}{4h^{2}}+a\cfrac{{f}_{i+1}-2{f}_{i}+{f}_{i-1}}{h^{2}} \end{align}\end{split}\]

where \({f}''_{i}\) represents the finite difference approximation to the second derivative at node \(i\). The relations between the coefficients \(a,b,c\) and \(\alpha,\beta\) are derived by matching the Taylor series coefficients of various orders. The first unmatched coefficient determines the formal truncation error of the approximation. These constraints are:

\[\begin{split}\begin{align} a+b+c & = 1+2\alpha+2\beta\quad\text{(second order)}\\ a+2^{2}b+3^{2}c & = \cfrac{4!}{2!}(\alpha+2^{2}\beta)\quad\text{(fourth order)}\\ a+2^{4}b+3^{4}c & = \cfrac{6!}{4!}(\alpha+2^{4}\beta)\quad\text{(sixth order)}\\ a+2^{6}b+3^{6}c & = \cfrac{8!}{6!}(\alpha+2^{6}\beta)\quad\text{(eighth order)}\\ a+2^{8}b+3^{8}c & = \cfrac{10!}{8!}(\alpha+2^{8}\beta)\quad\text{(tenth order)}\\ \end{align}\end{split}\]

By choosing \(\beta=0\) and \(c = 0\) a one-parameter family of fourth-order schemes is generated. This family has

\[\beta=0,\quad c=0,\quad a=\cfrac{4}{3}(1-\alpha),\quad b=\cfrac{1}{3}(-1+10\alpha)\]

\[\begin{split}\begin{array}{c} \beta=0,c=0\\ a+2^{2}b+3^{2}c = \cfrac{4!}{2!}(\alpha+2^{2}\beta)\quad\text{(fourth order)}\\ \Rightarrow a+2^{2}b = \cfrac{4!}{2!}(\alpha)\quad\text{(fourth order)}\\ \Rightarrow a+4b = 4\times 3(\alpha)\quad\text{(fourth order)}\\ \cfrac{4}{3}(1-\alpha)+4\times \cfrac{1}{3}(-1+10\alpha) =\cfrac{4}{3}-\cfrac{4}{3}\alpha-\cfrac{4}{3}+\cfrac{4}{3}(10\alpha)=\cfrac{4}{3}(9\alpha)=4\times (3\alpha) \end{array}\end{split}\]

For \(\alpha=\cfrac{1}{10}\) the classical Pade scheme is recovered.

then

\[\alpha=\cfrac{1}{10},\quad \beta=0,\quad b=0, \quad c=0,\quad a=\cfrac{4}{3}(1-\cfrac{1}{10})=\cfrac{12}{10}\]

\[\cfrac{1}{10}{f}''_{i-1}+{f}''_{i}+\cfrac{1}{10}{f}''_{i+1}=\cfrac{12}{10}\left(\cfrac{{f}_{i+1}-2{f}_{i}+{f}_{i-1}}{h^{2}}\right)\]

\[\cfrac{1}{12}{f}''_{i-1}+\cfrac{10}{12}{f}''_{i}+\cfrac{1}{12}{f}''_{i+1}=\cfrac{{f}_{i+1}-2{f}_{i}+{f}_{i-1}}{h^{2}}\]

\[\cfrac{\partial u}{\partial t}=\alpha\cfrac{\partial ^{2} u}{\partial x^{2}}\Rightarrow \cfrac{\partial ^{2} u}{\partial x^{2}}=\cfrac{1}{\alpha}\cfrac{\partial u}{\partial t}\]

We usually use the order of accuracy as an indicator of the ability of finite difference schemes to approximate the exact solution as it tells us how the discretization error will reduce with mesh refinement. Another way of measuring the order of accuracy is the modified wavenumber approach. In this approach, we see how much the modified wave number is different from the true wave number. The solution of a nonlinear partial differential equation usually contains several frequencies and the modified wavenumber approach provides a way to assess how well the different components of the solution are represented. The modified wavenumber varies for every finite difference scheme and can be found using Fourier analysis of the differencing scheme.

Compact finite difference schemes have very good resolution characteristics and can be used for capturing high-frequency waves. In compact formulation, we express the derivative of a function as a linear combination of values of function defined on a set of nodes. The compact method tries to mimic global dependence and hence has good resolution characteristics. Lele investigated the behavior of compact finite difference schemes for approximating a first and second derivative. The fourth-order accurate approximation to the second derivative is given by

\[\left.\frac{1}{12} \frac{\partial^{2} u}{\partial x^{2}}\right|_{i-1}+\left.\frac{10}{12} \frac{\partial^{2} u}{\partial x^{2}}\right|_{i}+\left.\frac{1}{12} \frac{\partial^{2} u}{\partial x^{2}}\right|_{i+1}=\frac{u_{i-1}-2 u_{i}+u_{i+1}}{\Delta x^{2}}\]

Taking the linear combination of heat equation and using the same coefficients as the above equation, we get

\[\left.\frac{1}{12} \frac{\partial u}{\partial t}\right|_{i-1}+\left.\frac{10}{12} \frac{\partial u}{\partial t}\right|_{i}+\left.\frac{1}{12} \frac{\partial u}{\partial t}\right|_{i+1}=\alpha\frac{u_{i-1}-2 u_{i}+u_{i+1}}{\Delta x^{2}}\]

We then use Crank-Nicolson numerical scheme for time discretization and we arrive to below equation

\[\begin{split}\begin{align} &\frac{1}{12} \frac{u^{n+1}_{i-1}-u^{n}_{i-1}}{\Delta t}+\frac{10}{12} \frac{u^{n+1}_{i}-u^{n}_{i}}{\Delta t}+\frac{1}{12} \frac{u^{n+1}_{i+1}-u^{n}_{i+1}}{\Delta t}\\ = &\frac{\alpha}{2\Delta x^{2}} (u^{n+1}_{i-1}-2 u^{n+1}_{i}+u^{n+1}_{i+1}+u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1}) \end{align}\end{split}\]

Simplifying above equation, we get the implicit compact Pade (ICP) scheme for heat equation

\[a_{i}u^{n+1}_{i-1}+b_{i}u^{n+1}_{i}+c_{i}u^{n+1}_{i+1}=d_{i}\]

where

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{12}{\Delta x^{2}}-\cfrac{2}{\alpha\Delta t}\\ b_{i}=\cfrac{-24}{\Delta x^{2}}-\cfrac{20}{\alpha\Delta t}\\ c_{i}=\cfrac{12}{\Delta x^{2}}-\cfrac{2}{\alpha\Delta t}\\ d_{i}=\cfrac{-2}{\alpha\Delta t}(u^{n+1}_{i-1}+10 u^{n+1}_{i}+u^{n+1}_{i+1}) -\cfrac{12}{\Delta x^{2}}(u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1})\\ \end{array}\end{split}\]

proof:

\[\begin{split}\begin{align} &\frac{1}{12\Delta t}u^{n+1}_{i-1}- \frac{1}{12\Delta t}u^{n}_{i-1}+\frac{10}{12\Delta t}u^{n+1}_{i} -\frac{10}{12\Delta t}u^{n}_{i}+\frac{1}{12\Delta t}u^{n+1}_{i+1}-\frac{1}{12\Delta t}u^{n}_{i+1}\\ = &\frac{\alpha}{2\Delta x^{2}} (u^{n+1}_{i-1}-2 u^{n+1}_{i}+u^{n+1}_{i+1}+u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1}) \end{align}\end{split}\]

\[\begin{split}\begin{align} &(\frac{1}{12\Delta t}-\frac{\alpha}{2\Delta x^{2}})u^{n+1}_{i-1}- \frac{1}{12\Delta t}u^{n}_{i-1}\\ +&(\frac{10}{12\Delta t}+\frac{2\alpha}{2\Delta x^{2}})u^{n+1}_{i} -\frac{10}{12\Delta t}u^{n}_{i}\\ +&(\frac{1}{12\Delta t}-\frac{\alpha}{2\Delta x^{2}})u^{n+1}_{i+1}-\frac{1}{12\Delta t}u^{n}_{i+1}\\ =&\frac{\alpha}{2\Delta x^{2}} (u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1}) \end{align}\end{split}\]

\[\begin{split}\begin{align} &(\frac{1}{12\Delta t}-\frac{\alpha}{2\Delta x^{2}})u^{n+1}_{i-1} +(\frac{10}{12\Delta t}+\frac{2\alpha}{2\Delta x^{2}})u^{n+1}_{i} +(\frac{1}{12\Delta t}-\frac{\alpha}{2\Delta x^{2}})u^{n+1}_{i+1}\\ =&\frac{\alpha}{2\Delta x^{2}} (u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1}) + \frac{1}{12\Delta t}u^{n}_{i-1} + \frac{10}{12\Delta t}u^{n}_{i} + \frac{1}{12\Delta t}u^{n}_{i+1} \end{align}\end{split}\]

\[\begin{split}\begin{align} &(\frac{1}{12}-\frac{\alpha\Delta t}{2\Delta x^{2}})u^{n+1}_{i-1} +(\frac{10}{12}+\frac{2\alpha\Delta t}{2\Delta x^{2}})u^{n+1}_{i} +(\frac{1}{12}-\frac{\alpha\Delta t}{2\Delta x^{2}})u^{n+1}_{i+1}\\ =&\frac{\alpha\Delta t}{2\Delta x^{2}} (u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1}) + \frac{1}{12}u^{n}_{i-1} + \frac{10}{12}u^{n}_{i} + \frac{1}{12}u^{n}_{i+1} \end{align}\end{split}\]

Simplifying above equation, we get the implicit compact Pade (ICP) scheme for heat equation

\[a_{i}u^{n+1}_{i-1}+b_{i}u^{n+1}_{i}+c_{i}u^{n+1}_{i+1}=d_{i}\]

where

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{1}{12}-\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ b_{i}=\cfrac{10}{12}+2\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ c_{i}=\cfrac{1}{12}-\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ d_{i}=\cfrac{\alpha\Delta t}{2\Delta x^{2}} (u^{n}_{i-1}-2 u^{n}_{i}+u^{n}_{i+1}) + \cfrac{1}{12}u^{n}_{i-1} + \cfrac{10}{12}u^{n}_{i} + \cfrac{1}{12}u^{n}_{i+1}\\ \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{1}{12}-\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ b_{i}=\cfrac{10}{12}+2\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ c_{i}=\cfrac{1}{12}-\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ d_{i}=(\cfrac{1}{12}+\cfrac{\alpha\Delta t}{2\Delta x^{2}})u^{n}_{i-1}+(\cfrac{10}{12}-2\cfrac{\alpha\Delta t}{2\Delta x^{2}})u^{n}_{i}+(\cfrac{1}{12}+\cfrac{\alpha\Delta t}{2\Delta x^{2}}) u^{n}_{i+1} \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{1}{12}-\cfrac{\alpha\Delta t}{2\Delta x^{2}}\quad \hat{a}_{i}=\cfrac{1}{12}+\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ b_{i}=\cfrac{10}{12}+2\cfrac{\alpha\Delta t}{2\Delta x^{2}}\quad \hat{b}_{i}=\cfrac{10}{12}-2\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ c_{i}=\cfrac{1}{12}-\cfrac{\alpha\Delta t}{2\Delta x^{2}}\quad \hat{c}_{i}=\cfrac{1}{12}+\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ d_{i}=\hat{a}_{i}u^{n}_{i-1}+ \hat{b}_{i}u^{n}_{i}+\hat{c}_{i}u^{n}_{i+1} \end{array}\end{split}\]

version 2

\[\begin{split}\begin{array}{l} a_{i}=12\cfrac{1}{12}-12\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ b_{i}=12\cfrac{10}{12}+12\times2\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ c_{i}=12\cfrac{1}{12}-12\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ d_{i}=(12\cfrac{1}{12}+12\cfrac{\alpha\Delta t}{2\Delta x^{2}})u^{n}_{i-1}+(12\cfrac{10}{12}-12\times2\cfrac{\alpha\Delta t}{2\Delta x^{2}})u^{n}_{i}+(12\cfrac{1}{12}+12\cfrac{\alpha\Delta t}{2\Delta x^{2}}) u^{n}_{i+1} \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{1}{1}-12\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ b_{i}=\cfrac{10}{1}+12\times2\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ c_{i}=\cfrac{1}{1}-12\cfrac{\alpha\Delta t}{2\Delta x^{2}}\\ d_{i}=(\cfrac{1}{1}+12\cfrac{\alpha\Delta t}{2\Delta x^{2}})u^{n}_{i-1}+(\cfrac{10}{1}-12\times2\cfrac{\alpha\Delta t}{2\Delta x^{2}})u^{n}_{i}+(\cfrac{1}{1}+12\cfrac{\alpha\Delta t}{2\Delta x^{2}}) u^{n}_{i+1} \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{1}{\alpha\Delta t}-\cfrac{12}{2\Delta x^{2}}\\ b_{i}=\cfrac{10}{\alpha\Delta t}+24\cfrac{1}{2\Delta x^{2}}\\ c_{i}=\cfrac{1}{\alpha\Delta t}-\cfrac{12}{2\Delta x^{2}}\\ d_{i}=(\cfrac{1}{\alpha\Delta t}+\cfrac{12}{2\Delta x^{2}})u^{n}_{i-1}+(\cfrac{10}{\alpha\Delta t}-\cfrac{24}{2\Delta x^{2}})u^{n}_{i}+(\cfrac{1}{\alpha\Delta t}+\cfrac{12}{2\Delta x^{2}}) u^{n}_{i+1} \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{2}{\alpha\Delta t}-\cfrac{12}{\Delta x^{2}}\\ b_{i}=\cfrac{20}{\alpha\Delta t}+\cfrac{24}{\Delta x^{2}}\\ c_{i}=\cfrac{2}{\alpha\Delta t}-\cfrac{12}{\Delta x^{2}}\\ d_{i}=(\cfrac{2}{\alpha\Delta t}+\cfrac{12}{\Delta x^{2}})u^{n}_{i-1}+(\cfrac{20}{\alpha\Delta t}-\cfrac{24}{\Delta x^{2}})u^{n}_{i}+(\cfrac{2}{\alpha\Delta t}+\cfrac{12}{\Delta x^{2}}) u^{n}_{i+1} \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{2}{\alpha\Delta t}-\cfrac{12}{\Delta x^{2}}\\ b_{i}=\cfrac{20}{\alpha\Delta t}+\cfrac{24}{\Delta x^{2}}\\ c_{i}=\cfrac{2}{\alpha\Delta t}-\cfrac{12}{\Delta x^{2}}\\ d_{i}=\cfrac{2}{\alpha\Delta t}(u^{n}_{i-1}+10u^{n}_{i}+u^{n}_{i+1}) +\cfrac{12}{\Delta x^{2}}(u^{n}_{i-1}-2u^{n}_{i}+u^{n}_{i+1})\\ \end{array}\end{split}\]

\[\begin{split}\begin{array}{l} a_{i}=\cfrac{12}{\Delta x^{2}}-\cfrac{2}{\alpha\Delta t}\\ b_{i}=\cfrac{-24}{\Delta x^{2}}+\cfrac{-20}{\alpha\Delta t}\\ c_{i}=\cfrac{12}{\Delta x^{2}}-\cfrac{2}{\alpha\Delta t}\\ d_{i}=\cfrac{-2}{\alpha\Delta t}(u^{n}_{i-1}+10u^{n}_{i}+u^{n}_{i+1}) +\cfrac{-12}{\Delta x^{2}}(u^{n}_{i-1}-2u^{n}_{i}+u^{n}_{i+1})\\ \end{array}\end{split}\]