一种基于球簇假设和Laguerre-Voronoi结构的生物质大颗粒热解仿真方法

一种基于球簇假设和laguerre-voronoi结构的生物质大颗粒热解仿真方法
技术领域
1.本发明属于生物质热解领域，具体涉及一种基于球簇假设和laguerre-voronoi结构的生物质大颗粒热解仿真方法，可用于探究单颗粒和多颗粒系统中生物质热解行为和现象。


背景技术：

2.由于原料的可再生性、产物的多样性等原因，生物质热解在新能源转换和利用中具有重要的应用前景。目前在该领域的研究，除了实验测量，利用数值模拟方法对其进行探究不仅能够提供热解实验中无法得到的细节，也可以进一步指导实验的设计和运行，减少研发和应用成本。
3.实际过程中，生物质热解多发生在流化床、固定床等反应器中。该尺度下的多颗粒模拟工作往往采取了过多简化假设，如假设颗粒为球形、忽略颗粒内部热梯度、假设颗粒等体积或等密度转化等。这些假设并不适用于各向异性的生物质颗粒，因此需要开发更加精细和有效的颗粒尺度热解模型，并将其用于多颗粒系统中，以提高反应器尺度模拟的精度和性能。
4.已有的考虑了颗粒内部流体输运的颗粒尺度模型，计算成本过大。而忽略颗粒内部流场的单颗粒热解模型大多只适用于球、方块等规则形状颗粒，且无法解析各向异性的颗粒内部热传输和收缩过程。


技术实现要素：

