化工过程中模型的建立与计算讲解

化工过程中模型的建立与计算
姓名：田保华
化工过程中模型的建立与计算
1. 概述
过程的状态监控或过程的在线监测都需要建立合适的数学模型。

化工过程的数学模
型主要有三大类方法，即机理模型，统计模型和混合模型。

描述过程的方程组由过程机理出发，经推导得到，并且由实验验证，这样建立起来的模型就是机理模型。

机理模型方法需要凭借可靠的规律及经验知识来建立原始微分方程式，这些规律和经验知识必须被表达为一般的形式。

机理模型是对实际过程直接的数学描述，是过程本质的反映，因此结果可以外推。

数学模型也可以根据实验装置、中型或者大型工业装置的实测数据，通过数据的回归分析得到的纯经验的数学关系式，这就是统计模型。

统计模型和过程机理无关，是根据实验从输出和输入变量之间的关系，经分析整理得到的。

它只是在实验范围内才是有效的，因而不宜外推或者以较大幅度外推。

由于实验条件的限制，统计模型的局限性很大，所以总是希望尽可能建立机理模型。

对于化工过程来说，由于经验模型受到实际条件的限制，应用范围有限，机理模型求解又十分困难，这样就产生了第三种数学模型，即混合模型。

混合模型是对实际过程进行抽象概括和合理简化，然后对简化的物理模型加以数学描述，混合模型主要是设法回避过程中一些不确定的和复杂的因素，代之以一些统计的结果和一定的当量关系，它是半经验半理论性质的。

在化工过程的数学模拟中，混合模型是应用最广的一种模型。

例如，混合模型用于粉仓中粉体流动数学模型分析等等。

近年来，人们将人工神经网络方法用于化工建模，并取得较好的效果。

化工过程中数学建模的建立一般是基于流体的性质。

流体的热力学性质主要是从状态方程（EOS）得到。

至今，文献报道的EOS已有一百五十种之多，有的从理论分析得到的、有的从实验数据分析归纳而来、还有一些是理论分析和实验数据相结合推出来。

比较经典的EOS有VDW方程，R-K方程，Soave方程，CS方程，33参数的MBWR方程。

这些经验、半经验的EOS只能在一定的温度和密度范围内对于某些流体适应，应用的范围比较窄，理论基础不强[34]。

随着计算机技术的迅速的发展，现代化学工程越来越要求EOS的精度高，应用范围广，可靠性好。

因此，EOS 的研究已逐步从经验、半经验型向理论型发展，它们的应用范围广而且具有较强的理论基础。

近年来，研究的热点是从径向分布函数
（Radial Distribution Function，简称RDF）来得到理论型EOS。

2. 径向分布函数理论
分布函数方法不以任何物理模型为依据，而从解决粒子间的相互作用势能
u(r)出发的统计热力学理论，该方法的特点是以求出少数几个 RDF 去解决整个体 系的构型积分，然后利用它计算流体的全部热力学量， 该方法是现代流体理论的 一个最为活跃的研究方向，分布函数理论是所有溶液理论中最为精确的部分。

2.1径向分布函数的应用
径向分布函数是研究流体性质的基本内容之一。

原则上 ，确定了径向分布函 数，就可以求出流体的热力学性质。

