CPLEX约束规划模型长时间运行无结果，请求问题排查

我用CPLEX的约束规划（Constraint Programming）模块运行模型时，耗时极久始终出不了结果，只能终止运行。哪怕只用少量样本数据，模型跑几个小时还是没结果，怀疑模型本身有问题。下面附上mod文件、dat文件及相关数据信息，麻烦帮忙排查问题，谢谢！

Mod文件内容 #

using CP;
int Z_val = 8;
int P_val = 5;
int T_val = 20;
int I_val = 4;
int J_val = 5;
int K_val = ...;
int N = ...;// Planning horizon (day);
float alpha = ...;// The average loading rate of delivery trucks;
float gamma = ...;// The maximum of containers waiting in the queue at each block.
int BigM = 1000;
int fact_K_minus_1 = ...;

range Z = 1 .. Z_val;// Number of vessels
range P = 1 .. P_val;// Number of appointment periods
range T = 1 .. T_val;// Number of time intervals
range I = 1 .. I_val;// Number of gate lanes
range J = 1 .. J_val;// Number of yard blocks
range K = 1 .. K_val;// Number of RTGCs

float T_z_A[Z] = ...;
float T_z_D[Z] = ...;
float V_z[Z] = ...;
float beta_zj[Z][J] = ...;

float C_s = ...;
range KSet = 0 .. K_val - 1;
float T_l = ...;
float T_k = ...;

float Y_j[J] = ...;
float c_o = ...;
float c_R = ...;
float factorial[KSet] = ...;

// Define sets
{int} Z_j[J];

execute {

  // Iterate over j and z to update Z sets
  for ( var j in J) {
    for ( var z in Z) {
      if (beta_zj[z][j] > 0) {
        Z_j[j].add(z); // Add z to the set Z_j[j]
      }
    }
  }

  // Print the result (optional)
  for ( var j in J) {
    writeln("Z_j[", j, "] = ", Z_j[j]);
  }
}


int scale = 1000;
dvar int pS[Z] ;// Starting appointment period for delivery trucks related to vessel z
dvar int pE[Z] ;


// Decision variables
dexpr float p_z_S[z in Z] = pS[z] / scale;// Starting appointment period for delivery trucks related to vessel z
dexpr float p_z_E[z in Z] = pE[z] / scale;// Ending appointment period for delivery trucks related to vessel z

// Deived variables
float d = ceil ( T_val / P_val );// Number of time intervals in one appointment period
int m = ftoi ( d );
execute {
  writeln(m);
}


dvar int lam_zp[Z][P]  ;// Appointment quota of export containers related to vessel z at appointment period p
dvar int lit_it[I][T]  ;// Average number of trucks waiting in queue at gate lane i at interval t
dvar int dit_it[I][T]  ;// Actual discharge rate of gate lane i at interval t (truck/min)
dvar int pit_it [I][T]  ;// Capacity utilization rate of gate lane i at interval t
dvar int wp_p [P]  ;// Average waiting time of trucks at terminal gate during appointment period p (min)
dvar int wg_g  ;// Average waiting time of trucks at terminal gate during the planning horizon (min)
dvar int djkt_jkt [J][K][T] ;// Discharge rate of RTGC k deployed to block j at interval t (natural container/min)
dvar int pjkt_jkt [J][K][T] ;// Average utilization rate of RTGC k deployed to block j at interval t
dvar int wjp_jp [J][P] ;// Average waiting time of trucks at block j in appointment period p (min)
dvar int wy_y ;// Average waiting time of trucks at yard in the planning horizon (min)


