Orthogonal Collocation on Finite Elements for Reaction-Diffusion Problems

While the orthogonal collocation method is a useful method for problems of diffusion and reaction, it has difficulty when the concentration profiles become steep. In those cases, the method of orthogonal collocation on finite elements is a better choice. In reaction-diffusion problems, the steepness of the profile is related to the Thiele modulus, f. For large f the solution is always steep. The method of orthogonal collocation on finite elements is applied here to a difficult nonlinear problem.

The first problem is

with b = 0.4, g = 30, and various f. In this case, multiple solutions are possible, that is, for the same parameters, there is more than one solution to the problem. This can be seen in Figure 1 of effectiveness factor versus Thiele modulus. For 0.21 < f < 0.56 there are three solutions. For f = 0.5, two of the solutions are shown in Figure 2.

Figure 1. Effectiveness Factor for b = 0.4, g = 30

Figure 2a. Two concentration solutions for b = 0.4, g = 30,f = 0.5

Figure 2b. The corresponding temperature solutions for b = 0.4, g = 30,f = 0.5

To apply orthogonal collocation on finite elements we divide the domain, 0 r 1, into finite elements, as shown in Figure 2. The elements are denoted by the superscript e (going from 1 to NE) and the points within an element are denoted by capital I or J (going from 1 to NP). It is convenient to expand the equation to the following form.

Figure 3. Collocation Points for OCFE

Within the e-th element, we make the transformation

Then the differential equation becomes

in each element, and the boundary conditions become

The orthogonal collocation method is applied at each collocation point interior to the e-th element.

In addition, continuity of the function and first derivative between elements requires the condition

The function is continuous by virtue of assigning it only one value, regardless of whether it is the last point on one element or the first point on the next element. The continuity of first derivative is enforced using

at the points between elements. The boundary conditions are

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

The problem is solved using the general form

In this form, we need

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

The code is available.

The solutions are shown in Figure 2. The third solution is difficult to find, since it is unstable. These two solutions were found by starting from an initial guess of c(r) = 1 or c(r) = 0.

For the second problem, take the case of b = 0.02, g = 30, and various f. This is a more physically realistic value of b. For such low values, though, multiple solutions are not possible, as shown in Figure xxx

Figure 4. Effectiveness Factors for b = 0.02, g = 30

The problem can be made more difficult, however, when there is heat and mass transfer resistance at the catalyst boundary. The equations are

For Bi_m = 250 and Bi = 5 multiple solutions are possible, as shown in Figure xxx. This problem is extremely difficult to solve, and only a powerful method can solve it.

Take Home Message: The method of orthogonal collocation on finite elements is ideally suited to problems whose solutions have extremely steep gradients.