(1) 由径向分布函数求系统的位能
文献⑵指出，在系统中，每一个分子对的平均位能应为：
u(r) = ° uWQ p ⑵(GDldAda
oO 2 1
(v )dr ! ° u(r)g(r)4二r dr 0 u(r)g(r)4：r 2dr
如设系统的位能为所有分子对的位能的总和，由于系统中各分子对的总合为
N(N-1)/2或N 2/2，因此，系统的位能的平均值为：
E p 二 N 2 u(r) 2 二少 ° u(r)g(r)4：r 2dr
(2) 径向分布函数导得状态方程
e -Ep/kT 少 1 dx N dy N dz N
u(r)二
V 2 (2-1) (2-2)
由方程: §ln ① P
"F T ,N
(2-3) 设流体装载在一边长为 a 的立方容器中，V=a 3，位能配分函数为:
(7 (：V )T ,N x'i = X j / a, y[ ’P (寸)T ,N = N J 不 N
点①'P
=NV N P 'P +V N (——)T N P eV
二 y i / a, z'^ z i / a
1 1 1 kT .|■严 E
P
、
e -E P / kTdx '1 …dz ' N
0 ( : V ) x'1 : ,X 'N e
(2-4) cE P "'P (一
)T ,N :V "p ( )T ,N :V NkT 3V i j _ ：-P 6kTV N , ——① V d ；P (r j )
r
j
dr j N
(V 1 )2 ::d ；P (r) dr
3
g(r)4「：r dr (2-5) N
P [1 (
6kT V )；p (r) 丿"0
)2 :: d ；p (r) 丿p
dr 1 N ( 0 V 6kT V 0 dr (3)径向分布函数计算真实流体的 P 、V 、T 性质 3
g (r )4 二 r dr ]
文献报道了用径向分布函数及 L-J 位能模型计算真实流体的PVT 性质。

文献
用顺序解析平衡法，结合L-J 位能函数模型计算径向分布函数， 然后直接代入统 计热力学所给出的理论状态方程预测真实流体的 PVT 性质。

L-J 位能模型参数有 临界温度T c 和临界压力P c 确定。

此方法可以能够比较成功的预测非极性和弱极 性纯流体的PVT 性
质。

当分子尺寸小于正已烷时，除临界温度附近之外，预测 的饱和液体体积的平均相对误差小于 5%,对于分子尺寸更大的流体，可以用象
Kihara 模型这样的3参数模型来提高预测精度。

具体公式如下: 在L-J 位能模型上应用分布函数理论的状态方程
其中M=6.02 x 1023, Vm 是摩尔体积，弓I 入对比半径
其中V b =N 0(二「3/6)，由于r ：：： 0.5二时，可认为g(r p 0 ,而时r • 3.5二，则可认 为g(r) 1，所以上式可以写成
上述方程中的积分项可由数值积分方法得出。

令温度T 等于某流体的临界温 度，调节二和；/k ，直到上式计算出的等温线出现水平拐点，而且水平拐点处压 力等于该流体的实测临界压力，此时对应的匚和；/k 值就是所需要确定的位能模 型的参数值。

当位能模型参数确定后，就可以预测该流体的 PVT 关系。

2.2 RDF 的获得
目前RDF 主要通过实验法和计算法来获得。

2.2.1实验法
RDF 的实验方法主要有X 射线衍射、电子衍射和中子衍射。

就国内的研究 情况来说，一般采用2二型液体X 射线衍射仪来获得径向分布函数。

2.2.2计算法
目前径向分布函数(RDF )的计算方法有计算机模拟法，积分法和微相平衡 法。

(1) 计算机模拟法：包括 MD 法和MC 法。

MD(Molecular Dynamics)法的关 键是解力学方程以便对这些微观状态对应的力学量作适当的平均。

而 MC(M onte Carlo)法借助计算机做随机抽样获取，要求取样的数目必须是大量的，且样品必N °kT
V m
3V ： dE P dr g (r)dr (2-6) -2 1 尹歹g(x)dx (2-7)
RT
3.5-2 1 0.5 (尹 F )g (x )dx 1 3 3.53 (2-8)
2二 2 0
须具有代表性，需要一种均匀分布的随机序列。

计算机模拟法求得的RDF 值精度 高，但是计算量大，需要很长的时间。

(2) 积分法：就是利用近似法解 Ornstein-
Zernike(OZ)积分方程来得到流
体的径向分布函数。

OZ 方程的表达式是：
h(r) =c(r)亠idr'c(r')h(r -r')
式中：r 表示分子1和分子2之间的距离，r 表示的它们之间的向量 r'表示
分子1和分子3之间的距离，r'表示它们之间的向量 h(r)—间接相关函数，c(r)—直接相关函数，T —流体的密度 h(r)和g(r)有如下公式：
h(r)二 g(r)—1
推导径向分布函数主要有以下近似方法：
1) Percuse-Yevick(PY)法：
h(r)f(r) =[1 f(r)]c(r) - f (r) f (r) =exp[- '-u(r)] -1
2) Mea n-Spherical Approximation (MSA)法：
c(r) - - - ；(r) 3) Hyper-Netted Chain (HNC)法:
c(r) = f(r)y(r)+y(r)-1-1n y(r)
y(r)二 g(r)/e(r)
e(r)二 exp[- ；(r)]
f (r) =e(r) -1
4) Expo nental (EXP)法：
g(r)二 g °(r)exp(g MSA (r) -g °(r))
式中：g MSA (r)表示MSA 法计算的径向分布函数 5) Tang 法:
~(k) =~(k)
「~(k)~(k) 其中：
〜 4二 c(k ^4T o rc(r)sin krdr
〜(k)和〜(k)分别表示h(r)和c(r)的三维傅利叶变换。

(3) 微相平衡法：陈光进等提出了一种计算径向分布函数的顺序解析法， 即微相平衡法。

在一个基准分子周围， 在半径方向上存在密度梯度，这就是径向 分布函数的物理〜 4 ■: h(k
"匚 oCi 0 rh (r ) sin krdr (2-9)
(2-10)
(2-11)
(2-12) (2-13)
(2-14) (2-15) (2-16)
意义。

对基准分子的周围空间沿半径方向微分分割，定义分割得
到的每一个微元体为一个微相。

不同的微相有不同的分子数密度和能量密度。

微相和一般的宏观相有区别，例如微相所含的分子数可以是分数；微相必须是以一个序列的形式出现，单一的微相是不存在的。

但对于一个平衡热力学体系，微相之间彼此应处于平衡，满足普通相平衡所需的条件。

利用相平衡的条件，就可以得到微相之间分子数密度的关系并进一步得到径向分布函数。

对于宏观的相平衡
问题，利用化学位准则，即可得到平衡相间的密度关系。

文献提出了描述相际分子转移行为的双阻力物理模型，并给出每项阻力的计算方法，导出一相中的分子向另一相转移的质量通量的计算通量准则公式，提出了相平衡的质量通量准则以代替化学位准则，并因此得到两相平衡时的密度关系方程。

该方法的特点是当已知前一位置的径向分布函数值就可以解析的求出。

此方法无须迭代，比分子模拟和积分方程的数值解法相比，节约了计算时间。

2.3位能函数模型
在计算位能配分函数时，需要计算系统的位能，所以根据位能函数的不同，计算径向分布函数也有几种模型。

(1) 硬球(Hard Sphere)模型
硬球模型是一种非常简单的模型，它粗略地反映了分子间极强的超短程的排斥力，它将分子看成没有吸引力的硬球，由于简单，常作为进一步处理的基础，其位能函数是：
‘0 r x
u(r)，(2-17)
oQ r <<y
(2) 方阱(Square Well)模型
方阱模型是一种理想模型，溶液理论和统计力学中得到广泛的应用。

它常用来作为实际流体的参考体系，相对于硬球模型来说，它是一种更好的参考体系，因为它以粗略的方式考虑了分子间吸引力和排斥力，而硬球模型认为只有斥力项没有引力项，所以在实际计算方阱模型为参考体系更恰当。

它将分子看成是直径为二并有吸引力的硬球，但是吸引力仅在分子间距为二处起作用，，称为对比阱宽，阱深为；，位能函数为:
u(r)=卜总
oO
(3) 萨日兰(Sutherland)模型萨日
兰模型将分子看作直径为
比，位能函数为:
(2-19)
(2-18) 二的有吸引力的硬球，但是吸引力与r6成反
这个模型比上面两个合理得多。

首先它对吸引力的考虑和色散作用一样。

另一方 面，实际的位能在斥力起作用时随分子间距缩小而上升很快， 这个模型和范德华 方程一致。

（4） Lennard-Jones （ LJ ）模型
Le nn ard-Jo nes （ LJ ）模型是一种很重要的模型，常用来研究简单流体的行为 以及汽液平衡、液液平衡。

对于复杂的分子，就用它来作它们的参考体系， LJ 模型对分子间的作用力作了很好的描述。

LJ 位能曲线能较好的描述位能与分子 间距的关系，应用比较广泛。

LJ 模型有两种形式：
LJ12-6 模型：
' 旳
r c R 灾）刊4订（二12一（与1 ,R
（2-20）
IL r r LJ9-6模型：
例如：引入LJ 位能函数对内压相进行修正的 van der Waals 模型，较仅对已占 体积相进行修正的van der Waals 模型明显更优。

（5） 基哈拉（Kihara ）模型
基哈拉模型将分子看作是一个由软的电子云包围着的一个硬核，
位能取决于
核间的最短距离，核为球形，半径为 a ，其位能模型为：
…2a 12 -…2a 、6i ；p=4[（ 「） -（ 2-
）6] （2-22）
r —2a r —2a （6） 斯托克迈尔（Stockmeyer ）模型
斯托克迈尔模型是在LJ 模型上增加一个静电作用项，其位能模型为：
.1.2
飞[2cos 寸j cos^j -sin 耳 sin 寸 j cos( ] _ ])] r 此式有3个特征参数；、；「、」，其中"■为偶极距可以单
独测定。

