CN115798625B - 一种耦合单步化学机理的反应流数值求解方法 - Google Patents

Abstract

本发明公开了一种耦合单步化学机理的反应流数值求解方法，属于反应流数值求解领域。本发明采用化学扩散模型计算化学反应速率、热释放速率和扩散系数，建立单步机理反应流数学模型，同时将化学扩散模型耦合到压力修正方法中，顺序求解反应流中的各个变量场。本发明能够处理非稳项、对流项、扩散项和源项的耦合问题，实现非稳态、可压、多维度的均匀反应流和非均匀反应流计算；同时将化学扩散模型的应用场景拓展到低速流领域。相比于采用详细机理或骨架机理化学反应模型，本发明无需求解中间组分，能够更快速解决反应流数值计算领域的相关问题。

Description

一种耦合单步化学机理的反应流数值求解方法

技术领域

本发明属于反应流数值求解领域，具体涉及一种耦合单步化学机理的反应流数值求解方法的构建。

背景技术

反应流是流体流动与传热传质相互耦合的复杂物理化学过程，涉及各种不同的时间和空间尺度，广泛存在于自然界与人造装置中，如燃烧、爆炸等。反应流数值求解过程通常使用详细化学机理或骨架化学机理来建立化学反应模型，这通常需要求解大量的组分方程；此外，不同中间组分和基元反应之间的差异会增加控制方程的刚性，使反应流的数值求解更加消耗计算资源和计算时间。而化学扩散模型是一种新型高效的单步化学机理，对化学反应动力学机理和扩散特性（热扩散、质量扩散和粘性效应）同时进行了简化，具有计算准确、效率高、节省计算资源等优势，能很好地复现火焰和燃爆的主要特性。

目前化学扩散模型的应用较多集中在预混燃烧的情景中，而对非均匀混合和低速反应流的应用非常少。

发明内容

为克服现有技术的缺点，本发明提供了一种耦合单步化学机理的反应流数值求解方法，其为将化学扩散模型和压力修正方法耦合的反应流数值求解方法，并将其应用扩展到低速、非均匀反应流。本发明使用单步化学机理来建立化学反应模型，能够以较高的精度和更快的计算速率再现均匀反应流和非均匀反应流的主要特性。本发明实现了反应流的数值求解，可处理多维度、可压缩、非定常反应流的数值模拟，包括预混燃烧、扩散燃烧以及部分预混燃烧。本发明使用单步化学机理来建立化学反应模型，能够以较高的精度和更快的计算速率再现均匀反应流和非均匀反应流的主要特性。

为达到上述目的，本发明采用的技术方案为：

一种耦合单步化学机理的反应流数值求解方法，包括如下步骤：

步骤一，建立非定常、有粘、可压的多维度反应流数学模型，化学反应速率、热释放速率、扩散系数使用化学扩散模型来确定；

步骤二，结合初始条件、边界条件和网格划分情况，对步骤一建立的多维度反应流数学模型进行数值离散，得到包括密度方程、速度方程、组分方程、温度方程和压力方程的代数方程组，并对代数方程组进行线性化；

步骤三，将化学扩散模型耦合到压力修正方法中，化学反应速率和热释放速率的计算位于组分方程和温度方程之前，且均在预估校正迭代内部，进行多次循环以保证收敛；在包括速度方程、组分方程和温度方程的代数方程组求解之前更新扩散系数；通过步骤三的上述步骤，建立耦合单步化学机理的反应流求解方法。

进一步地，所述步骤一中，采用化学扩散模型来建立反应机理模型和扩散模型，用来求解化学反应速率和热释放速率。

进一步地，所述步骤一中，所述化学扩散模型的参数的确定采用查表法。

进一步地，所述步骤二采用有限体积法，在物理网格划分的基础上，将步骤一建立的密度方程、速度方程、组分方程、温度方程以及压力方程在网格内体积积分，并在时间步长内进行时间积分，通过引入不同的空间离散格式和时间离散格式，最终得到定义在该网格中心点上的代数方程。

