考虑强度非线性的堆石料弹塑性本构模型研究

堆石料亚塑性本构模型及面板堆石坝数值分析堆石料是面板堆石坝的受力支撑体系,其受力变形特性在很大程度上决定了大坝的工作性态。

本文在堆石料工程特性分析研究的基础上,对其进行全新的亚塑性本构建模,结合理论分析和数值计算,指出该模型能够反映堆石料的主要力学特性,且又有自身特色,一些特有的研究成果对堆石坝进一步的应力变形分析研究具有一定的指导意义。

具体研究了以下主要内容: (1) 分析了堆石料的压缩性,抗剪特性及其变形机理,概括了堆石料的非线性、压硬性、剪胀性和各向异性等主要力学特性,为亚塑性本构建模奠定了力学分析基础。

此外,对常规堆石料的本构模型进行了概括比较分析。

(2) 系统地介绍了亚塑性本构理论的起源、发展和目前国内外研究现状。

给出了亚塑性理论的数学力学基础理论,分析了亚塑性模型的一般表达形式及其基本力学特性。

(3) 分析研究了散粒体材料的孔隙比变化规律。

结合堆石料的力学特性分析,将Gudehus-Bauer亚塑性模型作一些改进后引入到堆石料本构建模中。

理论分析和计算都表明,该模型能够反映堆石料的主要力学特性。

此外,针对改进模型的9个材料常数中部分参数确定较困难的现象,建议了一种基于部分试验数据和反演分析基础上的参数确定的复合识别法,可方便地确定堆石料亚塑性模型的本构参数。

(4)对几种常规室内土工试验进行了亚塑性数值模拟,结合计算结果对亚塑性模型做了进一步的力学特性分析研究。

(5)开发研制了改进的堆石料Gudehus-Bauer亚塑性模型的有限元法计算程序,将其应用于堆石坝数值分析计算中,得到了一些有意义的新的研究成果。

(6)总结了面板堆石坝三维非线性有限元法计算的基本理论。

(7)对堆石料的流变性采用粘性-亚塑性理论进行了初步建模分析。

建立了基于亚塑性理论的面板与垫层之间的新型接触面模型。

弹性-塑性材料本构模型与模拟方法研究弹性-塑性材料的本构模型与模拟方法是材料力学研究领域的重要内容之一。

本文将介绍弹性-塑性材料的本构模型和模拟方法，并探讨其在工程实践中的应用。

弹性-塑性材料是一类具有弹性和塑性行为的材料，其在受力作用下可以发生弹性变形和塑性变形。

弹性变形是指材料在受力后能够恢复到原始形状的能力，而塑性变形是指材料在受力后无法完全恢复到原始形状的能力。

弹性-塑性材料的本构模型是描述材料力学行为的数学模型。

常见的弹性-塑性本构模型包括线性弹性模型、非线性弹性模型和塑性模型等。

线性弹性模型假设材料的应力和应变之间存在线性关系，适用于小应变情况。

非线性弹性模型考虑材料的应力-应变曲线是非线性的情况，适用于大应变情况。

塑性模型描述材料的塑性行为，常用的塑性模型有屈服准则、硬化规律和流动规律等。

在弹性-塑性材料的模拟方法中，有限元法是最常用的方法之一。

有限元法将材料划分为许多小的有限元单元，并在每个单元内建立本构模型，通过求解有限元方程组得到材料的应力和应变分布。

有限元法具有较高的精度和灵活性，适用于各种复杂的材料和结构。

除了有限元法，还有其他一些模拟方法可以用于弹性-塑性材料的研究。

例如，离散元法可以用于描述材料的离散微观结构和颗粒之间的相互作用，适用于颗粒材料的模拟。

分子动力学方法可以用于模拟材料的原子尺度行为，适用于纳米材料的研究。

这些方法各有特点，可以根据研究对象的不同选择合适的方法。

弹性-塑性材料的本构模型和模拟方法在工程实践中有广泛的应用。

例如，在材料设计和优化中，可以通过模拟方法预测材料在不同载荷下的应力和应变分布，为工程设计提供参考。

在材料加工和成形过程中，可以通过模拟方法优化工艺参数，提高产品的质量和效率。

在材料疲劳和断裂分析中，可以通过模拟方法评估材料的寿命和安全性能。

总之，弹性-塑性材料的本构模型与模拟方法是材料力学研究的重要内容。

通过建立适合材料行为的本构模型，结合合适的模拟方法，可以更好地理解和预测材料的力学行为，为工程实践提供支持。

第 55 卷第 1 期2024 年 1 月中南大学学报(自然科学版)Journal of Central South University (Science and Technology)V ol.55 No.1Jan. 2024高心墙堆石坝材料本构模型计算的适用性研究程瑞林1, 2，汪泾周1，范钦煜1，湛正刚2，周伟1, 3，马刚1, 3(1. 武汉大学 水资源工程与调度全国重点实验室，湖北 武汉，430072；2. 中国电建集团 贵阳勘测设计研究院有限公司，贵州 贵阳，550081；3. 武汉大学 水工程科学研究院，湖北 武汉，430072)摘要：为研究不同材料本构模型对高心墙堆石坝应力、变形计算的适用性，以RM 300 m 级高心墙堆石坝堆石料的室内试验成果为基础，分别对邓肯E-B 模型、沈珠江双屈服面模型及MPZG 模型进行材料本构验证和坝体有限元分析。

研究结果表明：3种模型均能较好地表现堆石料三轴路径下的应力响应，而邓肯E-B 模型在低围压下体积响应表现效果欠佳；3种模型计算的坝体应力、变形分布均符合一般规律，应力相近，但变形计算中邓肯E-B 模型得到的沉降与顺河向位移均较大；在高心墙堆石坝建设初期，可采用邓肯E-B 模型进行早期分析，后期宜采用沈珠江双屈服面模型、MPZG 模型等进行决策评价。

