opt100

スケジューリング問題

スケジューリング問題に対する定式化とソルバー

スケジューリング問題

スケジューリング問題(scheduling problem)とは，稀少資源（機械）を諸活動（ジョブ，タスク，仕事，作業，オペレーションなどの総称）へ（時間軸を考慮して）割り振るための方法に対する理論体系である．スケジューリングの応用は，工場内での生産計画，計算機におけるジョブのコントロール，プロジェクトの遂行手順の決定など，様々である．

以下では，スケジューリングに関連する様々な定式化と専用ソルバーについて述べる．

1機械リリース時刻付き重み付き完了時刻和最小化問題

1機械リリース時刻付き重み付き完了時刻和最小化問題 $1|r_j|\sum w_j C_j$ は，以下に定義される問題である．

単一の機械で $n$ 個のジョブを処理する問題を考える．この機械は一度に1つのジョブしか処理できず，ジョブの処理を開始したら途中では中断できないものと仮定する．ジョブの集合の添え字を $j=1,2,\ldots,n$，ジョブの集合を $J=\{1,2,\ldots,n\}$ と書く．各ジョブ$j$ に対する処理時間 $p_j$，重要度を表す重み$w_j$，リリース時刻（ジョブの処理が開始できる最早時刻）$r_j$が与えられたとき，各ジョブ $j$ の処理完了時刻 $C_j$ の重み付きの和を最小にするようなジョブを機械にかける順番（スケジュール）を求めよ．

この問題は古くから多くの研究がなされているスケジューリング問題であり，$NP$-困難であることが知られている．ここでは，離接定式化(disjunctive formulation)と呼ばれる単純な定式化を示す．

ジョブ $j$ の開始時刻を表す連続変数 $s_{j}$ と，ジョブ $j$ が他のジョブ $k$ に先行する（前に処理される）とき $1$ になる $0$-$1$ 整数変数 $x_{jk}$ を用いる．

離接定式化は，非常に大きな数 $M$ を用いると，以下のように書ける．

$$\begin{array}{l l l} minimize & \displaystyle\sum_{j=1}^n w_j s_{j} + \displaystyle\sum_{j=1}^n w_j p_j & \\ s.t. & s_j + p_j -M (1-x_{jk}) \leq s_k & \forall j \neq k \\ & x_{jk} +x _{kj}= 1 & \forall j < k \\ & s_j \geq r_j & \forall j=1,2,\ldots,n \\ & x_{jk} \in \{ 0,1 \} & \forall j \neq k \end{array}$$

目的関数は重み付き完了時刻の和 $\sum w_j (s_j+p_j)$ を展開したものであり，それを最小化している．展開した式の第2項は定数であるので，実際には第1項のみを最小化すれば良い．最初の制約は，$x_{jk}$ が $1$ のとき，ジョブ$j$ の完了時刻 $s_j+p_j$ よりジョブ$k$ の開始時刻 $s_k$ が後になることを表す．（$x_{jk}$ が $0$ のときには，$M$ が大きな数なので制約にはならない．）2番目の制約は，ジョブ$j$ がジョブ$k$ に先行するか，ジョブ$k$ がジョブ$j$ に先行するかの何れかが成立することを表す．これを離接制約(disjunctive constraint)と呼ぶ．これが離接定式化の名前の由来である．3番目の制約は，ジョブ$j$ の開始時刻 $s_j$ がリリース時刻 $r_j$ 以上であることを表す．

大きな数 $M$ は，実際にはなるべく小さな値に設定すべきである．制約ごとに異なる $M$ を設定するならば，$r_j +\sum s_j -s_k$，すべての制約で同じ $M$ を用いるならば $\max (r_j) +\sum s_j$ とすれば十分である．

以下では，上の離接定式化を用いて5ジョブの例題を求解する．

defmake_data(n):"""    Data generator for the one machine scheduling problem.    """p,r,d,w={},{},{},{}J=range(1,n+1)forjinJ:p[j]=random.randint(1,4)w[j]=random.randint(1,3)T=sum(p)forjinJ:r[j]=random.randint(0,5)d[j]=r[j]+random.randint(0,5)returnJ,p,r,d,wdefscheduling_disjunctive(J,p,r,w):"""    scheduling_disjunctive: model for the one machine total weighted completion time problem    Disjunctive optimization model for the one machine total weighted    completion time problem with release times.    Parameters:        - J: set of jobs        - p[j]: processing time of job j        - r[j]: earliest start time of job j        - w[j]: weighted of job j,  the objective is the sum of the weighted completion time    Returns a model, ready to be solved.    """model=Model("scheduling: disjunctive")M=max(r.values())+sum(p.values())# big Ms,x=({},{},)# start time variable, x[j,k] = 1 if job j precedes job k, 0 otherwiseforjinJ:s[j]=model.addVar(lb=r[j],vtype="C",name="s(%s)"%j)forkinJ:ifj!=k:x[j,k]=model.addVar(vtype="B",name="x(%s,%s)"%(j,k))model.update()forjinJ:forkinJ:ifj!=k:model.addConstr(s[j]-s[k]+M*x[j,k]<=(M-p[j]),"Bound(%s,%s)"%(j,k))ifj<k:model.addConstr(x[j,k]+x[k,j]==1,"Disjunctive(%s,%s)"%(j,k))model.setObjective(quicksum(w[j]*s[j]forjinJ),GRB.MINIMIZE)model.update()model.__data=s,xreturnmodel

n=5# number of jobsJ,p,r,d,w=make_data(n)model=scheduling_disjunctive(J,p,r,w)model.optimize()s,x=model.__dataz=model.ObjVal+sum([w[j]*p[j]forjinJ])print("Opt.value by Disjunctive Formulation:",z)seq=[jfor(t,j)insorted([(int(s[j].X+0.5),j)forjins])]print("solition=",seq)

Academic license - for non-commercial use only - expires 2023-06-24Using license file /Users/mikiokubo/gurobi.licGurobi Optimizer version 9.1.1 build v9.1.1rc0 (mac64)Thread count: 8 physical cores, 16 logical processors, using up to 16 threadsOptimize a model with 30 rows, 25 columns and 80 nonzerosModel fingerprint: 0x7f60aa0eVariable types: 5 continuous, 20 integer (20 binary)Coefficient statistics:  Matrix range     [1e+00, 2e+01]  Objective range  [1e+00, 3e+00]  Bounds range     [1e+00, 4e+00]  RHS range        [1e+00, 2e+01]Found heuristic solution: objective 52.0000000Presolve removed 10 rows and 10 columnsPresolve time: 0.00sPresolved: 20 rows, 15 columns, 60 nonzerosVariable types: 5 continuous, 10 integer (10 binary)Root relaxation: objective 1.600000e+01, 6 iterations, 0.00 seconds    Nodes    |    Current Node    |     Objective Bounds      |     Work Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time     0     0   16.00000    0    6   52.00000   16.00000  69.2%     -    0sH    0     0                      42.0000000   16.00000  61.9%     -    0s     0     0   30.49035    0    6   42.00000   30.49035  27.4%     -    0sH    0     0                      37.0000000   30.49035  17.6%     -    0s     0     0   33.60360    0    4   37.00000   33.60360  9.18%     -    0sCutting planes:  Gomory: 3  Cover: 2  Implied bound: 2  MIR: 6  RLT: 1Explored 1 nodes (30 simplex iterations) in 0.02 secondsThread count was 16 (of 16 available processors)Solution count 3: 37 42 52 Optimal solution found (tolerance 1.00e-04)Best objective 3.700000000000e+01, best bound 3.700000000000e+01, gap 0.0000%Opt.value by Disjunctive Formulation: 60.0solition= [5, 1, 3, 4, 2]