进一步地，所述体积积分时，根据高斯定理将对流通量、扩散通量的体积分转化为面积分，然后根据不同的空间离散格式将定义在网格面的变量重构到网格点；在化学反应速率和热释放速率的体积积分过程引入变量在网格内均匀分布的假设；在所述时间积分的过程中，化学反应速率和热释放速率进行显式处理；所述时间离散格式包括显式、半隐式、隐式。

进一步地，所述步骤二中，所述压力方程的推导考虑密度的变化，实现可压缩流的数值求解。

进一步地，所述步骤三具体包括如下步骤：

步骤（1）确定全局时间步长，并进行时间推进；全局时间步长选取为化学时间步和对流时间步的较小值；采用全局隐式方法实现不同物理过程和化学过程之间的耦合，即非稳项、对流项、扩散项和源项同时进行离散化、线性化和求解；

步骤（2）求解离散化后的密度方程，得到密度预估值，用于预估校正迭代中的代数方程组的线性化；

步骤（3）求解离散后的速度方程，得到速度预测值；

步骤（4）根据化学扩散模型计算每个网格中的组分反应速率和热释放速率；

步骤（5）求解离散后的组分方程和温度方程，得到组分场和温度场的分布；

步骤（6）求解离散后的压力方程更新压力场，并将更新的压力场带入速度方程用于更新速度场，完成一次压力修正循环；达到指定的压力修正循环次数后跳出循环，进入步骤（7）；

步骤（7）重复步骤（3）到步骤（6），直到收敛或者达到最大的预估校正迭代次数；

步骤（8）求解理想气体状态方程，更新密度值；

步骤（9）判断是否到达指定的结束时间，否则重复步骤（1）到步骤（8）。

进一步地，所述步骤三中，在压力修正方法的预估校正迭代中，采用更新的密度、温度和组分质量分数，根据化学扩散模型的算法，计算各组分的化学反应速率和热释放速率，进而求解组分方程和温度方程，用以更新组分质量分数和温度，形成双向耦合。

进一步地，用于处理非稳项、对流项、扩散项和源项。

本发明与现有技术相比的优点为：

1. 将化学扩散模型耦合到纳维-斯托克斯方程组中，建立了单步机理反应流数学模型；

2. 采用单步反应机理来建立化学反应模型，无需求解中间组分，加快了计算速率；

3. 通过校准不同当量比下的化学扩散模型参数，可对富燃、贫燃等非均匀反应流进行计算，拓展了化学扩散模型的应用场景；

4. 实现了化学扩散模型和压力修正方法的耦合，解决了低流速下压力-密度弱耦合的问题，可实现低速反应流的数值求解。

附图说明

图 1 为本发明的耦合单步化学机理的反应流数值求解方法。

图 2 为本发明的实施例中二维氢/空气对向流扩散火焰的装置示意图。

图 3 为本发明的实施例中一阶欧拉隐式离散格式的基架点示意图。

图 4 为本发明的实施例工况一的氢/空气对向流扩散火焰的温度和组分廓线。

图 5 为本发明的实施例工况二的氢/空气对向流扩散火焰的温度和组分廓线。

具体实施方式

为了便于理解本发明，下文将结合实施例对本发明作更全面、细致的描述，但本发明的保护范围并不限于以下具体的实施例。

本发明在有限体积法的框架下将化学扩散模型与压力修正方法耦合，实现了反应流的数值求解，可处理非定常、可压缩、多维度反应流的数值模拟。

如图1所示，本发明的一种耦合单步化学机理的反应流数值求解方法具体包括如下步骤：

步骤一，建立非定常、有粘、可压的多维度反应流数学模型，其中化学反应速率、热释放速率、扩散系数使用化学扩散模型来确定；

