Chemkin 4.0.1 帮助 - 图文

更新时间:2023-09-25 08:08:01 阅读量: 综合文库 文档下载

说明:文章内容仅供预览,部分内容可能不全。下载后的文档,内容与下面显示的完全一致。下载之前请确认下面内容是否您想要的,是否完整无缺。

CHEMKIN 4.0.1入门指南

——《燃烧学1》辅助教程

文中如有任何错误,敬请指出,以便不断改进;如有任何问题,欢迎提出,共同探讨

助教博士生:卢智恒

联系方式:热能系系馆办公室201

(O)62782108

(H)62779574

luzhiheng@tsinghua.org.cn

2005.3

1

CHEMKIN 4.0.1入门指南 ............................................................................................................. 1 第一章 安装和简介 ...................................................................................................................... 3

1.1 安装.................................................................................................................................. 3 1.2 简介.................................................................................................................................. 3 第二章 简单使用入门................................................................................................................... 4

2.1 用户界面.......................................................................................................................... 4 2.2 求解过程.......................................................................................................................... 5

2.2.1 Gas Phase Kinetics(气相动力学)的处理 ........................................................ 5 2.2.2 Surface Kinetics(表面动力学)和Transport(传递过程)的处理 ................ 6 2.2.3 反应模型求解 ....................................................................................................... 6 2.2.4 后处理(Post-Process) ...................................................................................... 6

第三章 用户操作步骤................................................................................................................... 6

3.1 例1定压绝热燃烧 .......................................................................................................... 6

3.1.1 运行CHEMKIN ................................................................................................... 6 3.1.2 选择反应模型 ....................................................................................................... 6 3.1.3 前处理 ................................................................................................................... 7 3.1.4 边界条件 ............................................................................................................. 11 3.1.5 运行计算 ............................................................................................................. 12 3.2 例2全混燃烧器中的绝热燃烧 .................................................................................... 12

3.2.1 运行PSR ............................................................................................................. 12 3.2.2 前处理 ................................................................................................................. 13 3.2.3 对Cluster1中的设置 ......................................................................................... 13 3.2.4 运行计算 ............................................................................................................. 14

附录1 CHEMKIN热力学数据库中包含的组分 ...................................................................... 15 附录2 化学平衡和相平衡计算的模型输入文件常用关键字 .................................................. 16 附录3 全混反应器的模型输入文件常用关键字 ...................................................................... 17 附录4 与燃烧和热力学数据有关的网络资源 .......................................................................... 20

2

第一章 安装和简介

1.1 安装

Chemkin最早的版本始于1980,由美国Sandia实验室的Kee RJ等人编写,经过多年的不断发展日趋完善。后来由Reaction Design公司收购并继续开发,目前最新版为4.0.1。由于学习和科研需要,我们花费12000$向ReactionDesign公司订购了一套最新版本的CHEMKIN 4.0.1,其中包括可供20人同时在线计算的license,用于《燃烧学》课程的学习。

【安装】请登录ftp://combustion:combustion@166.111.56.155 下载相关文件,其中chemkin401_pc_setup.exe为CHEMKIN的安装程序,chemkin.lic为网络认证文件,详细的安装信息可以参看ftp上的“安装说明.txt”文件。安装完后会自动在桌面及开始菜单建立快捷方式。

【注意】1、本套教学用的CHEMKIN软件采用网络认证的方式,请确保电脑已经联网(校内),否则无法计算。

2、建议采用1024×768的分辨率,否则某些界面将无法完全显示。

1.2 简介

CHEMKIN是一种非常强大的求解复杂化学反应问题的软件包,常用于对燃烧过程、催化过程、化学气相沉积、等离子体及其他化学反应的模拟。CHEMKIN以气相动力学、表面动力学、传递过程这三个核心软件包为基础,提供了对21种常见化学反应模型及后处理程序。

三个核心程序模块为:

1) 气相动力学(Gas-Phase Kinetics):是所有程序计算的基础,包括气相成分组成、气相化学反应与相关的Arrhenius数据等信息。

2) 表面动力学(Surface Kinetics)。很多反应过程包括多相反应,如催化反应、化学气相沉积、固体腐蚀等。在这些反应里,Surface Kinetics提供两相反应所需的各种信息,如表面结构、表面和体内的成分组成及热力学数据、表面化学反应等。

3) 传递(Transport)。提供气相多组分粘度、热传导系数、扩散系数和热扩散系数等。 其中Surface Kinetics和Transport必须以Gas-Phase Kinetics为基础,因为它们中出现的成分都必须在Gas-Phase Kinetics中已定义。

Gas-Phase Kinetics、Surface Kinetics和Transport提供了化学反应的基本信息,生成动态链接库供后续程序调用。用户可以自己编写程序调用它们来实现指定的功能,但最方便的是使用CHEMKIN自带的反应模型,共分7类,包括21个反应模型和2个应用程序,将在下一节介绍。

3

第二章 简单使用入门

2.1 用户界面

CHEMKIN的用户界面如图1所示。除了传统的Windows菜单和按钮以外,CHEMKIN的工作界面分为任务栏、窗口栏、消息栏三部分。

图1 CHEMKIN 4.0.1的用户界面

其中任务栏包括Open Projects和Models两个选项卡,Models选项卡中列出了CHEMKIN提供的21种反应模型和2种应用程序,分7类排列,如表1所示。有关Open Projects选项卡的使用将在下文的例子中介绍。

