Page 323 - Numerical Methods for Chemical Engineering

312 6 Boundary value problems

can solve systems of equations of elliptic, hyperbolic, and parabolic type using FEM. The
software package FEMLAB TM , built upon MATLAB by the developers of the PDE toolkit
(www.comsol.com) can solve more general and complex BVPs in two and three dimensions.
Three dimensional BVPs do not pose any new conceptual issues, but elimination can no
longer be used to solve the resulting linear systems. Iterative methods are necessary, such
as pcg if the system is positive-definite and bicgstab or gmres if it is not. Preconditioners
improve significantly the efficiency of these methods, and one may either use cholinc or
luinc to perform an incomplete Cholesky or LU factorization respectively.

Problems

6.A.1. Use finite differences to discretize the following BVP in three dimensions:
2
2
2
2
−∇ ϕ = exp[− (x + y + z )/2] − 1 ≤ x, y, z ≤ 1
(6.231)
Dirichlet condition ϕ = 0 on all boundaries
Solve the linear system with pcg, for a grid of 50×50×50 points. How many iterations
are necessary with no preconditioner? Next, use an incomplete Cholesky preconditioner
with no fill-in, and report the number of iterations required for convergence. Finally, as a
function of droptol, plot the number of iterations and the number of nonzero elements in the
Cholesky factor. What value of drop tolerance do you recommend using? For this optimal
value of the drop tolerance, change the number of grid points, and report how the number
of iterations and CPU time (doc cputime) varies with the size of the grid.

6.A.2. Solve the system Ax = b that discretizes the BVP of problem 6.A.1 without storing
the matrix A in memory. Supply a routine that returns Av for input v.
6.A.3. Solve the 2-D version of problem 6.A.1 with z = 0 by the FEM using the PDE toolkit
or FEMLAB TM .
6.B.1. Consider the following heat transfer problem. A fluid with thermal conductivity λ,
ˆ
density ρ, and specific heat capacity C p flows at a Reynolds number Re < 100 through a
cylindrical pipe of radius R. The fluid viscosity is µ, and you may neglect viscous heating. If
z is the axial position, for z < 0 the wall temperature is T 0 and at z = 0 the wall temperature
jumps abruptly to T 1 for z > 0. This is known as the Graetz problem, and an analytical
solution exists if conduction in the axial direction is neglected. Transform this problem
into dimensionless variables and write a program to compute numerically the steady-state
temperature profile without neglecting axial conduction. Have your program take as input
the values of the dimensionless parameters that characterize the system. Report your results
2
−2
forvaluesoftheparametersintherange[10 ,10 ].Tryreducingthenumberofindependent
parameters through clever rescaling.
6.B.2. You are conducting the enzymatic conversion of a substrate S into a product P, with
the micromoles of S converted per minute per milligram of enzyme being described by
Michaelis Menten kinetics,
V m S µmol
−ˆ r S = V m = 200 K m = 0.2 M (6.232)
K m + S min mg E