步骤二，结合初始条件、边界条件和网格划分情况，对步骤一建立的多维度反应流数学模型进行数值离散，得到关于密度、速度、组分、温度和压力的代数方程组，并对代数方程组进行线性化；

步骤三，将化学扩散模型耦合到压力修正方法中；其中化学反应速率和热释放速率的计算位于组分方程和温度方程之前，且都在预估校正迭代内部，进行多次循环以保证收敛；此外，在速度、组分和温度代数方程组求解之前更新扩散系数；通过上述三个步骤，建立耦合单步化学机理的反应流求解方法。

具体地，如图1所示，所述步骤一中建立的多维度反应流数学模型包括以下子模型：

（一）非定常、有粘、可压缩纳维-斯托克斯方程组，其包含密度方程、速度方程、组分方程、温度方程、本构方程和理想气体状态方程。纳维-斯托克斯方程组反映了流动传热传质问题的基本规律。

（二）压力方程：压力方程由密度方程、速度方程和理想气体状态方程推导得来，是关于压力变量的输运方程，求解压力方程可得到最新的压力场。本发明在压力方程推导时考虑了密度的变化，因此能够实现可压缩流的数值求解。推导压力方程时首先将密度方程半离散，然后将半离散式中的密度和速度用压力代替，进而得到了压力方程。其中，密度和压力的关系可用理想气体状态方程来表达，速度和压力的关系可使用离散形式的速度方程来表达。

（三）化学扩散模型，用来求解化学反应速率和热释放速率。化学扩散模型是一种新型高效的简化机理，对化学反应动力学机理和扩散特性（热扩散、质量扩散和粘性效应）进行了简化，具有计算准确、效率高、节省计算资源等优势，能很好地复现火焰和燃爆的主要特性。所述化学扩散模型假设混合物共三种组分，即燃料、氧化剂和产物，分别对应三个组分守恒方程，燃料和氧化剂单步不可逆反应生成产物。化学反应速率由阿伦尼乌斯形式的公式进行表述。化学扩散模型包含指前因子、活化能、放热速率、扩散率、比热比和分子质量这六个参数，校准不同当量比下的化学扩散模型参数，可实现对非均匀反应流的计算。指前因子、活化能、放热速率、扩散率都采用查表法来根据当量比确定其数值。扩散率包括动量扩散、质量扩散和热扩散，是温度和密度的函数；路易斯数和普朗特数将三种扩散参数关联，是一组常量；此外，不同组分的比热比和分子质量也作为常数对待。所述查表法，即给定若干组已知的一一对应的当量比和目标值，组成一个查找表；然后随意给定一个当量比值，若其正好对应查找表中的一组数据，赋值即可；若没有对应的数据，则根据上下两组相近的数据组进行线性插值。

所述步骤二采用有限体积法，在物理网格划分的基础上，将步骤一建立的密度方程、速度方程、组分方程、温度方程以及压力方程在单个网格内体积积分，并在时间步长内进行时间积分，通过引入不同的空间离散格式和时间离散格式，最终得到了定义在该网格中心点上的代数方程；其中，体积积分时，根据高斯定理将对流通量、扩散通量的体积分转化为面积分，然后根据不同的空间离散格式将定义在网格面的变量重构到网格点；此外，化学反应速率项和热释放速率项的体积积分过程引入了变量在网格内均匀分布的假设；在时间积分过程中，化学反应速率和热释放速率进行显式处理；时间离散格式包括显式、半隐式、隐式。由于引入了时间离散格式，因此本发明可实现对非稳态问题的求解。完成了对单个网格内的体积积分和时间积分后，在所有网格内进行相同的运算，得到关于密度、速度、组分、温度和压力的代数方程组。将代数方程组线性化后进行迭代求解，得到相应的变量场。各个方程组之间的耦合求解将在步骤三中阐述。

如图1所示，所述步骤三将化学扩散模型耦合到压力修正方法中，用来耦合求解反应流中的各个变量场。其求解策略分为以下步骤：

