/*
   There are a total of 3 entries in the algebraic variable array.
   There are a total of 2 entries in each of the rate and state variable arrays.
   There are a total of 15 entries in the constant variable array.
 */
/*
 * VOI is time in component environment (hour).
 * CONSTANTS[0] is T_a in component model_parameters (celsius).
 * CONSTANTS[1] is T_b in component model_parameters (celsius).
 * CONSTANTS[2] is delta_T in component model_parameters (celsius).
 * CONSTANTS[3] is kinc in component model_parameters (W_per_kg_C2).
 * ALGEBRAIC[1] is M_c in component M_c (W_per_kg).
 * CONSTANTS[4] is t_day in component M_c (hour).
 * CONSTANTS[5] is t_night in component M_c (hour).
 * ALGEBRAIC[0] is tprime in component M_c (second).
 * CONSTANTS[6] is day_length in component M_c (second).
 * CONSTANTS[13] is M_day in component M_day (W_per_kg).
 * CONSTANTS[14] is M_night in component M_night (W_per_kg).
 * STATES[0] is M in component M (W_per_kg).
 * CONSTANTS[7] is km in component M (per_hour).
 * STATES[1] is T in component T (celsius).
 * CONSTANTS[8] is c in component T (kJ_per_kg_C).
 * ALGEBRAIC[2] is k in component k (W_per_kg_C).
 * CONSTANTS[12] is kb in component kb (W_per_kg_C).
 * CONSTANTS[10] is T_day in component T_day (celsius).
 * CONSTANTS[11] is T_night in component T_night (celsius).
 * CONSTANTS[9] is M_b in component kb (W_per_kg).
 * RATES[0] is d/dt M in component M (W_per_kg).
 * RATES[1] is d/dt T in component T (celsius).
 * There are a total of 2 condition variables.
 */
void
initConsts(double* CONSTANTS, double* RATES, double *STATES)
{
CONSTANTS[0] = 21.0;
CONSTANTS[1] = 38.0;
CONSTANTS[2] = 1.57;
CONSTANTS[3] = 0.0258;
CONSTANTS[4] = 17.5;
CONSTANTS[5] = 6.73;
CONSTANTS[6] = 86400;
STATES[0] = 3.5;
CONSTANTS[7] = 1.1375;
STATES[1] = 38.785;
CONSTANTS[8] = 3.47;
CONSTANTS[9] = 3.0;
CONSTANTS[10] = CONSTANTS[1]+CONSTANTS[2]/2.00000;
CONSTANTS[11] = CONSTANTS[1] - CONSTANTS[2]/2.00000;
CONSTANTS[12] = CONSTANTS[9]/(CONSTANTS[1] - CONSTANTS[0]);
CONSTANTS[13] =  (CONSTANTS[12]+ CONSTANTS[3]*(CONSTANTS[10] - CONSTANTS[1]))*(CONSTANTS[10] - CONSTANTS[0]);
CONSTANTS[14] =  (CONSTANTS[12]+ CONSTANTS[3]*(CONSTANTS[11] - CONSTANTS[1]))*(CONSTANTS[11] - CONSTANTS[0]);
RATES[0] = 0.1001;
RATES[1] = 0.1001;
}
void
computeResiduals(double VOI, double* CONSTANTS, double* RATES, double* OLDRATES, double* STATES,
                 double* OLDSTATES, double* ALGEBRAIC, double* CONDVARS)
{
resid[0] = RATES[0] -  - CONSTANTS[7]*(STATES[0] - ALGEBRAIC[1]);
resid[1] = RATES[1] -  pow(CONSTANTS[8], -1.00000)*(STATES[0] -  ALGEBRAIC[2]*(STATES[1] - CONSTANTS[0]));
}
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] =  (int)( VOI*3600.00*1.00000) % (int)(CONSTANTS[6]);
ALGEBRAIC[1] = (CONDVAR[0]>=0.00000&&CONDVAR[1]<0.00000 ? CONSTANTS[14] : CONSTANTS[13]);
ALGEBRAIC[2] = CONSTANTS[12]+ CONSTANTS[3]*(STATES[1] - CONSTANTS[1]);
}
void
getStateInformation(double* SI)
{
SI[0] = 1.0;
SI[1] = 1.0;
}
void
computeRoots(double VOI, double* CONSTANTS, double* RATES, double* OLDRATES, double* STATES,
             double* OLDSTATES, double* ALGEBRAIC, double* CONDVARS)
{
CONDVAR[0] = ALGEBRAIC[0]/3600.00 - CONSTANTS[5];
CONDVAR[1] = ALGEBRAIC[0]/3600.00 - CONSTANTS[4];
}