CN101840220B - 内部热耦合精馏塔平衡级节能控制系统及方法 - Google Patents

Abstract

一种内部热耦合精馏塔平衡级节能控制系统，包括与内部热耦合精馏塔连接的现场智能仪表、控制站、数据库以及上位机，上位机包括：信号采集模块，用以采集当前生产工况数据；节能控制模块，过程：设定塔的结构参数和操作参数，指定进料流量初值；假定各塔板液相组成；对每一个塔板，分别计算其平衡温度和汽相组成、汽液相的焓值和汽液相流量；判断条件是否成立，如果成立，则继续，否则，更新各塔板液相组成；判断产品纯度是否满足约束，如果不满足则结束迭代，如果满足则将进料流量增加一个迭代步长Δ，继续迭代。以及提出了一种内部热耦合精馏塔平衡级节能控制方法。本发明在当前生产工况条件下使得内部热耦合精馏塔单位能耗最小、节能性好。

Description

内部热耦合精馏塔平衡级节能控制系统及方法

技术领域

本发明涉及精馏节能领域，特别地，涉及一种内部热耦合精馏塔平衡级节能控制系统及方法。

背景技术

自二十世纪七十年代经历两次“石油危机”，能源成为全球关注的问题以来，节能作为解决能源危机问题的根本途径之一，逐渐受到各国的高度重视。

精馏过程是石油炼制、石油化工和其它化工过程中应用最为广泛的传质单元操作过程，也是石油化工领域中能耗最大的单元操作之一，其能耗约占化工厂总能耗的三分之一，有时甚至还更多。美国的统计数据表明，美国1976年精馏耗能约占全国总能耗的3％，如果精馏过程节能10％，相当于每天节省100000桶石油，若以当时的每桶石油40美元计算，相当于美国1976年的精馏过程的操作费用每天可以节省400万美元，全年仅精馏过程可以节省约15亿美元。

从二十世纪七十年代以来，研究者已经从热力学、操作原理等角度提出了大量的精馏过程的节能方法。内部热耦合精馏塔是这些方法中最吸引人的一种，目前精馏塔节能研究的一个前沿。研究结果表明，内部热耦合精馏塔与常规精馏塔最小回流比下的能耗和操作费用相比还可以节省30％以上。过程在线优化是指综合应用过程建模技术、优化技术、先进控制技术以及计算机技术，在满足生产安全要求以及产品质量约束等条件下，不断计算并改变过程的操作条件，使得生产过程始终运行在“最优状态”。

研究表明，在线优化所带来的经济效益，相当于DCS和各种先进控制手段所带来的经济效益的总和，其投资却只有先进控制的三分之一。可见在线优化的成功实施是企业获取最小经济效益的必要保证。内部热耦合精馏塔平衡级节能控制是指基于内部热耦合精馏塔平衡级模型，在满足产品纯度要求的前提下，不断计算并改变精馏塔的进料流量，使得精馏过程始终运行在产量最大即单位能耗最小的状态。

发明内容

为了克服目前内部热耦合精馏塔尚无节能控制系统、单位能耗较大的不足，本发明提供一种在当前生产工况条件下使得内部热耦合精馏塔单位能耗最小、节能性好的内部热耦合精馏塔平衡级节能控制系统及方法。

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

一种内部热耦合精馏塔平衡级节能控制系统，包括与内部热耦合精馏塔连接的现场智能仪表、控制站、数据库以及上位机，所述现场智能仪表与控制站、数据库和上位机连接，所述的上位机包括：

信号采集模块，用以采集当前生产工况数据；

节能控制模块，用以节能控制，采用以下过程来完成：

1)设定塔的结构参数和操作参数，指定进料流量初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，分别计算汽液相的焓值；

5)由式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，H^G和H^L分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算：

Q＝UAΔT        (3)

其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差；

6)判断下式是否成立，如果成立，则继续7)，否则，更新各塔板液相组成，返回3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j{z}_{i,j}-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

其中，x是液相组成，y是汽相组成，z是进料组成，下标i＝1，...M表示组分，M表示组分数；

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回2)继续迭代。

作为优选的一种方案：所述的上位机还包括：泡点法模块，用以由泡点法计算其平衡温度和汽相组成，其过程如下：

3.1)假定塔板平衡温度；

3.2)计算汽液平衡常数，采用以下过程完成：

$\ln {\mathrm{\&Phi;}}_i^L=\ln \frac{\mathrm{RT}}{P(v^L-b^L)}-\frac{b_i}{b^L}(1-Z^L)+{\mathrm{\&xi;}}^L{a}^L(\frac{b_i}{b^L}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^L})/b^L\mathrm{RT}---\left(5\right)$

$\ln {\mathrm{\&Phi;}}_i^G=\ln \frac{\mathrm{RT}}{P(v^G-b^G)}-\frac{b_i}{b^G}(1-Z^G)+{\mathrm{\&xi;}}^G{a}^G(\frac{b_i}{b^G}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^G})/b^G\mathrm{RT}---\left(6\right)$

$K_i={\mathrm{\&Phi;}}_i^L/{\mathrm{\&Phi;}}_i^G---\left(7\right)$

y_i＝K_ix_i           (8)

其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1，...M表示组分，M表示组分数，摩尔体积v、物性参数b^G、b^L、b_i、a^G、a^L、a_i，m、ξ^G、ξ^L、汽相压缩因子Z^G、液相压缩因子Z^L由物性模块计算；

3.3)检验 $|1-\underset{i}{\mathrm{\&Sigma;}}{y}_i|<0.0001$ 是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代。

作为优选的另一种方案：所述的上位机还包括：焓模块，用以计算汽液相混合焓，其过程如下：

$H_i^{*}=c_i+d_{i}T+e_i{T}^2+{f_{i}T}^3+h_i{T}^4---\left(9\right)$

$H^{*}=\underset{i}{\mathrm{\&Sigma;}}{y}_i{H}_i^{*}---\left(10\right)$

$H^G=H^{*}-\mathrm{RT}(1-Z^G)-{\mathrm{\&xi;}}^G(a^G-T\frac{d{a}^G}{\mathrm{dT}})/b^G---\left(11\right)$