1機械総納期遅れ最小化問題

1機械総納期遅れ最小化問題 $1||\sum T_j$ は，以下に定義される問題である．

単一の機械で $n$ 個のジョブを処理する問題を考える．この機械は一度に1つのジョブしか処理できず，ジョブの処理を開始したら途中では中断できないものと仮定する．ジョブの添え字を $1,2,\ldots,n$ と書く．各ジョブ$j$ に対する処理時間$p_j$，納期（ジョブが完了しなければならない時刻）$d_j$，が与えられたとき，各ジョブ $j$ の納期からの遅れの和（総納期遅れ）を最小にするようなジョブを機械にかける順番（スケジュール）を求めよ．

この問題は，$NP$-困難であることが知られている．ここでは，大きな数 $M$ を用いない定式化として線形順序付け定式化(linear ordering formulation)を考える．

ジョブ$j$ の納期遅れを表す連続変数 $T_{j}$と，ジョブ$j$ が他のジョブ$k$ に先行する（前に処理される）とき $1$ になる $0$-$1$ 整数変数 $x_{jk}$ を用いる．

上の記号を用いると，線形順序付け定式化は，以下のように書ける．

$$\begin{array}{l l l} minimize & \displaystyle\sum_{j} w_j T_j & \\s.t. & x_{jk} +x _{kj}= 1 & \forall j < k \\ & x_{jk}+x_{k \ell} +x_{\ell j} \leq 2 & \forall j \neq k \neq \ell \\ & \displaystyle\sum_{j=1}^n \displaystyle\sum_{k \neq j} p_k x_{kj} + p_j \leq d_j +T_j & \forall j=1,2,\ldots,n \\ & x_{jk} \in \{ 0,1 \} & \forall j \neq k \end{array}$$

上の定式化において，目的関数は，各ジョブ $j$ に対する納期遅れ $T_j$ を重み $w_j$ で乗じたものの和を最小化することを表す．最初の制約は，ジョブ$j$ がジョブ $k$ の前に処理されるか，その逆にジョブ$k$ がジョブ$j$ の前に処理されるかのいずれかが成立することを表す．（ここまでは離接定式化と同じである．）2番目の制約は，異なるジョブ$j,k,\ell$ に対して，ジョブ$j$ がジョブ$k$ に先行，ジョブ$k$ がジョブ$\ell$ に先行，ジョブ$\ell$ がジョブ$j$ に先行の3つが同時に成り立つことがないことを表す．この2つの制約によって，ジョブ上に線形順序（ジョブが一列に並べられること）が規定される．これが線形順序付け定式化の名前の由来である．3番目の制約は，変数 $x$ で規定される順序付けと納期遅れを表す変数 $T$ を繋ぐための制約であり，ジョブ$j$ の完了時刻が納期 $d_j$ より超えている時間が納期遅れ $T_j$ であることを規定する．

以下では，上の線形順序定式化を用いて5ジョブの例題を求解する．

defscheduling_linear_ordering(J,p,d,w):"""    scheduling_linear_ordering: model for the one machine total weighted tardiness problem    Model for the one machine total weighted tardiness problem    using the linear ordering formulation    Parameters:        - J: set of jobs        - p[j]: processing time of job j        - d[j]: latest non-tardy time for job j        - w[j]: weighted of job j,  the objective is the sum of the weighted completion time    Returns a model, ready to be solved.    """model=Model("scheduling: linear ordering")T,x={},{}# tardiness variable,  x[j,k] =1 if job j precedes job k, =0 otherwiseforjinJ:T[j]=model.addVar(vtype="C",name="T(%s)"%(j))forkinJ:ifj!=k:x[j,k]=model.addVar(vtype="B",name="x(%s,%s)"%(j,k))model.update()forjinJ:model.addConstr(quicksum(p[k]*x[k,j]forkinJifk!=j)-T[j]<=d[j]-p[j],"Tardiness(%r)"%(j),)forkinJ:ifk<=j:continuemodel.addConstr(x[j,k]+x[k,j]==1,"Disjunctive(%s,%s)"%(j,k))forellinJ:ifell>k:model.addConstr(x[j,k]+x[k,ell]+x[ell,j]<=2,"Triangle(%s,%s,%s)"%(j,k,ell),)model.setObjective(quicksum(w[j]*T[j]forjinJ),GRB.MINIMIZE)model.update()model.__data=x,Treturnmodel

n=5# number of jobsJ,p,r,d,w=make_data(n)model=scheduling_linear_ordering(J,p,d,w)model.optimize()z=model.ObjValx,T=model.__datafor(i,j)inx:ifx[i,j].X>0.5:print("x(%s) =%s"%((i,j),int(x[i,j].X+0.5)))foriinT:print("T(%s) =%s"%(i,int(T[i].X+0.5)))print("Opt.value by the linear ordering formulation=",z)

Gurobi Optimizer version 9.1.1 build v9.1.1rc0 (mac64)Thread count: 8 physical cores, 16 logical processors, using up to 16 threadsOptimize a model with 25 rows, 25 columns and 75 nonzerosModel fingerprint: 0xf6727fbfVariable types: 5 continuous, 20 integer (20 binary)Coefficient statistics:  Matrix range     [1e+00, 4e+00]  Objective range  [1e+00, 3e+00]  Bounds range     [1e+00, 1e+00]  RHS range        [1e+00, 4e+00]Found heuristic solution: objective 44.0000000Presolve removed 11 rows and 11 columnsPresolve time: 0.00sPresolved: 14 rows, 14 columns, 50 nonzerosVariable types: 0 continuous, 14 integer (10 binary)Root relaxation: objective 1.700000e+01, 4 iterations, 0.00 seconds    Nodes    |    Current Node    |     Objective Bounds      |     Work Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time     0     0   17.00000    0    1   44.00000   17.00000  61.4%     -    0sH    0     0                      19.0000000   17.00000  10.5%     -    0sH    0     0                      17.0000000   17.00000  0.00%     -    0s     0     0   17.00000    0    1   17.00000   17.00000  0.00%     -    0sExplored 1 nodes (4 simplex iterations) in 0.01 secondsThread count was 16 (of 16 available processors)Solution count 3: 17 19 44 Optimal solution found (tolerance 1.00e-04)Best objective 1.700000000000e+01, best bound 1.700000000000e+01, gap 0.0000%x((2, 1)) = 1x((2, 4)) = 1x((2, 5)) = 1x((3, 1)) = 1x((3, 2)) = 1x((3, 5)) = 1x((4, 1)) = 1x((4, 3)) = 1x((4, 5)) = 1x((5, 1)) = 1T(1) = 7T(2) = 0T(3) = 3T(4) = 0T(5) = 4Opt.value by the linear ordering formulation= 17.0

