留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

基于离散元的竹粉颗粒接触参数标定

刘东 刘子昕 王琦 李晓旭 闫承琳

张海涛, 宫渊波, 付万权, 等. 马尾松低效人工林不同改造模式下降雨及产流特征[J]. 浙江农林大学学报, 2018, 35(1): 29-34. DOI: 10.11833/j.issn.2095-0756.2018.01.004
引用本文: 刘东, 刘子昕, 王琦, 等. 基于离散元的竹粉颗粒接触参数标定[J]. 浙江农林大学学报, 2023, 40(4): 875-882. DOI: 10.11833/j.issn.2095-0756.20220662
ZHANG Haitao, GONG Yuanbo, FU Wanquan, et al. Rainfall and runoff in low efficiency Pinus massoniana forests with different silvicultural prescriptions[J]. Journal of Zhejiang A&F University, 2018, 35(1): 29-34. DOI: 10.11833/j.issn.2095-0756.2018.01.004
Citation: LIU Dong, LIU Zixin, WANG Qi, et al. Calibration of contact parameters for bamboo powder particles based on DEM[J]. Journal of Zhejiang A&F University, 2023, 40(4): 875-882. DOI: 10.11833/j.issn.2095-0756.20220662

基于离散元的竹粉颗粒接触参数标定

DOI: 10.11833/j.issn.2095-0756.20220662
基金项目: 中央级公益性科研院所基本科研业务费专项资金(CAFYBB2020SY040)
详细信息
    作者简介: 刘东(ORCID: 0009-0003-9644-1395),从事木质材料增材制造研究。E-mail: 1104496173@qq.com
    通信作者: 闫承琳(ORCID: 0009-0009-0014-8180),副研究员,博士,从事木质材料增材制造技术与人造板机械等研究。E-mail: yanchenglin@126.com
  • 中图分类号: S781.9