$H^L=H^{*}-\mathrm{RT}(1-Z^L)-{\mathrm{\&xi;}}^L(a^L-T\frac{d{a}^L}{\mathrm{dT}})/b^L---\left(12\right)$

其中H_i ^*表示第i个纯组分理想气体的焓值，H^*是混合物理想气体焓值，c、d、e、f、h为常数。

作为优选的再一种方案：所述的上位机还包括：物性模块，用以计算物性参数，其过程如下：

$a_{i,m}={\mathrm{\&Omega;}}_{\mathrm{ai},m}{R}^2{T}_{\mathrm{ci},m}^2/P_{\mathrm{ci},m}---\left(13\right)$

b_i＝Ω_bRT_ci/P_cia              (14)

$T_{\mathrm{ci},m}=\sqrt{T_{\mathrm{ci}}{T}_{\mathrm{cj}}}(1-k_{i,m})---\left(15\right)$

$V_{\mathrm{ci},m}=0.125{(V_{\mathrm{ci},m}^{1/3}+V_{\mathrm{ci},m}^{1/3})}^3---\left(16\right)$

Z_ci，m＝0.5(Z_ci+Z_cm)          (17)

P_ci，m＝RT_ci，mZ_ci，m/V_ci，m  (18)

Ω_ai，m＝0.5(Ω_ai+Ω_am)       (19)

对汽相：

$a^G=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{y}_i{y}_m{a}_{i,m}---\left(20\right)$

$b^G=\underset{i}{\mathrm{\&Sigma;}}{y}_i{b}_i---\left(21\right)$

令

A^G＝a^GP/R²T²                  (22)

B^G＝b^GP/RT                    (23)

α^G＝2B^G-1                    (24)

${\mathrm{\&beta;}}^G=A^G-3B^G-5B^{G^2}---\left(25\right)$

${\mathrm{\&gamma;}}^G=2(B^{G^3}+B^{G^2})-A^G{B}^G---\left(26\right)$

取初值为1-0.6P_r，用牛顿法解如下方程，即得到汽相压缩因子Z^G

$Z^{G^3}+{\mathrm{\&alpha;}}^G{Z}^{G^2}+\mathrm{\&beta;}{Z}^G+{\mathrm{\&gamma;}}^G=0---\left(27\right)$

则，

v^G＝RT/PZ^G                (28)

${\mathrm{\&xi;}}^G=0.242536\ln \frac{v^G+3.561553b^G}{v^G-0.561553b^G}---\left(29\right)$

对液相：

$a^L=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{x}_i{x}_m{a}_{i,m}---\left(30\right)$

$b^L=\underset{i}{\mathrm{\&Sigma;}}{x}_i{b}_i---\left(31\right)$

令

A^L＝a^LP/R²T²              (32)

B^L＝b^LP/RT                (33)

α^L＝2B^L-1                (34)

${\mathrm{\&beta;}}^L=A^L-3B^L-5B^{L^2}---\left(35\right)$

${\mathrm{\&gamma;}}^L=2(B^{L^3}+B^{L^2})-A^L{B}^L---\left(36\right)$

取初值为P_r(0.106+0.078P_r)，用牛顿法解如下方程，即得到液相压缩因子Z^L

$Z^{L^3}+{\mathrm{\&alpha;}}^L{Z}^{L^2}+\mathrm{\&beta;}{Z}^L+{\mathrm{\&gamma;}}^L=0---\left(37\right)$

则，

v^L＝RT/PZ^L     (38)

${\mathrm{\&xi;}}^L=0.242536\ln \frac{v^L+3.561553b^L}{v^L-0.561553b^L}---\left(39\right)$

Ω_ai＝C_i-D_iτ+E_iτ²-W_iτ³ (40)

Ω_b＝0.070721             (41)

τ＝0.01T                 (42)

其中，A、B、α、β、γ、τ、Ω_a、Ω_b是中间变量，C、D、E、W是常数，T_c、P_c、V_c、Z_c分别是临界温度、压力、体积和压缩因子，P_r是对比压力，R是气体常数，k_i，m表示第i组分和第m组分的二元交互系数，k_i，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物。

进一步，所述的上位机还包括：结果显示模块，用于将节能控制结果传给控制站进行显示，并通过现场总线将节能控制结果传递到现场操作站进行显示。

一种用所述的内部热耦合精馏塔平衡级节能控制系统实现的节能控制方法，所述的节能控制方法包括以下步骤：

1)设定塔的结构参数，采集塔的生产工况数据，采集当前进料流量作为初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，计算其汽液相的焓值；

5)联立式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，下标j-1、j、j+1分别表示第j-1、j、j+1块板，Q表示热耦合量，由下式计算：

Q＝UAΔT            (3)；

6)判断式(4)是否成立，如果成立，则继续步骤7)，否则，更新液相组成，返回步骤3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j{z}_{i,j}-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回步骤2)继续迭代。

作为优选的一种方案：所述步骤3)中，由泡点法计算平衡温度和汽相组成的过程如下：

3.1)假定塔板平衡温度；

3.2)计算汽液平衡常数，采用以下过程完成：

$\ln {\mathrm{\&Phi;}}_i^L=\ln \frac{\mathrm{RT}}{P(v^L-b^L)}-\frac{b_i}{b^L}(1-Z^L)+{\mathrm{\&xi;}}^L{a}^L(\frac{b_i}{b^L}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^L})/b^L\mathrm{RT}---\left(5\right)$

$\ln {\mathrm{\&Phi;}}_i^G=\ln \frac{\mathrm{RT}}{P(v^G-b^G)}-\frac{b_i}{b^G}(1-Z^G)+{\mathrm{\&xi;}}^G{a}^G(\frac{b_i}{b^G}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^G})/b^G\mathrm{RT}---\left(6\right)$

$K_i={\mathrm{\&Phi;}}_i^L/{\mathrm{\&Phi;}}_i^G---\left(7\right)$

y_i＝K_ix_i           (8)

