Tutorial 1 - 2010 - Solution

Due date(s):	15 September 2010
(PDF)	Tutorial questions
(PDF)	Solutions
Other instructions	Hand-in at class

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

.. rubric:: Tutorial objectives

- Three questions that will help you get comfortable (again?) with MATLAB or Python - Creating vectors and matrices - Mathematical manipulation of variables - Displaying simple plots of these vectors and matrices

.. rubric:: Recap of tutorial rules

• Tutorials can be completed by yourself, or with one other person.
• Tutorials must be handed in by the end of the tutorial session in the lab.
• Please give the TA a paper submission, no electronic submissions - thanks!

Question 1 [1]

===

.. note::

This question is a recap of creating vectors and generating plots.

A pellet that is shot vertically upward with initial speed :math:`u` will have vertical displacement, :math:`s`, after a time, :math:`t`, given by the formula:

.. math::

s = s_\text{up} + s_\text{down} = ut - \frac{gt^2}{2}

where :math:`g=9.81` m/:math:`\text{s}^2` is the acceleration due to gravity, and we ignore the effect of air resistance. We would like to plot a graph that shows the vertical displacement as a function of time. Also assume that the initial speed is 50 m/s.

a) Create a vector, called ``t``, that contains the time steps from :math:`t=0` to :math:`t=10.0` in increments of 0.2 seconds. b) How many elements in this vector? Write down the line of code used to calculate the vector's size. c) Create a vector, ``s_up`` that contains the first term in the above equation.

* How many elements in this vector? * What is the smallest value in the vector? * What is the largest value?

d) Also create a vector ``s_down``, the downward displacement of the pellet. Write down the line of code used to create this vector. e) Next create a vector ``s``, as the **sum** of ``s_up`` and ``s_down``. Note that we can add and subtract vectors as long as they have the same dimension. f) Create a plot of the upward displacement vs time. Your solution should contain the line of code you use. g) Create a plot of the the downward displacement vs time. h) Finally, create a plot of the total displacement vs time.

* label your x- and y-axes with labels that also contain the *units of measurement*, * give the plot a suitable title, * and make sure there is a grid in the background.

Your solution to this question must include the intermediate questions, but you will be graded on the lines of code and the plot from the last part.

Solution

---

.. twocolumncode::

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

Plots

---

• • MATLAB**

.. figure:: ../che3e4/Tutorials/Tutorial-1/code/tut-1-q1-8-MATLAB.jpg :alt: images/math :scale: 50 :align: center

       *Click on figures to enlarge them*

• • Python**

.. figure:: ../che3e4/Tutorials/Tutorial-1/code/tut-1-q1-8-Python.png :scale: 50 :align: center

Question 2 [2]

==

.. note::

This question is a recap of for-loops and while-loops.

The square root of :math:`y`, a positive number, can be approximated using Newton's method, a method we will learn about later in the course. The method requires an initial guess, and then iterates (repeats over and over) certain calculations, until the approximation does not improve any more. Here is some pseudo-code for the method:

1. Assign a value to :math:`y`

2. Set the initial guess of the final answer to :math:`x = y/2.0`

3. Repeat a few times, in this tutorial, repeat 5 times:

* Replace :math:`x` with :math:`\displaystyle \frac{x+y/x}{2.0}`. In the future we will write such statements as: :math:`x \leftarrow \displaystyle \frac{x+y/x}{2.0}`

* Print the value of :math:`x` to the screen

4. Stop and compare your answer with the built-in ``sqrt()`` command.

Answer these questions:

• Implement this pseudo code in either MATLAB or Python and show the printed output when calculating :math:`\sqrt{13}`. MATLAB users: please type ``format long`` so that the output shows 15 digits of precision; Python shows the full precision by default.
• Change your code so that the 3rd step is replaced with a ``while``-loop that will stop the code when the change in :math:`x` from the previous iteration to the current iteration is smaller than :math:`1 \times 10^{-8}`.
• How many iterations are required to calculate :math:`\sqrt{65536.0}` with this criterion?

