基于SQP-IPM交替求解的CHP机组变负荷过程动态优化方法

基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法
技术领域
1.本发明涉及chp机组变负荷动态过程控制技术领域，具体涉及一种基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法。


背景技术：

2.改造后的热电联产(chp)机组动态调节特性复杂，呈现出多变量、多过程耦合、多参数叠加和强非线性等特征，常采用微分代数方程组(differential algebraic equations,daes)描述其在运行区间的动态特性，往往导致机组变负荷响应能力难以精确建模和快速仿真计算。含chp机组daes模型约束的动态过程控制优化问题，常采用联立法进行求解。该方法虽然在处理复杂约束问题时具有优势，但其同时将状态变量和控制变量进行离散化的做法，常使离散化后的非线性规划(nonlinear programming,nlp)命题规模成倍扩大，导致优化计算难度加大、求解时间增加。
3.现有nlp问题求解研究主要在以下几个方面:对含大量等式约束且自由度有限的大规模优化问题，可采用简约空间方法通过空间分解策略降低序列2次规划(sequential quadratic programming,sqp)子问题的维数，提高算法的求解效率。但空间分解过程会额外增加计算成本。对于不等式约束的问题常用内点法(interior point method,ipm)处理。ipm是将不等式约束转化为障碍项加入目标函数，原问题求解转换成一系列等式约束的障碍子问题求解，从而避免求解过程对积极集的确定。虽对初始点的依赖性不强，但在求解初始点远离最优点的大规模复杂过程优化问题时，不能兼顾收敛性和收敛速度。
4.同时，收敛准则对于计算时间具有很大影响，计算延迟将导致控制器性能下降。
5.优化算法常用传统的终止准则进行全部或部分的收敛控制，传统准则只能给出“收敛”或“不收敛”的刚性结论。在考虑变量的改变程度，目标函数的减小(或增加)程度等时，常使优化计算指标在达到某一收敛水平后，就进入了计算代价远远大于收敛效果改进的阶段。基于传统准则的非线性优化求解算法在求解动态优化问题时，往往过度迭代消耗大量计算时间，使得控制作用不能及时作用于被控系统。


技术实现要素：