其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1，...M表示组分，M表示组分数，摩尔体积v、物性参数b^G、b^L、b_i、a^G、a^L、a_i，m、ξ^G、ξ^L、汽相压缩因子Z^G、液相压缩因子Z^L由物性参数计算单元计算；

3.3)检验 $|1-\underset{i}{\mathrm{\&Sigma;}}{y}_i|<0.0001$ 是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回步骤3.2)继续迭代。

作为优选的另一种方案：所述步骤4)中，计算汽液相混合焓，其过程如下：

$H_i^{*}=c_i+d_{i}T+e_i{T}^2+{f_{i}T}^3+h_i{T}^4---\left(9\right)$

$H^{*}=\underset{i}{\mathrm{\&Sigma;}}{y}_i{H}_i^{*}---\left(10\right)$

$H^G=H^{*}-\mathrm{RT}(1-Z^G)-{\mathrm{\&xi;}}^G(a^G-T\frac{d{a}^G}{\mathrm{dT}})/b^G---\left(11\right)$

$H^L=H^{*}-\mathrm{RT}(1-Z^L)-{\mathrm{\&xi;}}^L(a^L-T\frac{d{a}^L}{\mathrm{dT}})/b^L---\left(12\right)$

其中H_i ^*表示第i个纯组分理想气体的焓值，H^*是混合物理想气体焓值，c、d、e、f、h为常数。

作为优选的再一种方案：所述物性参数计算方法包括以下步骤：

$a_{i,m}={\mathrm{\&Omega;}}_{\mathrm{ai},m}{R}^2{T}_{\mathrm{ci},m}^2/P_{\mathrm{ci},m}---\left(13\right)$

b_i＝Ω_bRT_ci/P_cia            (14)

$T_{\mathrm{ci},m}=\sqrt{T_{\mathrm{ci}}{T}_{\mathrm{cj}}}(1-k_{i,m})---\left(15\right)$

$V_{\mathrm{ci},m}=0.125{(V_{\mathrm{ci},m}^{1/3}+V_{\mathrm{ci},m}^{1/3})}^3---\left(16\right)$

Z_ci，m＝0.5(Z_ci+Z_cm)         (17)

P_ci，m＝RT_ci，mZ_ci，m/V_ci，m (18)

Ω_ai，m＝0.5(Ω_ai+Ω_am)      (19)

对汽相：

$a^G=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{y}_i{y}_m{a}_{i,m}---\left(20\right)$

$b^G=\underset{i}{\mathrm{\&Sigma;}}{y}_i{b}_i---\left(21\right)$

令

A^G＝a^GP/R²T²                 (22)

B^G＝b^GP/RT                   (23)

α^G＝2B^G-1                   (24)

${\mathrm{\&beta;}}^G=A^G-3B^G-5B^{G^2}---\left(25\right)$

${\mathrm{\&gamma;}}^G=2(B^{G^3}+B^{G^2})-A^G{B}^G---\left(26\right)$

取初值为1-0.6P_r，用牛顿法解如下方程，即得到汽相压缩因子Z^G

$Z^{G^3}+{\mathrm{\&alpha;}}^G{Z}^{G^2}+\mathrm{\&beta;}{Z}^G+{\mathrm{\&gamma;}}^G=0---\left(27\right)$

则，

v^G＝RT/PZ^G                   (28)

${\mathrm{\&xi;}}^G=0.242536\ln \frac{v^G+3.561553b^G}{v^G-0.561553b^G}---\left(29\right)$

对液相：

$a^L=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{x}_i{x}_m{a}_{i,m}---\left(30\right)$

$b^L=\underset{i}{\mathrm{\&Sigma;}}{x}_i{b}_i---\left(31\right)$

令

A^L＝a^LP/R²T²              (32)

B^L＝b^LP/RT                (33)

α^L＝2B^L-1                (34)

${\mathrm{\&beta;}}^L=A^L-3B^L-5B^{L^2}---\left(35\right)$

${\mathrm{\&gamma;}}^L=2(B^{L^3}+B^{L^2})-A^L{B}^L---\left(36\right)$

取初值为P_r(0.106+0.078P_r)，用牛顿法解如下方程，即得到液相压缩因子Z^L

$Z^{L^3}+{\mathrm{\&alpha;}}^L{Z}^{L^2}+\mathrm{\&beta;}{Z}^L+{\mathrm{\&gamma;}}^L=0---\left(37\right)$

则，

v^L＝RT/PZ^L                (38)

${\mathrm{\&xi;}}^L=0.242536\ln \frac{v^L+3.561553b^L}{v^L-0.561553b^L}---\left(39\right)$

Ω_ai＝C_i-D_iτ+E_iτ²-W_iτ³ (40)

Ω_b＝0.070721             (41)

τ＝0.01T                 (42)

其中，A、B、α、β、γ、τ、Ω_a、Ω_b是中间变量，C、D、E、W是常数，T_c、P_c、V_c、Z_c分别是临界温度、压力、体积和压缩因子，P_r是对比压力，R是气体常数，k_i，m表示第i组分和第m组分的二元交互系数，k_i，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物。

进一步，在所述的步骤7)中，上位机将节能控制结果传给控制站进行显示，并通过现场总线将节能控制结果传递到现场操作站进行显示。

本发明的有益效果主要表现在：1、对内部热耦合精馏塔进行平衡级节能控制；2、可以用于指导生产，在满足当前生产工况要求的前提下提高产量；3、降低单位产品能耗，提高生产效益。

附图说明

图1是本发明所提出的平衡级节能控制系统的硬件结构图。

图2是本发明所述内部热耦合精馏塔结构示意图。

图3是本发明上位机的功能结构图。

具体实施方式

下面结合附图对本发明作进一步描述。

实施例1

参照图1、图2、图3，一种内部热耦合精馏塔平衡级节能控制系统，包括内部热耦合精馏塔1连接的现场智能仪表2和控制站4、数据库5以及上位机6，所述现场智能仪表2与数据接口3连接，所述数据接口3与控制站4、数据库5以及上位机6连接，所述的上位机6包括：

信号采集模块7，用以采集当前生产工况数据；

节能控制模块8，用以节能控制，采用以下过程来完成：