表1 CHEMKIN提供的7类21种反应器与2种应用程序 Miscellaneous (杂项) External Source of Inlet Gas,添加入口气体源,一般在有多路气体输入时才使用 Non-Reactive Gas Mixer,计算不反应气体组分的混合过程 Chemical and Phase Equilibrium Calculation,计算化学平衡和相平衡 Mechanism Analyser,分析气相和表面化学反应系统中的热化学、传递及动力学数据 Closed Internal Combustion Engine Simulator,模拟封闭的内燃机内的燃烧情况 Closed 0-D Reactors (封 4

闭0维反应器) Closed Homogenous Batch Reactor,模拟封闭的全混同性反应器,包括定压、定容反应器 Closed Partially Stirred Reactor(PaSR),模拟封闭的部分混合的反应器 Closed Plasma Reactor,模拟封闭的等离子体反应器 Perfectly Stirred Reactor(PSR),模拟稳态及瞬态的充分混合反应器 Plasma PSR,模拟充分混合等离子体反应器。 Partially Stirred Reactor(PaSR),模拟部分混合反应器 Open 0-D Reactors (开口0维反应器) Plug Flow Reactor,一维柱塞流反应器 Plasma Plug Flow Reactor,等离子体一维柱塞流反应器 Flow Reactors (流动反应器) Planar Shear Flow Reactor,平板层流反应器 Cylindrical Shear Flow Reactor,圆柱层流反应器 Honeycomb Monolith Reactor,蜂窝结构反应器 CVD Reactors (化学气相沉积反应器) Flame Simulators (火焰模拟) Premixed Laminar Burner-stabilized Flame,模拟层流预混的稳态火焰 Premixed Laminar Flame-speed Calculation,层流预混火焰传播速度的计算 Diffusion or Premixed Opposed-flow Flame,模拟扩散或预混的对冲火焰 Stagnation Flow CVD Reactor,模拟用于化学气相沉积的滞流反应器 Rotating Disk CVD Reactor,模拟用于化学气相沉积的转盘流反应器 Normal Incident Shock,模拟入射激波的化学动力学 Shock Tube Reactors (激波管道反应器) Normal Reflected Shock,模拟反射激波的化学动力学 2.2 求解过程

2.2.1 Gas Phase Kinetics(气相动力学)的处理

Gas Phase Kinetics的前处理器(Pre-processor)读取用户编写的气相动力学输入文件和自带的热力学数据库(therm.dat),生成包含元素、组分、热力学数据反应信息的Gas-Phase Kinetics连接文件。Gas-Phase Kinetics提供子程序库处理该连接文件。

5

2.2.2 Surface Kinetics(表面动力学)和Transport(传递过程)的处理

如果化学反应包含表面反应或传递过程,则需要相应地执行这两个核心程序块。 Surface Kinetics的前处理器读取用户编写的表面动力学输入文件,生成包含表面反应信息的Surface Kinetics连接文件,Surface Kinetics提供子程序库处理该连接文件。

Transport的前处理器根据Gas-phase Kinetics连接文件中的信息,自动从CHEMKIN自带的传递数据库(tran.dat)读取相应的数据,然后生成包含传递信息的Transport连接文件,Transport提供子程序库处理该连接文件。

2.2.3 反应模型求解

根据问题需要,CHEMKIN读取模型输入文件确定求解方法。在4.0及其以上版本的CHEMKIN中,新加入了通过相应的模型设置窗口中设置模型的参数的功能。用户设置完成后点击Create Input File即可生成模型输入文件,然后用户可以通过View Input File按钮查看输入文件的内容。而对于4.0以前的版本,用户需要通过手动编辑模型输入文件。但对于文件的内容,都是相同的,都是采用关键字的形式声明模型的功能调用和参数设置。有关常用的平衡计算模型和全混反应器模型的关键字列于附录2、3。

完成所有设置后,即可Run Model进行计算,CHEMKIN会自动调用上述Gas-Phase Kinetics、Surface Kinetics、Transport各自的子程序来读取反应信息,调用模型输入文件控制模型求解过程。程序计算结束后,会生成一数据文件xxxxxx.out供用户查阅数据,以及一动态连接文件XMLdata.zip供后处理(绘图)使用。

2.2.4 后处理(Post-Process)

CHEMKIN提供了统一的后处理器,用于对应用程序的结果进行分析和绘图。

第三章 用户操作步骤

下面结合实际例子介绍用户操作的步骤。

3.1 例1定压绝热燃烧

【例1】计算化学当量的H2与空气的定压绝热燃烧温度。

3.1.1 运行CHEMKIN

1.运行CHEMKIN,点击菜单Project->New,输入项目名称,这里我们定为H2-air。

3.1.2 选择反应模型

2. 决定问题的性质,选择适当的反应模型。由于绝热燃烧温度的计算是一个相平衡过程,不涉及具体的反应过程,于是我们很容易地想到用Chemical and Phase Equilibrium

6

Calculation模型。点击任务栏中的model选项卡,点击Chemical and Phase Equilibrium Calculation图标,此时在窗口栏的Diagram View窗口中将看到新加入一个相平衡计算的模型,最后点击窗口右下角黄色的Update Project按钮,如下图所示。

通常情况下,一般的化学反应问题通过适当的假设和简化,都可以对应到某一种CHEMKIN包含的反应模型,有时一个问题还可以有多种选择。选择恰当的应用程序是求解问题的第一步。但是如果实在不幸没有一种模型可以很好的解决你的问题,或者你想要更完美地解决一些问题,就只有自己编写程序,调用CHEMKIN里的子程序库进行计算了。

3.1.3 前处理

3.此时左侧任务栏会自动切换至Open Projects选项卡,双击Pre-Processing,窗口栏出现Pre-Processing的参数窗口。在窗口中的Working Dir一项中填入你希望的保存路径,或者通过右侧的Browse按钮点击选取。然后按New Chemistry Set按钮,点击Gas-Phase Kinetics Files项右端的编辑按钮,如下图所示。

7

在弹出的窗口中选择刚才的工作路径,输入文件名chem.inp,按Open/Create按钮后即可开始编辑气相动力学输入文件了。

↓--------------------------------------【CHEMKIN的文件规则】--------------------------------------↓

CHEMKIN的输入文件有它自己的规则,用户在编辑输入文件的时候应该遵守这些规则。在介绍Gas-Phase Kinetics输入文件之前,先介绍输入文件的一些通用规则:

① 注释符号“!”。符号“!”无论出现任何位置,此行后面的文本将作为注释文本而被忽略。

② 输入文件每行不应超过80个字符

③ 除了个别有极其严格规则的地方外(如热力学数据的定义等,均会特殊声明),空格作为分隔符,而且多个空格将被视为一个。

