具体实施方式

为了使本申请的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本申请进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本申请，并不用于限定本申请。

本申请提供一种高炉鼓风机的控制方法。需要说明的是，本申请提供的高炉鼓风机的控制方法的应用于高炉鼓风机。

此外，本申请提供的高炉鼓风机的控制方法不限制其执行主体。可选地，本申请提供的高炉鼓风机的控制方法的执行主体的可以为一个控制器。

如图1所示，在本申请的一实施例中，高炉鼓风机的控制方法包括如下S100至S900：

S100，获取高炉鼓风机的输入数据和输出数据，依据输入数据和输出数据建立鼓风机系统模型。

S300，基于所述鼓风机系统模型，对鼓风机系统中的各个变量进行解耦，得到多个解耦后的变量。

S500，基于PID控制算法，生成鼓风机的PID控制逻辑表达式。

S700，构建鼓风机三维模型。

S900，将多个解耦后的变量依次导入鼓风机三维模型中，以实时驱动鼓风机三维模型中各个设备的设备动作。

具体地，本实施例将数字孪生技术与系统建模控制方法相结合，提出一种高炉鼓风机的定风量定风压的控制策略与控制方法。S500中，使用PID控制器生成鼓风机的PID控制逻辑表达式。

本实施例中，本申请涉及一种高炉鼓风机的控制方法，首先从数据出发建立鼓风机系统模型，进而针对鼓风机系统存在的多变量强耦合性，采用了解耦设计，消除各变量之间的耦合作用，把多变量的转化为独立的单变量系统，以实现更高质量的控制性能。接下来基于PID控制算法，生成鼓风机的PID控制逻辑表达式，实现对高炉鼓风机的准确控制。最后根据高炉鼓风机实际设备的外观与内部结构，构建鼓风机三维模型，并将解耦后的变量导入鼓风机三维模型中，实现高炉鼓风机控制流程和模型结构两者同步动作，从而实现高炉鼓风机设备的自动化运作，以及实时优化设备当前的运行状态。

在本申请的一实施例中，所述S100包括如下S111至S190：

S111，设置采样条件，所述采样条件包括采样周期长度Ts，采样周期序号k，采样周期总数N和采样起点Ts。所述采样长度为每一个采样周期的时间长度。

具体地，例如，如果Ts取1秒，N取1000，那采样条件就是在1000秒内采样，每秒采样1次。

S112，获取开度数据和出口数据。所述开度数据包括放空阀开度、防喘阀开度、静叶开度和进气阀门开度中的一种或多种。所述出口数据包括鼓风机出口压力和鼓风机出口流量。

具体地，所述输入数据即开度数据，所述输出数据即出口数据。

在实际现场，放空阀开度、防喘阀开度、静叶开度、进气阀门开度都可能成为控制设备的运行状态的控制量，所以这里的四个变量采用一种或多种的“或”的关系。

而鼓风机出口压力和鼓风机出口流量，为高炉鼓风机系统当前状态的关键变量，缺一不可，因此是“和”的关系。本申请希望通过控制放空阀或防喘阀、静叶开度或进气阀门的开度，获得高炉鼓风机系统稳定的鼓风机出口压力与鼓风机出口流量。

S120，依据公式1-a、公式1-b和公式1-c构建鼓风机系统模型。

k×T_s≤t≤(k+1)×T_s 公式1-c

其中，为开度数据。y(t)为出口数据。x(t)为状态向量。为状态向量的一阶导数。w(t)为过程噪声。v(t)为测量噪声。A为第一矩阵。B为第二矩阵。C为第三矩阵。D为第四矩阵。k为采样周期序号。Ts为采样周期长度。所述采样周期长度为每一个采样周期的时间长度。

具体地，第一矩阵A可以为m行n列矩阵。第二矩阵B可以为m行m列矩阵。第三矩阵C可以为m行n列矩阵。第四矩阵D可以为单行m列矩阵。t的单位为秒。本申请中，为了不影响后续计算，过程噪声和测量噪声默认在一个采样周期内保持不变。

