技术领域本发明属于三维建模领域，涉及一种重力堆积结构和随机泡沫结构的构建方法，用于仿真多孔介质结构生成过程的三维建模。
背景技术：
燃料的燃烧是人类获取能源的主要途径，但目前全世界面临着燃烧污染严重和化石能源枯竭两大挑战。如何减少燃烧排放和提高燃烧效率已成为当前燃烧领域的两大热点课题，为此国内外学者提出了各种新型的燃烧技术。近二十年来，多孔介质燃烧技术[MujeebuMA,AbdullahMZ,AbuBakarMZ,etal,JournalofEnvironmentalManagemen,2009,90:2287-2312]以其超绝热、低污染等突出特点引起学术界的广泛关注，有可能成为近期、甚至相当长时间内解决这些问题的一项重要技术，为开发和设计新型高效清洁燃烧系统提供了一条广阔的新途径。多孔介质内的燃烧又称为过滤燃烧，是指采用多孔介质材料取代自由空间，使可燃混合气或燃油蒸汽经固体多孔介质或颗粒堆积床，在类似于过滤情况下发生的燃烧。过滤燃烧是自然界和工程中广泛存在的一种燃烧现象，其广泛性和重要性绝不亚于常见的自由空间燃烧。过滤燃烧与自由空间燃烧的差异在于多孔介质的存在，而多孔介质具有适用于燃烧的多个优良特性[M.M.Kamal,A.A.Mohamad,JournalofPowerandEnergy220(2006)487–508]，具体表现在：1)多孔介质独具的热能积累和反馈效应。2)多孔介质对气体的扰流和弥散作用。3)多孔介质内不存在宏观尺度上的燃烧火焰等。近年来，采用多孔介质燃烧技术燃烧液体燃料作为一项新型燃烧技术受到人们的重视[MAMujeebu,AbdullahMZ,AbuBakarMZ,etal.ProgressinEnergyandCombustionScience35(2009)216–230]，已在多孔介质发动[M.Weclas,J.Cypris,andT.M.A.Maksoud.Proceedings(CD)ofthe4thInternationalConferenceonPorousMediaanditsApplicationsinScienceandEngineeringICPM4June17-22,2012,Potsdam,Germany]、工业有毒废液焚烧[PastoreA,MastorakosE.Fuel,(2011)90:64-76]、微型燃烧器[JWLi,JHHuang,MiYan,DZhao.ExperimentalThermalandFluidScience,2014,52:47–58]等方面得到了较为广泛实际应用，在经济建设和社会发展中有着广阔的应用前景。但目前人们对过滤燃烧的认识还远远不够深入和全面，特别是许多基础性的科学问题，有待深入研究和探明。多孔介质是一种典型的多尺度结构的复杂几何系统，其孔隙尺度从毫米级跨越到纳米级，要深入理解多孔介质内的过滤燃烧机理是非常困难的。长期以来对过滤燃烧的研究，无论是理论分析还是数值模拟，一般都是把多孔介质视为均匀的连续介质[NijemeislandM.AICHEJournal,2006,50(5):906-921]，将多孔介质布置为对称结构或有序排列。这种模型与实际多孔介质的真实结构相差较远，多孔介质内的蒸发、传热传质及燃烧过程被掩盖，不利于揭示蒸发过程的微观机理和过滤燃烧的本质。因此建立准确的多孔介质几何模型是研究微孔流体流动与传热的前提[MullerGE.PowerTechnology,2010,23:626-633]，在此基础上得到流场与温度场的详细信息才更符合实际情况。迄今为止，各国研究者常采用的多孔介质的几何结构重构方法有：物理方法(如核磁共振成像)[SongYC,JiangLL,LiuY,YangMJ,AbudulaA.InternationalJournalofGreenhouseGasControl,2012,10:501-509]、数值方法[RoozbahaniMM,HuatBBK,AsadiA.AdvancedPowerTechnology,2013,24(1):26-35]、分形理论[TongjunMiao,BomingYu,YonggangDuan,QuantangFang,InternationalCommunicationsinHeatandMassTransfer58(2014)71-78]和基于网格微结构单元体的算法[LiouMF,KimH.ComputationalFluidDynamics,2009,8:363-374)等，由于多孔介质复杂的结构，因此生成的网格数量巨大，对计算要求较高。随电子运算能力日益发达，离散元方法的兴起推动了颗粒模拟软件的发展。近期，国际上出现一种开源的基于离散单元法的LIGGGHTS离散元软件[C.Kloss,C.Goniva,LIGGGHTS:anewopensourcediscreteelementsimulationsoftware,Proc.5thInt.Conf.onDiscreteElementMethods,London,UK,August2010,(website:www.liggghts.com)]，在充分考虑堆积过程中重力、球间碰撞及壁面碰撞引起的轨迹变化，构建重力作用下自由堆积的几何结构[JidongZhao,TongShan,PowderTechnology,2013,239:248-258]，但作为一种新型建模方式，该方法尚未受到广泛关注。
技术实现要素：
针对现有技术的不足，本发明提供一种重力堆积结构和随机泡沫结构的构建方法，本发明以离散元方法为基础，以随机多孔介质内过滤燃烧为研究背景，通过对LIGGGHTS离散元软件的改进与优化，建立多孔介质三维随机结构模型，为多孔介质燃烧技术的理论模型和实际堆积床建模提供技术指导。为达到上述目的，本发明的技术方案为：一种重力堆积结构和随机泡沫结构的构建方法，主要包括以下步骤：第一步，在离散元软件LIGGGHTS中，输入初始化参数，初始化参数包括堆积单元体几何特征、单元体材料力学性质、堆积床区域的几何特征、壁面硬度、下落高度和下落区面积。1.1堆积单元体几何特征利用rand随机函数定义单元体的大小，所述的单元体为长方体或椭球体。长方体三棱比例固定，椭球轴间比例固定，设定单元体的特征尺寸(长方体特征尺寸为三棱棱长，椭球体特征尺寸为长短轴轴长)，随机生成多个不同大小的单元体，所述的长方体最大棱长为0.8mm-1.2mm；椭球体长轴长度为0.8mm-1.2mm。1.2用离散元软件LIGGGHTS定义单元体材料力学性质，所述的力学性质包括硬度和弹性形变率。1.3堆积床区域的几何特征在离散元软件LIGGGHTS中，将堆积床区域设置为长方体或圆柱体，设定堆积床区域的几何尺寸。为了防止近壁面处堆积体的稀疏度对整体堆积床区域产生影响，堆积床区域的几何尺寸大于单元体最大特征尺寸的4倍；为了在后续的计算研究中得到有效结果，圆柱体堆积床区域的轴向长度大于单元体最大特征尺寸的8倍，长方体堆积床区域的最大边长大于单元体最大特征尺寸的8倍。1.4在离散元软件LIGGGHTS中，定义堆积床区域壁面硬度、单元体下落高度和下落区域面积。当堆积床区域为长方体时，单元体的下落高度高于长方体最大边长的1.5倍，当堆积床区域为圆柱体时，单元体的下落高度高于圆柱体轴向长度的1.5倍，保证单元体下落后能够随机分布。所述的下落区域面积小于堆积床区域的截面积，以保证小球均投入堆积床区域中，下落区域面积大于单元体的最大截面积，保证顺利生成堆积床结构。第二步，在离散元软件LIGGGHTS中设定重力、浮力、磁场力、范德瓦尔力的力场模型，建立力平衡方程；结合单元体材料力学性质，对每个时间步长中单元体受力变形进行分析，得到随时间变化的单元体形变值δ，这与实际结构更为接近；从模拟角度看，考虑单元体形变，使单元体之间具有接触面，可避免接触处为点接触所导致的网格无法划分的情况。第三步，通过离散元软件LIGGGHTS，确定单元体的虚构边界，单元体按初始化参数下落，在单元体下落过程中进行接触搜索与检测，确定单元体是否碰撞，所述的碰撞包括单元体之间的相互碰撞、单元体与壁面之间的碰撞，碰撞后单元体形变由第二步的单元体形变值确定，单元体碰撞停止后得到堆积床模型，形成堆积床数据信息。3.1确定单元体的虚构边界为了更好地检测非球体结构的碰撞情况，对非球体结构设置虚构边界，本发明中，将单元体的虚构边界设置为球体边界，使非球体结构的最大截面积与虚构边界的最大截面积相等。3.2接触搜索与检测采用格子模型搜索方式，检测单元体是否存在接触及接触点位置；所述的格子模型(图4)是将计算区域划分成正方体格子，以各个单元体的上顶点、中点或下顶点为关联点，将球体与格子关联起来，接触检测中只需对被检测单元体和其相邻格子内的单元体进行检测即可，这样既节省了计算资源又可以保证接触点检测的精度。而传统暴力搜索采用的格子搜索法(图3)，在检测某个单元体的接触点时，需要对该单元体及其他单元体一一进行检测，一次检测需要耗费庞大的计算资源与时间。第四步，将堆积床数据信息导入到Matlab软件中，提取出各单元体的几何中心坐标、特征尺寸和转角变量，得到各单元体对应的画图指令。所述的球形单元体的特征尺寸为半径，椭球单元体的特称尺寸为长轴和短轴，多面体单元体的特称尺寸为所有的边长，而且椭球及多面体单元体均需要确定多个转角变量以确定其最终堆积形态；将几何中心坐标、特征尺寸和转角变量建立为数组B(i，j)，其中i为单元体ID，j为数组ID，当j＝1时，B(i，j)表示第i个单元体的几何中心坐标；当j＝2，3...，1+n时，B(i，j)表示第i个单元体的特征尺寸,n为该单元体特征尺寸个数；当j＝2+n，3+n...，1+n+m时，B(i，j)表示第i个单元体的转角变量，m为该单元体转角变量个数。在Matlab软件中，通过for循环，将所有堆积床单元体数组B(i，j)的信息，通过fprintf语句生成gambit软件可识别的建模指令:{'volumecmove\volume.1\multiple1offset％d％d％d',B(i,1),B(i,2),B(i,3)