1)设定塔的结构参数和操作参数，指定进料流量初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，分别计算汽液相的焓值；

5)由式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，H^G和H^L分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算：

Q＝UAΔT                   (3)

其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差；

6)判断下式(4)是否成立，如果成立，则继续7)，否则，更新各塔板液相组成，返回3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j{z}_{i,j}-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

其中，x是液相组成，y是汽相组成，z是进料组成，下标i＝1，...M表示组分，M表示组分数；

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回2)继续迭代。

所述上位机6还包括：泡点法模块9，用以由泡点法计算其平衡温度和汽相组成，其过程如下：

3.1)假定塔板平衡温度；

3.2)计算汽液平衡常数，采用以下过程完成：

$\ln {\mathrm{\&Phi;}}_i^L=\ln \frac{\mathrm{RT}}{P(v^L-b^L)}-\frac{b_i}{b^L}(1-Z^L)+{\mathrm{\&xi;}}^L{a}^L(\frac{b_i}{b^L}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^L})/b^L\mathrm{RT}---\left(5\right)$

$\ln {\mathrm{\&Phi;}}_i^G=\ln \frac{\mathrm{RT}}{P(v^G-b^G)}-\frac{b_i}{b^G}(1-Z^G)+{\mathrm{\&xi;}}^G{a}^G(\frac{b_i}{b^G}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^G})/b^G\mathrm{RT}---\left(6\right)$

$K_i={\mathrm{\&Phi;}}_i^L/{\mathrm{\&Phi;}}_i^G---\left(7\right)$

y_i＝K_ix_i               (8)

其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1，...，M表示组分，M表示组分数，摩尔体积v、物性参数b^G、b^L、b_i、a^G、a^L、a_i，m、ξ^G、ξ^L、汽相压缩因子Z^G、液相压缩因子Z^L由物性模块计算；

3.3)检验 $|1-\underset{i}{\mathrm{\&Sigma;}}{y}_i|<0.0001$ 是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代；

所述上位机6还包括：焓模块10，用以计算汽液相混合焓，其过程如下：

$H_i^{*}=c_i+d_{i}T+e_i{T}^2+{f_{i}T}^3+h_i{T}^4---\left(9\right)$

$H^{*}=\underset{i}{\mathrm{\&Sigma;}}{y}_i{H}_i^{*}---\left(10\right)$

$H^G=H^{*}-\mathrm{RT}(1-Z^G)-{\mathrm{\&xi;}}^G(a^G-T\frac{d{a}^G}{\mathrm{dT}})/b^G---\left(11\right)$

$H^L=H^{*}-\mathrm{RT}(1-Z^L)-{\mathrm{\&xi;}}^L(a^L-T\frac{d{a}^L}{\mathrm{dT}})/b^L---\left(12\right)$

其中H_i ^*表示第i个纯组分理想气体的焓值，H^*是混合物理想气体焓值，c、d、e、f、h为常数。

所述上位机6还包括：物性模块11，用以计算物性参数，其过程如下：

$a_{i,m}={\mathrm{\&Omega;}}_{\mathrm{ai},m}{R}^2{T}_{\mathrm{ci},m}^2/P_{\mathrm{ci},m}---\left(13\right)$

b_i＝Ω_bRT_ci/P_cia               (14)

$T_{\mathrm{ci},m}=\sqrt{T_{\mathrm{ci}}{T}_{\mathrm{cj}}}(1-k_{i,m})---\left(15\right)$

$V_{\mathrm{ci},m}=0.125{(V_{\mathrm{ci},m}^{1/3}+V_{\mathrm{ci},m}^{1/3})}^3---\left(16\right)$

Z_ci，m＝0.5(Z_ci+Z_cm)           (17)

P_ci，m＝RT_ci，mZ_ci，m/V_ci，m   (18)

Ω_ai，m＝0.5(Ω_ai+Ω_am)        (19)

对汽相：

$a^G=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{y}_i{y}_m{a}_{i,m}---\left(20\right)$

$b^G=\underset{i}{\mathrm{\&Sigma;}}{y}_i{b}_i---\left(21\right)$

令

A^G＝a^GP/R²T²                   (22)

B^G＝b^GP/RT                     (23)

α^G＝2B^G-1                     (24)

${\mathrm{\&beta;}}^G=A^G-3B^G-5B^{G^2}---\left(25\right)$

${\mathrm{\&gamma;}}^G=2(B^{G^3}+B^{G^2})-A^G{B}^G---\left(26\right)$

取初值为1-0.6P_r，用牛顿法解如下方程，即得到汽相压缩因子Z^G

$Z^{G^3}+{\mathrm{\&alpha;}}^G{Z}^{G^2}+\mathrm{\&beta;}{Z}^G+{\mathrm{\&gamma;}}^G=0---\left(27\right)$

则，

v^G＝RT/PZ^G                       (28)

${\mathrm{\&xi;}}^G=0.242536\ln \frac{v^G+3.561553b^G}{v^G-0.561553b^G}---\left(29\right)$

对液相：

$a^L=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{x}_i{x}_m{a}_{i,m}---\left(30\right)$

$b^L=\underset{i}{\mathrm{\&Sigma;}}{x}_i{b}_i---\left(31\right)$

令

A^L＝a^LP/R²T²                     (32)

B^L＝b^LP/RT                       (33)

α^L＝2B^L-1                       (34)

${\mathrm{\&beta;}}^L=A^L-3B^L-5B^{L^2}---\left(35\right)$

${\mathrm{\&gamma;}}^L=2(B^{L^3}+B^{L^2})-A^L{B}^L---\left(36\right)$

取初值为P_r(0.106+0.078P_r)，用牛顿法解如下方程，即得到液相压缩因子Z^L

$Z^{L^3}+{\mathrm{\&alpha;}}^L{Z}^{L^2}+\mathrm{\&beta;}{Z}^L+{\mathrm{\&gamma;}}^L=0---\left(37\right)$

则，

v^L＝RT/PZ^L                       (38)

${\mathrm{\&xi;}}^L=0.242536\ln \frac{v^L+3.561553b^L}{v^L-0.561553b^L}---\left(39\right)$

