File last modified 29 October 1998
In this section we will explore some of the techniques that have been developed to handle parabolic PDEs. We won't, by any stretch of the imagination, be able to cover all of the techniques that have been applied (there are a lot of 'em), but these techniques have shown their usefulness in many applications and so we will present them to you here.
If we have a "one-dimensional" PDE of the form:
and if D is constant with respect to x, then:
and we can rewrite it in finite difference form as:
We've seen this method (and others such as the FUDM and SUDM) in our finite difference sections. And, as we did earlier in our finite difference lecture (section 11.4.2) we can define Fourier components that with a little algebra and "well known" trigonometric/exponential identities will lead us to the amplification factor and stability requirements.
Let
We can interpret this stability analysis as also defining a time scale governing the propagation (or diffusion) of information across the model grid. So we can define this diffusion time scale as:
Where 
is the corresponding spatial scale.
Normally we are interested in phenomena with 
, if we are limited to
that satisfy:
then we will be forced to take on the order of 
time steps before anything starts to
happen. This can be a prohibitively large number of time steps (computationally or just in terms
of our own impatience).
There is hope. One solution is to force all of those small-scale processes, that we are not
interested in and for which our large 
would be woefully inadequate (and inaccurate), to their
equilibrium state. We do this with a fully implicit form of the finite difference equations:
You will note that all of the dependent variables on the right hand side are at the n+1 time step. Regrouping and rearranging gives us an algorithm of the following form:
This may still look hopeless. There are, after all, a lot of un+1 terms on the lefthand side of the equation. But it turns out that this algorithm can be solved, with the appropriate boundary conditions, by the tridiagonal algorithm. This solution essentially uses the boundary conditions to solve backwards in time, a pretty cute trick. There's a good write-up of this algorithm in Appendix A of Roache (1976) and you will see this algorithm in Problem Set #10.
But what about stability? A stability analysis reveals that the coefficient of the Fourier components is:
which is unconditionally stable for any choice of 
, but it is only first order accurate in 
.
There is another possibility and this one is accurate to 2nd order in time. Allow those small scales to fluctuate in a sort of "background" mode of the scales we are interested in. To put it another way, average the fully explicit and fully implicit methods to form a type of hybrid. The finite difference algorithm looks like:
If you think about it you'll see that both the lefthand and right hand side of this equation is
centered on time step 
and hence it is 2nd order in time.
What about its stability? Again following the same stability analysis as before we find that the amplitude of the Fourier components comes out to be:
which is also unconditionally stable for any 
.
The algorithm laid out in equation 14.3.10 is known as the Crank-Nicholson scheme and is probably your best bet starting point for diffusive problems. You'll still need to solve the n+1 parts with the tridiagonal algorithm, but for most parabolic PDE problems involving diffusion only this is a good algorithm.
Now, of course, not all of our problems are one-dimensional (in fact most are not). So what would a two-dimensional Crank-Nicholson scheme look like? Consider the following two-dimensional form of the diffusion equation:
assuming that the diffusion coefficient is isotropic, here's what the 2-D Crank-Nicholson looks like:
and here I've used the j subscript and the l subscript for the spatial indexing and
and
and the same goes for the n+1 terms. It's messy, but doable.
There are other semi-implicit methods, such as time splitting methods, but they go beyond the scope of this course. If you are interested, I recommend Roach (1976, Computational Fluid Dynamics) as well as Press et al. (1990, Numerical Recipes).
GoTo Next Section