/*
   There are a total of 13 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 55 entries in the constant variable array.
 */
/*
 * ALGEBRAIC[0] is time in component environment (second).
 * STATES[0] is p in component p (dimensionless).
 * CONSTANTS[30] is v1 in component v1 (dimensionless).
 * ALGEBRAIC[1] is v2 in component v2 (dimensionless).
 * ALGEBRAIC[6] is v7 in component v7 (dimensionless).
 * VOI is tau in component normalised_constants (dimensionless).
 * STATES[1] is a in component a (dimensionless).
 * ALGEBRAIC[2] is v3 in component v3 (dimensionless).
 * CONSTANTS[48] is epsilon1 in component normalised_constants (dimensionless).
 * STATES[2] is c in component c (dimensionless).
 * ALGEBRAIC[3] is v4 in component v4 (dimensionless).
 * CONSTANTS[49] is epsilon2 in component normalised_constants (dimensionless).
 * STATES[3] is k in component k (dimensionless).
 * ALGEBRAIC[4] is v5 in component v5 (dimensionless).
 * ALGEBRAIC[5] is v6 in component v6 (dimensionless).
 * CONSTANTS[50] is epsilon3 in component normalised_constants (dimensionless).
 * STATES[4] is o in component o (dimensionless).
 * ALGEBRAIC[7] is v8 in component v8 (dimensionless).
 * CONSTANTS[51] is epsilon4 in component normalised_constants (dimensionless).
 * STATES[5] is n in component n (dimensionless).
 * ALGEBRAIC[10] is vresp in component vresp (dimensionless).
 * CONSTANTS[52] is epsilon5 in component normalised_constants (dimensionless).
 * STATES[6] is en in component en (dimensionless).
 * ALGEBRAIC[12] is vATP in component vATP (dimensionless).
 * ALGEBRAIC[8] is vANT in component vANT (dimensionless).
 * CONSTANTS[53] is epsilon6 in component normalised_constants (dimensionless).
 * STATES[7] is s in component s (dimensionless).
 * ALGEBRAIC[9] is vleak in component vleak (dimensionless).
 * CONSTANTS[54] is epsilon7 in component normalised_constants (dimensionless).
 * CONSTANTS[31] is beta2 in component normalised_constants (dimensionless).
 * CONSTANTS[32] is beta3 in component normalised_constants (dimensionless).
 * CONSTANTS[33] is beta4 in component normalised_constants (dimensionless).
 * CONSTANTS[34] is beta5 in component normalised_constants (dimensionless).
 * CONSTANTS[35] is delta_6 in component normalised_constants (dimensionless).
 * CONSTANTS[36] is beta6 in component normalised_constants (dimensionless).
 * CONSTANTS[37] is beta7 in component normalised_constants (dimensionless).
 * CONSTANTS[38] is beta8 in component normalised_constants (dimensionless).
 * CONSTANTS[39] is beta_ANT in component normalised_constants (dimensionless).
 * CONSTANTS[40] is beta_leak in component normalised_constants (dimensionless).
 * CONSTANTS[41] is beta_resp in component normalised_constants (dimensionless).
 * CONSTANTS[42] is delta_r1 in component normalised_constants (dimensionless).
 * CONSTANTS[43] is delta_r2 in component normalised_constants (dimensionless).
 * CONSTANTS[44] is beta_ATP in component normalised_constants (dimensionless).
 * CONSTANTS[45] is delta_atp in component normalised_constants (dimensionless).
 * ALGEBRAIC[11] is en_crit in component en_crit (dimensionless).
 * CONSTANTS[46] is Kapp_dash in component normalised_constants (dimensionless).
 * CONSTANTS[47] is delta_crit in component normalised_constants (dimensionless).
 * CONSTANTS[0] is At in component normalised_constants (millimolar).
 * CONSTANTS[1] is Nt in component normalised_constants (millimolar).
 * CONSTANTS[2] is Pyr_bar in component normalised_constants (dimensionless).
 * CONSTANTS[3] is Cit_bar in component normalised_constants (dimensionless).
 * CONSTANTS[4] is AcCoA_bar in component normalised_constants (dimensionless).
 * CONSTANTS[5] is KG_bar in component normalised_constants (dimensionless).
 * CONSTANTS[6] is OAA_bar in component normalised_constants (dimensionless).
 * CONSTANTS[7] is k1 in component normalised_constants (micromolar_per_second).
 * CONSTANTS[8] is k2 in component normalised_constants (second_order_rate_constant).
 * CONSTANTS[9] is k3 in component normalised_constants (second_order_rate_constant).
 * CONSTANTS[10] is k4 in component normalised_constants (second_order_rate_constant).
 * CONSTANTS[11] is k5 in component normalised_constants (third_order_rate_constant).
 * CONSTANTS[12] is k6 in component normalised_constants (first_order_rate_constant).
 * CONSTANTS[13] is k7 in component normalised_constants (second_order_rate_constant).
 * CONSTANTS[14] is k8 in component normalised_constants (first_order_rate_constant).
 * CONSTANTS[15] is kresp in component normalised_constants (millimolar_per_second).
 * CONSTANTS[16] is kATP in component normalised_constants (millimolar_per_second).
 * CONSTANTS[17] is kANT in component normalised_constants (dimensionless).
 * CONSTANTS[18] is kleak in component normalised_constants (molar_per_millivolt_per_second).
 * CONSTANTS[19] is Keq in component normalised_constants (dimensionless).
 * CONSTANTS[20] is K in component normalised_constants (millimolar).
 * CONSTANTS[21] is alpha in component normalised_constants (per_millivolt).
 * CONSTANTS[22] is b in component normalised_constants (per_millimolar).
 * CONSTANTS[23] is delta_psi_m in component normalised_constants (millivolt).
 * CONSTANTS[24] is R in component normalised_constants (joule_per_mole_kelvin).
 * CONSTANTS[25] is T in component normalised_constants (kelvin).
 * CONSTANTS[26] is F in component normalised_constants (coulomb_per_mole).
 * CONSTANTS[27] is C in component normalised_constants (millimolar_per_millivolt).
 * CONSTANTS[28] is Kapp in component normalised_constants (per_millimolar).
 * CONSTANTS[29] is Pi in component normalised_constants (millimolar).
 * RATES[0] is d/dt p in component p (dimensionless).
 * RATES[1] is d/dt a in component a (dimensionless).
 * RATES[2] is d/dt c in component c (dimensionless).
 * RATES[3] is d/dt k in component k (dimensionless).
 * RATES[4] is d/dt o in component o (dimensionless).
 * RATES[5] is d/dt n in component n (dimensionless).
 * RATES[6] is d/dt en in component en (dimensionless).
 * RATES[7] is d/dt s in component s (dimensionless).
 */