关键词：应力；变形；本构模型；心墙堆石坝；有限元分析；广义塑性模型中图分类号：TV311 文献标志码：A 开放科学(资源服务)标识码(OSID)文章编号：1672-7207（2024）01-0219-11Study on applicability of material constitutive model calculationfor high core wall rockfill damCHENG Ruilin 1, 2, WANG Jingzhou 1, FAN Qinyu 1, ZHAN Zhenggang 2, ZHOU Wei 1, 3, MA Gang 1, 3(1. State Key Laboratory of Water Resources Engineering and Management, Wuhan University,Wuhan 430072, China;2. Power China Corporation Guiyang Survey, Design and Research Institute Co. Ltd., Guiyang 550081, China;3. Institute of Water Engineering Science, Wuhan University, Wuhan 430072, China)Abstract: In order to study the applicability of different material constitutive models to the calculation of stress and deformation of high core rockfill dam, the material constitutive verification and the finite element analysis ofthe dam body of Duncan E-B model, Shen Zhujiang double yield surface model and MPZG model were carried收稿日期： 2023 −02 −20； 修回日期： 2023 −04 −16基金项目(Foundation item)：国家重点研发计划项目(2022YFC3005505)；国家自然科学基金资助项目(52179141)；大学生创新创业训练项目(202110486038) (Project(2022YFC3005505) supported by the National Key R&D Program of China; Project (52179141) supported by the National Natural Science Foundation of China; Project(202110486038) supported by Innovation and Entrepreneurship Training Program for College Students)通信作者：周伟，博士，教授，从事高坝结构数值仿真研究；E-mail ：**************.cnDOI: 10.11817/j.issn.1672-7207.2024.01.018引用格式： 程瑞林, 汪泾周, 范钦煜, 等.高心墙堆石坝材料本构模型计算的适用性研究[J].中南大学学报(自然科学版), 2024, 55(1): 219−229.Citation: CHENG Ruilin, WANG Jingzhou, FAN Qinyu, et al. Study on applicability of material constitutive model calculation for high core wall rockfill dam[J]. Journal of Central South University(Science and Technology), 2024, 55(1): 219−229.第 55 卷中南大学学报(自然科学版)out based on the laboratory test results of the rockfill material of RM 300 m grade high core rockfill dam. The results show that all the three models can better perform the stress response under the triaxial path of rockfill, while Duncan E-B model has poor performance in volume response under low confining pressure. The stress and deformation distribution of the dam body calculated by the three models are in accordance with the general law, and the stress is similar, but the settlement and along-river displacement obtained by Duncan E-B model are both larger. It is suggested that Duncan E-B model should be used for early analysis in the early stage of high core rockfill dam construction, and that Shen Zhujiang double yield surface model and MPZG model should be used for decision evaluation in the later stage.Key words: stress; deformation; constitutive model; core wall rockfill dam; finite element analysis; generalized plastic model近年来，随着堆石坝设计理论的发展与高坝填筑技术的革新，我国高堆石坝的建设水平逐渐从200 m级上升到300 m级。

加筋土静三轴试验和本构模型研究进展刘德顺;陈永杰【摘要】基于加筋土的加固机理,阐述了加筋土静力三轴试验和本构模型的研究进展,得出的结论如下:静三轴试验研究可获取加筋土复合体的静力本构模型参数,加筋层数、围压和填料性质是影响试验结果的主要因素,加筋后土体的强度得到较大提高;现有的加筋土本构模型难以反映加筋土的各向异性,而且参数难以获得,因此,目前尚未形成一种普遍适用的本构模型.【期刊名称】《黑龙江科技信息》【年(卷),期】2017(000)035【总页数】3页(P139-141)【关键词】静力三轴试验;加筋;本构模型【作者】刘德顺;陈永杰【作者单位】常州市轨道交通发展有限公司,江苏常州 213000;常州市轨道交通发展有限公司,江苏常州 213000【正文语种】中文加筋土工程是在土体中沿水平方向铺设多层加筋材料，使其与土体形成加筋土复合体。

筋材通过与土体之间的摩擦力将自身的抗拉强度与土体的抗压强度结合起来，从而提高了土体的强度，减小了土体的变形，提高了加筋土工程的稳定性。

法国工程师Herrividal在20世纪60年代初根据三轴实验结果提出了加筋土的设计理论，并于1965年在法国成功地建造了世界上第一座加筋土挡墙[1]，为加筋土技术在土木工程中的应用开创了先河。

1967-1968年法国至意大利高速公路建设上大量采用了加筋土结构，当时挡土墙的最大高度达到了23m。

20世纪70年代，美国工程师兵团采用了高强土工织物加固土堤地基，减少了软土地基上土堤的沉降[2]。

加筋土技术在我国的应用大约始于20世纪70年代末，在云南修建了一座高2m的储煤仓下挡土墙，这是我国修建的第一座加筋士结构物，从那以后经过20多年的探索和研究，加筋土技术在我国公路、铁路、水利、煤炭、城建等行业的工程建设中得到了迅速发展和推广应用，取得了良好的经济效益和环境效益。

同时，加筋土的加固研究在实验技术、设计理论、施工设备与工艺等方面也日臻成熟，取得了令人瞩目的成就。

材料塑性行为的本构模型与材料强度评估材料的塑性行为是指在外力作用下，材料会发生形变而不会恢复原状的性质。

研究材料的塑性行为对于工程设计和材料强度评估至关重要。

本文将探讨材料塑性行为的本构模型以及如何评估材料的强度。

一、材料塑性行为的本构模型材料塑性行为的本构模型是描述材料塑性变形的数学模型。

常见的本构模型有线性弹性模型、塑性流动模型和本构方程模型等。

1. 线性弹性模型线性弹性模型是最简单的本构模型之一。

它假设材料在小应变范围内具有线性的应力-应变关系。

这种模型适用于弹性变形较为明显的材料，如金属等。

线性弹性模型可以通过胡克定律来描述，即应力等于弹性模量乘以应变。

2. 塑性流动模型塑性流动模型是描述材料在塑性变形时的行为的模型。

它假设材料在塑性变形时具有一定的流动规律。

常见的塑性流动模型有屈服准则模型和流动规律模型。

屈服准则模型是通过确定材料的屈服准则来描述材料的塑性行为。

常见的屈服准则有屈服应力准则和屈服能准则等。

屈服应力准则是通过确定材料的屈服应力来描述材料的塑性行为，而屈服能准则是通过确定材料的屈服能来描述材料的塑性行为。

流动规律模型是通过确定材料的流动规律来描述材料的塑性行为。

常见的流动规律模型有希尔德布兰特模型和维斯科普模型等。

这些模型通过考虑应力、应变速率和温度等因素来描述材料的流动规律。

3. 本构方程模型本构方程模型是综合考虑材料的弹性和塑性行为的模型。

它通过建立材料的本构方程来描述材料的塑性行为。

