/*
   There are a total of 23 entries in the algebraic variable array.
   There are a total of 29 entries in each of the rate and state variable arrays.
   There are a total of 47 entries in the constant variable array.
 */
/*
 * VOI is time in component environment (millisecond).
 * STATES[0] is Gi in component Gi (millimolar).
 * ALGEBRAIC[2] is Jglut in component Jglut (millimolar_per_millisecond).
 * ALGEBRAIC[3] is Jgk in component Jgk (millimolar_per_millisecond).
 * STATES[1] is Ge in component Ge (millimolar).
 * ALGEBRAIC[0] is G_infinity in component Ge (millimolar).
 * CONSTANTS[0] is tau_G in component Ge (millisecond).
 * CONSTANTS[1] is Gmax in component Ge (millimolar).
 * CONSTANTS[2] is Gmin in component Ge (millimolar).
 * CONSTANTS[3] is SG in component Ge (dimensionless).
 * CONSTANTS[4] is I_ in component Ge (dimensionless).
 * STATES[2] is I in component I (dimensionless).
 * STATES[3] is G6P in component G6P (millimolar).
 * CONSTANTS[5] is kappa in component model_parameters (dimensionless).
 * STATES[10] is JPFK in component JPFK (micromolar_per_millisecond).
 * STATES[4] is FBP in component FBP (micromolar).
 * ALGEBRAIC[4] is JGPDH in component JGPDH (micromolar_per_millisecond).
 * ALGEBRAIC[1] is F6P in component F6P (millimolar).
 * CONSTANTS[6] is Kglut in component Jglut (millimolar).
 * CONSTANTS[7] is Vglut in component Jglut (millimolar_per_millisecond).
 * CONSTANTS[8] is Kgk in component Jgk (millimolar).
 * CONSTANTS[9] is Vgk in component Jgk (millimolar_per_millisecond).
 * CONSTANTS[10] is ngk in component Jgk (dimensionless).
 * CONSTANTS[11] is lambda in component JPFK (dimensionless).
 * CONSTANTS[12] is Vmax in component JPFK (micromolar_per_millisecond).
 * CONSTANTS[46] is w0000 in component w (dimensionless).
 * STATES[11] is w1000 in component w (dimensionless).
 * ALGEBRAIC[5] is w0100 in component w (dimensionless).
 * ALGEBRAIC[6] is w0010 in component w (dimensionless).
 * STATES[12] is w0001 in component w (dimensionless).
 * STATES[13] is w1100 in component w (dimensionless).
 * STATES[14] is w1010 in component w (dimensionless).
 * STATES[15] is w1001 in component w (dimensionless).
 * STATES[16] is w0101 in component w (dimensionless).
 * STATES[17] is w0011 in component w (dimensionless).
 * ALGEBRAIC[7] is w0110 in component w (dimensionless).
 * STATES[18] is w1110 in component w (dimensionless).
 * STATES[19] is w0111 in component w (dimensionless).
 * STATES[20] is w1011 in component w (dimensionless).
 * STATES[21] is w1101 in component w (dimensionless).
 * STATES[22] is w1111 in component w (dimensionless).
 * CONSTANTS[13] is f13 in component w (dimensionless).
 * CONSTANTS[14] is f23 in component w (dimensionless).
 * CONSTANTS[15] is f41 in component w (dimensionless).
 * CONSTANTS[16] is f42 in component w (dimensionless).
 * CONSTANTS[17] is f43 in component w (dimensionless).
 * CONSTANTS[18] is K1 in component w (micromolar).
 * CONSTANTS[19] is K2 in component w (micromolar).
 * CONSTANTS[20] is K3 in component w (micromolar).
 * CONSTANTS[21] is K4 in component w (micromolar).
 * STATES[23] is AMP in component AMP (micromolar).
 * STATES[24] is ATP in component ATP (micromolar).
 * CONSTANTS[22] is Atot in component ATP (micromolar).
 * STATES[5] is ADP in component ADP (micromolar).
 * CONSTANTS[23] is tau_a in component ADP (millisecond).
 * CONSTANTS[24] is r in component ADP (dimensionless).
 * CONSTANTS[25] is r1 in component ADP (micromolar).
 * ALGEBRAIC[8] is gamma in component gamma (dimensionless).
 * STATES[6] is Ca in component Ca (micromolar).
 * CONSTANTS[26] is v_gamma in component gamma (dimensionless).
 * CONSTANTS[27] is k_gamma in component gamma (micromolar_per_millisecond).
 * STATES[7] is v in component membrane (millivolt).
 * CONSTANTS[28] is cm in component membrane (femtofarad).
 * ALGEBRAIC[19] is I_Ca in component I_Ca (picoampere).
 * ALGEBRAIC[9] is I_K in component I_K (picoampere).
 * ALGEBRAIC[20] is I_K_Ca in component I_K_Ca (picoampere).
 * STATES[25] is I_K_ATP in component I_K_ATP (picoampere).
 * CONSTANTS[29] is gK_ in component I_K (picosiemens).
 * CONSTANTS[30] is vK in component model_parameters (millivolt).
 * STATES[8] is n in component n (dimensionless).
 * ALGEBRAIC[11] is n_infinity in component n (dimensionless).
 * ALGEBRAIC[10] is tau_n in component n (millisecond).
 * CONSTANTS[31] is gCa_ in component I_Ca (picosiemens).
 * CONSTANTS[32] is vCa in component model_parameters (millivolt).
 * ALGEBRAIC[12] is m_infinity in component m (dimensionless).
 * CONSTANTS[33] is gkCa_ in component I_K_Ca (picosiemens).
 * ALGEBRAIC[13] is gkCa in component I_K_Ca (picosiemens).
 * CONSTANTS[34] is KD in component I_K_Ca (micromolar).
 * CONSTANTS[35] is gkATP_ in component I_K_ATP (picosiemens).
 * STATES[26] is gkATP in component I_K_ATP (picosiemens).
 * STATES[27] is o_infinity in component I_K_ATP (dimensionless).
 * ALGEBRAIC[14] is MgADP in component I_K_ATP (micromolar).
 * ALGEBRAIC[15] is ADP3_ in component I_K_ATP (micromolar).
 * STATES[28] is ATP4_ in component I_K_ATP (micromolar).
 * CONSTANTS[36] is fcyt in component Ca (dimensionless).
 * ALGEBRAIC[21] is Jmem in component Jmem (micromolar_per_millisecond).
 * ALGEBRAIC[22] is Jer in component Jer (micromolar_per_millisecond).
 * STATES[9] is Caer in component Caer (micromolar).
 * CONSTANTS[37] is fer in component Caer (dimensionless).
 * CONSTANTS[38] is Vcyt_Ver in component Caer (dimensionless).
 * CONSTANTS[39] is kPMCA in component Jmem (first_order_rate_constant).
 * CONSTANTS[40] is alpha in component Jmem (micromolar_per_millisecond).
 * ALGEBRAIC[17] is Jleak in component Jleak (micromolar_per_millisecond).
 * ALGEBRAIC[16] is JSERCA in component JSERCA (micromolar_per_millisecond).
 * CONSTANTS[41] is kSERCA in component JSERCA (first_order_rate_constant).
 * CONSTANTS[42] is pleak in component Jleak (first_order_rate_constant).
 * ALGEBRAIC[18] is I_infinity in component I (dimensionless).
 * CONSTANTS[43] is tau_I in component I (millisecond).
 * CONSTANTS[44] is I_slope in component I (per_micromolar).
 * CONSTANTS[45] is Canull in component I (micromolar).
 * RATES[0] is d/dt Gi in component Gi (millimolar).
 * RATES[1] is d/dt Ge in component Ge (millimolar).
 * RATES[3] is d/dt G6P in component G6P (millimolar).
 * RATES[4] is d/dt FBP in component FBP (micromolar).
 * RATES[5] is d/dt ADP in component ADP (micromolar).
 * RATES[7] is d/dt v in component membrane (millivolt).
 * RATES[8] is d/dt n in component n (dimensionless).
 * RATES[6] is d/dt Ca in component Ca (micromolar).
 * RATES[9] is d/dt Caer in component Caer (micromolar).
 * RATES[2] is d/dt I in component I (dimensionless).
 * There are a total of 2 condition variables.
 */
