1、1CHEMKIN 4.0.1 入门指南燃烧学 1辅助教程文中如有任何错误,敬请指出,以便不断改进;如有任何问题,欢迎提出,共同探讨助教博士生:卢智恒联系方式:热能系系馆办公室 201(O)62782108 (H )2005.32一、CHEMKIN 的安装和简介11 安装 CHEMKINChemkin 最早的版本始于 1980,由美国 Sandia 实验室的 Kee RJ 等人编写,经过多年的不断发展日趋完善。后来由 Reaction Design 公司收购并继续开发,目前最新版为 4.0.1。由于学习和科研需要,我们花费 12000$向 ReactionDesign 公司订购了一套最新版本的C
2、HEMKIN 4.0.1,其中包括可供 20 人同时在线计算的 license,用于燃烧学课程的学习。【安装】请登录 ftp:/combustion:combustion166.111.56.155 下载相关文件,其中chemkin401_pc_setup.exe 为 CHEMKIN 的安装程序,chemkin.lic 为网络认证文件,详细的安装信息可以参看 ftp 上的“安装说明.txt”文件。安装完后会自动在桌面及开始菜单建立快捷方式。【注意】1、本套教学用的 CHEMKIN 软件采用网络认证的方式,请确保电脑已经联网(校内) ,否则无法计算。2、建议采用 1024768 的分辨率,否则某
3、些界面将无法完全显示。12 CHEMKIN 简介CHEMKIN 是一种非常强大的求解复杂化学反应问题的软件包,常用于对燃烧过程、催化过程、化学气相沉积、等离子体及其他化学反应的模拟。CHEMKIN 以气相动力学、表面动力学、传递过程这三个核心软件包为基础,提供了对 21 种常见化学反应模型及后处理程序。三个核心程序模块为:1) 气相动力学(Gas-Phase Kinetics):是所有程序计算的基础,包括气相成分组成、气相化学反应与相关的 Arrhenius 数据等信息。2) 表面动力学(Surface Kinetics) 。很多反应过程包括多相反应,如催化反应、化学气相沉积、固体腐蚀等。在这
4、些反应里,Surface Kinetics 提供两相反应所需的各种信息,如表面结构、表面和体内的成分组成及热力学数据、表面化学反应等。3) 传递(Transport) 。提供气相多组分粘度、热传导系数、扩散系数和热扩散系数等。其中 Surface Kinetics 和 Transport 必须以 Gas-Phase Kinetics 为基础,因为它们中出现的成分都必须在 Gas-Phase Kinetics 中已定义。Gas-Phase Kinetics、Surface Kinetics 和 Transport 提供了化学反应的基本信息,生成动态链接库供后续程序调用。用户可以自己编写程序调用它
5、们来实现指定的功能,但最方便的是使用 CHEMKIN 自带的反应模型,共分 7 类,包括 21 个反应模型和 2 个应用程序,将在下一节介绍。3二、CHEMKIN 的简单使用入门21 CHEMKIN 的用户界面CHEMKIN 的用户界面如图 1 所示。除了传统的 Windows 菜单和按钮以外,CHEMKIN 的工作界面分为任务栏、窗口栏、消息栏三部分。图 1 CHEMKIN 4.0.1 的用户界面其中任务栏包括 Open Projects 和 Models 两个选项卡,Models 选项卡中列出了 CHEMKIN提供的 21 种反应模型和 2 种应用程序,分 7 类排列,如表 1 所示。有关
6、 Open Projects 选项卡的使用将在下文的例子中介绍。表 1 CHEMKIN 提供的 7 类 21 种反应器与 2 种应用程序External Source of Inlet Gas,添加入口气体源,一般在有多路气体输入时才使用Non-Reactive Gas Mixer,计算不反应气体组分的混合过程Chemical and Phase Equilibrium Calculation,计算化学平衡和相平衡Miscellaneous (杂项)Mechanism Analyser,分析气相和表面化学反应系统中的热化学、传递及动力学数据Closed Internal Combustion
7、Engine Simulator,模拟封闭的内燃机内的燃烧情况Closed 0-D Reactors (封闭 0 维反应器) Closed Homogenous Batch Reactor,模拟封闭的全混同性反应器,包括定压、定容反应器4Closed Partially Stirred Reactor(PaSR),模拟封闭的部分混合的反应器Closed Plasma Reactor,模拟封闭的等离子体反应器Perfectly Stirred Reactor(PSR),模拟稳态及瞬态的充分混合反应器Plasma PSR,模拟充分混合等离子体反应器。Open 0-D Reactors (开口 0
8、维反应器)Partially Stirred Reactor(PaSR),模拟部分混合反应器Plug Flow Reactor,一维柱塞流反应器Plasma Plug Flow Reactor,等离子体一维柱塞流反应器Planar Shear Flow Reactor,平板层流反应器Cylindrical Shear Flow Reactor,圆柱层流反应器Flow Reactors (流动反应器)Honeycomb Monolith Reactor,蜂窝结构反应器Premixed Laminar Burner-stabilized Flame,模拟层流预混的稳态火焰Premixed Lami
9、nar Flame-speed Calculation,层流预混火焰传播速度的计算Flame Simulators (火焰模拟) Diffusion or Premixed Opposed-flow Flame,模拟扩散或预混的对冲火焰Stagnation Flow CVD Reactor,模拟用于化学气相沉积的滞流反应器CVD Reactors (化学气相沉积反应器) Rotating Disk CVD Reactor,模拟用于化学气相沉积的转盘流反应器Normal Incident Shock,模拟入射激波的化学动力学Shock Tube Reactors (激波管道反应器) Normal
10、 Reflected Shock,模拟反射激波的化学动力学22 CHEMKIN 的求解过程1、Gas Phase Kinetics(气相动力学)的处理Gas Phase Kinetics 的前处理器 (Pre-processor)读取用户编写的气相动力学输入文件和自带的热力学数据库(therm.dat),生成包含元素、组分、热力学数据反应信息的 Gas-Phase Kinetics 连接文件。Gas-Phase Kinetics 提供子程序库处理该连接文件。2、Surface Kinetics(表面动力学)和 Transport(传递过程)的处理如果化学反应包含表面反应或传递过程,则需要相应地
11、执行这两个核心程序块。Surface Kinetics 的前处理器读取用户编写的表面动力学输入文件,生成包含表面反应信息的 Surface Kinetics 连接文件,Surface Kinetics 提供子程序库处理该连接文件。5Transport 的前处理器根据 Gas-phase Kinetics 连接文件中的信息,自动从 CHEMKIN 自带的传递数据库(tran.dat) 读取相应的数据,然后生成包含传递信息的 Transport 连接文件,Transport 提供子程序库处理该连接文件。3、反应模型求解根据问题需要,CHEMKIN 读取模型输入文件确定求解方法。在 4.0 及其以上
12、版本的CHEMKIN 中,新加入了通过相应的模型设置窗口中设置模型的参数的功能。用户设置完成后点击 Create Input File 即可生成模型输入文件,然后用户可以通过 View Input File 按钮查看输入文件的内容。而对于 4.0 以前的版本,用户需要通过手动编辑模型输入文件。但对于文件的内容,都是相同的,都是采用关键字的形式声明模型的功能调用和参数设置。有关常用的平衡计算模型和全混反应器模型的关键字列于附录 2、3。完成所有设置后,即可 Run Model 进行计算,CHEMKIN 会自动调用上述 Gas-Phase Kinetics、Surface Kinetics、 Tr
13、ansport 各自的子程序来读取反应信息,调用模型输入文件控制模型求解过程。程序计算结束后,会生成一数据文件 xxxxxx.out 供用户查阅数据,以及一动态连接文件 XMLdata.zip 供后处理(绘图)使用。4、后处理(Post-Process )CHEMKIN 提供了统一的后处理器,用于对应用程序的结果进行分析和绘图。23 CHEMKIN 的用户操作步骤下面结合实际例子介绍用户操作的步骤。【例 1】计算化学当量的 H2 与空气的定压绝热燃烧温度。1运行 CHEMKIN,点击菜单 Project-New,输入项目名称,这里我们定为 H2-air。2. 决定问题的性质,选择适当的反应模型
14、。由于绝热燃烧温度的计算是一个相平衡过程,不涉及具体的反应过程,于是我们很容易地想到用 Chemical and Phase Equilibrium Calculation 模型。点击任务栏中的 model 选项卡,点击 Chemical and Phase Equilibrium Calculation 图标,此时在窗口栏的 Diagram View 窗口中将看到新加入一个相平衡计算的模型,最后点击窗口右下角黄色的 Update Project 按钮,如下图所示。6通常情况下,一般的化学反应问题通过适当的假设和简化,都可以对应到某一种CHEMKIN 包含的反应模型,有时一个问题还可以有多种选
15、择。选择恰当的应用程序是求解问题的第一步。但是如果实在不幸没有一种模型可以很好的解决你的问题,或者你想要更完美地解决一些问题,就只有自己编写程序,调用 CHEMKIN 里的子程序库进行计算了。3此时左侧任务栏会自动切换至 Open Projects 选项卡,双击 Pre-Processing,窗口栏出现 Pre-Processing 的参数窗口。在窗口中的 Working Dir 一项中填入你希望的保存路径,或者通过右侧的 Browse 按钮点击选取。然后按 New Chemistry Set 按钮,点击 Gas-Phase Kinetics Files 项右端的编辑按钮,如下图所示。在弹出的
16、窗口中选择刚才的工作路径,输入文件名 chem.inp,按 Open/Create 按钮后即可开始编辑气相动力学输入文件了。-【CHEMKIN 的文件规则 】-CHEMKIN 的输入文件有它自己的规则,用户在编辑输入文件的时候应该遵守这些规则。在介绍 Gas-Phase Kinetics 输入文件之前,先介绍输入文件的一些通用规则: 注释符号“!” 。符号“!”无论出现任何位置,此行后面的文本将作为注释文本而被忽略。 输入文件每行不应超过 80 个字符 除了个别有极其严格规则的地方外(如热力学数据的定义等,均会特殊声明) ,空格作为分隔符,而且多个空格将被视为一个。 数字格式:可以为整数(如
17、99) 、浮点数如(99.99 ) 、或 E 格式(如 9.999E2,E大小写均可) 。下面介绍如何编写 Gas-Phase Kinetics 的输入文件。该文件包括四部分的内容:元素、组分、热力学数据、化学反应,如下例所示:! 例: Gas-Phase Kinetics输入文件ELEMENTS H O END ! 元素定义SPECIES H2 H O2 O OH H2O END ! 组分定义THERMO ! 热力学数据( 本例中只重新定义了“OH ”的热力学数据)OH 121286O 1H 1 G 0300.00 5000.00 1000.00 10.02882730E+02 0.1013
18、9743E-02-0.02276877E-05 0.02174683E-09-0.05126305E-14 20.03886888E+05 0.05595712E+02 0.03637266E+02 0.01850910E-02-0.16761646E-05 370.02387202E-07-0.08431442E-11 0.03606781E+05 0.13588605E+01 4ENDREACTIONS ! 反应方程及Arrhenius系数H2+O2=2OH 0.170E+14 0.00 47780OH+H2=H20+H 0.117E+10 1.30 3626O+OH=O2+H 0.400
19、E+15 -0.50 0O+H2=OH+H 0.506E+05 2.67 62902OH=O+H2O 0.600E+09 1.30 0H+H+M=H2+M 0.100E+19 -1.00 0H2O/0.0/ H2/0.0/ ! 辅助数据H+H+H2=H2+H2 0.920E+17 -0.60 0H+H+H2O=H2+H2O 0.600E+20 -1.25 0END 元素(Elements)规则此部分以 ELEMENTS(或者 ELEM,两者等价)关键字开头;其后以空格为间隔符列出将在反应中出现的所有元素;最后以 END 关键字结束。例:ELEMENTS H O END ! 元素定义用户要按周期
20、表的元素名(两个字母均须大写)来定义元素,CHEMKIN 可以辨认的元素如表 2.1。如果用户定义同位素或新元素的话,可以以 12 个字母命名(不与元素表重名) ,并将其原子量以“/”括住跟随其后。如定义氢( H)的同位素氚(命名为 HH)如下:ELEM HH / 3.0 / END表 2 CHEMKIN 的元素表H, HE, LI, BE, B, C, N, O, F, NE, NA, MG, AL, SI, P, S, CL, AR, K, CA, SC, TI, V, CR, MN, FE, CO, NI, CU, ZN, GA, GE, AS, SE, BR, KR, RB, SR,
21、Y, ZR, NB, MO, TC, RU, RH, PD, AG, CD, IN, SN, SB, TE, I, XE, CS, BA, LA, CE, PR, ND, PM, SM, EU, GD, TB, DY, HO, ER, TM, YB, LU, HF, TA, W, RE, OS, IR, PT, AU, HG, TL, PB, BI, PO, AT, RN, FR, RA, AC, TH, PA, U, NP, PU, AM, CM, BK, CF, ES, FM, D, E其中,D 为氢(H)的同位素氘(D) ,E 为电子,如果化学反应中有离子参加,电子必须作为一个元素进行定
22、义。 组分(Species)规则此部分以 SPECIES(或 SPEC)开头;其后以空格为分隔符列出将在反 应中出现的所有组份;最后以关键字 END 结束。例:SPECIES H2 H O2 O OH H2O END ! 组分定义对于我们所要用到的组分,都可以从 CHEMKIN 的热力学数据库 them.dat 找到。CHEMKIN 数据库自带的组分名称(778 种)列入附录 1 中。我们要使用这些组分的数据,组分的命名就必须按照 CHEMKIN 的规则来,这样 CHEMKIN 才可以从数据库中自动获取该组分的元素组成及其热力学性质。当然,我们也可以用不超过 16 个字符而且以字母开头的任意字
23、符串来定义自己的组分,但组分的热力学数据同样需要自行定义。 热力学数据(Thermodynamic Data)规则此部分以 THERMO 开 头;其后以空格为分隔符列出定义组份的名称和热力学数据;最后以关键字 END 结束。例:THERMO ! 热力学数据( 本例中只重新定义了“OH ”的热力学数据)OH 121286O 1H 1 G 0300.00 5000.00 1000.00 10.02882730E+02 0.10139743E-02-0.02276877E-05 0.02174683E-09-0.05126305E-14 20.03886888E+05 0.05595712E+02
24、0.03637266E+02 0.01850910E-02-0.16761646E-05 30.02387202E-07-0.08431442E-11 0.03606781E+05 0.13588605E+01 4ENDCHEMKIN 的热力学数据是基于 的多项式模拟。这里不作介绍了,有兴00SHCP、趣的同学可参看帮助文件。 化学反应(Reaction)8规则起始行:关键字 REACTIONS (或 REAC),其后 为 Arrhenius 系数的单位(可选)。中间行:反应方程式,然后为该 方程的 Arrhenius 系数( 依次 为 Ai,i,Ei);有些反应需要辅助数据补充说明,出 现在
25、该反应方程式的下一行。结束行:关键字 END。例:REACTIONS ! 反应方程及Arrhenius系数H2+O2=2OH 0.170E+14 0.00 47780OH+H2=H20+H 0.117E+10 1.30 3626O+OH=O2+H 0.400E+15 -0.50 0O+H2=OH+H 0.506E+05 2.67 62902OH=O+H2O 0.600E+09 1.30 0H+H+M=H2+M 0.100E+19 -1.00 0H2O/0.0/ H2/0.0/ ! 辅助数据H+H+H2=H2+H2 0.920E+17 -0.60 0H+H+H2O=H2+H2O 0.600E+2
26、0 -1.25 0ENDArrhenius 定律: )exp(TREAKciifii1)在起始行内可以定义 Arrhenius 系数 Ai 和 Ei 的单位:Ei 的单位可以定义为 CAL/MOLE、 KCAL/MOLE、JOULES/MOLE、KJOULES/MOLE、KELVINS、EVOLTS;Ai 的单位可以用关键字 MOLES 或 MOLECULES 来定义,分别对应 cmmolesecK和 cmmoleculessecK。如果没有定义,Ei 和 Ai 的默认单位是:cal/mole 和 cmmolesecK。2)反应方程中,“=”和“”用于可逆反应(两者等价),“=”用于不可逆反应
27、,“+M”表示催化剂,“(+M) ”用于压力控制反应,“HV”表示光子,“E”表示电子。3)有些反应方程之后需要有补充的辅助数据说明。CHEMKIN 有很多的辅助数据说明,有兴趣的可参看帮助手册。这里只介绍两个比较常用的FORD 和 RORD,格式为FORD /组分 反应级数 /RORD /组分 反应级数 /前者是重新定义正向反应时某组分反应级数,后者是重新定义逆向反应时某组分的反应级数的。-对于本例,由于绝热燃烧温度不考虑具体的反应过程,所以不涉及具体的反应方程,故 Reaction 部分可以省略(即使写了也用不上) ,此时 CHEMKIN 将根据反应物和生成物的焓(定压时)或内能(定容时)
28、进行平衡计算。同时由于涉及的反应物和生成物的热力学数据都已在 CHEMKIN 的热力学数据库中有定义,故气相动力学输入文件只需要输入元素和组分两部分,如下所示:ELEM H O N ENDSPECIES O O2 H2 H OH HO2 H2O H2O2 N2 END保存文件后关闭编辑器,此时 Pre-Processing 窗口中的 Gas-Phase Kinetics Files 一项仍为空白,通过点击该项的浏览按钮选取刚才编辑的输入文件。然后需要在 Thermodynamic Data File 项中指定热力学数据库,点击该项的浏览按钮。在弹出的窗口中先点击右边的 Special Dire
29、ctory 中的 System Data,然后点击左边出现的therm.dat,按 select 按钮完成热力学数据库的选择,如下图所示。9回到 Pre-Processing 窗口后点击 Save As按钮,在弹出的窗口中直接点击 Save 按钮,以默认的文件名和路径保存。由于本例不涉及表面反应,所以不需要表面动力学输入和气体传递数据文件,故再次返回 Pre-Processing 窗口后即可点击 Run Pre-Processor 按钮运行预处理了。预处理的结果可以在下拉式按钮 View Results中查看。4. 双击任务栏中的 Cluster1,会弹出子菜单,双击其中的每项都会弹出相应的设
30、置窗口,其中: Cluster Properties:设置问题的计算方式,一般选择默认值 Normal Start 即可,即进行新的计算;此外还可以选择 Initialization of Reactor from Previous Solution 或 Initialization of Reactor from Solution File,对应从上次计算的结果开始继续进行计算和从已存在的解文件开始进行计算,在此不介绍,有兴趣的可参考帮助文件。 C1_R1 Equilibrium:包括 Reactor Physical Property 和 Species-specific Data 两个选
31、项卡,分别用来设置反应条件和反应物的组成,是重点的参数输入部分。对于本例要解决的化学当量的 H2 和空气的定压绝热燃烧温度即为定压、定焓问题,所以在反应器物性的选项卡中:a)Problem Type 中选择 Constant Pressure Ehthalpy;b)点击 Calculate Species Composition,即通过平衡计算确定组分的组成;c)点击 Temperature,输入 298,单位 K;d)点击 Pressure, 输入 1,单位 atm;在组分组成的选项卡中:a)点击 Reactant Fraction 选项卡;b)Unit Selection 选择 mole
32、fraction;c)Species 下拉式按钮中依次选择和输入: H2,2,Add 按钮;O2,1,Add 按钮;N2,3.76,Add 按钮。通过以上步骤完成反应条件和反应物的输入。 Solver:设置问题的求解方法。在本例中没有可供设置的参数。对于其它问题,稳态计算可能需要设置迭代的步长、迭代次数等,瞬态计算可能需要设置反应结束的时间、解的误差、灵敏度系数的误差等。可以参看下文中对本问题用不同模型求解时的例子。 Output Control:设置输出文件的大小限制、输出的灵敏度结果的个数等。对本例,我们不需做任何设置,均采用默认值即可。10 Continuations:设置重复计算的次数
33、,可以为 0、1、2、3。当为 0 时,表示不作重复计算。对于本例,我们不需要作重复计算,输入 0 即可;或者直接不点击Continuations 该选项,使用其默认值 0。5. 双击任务栏的 Run Model,在右侧窗口中点击 Create Input File 按钮,弹出提示输入文件名的窗口,直接按 OK 采用默认文件名即可。输入文件建立成功,此时可通过 View Input File 查看输入文件的内容,但请勿手动修改其中的内容,否则模型可能无法运行。如需要手动增加输入参数,请使用 Add Supplemental Input 按钮,在此不作介绍。点击 Run Model 按钮,消息栏
34、提示模型运行成功。如下图所示。点击 View Results 下拉式按钮中的 Output File 即可查看计算结果,本例中得到平衡状态的温度为 2.3892E+03(K),此即所求的定压绝热燃烧温度。我们可以尝试点击 Run Post Processor 按钮,在弹出窗口中选择要查看的的结果及其单位,按 OK 后发现弹出一个 Warning 的消息窗口,提示计算结果只有一组数值,不满足作图要求(至少两组数值),得到的图形也只是作为示例的 Cos(x)、Sin(x)图形。所以这里所要说明的是,使用后处理器(Post Processor)绘图需要有多个状态、多个时刻、多个计算次数的计算结果,一
35、般用于瞬态燃烧时各参数随时间的变化的绘图,这在下文可以看到。思考与实践:如果反应物的组分是化学当量比的 H2 和 O2,定压绝热燃烧温度和定容绝热燃烧温度又该如何设置求解?结果如何?【例 2】计算化学当量比的 H2 和空气在全混燃烧器中的绝热燃烧过程1. 建立新 Project,选择 Perfectly Stirred Reactor(PSR)全混反应器模型,Update Project。2. Pre-Processing 中,设定工作路径,点 New Chemistry Set,输入气相动力学文件chem.inp 如下:ELEMENTS O H N ENDSPECIES O O2 H H2 OH HO2 H2O H2O2 N2 ENDREACTIONSH+O2 = O+OH 5.1E16 -0.82 16510 !MILLER 81H2+O = H+OH 1.8E10 1.0 8830 !MILLER 77