void
initConsts(double* CONSTANTS, double* RATES, double *STATES)
{
STATES[0] = 0.154;
STATES[1] = 0.063;
STATES[2] = 0.44;
STATES[3] = 0.225;
STATES[4] = 0.005;
STATES[5] = 0.856;
STATES[6] = 3.536;
STATES[7] = 150.0;
CONSTANTS[0] = 4.160;
CONSTANTS[1] = 1.070;
CONSTANTS[2] = 0.161;
CONSTANTS[3] = 0.460;
CONSTANTS[4] = 0.105;
CONSTANTS[5] = 0.146;
CONSTANTS[6] = 0.004;
CONSTANTS[7] = 38;
CONSTANTS[8] = 152;
CONSTANTS[9] = 57142;
CONSTANTS[10] = 53;
CONSTANTS[11] = 40;
CONSTANTS[12] = 82361;
CONSTANTS[13] = 3.2e-3;
CONSTANTS[14] = 3.6;
CONSTANTS[15] = 2.5;
CONSTANTS[16] = 131.9;
CONSTANTS[17] = 0.1;
CONSTANTS[18] = 0.426;
CONSTANTS[19] = 0.3975;
CONSTANTS[20] = 2;
CONSTANTS[21] = 0.100;
CONSTANTS[22] = 0.004;
CONSTANTS[23] = 150.0;
CONSTANTS[24] = 8.314;
CONSTANTS[25] = 298;
CONSTANTS[26] = 96485;
CONSTANTS[27] = 6.75e-06;
CONSTANTS[28] = 4.4e-6;
CONSTANTS[29] = 2.440;
CONSTANTS[30] = 1.00000;
CONSTANTS[31] =  (CONSTANTS[8]/CONSTANTS[7])*CONSTANTS[1]*CONSTANTS[2];
CONSTANTS[32] =  (CONSTANTS[9]/CONSTANTS[7])*CONSTANTS[6]*CONSTANTS[4];
CONSTANTS[33] =  (CONSTANTS[10]/CONSTANTS[7])*CONSTANTS[1]*CONSTANTS[3];
CONSTANTS[34] =  (CONSTANTS[11]/CONSTANTS[7])*CONSTANTS[1]*CONSTANTS[0]*CONSTANTS[5];
CONSTANTS[35] = CONSTANTS[5]/( CONSTANTS[6]*CONSTANTS[19]);
CONSTANTS[36] =  (CONSTANTS[12]/CONSTANTS[7])*CONSTANTS[6];
CONSTANTS[37] =  (CONSTANTS[13]/CONSTANTS[7])*CONSTANTS[0]*CONSTANTS[2];
CONSTANTS[38] =  (CONSTANTS[14]/CONSTANTS[7])*CONSTANTS[6];
CONSTANTS[39] =  (CONSTANTS[17]/CONSTANTS[7])*CONSTANTS[0];
CONSTANTS[40] =  (CONSTANTS[18]/CONSTANTS[7])*CONSTANTS[23];
CONSTANTS[41] = CONSTANTS[15]/CONSTANTS[7];
CONSTANTS[42] = CONSTANTS[20]/CONSTANTS[1];
CONSTANTS[43] =  CONSTANTS[21]*CONSTANTS[23];
CONSTANTS[44] = CONSTANTS[16]/CONSTANTS[7];
CONSTANTS[45] =  CONSTANTS[22]*CONSTANTS[0];
CONSTANTS[46] =  CONSTANTS[28]*CONSTANTS[29];
CONSTANTS[47] =  3.00000*(( 1.20000*CONSTANTS[26]*CONSTANTS[23])/( CONSTANTS[24]*CONSTANTS[25]));
CONSTANTS[48] = CONSTANTS[4]/CONSTANTS[2];
CONSTANTS[49] = CONSTANTS[3]/CONSTANTS[2];
CONSTANTS[50] = CONSTANTS[5]/CONSTANTS[2];
CONSTANTS[51] = CONSTANTS[6]/CONSTANTS[2];
CONSTANTS[52] = CONSTANTS[1]/CONSTANTS[2];
CONSTANTS[53] = CONSTANTS[0]/CONSTANTS[2];
CONSTANTS[54] =  (CONSTANTS[23]/CONSTANTS[2])*CONSTANTS[27];
}
void
computeRates(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC)
{
ALGEBRAIC[1] =  CONSTANTS[31]*STATES[0]*STATES[5];
ALGEBRAIC[2] =  CONSTANTS[32]*STATES[4]*STATES[1];
RATES[1] = (ALGEBRAIC[1] - ALGEBRAIC[2])/CONSTANTS[48];
ALGEBRAIC[3] =  CONSTANTS[33]*STATES[2]*STATES[5];
RATES[2] = (ALGEBRAIC[2] - ALGEBRAIC[3])/CONSTANTS[49];
ALGEBRAIC[4] =  CONSTANTS[34]*STATES[3]*STATES[5]*(1.00000 - STATES[6]);
ALGEBRAIC[5] =  CONSTANTS[36]*(STATES[4] -  CONSTANTS[35]*STATES[3]);
RATES[3] = ((ALGEBRAIC[3]+ALGEBRAIC[5]) - ALGEBRAIC[4])/CONSTANTS[50];
ALGEBRAIC[6] =  CONSTANTS[37]*STATES[0]*STATES[6];
RATES[0] = CONSTANTS[30] - (ALGEBRAIC[1]+ALGEBRAIC[6]);
ALGEBRAIC[7] =  CONSTANTS[38]*STATES[4];
RATES[4] = ((ALGEBRAIC[4]+ALGEBRAIC[6]) - (ALGEBRAIC[2]+ALGEBRAIC[7]+ALGEBRAIC[5]))/CONSTANTS[51];
ALGEBRAIC[10] =  CONSTANTS[41]*((1.00000 - STATES[5])/((CONSTANTS[42]+1.00000) - STATES[5]))*(1.00000/(1.00000+exp( CONSTANTS[43]*(STATES[7] - 1.00000))));
RATES[5] = (ALGEBRAIC[10] - (ALGEBRAIC[1]+ALGEBRAIC[3]+ 2.00000*ALGEBRAIC[4]))/CONSTANTS[52];
ALGEBRAIC[11] = CONSTANTS[46]/(CONSTANTS[46]+exp( -1.00000*CONSTANTS[47]*STATES[7]));
ALGEBRAIC[12] =  CONSTANTS[44]*(2.00000/(1.00000+exp( CONSTANTS[45]*(STATES[6] -  ALGEBRAIC[11]*STATES[7]))) - 1.00000);
ALGEBRAIC[8] =  CONSTANTS[39]*STATES[6];
RATES[6] = ((ALGEBRAIC[12]+ALGEBRAIC[4]) - (ALGEBRAIC[8]+ALGEBRAIC[6]))/CONSTANTS[53];
ALGEBRAIC[9] =  CONSTANTS[40]*STATES[7];
RATES[7] = ( 10.0000*ALGEBRAIC[10] - ( 3.00000*ALGEBRAIC[12]+ALGEBRAIC[9]+ALGEBRAIC[8]))/CONSTANTS[54];
}
void
computeVariables(double VOI, double* CONSTANTS, double* RATES, double* STATES, double* ALGEBRAIC)
{
ALGEBRAIC[1] =  CONSTANTS[31]*STATES[0]*STATES[5];
ALGEBRAIC[2] =  CONSTANTS[32]*STATES[4]*STATES[1];
ALGEBRAIC[3] =  CONSTANTS[33]*STATES[2]*STATES[5];
ALGEBRAIC[4] =  CONSTANTS[34]*STATES[3]*STATES[5]*(1.00000 - STATES[6]);
ALGEBRAIC[5] =  CONSTANTS[36]*(STATES[4] -  CONSTANTS[35]*STATES[3]);
ALGEBRAIC[6] =  CONSTANTS[37]*STATES[0]*STATES[6];
ALGEBRAIC[7] =  CONSTANTS[38]*STATES[4];
ALGEBRAIC[10] =  CONSTANTS[41]*((1.00000 - STATES[5])/((CONSTANTS[42]+1.00000) - STATES[5]))*(1.00000/(1.00000+exp( CONSTANTS[43]*(STATES[7] - 1.00000))));
ALGEBRAIC[11] = CONSTANTS[46]/(CONSTANTS[46]+exp( -1.00000*CONSTANTS[47]*STATES[7]));
ALGEBRAIC[12] =  CONSTANTS[44]*(2.00000/(1.00000+exp( CONSTANTS[45]*(STATES[6] -  ALGEBRAIC[11]*STATES[7]))) - 1.00000);
ALGEBRAIC[8] =  CONSTANTS[39]*STATES[6];
ALGEBRAIC[9] =  CONSTANTS[40]*STATES[7];
rootfind_0(VOI, CONSTANTS, RATES, STATES, ALGEBRAIC, pret);
}

