Basic discretization techniques

Could maybe do an engineering approach and a mathematical approach... Two ways of presenting the same material

Introduction: discretization of derivatives

What are we trying to do?

Before discussing discretization, we must first discuss what we're discretizing

We start with some differential equation - either partial differential equation or ordinary differential equation

Let's say we have a function $ u $ of space and time $ x,t $

What is the general form of a partial differential equation for $ u $? (Introduction_to_partial_differential_equations)

$ F(x,t,u_x,u_xx,\dots,u_t,u_tt,\dots) = 0 $

Revisit example of heat equation, solved analytically using separation of variables (Analytical_solution_of_PDEs#Example_3:_Heat_Equation):

$ \Theta_t = \alpha^2 \Theta_{xx} $

How did we treat $ \Theta $? How did we solve the PDE? We assumed a form for the solution, plugged it in, and found out some additional information

For example, we assumed a solution of the form $ X(x) T(t) $, and then we found, for time:

$ T = c_1 \exp{( -\lambda^2 \alpha^2 t )} $

and something similar for $ X(x) $.

So the answer is: we solved the PDE by finding functions that satisfied certain conditions...

Solution approach:

• We assumed a form for the solution
  • Found two ordinary differential equations that the solution had to satisfy
• We then had a general solution (that worked for any boundary and initial conditions)
  • Found values for constants by plugging our boundary conditions into the general solution
• When then had a (less general) solution that worked for our boundary conditions, and worked for any initial conditions
  • Found values for the remaining constants by plugging our initial conditions into the solution
• Final result: we found a solution that worked for our boundary conditions and our initial conditions

Can we find analytical solutions using this approach with a computer?

Only when the partial differential equation is really simple... Only when an analytical solution actually exists...

But whether we can find analytical solutions is beside the point.

We want a different method, a more robust method, than trying every analytical solution technique

Functions are continuous. What is the opposite of continuous? Discrete - if we don't want an analytical (functional) solution, we need a discrete solution

Discrete - means, the solution only exists at discrete points, and is not an exact solution, but only an approximation

1st derivative approximations

(A note on notation: while the following definitions/formulae hold for either partial or full derivatives, the derivative will be written as partial derivatives to keep the presentation general)

Limit definition of a derivative:

$ \frac{\partial u(x)}{ \partial x} = \displaystyle{ \lim_{\Delta x \rightarrow 0} } \frac{ u(x_0 + \Delta x) - u(x_0) }{ \Delta x } $

This is saying, when $ \Delta x $ is "small enough", the algebraic difference $ \frac{ u(x_0 + \Delta x) - u(x_0) }{ \Delta x} $ is a good approximation for the derivative; if $ \Delta x $ is small enough, the representation becomes exact.

OK, so naturally we want to know: how small? How small does $ \Delta x $ have to be?

Using a Taylor series expansion for $ u(x_0 + \Delta x) $ in terms of $ u(x_0) $:

$ u(x_0 + \Delta x) = u(x_0) + \left. \frac{\partial u}{\partial x} \right|_{0} \Delta x + \left. \frac{\partial^2 u}{\partial x^2} \right|_{0} \frac{(\Delta x)^2}{2!} + \dots + \left. \frac{\partial^n u}{\partial x^n} \right|_{\xi} \frac{ (\Delta x)^n }{ n! } \qquad x_0 \leq \xi \leq x_0 + \Delta x $

What does the Taylor series tell us about how small $ \Delta x $ needs to be?

Rearrange the Taylor series to look like the limit definition of the derivative...

$ \left. \frac{\partial u}{\partial x} \right|_{x_0} = \frac{ u(x_0 + \Delta x) - u(x) }{ \Delta x } - \left. \frac{\partial^2 u}{\partial x^2} \right|_0 \frac{\Delta x}{2!} - \dots $

For a discrete representation, values of $ u $ can be indexed at each discrete point using some index $ i $, so that $ u_i = u(x_0) $, $ u_{i+1} = u(x_0 + \Delta x) $, etc...

$ \left. \frac{\partial u}{\partial x} \right|_{i} = \frac{ u_{i+1} - u_{i} }{ \Delta x } - \left. \frac{\partial^2 u}{\partial x^2} \right|_0 \frac{\Delta x}{2!} - \dots $

This is a forward difference.

Let's look at this last term... This is an "error" in our difference approximation. So we can replace it with:

$ - \left. \frac{ \partial^2 u}{\partial x^2} \right|_0 \frac{\Delta x}{2!} = O( \Delta ) $

where $ O() $ is the order of the error.

Tells us two things:

1. For the difference approximation of the derivative to be accurate, we have to have $ \Delta x $ become really small to make our error really small.

2. We ALSO need $ \frac{\partial^2 u}{\partial x^2} $ to be really small.

That means if we're trying to approximate a function with a large second derivative, we need an even smaller $ \Delta x $.

Applying this concept: the red line is a very high-order polynomial, and the blue line is a low-order polynomial.

Which one has a larger second derivative?

Which one will require a smaller $ \Delta x $?

We didn't get an exact answer to the question "how big should $ \Delta x $ be?", but we've got an idea now

Next, we can write another Taylor series, this time for the point $ u(x_0 - \Delta x) $:

$ u(x_0 - \Delta x) = u(x_0) - \left. \frac{ \partial u }{ \partial x } \right|_0 \Delta x + \left. \frac{ \partial^2 u}{ \partial x^2} \right|_0 \frac{ (\Delta x)^2 }{ 2! } - \left. \frac{ \partial^3 u}{ \partial x^3} \right|_0 \frac{ (\Delta x)^3 }{ 3! } + \dots $

This is a backward difference.

Rearranging it the same way we did the forward difference,

$ \left. \frac{ \partial u}{\partial x} \right|_{i} = - \frac{ u_i - u_{i-1} }{ \Delta x } + \frac{ \partial^2 u }{ \partial x^2 } \frac{ \Delta x }{2!} + \dots $

So this one has an extra error term, that is also $ O(\Delta x) $.

Can we represent the derivative some other way?

Yes!

In fact, we can represent the derivative a whole bunch of different ways.

We've shown approximations of the derivative of $ u $ at the point $ x $ using the value of $ u(x-\Delta x) $ and $ u(x + \Delta x) $ - but what about using both?

How to approximate the derivative $ u(x) $ at $ x $ using both the value of $ u(x-\Delta x) $ and $ u(x+\Delta x) $?

Subtracting the backward difference from the forward difference yields a central difference:

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

The order term is different, because the first-order error terms canceled out. What does that tell us?

1. To approximate a given function, the $ \Delta x $ for a central difference scheme can be larger than the $ \Delta x $ for a forward or backward difference scheme

2. If a function has large gradients (and thus a larger error term), a central difference scheme will approximate the derivative better than a forward or backward difference scheme for the same $ \Delta x $

2nd derivative approximations

As expected, one can also discretize second-order derivatives. Adding the backward and forward differencing schemes and rearranging yields a second-order derivative central difference scheme:

$ \displaystyle{ \left. \frac{\partial^2 u}{\partial x^2} \right|_{i} = \frac{ u_{i+1} - 2 u_{i} + u_{i-1} }{ (\Delta x)^2 } } + O( (\Delta x)^2 ) $

Like the first derivative, the second derivative may be approximated in a number of different ways. Alternative discretization methods (analogous forward and backward differences) can be derived by expanding various points neighboring $ x_i $ using Taylor series, and combining these Taylor series approximations to yield discrete expressions for second derivatives.

For example, the second derivative forward differencing scheme is:

$ \left. \frac{ \partial^2 u }{ \partial x^2 } \right|_i = \frac{ u_i - 2 u_{i+1} + u_{i+2} }{ (\Delta x)^2 } + O[ (\Delta x) ] $

and the second derivative backward differencing scheme is:

$ \left. \frac{ \partial^2 u }{ \partial x^2 } \right|_i = \frac{ u_i - 2 u_{i-1} + u_{i-2} }{ (\Delta x)^2 } + O[ (\Delta x) ] $

Alternative approaches

Integral Form

So far, we've been approximating derivatives using differences (via Taylor series).

An alternative approach, one that's more physically intuitive, is to approximate the integral of the derivative.

For example, if we are looking at a property of a fluid contained in a pipe, one way to approximate the change in the property is to use a macro balance:

$ \mbox{in} - \mbox{out} + \mbox{generated} - \mbox{consumed} = \mbox{accumulated} $

This can be applied to the volume of the entire pipe $ \Omega $ for a property $ \phi $.

First, each term can be represented using an integral quantity:

$ \begin{array}{rcl} \mbox{in} &=& \int_{S_{in}} \phi u_{\phi,in} d A_{in} \\ \mbox{out} &=& \int_{S_{out}} \phi u_{\phi,out} d A_{out} \\ \mbox{generated} &=& \int_{\Omega} G_{\phi} dV \\ \mbox{consumed} &=& \int_{\Omega} C_{\phi} dV \\ && \\ \mbox{accumulated} &=& \int_{\Omega} \frac{ \partial \phi }{ \partial t } dV \end{array} $

These can then be assembled into the integral form of the differential equation for $ \phi $:

$ \int_{\Omega} \frac{ \partial \phi }{ \partial t } dV = \int_{S_{in}} \phi u_{\phi,in} d A_{in} - \int_{S_{out}} \phi u_{\phi,out} d A_{out} + \int_{\Omega} G_{\phi} dV - \int_{\Omega} C_{\phi} dV $

When this approach is then applied to arbitrarily small control volumes in the fluid, and various assumptions made that simplify the above integral terms, this yields a large number of discrete equations, analogous to the set of algebraic equations that result in applying the finite difference method to a partial differential equation. This method is known as the Finite volume method.

Error Analysis

Tannehill Fletcher p. 52

References

Tannehill, John C.; Anderson, Dale A.; Pletcher, Richard H. (1997). Computational Fluid Mechanics and Heat Transfer (2nd ed.). Taylor and Francis. ISBN 978-1560320463. 

Anderson, John D. (1995). Computational Fluid Dynamics: The Basics with Applications. McGraw Hill. ISBN 0-07-001685-2.