污水处理过程递推双线性子空间建模及无模型自适应控制

作者：张帅周平来源：《自动化学报》日期：2022-08-26人气：655

活性污泥法是目前使用最为广泛的污水处理方法, 其工艺流程如图1所示. 由污水管网收集的污水首先经过格栅、洗砂间、初沉池等工序去除碎石、砂子等固体悬浮物; 然后, 污水在生化池内通过厌氧反硝化作用和好氧硝化作用, 分解有机氮和氨氮, 去除大部分污染物; 最后, 经过生化反应处理的污水在二沉池内进行二次沉淀, 进一步提高水质, 从而排入收纳河道. 图1所示污水处理过程中, 生化池厌氧区第二单元硝态氮浓度 $S_{N O, 2}$ 和好氧区第五单元溶解氧浓度 $D_{O, 5}$ 是关乎出水质量优劣的最重要过程变量[1-3]. 其中 $D_{O, 5}$ 受控于五区溶解氧转换系数 $K_{L a, 5}$ , $S_{N O, 2}$ 受控于由五区回流至二区的内回流流量 $Q_{a}$ , 因此它们之间有着较强的耦合性, 这在一定程度上增加了 $S_{N O, 2}$ 和 $D_{O, 5}$ 建模和控制的难度. 同时, 由于生化反应过程的时变性、不确定性和复杂性, 使得 $S_{N O, 2}$ 和 $D_{O, 5}$ 的机理建模和基于模型的控制难以实现.

图 1 污水处理工艺流程图

Fig. 1 Wastewater treatment process

目前, 针对难以建立机理和经验模型的复杂工业过程, 数据驱动子空间辨识方法得到了越来越广泛的应用[4-8]. 子空间辨识技术不需要系统的先验知识, 可直接利用输入输出数据, 通过统计、优化和线性代数的思想得到系统的动态模型, 可有效避免传统建模过程代价昂贵、计算复杂、应用范围较窄的缺点, 因而具有普遍的适用性和应用价值[9-10]. 目前, 子空间辨识在污水处理过程中的应用主要集中于离线或线性子空间辨识建模[6, 11], 而污水生化反应过程是一个典型的时变、非线性过程; 因此, 基于线性子空间辨识的污水生化过程建模和控制效果存在一定的局限性. 为此, 本文采用一种递推最小二乘双线性子空间辨识方法, 建立污水生化池 $S_{N O, 2}$ 和 $D_{O, 5}$ 的具有参数自适应能力的数据驱动递推双线性模型, 以更好地捕捉污水处理过程的非线性时变动态特性.

近年来, 针对污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 的控制问题, 基于数据驱动的控制策略引起了广泛关注, 如基于数据驱动的自适应神经网络控制[12]、模糊神经网络控制[13-14], 神经网络模型预测控制[15-16]等. 这些方法虽然可以在不使用精确数学模型的条件下实现 $S_{N O, 2}$ 和 $D_{O, 5}$ 的自适应控制, 但污水生化反应过程的模糊规则建立困难, 神经网络模型则需要大量的污水运行数据进行训练. 此外, 预测控制在每一时刻都要求解一个非线性约束优化问题, 计算成本较高. 同时, 上述算法控制器的设计依赖于模糊规则或神经网络训练模型, 因此基于模型的控制器设计中的未建模动态等问题依然存在. 其他在实际系统中应用较好的数据驱动控制和自适应控制方法, 如文献[17]事件驱动下的自学习最优调节控制和文献[18]自适应非线性鲁棒控制, 同样存在未建模动态等与上述方法中相同的问题. 在这种背景下, 针对污水处理过程研究直接数据驱动控制的无模型自适应控制(Model-free adaptive control, MFAC)具有十分重要的意义. MFAC于1994年由侯忠生教授首次提出[19], 并在后期不断发展和完善[20-25], 现已形成一套新的针对非线性动态系统的控制理论与方法. MFAC利用非线性系统的输入输出数据以及动态线性化技术[26], 在受控系统的每个工作点处, 建立虚拟动态线性化数据模型, 然后基于此虚拟模型设计加权一步向前控制器. 整个控制器设计过程仅依靠输入输出数据和动态线性化技术, 不需要过程先验知识和系统模型的参与, 且控制器结构简单、计算方便, 可有效解决未建模动态问题.

综上所述, 本文针对污水处理过程难以机理建模和直接数据驱动控制的问题, 采用一种改进的MFAC方法, 研究污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 的直接数据驱动控制. 由于MFAC控制器待调参数众多, 且参数的取值直接影响控制系统的性能指标, 为改善MFAC控制器参数的调整效率, 提高MFAC控制器性能, 采用基于蒙特卡罗(Monte Carlo)实验的多参数灵敏度分析(Multi-parameter sensitivity analysis, MPSA)和遗传粒子群优化(Genetic algorithm−particle swarm optimization, GA-PSO)方法整定控制器待调参数, 对常规MFAC方法进行改进, 完成的主要工作如下:

1) 由于在控制器参数的整定过程中, 需要基于污水处理过程做数千次Monte Carlo实验和数百次优化计算, 而基于实际污水处理系统的参数整定实验代价高昂, BSM1 (Benchmark simulation model No.1)标准模型[27]运行又较为复杂. 因此, 利用一种递推双线性子空间辨识(Recursive bilinear subspace identification, RBLSI)方法, 建立污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 的RBLSI模型, 与递推线性子空间辨识(Recursive linear subspace identification, RLSI)的模型测试效果对比后, 选择效果更优的RBLSI模型作为MFAC控制器参数整定的被控对象模型.

2) 采用紧格式动态线性化无模型自适应控制(Compact form dynamic linearization MFAC, CFDL-MFAC)的方法, 设计 $S_{N O, 2}$ 和 $D_{O, 5}$ 的CFDL-MFAC加权一步向前控制器; 基于 $S_{N O, 2}$ 和 $D_{O, 5}$ 的递推双线性子空间模型, 开展 $S_{N O, 2}$ 和 $D_{O, 5}$ 的CFDL-MFAC控制器待调参数灵敏度分析的Monte Carlo实验, 并采用MPSA方法, 依据实验结果选出控制器的敏感参数.

