/*
 * Electroblue, 2007
 * 
 * Model Status
 * 
 * Reproduces results as placed online by the iGEM07 Glasgow team,
 * as at February 2007. Unit formulation here is such that the
 * equations are dimensionally consistent.
 * 
 * Model Structure
 * 
 * This model is a Cellml version of the Glasgow team's award-winning
 * entry into iGEM07. Details of their formulation and model output
 * can be found here.
 * 
 * Here the Cellml aims to reproduce the original equations as
 * faithfully as possible. As such, this is not described using
 * 'standard parts' Cellml. At the time of writing, standard parts
 * Cellml is under development.
 * 
 * Some modularity considerations have been observed. The basic
 * architecture of the model has been broken up into four main
 * components. One for each 'device' or promoter-to-double terminator
 * DNA segment, and one for each 'bioenvironment' centered on a
 * specific reaction (TFS association in one case, PYO formation
 * from PCA in another). Degradation reactions are housed in the
 * appropriate BioEnvironment component.
 * 
 * Model Schematic
 * 
 * [[Image file: gilbert_2007.png]]
 * 
 * The component architecture of the Cellml formulation.
 * 
 * This Cellml formulation could be easily extended in at least
 * two ways: (1) Use Cellml 1.1 functionality and place the separate
 * components in different files via imports. (2) Tease apart the
 * Bioenvironment components into species-specific components.
 * At the moment, there is no reason to do this but it would make
 * the model more extensible.
 */

import nsrunit;
unit conversion on;
unit per_second=1 second^(-1);
unit uM=1E-3 meter^(-3)*mole^1;
unit uM_per_second=1E-3 meter^(-3)*second^(-1)*mole^1;
unit per_uM_per_second=1E3 meter^3*second^(-1)*mole^(-1);

math main {
	realDomain time second;
	time.min=0;
	extern time.max;
	extern time.delta;
	real alphaTF uM_per_second;
	alphaTF=0.05;
	real Device_TFgenerator.JTF uM_per_second;
	real TFS(time) uM;
	when(time=time.min) TFS=0;
	real Jgain_TF(time) uM_per_second;
	real Jgain_TFS uM_per_second;
	real JTF2(time) uM_per_second;
	real TF(time) uM;
	when(time=time.min) TF=8.6207;
	real s uM;
	s=5;
	real betaTFS per_uM_per_second;
	betaTFS=1e6;
	real kd per_second;
	kd=4e6;
	real deltaTFS per_second;
	deltaTFS=3.851e-4;
	real deltaTF per_second;
	deltaTF=5.8e-3;
	real Jassociation(time) uM_per_second;
	real gammaPhzMS uM;
	gammaPhzMS=5;
	real betaPhzMS uM_per_second;
	betaPhzMS=0.1;
	real JPhzMS(time) uM_per_second;
	real gammaTF uM;
	gammaTF=4;
	real betaTF uM_per_second;
	betaTF=0.07;
	real feedbackOn dimensionless;
	feedbackOn=1;
	real Jgain_PhzMS(time) uM_per_second;
	real Jgain_PYO uM_per_second;
	real PhzMS(time) uM;
	when(time=time.min) PhzMS=0;
	real PYO(time) uM;
	when(time=time.min) PYO=0;
	real deltaPhzMS per_second;
	deltaPhzMS=8.0225e-6;
	real alphaPYO per_second;
	alphaPYO=1.3;
	real deltaPYO per_second;
	deltaPYO=5.8e-1;
	real JPYOformation(time) uM_per_second;

	// <component name="temporalExistence">

	// <component name="Device_TFgenerator">
	Device_TFgenerator.JTF=alphaTF;

	// <component name="BioEnv_TFSAssociation_interface">
	Jgain_TFS=(0 uM_per_second);
	Jgain_TF=(Device_TFgenerator.JTF+JTF2);

	// <component name="BioEnv_TFSAssociation">
	Jassociation=(betaTFS*s*TF-kd*TFS);
	TF:time=((-1)*Jassociation+Jgain_TF-deltaTF*TF);
	TFS:time=(Jassociation-deltaTFS*TFS+Jgain_TFS);

	// <component name="Device_PhzMSgenerator">
	JTF2=(if (feedbackOn<>0) betaTF*TFS/(gammaTF+TFS) else (0 uM_per_second));
	JPhzMS=(betaPhzMS*TFS/(gammaPhzMS+TFS));

	// <component name="BioEnv_PhzMStoPYO_interface">
	Jgain_PYO=(0 uM_per_second);
	Jgain_PhzMS=JPhzMS;

	// <component name="BioEnv_PhzMStoPYO">
	JPYOformation=(alphaPYO*PhzMS);
	PhzMS:time=(Jgain_PhzMS-deltaPhzMS*PhzMS);
	PYO:time=(Jgain_PYO-deltaPYO*PYO+JPYOformation);
}