CN101884848B - 一种空分节能过程温度分布的非线性观测系统及方法 - Google Patents

Abstract

一种空分节能过程温度分布的非线性观测系统，包括与空分塔直接连接的现场智能仪表和DCS系统，所述DCS系统中包括存储装置、控制站和上位机，所述现场智能仪表与存储装置、控制站和上位机相连，所述的上位机包括用以预测空分塔将来时刻的温度分布的观测器，所述的观测器包括检测模块、温度分布函数拟合模块、将来时刻空分塔温度预测模块和温度预测误差补偿模块。本发明也提供了一种空分节能过程温度分布的非线性动态观测方法。本发明能够有效快速的在线预测空分节能过程的温度状态变化，准确地反应空分节能过程的强非线性性，耦合性等复杂非线性动态特性。

Description

一种空分节能过程温度分布的非线性观测系统及方法

技术领域

本发明涉及空分动态建模领域，特别地，涉及空分节能过程温度分布的非线性观测器系统设计及其方法。

背景技术

空分是对空气进行分离，得到氧、氮、氩等高纯工业气体的国民经济重要行业，其产品广泛用于石油、化工、冶金、电子、能源、航空航天、食品饮料、医疗保健等各种工业领域。

巨大的能量消耗一直为空分行业的瓶颈问题。在空分过程中，精馏塔是一个重要的传质单元操作过程，也是主要的耗能单元。由于空分精馏过程的强非线性，耦合性等复杂动态特性，传统的机理模型虽然具有较高的预测精度，但是在线运算效率太低，较难应用于现今较流行的先进控制方案，而目前商业使用的基于数据的线性辨识模型，往往具有较低的预测精度，对控制器的设计提出了较高的要求，从而建立空分精馏过程的高效率高精度非线性模型实现空分精馏过程中的温度分布等状态预测已经成为一项关键的空分节能技术。

发明内容

为了克服现有的空分节能过程温度分布的控制系统的预测精度低、精度较低、效率较低的不足，本发明的提供一种预测精度高、高效率和高精度的空分节能过程温度分布的非线性观测系统及方法。

本发明解决其技术问题所采用的技术方案是：

一种空分节能过程温度分布的非线性观测系统，包括与空分塔直接连接的现场智能仪表和DCS系统，所述DCS系统包括存储装置和上位机，所述现场智能仪表与存储装置、上位机相连，所述上位机包括用以预测空分塔将来时刻的温度分布的观测器，所述观测器包括：

温度分布函数拟合模块，用以采用历史数据库中的各块塔板温度T_i拟合温度分布函数，并将拟合参数存储到历史数据库当中，通过式(1)(2)得到：

${\hat{T}}_i=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r)}}$ i＝1，2，……，f-1       (1)

${\hat{T}}_i=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s)}}$ i＝f，f+1，……，n    (2)

其中，i为塔板编号，f为进料板编号，n为上塔塔底编号，

为第i块塔板处预估温度，T_min，r、T_max，r、k_r、T_min，s、T_max，s、k_s为拟合参数，S_r、S_s分别为空分塔下塔、上塔液相温度分布曲线的位置；

将来时刻空分塔温度预测模块，用以根据当前温度数据以及温度分布函数预测将来时刻的温度分布，该模块包括以下过程：

①当前时刻温度分布曲线位置变化速度计算，通过式(3)-(9)得到：

$X_{i,N}=\frac{P\×{\mathrm{\α}}_N\×{10}^{(\frac{T_i+c_N}{b_N}-a_N)}-1}{{\mathrm{\α}}_N-1}---\left(3\right)$

$X_{i,O}=\frac{P\×{\mathrm{\α}}_O\×{10}^{(\frac{T_i+c_O}{b_O}-a_O)}-1}{{\mathrm{\α}}_O-1}---\left(4\right)$

$y_{1,N}\left(k\right)=\frac{{\mathrm{\α}}_N{x}_{1,N}\left(k\right)}{({\mathrm{\α}}_N-1)x_{1,N}\left(k\right)+1}---\left(5\right)$

$y_{1,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{1,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{1,O}\left(k\right)+1}---\left(6\right)$

$y_{f,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{f,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{f,O}\left(k\right)+1}---\left(7\right)$

$\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)=\frac{[-V_1{y}_{1,N}\left(k\right)-L_n{x}_{n,N}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,N}\left(k\right)]({\mathrm{\α}}_N-1){(1+e^{-k_r(i-s_r)})}^2}{M(x_{n,N}\left(k\right)-x_{1,N}\left(k\right)){\mathrm{P\α}}_N(T_{\max ,r}-T_{\min ,r})k_r{e}^{-k_r(i-s_r)}\underset{i=1}{\overset{f-1}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_N-a_N{b}_N}}---\left(8\right)$

$\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)=\frac{[-V_f{y}_{f,O}\left(k\right)-L_n{x}_{n,O}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,O}\left(k\right)]({\mathrm{\α}}_O-1){(1+e^{-k_s(i-s_s)})}^2}{M(x_{n,O}\left(k\right)-x_{1,O}\left(k\right)){\mathrm{P\α}}_O(T_{\max ,s}-T_{\min ,s})k_s{e}^{-k_s(i-s_s)}\underset{i=f}{\overset{n}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_O-a_O{b}_O}}---\left(9\right)$

