Wyższa Szkoła Zarządzania i Bankowości
LECTURESNOTEBOOKSPACKAGES
3 Population Dynamics
3.1 Fibonacci Growth

1 Introduction to Simulation and Modeling
2 Discrete Medeling (L-Systems)
3 Population Dynamics
3.1 Fibonacci Growth
3.2 Malthusian Growth
3.3 Logistic Growth
3.4 Phase Trajectories and Limit Cycles
3.5 Lotka-Volterra Predator-Prey Populations
3.6 Gilpins model and transition Chaos
3.7 Additional Exercises
4 Number Representation and Error Propagation
5 Modeling with Random Numbers
6 Heat Transfer in a Rod (Connection Mathematica and C: MathLink)
7 Special Topics in Stochastic Finance
8 Appendix: Introduction to Mathematica
9 Population Dynamics in Vensim®PLE
     
 

III.1 Fibonacci Growth

Perhaps one of the first mathematical models of a population growth was due to Leonardo Pisano (Fibonacci) in 1202. Let us take a rabbit population model and we will  assume that we have a newborn pair of rabbits that will eventually mate.  We will  assume that rabbits are born in oppositely sexed pairs, which become sexually mature at age 1 month and that the gestational term is 1 month.  

Let [Graphics:Images/nb3_gr_2.gif] denote the number of pairs of rabbits at the end of each month.  Lets look at the different ways we can compute the number of pairs of rabbits at the end of the nth month.

III.1.1 Simple but Inefficient Mathematica Code for the Fibonacci Recurrence

The most natural way to describe the Fibonacci relation in Mathematica is as follows.

[Graphics:Images/nb3_gr_3.gif]

REQUIRED  
III.1.1a
Derive this function from the described rabbit population.

OPTIONAL
III.1.1b
Note that if we want to make a table of Fibonacci numbers, this formula for [Graphics:Images/nb3_gr_4.gif] is horribly inefficient. Explain why.  

REQUIRED
III.1.1c
Use the Timing function in Mathematica to measure how long a calculation of Fibonacci series takes for different sizes of the problem. Explain why there is such a difference.

REQUIRED
III.1.1d
Using the table command (see "Introduction to Mathamatica")  plot  and fit  the timings.  Any comments?
Hint : look at the function Fit to help you do this.

[Graphics:Images/nb3_gr_5.gif]
[Graphics:Images/nb3_gr_6.gif]

III.1.2 Improved Code

Mathematica can be convinced to cache its previously computed information about a sequence as follows.
ADVANCED
III.1.2a
What do we mean by 'caching the information'? (See "Introduction to Mathematica" for postponed and direct  definitions).

[Graphics:Images/nb3_gr_7.gif]
[Graphics:Images/nb3_gr_8.gif]

ADVANCED
III.1.2b
This approach works good so long as n not too large.  Exaplin  why?.

ADVANCED
III.1.2c
Experiment with this approach, and clear Mathematica's memory of the cached [Graphics:Images/nb3_gr_9.gif] values if you run Timing tests so you are not using previously stored information (explain this).

ADVANCED
III.1.2d
Compare the timings with the results found in III.1.1

III.1.3 Another Approach

The disadvantage in the above improvement is that in order to compute a Fibonacci number, all the intermediate numbers are stored as rules.  This will  eventually  require a large amount of storage for a term in the sequence if the index [Graphics:Images/nb3_gr_10.gif] is very large.  In fact, for sufficiently large [Graphics:Images/nb3_gr_11.gif], the process will fail due to exceeded capacity.  Also, the lookup procedure for rules is slow.  These problems can be overcome with an iterative procedure.

[Graphics:Images/nb3_gr_12.gif]
[Graphics:Images/nb3_gr_13.gif]
[Graphics:Images/nb3_gr_14.gif]
[Graphics:Images/nb3_gr_15.gif]
CLICK FOR BIG RESOLUTION

OPTIONAL
III.1.3a
Compare your timings with the previous results (plot in one graph for instance)

OPTIONAL
III.1.3b
This procedure really flies for the computation of a single value.  It does not have the recursion depth restriction of the previous procedure. Explain this in your own words.

III.1.4 Graphs

To determine the growth rate of the Fibonacci sequence, we can plot the points.

[Graphics:Images/nb3_gr_17.gif]

[Graphics:Images/nb3_gr_18.gif]

OPTIONAL
III.1.4a
Is this exponential?  (To see, you can plot the log of the data and let the data go further out) (These observations are made more precise in the next subsection)

III.1.5 A Closed Form Method

In this subsection, we will derive the standard closed form formula for the Fibonacci sequence to show that the terms  grow exponentially.

The method used to compute F3 can be rewritten in matrix form.

