邓肯-张模型参数求取

土体邓肯—张非线性弹性模型参数反演分析邓肯张非线性弹性模型是一种用于描述弹性材料反应的普遍通用的模型，它适用于许多矿物组成的各种材料，如金属、合金、复合材料、多孔介质、矿物质等。

土体及其它各种介质的相关研究，借助邓肯张非线性弹性模型表达材料性能，曾在国内外受到广泛关注。

反演分析作为实验室微观及宏观材料特性参数确定的技术工具，也被广泛应用到各类固体中。

由于土体及其它介质微观结构的复杂性，传统的正演分析法在反演其变形特性参数时具有相当的局限性。

由于邓肯张非线性弹性模型的普遍性，它被认为是反演土体与其它介质的理想模型。

然而，由于邓肯张非线性弹性模型的复杂性，无法通过简单的正演分析方法来确定相关模型参数。

为了能够准确描述土体及其它介质的变形特性，本文对邓肯张非线性弹性模型参数反演分析进行了一系列研究。

首先，根据实验测量所得的弹性模量和泊松比数据，建立了解析解的求解过程，以解决邓肯张非线性弹性模型参数反演的问题；其次，利用最小二乘法的概念，对解析解式做进一步优化处理，以期达到更高的拟合准确度；最后，利用数值拟合算法（二次函数拟合、三次样条拟合及支持向量机方法），对反演得到的数据进行多次拟合优化。

本文通过研究邓肯张非线性弹性模型参数的反演分析，实现了土体及其它介质变形特性参数的智能识别。

结果表明：经过解析解与数值拟合过程的叠加，该方法能够极大地提高参数的反演精度，从而大大降低了反演变形特性参数时的误差。

本文反演分析的结果表明，邓肯张非线性弹性模型尤其适合于反演土体及其它介质变形特性参数，这可以为后续研究工作提供很好的基础。

未来，可以继续研究其它的介质及反演方法，以进一步提高变形特性参数反演的准确性及效率。

总之，本文将邓肯张非线性弹性模型应用于反演土体及其它介质变形特性参数，研究结果证明了该模型在反演方面的可行性和有效性，为土体及其它介质物理性能参数确定提供了一种有效的方法。

可编辑修改精选全文完整版3.Duncan-Chang 模型参数的确定实验目的：Duncan 双曲线模型是一种建立在增量广义虎克定律基础上的非线性弹性模型，它在岩土工程界为人们所熟知和广泛应用。

这一类模型可以反映应力应变关系的非线性，参数的物理意义明确和易于确定， 本实验通过对不同围压的控制来模拟模型并确定其参数。

实验原理：点绘()a εσσ~31-曲线，如图3－1所示，Kondner 等人发现，可以用双曲线来拟和这些曲线。

对某一3σ，()a εσσ~31-关系可表示成：aab a εεσσ+=-31 （3－1）渐近线σ3＝常量E iE tσ1-σ3(σ1-σ3)uεa 0εa /(σ1-σ3)uεa ba图 3－1 ()a εσσ~31-关系曲线 图3－2 ()a a εσσε--31/关系曲线式中：a 和b 为试验常数。

上式也可以写成：a ab a εσσε+=-31 (3－2)以()31/σσε-a 为纵坐标，a ε为横坐标，构成新的坐标系，则双曲线转换成直线。

见图3－2。

其斜率为b ，截距为a 。

有增量广义虎克定律，如果只沿某一方向，譬如Z 方向，给土体施加应力增量ΔZσ，而保持其他方向的应力不变，可得：E zx σεΔΔ=(3－3) Ev zx σεΔΔ-= （3－4）则 xzE εσΔΔ= （3－5）zxv εεΔΔ-= （3－6）邓肯和张利用上述关系推导出弹性模量公式。

由式（3－5）得：()()aa E εσσεσσεσ∂-∂=-==313111ΔΔΔ （3－7）由此可见虎克定律中所用的弹性模量实际上是常规三轴试验()a εσσ~31-曲线的切线斜率。

这样的模量叫做切线弹性模量，可用t E 表示，见图3-1。

将式（3－1）代入式（3－7），得到：()2a tb a aE ε+= （3－8）由式（3－2）可得：ba a --=311σσε （3－9）式（3－9）代入式（3－8），得： ()[]23111σσ--=b a E t （3－10）由式（3－2）可得：当0→a ε时31→⎪⎪⎭⎫⎝⎛-=a aa εσσε（3－11）而双曲线的初始切线模量i E 为： 031→⎪⎪⎭⎫⎝⎛-=a a i E εεσσ （3－12） 见图3－1。