其中，k为当前采样时刻，下角标i为塔板编号，1为下塔塔顶编号，n为上塔塔底编号，下角标N、O分别代表氮和氧，上角标f代表进料，X_i，N、X_i，O分别为空分塔中第i块塔板处氮、氧的液相组分浓度，P为空分塔压强，T_i为第i块塔板处的温度，α_N、α_O分别为氮和氧相对于氩的相对挥发度，a_N、b_N、c_N、a_O、b_O、c_O为安东尼常数，F_i(k)为k采样时刻第i块塔板进料流量，L_n为上塔塔底液相流率，V₁为下塔塔顶汽相流率，x_n，N(k)、x_n，O(k)分别为k采样时刻上塔塔底氮、氧的液相组分浓度，x_l，N(k)、x_l，O(k)分别为k采样时刻下塔塔顶氮、氧的液相组分浓度，y_f，O(k)分别为k采样时刻进料板氧汽相组分浓度，y_1，N(k)、y_1，O(k)分别为k采样时刻下塔塔顶氮、氧汽相组分浓度，x^f _i，N(k)，x^f _i，O(k)分别为第i块塔板的进料液氮组分浓度和进料液氧组分浓度，M为塔板持液量，

分别为下塔上塔温度分布位置变化速度；

②下一时刻波形位置以及下一时刻温度分布，通过式(10)-(13)得到：

$S_r(k+1)=S_r\left(k\right)\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(10\right)$

$S_s(k+1)=S_s\left(k\right)\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(11\right)$

${\hat{T}}_i(k+1)=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r(k+1))}}$ i＝1，2，……，f-1    (12)

${\hat{T}}_i(k+1)=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s(k+1))}}$ i＝f，f+1，……，n    (13)

其中，δ为采样周期，S_r(k)、S_s(k)分别为k采样时刻下塔、上塔温度分布曲线位置，S_r(k+1)、S_s(k+1)分别为k+1采样时刻下塔、上塔温度分布曲线位置，为k+1时刻下塔或上塔第i块塔板温度的预测值；

温度预测误差补偿模块，用以从历史数据库中获取当前温度测量值，以及预测值，计算预测误差，修正温度预测值，通过式(14)(15)得到

$r_i\left(k\right)=T_i\left(k\right)-{\hat{T}}_i\left(k\right)$ i＝1，2，……，n    (14)

$T_i^{\%}(k+1)={\hat{T}}_i(k+1)+r_i\left(k\right)$ i＝1，2，……，n    (15)

其中，T_i(k)为k+1采样时刻第i块塔板的温度，r_i(k)为k采样时刻第i块塔板的温度预测误差，为k+1采样时刻第i块塔板的温度预测修正值。

作为优选的一种方案：所述的现场智能仪表与DCS系统中的存储装置、上位机通过现场总线相互连接。

再进一步，所述的上位机还包括人机界面模块，用以设定采样周期，设定预测时间长度，显示观测器输出的将来时刻温度预测值和当前时刻温度测量值。

一种空分节能过程温度分布的非线性观测方法，所述的观测方法包括以下步骤：

1)确定采样周期，并将氮和氧相对于氩的相对挥发度α_N、α_O，安东尼常数a_N、b_N、c_N、a_O、b_O、c_O保存在历史数据库当中；

2)从现场智能仪表获取k采样时刻压强P，各塔板温度T_i，计算第i块塔板处氮、氧的液相组分浓度值，通过式(3)-(4)得到：

$X_{i,N}=\frac{P\×{\mathrm{\α}}_N\×{10}^{(\frac{T_i+c_N}{b_N}-a_N)}-1}{{\mathrm{\α}}_N-1}---\left(3\right)$

$X_{i,O}=\frac{P\×{\mathrm{\α}}_O\×{10}^{(\frac{T_i+c_O}{b_O}-a_O)}-1}{{\mathrm{\α}}_O-1}---\left(4\right)$

其中，X_i，N为空分塔第i块塔板处氮的液相组分浓度，X_i，O为空分塔中第i块塔板处氧的液相组分浓度；

检测k采样时刻下塔塔顶汽相流量和上塔塔底液相流量值，同上塔压强数据，各塔板温度数据，组分浓度的测量值一起存储到观测器系统的历史数据库当中；

3)从历史数据库当中获取k采样时刻温度的测量值，拟合k采样时刻的温度分布函数如式(1)(2)：

${\hat{T}}_i=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r)}}$ i＝1，2，……，f-1    (1)

${\hat{T}}_i=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s)}}$ i＝f，f+1，……，n    (2)

其中，1为下塔塔顶编号，f为进料板编号，n为上塔塔底编号，

为第i块塔板处预估温度，T_min，r，T_max，r、k_r、T_min，s、T_max，s、k_s为拟合参数，S_r、S_s分别为空分塔下塔、上塔液相温度分布曲线的位置；

4)计算k采样时刻的下塔上塔温度分布位置变化速度，通过式(5)-(9)得到：

$y_{1,N}\left(k\right)=\frac{{\mathrm{\α}}_N{x}_{1,N}\left(k\right)}{({\mathrm{\α}}_N-1)x_{1,N}\left(k\right)+1}---\left(5\right)$