本构方程模型可以通过试验数据拟合得到，也可以通过数学模型推导得到。

二、材料强度评估材料强度评估是指通过实验或计算方法来评估材料的强度性能。

常见的材料强度评估方法有拉伸试验、压缩试验和扭转试验等。

1. 拉伸试验拉伸试验是最常用的评估材料强度的方法之一。

它通过施加拉力来测量材料的抗拉强度、屈服强度和延伸率等指标。

拉伸试验可以得到材料的应力-应变曲线，从而评估材料的强度性能。

2. 压缩试验压缩试验是评估材料强度的另一种常用方法。

岩土材料非线性本构模型的开发与验证ZHANG Liangyi;CHEN Tielin;ZHANG Dingli;AN Fei【摘要】在增量塑性理论框架下,开发一种适用于岩土材料的非线性数值本构模型,利用开发的模型模拟三轴压缩试验、单轴拉伸试验及三轴拉压组合试验,并与室内试验结果和常用理论模型计算结果进行对比分析.结果表明:本文模型能够更有效地反映岩土体材料压缩、拉伸、压缩到反向拉伸、应变软化及循环加卸载状态,并通过单参数调节实现曲线形式由线性到双曲线的转变,弥补了现有模型的不足.【期刊名称】《北京交通大学学报》【年(卷),期】2018(042)006【总页数】7页(P9-14,40)【关键词】岩土力学;本构模型;非线性;加卸载;拉压组合;应变软化【作者】ZHANG Liangyi;CHEN Tielin;ZHANG Dingli;AN Fei【作者单位】;;;【正文语种】中文【中图分类】TU45岩土体作为一种气液固三相混合材料，其应力应变关系复杂而多变.为了反映岩土体应力应变的非线性特质，Duncan等[1]提出了经典的双曲线型非线性弹性模型，得到了广泛的应用.Schanz等[2]在Duncan-Chang模型的基础上提出一种非线性硬化模型，以弥补Duncan-Chang模型在反映土体剪胀性等方面的不足.此外，相关学者提出抛物线型[3]、指数型[4]、幂数型[5]等非线性弹性模型，以克服双曲线拟合的不足.上述模型均采用外凸型、单调递增且存在渐近线的单一函数形式，然而单一固定函数形式难以描述复杂多变的土体本构关系，因此相继有学者提出更加灵活多变的复合型函数关系[6-7].随着岩土体工程复杂性的增加，仅仅考虑压缩状态下的应力应变非线性特质往往不够.在复杂的隧道、基坑等工程中，时常存在卸载、再加载甚至反向拉伸的情况，单纯以线弹性模型来分析岩土体时常存在较大的差异.故而针对复杂岩土体工程开发出简洁实用的本构模型具有重要的实践意义.在运用各种非线性弹性模型时，限制其曲线形式的往往是其渐近线，为了更好地贴合峰值强度，而导致整体曲线拟合效果变差.复合型函数虽然能有效适应各种多变的本构关系，但其复杂性及相关参数的确定难度大大增加，同时也难以适应峰后各式各样的应变软化曲线.上述不足，可以通过引入经典的弹塑性理论来克服.本文作者基于增量塑性理论框架，引入双参数变切线模量概念，编程开发一种改进的数值计算模型，通过与室内试验结果和现有理论模型进行对比，验证本文模型的正确性与适用性，以期获得能够反映土体从压缩至拉伸全过程、循环加卸载及应变软化特性的数值计算模型，非线性应力应变关系示意图如图1所示.图1 非线性应力应变关系示意图Fig.1 Diagram of nonlinear stress-strain relation1 理论描述为了便于理解和程序化，理论以模块化逻辑顺序进行描述，采用土力学符号约定，以压为正，以拉为负，其中主应力大小依次为σ1>σ2>σ3.1.1 压缩/拉伸曲线当前拉伸应力应变形态可归结为3类：应变软化型、应变硬化型和应变强化型[8].试验结果[8-10]表明拉伸曲线和常规三轴压缩曲线有着相似的形态，均具有典型的非线性特性.当土体处于压缩状态时，主应力大小均为正值，主应力最小值σ3≥0.根据Mohr-Coulomb准则，屈服函数fs可表示为(1)式中：c为黏聚力；φ为摩擦角；Nφ=(1+sin φ)/(1-sin φ).根据式(1)强度准则，引入应力水平相关变切线模量理念，可得压缩切线模量计算公式为(2)式中：Rc为形状系数，决定破坏应力下双曲线的形式，当Rc=0时，模型将退化为线弹性；Ei为初始模量，也是最大的切线模量.当土体处于拉伸状态时，其主应力最小值σ3<0.引入Rankine准则，屈服函数ft 表示为ft=-σt-σ3(3)式中：σt为抗拉强度，其值为正.计算中如果直接采用Mohr-Coulomb抗拉强度必然会高估土体的抗拉能力[11]，秉着简单易用原则，本文采用抗拉强度截断方式，仅对抗拉强度进行上限限制(4)参照式(2)，可获得拉伸切线模量的表达式(5)式中：Eit为拉伸初始模量；Rt为拉伸形状系数.为了简化，当可选参数Eit和Rt空缺的时候，默认令其与受压状态相同，即为Ei和Rc.1.2 加/卸载曲线相关研究表明[10, 12-13]，不论三轴拉伸试验还是单轴拉伸试验，其卸载再加载曲线与三轴压缩时的形式一致，均表现为滞回圈形式，这也说明拉伸情况下的卸载再加载可以与压缩情况下采用相同的形式.因此，考虑历史应力状态的影响，由应力水平s区分加/卸载情况.卸载情况下，切线模量表达式为Et=[1-Ru(sH-s)]2Eu(6)式中：Eu为卸载初始模量；sH为加载历史最大应力水平；Ru为形状系数，用于控制加卸载曲线的形式.在加载情况下，切线模量的计算方法为Et=[1-Rr(s-sL)]2Er(7)式中：sL为卸载历史最低应力水平，并设定再加载初始模量Er与形状系数Rr为可选参数，默认情况下与卸载情况相同.从式(6)和式(7)可以发现，当Ru=Rr=0时，加/卸载切线模量退化为线性，说明线性假设是本方法的一种特例.1.3 应力应变关系将式(2)、式(5)、式(6)和式(7)代入胡克增量定律即可获得应力增量Δσij与应变增量Δεij的关系(8)式中：K为体积模量；剪切模量G=3KEt/(9K-Et)；下标i,j,k=1,2,3.由于切线泊松比的范围限制0<υt<0.5，故需要对体积模量K进行限制，0.33Et≤K≤17Et.1.4 塑性修正根据塑性力学理论，引入非关联流动法则对达到剪切破坏的应力点进行修正，并使用关联流动法则修正张拉破坏的应力点.由Mohr-Coulomb屈服准则，确定剪切塑性势函数gs为gs=σ1-σ3Nψ(9)式中：ψ为剪胀角；Nψ=(1+sin ψ)/(1-sin ψ).当ψ=φ时，gs转变为关联流动法则.对于非关联流动法则，一般情况下ψ<φ，通常近似取ψ=0.剪胀角ψ的引入可以一定程度上反映土体的剪胀性.由Rankine准则，张拉塑性势函数gt为gt=-σ3(10)为了便于区分破坏类型，定义函数H(σ3,σ1)=0，表示fs=0与ft=0所围区域的外斜线，其示意图如图2所示，表达为H=σ3+σt+A(σ1-B)(11)式中：图2 破坏分区示意图Fig.2 Diagram of damage division当弹性试算应力点(σ3,σ1-σ3)位于区域1时，即fs>0且H>0，说明发生剪切破坏，由式(9)确定流动法则，将应力修正回fs=0的曲线上.当试算应力点位于区域2时，即ft>0且H<0，此时发生张拉破坏，应由式(10)确定流动法则，将应力修正回ft=0的曲线上.1.5 应变软化图3为应力应变关系示意图，可知总应变可表示为其中εe为弹性应变，与之和为塑性应变.发生于破坏点之前，可由非线性应力应变关系获得，而发生在破坏点之后，依据相应流动法则确定.图3 应力应变关系示意图Fig.3 Diagram of stress-strain relation在数学上，线性函数是最容易处理的，通过线性插值，使得计算复杂非线性函数的值更加简单、高效.对于前人提出的各式各样复杂的应变软化曲线[10, 14-15]，可将其简化为多段直线的形式进行近似，塑性应变与软化参数关系示意图如图4所示.图4 塑性应变与软化参数关系示意图Fig.4 Relation of softening parameter and plastic strain引入损伤变量来描述材料软化的思路[14-15]，由于概念明确、使用简单，并且适用于连续介质数值计算程序的使用.故将黏聚力内摩擦角剪胀角抗拉强度定义为相关塑性应变的分段线性函数.其中，表示塑性剪应变，表示塑性拉应变.以塑性偏应变增量第二不变量J2′定义塑性剪应变增量(12)式中：Δeij为偏应变增量，i,j=1,2,3.以最大塑性拉应变增量定义塑性拉应变增量(13)2 模型开发模型基于有限差分算法，由给定时刻t的应力状态与时步Δt总应变增量确定时刻t+Δt的应力状态.基于模块化思路，采用C++语言编写相应计算程序，程序计算流程见图5.图5 程序流程图Fig.5 Flow chart of program3 模型验证下述模型材料参数依照Duncan-Chang模型参数取法[1,16]选取，其中形状系数可按照破坏比Rf取值或进行适当微调以调节曲线形式.而应变软化部分根据1.5节的定义，依据室内试验结果选取关键节点进行拟合.3.1 常规三轴压缩试验针对某工程土坝含细粒土砂试样进行固结围压为50 kPa、100 kPa、150 kPa的常规三轴固结排水试验，试验结果[17]如图6所示.可以发现，当应力达到峰值后，有明显的应变软化现象，此时应当添加应变软化参数，由于砂土黏聚力为0，故仅需给出摩擦角φ随等效塑性剪应变关系，当为0、0.03、0.091、1时,φ为41.5°、41.24°、38.41°、38.41°.Ei取20.35 MPa,K取12.97 MPa,Rc取0.55,φ取41.5°,c取0.根据试验的流程，数值模拟采用位移加载的方式.首先在选定围压下达到初始应力平衡，保持围压不变，固定底部，顶部施加1.5×10-5m/时步的速度，计算10 000时步，记录相应的应力与应变.图6 常规三轴压缩试验结果图Fig.6 Results of triaxial compression test为了进一步验证程序计算的正确性和适用性，将本文模型与Duncan-Chang模型、Mohr-Coulomb模型进行对比.其中，Duncan-Chang模型(邓肯模型)参数取值为：k=147，n=0.8,kb=117,m=0.23,Rf=0.55,φ=41.5°,c=0.Mohr-Coulomb模型(摩尔模型)选择初始模量Ei=20.35 MPa作为弹性模量，其余所需参数与上文相同.选取围压150 kPa的试验数据进行对比，数值计算结果与室内试验结果对比见图7.图7 数值与试验的结果对比Fig.7 Comparison results of simulations andexperiment可以发现，本文模型与试验数据在应力应变全过程贴合度均很高.而邓肯模型在破坏应力后与试验曲线存在较大误差，这也说明邓肯模型较适合于荷载小于破坏应力条件下的情况，而对于具有明显峰值特征的应力应变关系具有一定的局限性.对于峰前非线性应力应变段，摩尔模型的拟合效果不佳，存在较大偏差.3.2 单轴拉伸试验张丙印等[13]研发了一种长方体试样的单轴拉伸试验，由于具有统一的横截面，其结果与单元测试应力应变状态相似，非常适合于本构模型的建立以及数值程序的验证.故本文选用直拉试验结果数据来验证数值程序拉伸部分及循环加/卸载部分的正确性.该试验材料选择糯扎渡高塑性黏土，通过位移控制式加载装置进行加载，当偏应力达到峰值应力70%时，通过施加反向位移以实现加/卸载的循环荷载效果.循环荷载下糯扎渡黏土拉应力应变曲线对比见图8.图8 循环荷载下糯扎渡黏土拉应力应变曲线对比Fig.8 Comparison of tensile stress-stain curves under cycle loading for Nuozadu clay由于是拉伸试验，并不发生剪切破坏，故黏聚力和摩擦角参数设置可以忽略.抗拉强度选择试验结果峰值70.45 kPa，卸载再加载初始模量取6个卸载再加载初始切线模量的平均值，单轴拉伸试验数值模型相关参数为Ei=158.82 MPa,K=52.3 MPa,Eu=300.07 MPa,Rt=0.77,Ru=0.82.当取0、1.5×10-4、2.8×10-4、6.5×10-4、10×10-4、11×10-4时，σt取70.45 kPa、69 kPa、67 kPa、21 kPa、3 kPa、0 kPa.依照试验流程，由于是无侧限单轴拉伸试验，无需施加围压.固定底部，在顶部施加正向与反向的速度来模拟试验中的加/卸载过程，循环加/卸载过程详见表1，其中正为拉，负为压.从图8可以发现，本文模型从起始到软化与试验数据均具有较高的吻合度，说明本文模型可以用于分析土体拉伸状态的问题.卸载和再加载曲线与试验数据吻合度基本一致，说明本文提出的近似方法可行，能够反映出土体循环加/卸载的特性.由于压缩情况与拉伸的一致，在此不再验证压缩状态下的循环加/卸载的情况.表1 循环加卸载流程Tab.1 Step of cyclic loading-unloading加载速率×10-7/(m/时步)时步加载速率×10-7/(m/时步)时步 24 000 24 600-22 000-22 10023 000 28 780-22 5003.3 三轴压缩-拉伸组合试验张琰等[10]应用自制新型卧式三轴仪，对糯扎渡心墙土料进行三轴压缩-拉伸试验.在室压为150 kPa时，分别在压缩应力水平s为22.5%、72.4%和94.6%时开始反向加载，直至试样破坏.依照试验流程，采用位移加载方式，固定底部，施加围压150 kPa达到初始应力平衡，然后在顶部施加1.5×10-5m/时步的速度，计算相应步数后，施加反向速度，模拟反向加载.三轴压缩-拉伸组合试验数值模型参数详见表2.由于拉伸软化部分并不是本小节的验证内容，故忽略软化部分.各组试验结果与数值模拟结果如图9所示.可以发现，本文模型可以较好地拟合三轴压缩-拉伸组合试验的结果.表2 三轴压缩-拉伸组合试验数值模型参数Tab.2 Numerical model parameters of triaxial compression-tension test参数数值参数数值Ei/MPa129.8Rc0.9Eit/MPa35Rt0.5Eu/MPa163Ru0.55K/MPa45σt/kPa60c/kPa 228.35φ/(°)27图9 三轴压缩-拉伸组合试验结果对比Fig.9 Comparison of triaxial compression-tension tests3.4 形状系数对应力应变曲线的影响为了分析形状系数对应力应变曲线的影响，采用单一变量法，其他所需材料参数见表2，不考虑应变软化，形状系数Rc由0至1.0变化.图10给出了不同形状系数情况下的应力应变曲线，可以发现，当Rc=0时，曲线形式与图7中Mohr-Coulomb曲线完全一致，说明Mohr-Coulomb模型是本文模型的一种特殊形式.随着Rc的增加，曲线形式从理想弹塑性逐渐向双曲线转变.相比于传统双曲线模型，塑性理论的引入，释放了破坏比对于求取峰值渐近线的约束，使得在相同渐近值的前提下，应力应变曲线呈现从曲线Rc=0到曲线Rc=1的多样性变化，更加符合岩土材料性质复杂多样性的特点，大大增加了数值程序的适用性.图10 不同形状系数情况下应力应变曲线Fig.10 Stress-strain curves under different valus of shape factor4 结论基于有限差分算法，编写开发了一种非线性本构模型.通过常规室内三轴试验、单轴循环拉伸试验及三轴压缩-拉伸组合试验，验证了所开发模型的正确性与适用性.对比结果表明，本文模型能够有效地反映土体压缩、拉伸状态下的非线性应力应变特性，能够反映循环加/卸载状况下及压缩到反向拉伸情况下的非线性特性，并且通过劣化强度参数有效地反映土体应变软化的特性.本文模型建立在变切线模量与增量塑性理论基础之上，材料参数个数少且均能从常规室内试验中获取，保持了本构模型简单易用的特性.同时，采用非关联流动法则及增量理论，在一定程度上可以反映岩土体剪胀性及应力路径的影响.另外，相较于常用的理论模型，可通过改变形状系数实现曲线形式由线性到双曲线型的转变，满足了岩土材料应力应变关系多样性的要求，具有更广泛的适用性.参考文献(References)：【相关文献】[1] DUNCAN J M. Nonlinear analysis of stress and strain in soils[J]. Journal of the Soil Mechanics and Foundation Division, 1970, 96(5): 1629-1653.[2] SCHANZ T, VERMEER P A, BONNIER P G. The hardening soil model: formulation and verification[J]. Beyond 2000 in Computational Geotechnics, 1999: 281-296.[3] 王伟，卢廷浩，周干武. 黏土非线性模型的改进切线模量[J]. 岩土工程学报, 2007, 29(3): 458-462.WANG Wei, LU Tinghao, ZHOU Ganwu. Improved tangent modulus of nonlinear soil model[J]. Chinese Journal of Geotechnical Engineering, 2007, 29(3):458-462.(in Chinese) [4] GITAU A N, GUMBE L O, BIAMAH E K. Influence of soil water on stress-strain behaviour of a compacting soil in semi-arid Kenya[J]. Soil & Tillage Research, 2006, 89(2): 144-154. [5] 殷德顺，王保田. 负乘幂本构模型的切线模量[J]. 岩土力学, 2009, 30(7): 2168-2172.YIN Deshun, WANG Baotian. Tangent modulus of negative-power constitutive model[J]. Rock and Soil Mechanics, 2009, 30(7):2168-2173.(in Chinese)[6] 王丽琴，鹿忠刚，邵生俊. 岩土体复合幂-指数非线性模型[J]. 岩石力学与工程学报, 2017, 36(5): 1269-1278.WANG Liqin, LU Zhonggang, SHAO Shengjun. A composite power exponential nonlinear model of rock and soil[J]. Chinese Journal of Geotechnical Engineering, 2017, 36(5):1269-1278.(in Chinese)[7] 王伟，宋新江，凌华，等. 滨海相软土应力-应变曲线复合指数-双曲线模型[J]. 岩土工程学报, 2010, 32(9): 1455-1459.WANG Wei, SONG Xinjiang, LING Hua, et al. Composite exponential-hyperbolic model for stress-strain curve of seashore soft soil[J]. Chinese Journal of Geotechnical Engineering, 2010, 32(9):1455-1459.(in Chinese)[8] 张云，王惠敏，鄢丽芬. 击实黏土单轴拉伸特性试验研究[J]. 岩土力学, 2013, 34(8): 2151-2157. ZHANG Yun, WANG Huimin, YAN Lifen. Test research on tensile properties of compacted clay[J]. Rock and Soil Mechanics, 2013, 34(8): 2151-2157. (in Chinese)[9] STIRLING R A, HUGHES P, DAVIE C T, et al. Tensile behaviour of unsaturated compacted clay soils-a direct assessment method[J]. Applied Clay Science, 2015, 112: 123-133.[10] 张琰，张丙印，李广信，等. 压实黏土拉压组合三轴试验和扩展邓肯张模型[J]. 岩土工程学报, 2010, 32(7): 999-1004.ZHANG Yan, ZHANG Bingyin, LI Guangxin, et al. Combined tension-compression triaxial tests and extended Duncan-Chang model of compacted clay[J]. Chinese Journal of Geotechnical Engineering, 2010, 32(7):999-1004. (in Chinese)[11] 宋焱勋，李荣建，刘军定，等. 结构性黄土的双曲线强度公式及其破坏应力修正[J]. 岩土力学, 2014, 35(6): 1534-1540.SONG Yanxun, LI Rongjian, LIU Junding, et al. A hyperbola strength formula of structuralloess and its modifying algorithm under failure stress state[J]. Rock and Soil Mechanics, 2014, 35(6):1534-1540. (in Chinese)[12] SUBRAMANIAM P, BANERJEE S. A correction to damping ratio for hyperbolichysteretic model for clayey soil[J]. International Journal of Geotechnical Engineering, 2013, 7(2): 124-129.[13] ZHANG B Y, LI Q M, YUAN H N, et al. Tensile fracture characteristics of compacted soils under uniaxial tension[J]. Journal of Materials in Civil Engineering, 2014, 27(10): 1-11.[14] 沈珠江. 结构性粘土的弹塑性损伤模型[J]. 岩土工程学报, 1993, 15(3): 21-28.SHEN Zhujiang. Elasto-plastic damage constitutive model for structured clay[J]. Chinese Journal of Geotechnical Engineering,1993, 15(3):21-28. (in Chinese)[15] 沈珠江. 结构性粘土的非线性损伤力学模型[J]. 水利水运工程学报, 1993(3): 247-255. SHEN Zhujiang. A nonlinear damage model for structured clay[J]. Hydro-Science and Engineering, 1993(3):247-255. (in Chinese)[16] 冯卫星，常绍东，胡万毅. 北京细砂土邓肯-张模型参数试验研究[J]. 岩石力学与工程学报, 1999, 18(3):327-330.FENG Weixing, CHANG Shaodong, HU Wanyi. Experimental study on parameters of Duncan-Chang model for Beijing fine sandy soil[J]. Chinese Journal of Rock Mechanics and Engineering, 1999, 18(3): 327-330. (in Chinese)[17] 唐洪祥，张丽娟. 平面应变试验与常规三轴试验条件下含细粒土砂力学特性比较研究[J]. 水利与建筑工程学报, 2015, 13(5): 223-227.TANG Hongxiang, ZHANG Lijuan. Comparative studies on mechanical characteristics of fine-grained sand with plane strain test and triaxial test[J]. Journal of Water Resources and Architectural Engineering, 2015,13(5):223-227. (in Chinese)。

