梯度遗传算法对一类非平滑优化问题的应用

王 科, 周中成, 颜文勇, 肖 翔

(1. 成都工业学院 信息与计算科学系, 四川 成都 611730; 2. 西南大学 数学与统计学院, 重庆 400715;3. 湖南师范大学 商学院, 湖南 长沙 410000)

脉冲系统可以应用于许多领域,如登陆月球的软着陆问题[1-2],具有脉冲控制影响的群体模型[3-6],生物系统中的农药的脉冲控制[7]和空气静力系统的脉冲响应[8]等.所以,研究脉冲系统的解有非常重要的实际意义.

已有的解决脉冲系统的论文,大多数都是通过时间缩放变换[9-12]来进行求解.文献[11]提出了一个等效优化问题,通过固定脉冲时间来进行脉冲系统控制.因为目标函数是关于状态变量的不连续函数,文献[3]引入了新的二进制变量进入目标函数以便将对象函数变换为平滑函数,并使用罚函数的方法来解决非线性平滑优化问题.

本文主要在文献[3,11]的基础上,进一步讨论了一类非平滑优化问题的解.通过引入约束的方法,将非平滑优化问题转化为标准的平滑优化问题,以便使用基于梯度的方法来进行求解.又因为基于梯度的算法易陷于局部最优解,所以使用GA提供的初始值来克服这个问题.在获得较好的初始解后,通过基于梯度的算法对解进行快速地优化，从而既充分利用GA提供的自由性,又利用了基于梯度法提供的高效性,简称为梯度遗传算法(简称GGA).

1 问题描述

首先考虑下面的脉冲系统：

i=1,2,…,p,

(1)

(2)

其中,x(t)∈Rn是状态变量,x0∈Rn是系统给定的初始状态;ζ∈Rr是需要优化的一组参数向量;fi:Rn×Rr→Rn(i=1,2,…,p),hi:Rn×Rr→Rn(i=1,2,…,p-1)是给定的函数;τ00,τpT和T> 0是给定的最终时间;τi(i=1,2,…,p-1)是转换时间并且满足约束

τi-τi - 1≥Δi,i=1,2,…,p,

(3)

其中Δi> 0是子系统的最小持续时间.显然

Δ1+Δ2+… +Δp≤T.

令T表示满足(3)式的全体元素

τ=[τ1,2,…,τp-1]∈Rp-1.

ζ为系统参数且满足约束：

aj≤ζj≤bj,j=1,2,…,r,

(4)

其中,aj和bj为给定实数,使得aj

假设1函数fi(i=1,2,…,p)和hi(i=1,2,…,p-1)是连续可微的.

假设2存在H1>0,使得

‖fi(u,w)‖≤H1(1+‖u‖),

(u,w)∈Rn×W,i=1,2,…,p,

其中‖·‖表示欧几里德范数.

由假设1和假设2知, 脉冲系统(1)～(2)对应于每一对(τ,ζ)∈T×W有唯一解x(·|τ,ζ).然后,定义一个成本函数J如下：

(5)

其中

Φi,l,k:R×Rr→R,i=1,2,…,p,

l=1,2,…,n,k=1,2,…,q

是单调非减函数,即

Φi,l,k(η,ζ)>Φi,l,k-1(η,ζ),

∀i,l,k,η,ζ,

(6)

其中Φi,l,0=0;αk(k=1,2,…,q)是给定的常数满足0=α0<α1<α2<…<αq,而且

将选择切换时间τ1,2,…,τp-1和系统参数ζ1,2,…,ζr以优化成本函数(5),同时满足约束(1)～(4).即得如下问题.

问题1存在(τ*,ζ*)∈T ×W,使得

2 等价问题

根据成本函数(5)式,由于问题1的不连续性,不能使用标准优化方法求解,以下将对问题1进行等价变换.

2.1问题陈述令(τ,ζ)∈T ×W和vilk∈R(i=1,2,…,p,l=1,2,…,n,k=1,2,…,q).考虑以下约束：

∀i,l,k,

(7)

∀i,l,k,

(8)

0≤vilk≤1,

(9)

以及

vilk(1-vilk)≤0.

(10)

注意到

其中x(·|τ,ζ)表示系统(1)～(2)的解.考虑如下优化问题.

