import nsrunit;
unit conversion on;
unit s=1 second^1;
unit um=1E-6 meter^1;
unit per_um=1E6 meter^(-1);
unit uM=1E-3 meter^(-3)*mole^1;
unit um3=1E-18 meter^3;
unit uM_per_s=1E-3 meter^(-3)*second^(-1)*mole^1;
unit per_uM_per_um3=1E21 mole^(-1);
unit s=1 second^1;
unit uM=1E-3 meter^(-3)*mole^1;
unit per_s=1 second^(-1);
unit s=1 second^1;
unit uM=1E-3 meter^(-3)*mole^1;
unit um2=1E-12 meter^2;
unit per_um2=1E12 meter^(-2);
unit per_um2_per_s=1E12 meter^(-2)*second^(-1);
unit per_s=1 second^(-1);
unit per_uM_per_s=1E3 meter^3*second^(-1)*mole^(-1);
unit um2_per_s=1E-12 meter^2*second^(-1);
unit s=1 second^1;
unit um2=1E-12 meter^2;
unit per_um2=1E12 meter^(-2);
unit per_um2_per_s=1E12 meter^(-2)*second^(-1);
unit per_s=1 second^(-1);
unit s=1 second^1;
unit uM_per_s=1E-3 meter^(-3)*second^(-1)*mole^1;
unit uM=1E-3 meter^(-3)*mole^1;
unit s=1 second^1;
unit um2=1E-12 meter^2;
unit per_um2=1E12 meter^(-2);
unit per_um2_per_s=1E12 meter^(-2)*second^(-1);
unit per_uM_per_s=1E3 meter^3*second^(-1)*mole^(-1);
unit um2_per_s=1E-12 meter^2*second^(-1);
unit per_s=1 second^(-1);
unit uM=1E-3 meter^(-3)*mole^1;
unit s=1 second^1;
unit um2=1E-12 meter^2;
unit per_um2=1E12 meter^(-2);
unit uM_um2=1E-15 meter^(-1)*mole^1;
unit per_um2_per_s=1E12 meter^(-2)*second^(-1);
unit uM=1E-3 meter^(-3)*mole^1;
unit per_s=1 second^(-1);
unit s=1 second^1;
unit uM=1E-3 meter^(-3)*mole^1;
unit uM_per_s=1E-3 meter^(-3)*second^(-1)*mole^1;
unit per_s=1 second^(-1);
unit um2=1E-12 meter^2;
unit uM_um2=1E-15 meter^(-1)*mole^1;
unit per_um2=1E12 meter^(-2);
unit per_um2_per_s=1E12 meter^(-2)*second^(-1);