堆石料Wang-Wu亚塑性本构模型参数求取的改进陈泽钦;刘国明;林伯江【摘要】研究Wang-Wu提出的亚塑性模型,并指出该模型应用在堆石料时存在的不足.推导了模型参数的求取方法,建立了参数与围压之间的关系,提出了新的模型参数.通过堆石料大型侧限压缩试验、大型三轴排水试验、邓肯E-B模型和改进的沈珠江双屈服面模型,反演其他围压的试验值并进行对比验证.所得成果表明,提出的7个新参数有效地解决该模型应用到堆石体材料时遇到的确定参数的难题,有利于将该模型推广到堆石体材料工程.%The Wang-Wu hypoplastic model was studied and its defect of the model used in rockfill materials was pointed out. According to the derivation of the model parameters equation, the relation-ship between the parameters and the confining pressure was established and the new model parameters were proposed. The method was verified by means of rockfill oedometric test, rockfill triaxial consoli-dated drained test, Duncan E-B model inversion test and modified Shen Zhujiang' s double-yield sur-face model inversion test. It was shown that the proposed seven new parameters could effectively solve the problem of determining the parameters when the model was introduced into the rockfill materials, and it is helpful for the model to be extended to the rockfill material engineerings.【期刊名称】《福州大学学报（自然科学版）》【年(卷),期】2017(045)002【总页数】7页(P199-205)【关键词】Wang-Wu模型;亚塑性;参数改进;堆石料【作者】陈泽钦;刘国明;林伯江【作者单位】福州大学土木工程学院,福建福州 350116;福州大学土木工程学院,福建福州 350116;福州大学土木工程学院,福建福州 350116【正文语种】中文【中图分类】TV32亚塑性本构模型最早由Kolymbas[1]借助张量函数提出，近十几年来在国内得到发展的一种新型本构理论.该理论以非线性连续介质力学为理论基础，以张量为计算工具，直接建立材料的应力率与应变率之间的关系.亚塑性理论与传统的弹塑性本构理论不同，它隐含了弹塑性理论中的屈服面、加工硬化定律、流动法则等基本概念，减少了大量的人为假设，具有极大的发展前景.目前已提出较多适用于无粘性土的本构方程，可用于堆石料应力变形特性分析.最早的亚塑性本构模型相当复杂，验证困难，适用性差.随着Wu等[2]在前人研究的基础上，针对无粘性土提出一个改进模型(简称Wu-Bauer亚塑性模型)后，该模型基本能够反映无粘性土的应力变形特性，随后亚塑性本构模型得到不断的发展.文献[3-4]在亚塑性理论中引入“向压因子”和“密度因子”形成了八参数亚塑性本构模型(简称Gudehus-Bauer亚塑性本构模型).该模型能够考虑无粘性土的物理状态量孔隙比的变化.Bauer等[5]用非线性张量值函数和相应的亚塑性模型分析了无粘性土的横观各向同性性质，熊保林等[6]在Wolffersdorff亚塑性模型中引入了考虑横观各向同性的矩阵，使得该模型能够考虑砂土的横观各向同性效应.张嘎等[7]在Wu-Bauer亚塑性模型引入临界状态参数，建立了粗粒土亚塑性模型.邵龙潭等[8]针对Wu-Bauer亚塑性模型参数确定问题，提出的单形调优法能够方便、准确地求取参数.文献[9-11]在亚塑性模型中考虑了颗粒破碎特征，使得改进的模型能够较好模拟颗粒材料力学特性.Wang等[12]针对Wu-Bauer亚塑性模型在三轴压缩模拟时具有过分的剪缩性缺点，提出新的无粘性土亚塑性本构方程，改善了这个缺点，并应用于隧洞开挖模拟.本研究在文献[12]提出的本构模型方程(简称为Wang-Wu亚塑性模型)基础上，针对堆石料求取参数时存在的问题，作了进一步改进，提出了新的模型参数.利用南京水利科学院福建仙游抽水蓄能电站面板堆石坝坝体填筑料静力特性试验[13]和邓肯的E-B模型[14]及改进的沈珠江模型[15]反演试验等试验方法验证其可行性.为了与岩土工程正负号一致，规定模型方程中以压应力、压应变及对应的率型为正，反之为负.文献[12]在亚塑性模型[2]基础上提出了新的亚塑性模型，如下式所示，并引入了项，该项在土体破坏状态时显然是不存在的，符合规律.式中：为Jaumann应力率张量；σ为应力张量；ci(i=1～4)为本构模型参数，无量纲；σ*为偏应力张量，σ*=σ-1/3(tr σ)；‖·‖为张量求Euclidean模；为应变率张量；tr为张量求迹.考虑常规三轴压缩试验边界条件和应变应力的对称性(见图1)，式(1)可以简化为式(2)和式(3).从式(1)可知，求取4个模型参数ci(i=1,2,3,4)，相应需要四个方程组成的求解关于ci(i=1,2,3,4)的线性方程组.为此，需从图1中选取合适的A、B两点(本研究建议A点尽量取在试验曲线中初始稳定阶段，B点尽量取值在试验曲线中段稳定阶段)，读取A、B两点的切线模量Et、倾角β、主应力差(σ1-σ3)，并引入应力比R=σ1/σ3、切线模量Et=d(σ1-σ3)/dε1、倾角).将相关的值代入式(2)和式(3)中即可得到线性方程组式(4)，并能求得模型的4个参数，求得的结果如表1所示.根据表1的模型参数进行模拟计算，如图2所示.从图2可知，Wang-Wu亚塑性模型可用于描述堆石料力学基本特性，如非线性性和剪胀性.但是从表1的模型参数可知：相同的堆石料下，不同围压确定的模型参数相差较大，模型参数对围压较为敏感.这样将该模型应用于堆石体材料,工程数值分析时就会遇到对同样的堆石料如何选用4个模型参数的问题.针对模型参数对围压较为敏感这一问题，对Wang-Wu亚塑性模型线性方程组(式(4))进行详细的推导，将A点取在常规三轴试验的初始状态即σ1=σ3，EtA=Ei，应力比RA=1，代入式(4)中.将前两式相减即可得到模型参数c1.式中：Ei为初始切线模量；为初始泊松比.根据文献[14]可知，三轴试验中初始切线模量Ei=Kpa(σ3/pa)n，其中pa=101.4 kPa，K和n为试验常数.vi是较为稳定的材料常数，假定vi与围压无关，并将Ei 代入式(5)中可得:通过式(6)建立了c1和围压σ3的关系，同时引进两个新参数K1和n1.将求得的c1代入式(4)可得其余参数c2、c3、c4与c1之间的关系,即ci=f(c1,RB,vi,βB)，i=2,3,4.为了简化c2、c3、c4与c1之间的关系，假定其余3个参数也具有c1与围压关系的表达式，详见下式.此时模型的4个模型参数ci(i=1，2，3，4)与围压建立了关系，并由4个模型参数扩展到8个模型参数,即Ki、ni，i=1，2，3，4.但是从表1参数c1、c2、c3、c4的结果可知，围压0.3 MPa堆石料模型参数c2的变化规律与其余3个参数不一致.从表1求出的模型参数c2变化规律看,随着围压的增大，模型参数c2趋于稳定.对此，本研究对c2处理的方法是取围压1.2 MPa求得的模型参数c2作为堆石料的最终模型参数，将c2视为常数.由此该模型由原来的4个参数变为7个参数. 分别对式(7)左右两边取对数，并作ln ci～ln(σ3/pa)，i=1，3，4关系图，其中ln Ki为截距，ni为斜率，如图3所示(限于篇幅仅展示垫层料和主堆料的结果).由此可得坝体填筑料的7个Wang-Wu模型新参数，见表2.由式(7)可知，根据表2给出的Wang-Wu模型7个新参数，即可求出不同围压下对应的Wang-Wu模型原来的4个参数c1～c4.与之前的模型参数相比，c1～c4不是常数，而是随围压改变而改变，解决了Wang-Wu模型在堆石体材料工程数值分析时参数确定的难题.为了验证求出的7个模型新参数的正确性，根据表2给出的模型参数，并用Wang-Wu亚塑性模型式(2)和式(3)进行常规三轴模拟计算，并与福建仙游抽水蓄能电站面板堆石坝坝体填筑料大型三轴试验成果对比，见图4(限于篇幅仅展示垫层料和主堆料的结果).同时根据试验成果求取邓肯E-B模型[14]和改进的沈珠江双屈服面模型[15]的模型参数，见表3，并将两个模型的模拟计算值与试验数据对比，为下一步反演验证数据做准备.由图4可知，Wang-Wu亚塑性模型模拟的常规三轴试验堆石料主应力差与轴向应变曲线和试验值均吻合程度很好，堆石料体积应变与轴向应变曲线和试验值总体吻合程度较好，能够反映堆石体材料的剪胀剪缩性.同时应当指出体变曲线达不到应力曲线的吻合程度，主要原因在于计算参数取值时与应力曲线有关的值有EtA、EtB、(σ1-σ3)A、(σ1-σ3)B，这四个值分别控制A、B两点的应力值和斜率，所以主应力差曲线吻合程度较好.只有倾角βA、βB与控制体变曲线相关，仅控制A、B 两点的斜率，所以体变曲线的吻合程度没有主应力差曲线好.为了进一步验证提出的7个模型新参数的正确性.根据表3给出的邓肯E-B模型和改进的沈珠江双屈服面模型相关参数，用邓肯E-B模型[14]反演不同围压的主应力差与轴向应变曲线，用罗刚等[15]改进的沈珠江双屈服面模型反演不同围压的体积应变与轴向应变曲线.同时也对Wang-Wu模型进行相应围压的模拟计算，并与反演值进行对比，如图5所示(限于篇幅仅展示垫层料和主堆料的结果).由图4和图5可知，邓肯E-B模型和改进的沈珠江双屈服面模型模拟值与试验数据的吻合程度较好，反演值也是合理的，具有可比性.Wang-Wu亚塑性模型模拟值与邓肯E-B模型和改进的沈珠江双屈服面模型反演值吻合程度总体上较好，从而较为充分地验证了模型新参数能够合理地反映出围压对原模型参数c1～c4的影响，说明提出的7个模型新参数是可行的.根据表2给出的Wang-Wu模型新参数对仙游抽水蓄能电站上库面板堆石坝坝体填筑料进行侧限压缩模拟计算，并将模拟值与试验值进行比较，见图6.从图6可知，Wang-Wu亚塑性本构模型模拟值与侧限压缩试验值吻合程度总体上较好，亚塑性本构模型模拟值总体上能反映堆石料孔隙比随压力增大而减小的趋势.从而进一步验证了7个模型新参数能够刻画出围压对原模型参数c1～c4的影响，也能够反映堆石料的力学性能，说明所提出的7个模型新参数是可行的.针对Wang-Wu亚塑性模型应用在堆石料时存在的模型参数问题，提出了新的模型参数.新的模型参数有效地解决了原模型参数c1～c4对围压较为敏感这一问题，使得新参数的Wang-Wu亚塑性模型能够描述堆石料力学性能如非线性性、剪胀剪缩性等，有效地解决了该模型应用于堆石料参数确定的难题，并有利于将该模型推广到堆石体材料工程.。

