Supermarket refrigeration system (JModelica)

This page contains the source code to solve the Supermarket refrigeration system problem with JModelica. The automatic differentiation tool CasADI and the solver IPOPT were used to solve the problem. As this problem is numerically challenging, the following implementation does not result in a fully converged optimal solution.

Model file (supermarket_opt.mop)

//-------------------------------------------------------------------------
//Supermarket refrigeration system with direct collocation using JModelica
//(c) Madeleine Schroter
//--------------------------------------------------------------------------
 
 
optimization supermarket_opt (objective = cost(finalTime),
                         startTime = 0,
                         finalTime=1)
 
    // Parameters
    parameter Real Q_airload = 3000.00;             
    parameter Real m_ref = 0.20;             
    parameter Real M_goods=200.00;
    parameter Real Cp_goods=1000.00;
    parameter Real UA_goods=300.00;
    parameter Real M_wall=260.00;
    parameter Real Cp_wall=385.00;
    parameter Real UA_air=500.00;
    parameter Real M_air=50.00;
    parameter Real Cp_air=1000.00;
    parameter Real UA_wall=4000.00;
    parameter Real t_fill=40.00;
    parameter Real T_sh=10.00;
    parameter Real M_ref=1.00;
    parameter Real V_suc=5.00;
    parameter Real V_sl=0.08;             
    parameter Real eta_vol=0.81;
    parameter Real tf(free=true,min=650,max=750)=650;
 
   // Differential states
    Real x0;
    Real x1;
    Real x2;
    Real x3;
    Real x4;
    Real x5;
    Real x6;
    Real x7;
    Real x8;
    Real cost(start=0);
 
    // The control signal
    input Real u0(free=true,min=0,max=1,initialGuess=0.5);
    input Real u1(free=true,min=0,max=1,initialGuess=0.5);
    input Real u2(free=true,min=0,max=1,initialGuess=1);
 
  equation
    der(x0) = ((1.0/(V_suc*(-0.0329*x0*x0*x0+0.2161*x0*x0-0.4742*x0+5.4817)))*((UA_wall/(M_ref*((0.0217*x0*x0-0.1704*x0+2.2988)*10^5)))*(x4*(x2-(-4.3544*x0*x0+29.224*x0-51.2005))+x8*(x6-(-4.3544*x0*x0+29.224*x0-51.2005)))+M_ref-eta_vol*V_sl*0.5*(u2)*(4.6073*x0+0.3798)))*tf;
    der(x1) =(-UA_goods*(x1-x3)/(M_goods*Cp_goods))*tf;
    der(x2)=((UA_air*(x3-x2)-(UA_wall/M_ref)*x4*(x2-(-4.3544*x0*x0+29.224*x0-51.2005)))/(M_wall*Cp_wall))*tf;
    der(x3)=((UA_goods*(x1-x3)+Q_airload-UA_air*(x3-x2))/(M_air*Cp_air))*tf;
    der(x4)=(((M_ref-x4)/t_fill)*u0-((UA_wall/(M_ref*((0.0217*x0*x0-0.1704*x0+2.2988)*10^5)))*x4*(x2-(-4.3544*x0*x0+29.224*x0-51.2005)))*(1-u0))*tf;
    der(x5)=(-UA_goods*(x5-x7)/(M_goods*Cp_goods))*tf;
    der(x6)=((UA_air*(x7-x6)-(UA_wall/M_ref)*x8*(x6-(-4.3544*x0*x0+29.224*x0-51.2005)))/(M_wall*Cp_wall))*tf;
    der(x7)=((UA_goods*(x5-x7)+Q_airload-UA_air*(x7-x6))/(M_air*Cp_air))*tf;
    der(x8)=(((M_ref-x8)/t_fill)*u1-((UA_wall/(M_ref*((0.0217*x0*x0-0.1704*x0+2.2988)*10^5)))*x8*(x6-(-4.3544*x0*x0+29.224*x0-51.2005)))*(1-u1))*tf;
    der(cost) = ((-1.0/(((1/*time*/)*tf)))*cost+(1.0/(((1/*time*/))*tf))*((u2)*0.5*eta_vol*V_sl*((0.0265*x0*x0*x0-0.4346*x0*x0+2.4923*x0+1.2189)*10^5)))*tf;
//    der(cost) = 1.0 / 650 * (((u2)*0.5*eta_vol*V_sl*((0.0265*x0*x0*x0-0.4346*x0*x0+2.4923*x0+1.2189)*10^5)))*tf;
 
  constraint
    x3>=2.0;
    x3<=5.0;
    x7>=2.0;
    x7<=5.0;
    x0<=1.7;
    cost(startTime)=0.0;
    x0(finalTime)=x0(startTime);
    x1(finalTime)=x1(startTime);
    x2(finalTime)=x2(startTime);
    x3(finalTime)=x3(startTime);
    x4(finalTime)=x4(startTime);
    x5(finalTime)=x5(startTime);
    x6(finalTime)=x6(startTime);
    x7(finalTime)=x7(startTime);
    x8(finalTime)=x8(startTime);
 
  end supermarket_opt;

Run file

#Import the function for transfering a model to CasADiInterface
from pyjmi import transfer_optimization_problem
 
op=transfer_optimization_problem("supermarket_opt", "supermarket_opt.mop")
 
res=op.optimize()