CN111414677B - 一种多相高含量颗粒材料的传输性能计算方法 - Google Patents

Abstract

本发明公开了一种多相高含量颗粒材料的传输性能计算方法，首先生成高含量的颗粒复合材料三相细观结构模型，在结构体中随机投入一定数量的行人，通过变步长和定步长结合的方式使行人在颗粒复合材料三相细观结构体中进行随机行走，以模拟现实中的布朗运动过程，当行人行走的时间满足一定的条件之后，行人停止行走，并统计每个行人的位移值，根据行人的位移统计平均值以及行走时间计算得出颗粒复合材料的传输系数，本发明方法克服了传统数值方法对高含量颗粒材料传输性能预测的精度差、效率低等缺陷，操作灵活简便，具有良好的推广价值。

Description

一种多相高含量颗粒材料的传输性能计算方法

技术领域

本发明公开了一种多相高含量颗粒材料的传输性能计算方法，属于材料和土木工程技术领域。

背景技术

颗粒材料如陶瓷、泡沫混凝土、纳米材料等作为隔热、导热材料广泛地运用于各种工程领域。传输性能(包括热导率、电导率、扩散系数、渗透率、介电常数等)是颗粒材料最重要的物理性能之一，准确地预测颗粒材料的传输性能对实际生产具有重要意义。颗粒材料是一类典型的多相复合材料体系，由于其微观结构的复杂性，导致其性能的预测工作一直是材料科学领域的重点与难点之一。因此建立颗粒材料多相组分-结构-性能之间的定量关联机制，这对材料的设计和优化有非常重要的参考价值。最近的微观试验研究表明，在颗粒材料内部，复杂形貌的颗粒表面周围存在一层具有低刚度、高传导性的中间相层(界面过渡区)。这些中间相层相互联结组成了颗粒材料中重要的传输通路，对颗粒材料传输性能的影响不容忽视。为了考虑中间相层带来的影响，一种由颗粒、颗粒周围软质中间相层及基体三种组分构成的壳-核模型被广泛应用。

复合材料传输性能预测工作的主要内容是通过已知的参数如颗粒面积分数、中间相层的厚度以及各组分的传输系数等来推测复合材料整体的有效传输系数。主要方法分为理论方法与数值方法。理论方法主要基于复合材料的有效均质化理论，目前已有大量的均质化理论模型，如Hashin-Strikman界限模型、Mori-Tanaka模型、双/多夹杂模型、微分有效介质理论等。这类的理论方法通常适用于球形或椭球形的颗粒含量较低、组分性能差异不显著的情况，对于高含量的颗粒材料显得力不从心。数值方法为解决这类局限提供了另一种解决思路。传统的数值方法如有限差分法、有限元法及边界元法等方法通过将结构体离散化的方式对传输过程的控制微分方程进行近似求解。然而，这类方法的结果稳定性强烈依赖于网格划分的精准程度。对于多相高含量颗粒材料而言，由于其微观结构过于复杂，往往不能够进行准确的网格划分，因此传统的数值方法对于多相高含量颗粒材料的传输性能的预测能力有限。基于以上背景，建立一种可适用于多相高含量颗粒材料的传输性能预测方法是非常有必要的。

发明内容

本发明目的：针对上述的背景技术的不足，提供了能够精确高效地计算多相高含量颗粒材料传输性能的方法。

技术方案：一种多相高含量颗粒材料的传输性能计算方法，包括以下步骤：

步骤1：构建高含量颗粒材料三相细观结构模型，该结构模型由基体、软质中间相层及颗粒这三种材料组分复合而成，颗粒随机堆积在基体内，每个颗粒表面构造等厚度的软质中间相层；结构模型满足周期性边界条件；

步骤2：在步骤1构建的模型中随机投放一行人，采用变步长和定步长结合的方式使行人在模型中进行随机行走，得到行人行走的时间，当行人行走的时间达到预设值T时，行人停止行走，通过比较行人在行走前后的位置变化，统计每个行人的位移值；