④ 数字格式:可以为整数(如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 1 0.02882730E+02 0.10139743E-02-0.02276877E-05 0.02174683E-09-0.05126305E-14 2 0.03886888E+05 0.05595712E+02 0.03637266E+02 0.01850910E-02-0.16761646E-05 3 0.02387202E-07-0.08431442E-11 0.03606781E+05 0.13588605E+01 4 END

REACTIONS ! 反应方程及Arrhenius系数 H2+O2=2OH 0.170E+14 0.00 47780 OH+H2=H20+H 0.117E+10 1.30 3626 O+OH=O2+H 0.400E+15 -0.50 0 O+H2=OH+H 0.506E+05 2.67 6290 2OH=O+H2O 0.600E+09 1.30 0 H+H+M=H2+M 0.100E+19 -1.00 0

H2O/0.0/ H2/0.0/ ! 辅助数据 H+H+H2=H2+H2 0.920E+17 -0.60 0 H+H+H2O=H2+H2O 0.600E+20 -1.25 0 END

8

① 元素(Elements)

[规则]此部分以ELEMENTS(或者ELEM,两者等价)关键字开头;其后以空格为间隔符列出将在反应中出现的所有元素;最后以END关键字结束。例:

ELEMENTS H O END

! 元素定义

用户要按周期表的元素名(两个字母均须大写)来定义元素,CHEMKIN可以辨认的元素如表2.1。如果用户定义同位素或新元素的话,可以以1~2个字母命名(不与元素表重名),并将其原子量以“/”括住跟随其后。如定义氢(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, 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为电子,如果化学反应中有离子参加,电子必须作为一个元素进行定义。

② 组分(Species)

[规则]此部分以SPECIES(或SPEC)开头;其后以空格为分隔符列出将在反应中出现的所有组份;最后以关键字END结束。例:

SPECIES H2 H O2 O OH H2O END ! 组分定义

对于我们所要用到的组分,都可以从CHEMKIN的热力学数据库therm.dat找到。CHEMKIN数据库自带的组分名称(778种)列入附录1中。我们要使用这些组分的数据,组分的命名就必须按照CHEMKIN的规则来,这样CHEMKIN才可以从数据库中自动获取该组分的元素组成及其热力学性质。

当然,我们也可以用不超过16个字符而且以字母开头的任意字符串来定义自己的组分,但组分的热力学数据同样需要自行定义。

③ 热力学数据(Thermodynamic Data)

[规则]此部分以THERMO开头;其后以空格为分隔符列出定义组份的名称和热力学数据;最后以关键字END结束。例:

THERMO ! 热力学数据(本例中只重新定义了“OH”的热力学数据)

OH 121286O 1H 1 G 0300.00 5000.00 1000.00 1 0.02882730E+02 0.10139743E-02-0.02276877E-05 0.02174683E-09-0.05126305E-14 2 0.03886888E+05 0.05595712E+02 0.03637266E+02 0.01850910E-02-0.16761646E-05 3 0.02387202E-07-0.08431442E-11 0.03606781E+05 0.13588605E+01 4 END

CHEMKIN的热力学数据是基于CP、H、S的多项式模拟。这里不作介绍了,有兴趣的同学可参看帮助文件。

000④ 化学反应(Reaction)

[规则]起始行:关键字REACTIONS (或REAC),其后为Arrhenius系数的单位(可选)。

中间行:反应方程式,然后为该方程的Arrhenius系数(依次为Ai,βi,Ei);

有些反应需要辅助数据补充说明,出现在该反应方程式的下一行。

9

结束行:关键字END。例:

REACTIONS H2+O2=2OH 0.170E+14 OH+H2=H20+H 0.117E+10 O+OH=O2+H 0.400E+15 O+H2=OH+H 0.506E+05 2OH=O+H2O 0.600E+09 H+H+M=H2+M 0.100E+19

H2O/0.0/ H2/0.0/ H+H+H2=H2+H2 0.920E+17 H+H+H2O=H2+H2O 0.600E+20 END

?Arrhenius定律:Kfi?AiTiexp(! 反应方程及Arrhenius系数 0.00 47780 1.30 3626 -0.50 0 2.67 6290 1.30 0 -1.00 0 ! 辅助数据 -0.60 0 -1.25 0

?Ei) RcT1)在起始行内可以定义Arrhenius系数Ai和Ei的单位: Ei的单位可以定义为CAL/MOLE、 KCAL/MOLE、JOULES/MOLE、KJOULES/MOLE、KELVINS、EVOLTS;

Ai的单位可以用关键字MOLES或MOLECULES来定义,分别对应cm?mole?sec?K和cm?molecules?sec?K。

如果没有定义,Ei和Ai的默认单位是:cal/mole和cm?mole?sec?K。

2)反应方程中,“=”和“<=>”用于可逆反应(两者等价),“=>”用于不可逆反应,“+M”表示催化剂,“(+M)”用于压力控制反应,“HV”表示光子,“E”表示电子。

3)有些反应方程之后需要有补充的辅助数据说明。CHEMKIN有很多的辅助数据说明,有兴趣的可参看帮助手册。这里只介绍两个比较常用的——FORD和RORD,格式为

FORD /组分 反应级数/ RORD /组分 反应级数/

前者是重新定义正向反应时某组分反应级数,后者是重新定义逆向反应时某组分的反应级数的。

↑----------------------------------------------------------------------------------------------------------------↑

对于本例,由于绝热燃烧温度不考虑具体的反应过程,所以不涉及具体的反应方程,故Reaction部分可以省略(即使写了也用不上),此时CHEMKIN将根据反应物和生成物的焓(定压时)或内能(定容时)进行平衡计算。同时由于涉及的反应物和生成物的热力学数据都已在CHEMKIN的热力学数据库中有定义,故气相动力学输入文件只需要输入元素和组分两部分,如下所示:

ELEM H O N END

SPECIES O O2 H2 H OH HO2 H2O H2O2 N2 END

保存文件后关闭编辑器,此时Pre-Processing窗口中的Gas-Phase Kinetics Files一项仍为空白,通过点击该项的浏览按钮选取刚才编辑的输入文件。