$y_{1,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{1,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{1,O}\left(k\right)+1}---\left(6\right)$

$y_{f,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{f,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{f,O}\left(k\right)+1}---\left(7\right)$

$\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)=\frac{[-V_1{y}_{1,N}\left(k\right)-L_n{x}_{n,N}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,N}\left(k\right)]({\mathrm{\α}}_N-1){(1+e^{-k_r(i-s_r)})}^2}{M(x_{n,N}\left(k\right)-x_{1,N}\left(k\right)){\mathrm{P\α}}_N(T_{\max ,r}-T_{\min ,r})k_r{e}^{-k_r(i-s_r)}\underset{i=1}{\overset{f-1}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_N-a_N{b}_N}}---\left(8\right)$

$\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)=\frac{[-V_f{y}_{f,O}\left(k\right)-L_n{x}_{n,O}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,O}\left(k\right)]({\mathrm{\α}}_O-1){(1+e^{-k_s(i-s_s)})}^2}{M(x_{n,O}\left(k\right)-x_{1,O}\left(k\right)){\mathrm{P\α}}_O(T_{\max ,s}-T_{\min ,s})k_s{e}^{-k_s(i-s_s)}\underset{i=f}{\overset{n}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_O-a_O{b}_O}}---\left(9\right)$

其中，k为当前采样时刻，下角标i为塔板编号，1为下塔塔顶编号，n为上塔塔底的编号，下角标N、O分别代表氮和氧，上角标f代表进料，X_i，N、X_i，O分别为空分塔中第i块塔板处氮、氧的液相组分浓度，P为空分塔压强，T_i为第i块塔板处的温度，α_N、α_O分别为氮和氧相对于氩的相对挥发度，a_N、b_N、c_N、a_O、b_O、c_O为安东尼常数，F_i(k)为k采样时刻第i块塔板进料流量，L_n为上塔塔底液相流率，V₁为下塔塔顶汽相流率，x_n，N(k)、x_n，O(k)分别为k采样时刻上塔塔底氮、氧的液相组分浓度，x_l，N(k)、x_l，O(k)分别为k采样时刻下塔塔顶氮、氧的液相组分浓度，y_f，O(k)分别为k采样时刻进料板氧汽相组分浓度，y_1，N(k)、y_1，O(k)分别为k采样时刻下塔塔顶氮、氧汽相组分浓度，x^f _i，N(k)，x^f _i，O(k)分别为第i块塔板的进料液氮组分浓度和进料液氧组分浓度，M为塔板持液量，

分别为下塔上塔温度分布位置变化速度；

5)下一时刻波形位置以及下一时刻温度分布，通过式(10)-(13)得到：

$S_r(k+1)=S_r\left(k\right)\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(10\right)$

$S_s(k+1)=S_s\left(k\right)\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(11\right)$

${\hat{T}}_i(k+1)=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r(k+1))}}$ i＝1，2，……，f-1    (12)

${\hat{T}}_i(k+1)=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s(k+1))}}$ i＝f，f+1，……，n    (13)

其中，δ为采样周期，S_r(k)、S_s(k)分别为k采样时刻下塔，上塔温度分布曲线位置，S_r(k+1)、S_s(k+1)分别为k+1采样时刻下塔，上塔温度分布曲线位置，

为k+1时刻下塔或上塔的第i块塔板温度的预测值；

6)从历史数据库中获取当前温度测量值，以及预测值，计算预测误差，修正温度预测值，通过式(14)(15)得到

$r_i\left(k\right)=T_i\left(k\right)-{\hat{T}}_i\left(k\right)$ i＝1，2，……，n    (14)

$T_i^{\%}(k+1)={\hat{T}}_i(k+1)+r_i\left(k\right)$ i＝1，2，……，n    (15)

其中，T_i(k)为k+1采样时刻第i块塔板的温度，r_i(k)为k采样时刻第i块塔板的温度预测误差，

为k+1采样时刻的温度预测修正值；

进一步，所述历史数据库为DCS系统的存储装置，控制站读取历史数据库，显示观测器工作过程状态。

本发明的有益效果主要表现在：1.基于空分机理建立高精度的简化预测模型，能够高效率地准确地跟踪空分节能过程温度的复杂非线性变化特征；2.在线运行效率高，温度状态变量预测比传统基于数据的辨识模型精确。

附图说明

图1是本发明所提出的空分节能过程温度非线性观测系统的结构图。

图2是上位机观测器实现方法的原理图。

具体实施方式

下面根据附图具体说明本发明。

实施例1

参照图1、图2，一种空分节能过程温度分布的非线性观测系统，包括与空分塔1直接连接的现场智能仪表2和DCS系统13，所述DCS系统13包括存储装置4、控制站5及上位机6，所述现场智能仪表2与存储装置4、控制站5和上位机6依次相连，上位机包括用以预测空分塔将来时刻的温度分布的观测器，所述观测器包括温度分布函数拟合模块9，将来时刻温度预测模块11，温度预测误差补偿模块10。

