October 2, 2002 1Fall 2002. 10.34. Numerical Methods Applied to Chemical EngineeringExam I. 10/4/2002There are two questions to this exam. You have the full hour to work on the exam, and may use all class material (notes, books, etc.). Submit your results for each question (1 and 2) on a separate blue book. You may keep your copy of the exam questions. Make sure that your name is on each blue book!Question 1. 60 points totalConsider the following set of three nonlinear algebraic equations,(EQ 1)1.A. (15 pts.) Calculate the analytical form of the Jacobian matrix of this system, expressed in terms of the unknowns .1.B. (10 pts.) Starting from an initial guess of , , , derive the set of linear algebraic equations that must be solved for the first iteration of Newton’s method.1.C. (20 pts.) Solve this set of equations for the full Newton update using Gaussian elimination with partial pivoting. Show all of your calculations by hand, and use the exact solution process followed by a computer (no human intuition, please!).You will not receive full credit unless you show all of your calculations as performed by hand. This is to avoid giving an unfair advantage to those with fancy calculators. You may use a calculator, however, to add or multiply simple numbers, e.g. 5 - 4*2.1.D. (15 pts.) We have seen that Newton’s method can be erratic when the initial guess is far away from the solution. If you were using a modification of Newton’s method that was more robust to the choice of initial guess, would you accept the Newton update that you have calculated above?If not, explain briefly how you would go about selecting the new estimate of the solution, .(QUESTION 2 ON NEXT PAGE)f1x1x2x3,,()4x143x1x22– x2x357–+0==f2x1x2x3,,()x245x1x3–9+0==f3x1x2x3,,()3x346x1x2x315–+0==x1x2x3,,x11=x21=x31=∆x0[]∆x0[]x1[]October 2, 2002 2Question 2. 40 points totalAs we discussed in class, eigenvalue analysis can be used to study the stability of a dynamic system. If we have a system governed by a system of differential equations(EQ 2)with a steady state at , that is,(EQ 3)then the condition that the steady state be stable is that ALL eigenvalues of the Jacobian,(EQ 4)must have real parts less than zero.(EQ 5)The full details of the derivation are repeated for your review (should not be needed to solve the problem) in the optional background section following the problem statement.Consider the case of a CSTR with the single chemical reaction,(EQ 6)For example, this model may describe the increase in cell concentration (species A) in a bioreactor containing a growth factor (species B). Since we see that the reaction is self-accelerating, we might be concerned about the stability of this reactor system.The mass balances for the reactor are(EQ 7)dx1dt--------x·1f1x1x2… xN,,,()==dx2dt--------x·2f2x1x2… xN,,,()==:dxNdt---------x·NfNx1x2… xN,,,()==xxˆ=fjxˆ1xˆ2…xˆN,,,()0= j 12…N,, ,=Jjkxˆ()xk∂∂fjxxˆ==Jwj[]λjwj[]= Re λj()0< j 12…N,, ,=AB+2A→rkcAcB=tddVcA[]QcAin()QcA– VkcAcB+=tddVcB[]QcBin()QcB– VkcAcB–=October 2, 2002 3 is the volumetric flow rate through the reactor, is the total reactor volume, and and are the inlet concentrations of A and B respectively. Defining the reactor mean residence time as , the two nonlinear equations for the steady state are(EQ 8)We can easily solve these equations by adding them to obtain the relation(EQ 9)Substitution into the balance for A yields(EQ 10)This quadratic equation is solved easily for the concentration of A.For the process conditions,(EQ 11)the concentrations of A and B in the reactor at steady state are(EQ 12)2.A. (20 pts.). Is the reactor stable under these process conditions?2.B. (15 pts.). This application of eigenvalue analysis assumes that the eigenvectors of the Jacobian are linearly independent, so that any vector may be written as a linear combina-tion of the eigenvectors. While this assumption is usually OK, we can only prove that it is valid for special cases of the Jacobian matrix and of its eigenvalues . What are the most general conditions on and that you can impose for the eigenvectors of to be linearly independent?2.C. (5 pts.). For the Jacobian matrix that you calculate in 2.A, are the eigenvectors expected to be orthogonal to each other?QVcAin()cBin()θVQ⁄=0cAin()cA– θkcAcB+=0cBin()cB– θkcAcB–=0cAin()cBin()cA– cB–+=cBcAin()cBin()cA–+=0 cAin()cA– θkcAcAin()cBin()cA–+[]+=0 cAin()cA– θkcAcAin()cBin()+[]θkcA2–+=θkcA21 θkcAin()cBin()+[]–[]cAcAin()–+0=θ 10= k 1= cAin()1= cBin()10=cA10.9092= cB0.0908=Jλ1λ2…λN,,,Jλ1λ2…λN,,,JOctober 2, 2002 4Background (should not need to solve problem)A common use of eigenvalue analysis is to study the stability of a dynamic system. Con-sider a system governed by the following set of N linear first-order differential equations,(EQ 13)Defining the vectors and matrix,(EQ 14)this system of differential equations takes the form,(EQ 15)We see that is a steady state since the time derivatives are zero, . But will it be a stable steady state? That is, if we change by a small value away from will it return to the steady state ( is stable) or will it move away ( is unstable)?To answer that question, we define the eigenvalues and eigenvectors of A as(EQ 16)We assume that the eigenvectors of A form a complete basis for all N-dimensional vectors, so that we may write any time dependent vector as the linear combination(EQ 17)If we substitute this form of into the differential equation, we obtain(EQ 18)dx1dt--------x·1a11x1a12x2… a1NxN+++==dx2dt--------x·2a21x1a22x2… a2NxN+++==:dxNdt---------x·NaN1x1aN2x2… aNNxN+++==xx1x2:xN= x·x·1x·2:x·N= Aa11a12… a1Na21a22… a2N:: :aN1aN2… aNN=x·Ax=x 0=x·0=xx 0=x 0=x 0=Awj[]λjwj[]= j 12…N,, ,=xt()xt() c1t()w1[]c2t()w2[]… cNt()wN[]+++=xt()tddc1w1[]c2w2[]… cNwN[]+++[]Ac1w1[]c2w2[]… cNwN[]+++[]=c·1w1[]c·2w2[]… c·NwN[]+++ c1λ1w1[]c2λ2w2[]… cNλNwN[]+++=October 2, 2002 5This yields the following set of N uncoupled differential equations,(EQ 19)If the state of the system at time is(EQ 20)the time evolution of the coefficients are(EQ 21)If we write each eigenvalue in terms of its real and imaginary parts,(EQ 22)and use Euler’s formula, , the time evolution of each coefficient becomes(EQ 23)We see then that for any possible values of the coefficients at , if the real parts of each eigenvalue are negative, as and the steady state at is