[发明专利]循环流化床生活垃圾焚烧锅炉汽包水位的预测系统及方法

权利要求书

1.一种循环流化床生活垃圾焚烧锅炉汽包水位的预测系统，该系统与循环流化床锅炉的集散控制系统相连，包括数据通讯接口和上位机，所述上位机包括：

第一信号采集模块；利用该模型采集CFB生活垃圾焚烧锅炉在正常运行时的运行工况状态参数和操作变量，并组成ANIFS汽包水位补偿模型输入变量的训练样本矩阵X(m×n)，m表示样本个数，n表示变量的个数；

数据预处理模块；对X(m×n)进行粗大误差处理和随机误差处理，剔除训练样本中的野值，排除异常工况，将训练样本输入变量经归一化处理后映射到[0，1]区间内，得到标准化后的训练样本X^*(m×n)；

专家知识库模块；将水位计测量出来的实际汽包水位值l(m×1)与机理模型计算得到的水位值l_d(m×1)相减，将两者之间的差值作为ANFIS汽包水位补偿模型训练样本的输出向量l^*(m×1)，X^*和l^*共同组成ANFIS模型的训练样本，并进行保存；

混合建模模块；该模块主要包括1)机理建模和2)ANFIS建模两个部分；1)机理建模的步骤如下：

1.1)对水冷壁中的工质进行质量守恒、动量守恒和能量守恒分析，可得：

其中，z为坐标轴，正向为工质的上升方向；ρ_s、ρ_l分别表示研究对象内饱和水蒸气的密度和饱和水的密度，密度的单位是kg/m³；α_s表示截面容积含汽率；v表示工质的流速，单位是m/s；τ_wm是壁面剪切应力，单位是Pa/m³；p_r是上升管中的压力，单位Pa；表示热流密度，单位是W/m³；g是重力加速度，单位是m/s²；h_l为饱和水的比焓；

饱和温度和工质压力之间存在着以下关系：

T＝T_b(1+a₁ln(P/P_b)+a₂ln²(P/P_b)+a₃ln³(P/P_b)+a₄ln⁴(P/P_b)) (5)

其中P是当前位置压力，T是当前位置的温度，T_b参考点的温度，P_b表示参考点的压力，a₁、a₂、a₃、a₄是常系数；

壁面的剪切应力可以通过下式求取：

f＝0.0014+0.125Re^-0.32(7)

ρ_m＝α_sρ_s+(1-α_s)ρ_l (8)

式中Re表示雷诺数，f为摩擦系数，d_t为水冷壁上升管的管径，ρ_m为汽水混合物的密度；

水冷壁的热流密度可以通过下式求取：

式中，A_t表示截面积，U_wm表示传热系数，T_w表示管壁的温度；

对(2)～(4)式进行展开，并将(5)～(9)式代入，可得到如下方程组：

其中，A₁₁＝(ρ_s-ρ_l)v，A₁₂＝α_sv，A₁₃＝α_sρ_s+(1-α_s)ρ_l，A₁₄＝0，A₁₅＝0；

A₂₅＝0；A₃₁＝(ρ_s-ρ_l)v²，A₃₂＝α_sv²，

A₃₃＝2[α_sρ_s+(1-α_s)ρ_l]v，A₃₄＝0，A₃₅＝1；A₄₁＝0，A₄₂＝0，A₄₃＝0，A₄₄＝1，

A₅₁＝0，A₅₂＝1，A₅₃＝0，

y₁＝α_s，y₂＝ρ_s，y₃＝v，y₄＝T，y₅＝P；B₁＝0，

B₄＝0，B₅＝0；式中，λ_V表示汽化潜热，表示饱和水的定压比热，表示饱和蒸汽的定压比热，T₀表示指定状态下工质的温度；最终可以算得：

1.2)通过对汽包-下降管-水冷壁中的所有工质进行质量守恒分析和能量守恒分析，可得：

其中，e₁＝ρ_l-ρ_s；

