Tutorial 5 - 2010 - Solution

From Process Model Formulation and Solution: 3E4
Revision as of 17:10, 29 October 2010 by Kevindunn (talk | contribs) (Created page with "{{OtherSidebar | due_dates = 25 October 2010 | dates_alt_text = Due date(s) | questions_PDF = Tutorial-5-2010.pdf | questions_text_alt = Tutorial questions | solutions_PDF = Tuto...")
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search
Due date(s): 25 October 2010
Nuvola mimetypes pdf.png (PDF) Tutorial questions
Nuvola mimetypes pdf.png (PDF) Solutions by Elliot Cameron
Other instructions Hand-in at class.

<rst> <rst-options: 'toc' = False/> <rst-options: 'reset-figures' = False/>

.. |m3| replace:: m\ :sup:`3` .. highlight:: python

.. rubric:: Solutions by Elliot Cameron.

.. note:: All questions below will consider the following problem.

The heat of reaction for a certain reaction is given by :math:`\Delta H_{r}^{0}(T)= -24097 -0.26 T+1.69\times 10^{-3}T^2 + {\displaystyle{\frac{1.5\times10^5}{T}}}\;` cal/mol. Compute the temperature at which :math:`\Delta H_{r}^{0}(T)= -23505` cal/mol.

Question 1 [1]

===

Show the first 4 iterations of the bisection method to solve for :math:`T`, justifying your choice for the initial bracket.

For iteration 2, 3, and 4, please also report these two relative errors:

.. math::

\varepsilon_\text{rel,x}^{(k)} = \displaystyle \left|\frac{x_m^{(k)} - x_m^{(k-1)}}{x_m^{(k)}} \right| \qquad\qquad \varepsilon_\text{rel,f}^{(k)} = \displaystyle \left|\frac{f(x_m^{(k)}) - f(x_m^{(k-1)})}{f(x_m^{(k)})} \right|

Solution


Since we are attempting to solve for the temperature associated with :math:`\Delta H_{r}^{0}(T) = -23505` cal/mol we are essentially asked to solve for the zero of the following function:

.. math::

f(T) = - 24097 - 0.26T + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T} - (-23505)

where

.. math::

\Delta H_{r}^{0}(T) = - 24097 - 0.26T + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T}

To apply the bisection method we first need to select a starting interval that the root is known to fall within. More specifically, we need to also be sure that the root is a "true root" in the sense that it crosses the y-axis (and therefore the function changes sign on either side of the root). **The bisection method is unable to locate roots that only touch the x-axis but do not cross it, e.g. a double root**. Given that this is a physically perceivable system of one variable we can select the lower and upper bounds graphically or based on prior engineering knowledge.

If we choose to attack the problem based on our previous knowledge then selecting the lower bound is easy since we are dealing with the Kelvin temperature scale. Therefore we know the lowest allowable temperature is 0 K. Selecting the upper bound based on prior knowledge is slightly trickier. If we knew what reaction this equation is describing then we might have some basic knowledge of how the heat of reaction changes with T. In this case we would be able to guess a suitable upper temperature. Other than this we are essentially stuck trying to overguess to ensure we have a root contained in our bounds. This technique, of course, carries with it the critical issue that if the function contains multiple roots (as it will be shown this one does) an over guess could easily contain more than one root and so would either confuse the solver or cause the solver to think that no roots exist (if the function is a quadratic shape, for example, then it is possible that the function could have the same sign on either end of a bounding interval despite containing both roots). For this reason I would recommend the graphical approach.

Plotting the function in MATLAB/Python is accomplished using the following scripts.

.. twocolumncode:: :code1: ../che3e4/Tutorials/Tutorial-5/code/Q1_Plot.m :language1: matlab :header1: MATLAB code :code2: ../che3e4/Tutorials/Tutorial-5/code/Q1_Plot.py :language2: python :header2: Python code

The result of this script is shown below. From the plot below it is clear that this function has two roots. One between the bounds (200,400) and one between the bounds (400,600).

.. figure:: ../che3e4/Tutorials/Tutorial-5/images/Q1_Figure.png :scale: 70 :align: center

As a final note before getting into the bisection solution, it is believed that some confusion has likely arisen over the function error tolerance equation. In this case f(x) refers to the original :math:`\Delta H_{r}^{0}(T)` function and not the vertically shifted function you are trying to find the roots of (i.e. :math:`\Delta H_{r}^{0}(T) + 23505`). Since in your adjusted function you are trying to solve for when the function is 0, it does not make sense to calculate a relative comparison to 0 (the numbers will explode and contain essentially no useful information). This is actually an error that can crop up when a root also happens to be zero, but naturally this is less easy to predict. So to quickly wrap up, when root finding, your best bets for stopping criteria are either the relative error on x from step to step (unless the root is zero), the value of the root finding function (this would be the most versatile, as you are no longer dividing by any value and you know your function will be going to zero), the relative error of the original function from step to step (unless this also is equal to zero at the root...).

For the bisection test we will consider the lower of the two bounds. Utilizing the procedure given on slide 8 of the *Section C: Nonlinear Algebraic Equations* slide set.

    • ITERATION 1**
  • STEP 0: INITIALIZATION*
  • :math:`T_{U}^{(0)} = 400`
  • :math:`T_{L}^{(0)} = 200`
  • :math:`\epsilon_{tol} = 10^{-6}` (for the sake of this question the selection of :math:`\epsilon_{tol}` was arbitrary)

Calculating the function value at upper and lower bounds.

.. math::