所述的温度分布函数拟合模块9，采用历史数据库中的各块塔板温度T_i拟合温度分布函数，并将拟合参数存储到历史数据库当中，通过式(1)(2)得到：

${\hat{T}}_i=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r)}}$ i＝1，2，……，f-1    (1)

${\hat{T}}_i=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s)}}$ i＝f，f+1，……，n    (2)

其中，1为下塔塔顶编号，f为进料板编号，n为上塔塔底编号，

为第i块塔板处预估温度，T_min，r、T_max，r、k_r、T_min，s、T_max，s、k_s为拟合参数，S_r、S_s分别为空分塔下塔、上塔液相温度分布曲线的位置。

所述的将来时刻空分塔温度预测模块11，根据当前温度数据以及温度分布函数预测将来时刻的温度分布，该模块包括以下几部分

①当前时刻温度分布曲线位置变化速度计算，通过式(3)-(9)得到：

$X_{i,N}=\frac{P\×{\mathrm{\α}}_N\×{10}^{(\frac{T_i+c_N}{b_N}-a_N)}-1}{{\mathrm{\α}}_N-1}---\left(3\right)$

$X_{i,O}=\frac{P\×{\mathrm{\α}}_O\×{10}^{(\frac{T_i+c_O}{b_O}-a_O)}-1}{{\mathrm{\α}}_O-1}---\left(4\right)$

$y_{1,N}\left(k\right)=\frac{{\mathrm{\α}}_N{x}_{1,N}\left(k\right)}{({\mathrm{\α}}_N-1)x_{1,N}\left(k\right)+1}---\left(5\right)$

$y_{1,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{1,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{1,O}\left(k\right)+1}---\left(6\right)$

$y_{f,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{f,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{f,O}\left(k\right)+1}---\left(7\right)$

$\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)=\frac{[-V_1{y}_{1,N}\left(k\right)-L_n{x}_{n,N}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,N}\left(k\right)]({\mathrm{\α}}_N-1){(1+e^{-k_r(i-s_r)})}^2}{M(x_{n,N}\left(k\right)-x_{1,N}\left(k\right)){\mathrm{P\α}}_N(T_{\max ,r}-T_{\min ,r})k_r{e}^{-k_r(i-s_r)}\underset{i=1}{\overset{f-1}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_N-a_N{b}_N}}---\left(8\right)$

$\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)=\frac{[-V_f{y}_{f,O}\left(k\right)-L_n{x}_{n,O}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,O}\left(k\right)]({\mathrm{\α}}_O-1){(1+e^{-k_s(i-s_s)})}^2}{M(x_{n,O}\left(k\right)-x_{1,O}\left(k\right)){\mathrm{P\α}}_O(T_{\max ,s}-T_{\min ,s})k_s{e}^{-k_s(i-s_s)}\underset{i=f}{\overset{n}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_O-a_O{b}_O}}---\left(9\right)$

其中，k为当前采样时刻，下角标i为塔板编号，1为下塔塔顶编号，n为上塔塔底的编号，下角标N、O分别代表氮和氧，上角标f代表进料，X_i，N、X_i，O分别为空分塔中第i块塔板处氮、氧的液相组分浓度，P为空分塔压强，T_i为第i块塔板处的温度，α_N、α_O分别为氮和氧相对于氩的相对挥发度，a_N、b_N、c_N、a_O、b_O、c_O为安东尼(Antonie)常数，F_i(k)为k采样时刻第i块塔板进料流量，L_n为上塔塔底液相流率，V₁为下塔塔顶汽相流率，x_n，N(k)、x_n，O(k)分别为k采样时刻上塔塔底氮、氧的液相组分浓度，x_l，N(k)、x_l，O(k)分别为k采样时刻下塔塔顶氮、氧的液相组分浓度，y_f，O(k)分别为k采样时刻进料板氧汽相组分浓度，y_1，N(k)、y_1，O(k)分别为k采样时刻下塔塔顶氮、氧汽相组分浓度，x^f _i，N(k)，x^f _i，O(k)分别为第i块塔板的进料液氮组分浓度和进料液氧组分浓度，M为塔板持液量，

分别为下塔、上塔温度分布位置变化速度；

②下一时刻波形位置以及下一时刻温度分布，通过式(10)-(13)得到：

$S_r(k+1)=S_r\left(k\right)\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(10\right)$

$S_s(k+1)=S_s\left(k\right)\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(11\right)$

${\hat{T}}_i(k+1)=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r(k+1))}}$ i＝1，2，……，f-1    (12)

${\hat{T}}_i(k+1)=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s(k+1))}}$ i＝f，f+1，……，n    (13)

其中，δ为采样周期，S_r(k)、S_s(k)分别为k采样时刻下塔、上塔温度分布曲线位置，S_r(k+1)、S_s(k+1)分别为k+1采样时刻下塔，上塔温度分布曲线位置，

为k+1时刻下塔或上塔的第i块塔板温度的预测值；

所述的温度预测误差补偿模块10，从历史数据库中获取当前温度测量值，以及预测值，计算预测误差，修正温度预测值，通过式(14)(15)得到：

$r_i\left(k\right)=T_i\left(k\right)-{\hat{T}}_i\left(k\right)$ i＝1，2，……，n    (14)