步骤3：根据行人的位移平均值以及行走时间，计算出多相高含量颗粒材料的整体传输系数。

进一步的，步骤2中，所述的采用变步长和定步长结合的方式使行人在模型中进行随机行走，得到行人行走的时间的具体步骤包括：

S21：判断行人是否落入颗粒内部，若行人落入颗粒内部，则直接判定此次行走结束，该行人发生的位移X＝0；若行人位于颗粒外部，则进行S22；

S22：判断行人当前位置到最近界面边界的距离d_min与判定参数δ的大小关系，若d_min>δ，则进入步骤S23，否则进入步骤S24；所述界面边界为颗粒与软质中间相层之间的边界、软质中间相层与基体之间的边界；

S23：采用变步长的方式进行行走，以行人当前位置为圆心，以r_i＝d_min为半径构造虚拟圈，在该虚拟圈上随机选取一点作为行人的新位置，采用式(2)计算得到这一步使用的时间τ(r_i)：

其中，K_i为行人第i步位置所在组分的传输系数；

进入S25；

S24：采用定步长r_j进行行走，以行人当前位置为圆心，以r_j为半径构造虚拟圈，界面边界将虚拟圈分为两弧段，在与行人当前位置同侧或对侧的弧段上随机选点作为行人的新位置，采用式(3)计算得到这一步的行走时间τ_s(r_j)：

进入S25；

S25：计算行走总时长：

式中，i表示以变步长的方式行走的步数，j表示以定步长的方式行走的步数；

进入S26；

S26：判断行走总时长t是否达到预设值T，若达到，则行人停止行走，记录当前行人位移X，n＝n+1，其中，n为当前行人数，判断n是否小于N，N为向模型中随机投放行人的总数，若小于，则进入步骤2；否则进入步骤3；若行走总时长t未达到预设值T，则转入S22。

进一步的，所述S24中，引入行人穿越界面边界的概率P：

其中，K₁表示行人所在组分的传输系数，K₂表示界面边界另一侧的组分的传输系数，K₁和K₂均为已知系数；

取随机数s∈[0,1]，判断s是否小于P，若是，则在与行人当前位置对侧的弧段上随机选点作为行人的新位置，否则在与行人当前位置同侧的弧段上随机选点作为行人的新位置。

进一步的，S21中，以椭圆颗粒的形函数方程作为判别依据，判断行人是否落入颗粒内部：

式中，(x_cm,y_cm)为第m个椭圆颗粒中心坐标，a_m、b_m分别为第m个椭圆颗粒的长半轴长与短半轴长，θ_m为第m个椭圆的偏转角，代入行人当前位置坐标(x_w,y_w)，若f(x_w,y_w)≤0，则判定行人落入颗粒内部；若f(x_w,y_w)>0，则判定行人位于颗粒外部。

进一步的，步骤3中，采用式(6)计算得出多相高含量颗粒材料的传输系数K_e：

其中，X表示单个行人行走前后的位移；<>表示对所有行人的平均值。

进一步的，S22中所述的行人当前位置到最近界面边界的距离d_min的计算方法为：

将行人当前位置到颗粒与中间相之间的边界的最小距离记为d_p ^min，行人到中间相与基体之间的边界的最小距离记为d_i ^min；

若d_p ^min大于等于中间相层的厚度h，则d_i ^min＝d_p ^min-h，此时d_min＝d_i ^min；若d_p ^min<h，则d_i ^min＝h-d_p ^min，此时d_min＝min{d_i ^min，d_p ^min}。

进一步的，采用黄金分割搜索的方法计算行人位置到单个颗粒表面的距离，遍历所有颗粒后，即可得到行人当前位置到颗粒与中间相之间的边界的最小距离d_p ^min。

进一步的，所述黄金分割搜索的方法为：

