一種用于二維穩(wěn)態(tài)固定床的通用計算方法

吳斯侃,曹志凱,周 華,陳秉輝

(廈門大學化學化工學院,福建廈門361005)

一種用于二維穩(wěn)態(tài)固定床的通用計算方法

吳斯侃,曹志凱*,周 華,陳秉輝

(廈門大學化學化工學院,福建廈門361005)

針對固定床反應器模擬過程難以形成其模型完整通用的計算模塊且編程工作量大的問題,提出了一種穩(wěn)態(tài)下固定床反應器中反應物和產(chǎn)物溫度及濃度分布的通用計算方法.該方法以關鍵組分的摩爾分數(shù)作為變量,推導出關鍵組分在物料與熱量衡算中的偏微分方程組中的通式,并利用正交配置法求解.同時考慮固定床反應器內(nèi)氣相流速的變化,將其描述為關鍵組分摩爾分數(shù)的函數(shù)進一步實現(xiàn)模型的通用性.對于特定的反應,只需提供目標反應的動力學方程及操作條件即可求解.最后利用所提出的方法實現(xiàn)乙苯催化脫氫制苯乙烯與乙炔法制氯乙烯反應過程的模擬,結果表明該方法具有通用性和高效性.

固定床反應器;計算模擬;正交配置法;二維擬均相模型

凡是反應床層由固定器壁所固定,流體通過固體催化劑或隔板進行反應的反應器,統(tǒng)稱為固定床反應器[1],固定床反應器在化工領域應用非常廣泛.為實現(xiàn)其優(yōu)化操作,常對固定床反應器進行模擬.早在20世紀60年代,Forment[2]提出了固定床反應器的二維數(shù)學模型,并給出熱量和物料衡算式.直至今日,在文獻[2]的基礎上已報道了不少有關二維模型的模擬研究工作[3-5].

隨著模擬軟件的不斷改進[6],結合反應動力學來模擬固定床反應器以研究其內(nèi)部的溫度場和濃度場的二維分布成為可能.利用單管或中試數(shù)據(jù)實現(xiàn)反應動力學的擬合[7-9]也有報道,但對于固定床二維模型的通用計算方法,近年來進展不大[10-11].文獻[12]通過假設反應器內(nèi)部氣相流速恒定來簡化反應器模型的求解,該方法較易實現(xiàn)單一反應體系固定床反應器的模擬.由于假設反應器內(nèi)部氣相流速恒定,對流速或溫度變化較大的體系可能導致模擬結果存在一定誤差.為此本文中將模型中的流速轉化為摩爾分數(shù)的函數(shù),提出通用的固定床反應器二維數(shù)學模型建立的方法,成功地解決了動力學方程通用模塊編程難的問題,并利用正交配置法在MATLAB平臺上實現(xiàn)了固定床內(nèi)部溫度場和組分濃度的二維分布問題的求解.

1 模型的建立

1.1模型的基本假設

確定的反應體系可通過具體的反應方程式、反應組分、反應動力學及反應裝置的特點來建立確實可靠的反應模型[13].但對反應組分或反應體系不確定的抽象對象建立通用的模擬模塊具有一定的難度.為了研究的方便,可對所研究的對象提出一些符合客觀條件的假設.針對二維擬均相固定床反應器假設如下:1)不考慮反應器內(nèi)的壓降;2)反應器為軸對稱的圓柱;3)反應器內(nèi)氣相沿軸向流動且不存在軸向返混和導熱.

1.2物料和熱量衡算

為建立固定床反應過程的二維擬均相模型,首先需確定其衡算方程、本構關系和定解條件.由于本構關系涉及具體反應體系的動力學速率關系,因此不在本研究的討論范圍之內(nèi),本節(jié)主要介紹物料和熱量的衡算方程.

如圖1所示,在床層中環(huán)繞床層且距床壁徑向距離R′的中心線上取徑向長度為r、徑向厚度為dr、高為lr的微元環(huán)柱體.在該微元內(nèi)對反應中的關鍵組分k作物料衡算,得到固定床反應器中關鍵組分k的二維擬均相濃度分布模型如下:

式中:ck為第k個關鍵組分的濃度,mol/m3;Er為徑向有效擴散系數(shù),m2/s;ρB為催化劑堆積密度, kg/m3;u為流速,m/s;l為軸向高度,m;r為徑向長度,m;vjk為第k個關鍵組分的第j個獨立反應的化學計量系數(shù)為第j個獨立反應的普遍反應速率, mol/(kg·s);m為關鍵組分數(shù).