步骤（1）确定全局时间步长，并进行时间推进。全局时间步长选取为化学时间步和对流时间步的较小值。采用全局隐式方法实现不同物理过程和化学过程之间的耦合，即非稳项、对流项、扩散项和源项将同时进行离散、线性化和求解。

步骤（2）求解离散化后的密度方程，得到密度预估值，用于预估校正迭代中的代数方程组的线性化。

步骤（3）求解离散后的速度方程，得到速度预测值。其中的压力值从上一时间步或初始条件中获取。

步骤（4）根据化学扩散模型计算每个网格中的组分反应速率和热释放速率。

步骤（5）求解离散后的组分方程和温度方程，得到组分场和温度场的分布。

步骤（6）求解离散后的压力方程更新压力场，并将更新的压力场带入速度方程用于更新速度场，完成一次压力修正循环。达到指定的压力修正循环次数后跳出循环，进入步骤（7）。

步骤（7）重复步骤（3）到步骤（6），直到收敛或者达到最大的预估校正迭代次数。

步骤（8）求解理想气体状态方程，更新密度值。

步骤（9）判断是否到达指定的结束时间，否则重复步骤（1）到步骤（8）。

在上述过程中，求解各代数方程组之前将代数方程组线性化，其中，代数方程中系数的值选取为当前可提供的最新值，以加快收敛速度；由于反应流问题中传热传质和流动过程相互耦合，因此化学反应速率和热释放速率的计算位于组分方程和温度方程之前，且都在预估校正迭代内部；在求解速度、组分和温度代数方程组之前进行扩散系数的更新；此外，所述步骤三的耦合机制在于：在压力修正方法的预估校正迭代中，采用最新的密度、温度和组分质量分数，根据化学扩散模型的算法，计算各组分的化学反应速率和热释放速率，进而求解组分方程和温度方程用以更新组分质量分数和温度，形成双向耦合。

实施例

以二维氢/空气对向流扩散火焰的数值计算为例，对本发明作出全面、详细的描述。该方法不局限于对向流场景，对所有扩散火焰以及预混火焰均适用。本发明实施例主要分为以下步骤：步骤（1）数学模型的描述；步骤（2）计算域、初始条件、边界条件的建立和描述；步骤（3）网格划分和计算流程描述；步骤（4）基于模型实现温度场、组分场等火焰特性的分析。

步骤（1）. 首先对数学模型进行描述，分为 3 个部分，具体包括：

第一部分，建立非定常、有粘、可压的多维度反应流数学模型，其中，根据质量守恒定律建立密度方程，根据牛顿运动定律在x方向和y方向上建立速度方程，根据能量守恒定律建立温度方程，此外还有关于燃料、氧化剂和生成物的三个组分方程，由此组成了制约反应流流动和传热传质的控制方程。控制方程的形式为：

，

其中，ρ为密度；φ为待求变量，包括密度、速度、组分、温度和压力；Γ为广义扩散系数，在速度方程、组分方程和温度方程中分别为运动粘度、质量扩散系数和热扩散系数；S为广义源项；u代表速度张量；t代表时间，

代表散度算子。

第二部分，根据化学扩散模型，本实施例选取氢、氧化剂和产物三种组分，且为单步不可逆反应，即氢和氧化剂单步不可逆反应生成产物。对不同当量比下氢/空气混合物的化学扩散模型参数进行校准，并写入查找表。热物理性质参数中，分子质量为25 g/mol，比热比为1.21，路易斯数和普朗特数分别为1和0.7。

第三部分，理想气体状态方程将压力和密度联系起来，在压力收敛后可用来更新对应的密度值。此外，本构方程将速度场和偏应力张量联系起来，用以封闭方程组。此外，由密度方程和速度方程推导出压力的输运方程，求解该方程可得到最新的压力场；压力方程推导时考虑了密度的变化，因此能够实现可压缩流的数值求解。

步骤（2） 计算域、初始条件、边界条件的建立和描述：