順列フローショップ問題

順列フローショップ問題(permutation flow shop problem) $F|prmu|C_{\max}$は，以下に定義される問題である．

$n$ 個のジョブを$m$台の機械で順番に処理することを考える．この機械は一度に1つのジョブしか処理できず，ジョブの処理を開始したら途中では中断できないものと仮定する．いま，各機械におけるジョブの処理順序が同一であるとしたとき，最後の機械 $m$ で最後に処理されるジョブの完了時刻を最小にする処理順を求めよ．

ここでは，この問題に対する位置データ定式化(positional data formulation)と呼ばれる定式化を示す．

いま，$n$個のジョブ（添え字は $j$）と $m$台の機械（添え字は $i$）があり，各ジョブは，機械 $1,2,\ldots,m$ の順番で処理されるものとする．（これがフローショップと呼ばれる所以である．）ジョブ$j$ の機械$i$ 上での処理時間を $p_{ij}$ とする．ジョブの追い越しがないので，解はジョブの投入順序，いいかえれば順列で表現できる．（これが順列フローショップと呼ばれる所以である．）ジョブを並べたときの順番を $\kappa$ で表すことにし，順番を表す変数 $x_{j \kappa}$ を用いる．

$x_{j \kappa}$は，ジョブ $j$ が $\kappa$ 番目に投入されるとき $1$ になる $0$-$1$ 変数であり，$s_{i \kappa}$は，機械 $i$ の $\kappa$ 番目に並べられてるジョブの開始時刻を表す実数変数であり，$f_{i \kappa}$は，機械 $i$ の $\kappa$ 番目に並べられてるジョブの完了（終了）時刻を表す実数変数である．

この定式化は，ジョブ自身の開始時刻ではなく，機械にかかるジョブの順番（位置）に対する開始時刻や終了時刻のデータを用いるので，位置データ定式化と呼ばれる．

定式化は以下のようになる．

$$\begin{array}{l l l} minimize & f_{m n} & \\s.t. & \displaystyle\sum_{\kappa} x_{j\kappa} = 1 & \forall j=1,2,\ldots,n \\ & \displaystyle\sum_{j} x_{j\kappa} = 1 & \forall \kappa=1,2,\ldots,n \\ & f_{i \kappa} \leq s_{i,\kappa +1} & \forall i=1,2,\ldots,m, \kappa=1,2,\ldots,n-1 \\ & s_{i \kappa}+\displaystyle\sum_{j} p_{ij} x_{j \kappa} \leq f_{i \kappa} & \forall i=1,2,\ldots,m, \kappa=1,2,\ldots,n \\ & f_{i \kappa} \leq s_{i+1,\kappa} & \forall i=1,2,\ldots,m-1, \kappa=1,2,\ldots,n \\ & x_{j \kappa} \in \{ 0,1 \} & \forall j=1,2,\ldots,n, \kappa=1,2,\ldots,n \\ & s_{i \kappa} \geq 0 & \forall i=1,2,\ldots,m, \kappa=1,2,\ldots,n \\ & f_{i \kappa} \geq 0 & \forall i=1,2,\ldots,m, \kappa=1,2,\ldots,n \end{array}$$

上の定式化において，目的関数は最後の機械($m$)の最後の（$n$番目の）ジョブの完了時刻を最小化している．最初の制約は，各ジョブがいずれかの位置（順番）に割り当てられることを表す．2番目の制約は，各位置（場所）にいずれかのジョブが割り当てられることを表す．3番目の制約は，$\kappa$番目のジョブの完了時刻より $\kappa+1$番目のジョブの開始時刻が遅いことを表す．4番目の制約は，機械$i$ の $\kappa$番目のジョブの開始時刻と完了時刻の関係を規定する制約であり，ジョブの位置を表す変数 $x$ と開始（終了）時刻を表す変数 $s$ ($f$) を繋ぐ制約である．これは，$x_{j\kappa}$ が $1$ のときに処理時間が $p_{ij}$ になることから導かれる．5番目の制約は，機械$i$ 上で $\kappa$番目に割り当てられたジョブの完了時刻より，機械$i+1$ 上で $\kappa$番目に割り当てられたジョブの開始時刻が遅いことを表す．

以下では，ジョブ数 $15$，機械数 $10$ の例題を求解する．

defpermutation_flow_shop(n,m,p):""" permutation_flow_shop problem    Parameters:        - n: number of jobs        - m: number of machines        - p[i,j]: processing time of job i on machine j    Returns a model, ready to be solved.    """model=Model("permutation flow shop")x,s,f={},{},{}forjinrange(1,n+1):forkinrange(1,n+1):x[j,k]=model.addVar(vtype="B",name="x(%s,%s)"%(j,k))foriinrange(1,m+1):forkinrange(1,n+1):s[i,k]=model.addVar(vtype="C",name="start(%s,%s)"%(i,k))f[i,k]=model.addVar(vtype="C",name="finish(%s,%s)"%(i,k))model.update()forjinrange(1,n+1):model.addConstr(quicksum(x[j,k]forkinrange(1,n+1))==1,"Assign1(%s)"%(j))model.addConstr(quicksum(x[k,j]forkinrange(1,n+1))==1,"Assign2(%s)"%(j))foriinrange(1,m+1):forkinrange(1,n+1):ifk!=n:model.addConstr(f[i,k]<=s[i,k+1],"FinishStart(%s,%s)"%(i,k))ifi!=m:model.addConstr(f[i,k]<=s[i+1,k],"Machine(%s,%s)"%(i,k))model.addConstr(s[i,k]+quicksum(p[i,j]*x[j,k]forjinrange(1,n+1))<=f[i,k],"StartFinish(%s,%s)"%(i,k),)model.setObjective(f[m,n],GRB.MINIMIZE)model.update()model.__data=x,s,freturnmodeldefmake_data_permutation_flow_shop(n,m):"""make_data: prepare matrix of m times n random processing times"""p={}foriinrange(1,m+1):forjinrange(1,n+1):p[i,j]=random.randint(1,10)returnp