Calibration of contact parameters for bamboo powder particles based on DEM

  • 摘要:   目的  建立竹粉颗粒材料的离散元模型并标定其相关接触参数。  方法  基于毛竹Phyllostachys edulis粉末物理堆积角实验结果,采用离散元法(DEM)模拟仿真和实验设计(DOE)相结合的方法,通过Plackett-Burman (P-BD)实验、爬坡实验和响应面实验获得竹粉堆积角和仿真参数之间的二次多项式回归模型,以物理堆积角为目标值预测得到接触参数的最优组合。  结果  实测竹粉物理堆积角为49.29°;P-BD实验结果表明了竹粉-竹粉滚动摩擦系数、竹粉-竹粉恢复系数和竹粉-不锈钢板静摩擦系数对堆积角影响显著(P<0.05);响应面实验结果进一步证实了P-BD实验的结果,同时表明竹粉-竹粉滚动摩擦系数自身交互作用、竹粉-竹粉滚动摩擦系数和竹粉-竹粉恢复系数的交互作用对堆积角有显著影响(P<0.05);最优因素组合为竹粉-竹粉恢复系数0.44,竹粉-竹粉滚动摩擦系数0.22,竹粉-不锈钢板静摩擦系数0.45。  结论  本研究在41.20°~55.27°堆积角范围内建立的竹粉颗粒堆积角定量模拟标定与预测平台具有较高的可靠性,DEM标定方法可为后续建立评估生物质颗粒运动行为及其与设备相互作用的材料模型提供参考。图5表6参27
  • 森林是陆地生态系统的主体,起着调节大气、涵养水源、为人类提供生产资料等作用。截至2009年,中国人工林栽植面积已达5 300万hm2,约占全世界人工林面积的40%[1],但受人为因素或诱导自然因素所致[2],中国人工林普遍存在地力衰退、生物多样性差、水土流失严重等问题,对林地生态安全造成了一定的隐患[3],因此,恢复和重建退化人工林生态系统势在必行。健康的森林生态系统功能,能够提高林地土壤肥力以及水土保持能力[4]。植被的存在能有效减少地表径流量,同时不同林地类型的减滞能力也不尽相同[5-6]。因此,降雨与径流及植被、土壤等因子之间的关系是目前研究的热点及难点问题[7-8]。近年来国内针对南方红壤区[9]和黄土高原[10]等区域坡地不同经营、利用方式下水土流失研究较多,针对川南地区的研究则鲜见报道。马尾松Pinus massoniana是广泛分布于中国南方的先锋树种,它具有耐贫瘠、速生、适应性强、经济价值高等特点,然而受不合理经营方式及人为活动的影响,长江流域低山丘陵区马尾松人工林普遍生长较差,生物多样性低,生态功能不强,已成为中国南方森林面积最大的退化类型之一。笔者针对长江上游低山丘陵区存在的生态安全以及生态工程建设中遇到的科学技术等问题,在前期研究工作的基础上,以川南低山丘陵区马尾松低效人工林为示范区,引入珍贵乡土树种,改造马尾松低效人工林,提升森林生态及经济功能,减少水土流失。本研究根据设立于四川省宜宾市高县来复镇的5个人工径流小区,分析马尾松低效人工林改造初期自然降雨与产流特征的关系及随植被恢复的演变规律,以期为川南马尾松低效人工林改造及经营管理提供参考。

    研究区位于四川省宜宾市高县来复镇毛颠坳(28°11′~28°47′N, 104°21′~104°48′E)。该地区地处四川盆地与云贵高原的过渡地带,宜宾市中南部,属乌蒙山余脉。地貌以平坝、丘陵、低山为主。土壤多为山地黄壤,间断分布有少量紫色土。全年平均日照时数为1 107.7 h,大于等于10 ℃以上年积温为6 523.1 ℃,年平均气温为18.0 ℃,1月平均气温为8.0 ℃,7月平均气温为27.0 ℃。年平均降水量为1 037.9 mm,年平均无霜期为346.2 d。

    2012年初在试验区选择坡度约22°,东西走向并栽植有马尾松纯林的坡地,建立5个相邻规格为20 m × 5 m的标准径流小区,各径流小区出口处均设有一个径流收集室。5个径流小区采用不同的改造更新方式,分别用Ⅰ,Ⅱ,Ⅲ,Ⅳ,Ⅴ标识。Ⅰ号径流小区为马尾松低效人工林下植被自然生长恢复,Ⅱ号径流小区为马尾松低效人工林下空隙更新樟树Cinnamomum camphora,Ⅲ和Ⅳ号径流小区为马尾松低效人工林皆伐后当年更新樟树,Ⅴ号径流小区为马尾松低效人工林皆伐后第2年更新樟树(由于有时会出现因为工期或者苗木原因不能当年更新造林的情况,因此,将Ⅴ号径流小区设计为皆伐后第2年更新樟树)。小区内樟树为当年生幼苗,株高约30 cm,按1.5 m × 1.5 m间距种植。2012年3月对Ⅲ,Ⅳ,Ⅴ号径流小区进行皆伐,随即在Ⅱ,Ⅲ,Ⅳ号内更新樟树;Ⅴ号径流小区则自然恢复后于2013年3月更新樟树。试验开始前5个径流小区内地表灌草及凋落物均清理完毕。试验分为3个时段,分别为2012年7-12月、2013年1-12月、2014年1-12月。在每次自然降雨后收集降水、计算24 h降雨量,并量取径流桶内径流水深,换算成径流量。雨量级别根据气象学上对于降雨量的定义判定,即24 h降雨量小于10 mm为小雨,10.1~24.9 mm为中雨,25.0~49.9 mm为大雨,50.0~99.9 mm为暴雨。

    径流系数是地表径流量与降雨量的比值,表示有比例的降水变成了径流,它能够一定程度反映该区域植被和土壤的水源涵养能力以及水土流失状况。可利用公式r=Rn/P×10-3An计算径流系数。其中:r为径流系数,Rn为径流量(m3),P为降雨量(mm),An为径流小区面积(m2[11]

    利用Excel 2007和SPSS 19.0软件进行数据处理及统计学分析。

    试验期年平均降雨总量约1 000 mm(表 1),与当地多年平均降雨量数据基本吻合。试验区降雨明显呈现夏季多、冬季少的特点。2013年和2014年全年多集中于每年6-9月,分别占全年降雨总量的73.23%和69.50%;1月、2月、11月、12月降雨极少,分别占全年降雨总量的3.72%和6.19%。

    表  1  马尾松低效人工林改造初期降雨量及分布特征
    Table  1.  Characteristics of rainfall capacity and distribution under the early of transformation on low eificiency Pinus massoniana forests
    年份 各月份降雨量 合计/mm
    1 2 3 4 5 6 7 8 9 10 11 12
    2012 250.3 114.0 277.2 49.4 14.1 7.6 712.6
    2013 3.2 8.1 23.0 85.9 92.1 130.8 279.4 231.7 137.8 44.4 12.7 15.6 1 064.7
    2014 5.3 14.5 88.2 17.7 45.9 137.9 104.2 236.3 168.9 74.6 20.8 17.0 931.3
    下载: 导出CSV 
    | 显示表格

    2012-2014年24 h最小降雨量分别为1.8,1.5和1.2 mm(表 2),最大降雨量分别为83.6,88.7和93.8 mm。2012年7-12月期间共出现5次暴雨天气,24 h最大降雨量为83.6 mm,累计降雨341.6 mm;大雨仅出现2次。2013年累计降雨73次,其中小雨47次,中雨16次,大雨6次,暴雨4次。2014年累计降雨78次,其中小雨47次,与2013年相同,降雨量减少28.8 mm;暴雨次数比2013年少2次,降雨量减少133.4 mm。纵观整个试验期,降雨次数以小雨和中雨居多,占85.56%,但雨量仅占49.05%;大雨和暴雨次数占14.44%,降雨量占比却达到50.95%。表明在1年中大到暴雨的次数比例虽然很低,但降雨量大且集中,是导致夏季地表产流较多的主要原因。

    表  2  马尾松低效人工林改造初期雨量特征
    Table  2.  Characteristics of rainfall under the early of transformation on low efficiency Pinus massoniana forests
    年份 总降雨量/mm 总降雨次数/次 小雨 中雨 大雨 暴雨 取大降雨量/mm 最小降雨量/mm
    雨量/mm 占总降雨量比例/% 次数/次 占总降雨次数比例/% 雨量/mm 占总降雨量比例/% 次数/次 占总降雨次数比例/% 雨量/mm 占总降雨量比例/% 次数/次 占总降雨次数比例/% 雨量/mm 占总降雨量比例/% 次数/次 占总降雨次数比例/%
    2012下半年 712.6 36 97.4 13.67 17 47.22 199.7 28.02 12 33.33 73.9 10.37 2 5.56 341.6 47.94 5 13.89 83.6 1.8
    2013年 1 064.7 73 238.2 22.37 47 64.38 275.0 25.83 16 21.92 240.7 22.61 6 8.22 310.8 29.19 4 5.48 88.7 1.5
    2014年 931.3 78 209.4 22.48 47 60.26 309.0 33.18 21 26.92 262.8 28.22 8 10.26 150.1 16.12 2 2.26 93.8 1.2
    合计 2 708.6 187 545.0 20.12 111 59.36 783.7 28.93 49 26.20 577.4 21.32 16 8.56 802.5 29.63 11 5.88 266.1 4.5
    下载: 导出CSV 
    | 显示表格

    表 3可以看出:并非每次降雨都能引起地表径流。一般情况下,在达到最小产流降雨量后,才能产生地表径流。2012年下半年降雨36次,产流20次,产流次数占降雨次数比例为55.60%;2013年下半年与2014年下半年降雨次数相同,均为42次,产流次数占降雨次数比例分别为33.33%和45.24%,可见下半年产流降雨比例在不同年份间有较大变化。从全年看,2013年和2014年分别降雨73次和78次,产流25次和29次,产流次数占降雨次数比例为34.25%和37.18%,比较接近。对比表 2表 3可以看出,2013年和2014年小雨都是47次,仅有1次产流,可见24 h降雨量小于10 mm时一般不会产生地表径流。而24 h降雨量达到中雨时,产流几率达到了85.00%以上,下半年甚至可能高于90.00%;在大雨及暴雨状态下,每次均有产流。

    表  3  马尾松低效人工林改造初期降雨及产流次数
    Table  3.  Times of rainfall and runoif under the early of transformation on low efficiency Pinus massoniana forests
    年份 降雨次数/次 产流次数/次 产流降雨比例/% 降雨量/mm 最小产流降雨量 小雨及产流次数/次 中雨及产流次数/次 大到暴雨及产流次数/次
    2012下半年 26 20 55.56 712.6 6.9 2(17) 11(12) 7(7)
    2012下半年 42 14 22.22 721.6 14.0 0(26) 7(8) 7(7)
    2014下半年 42 19 45.24 621.8 8.2 1(22) 11(12) 7(7)
    2012年 72 25 24.25 1 064.7 12.8 0(47) 15(16) 10(10)
    2014年 78 29 27.18 921.2 8.2 1(47) 18(21) 10(10)
    说明:表中括号内数字表示降雨次数,括号外表示产流次数。
    下载: 导出CSV 
    | 显示表格

    同时可以看出,观测时期内,24 h最小产流降雨量为6.9~14.0 mm。造成这种差异的原因可能为该次小雨前有降雨发生,地表凋落物及土壤持水能力趋近于饱和,小雨过后就能形成蓄满产流;若前期降雨较少,凋落物和表层土壤含水率低,在中雨条件下雨水填洼、入渗等较多,降雨过后土壤仍未达到饱和,因此,无地表径流产生。当24 h降雨量达到大雨或暴雨条件时,5个径流小区皆有地表径流产生。

    表 4看出:Ⅰ~Ⅴ号小区最大径流量均出现在下半年,2013年最大径流量比2012年分别增加18.21%,35.62%,15.42%,13.66%和33.13%。鉴于2013年和2012年下半年降雨量仅相差9 mm,说明改造初期地表径流量差异主要由植被恢复状况及地被物覆盖变化等因素导致。其原因在于进行不同模式低效人工林改造后,植被生长、更替和人为干扰直接影响林冠层和地表层,影响了雨水拦截、入渗能力,更多的降雨转变为地表径流。说明在低效人工林改造开始阶段,林地水土保持功能体现并不明显,水土流失现象反而有可能加剧,因此,该时期应是水土流失重点监测、保护期。该结果与廖承彬等[12]研究结论相似。5个径流小区2014年最大量与2013年同期相比分别降低32.25%,26.49%,22.26%,25.08%和25.22%;其中Ⅲ,Ⅳ,Ⅴ号小区最大径流量低于2012年同期水平,呈波状变化。造成这一现象的原因主要有:一是2014年7-12月降雨量比2013年同期减少近100 mm,更少降雨导致更少径流;二是改造后随着人为扰动减少,乔?鄄灌?鄄草体系逐渐恢复,凋落物增加,凋落物和土壤持水能力进一步增强。各径流小区间最小产流量差异不大,当降雨量为大雨到暴雨时,前期短时间内降水经过林冠截留、凋落物截留、地表填洼等阶段后,易形成超渗产流,这时降雨强度成为影响产流量大小的主导因子。

    表  4  马尾松低效人工林不同模式改造初期产流量特征
    Table  4.  Characteristics of runoff under different transformation patterns of low efficiency Pinus massoniana forests
    径流小区 产流量/m3
    2012年下半年最小 2012年下半年最大 2012年下半年最小 2012年下半年最大 2013年最小 2013年最大 2014年下半年最小 2014年下半年最大 2014年最小 2014年最大
    0.010 2.202 0.098 2.603 0.018 2.603 0.022 2.185 0.008 2.185
    0.008 2.125 0.102 2.882 0.004 2.882 0.003 2.200 0.003 2.200
    0.012 2.601 0.086 3.002 0.006 3.002 0.006 2.006 0.006 2.006
    0.005 2.547 0.075 2.895 0.012 2.895 0.005 2.102 0.005 2.102
    0.010 2.092 0.070 2.785 0.012 2.785 0.018 1.956 0.018 1.956
    下载: 导出CSV 
    | 显示表格

    表 5显示了2.5 a内5个径流小区累计产流量状况。由于Ⅲ号和Ⅳ号径流小区改造模式相同,其累计产流量十分接近,分别为36.339 m3和36.369 m3;而Ⅰ号由于郁闭度较高,林冠截留作用好于其他径流小区,且受人为扰动最小,因此,累计产流量最低,为31.315 m3,试验期内变化幅度最小。而平均单次产流量大小排序为Ⅲ>Ⅳ>Ⅴ>Ⅱ>Ⅰ,排序规律与径流小区所受扰动强度一致。

    表  5  马尾松低效人工林不同模式改造初期累计产流量特征
    Table  5.  Characteristics of cumulative runoff under different transformation patterns of low efficiency Pinus massoniana forests
    径流小区 不同年份产流量/m3 累计产流次数/次 累计产流量/m3 平均单次产流量/m3
    2012年下半年 2013年 2014年
    10.077 13.166 8.072 69 31.315 0.453 8
    10.165 14.355 8.562 72 33.082 0.459 5
    13.341 15.068 7.930 72 36.339 0.504 7
    13.045 14.768 8.556 74 36.369 0.491 5
    11.662 13.675 7.299 69 32.636 0.473 0
    下载: 导出CSV 
    | 显示表格

    对2012-2014年下半年各径流小区径流系数进行比较。由表 6可看出:除Ⅱ号径流小区呈现先增大后减小的变化以外,其他径流小区径流系数均呈逐年减小趋势。与2012年下半年相比,各径流小区2014年下半年径流系数分别减少29.36%,26.30%,47.37%,41.54%和43.97%。2012年,5个径流小区径流系数大小排序为Ⅲ>Ⅳ>Ⅴ>Ⅱ>Ⅰ,2014年则变化为Ⅳ>Ⅱ>Ⅰ>Ⅲ>Ⅴ。2014年与2013年相比,5个径流小区径流系数都有明显下降,均在25%以上,特别是Ⅲ和Ⅴ小区,径流系数下降更为明显,达36%,超出Ⅰ号小区10%。随着时间推移,效果会更加明显。从径流系数来看,Ⅰ号径流小区虽然在改造初期郁闭度最高,林冠截留能力最强,但由于林下缺乏凋落物,导致对于林内降雨的拦滞能力较弱;其他径流小区由于拥有更好的光照、水热条件,植被生长较快,生物多样性更为丰富,因此,随着改造进行Ⅰ号小区的径流系数可能会逐渐高于其他径流小区。而Ⅴ号径流小区由于在2012年未更新樟树,有1 a自然更新时间,土体扰动少于Ⅲ和Ⅳ号径流小区,因此,其径流系数在3个皆伐小区最小。

    表  6  马尾松低效人工林不同模式改造初期径流系数
    Table  6.  Runoff coefficient under different transformation patterns of low eificiency Pinus massoniana forests
    时间
    2012年下半年(r1) 0.139 3 0.140 7 0.1841 0.180 3 0.161 7
    2011年下半年 0.136 9 0.148 4 0.152 4 0.150 6 0.142 2
    2014年下半年(r2) 0.098 4 0.103 7 0.096 9 0.105 4 0.090 6
    2013年(r3) 0.117 1 0.127 7 0.1347 0.131 9 0.122 9
    2014年(r4) 0.086 7 0.091 9 0.085 2 0.091 9 0.078 4
    r1/r2/% 29.36 26.30 47.37 41.54 43.97
    r4/r3/% 25.99 28.03 36.75 30.33 36.21
    下载: 导出CSV 
    | 显示表格

    利用SPSS 19.0软件进行降雨量与各径流小区径流量相关性分析。由表 7可以看出:降雨量与径流量之间相关系数均大于0.905,径流量大小随降雨量变化而变化。与刘芝芹等[13]对30场降雨观测中得出的研究结论相吻合,说明降雨量是影响地表径流大小的重要因子。

    表  7  马尾松低效人工林改造初期降雨量与各径流小区径流量相关性分析
    Table  7.  Correlation coefficients of rainfall and runoff under different transformation patterns of low efficiency Pinus massoniana forests
    年份
    2012 0.971** 0.971** 0.988** 0.988** 0.986**
    2013 0.916** 0.909** 0.905** 0.913** 0.906**
    2014 0.947** 0.957** 0.959** 0.960** 0.951**
    说明:**表示在0.01水平(双侧)上极显著相关。
    下载: 导出CSV 
    | 显示表格

    试验区年均降水量约为1 000 mm,集中于6-9月,季节性分布不均。2013年与2014年6-9月降雨量分别占全年降雨总量的73.23%和69.50%;1,2,11,12月降雨量则仅占3.72%和6.19%。

    试验期内24 h最小降雨量为1.2 mm,最小产流降雨量为6.9 mm,最大降雨量为93.8 mm。全年产流次数约占降雨次数的30%~40%,雨量以中到暴雨为主,偶有小雨下产流和中雨下无产流状况发生。各径流小区最大径流量均呈现先增大后减小的情况,说明在低效人工林改造初期有水土流失加剧现象,是水土流失监测、治理的关键期。

    径流系数反映了降雨转变为地表径流的比例。改造后第1年下半年Ⅰ号和Ⅱ号径流小区由于较好的植被覆盖以及相对较少的人为干扰,径流系数低于Ⅲ,Ⅳ,Ⅴ号径流小区;但2013年和2014年Ⅰ号和Ⅱ号径流小区径流系数降低均超过30%,而Ⅲ,Ⅳ,Ⅴ号径流小区降低50%左右,基本与Ⅰ号和Ⅱ号径流小区持平或更低。表明各径流场水土保持功能均在恢复中。2012年与2013年Ⅰ号径流小区径流系数均为最小,2014年则逐渐接近甚至高于皆伐小区;Ⅲ,Ⅳ与Ⅴ号径流小区在试验期内产流变化过程相似。目前来看,皆伐径流小区水源涵养功能恢复速度快于Ⅰ号和Ⅱ号径流小区,减滞径流效果更佳。

    对降雨量和各径流小区径流量进行相关性分析发现,其相关系数均大于0.905。因此,减少林内降雨以及增强林地水源涵养能力是减少坡面径流的关键。由于试验尚处于改造初期,时间较短,相信随着改造的深入以及林地生态系统的恢复,樟树等阔叶林的水源涵养能力会显著优于马尾松低效纯林。

  • 图  1  堆积角测定装置

    Figure  1  Stacking angle measuring device

    图  2  竹粉颗粒模型(A)和堆积角模型(B)

    Figure  2  Bamboo powder particle model (A) and powder stacking angle model (B)

    图  3  堆积角图像处理

    Figure  3  Image processing of stacking angle

    图  4  有效堆积图像(A)和无效堆积图像(B)

    Figure  4  Effective (A) and invalid (B) stacking angles image

    图  5  交互项X1X3响应面

    Figure  5  Response surface of interaction term X1X3

    表  1  竹粉堆积角测定仿真参数

    Table  1.   Simulation parameters for measuring the stacking angle of bamboo powder

    参数符号参数名称取值范围参数水平
    −11
    μ 竹粉泊松比 0.20~0.50 0.20 0.40
    E 竹粉剪切模量/GPa 1~15 2 4
    X1 竹粉-竹粉恢复系数 0.15~0.30 0.25 0.50
    X2 竹粉-竹粉静摩擦系数 0.20~0.80 0.30 0.60
    X3 竹粉-竹粉滚动摩擦系数 0~0.40 0.15 0.30
    X4 竹粉-不锈钢板恢复系数 0.10~0.80 0.25 0.50
    X5 竹粉-不锈钢板静摩擦系数 0.30~0.70 0.30 0.60
    X6 竹粉-不锈钢板滚动摩擦系数 0~0.40 0.15 0.30
    J 表面能 /(J·m−2) 0.01~0.10 0.04 0.08
    下载: 导出CSV

    表  2  竹粉堆积角P-BD实验设计及结果

    Table  2.   Plackett-Burman test design and results of bamboo powder stacking angle

    实验号μEX1X2X3X4X5X6J模拟堆积角/(°)
    1−1−1111−111−151.54
    200000000047.17
    3−1−1−1111−11154.48
    4−1−1−1−1−1−1−1−1−142.63
    51−11−1−1−111141.68
    6−111−11−1−1−1146.47
    7111−111−11−137.22
    8−11−1−1−1111−140.56
    911−111−11−1−157.91
    101−111−11−1−1−123.71
    11−1111−111−1139.00
    1211−11−1−1−11137.79
    131−1−1−1111−1161.32
      说明:符号所表示参数名称见表1
    下载: 导出CSV

    表  3  竹粉堆积角爬坡实验设计及结果

    Table  3.   Climbing test design and results of bamboo powder stacking angle

    实验号一次爬坡二次爬坡
    X1X3X5模拟堆积角/(°)相对误差/%X1X3X5模拟堆积角/(°)相对误差/%
    10.500.150.3038.9221.040.500.150.3039.6719.52
    20.400.250.4048.49 4.500.450.200.3541.5315.75
    30.300.350.5057.9517.570.400.250.4049.325.20
    40.200.450.600.350.300.4552.827.17
    50.100.550.700.300.350.5057.2716.19
      说明:−表示数据无效。符号所表示参数名称见表1
    下载: 导出CSV

    表  4  竹粉堆积角响应面实验设计及结果

    Table  4.   Box-Behnken test design and results of bamboo powder stacking angle

    序号X1X3X5模拟堆积角/(°)
    10(0.40)1(0.30)−1(0.35)48.69
    200(0.25)0(0.40)48.89
    31(0.45)01(0.45)52.17
    410−143.18
    5−1(0.35)0−148.71
    600048.51
    7−11054.98
    80−1(0.20)−141.53
    9−10154.97
    100−1148.85
    111−1041.20
    1200050.54
    1311051.69
    1401155.27
    15−1−1048.11
      说明:括号内数值为参数值。参数X1X3X5所表示名称见表1
    下载: 导出CSV

    表  5  竹粉堆积角P-BD实验参数显著性分析

    Table  5.   Significance parameters analysis of bamboo powder stacking angle Plackett-Burman test

    参数效应系数P显著性排序
    μ−2.506−1.2530.2307
    E−2.735−1.3680.2036
    X1−9.179−4.5890.0252
    X2−0.904−0.4520.6019
    X313.9276.9630.0111
    X4−3.623−1.8110.1325
    X58.2864.1430.0303
    X6−1.297−0.6490.4708
    J4.5272.2640.0914
      说明:参数所表示名称见表1
    下载: 导出CSV

    表  6  Box-Behnken实验回归模型方差分析

    Table  6.   ANOVA of modified model of Box-Behnken

    方差来源均方自由度平方和FP
    模型 265.87 9 29.54 86.60 <0.000 1
    X1 38.22 1 38.22 112.04 <0.000 1
    X3 108.81 1 108.81 318.97 <0.000 1
    X5 103.35 1 103.35 302.96 <0.000 1
    X1X3 5.44 1 5.44 15.93 0.010 4
    X1X5 1.86 1 1.86 5.45 0.066 9
    X3X5 0.03 1 0.03 0.09 0.776 9
    X1X1 0.49 1 0.49 1.44 0.283 2
    X3X3 7.63 1 7.63 22.36 0.005 2
    X5X5 0.64 1 0.64 1.88 0.228 4
    残差 1.71 5 0.34
    纯误差 0.015 2 0.007
    合计 267.58 14
    R2=0.9936,Radj2=0.9822,CV=1.18,AP=30.541
      说明:参数X1X3X5所表示名称见表1
    下载: 导出CSV
  • [1] ZHANG Wenbo, HU Tao, CHANG Yanting, et al. Correlation between genetic characteristics, cell structure and material properties of moso bamboo (Phyllostachys edulis (Carriere) J. Houzeau) in different areas of China [J/OL]. Forests, 2022, 13(1): 107[2022-09-30]. doi:10.3390/F13010107.
    [2] 闫承琳, 刘东, 刘子昕, 等. 基于木塑基耗材的增材制造技术研究进展[J]. 林业工程学报, 2022, 7(4): 22 − 30.

    YAN Chenglin, LIU Dong, LIU Zixin, et al. Research progress of additive manufacturing technology and equipment based on wood-plastic consumables [J]. Journal of Forestry Engineering, 2022, 7(4): 22 − 30.
    [3] 陈铭, 郭琳, 郑笑, 等. 中国15个主产区毛竹纤维形态比较[J]. 南京林业大学学报(自然科学版), 2018, 42(6): 7 − 12.

    CHEN Ming, GUO Lin, ZHENG Xiao, et al. Comparison of cell morphology of moso bamboo fibers from fifteen main producing regions in China [J]. Journal of Nanjing Forestry University (Natural Sciences Edition), 2018, 42(6): 7 − 12.
    [4] RACKL M, TOP F, MOIHOEK C P, et al. Feeding system for wood chips: a DEM study to improve equipment performance [J]. Biomass and Bioenergy, 2017, 98: 43 − 52.
    [5] FRACZEK J, ZIOBECKI A, ZEMANEK J. Assessment of angle of repose of granular plant material using computer image analysis [J]. Journal of Food Engineering, 2007, 83(1): 17 − 22.
    [6] RAMIREZ-GOMEZ A, GALLEGO E, FUENTES J M, et al. Values for particle-scale properties of biomass briquettes made from agroforestry residues [J]. Particuology, 2014, 12: 100 − 106.
    [7] GRIMA A P, WYPYCH P W. Discrete element simulations of granular pile formation [J]. Engineering Computations, 2011, 28(3): 314 − 339.
    [8] KRUGGEL-EMDEN H, RICKELT S, WIRTZ S, et al. A study on the validity of the multi-sphere discrete element method [J]. Powder Technology, 2008, 188(2): 153 − 165.
    [9] IBARRA J, ESTAY D, PACHECO A, et al. Bond calibration method for macroparameters using the discrete element method framework [J/OL]. Engineering Fracture Mechanics, 2022, 262: 108223[2022-09-30]. doi:10.1016/j.engfracmech.2021.108223.
    [10] ADAJAR J B, AIFARO M, CHEN Y, et al. Calibration of discrete element parameters of crop residues and their interfaces with soil [J/OL]. Computers and Electronics in Agriculture, 2021, 188: 106349[2022-09-30]. doi: 10.1016/j.compag.2021.106349.
    [11] ZHAO Liang, ZHOU Hongping, XU Linyun, et al. Parameter calibration of coconut bran substrate simulation model based on discrete element and response surface methodology [J]. Powder Technology, 2022, 395: 183 − 194.
    [12] XIA Yidong, CHEN Feiyang, KLINGER J L, et al. Assessment of a tomography-informed polyhedral discrete element modelling approach for complex-shaped granular woody biomass in stress consolidation [J]. Biosystems Engineering, 2021, 205: 187 − 211.
    [13] PACHON-MORALES J, DO H, COLIN J, et al. DEM modelling for flow of cohesive lignocellulosic biomass powders: model calibration using bulk tests [J]. Advanced Powder Technology, 2019, 30(4): 732 − 750.
    [14] TAN Yuan, YU Yue, FOTTNER J, et al. Automated measurement of the numerical angle of repose (aMAoR) of biomass particles in EDEM with a novel algorithm [J]. Powder Technology, 2021, 388: 462 − 473.
    [15] 轻工业部日用化学工业科学研究所. 表面活性剂 粉体和颗粒休止角的测量: GB/T 11986—1989[S]. 北京: 中国标准出版社, 1989.

    Institute of Daily Chemical Industry Science, Ministry of Light Industry. Surface Active Agents- Powders and Granules-Measurement of the Angle of Repose: GB/T 11986−1989[S]. Beijing: Standards Press of China, 1989.
    [16] ALIZADEH M, ASACHI M, GHADIRI M, et al. A methodology for calibration of DEM input parameters in simulation of segregation of powder mixtures, a special focus on adhesion [J]. Powder Technology, 2018, 339: 789 − 800.
    [17] FENG Y T, HAN K, QWEN D R J, et al. Engineering computations on upscaling of discrete element models: similarity principles [J]. Engineering Computations, 2009, 26(6): 599 − 609.
    [18] OREFICE L, KHINAST J G. A novel framework for a rational, fully-automatised calibration routine for DEM models of cohesive powders [J]. Powder Technology, 2020, 361: 687 − 703.
    [19] COETZEE C. Calibration of the discrete element method: Strategies for spherical and non-spherical particles [J]. Powder Technology, 2020, 364: 851 − 878.
    [20] HOSHISHIMA C, OHSAKI S, NAKAMURA H, et al. Parameter calibration of discrete element method modelling for cohesive and non-spherical particles of powder [J]. Powder Technology, 2021, 386: 199 − 208.
    [21] 杨云芳, 刘志坤. 毛竹材抗拉弹性模量及抗拉强度[J]. 浙江林学院学报, 1996, 13(1): 21 − 27.

    YANG Yunfang, LIU Zhikun. Phyllostachys pubescens wood: tensile elastic modulus and tensile strength [J]. Journal of Zhejiang Forestry College, 1996, 13(1): 21 − 27.
    [22] 李荣荣, 贺楚君, 彭博, 等. 毛竹材不同部位纤维形态及部分物理性能差异[J]. 浙江农林大学学报, 2021, 38(4): 854 − 860.

    LI Rongrong, HE Chujun, PENG Bo, et al. Differences in fiber morphology and partial physical properties in different parts of Phyllostachys edulis [J]. Journal of Zhejiang A&F University, 2021, 38(4): 854 − 860.
    [23] 刘文政, 何进, 李洪文, 等. 基于离散元的微型马铃薯仿真参数标定[J]. 农业机械学报, 2018, 49(5): 125 − 135, 142.

    LIU Wenzheng, HE Jin, LI Hongwen, et al. Calibration of simulation parameters for potato minituber based on EDEM [J]. Transactions of the Chinese Society for Agricultural Machinery, 2018, 49(5): 125 − 135, 142.
    [24] GALLEGO E, FUENTES J, RUIZ Á, et al. Determination of mechanical properties for wood pellets used in DEM simulations [J]. International Agrophysics, 2020, 34(4): 485 − 494.
    [25] ZU E X, ZHOU P, JIANG Z H. Discrete element method of coke accumulation: calibration of the contact parameter [J]. IFAC-Papers on Line, 2018, 51(21): 241 − 245.
    [26] ZHU Jianzhong, ZOU Meng, LIU Yansong, et al. Measurement and calibration of DEM parameters of lunar soil simulant [J]. Acta Astronautica, 2022, 191: 169 − 177.
    [27] 周捍东, 徐长妍, 丁沪闽. 木材散碎物料基本堆积特性的研究[J]. 林业和草原机械, 2002(6): 9 − 12, 19.

    ZHOU Handong, XU Changyan, DING Humin. Studies on the basic characteristics of wooden bulk materials [J]. Forestry and Grassland Machinery, 2002(6): 9 − 12, 19.
  • [1] 叶金枝, 徐涓, 杨静, 杨海艳, 王大伟, 史正军.  没食子酸改性竹粉的制备及抑菌性能 . 浙江农林大学学报, 2023, 40(6): 1333-1340. doi: 10.11833/j.issn.2095-0756.20220718
    [2] 叶淑媛, 董雷鸣, 董昂, 喻卫武, 戴文圣, 曾燕如.  榧树半同胞子代幼苗生长性状的遗传参数估算 . 浙江农林大学学报, 2020, 37(4): 817-822. doi: 10.11833/j.issn.2095-0756.20190542
    [3] 金迪, 张明如, 王佳佳, 高磊, 侯平.  遮光与水分胁迫对盆栽芒萁光合与叶绿素荧光参数的影响 . 浙江农林大学学报, 2020, 37(6): 1054-1063. doi: 10.11833/j.issn.2095-0756.20190666
    [4] 周哲宇, 徐超, 胡策, 王海湘, 梁谢恩, 张汝民, 温国胜.  毛竹快速生长期的叶绿素荧光参数特征 . 浙江农林大学学报, 2018, 35(1): 75-80. doi: 10.11833/j.issn.2095-0756.2018.01.010
    [5] 武录义, 岳永杰, 刘果厚, 高润宏, 苏志成.  气候变化对元上都遗址区景观格局的影响 . 浙江农林大学学报, 2016, 33(2): 232-238. doi: 10.11833/j.issn.2095-0756.2016.02.007
    [6] 彭维亮, 顾蕾, 胡晨沛, 周鹏飞, 洪明慧, 李翠琴.  碳标签情景模拟下的消费者低碳竹(木)地板支付意愿 . 浙江农林大学学报, 2015, 32(5): 655-660. doi: 10.11833/j.issn.2095-0756.2015.05.001
    [7] 刘彬彬, 张绍勇, 周月英, 林家羽, 孙芳利.  浸提处理对竹粉霉变性能的影响 . 浙江农林大学学报, 2015, 32(1): 11-17. doi: 10.11833/j.issn.2095-0756.2015.01.002
    [8] 易武英, 苏维词, 周文龙, 唐金刚, 张凤太.  基于元胞自动机模型的贵阳市花溪区生态安全预警模拟研究 . 浙江农林大学学报, 2015, 32(3): 369-375. doi: 10.11833/j.issn.2095-0756.2015.03.006
    [9] 金明, 丁贵杰.  贵州马尾松单株木二元材种出材率表的编制 . 浙江农林大学学报, 2011, 28(4): 576-582. doi: 10.11833/j.issn.2095-0756.2011.04.009
    [10] 任斌斌, 冯久莹, 李树华.  模拟邯郸地区自然群落的植物景观设计 . 浙江农林大学学报, 2011, 28(6): 870-877. doi: 10.11833/j.issn.2095-0756.2011.06.006
    [11] 陈国良, 蒋晶, 乔桂荣, 金梁, 孙月华, 杨晔, 周婧, 卓仁英.  青杨双向电泳实验体系的建立 . 浙江农林大学学报, 2009, 26(5): 644-651.
    [12] 刘昊, 余树全, 江洪, 方江保.  模拟酸雨对山核桃叶绿素荧光参数、叶绿素和生长的影响 . 浙江农林大学学报, 2009, 26(1): 32-37.
    [13] 郭金剑, 池伟, 赵鑫珠, 王丽, 曾燕如, 张秋月.  山核桃AFLP实验技术体系的建立 . 浙江农林大学学报, 2008, 25(4): 532-537.
    [14] 马灵飞, 嵇伟兵.  一种实验用长形刨花的制备方法 . 浙江农林大学学报, 2007, 24(1): 96-99.
    [15] 叶玉珠, 王林伟, 王淑媛, 赵沛忠.  板栗粉锈病初步研究 . 浙江农林大学学报, 2003, 20(3): 328-330.
    [16] 田有圳, 黄金桃, 林照授, 涂育合, 叶功富.  凹叶厚朴一元立木材积方程的研究 . 浙江农林大学学报, 2002, 19(3): 255-258.
    [17] 高峰.  植物形态的模拟 . 浙江农林大学学报, 2001, 18(1): 76-79.
    [18] 吴延熊, 周国模, 郭仁鉴, .  区域森林资源系统的“三元论” . 浙江农林大学学报, 1999, 16(1): 28-33.
    [19] 刘志坤.  BX218型削片机三角带传动存在的问题及优化设计 . 浙江农林大学学报, 1995, 12(2): 166-171.
    [20] 童再康, 范义荣.  黄山松种群拓广遗传参数的研究 . 浙江农林大学学报, 1993, 10(1): 43-48.
  • 期刊类型引用(2)

    1. 刘江丽,吕倩,刘俊杰,陈小红,范川,李贤伟. 林窗尺度对马尾松人工林林下灌草层优势种生态位的早期影响. 西北植物学报. 2022(02): 312-325 . 百度学术
    2. 胡文杰,李阁,李冠喜. 马尾松松针挥发油化学成分及抗氧化活性研究. 中国粮油学报. 2018(12): 42-48 . 百度学术

    其他类型引用(3)

  • 加载中
  • 链接本文:

    https://zlxb.zafu.edu.cn/article/doi/10.11833/j.issn.2095-0756.20220662

    https://zlxb.zafu.edu.cn/article/zjnldxxb/2023/4/875

图(5) / 表(6)
计量
  • 文章访问数:  457
  • HTML全文浏览量:  82
  • PDF下载量:  29
  • 被引次数: 5
出版历程
  • 收稿日期:  2022-10-09
  • 修回日期:  2023-05-04
  • 网络出版日期:  2023-07-13
  • 刊出日期:  2023-08-20

基于离散元的竹粉颗粒接触参数标定

doi: 10.11833/j.issn.2095-0756.20220662
    基金项目:  中央级公益性科研院所基本科研业务费专项资金(CAFYBB2020SY040)
    作者简介:

    刘东(ORCID: 0009-0003-9644-1395),从事木质材料增材制造研究。E-mail: 1104496173@qq.com

    通信作者: 闫承琳(ORCID: 0009-0009-0014-8180),副研究员,博士,从事木质材料增材制造技术与人造板机械等研究。E-mail: yanchenglin@126.com
  • 中图分类号: S781.9

摘要:   目的  建立竹粉颗粒材料的离散元模型并标定其相关接触参数。  方法  基于毛竹Phyllostachys edulis粉末物理堆积角实验结果,采用离散元法(DEM)模拟仿真和实验设计(DOE)相结合的方法,通过Plackett-Burman (P-BD)实验、爬坡实验和响应面实验获得竹粉堆积角和仿真参数之间的二次多项式回归模型,以物理堆积角为目标值预测得到接触参数的最优组合。  结果  实测竹粉物理堆积角为49.29°;P-BD实验结果表明了竹粉-竹粉滚动摩擦系数、竹粉-竹粉恢复系数和竹粉-不锈钢板静摩擦系数对堆积角影响显著(P<0.05);响应面实验结果进一步证实了P-BD实验的结果,同时表明竹粉-竹粉滚动摩擦系数自身交互作用、竹粉-竹粉滚动摩擦系数和竹粉-竹粉恢复系数的交互作用对堆积角有显著影响(P<0.05);最优因素组合为竹粉-竹粉恢复系数0.44,竹粉-竹粉滚动摩擦系数0.22,竹粉-不锈钢板静摩擦系数0.45。  结论  本研究在41.20°~55.27°堆积角范围内建立的竹粉颗粒堆积角定量模拟标定与预测平台具有较高的可靠性,DEM标定方法可为后续建立评估生物质颗粒运动行为及其与设备相互作用的材料模型提供参考。图5表6参27

English Abstract

张海涛, 宫渊波, 付万权, 等. 马尾松低效人工林不同改造模式下降雨及产流特征[J]. 浙江农林大学学报, 2018, 35(1): 29-34. DOI: 10.11833/j.issn.2095-0756.2018.01.004
引用本文: 刘东, 刘子昕, 王琦, 等. 基于离散元的竹粉颗粒接触参数标定[J]. 浙江农林大学学报, 2023, 40(4): 875-882. DOI: 10.11833/j.issn.2095-0756.20220662
ZHANG Haitao, GONG Yuanbo, FU Wanquan, et al. Rainfall and runoff in low efficiency Pinus massoniana forests with different silvicultural prescriptions[J]. Journal of Zhejiang A&F University, 2018, 35(1): 29-34. DOI: 10.11833/j.issn.2095-0756.2018.01.004
Citation: LIU Dong, LIU Zixin, WANG Qi, et al. Calibration of contact parameters for bamboo powder particles based on DEM[J]. Journal of Zhejiang A&F University, 2023, 40(4): 875-882. DOI: 10.11833/j.issn.2095-0756.20220662
  • 竹材适应性强、生长期短、分布广,相比木材等其他农林生物质材料,竹材纤维长径比高,硬度和韧性更好,在家具、包装、建筑、人造板等行业占据越来越重要的市场地位[13]。生物质材料的处理、运输、存储和应用依赖于适当的加工设备,而准确的接触参数是实现生物质材料机械化加工的关键,也是相关设备研发和加工参数优化的需要[46]。但生物质材料的物理、力学和化学性质受与颗粒特性相关的形状、弹性、黏性、纤维性等影响,传统的研究方法难以准确获得加工设备与生物质材料粉末之间的接触参数[45]。离散元法(discrete element method,DEM)是设计、建模和研究颗粒系统常用的方法,能有效模拟生物质材料采集、筛选、搅拌、粉碎和运输等过程中的颗粒行为[79]。接触参数标定是颗粒系统离散元仿真计算研究的前提,国内外学者对接触参数的离散元研究主要涉及土壤、肥料、秸秆等材料。ADAJAR等[10]根据环剪实验和直接剪切测量结果对不同作物秸秆摩擦系数、颗粒的法向刚度和剪切刚度进行了标定。ZHAO等[11]以物理实验测定的椰糠实际堆积角为响应值,采用离散元方法获得了椰糠接触参数的最佳组合。XIA等[12]提出了一种层析成像的DEM方法,并建立了碾磨松木颗粒的近似模型,证明了DEM模型适用颗粒状林业生物质材料的可行性。 PACHON-MORALES等[13]分别用多球体法和Hertz-Mindlin with JKR接触模型模拟了杨木粉末颗粒的细长形状和内聚行为。但尚未发现竹材接触参数离散元研究的相关报道。

    本研究对竹粉颗粒的接触参数进行标定,以毛竹Phyllostachys edulis为研究对象,基于离散元模拟仿真和实验设计(design of experiment,DOE)方法,通过Plackett-Burman (P-BD)实验、爬坡实验和响应面实验建立了毛竹粉末接触参数与堆积角的回归模型,并预测最优接触参数组合,最后验证了模拟接触参数的可靠性和准确性,以期为竹粉加工过程的模拟实验研究和相关设备的研发优化提供有效的接触参数。

    • 毛竹采自安徽宣城,竹龄4~5年生,胸径为10.2 cm,密度为680 kg·m−3,静曲强度为156.97 MPa,抗压强度为66.24 MPa,取竹中部位,包含竹青和竹黄。将样品研磨成粉,干燥后,采用不同规格标准筛分得到100目竹粉样品。

    • 采用漏斗法测定竹粉颗粒堆积角,参照相关标准和文献[5, 1416],自制堆积角测定装置如图1所示,短颈漏斗颈口内径10 mm,底盘培养皿直径100 mm,漏斗颈口距底盘培养皿90 mm。测量时将竹粉缓慢倒入并轻微搅动,避免颗粒在漏斗壁上积聚;待料堆高度不再增加时停止添加竹粉。在相同位置使用相机拍摄料堆的正视图像,重复10次并记录。

      图  1  堆积角测定装置

      Figure 1.  Stacking angle measuring device

    • 采用离散元软件EDEM模拟毛竹粉末料堆形成过程,建立毛竹粉末颗粒模型和堆积角模型,并选择适当的接触模型:①颗粒模型。竹粉实际为纤维状,且不同竹粉颗粒长宽比差异较大,测定物理堆积角时采用统一粒径的竹粉是为了降低颗粒模型的复杂程度;在对比实际粒径和放大粒径形成的堆积角后,将竹粉颗粒的粒径(宽)放大至1 mm [17]。建立的竹粉颗粒模型如图2A所示,此颗粒模型的长宽比为2∶1,粒径分布设置为fixed。②堆积角模型。依照物理堆积角测定装置建立几何模型并导入EDEM软件,在漏斗上方建立总数为25 000的颗粒工厂,生成速率1 250个·s−1,静置4 s,总仿真时间设为24 s。待竹粉颗粒自由落下形成稳定的颗粒堆时,截取+X、−X、+Y、−Y共4个方向的料堆图像。建立的竹粉颗粒堆积角模型如图2B所示。③接触模型。毛竹粉末与加工设备之间的颗粒运动特征是密集和缓慢的,加上生物质材料的纤维性和联锁效应,干燥后的样品仍具有轻微的黏附性[1819]。选用EDEM软件提供的JKR Cohesion接触模型,该模型考虑了范德华力在接触区内的影响,也允许构建强黏合系统[16, 18]

      图  2  竹粉颗粒模型(A)和堆积角模型(B)

      Figure 2.  Bamboo powder particle model (A) and powder stacking angle model (B)

    • 利用图像处理软件对物理实验和仿真实验获得的堆积角图像进行处理,调用灰度函数rgb2gray、二值函数im2bw、形态学运算函数strel (square和disk)得到内部填充、边缘柔和的堆积角图像;在GUI界面利用canny算子提取料堆边缘,并利用最小二乘法进行一阶线性拟合得到斜率(精确到0.001);最后将斜率(正切值)转换为角度,堆积角处理过程如图3所示。实测物理堆积角为49.29°。

      图  3  堆积角图像处理

      Figure 3.  Image processing of stacking angle

    • 进行离散元模拟研究所需的仿真参数主要有颗粒密度、泊松比、恢复系数、剪切模量、静摩擦系数、动摩擦系数和表面能等[7, 20]。壁面材料采用不锈钢板,其密度为7 850 kg·m−3、泊松比为0.3、剪切模量为7×1010 Pa,而竹粉的接触参数难以通过物理实验测量,需通过模拟实验进行标定[7, 14]。结合预模拟实验、相关农林生物质材料参数标定文献[6, 11, 2124],本研究中竹粉和壁面材料设置仿真参数见表1。影响堆积效果的因素很多,先进行Plackett-Burman (P-BD)实验筛选显著因素,再进行爬坡实验缩小显著因素的参数范围,然后进行响应面实验获得堆积角和DEM参数之间的二次多项式回归模型,最后使用预测得到标定参数的最优组合进行验证。

      表 1  竹粉堆积角测定仿真参数

      Table 1.  Simulation parameters for measuring the stacking angle of bamboo powder

      参数符号参数名称取值范围参数水平
      −11
      μ 竹粉泊松比 0.20~0.50 0.20 0.40
      E 竹粉剪切模量/GPa 1~15 2 4
      X1 竹粉-竹粉恢复系数 0.15~0.30 0.25 0.50
      X2 竹粉-竹粉静摩擦系数 0.20~0.80 0.30 0.60
      X3 竹粉-竹粉滚动摩擦系数 0~0.40 0.15 0.30
      X4 竹粉-不锈钢板恢复系数 0.10~0.80 0.25 0.50
      X5 竹粉-不锈钢板静摩擦系数 0.30~0.70 0.30 0.60
      X6 竹粉-不锈钢板滚动摩擦系数 0~0.40 0.15 0.30
      J 表面能 /(J·m−2) 0.01~0.10 0.04 0.08

      ① P-BD实验。采用Minitab软件9因子表进行P-BD实验设计,总计13组(包括1个中心点),实验方案及结果见表2

      表 2  竹粉堆积角P-BD实验设计及结果

      Table 2.  Plackett-Burman test design and results of bamboo powder stacking angle

      实验号μEX1X2X3X4X5X6J模拟堆积角/(°)
      1−1−1111−111−151.54
      200000000047.17
      3−1−1−1111−11154.48
      4−1−1−1−1−1−1−1−1−142.63
      51−11−1−1−111141.68
      6−111−11−1−1−1146.47
      7111−111−11−137.22
      8−11−1−1−1111−140.56
      911−111−11−1−157.91
      101−111−11−1−1−123.71
      11−1111−111−1139.00
      1211−11−1−1−11137.79
      131−1−1−1111−1161.32
        说明:符号所表示参数名称见表1

      ② 爬坡实验。根据P-BD实验分析结果设计爬坡实验。显著因素效应值为正,则水平逐次增加;效应值为负,则水平逐渐减少[11, 2526]。一次爬坡实验设置步长为0.1,二次爬坡实验设置步长为0.05,分别进行5组实验,爬坡实验方案及结果见表3

      表 3  竹粉堆积角爬坡实验设计及结果

      Table 3.  Climbing test design and results of bamboo powder stacking angle

      实验号一次爬坡二次爬坡
      X1X3X5模拟堆积角/(°)相对误差/%X1X3X5模拟堆积角/(°)相对误差/%
      10.500.150.3038.9221.040.500.150.3039.6719.52
      20.400.250.4048.49 4.500.450.200.3541.5315.75
      30.300.350.5057.9517.570.400.250.4049.325.20
      40.200.450.600.350.300.4552.827.17
      50.100.550.700.300.350.5057.2716.19
        说明:−表示数据无效。符号所表示参数名称见表1

      ③ 响应面实验。采用Design Expert软件中的Box-Behnken方法(B-BD)进行响应面实验。选择与实际堆积角误差最小的一组参数作为响应面实验设计的水平值,根据参数个数在Design Expert软件中生成响应面实验顺序表(加3个中心点),响应面实验方案及结果见表4

      表 4  竹粉堆积角响应面实验设计及结果

      Table 4.  Box-Behnken test design and results of bamboo powder stacking angle

      序号X1X3X5模拟堆积角/(°)
      10(0.40)1(0.30)−1(0.35)48.69
      200(0.25)0(0.40)48.89
      31(0.45)01(0.45)52.17
      410−143.18
      5−1(0.35)0−148.71
      600048.51
      7−11054.98
      80−1(0.20)−141.53
      9−10154.97
      100−1148.85
      111−1041.20
      1200050.54
      1311051.69
      1401155.27
      15−1−1048.11
        说明:括号内数值为参数值。参数X1X3X5所表示名称见表1
    • 表5所示:模型的决定系数(R2)为0.989 3,表明该模型可用于解释因素对响应值的影响。X1X3X5对竹粉的堆积角有显著影响(P<0.05),影响堆积角的显著因素从大到小依次为X3X1X5,其他因素影响较小。原因为干燥后的竹粉黏附性较低,因此竹粉与竹粉之间运动状态主要取决于运动方式和碰撞力度,采用多球体方法构建的颗粒模型之间多为滚动;非球形颗粒的不规则形状减小了其流动性,因此影响竹粉与不锈钢板之间运动状态的主要原因是静摩擦。

      表 5  竹粉堆积角P-BD实验参数显著性分析

      Table 5.  Significance parameters analysis of bamboo powder stacking angle Plackett-Burman test

      参数效应系数P显著性排序
      μ−2.506−1.2530.2307
      E−2.735−1.3680.2036
      X1−9.179−4.5890.0252
      X2−0.904−0.4520.6019
      X313.9276.9630.0111
      X4−3.623−1.8110.1325
      X58.2864.1430.0303
      X6−1.297−0.6490.4708
      J4.5272.2640.0914
        说明:参数所表示名称见表1
    • 表3可以看出:一次爬坡实验前3组堆积角相对误差分别为21.04%、4.50%、17.57%,呈先减小后增大的趋势;后2组实验的堆积图像无效。原因为恢复系数过小,静摩擦系数过大,颗粒不易运动,造成无效堆叠,堆积效果如图4。缩短步长后二次爬坡实验5组实验堆积角均有效,且根据相对误差可知最优区间在3组附近,2、4组之间,则选择2、3、4组实验的参数值作为响应面实验的各水平值。

      图  4  有效堆积图像(A)和无效堆积图像(B)

      Figure 4.  Effective (A) and invalid (B) stacking angles image

    • 根据表4的数据,建立了竹粉堆积角和显著因素的二次回归模型:$ {\theta }_{1} $ =50.54−2.19 X1+3.69 X3+3.59 X5+1.17 X1X3+0.68 X1X5−0.087 X3X5−0.37 X12−1.44 X32−0.42 X52。从表6可以看出:模型P<0.000 1,表明该模型用于描述参数与响应值之间的关系极显著,可以根据模型预测堆积角;变异系数(CV)为1.18%,表明该测试具有很高的可靠性;R2为0.993 6,修正决定系数($R_{{\rm{adj}}}^2 $)为0.982 2,这些数值表明回归方程的拟合度较好(一般认为2个系数需大于80%),回归方程可以用来代替实际测试来分析结果;精度(AP)为30.541,表明该模型具有良好的精度。

      表 6  Box-Behnken实验回归模型方差分析

      Table 6.  ANOVA of modified model of Box-Behnken

      方差来源均方自由度平方和FP
      模型 265.87 9 29.54 86.60 <0.000 1
      X1 38.22 1 38.22 112.04 <0.000 1
      X3 108.81 1 108.81 318.97 <0.000 1
      X5 103.35 1 103.35 302.96 <0.000 1
      X1X3 5.44 1 5.44 15.93 0.010 4
      X1X5 1.86 1 1.86 5.45 0.066 9
      X3X5 0.03 1 0.03 0.09 0.776 9
      X1X1 0.49 1 0.49 1.44 0.283 2
      X3X3 7.63 1 7.63 22.36 0.005 2
      X5X5 0.64 1 0.64 1.88 0.228 4
      残差 1.71 5 0.34
      纯误差 0.015 2 0.007
      合计 267.58 14
      R2=0.9936,Radj2=0.9822,CV=1.18,AP=30.541
        说明:参数X1X3X5所表示名称见表1

      表6可以看出:线性项X1X3X5P均小于0.000 1,与P-BD实验结果完全一致,在模拟实验参考范围内对堆积角影响极显著;交互项X1X3P为0.010 4,对堆积角有显著影响,其余交互项X1X5X3X5对堆积角影响很小(P<0.05);二次项X3X3P为0.005 2,对堆积角有显著影响,X1X1X5X5对堆积角影响很小(P<0.05)。

      X5为0.40时,利用Design-Expert软件得到X1X3的响应曲面,如图5所示。在X1X3作用下竹粉堆积角整体表现为增大趋势,原因为竹粉-竹粉滚动摩擦系数的增大增加了竹粉颗粒间接触的摩擦阻力,使颗粒堆不易滑散、趋于稳定;而竹粉-竹粉恢复系数的增大增加了颗粒间的碰撞反弹程度,使颗粒堆易飞散,不易成形。但竹粉-竹粉滚动摩擦系数的响应曲面曲线较竹粉-竹粉恢复系数的响应曲面曲线更陡峭,这表明竹粉-竹粉滚动摩擦系数对堆积角的影响比竹粉-竹粉恢复系数的影响更大,这也与P-BD实验得到的结论一致。

      图  5  交互项X1X3响应面

      Figure 5.  Response surface of interaction term X1X3

    • 利用Design-Expert软件的预测模块,以物理堆积角49.29°作为目标值,选择一组与实际物理堆积角度数据相近的最优解,即竹粉-竹粉恢复系数X1=0.44,竹粉-竹粉滚动摩擦系数X3=0.22,竹粉-不锈钢板静摩擦系数X5=0.45,其余非显著参数取值与爬坡实验相同,即泊松比μ=0.33,剪切模量E=0.033 GPa,竹粉-竹粉静摩擦系数X2=0.4,竹粉-不锈钢板恢复系数X4=0.3,竹粉-不锈钢板滚动摩擦系数X6=0.1,表面能J=0.01 J·m−2

      采用上述标定的竹粉最优参数组合进行验证,测得堆积角为48.06°,相对误差为2.50%,误差在可接受的范围内;再选择2组测定实验目标值(47.95°和48.52°),预测得到最优解后进行验证,分别测得堆积角为46.58°和46.65°,相对误差分别为2.86%和3.85%,进一步证明了该模型在竹粉颗粒参数标定中的有效性。

    • 农林生物质材料的含水率影响其微观结构和物理性质,进而影响材料的生产、运输、存储和应用[24]。含水率受样品年龄、采样部位等多因素的影响,但在含水率较低的情况下,木材散碎物料的堆积特性变化不大[27]。本研究采用竹中部位的竹段,磨粉干燥后控制含水率低于2%,所以本研究未设计含水率对模拟堆积角的影响实验。

      模拟不同比例、不同粒径、不同长宽比的竹粉纤维会增加模拟实验的复杂性,因此构建4球面、长宽比2∶1的颗粒模型,并选择fixed平均分布,与物理堆积角采用统一粒径竹粉对应;当粒径>1.5 mm时会堵塞漏斗下端出口,而真实粒径0.15 mm需要的颗粒总量会按量级增加,相应的计算时间也会增加。当表面能超过0.12 J·m−2时,颗粒会黏附在漏斗下端出口造成堵塞;当摩擦系数超过0.85时,颗粒落到料堆上运动状态阻碍太大造成不规则团块堆叠。根据实际情况合理简化颗粒形状和缩放粒径,在正式实验前需要做预实验选择接触模型和参数范围。

      剪切模量量级、颗粒生成速率和计算机硬件配置都会影响仿真时间。在P-BD实验得出剪切模量对堆积角的影响很小的结论后,将剪切模量降低量级到107。在不同配置的中央处理器(CPU)和图像处理器(GPU)上仿真,并试验颗粒生成速率5 000、2 500、1 250个·s−1时的仿真时间,发现GPU比CPU的速度稳定,例如12核CPU仿真时间约46 h,而6核GPU仿真时间约22 h;减小颗粒生成速率初期仿真较慢,但是为了降低颗粒生成位置的计算时间,会缩短后期仿真时间。在颗粒数量更大时,选择GPU是更高效的。

    • 本研究所得主要结论如下:①通过P-BD、爬坡和B-BD响应面等方法,建立了检验因素与堆积角的回归模型,该模型能够准确描述堆积角与各参数的关系。证明了竹粉-竹粉滚动摩擦系数、竹粉-竹粉恢复系数、竹粉-不锈钢板静摩擦系数的线性项,竹粉-竹粉滚动摩擦系数的二次项,竹粉-竹粉滚动摩擦系数和竹粉-竹粉恢复系数的交互项对堆积角有显著影响。②以竹粉物理堆积角为目标,通过求解回归模型获得了最优因素组合,即竹粉-竹粉恢复系数0.44,竹粉-竹粉滚动摩擦系数0.22,竹粉-不锈钢板静摩擦系数0.45。将回归模型的最优解作为DEM参数进行验证,测得堆积角的相对误差在允许范围内。③本研究在41.20°~55.27°堆积角范围内建立的回归多项式方程可以对显著因素进行快速预测,结果具有较高的可靠性。当具体的工作环境和接触模型发生变化时,线性项和二次项的相互作用会增加标定接触参数的难度,最优解需要根据生物质物料的堆积角、料堆形状和动态行为确定。

      本研究得到的回归模型证明了离散元法应用于竹材粉末研究的可行性,生物质材料特性、不同响应指标检验方法、颗粒模型的准确性等是下一步的研究方向。

参考文献 (27)

目录

/

返回文章
返回