$T_i^{\%}(k+1)={\hat{T}}_i(k+1)+r_i\left(k\right)$ i＝1，2，……，n    (15)

其中，T_i(k)为k+1采样时刻第i块塔板的温度，r_i(k)为k采样时刻第i块塔板的温度预测误差，

为k+1采样时刻第i块塔板的温度预测修正值。

所述的现场智能仪表2、DCS系统中的存储装置4.、控制站5及上位机6通过现场总线和数据接口3依次相连；所述的智能仪表用以实现检测模块7的功能即检测得到空分塔的温度、压强数据，所述的检测模块7包括温度检测元件和压力检测元件，所述的数据接口用以实现I/O模块8的功能即数据的输入输出。

所述的上位机包括人机界面模块12，用于设定采样周期，设定预测时间长度，显示观测器输出的将来时刻温度预测值和当前时刻温度测量值。

实施例2

参照图1和图2，一种空分节能过程温度分布的非线性观测方法，包括以下步骤：

1)确定采样周期，并将氮和氧相对于氩的相对挥发度α_N、α_O，安东尼(Antonie)常数a_N、b_N、c_N、a_O、b_O、c_O保存在历史数据库当中；

2)从现场智能仪表获取k采样时刻压强P，各塔板温度T_i，计算第i块塔板处氮、氧的液相组分浓度值，通过式(3)-(4)得到：

$X_{i,N}=\frac{P\×{\mathrm{\α}}_N\×{10}^{(\frac{T_i+c_N}{b_N}-a_N)}-1}{{\mathrm{\α}}_N-1}---\left(3\right)$

$X_{i,O}=\frac{P\×{\mathrm{\α}}_O\×{10}^{(\frac{T_i+c_O}{b_O}-a_O)}-1}{{\mathrm{\α}}_O-1}---\left(4\right)$

其中，X_i，N为空分塔第i块塔板处氮的液相组分浓度，X_i，O为空分塔中第i块塔板处氧的液相组分浓度；

检测k采样时刻下塔塔顶汽相流量和上塔塔底液相流量值，同上塔压强数据，各塔板温度数据，组分浓度的测量值一起存储到观测器系统的历史数据库当中；

3)从历史数据库当中获取k采样时刻温度的测量值，拟合k采样时刻的温度分布函数如式(1)(2)：

${\hat{T}}_i=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r)}}$ i＝1，2，……，f-1    (1)

${\hat{T}}_i=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s)}}$ i＝f，f+1，……，n    (2)

其中，1为下塔塔顶编号，f为进料板编号，n为上塔塔底编号，

为第i块塔板处预估温度，T_min，r，T_max，r、k_r、T_min，s、T_max，s、k_s为拟合参数，S_r、S_s分别为空分塔下塔、上塔液相温度分布曲线的位置；

4)计算k采样时刻的下塔上塔温度分布位置变化速度，通过式(5)-(9)得到：

$y_{1,N}\left(k\right)=\frac{{\mathrm{\α}}_N{x}_{1,N}\left(k\right)}{({\mathrm{\α}}_N-1)x_{1,N}\left(k\right)+1}---\left(5\right)$

$y_{1,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{1,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{1,O}\left(k\right)+1}---\left(6\right)$

$y_{f,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{f,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{f,O}\left(k\right)+1}---\left(7\right)$

$\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)=\frac{[-V_1{y}_{1,N}\left(k\right)-L_n{x}_{n,N}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,N}\left(k\right)]({\mathrm{\α}}_N-1){(1+e^{-k_r(i-s_r)})}^2}{M(x_{n,N}\left(k\right)-x_{1,N}\left(k\right)){\mathrm{P\α}}_N(T_{\max ,r}-T_{\min ,r})k_r{e}^{-k_r(i-s_r)}\underset{i=1}{\overset{f-1}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_N-a_N{b}_N}}---\left(8\right)$

$\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)=\frac{[-V_f{y}_{f,O}\left(k\right)-L_n{x}_{n,O}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,O}\left(k\right)]({\mathrm{\α}}_O-1){(1+e^{-k_s(i-s_s)})}^2}{M(x_{n,O}\left(k\right)-x_{1,O}\left(k\right)){\mathrm{P\α}}_O(T_{\max ,s}-T_{\min ,s})k_s{e}^{-k_s(i-s_s)}\underset{i=f}{\overset{n}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_O-a_O{b}_O}}---\left(9\right)$

其中，k为当前采样时刻，下角标i为塔板编号，1为下塔塔顶编号，n为上塔塔底的编号，下角标N、O分别代表氮和氧，上角标f代表进料，X_i，N、X_i，O分别为空分塔中第i块塔板处氮、氧的液相组分浓度，P为空分塔压强，T_i为第i块塔板处的温度，α_N、α_O分别为氮和氧相对于氩的相对挥发度，a_N、b_N、c_N、a_O、b_O、c_O为安东尼常数，F_i(k)为k采样时刻第i块塔板进料流量，L_n为上塔塔底液相流率，V₁为下塔塔顶汽相流率，x_n，N(k)、x_n，O(k)分别为k采样时刻上塔塔底氮、氧的液相组分浓度，x_l，N(k)、x_l，O(k)分别为k采样时刻下塔塔顶氮、氧的液相组分浓度，y_f，O(k)分别为k采样时刻进料板氧汽相组分浓度，y_1，N(k)、y_1，O(k)分别为k采样时刻下塔塔顶氮、氧汽相组分浓度，，x^f _i，N(k)，x^f _i，O(k)分别为第i块塔板的进料液氮组分浓度和进料液氧组分浓度，M为塔板持液量，