n=15m=10p=make_data_permutation_flow_shop(n,m)model=permutation_flow_shop(n,m,p)model.optimize()x,s,f=model.__dataprint("Opt.value=",model.ObjVal)

Gurobi Optimizer version 9.1.1 build v9.1.1rc0 (mac64)Thread count: 8 physical cores, 16 logical processors, using up to 16 threadsOptimize a model with 455 rows, 525 columns and 3550 nonzerosModel fingerprint: 0xbc9b34f3Variable types: 300 continuous, 225 integer (225 binary)Coefficient statistics:  Matrix range     [1e+00, 1e+01]  Objective range  [1e+00, 1e+00]  Bounds range     [1e+00, 1e+00]  RHS range        [1e+00, 1e+00]Presolve removed 48 rows and 49 columnsPresolve time: 0.01sPresolved: 407 rows, 476 columns, 3170 nonzerosVariable types: 229 continuous, 247 integer (225 binary)Root relaxation: objective 1.236138e+02, 1043 iterations, 0.04 seconds    Nodes    |    Current Node    |     Objective Bounds      |     Work Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time     0     0  123.61383    0  112          -  123.61383      -     -    0sH    0     0                     158.0000000  123.61383  21.8%     -    0sH    0     0                     155.0000000  123.61383  20.2%     -    0s     0     0  123.95833    0  114  155.00000  123.95833  20.0%     -    0s     0     0  123.99812    0  114  155.00000  123.99812  20.0%     -    0s     0     0  124.00981    0  114  155.00000  124.00981  20.0%     -    0s     0     2  124.03088    0  114  155.00000  124.03088  20.0%     -    0sH   78    93                     151.0000000  124.47625  17.6%   179    0sH   85    93                     149.0000000  124.47625  16.5%   172    0s*  182   168              18     148.0000000  124.54493  15.8%   122    0sH  223   201                     146.0000000  124.54493  14.7%   118    0sH  377   285                     145.0000000  124.85476  13.9%   100    0sH  380   285                     142.0000000  124.85476  12.1%   100    0s  6177  2050  129.00000   28   56  142.00000  128.30382  9.65%  61.7    5s*35346 11391              65     141.0000000  136.27130  3.35%  43.9    9s 40739 12736  139.82068   59   34  141.00000  136.74569  3.02%  42.3   10s 42407 13177  138.92099   57   53  141.00000  136.89366  2.91%  41.9   15s*50385  9787              43     140.0000000  137.31864  1.92%  40.4   15s 85349  4252     cutoff   32       140.00000  139.00000  0.71%  37.1   20sCutting planes:  Gomory: 4  MIR: 16  Flow cover: 14Explored 94493 nodes (3386238 simplex iterations) in 20.96 secondsThread count was 16 (of 16 available processors)Solution count 10: 140 141 142 ... 158Optimal solution found (tolerance 1.00e-04)Best objective 1.400000000000e+02, best bound 1.400000000000e+02, gap 0.0000%Opt.value= 140.0

ジョブショップスケジューリング問題

ジョブショップスケジューリング問題 $J| | C_{\max}$ は，以下のように定義される古典的なスケジューリング問題である．

ジョブを $J_1, J_2,\cdots, J_n$，ジョブ $J_j$ に属するオペレーションを $O_{1j},O_{2j},\cdots,O_{m_j j}$，機械を $M_1,M_2, \cdots, M_m$ とする．オペレーションは $O_{1j},O_{2j},\cdots,O_{m_j j}$ の順で処理されなければならず，オペレーション $O_{ij}$ を処理するには機械 $\mu_{ij}$ を作業時間 $p_{ij}$ だけ占有する．オペレーションが途中で中断できないという仮定の下で，最後のオペレーションが完了する時刻を最小化する各機械上でのオペレーションの処理順序を求める問題．

この問題は中規模な問題例でさえ数理最適化ソルバーで解くことは難しい．実務的な付加条件のほとんどを考慮したスケジューリング最適化システムとしてOptSeqが開発されている． OptSeqについての詳細は，付録1を参照されたい．

ジョブショップスケジューリングのベンチマーク問題例も，OptSeqを使用すると比較的容易に良い解を得ることができる．

例として，以下のサイトから入手できるベンチマーク問題例を読み込んで解いてみよう．

OR-Lib.:http://people.brunel.ac.uk/~mastjjb/jeb/orlib/jobshopinfo.html
多くのジョブショップスケジューリング問題例を含んだサイト:http://jobshop.jjvh.nl/

例としてft06.txtを用いる．このデータは，以下のテキストファイルである．

 6 6 2  1  0  3  1  6  3  7  5  3  4  6 1  8  2  5  4 10  5 10  0 10  3  4 2  5  3  4  5  8  0  9  1  1  4  7 1  5  0  5  2  5  3  3  4  8  5  9  2  9  1  3  4  5  5  4  0  3  3  1 1  3  3  3  5  9  0 10  4  4  2  1

得られた最大完了時刻 $55$ は最適値である．

fromoptseqimport*# ジョブショップスケジューリング問題のベンチマーク問題例fname="../data/jssp/ft06.txt"#fname = "../data/jssp/ta01.txt"f=open(fname,"r")lines=f.readlines()f.close()n,m=map(int,lines[0].split())print("n,m=",n,m)model=Model()act,mode,res={},{},{}forjinrange(m):res[j]=model.addResource(f"machine[{j}]",capacity=1)# prepare data as dicmachine,proc_time={},{}foriinrange(n):L=list(map(int,lines[i+1].split()))forjinrange(m):machine[i,j]=L[2*j]proc_time[i,j]=L[2*j+1]fori,jinproc_time:act[i,j]=model.addActivity(f"Act[{i}_{j}]")mode[i,j]=Mode(f"Mode[{i}_{j}]",proc_time[i,j])mode[i,j].addResource(res[machine[i,j]],1)act[i,j].addModes(mode[i,j])foriinrange(n):forjinrange(m-1):model.addTemporal(act[i,j],act[i,j+1])model.Params.TimeLimit=1model.Params.Makespan=True#model.Params.OutputFlag = Truemodel.optimize()