3) 为不敏感参数取适当的值, 以敏感参数为优化变量, 采用GA-PSO优化算法, 做CFDL-MFAC控制器敏感参数的优化整定实验.

4) 将控制器参数整定的结果作用于CFDL-MFAC控制器, 开展污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 的设定值跟踪和干扰抑制的控制实验研究, 以充分验证所提方法的有效性和先进性.

本文基于多参数灵敏度分析(MPSA)和遗传粒子群优化(GA-PSO)算法的无模型自适应控制(MFAC)方法在硝态氮和溶解氧浓度的控制中具有以下优势. 首先, 控制器的设计仅需污水系统的 $Q_{a}$ 和 $K_{L a, 5}$ 数据以及 $S_{N O, 2}$ 和 $D_{O, 5}$ 数据, 不需要污水系统的模型信息, 可避免传统的未建模动态问题的存在, 非常适用于污水处理系统的实际应用. 其次, 本文算法控制器结构简单、计算量小、控制器设计不需精确的污水模型, 控制器计算成本较低. 最后, 本文方法可实现带有参数变化的非线性污水处理系统 $S_{N O, 2}$ 和 $D_{O, 5}$ 浓度的直接数据驱动控制. 因此, 本文方法无论在污水处理的理论研究方面还是实际应用方面都具有重要意义.

1. $S_{N O, 2}$ 和 $D_{O, 5}$ 的递推双线性子空间建模

为更加精确地描述污水处理过程的非线性时变动态特性, 为后期CFDL-MFAC控制器参数灵敏度分析的Monte Carlo实验和优化整定实验做好模型准备, 首先采用带有遗忘因子的递推最小二乘技术, 改进常规双线性子空间辨识方法[28-29], 建立 $S_{N O, 2}$ 和 $D_{O, 5}$ 的RBLSI模型, 并与文献[7]中递推线性子空间辨识(RLSI)建模方法做模型测试对比, 根据两种方法的模型测试结果及其对比分析, 选择效果更好的RBLSI模型作为控制器参数Monte Carlo实验和优化整定实验的被控对象模型, 具体过程如下.

1.1 基于最小二乘的递推双线性子空间辨识建模

1.1.1 递推双线性子空间辨识矩阵构造

用如下形式的MIMO双线性系统描述污水处理生化池内 $S_{N O, 2}$ 和 $D_{O, 5}$ 的变化过程:

_{}_{}_{}_{}_{}_{}

(1)

(2)

式中, $u_{k} \in R^{m}$ , $y_{k} \in R^{l}$ , $x_{k} \in R^{n}$ 分别表示系统的输入、输出和状态向量, $e_{k}$ 表示独立于输入 $u_{k}$ 的零均值白噪声序列, 其方差 $E [e_{p}, e_{q}] = S δ_{p q} > 0$ . $A$ , $B$ , $C$ , $D$ , $K$ 为具有适当维数的系统矩阵, $N \in$ $R^{n \times n m}$ 表示系统的双线性特性矩阵, $\otimes$ 表示克罗内克(Kronecker)积.

将式(2)中 $u_{k}$ 代入式(1), 通过化简, 可将式(1)中 $u_{k} \otimes x_{k}$ 项转化为输入和输入以及输入和输出向量之间的Kronecker积, 由此得到如下具有线性系统结构形式的双线性系统状态空间方程:

_{} {\overset{}{}}_{}

(3)

(4)

式中, $\tilde{B} = [B N P - N Q]$ , $\tilde{K} = [K - N u_{k} \otimes C^{†}],$ $P = d i a g {\underset{m 个分块}{\underset{⏟}{C^{†}, \dots, C^{†}}}}$ 和 $Q = d i a g {\underset{m 个分块}{\underset{⏟}{C^{†} D, \dots, C^{†} D}}}$ 表示两个分块对角阵, $C^{†}$ 代表 $C$ 的伪逆矩阵, ${\tilde{u}}_{k} =$ ${[u_{k}^{T} (u_{k} \otimes y_{k})^{T} (u_{k} \otimes u_{k})^{T}]}^{T}$ 表示广义输入向量.

由于式(3)和式(4)与线性状态空间模型有着相同的结构形式, 因此线性系统子空间辨识的一些原理也可应用于双线性系统(3)和(4)的辨识. 假定 ${A, \tilde{B}}$ 能观, ${A, C}$ 能控, 在不考虑噪声影响前提下, 将式(4)代入式(3), 通过对输出变量的前向迭代. 可得类似于线性系统子空间辨识中的矩阵输入输出方程

_{}_{}_{} {\overset{}{}}_{}^{} {\overset{}{}}_{}_{}^{}_{}

(5)

_{}_{}_{} {\overset{}{}}_{}^{} {\overset{}{}}_{}_{}^{}_{}

(6)

_{}^{}_{}_{} {\overset{}{}}_{}

(7)

式中, $U_{p}$ 和 $U_{f}$ 分别表示过去和当前时刻的输入汉克(Hankel)矩阵,

\begin{matrix} _{} & _{} & _{} \\ _{} & _{} & _{} \\ _{} & _{} & _{} \end{matrix}

\begin{matrix} _{} & _{} & _{} \\ _{} & _{} & _{} \\ _{} & _{} & _{} \end{matrix}

$Y_{p}$ 和 $Y_{f}$ 分别表示过去和当前时刻的输出Hankel矩阵, 其定义与上述 $U_{p}$ 和 $U_{f}$ 的定义相似. ${\tilde{U}}_{p}$ 和 ${\tilde{U}}_{f}$ 分别代表过去时刻和当前时刻的广义输入Hankel矩阵, 其定义类似于 $U_{p}$ 和 $U_{f}$ 的定义, 只需将 $U_{p}$ 和 $U_{f}$ 中的输入向量 $u_{k} (k = 0, \dots, 2 i + j - 2)$ 替换为式(3)中定义的广义输入向量 ${\tilde{u}}_{k} (k = 0, \dots, 2 i + j - 2)$ , 便可得到 ${\tilde{U}}_{p}$ 和 ${\tilde{U}}_{f}$ . $X_{p}$ , $X_{f}$ 分别表示过去和当前时刻的子空间状态矩阵, 其定义分别为