分别为下塔上塔温度分布位置变化速度；

5)下一时刻波形位置以及下一时刻温度分布，通过式(10)-(13)得到：

$S_r(k+1)=S_r\left(k\right)\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(10\right)$

$S_s(k+1)=S_s\left(k\right)\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(11\right)$

${\hat{T}}_i(k+1)=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r(k+1))}}$ i＝1，2，……，f-1        (12)

${\hat{T}}_i(k+1)=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s(k+1))}}$ i＝f，f+1，……，n    (13)

其中，δ为采样周期，S_r(k)、S_s(k)分别为k采样时刻下塔，上塔温度分布曲线位置，S_r(k+1)、S_s(k+1)分别为k+1采样时刻下塔，上塔温度分布曲线位置，

为k+1时刻下塔或上塔的第i块塔板温度的预测值；

6)从历史数据库中获取当前温度测量值，以及预测值，计算预测误差，修正温度预测值，通过式(14)(15)得到

$r_i\left(k\right)=T_i\left(k\right)-{\hat{T}}_i\left(k\right)$ i＝1，2，……，n    (14)

$T_i^{\%}(k+1)={\hat{T}}_i(k+1)+r_i\left(k\right)$ i＝1，2，……，n    (15)

其中，r_i(k)为k采样时刻第i块塔板的温度预测误差，

为k+1采样时刻的温度预测修正值。

所述的历史数据库为DCS系统中的存储装置，所述的DCS系统包括数据接口，存储装置，控制站，其中控制站可以读取历史数据库，显示观测器工作过程状态。

上述实施例用来解释说明本发明，而不是对本发明进行限制，在本发明的精神和权利要求的保护范围内，对本发明作出的任何修改和改变，都落入本发明的保护范围。

Claims (5)

1.一种空分节能过程温度分布的非线性观测系统，包括与空分塔直接连接的现场智能仪表和DCS系统，所述DCS系统包括存储装置、控制站和上位机，所述现场智能仪表与存储装置、控制站和上位机相连，其特征在于：所述上位机包括用以预测空分塔将来时刻的温度分布的观测器，所述观测器包括：

温度分布函数拟合模块，用以采用历史数据库中的各块塔板温度T_i拟合温度分布函数，并将拟合参数存储到历史数据库当中，通过式(1)(2)得到：

${\hat{T}}_i=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r)}}$ i＝1，2，……，f-1    (1)

${\hat{T}}_i=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s)}}$ i＝f，f+1，……，n    (2)

其中，i为塔板编号，f为进料板编号，n为上塔塔底编号，为第i块塔板处预估温度，T_min，r、T_max，r、k_r、T_min，s、T_max，s、k_s为拟合参数，S_r、S_s分别为空分塔下塔、上塔液相温度分布曲线的位置；

将来时刻空分塔温度预测模块，用以根据当前温度数据以及温度分布函数预测将来时刻的温度分布，该模块包括以下过程：

①当前时刻温度分布曲线位置变化速度计算，通过式(3)-(9)得到：

$X_{i,N}=\frac{P\×{\mathrm{\α}}_N\×{10}^{(\frac{T_i+c_N}{b_N}-a_N)}-1}{{\mathrm{\α}}_N-1}---\left(3\right)$

$X_{i,O}=\frac{P\×{\mathrm{\α}}_O\×{10}^{(\frac{T_i+c_O}{b_O}-a_O)}-1}{{\mathrm{\α}}_O-1}---\left(4\right)$

$y_{1,N}\left(k\right)=\frac{{\mathrm{\α}}_N{x}_{1,N}\left(k\right)}{({\mathrm{\α}}_N-1)x_{1,N}\left(k\right)+1}---\left(5\right)$

$y_{1,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{1,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{1,O}\left(k\right)+1}---\left(6\right)$

$y_{f,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{f,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{f,O}\left(k\right)+1}---\left(7\right)$

$\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)=\frac{[-V_1{y}_{1,N}\left(k\right)-L_n{x}_{n,N}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,N}\left(k\right)]({\mathrm{\α}}_N-1){(1+e^{-k_r(i-s_r)})}^2}{M(x_{n,N}\left(k\right)-x_{1,N}\left(k\right)){\mathrm{P\α}}_N(T_{\max ,r}-T_{\min ,r})k_r{e}^{-k_r(i-s_r)}\underset{i=1}{\overset{f-1}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_N-a_N{b}_N}}---\left(8\right)$

$\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)=\frac{[-V_f{y}_{f,O}\left(k\right)-L_n{x}_{n,O}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,O}\left(k\right)]({\mathrm{\α}}_O-1){(1+e^{-k_s(i-s_s)})}^2}{M(x_{n,O}\left(k\right)-x_{1,O}\left(k\right)){\mathrm{P\α}}_O(T_{\max ,s}-T_{\min ,s})k_s{e}^{-k_s(i-s_s)}\underset{i=f}{\overset{n}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_O-a_O{b}_O}}---\left(9\right)$