n,m= 6 6 ================ Now solving the problem ================ Solutions:    source   ---     0     0      sink   ---    55    55  Act[0_0]   ---     5     6  Act[0_1]   ---     6     9  Act[0_2]   ---    16    22  Act[0_3]   ---    30    37  Act[0_4]   ---    42    45  Act[0_5]   ---    49    55  Act[1_0]   ---     0     8  Act[1_1]   ---     8    13  Act[1_2]   ---    13    23  Act[1_3]   ---    28    38  Act[1_4]   ---    38    48  Act[1_5]   ---    48    52  Act[2_0]   ---     0     5  Act[2_1]   ---     5     9  Act[2_2]   ---     9    17  Act[2_3]   ---    18    27  Act[2_4]   ---    27    28  Act[2_5]   ---    42    49  Act[3_0]   ---     8    13  Act[3_1]   ---    13    18  Act[3_2]   ---    22    27  Act[3_3]   ---    27    30  Act[3_4]   ---    30    38  Act[3_5]   ---    45    54  Act[4_0]   ---    13    22  Act[4_1]   ---    22    25  Act[4_2]   ---    25    30  Act[4_3]   ---    38    42  Act[4_4]   ---    48    51  Act[4_5]   ---    52    53  Act[5_0]   ---    13    16  Act[5_1]   ---    16    19  Act[5_2]   ---    19    28  Act[5_3]   ---    28    38  Act[5_4]   ---    38    42  Act[5_5]   ---    42    43

importplotly.expressaspximportpandasaspdimportdatetimeasdtstart="2022/1/1"period="days"start=pd.to_datetime(start)deftime_convert_long(periods):ifperiod=="days":time_=start+dt.timedelta(days=float(periods))elifperiod=="hours":time_=start+dt.timedelta(hours=float(periods))elifperiod=="minutes":time_=start+dt.timedelta(minutes=float(periods))elifperiod=="seconds":time_=start+dt.timedelta(seconds=float(periods))else:raiseTypeError("pariod must be 'days' or 'seconds' or minutes' or 'days'")returntime_.strftime("%Y-%m-%d")L=[]fori,jinproc_time:a=act[i,j]res=Noneforrina.modes[0].requirement:res=r[0]breakst=time_convert_long(a.start)fi=time_convert_long(a.completion)L.append(dict(Task=a.name,Start=st,Finish=fi,Resource=res,Job=f"Job{i}"))df=pd.DataFrame(L)fig=px.timeline(df,x_start="Start",x_end="Finish",y="Resource",color="Job",opacity=0.5)plotly.offline.plot(fig);

OR-toolsを用いた求解

Googleが開発しているOR-toolsを用いると，制約伝播を用いて最適解を出すことが可能である．OR-toolsについての詳細は，以下を参照されたい．

https://developers.google.com/optimization

fname="../data/jssp/ft06.txt"#fname = "../data/jssp/ta01.txt"#fname = "../data/jssp/ta21.txt"f=open(fname,"r")lines=f.readlines()f.close()n,m=map(int,lines[0].split())print("n,m=",n,m)# prepare datamachine,proc_time={},{}foriinrange(n):L=list(map(int,lines[i+1].split()))forjinrange(m):machine[i,j]=L[2*j]proc_time[i,j]=L[2*j+1]jobs_data=[]foriinrange(n):row=[]forjinrange(m):row.append((machine[i,j],proc_time[i,j]))jobs_data.append(row)

n,m= 6 6

importcollectionsfromortools.sat.pythonimportcp_modelmodel=cp_model.CpModel()machines_count=1+max(task[0]forjobinjobs_datafortaskinjob)all_machines=range(machines_count)# Computes horizon dynamically as the sum of all durations.horizon=sum(task[1]forjobinjobs_datafortaskinjob)# Named tuple to store information about created variablestask_type=collections.namedtuple("task_type","start end interval")# Named tuple to manipulate solution informationassigned_task_type=collections.namedtuple("assigned_task_type","start job index duration")# Creates job intervals and add to the corresponding machine listsall_tasks={}machine_to_intervals=collections.defaultdict(list)forjob_id,jobinenumerate(jobs_data):fortask_id,taskinenumerate(job):machine=task[0]duration=task[1]suffix="_%i_%i"%(job_id,task_id)start_var=model.NewIntVar(0,horizon,"start"+suffix)end_var=model.NewIntVar(0,horizon,"end"+suffix)interval_var=model.NewIntervalVar(start_var,duration,end_var,"interval"+suffix)all_tasks[job_id,task_id]=task_type(start=start_var,end=end_var,interval=interval_var)machine_to_intervals[machine].append(interval_var)# Create and add disjunctive constraintsformachineinall_machines:model.AddNoOverlap(machine_to_intervals[machine])# Precedences inside a jobforjob_id,jobinenumerate(jobs_data):fortask_idinrange(len(job)-1):model.Add(all_tasks[job_id,task_id+1].start>=all_tasks[job_id,task_id].end)# Makespan objectiveobj_var=model.NewIntVar(0,horizon,"makespan")model.AddMaxEquality(obj_var,[all_tasks[job_id,len(job)-1].endforjob_id,jobinenumerate(jobs_data)],)model.Minimize(obj_var)# Solve modelsolver=cp_model.CpSolver()solver.parameters.max_time_in_seconds=360.0status=solver.Solve(model)# Create one list of assigned tasks per machine.assigned_jobs=collections.defaultdict(list)forjob_id,jobinenumerate(jobs_data):fortask_id,taskinenumerate(job):machine=task[0]assigned_jobs[machine].append(assigned_task_type(start=solver.Value(all_tasks[job_id,task_id].start),job=job_id,index=task_id,duration=task[1],))# Create per machine output lines.output=""formachineinall_machines:# Sort by starting time.assigned_jobs[machine].sort()sol_line_tasks="Machine "+str(machine)+": "sol_line="           "forassigned_taskinassigned_jobs[machine]:name="job_%i_%i"%(assigned_task.job,assigned_task.index)# Add spaces to output to align columns.sol_line_tasks+="%-10s"%namestart=assigned_task.startduration=assigned_task.durationsol_tmp="[%i,%i]"%(start,start+duration)# Add spaces to output to align columns.sol_line+="%-10s"%sol_tmpsol_line+="\n"sol_line_tasks+="\n"output+=sol_line_tasksoutput+=sol_line# Finally print the solution found.print("Optimal Schedule Length:%i"%solver.ObjectiveValue())print(output)

Optimal Schedule Length: 55Machine 0: job_0_1   job_3_1   job_2_3   job_5_3   job_1_4   job_4_4              [1,4]     [13,18]   [19,28]   [28,38]   [38,48]   [48,51]   Machine 1: job_1_0   job_3_0   job_5_0   job_0_2   job_4_1   job_2_4              [0,8]     [8,13]    [13,16]   [16,22]   [22,25]   [41,42]   Machine 2: job_0_0   job_2_0   job_1_1   job_4_0   job_3_2   job_5_5              [0,1]     [2,7]     [8,13]    [13,22]   [22,27]   [42,43]   Machine 3: job_2_1   job_5_1   job_3_3   job_0_3   job_1_5   job_4_5              [7,11]    [16,19]   [27,30]   [30,37]   [48,52]   [52,53]   Machine 4: job_1_2   job_4_2   job_3_4   job_5_4   job_0_5   job_2_5              [13,23]   [25,30]   [30,38]   [38,42]   [42,48]   [48,55]   Machine 5: job_2_2   job_5_2   job_1_3   job_0_4   job_4_3   job_3_5              [11,19]   [19,28]   [28,38]   [38,41]   [41,45]   [45,54]