如图2所示，氢/空气对向流扩散火焰的装置包括：两个直径均为10 mm的氢气喷嘴、氧化剂喷嘴上下对立排布，氢气喷嘴、氧化剂喷嘴出口相距10 mm。氢气作为燃料和空气作为氧化剂以相同的速度分别从氢气喷嘴、氧化剂喷嘴流出。本实施例计算两个工况，工况一和工况二对应的喷嘴速度分别为2 m/s和4 m/s。

本实施例的计算域为喷嘴出口中间的二维矩形区域。初始时计算域温度设定为2000 K用于引燃氢/空气混合物。上边界属于燃料入口边界，此处的温度设定为300 K，燃料质量分数为1，压力为零梯度条件；下边界属于氧化剂入口边界，温度也设定为300 K，氧化剂质量分数为1，压力采用零梯度条件；左边界和右边界是出口边界，压力被指定为从总压力中减去动压，其余变量都采用零梯度条件。

步骤（3）二维氢/空气对向流扩散火焰的计算流程如下：

首先将计算域进行物理离散，划分为二维均匀正交网格。经过网格分辨率的测试，工况一和工况二下的网格尺寸最终选取为0.025 mm。然后采用有限体积法将密度方程、速度方程、组分方程、温度方程、压力方程在所有的网格上数值离散；其中时间离散格式采用一阶欧拉隐式，即对流项和扩散项都归属在未知时间层中，一阶欧拉隐式的基架点如图3所示，其中φ为待求变量，包括密度、速度、组分、温度和压力；下标中，C代表主控制体，N、E、S、W分别代表主控制体的相邻节点；上标中的t代表已知时间层，

代表未知时间层。最终得到定义在网格中心点上的关于密度、速度、组分、温度和压力的代数方程组，并使用冻结系数法对代数方程组进行线性化。

随后根据化学扩散模型和压力修正方法相互耦合的反应流求解策略对已经线性化的代数方程组进行顺序求解，具体步骤为：

步骤1，求解密度代数方程组，得到密度值；

步骤2，求解速度代数方程组，得到速度预估值；

步骤3，根据化学扩散模型计算组分反应速率和热释放速率，将计算结果作为源项分别带入组分和温度代数方程组中并求解，解得到组分和温度值；

步骤4，求解压力代数方程组，将压力的计算结果带入速度方程中并求解，得到速度校正值；

步骤5，循环步骤4，不断校正速度值，直到达到指定的循环次数；

步骤6，循环步骤2到步骤5，直到收敛或者达到最大的循环次数；

步骤7，求解理想气体状态方程，更新密度值；

步骤8，判断是否到达指定的结束时间，否则进行时间步推进，并重复步骤1到步骤7。

上述步骤所提到的对密度、速度、组分、温度和压力的代数方程组的求解采用共轭梯度法。根据以上步骤得到定义在每个网格中心点上的变量值，实现耦合单步化学机理的氢/空气对向流扩散火焰的数值求解。

步骤（4）对二维氢/空气对向流扩散火焰的温度场和组分场等火焰特性进行分析：

基于本发明计算得到的温度和组分廓线，如图4和图5所示，分别对应工况一和工况二的结果，横坐标表示到燃料喷嘴的距离。燃料和氧化剂通过喷嘴进入计算域，并对流效应和扩散效应下相互混合，高温引燃后在喷嘴之间产生扩散火焰。燃料和氧化剂在通过火焰面时被完全消耗，释放出化学能。两个工况下的温度廓线均呈锥形分布，峰值温度位置更靠近燃料喷嘴。对于两个工况，峰值温度分别为2373 K和2369 K；火焰厚度（定义为半峰值温度所占据的宽度）分别为3.8mm和2.5 mm。可以发现随着入口速度增大，峰值温度和火焰厚度都随之减小。因为当入口速度增大时，系统的应变率增加，达姆科勒数（定义为反应物停留时间与反应时间之比）减小，导致燃烧不完全，峰值温度降低。本发明计算出了氢/空气对向流扩散火焰的结构和主要种类的分布，并发现火焰厚度和峰值温度会随入口速度的增加而降低。

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (8)