_{} \begin{matrix} _{} & _{} & _{} \end{matrix}

(8)

_{} \begin{matrix} _{} & _{} & _{} \end{matrix}

(9)

$Γ_{i}$ , $Δ_{i}$ 分别表示广义能观和广义能控矩阵, 其定义分别为

{\begin{matrix} ^{} & ^{} & {^{}}^{} \end{matrix}}^{}

\begin{matrix} ^{} \overset{}{} & ^{} \overset{}{} & \overset{}{} \end{matrix}

${\tilde{H}}_{i}^{d}$ , $H_{i}^{d}$ 分别代表分块下三角Toeplitz矩阵和一个分块对角阵, 其定义分别为

\overset{}{} \begin{matrix} \overset{}{} \\ \overset{}{} & \overset{}{} \\ ^{} \overset{}{} & ^{} \overset{}{} & ^{} \overset{}{} \end{matrix}

_{}^{} \begin{matrix}  \end{matrix}

1.1.2 递推参数矩阵求解

将式(5)和式(7)代入式(6), 在忽略噪声影响的情况下, 系统当前输出 $Y_{f}$ 的行空间等价于过去输入输出Hankel矩阵、当前输入Hankel矩阵以及过去和当前广义输入Hankel矩阵行空间的组合, 由此得到系统输出 $Y_{f}$ 的预测方程为

{\overset{}{}}_{}_{}_{}_{}_{} {\overset{}{}}_{} {\overset{}{}}_{}

(10)

式中, ${\hat{Y}}_{f}$ 表示 $Y_{f}$ 的预测, $L_{p}$ , $L_{u}$ , ${\tilde{L}}_{u}$ 为具有适当维数的输入输出参数矩阵, 其具体形式为

\begin{aligned} _{}^{} \begin{matrix} ^{}_{}^{}_{}^{} & ^{}_{}^{} \end{matrix} \\ _{}_{}^{} \end{aligned}

\overset{}{} \begin{matrix} _{}_{}_{}^{}_{}^{}_{} & {\overset{}{}}_{}^{} \end{matrix}

$W_{p} = {[U_{p}^{T} Y_{p}^{T}]}^{T}$ 表示由过去输入输出Hankel矩阵组成的分块矩阵, ${\tilde{W}}_{u} = [{\tilde{U}}_{p}^{T} {\tilde{U}}_{f}^{T}]^{T}$ 表示由过去和当前广义输入Hankel矩阵组成的分块矩阵.

为充分利用输入输出数据, 建立具有参数自适应能力的递推双线性子空间输入输出模型, 引入具有遗忘因子的递推最小二乘求解式(10)中的时变参数矩阵 $L_{p}$ , $L_{u}$ 和 ${\tilde{L}}_{u}$ 的估计值 ${\hat{L}}_{p}$ , ${\hat{L}}_{u}$ 和 ${\hat{\tilde{L}}}_{u}$ , 具体递推求解过程为

{\overset{}{}}_{} {\overset{}{}}_{}_{}_{}^{}_{}^{} {\overset{}{}}_{}

(11)

_{}_{}_{} {_{}^{}_{}_{}}^{}

(12)

_{} \frac{}{}_{}_{}^{}_{}

(13)

式中, ${\hat{L}}_{k} = [{\hat{L}}_{p} {\hat{L}}_{u} {\hat{\tilde{L}}}_{u}]_{k}^{T}$ 为 $k$ 时刻参数矩阵, $φ_{k} = [w_{p}^{T} u_{f}^{T} {\tilde{w}}_{u}^{T}]_{k}$ 为 $k$ 时刻数据向量, $g_{k}$ 代表增益向量, $P_{k}$ 表示协方差矩阵, $λ$ 为遗忘因子. 通过递推式(11) $\sim$ (13)即可得到如下具有参数自适应能力的RBLSI预测模型:

{\overset{}{}}_{} {\overset{}{}}_{}_{} {\overset{}{}}_{}_{} {\overset{}{\overset{}{}}}_{} {\overset{}{}}_{}

(14)

1.2 $S_{N O, 2}$ 和 $D_{O, 5}$ 的递推双线性子空间模型

将第1.1节所述RBLSI方法应用于污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 的建模, 以内回流流量 $Q_{a}$ 和五区溶解氧转换系数 $K_{L a, 5}$ 为控制变量, 即 $u_{1} = Q_{a},$ $u_{2} = K_{L a, 5}$ , $u = {[u_{1} u_{2}]}^{T}$ , 以 $S_{N O, 2}$ 和 $D_{O, 5}$ 为被控变量, 即 $y_{1} = S_{N O, 2}$ , $y_{2} = D_{O, 5}$ , $y = {[y_{1} y_{2}]}^{T}$ , 状态向量 $x$ 是输入数据 $u$ 和输出数据 $y$ 子空间矩阵的非线性组合, 可由式(5)和式(7)求得. 基于污水处理数据, 建立 $S_{N O, 2}$ 和 $D_{O, 5}$ 的递推双线性子空间模型, 并与文献[7]中递推线性子空间辨识(RLSI)建模方法做模型测试对比.

实验共使用实际污水处理数据550个, 其中400个为建模训练用数据, 150个为模型测试用数据. 由本文RBLSI方法得到最终的 $S_{N O, 2}$ 和 $D_{O, 5}$ 预测模型结构如式(14)所示, 基于最新一组数据得到的子空间参数矩阵为

_{} \begin{array}{r}  \end{array}

\overset{}{} {\begin{array}{r}  \end{array}}_{}

_{} \begin{array}{r}  \end{array}

{\overset{}{}}_{} {\overset{}{}}_{}_{} {\overset{}{}}_{}_{}

(15)

由文献[7]算法得到的 $S_{N O, 2}$ 和 $D_{O, 5}$ 预测模型结构如式(15)所示, 基于最新一组数据得到的子空间参数矩阵为

\begin{array}{r}  \end{array}

\begin{array}{r}  \end{array}

与式(14)描述的RBLSI模型相比, 式(15)描述的RLSI模型在结构上缺少了 ${\hat{\tilde{L}}}_{u} {\tilde{w}}_{u}$ 这一部分, 即未考虑到输入、输出数据向量的Kronecker积对输出预测的影响, 对输入、输出数据的利用不充分, 因此其预测效果相对于RBLSI模型来说具有一定的保守性.

