半无限minimax优化问题是一类非常重要的优化问题, 有着广泛的应用背景, 例如工程设计、最优控制、金融工程等领域的很多问题可以归结为求解这类优化问题, 很多学者对此进行了研究, 获得了丰富的研究成果(如文献[1-12]).由于其特殊的结构, 大多数传统的方法不再适用, 离散化方法是求解这类型问题的重要数值方法之一(如文献[1-8]).半无限minimax问题具有如下形式
其中指标集$Y = [ 1, 2 , \cdots m ]$, $ f : R ^ { n } \times Y \rightarrow R$, 关于$x, y$都连续可微关于$x$连续可微$g : R ^ { n } \times [ 0, 1 ] \rightarrow R$.为了方便起见, 记问题(1.1)的可行集$X$和水平集$l \left( x ^ { 0 } , \Omega \right)$:
再记:对于任意的$x \in R ^ { n } , F _ { Y } ( x ) = \max\limits _ { y \in Y } f ( x , y )$.对于$\hat { Y } \subset Y , F _ { \hat { Y } } ( x ) = \max\limits _ { y \in \hat { Y } } f ( x , y )$的方向导数
离散化方法的主要思想通过不断离散连续变量的连续区间来逼近约束函数, 将求解半无限规划问题转化为求解其离散后的一系列有限约束优化问题.将$\Omega$离散成有限集:$\Omega _ { \mathrm { E } } = \left\{ 0 , \frac { 1 } { q } , \frac { 2 } { q } , \cdots \frac { q - 1 } { q } , 1 \right\}$, 其中$q$反映了离散水平, $q$越大离散水平越好.定义$\Omega$和$\Omega_E$之间的Hausdorff距离为$\operatorname { dist } \left( \Omega _ { \mathrm { E } } , \Omega \right) = \max\limits _ { \omega \in \Omega } \min\limits _ { \omega \in \Omega _ { \mathcal { E } } } \left\| \omega _ { \mathrm { E } } - \omega \right\|$, 其中集列$\left\{ \Omega _ { \mathrm { E } } ^ { q _ { i } } \right\} _ { i \in N ^ { + } }$满足条件
基于离散化方法求解原问题(1.1)可归结为求解一系列具有如下形式的minimax离散化问题:
在一定的条件下, 当$\operatorname { dist } \left( \Omega _ { \mathrm { E } } , \Omega \right)\to 0$时式(1.3)的最优解趋向于原问题(1.1)的最优解.当$q$非常大的时候, 问题(1.3)的约束个数非常多, 求解的成本也会很高, 如何设计高效的算法求解问题(1.3)是解决半无限minimax问题的一个关键.文献[5]提出了一种求解半无限规划问题的超线性收敛的模松弛SQP算法, 每次迭代只需要求解一个QP子问题就可以获得搜索方向, 遗憾的是上述算法要求初始点可行, 而通常求解可行点的计算量很大.为了克服这一问题, 文献[6]提出了一种初始点任意的模松弛强次可行SQP算法.另外, 在模松弛强次可行方向法中常利用线搜索来确定步长, 传统的线搜索方法都要求目标函数值严格下降, 其缺点是当迭代点“陷入很窄的峡谷时”, 常常会导致小步长或出现锯齿现象, 而采用非单调线搜索技术可以克服这些缺点(如文献[13-15]).文献[9]提出一种约束积极集识别技术, 可以对积极集进行精确识别和估计来减少计算量.文献[16]提出了一种求解无约束优化问题的有限记忆BFGS修正规则(简称L-BFGS), L-BFGS修正方法无需存贮近似海森矩阵$H_k$, 从而大大减少了计算机存储量, 提高运行效率, 对大规模的优化问题更有效.受文献[5-9]启发, 本文结合离散化方法和模松弛强次可行SQP方法, 基于新型约束积极集识别技术, 采用非单调搜索和有限记忆L-BFGS修正方法更新$H_k$, 提出了一种新的求解半无限minimax问题的非单调SQP算法.
定义2.1 点$x$称为QP子问题(2.1)的稳定点(KKT点), 如果存在乘子向量$\lambda _ { \omega } , \left( \omega \in \Omega _ { \mathrm { E } } \right)$; $\gamma _ { y } , ( y \in Y )$使得
以下给出了式(1.3)的拉格朗日函数$L(x, \lambda, y)$, 可行集$X_E$和有效集$Y(x)$和$\Omega_E(x)$:
由式(2.1), (2.2), 定义如下的约束积极集识别函数
其中$e = ( 1, 1 , \cdots , 1 ) ^ { T } \in R ^ { n } , \quad 0 < \delta < 1$, 由文[5]可得到问题(1.3)的两个积极约束识别集
为了叙述方便, 对任给的$x ^ { k } \in R ^ { n } , \hat { Y } \subset Y ( x , \lambda , \gamma ) , \Omega _ { k } \subset \Omega _ { \mathrm { E } } ( x , \lambda , \gamma )$, 定义记号如下
构造如下的QP子问题
其中参数$r _ { 0 } , r , r _ { \omega } \left( \omega \in \Omega _ { \mathrm { E } } \right) , \theta$都是正的常数, $\sigma_k$是随迭代调整的正参数, $H_k$是$\nabla _ { x x } ^ { 2 } f \left( x ^ { k } , y \right)$的近似矩阵.我们引入一个重要的项$-\varepsilon_kz^2$是为了确保当$d^k\to 0$有$z_k\to 0$.以下引理说明QP子问题具有良好的性质.
假设2.1 对于任意的$y \in Y , \omega \in \Omega$, 函数$f ( \cdot , y )$和$g ( \cdot , \omega )$在可行集上至少一阶连续可微.
假设2.2 对于任意的$x\in X_E$, 弱MFCQ成立, 即存在$d\in R^n$使得$\nabla _ { x } g ( x , \omega ) ^ { T } d < 0$, 对于所有的$\omega\in\Omega(x)$成立.
引理2.1 设假设2.1和假设2.2成立, $x ^ { k } \in X , H _ { k }$是对称正定, 若$(z_k, d^k)$是子问题(2.5)的最优解, 则
(1) $r _ { 0 } z _ { k } + \frac { 1 } { 2 } \left( d ^ { k } \right) ^ { T } H _ { k } d ^ { k } \leq 0 , z _ { k } \leq 0 $;
(2) $ z _ { k } = 0 \Leftrightarrow d ^ { k } = 0 $;
(3) $z _ { k } = 0 \Rightarrow x ^ { k } \text { 是问题 } ( 1.3 ) \text { 的一个稳定点 } $;
(4) 如果$x ^ { k } \notin X \text { 则 } z _ { k } < 0 ;$
(5) 若$d^k\neq 0$, 则$z_k<0$, $d^k$为问题(1.3)在$x^k$处的可行下降方向.
证 (1), (2)的证明类似于文献[9](引理2.2, 2.3)其余证明可参考文献[7](引理2.1).
在文献[13]中Zhang和Hager提出了新的非单调搜索技术, 受文献[13]启发, 本文对其进行了改进, 具体形式如下
其中$\gamma \in [ 0, 2 ) , \delta \in \left( 0 , \frac { 1 } { 2 } \right) , s _ { k } = \frac { - \nabla _ { x } f ( x , y ) ^ { T } d ^ { k } } { L _ { k } \left\| d ^ { k } \right\| ^ { 2 } } , \quad t \in ( 0, 1 ) , \alpha _ { k }$为$\left\{ s _ { k } , t s _ { k } , t ^ { 2 } s _ { k } , \cdots \right\}$中满足上式的最大者, $L_k$的调整规则参见文献[15].
由于采用离散化技术后, 问题(1.3)的约束个数较多, 导致算法求解的计算量增加, 效率降低, 如何降低计算量成为本文的关键.在SQP算法中通常用BFGS公式更新$H_k$, 文[16]提出一种新的有限记忆L-BFGS更新规则, 它最显著的特点是不需要存储$H_k$, 对给定的$H_k$和非负整数$m$, 利用前$m$步的信息对$H_0$进行修正$m$次得到$H_k$, 仅存贮$m+1$个向量组就能计算$H_{k+1}$, 从而降低了算法对计算量和存储量的要求, 本文将其应用到算法设计中更新$H_k$.
算法A (求解离散化问题(1.3))
步1 初始化.取适当大正整数$q$, 将区间[0, 1]离散成有限集$\Omega_k$选取参数$r _ { 0 } , r , r _ { \omega }, \theta \in ( 0, 1 ), \sigma _ { k }, \eta \in ( 0, 1 )$.选取初始可行点$x^0\in X_E$, 初始对称正定矩阵$H_0$, $\left( \lambda ^ { 0 }, \gamma ^ { 0 } \right)^{T}=(1, 1, \cdots, 1)^{T}\in R^{m+q+1}$, 并令$k:=0$.
步2 由(2.3)式计算$\rho \left( x ^ { k } , \lambda ^ { k } , \gamma ^ { k } \right)$, 由(2.4)式生成积极约束识别集$Y \left( x ^ { k } , \lambda ^ { k } , \gamma ^ { k } \right)$和$\Omega_E \left( x ^ { k } , \lambda ^ { k } , \gamma ^ { k } \right)$.如果$\rho \left( x ^ { k } , \lambda ^ { k } , \gamma ^ { k } \right)=0$, 则$x^k$是问题(1.3)的一个稳定点, 停止.否则进入步3.
步3 计算搜索方向.对当前迭代点$x^k$求解QP子问题(2.5)得到一个KKT点对$(z_k, d^k)$.如果$d^k=0$, 则$x^k$是问题(1.3)的一个稳定点, 停止.否则进入步4.
步4 求搜索步长.由新的非单调搜索(2.6)获得步长$\alpha_k$.
步5 令$x ^ { k + 1 } = x ^ { k } + \alpha _ { k } d ^ { k }$, 对称正定矩阵$k _ { k+1 }$可按文献[16]中L-BFGS更新规则得到, $\sigma _ { k+1 }$由$\sigma _ { k + 1 } = \min \left\{ \bar { \sigma } , \sigma _ { 1 } \left\| d ^ { k } \right\| ^ { \xi } \right\} , \forall k \geq 1$获得, 其中$\bar { \sigma } , \sigma _ { 1 } , \xi$均为正常数.置$k=k+1$, 返回步2.
对于半无限minimax问题的离散化方法, 由文献[7]可知, 若下面的假设2.3成立, 且存在$x^0\in X$, 使得$l \left( x ^ { 0 } , \Omega _ { \mathrm { E } } ^ { q _ { 0 } } \right)$紧, 那么离散化问题(1.3)解序列的某个子列$\left\{ \tilde { \bf { x } } _ { q _ { i } } \right\} _ { i \in N ^ { k } } , N ^ { k } \subseteq N , k \in R$收敛于原问题(1.1)的解$\tilde { \bf { x }}$.
假设2.3 集列$\left\{ \Omega _ { \mathrm { E } } ^ { q _ { i } } \right\} _ { i \in N ^ { + } }$满足条件(1.2), 且$ \Omega _ { \mathrm { E } } ^ { q _ { i } } \subseteq \Omega _ { \mathrm { E } } ^ { q _ { i + 1 } } , i \in N ^ { + }$成立.
下面给出求解原问题(1.1)的算法
算法B
初始步 给定参数: $\left\{ \pi_i \right\} _ { i \in N ^ { + } }$满足$0 < \varepsilon < \pi _ { i + 1 } < \pi _ { i } , \forall i \in N ^ { + } $和$\lim\limits _ { i \rightarrow \infty } \pi _ { i } = 0$, $q _ { 0 } \in N ^ { + } \backslash \{ 0 \} , \varepsilon = 10 ^ { - 4 }$以及算法A的初始化参数.
步1 选取$x^0\in X$离散集$\Omega _ { \mathrm { E } } ^ {c_{0}} \subset \Omega$, 满足$\left| \Omega _ { \mathrm { E } } ^ { q _ { 0 } } \right| < \infty$.令$i=0$.
步2 利用算法A求解离散化问题(1.3), 其最优解记为$\tilde { x } _ { q _ { i } }$.
步3 如果$\operatorname { dist } \left( \Omega _ { \mathrm { E } } ^ { q _ { i } } , \Omega \right) \leq \varepsilon$, 则算法停止, 否则, 取更为精细的离散集$\Omega _ { \mathrm { E } } ^ { q _ { i + 1 } } \subset \Omega$, 使得$\left\{ \omega _ { q _ { i } } \right\} \cup \Omega _ { \mathrm { E } } ^ { q _ { i } } \subset \Omega _ { \mathrm { E } } ^ { q _ { i + 1 } }$, 其中$q _ { i } \in \Omega _ { \mathrm { E } } ^ { q _ { i + 1 } } , \left| \Omega _ { \mathrm { E } } ^ { q _ { i + 1 } } \right| < \infty$和$\operatorname { dist } \left( \Omega _ { \mathrm { E } } ^ { q _ { i + 1 } } , \Omega \right) \leq \operatorname { dist } \left( \Omega _ { \mathrm { E } } ^ { q _ { i } } \right) \leq \pi _ { i + 1 }$, 且满足$g \left( \widetilde { x } _ { q _ { i } } , \omega \right) = \max\limits _ { \omega \in \Omega _ { q _ { i + 1 } } } g \left( \widetilde { x } _ { q _ { i } } , \omega \right)$.
步4 令$x ^ { 0 } = \tilde { x } _ { q _ { i } } , i = i + 1$, 转步2.
本节对算法B进行收敛性分析, 首先给出如下的记号
对于算法B, 定义如下的水平集
假设3.1 水平集$l \left( x ^ { 0 } , \Omega _ { \mathrm { E } } ^ { q _ { 0 } } \right)$有界, 且存在Lipschitz常数$l_{g_x}$和常数$l_{g_w}$使得下式成立.
假设3.2 问题(1.1)在点$\tilde { x } \in X$时有$\left\{ \nabla _ { x } g ( \tilde { x } , \omega ) , \omega \in \Omega ( \tilde { x } ) \right\}$线性无关.
引理3.1[7]若$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$为算法B生成的迭代点列, 则序列$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$存在聚点$\tilde { x }$.
引理3.2 令$\left\{ \widetilde { x } _ { q _ { i } } \right\} _ { N ^ { k } } , N ^ { k } \subset N ^ { + } , \left| N ^ { k } \right| = \infty$是算法B生成的序列$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$存在聚点$\tilde { x }$的一个收敛于$\tilde { x }$的子列, 则$\left\{ \varphi _ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right) \right\} _ { i \in N ^ { k } }$收敛, 且$\lim\limits _ { i \rightarrow \infty , i \in N ^ { k } } \varphi _ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right) = \Phi ( \tilde { x } )$.
证 先证$0 \leq \Phi ( x ) - \varphi _ { q _ { i } } ( x ) \leq l _ { g _ { w } } \operatorname { dist } \left( \Omega _ { \mathrm { E } } ^ { q _ { i } } , \Omega \right) , \forall x \in l \left( x ^ { 0 } , \Omega _ { \mathrm { E } } ^ { q _ { 0 } } \right) , i \in N ^ { + }$.
若$\Phi(x)=0$, 则以上的关系式显然成立.若$\Phi(x)>0$, 取$\omega\in X$, 则存在$\omega _ { q _ { i } } \in \Omega _ { \mathrm { E } } ^ { q }$, 满足$\left\| \omega - \omega _ { q _ { i } } \right\| \leq l _ { g _ { o } } \operatorname { dist } \left( \Omega _ { \mathrm { E } } ^ { q _ { i } } , \Omega \right)$, 因而有
再证$\left\{ \varphi _ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right) \right\} _ { i \in N ^ { k } }$收敛.由式(3.3)可知$\varphi _ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right) \leq \Phi \left( \tilde { x } _ { q _ { i } } \right)$且有
由算法B可知$\lim\limits _ { \pi _ { i } \rightarrow \infty } \pi _ { i } = 0$.根据式(3.3)可知$\lim\limits _ { i \rightarrow \infty , i \in N ^ { k } } \varphi _ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right) = \Phi \left( \tilde { x } _ { q _ { i } } \right)$.
引理3.3 令$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \in N ^ { k } } , N ^ { k } \subseteq N ^ { + } , \left| N ^ { k } \right| = \infty$是算法B生成的序列$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$的一个收敛于$\tilde { x }$的子列, 则对于$\tilde { \omega } \in \Omega ( \tilde { x } )$, 存在点列$\left\{ \omega _ { q _ { i } } \right\} _ { i \in N ^ { k } }$使得$q _ { i } \in \Omega _ { \mathrm { E } } ^ { q _ { i } \left( \tilde { x } _ { q _ { i } } \right) }$对于充分大的满足$\omega _ { q _ { i } } \rightarrow \tilde { \omega }$.
证明可参考文献[7]的引理6.1.4.
定理3.1 若假设3.1, 3.2成立, 如果$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$是由算法B生成的迭代点列, 则存在$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$的一个聚点$\tilde { x }$是半无限minimax问题(1.1)的KKT点, 即算法B是收敛的.
证 类似于文献[8]中定理2.1的证明过程, 由引理3.1可知序列$\left\{ \tilde { x } _ { q _ { i } } \right\} _ { i \subset N ^ { + } }$存在聚点$\tilde { x }$, 取其收敛于$\tilde { x }$的子列, $\left\{ \tilde { x } _ { q _ { i } } \right\} , N ^ { k } \subseteq N ^ { + } , \left| N ^ { k } \right| = \infty$.由算法B可知$\tilde { x }_{qi}$是问题(1.3)的KKT点.因此, 由Caratheodory定理, 对于$m = n + 1 , i \in N ^ { k } \text { 和 } \omega _ { q _ { i } } \in \Omega _ { \mathrm { E } } ^ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right)$, 有
由文献[7]中引理2.3.2, 3.3.3可知, 乘子序列$\left\{ \xi _ { i } ^ { j } \right\} _ { i \in N ^ { k } } , j = 1, 2 , \cdots , m$有界, 故存在子列收敛于$\tilde { \xi } ^ { j } , j = 1, 2 , \cdots , m$, 不妨记为原数列.再记$J(\tilde { x })$为$\Omega_E^{qi}(\tilde { x })$的指标集, 由引理3.3知, 对于$\tilde { \omega } ^ { j } \in \Omega ( \tilde { x } ) , j \in J ( \tilde { x } )$, 存在点列$\left\{ \omega _ { q _ { i } } ^ { j } \right\} _ { i \in N ^ { k } }$使得$\omega _ { q _ { i } } ^ { j } \in \Omega _ { \mathrm { E } } ^ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right)$.对充分大的$i\in N^k$满足, 即若$i\in N^k$充分大, 则通过整理可使得$K(=1, 2, \cdots, m)$中前$l ( = | J ( \tilde { x } ) | )$个对应的$\Omega _ { \mathrm { E } } ^ { q _ { i } } \left( \tilde { x } _ { q _ { i } } \right)$中的指标$\omega _ { q _ { i } } ^ { j } \rightarrow \tilde { \omega } ^ { j } , j = l + 1 , \cdots , m$.此外由$|K|<\infty$和$\Omega$是紧集可知, 序列$\left\{ \omega _ { q _ { i } } ^ { j } \right\} _ { i \in N ^ { k } } , j = l + 1 , \cdots , m$收敛于$\tilde { \omega } ^ { j } \in \Omega , j = l + 1 , \cdots , m$.在式(3.4)-(3.7)中令$i \rightarrow \infty , i \in N ^ { k }$, 再结合引理3.2可得
由式(3.11)可知$\tilde { x }\in X$.结合文[7]中定义6.1.2, 并由式(3.8)-(3.11)可知, $\tilde { x }$是半无限minimax问题(1.1)的KKT点.
本节将对算法B在MATLAB2016a上编程序进行数值实验.其中P1, P2, P3三个算例来源于文献[4], 三个问题如下
将算法与文献[7]中的算法C (文献[7]第二章2.6节的算法)进行对比.参数选取$q = 100 , r _ { 0 } = 5 , r = 1 , r _ { \omega } = 0.01 , \theta = 0.1 , \gamma = 1 , \delta = 0.38 , L _ { k } = 1 , t = 0.87 , L _ { k } = 1$, $\bar { \sigma } = 0.01 , \sigma _ { 1 } = 100 , \zeta = 0.5 , H _ { 0 } = I$.算法的终止准则为$|z_k|<\leq 10^{-4}$, 其中Ni为迭代次数, $x^0$为初始点, $x^*$为算法终止时近似的最优解, $F(x^*)$为相应的最优值. P为所有迭代求解问题使用的约束个数, 数值结果见表 1.
数值实验表明, 算法B的迭代次数较算法C有明显减少, 近似最优值略有差距, 求解问题所使用的约束个数也有所减少.在精度要求不太高的情况下, 算法B采用的积极约束集识别技术和非单调技术可以降低求解计算量和迭代次数, 因而本文所提出的算法是有效的.