Team:Aberdeen Scotland/modeling/pde

Introduction
It was realised quite soon into the project that upon modelling the Inner Dynamics of the system, the individual modules could be collimated in to a single Partial Differential Equation – hereafter referred to as PDE – which would accurately describe the evolution of the system in space and time.

A PDE model was an excellent “next step” for our modelling effort – and although a physically meaningful model could not be generated; the following describes a concise modelling methodology which would merit implementation outside of iGEM – or indeed, in future iGEM competitions.

The Keller-Segel Chemotaxis system is a set of two coupled Nonlinear PDE’s; ostensibly Parabolic except for a mixed Hyperbolic-Elliptic term in the cross-diffusive flux. The terms “Parabolic”, “Hyperbolic” and “Elliptic” refer to the number of boundary conditions the system must satisfy – the dependence on which defines (or “sets”) the behaviour of the system.



Where is the cross-diffusive flux

The Keller-Segel system was perfect for modelling our PicoPlumber modules, inasmuch as Quorum Sensing could be effectively integrated into the formulation[1]. This allowed the advantage of analysing the dynamical evolution of the system in space & time - and in regard to the Internal Dynamics.

The following is an account of our researches into utilising these powerful Mathematical tools to describe our system.

One-Dimensional Diffusion Equation
As PDE’s are a challenging and actively-researched topic in Mathematics[2-3]; our modelling effort began by solving a much simpler problem by way of introduction to the field.

The one-dimensional Heat Equation is the classic introduction to Partial Differential Equations. It was modelled by utilising an Implicit Finite Difference method of our own composition. This allowed us to attain an aptitude for tackling the more difficult Keller-Segel system to come.

The Heat Equation and the Diffusion Equation of our chemo-attractant shared the same mathematical formulation – differing in fact, by the choice of the constant k – and so this introduction still retained a Biological pertinence.



Keller-Segel System for Chemotaxis
As mentioned above, the Keller-Segel System for Chemotaxis is the set of coupled Nonlinear PDE's which describe chemotactic motion. The mixed Hyperbolic-Elliptic term endows the system with rich behaviour – but also renders it insoluble by elementary Numerical Methods[4].

In searching for a more advanced numerical method, Discontinuous Galerkin Finite Elements were found in Systems Biology literature[4]. to have successfully modelled the system and identified parametric instabilities.

In their award-winning paper[5]. , Hillen and Painter describe augmentations to the standard Keller-Segel System, hereafter referred to as the "Minimal Model". These augmentations account for effects such as quorum sensing, nonlinear diffusion and cell-growth & multiplication. The following documents the theoretical knowledge necessary to create a numerical scheme to solve the system - which may be undertaken on a "Next steps" basis.

Comment is made on the extension of the Minimal Model to incorporate quorum sensing and cell-growth & multiplication - so as to be fully representative of the PicoPlumber system - in the final subsection of this page.

Minimal Model
The Chemotactic Minimal Model is the simplest "chassis" for describing a bacterial population chemotaxing towards a chemoattractant. It is of the form



where is the chemotactic sensitivity,  the bacterial density,  the concentration of the chemoattractant and,  their respective time derivatives.

In Ephsteyn & Kurganov[4], the aforementioned Discontinuous Galerkin Finite Element Method is applied to successfully solve the system. The following will provide an abridged summary of the method, covering all elements necessary to implement it.

Theoretical Considerations
The Finite Element Method is a numerical technique for the solution of PDE's by discretising the domain over which the PDE is solved into finite "elements". The discretised domain is then termed the mesh. Elements are interpolated by polynomial basis functions which serve to approximate the equation over the element without allowing errors in intermediate calculations to accumulate and render the final output meaningless.

Discontinuous Galerkin methods are an augmented family of Finite Element Methods developed by the Russian Engineer, Boris Galerkin. They are superior to other finite element methods for several reaons, including - but not restricted to - local, element-wise mass conservation and flexibility to choose non- and high order polynomial bases. Element-wise mass conservation is important, and naturally the bacteria must be able to "exist" computationally through all space & time.

To solve the Minimal Model we require a Mesh upon which to discretise the system and a reformulation of the Convection-Reaction-Diffusion to render it soluble by Discontinuous Galerkin methods.

The following documents the theoretical considerations necessary to obtain a robust solution. Ephsteyn & Kurganov[4] must be fully credited for their description of the numerical scheme outlined below - which is expanded upon concomitantly with the considerations of Hillen & Painter[1]- in a manner so as to provide clarity and insight to the reader.

Mesh Generation
The method begins by defining a Sobolev space over the domain  upon which the mesh is specified. A Sobolev space is a Linear Space of polynomial functions - strictly a Banach space of square-integrable functions - endowed with the Sobolev s-norm

,

For s a labelling-index and an element of the discrete space of discontinuous piecewise polynomials of degree  as a function of the nondegenerate quasi-uniform rectangular subdivisions



And for those w in the Sobolev space such that for w evaluated at an arbitary element ; w is in the space of polynomials of degree over, styled

While this may seem mathematically "verbose", the mesh itself and the non-degeneracy of the elements therein can be seen in the diagram of below - whilst the Linear space constructions remain abstract



This rigourously defines the mesh upon with the Discontinuous Galerkin method solves the Keller-Segel system - which is the key component of the algorithm.

Convection-Reaction-Diffusion Reformulation
Algebraically manipulating the Minimal Model by defining the vector, we obtain the variation



which is in a form soluble by Discontinous Galerkin methods. This form of the Minimal Model can in turn be considered as a Convection-Reaction-Diffusion system[4]



where the boldface 3-vectors are, , and.

Ephsteyn & Kurganov[4] however, report that such a re-formulation is prone to severe instabilities in the transition from a Hyperbolic region to an Elliptic one - upon implementation of a high-order Numerical Method. To implement the Discontinuous Galerkin method we therefore rewrite the system as



subject to the boundary conditions,. This may also be considered as a Convection-Reaction-Diffusion system of the form



once again where the boldface 3-vectors in this instance are where,  are the fluxes and  is the reaction term. It should similarly be noted that k=1 is a constant and the convective part of the system is requisitely Hyperbolic.

This renders the Minimal Model soluble by a Discontinuous Galerkin numerical scheme as proposed by Ephsteyn & Kurganov[4].

Together, the Mesh and the Convection-Reaction-Diffusion equation form the complete theoretical basis upon which the program is to be developed.

Analytic Form of Solution
In the nomenclature defined above, a solution to the Keller-Segel System must take the form



for some solution in the direct product space of the respective nondegenerate quasi-linear rectangular subdivisions, satisfying the weak formulation of the Minimal Model



where sums of the integrals are over and the set of interior edges, which in turn are split up into horizontal & vertical edges -  and  respectively - of the rectangular boundary.

Feasibility of developing Numerical Method
As seen in the precending section. even the simplest of analytic solutions to the Keller-Segel System using the Discontinuous Galerkin method is distinctly non-trivial - and would require a high degree of technical proficiency and time to implement. This was outwith the scope of our project, but were it to be attempted, would prove a robust and highly accurate means of collimating the individual PicoPlumber modelling modules.

Moreover, comment must be made on extensions to the basic Minimal Model "chassis" to incorporate effects such as quorum sensing and cell-growth & multiplication. These augmentations are discussed thoroughly in Hillen & Painter[1] but would require a re-derivation of some of the results given on this page - assuming compatibility with the Discontinuous Galerkin methods, that is.