然后需要在Thermodynamic Data File项中指定热力学数据库,点击该项的浏览按钮。在弹出的窗口中先点击右边的Special Directory中的System Data,然后点击左边出现的therm.dat,按select按钮完成热力学数据库的选择,如下页图所示。

回到Pre-Processing窗口后点击Save As…按钮,在弹出的窗口中直接点击Save按钮,以默认的文件名和路径保存。由于本例不涉及表面反应,所以不需要表面动力学输入和气体传递数据文件,故再次返回Pre-Processing窗口后即可点击Run Pre-Processor按钮运行预处理了。预处理的结果可以在下拉式按钮View Results…中查看。

10

3.1.4 边界条件

4. 双击任务栏中的Cluster1,会弹出子菜单,双击其中的每项都会弹出相应的设置窗口,其中:

① 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两个选项卡,分别用来设置反应条件和反应物的组成,是重点的参数输入部分。对于本例要解决的化学当量的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 fraction;

c)Species下拉式按钮中依次选择和输入:H2,2,Add按钮;O2,1,Add按钮;N2,3.76,Add按钮。

通过以上步骤完成反应条件和反应物的输入。

③ Solver:设置问题的求解方法。在本例中没有可供设置的参数。对于其它问题,稳态计算可能需要设置迭代的步长、迭代次数等,瞬态计算可能需要设置反应结束的时间、解的误差、灵敏度系数的误差等。可以参看下文中对本问题用不同模型求解时的例子。

④ Output Control:设置输出文件的大小限制、输出的灵敏度结果的个数等。对本例,我们不需做任何设置,均采用默认值即可。

⑤ Continuations:设置重复计算的次数,可以为0、1、2、3??。当为0时,表示不作重复计算。对于本例,我们不需要作重复计算,输入0即可;或者直接不点击Continuations该选项,使用其默认值0。

11

3.1.5 运行计算

5. 双击任务栏的Run Model,在右侧窗口中点击Create Input File按钮,弹出提示输入文件名的窗口,直接按OK采用默认文件名即可。输入文件建立成功,此时可通过View Input File查看输入文件的内容,但请勿手动修改其中的内容,否则模型可能无法运行。如需要手动增加输入参数,请使用Add Supplemental Input按钮,在此不作介绍。

点击Run Model按钮,消息栏提示模型运行成功。如下图所示。

点击View Results下拉式按钮中的Output File即可查看计算结果,本例中得到平衡状态的温度为2.3892E+03(K),此即所求的定压绝热燃烧温度。

我们可以尝试点击Run Post Processor按钮,在弹出窗口中选择要查看的的结果及其单位,按OK后发现弹出一个Warning的消息窗口,提示计算结果只有一组数值,不满足作图要求(至少两组数值),得到的图形也只是作为示例的Cos(x)、Sin(x)图形。所以这里所要说明的是,使用后处理器(Post Processor)绘图需要有多个状态、多个时刻、多个计算次数的计算结果,一般用于瞬态燃烧时各参数随时间的变化的绘图,这在下文可以看到。

思考与实践:如果反应物的组分是化学当量比的H2和O2,定压绝热燃烧温度和定容绝热燃烧温度又该如何设置求解?结果如何?

3.2 例2全混燃烧器中的绝热燃烧

【例2】计算化学当量比的H2和空气在全混燃烧器中的绝热燃烧过程

3.2.1 运行PSR

1. 建立新Project,选择Perfectly Stirred Reactor(PSR)全混反应器模型,Update Project。

12

3.2.2 前处理

2. Pre-Processing中设工作路径,点New Chemistry Set,输入气相动力学文件chem.inp如下

ELEMENTS O H N END

SPECIES O O2 H H2 OH HO2 H2O H2O2 N2 END REACTIONS

H+O2 = O+OH 5.1E16 -0.82 16510 !MILLER 81 H2+O = H+OH 1.8E10 1.0 8830 !MILLER 77

H2+OH = H2O+H 1.2E09 1.3 3630 !DIXON-LEWIS 77 OH+OH = H2O+O 6.0E08 1.3 0.0 !COHEN 79 H+OH+M = H2O+M 7.5E23 -2.6 0.0 !MILLER 77 H2O /20.0/

O2+M = O+O+M 1.9E11 0.5 95560 !MILLER 81 H2+M = H+H+M 2.2E12 0.5 92600 !MILLER 77 H2O /6.0/ H /2.0/ H2 /3.0/

H2+O2 = OH+OH 1.7E13 0.0 47780 !MILLER 77 H+O2+M = HO2+M 2.1E18 -1.0 0.0 !SLACK 77 H2O /21.0/ H2 /3.3/ O2 /0.0/ N2 /0.0/

H+O2+O2 = HO2+O2 6.7E19 -1.42 0.0 !SLACK 77 H+O2+N2 = HO2+N2 6.7E19 -1.42 0.0 !SLACK 77 HO2+H = H2+O2 2.5E13 0.0 700 !LLOYD 74 HO2+H = OH+OH 2.5E14 0.0 1900 !LLOYD 74 HO2+O = OH+O2 4.8E13 0.0 1000 !LLOYD 74 HO2+OH = H2O+O2 5.0E13 0.0 1000 !LLOYD 74 HO2+HO2 = H2O2+O2 2.0E12 0.0 0.0 !TROE 69 H2O2+M = OH+OH+M 1.2E17 0.0 45500 !BAULCH 72 H2O2+H = HO2 + H2 1.7E12 0.0 3750 !BAULCH 72 H2O2+OH = H2O+HO2 1.0E13 0.0 1800 !BAULCH 72 END

选择热力学数据库文件,由于本问题不涉及表面反应,所以表面动力学文件和传递数据库文件位置仍然留空。Run Pre-Processor进行预处理。

3.2.3 对Cluster1中的设置

3. 对Cluster1中的设置如下:

① Cluster Properties:选择默认值Normal Start。

② C1_R1 PSR:也包括Reactor Physical Property和Species-specific Data两个选项卡,但应注意这里的Species-specific Data是指反应器中预先充入的气体组分,而并不是反应物的组分,反应物的组分应在R1_IN1中设置,所以这里的Species-specific Data应留空。

对Reactor Physical Property选项卡:

a)Problem Type中选择Solve Gas Energy Equation,即通过求解气体能量方程获得生成物的温度;

b)点击Transient Solve,因为是求解燃烧过程各参数随时间的变化情况; c)点击Temperature,输入1000,单位K,指反应初始的温度为1000K; d)点击Pressure, 输入1,单位atm;

③ R1_IN1:包括Stream Property Data和Species-specific Property两个选项卡,设置入口物质的流量、温度以及反应物的组成。本例中,Stream Property Data留空,Species-specific Property可以通过以下两种方式的其中之一设置反应物的组成:

a)不勾选Equivalence Ratio,在Reactant Fraction选择组分并输入数据H2,2;O2,1;N2,3.76;

13

b)勾选Equivalence Ratio,输入数值1。在Fuel Fraction中输入H2,2;在Oxidizer Fraction中输入O2,1;N2,3.76;在Combustion Product中选择H2O,N2;

④ Solver:在End Time中输入计算结束的时刻0.001sec; ⑤ Output Control和Continuations均留空,不作设置

3.2.4 运行计算

4. 在Run Model中,依次点击Create Input File、Run Model按钮,提示运行成功,可点击View Results…中的Output File查看计算结果,这里我们点击Run Post Processor(后处理器)进行绘图观察,在弹出的窗口中Species/Variables按钮可以选择哪些参数要输出绘图,Units of Measure按钮可以选择各种物理量的计量单位。选择完毕后按OK,就会出现计算结果的曲线,此时选择菜单中的Plot->X-Y Plot,从弹出的窗口中选择曲线的X轴和Y轴(可多选)。选择菜单Plot->Contour可绘制等高线的变量,但要求解结果必须是二维的(如圆柱和平板的层流反应模型)。我们选择温度和时间为坐标轴绘制图形,得到如下曲线:

思考与实践:如果反应初始温度为1200K、1500K,计算结果将如何变化?如果反应初始温度过低,会出现什么情况?

本文作为一本入门教程,旨在使大家能够迅速了解和掌握CHEMKIN的基本功能和应用,所以很多复杂的高级技巧只做了一些简单介绍或者完全忽略。CHEMKIN的帮助系统提供对程序所有功能的详细介绍,请自行参阅。

CHEMKIN是一个非常好的计算工具,但只有深入了解化学反应和燃烧理论,才能真正发挥它的作用。

【几点使用中的小技巧】

1、注意看提示信息,包括将鼠标移动至相应项时弹出的提示信息、输出文件中的信息、弹出窗口的信息、消息栏的信息等;

2、注意加粗字体的项,这些一般都是模型正常运行所必须要求填写的,或者至少需要填写的;

3、灵活运用文件选择窗口中右侧Special Directory的选项快速定位目录,可以很大地提高工作效率。

14

附录1 CHEMKIN热力学数据库中包含的组分

(CH2O)3 (CH2O)4 AL AL2H6 AL2ME6 ALAS ALH ALH2 ALH3 ALME ALME2 ALME3 AR AR+ AS AS2 AS3 AS4 ASALME ASALME2 ASGAET ASGAET2 ASGAME ASGAME2 ASGAMEH ASH ASH2 ASH3 ASME ASME2 ASME3 B B(S) BCL BCL2 BCL3 BE BE(S) BE2SIO4(S) BE3B2O6(S) BE3N2(A) BE3N2(L) BEAL2O4(S) BEB2O4 BEBO2 BEBR BEBR2 BEBR2(S) BECL BECL2 BECL2(A) BECL2(B) BECLF BEF BEF2 BEF2(L) BEH BEH+ BEH2 BEH2O2 BEH2O2(A) BEH2O2(B) BEI BEI2 BEI2(S) BEN BEO(A) BEO(B) BEOH BES(S) BESO4(A) BESO4(B) BESO4(GAM) BN C C(S) C+ C- C2 C2- C2CL3 C2CL5 C2CL6 C2F6 C2H C2H2 C2H3 C2H4 C2H5 C2H6 C2HCL C2HCL5 C2N C2N2 C2O C3 C3H2 C3H2(S) C3H4 C3H4C C3H4P C3H6 C3H8 C3O2 C4 C4H C4H10 C4H2 C4H6 C4H8 C5 C5H C5H12 C5H2 C5H5 C5H6 C6H C6H10 C6H14 C6H2 C6H3 C6H4 C6H5 C6H5(L) C6H5O C6H5OH C6H6 C6H7 C8H C8H2 CA CA(A) CA(B) CA(L) CA+ CA2 CABR CABR2 CABR2(S) CACL CACL2 CACL2(S) CAF CAF2 CAH2O2 CAH2O2(S) CAI CAI2 CAO CAO(S) CAOH CAOH+ CAS CCL CCL2 CCL2CCLO CCL2CCLOH CCL2CH CCL2HOO CCL2OHCH2 CCL2OHCHCL CCL3 CCL3CCLH2 CCL3CCLO CCL3CH2 CCL3CHCL CCL3CHO CCL3OO CCL4 CCLH2OO CH CH+ CH2 CH2(S) CH2CCL CH2CCL2 CH2CCLOH CH2CHCCH CH2CHCCH2 CH2CHCH2 CH2CHCHCH CH2CHCHCH2 CH2CHCL CH2CL CH2CL2 CH2CLCCL2 CH2CLCCLO CH2CLCH2 CH2CLCH2CL CH2CLCHCL CH2CLCHCL2 CH2CLCHO CH2CO CH2F2 CH2HCO CH2O CH2OH CH2OHCCL2 CH2OHCHCL CH2SICL CH2SICL3 CH2SIH2CL CH2SIHCL2 CH3 CH3C(O)CL CH3CC CH3CCCH2 CH3CCCH3 CH3CCH2 CH3CCL CH3CCL2 CH3CCL3 CH3CCLO CH3CH2CCH CH3CH2CH2CH3 CH3CH2CL CH3CH2O CH3CHCH CH3CHCL CH3CHCL2 CH3CHOH CH3CL CH3CO CH3F CH3HCO CH3NO CH3NO2 CH3O CH3OCH3 CH3OCL CH3OH CH3ONO CH3ONO2 CH3SICL CH3SIH2SIH CH3SIH2SIH2CH3 CH3SIHCL2 CH4 CHCL CHCL2 CHCL2CCL2 CHCL2CCLO CHCL2CH2 CHCL2CHCL CHCL2CHCL2 CHCL3 CHCLCCL CHCLCCLOH CHCLCH CHCLCHCL CHCLCHOH CHCLOH CHCLOHCH2 CHCLOHCHCL CHF CHF3 CHOHCLCCL2 CHSICL CHSICL2 CHSICL3 CHSIH2CL CHSIHCL CHSIHCL2 CL CL(CH3)SICH2 CL2 CL2CCCL2 CL2CCHCL CL2CCHO CL2CHOH CL2CO CL2COH CL2HCO CL2SI(CH3)2 CL2SI(CH3)CH2 CL2SICH2 CL2SICH3 CL2SISI CL2SISICL CL2SISICL2 CL3CCO CL3CO CL3COH CL3SICH3 CL3SISI CL3SISICL CLCCCL CLCCO CLCH2OH CLCO CLCOH CLH2CO CLHCO CLO CLOCL CLOO CLSI(CH3)2 CLSI(CH3)2CH2 CLSI(CH3)3 CLSICH3 CLSISI CLSISICL CN CN+ CN- CN2 CNN CNO CO CO2 CO2- COS CS CS2 CSICL CSICL2 CSICL3 CSIH2CL CSIHCL CSIHCL2 D D2 DH DIOXANE E F F- F2 F2N2(C) F2O2 F2SINH F3SIN FNNF FNO3 FO FO2 FONO(C) FONO(T) FSIN GA GA2H6 GAAS \\\\GAET GAET2 GAET3 GAH GAH2 GAH3 GAME GAME2 GAME3 H H(CH3)SICH2 H+ H- H2 H2ALME H2ASME H2C4O H2CCC H2CCC(S) H2CCCCCH H2CCCCH H2CCCCH2 H2CCCH H2CCCLO H2CCH(SICL2H) H2CCH2OH H2CCHO H2CCHSI H2CCHSIH H2CCHSIH2 H2CCHSIH3 H2CLSICH3 H2CLSISICL3 H2CN H2CNCH2 H2CNCH2O H2CNCHO H2CNH H2CNNHO H2CNNO H2CNNO2 H2CNO H2CNO2 H2CONO H2GAET H2GAME H2NF H2NNO H2NO H2NOH H2O H2O(L) H2O(S) H2O2 H2S H2SI(CH3)2 H2SI(CH3)CH2 H2SI(NH2)2 H2SIC H2SICCH H2SICH H2SICH2 H2SICH3 H2SIN H2SINH H2SINH2 H2SINH3 H2SISIH2 H3ASGAET3 H3ASGAME3 H3CONHO H3SIC H3SICCH H3SICH H3SICH2 H3SICH3 H3SIN H3SINH H3SISIH H3SISIH3 HALME HALME2 HASALME HASGAET HASGAME HASME HASME2 HCCCHCCH HCCCL HCCHCCH HCCO HCCOH 15