Ω_ai＝C_i-D_iτ+E_iτ²-W_iτ³        (40)

Ω_b＝0.070721                    (41)

τ＝0.01T                      (42)

其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，T_c、P_c、V_c、Z_c分别是临界温度、压力、体积和压缩因子，P_r是对比压力，R是气体常数，k_i，m表示第i组分和第m组分的二元交互系数，k_i，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物，Ω_a、Ω_b是中间变量；

所述的上位机6还包括：结果显示模块12，用于将节能控制结果传给控制站进行显示，并通过现场总线将节能控制结果传递到现场操作站进行显示。

本实施例的内部热耦合精馏塔节能潜力优化系统的硬件结构图如附图1所示，所述的优化系统核心由包括信号采集模块7、节能控制模块8、泡点法模块9、焓模块10、物性模块11、结果显示模块12和人机界面的上位机6构成，此外还包括：现场智能仪表2，数据接口3、控制站4、数据库5和现场总线。内部热耦合精馏塔1、智能仪表2、数据接口3、控制站4、数据库5、上位机6通过现场总线依次相连，实现信息流的上传和下达。节能控制系统在上位机6上运行，可以方便地与底层系统进行信息交换。

本实施例的优化系统的功能结构图如附图3所示，主要包括信号采集模块7、节能控制模块8、泡点法模块9、焓模块10、物性模块11、结果显示模块12等。

所述的平衡级节能控制方法按照如下步骤进行实施：

1)设定塔的结构参数，采集塔的生产工况数据，采集当前进料流量作为初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，计算其汽液相的焓值；

5)联立式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，H^G和H^L分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算：

Q＝UAΔT           (3)；

6)判断式(4)是否成立，如果成立，则继续步骤7)，否则，更新液相组成，返回步骤3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j{z}_{i,j}-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回步骤2)继续迭代。

实施例2

参照图1、图2、图3，一种内部热耦合精馏塔平衡级节能控制方法，所述的平衡级节能控制方法包括以下步骤：

1)设定塔的结构参数，采集塔的生产工况数据，采集当前进料流量作为初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，计算其汽液相的焓值；

5)联立式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，H^G和H^L分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算：

Q＝UAΔT                   (3)；

6)判断式(4)是否成立，如果成立，则继续步骤7)，否则，更新液相组成，返回步骤3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j{z}_{i,j}-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回步骤2)继续迭代。

所述步骤3)中，由泡点法计算平衡温度和汽相组成的过程如下：

3.1)假定塔板平衡温度；

3.2)计算汽液平衡常数，采用以下过程完成：

$\ln {\mathrm{\&Phi;}}_i^L=\ln \frac{\mathrm{RT}}{P(v^L-b^L)}-\frac{b_i}{b^L}(1-Z^L)+{\mathrm{\&xi;}}^L{a}^L(\frac{b_i}{b^L}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^L})/b^L\mathrm{RT}---\left(5\right)$

$\ln {\mathrm{\&Phi;}}_i^G=\ln \frac{\mathrm{RT}}{P(v^G-b^G)}-\frac{b_i}{b^G}(1-Z^G)+{\mathrm{\&xi;}}^G{a}^G(\frac{b_i}{b^G}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^G})/b^G\mathrm{RT}---\left(6\right)$

$K_i={\mathrm{\&Phi;}}_i^L/{\mathrm{\&Phi;}}_i^G---\left(7\right)$

y_i＝K_ix_i               (8)

其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1，...M表示组分，M表示组分数，摩尔体积v、物性参数b^G、b^L、b_i、a^G、a^L、a_i，m、ξ^G、ξ^L、汽相压缩因子Z^G、液相压缩因子Z^L由物性参数计算单元计算；

3.3)检验 $|1-\underset{i}{\mathrm{\&Sigma;}}{y}_i|<0.0001$ 是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回步骤3.2)继续迭代。

所述步骤4)中，计算汽液相混合焓，其过程如下：

$H_i^{*}=c_i+d_{i}T+e_i{T}^2+{f_{i}T}^3+h_i{T}^4---\left(9\right)$

$H^{*}=\underset{i}{\mathrm{\&Sigma;}}{y}_i{H}_i^{*}---\left(10\right)$

$H^G=H^{*}-\mathrm{RT}(1-Z^G)-{\mathrm{\&xi;}}^G(a^G-T\frac{d{a}^G}{\mathrm{dT}})/b^G---\left(11\right)$

$H^L=H^{*}-\mathrm{RT}(1-Z^L)-{\mathrm{\&xi;}}^L(a^L-T\frac{d{a}^L}{\mathrm{dT}})/b^L---\left(12\right)$

其中H_i ^*表示第i个纯组分理想气体的焓值，H^*是混合物理想气体焓值，c、d、e、f、h为常数。

所述物性参数计算方法包括以下步骤：

$a_{i,m}={\mathrm{\&Omega;}}_{\mathrm{ai},m}{R}^2{T}_{\mathrm{ci},m}^2/P_{\mathrm{ci},m}---\left(13\right)$

b_i＝Ω_bRT_ci/P_cia              (14)

$T_{\mathrm{ci},m}=\sqrt{T_{\mathrm{ci}}{T}_{\mathrm{cj}}}(1-k_{i,m})---\left(15\right)$

$V_{\mathrm{ci},m}=0.125{(V_{\mathrm{ci},m}^{1/3}+V_{\mathrm{ci},m}^{1/3})}^3---\left(16\right)$

Z_ci，m＝0.5(Z_ci+Z_cm)          (17)

P_ci，m＝RT_ci，mZ_ci，m/V_ci，m  (18)

Ω_ai，m＝0.5(Ω_ai+Ω_am)       (19)

对汽相：

$a^G=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{y}_i{y}_m{a}_{i,m}---\left(20\right)$

$b^G=\underset{i}{\mathrm{\&Sigma;}}{y}_i{b}_i---\left(21\right)$

令

A^G＝a^GP/R²T²                  (22)

B^G＝b^GP/RT                    (23)

α^G＝2B^G-1                    (24)

${\mathrm{\&beta;}}^G=A^G-3B^G-5B^{G^2}---\left(25\right)$