由预测模型表达式(14)和式(15), 基于污水数据可得到两种算法的模型测试曲线, 如图2所示. 由两种算法实际的预测误差, 通过核心平滑密度估计技术, 可得到两种算法预测误差的概率密度函数(Probability density function, PDF)形状曲线, 如图3所示. 图3曲线服从高斯分布, 使用数据拟合工具, 通过核密度估计技术和高斯拟合算法, 可得到文献[7]的RLSI算法和本文RBLSI算法预测误差的PDF表达式, 分别如式(16), (17)和式(18), (19)所示.

图 2 不同算法

S_{N O, 2}

和

D_{O, 5}

的模型预测效果

Fig. 2 Model prediction effects of

S_{N O, 2}

and

D_{O, 5}

with different algorithms

图 3 不同算法

S_{N O, 2}

和

D_{O, 5}

的预测误差PDF曲线

Fig. 3 Error PDF shapes of

S_{N O, 2}

and

D_{O, 5}

prediction with different algorithms

^{{\frac{}{}}^{}}^{{\frac{}{}}^{}}

(16)

\begin{aligned} _{} & ^{{\frac{}{}}^{}} \\ ^{{\frac{}{}}^{}} \end{aligned}

(17)

\begin{aligned} _{} & ^{{\frac{}{}}^{}}^{{\frac{}{}}^{}} \\ ^{{\frac{}{}}^{}}^{{\frac{}{}}^{}} \end{aligned}

(18)

\begin{aligned} _{} & ^{{\frac{}{}}^{}} \\ ^{{\frac{}{}}^{}}^{{\frac{}{}}^{}} \end{aligned}

(19)

式中, $e x$ 代表相应变量的预测误差, $f_{p d f} (e x)$ 表示对应该误差的概率.

通过对比分析两种算法的模型测试效果(图2)、预测误差PDF形状(图3)及测试均方根误差(Root mean square error, RMSE) (表1)可知, 在对 $S_{N O, 2}$ 和 $D_{O, 5}$ 的预测过程中, 本文RBLSI算法具有更加接近于真实值的预测输出、形状更好的预测误差PDF曲线以及数值更小的测试均方根误差, 因而能够更加精确地描述 $S_{N O, 2}$ 和 $D_{O, 5}$ 的变化过程. 所以, 最终选择RBLSI模型作为CFDL-MFAC控制器参数灵敏度分析和优化整定实验的被控对象模型.

表 1 不同算法的RMSE比较

Table 1 Comparison of RMSE based on different algorithms

算法	RMSE $(S_{N O, 2})$	RMSE $(D_{O, 5})$
RBLSI	0.0165	0.0046
RLSI	0.0211	0.0085

2. $S_{N O, 2}$ 和 $D_{O, 5}$ 的CFDL-MFAC控制

本节首先给出了CFDL-MFAC控制器的设计方法, 设计了 $S_{N O, 2}$ 和 $D_{O, 5}$ 的CFDL-MFAC控制器. 其次, 为提高控制器参数调整效率, 基于1.2节建立的RBLSI模型, 详细介绍了采用MPSA方法对CFDL-MFAC控制器参数做Monte Carlo实验, 并依据实验结果选择灵敏参数的过程. 最后, 描述了采用GA-PSO方法对控制器灵敏参数做优化整定, 并利用整定好的参数对 $S_{N O, 2}$ 和 $D_{O, 5}$ 浓度进行跟踪控制的过程, 整个系统的控制框图如图4所示. 图4中, ① 表示输出信号流, ② 表示输入信号流, ③ 表示CFDL-MFAC控制器参数更新信号流, 其中包含伪雅可比矩阵 $Φ (k)$ 的初值 $Φ (0)$ , ④ 表示伪雅可比矩阵 $Φ (k)$ 更新信号流.

图 4 控制系统框图

Fig. 4 The control flow chart

2.1 CFDL-MFAC控制器设计

假设污水处理 $S_{N O, 2}$ 和 $D_{O, 5}$ 控制过程的离散非线性输入输出方程为

\begin{aligned} _{} \\ _{} \end{aligned}

(20)

式中, $y_{k} \in R^{l}$ , 代表被控输出 $S_{N O, 2}$ 和 $D_{O, 5}$ , $u_{k} \in R^{m}$ 代表控制输入 $Q_{a}$ 和 $K_{a, 5}$ , $n_{y}$ 、 $n_{u}$ 分别表示系统输出和输入阶次, $f (\cdot)$ 表示一个未知非线性向量函数.

2.1.1 紧格式动态线性化(CFDL)

由文献[20, 25]可知, $S_{N O, 2}$ 和 $D_{O, 5}$ 的CFDL-MFAC控制器的设计需要满足如下两个假设.

假设 1[20, 25]. 系统(20)中的非线性函数 $f (\cdot)$ 对输入 $u (k)$ 存在连续偏导数.

假设 2[20, 25]. 系统(20)满足广义Lipschitz条件, 即在任意时刻 $k$ , 若 $‖ Δ u (k) ‖ = ‖ u (k) - u (k - 1) ‖ \neq$ $0$ , 则 $‖ Δ y (k + 1) ‖ = ‖ y (k + 1) - y (k) ‖ \leq b ‖ Δ u (k) ‖$ , $b$ 是一个大于零的常数.

上述假设条件分别是对常规非线性控制系统设计和被控输出变化率的典型约束, 对实际污水处理系统来说, 由于它满足能量守恒定律, 有限范围内的输入变化必然引起有限范围内的输出变化. 因此, 任意两个时刻 $S_{N O, 2}$ 和 $D_{O, 5}$ 的输出变化量与 $Q_{a}$ 和 $K_{a, 5}$ 的控制变化量之间, 一定存在有界的系数关系. 故污水处理系统满足假设1和假设2, 由此可得如下引理.

