MXB261 – Modelling and Simulation Science
Assignment 2 – A Simulation Project
Due Date: end of Tuesday 22 October 2019 (Week 13)
This is a Group Project with 4 students per group, and is worth 40% of your final grade. There are four individual tasks available, and each student must choose a different task.
If circumstances dictate that your group has only 3 students, then you may omit Task 3. Note: Tasks 1, 2 and 4 are still to be completed.
The final mark for each student is made up of the score for their individual component (out of 20 marks), plus a contribution from the Group Report (out of 20 marks) – this contribution will depend on the percentage the student contributed to the report. For example, if all students contributed equally (25% each), then all will receive the same score for the report, otherwise the mark will be reduced according to the percentage effort.
The individual component includes an oral discussion (5 of the 20 marks) of your contribution, with your tutor. This discussion will cover your understanding of the programming (and any issues) associated with the individual task you have selected. You should find an opportunity to complete this oral discussion during the workshop in Week 11, 12 or 13.
The theme of this project is how parameter values can determine the system dynamics of a parasite model, in both a temporal and spatial setting.
The group is to submit a zip archive containing all the m-files, movie files (.avi), a pdf of the report, and a Statement of Contribution showing who has worked on which individual task, together with the percentage contributed to the Group Report; these percentages must total 100.
The Parasite Model
Consider the following differential equation system
dX1 = k1X1X2 − k2X1
dt
dX2 =k3 −k4X2 −k5X1
dt
where X1 represents a population of parasites that feeds off a population X2. Parameter k1 represents the birth rate of X1, k2 represents the death rate of X1, k3 represents the rate of food growth, k4 represents food decay, and k5 represents food consumption by the parasites. All these constants are positive. You will solve the differential equation system using MATLAB’s ODE45.
1
Task 1: Delay Dynamics and a Parameter Sweep
Consider first the one-dimensional logistic equation of parasite growth dX(t) X(t)
dt =rX(t) 1− K
where K is the carrying capacity and r is the growth rate. We can simulate this using the
difference scheme
Xn+1 =Xn(1+hr−hrXn/K)
Now introduce a delay term that could, for example, represent a delay in the food source
dX(t) X(t−τ) dt =rX(t) 1− K .
The initial value is X0, defined on [−τ,0], where τ represents the delay. The associated difference scheme is (with τ = s h)
Xn+1 =Xn(1+hr−hrXn−s/K).
For example, if s = 1, we will have an initial value for the system defined on [−h, 0] (with a delay of h), while s = 10 corresponds to the initial values defined on [−10 ∗ h, 0] with a delay of 10h.
(a) For K = 100,X0 = 50,r = 2,h = 0.01, investigate the dynamics of the numerical
solution for delays s = 50, 100, 150, 200, 250, 300, 350, 400 over 2000 steps. Repeat this
investigation for r = 1 , h = 0.05, with the same delays. You should see an oscillatory 2
solution in the delay case. What is the relationship between the oscillatory peaks, their period, and the delay term? What do you observe about the parasite population?
Now consider the 2D representation of the parasite/food system given on page 1 of the assignment specification.
Suppose the following parameters have values as indicated: k1 = 1, k2 = 2, k5 = 3
and that (after some scaling), the initial values are X1(0) = 1, X2(0) = 1. The time span for the investigation is [0, T ] where T = 20 units.
(b) In addition to these values, consider first the case that k4 = 4. You are to perform a parameter sweep on k3 ∈ [0, 50] to find k3 -values that result in X1 → 0 + Tol or X2 → 2 ± Tol, for a specified tolerance Tol (reasonable values for Tol are 10−1 or 10−2). It is also a requirement that both populations X1 and X2 do not drop below 0; if this happens, you must discard those parameter choices. Plot your successful k3-values, according to the characteristics you observe in the system dynamics. You should determine the boundary of the regions that characterise k3 with respect to the observed dynamics.
(c) As above, but now implement a parameter sweep on both k3 and k4 (over [0, 50]). Plot successful parameter pairs (k3,k4), showing the region(s) corresponding to different system dynamics.
2
with h being the stepsize. for the parasites:
(d) Now fix k3 = 10 and perform a similar parameter sweep on both k4 and k5 over [0, 50], plotting successful parameter pairs (k4, k5). Are your results consistent with the previous parameter sweep from (b)? Why/why not?
(e) Write a concise paragraph to explain your results.
Task 2: Latin Hypercube Sampling in 3D
In this task, we fix two of the parameters: k1 = 1, k2 = 2. You are to implement Latin Hypercube Sampling on the 3D space of parameters k3, k4, k5 (each in the range [0,50]). For an appropriate mesh size and number of trials, build a population of successful parameter 3-tuples, that reflect the system of equations exhibiting the characteristics that either X1 → 0 + Tol or X2 → 2 ± Tol, with the constraint that neither X1 nor X2 drops below zero in [0, T ]; if this happens, you must discard those parameter choices. Take T = 20 time units, and have X1(0) = 1, X2(0) = 1.
You are to write your own code to carry out the Latin Hypercube Sampling. You can check that your code is working by comparing results with the built-in MATLAB command lhsdesign.
In your exploration of the 3D parameter space, you are to group your successful parameter 3-tuples so that a region in the 3D space corresponds to particular characteristics of the system dynamics. In addition to this, you should find the boundary between these regions. Demonstrate your results in a 3D visualisation that is appropriately labelled (axis labels and title). Write a concise paragraph to explain your results.
Task 3: Latin Hypercube Sampling in 4D
This task requires a Latin Hypercube Sampling to be implemented as for Task 2, but over a 4D parameter space. Here k1 is fixed at k1 = 1. Again, the parameter space is explored to determine regions that correspond to certain system dynamics.
Visualisation of the successful parameter 4-tuples is required, with appropriate labelling. Discuss the approaches you take to represent the 4th dimension in your plots. Write a concise paragraph to explain your results.
The Latin Hypercube Sampling requirements of Task 3 are very similar to the coding required for Task 2. It is expected that the students involved in Tasks 2 and 3 may like to discuss their programming approach for the sampling, but then to develop their own visualisations.
3
Task 4: Spatial agent-based implementation
The spatial implementation allows exploration of the system dynamics when interactions take place at an individual level rather than at population level.
Construct a grid of cells (200 x 200) and populate that grid with F agents for food and P agents for parasites, according to various densities (see below). The parasites will be positioned randomly, while the positions of the food will be either random or localised according to the food-placement strategy being investigated. There is to be one agent only per cell. The rules for the simulation are as follows:
• Each parasite in turn will attempt to move to a neighbouring cell that is N, S, E or W of its current position; if the new cell is empty, the parasite moves; if the new cell is already occupied by a parasite, the move does not take place; if the new cell is occupied by food, the food is consumed (the food is replaced by the parasite), and a birth event takes place with the new parasite being placed in the original cell.
• A parasite dies after f1 iterations and is removed from the grid.
• For all food agents, each agent dies if a uniform random sample (u ∼ U(0,1)) u < f2,
and it is removed from the grid. Note u has to be different for each food agent.
• After all the parasites and food agents have been processed for the current step, create f3 new food agents and position them either randomly in empty cells, or localised (in empty cells). Localised means in a region such as a quadrant of the grid. Note that when choosing a cell at random, if it is already occupied, place the new food agent in a nearby unoccupied cell.
For a selection of parameter values, and for various initial F and P population densities of 10%, 20%, 30% and 40% of the grid (assume F and P have the same population density), simulate the evolution of the system and hence describe the relationship between initial population density, food-placement strategies, and the observed system dynamics. For the scale of this spatial simulation, it will be useful to consider a value of f1 ∈ [0, 15], a value of f2 ∈ [0, 0.1], and f3 values such as 100, 200, 300 or 400.
To demonstrate your simulation, create a movie of the spatial stochastic model above, capturing the output at regular time points as the simulation evolves.
4
The Group Report
The purpose of the Group Report is to combine the results from the individual tasks, and to compare and contrast the results obtained from the different approaches and strategies.
The Report (which should be approximately 10 pages) should have the following structure:
• an Introduction, setting the scene for the content;
• a Methods (by Task) section, where you describe the approaches taken, for each Task; include here an analysis of the equilibrium solutions of the parasite model;
• a Results and Discussion section, where you combine, compare and contrast the results for parameter regions and the characteristics of the system dynamics from Tasks 1, 2 and 3; include here plots that justify your answers; also include results from Task 4 (the spatial agent-based simulation) and the interpretation of the effect of initial population densities and food-placement strategies; include a screenshot taken from the movie;
• a Conclusions section, where you summarise the main results in no more than a few sentences.
Guide to the Marking Schedule
5
Marks Breakdown
Individual component
Structure 2 Solution 8
Methods (by task) 3
Results 2
Oral component 5
Subtotal 20 Group component
Structure and presentation 3
Introduction 3
Discussion and comparison 10
Conclusions 4
Subtotal 20 Total 40
2 Code is well structured and well documented.
0 Code is not structured and/or not documented. 8 Solution is correct, all aspects of the task have
been considered.
4 Parts of the solution may be incorrect, most as-
pects of the task have been considered.
0 Little/no solution.
3 Methods for the task have been clearly and cor-
rectly outlined for reproducibility.
1 Methods for the task have been outlined, but
may not be correct, clear or reproducible.
0 Little/no methods.
2 Relevant and concise figures that clearly demon-
strate the solution have been produced. Figures
are clear, correctly and completely labelled.
0 Irrelevant, excessive or no figures that do not
clearly demonstrate the solution.
5 Conveys a depth of understanding of the pro-
gramming and issues for the Individual Task. 2.5 Conveys some understanding of the program-
ming and issues for the Individual Task.
0 No oral discussion, or poor level of understand-
ing conveyed.
3 Report has a clear and logical structure. Presen- tation is of a professional standard. No major spelling/grammatical errors.
0 Report does not have a clear or logical structure. Presentation is not professional. Major spelling or grammatical errors.
3 Introduction provides a short relevant back- ground for the exercises performed in the report. Motivation for each task and an outline of the structure of the report is provided.
1.5 Introduction missing key components, but out- lines the main features of the report.
0 Little/no introduction.
10 Discussion shows concisely presented, but exten-
sive, detailed and relevant insight. Results from individual tasks have effectively been compared. Plots, figures and tables are clear, concise and relevant.
5 Discussion shows some relevant insight. Re- sults from individual tasks have been compared. Plots, figures and tables may not be clear, con- cise and relevant.
0 Little/no discussion or comparison.
4 Conclusions accurately and concisely summarise
main findings of the project.
2 Conclusions summarise some main findings of
the project.
0 Little/no conclusions.