HCCSICL2H HCL HCL2SICH3 HCL2SISICL2H HCLCCCLO HCLCCHO HCLSI(CH3)2 HCLSI(CH3)CH2 HCLSICH2 HCLSICH3 HCLSISI HCN HCNH HCNO HCO HCO+ HCOOH HE HE+ HF HG2BR2(S) HG2CL2(S) HG2F2(S) HG2I2(S) HGAET HGAET2 HGAME HGAME2 HGBR HGCL2 HGCL2(S) HGF2 HGF2(S) HGH HGI HGO HGO(S) HMEGAET HN(OH)2 HN3 HNC HNCN HNCNH HNCO HNF HNF2 HNNHO HNNO HNNONO HNO HNO2 HNO3 HNOH HO2 HOCH2OH HOCL HOCN HOCO HONO HONO2 HSI(CH3)2 HSI(CH3)2CH2 HSI(CH3)3 HSI(NH2)2 HSI(NH2)3 HSIC HSICCH HSICH2 HSICH3 HSICL HSIN HSINH HSINH2 HSISICL I*C3H7 K K(L) K+ K2 K2B4O7(S) K2B6O10(S) K2B8O13(S) K2CO3(S) K2H2O2 K2O(S) K2O2(S) K2SO4 K2SO4(A) K2SO4(B) K3CL6AL(S) K3CL9AL2(S) K3F6AL(S) KBF4 KBO2 KBO2(S) KBR KBR(L) KBR(S) KCL KCL(L) KCL(S) KCL4AL(S) KCLO4(S) KCN KCN(L) KCN(S) KF KF(L) KF(S) KH KH(S) KI KI(L) KI(S) KO KO- KO2(S) KOH KOH+ ME2GAET MEGAET MEGAET2 MG MG(L) MG(S) MG+ MG2 MG2BR4 MG2C3(S) MG2F4 MGAL2O4(S) MGB2(S) MGBR MGBR2 MGBR2+ MGC2(S) MGCL MGCL2 MGCL2(S) MGCO3(S) MGF MGF2 MGF2(S) MGF2+ MGH MGH2(S) MGH2O2 MGH2O2(S) MGN MGO(S) MGOH MGOH+ MGS MGS(S) MGSO4(S) N N*C3H7 N2 N2F2(C) N2F2(T) N2F4 N2H2 N2H3 N2H4 N2H4(L) N2O N2O+ N2O4 N3 NA NA(L) NA+ NA2 NA2B4O7(S) NA2B6O10(S) NA2C2N2 NA2F2 NA2H2O2 NA2O(L) NA2O(S) NA2O2 NA2SO4 NA2SO4(D) NA2SO4(i) NA2SO4(iii) NA2SO4(iv) NA2SO4(v) NA3CL6AL(S) NA3F6AL(S) NABH4(S) NABO2 NABO2(S) NABR NABR(L) NABR(S) NACH NACL NACL(S) NACL4AL(S) NACN NACN(S) NAF NAF2- NAF4AL NAH NAI(S) NAO NAO- NAO2(S) NAO2AL(S) NAOH NAOH(L) NAOH+ NCN NCO NF NF2 NF3 NFO NFO2 NH NH2 NH3 NNH NO NO+ NO2 NO2- NO2F NO3 NO3F NOF NOF3 O O+ O- O2 O2- O2F O3 OC(OH)2 OCHCHO OCHNNHO OF OF2 OH OH+ OH- ONHNHO ONHNOH P P2 P4 S S(L) S(S) S+ S-TRIAZINE S2 SH SI SI(CH3)2 SI(CH3)3 SI(CH3)3CH2 SI(CH3)4 SI(L) SI(NH2)3 SI(NH2)4 SI(S) SI2 SI2C SI2CL5 SI2CL5H SI2CL6 SI2F6 SI2H2 SI2H3 SI2H5 SI2H6 SI3 SI3H8 SI3N4(A) SIC SIC(B) SIC2 SICCH SICH SICH2 SICH3 SICL SICL2 SICL2H2 SICL3 SICL3CH2CH SICL3H SICL4 SICLH3 SIF SIF(NH2)2 SIF(NH2)3 SIF2 SIF2(NH2)2 SIF2N SIF2NH2 SIF3 SIF3NH SIF3NH2 SIF3NHSIH3 SIF3NSIH3 SIF4 SIFH3 SIFNH SIFNH2 SIH SIH2 SIH2CL SIH2F SIH2F2 SIH3 SIH3NH2 SIH3NHSIH3 SIH3NSIH3 SIH3SIH2CH3 SIH4 SIHCL2 SIHF SIHF2 SIHF3 SIN SINH SINH2 SN SO SO2 SO3 TI TICL TICL2 TICL3 TICL4 附录2 化学平衡和相平衡计算的模型输入文件常用关键字