式中，p表示汽包压力，单位Pa；Q表示整个蒸发区的单位时间内的吸热量，单位是W；V_lt表示汽包-下降管-水冷壁中所有饱和水的体积，单位是m³；V_st表示回路中所有水蒸气的体积，单位是m³；q_f、q_s分别表示给水流量和主汽流量和下降管中循环工质的流量，流量的单位是kg/s；h_f、h_l、h_s分别表示给水比焓、饱和水的比焓和饱和水蒸气的比焓，单位是kJ/kg；m_t表示水冷壁金属和汽包金属的总质量，单位是kg；C_p表示金属的比热容；t_s表示饱和温度；

1.3)对汽包中液面下的工质进行质能守恒分析，忽略液面上饱和蒸汽的凝结，可得：

其中，

式中，α_r表示上升管出口的截面质量含汽率；表示汽包内液面下工质的质量，单位是kg；q_r、q_cd、q_dc分别表示上升管中的流量、饱和蒸汽逸出液面的流量、饱和蒸汽凝结的流量和下降管中的流量，单位是kg/s；t_s表示饱和蒸汽的温度，单位是℃；m_d表示汽包金属的质量，单位kg；C_p表示水冷壁和汽包中金属的比热容，单位是kJ/(kg·℃)；A_d是汽包液位的截面面积，单位是m²；是饱和蒸汽逸出的速度，单位是m/s；σ是液面的表面张力，单位是N/m；V_l分别表示液面下饱和水蒸气的体积和饱和水的体积，单位是m³；

联立(12)～(15)式并进行展开，可以得到如下方程组：

其中，M₁₂＝0，M₁₃＝0，M₁₄＝0；

M₂₂＝ρ_l，M₂₃＝ρ_s，M₂₄＝0；M₃₂＝0，M₃₃＝ρ_s，M₃₄＝0；

M₄₂＝0，M₄₃＝0，M₄₄＝ρ_l-ρ_s；

N₁＝Q-q_f(h_l-h_f)-q_s(h_s-h_l)，

N₄＝q_f-q_s；o₁＝p，o₂＝V_l，o₃＝V_sd，o₄＝V_lt；

最终可得：

采用Runge-Kutta数值解法求解微分方程(11)和(16)；对仿真时间进行离散化处理，将[0，T₀]区间内的仿真时间，离散成n个时间间隔为ΔT的区间，离散后的区间表示为[0,ΔT]，[ΔT，2ΔT]，……，[(n-1)ΔT，nΔT]；ΔT为0.01s～0.1s；

1.4)汽包水位是由汽包中的饱和水和饱和水蒸气的体积共同决定的，忽略汽包体积和汽包水位之间的非线性，可得汽包水位的表达式：

式中，l_d表示机理模型计算得到的汽包水位，单位mm；l₀是汽包水位的基准水位线；

2)ANFIS建模的步骤如下：

2.1)初始化粒子群；以聚类半径r_α作为粒子，15个粒子作为一个种群，每个粒子随机赋予[0.20.9]区间内的随机值，其中第i个粒子的位置的向量标示为ri，i＝1，2，…，15；

2.2)以r_i为聚类半径，进行减法聚类分析，将每个数据点作为可能的聚类中心，并根据各个数据点周围的数据点密度来计算该点作为聚类中心的可能性；

每个数据点X_i作为聚类中心的可能性P_i由式(18)来定义：

式中m表示n维输入空间中全部的数据点数，X_i＝[X_i1,X_i2,...,X_in]、X_j＝[X_j1,X_j2,...,X_jn]是具体的数据点，r_i是一个正数，定义了该点的邻域半径，||·||符号表示欧式距离；被选为聚类中心的点具有最高的数据点密度，同时该该数据点周围的点被排除作为聚类中心的可能性；第一个聚类中心为X_C1，数据点密度为P_c1；选出第一个聚类中心后，继续采用类似的方法确定下一个聚类中心，但需消除已有聚类中心的影响，修改密度指标的山峰函数如下：

其中，r_β＝1.5r_i；循环重复上述过程直到所有剩余数据点作为聚类中心的可能性低于 某一阈值δ，即P_ck/P_c1＜δ；

2.3)ANFIS模型训练；根据减法聚类算法得到的聚类中心，按照ANFIS模型结构训练汽包水位补偿模型；对于模糊神经网络模型的所有参数，采用混合最小二乘法的梯度下降算法进行学习；

2.4)计算适应度值；利用训练得到的预测模型计算汽包水位补偿值l^*，并与实际补偿值l_e进行比较，并以误差平方和作为粒子的适应度值MSE，适应度计算公式如下：