问题2存在(τ*,ζ*,v*)∈F,使得

2.2问题1和问题2的等价对于每一对

(τ,ζ)∈T×W,

定义

(11)

而且

(12)

因此

J(τ,ζ).

(13)

∀i,l,k.

根据(6)式

从而

(14)

由(13)式得

(15)

生育期是评价品种能否作为直播再生稻安全种植的关键指标之一。从表1可以看出，10个品种中，两季全生育期比对照丰黄华占短3 d及以上的有金优38、两优33、天两优953、黄广油占共 4个品种，与对照相当（3 d以内）的有泰优2806、A优338、黄科香1号、黄科香2号共4个品种，比对照丰黄华占长3 d以上的只有甬优4949一个品种。甬优4949在11月5日才成熟，生育期偏长，作直播再生稻种植有风险。

另一方面,假设(τ*,ζ*,v*)是问题2的最优解.从而有

(16)

J(τ*,ζ*).

(17)

令(τ,ζ)是问题1的一个可行解.从而,根据(13)和(17)式

因为(τ,ζ)是任意选择的,上述不等式表明(τ*,ζ*)是问题1的最优解.因此问题1和问题2是等价的.

2.3时间缩放变换令

Θ={θ∈Rp:θi≥Δi,

i=1,2,…,p,θ1+…+θp=T}.

考虑下列脉冲系统

(18)

和

y(i)=y(i+)=

(19)

其中 (θ,ζ)∈Θ×W.令y(·|θ,ζ)表示系统(18)～(19)对于(θ,ζ)∈Θ×W的解.注意子系统的切换(18)～(19)式发生在固定的时间s=1,2,…,p-1.

系统参数vilk受以下约束:

(yl(i-|θ,ζ)-αk - 1)(1 -vilk)≤0,

∀i,l,k,

(20)

(yl(i-|θ,ζ)-αk - 1)vilk≥0,

∀i,l,k,

(21)

0≤vilk≤1 ,

(22)

vilk(1-vilk)≤0 .

(23)

Q(θ,ζ,v)=

根据文献[11],通过引入子系统的等价变换,将脉冲系统的切换时间变换到固定时间,即问题2和问题3是等价的,这里就不再证明了.

3 算法

在本节中,设计了一个用于计算问题3的梯度遗传算法(GGA).遗传算法(GA)是基于自然遗传学和生物进化的语言,可用于获得较好的初始值, 然后使用基于梯度的方法来改善结果.

3.1GGA的设计令

δu,i

和

对于每一个u=1,2,…,p,定义以下辅助系统：

δ(u,i)fi(y(s|θ,ζ),ζ),

s∈(i-1,i),i=1,2,…,p

(24)

和

ψu(i)=ψu(i+)=

(25)

其中(θ,ζ)∈Θ×W.令ψu(·|θ,ζ)表示系统(24)～(25)的解.

对于每一个j=1,2,…,r,定义另一个辅助系统如下：

s∈(i-1,i),i=1,2,…,p

(26)

和

φj(i)=φj(i+)=

(27)

其中(θ,ζ)∈Θ×W.令ψj(·|θ,ζ)表示系统(26)～(27)的解.

根据文献 [11],有以下重要的结果,它的证明与文献[11]非常相似的,这里就不再重复了.

定理1对于每一对(θ,ζ)∈Θ×W,

s∈[0,p],u= 1,2,…,p,

s∈[0,p],j= 1,2,…,r.

接下来,有Q的偏导数的公式,这对于GGA计算是非常重要的.由定理1和链式求导规则,

(28)

其中yy(·|θ,ζ),ψuψu(·|θ,ζ), 并且∂/∂x-分别表示第i次开关之前的状态的微分.

同理

(29)

从而,

[Φi,l,k(yl(i-),ζ)-Φi,l,k-1(yl(i-),ζ)] .

(30)

接下来,根据等式(28)～(30),计算Q及其偏导数.类似于文献[13-14],定义算法1和算法2如下所示.

算法1遗传算法，随机生成很多满足条件的个体 (θ,ζ):

(i) 通过(θ,ζ),求解脉冲系统(18)和(19),从而得到y(·|θ,ζ);

(ii) 用(θ,ζ)和y(·|θ,ζ)计算v,从而计算目标函数Q；

(iii) 选择一些好的个体遗传到下一代；

