搅拌机颗粒导热系数的测定

魏慧荣，李晓东

(1.成都理工大学 工程技术学院资源勘查与土木工程系， 四川 乐山 614000；2.华北电力大学 能源动力与机械工程学院， 北京 102206)

温度梯度引起的混凝土开裂位置、 裂纹扩展方向以及扩展速率与热传导系数密切相关[1]，如何用准确高效的方法对热传导系数开展研究，分析其影响因素，对相关设备与材料进行优化，从而提高混凝土使用的安全性与寿命一直是众多学者研究的课题，但目前已有的研究方法存在较大的局限性，且不能完全满足并符合实际搅拌混凝土的工作要求[2]。

离散元法具有计算大量颗粒在给定条件下如何运动的能力，被广泛应用于建筑科学、道路建设等行业[3-5]。余志祥等采用离散元法对Tsuchiya阶梯型屋面的积雪开展了研究，得到了耦合状态下雪颗粒的碰撞模型。[6]。王鹏程应用离散元法得到了不同抄板角度、长度与混同设备中颗粒混合程度之间的关系[7]。王家全等采用EDEM软件分析了土工格栅在新旧路堤的受力分布,以及加筋区域土体细观参数的变化[8]。从检索情况得知，目前尚未采用离散元法对建筑科学、道路建设中颗粒的导热系数开展研究。本文以“JS2250/1500”型号双卧轴搅拌机为研究载体，对其搅拌过程进行研究分析，在此基础上，以导热系数的评价指标，运用离散元软件EDEM作二次正交旋转组合试验，旨在利用一定的计算方法获得颗粒理想导热系数的最优加工参数组合，为搅拌机的优化提供理论依据，促进建筑与筑路等行业的发展。

1 模型的建立

1.1 离散元算法理论模型

在离散元法的研究中，若颗粒间表面没有粘附作用，设置颗粒与几何体、颗粒与颗粒之间的接触模型为Hertz-Mindlin (no slip) built-in[9]。该模型颗粒间法向力为

(1)

式(1)中，R*为等效粒子半径，α为接触半径，E*为等效弹性模量，计算公式如下：

(2)

式(2)中，E1、ν1为颗粒1的弹性模量和泊松比，E2、ν2为颗粒2的弹性模量和泊松比。

颗粒间纵向力为

(3)

式(3)中，δ为切向重叠量，G*为等效剪切模量，计算公式如下：

(4)

式(4)中，G1和G2分别是颗粒1和颗粒2的剪切模量。

1.2 搅拌过程模型

将搅拌过程中不同颗粒含量的变化视为随机的马尔柯夫过程，使得搅拌过程可以采用柯尔莫高洛夫微分-差分方程表达[10]，即：

(5)

式中，Nf为时间为t1，搅拌机内某成分结合数量的随机值；

PNb(t)为随机取值Nf取整数值Nb的概率，即PNb=p{Nf(t)=Nb}；γNb为常数。

圆柱壳体搅拌机(图1)的搅拌过程另可以采用柯尔莫高洛夫方程模拟搅拌动态变化过程的扩散过程，即：

(6)

式中，M为搅拌物中关键成分的含量(kg)；t为搅拌时间(s)；v为搅拌物的轴向速度(m/s)；c1、c2分别为搅拌物的轴向、径向扩散系数，l、r分别为搅拌机内腔的轴向与径向坐标，轴向、径向扩散系数的大小主要取决于颗粒尺寸、充盈率以及搅拌速度等。但除上述搅拌过程的分析之外，混凝土等属粘-塑混合物，且搅拌轴上的叶片运动破坏了颗粒行为的连续性，因此混凝土的导热系数受其内部多因素影响。本研究通过改变充盈率、搅拌轴转速、叶片安装角等参数，可以间接影响导热系数。

图1 搅拌机模型简图

2 基于EDEM的辅助参数设计和优化

通过前述对搅拌工作过程的分析，且为优化搅拌机的几何参数和运动参数，改变混凝土的导热系数，选择取搅拌轴转速、充盈率、叶片安装角3个因素进行仿真试验。试验采用二次回归正交旋转组合设计。

2.1 EDEM软件仿真

搅拌机中水泥、砂子、碎石的几何尺寸等为随机变量，同样是影响导热系数大小的因素。本研究以矿渣硅酸盐水泥、粗砂、石子为研究对象，采用塌落试验(图2)测得不同颗粒塌落后的休止角，通过离散元公司所提供相关网站输入不同颗粒休止角获得其部分物理特性与力学特性[11]。对所购买某厂家所产的双卧轴搅拌机(图3)进行测量，并建立其三维模型，由于在离散元仿真模拟中只需将与颗粒接触的几何部件导入，所以将双卧轴搅拌机模型进行简化，应用三维软件UG对其建模，并导入EDEM软件尝试生成单一颗粒如图4所示。

图2 塌落试验装置图3 双卧轴搅拌机

图4 双卧轴搅拌机仿真模型

根据文献、测试结果及材料库，设置水泥、砂子、碎石颗粒以及双卧轴搅拌机的物理特性和相互间的力学特性如表1、表2所示[12-14]。

表1 仿真物理特性

表2 仿真学特性