其中，k为当前采样时刻，下角标i为塔板编号，1为下塔塔顶编号，n为上塔塔底编号，下角标N、O分别代表氮和氧，上角标f代表进料，X_i，N、X_i，O分别为空分塔中第i块塔板处氮、氧的液相组分浓度，P为空分塔压强，T_i为第i块塔板处的温度，α_N、α_O分别为氮和氧相对于氩的相对挥发度，a_N、b_N、c_N、a_O、b_O、c_O为安东尼常数，F_i(k)为k采样时刻第i块塔板进料流量，L_n为上塔塔底液相流率，V₁为下塔塔顶汽相流率，x_n，N(k)、x_n，O(k)分别为k采样时刻上塔塔底氮、氧的液相组分浓度，x_l，N(k)、x_l，O(k)分别为k采样时刻下塔塔顶氮、氧的液相组分浓度，y_f，O(k)分别为k采样时刻进料板氧汽相组分浓度，y_1，N(k)、y_1，O(k)分别为k采样时刻下塔塔顶氮、氧汽相组分浓度，x^f _i，N(k)，x^f _i，O(k)分别为第i块塔板的进料液氮组分浓度和进料液氧组分浓度，M为塔板持液量，

分别为下塔、上塔温度分布位置变化速度；

②下一时刻波形位置以及下一时刻温度分布，通过式(10)-(13)得到：

$S_r(k+1)=S_r\left(k\right)\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(10\right)$

$S_s(k+1)=S_s\left(k\right)\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(11\right)$

${\hat{T}}_i(k+1)=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r(k+1))}}$ i＝1，2，……，f-1    (12)

${\hat{T}}_i(k+1)=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s(k+1))}}$ i＝f，f+1，……，n    (13)

其中，δ为采样周期，S_r(k)、S_s(k)分别为k采样时刻下塔、上塔温度分布曲线位置，S_r(k+1)、S_s(k+1)分别为k+1采样时刻下塔、上塔温度分布曲线位置，

为k+1时刻下塔或上塔第i块塔板温度的预测值；

温度预测误差补偿模块，用以从历史数据库中获取当前温度测量值，以及预测值，计算预测误差，修正温度预测值，通过式(14)(15)得到

$r_i\left(k\right)=T_i\left(k\right)-{\hat{T}}_i\left(k\right)$ i＝1，2，……，n    (14)

$X_i^{\%}(k+1)={\hat{T}}_i(k+1)+r_i\left(k\right)$ i＝1，2，……，n    (15)

其中，T_i(k)为k+1采样时刻第i块塔板的温度，r_i(k)为k采样时刻第i块塔板的温度预测误差，

为k+1采样时刻第i块塔板的温度预测修正值。

2.如权利要求1中所述的空分节能过程温度分布的非线性观测系统，其特征在于：所述的现场智能仪表与DCS系统中的存储装置、上位机通过现场总线相互连接。

3.如权利要求1或2所述的空分节能过程温度分布的非线性观测系统，其特征在于：所述的上位机还包括人机界面模块，用以设定采样周期，设定预测时间长度，显示观测器输出的将来时刻温度预测值和当前时刻温度测量值。

4.一种如权利要求1所述的空分节能过程温度分布的非线性观测系统实现的非线性观测方法，其特征在于：所述的观测方法包括以下步骤：

1)确定采样周期，并将氮和氧相对于氩的相对挥发度α_N、α_O，安东尼常数a_N、b_N、c_N、a_O、b_O、c_O保存在历史数据库当中；

2)从现场智能仪表获取k采样时刻压强P，各塔板温度T_i，计算第i块塔板处氮、氧的液相组分浓度值，通过式(3)-(4)得到：

$X_{i,N}=\frac{P\×{\mathrm{\α}}_N\×{10}^{(\frac{T_i+c_N}{b_N}-a_N)}-1}{{\mathrm{\α}}_N-1}---\left(3\right)$

$X_{i,O}=\frac{P\×{\mathrm{\α}}_O\×{10}^{(\frac{T_i+c_O}{b_O}-a_O)}-1}{{\mathrm{\α}}_O-1}---\left(4\right)$

其中，X_i，N为空分塔第i块塔板处氮的液相组分浓度，X_i，O为空分塔中第i块塔板处氧的液相组分浓度；

检测k采样时刻下塔塔顶汽相流量和上塔塔底液相流量值，同上塔压强数据，各塔板温度数据，组分浓度的测量值一起存储到观测器系统的历史数据库当中；

3)从历史数据库当中获取k采样时刻温度的测量值，拟合k采样时刻的温度分布函数如式(1)(2)：

${\hat{T}}_i=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r)}}$ i＝1，2，……，f-1        (1)

${\hat{T}}_i=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s)}}$ i＝f，f+1，……，n    (2)

其中，1为下塔塔顶编号，f为进料板编号，n为上塔塔底编号，

为第i块塔板处预估温度，T_min，r，T_max，r、k_r、T_min，s、T_max，s、k_s为拟合参数，S_r、S_s分别为空分塔下塔、上塔液相温度分布曲线的位置；