岩土类材料的弹塑性力学模型及本构方程摘要：本文主要结合岩土类材料的特性，开展研究其在受力变形过程中的弹性及塑性变形的特点，描述简化的力学模型特征及对应的适用条件，同时在分析研究其弹塑性力学模型的基础上，探究了关于岩土类介质材料的各种本构模型，如M-C、D-P、Cam、D-C、L—D及节理材料模型等,分析对应使用条件，特点及公式,从而推广到不同的材料本构模型的研究，为弹塑性理论更好的延伸发展做一定的参考性。

关键词：岩土类材料，弹塑性力学模型，本构方程不同的固体材料，力学性质各不相同.即便是同一种固体材料，在不同的物理环境和受力状态中，所测得的反映其力学性质的应力应变曲线也各不相同。

尽管材料力学性质复杂多变，但仍是有规律可循的，也就是说可将各种反映材料力学性质的应力应变曲线，进行分析归类并加以总结，从而提出相应的变形体力学模型.第一章岩土类材料地质工程或采掘工程中的岩土、煤炭、土壤,结构工程中的混凝土、石料，以及工业陶瓷等，将这些材料统称为岩土材料。

岩土塑性力学与传统塑性力学的区别在于岩土类材料和金属材料具有不同的力学特性。

岩土类材料是颗粒组成的多相体，而金属材料是人工形成的晶体材料。