引理 1[20, 25]. 对 $S_{N O, 2}$ 和 $D_{O, 5}$ 非线性控制系统(20), 若满足假设1和假设2, 且在任意时刻 $k,$ $‖ Δ u (k) ‖ \neq 0$ , 则存在一个伪雅克比矩阵(Pseudo Jacobian matrix, PJM) $Φ_{c} (k)$ , 可将系统(20)化为一种如下形式的虚拟CFDL模型:

_{}

(21)

式中, $Φ_{c} (k) \in R^{l \times m}$ 是一个时变、有界的参数矩阵.

2.1.2 控制器推导与伪雅可比矩阵(PJM)的估计

控制系统的首要目标是要控制 $S_{N O, 2}$ 和 $D_{O, 5}$ 的浓度跟随其设定值, 保证污水处理系统的顺利运行, 提高出水水质; 同时, 考虑到如果相邻两个时刻控制增量过大会引起系统输入剧烈跳变, 增加执行机构的负担, 因此, 以 $S_{N O, 2}$ 和 $D_{O, 5}$ 跟踪误差的二范数以及控制增量的加权二范数作为控制器推导的目标函数.

{_{}}^{}^{}

(22)

式中, $y_{r}$ 表示系统输出 $S_{N O, 2}$ 和 $D_{O, 5}$ 的目标值, $λ$ 为权重因子, $λ > 0$ .

将虚拟模型(21)代入目标函数(22), 并对输入 $u (k)$ 求偏导, 令其偏导数为零, 得到 $S_{N O, 2}$ 和 $D_{O, 5}$ 的加权一步向前控制器

\frac{_{}^{}_{}}{{_{}}^{}}

(23)

式中, $ρ$ 代表为避免求逆运算引入的步长因子, $0 < ρ \leq 1$ .

由控制器(23)可知, CFDL-MFAC算法的实现依赖于伪雅可比矩阵 $Φ_{c} (k),$ 而 $Φ_{c} (k)$ 是一个未知、时变的参数矩阵, 因此需要利用 $S_{N O, 2}$ 和 $D_{O, 5}$ 控制系统的输入输出数据对 $Φ_{c} (k)$ 进行实时估计. 考虑如下 $Φ_{c} (k)$ 估计的目标函数:

\begin{aligned} {\overset{}{}}_{} & {_{}}^{} \\ {_{} {\overset{}{}}_{}}^{} \end{aligned}

(24)

${\hat{Φ}}_{c} (k)$ 表示 $Φ_{c} (k)$ 的估计, ${\hat{Φ}}_{c} (k - 1)$ 表示 $Φ_{c} (k - 1)$ 的估计, $μ$ 表示权重因子, $μ > 0$ .

对目标函数(24), 关于 $Φ_{c} (k)$ 求极值, 得到 $Φ_{c} (k)$ 的估计式为

\begin{aligned} {\overset{}{}}_{} {\overset{}{}}_{} \\ \frac{{\overset{}{}}_{}^{}}{^{}} \end{aligned}

(25)

式中, $0 < η \leq 2$ 表示为避免求逆运算引入的步长因子.

由于 $S_{N O, 2}$ 和 $D_{O, 5}$ 模型是一个典型的时变系统, 参数变化较大, 为了使伪雅克比矩阵的估计式(25)能够更好地跟踪 $S_{N O, 2}$ 和 $D_{O, 5}$ 控制系统的时变参数, 引入以下算法重置机制:

\begin{aligned} {\overset{}{}}_{} {\overset{}{}}_{} & {\overset{}{}}_{}_{} {\overset{}{}}_{}_{} \\ {\overset{}{}}_{} {\overset{}{}}_{} \\ {\overset{}{}}_{} {\overset{}{}}_{} & {\overset{}{}}_{}_{} \\ {\overset{}{}}_{} {\overset{}{}}_{} \end{aligned}

(26)

式中, ${\hat{ϕ}}_{i j} (k), i = 1, \dots, l, j = 1, \dots, m,$ 表示 ${\hat{Φ}}_{c} (k)$ 的元素, ${\hat{ϕ}}_{i j} (1)$ 为 ${\hat{ϕ}}_{i j} (k)$ 的初值, $s i g n (\cdot)$ 表示符号函数, 常数 $α$ , $b_{1}$ , $b_{2}$ 满足 $α \geq 1$ , $b_{2} > b_{1} (2 α + 1) (m - 1)$ .

至此, 式(23), (25)和式(26)组成 $S_{N O, 2}$ 和 $D_{O, 5}$ 控制系统的CFDL-MFAC控制方法. 其稳定性的证明与文献[20, 25]中基本相同, 此处不再赘述. 假设CFDL-MFAC控制器求解时间为 $t_{0},$ 则该算法控制时间复杂度为 $O (t_{0})$ . 而对基于神经网络模型的预测控制而言, 若求解非线性约束优化问题时间为 $t_{1},$ 预测时域为 $P,$ 则其控制时间复杂度为 $O (P \times t_{1})$ . 此外, 预测控制中非线性约束优化问题一般采用SQP 方法求解, 需要一个基于变尺度法的迭代过程, 而本文算法控制量求解仅需简单的矩阵运算即可, 因此其控制量求解时间 $t_{0}$ 要小于非线性约束优化问题的求解时间 $t_{1}$ , 这在第3.2节实验结果所得的控制量平均求解时间上也得到了验证. 因此本文CFDL-MFAC算法控制环节的时间复杂度要低于基于神经网络的模型预测控制.

2.2 CFDL-MFAC控制器灵敏参数的选择和遗传−粒子群优化整定

由第2.1.2节可知, $S_{N O, 2}$ 和 $D_{O, 5}$ 控制系统的CFDL-MFAC控制算法共有包含PJM初值在内的11个待调参数, 而且并没有MFAC算法应用于实际污水处理过程的先例, 参数的选取毫无经验可循, 这给控制器参数调整带来极大的不便. 为此, 引入CFDL-MFAC控制器参数的灵敏性分析算法和参数优化整定算法, 对常规CFDL-MFAC算法进行改进. 首先, 采用文献[30]所描述的MPSA算法, 分析CFDL-MFAC控制器参数的灵敏度, 根据每个参数的累计频率曲线分离程度(Degree of separation, DS)选出对控制器性能较为敏感的参数, 以降低参数优化的复杂度; 然后, 采用GA-PSO算法对敏感参数做优化整定, 得到灵敏控制器参数的优化解; 有效改善了控制器参数的调整效率, 提升了控制器的性能.