类型 问题类型 (必选其一) 关键字 TP或PT TV或VT TS或ST PV或VP HP或PH PS或SP VU或UV VH或HV VS或SV CJ 参数 描述 定压定温条件 定容定温条件 定温定熵条件 定压定容条件 定压定焓条件 定压定熵条件 定容定内能条件 定容定焓条件 定容定熵条件 Chapman-Jouguet爆炸, 16

FREE FROZ FAZE CONX 初始条件 TEMP PRES VOL ENTH ENTR ENGY REAC 平衡估计 XML解文件选项 TSTR DIST USET SIZE CNTN T S N NSOL N TEST PEST XMLI A N N N N N N N A N N N 组分需要进行平衡计算 组分成分固定,与FREE相对 相成分恒定,如气相分子不会跑到表面去 指定组分A的摩尔分数固定为N,保持不变 起始温度为N(K) 起始压力为N(atm) 起始体积为N(cm3/g) 起始焓为N(erg/g) 1 erg=1 cm*g*s 起始熵为N(erg/gK) 起始内能N(erg/g) 反应物摩尔分数;组份A的摩尔分数为N 平衡温度估计值:N(K) 平衡压力估计值:N(atm) 从解文件XMLrestart.zip中读取初始条件。这一功能可用于处理其他应用程序(如CRESLAF、OPPDIF、SPIN等)的计算结果。 当解文件XMLrestart.zip中包含多个解时,NSOL决定采用哪个解来作为初始条件。 当解文件XMLrestart.zip中为瞬态解(如AURORA的瞬态解),TSTR规定采用T(sec)时刻的解作为初始条件 当解文件XMLrestart.zip中的解随轴向距离变化(如PLUG的解),DIST规定采用S(cm)距离处的解作为初始条件 当XMLI或restart(重启)时,压力和温度的值从输入文件中读入,取代从解文件XMLrestart.zip中读取的压力和温度 指定EQUIL产生的解文件XMLdata.zip的大小N(bytes) 继续计算关键字。在一次计算之后,可以改变某些条件进行下一次计算 约束条件 附录3 全混反应器的模型输入文件常用关键字

全混反应器指气体扩散非常快或经过充分搅拌,使反应器内气体分布处处均匀。全混反

应器是最简单的模型,在燃烧中有着广泛的应用,而且是测试化学反应的常用工具。其模型输入文件的常用关键字见下表。 类型 关键字 参数 含义 I、问题属性(Problem Type) 反应特性 是否求解能量方程 仅用于封闭系统 TRAN STST NPSR TGIV ENRG SOFT CONV CONT ICEN N 瞬态反应器(组分随时间变化) 稳态反应器 包含N个反应器的网络:最大个数为10 不求解能量方程,使用用户定义的温度(由TEMP定义) 可以只针对反应器网络中的某个反应器定义 求解能量方程来获取温度,TEMP定义的温度作为迭代初始温度。 求解能量方程式,忽略计算得到的中间温度,从用户给定的温度开始迭代 定容定焓瞬态系统 定温定容瞬态系统 内燃机瞬态模型。 17

II、瞬态反应器(Transient Solution) TIME TLIM 时间 STPT DELT DTSV SSTT TPRO 随时间变化量 PPRO VPRO QPRO WPRO CMPR VOLC RPM 内燃机模型 LOLR DEG0 DTDEG NREV ICHT ATOL 瞬态解容差 RTOL ATLS RTLS NNEG SDIR 重启和继续 TRST TRES USET N N N CNTT N T N N N N A B C A B C A B C A B C A B C N N N N N N N A,B,C,T N N N N 瞬态反应总时间N(sec) 点火时间(点火温度)。当温度达到T(K)时点火。 时间步长N(sec) 将结果输出为文本文件aurora.out的时间间隔N(sec) 将结果输出为解文件XMLdata.zip的时间间隔N(sec) 灵敏度分析的时间间隔N(sec) 第C个反应器在A(sec)时刻反应器温度为B(K),无C时应用于所有反应器 第C个反应器在a(sec)时刻反应器压力为B(atm) 第C个反应器在a(sec)时刻反应器容积为B(cm2) 第C个反应器在a(sec)时刻反应器放热量为B(cal/sec) 第C个反应器在A(sec)时刻等离子体反应器功率为B(watts) 压缩比为N,(最大体积与压缩后体积之比) 压缩后体积为N(cm3) 曲柄转速为N(rpm) 连接杆长度与曲柄半径之比为N 瞬态模型中曲柄起始角度N(degrees) 曲柄角度形式的时间步长N(degrees) 曲柄总转数N,可以取代TIME 瞬态内燃机模型中的对流换热关联式 绝对误差 相对误差 灵敏度系数的绝对误差 灵敏度系数的相对误差 保证所有解向量非负 使用参差直接法求解灵敏度方程 继续计算时,开始时间等于上次计算的结束时间。未定义此关键字,则继续计算时开始时间从0开始 重启计算时,初始值为Restart文件中N(sec)时刻的值。未定义此关键字,则初始值为Restart文件中最后时刻的值。 重启计算中定义新的初始时间N(sec) 重启计算中温度和压力从输入文件读取,而不是从Restart文件中读取 III、稳态反应器(steady solution) ATOL RTOL 容差分析 ATIM RTIM ABSL RELT 求解方法TIM1 N N N N N N A B 绝对容差 相对容差 时间步进分析时的绝对容差 时间步进分析时的相对容差 解变量的绝对摄动的值 解变量的相对摄动的值 如果牛顿法不收敛,将迭代A次,初始步长为B,用于固定温度情况(TGIV) 18