计算求解行人当前位置到某一颗粒的上下两条切线，其斜率分别记为k₁和k₂，对上下切线斜率构成的区间[k₂,k₁]进行黄金分割，形成两条黄金分割线，其斜率分别为g₁＝k₂+0.618×(k₁-k₂)、g₂＝k₁-0.618×(k₁-k₂)，四条直线与颗粒的交点与行人当前位置形成四条线段，记斜率g₁、g₂对应线段的长度分别为h_g1、h_g2，若h_g1≥h_g2，则舍弃区间(g₁,k₁]，对[k₂,g₁]区间进行新一轮分割；若h_g1<h_g2，则舍弃区间[k₂,g₂)，对[g₂,k₁]区间进行新一轮分割；不断迭代得到行人当前位置到该颗粒的距离大小。

有益效果：本发明具有以下优点：

1、本发明通过模拟布朗运动过程实现了对高含量多相颗粒复合材料传输系数包括热导率、电导率、磁导率、扩散系数、渗透率、介电常数等精确的数值预测。相比较于传统的数值模拟方法如有限差分法、有限元法和边界元法等需要对结构体的离散网格划分来求解控制方程而影响数值模拟结果的精度，本发明提出的方法无需对结构模型进行网格划分，其计算精度可控，操作更加简便，极具推广价值。

2、本发明对布朗运动模拟过程进行了优化，采取了变步长+定步长结合的方式进行随机行走，相较于以往单一地采用定步长行走的方式来说，计算效率及计算精度均得到了极大的提升。对于以往的定步长行走方式，行走步长严重影响计算精度和效率：若步长设置过小则计算效率低下；若步长设置过大则行人难以通过狭隘通路，影响计算精度。本发明提出的方法克服了这一缺点，根据结构体的实际情况自我调整行走步长，自适应能力极强。

附图说明

图1为高密度颗粒堆积模型构造过程；

图2为高含量三相颗粒复合材料模型；

图3为本发明方法的流程图；

图4为本发明方法中单个行人行走过程示意图；

图5为本发明提出的计算方法与理论结果及有限元结果的对比图，包括(a)、(b)、(c)三幅子图。

具体实施方式

现结合附图和实施例进一步阐述本发明的技术方案。

本实施例的多相高含量颗粒材料传输性能包括热导率、电导率、磁导率、扩散系数、渗透率、介电常数等传输系数，提出了一种精准高效计算传输性能数值的方法，包括以下步骤：

步骤1：生成高含量的颗粒复合材料三相细观结构模型，该高含量颗粒复合材料三相细观结构模型包括基体、软质中间相层及颗粒的三相复合材料，颗粒随机紧密地堆积在基体内，每个颗粒表面有等厚度的软质中间相层，将颗粒体的几何形状考虑为二维椭圆形，考虑颗粒为绝缘材料，中间相及基质为可导材料。

如图1所示，高含量的颗粒复合材料三相细观结构模型的构造方法为：首先在一个较大的方形容器中随机地投入椭圆形颗粒，保证颗粒之间互不重叠，以形成颗粒的松散随机堆积结构；再以一定的速率使容器边界向内收缩移动，椭圆形颗粒在边界的作用下向内运动，直至椭圆颗粒面积分数达到给定值；根据几何拓扑知识，在每个椭圆颗粒表面构造等厚度的中间相层。在颗粒堆积结构构造过程中满足周期性边界条件，即当颗粒与容器一条边相交时，在其对边处生成一个相应的复制颗粒；当颗粒与容器的角点相接触时，在另外三个角点处生成相应复制颗粒。最终形成的三相高含量颗粒复合材料结构模型如图2所示。

步骤2：基于步骤1生成的模型，使用优化的随机行走算法进行有效传导性测定，其步骤分为：

分步1：在模型中随机地投放数量为N(取N＝50000)的行人；以椭圆颗粒的形函数方程作为判别依据，判别行人是否落入颗粒内部：

