Description

Homework Assignment #4

Use Matlab to solve these problems. Generate script which is organized and commented. Please make one m file for the entire assignment. Turn in the m files AND pdf files of the output (using the publish function).

Boundary Value Problems

1) For the second order differential equation u u” 1, solve using finite differences using at least five internal nodes and plot the approximate solution as w xi ( )i for each of the following boundary conditions:
a) u(0)  0;u(1)  1
b) u(0)  0;u(1) u ‘(1)  1
c) u(0)  0;u ‘(1)  1

2) The equation below is the governing equation for heat conduction in an annulus with a uniform heat generation term. Use the method of finite differences to solve for ( ) for three different dimensionless heat generate rates, = [7.5, 8, 8.5] over the range 1with  = 0.8. Plot all three profiles on one plot.

d2 2  1 dd 2 1 with d(1)  0 d  
 0
d ( )

3) The expression below describes the steady state concentration, C, of a reactant inside a spherical catalyst pellet as a function of the radial position, r, which ranges from the center of the pellet (r = 0) to the outer edge of the pellet (r = R). The profile is a function of the diffusion coefficient, D, and the rate constant, k.

d C2 2 dC k
2  r dr  DC  0 dr

Rewrite this expression in terms of a dimensionless radius,  = r/R and dimensionless concentration,  = C/CR, where CR is the concentration at the outer edge of the pellet. One boundary condition is that the concentration at the outer edge is a constant (C (r=R) = CR) and the other is that the gradient at the center of the pellet must be zero. Once written in dimensionless form, there will be a Thiele modulus:
R2 / D
2
 
1/ k
Solve for the concentration profile in the pellet for the two cases of 2 = 0.01 and 10. Plot the two concentration profiles in a way that can be compared.

a) Use the shooting method first and comment on choices used for guessing and ensuring the correct solution
b) Use the method of finite differences and comment on the impact of step size on the solution.
c) Using your results from (a) or (b), determine the gradient of concentration at the outer edge of the pellet.

4) The governing equation for steady state one-dimensional transport of a species with both reaction and flow is given below. The concentration of the species, C(z), depends on the diffusivity, D, local fluid velocity, v, and reaction rate constant, k. This is valid over the range 0 z L.
Dd Cdz2 2 vdCdz kCn 0

For the case that the reaction is first order (n=1), plot the concentration profile C(z*) where z* = z/L using the method of finite differences. Use Dirichlet boundary conditions with the concentration at 1 M at the entrance and one tenth of that at the other boundary.

D 10 m2/sec
L 1 mm
 0.1 m/sec
k 5 x 10-3 sec-1

5) The governing equation for steady state one-dimensional transport of a species with both reaction and flow is given below. The concentration of the species, C(z), depends on the diffusivity, D, local fluid velocity, v, and reaction rate constant, k. This is valid over the range 0 z L.
Dd Cdz2 2  dCdz n 0 v kC 

For the case that the reaction is second order (n=2), plot the concentration profile C(z*) where z* = z/L using the method of finite differences. Use Dirichlet boundary conditions with the concentration at 1 M at the entrance and one tenth of that at the other boundary.

D 10 m2/sec
L 1 mm
 0.1 m/sec
k 5 x 10-5 M-1sec-1

6) One component of a model for a styrene monomer tubular reactor is the steady-state temperature profile of the solid phase catalyst. The governing boundary value problem is:

ddx2 2  1x dxd 2  0

d
 0, (1) 1.
dx x0
Here,  (T Tc )/(Tw T) is the nondimensional catalyst temperature, x is the nondimensional radial position, T is the temperature of the fluid in the reactor, Tw is the temperature at the wall of the reactor. The parameter 2 is given by

R2hA
2
  k(1)

Where R = 1.3cm is the radius of the reactor, h = 0.001 cal/ cm2 s oC is the heat transfer coefficient, =0.36 is the porosity of the packed bed reactor, A= 15 cm-1 is the surface area of the catalyst per unit volume and k = 0.0034 cal/cm s oC is the thermal diffusivity of the catalyst. Approximate (x) using x = 0.0025 using finite differences and plot (x).

7) The following BVP comes about when considering exothermic reactions in porous media,* where the steady state temperature profile (z) adheres to:

ddz2 2 1exp    0
 B
2  B    

(0)  0, (1)  0

The parameter B is the maximum possible temperature in the absence of natural convection, 2 is the ratio of the characteristic time for conduction to that for heat generation and g is the dimensionless activation energy. For B = 0.6, 2 = 0.25 and = 30.0, determine and plot (z) using finite differences.

*Reference if interested: Subramanian and Balakotaiah, Physics of Fluids 6(9): 2907-2922, 1994

Optimization – Single Variable
8) Being able to fit nonlinear models to data is important in many fields of engineering. The approach for nonlinear functions is to minimize the a cost function defined as:
Ndata
Fc ( )  [Spredict (tk ;) Sobs ( )]tk 2
k1
where the cost function is the sum of squared differences between the measured data (obs) and the predicted function of the fit (predict). As an example, the interfacial tension, EQ , which is a function of the bulk concentration of surfactant, Cbulk is described by a Langmuir isotherm:
  EQ  EQ Cbulk
 EQ  o  RT ln 1   where 
    Cbulk a

which really has two unknowns, a, and  . This is clearly a nonlinear function, so fitting the expression to data is best done using optimization. The table of data below gives the measured values of  (Cbulk). The interfacial tension of a clean interface (Cbulk = 0) is 52.2 mN/m and the molecular weight of the surfactant is 627 g/mole.
Cbulk (g/L) 5 x10-5 10-4 4×10-4 5×10-4 10-3 0.002 0.003
EQ (mN/m) 36.42 33.72 30.63 27.45 24.76 22.30 19.71

a) Write code that takes the input data and initial guesses for the parameter vector, , and determines a best fit.
b) Use your code to see if the results for  depend on the initial guess. How confident are you that you have the “right” answer?

9) Determine the length of the shortest ladder that can be made to lean against the side of a building without interfering with an adjoining room of height 10 ft and width 10 ft shown in the figure below.

State this problem as an optimization problem and solve this using an optimization routine in MatLAB .