资源描述：
《高等土力学课程-CamClay.doc》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。
1、基于修正剑桥模型模拟理想三轴不排水试验——两种积分算法的对比分析（CZQ-SpringGod）1、修正剑桥模型在塑性功中考虑体积塑性应变的影响，根据屈服面一致性原则，假定屈服函数对硬化参数的偏导为0，就获得了以理想三轴不排水试验为基础的修正剑桥模型屈服函数：（1）其中，，，，M为临界线斜率，为前期固结压力。硬化/软化法则：（2）式中为体积塑性应变，v为比体积，为正常固结线斜率，为回弹线斜率。由于不排水屈服面推导过程是基于硬化参数偏导为0，也就是说不排水试验中硬化参数同体积塑性应变无关，屈服面不变化，而若引入硬化法则就
2、同屈服面推导过程中的假定矛盾，因此计算时将模型处理为理想塑性模型。2、显式和隐式两种积分格式考虑应变增量驱动下，第n增量步到第n+1增量步之间的应力积分格式。显式积分格式的推导参考文献[1]，其中弹塑性矩阵中的塑性硬化模量H=0。隐式积分格式推导如下：（3）（4）（5）（6）（7）在这一组方程中没有硬化规律方程表明为理想塑性，并将式（3）-（7）合并化简得到：（8）式中求解（8）式方程组即可得到n+1增量步的各个增量。两种积分格式的matlab程序分别见邮件附件文件夹camclayexp和camclayimp，显式运
3、行主程序为camclayexp.m，而隐式运行主程序为camclayimp.m。3、数值分析（1）修正剑桥模型的参数设定：临界线斜率：M=1.1正常固结线斜率：=0.17回弹线斜率：=0.034初始比体积：v0=2.12前期固结压力：=100KPa剪切与体积模量的比值：GK=0.46155每个增量步体积模量的计算：剪切模量G=GK×K其中固结线方程为：。（2）计算结果：不排水有效应力路径：（a）显示算法(b)隐式算法图1不排水有效应力路径偏应力随轴向应变的变化：（a）显示算法（b）隐式算法图2偏应变随轴向应变的变化孔
4、隙水压力随轴向应变的变化：（a）显示算法（b）隐式算法图3孔压随轴向应变的变化两种算法的每个增量步同屈服面的偏移程度：（a）显式算法（b）隐式算法图4每个增量屈服面的偏移程度结论：两种算法在计算理想塑性修正剑桥模型时，数值解能很好地同理论屈服面符合。显示算法的误差是递增的，而隐式是收敛的。理想塑性模型的分析结果表明，经过屈服面修正后的显示算法在精度上要高于隐式算法，可能同收敛参数的设定有关，不过两者都是精确的。参考文献：[1]S.W.Sloan.A.J.Abbo.D.Sheng.Refinedexplicitinte
5、grationofelasoplasticmodelswithautomaticerrocontrol[J].EngineeringComputations.2001:18,121-19程序代码：显式积分算法：(ExplicitIntegrationAlgorithm)%functioncamclayexp%%Undrainedcondition(perfectplasticity)%%initializationofparameterek=0.034;%回弹斜率lam=0.17;%固结斜率M=1.1;%临界线斜率v
6、0=2.12;%初始比体积GK=0.46155;%剪切与体积模量的比值pc=100;%初始固结压力%%PreliminaryS=[pcpcpc000];[Pst,deviS]=deviT(S);[J2,J3,sJ2,q,lode]=invar(deviS);E=[000000];nstep=300;de1=0.0004;q1=0;dEpvol=0;devidEp=zeros(1,6);forn=1:nsteppcre(n)=pc;Sre(n,:)=S;pre(n)=Pst;qre(n)=q;q1re(n)=q1;%%
7、strainincrementdE=[de1-de1/2-de1/2000];v1=v0-lam*log(pc);%固结曲线K=v1*Pst/ek;%体积模量G=GK*K;%剪切模量m=[111000];De=K*m'*m+2*G*(eye(6)-m'*m/3);%弹性刚度矩阵dre(n)=det(De);var(n)=q/Pst;[meanE,devidE]=deviT(dE);dEvol=meanE*3;ddS=(De*dE')';%弹性应力增量pc=harden(pc,v1,lam,ek,0);%%px(1)=
8、Pst;py(1)=q;%%incrementofstrain:initializationYF1=ydfun(sJ2,Pst,pc,M);S1=S+ddS;[Pst1,deviS1]=deviT(S1);[J2p,J3p,sJ2p,qp,lodep]=invar(deviS1);YF2=ydfun(sJ2p,Pst1,pc,M);ifYF2<