其中，(x_cm,y_cm)为第m个椭圆颗粒中心坐标，a_m、b_m分别为第m个椭圆颗粒的长半轴长与短半轴长，θ_m为第m个椭圆的偏转角，代入行人位置坐标(x_w,y_w)，若f(x_w,y_w)≤0，则判定行人落入颗粒当中，由于考虑颗粒为绝缘相，落入颗粒中的行人无法进行行走，可直接判定此次行走结束，该行人发生的位移为0；若f(x_w,y_w)>0，则判定行人位于颗粒外部，允许进行下一步行走。

分步2：采用变步长和定步长结合的方式，使行人在颗粒复合材料三相细观模型体中进行随机行走，以模拟现实中的布朗运动过程，行走过程示意图如图4所示；

本分步中的变步长和定步长相结合的行走方式包括：

(a)计算行人到最近界面边界的距离d_min，具体计算方法叙述如下：

在本实施例涉及的问题中存在两种界面边界：颗粒与中间相之间的边界、中间相与基体之间的边界。将行人到前者边界的最小距离记为d_p ^min，行人到后者边界的距离记为d_i ^min，通过黄金分割搜索的方法搜寻行人当前位置到某个椭圆颗粒表面的距离，具体方法为：首先计算求解点到椭圆的上下两条切线，其斜率分别记为k₁和k₂，对上下切线斜率构成的区间进行黄金分割，形成两条黄金分割线，其斜率分别为g₁＝k₂+0.618×(k₁-k₂)、g₂＝k₁-0.618×(k₁-k₂)，记斜率g₁、g₂对应线段的长度分别为h_g1、h_g2，若h_g1≥h_g2，则舍弃区间(g₁,k₁]，对[k₂,g₁]区间进行新一轮分割；若h_g1<h_g2，则舍弃区间[k₂,g₂)，对[g₂,k₁]区间进行新一轮分割。如此不断迭代最终即可求得点到椭圆的距离大小。遍历所有椭圆颗粒并借助以上方法来搜寻到行人所在位置最近的距离，即可求得d_p ^min。

若d_p ^min大于等于中间相层的厚度h，即d_p ^min≥h，则表明行人在基质相内，d_i ^min＝d_p ^min-h，此时d_min＝d_i ^min；若d_p ^min<h则表明行人在中间相层内，则d_i ^min＝h-d_p ^min，此时d_min为d_p ^min与d_i ^min的较小值，即d_min＝min{d_i ^min，d_p ^min}。

(b)预先设置判定参数δ(如设置δ＝0.01r_eq，r_eq为颗粒的等效半径)作为判定依据，判断d_min与判定参数δ的大小关系，若d_min>δ，则进入步骤(c)，否则进入步骤(d)；

(c)当d_min>δ，采用变步长的方式进行行走，令行走步长大小r_i＝d_min，其中i表示以变步长的方式行走的步数，以行人的位置为中心以r_i为半径构建一个虚拟圆，在圆上随机选取一点作为行人新的位置，这一步使用的时间τ(r_i)可通过下式进行计算：

其中，K_i为行人第i步位置所在组分的传输系数；

返回步骤(b)进行下一步判断；

(d)当d_min≤δ，采用定步长r_j进行行走，其中j表示以定步长的方式行走的步数，以行人位置为中心r_j为半径构建虚拟圆，此时需要计算行人穿越界面边界的概率P：

其中，K₁表示行人所在组分的传输系数，K₂表示界面边界另一侧的组分的传输系数；

取随机数s∈[0,1]，判断s是否小于P，若是，则在与行人当前位置对侧的弧段上随机选点作为行人的新位置，否则在与行人当前位置同侧的弧段上随机选点作为行人的新位置。

这一步的行走时间τ_s(r_j)可表示为：

由于颗粒相传输系数为0，因此在行走过程中行人无法进入颗粒相中；

返回步骤(b)进行下一步判断；

