Finite Difference Method for Reaction-Diffusion Problems

The finite difference method is a useful method for problems of diffusion and reaction, especially when there are steep changes in the solution in a small region of space.

Consider first reaction and diffusion when the reaction rate is linear in concentration, the Biot number for mass is large, and the temperature is constant. The problem is

At the origin, the differential equation is derived using l'Hopital's rule (link)

The finite difference formulation at interior grid points is

while the boundary condition requires

A false boundary is used at the origin, so that the first equation is

Typically, the Newton-Raphson method is used to solve nonlinear problems. The nonlinear reaction rate term is expended to give

When the solution is expressed in the form

the terms are

After the solution for c is found, the effectiveness factor is obtained using the quadrature formula.

Simpson's rule is used (for an odd number of grid points), because the trapezoid rule does not give exact results when the reaction rate is constant (it doesn't integrate r² exactly). The code is avaiable. The solution is plotted in the figure; the values at the grid points are shown, and the interpolation is linear in between. For small f the concentration is relatively constant across

the catalyst pellet, whereas for larger f the concentration dips to small values.

A solution with 9 grid points works well for the smaller Thiele modulus but gives distinct broken lines for the larger Thiele modulus. Convergence of the concentration profiles with increasing number of grid points is shown in the next figure, and 17 points is sufficient for a good solution

The finite difference method is very effective for solutions with steep fronts, although many grid points may be necessary.

The effectiveness factor is plotted, too, and this figure also shows the convergence with increasing number of grid points.

Next consider nonisothermal problems. The first problem is

with b = 0.3, g = 18, and f = 0.5. In this case, the derivative of the rate expression is

which varies from -2.8 to 16 as c varies from zero to one. Thus, the problem is not highly nonlinear, despite the variation of temperature. The solution is shown in the figure, and only a few grid points are necessary for numerical convergence to the solution of the ordinary differential equation. Obvious errors occur when ngrid is 3 or 5, but 9 seems large enough for this problem. The Newton-Raphson method was used to solve the set of nonlinear equations.

Next consider a problem whose solution has steep gradients by taking b = 0.4 and g = 30. Now it is difficult to obtain convergence in the finite difference method for large values of Thiele modulus. When the reaction is almost complete throughout the catalyst, the initial guess is taken as c(r) = 0 to enhance convergence. In spite of that, however, solutions are not always obtained with the Newton Raphson method. For this problem, and harder problems, an initial value method is better.

Take Home Message: The finite difference method is ideally suited to problems whose solution has steep gradients. Even so, there are problems that are difficult to solve without a good initial guess.