2.5)更新极值；以适应度值为评价指标，比较当代粒子与上一代粒子之间的适应度值大小，如果当前粒子的适应度值优于上一代，则将当前粒子的位置设置为个体极值，否则个体极值保持不变；同时获取当代所有粒子适应度值最优的粒子，并与上一代最优粒子进行比较，如果当代最优粒子的适应度值优于上一代最优粒子的适应度值，则将当代粒子的最优适应度值设置为全局最优值，否则全局最优值保持不变；

2.6)更新粒子；根据最新的个体极值和全局极值，按照(29)式和(30)式更新粒子的速度v_id(t)和位置x_id(t)；

v_id(t+1)＝ωv_id(t)+c₁r₁(p_id-x_id(t))+c₂r₂(p_gd-x_id(t))(29)

x_id(t+1)＝x_id(t)+v_id(t+1)(30)

t是粒子群优化算法的寻优代数；更进一步，为了改善基本粒子群算法容易陷入局部极值和收敛速度慢的缺陷，在PSO算法的基础上引进了动态加速常数c₁、c₂和惯性权重ω：

其中，T_max为最大寻优代数，ω_max为最大惯性权重，ω_min为最小惯性权重，R₁、R₂、R₃、R₄为常数；

2.7)算法停止条件算法判定；判断是否达到最大迭代次数或者到达预测精度的要求，如果没有达到则返回步骤2.2)，利用更新的聚类半径继续搜索，否则退出搜索；

2.8)利用最终寻优得到的聚类半径，对样本进行聚类分析和ANFIS模型训练，得到 达到训练终止条件的ANFIS模型，即汽包水位的补偿模型；

第二信号采集模块；用于从数据库中选择需要预测汽包水位的运行工况，或者实时地采集当前锅炉的运行工况；

预测模块；该模块用于对指定的样本进行汽包水位的预测，或者对当前锅炉运行工况下的汽包水位进行实时预测；

结果显示模块；显示汽包水位的预测结果，或者对汽包水位的预测结果进行统计分析。

2.一种循环流化床生活垃圾焚烧入炉汽包水位的预测方法，其特征在于，该方法包括以下步骤：

1)选择混合模型的输入变量并采集训练样本；混合模型包括机理模型和ANFIS补偿模型；输入变量包括锅炉的设计参数、工质的物性参数和运行参数；设计参数包括水冷壁的当量管径、水冷壁的有效吸热长度、水冷壁的容积、循环水的流量、汽包的容积、汽包的平均截面面积、下降管的容积、汽包-下降管-水冷壁三者总共的金属质量；工质的物性参数包括压力、饱和温度、比焓、比容、汽化潜热、密度；运行参数包括给水流量、主蒸汽流量、垃圾给料量、辅助燃煤量、一次风量、二次风量以及炉膛中所有的温度测点；采集CFB生活垃圾焚烧锅炉在正常运行时的运行工况状态参数和操作变量，并组成混合模型输入变量的训练样本矩阵X(m×n)，m表示样本个数，n表示变量的个数；

2)训练样本预处理；对训练样本X(m×n)进行粗大误差处理和随机误差处理，剔除训练样本中的野值，排除异常工况，将训练样本输入变量经归一化处理后映射到[0，1]区间内，得到归一化后的训练样本X^*(m×n)；所述的异常工况包括锅炉停炉运行工况、炉膛压火运行状况和给料机堵塞工况；

3)利用机理建模的方法构建CFB生活垃圾焚烧锅炉汽包水位的机理模型；

3.1)对水冷壁中的工质进行质量守恒、动量守恒和能量守恒分析，可得：

其中，z为坐标轴，正向为工质的上升方向；ρ_s、ρ_l分别表示研究对象内饱和水蒸气的密度和饱和水的密度，密度的单位是kg/m³；α_s表示截面容积含汽率；v表示工质的流速， 单位是m/s；τ_wm是壁面剪切应力，单位是Pa/m³；p_r是上升管中的压力，单位Pa；表示热流密度，单位是W/m³；g是重力加速度，单位是m/s²；h_l为饱和水的比焓；

饱和温度和工质压力之间存在着以下关系：