由于砂子、碎石颗粒表面无粘附作用，加之前文对离散元模型的研究得知选择Hertz-MindLin(no slip)built-in为仿真接触模型。但水泥具有黏结的特性，且EDEM中的JKR模型适用于具有黏结作用的仿真，因此通过API接口对该接触模型进行二次开发[15]。命令水泥的接触模型为编译后的JKR模型，其余接触模型为Hertz-MindLin(no slip)built-in。本研究以混凝土强度等级为C25的配比为依据，同时为了保证仿真时双卧轴搅拌机工作有足够的颗粒和仿真时间，设置生成水泥颗粒数、砂子颗粒、碎石颗粒目分别为：90 000、32 000、17 000，半径平均值为1 mm、7 mm、23 mm，采用正态分布的方式生成颗粒尺寸，标准差分别为0.086、0.102、0.093，设置固定时间步长为Rayleigth时间步长17%，仿真时间总长22 s，前7 s为颗粒工厂生成各类颗粒的时间。同时为准确计算导热系数，对该双卧轴搅拌机进行网格划分如图5所示。

图5 网格划分

通过EDEM后处理器导出不同充盈率时，整个搅拌过程中颗粒的温度分布与运动轨迹如图6所示，其中不同颜色代表颗粒所具有的不同温度。

图6 不同充盈率颗粒的运动轨迹与温度

从图6可以看出搅拌机内不同充盈率会造成颗粒的运动轨迹与温度出现很大较大差别，为更加准确确定充盈率、搅拌轴转速、叶片安装角的变化对热导系数的影响趋势并确定了影响热导系数的主次因素，采用二次回归正交旋转组合试验对其开展研究。

2.2 仿真试验与指标

通过傅里叶定律，本研究的热导系数δy为[16]

(7)

根据前述理论分析及双卧轴搅拌机作业要求，合理的控制试验因素变化范围，每号试验重复10次，取其平均值作为试验结果，因素水平编码如表3所示，试验结果如表4所示。

表3 因素水平编码

表4 试验结果

运用SPSS数据分析软件对表3导热系数试验结果进行回归分析，其回归方程检验数据如表5所示。

表5 离散率回归方程检验

查F表得F0.01(10,13)=4.10，F=607>F0.01(10,13),因此回归方程高度显著，其二次回归方程模型为

A= 10.642-0.170X1-0.306X2-0.041X3+

0.002X1X2+0.001X2X3

运用Matlab绘制三维等值线图如图7所示。

分析图6可以得知，将充盈率固定在零水平时，随着搅拌轴转速的增大，导热系数呈现先下降后上升的趋势，将搅拌轴转速固定在零水平时，同样出现了随着充盈率的增加，导热系数先下降后上升的现象。另外，将叶片安装角固定在零水平时，随着搅拌轴转速的增加，导热系数出现了先下降后上升的趋势，而将搅拌轴转速固定在零水平时，随着叶片安装角的增加导热系数出现了缓慢下降的趋势。将叶片安装角固定在零水平时，随着充盈率的增加，导热系数出现了先缓慢下降后显著上升的现象，而将充盈率速固定在零水平时，叶片安装角的增加对导热系数影响不大。除此之外，由回归方程和三维等值线图可知，在试验范围内三因素对导热系数显著性影响大小顺序是：充盈率、搅拌轴转速、叶片安装角。

图7 导热系数等值线

根据优化数学模型和导热系数的回归方程，利用Matlab中非线性优化fmincon函数，以F(max)=δ即离散率最小为条件进行寻优处理。优化处理结果为：取搅拌轴转速为70.31 r/min，充盈率为29.42%，叶片安装角为54.41°。在此条件下进行仿真验证试验，试验表明导热系数为1.06 W/m·K。

3 验证

为验证离散元软件EDEM研究结果的可行性，开展试验验证。在充盈率为40%，叶片安装角为45°的情况下，使用变频器令搅拌轴转速分别为45 r/min、55 r/min、70 r/min、85 r/min、95 r/min，用导热系数测试仪XY-DRP-II(图8)对不同转速下颗粒的导热系数进行测定。试验与仿真结果如图9。

图8 XY-DRP-II型导热系数测试仪

图9 不同转速下试验与仿真的导热系数

从图8可以看出随着转速的增加，试验与仿真的导热系数变化趋势基本一致，试验所得导热系数略大。该现象主要原因有：试验中叶片倾角的加工存在误差，无法完全达到设定值，由于XY-DRP-II型导热系数测试仪所处室温变化，造成颗粒热流密度增大。当转速低于70 r/min，热导系数波动较小，试验所得导热系数与仿真所得导热系数差异较小，当转速高于70 r/min，导热系数波动较大，试验所得导热系数与仿真所得导热系数差异较大，最大差值为0.44 W/m·K。仿真分析结果与试验结果吻合，证明了采用离散元法研究搅拌机颗粒导热系数的可行性。

4 结论

运用EDEM作正交搅拌虚拟试验，通过二次回归正交旋转组合试验建立了以导热系数为考察指标的回归方程。运用Matlab绘制导热系数的等值线图，得到在试验范围内三因素对导热系数显著性影响大小顺序是：充盈率、搅拌轴转速、叶片安装角。

确定了该型双卧轴搅拌机的最佳参数组合，当拌轴转速为70.31 r/min，充盈率为29.42%，叶片安装角为54.41°。导数系数最小为1.06 W/m·K。

通过不同转速下颗粒导热系数的仿真与试验结果对比，得知两者结果变化趋势一致，结果基本相同。