/*
 * Bursting, Chaos and Birhythmicity Originating from Self-modulation
 * of the Inositol 1,4,5-triphosphate Signal in a Model for Intracellular
 * Ca2+ Oscillations
 * 
 * Model Status
 * 
 * This model has been built with the differential expressions
 * of Houart's 1999 paper for various types of intracellular Ca2+
 * oscillations bassed on self-modulation of the Inositol 1,4,5-triphosphate
 * signal. The parameter values for this model are different for
 * different types of intracellular Ca2+ oscillations. The parameter
 * values used for this particular model were taken from Figure
 * 2 on page 513 and reproduces the corresponding curves. Parameter
 * values from pages 515 and 520 may be changed to reflect curves
 * for Bursting, Chaos, Quasiperiodicity and Birhythmicity. The
 * initial conditions have been set after allowing the model to
 * run for approximately 5 minutes until steady state. This file
 * is known to run in COR and OpenCell.
 * 
 * Model Structure
 * 
 * We investigate the various types of complex Ca2+ oscillations
 * which can arise in a model based on the mechanism of Ca2+-induced
 * Ca2+ release (CICR), that takes into account the Ca2+-stimulated
 * degradation of inositol 1,4,5-trisphosphate (InsP3) by a 3-kinase.
 * This model was previously proposed in the course of an investigation
 * of plausible mechanisms capable of generating complex Ca2+ oscillations.
 * Besides simple periodic behavior, this model for cytosolic Ca2+
 * oscillations in nonexcitable cells shows complex oscillatory
 * phenomena like bursting or chaos. We show that the model also
 * admits a coexistence between two stable regimes of sustained
 * oscillations (birhythmicity). The occurrence of these various
 * modes of oscillatory behavior is analysed by means of bifurcation
 * diagrams. Complex oscillations are characterized by means of
 * Poincare sections, power spectra and Lyapounov exponents. The
 * results point to the role of self-modulation of the InsP3 signal
 * by 3-kinase as a possible source for complex temporal patterns
 * in Ca2+ signaling.
 * 
 * Bursting, Chaos and Birhythmicity Originating from Self-modulation
 * of the Inositol 1,4,5-triphosphate Signal in a Model for Intracellular
 * Ca2+ Oscillations, Houart G, Dupont G, Goldbeter A, 1999, Bulletin
 * of Mathematical Biology , 61, 507-530 PubMed ID: 17883229
 * 
 * Model Diagram
 * 
 * [[Image file: houart_1999.png]]
 * 
 * Schematic representation of the model based on the interplay
 * between CICR and the Ca2+-stimulated degradation of InsP3. Besides
 * simple periodic oscillations, this model can produce complex
 * Ca2+ oscillations including bursting, chaos, quasiperiodic behavior,
 * as well as birhythmicity.
 */

import nsrunit;
unit conversion on;
unit minute=60 second^1;
unit per_minute=.01666667 second^(-1);
unit per_litre=1E3 meter^(-3);
// unit micromolar predefined
unit micromolar_per_minute=1.6666667E-5 meter^(-3)*second^(-1)*mole^1;

math main {
	realDomain time minute;
	time.min=0;
	extern time.max;
	extern time.delta;
	real V_0 micromolar_per_minute;
	V_0=2;
	real V_1 micromolar_per_minute;
	V_1=2;
	real beta dimensionless;
	beta=0.6;
	real V_in micromolar_per_minute;
	real V_M2 micromolar_per_minute;
	V_M2=6;
	real Z(time) micromolar;
	when(time=time.min) Z=0.15;
	real K_2 micromolar;
	K_2=0.1;
	real V_2(time) micromolar_per_minute;
	real V_M3 micromolar_per_minute;
	V_M3=20;
	real K_Z micromolar;
	K_Z=0.5;
	real K_A micromolar;
	K_A=0.2;
	real K_Y micromolar;
	K_Y=0.2;
//	Var below replaced by constant in model eqns to satisfy unit correction
//	real m dimensionless;
//	m=2;
	real Y(time) micromolar;
	when(time=time.min) Y=1;
	real A(time) micromolar;
	when(time=time.min) A=0.42;
	real V_3(time) micromolar_per_minute;
	real V_M5 micromolar_per_minute;
	V_M5=5;
	real K_5 micromolar;
	K_5=1;
	real K_d micromolar;
	K_d=0.4;
//	Var below replaced by constant in model eqns to satisfy unit correction
//	real p dimensionless;
//	p=2;
//	Var below replaced by constant in model eqns to satisfy unit correction
//	real n dimensionless;
//	n=4;
	real V_5(time) micromolar_per_minute;
	real k per_minute;
	k=10;
	real k_f per_minute;
	k_f=1;
	real epsilon per_minute;
	epsilon=0.1;
	real V_4 micromolar_per_minute;
	V_4=2;

	// <component name="environment">

	// <component name="Vin">
	V_in=(V_0+V_1*beta);

	// <component name="V2">
	V_2=(V_M2*Z^2/(K_2^2+Z^2));

	// <component name="V3">
	V_3=(V_M3*Z^2/(K_Z^2+Z^2)*Y^2/(K_Y^2+Y^2)*A^4/(K_A^4+A^4));

	// <component name="V5">
	V_5=(V_M5*A^2/(K_5^2+A^2)*Z^4/(K_d^4+Z^4));

	// <component name="cytosol">
	Z:time=(V_in-V_2+V_3+k_f*Y-k*Z);

	// <component name="internal_pool">
	Y:time=(V_2-V_3-k_f*Y);

	// <component name="InsP3_conc">
	A:time=(beta*V_4-V_5-epsilon*A);
}