正是由于不同的材料特性决定了岩土类材料和金属材料的不同性质.归纳起来，岩土材料有3点基本特性：1.摩擦特性.2。

多相特性。

3.双强度特性。

另外岩土还有其特殊的力学性质:1。

岩土的压硬性，2.岩土材料的等压屈服特性与剪胀性,3.岩土材料的硬化与软化特性。

4.土体的塑性变形依赖于应力路径。

对于岩土类等固体材料往往在受力变形的过程中，产生的弹性及塑性变形具备相应的特点，物体本身的结构以及所加外力的荷载、环境和温度等因素作用,常使得固体物体在变形过程中具备如下的特点。

固体材料弹性变形具有以下特点：(1)弹性变形是可逆的.物体在变形过程中，外力所做的功以能量（应变能）的形式贮存在物体内，当卸载时，弹性应变能将全部释放出来，物体的变形得以完全恢复；（2)无论材料是处于单向应力状态，还是复杂应力状态，在线弹性变形阶段,应力和应变成线性比例关系；(3)对材料加载或卸载，其应力应变曲线路径相同。

堆石料弹塑性循环本构模型研究及应用过大的永久变形往往是导致堆石坝震害的主要原因。

准确预测堆石体的地震变形在堆石坝的抗震设计中尤为重要。