2.2.1 CFDL-MFAC控制器参数灵敏度分析

为选出CFDL-MFAC控制器的灵敏参数, 以第2.1节CFDL-MFAC控制算法控制第1.2节 $S_{N O, 2}$ 和 $D_{O, 5}$ 递推双线性系统, 以满足控制器参数范围的 $N$ 组服从均匀分布的随机参数集驱动CFDL-MFAC控制器, 将控制时长为 $T$ 的 $S_{N O, 2}$ 和 $D_{O, 5}$ 的累计跟踪误差平方和作为损失函数指标, 在该损失函数的基础上, 做Monte Carlo实验.

损失函数的计算式为

_{}_{}^{} \overset{}{}_{}^{}_{}

(27)

式中, $p a r a = (λ, μ, η, ρ, b_{1}, b_{2}, α, ϕ_{i i} (0), ϕ_{i j} (0))$ 表示控制器的待调参数, $\bar{y} (k)$ 表示在控制器(23)作用下 $S_{N O, 2}$ 和 $D_{O, 5}$ 递推双线性系统的输出, $y_{r} (k)$ 表示的设定值, $T$ 表示控制时长. 待调参数 $p a r a$ 中各控制器参数的约束范围依次为

\begin{aligned} _{}_{} \\ _{}_{}_{}_{}_{} \end{aligned}

(28)

依据式(27)计算得到的 $N$ 个损失函数值的数值大小, 将其对应的参数集分为可接受参数集和不可接受参数集. 考察每个参数在两个参数集中的分布情况, 若某个参数在两个参数集中分布情况相同, 则该参数为不灵敏参数, 反之, 该参数是灵敏的. 具体过程如下.

按式(29)计算某参数的累计频率 $c f_{j}$

_{}_{}^{} \frac{_{}}{}

(29)

其中, $v_{i}$ 表示第 $i$ 个参数值在其对应的参数集中出现的频数.

绘制每个参数的累计频率曲线, 按下式计算其分离程度 $D S$ :

\frac{_{}^{} {\overset{}{}}_{}_{}^{}}{_{}^{} {_{} {\overset{}{}}_{}}^{}}

(30)

式中, $a {\overset{˘}{f}}_{i}$ , $r f_{i}$ 分别代表可接受和不可接受参数集对应参数的累计频率值, $a {\bar{f}}_{i}$ 表示可接受参数集相应参数累计频率值的均值. 对某个参数而言, DS越大表示该参数越不灵敏, DS越小, 则该参数越灵敏, 由此可选出CFDL-MFAC控制器的灵敏参数.

2.2.2 CFDL-MFAC控制器灵敏参数的GA-PSO优化

以式(31)为适应度函数, 对控制器的灵敏参数进行GA-PSO优化整定

\begin{aligned} _{}^{} {\overset{}{}_{}}^{} \overset{}{}_{} \\ _{}_{}^{}_{}_{}_{} \\ _{}_{}_{}_{} \\ _{}_{}^{}_{}_{} \\ _{}_{}_{} \\ _{}^{}_{}_{} \end{aligned}

(31)

式中, $p_{i} (i = 1, \dots, 7)$ 分别对应 $p a r a$ 中的前7个控制器参数, $g (\cdot)$ 为分段函数, 满足 ${g (p_{i}) |}_{p_{i} > 0} = 0$ , ${g (p_{i}) |}_{p_{i} \leq 0} = 1$ ; $K_{1}$ , $K_{2}$ 和 $K_{3}$ 为控制器参数的约束惩罚系数, 其作用就是在任一参数不满足约束条件时, 在适应度函数值上增加一个较大的惩罚值, 使该参数对应的优化解不能进入最优解空间, 将式(27)和式(28)组合而成的灵敏控制器参数的约束优化问题转化为无约束优化问题, 使优化过程更加简便. 需要指出的是, 在优化之前需根据灵敏度分析结果, 将适应度函数(31)中不灵敏参数变量取为适当的固定值, 并去除其对应的约束惩罚项.

虽然常用智能优化算法粒子群和遗传算法均可应用于目标函数(31)的优化求解, 但由于单一粒子群优化(Particle swarm optimization, PSO)算法[31-32]和遗传算法(Genetic algorithm, GA)[33-34]易陷入早熟收敛和局部最优, 为提高控制器参数优化的收敛速度和优化精度, 将两种算法相结合, 采用收敛速度更快、优化精度更高的GA-PSO混合优化算法, 对敏感控制器参数做优化整定. GA-PSO算法首先使用GA对适应度函数进行遗传优化, 在达到设定精度或优化指定代数之后, 将得到的最优种群作为后续PSO算法的初始种群; 然后, 以PSO算法继续寻优, 直到满足终止条件, 此时得到的最优解为控制器参数最优值. 此算法继承了GA染色体变异的特点和PSO算法粒子具有记忆功能的特点, 可有效提升算法搜索能力和收敛速度. 单独GA和PSO算法的迭代过程可分别参考文献[33-34]和文献[31-32], GA-PSO混合优化算法的具体实现步骤如下.

算法1. GA-PSO混合优化算法

步骤 1. 初始化GA-PSO算法的参数.

步骤 2. 对灵敏控制器参数变量进行实数编码.

步骤 3. 根据式(31)选择染色体, 并对染色体进行交叉和变异操作.

步骤 4. 判断是否达到指定代数或优化精度, 若未达到, 返回步骤2; 若达到, 则将得到的最优种群作为PSO的初始化粒子种群.

步骤 5. 计算粒子群适应度值, 并作出评价.

步骤 6. 更新粒子的位置和速度.

步骤 7. 判断是否达到终止条件, 若未达到, 返回步骤5; 若达到, 则即得到最优控制器参数.

3. 数据实验分析