1.一种耦合单步化学机理的反应流数值求解方法，其特征在于，包括如下步骤：

步骤一，建立非定常、有粘、可压的多维度反应流数学模型，化学反应速率、热释放速率、扩散系数使用化学扩散模型来确定；

步骤二，结合初始条件、边界条件和网格划分情况，对步骤一建立的多维度反应流数学模型进行数值离散，得到包括密度方程、速度方程、组分方程、温度方程和压力方程的代数方程组，并对代数方程组进行线性化；

步骤三，将化学扩散模型耦合到压力修正方法中，化学反应速率和热释放速率的计算位于组分方程和温度方程之前，且均在预估校正迭代内部，进行多次循环以保证收敛；在包括速度方程、组分方程和温度方程的代数方程组求解之前更新扩散系数；通过步骤三的上述步骤，建立耦合单步化学机理的反应流求解方法，包括：

步骤（1）确定全局时间步长，并进行时间推进；全局时间步长选取为化学时间步和对流时间步的较小值；采用全局隐式方法实现不同物理过程和化学过程之间的耦合，即非稳项、对流项、扩散项和源项同时进行离散化、线性化和求解；

步骤（2）求解离散化后的密度方程，得到密度预估值，用于预估校正迭代中的代数方程组的线性化；

步骤（3）求解离散后的速度方程，得到速度预测值；

步骤（4）根据化学扩散模型计算每个网格中的组分反应速率和热释放速率；

步骤（5）求解离散后的组分方程和温度方程，得到组分场和温度场的分布；

步骤（6）求解离散后的压力方程更新压力场，并将更新的压力场带入速度方程用于更新速度场，完成一次压力修正循环；达到指定的压力修正循环次数后跳出循环，进入步骤（7）；

步骤（7）重复步骤（3）到步骤（6），直到收敛或者达到最大的预估校正迭代次数；

步骤（8）求解理想气体状态方程，更新密度值；

步骤（9）判断是否到达指定的结束时间，否则重复步骤（1）到步骤（8）。

2. 根据权利要求 1 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，所述步骤一中，采用化学扩散模型来建立反应机理模型和扩散模型，用来求解化学反应速率和热释放速率。

3. 根据权利要求 1 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，所述步骤一中，所述化学扩散模型的参数的确定采用查表法。

4. 根据权利要求 1 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，所述步骤二采用有限体积法，在物理网格划分的基础上，将步骤一建立的密度方程、速度方程、组分方程、温度方程以及压力方程在网格内体积积分，并在时间步长内进行时间积分，通过引入不同的空间离散格式和时间离散格式，最终得到定义在该网格中心点上的代数方程。

5. 根据权利要求4 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，所述体积积分时，根据高斯定理将对流通量、扩散通量的体积分转化为面积分，然后根据不同的空间离散格式将定义在网格面的变量重构到网格点；在化学反应速率和热释放速率的体积积分过程引入变量在网格内均匀分布的假设；在所述时间积分的过程中，化学反应速率和热释放速率进行显式处理；所述时间离散格式包括显式、半隐式、隐式。

6. 根据权利要求 1 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，所述步骤二中，所述压力方程的推导考虑密度的变化，实现可压缩流的数值求解。

7. 根据权利要求 1 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，所述步骤三中，在压力修正方法的预估校正迭代中，采用更新的密度、温度和组分质量分数，根据化学扩散模型的算法，计算各组分的化学反应速率和热释放速率，进而求解组分方程和温度方程，用以更新组分质量分数和温度，形成双向耦合。

8. 根据权利要求 1 所述的一种耦合单步化学机理的反应流数值求解方法，其特征在于，用于处理非稳项、对流项、扩散项和源项。

Images