T＝T_b(1+a₁ln(P/P_b)+a₂ln²(P/P_b)+a₃ln³(P/P_b)+a₄ln⁴(P/P_b)) (5)

其中P是当前位置压力，T是当前位置的温度，T_b参考点的温度，P_b表示参考点的压力，a₁、a₂、a₃、a₄是常系数；

壁面的剪切应力可以通过下式求取：

f＝0.0014+0.125Re^-0.32 (7)

ρ_m＝α_sρ_s+(1-α_s)ρ_l (8)

式中Re表示雷诺数，f为摩擦系数，d_t为水冷壁上升管的管径，ρ_m为汽水混合物的密度；

水冷壁的热流密度可以通过下式求取：

式中，A_t表示截面积，U_wm表示传热系数，T_w表示管壁的温度；

对(2)～(4)式进行展开，并将(5)～(9)式代入，可得到如下方程组：

其中，A₁₁＝(ρ_s-ρ_l)v，A₁₂＝α_sv，A₁₃＝α_sρ_s+(1-α_s)ρ_l，A₁₄＝0，A₁₅＝0；

A₂₅＝0；A₃₁＝(ρ_s-ρ_l)v²，A₃₂＝α_sv²，

A₃₃＝2[α_sρ_s+(1-α_s)ρ_l]v，A₃₄＝0，A₃₅＝1；A₄₁＝0，A₄₂＝0，A₄₃＝0，A₄₄＝1，

A₅₁＝0，A₅₂＝1，A₅₃＝0，

y₁＝α_s，y₂＝ρ_s，y₃＝v，y₄＝T，y₅＝P；B₁＝0，

B₄＝0，B₅＝0；式中，λ_V表示汽化潜热，表示饱和水的定压比热，表示饱和蒸汽的定压比热，T₀表示指定状态下工质的温度；最终可以算得：

3.2)通过对汽包-下降管-水冷壁中的所有工质进行质量守恒分析和能量守恒分析，可得：

其中，e₁＝ρ_l-ρ_s；

式中，p表示汽包压力，单位Pa；Q表示整个蒸发区的单位时间内的吸热量，单位是W；V_lt表示汽包-下降管-水冷壁中所有饱和水的体积，单位是m³；V_st表示回路中所有水蒸气的体积，单位是m³；q_f、q_s分别表示给水流量和主汽流量和下降管中循环工质的流量，流量的单位是kg/s；h_f、h_l、h_s分别表示给水比焓、饱和水的比焓和饱和水蒸气的比焓，单位是kJ/kg；m_t表示水冷壁金属和汽包金属的总质量，单位是kg；C_p表示金属的比热容；t_s表示饱和温度；

3.3)对汽包中液面下的工质进行质能守恒分析，忽略液面上饱和蒸汽的凝结，可得：

其中，

式中，α_r表示上升管出口的截面质量含汽率；表示汽包内液面下工质的质量，单位是kg；q_r、q_cd、q_dc分别表示上升管中的流量、饱和蒸汽逸出液面的流量、饱和蒸汽凝结的流量和下降管中的流量，单位是kg/s；t_s表示饱和蒸汽的温度，单位是℃；m_d表示汽包金属的质量，单位kg；C_p表示水冷壁和汽包中金属的比热容，单位是kJ/(kg·℃)；A_d是汽包液位的截面面积，单位是m²；是饱和蒸汽逸出的速度，单位是m/s；σ是液面的表面张力，单位是N/m；V_l分别表示液面下饱和水蒸气的体积和饱和水的体积，单位是m³；

联立(12)～(15)式并进行展开，可以得到如下方程组：

其中，M₁₂＝0，M₁₃＝0，M₁₄＝0；

M₂₂＝ρ_l，M₂₃＝ρ_s，M₂₄＝0；M₃₂＝0，M₃₃＝ρ_s，M₃₄＝0；

M₄₂＝0，M₄₃＝0，M₄₄＝ρ_l-ρ_s；

N₁＝Q-q_f(h_l-h_f)-q_s(h_s-h_l)，

N₄＝q_f-q_s；o₁＝p，o₂＝V_l，o₃＝V_sd，o₄＝V_lt；

最终可得：