分步3：当行人行走总时间达到预设值T(取T＝400)时则判定行人完成行走，行人停止行走，此时需要比较行人在行走前后的位置变化，统计每个行人的位移值；

行人的总行走时间t为定步长+变步长两种情况所用时间的叠加：

分步4：根据行人的位移统计平均值以及行走时间，根据如下公式计算得出颗粒复合材料的传输系数K_e：

其中，X表示单个行人行走前后的位移；<>表示对所有行人的平均值。

为验证本方法的准确性，使用本方法分别与有限元方法结果与多夹杂模型理论结果相对比。图5中的(a)为在颗粒长径比κ＝2，中间相层厚度h与颗粒的等效半径r_eq之比h/r_eq＝0.1以及中间相层与基体之间的传输系数比K_i/K_m＝2情况下得到的有效传输系数比K_e/K_m与颗粒面积分数的关系；图5中的(b)为在颗粒面积分数f_p＝0.5，h/r_eq＝0.1以及K_i/K_m＝2情况下得到的有效传输系数比K_e/K_m与颗粒长径比κ的关系；图5中的(c)为在颗粒面积分数f_p＝0.6，颗粒长径比κ＝2以及K_i/K_m＝2情况下得到的有效传输系数比K_e/K_m与h/r_eq的关系。图5的结果说明本方法的计算结果是可靠的。表1对比了在f_p＝0.4，K_m＝1，κ＝2以及颗粒等效半径r_eq＝3.568情况下，多夹杂理论模型、有限元法、传统随机游走算法及本专利方法这四种方法的结果。

表1

方法	K<sub>e</sub>/K<sub>m</sub>	相对误差	计算时间(小时)
				理论结果	0.4	基准	-
有限元法	0.381	4.6％	-
				传统随机游走算法	0.384	3.80％	11.04
本专利方法	0.386	3.50％	1.96

由此可知本发明提出的方法在保证精度的前提下，极大地提升了计算效率，具有良好的推广价值。

Claims (7)

1.一种多相高含量颗粒材料的传输性能计算方法，其特征在于：包括以下步骤：

步骤1：构建高含量颗粒材料三相细观结构模型，该结构模型由基体、软质中间相层及颗粒这三种材料组分复合而成，颗粒随机堆积在基体内，每个颗粒表面构造等厚度的软质中间相层；结构模型满足周期性边界条件；

步骤2：在步骤1构建的模型中随机投放一行人，采用变步长和定步长结合的方式使行人在模型中进行随机行走，得到行人行走的时间，当行人行走的时间达到预设值T时，行人停止行走，通过比较行人在行走前后的位置变化，统计每个行人的位移值；

步骤3：根据行人的位移平均值以及行走时间，计算出多相高含量颗粒材料的整体传输系数；

步骤2中，所述的采用变步长和定步长结合的方式使行人在模型中进行随机行走，得到行人行走的时间的具体步骤包括：

S21：判断行人是否落入颗粒内部，若行人落入颗粒内部，则直接判定此次行走结束，该行人发生的位移X＝0；若行人位于颗粒外部，则进行S22；

S22：判断行人当前位置到最近界面边界的距离d_min与判定参数δ的大小关系，若d_min>δ，则进入步骤S23，否则进入步骤S24；所述界面边界为颗粒与软质中间相层之间的边界、软质中间相层与基体之间的边界；

S23：采用变步长的方式进行行走，以行人当前位置为圆心，以r_i＝d_min为半径构造虚拟圈，在该虚拟圈上随机选取一点作为行人的新位置，采用式(2)计算得到这一步使用的时间τ(r_i)：

其中，K_i为行人第i步位置所在组分的传输系数；

进入S25；

S24：采用定步长r_j进行行走，以行人当前位置为圆心，以r_j为半径构造虚拟圈，界面边界将虚拟圈分为两弧段，在与行人当前位置同侧或对侧的弧段上随机选点作为行人的新位置，采用式(3)计算得到这一步的行走时间τ_s(r_j)：