dexpr float lambda_zp[z in Z][p in P] = lam_zp[z][p] / scale;// Appointment quota of export containers related to vessel z at appointment period p
dexpr float l_it_g[i in I][t in T] = lit_it[i][t] / scale;// Average number of trucks waiting in queue at gate lane i at interval t
dexpr float d_it_g[i in I][t in T] = dit_it[i][t] / scale;// Actual discharge rate of gate lane i at interval t (truck/min)
dexpr float p_it_g[i in I][t in T] = pit_it [i][t] / scale;// Capacity utilization rate of gate lane i at interval t
dexpr float w_p_g[p in P] = wp_p[p] / scale;// Average waiting time of trucks at terminal gate during appointment period p (min)
dexpr float w_g = wg_g/ scale;// Average waiting time of trucks at terminal gate during the planning horizon (min)
dexpr float d_jkt_y[j in J][k in K][t in T] = djkt_jkt[j][k][t] / scale;// Discharge rate of RTGC k deployed to block j at interval t (natural container/min)
dexpr float p_jkt_y[j in J][k in K][t in T] = pjkt_jkt [j][k][t] / scale;// Average utilization rate of RTGC k deployed to block j at interval t
dexpr float w_jp_y[j in J][p in P] = wjp_jp [j][p] / scale;// Average waiting time of trucks at block j in appointment period p (min)
dexpr float w_y = wy_y / scale;// Average waiting time of trucks at yard in the planning horizon (min)
dvar boolean d_zp[Z][P];// Binary variable, 1 if vessel z has departed at appointment period p, 0 otherwise

dvar int  lambda_zt_g[Z][T] ;// Number of trucks related to vessel z arriving at terminal gate at interval t
dvar int lambda_it_g[I][T] ;// Number of trucks arriving at gate lane i at interval t

dvar int lambda_t_y[ T] ;// Number of export containers arriving at yard at interval t
dvar int lambda_jt_y[j in J][T] ;// Number of export containers arriving at block j at interval t
dvar int l_jt_y[j in J][T] ;// Number of export containers waiting at block j at interval t

// Objective function
minimize
  ( sum ( t in T ) sum ( i in I ) l_it_g[i][t] + sum ( t in T ) sum ( j in J )
     l_jt_y[j][t] ) * ( 24 * N * c_o / T_val ) + sum ( j in J ) sum ( k in K )
     sum ( t in T ) ( 1 - p_jkt_y[j][k][t] ) * 24 * N * c_R / T_val;
