内部热耦合空分塔动态流程模拟系统及方法专利检索- .数学公式建模图表或曲线;设计方法专利检索查询-专利查询网

内部热耦合空分塔动态流程模拟系统及方法
申请号	CN200910157179.0	申请日	2009-12-23	公开(公告)号	CN101716427A	公开(公告)日	2010-06-02
申请人	浙江大学;			发明人	刘兴高; 闫正兵;
摘要	一种内部热耦合空分塔动态流程模拟系统，包括与内部热耦合空分塔连接的现场智能仪表、以及控制站、数据库和上位机，上位机包括信号采集模块，用以采集当前生产工况数据；求解计算主模块，过程：设定塔的结构参数和操作参数，设定起始时刻tstart，终止时刻tend；指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart；对每一个塔板，分别由泡点法计算其平衡温度和汽相组成、汽液相的焓值；计算(t+Δt)时刻的各塔板液相组成和液相流量；令t＝t+Δt，用新的各塔板液相组成和液相流量返回迭代，直到t≥tend，结束迭代，输出结果。以及提出了一种内部热耦合空分塔动态流程模拟方法。本发明提供一种能够准确模拟内部热耦合空分塔动态流程的系统及方法。
权利要求	1.一种内部热耦合空分塔动态流程模拟系统，包括与内部热耦合空分塔连接的现场智能仪表、以及控制站、数据库和上位机，智能仪表与控制站、数据库和上位机连接，其特征在于：所述的上位机包括：信号采集模块，用以采集当前生产工况数据；求解计算主模块，用以求解计算，采用以下过程来完成： 1)设定塔的结构参数和操作参数，设定起始时刻tstart，终止时刻tend； 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3) 其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{{dx}_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{{dx}_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{{dU}_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回3)迭代，直到t≥tend，结束迭代，输出结果。 2.如权利要求1所述的内部热耦合空分塔动态流程模拟系统，其特征在于：所述上位机还包括：泡点法模块，用以山泡点法计算其平衡温度和汽相组成，其过程如下： 3.1)假定塔板平衡温度； 3.2)计算汽液平衡常数，采用以下过程完成： $\ln Φ_{i}^{L} = \ln \frac{RT}{P (v^{L} - b^{L})} - \frac{b_{i}}{b^{L}} (1 - z^{L}) + ξ^{L} a^{L} (\frac{b_{i}}{b^{L}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{L}}) / b^{L} RT - - - (9)$ $\ln Φ_{i}^{G} = \ln \frac{RT}{P (v^{G} - b^{G})} - \frac{b_{i}}{b^{G}} (1 - z^{G}) + ξ^{G} a^{G} (\frac{b_{i}}{b^{G}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{G}}) / b^{G} RT - - - (10)$ $K_{i} = Φ_{i}^{L} / Φ_{i}^{G} - - - (11)$ yi＝Kixi (12) 其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1、2、3表示组分，依次对应氮、氩、氧，摩尔体积v、参数bG、bL、bi、aG、aL、ai，m、ξG、ξL、汽相压缩因子ZG、液相压缩因子ZL由物性模块计算； 3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代。 3.如权利要求1或2所述的内部热耦合空分塔动态流程模拟系统，其特征在于：所述上位机还包括：焓模块，用以计算汽液相混合焓，其过程如下： $H_{i}^{} = c_{i} + d_{i} T + e_{i} T^{2} + f_{i} T^{3} + h_{i} T^{4} - - - (13)$ $H^{} = \underset{i}{Σ} y_{i} H_{i}^{} - - - (14)$ $H^{G} = H^{} - RT (1 - Z^{G}) - ξ^{G} (a^{G} - T \frac{d a^{G}}{dT}) / b^{G} - - - (15)$ $H^{L} = H^{} - RT (1 - Z^{L}) - ξ^{L} (a^{L} - T \frac{d a^{L}}{dT}) / b^{L} - - - (16)$ 其中Hi表示第i个纯组分理想气体的焓值，H是混合物理想气体焓值，c、d、e、f、h为常数。 4.如权利要求2所述的内部热耦合空分塔动态流程模拟系统，其特征在于：所述上位机还包括：物性模块，用以计算物性参数，其过程如下： $a_{i, m} = Ω_{ai, m} R^{2} T_{ci, m}^{2} / P_{ci, m} - - - (17)$ bi＝ΩbRTci/Pcia (18) $T_{ci, m} = \sqrt{T_{ci} T_{cj}} (1 - k_{i, m}) - - - (19)$ $V_{ci, m} = 0.125 {(V_{ci, m}^{1 / 3} + V_{ci, m}^{1 / 3})}^{3} - - - (20)$ Zci，m＝0.5(Zci+Zcm) (21) Pci，m＝RTci，mZci，m/Vci，m (22) Ωai，m＝0.5(Ωai+Ωam) (23) 对汽相： $a^{G} = \underset{i}{Σ} \underset{m}{Σ} y_{i} y_{m} a_{i, m} - - - (24)$ $b^{G} = \underset{i}{Σ} y_{i} b_{i} - - - (25)$ 令 AG＝aGP/R2T2 (26) BG＝bGP/RT (27) αG＝2BG-1 (28) $β^{G} = A^{G} - 3 B^{G} - 5 B^{G^{2}} - - - (29)$ $γ^{G} = 2 (B^{G^{3}} + B^{G^{2}}) - A^{G} B^{G} - - - (30)$ 取初值为1-0.6Pr，用牛顿法解如下方程，即得到汽相压缩因子ZG $Z^{G^{3}} + α^{G} Z^{G^{2}} + β Z^{G} + γ^{G} = 0 - - - (31)$ 则， vG＝RT/PZG (32) $ξ^{G} = 0.242536 \ln \frac{v^{G} + {3.561553 b}^{G}}{v^{G} - {0.561553 b}^{G}} - - - (33)$ 对液相： $a^{L} = \underset{i}{Σ} \underset{m}{Σ} x_{i} x_{m} a_{i, m} - - - (34)$ $b^{L} = \underset{i}{Σ} x_{i} b_{i} - - - (35)$ 令 AL＝aLP/R2T2 (36) BL＝bLP/RT (37) αL＝2BL-1 (38) $β^{L} = A^{L} - 3 B^{L} - 5 B^{L^{2}} - - - (39)$ $γ^{L} = 2 (B^{L^{3}} + B^{L^{2}}) - A^{L} B^{L} - - - (40)$ 取初值为Pr(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子ZL $Z^{L^{3}} + α^{L} Z^{L^{2}} + β Z^{L} + γ^{L} = 0 - - - (41)$ 则， vL＝RT/PZL (42) $ξ^{L} = 0.242536 \ln \frac{v^{L} + {3.561553 b}^{L}}{v^{L} - {0.561553 b}^{L}} - - - (43)$ Ωai＝Ci-Diτ+Eiτ2-Wiτ3 (44) Ωb＝0.070721 (45) τ＝0.01T (46) 其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，Tc、Pc、Vc、Zc分别是临界温度、压力、体积和压缩因子，Pr是对比压力，R是气体常数，ki，m表示第i组分和第m组分的二元交互系数，ki，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物，Ωa、Ωb是中间变量。 5.如权利要求1或2所述的内部热耦合空分塔动态流程模拟系统，其特征在于：所述的上位机还包括：结果显示模块，用于将计算结果传给控制站进行显示，并通过现场总线将计算结果传递到现场操作站进行显示。 6.一种用如权利要求1所述的内部热耦合空分塔动态流程模拟系统实现的流程模拟方法，其特征在于：所述的流程模拟方法包括以下步骤： 1)设定塔的结构参数，采集生产工况数据，设定起始时刻tstart，终止时刻tend； 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3)；其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{d x_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{{dx}_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{{dU}_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回步骤3)迭代，直到t≥tend，结束迭代，输出结果。 7.如权利要求6所述的流程模拟方法，其特征在于：所述的步骤3)中，泡点法计算其平衡温度和汽相，采用以下过程完成： 3.1)假定塔板平衡温度； 3.2)计算汽液平衡常数，采用以下过程完成： $\ln Φ_{i}^{L} = \ln \frac{RT}{P (v^{L} - b^{L})} - \frac{b_{i}}{b^{L}} (1 - z^{L}) + ξ^{L} a^{L} (\frac{b_{i}}{b^{L}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{L}}) / b^{L} RT - - - (9)$ $\ln Φ_{i}^{G} = \ln \frac{RT}{P (v^{G} - b^{G})} - \frac{b_{i}}{b^{G}} (1 - z^{G}) + ξ^{G} a^{G} (\frac{b_{i}}{b^{G}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{G}}) / b^{G} RT - - - (10)$ $K_{i} = Φ_{i}^{L} / Φ_{i}^{G} - - - (11)$ yi＝Kixi (12) 其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1、2、3表示组分，依次对应氮、氩、氧，摩尔体积v、参数bG、bL、bi、aG、aL、ai，m、ξG、ξL、汽相压缩因子ZG、液相压缩因子ZL由物性参数计算方法计算； 3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代。 8.如权利要求6或7所述的流程模拟方法，其特征在于：所述步骤4)中，所述的焓计算方法过程如下： $H_{i}^{} = c_{i} + d_{i} T + e_{i} T^{2} + f_{i} T^{3} + h_{i} T^{4} - - - (13)$ $H^{} = \underset{i}{Σ} y_{i} H_{i}^{} - - - (14)$ $H^{G} = H^{} - RT (1 - Z^{G}) - ξ^{G} (a^{G} - T \frac{d a^{G}}{dT}) / b^{G} - - - (15)$ $H^{L} = H^{} - RT (1 - Z^{L}) - ξ^{L} (a^{L} - T \frac{d a^{L}}{dT}) / b^{L} - - - (16)$ 其中Hi表示第i个纯组分理想气体的焓值，H是混合物理想气体焓值，c、d、e、f、h为常数。 9.如权利要求7所述的流程模拟方法，其特征在于：所述物性参数计算方法的过程如下： $a_{i, m} = Ω_{ai, m} R^{2} T_{ci, m}^{2} / P_{ci, m} - - - (17)$ bi＝ΩbRTci/Pcia (18) $T_{ci, m} = \sqrt{T_{ci} T_{cj}} (1 - k_{i, m}) - - - (19)$ $V_{ci, m} = 0.125 {(V_{ci, m}^{1 / 3} + V_{ci, m}^{1 / 3})}^{3} - - - (20)$ Zci，m＝0.5(Zci+Zcm) (21) Pci，m＝RTci，mZci，m/Vci，m (22) Ωai，m＝0.5(Ωai+Ωam) (23) 对汽相： $a^{G} = \underset{i}{Σ} \underset{m}{Σ} y_{i} y_{m} a_{i, m} - - - (24)$ $b^{G} = \underset{i}{Σ} y_{i} b_{i} - - - (25)$ 令 AG＝aGP/R2T2 (26) BG＝bGP/RT (27) αG＝2BG-1 (28) $β^{G} = A^{G} - 3 B^{G} - 5 B^{G^{2}} - - - (29)$ $γ^{G} = 2 (B^{G^{3}} + B^{G^{2}}) - A^{G} B^{G} - - - (30)$ 取初值为1-0.6Pr，用牛顿法解如下方程，即得到汽相压缩因子ZG $Z^{G^{3}} + α^{G} Z^{G^{2}} + β Z^{G} + γ^{G} = 0 - - - (31)$ 则， vG＝RT/PZG (32) $ξ^{G} = 0.242536 \ln \frac{v^{G} + {3.561553 b}^{G}}{v^{G} - {0.561553 b}^{G}} - - - (33)$ 对液相： $a^{L} = \underset{i}{Σ} \underset{m}{Σ} x_{i} x_{m} a_{i, m} - - - (34)$ $b^{L} = \underset{i}{Σ} x_{i} b_{i} - - - (35)$ 令 AL＝aLP/R2T2 (36) BL＝bLP/RTv (37) αL＝2BL-1 (38) $β^{L} = A^{L} - 3 B^{L} - 5 B^{L^{2}} - - - (39)$ $γ^{L} = 2 (B^{L^{3}} + B^{L^{2}}) - A^{L} B^{L} - - - (40)$ 取初值为Pr(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子ZL $Z^{L^{3}} + α^{L} Z^{L^{2}} + β Z^{L} + γ^{L} = 0 - - - (41)$ 则， vL＝RT/PZL (42) $ξ^{L} = 0.242536 \ln \frac{v^{L} + {3.561553 b}^{L}}{v^{L} - {0.561553 b}^{L}} - - - (43)$ Ωai＝Ci-Diτ+Eiτ2-Wiτ3 (44) Ωb＝0.070721 (45) τ＝0.01T (46) 其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，Tc、Pc、Vc、Zc分别是临界温度、压力、体积和压缩因子，Pr是对比压力，R是气体常数，ki，m表示第i组分和第m组分的二元交互系数，ki，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物，Ωa、Ωb是中间变量。 10.如权利要求6或7所述的流程模拟方法，其特征在于：在所述的步骤7)中，上位机将计算结果传给控制站进行显示，并通过现场总线将计算结果传递到现场操作站进行显示。
说明书全文	技术领域本发明涉及空分领域，特别地，涉及一种内部热耦合空分塔动态流程模拟系统及方法。背景技术氧气、氮气和氩气的应用十分广泛。低温空气分离法是当前国内外空气分离行业应用最为广泛的方法。在空气分离行业中，能源成本占了空气产品价格的75％。在能源危机不断加深的形势下，提高空气分离技术的能量效率具有重要的社会和经济意义。内部热耦合精馏技术最初是至今为止的四大节能精馏技术中节能效能最高的一项节能技术，在世界范围内得到了巨大的重视。在低温空分过程中，空分塔是一个重要的操作单元，也是最重要的耗能单元。通过对内部热耦合技术的研究发现：该技术的实现，需要过程装置中具有不同高低压强的区域，而常规空分装置本身的工艺流程就给出了这样区域；研究表明内部热耦合技术在高纯度产品生产过程上将具有更大的经济效益，而空分过程恰好是一个产品纯度要求极高的过程；内部热耦合精馏技术在低操作温度情况下更加有利于发挥它的节能潜力，而常规空分装置本身就是低温精馏；内部热耦合精馏技术在处理沸点接近的物系时，能体现出很高的节能效率，而空气中氮--氩--氧三元物系又恰好满足了这个要求，也为热耦合技术的应用开发提高了非常有利条件。因此，内部热耦合技术可以被应用到空分过程，改变传统空分塔结构，达到良好的节能效果。化学工程的传统研究方法是以经验归纳为主，即运用实验手段探索系统变量之间的关系。随着生产过程的大型化和自动化水平的提高，对于复杂化工过程的研究，因次分析和相似方法往往不能满足研究的需要。内部热耦合空分塔过程是包含有传热、传质、流体流动的复杂过程，许多变量之间呈现强烈的非线性关系，要了解过程中浓度分布、速度分布等随时间的变化规律，一个有效的方法就是利用计算机进行动态流程模拟。通过内部热耦合空分塔动态流程模拟分析空分塔的动态特性和响应，是流程设计和控制的先决条件。发明内容为了克服现有的空分塔的内部热耦合精馏过程的不能准确模拟内部热耦合空分塔动态流程的不足，本发明提供一种能够准确模拟内部热耦合空分塔动态流程的内部热耦合空分塔动态流程模拟系统及方法。本发明解决其技术问题所采用的技术方案是：一种内部热耦合空分塔动态流程模拟系统，包括与内部热耦合空分塔连接的现场智能仪表、以及控制站、数据库和上位机，智能仪表与控制站、数据库和上位机连接，所述的上位机包括：信号采集模块，用以采集当前生产工况数据；求解计算主模块，用以求解计算，采用以下过程来完成： 1)设定塔的结构参数和操作参数，设定起始时刻tstart，终止时刻tend； 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3) 其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{{dx}_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{d x_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{d U_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回3)迭代，直到t≥tend，结束迭代，输出结果。作为优选的一种方案：所述上位机还包括：泡点法模块，用以由泡点法计算其平衡温度和汽相组成，其过程如下： 3.1)假定塔板平衡温度； 3.2)计算汽液平衡常数，采用以下过程完成： $\ln Φ_{i}^{L} = \ln \frac{RT}{P (v^{L} - b^{L})} - \frac{b_{i}}{b^{L}} (1 - Z^{L}) + ξ^{L} a^{L} (\frac{b_{i}}{b^{L}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{L}}) / b^{L} RT - - - (9)$ $\ln Φ_{i}^{G} = \ln \frac{RT}{P (v^{G} - b^{G})} - \frac{b_{i}}{b^{G}} (1 - Z^{G}) + ξ^{G} a^{G} (\frac{b_{i}}{b^{G}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{G}}) / b^{G} RT - - - (10)$ $K_{i} = Φ_{i}^{L} / Φ_{i}^{G} - - - (11)$ yi＝Kixi (12) 其中，Φ表示逸度系数，上标G表示汽相，Q表示热耦合量，R是气体常数，T是温度，P是塔板压强，下标m＝1、2、3表示组分，依次对应氮、氩、氧，摩尔体积v、参数bG、bL、bi、aG、aL、ai，m、ξG、ξL、汽相压缩因子ZG、液相压缩因子ZL由物性模块计算； 3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代。作为优选的另一种方案：所述上位机还包括：焓模块，用以计算汽液相混合焓，其过程如下： $H_{i}^{} = c_{i} + d_{i} T + e_{i} T^{2} + f_{i} T^{3} + h_{i} T^{4} - - - (13)$ $H^{} = \underset{i}{Σ} y_{i} H_{i}^{} - - - (14)$ $H^{G} = H^{} - RT (1 - Z^{G}) - ξ^{G} (a^{G} - T \frac{d a^{G}}{dT}) / b^{G} - - - (15)$ $H^{L} = H^{} - RT (1 - Z^{L}) - ξ^{L} (a^{L} - T \frac{d a^{L}}{dT}) / b^{L} - - - (16)$ 其中Hi表示第i个纯组分理想气体的焓值，H是混合物理想气体焓值，c、d、e、f、h为常数。作为优选的再一种方案：所述上位机还包括：物性模块，用以计算物性参数，其过程如下： $a_{i, m} = Ω_{ai, m} R^{2} T_{ci, m}^{2} / P_{ci, m} - - - (17)$ bi＝ΩbRTci/Pcia (18) $T_{ci, m} = \sqrt{T_{ci} T_{cj}} (1 - k_{i, m}) - - - (19)$ $V_{ci, m} = 0.125 {(V_{ci, m}^{1 / 3} + V_{ci, m}^{1 / 3})}^{3} - - - (20)$ Zci，m＝0.5(Zci+Zcm) (21) Pci，m＝RTci，mZci，m/Vci，m (22) Ωai，m＝0.5(Ωai+Ωam) (23) 对汽相： $a^{G} = \underset{i}{Σ} \underset{m}{Σ} y_{i} y_{m} a_{i, m} - - - (24)$ $b^{G} = \underset{i}{Σ} y_{i} b_{i} - - - (25)$ 令 AG＝aGP/R2T2 (26) BG＝bGP/RT (27) αG＝2BG-1 (28) $β^{G} = A^{G} - 3 B^{G} - 5 B^{G^{2}} - - - (29)$ $γ^{G} = 2 (B^{G^{3}} + B^{G^{2}}) - A^{G} B^{G} - - - (30)$ 取初值为1-0.6Pr，用牛顿法解如下方程，即得到汽相压缩因子ZG $Z^{G^{3}} + α^{G} Z^{G^{2}} + β Z^{G} + γ^{G} = 0 - - - (31)$ 则， vG＝RT/PZG (32) $ξ^{G} = 0.242536 \ln \frac{v^{G} + 3.561553 b^{G}}{v^{G} - 0.561553 b^{G}} - - - (33)$ 对液相： $a^{L} = \underset{i}{Σ} \underset{m}{Σ} x_{i} x_{m} a_{i, m} - - - (34)$ $b^{L} = \underset{i}{Σ} x_{i} b_{i} - - - (35)$ 令 AL＝aLP/R2T2 (36) BL＝bLP/RT (37) αL＝2BL-1 (38) $β^{L} = A^{L} - 3 B^{L} - 5 B^{L^{2}} - - - (39)$ $γ^{L} = 2 (B^{L^{3}} + B^{L^{2}}) - A^{L} B^{L} - - - (40)$ 取初值为Pr(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子ZL $Z^{L^{3}} + α^{L} Z^{L^{2}} + β Z^{L} + γ^{L} = 0 - - - (41)$ 则， vL＝RT/PZL (42) $ξ^{L} = 0.242536 \ln \frac{v^{L} + 3.561553 b^{L}}{v^{L} - 0.561553 b^{L}} - - - (43)$ Ωai＝Ci-Diτ+Eiτ2-Wiτ3 (44) Ωb＝0.070721 (45) τ＝0.01T (46) 其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，Tc、Pc、Vc、Zc分别是临界温度、压力、体积和压缩因子，Pr是对比压力，R是气体常数，ki，m表示第i组分和第m组分的二元交互系数，ki，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物，Ωa、Ωb是中间变量。进一步，所述的上位机还包括：结果显示模块，用于将计算结果传给控制站进行显示，并通过现场总线将计算结果传递到现场操作站进行显示。一种内部热耦合空分塔动态流程模拟方法，所述的流程模拟方法包括以下步骤： 1)设定塔的结构参数，采集生产工况数据，设定起始时刻tstart，终止时刻tend； 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3)；其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{{dx}_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{d x_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{d U_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回步骤3)迭代，直到t≥tend，结束迭代，输出结果。作为优选的一种方案：所述的步骤3)中，泡点法计算其平衡温度和汽相，采用以下过程完成： 3.1)假定塔板平衡温度； 3.2)计算汽液平衡常数，采用以下过程完成： $\ln Φ_{i}^{L} = \ln \frac{RT}{P (v^{L} - b^{L})} - \frac{b_{i}}{b^{L}} (1 - Z^{L}) + ξ^{L} a^{L} (\frac{b_{i}}{b^{L}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{L}}) / b^{L} RT - - - (9)$ $\ln Φ_{i}^{G} = \ln \frac{RT}{P (v^{G} - b^{G})} - \frac{b_{i}}{b^{G}} (1 - Z^{G}) + ξ^{G} a^{G} (\frac{b_{i}}{b^{G}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{G}}) / b^{G} RT - - - (10)$ $K_{i} = Φ_{i}^{L} / Φ_{i}^{G} - - - (11)$ yi＝Kixi (12) 其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1、2、3表示组分，依次对应氮、氩、氧，摩尔体积v、参数bG、bL、bi、aG、aL、ai，m、ξG、ξL、汽相压缩因子ZG、液相压缩因子ZL由物性参数计算方法计算； 3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代。作为优选的另一种方案：所述步骤4)中，所述的焓计算方法过程如下： $H_{i}^{} = c_{i} + d_{i} T + e_{i} T^{2} + f_{i} T^{3} + h_{i} T^{4} - - - (13)$ $H^{} = \underset{i}{Σ} y_{i} H_{i}^{} - - - (14)$ $H^{G} = H^{} - RT (1 - Z^{G}) - ξ^{G} (a^{G} - T \frac{d a^{G}}{dT}) / b^{G} - - - (15)$ $H^{L} = H^{} - RT (1 - Z^{L}) - ξ^{L} (a^{L} - T \frac{d a^{L}}{dT}) / b^{L} - - - (16)$ 其中Hi表示第i个纯组分理想气体的焓值，H是混合物理想气体焓值，c、d、e、f、h为常数。作为优选的再一种方案：所述物性参数计算方法的过程如下： $a_{i, m} = Ω_{ai, m} R^{2} T_{ci, m}^{2} / P_{ci, m} - - - (17)$ bi＝ΩbRTci/Pcia (18) $T_{ci, m} = \sqrt{T_{ci} T_{cj}} (1 - k_{i, m}) - - - (19)$ $V_{ci, m} = 0.125 {(V_{ci, m}^{1 / 3} + V_{ci, m}^{1 / 3})}^{3} - - - (20)$ Zci，m＝0.5(Zci+Zcm) (21) Pci，m＝RTci，mZci，m/Vci，m (22) Ωai，m＝0.5(Ωai+Ωam) (23) 对汽相： $a^{G} = \underset{i}{Σ} \underset{m}{Σ} y_{i} y_{m} a_{i, m} - - - (24)$ $b^{G} = \underset{i}{Σ} y_{i} b_{i} - - - (25)$ 令 AG＝aGP/R2T2 (26) BG＝bGP/RT (27) αG＝2BG-1 (28) $β^{G} = A^{G} - 3 B^{G} - 5 B^{G^{2}} - - - (29)$ $γ^{G} = 2 (B^{G^{3}} + B^{G^{2}}) - A^{G} B^{G} - - - (30)$ 取初值为1-0.6Pr，用牛顿法解如下方程，即得到汽相压缩因子ZG $Z^{G^{3}} + α^{G} Z^{G^{2}} + β Z^{G} + γ^{G} = 0 - - - (31)$ 则， vG＝RT/PZG (32) $ξ^{G} = 0.242536 \ln \frac{v^{G} + 3.561553 b^{G}}{v^{G} - 0.561553 b^{G}} - - - (33)$ 对液相： $a^{L} = \underset{i}{Σ} \underset{m}{Σ} x_{i} x_{m} a_{i, m} - - - (34)$ $b^{L} = \underset{i}{Σ} x_{i} b_{i} - - - (35)$ 令 AL＝aLP/R2T2 (36) BL＝bLP/RT (37) αL＝2BL-1 (38) $β^{L} = A^{L} - 3 B^{L} - 5 B^{L^{2}} - - - (39)$ $γ^{L} = 2 (B^{L^{3}} + B^{L^{2}}) - A^{L} B^{L} - - - (40)$ 取初值为Pr(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子ZL $Z^{L^{3}} + α^{L} Z^{L^{2}} + β Z^{L} + γ^{L} = 0 - - - (41)$ 则， vL＝RT/PZL (42) $ξ^{L} = 0.242536 \ln \frac{v^{L} + 3.561553 b^{L}}{v^{L} - 0.561553 b^{L}} - - - (43)$ Ωai＝Ci-Diτ+Eiτ2-Wiτ3 (44) Ωb＝0.070721 (45) τ＝0.01T (46) 其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，Tc、Pc、Vc、Zc分别是临界温度、压力、体积和压缩因子，Pr是对比压力，R是气体常数，ki，m表示第i组分和第m组分的二元交互系数，ki，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物，Ωa、Ωb是中间变量。进一步，在所述的步骤7)中，上位机将计算结果传给控制站进行显示，并通过现场总线将计算结果传递到现场操作站进行显示。本发明的有益效果主要表现在：对内部热耦合空分塔进行动态流程模拟，计算速度快，模拟结果准确，可以用于指导生产和进一步的优化、控制研究，从而提高生产效益。附图说明图1是本发明所提出的动态流程模拟系统的硬件结构图。图2是本发明所述内部热耦合空分塔结构示意图。图3是本发明上位机的功能模块图。具体实施方式下面结合附图对本发明作进一步描述。实施例1 参照图1、图2、图3，一种内部热耦合空分塔动态流程模拟系统，包括内部热耦合空分塔1连接的现场智能仪表2、数据接口3、控制站4、数据库5以及上位机6，智能仪表2连接现场总线，所述现场总线连接数据接口3，所述数据接口3与控制站4、数据库5和上位机6连接，所述的上位机6包括：信号采集模块7，用以采集当前生产工况数据；求解计算主模块8，用以求解计算，采用以下过程来完成： 1)设定塔的结构参数和操作参数，设定起始时刻tstart，终止时刻tend； 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3) 其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{{dx}_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{d x_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{d U_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回3)迭代，直到t≥tend，结束迭代，输出结果。所述上位机6还包括：泡点法模块9，用以由泡点法计算其平衡温度和汽相组成，其过程如下： 3.1)假定塔板平衡温度； 3.2)计算汽液平衡常数，采用以下过程完成： $\ln Φ_{i}^{L} = \ln \frac{RT}{P (v^{L} - b^{L})} - \frac{b_{i}}{b^{L}} (1 - Z^{L}) + ξ^{L} a^{L} (\frac{b_{i}}{b^{L}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{L}}) / b^{L} RT - - - (9)$ $\ln Φ_{i}^{G} = \ln \frac{RT}{P (v^{G} - b^{G})} - \frac{b_{i}}{b^{G}} (1 - Z^{G}) + ξ^{G} a^{G} (\frac{b_{i}}{b^{G}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{G}}) / b^{G} RT - - - (10)$ $K_{i} = Φ_{i}^{L} / Φ_{i}^{G} - - - (11)$ yi＝Kixi (12) 其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1、2、3表示组分，依次对应氮、氩、氧，摩尔体积v、物性参数bG、bL、bi、aG、aL、ai，m、ξG、ξL、汽相压缩因子ZG、液相压缩因子ZL由物性模块计算； 3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代；所述上位机6还包括：焓模块10，用以计算汽液相混合焓，其过程如下： $H_{i}^{} = c_{i} + d_{i} T + e_{i} T^{2} + f_{i} T^{3} + h_{i} T^{4} - - - (13)$ $H^{} = \underset{i}{Σ} y_{i} H_{i}^{} - - - (14)$ $H^{G} = H^{} - RT (1 - Z^{G}) - ξ^{G} (a^{G} - T \frac{d a^{G}}{dT}) / b^{G} - - - (15)$ $H^{L} = H^{} - RT (1 - Z^{L}) - ξ^{L} (a^{L} - T \frac{d a^{L}}{dT}) / b^{L} - - - (16)$ 其中Hi表示第i个纯组分理想气体的焓值，H是混合物理想气体焓值，c、d、e、f、h为常数。所述上位机6还包括：物性模块11，用以计算物性参数，其过程如下： $a_{i, m} = Ω_{ai, m} R^{2} T_{ci, m}^{2} / P_{ci, m} - - - (17)$ bi＝ΩbRTci/Pcia (18) $T_{ci, m} = \sqrt{T_{ci} T_{cj}} (1 - k_{i, m}) - - - (19)$ $V_{ci, m} = 0.125 {(V_{ci, m}^{1 / 3} + V_{ci, m}^{1 / 3})}^{3} - - - (20)$ Zci，m＝0.5(Zci+Zcm) (21) Pci，m＝RTci，mZci，m/Vci，m (22) Ωai，m＝0.5(Ωai+Ωam) (23) 对汽相： $a^{G} = \underset{i}{Σ} \underset{m}{Σ} y_{i} y_{m} a_{i, m} - - - (24)$ $b^{G} = \underset{i}{Σ} y_{i} b_{i} - - - (25)$ 令 AG＝aGP/R2T2 (26) BG＝bGP/RT (27) αG＝2BG-1 (28) $β^{G} = A^{G} - 3 B^{G} - 5 B^{G^{2}} - - - (29)$ $γ^{G} = 2 (B^{G^{3}} + B^{G^{2}}) - A^{G} B^{G} - - - (30)$ 取初值为1-0.6Pr，用牛顿法解如下方程，即得到汽相压缩因子ZG $Z^{G^{3}} + α^{G} Z^{G^{2}} + β Z^{G} + γ^{G} = 0 - - - (31)$ 则， vG＝RT/PZG (32) $ξ^{G} = 0.242536 \ln \frac{v^{G} + 3.561553 b^{G}}{v^{G} - 0.561553 b^{G}} - - - (33)$ 对液相： $a^{L} = \underset{i}{Σ} \underset{m}{Σ} x_{i} x_{m} a_{i, m} - - - (34)$ $b^{L} = \underset{i}{Σ} x_{i} b_{i} - - - (35)$ 令 AL＝aLP/R2T2 (36) BL＝bLP/RT (37) αL＝2BL-1 (38) $β^{L} = A^{L} - 3 B^{L} - 5 B^{L^{2}} - - - (39)$ $γ^{L} = 2 (B^{L^{3}} + B^{L^{2}}) - A^{L} B^{L} - - - (40)$ 取初值为Pr(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子ZL $Z^{L^{3}} + α^{L} Z^{L^{2}} + β Z^{L} + γ^{L} = 0 - - - (41)$ 则， vL＝RT/PZL (42) $ξ^{L} = 0.242536 \ln \frac{v^{L} + 3.561553 b^{L}}{v^{L} - 0.561553 b^{L}} - - - (43)$ Ωai＝Ci-Diτ+Eiτ2-Wiτ3 (44) Ωb＝0.070721 (45) τ＝0.01T (46) 其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，Tc、Pc、Vc、Zc分别是临界温度、压力、体积和压缩因子，Pr是对比压力，R是气体常数，ki，m表示第i组分和第m组分的二元交互系数，ki，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)设定塔的结构参数，采集生产工况数据，设定起始时刻tstart，终止时刻tend： 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3)；其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{{dx}_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{d x_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{d U_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回3)迭代，直到t≥tend，结束迭代，输出结果. 实施例2 参照图1、图2、图3，一种内部热耦合空分塔动态流程模拟方法，所述的流程模拟方法包括以下步骤： 1)设定塔的结构参数，采集生产工况数据，设定起始时刻tstart，终止时刻tend； 2)指定初始时刻的各塔板液相组成和液相流量，令当前迭代时间t＝tstart； 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 - - - (1)$ $V_{j + 1} + U_{j - 1} + F_{j}^{G} + F_{j}^{L} - (V_{j} + S_{j}^{G}) - (U_{j} + S_{j}^{L}) = 0 - - - (2)$ 其中，V表示汽相流量，U表示液相流量，F表示进料流量，HF表示进料焓值，S表示侧提流量，HG和HL分别是汽液相焓值，下标j-1、j、j+1分别表示第j-1、j、j+1块板，上标L表示液相，上标G表示汽相，Q表示热耦合量，由下式计算： Q＝UAΔT (3)；其中，UA表示热耦合系数，ΔT表示耦合塔板间的温差； 6)计算(t+Δt)时刻的各塔板液相组成和液相流量： $\frac{{dx}_{i, j}}{dt} M_{j} = U_{j - 1} x_{i, j - 1} + V_{j + 1} y_{i, j + 1} + F_{j} z_{i, j} - (U_{j} + S_{j}^{L}) x_{i, j} - (V_{j} + S_{j}^{G}) y_{i, j} - - - - (4)$ $[U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] x_{i, j}$ $\frac{d U_{j}}{dt} = 518.13 ρ_{j}^{- 1 / 3} A_{a}^{- 1} l^{2 / 3} U_{j}^{1 / 3} [U_{j - 1} + V_{j + 1} + F_{j} - (U_{j} + S_{j}^{L}) - (V_{j} + S_{j}^{G})] - - - (5)$ $x_{i, j}^{(t + Δt)} = x_{i, j}^{(t)} + \frac{d x_{i, j}}{dt} Δt - - - (6)$ $U_{j}^{(t + Δt)} = U_{j}^{(t)} + \frac{d U_{j}}{dt} Δt - - - (7)$ 其中，Δt是迭代步长，x是液相组成，y是汽相组成，z是进料组成，ρ是液相密度，Aa是塔板有效面积，l是塔板堰长，下标i＝1、2、3表示组分，依次对应氮、氩、氧，上标(t)和(t+Δt)分别表示时刻t和t+Δt，M表示塔板持液量，由下式计算： $M_{j} = ρ_{j} A_{aj} [h_{wj} + 0.00284 {(\frac{U_{j}}{ρ_{j} l_{j}})}^{2 / 3}] - - - (8)$ 其中hw是溢流堰高； 7)令t＝t+Δt，用新的各塔板液相组成和液相流量返回步骤3)迭代，直到t≥tend，结束迭代，输出结果。所述的步骤3)中，泡点法计算其平衡温度和汽相，采用以下过程完成： 3.1)假定塔板平衡温度； 3.2)计算汽液平衡常数，采用以下过程完成： $\ln Φ_{i}^{L} = \ln \frac{RT}{P (v^{L} - b^{L})} - \frac{b_{i}}{b^{L}} (1 - Z^{L}) + ξ^{L} a^{L} (\frac{b_{i}}{b^{L}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{L}}) / b^{L} RT - - - (9)$ $\ln Φ_{i}^{G} = \ln \frac{RT}{P (v^{G} - b^{G})} - \frac{b_{i}}{b^{G}} (1 - Z^{G}) + ξ^{G} a^{G} (\frac{b_{i}}{b^{G}} - \frac{2 \underset{m}{Σ} x_{m} a_{i, m}}{a^{G}}) / b^{G} RT - - - (10)$ $K_{i} = Φ_{i}^{L} / Φ_{i}^{G} - - - (11)$ yi＝Kixi (12) 其中，Φ表示逸度系数，上标L表示液相，上标G表示汽相，R是气体常数，T是温度，P是塔板压强，下标m＝1、2、3表示组分，依次对应氮、氩、氧，摩尔体积v、参数bG、bL、bi、aG、aL、ai，m、ξG、ξL、汽相压缩因子ZG、液相压缩因子ZL由物性参数计算方法计算； 3.3)检验是否成立，成立则结束迭代，返回计算结果，否则，更新塔板平衡温度，返回3.2)继续迭代。所述步骤4)中，所述的焓计算方法过程如下： $H_{i}^{} = c_{i} + d_{i} T + e_{i} T^{2} + f_{i} T^{3} + h_{i} T^{4} - - - (13)$ $H^{} = \underset{i}{Σ} y_{i} H_{i}^{} - - - (14)$ $H^{G} = H^{} - RT (1 - Z^{G}) - ξ^{G} (a^{G} - T \frac{d a^{G}}{dT}) / b^{G} - - - (15)$ $H^{L} = H^{} - RT (1 - Z^{L}) - ξ^{L} (a^{L} - T \frac{d a^{L}}{dT}) / b^{L} - - - (16)$ 其中Hi表示第i个纯组分理想气体的焓值，H是混合物理想气体焓值，c、d、e、f、h为常数。所述物性参数计算方法的过程如下： $a_{i, m} = Ω_{ai, m} R^{2} T_{ci, m}^{2} / P_{ci, m} - - - (17)$ bi＝ΩbRTci/Pcia (18) $T_{ci, m} = \sqrt{T_{ci} T_{cj}} (1 - k_{i, m}) - - - (19)$ $V_{ci, m} = 0.125 {(V_{ci, m}^{1 / 3} + V_{ci, m}^{1 / 3})}^{3} - - - (20)$ Zci，m＝0.5(Zci+Zcm) (21) Pci，m＝RTci，mZci，m/Vci，m (22) Ωai，m＝0.5(Ωai+Ωam) (23) 对汽相： $a^{G} = \underset{i}{Σ} \underset{m}{Σ} y_{i} y_{m} a_{i, m} - - - (24)$ $b^{G} = \underset{i}{Σ} y_{i} b_{i} - - - (25)$ 令 AG＝aGP/R2T2 (26) BG＝bGP/RT (27) αG＝2BG-1 (28) $β^{G} = A^{G} - 3 B^{G} - 5 B^{G^{2}} - - - (29)$ $γ^{G} = 2 (B^{G^{3}} + B^{G^{2}}) - A^{G} B^{G} - - - (30)$ 取初值为1-0.6Pr，用牛顿法解如下方程，即得到汽相压缩因子ZG $Z^{G^{3}} + α^{G} Z^{G^{2}} + β Z^{G} + γ^{G} = 0 - - - (31)$ 则， vG＝RT/PZG (32) $ξ^{G} = 0.242536 \ln \frac{v^{G} + 3.561553 b^{G}}{v^{G} - 0.561553 b^{G}} - - - (33)$ 对液相： $a^{L} = \underset{i}{Σ} \underset{m}{Σ} x_{i} x_{m} a_{i, m} - - - (34)$ $b^{L} = \underset{i}{Σ} x_{i} b_{i} - - - (35)$ 令 AL＝aLP/R2T2 (36) BL＝bLP/RT (37) αL＝2BL-1 (38) $β^{L} = A^{L} - 3 B^{L} - 5 B^{L^{2}} - - - (39)$ $γ^{L} = 2 (B^{L^{3}} + B^{L^{2}}) - A^{L} B^{L} - - - (40)$ 取初值为Pr(0.106+0.078Pr)，用牛顿法解如下方程，即得到液相压缩因子ZL $Z^{L^{3}} + α^{L} Z^{L^{2}} + β Z^{L} + γ^{L} = 0 - - - (41)$ 则， vL＝RT/PZL (42) $ξ^{L} = 0.242536 \ln \frac{v^{L} + 3.561553 b^{L}}{v^{L} - 0.561553 b^{L}} - - - (43)$ Ωai＝Ci-Diτ+Eiτ2-Wiτ3 (44) Ωb＝0.070721 (45) τ＝0.01T (46) 其中，A、B、α、β、γ、τ是中间变量，C、D、E、W是常数，Tc、Pc、Vc、Zc分别是临界温度、压力、体积和压缩因子，Pr是对比压力，R是气体常数，ki，m表示第i组分和第m组分的二元交互系数，ki，m是常数，下标c表示临界点的性质，下标r表示对比态，下标i，m表示第i组分和第m组分的二元混合物，Ωa、Ωb是中间变量。在所述的步骤7)中，上位机将计算结果传给控制站进行显示，并通过现场总线将计算结果传递到现场操作站进行显示。本发明所提出的内部热耦合空分塔动态流程模拟系统及方法，已通过上述具体实施步骤进行了描述，相关技术人员明显能在不脱离本发明内容、精神和范围内对本文所述的装置和操作方法进行改动或适当变更与组合，来实现本发明技术。特别需要指出的是，所有相类似的替换和改动对本领域的技术人员是显而易见的，它们都会被视为包括在本发明精神、范围和内容中。

意见反馈