采用Runge-Kutta数值解法求解微分方程(11)和(16)；对仿真时间进行离散化处理，将[0，T₀]区间内的仿真时间，离散成n个时间间隔为ΔT的区间，离散后的区间表示为[0,ΔT]，[ΔT，2ΔT]，……，[(n-1)ΔT，nΔT]；ΔT为0.01s～0.1s；

3.4)汽包水位是由汽包中的饱和水和饱和水蒸气的体积共同决定的，忽略汽包体积和汽包水位之间的非线性，可得汽包水位的表达式：

式中，l_d表示机理模型计算得到的汽包水位，单位mm；l₀是汽包水位的基准水位线；

4)建立ANFIS补偿模型，步骤如下：

4.1)初始化粒子群；以聚类半径r_α作为粒子，15个粒子作为一个种群，每个粒子随机赋予[0.20.9]区间内的随机值，其中第i个粒子的位置的向量标示为ri，i＝1，2，…，15；

4.2)以r_i为聚类半径，进行减法聚类分析，将每个数据点作为可能的聚类中心，并根据各个数据点周围的数据点密度来计算该点作为聚类中心的可能性；

每个数据点X_i作为聚类中心的可能性P_i由式(18)来定义：

式中m表示n维输入空间中全部的数据点数，X_i＝[X_i1,X_i2,...,X_in]、X_j＝[X_j1,X_j2,...,X_jn]是具体的数据点，r_i是一个正数，定义了该点的邻域半径，||·||符号表示欧式距离；被选为聚类中心的点具有最高的数据点密度，同时该该数据点周围的点被排除作为聚类中心的可能性；第一个聚类中心为X_C1，数据点密度为P_c1；选出第一个聚类中心后，继续采用类似的方法确定下一个聚类中心，但需消除已有聚类中心的影响，修改密度指标的山峰函数如下：

其中，r_β＝1.5r_i；循环重复上述过程直到所有剩余数据点作为聚类中心的可能性低于某一阈值δ，即P_ck/P_c1＜δ；

4.3)ANFIS模型训练；根据减法聚类算法得到的聚类中心，按照ANFIS模型结构训练汽包水位补偿模型；对于模糊神经网络模型的所有参数，采用混合最小二乘法的梯度下降算法进行学习；

4.4)计算适应度值；利用训练得到的预测模型计算汽包水位补偿值l^*，并与实际补偿值l_e进行比较，并以误差平方和作为粒子的适应度值MSE，适应度计算公式如下：

4.5)更新极值；以适应度值为评价指标，比较当代粒子与上一代粒子之间的适应度值大小，如果当前粒子的适应度值优于上一代，则将当前粒子的位置设置为个体极值，否则个体极值保持不变；同时获取当代所有粒子适应度值最优的粒子，并与上一代最优粒子进行比较，如果当代最优粒子的适应度值优于上一代最优粒子的适应度值，则将当代粒子的最优适应度值设置为全局最优值，否则全局最优值保持不变；

4.6)更新粒子；根据最新的个体极值和全局极值，按照(29)式和(30)式更新粒子的速度v_id(t)和位置x_id(t)；

v_id(t+1)＝ωv_id(t)+c₁r₁(p_id-x_id(t))+c₂r₂(p_gd-x_id(t))(29)

x_id(t+1)＝x_id(t)+v_id(t+1)(30)

t是粒子群优化算法的寻优代数；更进一步，为了改善基本粒子群算法容易陷入局部极值和收敛速度慢的缺陷，在PSO算法的基础上引进了动态加速常数c₁、c₂和惯性权重ω：

其中，T_max为最大寻优代数，ω_max为最大惯性权重，ω_min为最小惯性权重，R₁、R₂、R₃、 R₄为常数；

4.7)算法停止条件算法判定；判断是否达到最大迭代次数或者到达预测精度的要求，如果没有达到则返回步骤2.b)，利用更新的聚类半径继续搜索，否则退出搜索；

4.8)利用最终寻优得到的聚类半径，对样本进行聚类分析和ANFIS模型训练，得到达到训练终止条件的ANFIS模型，即汽包水位的补偿模型；

5)汽包水位的预测；选择需要预测的工况，通过机理模型得到汽包水位的预测值l_d，通过ANFIS模型得到汽包水位的补偿值l^*，则混合模型的输出值即最终预测的汽包水位为