多区互联电力系统负荷频率的约束GPC优化控制方法

摘要

本发明公开一种多区互联电力系统负荷频率的约束GPC优化控制方法，通过小信号机理分析建模方法和Z变换建立多区互联电力系统负荷频率的受控自回归积分滑动平均模型，采用约束广义预测控制作为控制方法，将系统预测输出与参考轨迹之间误差的平方值、控制器输出增量的平方值和控制输出的平方值这三者的加权值作为优化目标函数，并设计一种高效的自适应进化优化求解器实现约束GPC的在线滚动优化，实时更新系统各个变量从而实现负荷频率优化控制。采用本发明可实现多区互联电力系统负荷频率优化控制效果，相比现有技术具有更优良的稳态性能和动态性能以及更健壮的鲁棒性。

主权项

一种多区互联电力系统负荷频率的约束GPC优化控制方法，其特征在于，该方法包括以下步骤：(1)通过小信号机理分析建模方法建立多区互联电力系统负荷频率的状态空间模型：$\stackrel{·}{X}=AX+BU+EW---\left(1\right)$Y＝CX    (2)其中，状态变量X＝[X₁,X₂,…,X_n]，系统输入U＝[U₁,U₂,…,U_n]，负荷扰动W＝[W₁,W₂,…,W_n]，系统输出Y＝[Y₁,Y₂,…,Y_n]，n表示分区的数量，第i个状态向量为X_i＝[Δf_i,ΔP_mi,ΔP_vi,ΔP_ti]^T，Δf_i表示第i个分区的频率偏差，ΔP_mi表示第i个分区的发电功率偏差，ΔP_vi表示第i个分区的调节阀位置偏差，ΔP_ti表示第i个分区的联络线功率偏差；第i个系统输入为U_i＝ΔP_ri，ΔP_ri表示第i个分区的负荷参考设定值；第i个负荷扰动为W_i＝ΔP_Li，ΔP_Li表示第i个分区的负荷扰动偏差；第i个系统输出为Y_i＝ACE_i＝β_1iΔf_i+β_2iΔP_ti，ACE_i表示第i个区域控制偏差，β_1i和β_2i分别表示第i个分区频率偏差和联络线功率偏差对应的权重系数；A、B、C、E分别表示多区互联电力系统的状态矩阵、输入矩阵、输出矩阵和扰动矩阵，具体描述如下：其中，D_i表示第i个分区发电机的负荷频率变化量比值系数，M_ai表示第i个分区发电机组的总惯性时间系数，T_Gi表示第i个分区的调速器时间常数，T_CHi表示第i个分区的汽轮机时间常数，T_ij表示第i个分区和第j个分区之间联络线功率同步系数，R_fi表示第i个分区的发电机组的下垂系数。(2)通过Z变换建立多区互联电力系统负荷频率的受控自回归积分滑动平均模型：$\tilde{A}\left(z^{-1}\right)y\left(k\right)=\tilde{B}\left(z^{-1}\right)u(k-1)+\tilde{E}\left(z^{-1}\right)\omega (k-1)+\frac{ϵ\left(k\right)}{{\Delta }_d}---\left(3\right)$其中y(k)表示系统在k时刻的输出信号，u(k‑1)表示在k‑1时刻的控制输出信号，ω(k‑1)表示在k‑1时刻的负荷扰动信号，ε(k)表示均值为0且方差为δ²的n维白噪声，Δ_d表示差分因子，Δ_d＝1‑z^‑1；表示关于z^‑1的n×n维对角矩阵，表示关于z^‑1的n×m维矩阵，表示关于z^‑1的n×p维矩阵，其中的第i个对角元素的第i行第j列元素的第i行第j列元素的表达式如下：${\tilde{A}}_i\left(z^{-1}\right)=1+a_{i1}{z}^{-1}+...+a_{{\mathrm{in}}_{a_i}}{z}^{-n_{a_i}}$${\tilde{B}}_{ij}\left(z^{-1}\right)=b_{ij0}+b_{ij1}{z}^{-1}+...+b_{{\mathrm{ijn}}_{b_{ij}}}{z}^{-n_{b_{ij}}}$${\tilde{E}}_{ij}\left(z^{-1}\right)=e_{ij0}+e_{ij1}{z}^{-1}+...+e_{{\mathrm{ijn}}_{b_{ij}}}{z}^{-n_{b_{ij}}}$其中，表示各阶次的系数，表示各阶次的系数，表示各阶次的系数，n_ai和n_bij分别表示和的最高阶次。(3)设置参数值：预测时域长度N，控制域长度N_u，权重系数矩阵Q、R₁和R₂,柔化系数矩阵α^j，种群规模p，变异系数b，最大迭代次数I_max。(4)读取当前时刻k(k＝1,2,…,t_max)的历史信息，包括系统输出信号y(k‑1),…,y(k‑n_a)和控制器输出信号u(k‑1),…,u(k‑n_b)，t_max表示系统运行时窗的最大值。(5)随机产生满足约束条件的实数编码的种群P＝{S₁,S₂,…,S_p}，其中第i个体S_i表示待优化的控制增量序列{Δu(k),Δu(k+1),…,Δu(k+N_u‑1)}，具体产生过程如下：S_i＝Δu_min+r_i(Δu_max‑Δu_min),i＝1,2,...,p其中，Δu_min和Δu_max分别表示控制增量序列的下限和上限，r_i表示一组在0和1之间产生的随机数。(6)按照公式(4)～(10)所示的约束GPC优化目标对种群P进行适应度函数评价，获得当前时刻k种群P的最好适应度函数F_best＝min{F(P,k)}，将对应的个体设置为最好个体S_best，其中F(P,k)表示当前时刻k种群P的适应度集合；$\begin{array}{c}\underset{P}{\min }F\left(k\right)=\underset{P}{\min }[\underset{j=1}{\overset{N}{\Sigma }}{\left(y\right(k+j)-\omega (k+j\left)\right)}^{T}Q\left(y\right(k+j)-\omega (k+j\left)\right)+\\ \underset{i=1}{\overset{N_u}{\Sigma }}{\left(\Delta u\left(k+i-1\right)\right)}^T{R}_1\Delta u\left(k+i-1\right)+\underset{i=1}{\overset{N_u}{\Sigma }}{\left(u\left(k+i-1\right)\right)}^T{R}_{2}u(k+i-1)]\end{array}---\left(4\right)$s.t.Δu_min≤Δu(k+i‑1)≤Δu_max,i＝1,2,...,N_u   (5)u_min≤u(k+i‑1)≤u_max,i＝1,2,...,N_u  (6)y_min≤y(k+j)≤y_max,j＝1,2,...,N  (7)Δu(k+i‑1)＝u(k+i‑1)‑u(k+i‑2)  (8)ω(k+j)＝α^jω(k+j‑1)+(1‑α^j)y_r(k),ω(k)＝y_p(k),j＝1,2,...,N  (9)y(k+j)＝G_j(z^‑1)Δu(k+j‑1)+H_j(z^‑1)Δu(k‑1)+F_j(z^‑1)y(k),j＝1,2,...,N  (10)其中，y(k+j)表示系统的预测输出信号，ω(k+j)表示参考轨迹，y_r(k)表示设定输入信号，y_p(k)为系统实际输出信号，Q、R₁、R₂表示权重系数矩阵，u_min和u_max分别表示控制序列的下限和上限，y_min和y_max分别表示预测输出信号的下限和上限，公式(10)所示的y(k+j)的计算方法为求取如下丢番图方程(11)～(12)获得：$\begin{array}{c}1=L_j\left(z^{-1}\right)\tilde{A}\left(z^{-1}\right)\Delta +z^{-j}{F}_j\left(z^{-1}\right)\\ \left\{\begin{array}{c}{L}_j\left(z^{-1}\right)=l_0+l_1{z}^{-1}+\mathrm{...}+l_{j-1}{z}^{-j+1},\\ F_j\left(z^{-1}\right)=f_{j0}+f_{j1}{z}^{-1}+\mathrm{...}+f_{{\mathrm{jn}}_a}{z}^{-n_a}.\end{array}\right.\end{array}---\left(11\right)$$\begin{array}{c}{\tilde{E}}_j\left(z^{-1}\right)\tilde{B}\left(z^{-1}\right)=G_j\left(z^{-1}\right)+z^{-j}{H}_j\left(z^{-1}\right)\\ \left\{\begin{array}{c}{G}_j\left(z^{-1}\right)=g_0+g_1{z}^{-1}+\mathrm{...}+g_{j-1}{z}^{-j+1},\\ H_j\left(z^{-1}\right)=h_{j0}+h_{j1}{z}^{-1}++h_{j(n_b-1)}{z}^{-n_b+1}.\end{array}\right.\end{array}---\left(12\right)$其中，l₀,l₁,…,l_j‑1表示多项式L_j(z^‑1)各阶次的系数，表示多项式F_j(z^‑1)各阶次的系数，g₀,g₁,…,g_j‑1表示多项式G_j(z^‑1)各阶次的系数，表示多项式H_j(z^‑1)各阶次的系数，n_a和n_b分别表示F_j(z^‑1)和H_j(z^‑1)的最高阶次。(7)按照公式(13)～(14)所示的实数变异操作因子产生新的种群P_m＝{S_m1,S_m2,…,S_mp}，并无条件地接受P＝P_m；$S_{mi}=\left\{\begin{array}{c}{S}_i+\gamma (t,{\mathrm{\Delta u}}_{\max }-S_i),ifr<0.5\\ S_i+\gamma (t,S_i-{\mathrm{\Delta u}}_{\min }),ifr\ge 0.5\end{array}\right.---\left(13\right)$$\gamma (t,\alpha )=\alpha (1-{r_1}^{{(1-(t-I_{\max }\left)\right)}^b})---\left(14\right)$其中r、r₁是[0,1]范围内产生的随机数，t表示当前迭代次数，b为变异操作的可调系数，I_max为用户设定的最大迭代次数，γ(t,α)表示关于t和α的函数，α的具体描述如下：当r<0.5时，α＝Δu_max‑S_i，否则α＝Δu_max‑S_i。(8)重复步骤6～7直到优化求解器的迭代次数达到I_max，并保存S_best＝{Δu*(k+j),j＝0,1,…,N_u‑1}。(9)计算在最优控制增量Δu*(k)作用下的当前时刻k对应的系统输出y*(k)和u*(k)，并设置为Δu(k)＝Δu*(k)，y(k)＝y*(k)，u(k)＝u*(k)。(10)重复步骤4～9直到系统运行时窗达到t_max，输出多区互联电力系统最优系统输出曲线、最优控制增量信号曲线和最优控制信号曲线。