${\mathrm{\&gamma;}}^G=2(B^{G^3}+B^{G^2})-A^G{B}^G---\left(26\right)$

取初值为1-0.6P_r，用牛顿法解如下方程，即得到汽相压缩因子Z^G

$Z^{G^3}+{\mathrm{\&alpha;}}^G{Z}^{G^2}+\mathrm{\&beta;}{Z}^G+{\mathrm{\&gamma;}}^G=0---\left(27\right)$

则，

v^G＝RT/PZ^G                 (28)

${\mathrm{\&xi;}}^G=0.242536\ln \frac{v^G+3.561553b^G}{v^G-0.561553b^G}---\left(29\right)$

对液相：

$a^L=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{x}_i{x}_m{a}_{i,m}---\left(30\right)$

$b^L=\underset{i}{\mathrm{\&Sigma;}}{x}_i{b}_i---\left(31\right)$

令

A^L＝a^LP/R²T²               (32)

B^L＝b^LP/RT                 (33)

α^L＝2B^L-1                 (34)

${\mathrm{\&beta;}}^L=A^L-3B^L-5B^{L^2}---\left(35\right)$

${\mathrm{\&gamma;}}^L=2(B^{L^3}+B^{L^2})-A^L{B}^L---\left(36\right)$

取初值为P_r(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子Z^L

$Z^{L^3}+{\mathrm{\&alpha;}}^L{Z}^{L^2}+\mathrm{\&beta;}{Z}^L+{\mathrm{\&gamma;}}^L=0---\left(37\right)$

则，

v^L＝RT/PZ^L                 (38)

${\mathrm{\&xi;}}^L=0.242536\ln \frac{v^L+3.561553b^L}{v^L-0.561553b^L}---\left(39\right)$

Ω_ai＝C_i-D_iτ+E_iτ²-W_iτ³  (40)

Ω_b＝0.070721              (41)

τ＝0.01T                  (42)

其中，A、B、α、β、γ、τ、Ω_a、Ω_b是中间变量，C、D、E、W是常数，T_c、P_c、V_c、Z_c分别是临界温度、压力、体积和压缩因子，P_r是对比压力，R是气体常数，k_i，m表示第i组分和第m组分的二元交互系数，k_i，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物。

在所述的步骤7)中，上位机将节能控制结果传给控制站进行显示，并通过现场总线将计算结果传递到现场操作站进行显示。

本发明所提出的内部热耦合精馏塔平衡级节能控制系统及方法，已通过上述具体实施步骤进行了描述，相关技术人员明显能在不脱离本发明内容、精神和范围内对本文所述的装置和操作方法进行改动或适当变更与组合，来实现本发明技术。特别需要指出的是，所有相类似的替换和改动对本领域的技术人员是显而易见的，它们都会被视为包括在本发明精神、范围和内容中。

Claims (4)

1.一种内部热耦合精馏塔平衡级节能控制系统，包括与内部热耦合精馏塔连接的现场智能仪表、控制站、数据库以及上位机，所述现场智能仪表与控制站、数据库和上位机连接，其特征在于：所述的上位机包括：

信号采集模块，用以采集当前生产工况数据；

节能控制模块，用以节能控制，采用以下过程来完成：

1)设定塔的结构参数和操作参数，指定进料流量初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，分别计算汽液相的焓值；

5)由式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，H^G和H^L分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算：

Q＝UAΔT                                           (3)

其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差；

6)判断下式(4)是否成立，如果成立，则继续7)，否则，更新各塔板液相组成，返回3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j{z}_{i,j}-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

其中，x是液相组成，y是汽相组成，z是进料组成，下标i＝1，...M表示组分，M表示组分数；

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回2)继续迭代；

所述的上位机还包括：泡点法模块，用以由泡点法计算其平衡温度和汽相组成，其过程如下：

3.1)假定塔板平衡温度；

3.2)计算汽液平衡常数，采用以下过程完成：

$\ln {\mathrm{\&Phi;}}_i^L=\ln \frac{\mathrm{RT}}{P(v^L-b^L)}-\frac{b_i}{b^L}(1-Z^L)+{\mathrm{\&xi;}}^L{a}^L(\frac{b_i}{b^L}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^L})/b^L\mathrm{RT}---\left(5\right)$

$\ln {\mathrm{\&Phi;}}_i^G=\ln \frac{\mathrm{RT}}{P(v^G-b^G)}-\frac{b_i}{b^G}(1-Z^G)+{\mathrm{\&xi;}}^G{a}^G(\frac{b_i}{b^G}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^G})/b^G\mathrm{RT}---\left(6\right)$

$K_i={\mathrm{\&Phi;}}_i^L/{\mathrm{\&Phi;}}_i^G---\left(7\right)$

y_i＝K_ix_i                                               (8)

其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1，...M表示组分，M表示组分数，摩尔体积v、物性参数b^G、b^L、b_i、a^G、a^L、a_i，m、ξ^G、ξ^L、汽相压缩因子Z^G、液相压缩因子Z^L由物性模块计算；

3.3)检验

是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代；

所述的上位机还包括：焓模块，用以计算汽液相混合焓，其过程如下：

$H_i^{*}=c_i+d_{i}T+e_i{T}^2+f_i{T}^3+h_i{T}^4---\left(9\right)$

$H^{*}=\underset{i}{\mathrm{\&Sigma;}}{y}_i{H}_i^{*}---\left(10\right)$

$H^G=H^{*}-\mathrm{RT}(1-Z^G)-{\mathrm{\&xi;}}^G(a^G-T\frac{{\mathrm{da}}^G}{\mathrm{dT}}){/b}^G---\left(11\right)$

$H^L=H^{*}-\mathrm{RT}(1-Z^L)-{\mathrm{\&xi;}}^L(a^L-T\frac{{\mathrm{da}}^L}{\mathrm{dT}})/b^L---\left(12\right)$

其中

表示第i个纯组分理想气体的焓值，H^*是混合物理想气体焓值，c、d、e、f、h为常数；

所述的上位机还包括：物性模块，用以计算物性参数，其过程如下：