数理最適化ソルバーを用いた求解

数理最適化による定式化も可能である．以下では，次の資源制約付きプロジェクトスケジューリング問題で導入する時刻添え字定式化と，1機械リリース時刻付き完了時刻和最小化問題で用いた離接定式化による求解を試みる．

ここでは，各ジョブに含まれるオペレーション数が機械の数 $m$ に一致するものと仮定する．オペレーションをジョブと順序の組として $(i,j)$ もしくは $(k,l)$ と書く．

ジョブショップスケジューリング問題は資源制約付きプロジェクトスケジューリング問題の特殊形であるので，時刻添え字定式化は省略する．

離接定式化では，オペレーション $(i,j)$ の開始時刻を表す連続変数 $s_{ij}$ と，オペレーション $(i,j)$ が他のオペレーション $(k,l)$ に先行する（前に処理される）とき $1$ になる $0$-$1$ 整数変数 $x_{ijkl}$ を用いる．また，最大完了時刻を表す実数変数を $z$ とする．

非常に大きな数 $M$ を用いると，離散定式化は以下のように書ける．

$$\begin{array}{l l l} minimize & z & \\ s.t. & s_{ij} + p_{ij} -M (1-x_{ijkl}) \leq s_{kl} & \forall j \neq k \\ & x_{ijkl} +x _{klij}= 1 & \forall (i,j) < (k,l) \\ & s_{ij} + p_{ij} \leq s_{i,j+1} & \forall i, j=1,\ldots,m-1 \\ & s_{im} \leq z & \forall i \\ & s_{i1} \geq 0 & \forall i \\ & x_{ijkl} \in \{ 0,1 \} & \forall (i,j) \neq (k,l) \end{array}$$

fname="../data/jssp/ft06.txt"# fname = "../data/jssp/ta21.txt"f=open(fname,"r")lines=f.readlines()f.close()n,m=map(int,lines[0].split())print("n,m=",n,m)machine,proc_time={},{}foriinrange(n):L=list(map(int,lines[i+1].split()))forjinrange(m):machine[i,j]=L[2*j]proc_time[i,j]=L[2*j+1]

n,m= 6 6

fromgurobipyimportModel,GRB,quicksummodel=Model("scheduling: time-index")T=55# メイクスパン以上の計画期間s,x={},{}# s - start time variable,   x=1 if job (i,j) starts at period tfori,jinproc_time:s[i,j]=model.addVar(vtype="C",name=f"s({i},{j})")fortinrange(1,T-proc_time[i,j]+2):x[i,j,t]=model.addVar(vtype="B",name=f"x({i},{j},{t})")z=model.addVar(vtype="C",name="z")# completion timemodel.update()fori,jinproc_time:# job execution constraintsmodel.addConstr(quicksum(x[i,j,t]fortinrange(1,T-proc_time[i,j]+2))==1)# start time constraintsmodel.addConstr(quicksum((t-1)*x[i,j,t]fortinrange(2,T-proc_time[i,j]+2))==s[i,j])# resource upper bound constraintsfortinrange(1,T+1):forkinrange(m):# machinemodel.addConstr(quicksum(x[i,j,t_]fori,jinproc_timeifmachine[i,j]==kfort_inrange(max(t-proc_time[i,j]+1,1),min(t+1,T-proc_time[i,j]+2)))<=1)# time (precedence) constraints, i.e., s[k]-s[j] >= p[j]foriinrange(n):forjinrange(m-1):model.addConstr(s[i,j]+proc_time[i,j]<=s[i,j+1])model.addConstr(s[i,m-1]+proc_time[i,m-1]<=z)model.setObjective(z,GRB.MINIMIZE)model.optimize()

Set parameter UsernameAcademic license - for non-commercial use only - expires 2023-11-07Gurobi Optimizer version 9.5.2 build v9.5.2rc0 (mac64[x86])Thread count: 8 physical cores, 16 logical processors, using up to 16 threadsOptimize a model with 438 rows, 1856 columns and 13375 nonzerosModel fingerprint: 0x81d78eefVariable types: 37 continuous, 1819 integer (1819 binary)Coefficient statistics:  Matrix range     [1e+00, 5e+01]  Objective range  [1e+00, 1e+00]  Bounds range     [1e+00, 1e+00]  RHS range        [1e+00, 1e+01]Presolve removed 108 rows and 1160 columnsPresolve time: 0.03sPresolved: 330 rows, 696 columns, 4424 nonzerosVariable types: 0 continuous, 696 integer (665 binary)Root relaxation: objective 4.955932e+01, 325 iterations, 0.01 seconds (0.00 work units)    Nodes    |    Current Node    |     Objective Bounds      |     Work Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time     0     0   49.55932    0  105          -   49.55932      -     -    0s     0     0   50.37220    0  146          -   50.37220      -     -    0s     0     0   50.37507    0  146          -   50.37507      -     -    0s     0     0   50.37507    0  152          -   50.37507      -     -    0s     0     0   51.83488    0  168          -   51.83488      -     -    0s     0     0   52.09226    0  170          -   52.09226      -     -    0s     0     0   52.09238    0  171          -   52.09238      -     -    0s     0     0   53.31236    0  169          -   53.31236      -     -    0s     0     0   53.39541    0  160          -   53.39541      -     -    0s     0     0   53.42883    0  155          -   53.42883      -     -    0sH    0     0                      55.0000000   53.84137  2.11%     -    0s     0     0     cutoff    0        55.00000   55.00000  0.00%     -    0sCutting planes:  Cover: 7  Implied bound: 2  Clique: 41  MIR: 20  StrongCG: 4  GUB cover: 37  Zero half: 21  RLT: 5  Relax-and-lift: 6Explored 1 nodes (3003 simplex iterations) in 0.34 seconds (0.23 work units)Thread count was 16 (of 16 available processors)Solution count 1: 55 Optimal solution found (tolerance 1.00e-04)Best objective 5.500000000000e+01, best bound 5.500000000000e+01, gap 0.0000%