\begin{array}{rl} f\left(T_{U}^{(0)}\right) &= - 24097 - 0.26\left(T_{U}^{(0)}\right) + 1.69\text{x}10^{-3}\left(T_{U}^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{U}^{(0)}\right)} - \left(-23505\right)\\ f\left(T_{U}^{(0)}\right) &= - 24097 - 0.26\left(400\right) + 1.69\text{x}10^{-3}\left(400\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(400\right)} - \left(-23505\right)\\ f\left(T_{U}^{(0)}\right) &= -50.6000\\ & \\ f\left(T_{L}^{(0)}\right) &= - 24097 - 0.26\left(T_{L}^{(0)}\right) + 1.69\text{x}10^{-3}\left(T_{L}^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{L}^{(0)}\right)} - \left(-23505\right)\\ f\left(T_{L}^{(0)}\right) &= - 24097 - 0.26\left(200\right) + 1.69\text{x}10^{-3}\left(200\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(200\right)} - \left(-23505\right)\\ f\left(T_{L}^{(0)}\right) &= 173.6000\\ & \\ \Delta H_{r}^{0}\left(T_{U}^{(0)}\right) &= - 24097 - 0.26\left(T_{U}^{(0)}\right) + 1.69\text{x}10^{-3}\left(T_{U}^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{U}^{(0)}\right)}\\ \Delta H_{r}^{0}\left(T_{U}^{(0)}\right) &= - 24097 - 0.26\left(400\right) + 1.69\text{x}10^{-3}\left(400\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(400\right)}\\ \Delta H_{r}^{0}\left(T_{U}^{(0)}\right) &= -23555.6\\ & \\ \Delta H_{r}^{0}\left(T_{L}^{(0)}\right) &= - 24097 - 0.26\left(T_{L}^{(0)}\right) + 1.69\text{x}10^{-3}\left(T_{L}^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{L}^{(0)}\right)}\\ \Delta H_{r}^{0}\left(T_{L}^{(0)}\right) &= - 24097 - 0.26\left(200\right) + 1.69\text{x}10^{-3}\left(200\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(200\right)}\\ \Delta H_{r}^{0}\left(T_{L}^{(0)}\right) &= -23331.4 \end{array}

Testing if the initial bounds contain a root.

.. math::

\begin{array}{rl} Check &= f\left(T_{U}^{(0)}\right)\times f\left(T_{L}^{(0)}\right) \geq 0 ?\\ Check &= \left(-50.6000\right)\times\left(173.6000\right) \geq 0 ?\\ Check &= \left(-8784.16\right) \geq 0 ? \;\; \rightarrow \;\; \text{No, therefore a root is bounded} \end{array}

  • STEP 1: CALCULATE MID-POINT* :math:`T_{M}^{(1)}`

.. math::

\begin{array}{rl} T_{M}^{(1)} &= \frac{T_{U}^{(0)} + T_{L}^{(0)}}{2}\\ T_{M}^{(1)} &= \frac{\left(400\right) + \left(200\right)}{2}\\ T_{M}^{(1)} &= 300 \end{array}

.. math::

\begin{array}{rl} f\left(T_{M}^{(1)}\right) &= - 24097 - 0.26\left(T_{M}^{(1)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(1)}\right)^{2} + \frac{1.5\text{x}10^{5}}{T_{M}^{(1)}} - \left(-23505\right)\\ f\left(T_{M}^{(1)}\right) &= - 24097 - 0.26\left(300\right) + 1.69\text{x}10^{-3}\left(300\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(300\right)} - \left(-23505\right)\\ f\left(T_{M}^{(1)}\right) &= -17.9000\\ & \\ \Delta H_{r}^{0}\left(T_{M}^{(1)}\right) &= - 24097 - 0.26\left(T_{M}^{(1)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(1)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{M}^{(1)}\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(1)}\right) &= - 24097 - 0.26\left(300\right) + 1.69\text{x}10^{-3}\left(300\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(300\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(1)}\right) &= -23522.9 \end{array}

  • STEP 2: CALCULATE* :math:`f(T_{L}^{(0)}) \times f(T_{M}^{(1)})`

.. math::

\begin{array}{rl} Check &= f\left(T_{L}^{(0)}\right)\times f\left(T_{M}^{(1)}\right) < 0 ?\\ Check &= \left(173.6000\right)\times\left(-17.9000\right) < 0 ?\\ Check &= \left(-3107.44\right) < 0 ? \;\; \rightarrow \;\; \text{Yes, therefore the root falls between }T_{L}^{(0)}\text{ and }T_{M}^{(1)}. \end{array}

  • :math:`T_{U}^{(1)} = 300`
  • :math:`T_{L}^{(1)} = 200`
  • :math:`f(T_{U}^{(1)}) = -17.9000`
  • :math:`f(T_{L}^{(1)}) = 173.6000`
  • :math:`\Delta H_{r}^{0}(T_{U}^{(1)}) = -23522.9`
  • :math:`\Delta H_{r}^{0}(T_{L}^{(1)}) = -23331.4`
    • ITERATION 2**
  • STEP 1: CALCULATE MID-POINT* :math:`T_{M}^{(2)}`

.. math::

\begin{array}{rl} T_{M}^{(2)} &= \frac{T_{U}^{(1)} + T_{L}^{(1)}}{2}\\ T_{M}^{(2)} &= \frac{\left(300\right) + \left(200\right)}{2}\\ T_{M}^{(2)} &= 250 \end{array}


.. math::