S130，依据公式2-a、公式2-b和公式2-c定义第一展开矩阵、第二展开矩阵和第三展开矩阵。

Γ_i＝[C^T (CA)^T ... (CAⁱ)^T]^T 公式2-a

其中，Γ_i为第一展开矩阵。H_i为第二展开矩阵。G_i为第三展开矩阵。A为第一矩阵。B为第二矩阵。C为第三矩阵。D为第四矩阵。i为预设阶数。T为矩阵的转置符号。

具体地，第一展开矩阵可以为增广观测矩阵。第二展开矩阵和第三展开矩阵可以为Toeplitz矩阵。

S140，设置预设阶数i，对公式1-b求i阶导数，生成出口数据的i阶导数表达式，所述出口数据的i阶导数表达式如公式3所示。

其中，y_i(t)为出口数据y(t)的i阶导数。为开度数据。x(t)为状态向量。Γ_i为第一展开矩阵。H_i为第二展开矩阵。G_i为第三展开矩阵。w(t)为过程噪声。v(t)为测量噪声。i为出口数据y(t)的导数的阶数。

S150，分别构造出口数据，开度数据，过程噪声和测量噪声的栈向量的表达式。

注意，这里并没有对状态向量的栈向量的表达式的构建。

S160，将出口数据，开度数据，过程噪声和测量噪声的栈向量的表达式分别进行矩阵式扩展，得到出口数据，开度数据，过程噪声和测量噪声各自对应的矩阵扩展表达式。

S170，将出口数据，开度数据，过程噪声和测量噪声各自对应的矩阵扩展表达式代入至公式2，得到鼓风机系统模型的表达式。

S180，求取鼓风机系统模型的表达式中的第一矩阵，第二矩阵，第三矩阵和第四矩阵。

S190，将第一矩阵，第二矩阵，第三矩阵和第四矩阵代入至鼓风机系统模型的表达式。

具体地，本实施例就是整个鼓风机系统模型的建立过程，先建立总表达式(包括公式1-a、公式1-b和公式1-c)，然后求取总表达式中核心的四个参数，也就是第一矩阵，第二矩阵，第三矩阵和第四矩阵的具体表达式，再代入总表达式中，得到完整的总表达式。

在本申请的一实施例中，所述S150包括如下步骤：

S151，依据公式4，将出口数据转化为出口数据的栈向量的表达式。

其中，y_i(t)为出口数据的栈向量。y(t)为出口数据。为出口数据y(t)的一阶导数。y⁽ⁱ⁾为出口数据y(t)的i阶导数。T为矩阵的转置符号。i为预设阶数。

具体地，开度数据过程噪声w(t)和测量噪声v(t)也使用与公式4同样的方式转化为栈向量的表达式。由于噪声都是离散信号，并不能够准确的得到其导数，因此这样表示为了算法表达的方便，因为算法是基于采用数据的，这并不会影响算法的计算问题。

在本申请的一实施例中，S160包括如下步骤：

S161，基于公式4对公式3进行矩阵式扩展，并将采样数据扩展至N个，生成出口数据，开度数据，过程噪声和测量噪声各自对应的矩阵扩展表达式。出口数据，开度数据，过程噪声和测量噪声各自对应的矩阵扩展表达式如公式5所示。

其中，y_i(t_N)为t_N时刻下的出口数据的栈向量。为t_N时刻下的开度数据的栈向量。w(t_N)为t_N时刻下的过程噪声的栈向量。v(t_N)为t_N时刻下的测量噪声的栈向量。N为采样周期总数。

具体地，可以看出，公式5中的四个参数转化的过程是一致的。

在本申请的一实施例中，S170包括如下步骤如下S171至S172：

S171，将状态向量的采样数据拓展至N个，得到状态向量对应的扩展表达式。状态向量对应的扩展表达式为X_1,N＝[x(t₁)x(t₂)x(t₃)...x(t_N)]。

