/*
 * Modelling Actin-Myosin Binding in Airway Smooth Muscle
 * 
 * Model Status
 * 
 * This CellML model runs in both OpenCell and COR but to replicate
 * the published results (the first part of figure 3b). The units
 * have been checked and they are consistent.
 * 
 * Model Structure
 * 
 * ABSTRACT: We carried out a detailed mathematical analysis of
 * the effects of length fluctuations on the dynamically evolving
 * cross-bridge distributions, simulating those that occur in airway
 * smooth muscle during breathing. We used the latch regulation
 * scheme of Hai and Murphy (Am. J. Physiol. Cell Physiol. 255:C86-C94,
 * 1988) integrated with Huxley's sliding filament theory of muscle
 * contraction. This analysis showed that imposed length fluctuations
 * decrease the mean number of attached bridges, depress muscle
 * force and stiffness, and increase force-length hysteresis. At
 * frequencies >0.1 Hz, the bond-length distribution of slowly
 * cycling latch bridges changed little over the stretch cycle
 * and contributed almost elastically to muscle force, but the
 * rapidly cycling cross-bridge distribution changed substantially
 * and dominated the hysteresis. By contrast, at frequencies greater
 * than 0.033 Hz this behavior was reversed: the rapid cycling
 * cross-bridge distribution changed little, effectively functioning
 * as a constant force generator, while the latch bridge bond distribution
 * changed substantially and dominated the stiffness and hysteresis.
 * The analysis showed the dissociation of force/length hysteresis
 * and cross-bridge cycling rates when strain amplitude exceeds
 * 3%; that is, there is only a weak coupling between net external
 * mechanical work and the ATP consumption required for cycling
 * cross-bridges during the oscillatory steady state. Although
 * these results are specific to airway smooth muscle, the approach
 * generalizes to other smooth muscles subjected to cyclic length
 * fluctuations.
 * 
 * The original paper reference is cited below:
 * 
 * Perturbed Equilibria of Myosin Binding in Airway Smooth Muscle:
 * Bond-Length Distributions, Mechanics, and ATP Metabolism, Srboljub
 * M. Mijailovich, James P. Butler, and Jeffrey J. Fredberg, 2000,
 * Biophysical Journal, 79, 2667-2681. PubMed ID: 11053139
 * 
 * reaction diagram
 * 
 * [[Image file: mijailovich_2000.png]]
 * 
 * Hai and Murphy's four-state model: the latch regulatory scheme
 * for Ca2+-dependent smooth muscle activation and Huxley's slidin
 * filament model. A represents the actin filament, M represents
 * detached myosin, Mp is detached, phosphorylated myosin, AM is
 * the actin-myosin complex, and AMp is the phosphorylated actin-myosin
 * complex.
 */

import nsrunit;
unit conversion on;
// unit nanometre predefined
unit first_order_rate_constant=1 second^(-1);

math main {
	realDomain time second;
	time.min=0;
	extern time.max;
	extern time.delta;
	real M(time) dimensionless;
	when(time=time.min) M=1.0;
	real AM(time) dimensionless;
	when(time=time.min) AM=0.0;
	real Mp(time) dimensionless;
	when(time=time.min) Mp=0.0;
	real k1(time) first_order_rate_constant;
	real k2 first_order_rate_constant;
	k2=0.1;
	real gx first_order_rate_constant;
	gx=0.11;
	real AMp(time) dimensionless;
	when(time=time.min) AMp=0.0;
	real fp first_order_rate_constant;
	fp=0.44;
	real gp first_order_rate_constant;
	gp=0.11;
	real k5 first_order_rate_constant;
	k5=0.1;
	real k6(time) first_order_rate_constant;
	real fp1 first_order_rate_constant;
	fp1=0.88;
	real gp1 first_order_rate_constant;
	gp1=0.22;
	real g1 first_order_rate_constant;
	g1=0.01;
	real gp3 first_order_rate_constant;
	real g3 first_order_rate_constant;
	real g2 first_order_rate_constant;
	real gp2 first_order_rate_constant;

	// <component name="environment">

	// <component name="M">
	M:time=((-1)*(k1*M)+k2*Mp+gx*AM);

	// <component name="Mp">
	Mp:time=(gp*AMp+k1*M-(k2+fp)*Mp);

	// <component name="AMp">
	AMp:time=(fp*Mp+k6*AM-(k5+gp)*AMp);

	// <component name="AM">
	AM:time=(k5*AMp-(k6+gx)*AM);

	// <component name="model_parameters">
	k1=(if ((time>(0 second)) and (time<(5 second))) (.35 first_order_rate_constant) else (.06 first_order_rate_constant));
	k6=k1;
	gp2=(4*(fp1+gp1));
	g2=(20*g1);
	g3=(3*g1);
	gp3=(3*gp1);
}