6.为解决上述技术问题，本发明提供一种基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，能够有效解决chp机组变负荷动态过程最优控制问题，且连续控制过程能够快速收敛。
7.本发明采取的技术方案为：
8.基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，包括以下步骤：
9.步骤1：将发电厂热电联产chp机组变负荷优化控制问题表述成命题；
10.步骤2：求解步骤1中的动态优化命题，利用基于legendre多项式离散格式的配置法同时离散状态变量和控制变量，形成非线性规划nlp问题；
11.步骤3：设计基于sqp-ipm交替求解方法进行nlp问题求解，并在交替求解方法中引
入了收敛深度控制策略，以提高方法实时性。
12.所述步骤1中，动态优化命题的形式为：
[0013][0014]
式(1)中：min j(x(t0),t0,x(tf),tf)表示最小化动态优化命题关于状态变量的目标函数；t0、tf分别为机组优化控制的起始、终止时间点；t表示控制时域时间变量；
[0015]
x(t)、y(t)、u(t)分别表示t时刻chp机组微分状态变量、代数输出变量、控制变量；f(x(t),y(t),u(t),t)表示状态变量、输出变量和控制变量之间的关系函数模型；
[0016]
ge、g
ie
为等式与不等式路径约束；he、h
ie
为在时刻tf的终值约束；
[0017]
x(t0)、y(t0)、x(tf)、y(tf)分别表示状态变量初值、输出变量初值、状态变量终值、输出变量终值；
[0018]
最优控制性能目标函数描述为：
[0019][0020]
式(2)中：t0和tf分别为机组优化控制的起始和终止时间点；u为chp机组控制变量组成的向量；y为chp机组代数输出变量组成的向量，y
*
为稳态优化得到的输出变量的设定值；r和s是正定权矩阵，能够通过对变量做归一化来确定初始值；j为目标函数，u(t)为t时刻控制变量值、u(t-1)为t-1时刻控制变量值；y(t)为t时刻输出变量值。
[0021]
热电联产机组变负荷优化控制问题中，u为控制变量，y为输出变量，x为状态变量，其中：
[0022][0023]
式(3)中：x为状态变量，对应ψd为汽包压力，ph为机组发电功率，ψz为汽轮机供热抽汽压力，vf为锅炉燃烧率；
[0024]
u为控制变量，对应v
t
为汽轮机高压缸进汽调节阀开度，vb为机组燃料量，vh为供热抽汽调节蝶阀开度；
[0025]
y为输出变量，对应ψ
t
为汽轮机前压力，ph为机组发电功率，mh为供热抽汽流量。
[0026]
所述步骤1中，
[0027]
①
、chp机组机理模型如下：
[0028]vm
(t)＝vb(t-tb)
ꢀꢀꢀ
(4)；
[0029][0030]
式(4)～式(5)表征燃煤制粉系统动态关系，vb(t)为机组给煤质量流量，vm(t)为制粉系统中实际进入磨煤机的给煤质量流量，tb为制粉过程迟延时间常数，vf(t)为锅炉燃烧率，tf为制粉惯性时间常数，t表示控制过程时域值，vb(t-tb)表示t-tb时刻的机组给煤质量流量；
[0031][0032]
ψd(t)-ψ
t
(t)＝k2(k
1vf
(t))
1.5
ꢀꢀꢀ
(7)；
[0033]
式(6))～式(7)表征汽包锅炉系统动态关系，v
t
(t)为汽轮机调门开度，ψd(t)为汽包压力，ψ
t
(t)为汽轮机前压力，cd为锅炉蓄热系数，δ为常系数，k1为燃料增益静态参数，k2为过热器阻力系数，k3为汽轮机增益静态参数；
[0034][0035]
ψ1(t)＝0.01ψ
t
(t)v
t
(t)
ꢀꢀꢀ
(9)；
[0036][0037]
式(8)～式(10)表征汽轮机系统动态关系，vh(t)为抽汽调节蝶阀开度，ph(t)为机组发电功率，ψz(t)为中压缸排汽压力，ψ1(t)为汽轮机一级压力，mr(t)为循环水质量流量，εr(t)为循环水回水温度，t
t
为汽轮机惯性时间常数，cz为热网加热器蓄热系数，k4为高中压缸占汽轮机做功比例，k5为低压缸增益静态参数，k6为热网循环水的有效比热容；
[0038]
mh(t)＝k7k6mr(t)(96ψz(t)-εr(t)+103)
ꢀꢀꢀ
(11)；
[0039]
式(11)表征供热系统动态关系，mh(t)为供热抽汽流量；k7为供热抽汽有效热量折合蒸汽流量系数。
[0040]
②
、chp机组控制模型：
[0041]
设计一种比例-积分控制算法协同调节vb(t)、v
t
(t)、vh(t)的方向与大小，使输出变量ψ
t
(t)、mh(t)、ph(t)变化，其数学模型如下：
[0042][0043]
其中，k
pt
、k
it
分别为控制v
t
(t)变量的比例系数、积分系数；k
pb
、k
ib
、分别为控制vb(t)变量的比例系数、积分系数；k
ph
、k
ih
、分别为控制vh(t)变量的比例系数、积分系数；sp
ψ
、sp
p
、spm分别为主蒸汽、机组发电功率、供热抽汽流量的设定值；ψ
t
(t)为主蒸汽压力控制过程当前值，ph(t)为机组发电功率控制过程当前值，mh(t)为供热抽汽流量控制过程当前值。
[0044]
考虑到chp机组运行的安全稳定，在控制过程中输出变量过程值的波动范围、稳态值的误差范围不能超出允许的限值，应满足如下约束：
[0045][0046][0047]
其中，
△
ψ
、
△
p
、
△m分别为主蒸汽压力、电功率、抽汽蒸汽流量的波动范围；δ
ψ
、δ
p
、δm分别为主蒸汽压力、电功率、抽汽蒸汽流量的误差范围；优化控制周期内，ψ
t
(tf)、ph(tf)、mh(tf)分别为终点时刻tf主蒸汽压力、电功率、抽汽蒸汽流量的终值。
[0048]
③
、边界约束：
[0049]
chp机组在[t0,tf]时段内完成变负荷调节，各控制变量、输出变量的初始状态和终止状态如下：
[0050][0051][0052]
式中，t0和tf分别为机组优化控制的起始和终止时间点；
[0053]
分别为起始时刻t0输出量汽轮机前压力、机组发电功率、供热抽汽流量的取值；
[0054]
分别为终止时刻tf输出量汽轮机前压力、机组发电功率、供热抽汽流量的取值；
[0055]
分别为起始时刻t0控制量高压缸进汽调节阀开度、机组燃料量、抽汽调节蝶阀开度的取值；
[0056]
分别为终止时刻tf控制量高压缸进汽调节阀开度、机组燃料量、抽汽调节蝶阀开度的取值。
[0057]
④
、路径约束：
[0058]
在变负荷调节过程中，控制变量v
t
(t)、vb(t)、vh(t)需要满足以下约束：
[0059][0060]
式中，和分别为汽轮机调门开度变化的最小值和最大值；和分别为机组燃料量变化的最小值和最大值；和分别为抽汽调节蝶阀开度变化的最小值和最大值。所述步骤2中，包括如下步骤：
[0061]
1)时域变换：
[0062]
将最优控制时域t∈[t0,tf]划分为p段有限元，即0＝t0《t1《
…
t
q-1
《tq《
…
t
p
＝tf，由此可知，每段时间域的步长为：
[0063]
△
t＝t
q-t
q-1
,q∈[1,p]
ꢀꢀꢀ
(18)；
[0064]
其中：
△
t为每段时间域的步长，q为有限元个数，p为有限元总段数，tq为第q个有限元的时间变量，t
q-1
为第q-1个有限元的时间变量；
[0065]
在高斯伪谱法中配置点区间被选取为τ∈[-1,1]，因此通过时域变换将第q段上时间变量转换为变量τ，其变换如下：
[0066][0067]
其中：τ为第q段上时间变量通过时域变换转换的变量；
[0068]
则可得每段时间区域内的最优控制问题的描述如下：
[0069][0070]
式(20)中：j
(q)
为第q段时间区域上的目标函数，τ0为第q段上时间变量时域变换的初始配置点，τf第q段上时间变量时域变换的终端配置点；
[0071]
x
(q)
(τ)为第q段时间上第τ个配置点的微分状态变量，x
(q)
(τ0)为第q段时间上第τ0个配置点的微分状态变量，x
(q)
(τf)为第q段时间上第τf个配置点的微分状态变量；
[0072]y(q)
(τ)为第q段时间上第τ个配置点的代数输出变量，y
(q)
(τ0)为第q段时间上第τ0个配置点的代数输出变量，y
(q)
(τf)为第q段时间上第τf个配置点的代数输出变量；
[0073]u(q)
(τ)为第q段时间上第τ个配置点的控制变量，f为daes动态模型，
△
t为每段时间域的步长。
[0074]
2)变量离散化：
[0075]
高斯伪谱法是同时将控制变量和状态输出变量在legendre-gauss(lg)配置点上进行离散化，再通过构建拉格朗日插值函数逼近控制变量和状态变量，在区间[τ0,τf]上，n阶的legendre多项式为：
[0076][0077]
其中：pn(τ)为n阶legendre多项式，n为阶次，为函数对τ求n阶导数；
[0078]
以式(21)的零点和τ0组成n+1个配置点，在每个t
(q)
时间域中，微分状态变量、控制变量由多项式近似表示为：
[0079][0080]
其中：x
(q)
(τ)为第q段时间上第τ个配置点的微分状态变量，为状态变量的插值函数，x
(q)
(τi)为第q个有限元第i个配置点上状态变量的值；
[0081]u(q)
(τ)为第q段时间上第τ个配置点的控制变量，为控制变量的插值函数，u
(q)
(τi)为第q个有限元第i个配置点上控制变量的值。
[0082]
式(22)中，lagrange插值基函数为：
[0083][0084]
其中：n为插值阶次；为状态变量的插值函数，为控制变量的插值函数，第j个配置点上插值函数自变量，τi为第i个配置点上插值函数自变量；τ为第q段上时间变量通过时域变换转换的插值函数自变量；
[0085]
由式(22)状态变量求导，可将插值函数微分逼近状态变量的微分形式，得到在τk点的状态变量导数，可近似表示为：
[0086][0087]
其中：为第p段有限元第k个配置点上微分状态变量的导数，为第p段有限元第k个配置点上状态变量的插值函数，x
(q)
(τi)为第q个有限元第i个配置点上状态变量；
[0088]
结合式(20)和式(24)可得微分代数等式约束：
[0089][0090]
其中：
△
t为每段时间域的步长；f为daes动态模型，x(τk)为第τk个配置点上状态变量，m(τk)第τk个配置点上其他变量，τk为第q段上时间变量时域变换的第k个配置点；
[0091]
在机组变负荷控制过程中要保证各个状态输出变化轨迹连续性，因此状态变量也必须为连续变化，即每段终端时刻状态变量作为下一段初始时刻状态变量，其连接处约束等式为
[0092]
x
(q-1)
(τf)＝x
(q)
(τ0)
ꢀꢀꢀ
(26)；
[0093]
其中：x
(q-1)
(τf)为为第q-1个有限元第τf个配置点上状态变量，x
(q)
(τ0)为第q个有限元第τ0个配置点上状态变量；
[0094]
由于lg配置点中不包含终端时刻状态值，由式(26)，根据legendre-gauss型求积公式计算终端时刻状态变量：
[0095][0096]
式中，ak和f(τk)分别为高斯积分权重系数和求积函数。
[0097]
经过上述离散过程，能够将chp机组变负荷动态过程的最优控制问题转化为非线性规划问题进行求解。
[0098]
所述步骤3中，收敛深度控制策略具体如下：
[0099]
原命题式(20)所示的非线性规划问题可简写为：
[0100][0101]
其中，为标量目标函数，x∈rn为n维实数向量，c∈rm为m维等式方程。
[0102]
问题(28)的解必须满足如下一阶最优性条件:
[0103][0104]
其中，为目标函数的梯度向量，为约束方程的雅可比矩阵。λ为对应于约束方程的拉格朗日乘子，ν为对应于x的对偶变量，x＝diag(x)。
[0105]
假定由优化算法产生的迭代序列为{xk}，定义当前迭代点xk可行性误差目标函数预测改进量
[0106][0107][0108]
其中，为状态向量的第i个分量。c(xk)为系统的约束条件，xk为当前迭代点，x
k-1
为前一个迭代点，为当前迭代点梯度向量，dk为当前迭代点搜索方向。
[0109]
引入s为变形的sigmoid函数柔化传统的刚性收敛准则，其定义为：
[0110]
s(δk,ε0)＝[tanh(ξ
*
(logδk/logε0))]/[tanh(ξ)]
ꢀꢀꢀ
(32)
[0111]
其中，函数s(δk,ε0)将区间[ε0,1/ε0]光滑地连接，ξ用于衡量函数s在该区间内的平缓程度，ξ取1.5。
[0112]
则收敛深度
[0113][0114]
其中，ε0为给定的容差，s为变形的sigmoid函数，为可行性误差，为目标
函数预测改进量。
[0115]
所述步骤3中，sqp-ipm交替规划算法：
[0116]
基于cdc的sqp-ipm交替规划算法是采用全局收敛性强的ipm算法将初始点向最优点推进，利用收敛深度来衡量当前迭代点与最优点的距离，控制推进程度，进而切换到局部收敛性较好的sqp算法使迭代点连续控制过程快速收敛到最优点。本发明引入进展程度指标用于sqp-ipm交替规划算法之间的切换；
[0117]
可行性变化量以及目标函数变化量为：
[0118][0119][0120]
其中，分别为当前迭代点、前一个迭代点的可行性误差，f(xk)、f(x
k-1
)分别为xk、x
k-1
迭代点目标函数值。
[0121]
由上述公式可将和进展程度定义为：
[0122][0123]
其中，为可行性变化量，为目标函数变化量，ε0为给定的容差。
[0124]
当进展程度达到1.0时，则算法在可行性变化量和目标函数改进量上已经达到ε0，也就是说算法在目标函数和可行性改进上几乎已经没有进展了。
[0125]
基于上述引入的收敛深度和进展程度对基于cdc的sqp-ipm交替规划算法框架进行设计，具体步骤如下：
[0126]
步骤(1)：令k＝0，给定初始点x0，指定收敛深度阈值θc，进展程度阈值θ
p
，收敛深度推进步长d
θ
，最大收敛深度值θ
max
，并将ipm的标志位设为0；
[0127]
步骤(2)：选取ipm，将其标志位设为1，用于求解优化问题；
[0128]
步骤(3)：迭代一步，k＝k+1，若ipm的终止准则满足，则求解成功，整个算法终止，否则转步骤(4)；
[0129]
步骤(4)：计算收敛深度，检查是否满足，若满足，则将当前迭代点xk赋给x0，并将sqp的标准位都清零，转步骤(6)，否则转步骤(5)；
[0130]
步骤(5)：若ipm的终止准则求解失败未被满足，则转步骤(3)，否则求解失败整个算法终止；
[0131]
步骤(6)选取sqp并将其标志位设为1，以x0为初始点求解优化问题；
[0132]
步骤(7)迭代一步，k＝k+1，若sqp的终止准则满足，则求解成功，整个算法终止，否则转步骤(8)；
[0133]
步骤(8)计算进展程度，若和sqp的终止准则都不满足，则转步骤(7)，否则转步骤(9)；
[0134]
步骤(9)将当前迭代点xk赋给x0，转步骤(10)；
[0135]
步骤(10)令θc＝θc+dθ，若θc≥θ
max
则求解失败整个算法终止，否则将ipm的标志位清零，转步骤(2)。
[0136]
通过上述步骤，实现chp机组变负荷过程动态优化。
[0137]
步骤4：以某300mw抽汽式chp机组变负荷控制需求为算例进行仿真，验证了算法的可行性和连续控制过程快速收敛性。
[0138]
采用matlab整理模型输出数据，绘制chp机组连续变负荷控制变量变化轨迹、输出变量曲线，包括图给煤质量流量变化轨迹、汽轮机高调门开度变化轨迹、抽汽蝶阀开度变化轨迹、电出力输出曲线、主蒸汽压力输出曲线、供热蒸汽流量输出曲线。
[0139]
本发明一种基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，技术效果如下：
[0140]
1)本发明方法提出的收敛深度控制准则，可对chp机组变负荷迭代求解过程的动态变化进行衡量。机组优化求解中，对于过度收敛的求解过程可快速终值迭代，得到较好的解；对于求解失败的情况也可提前结束求解进程，减少计算代价的消耗。该方法可以实现显著减少chp机组变负荷过程优化计算时间。
[0141]
2)本发明方法按照全局和局部收敛性选择ipm和sqp，来求解chp机组变负荷动态优化问题，利用收敛深度来控制算法的相互协作、交替求解。充分发挥不同算法的优势，方法在处理chp机组变负荷过程优化的大规模复杂优化问题时，在计算效率和收敛性方面有明显提升。
附图说明
[0142]
图1为刚性收敛准则曲线柔性收敛准则曲线图。
[0143]
图2为连续控制过程动态优化问题的快速收敛sqp-ipm交替求解流程图。
[0144]
图3(1)为chp机组连续变负荷过程给煤质量流量变化轨迹图；
[0145]
图3(2)为chp机组连续变负荷过程汽轮机高调门开度变化轨迹图；
[0146]
图3(3)为chp机组连续变负荷过程抽汽蝶阀开度变化轨迹图；
[0147]
图3(4)为chp机组连续变负荷过程电出力输出曲线图；
[0148]
图3(5)为chp机组连续变负荷过程主蒸汽压力输出曲线图；
[0149]
图3(6)为chp机组连续变负荷过程供热蒸汽流量输出曲线图。
具体实施方式
[0150]
基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，包括以下步骤：
[0151]
步骤一、给出含微分代数方程的等式与不等式约束动态优化问题具体描述，确定连续控制过程动态优化命题：
[0152]
将发电厂热电联产机组变负荷优化控制问题表述成命题，动态优化命题的一般形式为：
[0153][0154]
式(1)中：f为chp机组daes动态模型，t0和tf分别为机组优化控制的起始和终止时间点；ge、g
ie
为等式与不等式路径约束；he、h
ie
为在时刻tf的终值约束；x(t),y(t),u(t)分别表示t时刻chp机组微分状态变量、代数输出变量、控制变量；
[0155]
1)优化命题的目标函数描述：
[0156]
为保证chp机组控制过程准确的、快速平稳的实现变负荷调节，选用控制变量与输出变量结合的数学形式作为最优控制性能目标函数，其描述为
[0157][0158]
式(2)中：t0和tf分别为机组优化控制的起始和终止时间点；u为chp机组控制变量组成的向量，y为chp机组代数输出变量组成的向量，y
*
为稳态优化得到的输出变量的设定值，r和s是正定权矩阵，能够通过对变量做归一化来确定初始值；j为目标函数，u(t)为t时刻控制变量值、u(t-1)为t-1时刻控制变量值、y(t)为t时刻输出变量值。
[0159]
热电联产机组变负荷优化控制问题中u为控制变量，y为输出变量，x为状态变量，其中：
[0160][0161]
式(3)中：x为状态变量，对应ψd为汽包压力，ph为机组发电功率，ψz为汽轮机供热抽汽压力，vf为锅炉燃烧率；u为控制变量，对应v
t
为汽轮机高压缸进汽调节阀开度，vb为机组燃料量，vh为供热抽汽调节蝶阀开度；y为输出变量，对应ψ
t
为汽轮机前压力，ph为机组发电功率，mh为供热抽汽流量。
[0162]
2)优化命题的等式与不等式约束描述：
[0163]
①
chp机组机理模型：
[0164]
通过对供热机组热力组件特性及蒸汽在汽轮机内做功过程分析，可得到简化的chp机组非线性机理模型描述：
[0165]vm
(t)＝vb(t-tb)
ꢀꢀꢀ
(4)
[0166][0167]
[0168]
ψd(t)-ψ
t
(t)＝k2(k
1vf
(t))
1.5
ꢀꢀꢀ
(7)
[0169][0170]
ψ1(t)＝0.01ψ
t
(t)v
t
(t)
ꢀꢀꢀ
(9)
[0171][0172]
mh(t)＝k7k6mr(t)(96ψz(t)-εr(t)+103)
ꢀꢀꢀ
(11)
[0173]
式(4)-(5)表征燃煤制粉系统动态关系，vb(t)为机组给煤质量流量，vm(t)为制粉系统中实际进入磨煤机的给煤质量流量，tb为制粉过程迟延时间常数，vf(t)为锅炉燃烧率，tf为制粉惯性时间常数，t表示控制过程时域值，vb(t-tb)表示t-tb时刻的机组给煤质量流量；
[0174]
式(6)-(7)表征汽包锅炉系统动态关系，v
t
(t)为汽轮机调门开度，ψd(t)为汽包压力，ψ
t
(t)为汽轮机前压力，cd为锅炉蓄热系数，δ为常系数，k1为燃料增益静态参数，k2为过热器阻力系数，k3为汽轮机增益静态参数；
[0175]
式(8)-(10)表征汽轮机系统动态关系，vh(t)为抽汽调节蝶阀开度，ph(t)为机组发电功率，ψz(t)为中压缸排汽压力，ψ1(t)为汽轮机一级压力，mr(t)为循环水质量流量，εr(t)为循环水回水温度，t
t
为汽轮机惯性时间常数，cz为热网加热器蓄热系数，k4为高中压缸占汽轮机做功比例，k5为低压缸增益静态参数，k6为热网循环水的有效比热容；
[0176]
式(11)表征供热系统动态关系，mh(t)为供热抽汽流量；k7为供热抽汽有效热量折合蒸汽流量系数。
[0177]
②
chp机组控制模型：
[0178]
针对chp机组机理模型，控制算法要实现电功率准确跟踪、热功率快速恢复和系统压力安全三个关键功能，设计一种比例-积分控制算法协同调节vb(t)、v
t
(t)、vh(t)的方向与大小，使输出变量ψ
t
(t)、mh(t)、ph(t)变化。其数学模型：
[0179][0180]
其中，k
pt
、k
it
分别为控制v
t
(t)变量的比例系数、积分系数；k
pb
、k
ib
、分别为控制vb(t)变量的比例系数、积分系数；k
ph
、k
ih
、分别为控制vh(t)变量的比例系数、积分系数；sp
ψ
、sp
p
、spm分别为主蒸汽、机组发电功率、供热抽汽流量的设定值；ψ
t
(t)为主蒸汽压力控制过程当前值，ph(t)为机组发电功率控制过程当前值，mh(t)为供热抽汽流量控制过程当前值。
[0181]
考虑到chp机组运行的安全稳定，在控制过程中输出变量过程值的波动范围、稳态值的误差范围不能超出允许的限值，应满足如下约束：
[0182]
[0183][0184]
其中，
△
ψ
、
△
p
、
△m分别为主蒸汽压力、电功率、抽汽蒸汽流量的波动范围；δ
ψ
、δ
p
、δm分别为主蒸汽压力、电功率、抽汽蒸汽流量的误差范围；优化控制周期内，ψ
t
(tf)、ph(tf)、mh(tf)分别为终点时刻tf主蒸汽压力、电功率、抽汽蒸汽流量的终值。
[0185]
③
边界约束：
[0186]
chp机组在[t0,tf]时段内完成变负荷调节，各控制变量、输出变量的初始状态和终止状态如下：
[0187][0188][0189]
式中，t0和tf分别为机组优化控制的起始和终止时间点；
[0190]
分别为起始时刻t0输出量汽轮机前压力、机组发电功率、供热抽汽流量的取值；
[0191]
分别为终止时刻tf输出量汽轮机前压力、机组发电功率、供热抽汽流量的取值；
[0192]
分别为起始时刻t0控制量高压缸进汽调节阀开度、机组燃料量、抽汽调节蝶阀开度的取值；
[0193]
分别为终止时刻tf控制量高压缸进汽调节阀开度、机组燃料量、抽汽调节蝶阀开度的取值。
[0194]
④
路径约束：
[0195]
在变负荷调节过程中，控制变量v
t
(t)、vb(t)、vh(t)需要满足以下约束：
[0196][0197]
式中，和分别为汽轮机调门开度变化的最小值和最大值，和分别
为机组燃料量变化的最小值和最大值，和分别为抽汽调节蝶阀开度变化的最小值和最大值。
[0198]
步骤二、求解步骤一中的动态优化命题，利用基于legendre多项式离散格式的配置法，同时离散状态变量和控制变量，形成的大规模nlp问题：
[0199]
1)时域变换：
[0200]
将最优控制时域t∈[t0,tf]划分为p段有限元，即0＝t0《t1《
…
t
q-1
《tq《
…
t
p
＝tf，由此可知，每段时间域的步长为
[0201]
△
t＝t
q-t
q-1
,q∈[1,p]
ꢀꢀꢀ
(18)
[0202]
其中：
△
t为每段时间域的步长，q为有限元个数，p为有限元总段数，tq为第q个有限元的时间变量，t
q-1
为第q-1个有限元的时间变量；
[0203]
在高斯伪谱法中配置点区间被选取为τ∈[-1,1]，因此通过时域变换将第q段上时间变量转换为变量τ，其变换如下：
[0204][0205]
其中：τ为第q段上时间变量通过时域变换转换的变量；
[0206]
则可得每段时间区域内的最优控制问题的描述如下：
[0207][0208]
式(20)中：j
(q)
为第q段时间区域上的目标函数，τ0为第q段上时间变量时域变换的初始配置点，τf第q段上时间变量时域变换的终端配置点；
[0209]
x
(q)
(τ)为第q段时间上第τ个配置点的微分状态变量，x
(q)
(τ0)为第q段时间上第τ0个配置点的微分状态变量，x
(q)
(τf)为第q段时间上第τf个配置点的微分状态变量；
[0210]y(q)
(τ)为第q段时间上第τ个配置点的代数输出变量，y
(q)
(τ0)为第q段时间上第τ0个配置点的代数输出变量，y
(q)
(τf)为第q段时间上第τf个配置点的代数输出变量；
[0211]u(q)
(τ)为第q段时间上第τ个配置点的控制变量，f为daes动态模型，
△
t为每段时间域的步长；
[0212]
2)变量离散化：
[0213]
高斯伪谱法是同时将控制变量和状态输出变量在legendre-gauss(lg)配置点上进行离散化，再通过构建拉格朗日插值函数逼近控制变量和状态变量，在区间[τ0,τf]上，n阶的legendre多项式为
[0214]
[0215]
其中：pn(τ)为n阶legendre多项式，n为阶次，为函数对τ求n阶导数；
[0216]
以式(21)的零点和τ0组成n+1个配置点，在每个t
(q)
时间域中，微分状态变量、控制变量由多项式近似表示为
[0217][0218]
其中：x
(q)
(τ)为第q段时间上第τ个配置点的微分状态变量，为状态变量的插值函数，x
(q)
(τi)为第q个有限元第i个配置点上状态变量的值，u
(q)
(τ)为第q段时间上第τ个配置点的控制变量，为控制变量的插值函数，u
(q)
(τi)为第q个有限元第i个配置点上控制变量的值。
[0219]
式(22)中，lagrange插值基函数为
[0220][0221]
其中：n为插值阶次；为状态变量的插值函数，为控制变量的插值函数，第j个配置点上插值函数自变量，τi为第i个配置点上插值函数自变量；τ为第q段上时间变量通过时域变换转换的插值函数自变量；
[0222]
由式(22)状态变量求导，可将插值函数微分逼近状态变量的微分形式，得到在τk点的状态变量导数，可近似表示为
[0223][0224]
其中：为第p段有限元第k个配置点上微分状态变量的导数，为第p段有限元第k个配置点上状态变量的插值函数，x
(q)
(τi)为第q个有限元第i个配置点上状态变量；
[0225]
结合式(20)和式(24)可得微分代数等式约束：
[0226][0227]
其中：
△
t
为每段时间域的步长；f为daes动态模型，x(τk)为第τk个配置点上状态变量，m(τk)第τk个配置点上其他变量，τk为第q段上时间变量时域变换的第k个配置点；
[0228]
在机组变负荷控制过程中要保证各个状态输出变化轨迹连续性，因此状态变量也必须为连续变化，即每段终端时刻状态变量作为下一段初始时刻状态变量，其连接处约束等式为
[0229]
x
(q-1)
(τf)＝x
(q)
(τ0)
ꢀꢀꢀ
(26)
[0230]
其中：x
(q-1)
(τf)为为第q-1个有限元第τf个配置点上状态变量，x
(q)
(τ0)为第q个有
限元第τ0个配置点上状态变量；
[0231]
由于lg配置点中不包含终端时刻状态值，由式(26)，根据legendre-gauss型求积公式计算终端时刻状态变量：
[0232][0233]
式中，ak和f(τk)分别为高斯积分权重系数和求积函数。
[0234]
经过上述离散过程，可将chp机组变负荷动态过程的最优控制问题转化为非线性规划问题(nlp)进行求解。
[0235]
步骤三、设计基于sqp-ipm交替求解方法进行问题求解，并在交替求解方法中引入了收敛深度控制策略提高方法实时性；
[0236]
1)收敛深度控制策略：
[0237]
原命题式(20)所示的非线性规划问题可简写为
[0238][0239]
其中，为标量目标函数，x∈rn为n维实数向量，c∈rm为m维等式方程。
[0240]
问题(28)的解必须满足如下一阶最优性条件:
[0241][0242]
其中，为目标函数的梯度向量，为约束方程的雅可比矩阵。λ为对应于约束方程的拉格朗日乘子，ν为对应于x的对偶变量，x＝diag(x)。
[0243]
假定由优化算法产生的迭代序列为{xk}，定义当前迭代点xk可行性误差目标函数预测改进量
[0244][0245][0246]
其中，为状态向量的第i个分量。c(xk)为系统的约束条件，xk为当前迭代点，x
k-1
为前一个迭代点，为当前迭代点梯度向量，dk为当前迭代点搜索方向。
[0247]
引入s为变形的sigmoid函数柔化传统的刚性收敛准则，其定义为：
[0248]
s(δk,ε0)＝[tanh(ξ
*
(logδk/logε0))]/[tanh(ξ)]
ꢀꢀꢀ
(32)
[0249]
其中，函数s(δk,ε0)将区间[ε0,1/ε0]光滑地连接，ξ衡量函数s在该区间内的平缓程度，ξ取1.5。
[0250]
采用sigmoid函数对传统刚性收敛准则的柔化，给定不同容差ε0，如图1所示，曲线1为1
×
10-10
，曲线2为1
×
10-5
，曲线3为1
×
10-3
，曲线4为刚性准则。收敛深度描述为：
[0251][0252]
其中，ε0为给定的容差，s为变形的sigmoid函数，为可行性误差，为目标
函数预测改进量。
[0253]
收敛深度相比于一阶最优性条件更加直观和统一地表征了优化算法产生的每一个迭代点的质量，包括可行性和最优性，并能定量地描述其与最优点之间的距离程度。当收敛深度达到1.0时，可行性误差和目标函数预测改进量也就达到了ε0，此时迭代点已经是最优点的一个很好的近似值。
[0254]
2)sqp-ipm交替规划算法：
[0255]
基于cdc的sqp-ipm交替规划算法是采用全局收敛性强的ipm算法将初始点向最优点推进，利用收敛深度来衡量当前迭代点与最优点的距离，控制推进程度，进而切换到局部收敛性较好的sqp算法使迭代点连续控制过程快速收敛到最优点。本发明引入进展程度指标用于sqp-ipm交替规划算法之间的切换。
[0256]
可行性变化量以及目标函数变化量为：
[0257][0258][0259]
其中，分别为当前迭代点、前一个迭代点的可行性误差，f(xk)、f(x
k-1
)分别为xk、x
k-1
迭代点目标函数值。
[0260]
由上述公式可将和进展程度定义为：
[0261][0262]
其中，为可行性变化量，为目标函数变化量，ε0为给定的容差。
[0263]
当进展程度达到1.0时，则算法在可行性变化量和目标函数改进量上已经达到ε0，也就是说算法在目标函数和可行性改进上几乎已经没有进展了。
[0264]
基于上述引入的收敛深度和进展程度对基于cdc的sqp-ipm交替规划算法框架进行设计。首先选取ipm用于求解nlp问题，在迭代过程中采用式(33)计算每一步的收敛深度，当收敛深度满足指定的阈值θc时，保存迭代点的信息，终止当前算法并切换到sqp，使用上一次保存的迭代点作为初始点继续迭代。由于不同的算法，其具有的局部收敛半径也不一样，故本发明先给定阈值θc一个初值若sqp不成功，则从ipm所得到的迭代点未进入这些算法的局部收敛半径，将收敛深度阈值增加d
θ
∈(0,1)，回到ipm进行迭代。交替规划算法的具体步骤如下：
[0265]
步骤(1)：令k＝0，给定初始点x0，指定收敛深度阈值θc，进展程度阈值θ
p
，收敛深度推进步长d
θ
，最大收敛深度值θ
max
，并将ipm的标志位设为0；
[0266]
步骤(2)：选取ipm，将其标志位设为1，用于求解优化问题；
[0267]
步骤(3)：迭代一步，k＝k+1，若ipm的终止准则满足，则求解成功，整个算法终止，否则转步骤(4)；
[0268]
步骤(4)：计算收敛深度，检查是否满足，若满足，则将当前迭代点xk赋给x0，并将sqp的标准位都清零，转步骤(6)，否则转步骤(5)；
[0269]
步骤(5)：若ipm的终止准则(求解失败)未被满足，则转步骤(3)，否则求解失败整个算法终止；
[0270]
步骤(6)：选取sqp并将其标志位设为1，以x0为初始点求解优化问题；
[0271]
步骤(7)：迭代一步，k＝k+1，若sqp的终止准则满足，则求解成功，整个算法终止，否则转步骤(8)；
[0272]
步骤(8)：计算进展程度，若和sqp的终止准则都不满足，则转步骤(7)，否则转步骤(9)；
[0273]
步骤(9)：将当前迭代点xk赋给x0，转步骤(10)；
[0274]
步骤(10)：令θc＝θc+dθ，若θc≥θ
max
则求解失败整个算法终止，否则将ipm的标志位清零，转步骤(2)。
[0275]
本发明基于cdc的sqp-ipm交替求解动态优化问题的具体实现框图，如图2所示。
[0276]
步骤四、以某300mw抽汽式chp机组变负荷控制需求为算例进行仿真，验证了算法的可行性和连续控制过程快速收敛性。采用matlab整理模型输出数据，绘制chp机组连续变负荷控制变量变化轨迹、输出变量曲线，包括图给煤质量流量变化轨迹、汽轮机高调门开度变化轨迹、抽汽蝶阀开度变化轨迹、电出力输出曲线、主蒸汽压力输出曲线、供热蒸汽流量输出曲线。
[0277]
实施例：
[0278]
为验证本发明提出方法的有效性，在matlab平台搭建300mw抽汽式chp机组模型进行仿真，机组相关参数如表1、表2、表3所示，表4为应用场景的连续变负荷需求算例数据。
[0279]
表1 chp机组工作点参数
[0280][0281]
表2 chp机组工况数据
[0282][0283][0284]
表3 chp机组控制系统参数
[0285]vt
控制参数vb控制参数vh控制参数k
pt
＝-0.9k
pb
＝0.1k
ph
＝0.01kit
＝-1k
ib
＝0.01k
ih
＝10
[0286]
表4 供热工况下chp机组连续变负荷需求
[0287]
时段初始值终端值300s-600s(250,16.67,390)(220,16.67,380)600s-900s(220,16.67,380)(200,16.67,390)900s-1200s(200,16.67,390)(230,16.67,395)1200s-1500s(230,16.67,395)(210,16.67,385)1500s-1800s(210,16.67,385)(240,16.67,395)1800s-2100s(240,16.67,395)(250,16.67,385)2100s-2400s(250,16.67,385)(220,16.67,380)2400s-2700s(220,16.67,380)(200,16.67,390)2700s-3000s(200,16.67,390)(230,16.67,395)3000s-3300s(230,16.67,395)(210,16.67,385)3300s-3600s(210,16.67,385)(240,16.67,395)3600s-3900s(240,16.67,395)(260,16.67,385)
[0288]
1.仿真计算性能分析：
[0289]
采用基于cdc的ipm-sqp交替求解方法、基于cdc的ipm以及单独的ipm，分别对chp机组连续变负荷过程最优控制问题进行计算，其性能指标如表5所示。由表5可知，在不同工况与变负荷需求下，3种方法的最优控制终端值均可达到设定值，且可行性误差也符合实际要求，各时段均可收敛成功，但离散模型的迭代次数和求解耗时有较大差异。
[0290]
表5 不同算法计算性能比较
[0291]
性能指标cdc-ipm-sqpcdc-ipmipm求解耗时10.17s19.92s38.68s收敛性各时段均成功各时段均成功各时段均成功终端值各时段均达到终值各时段均达到终值各时段均达到终值
[0292]
采用cdc方法改进的方法求解时，可使两次计算之间的时间间隔变短，初值与最优解之间的偏差也变小，因而迭代次数减少，这一特点可以有效地改进算法的计算效率，可大大缩短运算时间。算例计算过程中，与单独使用ipm求解相比，cdc-ipm求解采用的cdc准则，使约束违反略大于传统的h准则的结果，同时目标函数值略小于最优值，但计算时间缩短了约50％。采用cdc准则改进及ipm-sqp交替求解的方法使算例的12个时段的变负荷计算时间减少了70％左右，有效地提升了场景应用的实时性。
[0293]
2.仿真变负荷效果分析：
[0294]
针对综合能源系统应用场景中连续变负荷的需求，选取chp机组供热工况下12个时段连续改变电、热负荷的控制过程进行仿真，每个时段300s。采用本发明提出的方法(cdc-ipm-sqp)进行连续的动态优化求解，该过程中控制量轨迹如图3(1)、图(2)、图(3)所示与输出量曲线如图3(4)、图3(5)、图3(6)所示。
[0295]
结果表明，控制变量v
t
、vb、vh轨迹使chp机组变负荷动态过程优化控制实现连贯的规划性，为机组运行时反复调整流量与阀门开度提供精确量化指标，可指导实际chp机组连续变负荷调节。chp机组电功率ph、热功率(相关变量mh)的终端值均可准确的实现变负荷需
求，机组主蒸汽压力ψ
t
、供热蒸汽流量mh曲线连续的刻画了变负荷过程引起的压力、热量波动变化过程，有利于综合能源系统计算场景对安全、稳定运行的分析。

