/*
   There are a total of 10 entries in the algebraic variable array.
   There are a total of 8 entries in each of the rate and state variable arrays.
   There are a total of 29 entries in the constant variable array.
 */
/*
 * VOI is time in component environment (second).
 * STATES[0] is c in component c (micromolar).
 * CONSTANTS[0] is delta in component c (dimensionless).
 * ALGEBRAIC[2] is J_serca in component J_serca (flux).
 * CONSTANTS[28] is J_in in component J_in (flux).
 * ALGEBRAIC[4] is J_pm in component J_pm (flux).
 * ALGEBRAIC[0] is J_IPR in component J_IPR (flux).
 * STATES[1] is ce in component ce (micromolar).
 * CONSTANTS[1] is gamma in component ce (dimensionless).
 * CONSTANTS[2] is p in component p (micromolar).
 * CONSTANTS[3] is kf in component J_IPR (first_order_rate_constant).
 * CONSTANTS[4] is g1 in component J_IPR (first_order_rate_constant).
 * STATES[2] is A in component A (dimensionless).
 * STATES[3] is O in component O (dimensionless).
 * STATES[4] is R in component R (dimensionless).
 * ALGEBRAIC[1] is phi_1 in component IPR_parameters (first_order_rate_constant).
 * ALGEBRAIC[3] is phi_2 in component IPR_parameters (second_order_rate_constant).
 * ALGEBRAIC[5] is phi_2_ in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[5] is k_1_ in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[6] is l_2_ in component IPR_parameters (first_order_rate_constant).
 * STATES[5] is I_1 in component I_1 (dimensionless).
 * STATES[6] is S in component S (dimensionless).
 * ALGEBRAIC[6] is phi_3 in component IPR_parameters (second_order_rate_constant).
 * ALGEBRAIC[7] is phi_4 in component IPR_parameters (first_order_rate_constant).
 * ALGEBRAIC[8] is phi_4_ in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[7] is k_3_ in component IPR_parameters (first_order_rate_constant).
 * STATES[7] is I_2 in component I_2 (dimensionless).
 * ALGEBRAIC[9] is phi_5 in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[8] is k_1 in component IPR_parameters (second_order_rate_constant).
 * CONSTANTS[9] is k_2 in component IPR_parameters (second_order_rate_constant).
 * CONSTANTS[10] is k_2_ in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[11] is k_3 in component IPR_parameters (second_order_rate_constant).
 * CONSTANTS[12] is k_4 in component IPR_parameters (second_order_rate_constant).
 * CONSTANTS[13] is k_4_ in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[14] is l_2 in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[15] is l_4 in component IPR_parameters (second_order_rate_constant).
 * CONSTANTS[16] is l_4_ in component IPR_parameters (second_order_rate_constant).
 * CONSTANTS[17] is l_6 in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[18] is l_6_ in component IPR_parameters (first_order_rate_constant).
 * CONSTANTS[19] is L_1 in component IPR_parameters (micromolar).
 * CONSTANTS[20] is L_3 in component IPR_parameters (micromolar).
 * CONSTANTS[21] is L_5 in component IPR_parameters (micromolar).
 * CONSTANTS[22] is Vs in component J_serca (micromolar2_per_second).
 * CONSTANTS[23] is Ks in component J_serca (micromolar).
 * CONSTANTS[24] is Vp in component J_pm (flux).
 * CONSTANTS[25] is Kp in component J_pm (micromolar).
 * CONSTANTS[26] is alpha1 in component J_in (flux).
 * CONSTANTS[27] is alpha2 in component J_in (first_order_rate_constant).
 * RATES[0] is d/dt c in component c (micromolar).
 * RATES[1] is d/dt ce in component ce (micromolar).
 * RATES[4] is d/dt R in component R (dimensionless).
 * RATES[3] is d/dt O in component O (dimensionless).
 * RATES[5] is d/dt I_1 in component I_1 (dimensionless).
 * RATES[7] is d/dt I_2 in component I_2 (dimensionless).
 * RATES[6] is d/dt S in component S (dimensionless).
 * RATES[2] is d/dt A in component A (dimensionless).
 */