由于第1.2节中已经完成 $S_{N O, 2}$ 和 $D_{O, 5}$ 的递推双线性子空间建模实验, 这里只需在第1.2节所建立的 $S_{N O, 2}$ 和 $D_{O, 5}$ 递推双线性子空间模型的基础上, 完成 $S_{N O, 2}$ 和 $D_{O, 5}$ 控制过程的CFDL-MFAC控制器参数敏感性分析实验和敏感参数优化整定实验, 然后对参数优化整定后的CFDL-MFAC控制器做性能分析实验即可. 本文采用戴尔台式电脑XPS 8930进行实验, 其CPU为: i7-8700 (6核), 内存: 16 GB, 显卡: NVIDIA GTX1070 (独立8 GB), 操作系统为win10 64位家庭中文版, 实验时采用软件为MATLAB 2016b.

3.1 CFDL-MFAC控制器参数灵敏度分析和优化整定实验

基于污水处理数据, 按照第2.2.1节中所述的参数灵敏度分析方法, 做 $S_{N O, 2}$ 和 $D_{O, 5}$ 的CFDL-MFAC控制器参数灵敏度分析实验, Monte Carlo实验的次数 $N$ 取5000, 得到的参数累计频率分布情况如图5所示, 其中实线对应某参数在可接受参数集内的累计频率分布情况, 虚线对应不可接受参数集内的累计频率分布情况. 由图5可知, 控制器参数(8) $\sim$ (11)在两个参数集内的累计频率分布曲线差别较大, 其余控制器参数在两个参数集内的分布曲线形状较为接近, 由此可初步判断控制器参数(8) $\sim$ (11)比较灵敏.

图 5 控制器参数累计频率分布曲线

Fig. 5 Cumulative frequency distribution curve of controller parameters

用式(30)计算图5中各参数累计频率曲线的分离程度DS, 依据计算结果分析各参数的灵敏情况, 并为不灵敏的参数取值. 灵敏度分析结果如表2所示, 由表中第4列可选出对应DS值较小的参数作为CFDL-MFAC控制器的灵敏参数, 即 $ϕ_{11} (0)$ , $ϕ_{12} (0)$ , $ϕ_{21} (0)$ , $ϕ_{22} (0)$ , 同时对不灵敏的参数给出固定值, 如表2第5列. 以式(31)为适应度函数, 分别按第2.2.2节的GA-PSO混合优化算法、文献[33]中的GA和文献[31]中的惯性权重线性递减PSO算法对灵敏参数做优化整定实验. GA-PSO算法的参数设置如表3所示, 终止条件为连续180次迭代适应度值误差不超过0.002, 或达到最大迭代次数300, 则终止迭代. GA和惯性权重线性递减PSO算法的最大迭代次数均为300, 其他参数值与表3中对应参数相同.

表 2 CFDL-MFAC控制器参数灵敏度分析结果

Table 2 Sensitivity analysis results of CFDL-MFAC controller parameters

序号	参数	含义	DS	不灵敏参数固定值
1	$λ$	输入权重因子	0.9981	0.5
2	$μ$	PJM权重因子	0.9985	0.6
3	$η$	PJM步长因子	0.9996	0.5
4	$ρ$	输入步长因子	0.9378	1.0
5	$α$	PJM重置限定参数	0.998	1.5
6	$b_{1}$	PJM重置限定参数	0.9991	0.55
7	$b_{2}$	PJM重置限定参数	0.9995	0.8
8	$ϕ_{11} (0)$	PJM初值	0.5551	−
9	$ϕ_{12} (0)$	PJM初值	−0.6339	−
10	$ϕ_{21} (0)$	PJM初值	0.7289	−
11	$ϕ_{22} (0)$	PJM初值	0.5779	−

表 3 GA-PSO算法参数

Table 3 GA-PSO algorithm parameters

前期 GA 参数		后期 PSO 参数
种群规模	$M = 40$	种群规模	$M = 40$
指定遗传代数	60	最大迭代次数	240
基因重组概率	0.7	加速系数	$c_{1} = 4, c_{2} = 2$
变异概率	0.25	最大惯性权重系数	$W_{max} = 0.9$
子代选择系数	T = 1, $α_{T} = 0.5$	最小惯性权重系数	$W_{min} = 0.1$

图6给出了混合GA-PSO优化过程、GA优化过程和PSO优化过程的收敛曲线. 由图6可知, 所提CFDL-MFAC控制器灵敏参数的GA-PSO优化终止于第281次迭代, 最优适应度值为 $f i t_{o p t} =$ $2.100$ . 而GA优化至最大迭代次数300时, 最优适应度值为 $f i t_{o p t} = 2.181$ , 且算法仍未收敛. PSO优化终止于第300次迭代, 最优适应度值为 $f i t_{o p t} =$ $2.114$ . 由此可见, 在CFDL-MFAC控制器灵敏参数的优化过程中, GA-PSO算法无论在收敛速度还是优化精度上都要优于单独的GA算法和PSO算法, 因此, 最终CFDL-MFAC的敏感控制器参数的最优解为GA-PSO算法的优化结果: $ϕ_{11} (0) = 0.0835$ , $ϕ_{12} (0) = - 0.3469,$ $ϕ_{21} (0) = - 0.6992,$ $ϕ_{22} (0) = - 0.0467$ .

图 6 CFDL-MFAC控制器灵敏参数优化收敛过程

Fig. 6 The optimization convergence process of sensitive parameter based on CFDL-MFAC controller

3.2 CFDL-MFAC控制器性能分析实验

采用控制器参数优化整定后的CFDL-MFAC控制器对污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 的设定值做跟踪控制实验, 实验采样间隔为0.25 h, 实验控制时长为62.5 h. 实验中, 以跳变的阶跃信号作为 $S_{N O, 2}$ 和 $D_{O, 5}$ 的跟踪目标 $y_{r} (k)$ , 以幅值为0.02的随机白噪声序列描述污水处理过程存在的噪声. 同时, 为考察改进CFDL-MFAC控制器针对不同类型扰动的鲁棒性, 在系统输出端分别施加幅值为0.1、周期为30 h的方波扰动和幅值为0.15、周期为60 h的正弦扰动, 做CFDL-MFAC控制器的抗扰性分析实验. 实验时, 方波扰动和正弦扰动在初始时刻即施加于系统 $S_{N O, 2}$ 和 $D_{O, 5}$ 的输出值. 跟踪目标 $y_{r} (k)$ 的具体设置为