$a_{i,m}={\mathrm{\&Omega;}}_{\mathrm{ai},m}{R}^2{T}_{\mathrm{ci},m}^2/P_{\mathrm{ci},m}---\left(13\right)$

b_i＝Ω_bRT_ci/P_cia                                       (14)

$T_{\mathrm{ci},m}=\sqrt{T_{\mathrm{ci}}{T}_{\mathrm{cj}}}(1-k_{i,m})---\left(15\right)$

$V_{\mathrm{ci},m}=0.125{(V_{\mathrm{ci},m}^{1/3}+V_{\mathrm{ci},m}^{1/3})}^3---\left(16\right)$

Z_ci，m＝0.5(Z_ci+Z_cm)                                   (17)

P_ci，m＝RT_ci，mZ_ci，m/V_ci，m                           (18)

Ω_ai，m＝0.5(Ω_ai+Ω_am)                                (19)

对汽相：

$a^G=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{y}_i{y}_m{a}_{i,m}---\left(20\right)$

$b^G=\underset{i}{\mathrm{\&Sigma;}}{y}_i{b}_i---\left(21\right)$

令

A^G＝a^GP/R²T²                                           (22)

B^G＝b^GP/RT                                             (23)

α^G＝2B^G-1                                             (24)

${\mathrm{\&beta;}}^G=A^G-3B^G-5B^{G^2}---\left(25\right)$

${\mathrm{\&gamma;}}^G=2(B^{G^3}+B^{G^2})-A^G{B}^G---\left(26\right)$

取初值为1-0.6P_r，用牛顿法解如下方程，即得到汽相压缩因子Z^G

$Z^{G^3}+{\mathrm{\&alpha;}}^G{Z}^{G^2}+{\mathrm{\&beta;Z}}^G+{\mathrm{\&gamma;}}^G=0---\left(27\right)$

则，

v^G＝RT/PZ^G                                             (28)

${\mathrm{\&xi;}}^G=0.242536\ln \frac{v^G+3.561553b^G}{v^G-0.561553b^G}---\left(29\right)$

对液相：

$a^L=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{x}_i{x}_m{a}_{i,m}---\left(30\right)$

$b^L=\underset{i}{\mathrm{\&Sigma;}}{x}_i{b}_i---\left(31\right)$

令

A^L＝a^LP/R²T²                                           (32)

B^L＝b^LP/RT                                             (33)

α^L＝2B^L-1                                             (34)

${\mathrm{\&beta;}}^L=A^L-3B^L-5B^{L^2}---\left(35\right)$

${\mathrm{\&gamma;}}^L=2(B^{L^3}+B^{L^2})-A^L{B}^L---\left(36\right)$

取初值为P_r(0.106+0.078P_r)，用牛顿法解如下方程，即得到液相压缩因子Z^L

$Z^{L^3}+{\mathrm{\&alpha;}}^L{Z}^{L^2}+{\mathrm{\&beta;Z}}^L+{\mathrm{\&gamma;}}^L=0---\left(37\right)$

则，

v^L＝RT/PZ^L                                             (38)

${\mathrm{\&xi;}}^L=0.242536\ln \frac{v^L+3.561553b^L}{v^L-0.561553b^L}---\left(39\right)$

Ω_ai＝C_i-D_iτ+E_iτ²-W_iτ³                              (40)

Ω_b＝0.070721                                          (41)

τ＝0.01T                                              (42)

其中，A、B、α、β、γ、τ、Ω_a、Ω_b是中间变量，C、D、E、W是常数，T_c、P_c、V_c、Z_c分别是临界温度、压力、体积和压缩因子，P_r是对比压力，R是气体常数，k_i，m表示第i组分和第m组分的二元交互系数，k_i，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物。

2.如权利要求1所述的内部热耦合精馏塔平衡级节能控制系统，其特征在于：所述的上位机还包括：

结果显示模块，用于将节能控制结果传给控制站进行显示，并通过现场总线将节能控制结果传递到现场操作站进行显示。

3.一种用如权利要求1所述的内部热耦合精馏塔平衡级节能控制系统实现的节能控制方法，其特征在于：所述的节能控制方法包括以下步骤：

1)设定塔的结构参数，采集塔的生产工况数据，采集当前进料流量作为初值；

2)假定各塔板液相组成；

3)对每一个塔板，分别由泡点法计算其平衡温度和汽相组成；

4)对每一个塔板，计算其汽液相的焓值；

5)联立式(1)(2)计算各塔板的汽液相流量：

$V_{j+1}{H}_{j+1}^G+U_{j-1}{H}_{j-1}^L+F_j{H}_j^F-(V_j+S_j^G)H_j^G-(U_j+S_j^L)H_j^L-Q_j=0---\left(1\right)$

$V_{j+1}+U_{j-1}+F_j^G+F_j^L-(V_j+S_j^G)-(U_j+S_j^L)=0---\left(2\right)$

其中，V表示汽相流量，U表示液相流量，F表示进料流量，H^F表示进料焓值，S表示侧提流量，下标j-1、j、j+1分别表示第j-1、j、j+1块板，Q表示热耦合量，由下式计算：

Q＝UAΔT                                               (3)；

6)判断式(4)是否成立，如果成立，则继续步骤7)，否则，更新液相组成，返回步骤3)迭代；

$V_{j+1}{y}_{i,j+1}+U_{j-1}{x}_{i,j-1}+F_j^G{z}_{i,j}^G+F_j^L{z}_{i,j}^L-(V_j+S_j^G)y_{i,j}-(U_j+S_j^L)x_{i,j}<0.0001---\left(4\right)$

7)判断产品纯度是否满足约束，如果不满足则结束迭代，输出结果，前一步的进料流量即为最大进料量，如果满足则将进料流量增加一个迭代步长Δ，返回步骤2)继续迭代；

所述步骤3)中，由泡点法计算平衡温度和汽相组成的过程如下：

3.1)假定塔板平衡温度；

3.2)计算汽液平衡常数，采用以下过程完成：