目前广泛使用的等价粘弹塑性模型对堆石体在地震荷载作用下变形的描述具有很大经验性,只能粗略估算地震荷载引起的残余变形,采用循环弹塑性模型可弥补这一缺陷。

本文在前人研究的基础上,着重对堆石料在循环荷载作用下的应力应变规律进行了总结,建立了一个实用化的三维弹塑性循环本构模型,提出了基于该模型的数值化方法,开发了基于开源动力有限元软件OpenSees的面板堆石坝动力计算程序。

本文取得的研究成果主要包括以下几个方面:第一,建立了一个可合理地反映堆石料在循环荷载作用下应力-应变特性的三维弹塑性本构模型。

与目前已有的弹塑性循环本构模型相比,该模型具有以下几个明显优势:1)提出了一个基于修正空间滑动面强度理论的双参数三维非线性强度准则,使复杂循环加载下堆石料强度的描述更为简便合理;2)可反映偏应力循环、球应力循环及二者共同作用下不可逆体应变的累积规律,使循环荷载作用下堆石料体应变的描述更为准确合理;3)可反映循环剪切荷载作用下剪应力-剪应变关系的滞回特性、硬化特性和应变累积特性,使循环荷载作用下堆石料剪应变的描述更为精细合理;4)共有11个物理意义明确的模型参数,均可由振动三轴试验测得的最大剪切模量、模量衰减曲线和残余变形累积曲线等常规试验结果来确定,使该模型更加实用化。

第二,提出了该本构模型的数值实现方法,开发了基于开源动力有限元软件OpenSees的用户材料子程序以及针对面板堆石坝动力计算的辅助性前后处理程序,实现了较为高效的数值建模。

通过对二维理想坝模型、三维均质坝模型的系列计算和对比分析,完成了对程序有效性的测试。

第三,通过三维动力弹塑性有限元分析,较为合理地再现和揭示了紫坪铺高面板堆石坝的震动响应和变形规律,模拟结果与汶川地震后实测残余变形的分布基本上吻合,进一步检验了本文所建立模型及相应数值方法的适用性。