一种求解炼油厂连续时间调度模型的Lagrange分解算法
施磊, 江永亨 , 王凌, 黄德先    
清华大学 自动化系, 北京 100084
摘要:在炼油厂连续时间调度模型中,随着调度问题规模的增大,求解耗时会显著增长。该文提出了一种基于Lagrange分解的求解算法。根据炼油厂生产流程特点,将调度模型分解成9个子问题,并在子问题中加入辅助约束加快Lagrange乘子收敛。针对问题特点设计了乘子初始化方案、乘子迭代方案和对偶解可行化方法。案例仿真选用了3个具有不同调度周期和订单数量的案例进行仿真,结果表明:采用该文提出的算法能够显著提高模型的求解效率,算法求解时间与直接求解和普通Lagrange分解算法相比都要少,且随着问题规模的增大优势会更明显。从求解结果上看,算法能够得到原问题的最优解或者近似最优解。
关键词炼油厂调度    连续时间表达    Lagrange分解    
Lagrangian decomposition approach for solving continuous-time scheduling models of refinery production problems
SHI Lei, JIANG Yongheng , WANG Ling, HUANG Dexian    
Department of Automation, Tsinghua University, Beijing 100084, China
Abstract: Continuous-time models need more computational effort to solve refinery production scheduling problems as the scheduling problem size increases. A new Lagrangian decomposition approach was used which divides the whole scheduling problem into nine subproblems. The convergence of Lagrange multipliers is accelerated by adding auxiliary constraints to the subproblems. This paper gives an initialization scheme for the Lagrange multipliers, a hybrid method to update the Lagrange multipliers and a heuristic algorithm to find feasible solutions. Computational results for three cases with different time horizons and different numbers of orders show that the Lagrangian scheme improves the computational efficiency and obtains optimal or near-optimal solutions.
Key words: refinery scheduling    continuous-time representation    Lagrangian decomposition    

石油工业是中国国民经济的支柱产业之一。在炼油厂中,尤其是在大型的现代化炼油厂中,生产调度建模是一个在生产和管理方面的有效工具。通过调度系统,企业可以获得更好的经济效益。由于生产过程的复杂性,炼油厂短期调度是一个极具挑战性的问题。在过去的20 a里,研究人员提出了很多相关模型和算法来解决此问题。Pinto等[1]把炼油厂调度问题归纳成混合整数优化模型,模型中分别采用了连续时间表达和离散时间表达的描述方式。Göthe-Lundgren等[2]针对一个原油炼制企业归纳出生产计划和调度问题,里面考虑了生产过程的计划和利用。生产过程包括了1套精馏装置和2套加氢装置。调度问题的目标是确定在每个时间点,每套装置所处的生产模式。Jia和 Ierapetritou[3]针对原油短期调度和炼油厂生产调度建立了一个综合的数学规划模型。其中在炼油厂生产模型中,整个炼油厂生产流程分成了3个阶段:原油卸货和调合、生产装置加工、成品油调合和交货。Luo和 Rong[4]对炼油厂短期调度问题提出了一个具有2层结构的分层算法。上层是优化模型,下层是用在仿真系统里的启发式规则。但是对于装置操作模式是如何实现的没有讨论。Mouret等[5]把炼油厂计划问题和原油调度问题结合起来提出了一个综合问题。他们把炼油厂计划问题归纳成一个多阶段流程优化问题。并且假设在每个阶段炼油厂生产装置都处在一个稳定生产状态。Cao等[6]基于数据驱动模型,开发了一种滚动时间优化控制策略,来解决一个实际炼油厂的柴油生产在线调度问题。他们提出了一种混合整数非线性规划调度模型,在模型中考虑了调合时质量指标非线性关系和质量守恒原理。Shah等[7]对炼油厂调度、计划和供应链管理方法论做了综述。Joly[8]说明了在巴西石油公司的业务中炼油厂计划调度活动的核心地位,概述了过去和现在取得的重要成果,并介绍了应对现在和未来所面临挑战的长期战略。这些研究成果对炼油厂调度水平的发展具有重要作用。

对于炼油厂调度问题来说,调度方案能否顺利执行非常重要。调度方案能否执行与调度模型是否准确描述生产过程相关。对于炼油厂调度模型来说,需要考虑的就是生产装置的生产模式如何描述。Shi等[9]基于由装置级先进控制得到的生产装置操作模式提出了一种炼油厂全厂调度模型。该模型考虑了装置模式发生切换的过渡过程,因此能更好地反映生产过程动态性。

但是随着调度问题的规模进一步增大,相应的调度周期增长、订单数量增多,连续时间模型中所需要的事件点个数也会相应增多。同时由于连续时间模型中的事件点对应的时间是通过模型求解优化得到的,当调度周期增长之后,事件点发生时间的可行域变大。这两方面原因都会导致连续时间调度模型的求解时间变长。Lagrange分解方法是一种能求解大规模数学规划问题的算法[10, 11, 12],曾成功应用于炼钢生产过程[13]、生产计划与调度[14]、电动机供应链问题[15]和页岩气系统[16]等。

为了提高连续时间调度模型的求解效率,本文提出一种针对连续时间调度模型的Lagrange分解算法,根据连续时间模型的特点,利用供需关系,通过计算中间组分油的等价价格来设置 Lagrange 乘子初始值;结合次梯度算法、割平面方法、box-step方法和指数加权次梯度设计乘子迭代方法来加快算法的收敛性;在成品油调合交货子问题中增加辅助约束提高解的质量。

1 炼油厂连续时间调度模型简介

本文研究的某炼油厂生产流程示意图如图 1所示。整个炼油厂生产可以分成3个阶段:原油供应、生产加工和成品油调合交货。原油从原油罐中向常压蒸馏装置输送,并假设原油储量足够生产需要。生产加工装置包括了常压蒸馏装置、减压蒸馏装置、催化裂化装置、加氢脱硫装置、醚化装置、甲基叔丁基醚(MTBE)生产装置、柴油加氢精制装置1号和装置2号以及催化重整装置。在成品油调合交货阶段,一共涉及8种调合组分油和8种成品油。调合组分油和成品油储存在油罐中。所建立的模型是混合整数线性规划模型(MILP)。

图 1 炼油厂生产流程示意图

在本节中,主要介绍模型中的两大类约束和目标函数。第一类约束是和装置生产模式相关的约束,第二类约束是装置生产相关约束,

1.1 操作模式相关约束