邓肯模型参数k、n、R f整理方法研究张波（西安地质矿产勘查开发院，陕西，西安710100.）摘要：对四种石渣料进行三轴CD试验，采用全部法（全部试验点）、70~95法（应力水平70~95%的试验点）及分段法（小应变试验点求参数k、n，应力水平70~95%的试验点求参数R f）三种方法来求取Duncan-Chang模型参数k、n、R f，并采用三种方法求取的参数拟合应变应变曲线与试验曲线对比，研究了不同求取方法的优劣。

结果表明：三轴试验应力应变曲线并不完全符合双曲线关系，ε1/（σ1-σ3）~ε1关系曲线表现出明显的双线性，转折点应变约为1%；石渣料虽然是无粘性土，但其粘聚力并不为零，在工程中应予适当考虑；全部法和70~95法求取的参数拟合曲线与试验曲线较为吻合，拟合效果较好，但考虑到应力应变软化问题，70~95法求取的参数能较好的模拟试验曲线，分段法求取的参数虽然符合物理力学意义，但由于求取参数k、n仅采用小应变的试验点，拟合曲线与试验曲线差别较大，建议采用应力水平70~95%的试验点求取Duncan-Chang 模型参数k、n、R f。

关键词：三轴试验；石渣料；Duncan-Chang模型；应力-应变Duncan model parameters k, n, Rf sorting method/Zhang Bo,(Xi’an Institute of Geological And Mineral exploration,Xi’an Shannxi710100,China)Abstract: The four kinds of gravel material for triaxial CD tests, using all methods (all test points), 70 - 95 method (stress level 70 ~ 95% of test points) and section method (small strain test point demand parameter k , n, the stress level of 70 to 95% of the test points request parameter Rf) are three ways to strike a Duncan-Chang model parameters k, n, Rf, and the use of three methods of parameter fitting to strike a strain-strain curves and test curves comparative study of the pros and cons of different methods to strike.The results showed that: Triaxial stress-strain curve does not fully comply with the h yperbolic relationship, ε1 / (σ1-σ3) ~ ε1 curve showed a clear bilinear, the turning point is about 1% strain; carbide material although no viscous soil, but its cohesion is not zero, the project should be properly taken into account; all law and to strike the law from 70 to 95 parameter fitting curve is more consistent with the experimental curve, fit better, but considering the stress strain-softening problems, 70 to 95 law parameters can strike a better simulation curve, although the parameters of sub-strike law in line with the physical and mechanical sense, but to strike a parameter k, n using only a small strain of test points, the proposed curve together with the experimental curves vary greatly, the stress level is recommended 70 to 95% of the test points to strike Duncan-Chang model parameters k, n, Rf.Keywords: triaxial tests; carbide material; Duncan-Chang model; stress - strain1 引言土的本构模型研究是目前岩土工程的重要课题之一，目前主要有非线性弹性模型及弹塑性模型[1]。

在ansys中导入自定义本构模型---邓肯-张模型(转载)邓肯-张模型的关键点是材料的弹性模量随大小主应力差及小主应力(围压)的变化而变化，用APDL实现之的基本思路是：给每个单元定义一个材料号，分级施加荷载，在每个荷载步结束时提取出各单元的大小主应力，据此计算出下个荷载步的弹性模量Et,修改各单元之MP，用于下一步计算。

以下是一个简单算例，copy出去可直接运行。