具体地，状态向量的对应表达式包括在不同时间周期内的状态向量的采样数据。不同时间周期内的状态向量的采样数据是在初始采样步骤里可以直接获取的，不需要进行栈向量的转化，更不需要进行矩阵扩展。

S172，将公式5中出口数据，开度数据，过程噪声和测量噪声各自对应的矩阵扩展表达式代入至公式2，得到鼓风机系统模型的表达式，鼓风机系统模型的表达式如公式6所示

具体地，X_1,N为状态向量对应的扩展表达式。N为采样周期总数。x(t₁)为t_N时刻下的状态向量。

本实施例一方面将公式5拓展到矩阵层面，另一方面，将公式5扩展至N个输入数据和输出数据的层面。这样生成的公式6可以呈现高炉鼓风机中所有设备的状态。

在本申请的一实施例中，S180包括如下步骤如下S181至S189：

S181，在鼓风机系统模型的表达式加入EIV噪声，并对鼓风机系统模型的表达式进行j阶GPMF变换，生成公式7。

其中，F^j(Y_i,1,N)为对Y_i,1,N的j阶变换。F^j(X_i,1,N)为对X_i,1,N的j阶变换，为对的j阶变换。F^j(W_i,1,N)为对W_i,1,N的j阶变换。F^j(V_i,1,N)为对V_i,1,N的j阶变换。L_i,1,N为EIV噪声。F^j(L_i,1,N)为对L_i,1,N的j阶变换。

具体地，直接从离散的输入输出信号中计算其时间导数是非常困难的。前述形成鼓风机系统模型的表达式，即公式6属于一个理想型的鼓风机系统模型，忽略了EIV噪声，因此在本实施例需要对鼓风机系统模型的表达式进行修正，补充EIV噪声的部分。

这里我们考虑n阶EIV连续系统的辨识问题，因此需要等式中的输入信号含有噪声。设EIV噪声数据为l(t)，那么我们就在开度数据中加入EIV噪声数据，形成复合EIV噪声后的开度数据。

其中，u(t)为复合EIV噪声后的开度数据u为开度控制样本，l(t)为EIV噪声数据，为开度数据。

假设l(t)，w(t)，v(t)具有同样的性质，并满足如公式7-b所示的相关性质。

其中，E为数学期望算子，δ_qj为Kronecker符号。

GPMF变换相对与其他线性动态运算方法具有以下优点：对采样周期不敏感，对噪声不敏感。

信号s(t)在[0,t]时间内的j阶GPMF变换可由下式7-c卷积得到

其中，g_j(τ)是j阶Poisson冲激函数。g_j*s(t)表示对s(t)的卷积。λ为滤波器的截断频率。λ满足λ∈(0,+∞)。β为Poisson滤波器常数。β满足β∈(0,+∞)。λ＝β。

下面对信号s(t)的导数向量进行定义，如公式(7-d)所示。

其中，s_i(t)为s(t)的导数向量。T为矩阵的转置符号。为s(t)的一阶导数，s^(i-1)为s(t)的i-1阶导数。

对F[s(t)]进行卷积计算可得：

因此将对s(t)进行导数计算转换成对g_j(t)的导数计算。而从g_j(t)的定义中很容易的计算其导数。

GPMF变换可以认为是一种线性滤波方法，将信号s(t)通过如下滤波器可以得到其j阶GPMF变换

为了方便数字计算机的处理，需要将上式转换成数字滤波器。比较简单的离散方法有零阶保持器(ZOH)和Tustin方法。通过选取较大的观测时间，GPMF方法中的初始值条件在一定时间后可以忽略。截断频率λ的选择,一般来讲，截断频率λ应该取为比辨识对象的频带宽度稍大。