技术特征：
1.基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于包括以下步骤：步骤1：将发电厂热电联产chp机组变负荷优化控制问题表述成命题；步骤2：求解步骤1中的动态优化命题，利用基于legendre多项式离散格式的配置法同时离散状态变量和控制变量，形成nlp问题；步骤3：设计基于sqp-ipm交替求解方法进行nlp问题求解，并在交替求解方法中引入了收敛深度控制策略；通过上述步骤，实现chp机组变负荷过程动态优化。2.根据权利要求1所述基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于：所述步骤1中，动态优化命题的形式为：式(1)中：min j(x(t0),t0,x(t
f
),t
f
)表示最小化动态优化命题关于状态变量的目标函数；t0、t
f
分别为机组优化控制的起始、终止时间点；t表示控制时域时间变量；x(t)、y(t)、u(t)分别表示t时刻chp机组微分状态变量、代数输出变量、控制变量；f(x(t),y(t),u(t),t)表示状态变量、输出变量和控制变量之间的关系函数模型；g
e
、g
ie
为等式与不等式路径约束；h
e
、h
ie
为在时刻t
f
的终值约束；x(t0)、y(t0)、x(t
f
)、y(t
f
)分别表示状态变量初值、输出变量初值、状态变量终值、输出变量终值；最优控制性能目标函数描述为：式(2)中：t0和t
f
分别为机组优化控制的起始和终止时间点；u为chp机组控制变量组成的向量；y为chp机组代数输出变量组成的向量，y
*
为稳态优化得到的输出变量的设定值；r和s是正定权矩阵，能够通过对变量做归一化来确定初始值；j为目标函数，u(t)为t时刻控制变量值、u(t-1)为t-1时刻控制变量值；y(t)为t时刻输出变量值；热电联产机组变负荷优化控制问题中，u为控制变量，y为输出变量，x为状态变量，其中：式(3)中：x为状态变量，对应ψ
d
为汽包压力，p
h
为机组发电功率，ψ
z
为汽轮机供热抽汽压
力，v
f
为锅炉燃烧率；u为控制变量，对应v
t
为汽轮机高压缸进汽调节阀开度，v
b
为机组燃料量，v
h
为供热抽汽调节蝶阀开度；y为输出变量，对应ψ
t
为汽轮机前压力，p
h
为机组发电功率，m
h
为供热抽汽流量。3.根据权利要求1所述基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于：所述步骤1中，
①
、chp机组机理模型如下：v
m
(t)＝v
b
(t-t
b
)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)；式(4)～式(5)表征燃煤制粉系统动态关系，v
b
(t)为机组给煤质量流量，v
m
(t)为制粉系统中实际进入磨煤机的给煤质量流量，t
b
为制粉过程迟延时间常数，v
f
(t)为锅炉燃烧率，t
f
为制粉惯性时间常数，t表示控制过程时域值，v
b
(t-t
b
)表示t-t
b
时刻的机组给煤质量流量；ψ
d
(t)-ψ
t
(t)＝k2(k
1vf
(t))
1.5
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)；式(6))～式(7)表征汽包锅炉系统动态关系，v
t
(t)为汽轮机调门开度，ψ
d
(t)为汽包压力，ψ
t
(t)为汽轮机前压力，c
d
为锅炉蓄热系数，δ为常系数，k1为燃料增益静态参数，k2为过热器阻力系数，k3为汽轮机增益静态参数；ψ1(t)＝0.01ψ
t
(t)v
t
(t)
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)；式(8)～式(10)表征汽轮机系统动态关系，v
h
(t)为抽汽调节蝶阀开度，p
h
(t)为机组发电功率，ψ
z
(t)为中压缸排汽压力，ψ1(t)为汽轮机一级压力，m
r
(t)为循环水质量流量，ε
r
(t)为循环水回水温度，t
t
为汽轮机惯性时间常数，c
z
为热网加热器蓄热系数，k4为高中压缸占汽轮机做功比例，k5为低压缸增益静态参数，k6为热网循环水的有效比热容；m
h
(t)＝k7k6m
r
(t)(96ψ
z
(t)-ε
r
(t)+103)
ꢀꢀꢀꢀ
(11)；式(11)表征供热系统动态关系，m
h
(t)为供热抽汽流量；k7为供热抽汽有效热量折合蒸汽流量系数；
②
、chp机组控制模型：设计一种比例-积分控制算法协同调节v
b
(t)、v
t
(t)、v
h
(t)的方向与大小，使输出变量ψ
t
(t)、m
h
(t)、p
h
(t)变化，其数学模型如下：其中，k
pt
、k
it
分别为控制v
t
(t)变量的比例系数、积分系数；k
pb
、k
ib
、分别为控制v
b
(t)变
量的比例系数、积分系数；k
ph
、k
ih
、分别为控制v
h
(t)变量的比例系数、积分系数；sp
ψ
、sp
p
、sp
m
分别为主蒸汽、机组发电功率、供热抽汽流量的设定值；ψ
t
(t)为主蒸汽压力控制过程当前值，p
h
(t)为机组发电功率控制过程当前值，m
h
(t)为供热抽汽流量控制过程当前值；考虑到chp机组运行的安全稳定，在控制过程中输出变量过程值的波动范围、稳态值的误差范围不能超出允许的限值，应满足如下约束：误差范围不能超出允许的限值，应满足如下约束：其中，
△
ψ
、
△
p
、
△
m
分别为主蒸汽压力、电功率、抽汽蒸汽流量的波动范围；δ
ψ
、δ
p
、δ
m
分别为主蒸汽压力、电功率、抽汽蒸汽流量的误差范围；优化控制周期内，ψ
t
(t
f
)、p
h
(t
f
)、m
h
(t
f
)分别为终点时刻t
f
主蒸汽压力、电功率、抽汽蒸汽流量的终值；
③
、边界约束：chp机组在[t0,t
f
]时段内完成变负荷调节，各控制变量、输出变量的初始状态和终止状态如下：态如下：式中，t0和t
f
分别为机组优化控制的起始和终止时间点；分别为起始时刻t0输出量汽轮机前压力、机组发电功率、供热抽汽流量的取值；分别为终止时刻t
f
输出量汽轮机前压力、机组发电功率、供热抽汽流量的取值；分别为起始时刻t0控制量高压缸进汽调节阀开度、机组燃料量、抽汽调节蝶阀开度的取值；
分别为终止时刻t
f
控制量高压缸进汽调节阀开度、机组燃料量、抽汽调节蝶阀开度的取值；
④
、路径约束：在变负荷调节过程中，控制变量v
t
(t)、v
b
(t)、v
h
(t)需要满足以下约束：式中，和分别为汽轮机调门开度变化的最小值和最大值；和分别为机组燃料量变化的最小值和最大值；和分别为抽汽调节蝶阀开度变化的最小值和最大值。4.根据权利要求1所述基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于：所述步骤2包括如下步骤：1)时域变换：将最优控制时域t∈[t0,t
f
]划分为p段有限元，即0＝t0<t1<
…
t
q-1
<t
q
<
…
t
p
＝t
f
，由此可知，每段时间域的步长为：
△
t＝t
q-t
q-1
,q∈[1,p]
ꢀꢀꢀꢀꢀꢀꢀꢀ
(18)；其中：
△
t为每段时间域的步长，q为有限元个数，p为有限元总段数，t
q
为第q个有限元的时间变量，t
q-1
为第q-1个有限元的时间变量；在高斯伪谱法中配置点区间被选取为τ∈[-1,1]，因此通过时域变换将第q段上时间变量转换为变量τ，其变换如下：其中：τ为第q段上时间变量通过时域变换转换的变量；则每段时间区域内的最优控制问题的描述如下：式(20)中：j
(q)
为第q段时间区域上的目标函数，τ0为第q段上时间变量时域变换的初始配置点，τ
f
第q段上时间变量时域变换的终端配置点；x
(q)
(τ)为第q段时间上第τ个配置点的微分状态变量，x
(q)
(τ0)为第q段时间上第τ0个配置点的微分状态变量，x
(q)
(τ
f
)为第q段时间上第τ
f
个配置点的微分状态变量；y
(q)
(τ)为第q段时间上第τ个配置点的代数输出变量，y
(q)
(τ0)为第q段时间上第τ0个配置点的代数输出变量，y
(q)
(τ
f
)为第q段时间上第τ
f
个配置点的代数输出变量；u
(q)
(τ)为第q
段时间上第τ个配置点的控制变量，f为daes动态模型，
△
t为每段时间域的步长；2)变量离散化：高斯伪谱法是同时将控制变量和状态输出变量在legendre-gauss(lg)配置点上进行离散化，再通过构建拉格朗日插值函数逼近控制变量和状态变量，在区间[τ0,τ
f
]上，n阶的legendre多项式为：其中：p
n
(τ)为n阶legendre多项式，n为阶次，为函数对τ求n阶导数；以式(21)的零点和τ0组成n+1个配置点，在每个t
(q)
时间域中，微分状态变量、控制变量由多项式近似表示为：其中：x
(q)
(τ)为第q段时间上第τ个配置点的微分状态变量，为状态变量的插值函数，x
(q)
(τ
i
)为第q个有限元第i个配置点上状态变量的值；u
(q)
(τ)为第q段时间上第τ个配置点的控制变量，为控制变量的插值函数，u
(q)
(τ
i
)为第q个有限元第i个配置点上控制变量的值；式(22)中，lagrange插值基函数为：其中：n为插值阶次；为状态变量的插值函数，为控制变量的插值函数，第j个配置点上插值函数自变量，τ
i
为第i个配置点上插值函数自变量；τ为第q段上时间变量通过时域变换转换的插值函数自变量；由式(22)状态变量求导，可将插值函数微分逼近状态变量的微分形式，得到在τ
k
点的状态变量导数，可近似表示为：其中：为第p段有限元第k个配置点上微分状态变量的导数，为第p段有限元第k个配置点上状态变量的插值函数，x
(q)
(τ
i
)为第q个有限元第i个配置点上状态变量；结合式(20)和式(24)可得微分代数等式约束：其中：
△
t为每段时间域的步长；f为daes动态模型，x(τ
k
)为第τ
k
个配置点上状态变量，m(τ
k
)第τ
k
个配置点上其他变量，τ
k
为第q段上时间变量时域变换的第k个配置点；
每段终端时刻状态变量作为下一段初始时刻状态变量，其连接处约束等式为：x
(q-1)
(τ
f
)＝x
(q)
(τ0)
ꢀꢀꢀꢀꢀꢀ
(26)；其中：x
(q-1)
(τ
f
)为为第q-1个有限元第τ
f
个配置点上状态变量，x
(q)
(τ0)为第q个有限元第τ0个配置点上状态变量；由于lg配置点中不包含终端时刻状态值，由式(26)，根据legendre-gauss型求积公式计算终端时刻状态变量：式中，a
k
和f(τ
k
)分别为高斯积分权重系数和求积函数。5.根据权利要求1所述基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于：所述步骤3中，收敛深度控制策略具体如下：原命题式(20)所示的非线性规划问题可简写为：其中，为标量目标函数，x∈r
n
为n维实数向量，c∈r
m
为m维等式方程；问题(28)的解必须满足如下一阶最优性条件:其中，为目标函数的梯度向量，为约束方程的雅可比矩阵；λ为对应于约束方程的拉格朗日乘子，ν为对应于x的对偶变量，x＝diag(x)；假定由优化算法产生的迭代序列为{x
k
}，定义当前迭代点x
k
可行性误差目标函数预测改进量数预测改进量数预测改进量其中，为状态向量的第i个分量；c(x
k
)为系统的约束条件，x
k
为当前迭代点，x
k-1
为前一个迭代点，为当前迭代点梯度向量，d
k
为当前迭代点搜索方向；引入s为变形的sigmoid函数柔化传统的刚性收敛准则，其定义为：s(δ
k
,ε0)＝[tanh(ξ
*
(logδ
k
/logε0))]/[tanh(ξ)]
ꢀꢀꢀꢀꢀꢀ
(32)其中，函数s(δ
k
,ε0)将区间[ε0,1/ε0]光滑地连接，ξ用于衡量函数s在该区间内的平缓程度；则收敛深度则收敛深度其中，ε0为给定的容差，s为变形的sigmoid函数，为可行性误差，为目标函数
预测改进量。6.根据权利要求1所述基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于：所述步骤3中，sqp-ipm交替规划算法具体如下：引入进展程度指标用于sqp-ipm交替规划算法之间的切换；可行性变化量以及目标函数变化量为：为：其中，分别为当前迭代点、前一个迭代点的可行性误差，f(x
k
)、f(x
k-1
)分别为x
k
、x
k-1
迭代点目标函数值；由上述公式可将和进展程度定义为：其中，为可行性变化量，为目标函数变化量，ε0为给定的容差。7.根据权利要求6所述基于sqp-ipm交替求解的chp机组变负荷过程动态优化方法，其特征在于：基于引入的收敛深度和进展程度对基于cdc的sqp-ipm交替规划算法框架进行设计，具体步骤如下：步骤(1)：令k＝0，给定初始点x0，指定收敛深度阈值θ
c
，进展程度阈值θ
p
，收敛深度推进步长d
θ
，最大收敛深度值θ
max
，并将ipm的标志位设为0；步骤(2)：选取ipm，将其标志位设为1，用于求解优化问题；步骤(3)：迭代一步，k＝k+1，若ipm的终止准则满足，则求解成功，整个算法终止，否则转步骤(4)；步骤(4)：计算收敛深度，检查是否满足，若满足，则将当前迭代点x
k
赋给x0，并将sqp的标准位都清零，转步骤(6)，否则转步骤(5)；步骤(5)：若ipm的终止准则求解失败未被满足，则转步骤(3)，否则求解失败整个算法终止；步骤(6)选取sqp并将其标志位设为1，以x0为初始点求解优化问题；步骤(7)迭代一步，k＝k+1，若sqp的终止准则满足，则求解成功，整个算法终止，否则转步骤(8)；步骤(8)计算进展程度，若和sqp的终止准则都不满足，则转步骤(7)，否则转步骤(9)；步骤(9)将当前迭代点x
k
赋给x0，转步骤(10)；步骤(10)令θ
c
＝θ
c
+dθ，若θ
c
≥θ
max
则求解失败整个算法终止，否则将ipm的标志位清零，转步骤(2)。

技术总结
基于SQP-IPM交替求解的CHP机组变负荷过程动态优化方法，包括以下步骤：将发电厂热电联产CHP机组变负荷优化控制问题表述成命题；求解中的动态优化命题，利用基于Legendre多项式离散格式的配置法同时离散状态变量和控制变量，形成的非线性规划问题；设计基于SQP-IPM交替求解方法进行非线性规划问题求解，并在交替求解方法中引入了收敛深度控制策略提高方法实时性；以某300MW抽汽式CHP机组变负荷控制需求为算例进行仿真，验证了算法的可行性和连续控制过程快速收敛性。本发明能够有效解决CHP机组变负荷动态过程最优控制问题，且连续控制过程能够快速收敛。控制过程能够快速收敛。控制过程能够快速收敛。


技术研发人员：陈庆 黄悦华 张磊 叶婧 卢天林 刘兴韬 张子豪
受保护的技术使用者：三峡大学
技术研发日：2022.07.22
技术公布日：2022/11/1