常规三轴试验模拟by taomingxing,NWPU2003.7.16FINISH/CLEAR/TITLE,Numerical Simulation of three axes testing of soils /PREP7*dim,SUy,array,50 !Settlement records*dim,MaxPs,array,120 !Max history p1-p3*dim,MaxDs,array,120 !Max history Ds!*dim,EEt,array,50 !Et of elememtDuncan-Chang ModelSymbols：c-粘滞力，Fai-内摩擦角，Sf-破坏强度(p1-p3)f，Ds-应力水平，Pa-大气压,P3-围压*CREATE,Duncan-Chang !Creat Macro file*afun,deg !Unit of angle*set,Pa,1e5*set,P1,-ArrS3(i) !注意：岩土工程中应力为拉负压正*set,P3,-ArrS1(i)*if,P3,LT,0.1*Pa,thenP3=0.1*Pa !围压最小取值*endifSf=2*(c*cos(Fai)+P3*sin(Fai))/(1-sin(Fai)) !Mohr-Coulomb破坏强度(p1-p3)fDs=(P1-P3)/Sf !应力水平，*if,Ds,GT,0.95,thenDs=0.95 !应力水平最大取值*endif!判断加卸荷，如果(P1-P3)小于历史最大值视为卸荷－再加荷过程 *if,MaxPs(i),LT,P1-P3,thenEi=k*Pa*(P3/Pa)**nEt=Ei*(1-Rf*Ds)**2 !加荷情况的切线模量MaxPs(i)=P1-P3 !保存历史最大应力*elseif,MaxPs(i),GE,P1-P3Et=Kur*Pa*(P3/Pa)**n !卸荷模量*endifmp,ex,i,Et !修改单元i的Etmp,nuxy,i,Mu*END单元类型et,1,42 !平面四节点单元KEYOPT,1,3,2 !平面应变以下定义材料初始模量mp,ex,1,3.728e7 !砂土的弹性模量mp,nuxy,1,0.33mp,dens,1,1800建立几何模型blc4,0,0,0.08,0.15 !8cm X 15cm /PNUM,AREA,1/REPLOT网格划分aesize,all,0.01mat,1amesh,all边界条件nsel,s,loc,y,0d,all,Uy !底边界竖向约束nsel,s,loc,x,0d,all,Ux !左侧边界水平向约束nsel,all/replotfini/SOLUtime,0.01 !施加围压sfl,all,pres,2e5 !200kPasolve分级施加荷载，实现非线性计算荷载增量10kPa,共50级*DO,ti,1,50取出计算结果，修改弹性模量 /POST1*get,SUy(ti),node,29,u,y !Settlement record of time ti ETABLE,EtabS1,S,1 !取各单元第一主应力ETABLE,EtabS3,S,3 !取各单元第三主应力*dim,ArrS1,array,120*dim,ArrS3,array,120*do,Num,1,120 !Num为单元编号*get,ArrS1(Num),elem,Num,etab,EtabS1 !将单元结果存入数组*get,ArrS3(Num),elem,Num,etab,EtabS3*enddo/PREP7!^^修改砂土单元的Et,单元号1-120c=0 $Fai=35 $Rf=0.7 $k=400 $n=0.6 $Mu=0.33 $Kur=326.7 *do,i,1,120 !各单元循环计算*use,Duncan-Chang,c,Fai,Rf,k,n,Mu,Kur !调用Duncan-Chang宏文件*enddo!EEt(ti)=ET !保存第120单元之ET/SOLUtime,tisfl,3,pres,2e5+1e4*ti !施加荷载，增量1e4solve !对ti级荷载情况求解*ENDDO。