void
initConsts(double* CONSTANTS, double* RATES, double *STATES)
{
STATES[0] = 0;
STATES[1] = 7;
CONSTANTS[0] = 450000;
CONSTANTS[1] = 15;
CONSTANTS[2] = 1;
CONSTANTS[3] = 1;
CONSTANTS[4] = 5;
STATES[2] = 0;
STATES[3] = 200;
CONSTANTS[5] = 0.005;
STATES[4] = 40;
CONSTANTS[6] = 7;
CONSTANTS[7] = 8;
CONSTANTS[8] = 7;
CONSTANTS[9] = 0.8;
CONSTANTS[10] = 4;
CONSTANTS[11] = 0.06;
CONSTANTS[12] = 2;
CONSTANTS[13] = 0.02;
CONSTANTS[14] = 0.2;
CONSTANTS[15] = 20;
CONSTANTS[16] = 20;
CONSTANTS[17] = 20;
CONSTANTS[18] = 30;
CONSTANTS[19] = 1;
CONSTANTS[20] = 50000;
CONSTANTS[21] = 1000;
CONSTANTS[22] = 3000;
STATES[5] = 780;
CONSTANTS[23] = 300000;
CONSTANTS[24] = 0.5;
CONSTANTS[25] = 0.35;
STATES[6] = 0.1;
CONSTANTS[26] = 2.2;
CONSTANTS[27] = 0.1;
STATES[7] = -60;
CONSTANTS[28] = 5300;
CONSTANTS[29] = 2700;
CONSTANTS[30] = -75;
STATES[8] = 0;
CONSTANTS[31] = 1000;
CONSTANTS[32] = 25;
CONSTANTS[33] = 400;
CONSTANTS[34] = 0.5;
CONSTANTS[35] = 40000;
CONSTANTS[36] = 0.01;
STATES[9] = 185;
CONSTANTS[37] = 0.01;
CONSTANTS[38] = 31;
CONSTANTS[39] = 0.18;
CONSTANTS[40] = 4.5e-6;
CONSTANTS[41] = 0.4;
CONSTANTS[42] = 0.0002;
CONSTANTS[43] = 10000;
CONSTANTS[44] = 210;
CONSTANTS[45] = 0.055;
CONSTANTS[46] = 1.00000;
STATES[10] = 0.1001;
STATES[11] = 0.1001;
STATES[12] = 0.1001;
STATES[13] = 0.1001;
STATES[14] = 0.1001;
STATES[15] = 0.1001;
STATES[16] = 0.1001;
STATES[17] = 0.1001;
STATES[18] = 0.1001;
STATES[19] = 0.1001;
STATES[20] = 0.1001;
STATES[21] = 0.1001;
STATES[22] = 0.1001;
STATES[23] = 0.1001;
STATES[24] = 0.1001;
STATES[25] = 0.1001;
STATES[26] = 0.1001;
STATES[27] = 0.1001;
STATES[28] = 0.1001;
RATES[0] = 0.1001;
RATES[1] = 0.1001;
RATES[3] = 0.1001;
RATES[4] = 0.1001;
RATES[5] = 0.1001;
RATES[7] = 0.1001;
RATES[8] = 0.1001;
RATES[6] = 0.1001;
RATES[9] = 0.1001;
RATES[2] = 0.1001;
}
void
computeResiduals(double VOI, double* CONSTANTS, double* RATES, double* OLDRATES, double* STATES,
                 double* OLDSTATES, double* ALGEBRAIC, double* CONDVARS)
{
resid[0] = RATES[0] - ALGEBRAIC[2] - ALGEBRAIC[3];
resid[1] = RATES[1] - (ALGEBRAIC[0] - STATES[1])/CONSTANTS[0];
resid[2] = RATES[3] -  CONSTANTS[5]*(ALGEBRAIC[3] - STATES[10]);
resid[3] = RATES[4] -  CONSTANTS[5]*(STATES[10] -  0.500000*ALGEBRAIC[4]);
resid[4] = STATES[10] - ( CONSTANTS[12]*( (1.00000 - CONSTANTS[11])*STATES[18]+ CONSTANTS[11]*(ALGEBRAIC[6]+STATES[14]+STATES[17]+STATES[18]+STATES[19]+STATES[20]+ALGEBRAIC[7]+STATES[22])))/(CONSTANTS[46]+STATES[11]+ALGEBRAIC[5]+ALGEBRAIC[6]+STATES[12]+STATES[13]+STATES[14]+STATES[15]+STATES[16]+STATES[17]+ALGEBRAIC[7]+STATES[18]+STATES[19]+STATES[20]+STATES[21]+STATES[22]);
resid[5] = STATES[11] - pow(STATES[23]/CONSTANTS[18], 1.00000);
resid[6] = STATES[12] - pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[7] = STATES[13] -  pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(STATES[4]/CONSTANTS[19], 1.00000);
resid[8] = STATES[14] -  (1.00000/pow(CONSTANTS[13],  1.00000*1.00000))*pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000);
resid[9] = STATES[15] -  (1.00000/pow(CONSTANTS[15],  1.00000*1.00000))*pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[10] = STATES[16] -  (1.00000/pow(CONSTANTS[16],  1.00000*1.00000))*pow(STATES[4]/CONSTANTS[19], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[11] = STATES[17] -  (1.00000/pow(CONSTANTS[17],  1.00000*1.00000))*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[12] = STATES[18] -  (1.00000/( pow(CONSTANTS[13],  1.00000*1.00000)*pow(CONSTANTS[14],  1.00000*1.00000)))*pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(STATES[4]/CONSTANTS[19], 1.00000)*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000);
resid[13] = STATES[19] -  (1.00000/( pow(CONSTANTS[14],  1.00000*1.00000)*pow(CONSTANTS[16],  1.00000*1.00000)*pow(CONSTANTS[17],  1.00000*1.00000)))*pow(STATES[4]/CONSTANTS[19], 1.00000)*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[14] = STATES[20] -  (1.00000/( pow(CONSTANTS[13],  1.00000*1.00000)*pow(CONSTANTS[15],  1.00000*1.00000)*pow(CONSTANTS[17],  1.00000*1.00000)))*pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[15] = STATES[21] -  (1.00000/( pow(CONSTANTS[15],  1.00000*1.00000)*pow(CONSTANTS[16],  1.00000*1.00000)))*pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(STATES[4]/CONSTANTS[19], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[16] = STATES[22] -  (1.00000/( pow(CONSTANTS[13],  1.00000*1.00000)*pow(CONSTANTS[14],  1.00000*1.00000)*pow(CONSTANTS[15],  1.00000*1.00000)*pow(CONSTANTS[16],  1.00000*1.00000)*pow(CONSTANTS[17],  1.00000*1.00000)))*pow(STATES[23]/CONSTANTS[18], 1.00000)*pow(STATES[4]/CONSTANTS[19], 1.00000)*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000)*pow(pow(STATES[24], 2.00000)/CONSTANTS[21], 1.00000);
resid[17] = STATES[24] - CONSTANTS[22] - (STATES[5]+STATES[23]);
resid[18] = STATES[23] - pow(STATES[5], 2.00000)/STATES[24];
resid[19] = RATES[5] -  (1.00000/CONSTANTS[23])*(STATES[24] -  STATES[5]*exp( (CONSTANTS[24]+ALGEBRAIC[8])*(1.00000 - STATES[6]/CONSTANTS[25])));
resid[20] = RATES[7] - - (ALGEBRAIC[9]+ALGEBRAIC[19]+ALGEBRAIC[20]+STATES[25])/CONSTANTS[28];
resid[21] = RATES[8] - (ALGEBRAIC[11] - STATES[8])/ALGEBRAIC[10];
resid[22] = STATES[25] -  STATES[26]*(STATES[7] - CONSTANTS[30]);
resid[23] = STATES[26] -  CONSTANTS[35]*STATES[27];
resid[24] = STATES[27] - ( 0.0800000*(1.00000+( 2.00000*ALGEBRAIC[14])/17.0000)+ 0.890000*pow(ALGEBRAIC[14]/17.0000, 2.00000))/( pow(1.00000+ALGEBRAIC[14]/17.0000, 2.00000)*(1.00000+ALGEBRAIC[15]/26.0000+STATES[28]/1.00000));
resid[25] = STATES[28] -  0.00500000*STATES[24];
resid[26] = RATES[6] -  CONSTANTS[36]*(ALGEBRAIC[21]+ALGEBRAIC[22]);
resid[27] = RATES[9] -  - CONSTANTS[37]*CONSTANTS[38]*ALGEBRAIC[22];
resid[28] = RATES[2] - (ALGEBRAIC[18] - STATES[2])/CONSTANTS[43];
}
void
computeVariables(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC)
{
}
void
computeEssentialVariables(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC)
{
ALGEBRAIC[0] = CONSTANTS[2]+(CONSTANTS[1] - CONSTANTS[2])/(1.00000+exp((STATES[2] - CONSTANTS[4])/CONSTANTS[3]));
ALGEBRAIC[1] =  0.300000*STATES[3];
ALGEBRAIC[2] = ( CONSTANTS[7]*(STATES[1] - STATES[0])*CONSTANTS[6])/( (CONSTANTS[6]+STATES[1])*(CONSTANTS[6]+STATES[0]));
ALGEBRAIC[3] = ( CONSTANTS[9]*pow(STATES[0], CONSTANTS[10]))/(pow(CONSTANTS[8], CONSTANTS[10])+pow(STATES[0], CONSTANTS[10]));
ALGEBRAIC[4] =  0.200000* pow((STATES[4]/1.00000), 1.0 / 2);
ALGEBRAIC[5] = pow(STATES[4]/CONSTANTS[19], 1.00000);
ALGEBRAIC[6] = pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000);
ALGEBRAIC[7] =  (1.00000/pow(CONSTANTS[14],  1.00000*1.00000))*pow(STATES[4]/CONSTANTS[19], 1.00000)*pow(pow(ALGEBRAIC[1], 2.00000)/CONSTANTS[20], 1.00000);
ALGEBRAIC[8] = ( CONSTANTS[26]*ALGEBRAIC[4])/(CONSTANTS[27]+ALGEBRAIC[4]);
ALGEBRAIC[9] =  CONSTANTS[29]*STATES[8]*(STATES[7] - CONSTANTS[30]);
ALGEBRAIC[10] = 1.00000/( 0.0350000*cosh((STATES[7]+16.0000)/22.4000));
ALGEBRAIC[11] = 1.00000/(1.00000+exp(- (STATES[7]+16.0000)/5.60000));
ALGEBRAIC[14] =  0.165000*STATES[5];
ALGEBRAIC[15] =  0.135000*STATES[5];
ALGEBRAIC[18] = (CONDVAR[0]>=0.00000 ?  CONSTANTS[44]*(STATES[6] - CONSTANTS[45]) : CONDVAR[1]<0.00000 ? 0.00000 : 0.0/0.0);
ALGEBRAIC[12] = 1.00000/(1.00000+exp(- (STATES[7]+20.0000)/12.0000));
ALGEBRAIC[19] =  CONSTANTS[31]*ALGEBRAIC[12]*(STATES[7] - CONSTANTS[32]);
ALGEBRAIC[13] = ( CONSTANTS[33]*pow(STATES[6], 2.00000))/(pow(CONSTANTS[34], 2.00000)+pow(STATES[6], 2.00000));
ALGEBRAIC[20] =  ALGEBRAIC[13]*(STATES[7] - CONSTANTS[30]);
ALGEBRAIC[21] = - ( CONSTANTS[40]*ALGEBRAIC[19]+ CONSTANTS[39]*STATES[6]);
ALGEBRAIC[17] =  CONSTANTS[42]*(STATES[9] - STATES[6]);
ALGEBRAIC[16] =  CONSTANTS[41]*STATES[6];
ALGEBRAIC[22] = ALGEBRAIC[17] - ALGEBRAIC[16];
}
void
getStateInformation(double* SI)
{
SI[0] = 1.0;
SI[1] = 1.0;
SI[2] = 1.0;
SI[3] = 1.0;
SI[10] = 0.0;
SI[4] = 1.0;
SI[11] = 0.0;
SI[12] = 0.0;
SI[13] = 0.0;
SI[14] = 0.0;
SI[15] = 0.0;
SI[16] = 0.0;
SI[17] = 0.0;
SI[18] = 0.0;
SI[19] = 0.0;
SI[20] = 0.0;
SI[21] = 0.0;
SI[22] = 0.0;
SI[23] = 0.0;
SI[24] = 0.0;
SI[5] = 1.0;
SI[6] = 1.0;
SI[7] = 1.0;
SI[25] = 0.0;
SI[8] = 1.0;
SI[26] = 0.0;
SI[27] = 0.0;
SI[28] = 0.0;
SI[9] = 1.0;
}
void
computeRoots(double VOI, double* CONSTANTS, double* RATES, double* OLDRATES, double* STATES,
             double* OLDSTATES, double* ALGEBRAIC, double* CONDVARS)
{
CONDVAR[0] = STATES[6] - CONSTANTS[45];
CONDVAR[1] = STATES[6] - CONSTANTS[45];
}