CN113343474B
有效
一种碎石增强高聚物复合材料的细观分析方法
技术领域
[0001]本发明属于高聚物应用工程领域，具体涉及一种碎石增强高聚物复合材料的细观分析方法。
背景技术
[0002]碎石增强高聚物复合材料(Crushed-stone-reinforced Polymer，SRP)应用在复合地基中，负责承受重量和抵抗压力、冲击等荷载作用，其弹性模量和强度等力学性能是考察该新型材料特性和应用性的重要指标。SRP力学性能的研究是该新型材料应用和发展的前提，已有研究通过力学试验、理论模型和有限元模型对SRP的基本力学性能以及细观组份和宏观力学性能的关系进行了初步探索。
[0003]SRP和高聚物的强度和弹性模量均与基体高聚物的密度呈线性关系。SRP的抗压强度是等标称密度高聚物的2.5～3倍，弹性模量约为等标称密度高聚物的10倍。现有的关于SRP的研究，主要集中在SRP的宏观力学试验和高聚物碎石桩的模型试验及数值分析。关于
[0004]该材料的组成、结构和力学性能等属性值得进一步研究。因此，需要采用细观力学方法，建立细观结构与宏观性能之间的联系，探究SRP的力学性能和破坏机理，为SRP的工程应用提供依据，具有重要的科研意义和工程价值。
发明内容
[0005]本发明提出了一种碎石增强高聚物复合材料的细观分析方法，通过SRP中的碎石的结构特征，建立SRP数值模型和结构模型，进一步的通过对数值模型的处理，根据SRP内组成单元的本构关系、破坏准则和材料属性，计算SRP的力学性能。
[0006]为实现上述目的，本发明提供了如下方案：
[0007]一种碎石增强高聚物复合材料的细观分析方法，包括如下步骤：
[0008]获取碎石增强高聚物复合材料中碎石的结构特征的参数值和统计模型；
[0009]基于所述参数值和所述统计模型，建立基于统计相似性的随机碎石数值模型；
[0010]对所述随机碎石数值模型中的所述碎石进行扩展，建立所述碎石增强高聚物复合材料的组成单元的坐标位置和组分，生成所述碎石增强高聚物复合材料的结构模型；
[0011]基于所述结构模型和给定的所述组成单元的材料属性，对所述随机碎石数值模型进行网格单元剖分，根据所述组成单元的受拉状态和/或受压状态，建立所述组成单元的本构关系和破坏准则；
[0012]建立所述随机碎石数值模型的边界条件，根据所述本构关系、所述破坏准则和所述组成单元的所述材料属性，计算所述碎石增强高聚物复合材料的力学性能。
[0013]优选的，获取所述碎石的结构特征的参数值和统计模型的方法包括：
[0014]沿所述碎石增强高聚物复合材料的底面直径截开，获取所述碎石增强高聚物复合材料的截面图像；
[0015]对所述截面图像中的所述碎石界面进行等效变换，获取所述参数值；
[0016]获取所述截面图像中所述碎石的所述结构特征的变量值，基于所述变量值建立所述统计模型。
[0017]优选的，建立所述统计模型的方法包括：
[0018]对所述截面图像进行图片重建，得到重建图像；
[0019]对所述重建图像内相互接触的所述碎石图像进行分割，得到散列图像；
[0020]对所述散列图像进行像素与单位转换，得到所述碎石的所述结构特征的变量值；
[0021]基于所述结构特征的所述变量值，建立所述结构特征的所述统计模型。
[0022]优选的，建立所述随机碎石数值模型的方法包括：
[0023]基于所述参数值和所述统计模型，按照预设的碎石形心距离，采用蒙特卡罗法随机生成碎石形心，构建基于统计相似性的所述随机碎石数值模型。
[0024]优选的，生成所述碎石增强高聚物复合材料的所述结构模型的方法包括：
[0025]对所述随机碎石数值模型中的所述碎石进行扩展，得到所述组成单元的位置坐标；
[0026]基于所述组成单元的所述位置坐标，建立所述组成单元的组分；
[0027]基于所述位置坐标和所述组分，生成所述碎石增强高聚物复合材料的所述结构模型。
[0028]优选的，所述网格单元剖分的方法包括：
[0029]基于所述结构模型和给定的所述组成单元的所述材料属性，对所述随机碎石数值模型按照预设尺寸进行网格剖分，并对所述组成单元进行等分剖分。
[0030]优选的，所述边界条件的建立方法包括：
[0031]采用位移加载法，对所述随机碎石数值模型的加载端位置施加竖向位移荷载，并对所述随机碎石数值模型的非加载端位置施加位移固定约束。
[0032]优选的，所述力学性能的计算方法包括：
[0033]502.按照预设的加载位移步长，对所述加载端位置施加所述竖向位移荷载；
[0034]504.根据所述本构关系和所述破坏准则，判断所述组成单元是否破坏，如果所述组成单元没有破坏，转入506；如果所述组成单元有破坏，改变所述组成单元的所述材料属性，重复502-504；
[0035]506.计算所述碎石增强高聚物复合材料的力学性能；如果没有达到预设的加载步数，改变所述加载位移步长，重复502-506。
[0036]本发明的有益效果为：
[0037]本发明公开了一种碎石增强高聚物复合材料的细观分析方法，通过获取SRP中碎石的结构特征的参数值和统计模型，建立了基于统计相似性的随机碎石数值模型，该数值模型能够反映真实界面中的碎石形貌。进一步的，根据SRP结构模型、SRP组成单元的材料属性和受力状态，建立组成单元的本构关系和破坏准则，完成SRP的力学性能计算。本发明的细观分析方法，具有较强的适用性和稳定性，探究了SRP的力学性能和破坏机理，为建立更加完善的SRP细观数值模型，找到微观结构与宏观性能之间的定量关系，提供了技术基础，为SRP的工程应用提供了技术依据，具有重要的科研意义和工程价值，具有广阔的推广空间和使用价值。
附图说明
[0038]为了更清楚地说明本发明的技术方案，下面对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
[0039]图1为本发明实施例一种碎石增强高聚物复合材料的细观分析方法流程示意图；
[0040]图2为本发明实施例建立的碎石的结构特征统计模型的步骤示意图；
[0041]图3为本发明实施例建立的碎石椭圆度统计模型；
[0042]图4为本发明实施例建立的碎石粒径统计模型；
[0043]图5为本发明实施例建立的Matlab模型与真实截面图像的对比；
[0044]图6为本发明实施例的SRP整体数值模型网格剖分示意图；
[0045]图7为本发明实施例建立的碎石的本构关系示意图；
[0046]图8为本发明实施例建立的界面和基体高聚物的本构关系示意图；
[0047]图9为本发明实施例中SRP试件数值模拟加载流程示意图；
[0048]图10为本发明实施例建立的SRP圆柱体试件截面轴心抗压模型中边界条件设置的示意图；
[0049]图11为本发明实施例中验证网格剖分尺寸敏感性时的四种网格尺寸下SRP试件轴心受压破坏示意图；
[0050]图12为本发明实施例中验证本数值模型合理性时的四种不同碎石分别的SRP试件和轴心受压破坏示意图。
具体实施方式
[0051]下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0052]为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对本发明作进一步详细的说明。
[0053]本实施例提供了一种碎石增强高聚物复合材料的细观分析方法，如图1所示。
[0054]S101.获取SRP中碎石的结构特征的参数值和统计模型；
[0055]由于SRP为碎石夹杂复合材料，决定其内部碎石结构特征的统计变量和相应的特征参数可以根据实际需要进行选择，例如，采用碎石形状、碎石尺寸、碎石含量。在本实施例中，采用截面碎石的椭圆度、粒径和含量来表征上述结构特征的统计变量。根据物理量的意义，各变量的定义及计算式如下：
[0056](1)截面碎石椭圆度(e)：将各碎石截面等效为椭圆，定义截面碎石椭圆度为椭圆短轴与长轴的比值，由式(1)计算。
[0057]
[0058]式中，dmin为单个碎石截面形心的最小直径，dmax为单个碎石截面形心的最大直径。若e＝1时，则碎石截面为圆形。
[0059](2)截面碎石粒径(de)：统计截面内碎石的尺寸，将每个碎石等效为圆形，截面碎石粒径由圆的等效直径表示，即式(2)所示
[0060]
[0061]式中，Ae为截面内单个碎石面积。
[0062](3)截面碎石含量(fe)：SRP试件截面内各碎石总面积占截面总面积的百分比，如式(3)所示：
[0063]
[0064]其中，AE为截面内各碎石总面积，Af为SRP截面的总面积。
[0065]如图2所示，在本实施例中，将SRP试件沿底面直径截开，根据上述碎石结构变量的定义，提取截面图像中各结构变量的几何参数，建立碎石的结构特征的统计模型，主要步骤如下：
[0066](1)图像重建，例如，利用ImageJ和MATLAB软件处理图像，为了便于区分不同结构单元，可以填充不同颜色，例如，截面碎石用黑色填充表示，基体高聚物部分用白色填充表示；
[0067](2)碎石分割，将相互接触的碎石进行分割，例如，使用Image-ProPlus软件进行分割操作；
[0068](3)像素分析与单位转换，继续使用Image-ProPlus软件，对图像进行像素分析与单位转换；
[0069](4)碎石结构变量统计，由图像重建步骤处理后的图像可求得碎石含量。根据截面碎石椭圆度和粒径的定义，由碎石分割后的图像求得相应的结构变量值。
[0070]统计所有已测样本，建立碎石结构特征的统计模型，如图3-4所示。
[0071]S102.建立SRP数值模型
[0072]基于图形处理软件在二维平面上进行碎石投放，例如Matlab，以真实截面碎石的椭圆度、粒径和含量三个方面为依据，将截面碎石真实形貌简化考虑为多边形，构建数值模型。在本实施例中，为了使SRP数值模型与真实细观结构更加接近，在确定了截面碎石的椭圆度、粒径和含量后，例如，假设碎石在SRP数值模型截面上位置的分布具有均匀随机性，可以借助蒙特卡罗法，考虑数值模型中碎石的随机分布，构建基于统计相似性的随机碎石数值模型。本实施例的构造过程如下：
[0073](1)构造长为30mm、宽为15mm的平面空间，空间内碎石的形心采用蒙特卡罗法随机生成，每个形心，均验证其与已有形心的距离，若不满足两形心距离大于对应两多边形的半径和，则重新生成，直至所有形心满足要求；
[0074](2)模型截面内碎石椭圆度取值服从图3的频率分布。
[0075](3)模型截面内碎石粒径取值服从图4的频率分布。
[0076](4)模型截面内碎石含量取50％。
[0077]通过上述步骤建立的数值模型，与真实截面图像相比较，如图5所示，截面中碎石形貌基本相符。
[0078]S103.建立SRP细观结构模型
[0079]进一步的，在本实施例中，对已生成的碎石颗粒向外扩展一定宽度，生成高聚物-碎石界面，提取碎石及界面多边形顶点参数。在获取各相材料位置坐标后，进一步完成基体高聚物、碎石和界面各相组分的建立，生成SRP细观结构模型，例如，使用ANSYS软件参数化
[0080]设计语言(APDL)。需要特别提出的是，SRP的界面过渡区，存在两种类型，一种类型的密度低于基体高聚物，会使宏观力学性能低于基体高聚物；另一种类型的密度高于基体高聚物，会使宏观力学性能高于基体高聚物。
[0081]S104.对SRP数值模型进行网格单元剖分，建立各组成单元的本构关系和破坏准则
[0082]在使用ANSYS软件软件对SRP进行细观数值模拟时，发现网格的尺寸大小不但影响计算的精度，而且在极大程度上决定了计算规模和计算的时间长短。对此，在本实施例中，基于SRP细观结构及材料参数，对SRP整体数值模型进行网格剖分，网格单元尺寸均选取为
[0083]2mm(经验证，2mm为最优尺寸，后续有详细对比)，网格剖分单元采用平面182单元，将碎石颗粒、基体高聚物采用三角形单元进行等分剖分，高聚物-碎石界面采用四边形单元进行等分剖分，即根据碎石颗粒粒径大小将其边界划分成数个四边形或三角形，以保证网格划分的均匀性，如图6所示。进一步的，对于不同单元赋予相应的材料属性，同时，为了使具有不同材料属性的单元易于在细观数值计算中进行区分，在本实施例中，将不同材料属性的单元赋予不同的代表颜色予以区分。
[0084]其中，材料参数包括力学性能参数，拟合公式如下：
[0085]E＝763.88ρ2+130.18ρ
[0086]fc＝34.78ρ2+10.16ρ
[0087]ft＝34.38ρ2+2.22ρ
[0088]式中，E为高聚物的弹性模量，单位Mpa；fc为高聚物的抗压强度，单位Mpa；tt为高聚物的抗拉强度，单位Mpa；ρ为高聚物的密度，单位g/cm3。
[0089]通过上述过程可以看出，本实施例建立的SRP细观数值模型，由碎石单元、高聚物单元和界面单元组成。在进行SRP试件力学测试加载过程中，各单元处于受拉状态或受压状态。为了综合考虑单元的两种受力形式，各单元需给出受压本构关系与受拉本构关系，以及
[0090]对应的破坏准则。
[0091]为了便于模拟和计算，在本实施例中，将碎石受压时看作线弹性均质材料，基本的线弹性本构模型符合其特性，受拉时看作弹脆性均质材料，理想的弹脆性本构模型符合其特性，这样，碎石的本构关系如图7所示，碎石单元本构方程如下所示：
[0092]
[0093]式中，ε为单元的第一主应变；εt为单元的极限应变；σ为单元的第一主应力，单位MPa；σm为单元的残余强度，单位MPa；E为单元的初始弹性模量，单位MPa，图7中ft为单元的抗拉强度，单位MPa。
[0094]相应的，将高聚物材料受压时看作是弹塑性均质材料，理想的弹塑性本构模型符合其特性，受拉时界面和基体的高聚物是SRP中裂缝演化发展的主要空间，适合采用线弹性本构关系，这样，截面和基体高聚物的本构关系如图8所示，界面和基体高聚物单元本构方程如下所示：
[0095]
[0096]式中：ε为单元的第一主应变；εt为单元的极限应变；εc为单元的屈服应变；σ为单元的第一主应力，单位MPa；fc为单元的屈服强度，单位MPa；E为单元的初始弹性模量，单位MPa，σm为单元的残余强度，单位MPa；图8中ft为单元的抗拉强度，单位MPa。
[0097]通过以上过程，建立了SRP细观数值模型，该模型能够真实反映界面中的碎石形貌。进一步的，建立了各组成单元的本构关系和破坏准则，为接下来SRP的力学性能计算奠定了基础。
[0098]S105.计算SRP力学性能
[0099]在建立了SRP数值模型和本构关系、破坏准则后，开始测试SRP的力学性能。
[0100]在本实施例中，利用ANSYS有限元技术建立了SRP试件模型，SRP试件采用圆柱体状，实施轴心抗压试验，采用位移加载方式，如图9所示，测试步骤如下：
[0101]S201.施加位移荷载和固定约束，即在加载端位置施加Y轴负向竖向位移荷载，非加载端位置施加位移全固定约束，设定加载步数和每步加载位移大小，即步长；如图10所示；
[0102]S202.对SRP试件顶部施加位移荷载，根据前述的破坏准则和本构关系，判断每步加载后的试件中是否有单元被破坏；
[0103]对有破坏单元产生的加载步，改变破坏单元的材料属性，重复S201-S202，直到没有新的破坏单元再出现后，转入下一荷载步；对于没有破坏单元产生的加载步则直接转入下一荷载步；
[0104]S203.输出应力应变数据；
[0105]S204.判断是否达到加载步数；如果达到加载步数，求解终止；如果没有达到加载步数，重复S201-S204。
[0106]通过以上过程，能够较好的预测和揭示SRP承受外部荷载作用时的弹性模量、强度变化和破坏过程，也证明了本数值模型的准确性和适用性，为SRP的工程应用提供了技术依据。
[0107]SRP的弹性模量和强度特性是该新型材料最重要的性质，也是实际应用中最关键的因素。SRP力学性能的细观数值模拟能够较好的预测和揭示其承受外部荷载作用时的弹性模量、强度变化和破坏过程，对SRP特性的探索有十分重要的意义。如前所述，在使用ANSYS软件对SRP进行细观数值模拟时，网格的尺寸大小能够影响计算的精确程度，也在很大程度上决定了计算的时间长短和计算规模，因此网格单元剖分对计算结果有着直接影响。
[0108]为了证明网格尺寸的敏感性，在本实施例中，分析了模拟基体高聚物密度为0.45g/cm3，尺寸为Φ150mm×300mm的圆柱体SRP试件的截面受压破坏，SRP试件尺寸和参数，如表1所示：
[0109]表1
[0110]
[0111]对此试件分别进行了四种网格尺寸的数值模拟计算，网格单元尺寸分别为1、2、3、4mm，得到单元剖分情况、计算耗时以及SRP力学性能的变化情况，结果如表2所示：
[0112]表2
[0113]
[0114]可见，随着网格尺寸增大，有限元模型单元数量减少，弹性模量和抗压强度有所上升，力学参数最大变化率均在5％以内。
[0115]图11为四种网格尺寸下试件最终破坏图，可以看出，网格尺寸对SRP轴心压缩裂缝发展没有太大的影响。裂缝主要集中在界面周围，由于中密度基体高聚物中的碎石不会发生破坏，因此，碎石对裂缝的发展有阻断的作用。综合考虑模拟结果的精确性和计算耗时，取网格单元尺寸为2mm最优。
[0116]进一步的，为了验证本数值模型的合理性，在本实施例中，分析了模拟基体高聚物密度为0.45g/cm3，尺寸为Φ150mm×300mm的圆柱体SRP试件的截面轴心受压破坏过程，SRP试件尺寸和参数采用表1所示，利用Matlab随机碎石程序，随机生成四组不同碎石分布的样本，计算其在相同荷载下的SRP轴心受压破坏过程，模拟结果如表3所示：
[0117]表3
[0118]
[0119]可以看出，四组模拟结果离散性较小。最大误差不超过3.29％，最大变异系数不超过1.48％，可以证明本模型具有较高的稳定性。
[0120]图12给出了四组SRP样本的最终破坏图，可以看出，不同碎石分布情况下，裂缝的位置不一定相同，且高聚物-碎石界面过渡区域容易出现裂缝。由于碎石的强度大，单个颗粒所占空间较大，裂缝的破坏位置与碎石的分布集中程度有一定关系：碎石集中的地方，裂缝的扩展区域比较大；该基体密度下，裂缝不会穿过碎石发展，因此中密度试件中碎石对裂缝的发展方向有阻断作用。
[0121]经过上述验证，本数值模型能够较好地对SRP试件的破坏进行数值仿真计算，具有较强的适用性和稳定性。
[0122]以上所述的实施例仅是对本发明优选方式进行的描述，并非对本发明的范围进行限定，在不脱离本发明设计精神的前提下，本领域普通技术人员对本发明的技术方案做出的各种变形和改进，均应落入本发明权利要求书确定的保护范围内。
[0058]式中，dmin为单个碎石截面形心的最小直径，dmax为单个碎石截面形心的最大直径。若e＝1时，则碎石截面为圆形。
[0059](2)截面碎石粒径(de)：统计截面内碎石的尺寸，将每个碎石等效为圆形，截面碎石粒径由圆的等效直径表示，即式(2)所示
[0060]
[0061]式中，Ae为截面内单个碎石面积。
[0062](3)截面碎石含量(fe)：SRP试件截面内各碎石总面积占截面总面积的百分比，如式(3)所示：
[0063]
[0064]其中，AE为截面内各碎石总面积，Af为SRP截面的总面积。
[0065]如图2所示，在本实施例中，将SRP试件沿底面直径截开，根据上述碎石结构变量的定义，提取截面图像中各结构变量的几何参数，建立碎石的结构特征的统计模型，主要步骤如下：
[0066](1)图像重建，例如，利用ImageJ和MATLAB软件处理图像，为了便于区分不同结构单元，可以填充不同颜色，例如，截面碎石用黑色填充表示，基体高聚物部分用白色填充表示；
[0067](2)碎石分割，将相互接触的碎石进行分割，例如，使用Image-ProPlus软件进行分割操作；
[0068](3)像素分析与单位转换，继续使用Image-ProPlus软件，对图像进行像素分析与单位转换；
[0069](4)碎石结构变量统计，由图像重建步骤处理后的图像可求得碎石含量。根据截面碎石椭圆度和粒径的定义，由碎石分割后的图像求得相应的结构变量值。
[0070]统计所有已测样本，建立碎石结构特征的统计模型，如图3-4所示。
[0071]S102.建立SRP数值模型
[0072]基于图形处理软件在二维平面上进行碎石投放，例如Matlab，以真实截面碎石的椭圆度、粒径和含量三个方面为依据，将截面碎石真实形貌简化考虑为多边形，构建数值模型。在本实施例中，为了使SRP数值模型与真实细观结构更加接近，在确定了截面碎石的椭圆度、粒径和含量后，例如，假设碎石在SRP数值模型截面上位置的分布具有均匀随机性，可以借助蒙特卡罗法，考虑数值模型中碎石的随机分布，构建基于统计相似性的随机碎石数值模型。本实施例的构造过程如下：
[0073](1)构造长为30mm、宽为15mm的平面空间，空间内碎石的形心采用蒙特卡罗法随机生成，每个形心，均验证其与已有形心的距离，若不满足两形心距离大于对应两多边形的半径和，则重新生成，直至所有形心满足要求；
[0074](2)模型截面内碎石椭圆度取值服从图3的频率分布。
[0075](3)模型截面内碎石粒径取值服从图4的频率分布。
[0076](4)模型截面内碎石含量取50％。
[0077]通过上述步骤建立的数值模型，与真实截面图像相比较，如图5所示，截面中碎石形貌基本相符。
[0078]S103.建立SRP细观结构模型
[0079]进一步的，在本实施例中，对已生成的碎石颗粒向外扩展一定宽度，生成高聚物-碎石界面，提取碎石及界面多边形顶点参数。在获取各相材料位置坐标后，进一步完成基体高聚物、碎石和界面各相组分的建立，生成SRP细观结构模型，例如，使用ANSYS软件参数化
[0080]设计语言(APDL)。需要特别提出的是，SRP的界面过渡区，存在两种类型，一种类型的密度低于基体高聚物，会使宏观力学性能低于基体高聚物；另一种类型的密度高于基体高聚物，会使宏观力学性能高于基体高聚物。
[0081]S104.对SRP数值模型进行网格单元剖分，建立各组成单元的本构关系和破坏准则
[0082]在使用ANSYS软件软件对SRP进行细观数值模拟时，发现网格的尺寸大小不但影响计算的精度，而且在极大程度上决定了计算规模和计算的时间长短。对此，在本实施例中，基于SRP细观结构及材料参数，对SRP整体数值模型进行网格剖分，网格单元尺寸均选取为
[0083]2mm(经验证，2mm为最优尺寸，后续有详细对比)，网格剖分单元采用平面182单元，将碎石颗粒、基体高聚物采用三角形单元进行等分剖分，高聚物-碎石界面采用四边形单元进行等分剖分，即根据碎石颗粒粒径大小将其边界划分成数个四边形或三角形，以保证网格划分的均匀性，如图6所示。进一步的，对于不同单元赋予相应的材料属性，同时，为了使具有不同材料属性的单元易于在细观数值计算中进行区分，在本实施例中，将不同材料属性的单元赋予不同的代表颜色予以区分。
[0084]其中，材料参数包括力学性能参数，拟合公式如下：
[0085]E＝763.88ρ2+130.18ρ
[0086]fc＝34.78ρ2+10.16ρ
[0087]ft＝34.38ρ2+2.22ρ
[0088]式中，E为高聚物的弹性模量，单位Mpa；fc为高聚物的抗压强度，单位Mpa；tt为高聚物的抗拉强度，单位Mpa；ρ为高聚物的密度，单位g/cm3。
[0089]通过上述过程可以看出，本实施例建立的SRP细观数值模型，由碎石单元、高聚物单元和界面单元组成。在进行SRP试件力学测试加载过程中，各单元处于受拉状态或受压状态。为了综合考虑单元的两种受力形式，各单元需给出受压本构关系与受拉本构关系，以及
[0090]对应的破坏准则。
[0091]为了便于模拟和计算，在本实施例中，将碎石受压时看作线弹性均质材料，基本的线弹性本构模型符合其特性，受拉时看作弹脆性均质材料，理想的弹脆性本构模型符合其特性，这样，碎石的本构关系如图7所示，碎石单元本构方程如下所示：
[0092]
[0093]式中，ε为单元的第一主应变；εt为单元的极限应变；σ为单元的第一主应力，单位MPa；σm为单元的残余强度，单位MPa；E为单元的初始弹性模量，单位MPa，图7中ft为单元的抗拉强度，单位MPa。
[0094]相应的，将高聚物材料受压时看作是弹塑性均质材料，理想的弹塑性本构模型符合其特性，受拉时界面和基体的高聚物是SRP中裂缝演化发展的主要空间，适合采用线弹性本构关系，这样，截面和基体高聚物的本构关系如图8所示，界面和基体高聚物单元本构方程如下所示：
[0095]
[0096]式中：ε为单元的第一主应变；εt为单元的极限应变；εc为单元的屈服应变；σ为单元的第一主应力，单位MPa；fc为单元的屈服强度，单位MPa；E为单元的初始弹性模量，单位MPa，σm为单元的残余强度，单位MPa；图8中ft为单元的抗拉强度，单位MPa。
[0097]通过以上过程，建立了SRP细观数值模型，该模型能够真实反映界面中的碎石形貌。进一步的，建立了各组成单元的本构关系和破坏准则，为接下来SRP的力学性能计算奠定了基础。
[0098]S105.计算SRP力学性能
[0099]在建立了SRP数值模型和本构关系、破坏准则后，开始测试SRP的力学性能。
[0100]在本实施例中，利用ANSYS有限元技术建立了SRP试件模型，SRP试件采用圆柱体状，实施轴心抗压试验，采用位移加载方式，如图9所示，测试步骤如下：
[0101]S201.施加位移荷载和固定约束，即在加载端位置施加Y轴负向竖向位移荷载，非加载端位置施加位移全固定约束，设定加载步数和每步加载位移大小，即步长；如图10所示；
[0102]S202.对SRP试件顶部施加位移荷载，根据前述的破坏准则和本构关系，判断每步加载后的试件中是否有单元被破坏；
[0103]对有破坏单元产生的加载步，改变破坏单元的材料属性，重复S201-S202，直到没有新的破坏单元再出现后，转入下一荷载步；对于没有破坏单元产生的加载步则直接转入下一荷载步；
[0104]S203.输出应力应变数据；
[0105]S204.判断是否达到加载步数；如果达到加载步数，求解终止；如果没有达到加载步数，重复S201-S204。
[0106]通过以上过程，能够较好的预测和揭示SRP承受外部荷载作用时的弹性模量、强度变化和破坏过程，也证明了本数值模型的准确性和适用性，为SRP的工程应用提供了技术依据。
[0107]SRP的弹性模量和强度特性是该新型材料最重要的性质，也是实际应用中最关键的因素。SRP力学性能的细观数值模拟能够较好的预测和揭示其承受外部荷载作用时的弹性模量、强度变化和破坏过程，对SRP特性的探索有十分重要的意义。如前所述，在使用ANSYS软件对SRP进行细观数值模拟时，网格的尺寸大小能够影响计算的精确程度，也在很大程度上决定了计算的时间长短和计算规模，因此网格单元剖分对计算结果有着直接影响。
[0108]为了证明网格尺寸的敏感性，在本实施例中，分析了模拟基体高聚物密度为0.45g/cm3，尺寸为Φ150mm×300mm的圆柱体SRP试件的截面受压破坏，SRP试件尺寸和参数，如表1所示：
[0109]表1
[0110]
[0111]对此试件分别进行了四种网格尺寸的数值模拟计算，网格单元尺寸分别为1、2、3、4mm，得到单元剖分情况、计算耗时以及SRP力学性能的变化情况，结果如表2所示：
[0112]表2
[0113]
[0114]可见，随着网格尺寸增大，有限元模型单元数量减少，弹性模量和抗压强度有所上升，力学参数最大变化率均在5％以内。
[0115]图11为四种网格尺寸下试件最终破坏图，可以看出，网格尺寸对SRP轴心压缩裂缝发展没有太大的影响。裂缝主要集中在界面周围，由于中密度基体高聚物中的碎石不会发生破坏，因此，碎石对裂缝的发展有阻断的作用。综合考虑模拟结果的精确性和计算耗时，取网格单元尺寸为2mm最优。
[0116]进一步的，为了验证本数值模型的合理性，在本实施例中，分析了模拟基体高聚物密度为0.45g/cm3，尺寸为Φ150mm×300mm的圆柱体SRP试件的截面轴心受压破坏过程，SRP试件尺寸和参数采用表1所示，利用Matlab随机碎石程序，随机生成四组不同碎石分布的样本，计算其在相同荷载下的SRP轴心受压破坏过程，模拟结果如表3所示：
[0117]表3
[0118]
[0119]可以看出，四组模拟结果离散性较小。最大误差不超过3.29％，最大变异系数不超过1.48％，可以证明本模型具有较高的稳定性。
[0120]图12给出了四组SRP样本的最终破坏图，可以看出，不同碎石分布情况下，裂缝的位置不一定相同，且高聚物-碎石界面过渡区域容易出现裂缝。由于碎石的强度大，单个颗粒所占空间较大，裂缝的破坏位置与碎石的分布集中程度有一定关系：碎石集中的地方，裂缝的扩展区域比较大；该基体密度下，裂缝不会穿过碎石发展，因此中密度试件中碎石对裂缝的发展方向有阻断作用。
[0121]经过上述验证，本数值模型能够较好地对SRP试件的破坏进行数值仿真计算，具有较强的适用性和稳定性。
[0122]以上所述的实施例仅是对本发明优选方式进行的描述，并非对本发明的范围进行限定，在不脱离本发明设计精神的前提下，本领域普通技术人员对本发明的技术方案做出的各种变形和改进，均应落入本发明权利要求书确定的保护范围内。