與物料衡算不同的是,熱量衡算不對反應中的關鍵組分k作熱量衡算,而是著眼于整個固定床床層溫度的變化.由圖1可得固定床反應器中二維擬均相溫度分布的微分方程:

式中:G為總質量流量,kg/(m2·s)；為混合物平均熱容,J/(kg·K);T為溫度,K;λer為擬均相模擬導熱系數(shù),W/(m·K);ΔHj為第j個獨立反應的反應熱,J/mol.

圖1　固定床中的環(huán)狀微元體積示意[14]Fig.1 Cyclic annular micro-unit volume schematic[14]

1.3邊界條件的建立

數(shù)學模型的定解條件為其解是否存在、是否唯一且穩(wěn)定的依據(jù).定解條件包括初始條件與邊界條件,而穩(wěn)態(tài)過程則只需考慮邊界條件.據(jù)1.1的假設,邊界條件可確定如下:

式中:r0為反應器半徑,m;T0為初始溫度,K;L為管長,m;hw為壁面?zhèn)鳠嵯禂?shù),W/(m2·K);Tw為壁面溫度,K.

2 模型參數(shù)的計算

在1.2中,利用衡算方程僅僅得到一般性的二維擬均相固定床反應器模型,而對于不同的反應體系,反應動力學方程所涉及的自變量不盡相同.為使模型具有更好的通用性,對原有模型中的參數(shù)和自變量進行替換以簡化計算和編程難度.

2.1模型通用化衡算方程

由式(1)和(2)可知,反應體系的物料衡算與熱量衡算的對象分別為關鍵組分的濃度和反應體系的溫度.反應體系不同其動力學的表達形式也各異[15],反應速率實際上可以描述為相關參數(shù)(如濃度、分壓等變量)的函數(shù).若仍以式(1)中關鍵組分濃度作自變量,將對模型的通用化造成困難,因此接下來討論如何對衡算方程中的變量進行替換.

摩爾分數(shù)表示混合物中某組分物質的量與混合物總物質的量之比,故采用關鍵組分的摩爾分數(shù)轉化反應動力學中相關參數(shù),將使得模型具有通用性.由氣體狀態(tài)方程可知

式中:p為總壓力,Pa;Z為氣體壓縮系數(shù);Rg為摩爾氣體常數(shù),J/(mol·K);yk為組分k的摩爾分數(shù).

氣體混合物體積流率及流速u隨反應溫度的變化和反應進行的程度而變化(對于非等摩爾反應),但直接求解流速u較難,所以需對方程進行一些代數(shù)變換,以簡化其計算難度.由式(7)可知

式中:F為混合物的摩爾流率,mol/s;S為反應器橫截面積,m2.將式(8)對l求導得

由于

將式(9)和(10)移項整理后代入式(1)得

需要說明的是,模型求解過程中物料衡算與熱量衡算方程須聯(lián)立求解.若反應體系中獨立反應數(shù)為m,則需選取m個關鍵組分并對其進行物料衡算.即式(11)實際涵蓋m個關鍵組分的物料衡算方程,將它與式(2)聯(lián)立后,即可得到變量+1元一次方程組.具體如下:

式(12)和(13)即為固定床反應器通用的二維擬均相穩(wěn)態(tài)物料衡算與熱量衡算的微分方程.

2.2混合物摩爾分數(shù)的轉化

根據(jù)2.1可知,式(12)中F為混合物的摩爾流率,它實際上可描述為關于m個關鍵組分的摩爾分數(shù)y的函數(shù),即

若反應體系中組分數(shù)為n,為方便公式推導,先將化學計量系數(shù)矩陣進行初等行變換,得到獨立反應的化學計量系數(shù)矩陣v,使得m個關鍵組分的化學計量系數(shù)為1,變換后的系數(shù)矩陣中不為1的系數(shù)用vij表示,即:

對所有組分進行物料衡算分析,以探求混合物的摩爾流率與關鍵組分的摩爾分數(shù)之間的關系,相關表達式如表1所示.表中列出了每個組分初始狀態(tài)和瞬時狀態(tài)下的摩爾分數(shù)以及摩爾流率的數(shù)學表達式,根據(jù)質量守恒原理可得第i個組分瞬時狀態(tài)下的摩爾流率表達式為:

式中:F0為混合物的初始摩爾流率,mol/s;yi,in為組分i的進料摩爾分數(shù).由于各組分的摩爾分數(shù)之和等于1,即:

因此將方程聯(lián)立后即可求得混合物摩爾流率與關鍵組分的摩爾分數(shù)之間的關系:

最后將式(21)代入式(12)即可.