(iv) 随机交叉一些选定的个体；

(v) 随机变异一些个体；

(vi) 重复(iii)～(v)以m次,其中m是常数,通常m=200 .

由算法1,获得一个近似最优值(θ*,ζ*,v*).

算法2梯度算法:用遗传算法的得到的近似最优值(θ*,ζ*,v*)作为初值进行迭代.

(i) 求解脉冲系统(18)～(19)、(24)～(25)和(26)～(27),得到y(·|θ,ζ)、ψu(·|θ,ζ)和φj(·|θ,ζ)；

(ii) 用y(·|θ,ζ)和v计算目标函数Q；

3.2算例使用Yu等[15]提出的动态模型作为计算实例,该模型通过以下微分方程,描述了对虾的生物死亡率和生长过程：

(31)

(32)

其中,t是以周计的时间,x1(t)是时间t的剩余虾数,x2(t)是在时间t时以克为单位的单个虾的平均重量,m=-0.03 是对虾的自然死亡率的给定常数,a=3.5 是给定的常数,b=10- 5是给定的常数;x10=40 000和x20=1是给定的初始条件.令[0,T] 表示单个的生产周期(最终收获发生在时间T=13.2 ).

假设在生产周期中进行p次收获(p-1次中间收获和1次最终收获).令τi表示第i次收获的时间,其中τp是最后的收获时间.然后有以下约束：

0=τ0≤τ1≤τ2≤…≤τp=T.

令ζj表示在时间τi收获虾的比例数,很显然,

0≤ζj≤1 , i=1,2,…,p.

状态变量x1和x2是在t=τi收获后,约束于下面的脉冲条件

(33)

(34)

其中,对于一般函数h(t),采用符号

表示.

另外,还采用Blanchard等[3]提供的模型作为生产周期[0,T] 的总收入：

(35)

其中ξ[·] 表示相对于不同权重的价格函数,h=50,ci是常数且ci=0(i= 0,1,2,…,p),βk和αk(k=1,2,…,q)也是常数且在表1中给出.

表1 价格函数参数

通过本文第2节的等价变换,将问题转换为以下：

(36)

这里解决了这个问题在p=3和p=4的情况.当p=3时,最优解是3 683.1.当p=4时,最优解为3 648.7.与文献[3]中的结果相比,这些结果要好得多.

3.3结果虽然GA是一种随机算法,但得到一组近似最优解后,梯度算法很容易的将其收敛到最优解.如下系统有很多解决方案.当收获次数p=3时,得到如图1所示的结果.

图 1 当p=3时GA的收敛过程

图1中左图是遗传算法得到初值的过程,右图是梯度算法的收敛过程.

它表明目标函数通过GA收敛到3 672.7,最优解也稳定在约第50代.这意味着可以通过GA获得最大收入3 672.7.但GA只能获得近似的最优解,需要使用基于梯度的方法(算法2)进行改进.通过梯度算法,得到的最优解是3 683.1.决策变量和状态变量取值如下:

ζ1=0.193 3 ,ζ2= 0.433 9,ζ3=1;

τ1= 2.514 0,τ2= 6.058 3,τ3=13.2;

x2(τ2)=10.102 1,x2(τ3)=20.009 3.

同理,当收获次数p=4时,得到如图2所示的结果. 图2左图是遗传算法得到初值的过程,右图是梯度算法的收敛过程.它表明目标函数通过GA收敛到3 621.6,最优解也稳定在约第110代.通过梯度算法改进后,最优解是3 648.7.决策变量和状态变量取值如下:

ζ1=0.220 7,ζ2=0.337 4,

ζ3=0.157 4,ζ4=1;

τ1=3.324 4,τ2=6.550 8,

τ3=9.596 0,τ4=13.2;

x2(τ1)=6.293 3,x2(τ2)=10.247 2,

x2(τ3)=15.000 0,x2(τ4)=20.000 0.

图 2 当p=4时,GGA的收敛过程

4 结论

在本文中,引入了一种新的方法来解决一类非平滑动态优化问题,其中目标函数是一个不连续的函数.的新方法是通过引入不等式约束将非平滑最优问题转化为具有平滑目标函数的优化问题.特别的,梯度遗传算法对于解决一类非平滑优化问题是行之有效的算法.