$\ln {\mathrm{\&Phi;}}_i^L=\ln \frac{\mathrm{RT}}{P(v^L-b^L)}-\frac{b_i}{b^L}(1-Z^L)+{\mathrm{\&xi;}}^L{a}^L(\frac{b_i}{b^L}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^L})/b^L\mathrm{RT}---\left(5\right)$

$\ln {\mathrm{\&Phi;}}_i^G=\ln \frac{\mathrm{RT}}{P(v^G-b^G)}-\frac{b_i}{b^G}(1-Z^G)+{\mathrm{\&xi;}}^G{a}^G(\frac{b_i}{b^G}-\frac{2\underset{m}{\mathrm{\&Sigma;}}{x}_m{a}_{i,m}}{a^G})/b^G\mathrm{RT}---\left(6\right)$

$K_i={\mathrm{\&Phi;}}_i^L/{\mathrm{\&Phi;}}_i^G---\left(7\right)$

y_i＝K_ix_i                                               (8)

其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1，...M表示组分，M表示组分数，摩尔体积v、物性参数b^G、b^L、b_i、a^G、a^L、a_i，m、ξ^G、ξ^L、汽相压缩因子Z^G、液相压缩因子Z^L由物性参数计算单元计算；

3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回步骤3.2)继续迭代；

所述步骤4)中，计算汽液相混合焓，其过程如下：

$H_i^{*}=c_i+d_{i}T+e_i{T}^2+f_i{T}^3+h_i{T}^4---\left(9\right)$

$H^{*}=\underset{i}{\mathrm{\&Sigma;}}{y}_i{H}_i^{*}---\left(10\right)$

$H^G=H^{*}-\mathrm{RT}(1-Z^G)-{\mathrm{\&xi;}}^G(a^G-T\frac{{\mathrm{da}}^G}{\mathrm{dT}}){/b}^G---\left(11\right)$

$H^L=H^{*}-\mathrm{RT}(1-Z^L)-{\mathrm{\&xi;}}^L(a^L-T\frac{{\mathrm{da}}^L}{\mathrm{dT}})/b^L---\left(12\right)$

其中

表示第i个纯组分理想气体的焓值，H^*是混合物理想气体焓值，c、d、e、f、h为常数；

所述物性参数计算方法包括以下步骤：

$a_{i,m}={\mathrm{\&Omega;}}_{\mathrm{ai},m}{R}^2{T}_{\mathrm{ci},m}^2/P_{\mathrm{ci},m}---\left(13\right)$

b_i＝Ω_bRT_ci/P_cia                                       (14)

$T_{\mathrm{ci},m}=\sqrt{T_{\mathrm{ci}}{T}_{\mathrm{cj}}}(1-k_{i,m})---\left(15\right)$

$V_{\mathrm{ci},m}=0.125{(V_{\mathrm{ci},m}^{1/3}+V_{\mathrm{ci},m}^{1/3})}^3---\left(16\right)$

Z_ci，m＝0.5(Z_ci+Z_cm)                                   (17)

P_ci，m＝RT_ci，mZ_ci，m/V_ci，m                           (18)

Ω_ai，m＝0.5(Ω_ai+Ω_am)                                (19)

对汽相： $a^G=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{y}_i{y}_m{a}_{i,m}---\left(20\right)$

$b^G=\underset{i}{\mathrm{\&Sigma;}}{y}_i{b}_i---\left(21\right)$

令

A^G＝a^GP/R²T²                                           (22)

B^G＝b^GP/RT                                             (23)

α^G＝2B^G-1                                             (24)

${\mathrm{\&beta;}}^G=A^G-3B^G-5B^{G^2}---\left(25\right)$

${\mathrm{\&gamma;}}^G=2(B^{G^3}+B^{G^2})-A^G{B}^G---\left(26\right)$

取初值为1-0.6P_r，用牛顿法解如下方程，即得到汽相压缩因子Z^G

$Z^{G^3}+{\mathrm{\&alpha;}}^G{Z}^{G^2}+{\mathrm{\&beta;Z}}^G+{\mathrm{\&gamma;}}^G=0---\left(27\right)$

则，

v^G＝RT/PZ^G                                             (28)

${\mathrm{\&xi;}}^G=0.242536\ln \frac{v^G+3.561553b^G}{v^G-0.561553b^G}---\left(29\right)$

对液相： $a^L=\underset{i}{\mathrm{\&Sigma;}}\underset{m}{\mathrm{\&Sigma;}}{x}_i{x}_m{a}_{i,m}---\left(30\right)$

$b^L=\underset{i}{\mathrm{\&Sigma;}}{x}_i{b}_i---\left(31\right)$

令

A^L＝a^LP/R²T²                                           (32)

B^L＝b^LP/RT                                             (33)

α^L＝2B^L-1                                             (34)

${\mathrm{\&beta;}}^L=A^L-3B^L-5B^{L^2}---\left(35\right)$

${\mathrm{\&gamma;}}^L=2(B^{L^3}+B^{L^2})-A^L{B}^L---\left(36\right)$

取初值为P_r(0.106+0.078P_r)，用牛顿法解如下方程，即得到液相压缩因子Z^L

$Z^{L^3}+{\mathrm{\&alpha;}}^L{Z}^{L^2}+{\mathrm{\&beta;Z}}^L+{\mathrm{\&gamma;}}^L=0---\left(37\right)$

则，

v^L＝RT/PZ^L                                             (38)

${\mathrm{\&xi;}}^L=0.242536\ln \frac{v^L+3.561553b^L}{v^L-0.561553b^L}---\left(39\right)$

Ω_ai＝C_i-D_iτ+E_iτ²-W_iτ³                              (40)

Ω_b＝0.070721                                          (41)

τ＝0.01T                                              (42)

其中，A、B、α、β、γ、τ、Ω_a、Ω_b是中间变量，C、D、E、W是常数，T_c、P_c、V_c、Z_c分别是临界温度、压力、体积和压缩因子，P_r是对比压力，R是气体常数，k_i，m表示第i组分和第m组分的二元交互系数，k_i，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物。

4.如权利要求3所述的节能控制方法，其特征在于：在所述的步骤7)中，上位机将节能控制结果传给控制站进行显示，并通过现场总线将节能控制结果传递到现场操作站进行显示。

Images