为了更好地表达连续系统的子空间辨识方法，记F^j[s(t)]为s(t)的j阶变换。通过GPMF变换，将输入输出的导数问题转换成计算其GPMF变换。而GPMF变换的性质，可以很容易地使用低阶的变换进行线性变换得到高阶的变换。

因此，考虑EIV情况在鼓风机系统模型的表达式(公式6)中加入EIV噪声，并对鼓风机系统模型的表达式(公式6)进行j阶GPMF变换，就可以生成

S182，在公式7中掺入系统动态常数，生成公式8。

其中，f为系统动态常数。

S183，将公式8进行形式调整，得到公式9。

其中，I为单位矩阵。

具体地，将公式8中的移动至等式左边，并通过矩阵运算进行整合获得公式9。

S184，在公式9的两边各乘以第一展开矩阵的正交补矩阵，生成公式10。

其中，(Γ_i ^⊥)^T为第一展开矩阵的正交补矩阵。

具体地，公式9中的在公式9的两边各乘以第一展开矩阵的正交补矩阵的过程中，被消掉，得到公式10。

S185，设定辅助变量以及辅助变量的条件，辅助变量如公式11所示，辅助变量满足的条件如公式12所示。

其中，ζ为辅助变量。q为阶次。n为矩阵F^q(U_i,1,N-f)的行数。

其中，N为采样周期总数。f为系统动态常数。

具体地，为了消除等式右边的噪声项，最简洁的方法则是采用辅助变量。这在很多离散子空间辨识中都有应用。对于离散子空间辨识，都使用过去的输入输出向量作为辅助变量，并且在不会破坏一致性估计的情况下，过去的输入输出向量的长度越长，则越能提高其辨识精度。

为了能够得到一致性估计，使用辅助变量如选择过去时刻的联合输入输出的GPMF线性变换，其阶次的高低，根据所选择的导数的阶次变换。为了能够应用随机噪声项，假设N趋向于无穷，那么辅助变量需要满足公式12。

S186，将公式9乘以辅助变量，生成公式13。

S187，对公式13中的进行奇异值分解SVD，得到公式14。

其中，K₁为主左奇异矩阵，K₂为残差左奇异矩阵，Λ₁为奇异值矩阵，O₁为主右奇异矩阵，O₂为残差右奇异矩阵。

具体地，等于的正交列空间等于的列空间，因此对进行奇异值分解SVD得到公式14。

S188，依据公式14推导得出状态空间方程,状态空间方程如公式15所示。

其中，Q为非奇异矩阵。

S189，对状态空间方程中的矩阵K₂进行分解，得到Γ_i和H_i的估计值，从状态空间方程中提取出第一矩阵、第二矩阵、第三矩阵和第四矩阵，如公式16所示。

其中，I为单位矩阵。

具体地，本实施例最后可以提取出第一矩阵、第二矩阵、第三矩阵和第四矩阵的具体表达式。

在本申请的一实施例中，S300包括如下S310至S380：

S310，由鼓风机系统模型推导控制律，所述控制律如公式17所示。

u＝-Kx+Fr 公式17

其中，u为开度控制样本。x为实际状态样本。K为第一控制矩阵。F为第二控制矩阵。r为出口数据设定值。

具体地，r为出口数据设定值，出口数据设定值包括出口压力的设定值与流量的设定值。要想对鼓风机系统中的各个变量进行解耦，需要推导控制律。

S320，将控制律代入状态空间方程中，得到公式18；

其中，A为第一矩阵。B为第二矩阵。C为第三矩阵。K为第一控制矩阵。F为第二控制矩阵。r为出口压力与流量的设定值。x为实际状态样本。y为出口估计样本。为x的一阶导数。

S330，基于公式18得出目标矩阵，目标矩阵如公式19所示；

其中，H(s)为目标矩阵。G(s)为传递函数矩阵。I为单位矩阵。K为第一控制矩阵。F为第二控制矩阵。A为第一矩阵。B为第二矩阵。s为目标矩阵中的元素。