5.本发明针对以上技术难点，提出了一种基于球簇假设和laguerre-voronoi结构的生物质大颗粒热解仿真方法。所述方法属于生物质热解领域，可用于探究单颗粒和多颗粒系统中生物质热解的行为和现象，其主要包括初始球簇生成、时间步离散及热解迭代、结果输出三个步骤。本发明利用球簇假设实现对任意不规则形状颗粒在任意热边界条件下的热解表征；利用组成子球的体积减小和移动共同实现母颗粒的收缩；利用laguerre-voronoi结构对颗粒内部区域进行划分，不仅实现了颗粒内部(子球间)热传输过程，而且解析了颗粒内部孔隙结构，同时也为子球移动计算提供支持；利用方向相关的热导率和收缩因子考虑了颗粒各向异性的热传导和收缩过程。本方法在使用球簇假设和laguerre-voronoi结构的基础上，提供了一个快速有效的生物质大颗粒热解仿真方法，有助于获取生物质热解过程中实验难以观察得到的内部细节，也可以扩展用于含多生物质颗粒的反应器尺度综合模型中，进一步提高大尺度模拟的预测精度和准确度。
6.具体地，本发明提供了一种基于球簇假设和laguerre-voronoi结构的生物质大颗粒热解仿真方法，其包括如下步骤：
7.s1，初始球簇生成
8.采用具有大小、密度、温度、物质组成属性的子球体模拟真实生物质材料，由若干
子球体共同构成球簇，以球簇等效表示原始生物质大颗粒；
9.s2，在时间尺度上将热解过程划分为若干时间区间，并在每一时间区间上进行颗粒热解迭代；颗粒热解迭代过程包括颗粒结构更新、传热计算、热解反应计算和颗粒收缩计算四个子步骤；
10.s3，完成迭代并输出计算结果。
11.作为本发明的优选方案，所述的初始球簇生成的方法为：利用网格生成技术对原始真实生物质进行网格划分；根据原始真实生物质的孔隙率以及网格位置和体积，依次对应生成特定位置和半径的子球体，子球体具有与原始真实生物质一致的本征密度；子球体间的空间对应真实生物质材料内部的孔隙；子球体的温度、物质组成属性与原始真实生物质保持一致。
12.作为本发明的优选方案，所述的颗粒结构更新包括外包络面更新和内laguerre-voronoi结构更新；
13.其中，外包络面更新具体方法为：采用若干与球簇表面子球体相切的平面将球簇整体包络在一凸多面体内，所形成的包络面集合构成所表征生物质颗粒的边界；
14.内laguerre-voronoi结构更新具体方法为：利用laguerre-voronoi结构对任意子球体所属的空间区域进行划分：中心位置为xi，半径为ri的子球i所属的空间区域定义为laguerre-voronoi单元，其中任意一点x满足
[0015][0016]
对任意子球j(j≠i)恒成立，其中xj，rj分别为子球j的中心位置和半径，d(x,xi)和d(x,xj)表示点x与子球i和j中心的距离。
[0017]
作为本发明的优选方案，所述的传热计算包括外部对流辐射计算和内部导热计算；
[0018]
外部对流辐射计算仅影响外表面子球体的温度变化，其热通量根据包络面面积以及当地环境对流和辐射相关参数计算得到；
[0019]
内部导热计算影响任意子球体的温度变化，其热通量通过导热系数、相邻子球的温度差和相对位置以及laguerre-voronoi单元交界面面积计算得到；若将导热系数与相邻子球的相对位置关联，则可以模拟生物质颗粒内各向异性的导热过程。
[0020]
作为本发明的优选方案，所述的热解反应计算采用热解机理和动力学模型依次对子球体进行热解计算，根据反应后的产物得到子球质量和密度变化。
[0021]
进一步地，所述的热解机理和动力学模型可选择arrhenius模型。
[0022]
作为本发明的优选方案，所述的颗粒收缩计算包括子球体积更新和子球平移计算；
[0023]
其中，子球体积更新由热解过程中的质量变化和密度变化实现，其公式为
[0024][0025]
其中vi，mi，ρi分别表示子球i的体积、质量和密度；
[0026]
子球平移计算包括主动平移子球计算和被动平移子球计算。
[0027]
所述的主动平移子球计算方法为：首先定义收缩阈值γ
threshold
，其值在[0,1]范围
内选取；随后对子球体进行穷举，满足
[0028]
ωi＞γ
threshold
·
max(ω)
[0029]
关系式的子球均发生主动移动，相反则不具备主动移速，其中ωi表示子球i的单位体积固体消耗速率，max(ω)则表示当前时刻最大的单位体积固体消耗速率；
[0030]
子球主动移速基于在时间步内对应laguerre-voronoi单元孔隙率不变的假设，其计算公式为：
[0031][0032]
其中ui为子球i的主动移速，s
celli
和v
celli
则表示子球i所属的laguerre-voronoi单元的表面积和体积，(xc,yc,zc)和(xi,yi,zi)分别是生物质颗粒和其第i个组成子球的中心坐标，χ
x
，χy和χz分别是x，y，z方向上的收缩因子，用于实现各向异性的生物质颗粒收缩。
[0033]
所述的被动平移子球计算方法为：针对任意主动移动的子球i，首先在边界上搜索得到子球j，使得其与子球i的连线向量与子球i的主动移速方向最接近；随后在子球i与j之间搜索得到最短路径；最后所有该路径上除子球i外的子球体均施加被动平移，相应的被动移速为2ui。
[0034]
进一步地，所述的输出计算结果可以包括颗粒内部温度分布及演化、密度及孔隙率分布及演化、颗粒形态演化、颗粒质量演化等生物质大颗粒热解细节，其中，颗粒内部的相关信息和演化过程是实验方法难以观察和获得的。
[0035]
本发明与现有技术相比，具有如下有益效果：
[0036]
(1)采用的球簇假设能够实现对任意不规则形状生物质颗粒的表征，以及解析其在热解过程中的形态变化；
[0037]
(2)通过laguerre-voronoi结构对颗粒内部区域进行划分，不仅能够实现颗粒内部(子球间)热传输过程，而且能够解析颗粒内部孔隙结构，并进一步为颗粒收缩计算提供支持；
[0038]
(3)利用方向相关的热导率和收缩因子可以求解颗粒各向异性的热传导和收缩过程；
[0039]
(4)不仅能得到单生物质颗粒热解过程中实验难以观察到的颗粒内部细节信息，也可以与多颗粒反应器尺度计算流体力学模型集成，进而显著提高大尺度模拟的预测精度和准确度。
附图说明
[0040]
图1为本发明方法的流程图。
[0041]
图2为典型的采用球簇假设的颗粒表示示意图。
[0042]
图3为某球颗粒结构。
[0043]
图4为传热计算示意图。
[0044]
图5为子球被动收缩计算示意图。
[0045]
图6为某球形生物质颗粒在某工况下颗粒中心温度及剩余质量分数演化图。
[0046]
图7为某立方体生物质颗粒在底部辐射加热状态下热解过程形态演化。
具体实施方式
[0047]
为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方案作进一步详细说明。
[0048]
如图1所示，为本发明基于球簇假设和laguerre-voronoi结构的生物质大颗粒热解过程内部细节参数获取方法的流程图。所述方法主要包括以下步骤：
[0049]
步骤一，初始球簇生成。如图2所示，采用若干具有大小、密度、温度、物质组成等属性的子球体共同构成初始球簇，以等效表示原始生物质大颗粒；子球体具有当地真实的生物质材料的本征密度；子球间的空间抽象为材料内部的孔隙结构；图2示出了球形和不规则形两种原始生物质大颗粒的子球体表示方法。
[0050]
具体地，初始球簇生成的可选方法为：首先利用网格生成技术对原始生物质进行网格划分，如六面体结构化网格；其次根据原始生物质的孔隙率以及网格位置和体积，依次对应生成特定位置和半径的子球体，子球的温度、物质组成等其他属性与原始生物质保持一致。
[0051]
步骤二，在时间尺度上对热解过程划分为若干区间，一般可采用将时间区间等分的形式进行区间的划分，并在每一区间上进行颗粒热解迭代，包括颗粒结构更新、传热计算、热解反应计算和颗粒收缩计算四个子步骤；
[0052]
颗粒结构更新子步骤包括外包络面更新和内laguerre-voronoi结构更新。外包络面更新方法为采用若干与表面子球体相切的平面将球簇整体包络在一凸多面体内，所形成的包络面集合构成所表征生物质颗粒的边界。内laguerre-voronoi结构更新方法为利用laguerre-voronoi结构对任意子球体所属的空间区域进行划分；中心位置为xi，半径为ri的子球i所属的空间区域(下称laguerre-voronoi单元)中任意一点x满足
[0053][0054]
对任意子球j(j≠i)恒成立，其中xj，rj分别为子球j的中心位置和半径，d(x,xi)和d(x,xj)表示点x与子球i和j中心的距离。某典型的球形颗粒结构如图3所示。
[0055]
传热计算子步骤目的在于更新子球体的温度，包括外部对流辐射计算和内部导热计算，如图4所示。外部对流辐射计算仅对外表面子球体有效，其热通量根据包络面面积以及当地环境对流和辐射相关参数计算得到。内部导热热通量通过导热系数、相邻子球的温度差和相对位置以及laguerre-voronoi单元交界面面积计算得到；导热系数与相邻子球的相对位置有关。
[0056]
热解反应子步骤计算采用特定的热解机理和动力学模型(如arrhenius模型)依次对子球体进行热解计算，根据反应后的产物得到子球质量和密度变化。
[0057]
颗粒收缩计算子步骤包括子球体积更新和子球平移计算。子球体积更新由下式计算得到
[0058][0059]
其中vi，mi，ρi分别表示子球i的体积、质量和密度。子球平移计算包括主动平移子球计算和被动平移子球计算。
[0060]
主动平移子球判定方法为：首先定义收缩阈值γ
threshold
，其值在[0,1]范围内选
取；随后对子球体进行穷举，满足
[0061]
ωi＞γ
threshold
·
max(ω)
[0062]
关系式的子球均发生主动移动，相反则不具备主动移速，其中ωi表示子球i的单位体积固体消耗速率，max(ω)则表示当前时刻最大的单位体积固体消耗速率；子球主动移速计算公式为：
[0063][0064]
其中ui为子球i的主动移速，s
celli
和v
celli
则表示子球i所属的laguerre-voronoi单元的表面积和体积，(xc,yc,zc)和(xi,yi,zi)分别是生物质颗粒和其第i个组成子球的中心坐标，χ
x
，χy和χz分别是x，y，z方向上的收缩因子，用于实现各向异性的生物质颗粒收缩；
[0065]
如图5所示，被动平移子球判定和求解方法为：针对任意主动移动的子球i，首先在边界上搜索得到子球j，使得其与子球i的连线向量与子球i的主动移速方向最接近；随后在子球i与j之间搜索得到最短路径；最后所有该路径上(除子球i外)的子球体均施加被动平移，相应的被动移速为2ui。
[0066]
步骤三，完成迭代并输出计算结果。结果可以包括颗粒内部温度分布及演化、密度及孔隙率分布及演化、颗粒形态演化、颗粒质量演化等生物质大颗粒热解细节等。
[0067]
本发明的方法可在任意可编程平台中进行植入、运行和仿真，仿真时，需自行指定颗粒表面的对流和辐射状态。例如，采用c 语言实现该仿真方法，并对某球形和某立方体颗粒热解过程求解得到图6和图7，其中，图6为某球形颗粒在某工况下颗粒中心温度及剩余质量分数演化图，图7为某立方体颗粒在底部辐射加热状态下的热解过程形态演化；从结果可见，球形颗粒中心温度和质量损失与实验对比良好，立方体颗粒的形态变化也与实际过程较为符合。本发明能得到单生物质颗粒热解过程中实验难以观察到的颗粒内部细节信息，也可以与多颗粒反应器尺度计算流体力学模型集成，进而显著提高大尺度模拟的预测精度和准确度。