Solution

---

.. twocolumncode::

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

• • Selected MATLAB output**::

sqrt(13) = 3.605551e+000

while loop results: 6.500000000000000e+000 4.250000000000000e+000 3.654411764705882e+000 3.605877914546100e+000 3.605551290258318e+000 3.605551275463990e+000 3.605551275463990e+000

Number of while loop iterations: 6

sqrt(65536) = 256

while loop results: 3.276800000000000e+004 1.638500000000000e+004 8.194499877937138e+003 4.101248718697424e+003 2.058614121064604e+003 1.045224565257789e+003 5.539624834270534e+002 3.361332618959818e+002 2.655517766166163e+002 2.561717865301000e+002 2.560000575992625e+002 2.560000000000065e+002 2.560000000000000e+002

Number of while loop iterations: 12

• • Selected Python output**::

Used 6 iterations to calculate sqrt(13.000000) = 3.605551275463990; true value = 3.605551275463989

16385.0 8194.49987794 4101.2487187 2058.61412106 1045.22456526 553.962483427 336.133261896 265.551776617 256.17178653 256.000057599 256.0 256.0 Used 12 iterations to calculate sqrt(65536.000000) = 256.000000; true value = 256.000000

Bonus question [1.5]

========

The population count of the United States can be modelled using the formula:

.. math::

P(t) = \frac{197\, 252\, 000}{1 + e^{-0.03134(t-1913.25)}}

where the time, :math:`t`, is given in years. Write some code that will calculate the population count every **ten** years, from 1800 to 2010. Plot the results and label the plot appropriately.

Does the population ever reach steady-state? If so, roughly when, and what is the steady-state value?

Comment on whether or not this model is realistic. Is it useful?

Solution

---

.. note:: Please note that Elliot's solution is much more than was required to get full grade for this question. But please read his solution to gain some extra insight.

By the very definition of the model the population will approach but never reach the steady state value of 197252000 (the carrying capacity). Extrapolating the model, however, would seem to indicate that the steady state will be reached around **2010 - 2050** (it's hard to get a firm number).

In terms of the model's "usefulness" and "realism", as with anything in modelling, the answer really depends on what you are trying to achieve. The model is the classic Verhulst-Pearl population equation or simply the "Logistic equation" (for more information consult http://en.wikipedia.org/wiki/Logistic_growth). The basic form of the equation is :math:`\frac{dP}{dt} = rP\left(1 - \frac{P}{K}\right)`, which has an analytical solution of:

.. math:: P(t) = \frac{K}{1 + e^{-rt - \ln{\frac{P_0}{K - P_0}}}}

where:

• :math:`K` = carrying capacity (in our case 197252000)
• :math:`r` = growth rate (in our case 0.03134)
• :math:`P_0` = the initial population

The logistic equation is an ideal case where growth is assumed proportional only to the current population as well as the level of renewable resources available. The steady state achieved by the model represents the hypothetical population that could be sustained indefinitely given a constant resource level. Of course this is a gross over simplification and so the model is not very realistic.

However, depending on the accuracy of prediction one desires, this model could be useful in making short term predictions of population growth. Given that the population and infrastructure of a nation is a dynamic thing one would expect :math:`K` and :math:`r` to change over time (as various resources are discovered, exhausted, improved, etc.) but in the short term (assuming up to date estimates of :math:`K`, :math:`r`, and :math:`P_0`) this could be a useful tool for mapping population growth. As with any model though, always remember that the wider you make your extrapolation window, the more erroneous your approximation is likely to be.

As a famous statistician once said, *Essentially, all models are wrong, but some are useful*. (George E .P. Box)

.. twocolumncode::

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

.. raw:: latex

\newpage

• • Python figure**

.. figure:: ../che3e4/Tutorials/Tutorial-1/code/tut-1-qb1-Python.jpg

   :scale: 50 
   :align: center

.. raw:: latex

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

</rst>