math main {
	realDomain t s;
	t.min=0;
	extern t.max;
	extern t.delta;
	real Cpc uM_um2;
	real Cc uM;
	real Cp per_um2;
	real Vc um3;
	Vc=2550;
	real Rpc per_um;
	Rpc=4.61;
	real gpcrJ2(t) per_um2_per_s;
	real gpcrJ3(t) per_um2_per_s;
	real plcJ5(t) per_um2_per_s;
	real plcJ6(t) per_um2_per_s;
	real J_gain_Gd(t) per_um2_per_s;
	real gpcrJ6(t) per_um2_per_s;
	real plcJ2(t) per_um2_per_s;
	real plcJ3(t) per_um2_per_s;
	real J_gain_Gt(t) per_um2_per_s;
	real J_gain_Ca(t) uM_per_s;
	real plcJ1(t) per_um2_per_s;
	real plcJ4(t) per_um2_per_s;
	real J_gain_IP3(t) uM_per_s;
	real pip2J1(t) per_um2_per_s;
	real pip2J2(t) per_um2_per_s;
	real IP3degJ(t) uM_per_s;
	real L(t) uM;
	real Ls uM;
	Ls=0.100;
	real ts s;
	ts=10;
	real Gd(t) per_um2;
	when(t=t.min) Gd=10000;
	real GPCR_Cycle.J1(t) per_um2_per_s;
	real GPCR_Cycle.kf1 per_uM_per_s;
	GPCR_Cycle.kf1=0.0003;
	real GPCR_Cycle.kr1 per_s;
	real Kd1 uM;
	Kd1=3.00E-5;
	real GPCR_Cycle.kf2 um2_per_s;
	GPCR_Cycle.kf2=2.75E-4;
	real GPCR_Cycle.kr2 per_s;
	real Kd2 per_um2;
	Kd2=27500;
	real GPCR_Cycle.kf3 um2_per_s;
	GPCR_Cycle.kf3=1.00;
	real GPCR_Cycle.kr3 per_s;
	GPCR_Cycle.kr3=0.00100;
	real GPCR_Cycle.J4(t) per_um2_per_s;
	real GPCR_Cycle.kf4 per_uM_per_s;
	GPCR_Cycle.kf4=0.3;
	real GPCR_Cycle.kr4 per_s;
	real GPCR_Cycle.Kd4 uM;
	GPCR_Cycle.Kd4=3.00E-5;
	real GPCR_Cycle.J5(t) per_um2_per_s;
	real GPCR_Cycle.kf5 per_s;
	GPCR_Cycle.kf5=0.0004;
	real GPCR_Cycle.kf6 per_s;
	GPCR_Cycle.kf6=1.00;
	real R(t) per_um2;
	when(t=t.min) R=13.9;
	real Rl(t) per_um2;
	when(t=t.min) Rl=0.00;
	real Rg(t) per_um2;
	when(t=t.min) Rg=5.06;
	real Rlg(t) per_um2;
	when(t=t.min) Rlg=0.00;
	real Rlgp(t) per_um2;
	when(t=t.min) Rlgp=0.00;
	real Gt(t) per_um2;
	when(t=t.min) Gt=0.00;
	real Galpha.J1(t) per_um2_per_s;
	real Galpha.kf1 per_s;
	Galpha.kf1=0.150;
	real Ca(t) uM;
	when(t=t.min) Ca=0.100;
	real Pc(t) per_um2;
	when(t=t.min) Pc=9.09;
	real Pcg(t) per_um2;
	when(t=t.min) Pcg=0.00;
	real PLC_Cycle.kf1 per_uM_per_s;
	PLC_Cycle.kf1=0.0167;
	real PLC_Cycle.kr1 per_s;
	PLC_Cycle.kr1=0.0167;
	real PLC_Cycle.kf2 um2_per_s;
	PLC_Cycle.kf2=0.00420;
	real PLC_Cycle.kr2 per_s;
	PLC_Cycle.kr2=1.00;
	real PLC_Cycle.kf3 um2_per_s;
	PLC_Cycle.kf3=0.0420;
	real PLC_Cycle.kr3 per_s;
	PLC_Cycle.kr3=1.00;
	real PLC_Cycle.kf4 per_uM_per_s;
	PLC_Cycle.kf4=0.0334;
	real PLC_Cycle.kr4 per_s;
	real PLC_Cycle.Kd4 uM;
	PLC_Cycle.Kd4=0.1;
	real PLC_Cycle.kf5 per_s;
	PLC_Cycle.kf5=6.00;
	real PLC_Cycle.kf6 per_s;
	PLC_Cycle.kf6=6.00;
	real P(t) per_um2;
	when(t=t.min) P=90.9;
	real Pg(t) per_um2;
	when(t=t.min) Pg=0.00;
	real PIP2 per_um2;
	PIP2=4000;
	real PIP2.kf1 per_s;
	PIP2.kf1=0.444;
	real Km1 uM;
	Km1=19.8;
	real PIP2.kf2 per_s;
	PIP2.kf2=3.8;
	real Km2 uM;
	Km2=5.00;
	real IP3(t) uM;
	when(t=t.min) IP3=0.0150;
	real IP3_degradation.kf1 per_s;
	IP3_degradation.kf1=1.25;

	// <component name="environment">

	// <component name="geometry">
	Cc=(1/(Vc*(602.2 per_uM_per_um3)));
	Cp=(1/(Vc*Rpc));
	Cpc=(Cc/Cp);

	// <component name="Galpha_interface">
	J_gain_Gd=(plcJ5+plcJ6-(gpcrJ2+gpcrJ3));
	J_gain_Gt=(gpcrJ6-(plcJ2+plcJ3));

	// <component name="calcium_interface">
	J_gain_Ca=(Cpc*((-1)*(plcJ1+plcJ4)));

	// <component name="IP3_interface">
	J_gain_IP3=(Cpc*(pip2J1+pip2J2)-IP3degJ);

	// <component name="ligand">
	L=(if ((t<(ts+(.15 s))) and (t>=ts)) Ls/(1+exp(((-80 per_s))*(t-ts-(.05 s)))) else if (t>=(ts+(.15 s))) Ls else (0 uM));

	// <component name="GPCR_Cycle">
	GPCR_Cycle.kr1=(GPCR_Cycle.kf1*Kd1);
	GPCR_Cycle.J1=(GPCR_Cycle.kf1*R*L-GPCR_Cycle.kr1*Rl);
	GPCR_Cycle.kr2=(GPCR_Cycle.kf2*Kd2);
	gpcrJ2=(GPCR_Cycle.kf2*R*Gd-GPCR_Cycle.kr2*Rg);
	R:t=((-1)*GPCR_Cycle.J1-gpcrJ2);
	gpcrJ3=(GPCR_Cycle.kf3*Rl*Gd-GPCR_Cycle.kr3*Rlg);
	Rl:t=(GPCR_Cycle.J1+gpcrJ6-gpcrJ3);
	GPCR_Cycle.kr4=(GPCR_Cycle.kf4*GPCR_Cycle.Kd4);
	GPCR_Cycle.J4=(GPCR_Cycle.kf4*L*Rg-GPCR_Cycle.kr4*Rlg);
	Rg:t=(gpcrJ2-GPCR_Cycle.J4);
	GPCR_Cycle.J5=(GPCR_Cycle.kf5*Rlg);
	Rlgp:t=GPCR_Cycle.J5;
	gpcrJ6=(GPCR_Cycle.kf6*Rlg);
	Rlg:t=(gpcrJ3-GPCR_Cycle.J5+(GPCR_Cycle.J4-gpcrJ6));

	// <component name="Galpha">
	Gd:t=(Galpha.J1+J_gain_Gd);
	Gt:t=(J_gain_Gt-Galpha.J1);
	Galpha.J1=(Galpha.kf1*Gt);

	// <component name="calcium">
	Ca:t=J_gain_Ca;

	// <component name="PLC_Cycle">
	plcJ2=(PLC_Cycle.kf2*P*Gt-PLC_Cycle.kr2*Pg);
	plcJ1=(PLC_Cycle.kf1*P*Ca-PLC_Cycle.kr1*Pc);
	plcJ3=(PLC_Cycle.kf3*Pc*Gt-PLC_Cycle.kr3*Pcg);
	PLC_Cycle.kr4=(PLC_Cycle.kf4*PLC_Cycle.Kd4);
	plcJ4=(PLC_Cycle.kf4*Pg*Ca-PLC_Cycle.kr4*Pcg);
	plcJ5=(PLC_Cycle.kf5*Pcg);
	plcJ6=(PLC_Cycle.kf6*Pg);
	P:t=(plcJ6-(plcJ2+plcJ1));
	Pg:t=(plcJ2-(plcJ4+plcJ6));
	Pc:t=(plcJ1+plcJ5-plcJ3);
	Pcg:t=(plcJ3+plcJ4-plcJ5);

	// <component name="PIP2">
	pip2J1=(PIP2.kf1*Pc*PIP2/(Km1/Cpc+PIP2));
	pip2J2=(PIP2.kf2*Pcg*PIP2/(Km2/Cpc+PIP2));

	// <component name="IP3">
	IP3:t=J_gain_IP3;

	// <component name="IP3_degradation">
	IP3degJ=(IP3_degradation.kf1*IP3);
}