void
initConsts(double* CONSTANTS, double* RATES, double *STATES)
{
STATES[0] = 0.01;
CONSTANTS[0] = 0.1;
STATES[1] = 0.01;
CONSTANTS[1] = 5.4;
CONSTANTS[2] = 10.0;
CONSTANTS[3] = 0.96;
CONSTANTS[4] = 0.002;
STATES[2] = 0.16;
STATES[3] = 0.01;
STATES[4] = 0.16;
CONSTANTS[5] = 0.04;
CONSTANTS[6] = 0.8;
STATES[5] = 0.16;
STATES[6] = 0.16;
CONSTANTS[7] = 29.8;
STATES[7] = 0.16;
CONSTANTS[8] = 0.64;
CONSTANTS[9] = 37.4;
CONSTANTS[10] = 1.4;
CONSTANTS[11] = 0.11;
CONSTANTS[12] = 4.0;
CONSTANTS[13] = 0.54;
CONSTANTS[14] = 1.7;
CONSTANTS[15] = 1.7;
CONSTANTS[16] = 2.5;
CONSTANTS[17] = 4707.0;
CONSTANTS[18] = 11.4;
CONSTANTS[19] = 0.12;
CONSTANTS[20] = 0.025;
CONSTANTS[21] = 54.7;
CONSTANTS[22] = 120.0;
CONSTANTS[23] = 0.18;
CONSTANTS[24] = 28.0;
CONSTANTS[25] = 0.42;
CONSTANTS[26] = 0.03;
CONSTANTS[27] = 0.2;
CONSTANTS[28] = CONSTANTS[26]+ CONSTANTS[27]*CONSTANTS[2];
}
void
computeRates(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC)
{
ALGEBRAIC[1] = ( ( CONSTANTS[8]*CONSTANTS[19]+CONSTANTS[14])*STATES[0])/(CONSTANTS[19]+ STATES[0]*(1.00000+CONSTANTS[19]/CONSTANTS[20]));
RATES[5] =  ALGEBRAIC[1]*STATES[4] -  (CONSTANTS[5]+CONSTANTS[6])*STATES[5];
ALGEBRAIC[2] =  (( CONSTANTS[22]*STATES[0])/(CONSTANTS[23]+STATES[0]))*(1.00000/STATES[1]);
ALGEBRAIC[0] =  ( CONSTANTS[3]*pow( 0.100000*STATES[3]+ 0.900000*STATES[2], 4.00000)+CONSTANTS[4])*(STATES[1] - STATES[0]);
RATES[1] =  CONSTANTS[1]*(ALGEBRAIC[2] - ALGEBRAIC[0]);
ALGEBRAIC[4] = ( CONSTANTS[24]*pow(STATES[0], 2.00000))/(pow(CONSTANTS[25], 2.00000)+pow(STATES[0], 2.00000));
RATES[0] = (ALGEBRAIC[0] - ALGEBRAIC[2])+ CONSTANTS[0]*(CONSTANTS[28] - ALGEBRAIC[4]);
ALGEBRAIC[3] = ( CONSTANTS[9]*CONSTANTS[20]+ CONSTANTS[15]*STATES[0])/(CONSTANTS[20]+ STATES[0]*(1.00000+CONSTANTS[20]/CONSTANTS[19]));
ALGEBRAIC[5] = (CONSTANTS[10]+ CONSTANTS[16]*STATES[0])/(1.00000+STATES[0]/CONSTANTS[21]);
RATES[4] = ( ALGEBRAIC[5]*STATES[3] - ( ALGEBRAIC[3]*CONSTANTS[2]*STATES[4]+ ALGEBRAIC[1]*STATES[4]))+ (CONSTANTS[6]+CONSTANTS[5])*STATES[5];
ALGEBRAIC[6] = ( CONSTANTS[11]*CONSTANTS[21])/(STATES[0]+CONSTANTS[21]);
RATES[6] =  1.00000*ALGEBRAIC[6]*STATES[3] -  CONSTANTS[7]*STATES[6];
ALGEBRAIC[7] = ( ( CONSTANTS[12]*CONSTANTS[21]+CONSTANTS[17])*STATES[0])/(STATES[0]+CONSTANTS[21]);
ALGEBRAIC[8] = ( CONSTANTS[19]*(CONSTANTS[13]+CONSTANTS[18]))/(STATES[0]+CONSTANTS[19]);
RATES[3] = ( ALGEBRAIC[3]*CONSTANTS[2]*STATES[4] -  (ALGEBRAIC[5]+ALGEBRAIC[7]+ 1.00000*ALGEBRAIC[6])*STATES[3])+ ALGEBRAIC[8]*STATES[2]+ CONSTANTS[7]*STATES[6];
ALGEBRAIC[9] = ( ( CONSTANTS[8]*CONSTANTS[19]+CONSTANTS[14])*STATES[0])/(STATES[0]+CONSTANTS[19]);
RATES[7] =  ALGEBRAIC[9]*STATES[2] -  (CONSTANTS[5]+CONSTANTS[6])*STATES[7];
RATES[2] = ( ALGEBRAIC[7]*STATES[3] - ( ALGEBRAIC[8]*STATES[2]+ ALGEBRAIC[9]*STATES[2]))+ (CONSTANTS[5]+CONSTANTS[6])*STATES[7];
}
void
computeVariables(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC)
{
ALGEBRAIC[1] = ( ( CONSTANTS[8]*CONSTANTS[19]+CONSTANTS[14])*STATES[0])/(CONSTANTS[19]+ STATES[0]*(1.00000+CONSTANTS[19]/CONSTANTS[20]));
ALGEBRAIC[2] =  (( CONSTANTS[22]*STATES[0])/(CONSTANTS[23]+STATES[0]))*(1.00000/STATES[1]);
ALGEBRAIC[0] =  ( CONSTANTS[3]*pow( 0.100000*STATES[3]+ 0.900000*STATES[2], 4.00000)+CONSTANTS[4])*(STATES[1] - STATES[0]);
ALGEBRAIC[4] = ( CONSTANTS[24]*pow(STATES[0], 2.00000))/(pow(CONSTANTS[25], 2.00000)+pow(STATES[0], 2.00000));
ALGEBRAIC[3] = ( CONSTANTS[9]*CONSTANTS[20]+ CONSTANTS[15]*STATES[0])/(CONSTANTS[20]+ STATES[0]*(1.00000+CONSTANTS[20]/CONSTANTS[19]));
ALGEBRAIC[5] = (CONSTANTS[10]+ CONSTANTS[16]*STATES[0])/(1.00000+STATES[0]/CONSTANTS[21]);
ALGEBRAIC[6] = ( CONSTANTS[11]*CONSTANTS[21])/(STATES[0]+CONSTANTS[21]);
ALGEBRAIC[7] = ( ( CONSTANTS[12]*CONSTANTS[21]+CONSTANTS[17])*STATES[0])/(STATES[0]+CONSTANTS[21]);
ALGEBRAIC[8] = ( CONSTANTS[19]*(CONSTANTS[13]+CONSTANTS[18]))/(STATES[0]+CONSTANTS[19]);
ALGEBRAIC[9] = ( ( CONSTANTS[8]*CONSTANTS[19]+CONSTANTS[14])*STATES[0])/(STATES[0]+CONSTANTS[19]);
}