操作模式相关约束包括时间顺序约束、过渡过程变量和稳态变量约束、过渡持续时间约束。

1) 时间顺序约束

${T_{n - 1}} + {T_{\min }} \le {T_n},\;\;2 \le n \le {n_{\max }},$ (1)
${T_1} \ge {T_{\min }},$ (2)
${T_{{n_{\max }}}} = {\rm{TH}}{\rm{.}}$ (3)

2) 过渡过程变量约束

$\begin{array}{l} \sum\limits_m {\sum\limits_{m'} {y_{u,m,m',n}^s} + \sum\limits_m {\sum\limits_{m'} {y_{u,m,m',n}^e} } } \le 1,\\ \;\;\;\;\;\;\forall u \in U,m \in {M_u},\\ \;\;\;m' \in {M_u},m \ne m',n \ge 2. \end{array}$ (4)
$\begin{array}{l} 0 \le \sum\limits_{k \le n} {\sum\limits_m {\sum\limits_{m'} {y_{u,m,m',k}^s} } } - \sum\limits_{k \le n} {\sum\limits_m {\sum\limits_{m'} {y_{u,m,m',k}^e} } } \le 1,\\ \;\;\forall u \in U,m \in {M_u},m' \in {M_u},m \ne m',n \ge 2. \end{array}$ (5)
$\begin{array}{l} \sum\limits_{k \le n} {\sum\limits_m {\sum\limits_{m'} {y_{u,m,m',k}^s} } } - \sum\limits_{k \le n} {\sum\limits_m {\sum\limits_{m'} {y_{u,m,m',k}^e} } } = 0,\\ \;\;\;\;\;u \in U,m \in {M_u},m' \in {M_u},m \ne m'. \end{array}$ (6)

3) 稳态变量约束

$\begin{array}{l} {w_{u,m,n}} + {w_{u,m',n + 1}} \le 1,\\ \;\;\forall u \in U,m \in {M_u},\\ m' \in {M_u},m \ne m',n < {n_{\max }}. \end{array}$ (7)
$\forall u \in U,m \in {M_u},1 < n < {n_{\max }}.$ (8)
$\begin{array}{l} \;\;\;y_{u,m,m',n}^s \le {w_{u,m,n}},\\ \;\;\;\forall u \in U,m \in {M_u},\\ m' \in {M_u},m \ne m',n \ge 2. \end{array}$ (9)
$m \in {M_u},m' \in {M_u},m \ne m',n < {n_{\max }}.$ (10)
$\begin{array}{l} \sum\limits_{k \le n} {\sum\limits_m {\sum\limits_{m'} {y_{u,m,m',k}^s} } } - \sum\limits_{k \le n} {\sum\limits_{}^m {\sum\limits_{m'} {y_{u,m,m',k}^e} } } = \\ \;\;\;\;\;1 - \sum\limits_{m'} {{w_{u,m,n + 1}}} ,\forall u \in U,\\ m \in {M_u},m' \in {M_u},m \ne m',n < {n_{\max }}. \end{array}$ (11)
$\begin{array}{l} \;\;\;\;\;\;\;\;\;\;\sum\limits_{m'} {{w_{u,m,1}}} = 1,\\ \forall u \in U,m \in {M_u},m' \in {M_u},n = 1. \end{array}$ (12)
$\begin{array}{l} \;\;\;\;\;\;{w_{{\rm{FCCU,}}m,n}} = {w_{{\rm{HDS,}}m,n}},\\ m \in {M_{{\rm{FCCU}}}},{M_{{\rm{FCCU}}}} = {M_{{\rm{HDS}}}}. \end{array}$ (13)
$\begin{array}{l} \;\;\;\;\;\;{w_{{\rm{HDS,}}m,n}} = {w_{{\rm{ETH,}}m,n}},\\ m \in {M_{{\rm{HDS}},m,n}},{M_{{\rm{HDS}}}} = {M_{{\rm{ETH}}}}. \end{array}$ (14)

4) 过渡过程持续时间约束

$\begin{array}{l} - T{T_{u,m,m'}}\left[ {2 - y_{u,m,m',n}^s - y_{u,m,m',n'}^e + \sum\limits_{k = n + 1}^{n' - 1} {\left( {y_{u,m,m',k}^s + y_{u,m,m',k}^e} \right)} } \right]\\ \;\;\;\;\;\;\;\;\;\;\;\; \le {T_{n'}} - {T_n} - T{T_{u,m,m'}} \le \left( {TH - T{T_{u,m,m'}}} \right)\\ \;\;\;\;\;\;\;\left[ {2 - y_{u,m,m',n}^s - y_{u,m,m',n'}^e + \sum\limits_{k = n + 1}^{n' - 1} {\left( {y_{u,m,m',k}^s + y_{u,m,m',k}^e} \right)} } \right],\\ \;\;\;\;\;\;\;\;\;\forall u \in U,m \in {M_u},m' \in {M_u},m \ne m',n' \ge n. \end{array}$ (15)
1.2 装置生产相关约束

装置生产相关约束包括质量平衡约束、装置容量约束、成品油调合及交货约束。

1) 质量平衡约束