\begin{array}{rl} f\left(T_{M}^{(2)}\right) &= - 24097 - 0.26\left(T_{M}^{(2)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(2)}\right)^{2} + \frac{1.5\text{x}10^{5}}{T_{M}^{(2)}} - \left(-23505\right)\\ f\left(T_{M}^{(2)}\right) &= - 24097 - 0.26\left(250\right) + 1.69\text{x}10^{-3}\left(250\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(250\right)} - \left(-23505\right)\\ f\left(T_{M}^{(2)}\right) &= 48.6250\\ & \\ \Delta H_{r}^{0}\left(T_{M}^{(2)}\right) &= - 24097 - 0.26\left(T_{M}^{(2)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(2)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{M}^{(2)}\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(2)}\right) &= - 24097 - 0.26\left(250\right) + 1.69\text{x}10^{-3}\left(250\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(250\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(2)}\right) &= -23456.375 \end{array}


  • STEP 2: CALCULATE* :math:`f(T_{L}^{(1)}) \times f(T_{M}^{(2)})`

.. math::

\begin{array}{rl} Check &= f\left(T_{L}^{(1)}\right)\times f\left(T_{M}^{(2)}\right) < 0 ?\\ Check &= \left(173.6000\right)\times\left(48.6250\right) < 0 ?\\ Check &= \left(8441.3\right) < 0 ? \;\; \rightarrow \;\; \text{No, therefore the root falls between }T_{M}^{(2)}\text{ and }T_{U}^{(1)}. \end{array}

  • :math:`T_{U}^{(2)} = 300`
  • :math:`T_{L}^{(2)} = 250`
  • :math:`f(T_{U}^{(2)}) = -17.9000`
  • :math:`f(T_{L}^{(2)}) = 48.6250`
  • :math:`\Delta H_{r}^{0}(T_{U}^{(2)}) = -23522.9`
  • :math:`\Delta H_{r}^{0}(T_{L}^{(2)}) = -23456.375`
  • STEP 3: CHECK STOPPING CRITERION*

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(2)} &= \left| \frac{T_{M}^{(2)}-T_{M}^{(1)}}{T_{M}^{(1)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(2)} &= \left| \frac{\left(250\right)-\left(300\right)}{\left(300\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(2)} &= 0.1\bar{6} < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\Delta H_{r}^{0}\left(T_{M}^{(2)}\right)-\Delta H_{r}^{0}\left(T_{M}^{(1)}\right)}{\Delta H_{r}^{0}\left(T_{M}^{(1)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\left(-23456.375\right)-\left(-23522.9\right)}{\left(-23522.9\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(2)} &= 0.0028 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

    • ITERATION 3**
  • STEP 1: CALCULATE MID-POINT* :math:`T_{M}^{(3)}`

.. math::

\begin{array}{rl} T_{M}^{(3)} &= \frac{T_{U}^{(2)} + T_{L}^{(2)}}{2}\\ T_{M}^{(3)} &= \frac{\left(300\right) + \left(250\right)}{2}\\ T_{M}^{(3)} &= 275 \end{array}

.. math::

\begin{array}{rl} f\left(T_{M}^{(3)}\right) &= - 24097 - 0.26\left(T_{M}^{(3)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(3)}\right)^{2} + \frac{1.5\text{x}10^{5}}{T_{M}^{(3)}} - \left(-23505\right)\\ f\left(T_{M}^{(3)}\right) &= - 24097 - 0.26\left(275\right) + 1.69\text{x}10^{-3}\left(275\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(275\right)} - \left(-23505\right)\\ f\left(T_{M}^{(3)}\right) &= 9.7608\\ & \\ \Delta H_{r}^{0}\left(T_{M}^{(3)}\right) &= - 24097 - 0.26\left(T_{M}^{(3)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(3)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{M}^{(3)}\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(3)}\right) &= - 24097 - 0.26\left(275\right) + 1.69\text{x}10^{-3}\left(275\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(275\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(3)}\right) &= -23495.2392 \end{array}

  • STEP 2: CALCULATE* :math:`f(T_{L}^{(2)})\times f(T_{M}^{(3)})`

.. math::

\begin{array}{rl} Check &= f\left(T_{L}^{(2)}\right)\times f\left(T_{M}^{(3)}\right) < 0 ?\\ Check &= \left(48.6250\right)\times\left(9.7608\right) < 0 ?\\ Check &= \left(474.6189\right) < 0 ? \;\; \rightarrow \;\; \text{No, therefore the root falls between }T_{M}^{(3)}\text{ and }T_{U}^{(2)}. \end{array}

  • :math:`T_{U}^{(3)} = 300`
  • :math:`T_{L}^{(3)} = 275`
  • :math:`f(T_{U}^{(3)}) = -17.9000`
  • :math:`f(T_{L}^{(3)}) = 9.7608`
  • :math:`\Delta H_{r}^{0}(T_{U}^{(3)}) = -23522.9`
  • :math:`\Delta H_{r}^{0}(T_{L}^{(3)}) = -23495.2392`
  • STEP 3: CHECK STOPPING CRITERION*

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(3)} &= \left| \frac{T_{M}^{(3)}-T_{M}^{(2)}}{T_{M}^{(2)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(3)} &= \left| \frac{\left(275\right)-\left(250\right)}{\left(250\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(3)} &= 0.1 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\Delta H_{r}^{0}\left(T_{M}^{(3)}\right)-\Delta H_{r}^{0}\left(T_{M}^{(2)}\right)}{\Delta H_{r}^{0}\left(T_{M}^{(2)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\left(-23495.2392\right)-\left(-23456.375\right)}{\left(-23456.375\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(3)} &= 0.0017 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

    • ITERATION 4**
  • STEP 1: CALCULATE MID-POINT* :math:`T_{M}^{(4)}`

.. math::

\begin{array}{rl} T_{M}^{(4)} &= \frac{T_{U}^{(3)} + T_{L}^{(3)}}{2}\\ T_{M}^{(4)} &= \frac{\left(300\right) + \left(275\right)}{2}\\ T_{M}^{(4)} &= 287.5 \end{array}

.. math::

\begin{array}{rl} f\left(T_{M}^{(4)}\right) &= - 24097 - 0.26\left(T_{M}^{(4)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(4)}\right)^{2} + \frac{1.5\text{x}10^{5}}{T_{M}^{(4)}} - \left(-23505\right)\\ f\left(T_{M}^{(4)}\right) &= - 24097 - 0.26\left(287.5\right) + 1.69\text{x}10^{-3}\left(287.5\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(287.5\right)} - \left(-23505\right)\\ f\left(T_{M}^{(4)}\right) &= -5.3218\\ & \\ \Delta H_{r}^{0}\left(T_{M}^{(4)}\right) &= - 24097 - 0.26\left(T_{M}^{(4)}\right) + 1.69\text{x}10^{-3}\left(T_{M}^{(4)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T_{M}^{(4)}\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(4)}\right) &= - 24097 - 0.26\left(287.5\right) + 1.69\text{x}10^{-3}\left(287.5\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(287.5\right)}\\ \Delta H_{r}^{0}\left(T_{M}^{(4)}\right) &= -23510.3218 \end{array}

  • STEP 2: CALCULATE* :math:`f(T_{L}^{(3)})\times f(T_{M}^{(4)})`

.. math::

\begin{array}{rl} Check &= f\left(T_{L}^{(3)}\right)\times f\left(T_{M}^{(4)}\right) < 0 ?\\ Check &= \left(9.7608\right)\times\left(-5.3218\right) < 0 ?\\ Check &= \left(-51.945\right) < 0 ? \;\; \rightarrow \;\; \text{Yes, therefore the root falls between }T_{L}^{(3)}\text{ and }T_{M}^{(4)}. \end{array}

  • :math:`T_{U}^{(4)} = 287.5`
  • :math:`T_{L}^{(4)} = 275`
  • :math:`f(T_{U}^{(4)}) = -5.3218`
  • :math:`f(T_{L}^{(4)}) = 9.7608`
  • :math:`\Delta H_{r}^{0}(T_{U}^{(4)}) = -23510.3218`
  • :math:`\Delta H_{r}^{0}(T_{L}^{(4)}) = -23495.2392`
  • STEP 3: CHECK STOPPING CRITERION*

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(4)} &= \left| \frac{T_{M}^{(4)}-T_{M}^{(3)}}{T_{M}^{(3)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(4)} &= \left| \frac{\left(287.5\right)-\left(275\right)}{\left(275\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(4)} &= 0.0\bar{45} < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(4)} &= \left| \frac{\Delta H_{r}^{0}\left(T_{M}^{(4)}\right)-\Delta H_{r}^{0}\left(T_{M}^{(3)}\right)}{\Delta H_{r}^{0}\left(T_{M}^{(3)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(4)} &= \left| \frac{\left(-23510.3218\right)-\left(-23495.2392\right)}{\left(-23495.2392\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(4)} &= 0.00064 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}


Question 2 [1]

===
  1. . Derive a :math:`g(x) = x` function to use in the fixed-point algorithm.
  2. . Show the first 3 iterations of using the fixed-point algorithm, starting with an initial guess of :math:`T = 380` K.
  3. . Will the fixed-point method converge for this problem, using your :math:`g(x)`?


Solution


  1. . The goal of fixed-point iterative methods is to locate "fixed-point(s)" associated with a given function. A fixed-point is special case of a function where in it maps a value :math:`x` to itself (what this means in layman's terms is if you put a value :math:`x` into a function you get :math:`x` out, i.e. :math:`x = g(x)`). Knowing that our objective is to find the root(s) of a known objective function (we could have just had a black box function...) we may then construct :math:`g(x)` in way that gives it the greatest chance of converging. Let us start by looking at our objective function:

.. math::

f(T) = 0 = - 24097 - 0.26T + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T} - (-23505)

Likely the first choice of g(x) that popped to everyone's mind was

.. math::

T = \frac{- 24097 + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T} - (-23505)}{0.26} = g(T)_{1}

Given that we have a function :math:`f(T) = 0`, another good choice would be one of :math:`g(T) = T + f(T)` or :math:`g(T) = T - f(T)`. Let's say we were trying to identify the lower of the two roots. In this case we would want :math:`f(T)` to yield a negative value when :math:`T > T_{root}` and a positive value when :math:`T < T_{root}`. Looking at the plot generated in Question 1 we see that the function is positive when :math:`T` is smaller than the lower root and greater than the upper root and is negative when :math:`T` is between the two roots. Therefore if we wished to locate the lower root :math:`g(T) = T - f(T)` would locally satisfy the desired properties and if we wished to find the upper root :math:`g(T) = T + f(T)` would locally satisfy the desired properties (Note that these are gross generalizations and do not prove stability or ensure convergence. They are merely educated guesses as to potentially good functions). Therefore we will also investigate both:

.. math::

g(T)_{2} = T + (- 24097 - 0.26T + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T} - (-23505))

.. math::

g(T)_{3} = T - (- 24097 - 0.26T + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T} - (-23505))


  1. . Starting with a value of :math:`T^{(0)} = 380 K`:

**FUNCTION** :math:`g(T)_{1}`

.. math::

\begin{array}{rl} \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(T^{(0)}\right) + 1.69\text{x}10^{-3}\left(T^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(0)}\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(380\right) + 1.69\text{x}10^{-3}\left(380\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(380\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= -23557.02715 \end{array}

*ITERATION 1*

.. math::

\begin{array}{rl} T^{(1)} &= \frac{- 24097 + 1.69\text{x}10^{-3}(T^{(0)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(0)}} - (-23505)}{0.26}\\ T^{(1)} &= \frac{- 24097 + 1.69\text{x}10^{-3}(380)^{2} + \frac{1.5\text{x}10^{5}}{(380)} - (-23505)}{0.26}\\ T^{(1)} &= 179.8955\\ & \\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(T^{(1)}\right) + 1.69\text{x}10^{-3}\left(T^{(1)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(1)}\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(179.8955\right) + 1.69\text{x}10^{-3}\left(179.8955\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(179.8955\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= -23255 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(1)} &= \left| \frac{T^{(1)}-T^{(0)}}{T^{(0)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(1)} &= \left| \frac{\left(179.8955\right)-\left(380\right)}{\left(380\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(1)} &= 0.5266 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(1)}\right)-\Delta H_{r}^{0}\left(T^{(0)}\right)}{\Delta H_{r}^{0}\left(T^{(0)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\left(-23255\right)-\left(-23557.02715\right)}{\left(-23557.02715\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(1)} &= 0.0128 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 2*

.. math::

\begin{array}{rl} T^{(2)} &= \frac{- 24097 + 1.69\text{x}10^{-3}(T^{(1)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(1)}} - (-23505)}{0.26}\\ T^{(2)} &= \frac{- 24097 + 1.69\text{x}10^{-3}(179.8955)^{2} + \frac{1.5\text{x}10^{5}}{(179.8955)} - (-23505)}{0.26}\\ T^{(2)} &= 1140.4\\ & \\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(T^{(2)}\right) + 1.69\text{x}10^{-3}\left(T^{(2)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(2)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(1140.4\right) + 1.69\text{x}10^{-3}\left(1140.4\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(1140.4\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= -22064 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(2)} &= \left| \frac{T^{(2)}-T^{(1)}}{T^{(1)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(2)} &= \left| \frac{\left(1140.4\right)-\left(179.8955\right)}{\left(179.8955\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(2)} &= 5.3394 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(2)}\right)-\Delta H_{r}^{0}\left(T^{(1)}\right)}{\Delta H_{r}^{0}\left(T^{(1)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\left(-22064\right)-\left(-23255\right)}{\left(-23255\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(2)} &= 0.0512 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 3*

.. math::

\begin{array}{rl} T^{(3)} &= \frac{- 24097 + 1.69\text{x}10^{-3}(T^{(2)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(2)}} - (-23505)}{0.26}\\ T^{(3)} &= \frac{- 24097 + 1.69\text{x}10^{-3}(1140.4)^{2} + \frac{1.5\text{x}10^{5}}{(1140.4)} - (-23505)}{0.26}\\ T^{(2)} &= 6682.6\\ & \\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= - 24097 - 0.26\left(T^{(3)}\right) + 1.69\text{x}10^{-3}\left(T^{(3)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(3)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(6682.6\right) + 1.69\text{x}10^{-3}\left(6682.6\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(6682.6\right)}\\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= 49659 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(3)} &= \left| \frac{T^{(3)}-T^{(1)}}{T^{(2)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(3)} &= \left| \frac{\left(6682.6\right)-\left(1140.4\right)}{\left(1140.4\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(3)} &= 4.8598 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(3)}\right)-\Delta H_{r}^{0}\left(T^{(2)}\right)}{\Delta H_{r}^{0}\left(T^{(2)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\left(49659\right)-\left(-22064\right)}{\left(-22064\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(3)} &= 3.2507 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}


**FUNCTION** :math:`g(T)_{2}`

.. math::

\begin{array}{rl} \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(T^{(0)}\right) + 1.69\text{x}10^{-3}\left(T^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(0)}\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(380\right) + 1.69\text{x}10^{-3}\left(380\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(380\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= -23557.02715 \end{array}

*ITERATION 1*

.. math::

\begin{array}{rl} T^{(1)} &= T^{(0)} - 24097 - 0.26T^{(0)} + 1.69\text{x}10^{-3}(T^{(0)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(0)}} - (-23505)\\ T^{(1)} &= (380) - 24097 - 0.26(380) + 1.69\text{x}10^{-3}(380)^{2} + \frac{1.5\text{x}10^{5}}{(380)} - (-23505)\\ T^{(1)} &= 283.2531\\ & \\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(T^{(1)}\right) + 1.69\text{x}10^{-3}\left(T^{(1)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(1)}\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(283.2531\right) + 1.69\text{x}10^{-3}\left(283.2531\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(283.2531\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= -23505.49146 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(1)} &= \left| \frac{T^{(1)}-T^{(0)}}{T^{(0)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(1)} &= \left| \frac{\left(283.2531\right)-\left(380\right)}{\left(380\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(1)} &= 0.080 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(1)}\right)-\Delta H_{r}^{0}\left(T^{(0)}\right)}{\Delta H_{r}^{0}\left(T^{(0)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\left(-23505.49146\right)-\left(-23557.02715\right)}{\left(-23557.02715\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(1)} &= 0.0022 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 2*

.. math::

\begin{array}{rl} T^{(2)} &= T^{(1)} - 24097 - 0.26T^{(1)} + 1.69\text{x}10^{-3}(T^{(1)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(1)}} - (-23505)\\ T^{(2)} &= (283.2531) - 24097 - 0.26(283.2531) + 1.69\text{x}10^{-3}(283.2531)^{2} + \frac{1.5\text{x}10^{5}}{(283.2531)} - (-23505)\\ T^{(2)} &= 282.7616\\ & \\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(T^{(2)}\right) + 1.69\text{x}10^{-3}\left(T^{(2)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(2)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(282.7616\right) + 1.69\text{x}10^{-3}\left(282.7616\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(282.7616\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= -23504.913333 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(2)} &= \left| \frac{T^{(2)}-T^{(1)}}{T^{(1)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(2)} &= \left| \frac{\left(282.7616\right)-\left(283.2531\right)}{\left(283.2531\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(2)} &= 0.0017 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(2)}\right)-\Delta H_{r}^{0}\left(T^{(1)}\right)}{\Delta H_{r}^{0}\left(T^{(1)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\left(-23504.913333\right)-\left(-23505.49146\right)}{\left(-23505.49146\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(2)} &= 0.000026 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 3*

.. math::

\begin{array}{rl} T^{(3)} &= T^{(2)} - 24097 - 0.26T^{(2)} + 1.69\text{x}10^{-3}(T^{(2)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(2)}} - (-23505)\\ T^{(3)} &= (282.7616) - 24097 - 0.26(282.7616) + 1.69\text{x}10^{-3}(282.7616)^{2} + \frac{1.5\text{x}10^{5}}{(282.7616)} - (-23505)\\ T^{(2)} &= 282.8483\\ & \\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= - 24097 - 0.26\left(T^{(3)}\right) + 1.69\text{x}10^{-3}\left(T^{(3)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(3)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(282.8483\right) + 1.69\text{x}10^{-3}\left(282.8483\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(282.8483\right)}\\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= -23505.0156 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(3)} &= \left| \frac{T^{(3)}-T^{(1)}}{T^{(2)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(3)} &= \left| \frac{\left(282.8483\right)-\left(282.7616\right)}{\left(282.7616\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(3)} &= 0.00031 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(3)}\right)-\Delta H_{r}^{0}\left(T^{(2)}\right)}{\Delta H_{r}^{0}\left(T^{(2)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\left(-23505.0156\right)-\left(-23504.913333\right)}{\left(-23504.913333\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(3)} &= 0.0000044 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}


**FUNCTION** :math:`g(T)_{3}`

.. math::

\begin{array}{rl} \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(T^{(0)}\right) + 1.69\text{x}10^{-3}\left(T^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(0)}\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(380\right) + 1.69\text{x}10^{-3}\left(380\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(380\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= -23557.02715 \end{array}

*ITERATION 1*

.. math::

\begin{array}{rl} T^{(1)} &= T^{(0)} - \left(- 24097 - 0.26T^{(0)} + 1.69\text{x}10^{-3}(T^{(0)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(0)}} - (-23505) \right)\\ T^{(1)} &= (380) - \left(- 24097 - 0.26(380) + 1.69\text{x}10^{-3}(380)^{2} + \frac{1.5\text{x}10^{5}}{(380)} - (-23505) \right)\\ T^{(1)} &= 432.0272\\ & \\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(T^{(1)}\right) + 1.69\text{x}10^{-3}\left(T^{(1)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(1)}\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(432.0272\right) + 1.69\text{x}10^{-3}\left(432.0272\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(432.0272\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= -23547 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(1)} &= \left| \frac{T^{(1)}-T^{(0)}}{T^{(0)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(1)} &= \left| \frac{\left(432.0272\right)-\left(380\right)}{\left(380\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(1)} &= 0.1369 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(1)}\right)-\Delta H_{r}^{0}\left(T^{(0)}\right)}{\Delta H_{r}^{0}\left(T^{(0)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\left(-23547\right)-\left(-23557.02715\right)}{\left(-23557.02715\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(1)} &= 0.00044 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 2*

.. math::

\begin{array}{rl} T^{(2)} &= T^{(1)} - \left(- 24097 - 0.26T^{(1)} + 1.69\text{x}10^{-3}(T^{(1)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(1)}} - (-23505) \right)\\ T^{(2)} &= (432.0272) - \left(- 24097 - 0.26(432.0272) + 1.69\text{x}10^{-3}(432.0272)^{2} + \frac{1.5\text{x}10^{5}}{(432.0272)} - (-23505) \right)\\ T^{(2)} &= 473.7196\\ & \\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(T^{(2)}\right) + 1.69\text{x}10^{-3}\left(T^{(2)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(2)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(473.7196\right) + 1.69\text{x}10^{-3}\left(473.7196\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(473.7196\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= -23524 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(2)} &= \left| \frac{T^{(2)}-T^{(1)}}{T^{(1)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(2)} &= \left| \frac{\left(473.7196\right)-\left(432.0272\right)}{\left(432.0272\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(2)} &= 0.0965 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(2)}\right)-\Delta H_{r}^{0}\left(T^{(1)}\right)}{\Delta H_{r}^{0}\left(T^{(1)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\left(-23524\right)-\left(-23547\right)}{\left(-23547\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(2)} &= 0.00095 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 3*

.. math::

\begin{array}{rl} T^{(3)} &= T^{(2)} - \left(- 24097 - 0.26T^{(2)} + 1.69\text{x}10^{-3}(T^{(2)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(2)}} - (-23505) \right)\\ T^{(3)} &= (473.7196) - 24097 - 0.26(473.7196) + 1.69\text{x}10^{-3}(473.7196)^{2} + \frac{1.5\text{x}10^{5}}{(473.7196)} - (-23505)\\ T^{(2)} &= 492.9904\\ & \\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= - 24097 - 0.26\left(T^{(3)}\right) + 1.69\text{x}10^{-3}\left(T^{(3)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(3)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(492.9904\right) + 1.69\text{x}10^{-3}\left(492.9904\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(492.9904\right)}\\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= -23510 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(3)} &= \left| \frac{T^{(3)}-T^{(1)}}{T^{(2)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(3)} &= \left| \frac{\left(492.9904\right)-\left(473.7196\right)}{\left(473.7196\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(3)} &= 0.0407 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(3)}\right)-\Delta H_{r}^{0}\left(T^{(2)}\right)}{\Delta H_{r}^{0}\left(T^{(2)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\left(-23510\right)-\left(-23524\right)}{\left(-23524\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(3)} &= 0.00060 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}


  1. . It would appear that :math:`g(T)_{1}` is rapidly diverging while :math:`g(T)_{2}` is converging to the lower root and :math:`g(T)_{3}` is converging to the upper root.

Question 3 [1]

==
  1. . Write the Newton-Raphson iteration formula that you would use to solve this nonlinear equation.
  2. . Apply 3 iterations of this formula, also starting from :math:`T = 380` K, and calculate the error tolerances.

Solution


  1. . The Newton-Raphson algorithm is given on slide 15-17 of the *Section C: Nonlinear Algebraic Equations* slide set.

To apply the Newton-Raphson method we must first calculate the derivative of our function:

.. math::

f'(T) = - 0.26 + 3.38\text{x}10^{-3}T - \frac{1.5\text{x}10^{5}}{T^{2}}

Therefore the Newton-Raphson iteration formula we use is the following:

.. math::

T^{(k+1)} = T^{(k)} - \frac{f(T^{(k)})}{f'(T^{(k)})} = T^{(k)} - \frac{- 24097 - 0.26T + 1.69\text{x}10^{-3}T^{2} + \frac{1.5\text{x}10^{5}}{T} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}T - \frac{1.5\text{x}10^{5}}{T^{2}}}

  1. .

We start with :math:`T^{(0)} = 380 K`:

.. math::

\begin{array}{rl} \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(T^{(0)}\right) + 1.69\text{x}10^{-3}\left(T^{(0)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(0)}\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= - 24097 - 0.26\left(380\right) + 1.69\text{x}10^{-3}\left(380\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(380\right)}\\ \Delta H_{r}^{0}\left(T^{(0)}\right) &= -23557.02715 \end{array}

*ITERATION 1*

.. math::

\begin{array}{rl} T^{(1)} &= T^{(0)} - \frac{- 24097 - 0.26T^{(0)} + 1.69\text{x}10^{-3}(T^{(0)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(0)}} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}T^{(0)} - \frac{1.5\text{x}10^{5}}{(T^{(0)})^{2}}}\\ T^{(1)} &= (380) - \frac{- 24097 - 0.26(380) + 1.69\text{x}10^{-3}(380)^{2} + \frac{1.5\text{x}10^{5}}{(380)} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}(380) - \frac{1.5\text{x}10^{5}}{(380)^{2}}}\\ T^{(1)} &= -3237.72940\\ & \\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(T^{(1)}\right) + 1.69\text{x}10^{-3}\left(T^{(1)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(1)}\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= - 24097 - 0.26\left(-3237.72940\right) + 1.69\text{x}10^{-3}\left(-3237.72940\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(-3237.72940\right)}\\ \Delta H_{r}^{0}\left(T^{(1)}\right) &= -5585.4322 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(1)} &= \left| \frac{T^{(1)}-T^{(0)}}{T^{(0)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(1)} &= \left| \frac{\left(-3237.72940\right)-\left(380\right)}{\left(380\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(1)} &= 9.52 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(1)}\right)-\Delta H_{r}^{0}\left(T^{(0)}\right)}{\Delta H_{r}^{0}\left(T^{(0)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(1)} &= \left| \frac{\left(-5585.4322\right)-\left(-23557.02715\right)}{\left(-23557.02715\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(1)} &= 0.763 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 2*

.. math::

\begin{array}{rl} T^{(2)} &= T^{(1)} - \frac{- 24097 - 0.26T^{(1)} + 1.69\text{x}10^{-3}(T^{(1)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(1)}} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}T^{(1)} - \frac{1.5\text{x}10^{5}}{(T^{(1)})^{2}}}\\ T^{(2)} &= (-3237.72940) - \frac{- 24097 - 0.26(-3237.72940) + 1.69\text{x}10^{-3}(-3237.72940)^{2} + \frac{1.5\text{x}10^{5}}{(-3237.72940)} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}(-3237.72940) - \frac{1.5\text{x}10^{5}}{(-3237.72940)^{2}}}\\ T^{(2)} &= -1640.311\\ & \\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(T^{(2)}\right) + 1.69\text{x}10^{-3}\left(T^{(2)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(2)}\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= - 24097 - 0.26\left(-1640.311\right) + 1.69\text{x}10^{-3}\left(-1640.311\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(-1640.311\right)}\\ \Delta H_{r}^{0}\left(T^{(2)}\right) &= -19214.81711 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(2)} &= \left| \frac{T^{(2)}-T^{(1)}}{T^{(1)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(2)} &= \left| \frac{\left(-1640.311\right)-\left(-3237.72940\right)}{\left(-3237.72940\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(2)} &= 0.493 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(2)}\right)-\Delta H_{r}^{0}\left(T^{(1)}\right)}{\Delta H_{r}^{0}\left(T^{(1)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(2)} &= \left| \frac{\left(-19214.81711\right)-\left(-5585.4322\right)}{\left(-5585.4322\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(1)} &= 2.44 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

*ITERATION 3*

.. math::

\begin{array}{rl} T^{(3)} &= T^{(2)} - \frac{- 24097 - 0.26T^{(2)} + 1.69\text{x}10^{-3}(T^{(2)})^{2} + \frac{1.5\text{x}10^{5}}{T^{(2)}} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}T^{(2)} - \frac{1.5\text{x}10^{5}}{(T^{(2)})^{2}}}\\ T^{(3)} &= (-1640.311) - \frac{- 24097 - 0.26(-1640.311) + 1.69\text{x}10^{-3}(-1640.311)^{2} + \frac{1.5\text{x}10^{5}}{(-1640.311)} - (-23505)}{- 0.26 + 3.38\text{x}10^{-3}(-1640.311) - \frac{1.5\text{x}10^{5}}{(-1640.311)^{2}}}\\ T^{(3)} &= -908.1982\\ & \\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= - 24097 - 0.26\left(T^{(3)}\right) + 1.69\text{x}10^{-3}\left(T^{(3)}\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(T^{(3)}\right)}\\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= - 24097 - 0.26\left(-908.1982\right) + 1.69\text{x}10^{-3}\left(-908.1982\right)^{2} + \frac{1.5\text{x}10^{5}}{\left(-908.1982\right)}\\ \Delta H_{r}^{0}\left(T^{(3)}\right) &= -22632.0781 \end{array}

.. math::

\begin{array}{rl} \epsilon_{tol,x}^{(3)} &= \left| \frac{T^{(3)}-T^{(2)}}{T^{(2)}} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,x}^{(3)} &= \left| \frac{\left(-908.1982\right)-\left(-1640.311\right)}{\left(-1640.311\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,x}^{(3)} &= 0.45 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \\ & \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\Delta H_{r}^{0}\left(T^{(3)}\right)-\Delta H_{r}^{0}\left(T^{(2)}\right)}{\Delta H_{r}^{0}\left(T^{(2)}\right)} \right| < \epsilon_{tol} ? \\ \epsilon_{tol,f}^{(3)} &= \left| \frac{\left(-22632.0781\right)-\left(-19214.81711\right)}{\left(-19214.81711\right)} \right| < \left(10^{-6}\right) ? \\ \epsilon_{tol,f}^{(3)} &= 0.178 < \left(10^{-6}\right) ? \;\; \rightarrow \;\; \text{No, therefore keep going}. \end{array}

Question 4 [1]

==

Comment on the 3 approaches used so far. Are your calculations what you would expect from each method?

Solution


From question 1 we observed that the bisection method is a robust yet slow root finding tool. As long as you know that your bounds contain a single root and that the function changes sign on either side of this root, the bisection method is guaranteed to converge to within any pre-specified tolerance given enough time. This was, of course, expected.

In question 2 we observed that fixed-point iteration methods converge much more rapidly than the bisection method but performance is highly dependent on the iterative function selected. Again this was to be expected.

Finally, we observed that the Newton-Raphson method wildly diverged on the first iteration but appeared to be rapidly coming back to the "root region" in the second and third iterations. It was expected that the Newton-Raphson method would be the fastest to converge and arguably it had the fastest iteration to iteration convergence. The problem in this case was the poor initial guess of 380 K. If we go back to the plot presented in the solution to question 1 we see that the function is near its minimum at 380 K. Thus it is no surprise that the function's derivative was near zero and so this caused the first Newton-Raphson iterate to explode (you are dividing by the derivative after all). A better initial guess would likely greatly improve convergence.

Bonus question [1]

=======

A naive code for the bisection algorithm would evaluate the function :math:`f(x)` at the three points, :math:`[x_\ell, x_m, x_u]` in every iteration. Fewer function evaluations can be obtained though.

Write a function, in either MATLAB or Python, that implements the bisection method, that evaluates :math:`f(x)` as few times as possible. You should report the following 8 outputs in each iteration: :math:`[x_\ell^{(k)},\, x_m^{(k)},\, x_u^{(k)},\, f(x_\ell^{(k)}),\, f(x_m^{(k)}),\, f(x_u^{(k)}),\, \varepsilon_\text{rel,x}^{(k)},\, \varepsilon_\text{rel,f}^{(k)}]`.

Use this code to find the solution to the above problem, within a tolerance of :math:`\sqrt{\text{eps}}` based on :math:`\varepsilon_\text{rel,x}^{(k)}`.

Solution


.. twocolumncode:: :code1: ../che3e4/Tutorials/Tutorial-5/code/Bisection.m :language1: matlab :header1: MATLAB code :code2: ../che3e4/Tutorials/Tutorial-5/code/bisection.py :language2: python :header2: Python code


You should run the MATLAB code as follows:

>> func = @(T)(-24097-(0.26.*T)+(0.00169.*T.^2)+(150000./T)-(-23505)) >> root = Bisection(200,400,func,1000,sqrt(eps))

or just run the Python code as shown above.


.. raw:: latex

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

</rst>