2.3小结：RDF 理论在流体的热力学性质的研究中有着重要的作用， 而其应用的 关键是RDF 的求取。

微相平衡法计算有一定的简化。

目前也有人采用人工神经 网络方法来研究RDF 。

3. 传统的汽液平衡热力学计算模型
目前，汽液平衡计算方法主要有状态方程法和状态方程加活度系数法
（简称 活度系数法）两种 Q Q 心,（；、；）6 r ：
： R (2-21)
(2-23)
12 u(r) =4；[(—)
3.1状态方程法（EOS 法）
采用状态方程计算汽液平衡，简称状态方程法（EOS 法） 由公式
已经建立的状态方程在一定条件下都可以同时用来描述汽、液两相的逸度行为, 这些方程有 R-K 方程，SRK （ Soave ）方程，Peng-Robinson （ PR ）方程，BWR 、 BWRS 方程等
BWR 和BWRS 方程早已被应用来计算汽液平衡，也比较精确。

但是由于方 程本身及其混合规则的复杂性限制了其应用。

现在应用最多的还是简单的立方型 方程，主要是Soave 和PR 方程，后一方程对于液相体积的再现优于前一方程， 但对于汽液平衡，总的来说，其精度并不比前者好。

对于烷烃混合物或烷烃和非 烷烃气体（如 N 2,O 2,CO 和 CO 2等）混合物体系，Daubert, Grabosk 与 Danner （1978） 经比较和评价，认为最为可信的还是 Soave 方程。

状态方程法作汽液平衡的具体计算时往往比较冗长，许多量需要迭代计算， 因此常用计算机来完成，下面示出已知 T/用R-K 方程求y i ,P 的计算框图，如
图1-1所示。

3.2活度系数法
根据溶液热力学理论，将液相中组分的逸度与组分活度系数相联系，
简称活
度系数法。

活度系数的关联，常用的有伍尔（ Whol ）型方程，局部组成型方程以及基 团贡献型模型。

Whol 型方程（包括Margules ，Van Laar 等方程）迄今仍在应用，优点是计 算比较简单，缺点是不能用二元系数据直接推算多元系的汽液平衡，
而必须要用 多元系参数。

此外，对于含有强极性组分，非理想性很高的体系，
Whol 型方程
往往难以发挥作用。

/ A L /A 7 K 汀.Xi -I
< 可求出给定条件下组分的 ,、V ,、L
VLE ，式中「、分别表示汽、液相中组分i 的逸度 系数。

状态方程法计算汽液平衡的关键是求出两相的分逸度系数 V L
：i 和i ，目前
目前，广泛使用的局部组成型方程主要有著名的Wils on方程和UNIQUAC 方程。

基团贡献模型主要包括 ASOG 法（Analytical Solutions of Groups ，基团解析 法）和 UNIFAC 模型（Universal Quasichemical Functional Group Activity Coefficient ，通用基团活度系数）。

由以上活度系数法计算出组分的活度系数，由下式即可计算出体系的汽液平 衡数据：
斗 ，0
K i 式中y 「Xi ——汽、液相中组分
V
i ――汽相混合物中组分i 在体系温度T 和体系压力P 的逸度系数；
i
――组分i 的活度系数； 图1-2给出以活度系数法计算泡点汽液平衡的计算框图。

4.编程
目前，用于编程的语言很多，如C 语言，VB 语言，MATLAB 语言、labVIEW 及ActiveOberon 等等。

在计算中为了取得更高的精度，可以采用几种语言混合 编程，例如VB 与MATLAB6.5相结合的混合编程，综合了二者的优点，运行 结果明显比单程序运行结果好得多。

其中：
VB 具有良好的可视化界面、程序集成化程度高、易学易用等优点。

但在数 值计算和图形处理方面显得力不从心，而 MATLAB 正好弥补了这一缺点。

MATLAB 是一种功能强大的高级计算语言，它是由美国 Mathworks 公司80 年代逐步开发完善的一种新型语言，可以运行于 MS-Windows3.1、Windows95、 Windows NT 、OS/2、Unix 、Mac 及VMS 等操作系统。

它是一种主要用于数值计 算及可视 化图形处理的高级计算语言。

在矩阵运算方面更显示 其优势。

MATLAB [31]具有一系列专题应用软件包，称之为工具箱（Toolbox ）。

工具箱是极 为广泛的集合了 MATLAB 函数（M-files ），而这些函数大大地扩展了 MATLAB 环境，使之能够解决某些专门的问题。

工具箱所涉及的领域，包括信号处理、控 制系统设计、鲁棒控制、动态系统仿真、系统辨识、神经网络、图像处理、统计 学、卩分析、模型逻辑、优
化、样条函数、符号运算、图形设计、符号数学处理、 偏微分方程处理、高阶谱分析等。

y i
X i
i 的摩尔分率;
图1-1 已知T , 1,
勺用R-K 方程求P, y i 的框图
参考文献
[1] 黄忠霜，控制系统MATLAB 计算与仿真.国防工业大学出版社.2004.
[2] 胡英，流体的分子热力学[M]. 北京：高等教育出版社. 1983，164
[3] WangYu,Computer Simulation of Martensitic Transformation in Ni-Al Alloy[J].Journal of Hu
Zhou Teacher College.2004,26(2):33-36
[4] 陈光进、郭天民. 用径向分布函数及Lennard-Jones 位能模型计算真实流体的PVT 性质[J].
高校化学工程学报. 1999, 13(2): 101~106
[5] 周世琦，LJ 模型的热力学特性、相平衡及用于真实分子的研究. 化学物理学报.
2005,18(4): 487-488
⑹李志义等，超临界C02二元系统相平衡的数值模拟•化学工程• 2005, 33 ( 2)
[7] 窦晓冬等，修正的van der Waals 模型对其液相流体的应用. 甘肃化工. 2005,(2):18-20。