$\begin{array}{l} Q{O_{u,oi,n}} = \sum\limits_m {\sum\limits_{m'} {{\rm{tQ}}{{\rm{I}}_{u,m,m',n}}{\rm{tYiel}}{{\rm{d}}_{u,oi,m,m'}}} } + \\ \;\;\;\;\;\;\;\;\;\;\sum\limits_m {{\rm{sQ}}{{\rm{I}}_{u,m,n}}{\rm{Yiel}}{{\rm{d}}_{u,oi,m}}} ,\\ \;\;\;\;\;\;\;\;\;\forall u \in U,{\rm{oi}} \in {\rm{oi}},n \in N. \end{array}$ (16)
$\begin{array}{l} {\rm{Q}}{{\rm{I}}_{u,n}} = \sum\limits_m {\sum\limits_{m'} {{\rm{tQ}}{{\rm{I}}_{u,m,m',n}}} + } \sum\limits_m {{\rm{sQ}}{{\rm{I}}_{u,m,n}}} ,\\ \;\;\;\;\;\;\;\;\;\;\;\forall u \in U,n \in N. \end{array}$ (17)
$\begin{array}{l} \sum\limits_{{u_1}}^{} {Q{O_{{u_1},oi,n}} = } \sum\limits_{{u_2}}^{} {Q{O_{{u_2},oi,n}}} ,\forall {\rm{om}} \in {\rm{OM,}}\\ \;\;\;\;\;\;n \in N,{u_1} \in U_{{\rm{om}}}^{{\rm{out}}},{u_2} \in U_{{\rm{om}}}^{{\rm{in}}}. \end{array}$ (18)
$\begin{array}{l} \;\;\;\;\;\;\;\;\;\;\;\;\;\sum\limits_u^{} {Q{O_{u,oi,n}} = } \\ Q{I_{{\rm{oc,n}}}},\forall {\rm{oc}} \in {\rm{OC,}}n \in N,u \in U_{{\rm{oc}}}^{{\rm{out}}}. \end{array}$ (19)
$\begin{array}{l} \;\;\;IN{V_{{\rm{oc}},1}} = {\rm{IN}}{{\rm{V}}_{{\rm{oc}},{\rm{ini}}}} + \\ Q{I_{{\rm{oc}},1}} - Q{O_{{\rm{oc}},1}},\forall {\rm{oc}} \in {\rm{OC}}. \end{array}$ (20)
$\begin{array}{l} \;\;\;IN{V_{o,1}} = {\rm{IN}}{{\rm{V}}_{o,{\rm{ini}}}} + \\ Q{I_{o,1}} - \sum\limits_l {{D_{l,o,1}}} ,\forall o \in O. \end{array}$ (21)
$\begin{array}{l} IN{V_{{\rm{oc}},n}} = {\rm{IN}}{{\rm{V}}_{{\rm{oc}},n - 1}} + Q{I_{{\rm{oc}},n}} - {\rm{Q}}{{\rm{O}}_{{\rm{oc,}}n}},\\ \;\;\;\;\;\;\;\forall {\rm{oc}} \in {\rm{OC}},n > 1. \end{array}$ (22)
$\begin{array}{l} IN{V_{o,n}} = {\rm{IN}}{{\rm{V}}_{o,n - 1}} + Q{I_{o,n}} - \sum\limits_l {{D_{l,o,n}}} ,\\ \;\;\;\;\;\;\;\;\forall o \in O,n > 1. \end{array}$ (23)
$\begin{array}{l} \sum\limits_{{\rm{oc}}} {{Q_{{\rm{oc}},o,n}}} = {\rm{Q}}{{\rm{I}}_{o,n}},\\ \forall o \in O,n \in N. \end{array}$ (24)
$\begin{array}{l} \sum\limits_o {{Q_{{\rm{oc}},o,n}}} = {\rm{Q}}{{\rm{I}}_{oc,n}},\\ \forall {\rm{oc}} \in {\rm{OC}},n \in N. \end{array}$ (25)

2) 容量约束