fromgurobipyimportModel,GRB,quicksummodel=Model("scheduling: disjunctive")M=3000# big Ms,x={},{}# start time variable, x[i,j,k,l] = 1 if job (i,j) precedes job (k,l), 0 otherwisefori,jinproc_time:s[i,j]=model.addVar(vtype="C",name=f"s({i},{j})")fork,linproc_time:ifmachine[i,j]==machine[k,l]and(i,j)!=(k,l):x[i,j,k,l]=model.addVar(vtype="B",ub=1.,name=f"x({i},{j},{k},{l})")z=model.addVar(vtype="C",name="z")#completion timemodel.update()fori,jinproc_time:fork,linproc_time:ifmachine[i,j]==machine[k,l]and(i,j)!=(k,l):model.addConstr(s[i,j]-s[k,l]+M*x[i,j,k,l]<=M-proc_time[i,j])ifmachine[i,j]==machine[k,l]and(i,j)<(k,l):model.addConstr(x[i,j,k,l]+x[k,l,i,j]==1)foriinrange(n):forjinrange(m-1):model.addConstr(s[i,j]+proc_time[i,j]<=s[i,j+1])model.addConstr(s[i,m-1]+proc_time[i,m-1]<=z)model.setObjective(z,GRB.MINIMIZE)model.Params.LogFile="gurobi_jssp.log"model.Params.TimeLimit=60model.Params.OutputFlag=1model.optimize()#パラメータチューニング# model.Params.TuneTimeLimit = 3600# model.Params.TuneCriterion = 2 # best feasible solution found.# model.tune()

Set parameter LogFile to value "gurobi_jssp.log"Set parameter TimeLimit to value 60Gurobi Optimizer version 9.5.2 build v9.5.2rc0 (mac64[x86])Thread count: 8 physical cores, 16 logical processors, using up to 16 threadsOptimize a model with 306 rows, 217 columns and 792 nonzerosModel fingerprint: 0xe048963cVariable types: 37 continuous, 180 integer (180 binary)Coefficient statistics:  Matrix range     [1e+00, 3e+03]  Objective range  [1e+00, 1e+00]  Bounds range     [1e+00, 1e+00]  RHS range        [1e+00, 3e+03]Presolve removed 90 rows and 90 columnsPresolve time: 0.00sPresolved: 216 rows, 127 columns, 612 nonzerosVariable types: 37 continuous, 90 integer (90 binary)Found heuristic solution: objective 170.0000000Root relaxation: objective 4.700000e+01, 61 iterations, 0.00 seconds (0.00 work units)    Nodes    |    Current Node    |     Objective Bounds      |     Work Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time     0     0   47.00000    0   10  170.00000   47.00000  72.4%     -    0sH    0     0                      73.0000000   47.00000  35.6%     -    0sH    0     0                      60.0000000   47.00000  21.7%     -    0s     0     0   47.00000    0   11   60.00000   47.00000  21.7%     -    0sH    0     0                      58.0000000   47.00000  19.0%     -    0s     0     0   47.00000    0   12   58.00000   47.00000  19.0%     -    0s     0     0   47.00000    0   12   58.00000   47.00000  19.0%     -    0s     0     0   47.00299    0   18   58.00000   47.00299  19.0%     -    0sH    0     0                      57.0000000   47.00299  17.5%     -    0sH    0     0                      56.0000000   50.00000  10.7%     -    0s     0     0   50.00000    0    5   56.00000   50.00000  10.7%     -    0sH    0     0                      55.0000000   50.00000  9.09%     -    0s     0     0   50.00000    0    7   55.00000   50.00000  9.09%     -    0s     0     0   50.00000    0    8   55.00000   50.00000  9.09%     -    0s     0     0 infeasible    0        55.00000   55.00000  0.00%     -    0sExplored 1 nodes (568 simplex iterations) in 0.07 seconds (0.02 work units)Thread count was 16 (of 16 available processors)Solution count 9: 55 55 55 ... 170Optimal solution found (tolerance 1.00e-04)Best objective 5.499999999440e+01, best bound 5.499999999440e+01, gap 0.0000%

importgrblogtoolsasgltresults=glt.parse(model.Params.LogFile)summary=results.summary()nodelog_progress=results.progress("nodelog")nodelog_progress.Incumbent.plot()nodelog_progress.BestBd.plot();#summary.T

資源制約付きプロジェクトスケジューリング問題

様々なスケジューリング問題を一般化した，以下の資源制約付き（プロジェクト）スケジューリング問題を考える．

ジョブの集合 $Job$ ，各時間ごとの使用可能量の上限をもつ資源 $Res$ が与えられている．資源は，ジョブごとに決められた作業時間の間はジョブによって使用されるが，作業完了後は，再び使用可能になる．ジョブ間に与えられた時間制約を満たした上で，ジョブの作業開始時刻に対する任意の関数の和を最小化するような，資源のジョブへの割り振りおよびジョブの作業開始時刻を求める．

時刻を離散化した期の概念を用いた定式化（時刻添え字定式化）を示す．

まず，定式化で用いる集合，入力データ，変数を示す．

集合:

$Job$: ジョブの集合．添え字は $j,k$
$Res$: 資源の集合．添え字は $r$
$Prec$: ジョブ間の時間制約を表す集合． $Job \times Job$ の部分集合で， $(j,k) \in Prec$ のとき，ジョブ $j$ とジョブ $k$ の開始（完了）時刻間に何らかの制約が定義されるものとする

入力データ:

$T$: 最大の期数．期の添え字は $t=1,2,\ldots,T$．連続時間との対応づけは，時刻 $t-1$ から $t$ までを期 $t$ と定義する．なお，期を表す添え字は $t$ の他に $s$ も用いる．
$p_j$: ジョブ $j$ の処理時間．非負の整数を仮定
$Cost_{jt}$: ジョブ $j$ を期$t$ に開始したときの費用
$a_{jrt}$: ジョブ $j$ の開始後 $t(=0,1,\ldots,p_j-1)$期経過時の処理に要する資源$r$ の量
$RUB_{rt}$: 時刻 $t$ における資源 $r$ の使用可能量の上限

変数:

$x_{jt}$: ジョブ $j$ を時刻 $t$ に開始するとき $1$，それ以外のとき $0$ を表す $0$-$1$変数

以下に資源制約付きスケジューリング問題の定式化を示す．

目的関数:$$minimize \sum_{j \in Job} \sum_{t=1}^{T-p_j+1} Cost_{jt} x_{jt}$$
ジョブ遂行制約：$$\sum_{t=1}^{T-p_j+1} x_{jt} =1　 \forall j \in Job$$すべてのジョブは必ず1度処理されなければならないことを表す．

資源制約：

$$ \sum_{j \in Job} \sum_{s =\max\{t-p_j+1,1\}}^{\min\{t,T-p_j+1\}} a_{jr,t-s} x_{js} \leq RUB_{rt} \ \ \ \forall r \in Res , t \in \{1,2,\ldots,T \}$$

ある時刻 $t$ に作業中であるジョブの資源使用量の合計が，資源使用量の上限を超えないことを規定する．時刻 $t$ に作業中であるジョブは，時刻 $t-p_j+1$ から $t$ の間に作業を開始したジョブであるので，上式を得る．

時間制約：