表1　物料衡算表Tab.1 Material balance sheet

2.3正交配置法的應用

正交配置法是一種特殊的點配置方法,其具有接近于Galerkin方法的高精度,又保留了配置方法的簡單形式[16-17].根據(jù)1.1的假設,利用圓柱對稱性正交配置常數(shù),選取適當?shù)恼慌渲命c個數(shù)N后,即可將原偏微分方程組變?yōu)槌Ｎ⒎址匠探M,而后利用MATLAB工具箱函數(shù)求解.

r方向的正交配置即對式(12)和(13)的二次偏導項和一次偏導項展開,可得到內(nèi)配置點上的微分方程:

式中Aji、Bji、Rj是正交配置常數(shù).式(22)和(23)中所包含的方程數(shù)為2N,而變量數(shù)為2(N+1),為此還需增加濃度和溫度的邊界條件,即:

至此,式(22)～(25)組成常微分方程組封閉,然后利用MATLAB工具箱函數(shù)求解即可.

圖3　實例1中固定床反應器軸徑向溫度分布Fig.3 Axial and radial temperature distribution of fixed-bed reactor of example 1

3 程序設計與驗證

3.1程序設計框圖

圖2給出了所提出方法的程序設計框圖,并求解乙苯催化脫氫制苯乙烯與乙炔法制氯乙烯過程的模型來驗證所提出方法的有效性和通用性.

圖2　程序設計流程圖Fig.2 Program flow diagram

3.2實例分析

實例1為乙苯催化脫氫制苯乙烯的反應,利用所提出的方法結合文獻[12]所提供的動力學模型及操作條件對該過程進行模擬,模擬所得到的固定床軸徑向溫度分布和乙苯平均轉化率沿管長的分布如圖3～4.

文獻[12]模擬得到反應器出口中心處溫度為827.11 K,乙苯最終轉化率為46.83%.采用本文中提出的方法模擬得到反應器出口中心處溫度為819.57K,乙苯最終轉化率為43.15%.

圖4　乙苯平均轉化率沿管長的分布Fig.4 Axial average percent conversion distribution of ethyl benzene

圖5　實例2中固定床反應器軸徑向溫度分布Fig.5 Axial and radial temperature distribution of fixed-bed reactor of example 2

圖6　乙炔轉化率軸徑向分布Fig.6 Axial and radial percent conversion distribution of acetylene

為不失一般性,將文獻[13]中的乙炔法制氯乙烯體系作為實例2進行研究,將所得到的固定床軸徑向溫度分布和乙炔轉化率軸徑向分布的模擬結果列于圖5～6中.文獻[13]模擬結果為:固定床反應器內(nèi)最高溫度391.7 K,出口溫度374.3 K,反應器出口處乙炔平均轉化率81.52%.本文方法模擬結果為:固定床反應器內(nèi)最高溫度387.14 K,出口溫度372.65 K,反應器出口處乙炔平均轉化率88.09%.

通過2個實例的模擬結果分析,可知本文中所提出的方法解決了二維穩(wěn)態(tài)固定床反應器模型的通用求解問題,所得到的模擬結果與文獻變化趨勢吻合,但具體數(shù)據(jù)略有差異,這是由于原文獻均未考慮流速u的變化.

4 結 論

本文中提出了一種固定床反應器穩(wěn)態(tài)模型求解的通用計算方法,該方法有效處理了床層反應中氣相流速隨著反應進程的變化而不斷變化的問題,且該方法適用于反應動力學函數(shù)的不同形式.為滿足模型的通用性,首先將模型中的氣相流速轉化為關鍵組分摩爾分數(shù)的函數(shù),避免了將模型中流速視為常數(shù)而引起求解誤差的問題,從而提高了模擬結果的準確度,又降低了編程工作量.然后在此基礎上確立程序設計的基本步驟和算法框圖,并在MATLAB上實現(xiàn)了相應的功能.最后運用該方法模擬了固定床反應器中的2個具體實例,模擬結果與文獻結果趨勢相符,驗證了該方法的通用性和有效性.

[1] 朱炳辰.化學反應工程[M].北京:化學工業(yè)出版社,2009: 140-141.

[2] FROMENT G F.Fixed bed catalytic reactors-current design status[J].Industrial&Engineering Chemistry, 1967,59(2):18-27.

[3] AKPAN E,AKANDE A,ABOUDHEIR A,et al.Experimental,kinetic and 2-D reactor modeling for simulation of the production of hydrogen by the catalytic reforming of concentrated crude ethanol(CRCCE)over a Ni-based commercial catalyst in a packed-bed tubular reactor[J]. Chemical Engineering Science,2007,62(12):3112-3126.