[Graphics:Images/nb3_gr_19.gif]

Now the Fibonacci sequence can be described by a matrix equation with appropriate initial conditions (i.e. closed form).

[Graphics:Images/nb3_gr_20.gif]

Since the eigenvalues of the matrix A are distinct, the basis can be changed to an eigenbasis (a basis consisting of independent eigenvectors) in order to solve the system. In particular, we can write the start vector X[1] as a linear combination of the eigenvectors of A. If we denote the eigenvalues of A with [Graphics:Images/nb3_gr_21.gif] and [Graphics:Images/nb3_gr_22.gif] and the associated eigenvectors with [Graphics:Images/nb3_gr_23.gif] and [Graphics:Images/nb3_gr_24.gif] respectively, the start vector [Graphics:Images/nb3_gr_25.gif], where [Graphics:Images/nb3_gr_26.gif]denotes the transposed form of a vector, can be written as

[Graphics:Images/nb3_gr_27.gif]

for some [Graphics:Images/nb3_gr_28.gif] and [Graphics:Images/nb3_gr_29.gif]. The pair [Graphics:Images/nb3_gr_30.gif] forms the coordinates of our original start vector with respect to the basis of eigenvectors of A.

Returning to the matrix formulation of the Fibonacci sequence, it is not difficult to see that for an arbitrary [Graphics:Images/nb3_gr_31.gif]

[Graphics:Images/nb3_gr_32.gif]

Combining this result with the previous expression for [Graphics:Images/nb3_gr_33.gif], we find

[Graphics:Images/nb3_gr_34.gif] = [Graphics:Images/nb3_gr_35.gif]
= [Graphics:Images/nb3_gr_36.gif]
= [Graphics:Images/nb3_gr_37.gif]
= [Graphics:Images/nb3_gr_39.gif] [Graphics:Images/nb3_gr_40.gif]

From the definition of the matrix A the eigenvectors [Graphics:Images/nb3_gr_41.gif] and [Graphics:Images/nb3_gr_42.gif] can be computed. It appears that both are of the form [Graphics:Images/nb3_gr_43.gif]. Substituting this into the expression for [Graphics:Images/nb3_gr_44.gif] gives us

[Graphics:Images/nb3_gr_45.gif] = [Graphics:Images/nb3_gr_46.gif] [Graphics:Images/nb3_gr_47.gif]

This brings us exactly where we want to be. Remember that in the vector [Graphics:Images/nb3_gr_48.gif] the first element denotes the n-th number in the Fibonacci sequence, and the second element the (n-1)-th number. Thus, by concentrating on the first element of [Graphics:Images/nb3_gr_49.gif] the closed form emerges immediately. We only need to know the eigenvalues  [Graphics:Images/nb3_gr_50.gif] and [Graphics:Images/nb3_gr_51.gif] of A, which can easily be determined, and the values for [Graphics:Images/nb3_gr_52.gif] and [Graphics:Images/nb3_gr_53.gif]. These last values easily follow from the equation [Graphics:Images/nb3_gr_54.gif] and the specific form of the eigenvectors. In vector form this equation reduces to

[Graphics:Images/nb3_gr_55.gif] [Graphics:Images/nb3_gr_56.gif]

The computation of p and q requires the solution of this linear system.
A formula for [Graphics:Images/nb3_gr_57.gif] can now be determined as follows.

[Graphics:Images/nb3_gr_58.gif]
[Graphics:Images/nb3_gr_59.gif]
[Graphics:Images/nb3_gr_60.gif]
[Graphics:Images/nb3_gr_61.gif]

REQUIRED
III.1.5a
Determine the complete closed form formula for the n-th term of the Fibonacci sequence given the values for [Graphics:Images/nb3_gr_62.gif] and [Graphics:Images/nb3_gr_63.gif] and the values for [Graphics:Images/nb3_gr_64.gif] and [Graphics:Images/nb3_gr_65.gif]. Describe the behaviour of the Fibonacci sequence for large [Graphics:Images/nb3_gr_66.gif].

To implement this closed form description as a Mathematica definition is easy.

[Graphics:Images/nb3_gr_67.gif]

OPTIONAL
III.1.5b
Notice that there is a considerable amount of algebra being done by Mathematica when using this formula to compute large terms in the sequence. Time this function and compare with previously obtained results. Do you have any idea on the differences in memory requirements?

OPTIONAL
III.1.5c
In comparing the execution times of the different implementations, you should have discovered large differences. Can you account for the differences by counting the number of operations involved?  


 
     
  Lectures | Notebooks | Packages

 
  Copyright © 2003 Wyższa Szkoła Zarządzania i Bankowości. Wszystkie prawa zastrzeżone
webmaster@wszib.edu.pl