$\begin{array}{l} w{T_{u,m,1}}{\rm{QI}}_u^{\min } \le s{\rm{Q}}{{\rm{I}}_{u,m,1}} \le w{T_{u,m,1}}{\rm{QI}}_u^{\max },\\ \;\;\;\;\;\forall u \in U,m \in {M_u},n = 1. \end{array}$ (26)
$\begin{array}{l} w{T_{u,m,n}}{\rm{QI}}_u^{\min } \le s{\rm{Q}}{{\rm{I}}_{u,m,n}} \le w{T_{u,m,n}}{\rm{QI}}_u^{\max },\\ \;\;\;\;\;\forall u \in U,m \in {M_u},m \ne m',n > 1. \end{array}$ (27)
$\begin{array}{l} y{T_{u,m,m',n}}{\rm{QI}}_u^{\min } \le t{\rm{Q}}{{\rm{I}}_{u,m,m',n}} \le y{T_{u,m,m',n}}{\rm{QI}}_u^{\max },\\ \;\;\;\;\;\;\;\;\;\;\;\;\forall u \in U,m \in {M_u},\\ \;\;\;\;\;\;\;\;\;\;\;m' \in {M_u},m \ne m',n > 1. \end{array}$ (28)
$\begin{array}{l} w{T_{u,m,1}} + wT{1_{u,m,1}} = {T_1},\\ \forall u \in U,m \in {M_u},n = 1. \end{array}$ (29)
$\begin{array}{l} w{T_{u,m,n}} + wT{1_{u,m,n}} = {T_n} - {T_{n - 1}},\\ \;\;\;\forall u \in U,m \in {M_u},n > 1. \end{array}$ (30)
$\begin{array}{l} w{T_{u,m,n}} \le {w_{u,m,n}} \cdot {\rm{TH,}}\\ \forall u \in U,m \in {M_u},n \in N. \end{array}$ (31)
$\begin{array}{l} wT{1_{u,m,n}} \le \left( {1 - {w_{u,m,n}}} \right) \cdot {\rm{TH,}}\\ \forall u \in U,m \in {M_u},n \in N. \end{array}$ (32)
$w{T_{u,m,n}} \ge 0,\forall u \in U,m \in {M_u},n \in N.$ (33)
$wT{1_{u,m,n}} \ge 0,\forall u \in U,m \in {M_u},n \in N.$ (34)
$\begin{array}{l} y{T_{u,m,m',n}} + yT{1_{u,m,m',n}} = {T_n} - {T_{n - 1}},\forall u \in U,\\ \;\;\;\;m \in {M_u},m' \in {M_u},m \ne m',n > 1. \end{array}$ (35)
$\begin{array}{l} y{T_{u,m,m',n}} \le \left( {\sum\limits_{k \le n - 1} {y_{u,m,m',k}^s} - \sum\limits_{k \le n - 1} {y_{u,m,m',k}^e} } \right) \cdot {\rm{TH}}\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\forall u \in U,m \in {M_u},\\ \;\;\;\;\;\;\;\;\;\;\;\;\;m' \in {M_u},m \ne m',n > 1. \end{array}$ (36)
$\begin{array}{l} yT{1_{u,m,m',n}} \le \left[ {1 - \left( {\sum\limits_{k \le n - 1} {y_{u,m,m',k}^s} - \sum\limits_{k \le n - 1} {y_{u,m,m',k}^e} } \right)} \right] \cdot {\rm{TH}}\\ \;\;\;\;\forall u \in U,m \in {M_u},m' \in {M_u},m \ne m',n > 1. \end{array}$ (37)
$\begin{array}{l} y{T_{u,m,m',n}} \ge 0,\forall u \in U,m \in {M_u},\\ \;\;\;\;m' \in {M_u},m \ne m',n > 1. \end{array}$ (38)
$\begin{array}{l} yT{1_{u,m,m',n}} \ge 0,\forall u \in U,m \in {M_u},\\ \;\;\;\;m' \in {M_u},m \ne m',n > 1. \end{array}$ (39)
$\begin{array}{l} {\rm{INV}}_{{\rm{oc}}}^{{\rm{min}}} \le {\rm{IN}}{{\rm{V}}_{{\rm{oc,n}}}} \le {\rm{INV}}_{{\rm{oc}}}^{{\rm{max}}},\\ \;\;\;\;\;\forall {\rm{oc}} \in {\rm{OC}},n \in N. \end{array}$ (40)
$\begin{array}{l} {\rm{INV}}_o^{{\rm{min}}} \le {\rm{IN}}{{\rm{V}}_{o,n}} \le {\rm{INV}}_o^{{\rm{max}}},\\ \;\;\;\;\;\forall o \in O,n \in N. \end{array}$ (41)

3) 成品油调合约束

$\begin{array}{l} PRO_{o,p}^{\min } \cdot \sum\limits_{{\rm{oc}}} {{Q_{oc,o,n}}} \le \sum\limits_{{\rm{oc}}} {{\rm{PR}}{{\rm{O}}_{oc,p}}} \cdot {Q_{oc,o,n}} \le \\ \;\;\;\;\;\;\;\;\;\;\;PRO_{o,p}^{\min } \cdot \sum\limits_{{\rm{oc}}} {{Q_{oc,o,n}}} ,\\ \;\;\;\;\;\;\;\;\;\forall o \in O,p \in {P_o},n \in N. \end{array}$ (42)
$\begin{array}{l} r_{{\rm{oc,o}}}^{{\rm{min}}}\sum\limits_{{\rm{oc'}}} {{Q_{oc',o,n}}} \le {Q_{oc,o,n}} \le r_{{\rm{oc,o}}}^{{\rm{max}}}{Q_{oc',o,n}},\\ \;\;\forall {\rm{oc}} \in {\rm{OC}},o \in O,n \in N. \end{array}$ (43)

4) 交货约束

$0 \ge {T_{l1}}\;y{d_{l,o,1}},\forall l \in L,o \in O,n = 1.$ (44)
${T_{n - 1}} \ge {T_{l1}}\;y{d_{l,o,1}},\forall l \in L,o \in O,n \ge 2.$ (45)
$\begin{array}{l} {T_n} \le {T_{l2}} + TH\left( {1 - y{d_{l,o,n}}} \right),\\ \;\;\forall l \in L,o \in O,n \in N. \end{array}$ (46)
$\begin{array}{l} {D_{\min }}y{d_{l,o,n}} \le {D_{l,o,n}} \le {D_{\max }}y{d_{l,o,n}},\\ \;\;\;\;\forall l \in L,o \in O,n \in N. \end{array}$ (47)
$\sum\limits_n {{D_{l,o,n}} \le {R_{l,o}}} ,\forall l \in L,o \in O.$ (48)
1.3 目标函数

最小化成本总和是所建立的连续时间调度模型的目标函数。成本包括4个部分:原油成本、生产装置操作成本、库存成本和延迟交货惩罚,即有

$\begin{array}{l} \;\;\;\sum\limits_n {\left( {\sum\limits_u {\sum\limits_m {s{\rm{Q}}{{\rm{I}}_{u,m,n}}} } {\rm{OpCos}}{{\rm{t}}_{u,m}} + } \right.} \\ \left. {\;\sum\limits_u {\sum\limits_m {\sum\limits_{m'} t } } {\rm{Q}}{{\rm{I}}_{u,m,m',n}}{\rm{tOpCos}}{{\rm{t}}_{u,m,m'}}} \right) + \\ \alpha \frac{{\sum\limits_n {\left( {\sum\limits_{{\rm{oc}}} {{\rm{IN}}{{\rm{V}}_{{\rm{oc}},n}}} + \sum\limits_o {{\rm{IN}}{{\rm{V}}_{o,n}}} } \right)} }}{{{n_{\max }}}} \cdot {\rm{TH + }}\\ \;\;\;\sum\limits_l {\sum\limits_o {{\beta _l}} } \cdot \left( {{R_{l,o}} - \sum\limits_n {{D_{l,o,n}}} } \right). \end{array}$ (49)

其中:OPC是常数,表示单位质量原油价格;OpCostu,m表示装置u在操作模式m稳态下的单位进料量操作成本;tOpCostu,m,m表示装置u在操作模式从m切换成m′的过渡过程中单位进料量操作成本;α是常数,表示单位质量单位时间库存成本;βl表示订单l单位质量成品油缺货惩罚值。

完整的连续时间表达炼油厂调度模型为:

$\begin{array}{l} \;\;\;\;\;\;\;\;\;\;P:\min \;f,\\ {\rm{s}}{\rm{.t}}{\rm{.式}}\left( 1 \right) - 式\left( {48} \right)。 \end{array}$
2 Lagrange分解算法 2.1 Lagrange分解方案

本文采用空间分解的Lagrange分解算法。在本文中,假设原油从原油储罐中供应,并且量足够多。因此,将后面的2个生产阶段划分成若干个部分来得到对偶子问题。装置生产部分,分解的依据是按照生产装置分解,除了常压蒸馏装置和减压蒸馏装置之外,其余每个装置作为单独的子问题。常压蒸馏装置和减压蒸馏装置共同构成一个子问题。成品油调合交货阶段作为一个单独的子问题。这样的划分能够使得子问题的规模尽可能小,有利于在后续分解算法求解时提高求解子问题的计算效率。划分方案如图 2所示。

图 2 炼油厂生产流程划分示意图

图 2中,整个炼油厂全厂生产流程被分成了9个部分:常压-减压蒸馏、催化裂化、加氢脱硫、醚化、MTBE生产、柴油加氢精制装置1号和装置2号、催化重整以及成品油调合交货等。其中前8个部分是生产装置相关,第1个部分包含了2个生产装置,其余每个部分包含了1个生产装置;第9个部分是成品油相关。

为了将完整的调度问题分解成子问题,需要松弛与多个子问题有关联的约束。在本文中,这些约束包括中间组分油质量平衡约束、生产装置模式相关约束和事件点时间约束。

中间组分油包括了生产装置(除了常压蒸馏装置)原料油和调合组分油。在连续时间完整调度模型中,生产装置组分油质量平衡约束是式(18)。

将不同装置的进料油按照对应的装置命名,如减压装置进料油为iVDU、催化裂化装置进料油为iFCCU。式(18)可以拆分成如式(50)—式(57)形式,即

${\rm{Q}}{{\rm{O}}_{{\rm{ATM}},i{\rm{VDU}},n}} = {\rm{Q}}{{\rm{I}}_{{\rm{VDU}},n}},\forall n \in N.$ (50)
${\rm{Q}}{{\rm{O}}_{{\rm{VDU}},i{\rm{FCCU}},n}} = {\rm{Q}}{{\rm{I}}_{{\rm{FCCU}},n}},\forall n \in N.$ (51)
${\rm{Q}}{{\rm{O}}_{{\rm{FCCU}},i{\rm{HDS}},n}} = {\rm{Q}}{{\rm{I}}_{{\rm{HDS}},n}},\forall n \in N.$ (52)
${\rm{Q}}{{\rm{O}}_{{\rm{HDS}},i{\rm{ETH}},n}} = {\rm{Q}}{{\rm{I}}_{{\rm{ETH}},n}},\forall n \in N.$ (53)
${\rm{Q}}{{\rm{O}}_{{\rm{FCCU}},i{\rm{MTBE}},n}} = {\rm{Q}}{{\rm{I}}_{{\rm{MTBE}},n}},\forall n \in N.$ (54)
${\rm{Q}}{{\rm{O}}_{{\rm{ATM}},i{\rm{HTU1}},n}} = {\rm{Q}}{{\rm{I}}_{{\rm{HTU1}},n}},\forall n \in N.$ (55)
$\begin{array}{l} {\rm{Q}}{{\rm{O}}_{{\rm{ATM}},i{\rm{HTU2}},n}} + {\rm{Q}}{{\rm{O}}_{{\rm{VDU}},i{\rm{FCCU2}},n}} + {\rm{Q}}{{\rm{O}}_{{\rm{FCCU}},i{\rm{HTU2}},n}} = \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\rm{Q}}{{\rm{I}}_{{\rm{HTU2}},n}},\forall n \in N. \end{array}$ (56)
$\begin{array}{l} {\rm{Q}}{{\rm{O}}_{{\rm{ATM}},i{\rm{RF}},n}} + {\rm{Q}}{{\rm{O}}_{{\rm{HTU2,}}i{\rm{RF}},n}} = \\ \;\;\;\;\;\;{\rm{Q}}{{\rm{I}}_{{\rm{RF}},n}},\forall n \in N. \end{array}$ (57)

调合组分油质量平衡约束是式(19)。

在生产过程中,催化裂化装置、加氢脱硫装置和醚化装置的生产模式保持一致。对应的是式(13)和式(14)。

连续时间调度模型中采用的是全局事件点描述,所有的装置相同的事件点均对应同一时刻。当分解成子问题时,不同的子问题之间需要相互独立,因此在每个子问题中需要有单独的表示事件点对应时间的变量。在连续时间调度模型中和时间相关的约束有:式(1)—式(3)、式(15)、式(29)—式(30)、式(35)、式(45)—式(46)。

Tnu表示子问题中的事件点对应时刻变量。变量Tnu包括了TnAVTnFCCUTnHDSTnETHTnMTBETnHTU1TnHTU2TnRFTnBD,分别表示常压-减压蒸馏子问题、催化裂化子问题、加氢脱硫子问题、醚化子问题、MTBE生产子问题、柴油加氢精制子问题1和子问题2、催化重整子问题以及成品油调合交货子问题中的事件点对应时刻。式(1)—式(3)、式(15)、式(29)—式(30)、式(35)、式(45)—式(46)可以改写成如式(58)—式(66)的等价形式,即

$T_{n - 1}^u + {T_{\min }} \le T_n^u,2 \le n < {n_{\max }}.$ (58)
$T_1^u \ge {T_{\min }}.$ (59)
$T_{{n_{\max }}}^u = {\rm{TH}}{\rm{.}}$ (60)
$\begin{array}{l} - {\rm{T}}{{\rm{T}}_{u,m,m'}}\left[ {2 - y_{u,m,m',n'}^s - y_{u,m,m',n'}^e + \sum\limits_{k = n + 1}^{n' - 1} {\left( {y_{u,m,m',k}^s + y_{u,m,m',k}^e} \right)} } \right]\\ \;\;\;\;\;\;\;\;\;\;\;\; \le T_{n'}^u - T_n^u - {\rm{T}}{{\rm{T}}_{u,m,m'}} \le \left( {{\rm{TH}} - {\rm{T}}{{\rm{T}}_{u,m,m'}}} \right)\\ \;\;\;\;\;\;\;\left[ {2 - y_{u,m,m',n'}^s - y_{u,m,m',n'}^e + \sum\limits_{k = n + 1}^{n' - 1} {\left( {y_{u,m,m',k}^s + y_{u,m,m',k}^e} \right)} } \right],\\ \;\;\;\;\;\;\;\;\;\;\forall u \in U,m \in {M_u},m' \in {M_u},m \ne m',n' \ge n. \end{array}$ (61)
$\begin{array}{l} w{T_{u,m,1}} + wT{1_{u,m,1}} = T_1^u,\\ \forall u \in U,m \in {M_u},n = 1. \end{array}$ (62)
$\begin{array}{l} w{T_{u,m,n}} + wT{1_{u,m,n}} = T_n^u - T_{n - 1}^u,\\ \;\;\forall u \in U,m \in {M_u},n > 1. \end{array}$ (63)
$\begin{array}{l} \;\;\;\;y{T_{u,m,m',n}} + yT{1_{u,m,m',n}} = T_n^u - T_{n - 1}^u,\\ \forall u \in U,m \in {M_u},m' \in {M_u},m \ne m',n > 1. \end{array}$ (64)
$T_{n - 1}^{{\rm{BD}}} \ge {T_{l1}}\;y{d_{l,o,n}},\forall l \in L,o \in O,n \ge 2.$ (65)
$\begin{array}{l} T_n^{{\rm{BD}}} \le {T_{l2}} + TH\left( {1 - y{d_{l,o,n}}} \right),\\ \;\;\;\forall l \in L,o \in O,n \in N. \end{array}$ (66)

同时需要式(67)来确保所有子问题中相同事件点对应的时刻相同,即

$T_n^{{u_1}} = T_n^{{u_2}}.$ (67)

这样式(58)—式(67)就等价于式(1)—式(3)、式(15)、式(29)—式(30)、式(35)、式(45)—式(46)。

根据图 2所示的子问题划分方案,松弛式(13)—式(14)、式(19)、式(51)—式(57)、式(67),可以将完整调度问题分解成9个子问题。经过松弛约束的问题定义成PR(λ),生产子问题定义成PAVPFCCUPHDSPETHPMTBEPHTU1PHTU2PRF,成品油调合交货子问题定义成PBD

Lagrange对偶问题定义成PD

${P_D}:{\min _\lambda }v\left( {{P_R}\left( \lambda \right)} \right)$

整个Lagrange分解算法的框架如图 3所示。算法中的难点是初始化和更新Lagrange乘子,求解子问题和对偶解可行化。

图 3 Lagrange分解算法流程图
2.2 Lagrange乘子初始化

在此分解方案中,中间组分油质量平衡约束被松弛。中间组分油既是上游生产装置的产品又是下游进料对应装置的原料。这部分松弛约束对应的Lagrange乘子可以看做是中间组分油的价格,可以通过计算等价价格来初始化。初始化相应约束形式如式(68)所示,即

${\rm{outCos}}{{\rm{t}}_u} = \frac{{{\rm{inCos}}{{\rm{t}}_u} + {\rm{OpCos}}{{\rm{t}}_u}}}{{{\rm{Yiel}}{{\rm{d}}_u}}},\;\;\forall u \in U.$ (68)

其中:inCostu表示装置u的进料等效价格;outCostu表示装置u出料的等效价格;OpCostu表示装置u的单位操作费用;Yieldu表示装置u的收率。在炼油厂调度模型中,由于已知原油价格和各生产装置的操作费用和收率,因此每种中间组分油的等效价格都可以通过式(68)计算。将计算得到的值作为相应Lagrange乘子的初始值。

2.3 Lagrange乘子的更新

在已有的Lagrange分解算法的研究中,有很多经典的乘子更新方法,例如次梯度方法[17, 18]、割平面方法[19, 20]和boxstep方法[21]等。Mouret等[5]提出了一种结合这3种方法的混合方法来更新乘子。在本文中,在Mouret等提出的方法基础上,结合指数加权次梯度[22]方法来加快乘子的收敛速度。在Lagrange分解算法第(IT+1)次迭代时,乘子更新为优化问题$\hat P_{\rm{D}}^{{\rm{IT + 1}}}$的解,即

$\begin{array}{l} \;\;\;\;\;\;\;\;\;\;\hat P_{\rm{D}}^{{\rm{IT + 1}}}:\min \;\eta ,\\ {\rm{s}}{\rm{.}}t.\\ \;\eta \le {f^{{\rm{it}}}} + \sum\limits_{{\lambda _n} \in \Lambda } {\sum\limits_n {{\lambda _n}{g_{{\rm{it}}}},{\rm{it = 1,}} \cdots {\rm{,IT}}} .} \end{array}$ (69)
$\begin{array}{l} {\lambda _n} = \lambda _n^{{\rm{IT}}} + \mu \left( {{f^{{\rm{up}}}} - {f^{{\rm{low}}}}} \right)\frac{{{g^{{\rm{IT}}}}}}{{{{\left\| {{g^{{\rm{IT}}}}} \right\|}^2}}} + \delta ,{\lambda _n} \in \Lambda ,\\ \;\;\;\;\;\;\;\;{g^{{\rm{IT + 1}}}} = \left( {1 - \gamma } \right){g^{{\rm{IT}}}} + \gamma {g^{{\rm{IT + 1}}}}'. \end{array}$ (70)
$\begin{array}{l} \;\;n \in N,\mu \in \left[ { - \infty ,\bar \mu } \right],\\ \delta \in \left[ { - \bar \delta ,\bar \delta } \right],\gamma \in \left[ {0,1} \right]. \end{array}$ (71)

其中:gIT+1′表示由子问题的解计算得到的次梯度;gIT+1表示经过加权得到的次梯度,用于迭代乘子;上标IT表示算法迭代次数;flowfup表示分解前原问题目标函数值的下界值和上界值;${\bar \mu }$和 ${\bar \delta }$是给定的参数;Λ是Lagrange乘子的集合。

2.4 子问题的解

当通过分解得到子问题之后,子问题的求解效率关系到整个Lagrange分解算法的效率。在完整调度模型中,对于调合组分油是没有上下限的限制,仅仅有相应质量平衡约束。将完整模型分解成子问题时,这部分约束被松弛。所得到的成品油调合交货子问题是线性规划问题,如果直接求解,这些调合组分油的解很可能取0或者是最大值。这样的解和最优解相比有很大差距。从生产流程分析,由于原油进料有上下限限制,每个生产装置的产品收率都是已知的。这样每种调合组分油的流量上下限可以确定如式(72)所示,即

$\begin{array}{l} QI_{{\rm{oc}}}^{\min } \le \sum\limits_u {Q{I_{u,{\rm{oc,t}}}}} \le {\rm{QI}}_{{\rm{oc}}}^{\max },\\ \;\;\;\;\forall {\rm{oc}} \in {\rm{OC}},n \in N. \end{array}$ (72)

增加式(72)并不会改变模型的解,同时可以能够提高成品油调合交货子问题的解的质量。

2.5 可行化方案

通常情况下,分解得到的子问题解不是完整问题的解,被松弛的约束不能满足。对于好的 Lagrange 乘子而言,可以认为子问题解中的0-1变量的取值接近于最优解。在此基础上,本节设计一种启发式算法来获得原问题的可行解,使得能够满足原问题的所有约束,包括所松弛的中间组分油的质量平衡约束在内。在启发式算法中,一部分0-1变量根据对偶子问题的解固定下来,装置分为若干组,每次迭代都会固定其中一组装置的操作状态,经过几次迭代完成所有装置的遍历。具体由式(73)—式(75)表示,即

$\begin{array}{l} \;\;\;\;\;{w_{u,m,n}} = {\left( {{w_{u,m,n}}} \right)^{{\rm{IT}}}},\\ \forall u \in {U^{{\rm{IT}}}},m \in {M_u},n \in N. \end{array}$ (73)
$\begin{array}{l} \;\;y_{u,m,m',n}^s = {\left( {y_{u,m,m',n}^s} \right)^{{\rm{IT}}}},\\ \forall u \in {U^{{\rm{IT}}}},m \ne m',n \in N. \end{array}$ (74)
$\begin{array}{l} \;\;y_{u,m,m',n}^e = {\left( {y_{u,m,m',n}^e} \right)^{{\rm{IT}}}},\\ \forall u \in {U^{{\rm{IT}}}},m \ne m',n \in N. \end{array}$ (75)

其中:UIT表示在第IT次迭代中固定离散变量的生产装置集合,经过几次迭代所有生产装置对应的离散变量会轮流被固定;(wu,m,n)IT、(yu,m,m′,ns)IT和(yu,m,m′,ne)IT的值是由第IT次迭代时对应子问题的解确定的,其中上标中的s和e分别表示过渡过程开始和结束。

这样可行化问题的形式为

${P_{{\rm{fea}}}}:\min \;f,$

s.t. 式(1)—式(48),式(73)—式(75).

求解Pfea就能得到完整问题P的可行解。

2.6 算法停止准则

Lagrange分解算法的停止准则是基于对偶间隙。对偶间隙必须小于给定值$\varepsilon $,即

${\rm{gap}} \le \varepsilon {\rm{.}}$ (76)

用完整问题P的下界值flow和上界值fup来计算对偶间隙,即

${\rm{gap}} = \frac{{{f^{{\rm{up}}}} - {f^{{\rm{low}}}}}}{{{f^{{\rm{up}}}}}}.$ (77)
2.7 完整算法流程

综合上述内容,完整的针对炼油厂全场调度离散时间模型的多阶段Lagrange分解算法步骤如下:

Step 1 将完整模型P分解成9个子问题:PAVPFCCUPHDSPETHPMTBEPHTU1PHTU2PRFPB&D

Step 2 定义可行化问题Pfea

Step 3 根据式(68)初始化Lagrange乘子;

Step 4 求解子问题:将子问题的目标函数值求和更新完整问题P的下界值flow

Step 5 求解可行化问题Pfea:通过可行化问题的解计算原问题的目标函数值并更新完整问题P的目标函数值的上界fup

Step 6 根据式(77)计算对偶间隙:如果停止准则满足,则算法终止。否则,求解问题$\hat P_{\rm{D}}^{{\rm{IT + 1}}}$,更新Lagrange乘子并返回Step 4。

3 案例仿真

为了测试本章中所提出的多阶段Lagrange分解算法效果,在本节中用3个具有不同调度周期长度和不同订单数量的案例进行仿真。3个案例分别用直接求解、普通Lagrange分解方法和本文提出的Lagrange分解方法求解。在普通Lagrange分解方法中:乘子初始化取值全为0;乘子迭代采用的是次梯度方法;成品油调合交货子问题中不加入辅助约束;对偶可行化方法采用本文提出的启发式算法。这3个案例是基于图 1所表示的炼油厂生产流程。3个案例中的调度周期时长被均分成长度为8 h的时间间隔。在乘子迭代中,μ取值为1,δ取值为100,γ取值为0.6。算法停止准则中对偶间隙定为8%。求解模型所用软件是GAMS版本 24.2.2 ,求解器是CPLEX12.6.0.0 。电脑配置信息是CPU:Intel(R) Xeon(R) CPU E5-2609 v2 @2.5 GHz,内存:32.0 GB。

3个仿真案例的规模如表 1所示。从表 1中可以看出3个案例中调度时间长度依次增长;案例2和案例3具有相同的订单数量。模型中所需的事件点个数、模型中的约束数和变量数如表 2所示。

表 1 仿真案例规模
案例调度时间长度/h订单数/个
11283
21925
32405
表 2 调度模型规模
案例事件点个数约束数 变量数 离散变量数
177 3605 140836
2910 7306 9161 244
3910 7306 9161 244
3.1 案例1

案例1的成品油订单信息见表 3,案例求解结果见表 4。Lagrange分解算法每次迭代结果如图 4所示。

表 3 例1成品油订单信息(交货时间和交货量)
订单开始时间/h结束时间/h 汽油/t 柴油/t
JIV93JIV97GIII90GIII93GIII97 GIII0GIIIM10GIV0
10482707302702701 2509301 3302 700
2489613073001301 0709301 6001 960
3801284004501301309306701 3302 330
图 4 案例1连续时间调度模型Lagrange分解算法原问题目标函数上下界值收敛趋势

通过表 4的仿真结果可知,本文提出的 Lagrange 分解算法能够在9次迭代之后得到最优解,求解时间是完整模型求解时间的51.9%。

表 4 案例1仿真结果统计
算法目标函数对偶间隙求解时间迭代次数
%s
直接求解165 516 922.3-264.5-
普通Lagrange分解算法165 516 922.334.33497.230
本文Lagrange分解算法165 516 922.37.52137.29
3.2 案例2

案例2的成品油订单信息见表 5。案例求解结果见表 6。Lagrange分解算法每次迭代结果如图 5所示。

表 5 案例2成品油订单信息(交货时间和交货量)
订单开始时间/h结束时间/h 汽油/t 柴油/t
JIV93JIV97GIII90GIII93GIII97 GIII0GIIIM10GIV0
104820050004001 0006001 0001 500
2328040030010001 1504001 6501 500
3721204005002001505003008003 000
411216020060015008505001 4002 000
514419210080002501 2006001 5002 000
表 6 案例2仿真结果统计
算法目标函数对偶间隙求解时间迭代次数
%s
直接求解230 703 098.4-1 726.7-
普通Lagrange分解算法230 703 098.429.331 557.230
本文Lagrange分解算法230 703 098.47.64348.515
图 5 案例2连续时间调度模型Lagrange分解算法原问题目标函数上下界值收敛趋势

与案例1相比,案例2的订单个数增长了66.7%,所需事件点个数增加了28.6%。相应的模型中约束数增加了45.8%,变量数增加了34.6%,离散变量数增加了48.8%。由于约束和变量个数的增加,案例2求解完整模型耗时是案例1的6.5倍。而采用Lagrange分解算法能够在15次迭代之后得到最优解,采用Lagrange分解算法求解时间是完整模型求解时间的20.2%。案例2采用 Lagrange 分解算法求解时间是案例1所需时间的2.5 倍。

3.3 案例3

案例3的成品油订单信息见表 7。案例求解结果见表 8。Lagrange分解算法每次迭代结果如图 6所示。

表 7 案例3成品油订单信息(交货时间和交货量)
订单开始时间结束时间 汽油/t 柴油/t
hhJIV93JIV97GIII90GIII93GIII97 GIII0GIIIM10GIV0
10643004002002001 3006001 8002 200
2481125005002502001 2009002 0002 400
3961603007001004001 0006002 2002 800
41442003004002003008008001 7002 400
51922404005003003001 0007002 2002 500
表 8 案例3仿真结果统计
算法目标函数对偶间隙求解时间迭代次数
% s
直接求解298 545 020.2-3 473.9-
普通Lagrange分解算法298 545 020.231.732 785.230
本文Lagrange分解算法299 393 430.07.64313.512
图 6 案例3连续时间调度模型Lagrange分解算法原问题目标函数上下界值收敛趋势

与案例2相比,案例3的订单个数和所需事件点个数都不变。相应的模型中约束数、变量数和离散变量数也都相同。由于调度时间长度的增加,案例3求解完整模型耗时是案例2的2.0倍。而采用Lagrange分解算法能够在12次迭代之后得到解,该解和用完整模型得到的最优解相比有0.28%的差距,可以认为是近似最优解。采用Lagrange分解算法求解时间是完整模型求解时间的9.0%。由于迭代次数较少,案例3采用Lagrange分解算法求解时间是案例2所需时间的90.0%。

4 结 论

本文针对炼油厂连续时间调度模型的求解时间随着问题规模的增大显著增长的问题,设计了一种Lagrange分解算法。通过松弛中间组分油质量平衡约束、生产装置模式相关约束和事件点时间约束将原完整问题分解成9个子问题。为了加快算法收敛效率,通过经济价格等价计算设定Lagrange乘子初始值;采用一种混合方法更新乘子。为了得到原问题的可行解,设计了一个启发式算法。炼油厂连续时间调度模型的求解时间受到调度问题中订单数、事件点个数和调度时间长度的影响。随着问题规模的增大,调度模型的求解时间显著增长。采用本文所提出的Lagrange分解算法,能够有效地得到调度问题的最优解或近似最优解。Lagrange分解算法的求解时间比完整模型求解所需时间短,并且随着问题规模的增大优势会更明显。从3个案例的求解结果可以看出,原问题目标函数值上下界的收敛规律相似,上界值收敛很快,下界值要慢很多。从仿真结果上看,Lagrange分解算法中停止准则设为对偶间隙小于8%足以让算法找到调度问题的最优解或者近似最优解。案例仿真结果表明所提出的Lagrange分解算法的有效性,尤其是在大规模问题中算法优势更加明显。

参考文献
[1] Pinto J M, Joly M, Moro L F L. Planning and scheduling models for refinery operations[J]. Computers & Chemical Engineering, 2000, 24(9-10), 2259-2276.
[2] Göthe-Lundgren M, Lundgren J T, Persson JA. An optimization model for refinery production scheduling[J]. International Journal of Production Economics, 2002, 78(3), 255-270.
[3] JIA Zhenya, Ierapetritou M. Efficient short-term scheduling of refinery operations based on a continuous time formulation[J]. Computers & Chemical Engineering, 2004, 28(6-7), 1001-1019.
[4] LUO Chunpeng, RONG Gang. Hierarchical approach for short-term scheduling in refineries[J]. Industrial & Engineering Chemistry Research, 2007, 46(11), 3656-3668.
[5] Mouret S, Grossmann I E, Pestiaux P. A new Lagrangian decomposition approach applied to the integration of refinery planning and crude-oil scheduling[J]. Computers & Chemical Engineering, 2011, 35(12), 2750-2766.
[6] CAO Cuiwen, GU Xingsheng, XIN Zhong. A data-driven rolling-horizon online scheduling model for diesel production of a real-world refinery[J]. AIChE Journal, 2013, 59(4), 1160-1174.
[7] Shah N K, LI Zukui, Ierapetritou M G. Petroleum refining operations:Key issues, advances, and opportunities[J]. Industrial & Engineering Chemistry Research, 2011, 50(3), 1161-1170.
[8] Joly M. Refinery production planning and scheduling:The refining core business[J]. Brazilian Journal of Chemical Engineering, 2012, 29(2), 371-384.
[9] SHI Lei, JIANG Yongheng, WANG Ling, et al. Refinery production scheduling involving operational transitions of mode switching under predictive control system[J]. Industrial & Engineering Chemistry Research, 2014, 53(19), 8155-8170.
[10] Terrazas-Moreno S, Trotter P A, Grossmann I E. Temporal and spatial Lagrange an decompositions in multi-site, multi-period production planning problems with sequence-dependent changeovers[J]. Computers & Chemical Engineering, 2011, 35, 2913-2928.
[11] Neiro S M, Pinto J M. Langrange an decomposition applied to multiperiod planning of petroleum refineries under uncertainty[J]. Latin American Applied Research, 2006, 36(4), 213-220.
[12] Shah N, Saharidis G, JIA Zhenya, et al. Centralized-decentralized optimization for refinery scheduling[J]. Computers & Chemical Engineering, 2009, 33(12):2091-2105.
[13] TANG Lixin, Luh P B, LIU Jiyin, et al. Steel-making process scheduling using Lagrangian relaxation[J]. International Journal of Production Research, 2002, 40(1), 55-70.
[14] LI Zukui, Ierapetritou M. Production planning and scheduling integration through augmented Lagrangian optimization[J]. Computers & Chemical Engineering, 2010, 34, 996-1006.
[15] JIANG Yongheng, Rodriguez M A, Harjunkoski I, et al. Optimal supply chain design and management over a multi-period horizon under demand uncertainty. Part Ⅱ:A Lagrangean decomposition algorithm[J]. Computers & Chemical Engineering, 2014, 62, 211-224.
[16] Knudsen B R, Grossmann I E, Foss B, et al. Lagrangian relaxation based decomposition for well scheduling in shale-gas systems[J]. Computers & Chemical Engineering, 2014, 63, 234-249.
[17] Held M, Karp R M. The traveling-salesman problem and minimum spanning trees:Part Ⅱ[J]. Mathematical Programming, 1971, 1(1):6-25.
[18] Held M, Wolfe P, Crowder H P. Validation of subgradient optimization[J]. Mathematical Programming, 1974, 6(1), 62-88.
[19] Cheney E W, Goldstein A A. Newton's method for convex programming and Tchebycheff approximation[J]. Numerische Mathematik, 1959, 1(1), 253-268.
[20] Kelley J J E. The cutting-plane method for solving convex programs[J]. Journal of the Society for Industrial & Applied Mathematics, 1960, 8(4), 703-712.
[21] Marsten R E, Hogan W W, Blankenship J W. The boxstep method for large-scale optimization[J]. Operations Research, 1975, 23(3), 389-405.
[22] Baker B M, Sheasby, J. Accelerating the convergence of subgradient optimization[J]. European Journal of Operational Research, 1999, 117, 136-144.