Supermarket refrigeration system (Muscod)

The printable version is no longer supported and may have rendering errors. Please update your browser bookmarks and please use the default browser print function instead.

The differential equations for the Supermarket refrigeration system in C code read as follows

// number of compressors
#define NVALVES 2

// constants
#define M_GOODS 200.0
#define C_P_GOODS 1000.0
#define UA_GOODS_AIR 300.0
#define M_WALL 260.0
#define C_P_WALL 385.0
#define UA_AIR_WALL 500.0
#define M_AIR 50.0
#define C_P_AIR 1000.0
#define UA_WALL_REF_MAX 4000.0
#define M_REF_MAX 1.0
#define TAU_FILL 40.0
#define T_SH 10.0
#define V_SUC 5.0
#define V_SL 0.08 // 2 display cases - 2 compressors
// #define V_SL 0.095 // 3 display cases - 3 compressors
#define ETA_VOL 0.81

// disturbances - day scenario
#define Q_AIRLOAD 3000.0
#define M_REF_CONST 0.2

// disturbances - night scenario
// #define Q_AIRLOAD 1800.0
// #define M_REF_CONST 0.0

double delta_h = (0.0217*xd[0]*xd[0] - 0.1704*xd[0] + 2.2988)*1e5;
double T_e = -4.3544*xd[0]*xd[0] + 29.224*xd[0] - 51.2005;
double rho_suc = 4.6073*xd[0] + 0.3798;
double rho_suc__P_suc = -0.0329*xd[0]*xd[0]*xd[0] + 0.2161*xd[0]*xd[0] - 0.4742*xd[0] + 5.4817;
double f = (0.0265*xd[0]*xd[0]*xd[0] -0.4346*xd[0]*xd[0] + 2.4923*xd[0] + 1.2189)*1e5;

double Q_goods_air[NVALVES];
double Q_air_wall[NVALVES];
double UA_wall_ref[NVALVES];
double Q_e[NVALVES];
double m[NVALVES];

double m_in_suc = 0.0;

int i;

for (i=0; i<NVALVES; i++){
    Q_goods_air[i] = UA_GOODS_AIR*(xd[1 + i*4] - xd[3 + i*4]);
    Q_air_wall[i] = UA_AIR_WALL*(xd[3 + i*4] - xd[2 + i*4]);
    UA_wall_ref[i] = UA_WALL_REF_MAX * xd[4 + 4*i]/M_REF_MAX;
    Q_e[i] = UA_wall_ref[i]*(xd[2 + 4*i] - T_e);

    m[i] = Q_e[i]/delta_h;
    m_in_suc += m[i];
}

double V_comp = 0.0;
double comp_scale = (double) 1.0/NCOMPS;
V_comp = comp_scale*u[NVALVES]*ETA_VOL*V_SL;

// suction pressure
rhs[0] = (m_in_suc + M_REF_CONST - V_comp*rho_suc) / (V_SUC*rho_suc__P_suc);

// for each display/valve
for (i=0; i<NVALVES; i++){

    // temperatures:

    // goods
    rhs[1 + i*4] = - Q_goods_air[i]/(M_GOODS*C_P_GOODS);

    // wall
    rhs[2 + i*4] = (Q_air_wall[i] - Q_e[i])/(M_WALL*C_P_WALL);

    // air
    rhs[3 + i*4] = ((Q_goods_air[i] + Q_AIRLOAD - Q_air_wall[i]) /(M_AIR*C_P_AIR));

    // mass of liquefied refrigerant:
    rhs[4 + i*4] = ((M_REF_MAX - xd[4 + 4*i])/TAU_FILL * u[i] - m[i] * (1 - u[i]));
}