Page 320 - Numerical Methods for Chemical Engineering

FEM in MATLAB 309

[d]
Therefore for each node p /∈ ∂ , we have an algebraic equation

0 = A pq ϕ q − b p (ϕ) (6.228)
q /∈∂ [d]
where

[q] [q] [q] [q]
b p (ϕ) = S pq f r ,ϕ q + S pq f r ,ϕ B r − A pq ϕ B r
q /∈∂ [d] q∈∂ [d] q∈∂ [d]
(6.229)
For a source term with no field-dependence, (6.228) is a sparse linear system. If not, it must
be solved with Newton’s method; however, the Jacobian matrix is also sparse.

Convection terms in FEM

Above, we have considered only diffusive transport; however, it is not difficult to extend
the method to include convective transport as well. The choice of linear shape functions
used above results in a system of algebraic equations with similar accuracy to central
finite difference approximations and shares the tendency of CDS to exhibit unphysical
oscillations when convection is strong. The oscillations are reduced by mesh refinement (to
reduce the local Peclet number) or by adding additional numerical diffusion, especially in the
streamline direction. Upwind versions of FEM are obtained by biasing the shape functions
to weight more heavily the upstream direction. Such subjects are beyond the scope of
this chapter, and the reader is referred to a dedicated text on FEM such as Akin (1994).
The relationship between FEM with linear shape functions and CDS finite differences is
examined in the supplemental material in the accompanying website for the case of the 1-D
convection/diffusion equation.

FEM in MATLAB

In MATLAB, an optional pdetool toolbox allows one to solve 2-D BVPs and to perform var-
ious low-level mesh generation and assembly operations. There exists also a more advanced
software package, FEMLAB TM ,(www.comsol.com), from the developers of pdetool, that
can solve multiple PDEs of various types in two and three dimensions. The supplemen-
tal material in the accompanying website contains an example of the use of FEMLAB TM
to model a microfluidic H-filer, requiring the simultaneous solution of the Navier–Stokes
equations and a mass balance. A second FEMLAB TM example in the supplemental mate-
rial models natural convection between two flat, vertical plates at different temperatures
(Figure 6.30).
In general, if you have a simple geometry and are solving a system of transport-type
PDEs of the form of (6.9), it is not too difficult to write your own finite difference program.
For more complicated systems, it is not worth writing an FEM program from scratch, as
the FEM is more tedious to implement than finite differences, and canned software libraries
and packages already exist that are suitable for these problems (unless your system has
PDEs that are not of standard type, e.g. when solving a fluid mechanics problem with a