void objfunc_0(double *p, double *hx, int m, int n, void *adata)
{
  struct rootfind_info* rfi = (struct rootfind_info*)adata;
#define VOI rfi->aVOI
#define CONSTANTS rfi->aCONSTANTS
#define RATES rfi->aRATES
#define STATES rfi->aSTATES
#define ALGEBRAIC rfi->aALGEBRAIC
#define pret rfi->aPRET
  ALGEBRAIC[0] = *p;
  *hx = (VOI) - ( (CONSTANTS[7]/CONSTANTS[2])*ALGEBRAIC[0]);
#undef VOI
#undef CONSTANTS
#undef RATES
#undef STATES
#undef ALGEBRAIC
#undef pret
}
void rootfind_0(double VOI, double* CONSTANTS, double* RATES, 
double* STATES, double* ALGEBRAIC, int* pret)
{
  static double val = 0.1;
  double bp, work[LM_DIF_WORKSZ(1, 1)];
  struct rootfind_info rfi;
  rfi.aVOI = VOI;
  rfi.aCONSTANTS = CONSTANTS;
  rfi.aRATES = RATES;
  rfi.aSTATES = STATES;
  rfi.aALGEBRAIC = ALGEBRAIC;
  rfi.aPRET = pret;
  do_levmar(objfunc_0, &val, &bp, work, pret, 1, &rfi);
  ALGEBRAIC[0] = val;
}