控制 TIM2 DTMN DTMX IRET UFAC DFAC NFAC TFAC ISTP SFLR A B N N N N N N N N N 同TIM1,应用于求解能量方程情况(ENRG) 最小时间步长 最大时间步长 TWOPNT每当时间步长数达到N时增加时间步长 TWOPNT时间步长达到IRET时,步长为原来的N倍 当当前步长不收敛时,步长为原来的1/N TWOPNT牛顿法中新雅可比矩阵计算前迭代的次数 TWOPNT时间步长中新雅可比矩阵计算前迭代的次数 牛顿积分前初始时间步长数目 组分分数变量中的最小值 IV、反应器条件和入口条件 循环 RECY INLET TINL REAC EQUI 入口条件 FUEL OXID PROD ADD TEMP PRES TAU FLRT SCCM VOL AREA 反应器条件 QLOS HTRN AREA AFRA TSRF XEST A B C Name A Name A Name Spec A Name A Name Spec A Name Spec A Name Spec A Name Spec A A B A B A B A B A B A B A B Material A B Material H T B A Material A B Material A B A B 第A个反应器流出量的C倍(C<=1)流回第B个反应器 定义第A个反应器的入口气流名字为Name 定义名为Name的入口气流温度为A(K) 定义名为Name的入口气流中Spec组分的含量为A 入口气流Name的空燃比为A; 入口气流Name的燃料中Spec组分的含量为A 入口气流Name的氧化剂中Spec组分的含量为A 入口气流Name的产物中Spec组分的含量为A 将不参加化学反应的组分加入到入口气流中 入口气流Name中Spec组分的含量为A 第B个反应器气体温度A(K);无B时针对所有反应器 第B个反应器内压力A(atm);无B针对所有反应器 第B个反应器气体停留时间A(sec);无B针对所有 质量流率A(g/sec),其余同上 体积流率A(cm3/min),其余同上 体积A(cm3) 总表面积A(cm2) 第B个反应器的Material(代指material的名字,下同)上放热为A(cal/sec) 第B个反应器的Material上的传热系数为H (cal/ (cm2*K*sec)),环境温度为T(K) 总表面积为A(cm2) 第B个反应器名的Material上表面积占总表面积的B倍。 第B个反应器的Material表面温度为A(K) 初始或起始气相摩尔分数 19

SURF BULK SDEN A B A B A B 初始或起始表面摩尔分数 初始或起始体内摩尔分数 表面相中表面site A的密度 V、灵敏度和生成率分析 ASEN SEN 灵敏度 SENT SENG EPSS EPST EPSG 生成率分析 AROP ROP EPSR N A?? N N N N N N A?? N 计算所有组分和气体温度的一次灵敏度系数 计算A等组分的一次灵敏度系数 计算气相温度的一次灵敏度系数 计算体内相增长率的一次灵敏度系数 组分一次灵敏度系数的极小值 气体温度一次灵敏度系数的极小值 体内相增长率一次灵敏度系数的极小值 所有组分计算生成率 计算A等组分的生成率 生成率系数极小值 VI、解文件选项 XMLI NSOL TSTR DIST USET USEV N T S 从解文件XMLrestart.zip中读取数据作为初始值。常用于不同应用程序之间的交互 读取解文件中的第几组解 读取瞬态解文件中时间T时刻的值 如果解文件中的解是轴位移的函数,此关键字决定读取解文件中位移S处的值 当XMLI或RSTR时,温度和压力仍从输入文件中读取 当XMLI或RSTR时,速度和流率仍从输入文件中读取 VII、其它控制 SIZE 输出控制 PRNT RSTR CNTN END N N 限定解文件XMLdata.zip不能超过N bytes 文本输出控制。N=0仅输出解,N=1还提供积分路径,N=3还包括每一段牛顿积分及时间步长 结束关键字 从重启解文件XMLrestart.zip中读取数据然后按当前反应器条件求解。 结束之前,有此关键字表明程序还须继续计算 一个条件的输入完毕关键字 附录4 与燃烧和热力学数据有关的网络资源

通常我们并不需要自己编写反应机理,很多反应机理都有人研究过,常用网络资源如下 1. Lawrence Livermore National Laboratory,网址为:

http://www-cms.llnl.gov/combustion/combustion2.html

提供许多燃料的燃烧机理,直接为CHEMKIN格式,方便易用。

2. 美国NIST(National Institute of Standard and Technology)的化学动力学数据库,网址为:

http://kinetics.nist.gov/index.php。

可以查询任意反应的反应速率,而且可以查到该数据出自那篇文献。 3. GRI (Gas Research Institute of University of California at Berkeley),网址为:

http://www.me.berkeley.edu/gri-mech/data/frames.html。

可以查到大部分反应的Arrhenius系数。该研究所对天然气燃烧的研究尤为深入。

20

本文来源:https://www.bwwdw.com/article/xrkd.html

Top