其中，K₁表示行人所在组分的传输系数，K₂表示界面边界另一侧的组分的传输系数，K₁和K₂均为已知系数；

进入S25；

S25：计算行走总时长：

式中，i表示以变步长的方式行走的步数，j表示以定步长的方式行走的步数；

进入S26；

S26：判断行走总时长t是否达到预设值T，若达到，则行人停止行走，记录当前行人位移X，n＝n+1，其中，n为当前行人数，判断n是否小于N，N为向模型中随机投放行人的总数，若小于，则进入步骤2；否则进入步骤3；若行走总时长t未达到预设值T，则转入S22。

2.根据权利要求1所述的一种多相高含量颗粒材料的传输性能计算方法，其特征在于：所述S24中，引入行人穿越界面边界的概率P：

取随机数s∈[0,1]，判断s是否小于P，若是，则在与行人当前位置对侧的弧段上随机选点作为行人的新位置，否则在与行人当前位置同侧的弧段上随机选点作为行人的新位置。

3.根据权利要求1所述的一种多相高含量颗粒材料的传输性能计算方法，其特征在于：S21中，以椭圆颗粒的形函数方程作为判别依据，判断行人是否落入颗粒内部：

式中，(x_cm,y_cm)为第m个椭圆颗粒中心坐标，a_m、b_m分别为第m个椭圆颗粒的长半轴长与短半轴长，θ_m为第m个椭圆的偏转角，代入行人当前位置坐标(x_w,y_w)，若f(x_w,y_w)≤0，则判定行人落入颗粒内部；若f(x_w,y_w)>0，则判定行人位于颗粒外部。

4.根据权利要求1所述的一种多相高含量颗粒材料的传输性能计算方法，其特征在于：步骤3中，采用式(6)计算得出多相高含量颗粒材料的传输系数K_e：

其中，X表示单个行人行走前后的位移；<>表示对所有行人的平均值，T为预设值。

5.根据权利要求1所述的一种多相高含量颗粒材料的传输性能计算方法，其特征在于：S22中所述的行人当前位置到最近界面边界的距离d_min的计算方法为：

将行人当前位置到颗粒与中间相之间的边界的最小距离记为d_p ^min，行人到中间相与基体之间的边界的最小距离记为d_i ^min；

若d_p ^min大于等于中间相层的厚度h，则d_i ^min＝d_p ^min-h，此时d_min＝d_i ^min；若d_p ^min<h，则d_i ^min＝h-d_p ^min，此时d_min＝min{d_i ^min，d_p ^min}。

6.根据权利要求5所述的一种多相高含量颗粒材料的传输性能计算方法，其特征在于：采用黄金分割搜索的方法计算行人位置到单个颗粒表面的距离，遍历所有颗粒后，即可得到行人当前位置到颗粒与中间相之间的边界的最小距离d_p ^min。

7.根据权利要求6所述的一种多相高含量颗粒材料的传输性能计算方法，其特征在于：所述黄金分割搜索的方法为：

计算求解行人当前位置到某一颗粒的上下两条切线，其斜率分别记为k₁和k₂，对上下切线斜率构成的区间[k₂,k₁]进行黄金分割，形成两条黄金分割线，其斜率分别为g₁＝k₂+0.618×(k₁-k₂)、g₂＝k₁-0.618×(k₁-k₂)，四条直线与颗粒的交点与行人当前位置形成四条线段，记斜率g₁、g₂对应线段的长度分别为h_g1、h_g2，若h_g1≥h_g2，则舍弃区间(g₁,k₁]，对[k₂,g₁]区间进行新一轮分割；若h_g1<h_g2，则舍弃区间[k₂,g₂)，对[g₂,k₁]区间进行新一轮分割；不断迭代得到行人当前位置到该颗粒的距离大小。

Images