subject to {
  // Constraints 1


    forall ( z in Z )//(1)
    ( p_z_E[z] - p_z_S[z] + 1 ) * 24 * N / P_val >= T_l;

    forall ( z in Z )//(2)
    ( p_z_E[z] - p_z_S[z] + 1 ) * 24 * N / P_val <= T_k;

    forall ( z in Z )//(3)
    p_z_E[z] * 24 * N / P_val <= T_z_A[z];

    forall ( z in Z, p in P )//(4)
    lambda_zp[z][p] * ( p_z_E[z] - p_z_S[z] + 1 ) - V_z[z]<=0.1;

  forall (z in Z, p in P) {
  if ((p - 1) * 24 * N / P_val >= T_z_D[z]) {
    d_zp[z][p] == 1;
  } else {
    d_zp[z][p] == 0;
  }
}

 
    forall ( j in J, p in P )//(6)
    sum ( z in Z_j[j] ) ( sum ( a in 1 .. p ) 
      ( lambda_zp[z][a] - V_z[z] * d_zp[z][a] ) ) * beta_zj[z][j] <= Y_j[j];

  // Constraints at gate
    forall ( z in Z, p in P, t in ( ( p - 1 ) * m + 1 ) .. p * m )//(7)
    lambda_zt_g[z][t] -(lambda_zp[z][p] / ( m * alpha ))<=0.1;

    forall ( i in I, t in T )//(8)
    lambda_it_g[i][t] -( sum ( z in Z ) lambda_zt_g[z][t] / I_val)<=0.1;

    forall ( i in I, t in 1 .. T_val - 1 )//(9)
    l_it_g[i][t + 1] - (maxl ( l_it_g[i][t] + lambda_it_g[i][t] - 
      ( 24 * 60 * N * d_it_g[i][t] / T_val ), 0 ))<=0.1;

    forall ( i in I, t in T )//(10)
    d_it_g[i][t] - ( 19 / 60 ) * p_it_g[i][t]<=0.1;

    forall ( i in I, t in T )//(11)
    //l_it_g[i][t] == p_it_g[i][t] / ( 1 - p_it_g[i][t] );
    l_it_g[i][t] * ( 1 - p_it_g[i][t] ) - p_it_g[i][t] <=0.1;

    forall ( p in P )//(12)
    w_p_g[p] * sum ( t in ( ( p - 1 ) * m + 1 ) .. p * m ) sum ( i in I )
       d_it_g[i][t] - (sum ( t in ( ( p - 1 ) * m + 1 ) .. p * m ) sum 
      ( i in I ) l_it_g[i][t])<=0.1;

  w_g * sum ( t in T ) sum ( i in I ) d_it_g[i][t] - (sum ( t in T ) sum 
  ( i in I ) l_it_g[i][t])<=10;//(13)

  // Constraints at yard
    forall ( t in T )//(14)
    lambda_t_y[t] -( ( 24 * 60 * N / T_val ) * sum ( i in I ) d_it_g[i][t])<=0.1;


    forall ( j in J, t in T )//(15)
    lambda_jt_y[j][t] - (lambda_t_y[t] * alpha * sum ( z in Z_j[j] )
       beta_zj[z][j] * ( lambda_zt_g[z][t] / sum ( i in I ) lambda_it_g[i][t] ))<=0.1;
      

    forall ( j in J, t in 1 .. T_val - 1 )//(16)
    l_jt_y[j][t + 1] - ( maxl ( l_jt_y[j][t] + lambda_jt_y[j][t] - 
      ( 24 * 60 * N / T_val ) * sum ( k in K ) d_jkt_y[j][k][t], 0 ))<=0.1;

  forall ( j in J, t in T, k in K //(17)
  )
    d_jkt_y[j][k][t] - ( 9 / 60 ) * p_jkt_y[j][k][t] <= 0.1;

  forall ( j in J, t in T, k in K //(18)
  )
    l_jt_y[j][t] - (( ( p_jkt_y[j][k][t] * ( 1 + C_s * C_s ) ) / ( 2 * 
      ( K_val - p_jkt_y[j][k][t] ) ) ) * ( 1 + ( fact_K_minus_1 * 
      ( K_val - p_jkt_y[j][k][t] ) ) * ( sum ( n in KSet ) ( 1 / 
      ( factorial[n] * ( p_jkt_y[j][k][t] ^ ( K_val - n ) ) ) ) ) ) ^ ( - 1 )
       + p_jkt_y[j][k][t] )<= 0.1;

    forall ( j in J, t in T )//(19)
    l_jt_y[j][t] <= gamma;

    forall ( j in J, p in P )//(20)
    w_jp_y[j][p] * sum ( t in ( ( p - 1 ) * m + 1 ) .. p * m, k in K )
       d_jkt_y[j][k][t] -( sum ( t in ( ( p - 1 ) * m + 1 ) .. p * m )
       l_jt_y[j][t])<=0.1;

  w_y * sum ( t in T, j in J, k in K ) d_jkt_y[j][k][t] -( sum ( t in T ) sum 
  ( j in J ) l_jt_y[j][t])<=0.1;//(21)



}

Dat文件内容 #

SheetConnection my_sheet("eee.xlsx");
N = 7;
alpha = 1.4;
gamma = 8;
C_s = 0.42687;
T_l= 0.6;
T_k= 2.4;
c_o = 5.728;
c_R = 15.48;

K_val from SheetRead (my_sheet,"Sheet1!B4");
fact_K_minus_1 from SheetRead (my_sheet,"Sheet1!D4");

T_z_A from SheetRead (my_sheet,"Sheet1!B7:B14");
T_z_D from SheetRead (my_sheet,"Sheet1!C7:C14");
V_z from SheetRead (my_sheet,"Sheet1!D7:D14");
Y_j from SheetRead (my_sheet,"Sheet1!L7:L11");
factorial from SheetRead(my_sheet,"Sheet1!I7:I10");

beta_zj from SheetRead (my_sheet,"Sheet1!:B22:I26");

数据说明 #

模型输入数据取自Excel表格对应区域：

• K_val：Sheet1!B4
• fact_K_minus_1：Sheet1!D4
• T_z_A：Sheet1!B7:B14
• T_z_D：Sheet1!C7:C14
• V_z：Sheet1!D7:D14
• Y_j：Sheet1!L7:L11
• factorial：Sheet1!I7:I10
• beta_zj：Sheet1!B22:I26

内容的提问来源于stack exchange，提问作者Diem Nguyen