/*
 * A delay-differential equation model of the feedback-controlled
 * hypothalamus-pituitary-adrenal axis in humans
 * 
 * Model Status
 * 
 * This CellML model runs in OpenCell but it does not replicate
 * the published results (paramteter values for delta1 2 and 3
 * have been taken from the legend of figure 3 in the published
 * paper, and initial conditions for x, y and z have been taken
 * from the grpahs in figure 3). The inability of the CellML model
 * to replicate the published results is most likely to be due
 * to the lack of time delays in the CellML model (which currently
 * can not be described). The units are all dimensionless and are
 * therefore consistent, however expressing time as "dimensionless"
 * means that the model cannot be run in COR.
 * 
 * Model Structure
 * 
 * ABSTRACT: The present work develops and analyses a model system
 * of delay-differential equations which describes the core dynamics
 * of the stress-responsive hypothalamus-pituitary-adrenal axis.
 * This neuroendocrine ensemble exhibits prominent pulsatile secretory
 * patterns governed by nonlinear and time-delayed feedforward
 * and feedback signal interchanges. Formulation and subsequent
 * bifurcation analysis of the model provide a qualitative and
 * mathematical frame work for a better understanding of the delayed
 * responsive mechanisms as well as the dynamic variations in different
 * pathological situations.
 * 
 * model diagram
 * 
 * [[Image file: lenbury_2005.png]]
 * 
 * Schematic diagram of the hypothalamus-pituitary-adrenal (HPA)
 * axis. Stimulatory and inhibitory paths are indicated by the
 * arrows and + or - signs respectively. CRH represents corticotropin-releasing
 * hormone and ACTH represents corticotropin.
 * 
 * The original paper reference is cited below:
 * 
 * A delay-differential equation model of the feedback-controlled
 * hypothalamus-pituitary-adrenal axis in humans, Yongwimon Lenbury
 * and Pornsarp Pornsawad, 2005, Mathematical Medicine and Biology,
 * 22, 15-33. PubMed ID: 15716298
 * 
 * It should be noted that in its current form, the CellML description
 * of the this model is unable to perfectly capture the simulation
 * results of the published model, this is due to the time delays
 * which are difficult to describe in the CellML code.
 */

import nsrunit;
unit conversion on;

math main {
	realDomain time dimensionless;
	time.min=0;
	extern time.max;
	extern time.delta;
	real x(time) dimensionless;
	when(time=time.min) x=1.0;
	real y(time) dimensionless;
	when(time=time.min) y=0.01;
	real beta dimensionless;
	beta=1.091;
	real delta1 dimensionless;
	delta1=0.5;
	real z(time) dimensionless;
	when(time=time.min) z=0.01;
	real delta2 dimensionless;
	delta2=0.38;
	real delta3 dimensionless;
	delta3=0.6;

	// <component name="environment">

	// <component name="x">
	x:time=(delta1*exp(beta*(1-y^2))-delta1*x);

	// <component name="y">
	y:time=(delta2*x*exp(beta*(1-z^2))-delta2*y);

	// <component name="z">
	z:time=(delta3*y-delta3*z);

	// <component name="model_parameters">
}