 |
|
| |
Three-dimensional Simulation of Phytoplankton Dynamics
|
Exacte Wetenschappen
Programma Computational Science
The project is entitled: Three-dimensional Simulation of Phytoplankton Dynamics
The title is sufficiently short to avoid the use of an acronym.
Dr. B.P. Sommeijer
Centre for Mathematics and Computer Science (CWI)
Department Modelling, Analysis and Simulation
Kruislaan 413, 1098 SJ Amsterdam
P.O.Box 94079, 1090 GB Amsterdam
tel. 020 - 592 4192, fax 020 - 592 4199
e-mail: bsom@cwi.nl
In collaboration with Dr. J. Huisman, Aquatic Microbiology, UvA
Dr. B.P. Sommeijer will act as contact person.
Phytoplankton takes up CO2 from the atmosphere and provides the basis
of the food chains in the oceans. Modelling of phytoplankton is vital to
understanding the dynamics of phytoplankton populations and the feedback of
marine phytoplankton on the greenhouse effect.
For many years, phytoplankton models have been studied in a well-mixed environment,
leading to the common paradigm that
phytoplankton blooms cannot develop if the upper water layer of a
stratified water column is deeper
than some critical depth [25]. Huisman et al. [16] discovered that
also the vertical mixing rate is a critical parameter for the development of blooms.
Using a diffusion-reaction model, they found that if the
vertical diffusion is less than a critical value, a bloom can develop,
irrespective of the depth of the water column. Recently, these
one-dimensional models were extended by including vertical transport processes
[5, 14]. This is motivated by the fact that phytoplankton species may have
a specific weight larger/smaller than water,
causing sinking/buoyancy of these species. These studies
revealed many interesting and surprising results.
In all the above models, only variation in the vertical direction was
considered. An important aim of this project is to incorporate horizontal flow
of the water. The resulting three-dimensional PDEs will be used to study
the so-called `biological pump'. By this we mean the
sequestration of CO2 from the atmosphere which is then transported to the bottom
of the ocean by sinking phytoplankton species. This mechanism significantly contributes
to the reduction of the greenhouse effect. First results obtained by the 1D model provide
new insight in the relevant parameters [18].
Furthermore, we will include a module for the nutrients (which, up to now, were assumed
to be amply available), and we will take into account the influence
of water column stratification on phytoplankton dynamics.
The combined model is a system of 3D integro-partial differential equations, requiring
novel and tailor-made numerical algorithms.
The research topic can be classified in the following disciplines:
- Numerical mathematics;
- Population dynamics.
The research team will be composed as follows:
|
Researcher | Specialization | Financed by | Fte/year |
|
dr. B.P. Sommeijer | numerical analysis | CWI | 0.3 |
|
dr. J. Huisman | biology/model | UvA | 0.1 |
|
prof. dr. J.G. Verwer | numerical analysis | CWI | pm. |
|
drs. J. Kok | software/HPC/visualization | CWI | 0.1 |
|
drs. J. Passarge | biology/model | UvA | 0.1 |
|
OIO, NN | algorithms/software | EW, Comp. Sc. | 1.0 |
We request funding for a Ph.D. student (OIO).
He/she should have a strong background in numerical methods for
partial differential equations and ordinary
differential equations and, preferably, in parallel computing.
Also some expertise in biological modelling is very welcome.
The project should result in a Ph.D. thesis. As promotors will act prof. dr.
J.G. Verwer (for the numerical aspects) and dr. J. Huisman (for the
biological/modelling aspects). It is to be expected that dr. Huisman will be shortly
appointed as professor of Aquatic Microbiology at the University of Amsterdam (UvA).
Prof. Verwer is head of the department Modelling, Analysis and
Simulation at CWI
and professor of numerical analysis at the University of Amsterdam.
Dr. Huisman is the new leader of the
department of Aquatic Microbiology at the Institute for
Biodiversity and Ecosystem Dynamics of the UvA and expert in the field of
phytoplankton dynamics.
The integration of the specific expertise of numerical mathematicians (from the CWI-MAS group)
and biologists (from the UvA-IBED group) can lead to a major step forward in a better
understanding of the dynamics of aquatic ecosystems. Therefore, a close co-operation is of
paramount relevance for the success of the project.
The members of
the research team actively participate in two research schools:
Thomas Stieltjes School for Mathematics, for MAS-CWI, and
Research School SENSE (Socio-Economic
and Natural Sciences of the Environment), for Aquatic Microbiology, UvA.
Modelling the dynamics of phytoplankton is of great importance to
many aspects of human interest.
The most important aspect is the role that sinking phytoplankton species play in
the downward export of biological carbon to the bottom of the ocean [24, 23], the so-called `biological pump'.
By this mechanism, several gigatons per year of carbon dioxide are removed from the
atmosphere [7], thus making a significant contribution to the
reduction of the greenhouse problem on earth. In fact, it has been suggested to increase
the phytoplankton concentration in the upper water layer by nutrient enrichment of
the ocean. In this way it is hoped to obtain an enhanced sequestration of carbon from
the atmosphere. Some first experiments in this direction, by means of iron-fertilization
of the Southern Ocean, have taken place. These experiments indicate, however, that the
measured amount of
downward carbon export hardly increases, in spite of the enhanced phytoplankton
concentration in the upper layer [3].
Based on predictions resulting from our model analysis, we describe in [18] a
possible explanation for this failure. It is made plausible that the increased turbidity
of the water plays a crucial role. However, the models that we are using are still quite
limited and require significant refinement in order to produce predictions that are
sufficiently reliable for policy makers to base decisions on.
Along with the two classical scientific lines of experimentation
and theoretical analysis, the value of computer simulation of
phytoplankton dynamics has now been established.
The ever-increasing capacity of computers
has given an enormous impulse to the development of codes for
performing real-life simulations.
The complexity of the full system, however, is so large,
that a number of simplifications still has to be incorporated in the
current models, simply to keep the amount of storage and CPU at
a realistic level.
Nevertheless, a possibility to bridge part of the gap between
present-day practice and the ultimate goal is to use
efficient, tailor-made numerical algorithms and innovative
computer science techniques. In this way, satisfactory results are
feasible to make a significant step towards
understanding the full complexity of the system.
Initially, the development and testing of the algorithms will be performed
on workstations (which are available at CWI). At a later stage the simulation
of real-life situations will be carried out on a parallel machine, like the
new SGI Origin 3000 (`Teras') at SARA.
Critical parameters.
In recent years, research on phytoplankton models has received much attention.
Initially, only models in a well-mixed environment were studied. This led to
the conclusion that phytoplankton blooms can develop only in waters where the depth
of water-column stratification is
shallower than some critical depth (see e.g. [25]).
This can be understood by noticing that only in shallow water the light
conditions are sufficient to let (depth-integrated) photosynthesis exceed
the (depth-integrated) natural losses
of phytoplankton. However, field observations did not
always confirm this critical depth concept
[27, 6], but a firm theory to explain
these observations was lacking. A crucial step forward was to release the
assumption of vertical homogeneity and to study the phytoplankton distribution
over depth [16]. As a result of this approach, the model changes to an
integro-partial differential equation, with variation in time and in space.
A surprising result of this study is that there exist two different mechanisms
that may cause a phytoplankton bloom. Apart from the critical depth, it was found
that the diffusive mixing rate plays a crucial role: below a certain critical
mixing rate (i.e., in quiet waters) blooms can always develop, independent of
the water depth.
Transport phenomena. A next step is to include vertical transport
processes (other than by diffusion) into the models. This is motivated by the fact
that phytoplankton species often have a specific weight different from that of
water. Species with a specific weight larger than water
will sink downwards. On the other hand, species with a smaller
weight give rise to buoyancy. Since phytoplankton
requires light for its photosynthesis, buoyant species have an
advantage over sinking species,
because buoyant species can easily reside in the well-lit
zone (i.e., close to the water surface). Therefore, including the
vertical velocity of specific species in the model has a significant
impact. In mathematical terms we have now changed
to a advection-diffusion-reaction model. First results have been
reported in [5, 14, 19].
Competition models. In each of the aforementioned
studies only one particular phytoplankton
species has been considered. A natural step is to study
systems in which various species are considered in a
competition context. An important factor that will play a role is that species that
stay in the upper layer (and hence, capture a lot of light) have the additional
advantage to shade other species at deeper levels. Nevertheless, as
it turns out, the balance between survival and depth for the various species is
very delicate and depends on values of other parameters as well, such as background turbidity,
incident light intensity, and specific light attenuation coefficients of the
different species [17].
Competition models have been studied before by many authors, but always using
simpler models (see [26, 9, 15]). Recent work on horizontal
turbulence in a competition context seems very promising [4, 20].
All topics described in the preceding subsection are very promising, in the sense that
new, and sometimes surprising, results were obtained with the model. However, a common
restriction in all applications so far is that either only variation in the
vertical direction
has been considered or only variation in the horizontal direction.
The main aim of the current project is to combine the effects due to the various
spatial directions. Hence, we will study the full three-dimensional growth
and transport of phytoplankton,
for example by imposing a horizontal
wind field or taking into account differences in barometric pressure.
A priori, the influence of such effects is hard to predict, but it will certainly
play a significant role. For example,
it is already well-known that horizontal turbulence may disrupt bloom formation if
horizontal dispersal rates exceed local growth rates (see e.g., [21]).
Little is known, however, of the interactive effects of horizontal and vertical
turbulent diffusion on the dynamics of sinking or buoyant phytoplankton, i.e. the
three-dimensional structure of phytoplankton blooms. Generally speaking,
simple flows in three dimensions can potentially generate very complicated
patterns [8].
Another important feature for the dynamics is the influence of nutrients. In all above
studies it was tacitly assumed that an ample supply of nutrients was available. It is
of course more realistic to release this ideal nutrient condition and to study
combined nutrient-phytoplankton models. Such models have been studied before (see
e.g., [26, 9]), but usually in a well-mixed, homogeneous environment. Hence, the
spatial variation that is incorporated in our models gives (literally) `an extra dimension'.
This topic nicely fits in the 3-dimensional model extension proposed above since
nutrient transport will be largely determined by the horizontal flow as well as
by upwelling from below.
A third relevant issue, as already pointed out in Section 6.1, is the role of the
turbidity of the water column. In [18] it was shown that the maximal sinking
velocity that can be sustained
by a sinking phytoplankton population is inversely proportional to the background
turbidity. Consequently, an increased turbidity (for example by iron-fertilization)
will lead to a decreased maximal vertical velocity, which may have a negative
influence on the export of carbon to the deep ocean.
Since this export is directly related to the greenhouse problem, it is of high
relevance for society.
Finally, we want to develop models to study the competition between many
phytoplankton species. In [17], some results have been published based on a
model that did not take into account the effects of sinking and buoyancy. Moreover, only
one spatial dimension (i.e., the vertical) was considered
and only a few resources. This model predicts a low species
diversity: blooms are dominated by one or a few species only. It is highly interesting
to study whether this prediction remains valid in the more sophisticated models where
three spatial dimensions and vertical transport have been incorporated,
and where both the number of species
and the number of resources will be increased.
To summarize:
the new modelling issues of this proposal are:
- the combination of horizontal and vertical variability;
- a combined three-dimensional nutrient-phytoplankton model;
- the role of turbidity on carbon export by sinking phytoplankton;
- three-dimensional competition models.
Numerical algorithms.
As said before, the mathematical model is of
advection-diffusion-reaction type with a highly
nonlinear `reaction' term originating from the space-dependent light gradient.
This reaction term, which is in fact the production term generated by phytoplankton
photosynthesis, appears in the model as
an integral (over depth) of the unknown concentrations. Hence, the model is of so-called
`integro-partial differential equation of Fredholm type' in 3D and for its solution
in the space-time domain we need a numerical approach. We will follow the
Method of Lines approach, that is, first the spatial differential
operators, as well as the integral term, will be replaced by discrete
approximations and subsequently the huge resulting system of ordinary differential
equations (ODEs) -- which is still continuous in time -- will be
integrated numerically.
For the spatial discretization we will use a (conservative) finite volume technique based
on so-called upwind approximations
for the advection terms and central approximations for the diffusion terms
(see e.g., [11] for details). The integral in the reaction term will
be approximated by a suitable quadrature formula. The population density over space
is often far from uniform. This means that efficient use can be made of a non-uniform
spatial grid: in regions with a high spatial activity a fine mesh will be used and
a coarse mesh elsewhere. In this way the total number of grid points, and hence
the CPU time, can be drastically reduced.
Next, the resulting system of ODEs is integrated in time.
Because of the wide range in temporal solution scales, this ODE system is stiff and,
to avoid a severe timestep restriction for stability reasons, an implicit integration
technique seems to be an appropriate choice. In the literature (see e.g. [10])
a wide variety of implicit integration techniques has been described. However, in
using this type of methods we are faced with the task to solve, in each time step,
a system of implicit relations to obtain the solution at the next point in time.
Due to the three spatial dimensions in combination with the integral term, the
resulting linear algebra problem is extremely expensive. Therefore, standard
implicit techniques such as Crank-Nicolson or BDF-type methods are not feasible and
we have to change to a time integration process that is adapted to this particular
application.
A possibility is to consider special time-stepping techniques, such as
- Splitting methods; see e.g., [22, 12]
- IMEX (IMplicit-EXplicit) methods; see e.g., [1, 28]
- Implicit methods using approximate matrix factorization, see e.g.,
[13].
A common feature of these methods is that they circumvent the coupled three-dimensional
computations. In fact, this huge problem is replaced by a sequence of simpler ones.
However, all these techniques are rather `special purpose' in nature and a successful
use requires adaptation to the particular application.
In our case of phytoplankton dynamics a
thorough analysis should reveal whether the above techniques lead to a feasible
solution process and, if so, which is the most efficient one. If
none of these techniques is satisfactory, then new algorithms have
to be developed.
Another feature which renders these models numerically interesting is the coupling
of the dynamical model with the equations for the nutrients and for dead organic matter.
It may turn out that here an operator-splitting technique (of Strang type) will be
most appropriate.
To the best of our knowledge, the approach to model the dynamics of phytoplankton
by means of
three-dimensional advection-diffusion-reaction equations, in the form of
integro-PDEs is new and has not been done
before. However, our own recent work on one-dimensional vertical plankton models
[14, 16, 18] and
recent work by others on two-dimensional horizontal plankton
models [4, 20] suggest that a three-dimensional approach will be
very promising, and should be given high priority if we wish to keep ahead in this
exciting field of research.
Inside the Netherlands
- prof. W. Admiraal
from the University of Amsterdam, Department of Aquatic Ecology and Ecotoxicology,
with interest for the interaction of toxicants/nutrients with the functioning of
aquatic ecosystems.
- dr. F.J. Weissing from the University of Groningen, Department of Genetics.
- collaboration has been sought with prof. H.J.W. de Baar, from the NIOZ, Texel, who
participated in several experiments with
iron-fertilization in the Southern Ocean [2].
Outside the Netherlands
- prof. U. Sommer from the Marine Institute, University of Kiel, Germany, who is
one of the leading experts on phytoplankton competition.
- we have been in touch with prof. A. Ménesguen,
director of the Département Ecologie Côtière from IFREMER (the French
institute for research and exploitation of the sea) in Plouzané, France.
In the department MAS at CWI, `Applications to the Life Sciences' is one of
the major research topics. This new topic (started January 2000) focuses on
mathematical problems from biology and medicine. Co-operation has been started
with groups working in cell-, neuro- and microbiology. All models have been
formulated in terms of PDEs, the numerical solution of which has been a major
research topic at CWI for the last 25 years.
Furthermore, during the last 10 years,
expertise has been acquired in constructing parallel methods for PDEs.
Hence, the proposed research nicely fits in the field of interest of the department MAS.
For the department of Aquatic Microbiology, which is embedded in the Institute for
Biodiversity and Ecosystem Dynamics (IBED) of the University of Amsterdam,
the study of the dynamics of phytoplankton species is the main activity. This has been the case
for many years, and certainly will be continued in the years to come since dr. Huisman
recently received a `NWO-PIONIER' grant to
further stimulate the research in this direction.
The Aquatic Microbiology group combines mathematical modelling, laboratory experiments
and field research on phytoplankton, spanning the entire range from molecular biology
and physiology to the population dynamics and ecosystem impacts of phytoplankton.
Hence, it is obvious that the proposed project is a natural part of the
current activities of this group.
-
The results, including both the analysis of the appropriate algorithms and the biological
implications of the model predictions, will be
submitted to international scientific journals.
-
In collaboration with dr. Huisman we want to define several case studies. These
problems will act as a benchmark and should contain sufficiently many
`real life' aspects to build conclusions on.
The simulation results based on these test examples will be submitted to
the open literature.
-
An important instrument for biologists working in the field
(as well as for students to get some feeling for these matters) is
a graphical tool for visualizing the computed results.
We will extend the numerical solvers with such a graphical tool,
in order to
easily visualize all kinds of three-dimensional scenarios and to
facilitate a better understanding of the aquatic ecosystem.
-
At the end of the project a symposium is planned to
present the obtained results and to exchange ideas with other
groups working on similar or related topics.
The 4-year research programme comprises six activities:
-
The first 6 months of the project are reserved for the
OIO to get acquainted with the relevant literature on this subject,
including the recent papers by Huisman et al.
-
The second step will be to extend the one-dimensional model to a
three-dimensional model. This step requires the construction and analysis of
a suitable numerical solver, as well as its implementation. Since this
extension is an important task in the project, it is expected to take a full year.
-
Then we will investigate the influence of nutrients, when incorporated in this
three-dimensional model. In this way, our usual assumption of an ample supply
of nutrients is made more realistic. This activity is planned to require 6 months.
-
Next, in a stratified version of the model, carbon export will be studied.
Special attention will be paid to the role of the turbidity since it turned out to
be a critical parameter in the one-dimensional model. This part is estimated to
be completed within 9 months.
-
The last topic concerns competition models. Here the models will be modified as
to allow the study of multi-species dynamics. It is to be expected that the
computing time and memory requirements will drastically increase in these kind
of simulations. Therefore, the use
of a parallel computer like the `Teras' is necessary in this part
of the project. Also the numerical algorithms will probably need some
modifications. This topic will take another year.
-
The remaining time will be reserved for documentation and preparation of
the thesis.
-
This project does not require special equipment.
-
In the second part of the project the use of
a supercomputer will be necessary (see also point 5 of Section 7). A total of
250 hours CPU time is a good estimate.
The costs for the full 4-year period are specified in the
following table.
|
Personnel | Fte/year | EW, Comp. Sc. costs |
| | fl |
|
dr. B.P. Sommeijer | 0.3 |
0 |
|
dr. J. Huisman | 0.1 |
0 |
|
prof. dr. J. G. Verwer | pm. |
0 |
|
drs. J. Kok | 0.1 |
0 |
|
drs. J. Passarge | 0.1 |
0 |
|
OIO, NN | 1.0 |
223,000 |
|
Travel budget | |
7,350 |
|
Consumables (software, etc.) | |
10,000 |
| |
|
Totals | |
240,350 |
Five most relevant, recent publications of the research team
- P.J. van der Houwen and B.P. Sommeijer, Approximate factorization for time-dependent partial differential equations,
J. Comput. Appl. Math. 128, 447 - 466 (2001).
- J.G. Verwer and B.P. Sommeijer, A numerical study of mixed parabolic-gradient systems,
J. Comput. Appl. Math. 132, 191 - 210 (2001).
-
J. Huisman and B.P. Sommeijer, Population dynamics and export production of sinking phytoplankton in stratified waters, final version in preparation,
CWI Report, MAS-R01xx, submitted to Limnol. Oceanogr. (2001).
-
J. Huisman, P. van Oostveen and F.J. Weissing,
Critical depth and critical turbulence: two different mechanisms for the development of phytoplankton blooms, Limnology and Oceanography 44,
1781-1788 (1999).
-
J. Huisman and F.J. Weissing, Biodiversity of plankton by species oscillations and chaos, Nature 402, 407-410 (1999).
References
-
U.M. Ascher, S.J. Ruuth and B. Wetton, Implicit-explicit methods for time-dependent PDEs, SIAM J. Numer. Anal. 32, 797-823 (1997).
-
H.J.W. de Baar, et al., Importance of iron for plankton blooms and carbon-dioxide drawdown in the Southern Ocean, Nature 373, 412-415 (1995).
-
P.W. Boyd, et al., A mesoscale phytoplankton bloom in the polar Southern Ocean stimulated by iron fertilization, Nature 407, 695-702 (2000).
-
A. Bracco, A. Provenzale and I. Scheuring, Mesoscale vortices and the paradox of the plankton, Proc. R. Soc. Lond. B 267, 1795-1800 (2000).
-
U. Ebert, M. Arrayás, N.M. Temme, B.P. Sommeijer and J. Huisman, Critical conditions for phytoplankton bloom, final version in preparation, CWI Report, MAS-R01xx,
submitted to J. of Math. Biology (2001).
-
H.C. Eilertsen, Spring blooms and stratification, Nature 363, 24 (1993).
-
P.G. Falkowski, R.T. Barber and V. Smetacek, Biogeochemical controls and feedbacks on ocean primary production, Science 281, 200-206 (1998).
-
G.O. Fountain, D.V. Khakhar, I. Mezic and J.M. Ottino, Chaotic mixing in a bounded three-dimensional flow, J. Fluid. Mech. 417, 265-301 (2000).
-
J.P. Grover, Resource Competition, Chapman & Hall, London (1997).
-
E. Hairer and G. Wanner, Solving Ordinary Differential Equations II; Stiff and Differential-algebraic Problems,
Springer Series in Computational Mathematics 14,
Springer-Verlag, Berlin, 2nd edition (1996).
-
C. Hirsch, Numerical Computation of Internal and External Flows, I: Fundamentals of Numerical Discretization, Wiley and Sons, Chichester (1988).
-
P.J. van der Houwen and B.P. Sommeijer, Splitting methods for three-dimensional transport models with interaction terms, J. Scient. Comput. 12, 215-231 (1997).
-
P.J. van der Houwen and B.P. Sommeijer, Approximate factorization for time-dependent partial differential equations,
J. Comput. Appl. Math. 128, 447 - 466 (2001).
-
J. Huisman, M. Arrayás, U. Ebert and B.P. Sommeijer, How do sinking phytoplankton species manage to persist in unstratified waters?,
CWI Report, MAS-R0105, submitted to Am. Nat. (2001).
-
J. Huisman, R.R. Jonker, C. Zonneveld and F.J. Weissing, Competition for light between phytoplankton species: experimental tests of mechanistic theory,
Ecology 80, 211-222 (1999).
-
J. Huisman, P. van Oostveen and F.J. Weissing, Critical depth and critical turbulence: two different mechanisms for the development of phytoplankton blooms, Limnol. Oceanogr. 44, 1781-1788 (1999).
-
J. Huisman, P. van Oostveen and F.J. Weissing, Species dynamics in phytoplankton blooms: incomplete mixing and competition for light,
Am. Nat. 154, 46-68 (1999).
-
J. Huisman and B.P. Sommeijer, Carbon export by marine phytoplankton requires clear water, Internal Report.
-
J. Huisman and B.P. Sommeijer, Population dynamics and export production of sinking phytoplankton in stratified waters, final version in preparation,
CWI Report, MAS-R01xx, submitted to Limnol. Oceanogr. (2001).
-
G. Karolyi, A. Pentek, I. Scheuring, T. Tel and Z. Toroczkai, Chaotic flow: the physics of species coexistence, Proc. Natl. Acad. Sci. 97, 13661-13665 (2000).
-
H. Kierstead and L. Slobodkin, The size of water masses containing plankton blooms,
J. Mar. Res. 12, 141-147 (1953).
-
G.I. Marchuk, Splitting and alternating direction methods, In: Handbook
of Numerical Analysis I (Eds: P.G. Ciarlet and J.L. Lions), North-Holland,
Amsterdam (1990).
-
J.L. Sarmiento, T.M.C. Hughes, R.J. Stouffer and S. Manabe, Simulated response of the ocean carbon cycle to anthropogenic climate warming, Nature
393, 245-249 (1998).
-
J.L. Sarmiento and C. le Quéré, Ocean carbon dioxide uptake in a model of century-scale global warming, Science 274, 1346-1350 (1996).
-
H.U. Sverdrup, On conditions for the vernal blooming of phytoplankton,
J. Cons. Perm. Int. Explor. Mer 18, 287-295 (1953).
-
D. Tilman, Resource Competition and Community Structure, Princeton Univ. Press,
Princeton (1982).
-
D.W. Townsend, M.D. Keller, M.E. Sieracki and S.G. Ackleson, Spring phytoplankton blooms in the absence of vertical water column stratification,
Nature 360, 59-62 (1992).
-
J.G. Verwer, W.H. Hundsdorfer and J.G. Blom, Numerical time integration for air pollution models, to appear in Surveys on Math. for Industry (2001).

Document
http://www.cwi.nl ,
last modified: Jul 10 11:58