Assignment 3 - 2010

Due date(s):	01 November 2010
(PDF)	Assignment questions
Other instructions	Hand-in at class.

<rst> <rst-options: 'toc' = False/> <rst-options: 'reset-figures' = False/> .. |m3| replace:: m\ :sup:`3`

Question 1 [2]

====

The Secant method is the most widely used algorithm for solving a nonlinear equation, in chemical engineering, and other areas. Write a MATLAB or Python function that implements the Secant method, using a minimum number of function evaluations. (Python users: please do not just copy the built-in Secant method).

Your function should accept the following inputs:

• a function handle to :math:`f(x)` [see the `software tutorial <http://modelling3e4.connectmv.com/wiki/Software_tutorial/Functions_as_objects>`_ on the course website for a description of function handles]
• a relaxation coefficient, :math:`0 < \alpha \leq 1.0`
• the stopping tolerance, :math:`\varepsilon_\text{tol}`
• the maximum number of iterations, ``maxiter`` (default = 200)
• a single starting guess, :math:`x^{(0)}`

Your code should automatically calculate :math:`x^{(-1)} = 1.001x^{(0)}`, and should use a minimum number of function evaluations. The stopping criterion should be as shown in the course notes:

.. math::

\left|f(x^{(k)})\right| \leq \epsilon_{\rm tol} \qquad \text{and/or} \qquad \left|\frac{x^{(k)}-x^{(k-1)}}{x^{(k)}}\right| \leq \epsilon_{\rm tol}

Please print out a table at each iteration that shows these 5 outputs:

• Iteration :math:`k`
• :math:`x^{(k)}`
• :math:`f(x^{(k)})`
• :math:`\displaystyle \left|\frac{x^{(k)}-x^{(k-1)}}{x^{(k)}}\right|`
• Relaxed :math:`\Delta x`

Question 2 [2]

===

In your reactor design course you will learn about adiabatic CSTR's exhibiting multiple steady states, depending on the reactor inlet temperature. In practice reactors are sometimes operated at an unstable steady point, due to process economics or improved conversion. The equation below was derived for such a system:

.. math::

f(T) = (T-T_\text{feed}) - 150 \times \frac{1.34\times 10^9 e^{-\frac{7554}{T}}}{1+1.34\times 10^9 e^{-\frac{7554}{T}}}

There is more than one temperature at which :math:`f(T) = 0`; each :math:`f(T) = 0` corresponds to a reactor steady state. Use the secant method, with relaxation factor :math:`\alpha = 0.9`, and :math:`T_\text{feed} = 300` K to find one of these steady state temperatures, using starting values of :math:`T^{(0)}=320` K and :math:`T^{(-1)}= 1.001 \times 320` K.

1. . Show the first 3 iterations by hand. Explains what happens in the 3rd iteration, by means of a plot of :math:`f(x)`, and marking the iterations used so far on this plot.
2. . What could you change to avoid the problem in iteration 3?
3. . Compare your manual calculations in step 1 with the computer output from the previous question. Also report the temperature your computer code finds for :math:`f(T)=0` at this starting guess. Select a reasonable tolerance value.
4. . Use the computer code from question 1, but with a starting value of :math:`T = 400` K, and report the steady state operating temperature found.
5. . Comment on the solutions found from the two different starting points.

.. Solution .. --------

.. Part 2: change initial guess, or change \alpha value

Question 3 [2]

==

In assignment 1 we derived the dynamic balance for a bioreactor system. The steady state nutrient concentration, :math:`{\sf \overline{N}}` and the steady state biomass concentration, :math:`{\sf \overline{B}}`, can be shown to be:

.. math::

\frac{dN}{dt} &= 0 = {\sf \overline{N}}_\text{in}\frac{F}{V} - {\sf \overline{N}}\frac{F}{V} - \left(\frac{1}{Y_{B}}\right)\left(\mu_{\text{max}}\frac{{\sf \overline{N}}}{K+{\sf \overline{N}}}\right){\sf \overline{B}}\\ \frac{dB}{dt} &= 0 = - {\sf \overline{B}}\frac{F}{V} + \left(\mu_{\text{max}}\frac{{\sf \overline{N}}}{K+{\sf \overline{N}}}\right){\sf \overline{B}}

Although we showed in assignment 1 that we can solve one equation and substitute it into the other, you should not do that for this question.

Rather, let :math:`\displaystyle \frac{dN}{dt} = f_1({\sf \overline{N}}, {\sf \overline{B}}) = 0` and let :math:`\displaystyle \frac{dB}{dt} = f_2({\sf \overline{N}}, {\sf \overline{B}}) = 0`. Solve this coupled set of nonlinear equations for :math:`{\bf x} = ({\sf \overline{N}}, {\sf \overline{B}})` using Newton's method.

Use the following approach in your answer:

#. Construct the Jacobian matrix, :math:`{\bf J}`, symbolically in terms of the variables used in the equations. #. Simplify the Jacobian by substituting in these known constants:

* :math:`{\sf \overline{N}}_\text{in} = 100` g/|m3| * :math:`F = 5000` |m3|/day * :math:`V = 1600` |m3| * :math:`Y_B = 0.80` conversion efficiency * :math:`\mu_\text{max} = 5` g/|m3| * :math:`K = 20` g/|m3|

#. Let your initial guess be :math:`{\bf x}^{(0)} = (40 \text{g/m$^3$}, 20 \text{g/m$^3$})`. Calculate 2 iterations of the multivariate Newton-Raphson method, to show :math:`{\bf x}^{(1)}` and :math:`{\bf x}^{(2)}`. Take a full Newton step at every iteration. #. Calculate and show the convergence tolerance, using the 2-norm for each iteration. #. Compare these iterations to the true solution presented by Elliot in assignment 1. Is Newton's method converging to this answer?

Feel free to code up the Jacobian matrix, and the functions :math:`f_1` and :math:`f_2` in MATLAB or Python and use your code to evaluate :math:`{\bf J}({\sf \overline{N}}, {\sf \overline{B}})`, to avoid manual calculation errors.

.. raw:: latex

\vspace{0.5cm} \hrule \begin{center}END\end{center}

</rst>