（1）根据邓肯等人总结的经验公式计算参数a 、b ：总r] +(泊~[^v\[(勺「+(勺人订5 丿95% k a i ^3 J TO % k 0!5 丿』计算得到表一如下。

围压(kpa)（5一。

3开（W_。

3）95%(£1)95%（5_。

3）70%(^1)70%b (5 一a3)ulta Rf100289.4274.9300233 8 2015800103 1 0 00 261 26 382.765 09E -05 0756 079300 805.8 765 51 00317 6 564 060011960 00 102 45 976 082.12E -05 0 825 549500 1323 9 1257 71 00339 7926 730 0129 20.00 062 081610.821 39E -05 0 821882对Rf 取平均值可得:又因为a 为起始变形模量§的倒数，即1r可得表二，并绘制lg （Ei/Pa ）与Ig （o3/Pa ）的试验关系图如图一所示。

表二围压(kpa)a Eilg (Ei/Pa) Ig(o3/Pa) 100 5.09E-05 19648.88458 1287299947 -0.006037955 3002.12E-054716120736 2.667556168 047108330 500 1.39E-0571728.328172.8496527540 692932049(5 -二 b(q - s )f1 E 1 aPa Pa PaR^ + Rf 2 + Rf 33 = 0.80117（6 —。

3）1讥 （^1）95% 一（“）70%2At对图一中的试验点进行拟合，得到lg (Ei/Pa)与Ig(o3/Pa)的直线关系: 尸0.8033X+2.2914.根据公式：E=5③可知K、n分别代表lg (Ei/Pa)与lg(a 3/Pa)直线的截距和斜率，故可得K=2.2914： n=0.8033oE-v 法在常规三轴试验中，轴向应变£ 1与侧向应变一£ 3之间也存在双曲线关系，经 变换之后可得如下公式：由上式知一£ 3/8 1与一£ 3为直线关系，但实际上，二者并不是严格的直线关 系，需先对试验结果进行収舍，然后选取某一区间进行拟合。

应用MATLAB确定邓肯－张双曲线模型中的K,n参数简介：接合承德中密砂常规三轴试验数据，介绍应用Matlab语言编写计算及绘图程序来处理试验数据的方法，可显著提高试验研究的数据处理效率和结果的可视化程度。

关键字：Matlab 三轴试验邓肯－张模型1 前言基于广义胡克定律的线弹性理论形式简单，参数少，物理意义明确，而且在工程界有广泛深厚的基础，得以应用于许多工程领域中。

早期土力学中的变形计算主要是基于线弹性理论的，只有在计算机得到迅速发展之后，非线性理论模型才得到较广泛的应用。

邓肯－张模型是建立在增量广义胡克定律基础之上的变模量的弹性模型，可以反映土变形的非线性，并在一定程度上反映土变形的弹塑性，很容易为工程界所接受，加之所用参数和材料参数不多，物理意义明确，只需用常规三轴压缩试验即可确定这些参数及材料常数适应的土类比较广，所以该模型为岩土工程界所熟知，并得到了广泛的应用，成为土的最为普及的本构模型之一。

本文主要是应用MATLAB编写计算及绘图程序来处理承德中密砂常规三轴试验数据。

2 基于MATLAB的计算过程实现现场的观测数据经过采集和整理后，按照一定的格式把数据存储在数据文件中，然后可以使用MATLAB丰富的数值运算功能可以非常容易地编制出数据处理程序，先用函数fope n()打开数据文件，fid=fopen(‘filename’,’r’)再用fscanf 函数依次从文件中读取格式化数据来完成对各变量地赋值，其使用语法为：matrix=fscanf(fid,format)。

本文由于数据不是太多，所以在计算过程中没有采取调用存储文件地形式。

直接在计算过程中输入试验数据计算。

2.1 数据的处理对第一组数据，通过编写Matlab语言，由轴向应变和应力差的试验数据可以作出～（）和～双曲线关系图形，主要用到的MATLAB命令为：plot(x1,y1)；axis([0 0.04 0 3]) ；hold on%（1）plot(x1, x1./y1)；a=polyfit(x1, x1./y1,1)；t1=0:0.001:0.07;plot(x1, x1./y1,'.',t1,a(1)*t1 +a(2))%（2）其中x1代表第一组轴向应变，x2代表第一组应力差。

研究生课程作业邓肯张模型参数计算学生姓名李俊学科专业岩土工程学号201420105614任课教师周小文教授作业提交日期2014年12月1．计算轴向应变ch h∆∑=1ε式中 1ε－轴向应变；h ∆∑－固结下沉量，由轴向位移计测得0h －土样初始高度c h —按实测固结下沉的试样高度c h ∆—试样固结下沉量2．计算按实测固结下沉的试样高度，面积：式中 Ac －按实测固结下沉的试样面积0V －土样初始体积3．计算剪切过程中试样的平均面积：式中 a A －剪切过程中平均断面积c V －按实测固结下沉的试样的体积i V ∆－排水剪中剪切时的试样体积变化 按体变管或排水管读数求得1h ∆－固结下沉量，由轴向位移计测得 3. 计算主应力差cic h V V A ∆-=01h h V V A c i c a ∆-∆-=Cc c A h V ⨯=1031⨯=-aA CR σσ 式中 31σσ- － 主应力差 1σ―大主应力 3σ－小主应力 C －测力计率定系数 R －测力计读数2 数据处理2.1 3σ=100kPa 数据初步计算当3σ=100kPa 时，各数据初步计算如表1所示。

围压100kPa 数据初步计算表 表12.1.1 由切线模量计算数据 对公式)(311σσε-=a +b 1ε进行直线拟合，如图1所示。

图11131/()~εσσε-拟合曲线 a =0.0002，1i E a==5000kPa b ==0.0028，()131ult bσσ-==263.16kPa ()13f σσ-=204.26kPa ，()()1313f fultR σσσσ-=-=0.77622.1.2 由泊松比计算数据对公式()313/f D εεε-=+-进行直线拟合，如图2所示。

图2 313/~εεε--拟合曲线f=i ν=0.2122 D=2.72972.2 3σ=200kPa 数据初步计算当3σ=200kPa 时，各数据初步计算如表2所示。