ジョブ $j$ の開始時刻 $S_j$ は $\sum_{t=2}^{T-p_j+1} (t-1)x_{jt}$ であるので，完了時刻 $C_j$ はそれに $p_j$ を加えたものになる．ジョブの対 $(j,k) \in Prec$ に対しては，ジョブ $j,k$ の開始（完了）時刻間に何らかの制約が付加される．たとえば，ジョブ $j$ の処理が終了するまで，ジョブ $k$ の処理が開始できないことを表す先行制約は，

$$\sum_{t=2}^{T-p_j+1} (t-1)x_{jt} +p_j \leq \sum_{t=2}^{T-p_k+1} (t-1) x_{kt} \ \ \ \forall (j,k) \in Prec$$

と書くことができる．ここで，左辺の式は，ジョブ $j$の完了時刻を表し，右辺の式はジョブ $k$ の開始時刻を表す．

以下に，数理最適化ソルバー Gurobiを用いたモデルを示す．実際には，Gurobiでは中規模問題例までしか解くことができない．ジョブショップスケジューリング問題と同様に，OptSeqやOR-toolsを用いることによって求解できる．

ただし，実務の問題はさらなる付加条件がつくことが多い．OptSeqは，より一般的なスケジューリング問題をモデル化できるように設計されている．

fromgurobipyimportModel,quicksum,GRB,multidict# from mypulp import Model, quicksum, GRB, multidictdefrcs(J,P,R,T,p,c,a,RUB):"""rcs -- model for the resource constrained scheduling problem    Parameters:        - J: set of jobs        - P: set of precedence constraints between jobs        - R: set of resources        - T: number of periods        - p[j]: processing time of job j        - c[j,t]: cost incurred when job j starts processing on period t.        - a[j,r,t]: resource r usage for job j on period t (after job starts)        - RUB[r,t]: upper bound for resource r on period t    Returns a model, ready to be solved.    """model=Model("resource constrained scheduling")s,x={},{}# s - start time variable,   x=1 if job j starts on period tforjinJ:s[j]=model.addVar(vtype="C",name="s(%s)"%j)fortinrange(1,T-p[j]+2):x[j,t]=model.addVar(vtype="B",name="x(%s,%s)"%(j,t))model.update()forjinJ:# job execution constraintsmodel.addConstr(quicksum(x[j,t]fortinrange(1,T-p[j]+2))==1,"ConstrJob(%s,%s)"%(j,t),)# start time constraintsmodel.addConstr(quicksum((t-1)*x[j,t]fortinrange(2,T-p[j]+2))==s[j],"ConstrJob(%s,%s)"%(j,t),)# resource upper bound constraintsfortinrange(1,T+1):forrinR:model.addConstr(quicksum(a[j,r,t-t_]*x[j,t_]forjinJfort_inrange(max(t-p[j]+1,1),min(t+1,T-p[j]+2)))<=RUB[r,t],"ResourceUB(%s,%s)"%(r,t))# time (precedence) constraints, i.e., s[k]-s[j] >= p[j]for(j,k)inP:model.addConstr(s[k]-s[j]>=p[j],"Precedence(%s,%s)"%(j,k))model.setObjective(quicksum(c[j,t]*x[j,t]for(j,t)inx),GRB.MINIMIZE)model.update()model.__data=x,sreturnmodeldefmake_1r():J,p=multidict({# jobs, processing times1:1,2:3,3:2,4:2,})P=[(1,2),(1,3),(2,4)]R=[1]T=6c={}forjinJ:fortinrange(1,T-p[j]+2):c[j,t]=1*(t-1+p[j])a={(1,1,0):2,(2,1,0):2,(2,1,1):1,(2,1,2):1,(3,1,0):1,(3,1,1):1,(4,1,0):1,(4,1,1):2,}RUB={(1,1):2,(1,2):2,(1,3):1,(1,4):2,(1,5):2,(1,6):2}return(J,P,R,T,p,c,a,RUB)defmake_2r():J,p=multidict({# jobs, processing times1:2,2:2,3:3,4:2,5:5,})P=[(1,2),(1,3),(2,4)]R=[1,2]T=6c={}forjinJ:fortinrange(1,T-p[j]+2):c[j,t]=1*(t-1+p[j])a={# resource 1:(1,1,0):2,(1,1,1):2,(2,1,0):1,(2,1,1):1,(3,1,0):1,(3,1,1):1,(3,1,2):1,(4,1,0):1,(4,1,1):1,(5,1,0):0,(5,1,1):0,(5,1,2):1,(5,1,3):0,(5,1,4):0,# resource 2:(1,2,0):1,(1,2,1):0,(2,2,0):1,(2,2,1):1,(3,2,0):0,(3,2,1):0,(3,2,2):0,(4,2,0):1,(4,2,1):2,(5,2,0):1,(5,2,1):2,(5,2,2):1,(5,2,3):1,(5,2,4):1,}RUB={(1,1):2,(1,2):2,(1,3):2,(1,4):2,(1,5):2,(1,6):2,(1,7):2,(2,1):2,(2,2):2,(2,3):2,(2,4):2,(2,5):2,(2,6):2,(2,7):2,}return(J,P,R,T,p,c,a,RUB)

(J,P,R,T,p,c,a,RUB)=make_2r()model=rcs(J,P,R,T,p,c,a,RUB)model.optimize()x,s=model.__dataprint("Opt.value=",model.ObjVal)for(j,t)inx:ifx[j,t].X>0.5:print(x[j,t].VarName,x[j,t].X)forjins:ifs[j].X>0.:print(s[j].VarName,s[j].X)

Gurobi Optimizer version 9.5.2 build v9.5.2rc0 (mac64[x86])Thread count: 8 physical cores, 16 logical processors, using up to 16 threadsOptimize a model with 25 rows, 26 columns and 127 nonzerosModel fingerprint: 0x164fed55Variable types: 5 continuous, 21 integer (21 binary)Coefficient statistics:  Matrix range     [1e+00, 4e+00]  Objective range  [2e+00, 6e+00]  Bounds range     [1e+00, 1e+00]  RHS range        [1e+00, 2e+00]Found heuristic solution: objective 23.0000000Presolve removed 25 rows and 26 columnsPresolve time: 0.00sPresolve: All rows and columns removedExplored 0 nodes (0 simplex iterations) in 0.01 seconds (0.00 work units)Thread count was 1 (of 16 available processors)Solution count 1: 23 Optimal solution found (tolerance 1.00e-04)Best objective 2.300000000000e+01, best bound 2.300000000000e+01, gap 0.0000%Opt.value= 23.0x(1,1) 1.0x(2,3) 1.0x(3,4) 1.0x(4,5) 1.0x(5,1) 1.0s(2) 2.0s(3) 3.0s(4) 4.0

Movatterモバイル変換

スケジューリング問題