Assignment 2 - 2010

Due date(s):	13 October, at class
(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 [1]

====

The Gauss elimination method is an effective algorithm for solving a set of linear equations, :math:`Ax = b`. Solve this set of equations by hand, showing the steps taken.

.. math::

\left\{\begin{array}{r} x_1 + x_2 + x_3 + x_4 = 0\hphantom{-}\\ x_1 - 2x_2 - 2x_3 - 2x_4 = 4\hphantom{-}\\ x_1 + 4x_2 - 4x_3 + x_4 = 2\hphantom{-}\\ x_1 - 5x_2 - 5x_3 - 3x_4 = -4\\ \end{array}\right.\\

Question 2 [2]

====

1. . Derive, by hand, the lower-triangular matrix :math:`L` (with ones on the diagonal), and upper triangular matrix :math:`U` so that :math:`A = LU`, using the matrix :math:`A` from question 1. Afterwards, show that :math:`A = LU`.

1. . Use MATLAB or Python's built-in functions to solve for :math:`L` and :math:`U`. Show your code and show the matrices from the software. Explain any disagreement. See `the software tutorial on matrix operations <http://modelling3e4.connectmv.com/wiki/Software_tutorial/Matrix_operations>`_ for help.

1. . Use **your** matrices for :math:`L` and :math:`U` from the first part of this question, and calculate the inverse of :math:`A`. Compare your inverse with the inverse calculated by MATLAB or Python - does it agree with the software's values? Also, does :math:`A A^{-1} = I` for your inverse?

Question 3 [4]

====

.. note:: This is an example of the type of question you will get in the take-home midterm.

An impure gas stream can be cleaned, by absorbing the impurity into the liquid phase of an appropriate liquid solvent. For example, carbon dioxide can be removed from flue gas by scrubbing it with monoethanolamine; or NO\ :sub:`2` can be absorbed into water to produce nitric acid.

In this question we will derive the model equations for a sequence of absorber stages (trays), in which a component in the gas phase is absorbed into the liquid phase.

First consider the :math:`n^\text{th}` stage in an absorber:

.. figure:: images/absorber-tray.png :alt: images/math :scale: 16 :align: center

where the molar gas flow is :math:`G` and the liquid molar flow is :math:`L`, and assume these are constant throughout the column. Let the number of moles held-up in the liquid phase of the tray be :math:`H_n`, the hold-up. And let's only consider a binary system, where the species we are interested in is called A, and the absorbing liquid is species B.

The mol **fraction** of species A, in the gas phase, leaving the tray is :math:`y_{A,n}`, and we also have :math:`y_{B,n}`, the mol fraction of B leaving: this gives of course :math:`y_{A,n} + y_{B,n} = 1`. But since we are interested only in species A in this question, and because we can always calculate B by subtraction from 1.0, we will drop the A and B subscripts, and only write -- for example -- :math:`y_n` to denote the gas mol fraction of A leaving the :math:`n^\text{th}` stage of the absorber.

Similarly, the mol fraction of A, in the liquid phase, leaving the stage is :math:`x_n`. And we will denote the liquid mol fraction of A entering as :math:`x_{n+1}`, while the gas-phase mol fraction of A entering is defined as :math:`y_{n-1}`.

Using Raoult's Law, and assuming fairly ideal conditions, and that the liquid and vapour phase are well-mixed and in equilibrium, we may relate the concentrations in the gas and liquid phase as:

.. math::

y_n P_\text{total} &= x_n P_\text{A} \gamma_A \\ y_n &= \beta x_n

where :math:`P_\text{total}` is the total pressure in the absorber, and :math:`P_\text{A}` is the pure-component vapour pressure of A, at the given conditions, and :math:`\gamma_A` is an activity coefficient that corrects for non-idealities (cross-reference your thermodynamics course for more background on this equation). In this question you may assume that :math:`\beta` is constant throughout the absorber.

1. . Write down a **dynamic** mass balance for species A, in terms of the nomenclature above, for this :math:`n^\text{th}` stage in the absorber tower.
2. . Now write down the **steady-state** mass-balance equation for the above, but only in terms of liquid-phase mol fractions. Also divide the steady-state tray balance equations through by :math:`\beta G` and simplify them using the term :math:`\delta = \displaystyle \frac{L}{\beta G}`
3. . Now repeat this previous step, and write the steady-state mass-balance equations for the five-stage system shown here.

.. figure:: images/absorber-column.png :scale: 25

		:align: center

1. . Rewrite these mass-balance equations in matrix form, :math:`Ax = b`, where:

.. math::

x = \left[\begin{array}{c} y_F \\ x_1 \\ x_2 \\ x_3 \\ x_4 \\ x_5 \\ x_F \\ y_5 \end{array}\right]

1. . Write a function that will generate the :math:`A` matrix and :math:`b` vector, given the variables in the system. Use this function and solve the system of equations for the two cases below:

======= ========= ========= ============= =========== =========== Case :math:`G` :math:`L` :math:`\beta` :math:`x_F` :math:`y_F` ======= ========= ========= ============= =========== =========== Case I 100 150 0.6 0.05 0.35 ------- --------- --------- ------------- ----------- -----------

	Case II 100       50        0.6           0.05        0.35

======= ========= ========= ============= =========== ===========

In each case, we are obviously interested in the liquid-phase concentration profile along the column, :math:`x_F` down to :math:`x_1`. Report these values and explain whether your answers make physical sense, or not, for the two cases.

Question 4 [3]

===

Write a MATLAB or Python function that implements the Jacobi (simultaneous update) method for solving a general system of equations, :math:`Ax = b`. The function should accept the following inputs:

• :math:`A`: any :math:`N \times N` matrix.
• :math:`b`: any :math:`N \times 1` vector.
• :math:`x_0`: any :math:`N \times 1` vector that is the initial guess for :math:`x`.
• :math:`\text{tol}`: a small number, so that we will stop the algorithm when the relative difference between successive iterations is small: i.e. stop when :math:`\displaystyle \frac{||x^{(k)} - x^{(k-1)}|| }{||x^{(k)}||} < \text{tol}`, where the norms are assumed to be the Euclidean norm.
• :math:`\text{maxiter}`: an integer, for example, 200, representing an upper bound on the number of iterations we are willing to perform.

The function does not have to perform pivoting to achieve diagonal dominance. The function should calculate and print out, at every iteration, the relative difference between the :math:`k^\text{th}` and :math:`k-1^\text{th}` iteration.

Consider this system of equations:

.. math::

y_F &= 0.5\\ G/L y_F - (\delta+1)x_1 + x_2 &= 0 \\ \delta x_1 - (\delta+1)x_2 + x_3 &= 0 \\ \delta x_2 - (\delta+1)x_3 + x_F &= 0 \\ x_F &= 0.1

where :math:`G/L = 1.6`, and :math:`\delta = 0.8`

Your answer should include:

1. . A print out of the function
2. . The relative difference between the estimates for the first 5 iterations of solving the system of **five** equations, start from :math:`y^{(0)}_F = x^{(0)}_1 = x^{(0)}_2 = x^{(0)}_3 = x^{(0)}_F = 0`.
3. . The number of iterations required, and the final solution, using a termination tolerance of :math:`\epsilon = 1 \times 10^{-9}`.
4. . What is the condition number of :math:`A`?

Question 5 [2]

=====

.. note:: From a previous exam, 5 marks out of 50, 3 hours. This question should be solved by hand, and is a typical question you can expect in the midterm and final.

For the following system of linear algebraic equations:

.. math:: -5x_1 + 0x_2 + 12x_3 =&\ 80 \\ 4x_1 -x_2 - x_3 =&\ -2 \\ 6x_1 + 8x_2 +0x_3 =&\ 45

1. . Use the Jacobi method to compute :math:`x^{(1)}` starting from :math:`x^{(0)} = \bf 0`.
2. . Use the Gauss-Seidel method to compute :math:`x^{(1)}` starting from :math:`x^{(0)} = \bf 0`.
3. . Use the Gauss-Seidel method with relaxation, :math:`\omega=0.9`, to compute :math:`x^{(1)}` starting from :math:`x^{(0)} = \bf 0`.

Bonus question [1]

=======

Use the pseudo code developed in the course notes to write a MATLAB or Python function that implements Gauss elimination, without pivoting. The function should take :math:`A` and :math:`b` as inputs, and return vector :math:`x`.

Use your function to solve this set of equations from question 1, and print the updated :math:`A` matrix and :math:`b` vector after every step.

.. raw:: latex

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

</rst>