4)计算k采样时刻的下塔上塔温度分布位置变化速度，通过式(5)-(9)得到：

$y_{1,N}\left(k\right)=\frac{{\mathrm{\α}}_N{x}_{1,N}\left(k\right)}{({\mathrm{\α}}_N-1)x_{1,N}\left(k\right)+1}---\left(5\right)$

$y_{1,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{1,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{1,O}\left(k\right)+1}---\left(6\right)$

$y_{f,O}\left(k\right)=\frac{{\mathrm{\α}}_O{x}_{f,O}\left(k\right)}{({\mathrm{\α}}_O-1)x_{f,O}\left(k\right)+1}---\left(7\right)$

$\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)=\frac{[-V_1{y}_{1,N}\left(k\right)-L_n{x}_{n,N}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,N}\left(k\right)]({\mathrm{\α}}_N-1){(1+e^{-k_r(i-s_r)})}^2}{M(x_{n,N}\left(k\right)-x_{1,N}\left(k\right)){\mathrm{P\α}}_N(T_{\max ,r}-T_{\min ,r})k_r{e}^{-k_r(i-s_r)}\underset{i=1}{\overset{f-1}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_N-a_N{b}_N}}---\left(8\right)$

$\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)=\frac{[-V_f{y}_{f,O}\left(k\right)-L_n{x}_{n,O}\left(k\right)+\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{F}_i\left(k\right){x^f}_{i,O}\left(k\right)]({\mathrm{\α}}_O-1){(1+e^{-k_s(i-s_s)})}^2}{M(x_{n,O}\left(k\right)-x_{1,O}\left(k\right)){\mathrm{P\α}}_O(T_{\max ,s}-T_{\min ,s})k_s{e}^{-k_s(i-s_s)}\underset{i=f}{\overset{n}{\mathrm{\Σ}}}\frac{\ln 10}{T_i+c_O-a_O{b}_O}}---\left(9\right)$

其中，k为当前采样时刻，下角标i为塔板编号，1为下塔塔顶编号，n为上塔塔底的编号，下角标N、O分别代表氮和氧，上角标f代表进料，X_i，N、X_i，O分别为空分塔中第i块塔板处氮、氧的液相组分浓度，P为空分塔压强，T_i为第i块塔板处的温度，α_N、α_O分别为氮和氧相对于氩的相对挥发度，a_N、b_N、c_N、a_O、b_O、c_O为安东尼常数，F_i(k)为k采样时刻第i块塔板进料流量，L_n为上塔塔底液相流率，V₁为下塔塔顶汽相流率，x_n，N(k)、x_n，O(k)分别为k采样时刻上塔塔底氮、氧的液相组分浓度，x_l，N(k)、x_l，O(k)分别为k采样时刻下塔塔顶氮、氧的液相组分浓度，y_f，O(k)分别为k采样时刻进料板氧汽相组分浓度，y_1，N(k)、y_1，O(k)分别为k采样时刻下塔塔顶氮、氧汽相组分浓度，x^f _i，N(k)，x^f _i，O(k)分别为第i块塔板的进料液氮组分浓度和进料液氧组分浓度，M为塔板持液量，

分别为下塔上塔温度分布位置变化速度；

5)下一时刻波形位置以及下一时刻温度分布，通过式(10)-(13)得到：

$S_r(k+1)=S_r\left(k\right)+\frac{{\mathrm{dS}}_r}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(10\right)$

$S_s(k+1)=S_s\left(k\right)+\frac{{\mathrm{dS}}_s}{\mathrm{dt}}\left(k\right)\×\mathrm{\δ}---\left(11\right)$

${\hat{T}}_i(k+1)=T_{\min ,r}+\frac{T_{\max ,r}-T_{\min ,r}}{{1+e}^{-k_r(i-S_r(k+1))}}$ i＝1，2，……，f-1    (12)

${\hat{T}}_i(k+1)=T_{\min ,s}+\frac{T_{\max ,s}-T_{\min ,s}}{{1+e}^{-k_s(i-S_s(k+1))}}$ i＝f，f+1，……，n    (13)

其中，δ为采样周期，S_r(k)、S_s(k)分别为k采样时刻下塔，上塔温度分布曲线位置，S_r(k+1)、S_s(k+1)分别为k+1采样时刻下塔，上塔温度分布曲线位置，

为k+1时刻下塔或上塔的第i块塔板温度的预测值；

6)从历史数据库中获取当前温度测量值，以及预测值，计算预测误差，修正温度预测值，通过式(14)(15)得到

$r_i\left(k\right)=T_i\left(k\right)-{\hat{T}}_i\left(k\right)$ i＝1，2，……，n    (14)

$T_i^{\%}(k+1)={\hat{T}}_i(k+1)+r_i\left(k\right)$ i＝1，2，……，n    (15)

其中，，T_i(k)为k+1采样时刻第i块塔板的温度，r_i(k)为k采样时刻第i块塔板的温度预测误差，

为k+1采样时刻的温度预测修正值。

5.如权利要求4所述的非线性观测方法，其特征在于：所述历史数据库为DCS系统的存储装置，控制站读取历史数据库，显示观测器工作过程状态。

Images