具体地，解耦控制方法的主要问题就转化为：寻找两个控制矩阵K和F，使得H(s)成为一个对角的非奇异的矩阵。相当于采用H(s)＝G(s)×[I+K(sI-A)^-1B]^-1对原有的对象进行补偿，从而获得我们所期望的系统响应。

S340，将第三矩阵转化为公式20的形式。

其中，m为元素个数。T为矩阵的转置符号。

S350，依据公式21定义解耦参数；

其中，σ_i为解耦参数。j为序号。

S360，依据公式22求取第五矩阵和第六矩阵；

其中，B^*为第五矩阵。C^*为第六矩阵。σ_m为解耦参数。m为矩阵内的元素个数。

S370，依据公式19、公式21和公式22导出第一控制矩阵和第二控制矩阵的表达式。第一控制矩阵和第二控制矩阵的表达式如公式22所示。

其中，F为第一控制矩阵。K为第二控制矩阵。B^*为第五矩阵。C^*为第六矩阵。

S380，依据公式23计算解耦后的闭环传递函数；

在所述S380之后，所述S300还包括：

S390，基于解耦后的闭环传递函数对鼓风机系统中的各个变量进行解耦，得到多个解耦后的变量。

在本申请一实施例中，所述S500包括：

S510，依据公式24，生成鼓风机的PID控制逻辑表达式；

x(k)＝K_p{err(k)+K_i∑err(k)+K_d[err(k)-err(k-1)]} 公式24

其中，x(k)为第k时刻控制器的输入值。K_p为比例参数。K_i为积分参数。K_d为微分参数。err(k)为当前时刻控制器的输出值与预设值的偏差。err(k-1)为前一时刻控制器的输出值与预设值的偏差。

在本申请一实施例中，在所述S500之后，所述高炉鼓风机的控制方法还包括：

S610，构建Redis数据库。

S620，将所有解耦后的变量存储入Redis数据库。

在执行S500根据计算出大量离散的，解耦后的变量后，将这些所有解耦后的变量存在Redis数据库中，然后执行S700。Redis数据库可以搭载在一个Redis服务器中。

具体地，在S700中，高炉、鼓风机及TRT设备的三维建模均采用Unity3D完成，所有模型均来源于实际现场设备，在外观上做到和现场设备一致，并拥有透视功能，实现对设备的全面解析。

在本申请一实施例中，所述S900包括如下S910至S930：

S910，向Redis数据库发送数据获取请求。

S920，从Redis数据库中获取所有解耦后的变量。

S930，将所有解耦后的变量依次导入鼓风机三维模型中。

具体地，通过Redis数据库与S700中构建好的鼓风机三维模型建立通信连接。通过连接，可以向入鼓风机三维模型导入解耦后的变量，可以实现实时控制鼓风机三维模型中各设备动作，例如：管道中气体流速、鼓风机叶片转速、高炉出口气体温度等。

本申请不但通过将解耦后的变量导入鼓风机三维模型中，实现高炉鼓风机控制流程和模型结构两者同步动作，从而实现高炉鼓风机设备的自动化运作，以及实时优化设备当前的运行状态。而且，结合数字孪生技术对于数据显示的优越性，能够更为清晰地显示过程中各个变量的数值，并且结合高炉鼓风机定风量定风压的控制方法，使得鼓风机三维模型可以更好的展现相关工艺的控制流程。

本申请还可以通过将输入PID控制器的控制参数和鼓风机三维模型进行交互，实现对高炉鼓风机在不同控制参数情况下的控制仿真，鼓风机三维模型也会依据相应变量的数值变化，实现了对三维资源的充分利用。

以上所述实施例的各技术特征可以进行任意的组合，各方法步骤也并不做执行顺序的限制，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。

以上所述实施例仅表达了本申请的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对本申请专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本申请构思的前提下，还可以做出若干变形和改进，这些都属于本申请的保护范围。因此，本申请的保护范围应以所附权利要求为准。