_{} \begin{aligned} _{} \begin{array}{l}  \end{array} \\ _{} \begin{array}{l}  \end{array} \end{aligned}

为进一步展现本文改进CFDL-MFAC算法对 $S_{N O, 2}$ 和 $D_{O, 5}$ 的控制性能, 在相同噪声和扰动下, 与文献[7]中基于数据驱动的预测控制(Recursive bilinear subspace identification based model predictive control, RBL-MPC)算法进行比较. 两种算法在方波扰动和正弦扰动下的抗扰性控制结果对比如图7和图8所示, 控制性能指标对比如表4所示. 由于实验输出与扰动变量数值不在同一数量级, 为清晰地展现扰动变量对输出值的影响, 图中扰动变量曲线均为加上输出设定值之后的数值.

图 7 方波扰动抗扰性实验

Fig. 7 Immunity test with square wave disturbance

图 8 正弦扰动抗扰性实验

Fig. 8 Immunity test with sine disturbance

表 4 不同算法控制性能对比

Table 4 Comparison of control performance based on different algorithms

算法		RBL-MPC	CFDL-MFAC
控制量平均求解时间 (s)		0.0987	0.00003528
方波扰动 RMSE	$S_{N O, 2}$	0.0865	0.0297
方波扰动 RMSE	$D_{O, 5}$	0.0158	0.0107
正弦扰动 RMSE	$S_{N O, 2}$	0.0859	0.0249
正弦扰动 RMSE	$D_{O, 5}$	0.0099	0.0089

图7描述了在随机噪声和方波扰动的影响下, 两种算法对 $S_{N O, 2}$ 和 $D_{O, 5}$ 设定值的跟踪效果. 由图7(a)可知, 本文提出的改进CFDL-MFAC算法虽然在方波扰动的切换处和跟踪目标的跳变处有较大波动, 但能在控制输入的实际可操作范围内, 自适应地调整控制量的输入值(如图7(b)所示), 并快速收敛于设定目标, 完成对 $S_{N O, 2}$ 和 $D_{O, 5}$ 的稳定跟踪, 可有效抑制方波扰动信号的影响, 而RBL-MPC算法则无法在方波扰动下平稳地跟踪 $S_{N O, 2}$ 和 $D_{O, 5}$ 的设定目标, 抑制方波扰动的能力较差. 此外, 由两种算法在方波扰动下的跟踪误差绝对值曲线图7(d)可知, 除部分设定值的跳变时刻和扰动信号幅值的切换处之外, RBL-MPC方法的跟踪误差绝对值均大于本文CFDL-MFAC算法, 由此可进一步表明本文CFDL-MFAC算法在 $S_{N O, 2}$ 和 $D_{O, 5}$ 跟踪控制中的有效性和优越性.

由图8正弦扰动信号下的抗扰性实验结果可知. 在对 $S_{N O, 2}$ 和 $D_{O, 5}$ 的控制过程中, 改进CFDL-MFAC算法能有效抑制外部正弦扰动的影响, 实现对 $S_{N O, 2}$ 和 $D_{O, 5}$ 的精确跟踪, 且控制输入满足实际可操作范围(如图8(b)所示), 而RBL-MPC算法在正弦扰动下, 虽然能实现对 $D_{O, 5}$ 设定值的稳定跟踪, 但对 $S_{N O, 2}$ 的跟踪过程一直存在较大跟踪误差, 无法实现对 $S_{N O, 2}$ 和 $D_{O, 5}$ 的同时跟踪. 两种算法在正弦扰动下的跟踪误差绝对值曲线如图8(d)所示, 本文CFDL-MFAC算法仅在个别设定值的切换处有较大波动, 其余时刻的误差绝对值均趋于零, 而RBL-MPC方法的跟踪误差则一直存在较大波动, 且除少数设定值的切换时刻外, 其跟踪误差绝对值均大于本文CFDL-MFAC算法的跟踪误差绝对值. 由此可见, 本文CFDL-MFAC算法在污水处理 $S_{N O, 2}$ 和 $D_{O, 5}$ 控制中的跟踪效果要优于RBL-MPC算法.

由上述 $S_{N O, 2}$ 和 $D_{O, 5}$ 的设定值跟踪和抗扰性控制实验的对比分析可知, 与RBL-MPC方法相比, 本文提出的改进CFDL-MFAC算法不仅控制性能更加精确(见表4第3行和第4行), 控制量求解速度更加迅速(RBL-MPC方法控制量求解时间数量级为 $10^{- 2} s$ , 而CFDL-MFAC算法控制量求解时间数量级为 $10^{- 5} s$ , 见表4第2行), 且具有较强的抑制外部干扰的能力. 两种算法具体的控制性能指标可由表4所知.

4. 结束语

针对污水处理过程 $S_{N O, 2}$ 和 $D_{O, 5}$ 难以机理建模和直接数据驱动控制的问题, 本文提出一种改进CFDL-MFAC污水处理过程控制方法, 即在传统CFDL-MFAC算法的基础上, 引进控制器参数灵敏度分析和GA-PSO优化整定机制, 采用MPSA方法分析控制器参数的灵敏度, 并使用GA-PSO混合优化算法对敏感控制器参数做优化整定, 以改进传统CFDL-MFAC算法控制器参数的调整效率, 提高控制器的控制性能. 同时, 为便于进行灵敏度分析的Monte Carlo实验和灵敏参数的GA-PSO优化整定实验, 采用一种RBLSI建模方法, 建立 $S_{N O, 2}$ 和 $D_{O, 5}$ 的递推双线性子空间模型, 为参数优化整定实验和控制器性能分析实验提供较为精确的被控对象模型. 基于污水处理数据的控制器性能分析实验表明, 本文基于递推双线性子空间数据驱动建模的改进CFDL-MFAC算法不仅可以精确跟踪 $S_{N O, 2}$ 和 $D_{O, 5}$ 的目标值, 而且能够有效抑制外部扰动的影响, 具有优良的控制性能和良好的鲁棒性.

关键字：优秀论文

上一篇：基于移动边缘计算的GI/GI/1排队建模与调度算法
下一篇：基于低秩约束的熵加权多视角模糊聚类算法

栏目分类

热门排行

推荐信息

期刊知识