[4] WANG Y N,XU Y Y,LI Y W,et al.Heterogeneous modeling for fixed-bed Fischer-Tropsch synthesis:reactor model and its applications[J].Chemical Engineering Science, 2003,58(3/4/5/6):867-875.

[5] ODENDAAL W,GOBIE W,CARBERRY J.Thermal parameter sensitivity in the simulation of the nonisothermal,non-adiabatic fixed bed catalytic reactor:the two-dimensional heterogeneous model[J].Chemical Engineering Communications,1987,58(1/2/3/4/5/6):37-62.

[6] ALI N B,PETERSEN K,WOHLIN C.A systematic literature review on the industrial use of software process simulation[J].Journal of Systems and Software,2014,97:65-85.

[7] YANG Y B,YAMAUCHI H,NASSERZADEH V,et al. Effects of fuel devolatilisation on the combustion of wood chips and incineration of simulated municipal solid wastes in a packed bed[J].Fuel,2003,82(18):2205-2221.

[8] SHIN D,CHOI S.The combustion of simulated waste particles in a fixed bed[J].Combustion and Flame,2000, 121(1/2):167-180.

[9] NIJEMEISLAND M,DIXON A G.CFD study of fluid flow and wall heat transfer in a fixed bed of spheres[J]. AICh E Journal,2004,50(5):906-921.

[10] BERNINGER J A,WHITLEY R D,ZHANG X,et al.A versatile model for simulation of reaction and nonequilibrium dynamics in multicomponent fixed-bed adsorption processes[J].Computers&Chemical Engineering,1991,15(11):749-768.

[11] LIAPIS A I,RIPPINT D W T.A general model for the simulaiton of multi-component adsorption from a finite bath [J].Chemical Engineering Science,1977,32(6):619-627.

[12] 黃華江.實用化工計算機模擬:MATLAB在化學工程中的應用[M].北京:化學工業(yè)出版社,2004:144-148.

[13] 陳靜.乙炔法氯乙烯反應器模擬與優(yōu)化研究[D].杭州:浙江大學,2011:33-35.

[14] 陳甘棠.化學反應工程[M].3版.北京:化學工業(yè)出版社, 2007:203-205.

[15] KISSINGER H E.Reaction kinetics in differential thermal analysis[J].Analytical Chemistry,1957,29 (11):1702-1706.

[16] 李希.化工問題的建模與數(shù)學分析方法[M].北京:化學工業(yè)出版社,2006:288-291.

[17] FINLAYSON B A.Packed bed reactor analysis by orthogonal collocation[J].Chemical Engineering Science, 1971,26(7):1081-1091.

General Calculation Method for Two-dimensional Fixed-bed at Steady State

WU Sikan,CAO Zhikai*,ZHOU Hua,CHEN Binghui

(College of Chemistry and Chemical Engineering,Xiamen University,Xiamen 361005,China)

Fixed-bed reactor is the most widely used reactor in the chemical process.However,the programming for solving its mathematical model has been case by case via solving partial differential equations till now.In such a situation,it is impossible to build up a complete calculation module for the fixed-bed reactor model and incur the workload and burden in simulation.In this work, a general method to obtain the distribution of temperature and concentration of fixed-bed reactor at steady state is proposed.The method makes use of mole fractions of key components as variables and then applies orthogonal collocation to solve the formula of partial differential equations,which is deduced in the material and heat balance.Molar flow rate of reaction mixture is introduced at the same time and the relationship between the instantaneous molar flow rate,the initial value and mole fraction of key components are also built with a variable transfer method,to make the solution have versatility.Therefore,the computing module,which is used to solve partial differential equations,becomes programmable independent of a specific reaction,i.e.,a common program module.For a specific reaction,only a kinetic equation program is required to be added into the module.

fixed-bed reactor;computational simulation;orthogonal collocation method;two-dimensional pseudo-homogeneous model

TQ 015

A

0438-0479(2016)06-0830-06

10.6043/j.issn.0438-0479.201603015

2016-03-08 錄用日期:2016-04-08

國家自然科學基金(21576228)

zhikai@xmu.edu.cn

吳斯侃,曹志凱,周華,等.一種用于二維穩(wěn)態(tài)固定床的通用計算方法[J].廈門大學學報(自然科學版),2016,55(6): 830-835.

WU S K,CAO Z K,ZHOU H,et al.General calculation method for two-dimensional fixed-bed at steady state[J].Journal of Xiamen University(Natural Science),2016,55(6):830-835.(in Chinese)