FLAC3D命令流

命令流1genzonebrip0000p11000p20100p30010&p415150p501515p615010p7202020&size101010rat1.00.91.1groupbrick_1genzonebrip02000p1add1000p2add0200p3add0015&size101010rat1.00.91.1groupbrick_2genzonebrip04000edge10size101010rat1.00.91.1groupbrick_3plotsur2-1newgenzonbrisize333modelelaspropbulk3e8shear1e8inidens2000fixzranz-.1.1fixxranx-.1.1fixxranx2.93.1fixyrany-.1.1fixyrany2.93.1setgrav00-10solveappnstress-10e4ranz3x12y12histgpvel003histgpvel033plohist1redploaddhist2bluesolve3-1;-------------------------------------信息;ProjectRecordTreeexport;Title:Simpletest;---------------------------------计算第一步;...STATE:STATE1....configgrid10,10modelelasticgroup'User:Soil'notnullmodelelasticnotnullgroup'User:Soil'propdensity=1500.0bulk=3E6shear=1E6notnullgroup'User:Soil'fixxyj1fixxi1fixxi11setgravity=9.81history999unbalancedsolvesavestate1.sav;----------------------------------计算第二步;...STATE:STATE2....initialxdisp0ydisp0initialxvel0yvel0modelnulli47j810group'null'i47j810groupdelete'null'history1xdispi=4,j=11solvesavestate2.sav;--------------------------------绘图命令;***plotcommands****;plotname:syyplotholdgridsyyfill;plotname:Unbalancedforceplotholdhistory999;plotname:gridplotholdgridmagnify20.0lredgriddisplacement;plotname:Xdis-Aplotholdhistory1line5-1new;===============================;定义球体半径和半径方向上单元网格数;===============================defparmrad=10.0rad_size=5endparm;===============================;建立八分之一球体外接立方体网格;===============================genzonepyramidp0rad00p1rad0radp2radrad0p3000&p4radradradsizerad_sizerad_sizerad_sizegroup1genzonepyramidp00rad0p1radrad0p20radradp3000&p4radradradsizerad_sizerad_sizerad_sizegroup2genzonepyramidp000radp10radradp2rad0radp3000&p4radradradsizerad_sizerad_sizerad_sizegroup3;==================================;利用FISH语言将内部立方体节点调整到球面;==================================defmake_spherep_gp=gp_headloopwhilep_gp#null;获取节点点坐标值：P=(px,py,pz)px=gp_xpos(p_gp)py=gp_ypos(p_gp)pz=gp_zpos(p_gp)dist=sqrt(px*px+py*py+pz*pz)ifdist>0then;节点位置调整maxp=max(px,max(py,pz))k=(maxp/rad)*(rad/dist)gp_xpos(p_gp)=k*pxgp_ypos(p_gp)=k*pygp_zpos(p_gp)=k*pzend_ifp_gp=gp_next(p_gp)end_loopendmake_sphere;===============================;利用镜像生成完整球体网格;===============================genzonerefgenzonerefdip90genzonerefdip90dd90;===============================;显示球体网格;===============================plotsurfplsetbackwhplblgr5-3ngenzonbrisize112groupsoilranz11grouprockranz01expgrid1.flac3d6-1genzonbrisize333modelmohrpropbu3e6sh1e6coh10e3fric15fixzranz-.1.1fixxranx-.1.1fixxranx2.93.1fixyrany-.1.1fixyrany2.93.1inidens2000histunbalsetgrav10solveelasticsave6-1.sav6-2rest6-1.savinixd0yd0zd0xv0yv0zv0appnstress-100e3ranz2.93.1x12y12solvesave6-2.sav6-3rest6-1.savinixd0yd0zd0xv0yv0zv0appnstress-100e3ranz2.93.1x12y12histid=2gpzdis113histid=3gpzdis112histid=4gpxdis113histid=5gpxdis113histid=6zoneszz113histid=7zoneszz1.51.52.5histid=8zonesxz1.51.52.5solvesave6-3.sav6-4rest6-3.savsetlogonsetlogfile6-2.logprintzonestressprintgpdissetlogoff6-5rest6-1.savinixd0yd0zd0xv0yv0zv0appnstress-100e3ranz2.93.1x12y12plotsetrot20030plotconszzouonmagf10plotaddhist1setmovieavistep1file6-5.avimoviestartsolvemoviefinish7-1ngenzonbrisize112modelelaspropbulk3e7shear1e7fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solveploconsz7-2ngenzonbrisize112modelmohrpropbulk3e7shear1e7c1e10f15tension1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solvepropbulk3e7shear1e7c10e3f15ten0solveploconsz7-3ngenzonebricksize112modelmohrpropbulk3e7shear1e7coh10e3fri15ten0fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000setgrav00-10solveelasploconsz7-4newgenzonebricksize112modelmohrpropbulk3e7shear1e7coh10e3fri15ten0fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000iniszz-40e3grad0020e3ranz02inisyy-20e3grad0010e3ranz02inisxx-20e3grad0010e3ranz02setgrav00-10solveploconsz7-5ngenzonbrisize112modelmpropbulk3e7shear1e7c10e10f15ten1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000ranz01inidens1500ranz12iniszz-35e3grad0020e3ranz01inisyy-17.5e3grad0010e3ranz01inisxx-17.5e3grad0010e3ranz01iniszz-15e3grad0015e3ranz12inisyy-7.5e3grad007.5e3ranz12inisxx-7.5e3grad007.5e3ranz12inipp10e3grad00-10e3ranz01setgrav00-10solveploconsz7-6ngenzonbrisize112modelmpropbulk3e7shear1e7c10e10f15ten1e10fixzranz0fixxranx0fixxranx1fixyrany0fixyrany1inidens2000ranz02iniszz-50e3grad0020e3ranz01inisyy-30e3grad0010e3ranz01inisxx-30e3grad0010e3ranz01inipp30e3grad00-10e3ranz02appnstress-10e3ranz2setgrav00-10solveploconsz7-7newgenzonebrickp0000p16000p20600p30090&p460600p506090p6600150p76060150&size6610modelelasprobulk10e10she10e10iniden2500applysxx-1e9grad001.1111111e7rangex-.1.1applysxx-1e9grad006.6666666e6rangex59.960.1applysyy-1e9grad008.3333333e6rangey-.1.1applysyy-1e9grad008.3333333e6rangey59.960.1applyszz-1e8grad008.3333333e5ranz0120setgrav00-10step30000inixdisp0ydisp0zdisp0inixvel0yvel0zvel0plocontszz7-8newgenzonebrickp0000p16000p20600p30090&p460600p506090p6600150p76060150&size6610modelelasprobulk10e10she10e10iniden2500inisxx-1e9grad001.1111111e7rangex-.1.1inisxx-1e9grad006.6666666e6rangex59.960.1inisyy-1e9grad008.3333333e6rangey-.1.1inisyy-1e9grad008.3333333e6rangey59.960.1iniszz-1e8ranz-.1.1fixxyzranz-.1.1setgrav00-10solveinixdisp0ydisp0zdisp0inixvel0yvel0zvel0plocontszz7-9newgenzonebrickp000-50p127.50-50p205-50p300-10size8110groupclaygenzonebrickp027.50-50p11000-50p227.55-50p327.50-10ratio1.111size12110groupclaygenzonebrickp000-10p127.50-10p205-10p3000ratio110.8size814groupsoilgenzonebrickp027.50-10p11000-10p227.55-10p327.500ratio1.110.8size1214groupsoilgenzonebrickp0000p127.500p2050p3005p427.550&p5055p62005p72055size815groupdamfixxyzranz-49.9-50.1fixxranx-.1.1fixxranx99.9100.1fixymodelmohrranz-500modelnullranz05propbulk7.8e6shear3.0e6coh10e10tension1e10rangroupsoilinidens1500rangroupsoilpropbulk3.91e6shear1.5e6coh10e10tension1e10rangroupclayinidens1800rangroupclaysetgrav00-9.8histid=1unbalsolvepropbulk7.8e6shear3.0e6coh10e3fric15rangroupsoilpropbulk3.91e6shear1.5e6coh20e3fric20rangroupclaysolvesaveelastic.savinixdis0ydis0zdis0;将节点位移清零inixvel0yvel0zvel0;将节点速度清零histid=2gpzdis000;记录地基顶部中心点的沉降histid=3gpzdis27.500;记录路基坡脚处的沉降histid=4gpxdis27.500;记录路基坡脚处的水平位移modelelasticranz01;;激活0m~1m的单元propbulk7.8e6shear3.0e6ranz01inidens1500ranz01solve;按软件默认精度求解savefill-1.savmodelelasticranz12propbulk7.8e6shear3.0e6ranz12inidens1500ranz12solvesavefill-2.savmodelelasticranz23propbulk7.8e6shear3.0e6ranz23inidens1500ranz23solvesavefill-3.savmodelelasticranz34propbulk7.8e6shear3.0e6ranz34inidens1500ranz34solvesavefill-4.savmodelelasticranz45propbulk7.8e6shear3.0e6ranz45inidens1500ranz45solvesavefill-5.savpau;ploblgr;genzonebrickp0000p110000p2050p3005sizegenzonebrickp000-50p127.50-50p205-50p300-10size8110groupclaygenzonebrickp027.50-50p11000-50p227.55-50p327.50-10ratio1.111size12110groupclaygenzonebrickp000-10p127.50-10p205-10p3000ratio110.8size814groupsoilgenzonebrickp027.50-10p11000-10p227.55-10p327.500ratio1.110.8size1214groupsoilgenzonebrickp0000p127.500p2050p3005p427.550&p5055p62005p72055size815groupdam7-10setlogon;打开log记录setlogfile1.log;设置记录文件名为：1.logrestorefill-1.sav;调用保存的文件printgpdisrangeid517anyid533any;输出两个节点的变形值restorefill-2.savprintgpdisrangeid517anyid533anyrestorefill-3.savprintgpdisrangeid517anyid533anyrestorefill-4.savprintgpdisrangeid517anyid533anyrestorefill-5.savprintgpdisrangeid517anyid533anysetlogoff;关闭log记录8-1defabcabc=1+2*3abcd=1.0/2.0endabcprintfish8-2defabcifaa<0thenabc=0.0elseabc=2.0*aaendifendabc8-3defabcloopaa(1,2.5)commandprintaaendcommandendloopend8-4newgenzonbrisize333modelelasticpropbu3e7sh1e7inidens2000fixxyzranz-.1.1fixxranx-.1.1fixxranx2.93.1fixyrany-.1.1fixyrany2.93.1setgrav10solveinixd0yd0zd0xv0yv0zv0save8-4.sav8-5rest8-4.savdefE_modifyp_z=zone_headd_k=704d_n=0.38d_pa=101325.0;//大气压loopwhilep_z#nullsigma_3=-1.0*z_sig1(p_z)E_new=d_k*d_pa*(sigma_3/d_pa)^d_nz_prop(p_z,'young')=E_newp_z=z_next(p_z)endloopendE_modify8-6rest8-5.savtable1nameload_settlementdefadd_loadp_gp=gp_near(2,1,3)loopn(1,5)app_load=n*(-1000e3)file_name='7-6_add_step'+string(n)+'.sav'commandappnstressapp_loadranz2.93.1x12y12solvesavefile_nameendcommandxtable(1,n)=-1.0*app_loadytable(1,n)=gp_zdisp(p_gp)endloopendadd_loadsave8-6.sav8-7rest8-6.savdeffind_max_dispp_gp=gp_headmaxdisp_value=0.0maxdisp_gpid=0loopwhilep_gp#nulldisp_gp=sqrt(gp_xdisp(p_gp)^2+gp_ydisp(p_gp)^2+gp_zdisp(p_gp)^2)ifdisp_gp>maxdisp_valuemaxdisp_value=disp_gpmaxdisp_gpid=gp_id(p_gp)endifp_gp=gp_next(p_gp)endloopendfind_max_dispprintmaxdisp_valuemaxdisp_gpidrest8-6.savconfigzextra1defget_sigma_difp_z=zone_headloopwhilep_z#nullsigma_dif=z_sig3(p_z)-z_sig1(p_z)z_extra(p_z,1)=sigma_difp_z=z_next(p_z)endloopendget_sigma_difplotconzextra19-1;---------------------------------------------------;;移来移去法接触面的建立;---------------------------------------------------ngenzoneradcylp0(0,0,0)p1(8,0,0)p2(0,0,-5)p3(0,8,0)&p4(8,0,-5)p5(0,8,-5)p6(8,8,0)p7(8,8,-5)&p8(.3,0,0)p9(0,.3,0)p10(.3,0,-5)p11(0,.3,-5)&size310615ratio1111.15genzoneradcylp0(0,0,-5)p1(8,0,-5)p2(0,0,-8)p3(0,8,-5)&p4(8,0,-8)p5(0,8,-8)p6(8,8,-5)p7(8,8,-8)&p8(.3,0,-5)p9(0,.3,-5)p10(.3,0,-8)p11(0,.3,-8)&size36615ratio1111.15fillgenzonereflectdd270dip90groupclay;interface1facerangecylinderend1(0,0,0)end2(0,0,-5.1)radius.31&cylinderend1(0,0,0)end2(0,0,-5.1)radius.29notinterface2facerangecylinderend1(0,0,-4.9)end2(0,0,-5.1)radius.31;pausegenzonecylp0(0,0,6)p1(.3,0,6)p2(0,0,1)p3(0,.3,6)&p4(.3,0,1)p5(0,.3,1)&size3106genzonecylp0(0,0,6.1)p1(.3,0,6.1)p2(0,0,6)p3(0,.3,6.1)&p4(.3,0,6)p5(0,.3,6)&size316genzonereflectdd270dip90rangez16.1grouppilerangez16.1pauseinizadd-6.0rangegrouppilesavepile_geom.sav9-2;导来导去法;---------------------------------------------ngenzonebricksize333group2rangex12y12z12group1rangegr2notsave1.savdelrangroup2notinterface1facerangex1y12z12interface1facerangex2y12z12interface1facerangex12y1z12interface1facerangex12y12z1interface1facerangex12y12z2rest1.savdelrangroup2expgrid1.fac3dimpgrid1.flac3dmodelelafixxyzranz0iniden2000setgrav00-10interface1propkn20e6ks20e6coh10e3fri15appnstr-200e3ranx01y12z3solve9-3;切割模型法;------------------------------ngenzonebricksize333group1rangex12y12z23group2rangegroup1notgenseparate1int1wrap12int1maxedge0.5plointred9-4;CreateMaterialZonesgenzonebricksize555&p0(0,0,0)p1(3,0,0)p2(0,3,0)p3(0,0,5)&p4(3,3,0)p5(0,5,5)p6(5,0,5)p7(5,5,5)genzonebricksize555p0(0,0,5)edge5.0groupMaterial;CreateBinZonesgenzonebricksize155&p0(3,0,0)p1add(3,0,0)p2add(0,3,0)&p3add(2,0,5)p4add(3,6,0)p5add(2,5,5)&p6add(3,0,5)p7add(3,6,5)genzonebricksize155&p0(5,0,5)p1add(1,0,0)p2add(0,5,0)&p3add(0,0,5)p4add(1,6,0)p5add(0,5,5)&p6add(1,0,5)p7add(1,6,5)genzonebricksize515&p0(0,3,0)p1add(3,0,0)p2add(0,3,0)&p3add(0,2,5)p4add(6,3,0)p5add(0,3,5)&p6add(5,2,5)p7add(6,3,5)genzonebricksize515&p0(0,5,5)p1add(5,0,0)p2add(0,1,0)&p3add(0,0,5)p4add(6,1,0)p5add(0,1,5)&p6add(5,0,5)p7add(6,1,5)groupBinrangegroupMaterialnot;Createnamedrangesynonymsrangename=BingroupBinrangename=MaterialgroupMaterial;AssignmodelstogroupsmodelmohrrangeMaterialmodelelasrangeBingenseparateMaterialinterface1wrapMaterialBinrangeplaneori000normal1-10aboveinterface2wrapMaterialBinrangeplaneori000normal1-10belowint1maxedge0.55int2maxedge0.55;Assignpropertiespropshear1e8bulk2e8fric30rangeMaterialpropshear1e8bulk2e8rangeBininiden2000int1propks2e9kn2e9fric15int2propks2e9kn2e9fric15;AssignBoundaryConditionsfixxrangex-0.10.1anyx5.96.1anyfixyrangey-0.10.1anyy5.96.1anyfixzrangez-0.10.1Bin;Monitorhistorieshistunbalhistgpzdisp(6,6,10)histgpzdisp(0,0,10)histgpzdisp(0,0,0);Settingssetlargesetgrav0,0,-10;Cyclingstep4000savebin.sav9-5;;------------------------------------------------------------------restpile_geom.savmodelmohrrangegroupclaypropbulk8.333e7shear3.846e7coh30000fric0rangegroupclaymodelelasrangegrouppilepropbulk8.333e7shear3.846e7rangegrouppileinterface1propkn1e8ks1e8fric20coh30000interface2propkn1e8ks1e8fric20coh30000;inidens1230rangegroupclayinidens1230rangegrouppilemodelnullrangez-0.10.15;fixzrangez-8.1-7.9fixxrangex-8.1-7.9fixxrangex7.98.1fixyrangey-.1.1fixyrangey7.98.1setgrav00-10iniszz0.grad0012300.rangez-5.50.iniszz17600grad0015500rangez-8-5.5inisxx0.grad005271.4rangez-5.50.inisxx7542.86grad006642.86rangez-8-5.5inisxxadd31428.6grad005714.3rangez-8-5.5inisyy0.grad005271.4rangez-5.50.inisyy7542.86grad006642.86rangez-8-5.5inisyyadd31428.6grad005714.3rangez-8-5.5;waterdensity1000watertableorigin0,0,-5.5normal00-1inidens1550rangez-8-5.5histunbal;solverat1.e-6savepile0.sav;modelelasrangegrouppilepropbulk13.9e9shear10.4e9rangegrouppileinidens2500rangegrouppilecallfind_add.fissolverat1.e-6savepile1.savrestpile1.sav;调用保存文件inistate0inixdis0.0ydis0.0zdis0.0;位移清零applyszz-0.4e6rangez0.050.15grouppile;桩顶加第一级荷载solvesaveapp0.4.savprintgpdisprangeid1;输出第一级荷载下的桩顶位移，假定桩顶中心的id号为1applyszz-0.6e6rangez0.050.15grouppile;桩顶加第二级荷载solvesaveapp0.6.savprintgpdisprangeid1;输出第二级荷载下的桩顶位移…………………………………………;依次加载，直到桩破坏;--------------------------------------------------------------------;速度加载法restpile1.savinistate0inixdis0ydis0zdis0defzs_top;检测桩顶竖向荷载ad=top_headzftot=0.0loopwhilead#nullgp_pnt=mem(ad+1)zf=gp_zfunbal(gp_pnt)zftot=zftot+zfad=mem(ad)endloopzs_top=zftot/0.1414endfixzrangez0.05.15grouppile;固定桩顶速度，用速度来确定位移deframpwhile_steppingifstep<ncutthenudapp=float(step)*udmax/float(ncut)elseudapp=udmaxendifad=top_headloopwhilead#nullgp_pnt=mem(ad+1)gp_zvel(gp_pnt)=udappad=mem(ad)endloopendhistgpzdis0,0,0histgpzvel0,0,0histzs_tophistzoneszz0,0,-.1setmechdampcombsetudmax=-1e-8ncut30000step225000savepile2.sav;;--------------------------------------------------------------------;位移控制法defsolve_stepsloopn(1,21)save_file=string(n)+'-step.sav'commandstep40000savesave_fileprizonestressranid2381aid2361aid2341a;输出桩顶网格单元的应力endcommandendloopendsolve_steps;-----------------------------------------------------------------------------10-1selbeambeg000end200nseg2selbeambeg200end40-1nseg3selbeamid=2beg40-1end50-2nseg2plotselgeoidonnodonscale0.04plotadax10-2selnodeid=1000selnodeid=2200selnodeid=340-1selnodeid=450-2selbeamselid=1cid=1node12selbeamselid=1cid=2node23selbeamselid=1cid=3node34plotselgeoidonnodonscale0.04plotadax10-3defset_valsptA=25.0*sin(40.0*degrad)ptB=25.0*cos(40.0*degrad)endset_valsgenzonecylinderp0=(0.0,0.0,0.0)&p1=(ptA,0.0,ptB)&p2=(0.0,25.0,0.0)&p3=(0.0,0.0,25.0)&p4=(ptA,25.0,ptB)&p5=(0.0,25.0,25.0)&size=(1,2,2)selshellid=5rangecylinderend1=(0.0,0.0,0.0)&end2=(0.0,25.0,0.0)radius=24.5notplotblogroplotadselgeomblackblackcidonscale=0.03plotadaxpaudelete;deleteallzonesselnodeinitzposadd-25.010-4newgenzonebricksize688modelmohrpropbulk1e8shear0.3e8fric35propcoh1e10tens1e10setgrav00-9.81inidens1000fixxrangex-0.10.1fixxrangex5.96.1fixyrangey-0.10.1fixyrangey7.98.1fixzrangez-0.10.1histn5histunbalsetmechforce50solvesavebeam-brace0.sav;propcoh1e3tens1e3modelnullrangex24y26z510setlargeinixdis0ydis0zdis0selbeambegin=(2,4,8)end=(4,4,8)nseg=2selbeampropemod=2.0e11nu=0.30selbeampropXCArea=6e-3XCIz=200e-6XCIy=200e-6XCJ=0.0histgpzdisp448solvesavebeam-brace1.sav;plotcreateGravVplotsetplanedip90dd0origin340plotsetrot15020plotsetcenter2.54.24.0plotsetcapsize25plotaddcontdispplanebehindshadeonplotaddselbeamforcefxplotaddselgeomblackblacknode=offshrinkfac=0.0plotaddaxesplotshow10-5;非全长锚固、预紧力锚杆（锚索）模拟;方法1、通过删除－建立link链接来模拟托盘genzoneradtunp00,0,0p125,0,0p20,50,0p30,0,25size425410dim4444ratio1111.1fillgenzonereflectnormal100ori000genzonereflectnormal001ori000momohrprobulk2.2e9she1.3e9fric30coh1.3e6ten1.5e5inidens2000fixxrangex-25.1,-24.9fixxrangex24.925.1fixyrangey49.950.1fixzrangez-25.1-24.9fixzrangez24.925.1selcableid=1beg0,0,0end0,29,0nseg10selcableid=1beg0,29,0end0,35,0nseg6selcableid=1propemod2e10ytension310e3xcarea0.0004906&gr_coh1gr_k1gr_per0.0785rangecid1,10selcableid=1propemod2e10ytension310e3xcarea0.0004906&gr_coh10e5gr_k2e7rangecid11,17seldeletelinkrangeid1sellinkid=1001targetzonesellinkattachxdir=rigidydir=rigidzdir=rigidxrdir=rigidyrdir=rigidzrdir=rigidrangeid100selcableid=1pretension60e3rangecid1,10step2000sav10-5.sav10-6;非全长锚固、预紧力锚杆（锚索）模拟;方法2、通过设置极大锚固剂参数模拟托盘genzoneradtunp00,0,0p125,0,0p20,50,0p30,0,25size425410dim4444ratio1111.1fillgenzonereflectnormal100ori000genzonereflectnormal001ori000momohrprobulk2.2e9she1.3e9fric30coh1.3e6ten1.5e5inidens2000fixxrangex-25.1,-24.9fixxrangex24.925.1fixyrangey49.950.1fixzrangez-25.1-24.9fixzrangez24.925.1selcableid=1beg0,0,0end0,29,0nseg10selcableid=1beg0,29,0end0,35,0nseg6selcablepropemod2e10ytension310e3xcarea0.0004906&gr_coh1gr_k1gr_per0.0785rangecid2,10selcablepropemod2e10ytension310e3xcarea0.0004906&gr_coh10e5gr_k2e7rangecid11,17selcablepropemod2e10ytension310e3xcarea0.0004906&gr_coh10e8gr_k2e10rangecid1,1selcableid=1pretension60e3rangecid1,10step2000sav10-6.sav10-7;非全长锚固、预紧力锚杆（锚索）模拟;方法3：借助别的结构单元（如liner单元）来模拟托盘genzoneradtunp00,0,0p125,0,0p20,50,0p30,0,25size425410dim4444ratio1111.1fillgenzonereflectnormal100ori000genzonereflectnormal001ori000momohrprobulk2.2e9she1.3e9fric30coh1.3e6ten1.5e5inidens2000fixxrangex-25.1,-24.9fixxrangex24.925.1fixyrangey49.950.1fixzrangez-25.1-24.9fixzrangez24.925.1selcableid=1beg0,0,0end0,29,0nseg10selcableid=1beg0,29,0end0,35,0nseg6selcableid=1propemod2e10ytension310e3xcarea0.0004906&gr_coh1gr_k1gr_per0.0785rangecid1,10selcableid=1propemod2e10ytension310e3xcarea0.0004906&gr_coh10e5gr_k2e7rangecid11,17sellinerrangey=-.1,.1x=-1,1z=-1,1sellinerPROPiso=(25e9,0.15)thick=0.1;concretesellinerPROPcs_nk=8e8cs_sk=8e8&cs_ncut=0.0cs_scoh=0.0cs_scohres=0.0cs_sfric=0.0seldeletelinkrangeid1sellinkid=1001targetnodetgt_num18sellinkattachxdir=rigidydir=rigidzdir=rigidxrdir=rigidyrdir=rigidzrdir=rigidrangeid100selcableid=1pretension60e3rangecid1,10step2000sav10-7.sav10-8ntitleStructure_dynamic_analysis_lakewaterconfigdynselpileid=1beg000end001selpilepropdens2400&Emod1.0e10Nu0.3XCArea0.3&XCJ0.16375XCIy0.00625XCIz0.01575&Per2.8&CS_sK1.3e11CS_sCoh0.0CS_sFric10.0&CS_nK1.3e11CS_nCoh0.0CS_nFric0.0CS_nGapoffdeff1whilesteppingf0=10000*sin(10*dytime)np=nd_headloopwhilenp#nullifnd_pos(np,1,3)=1nd_apply(np,1)=f0endifnp=nd_next(np)endloopendselnodefixxyzxryrzrranid=1selsetdampcombinedplocrepileplocurrentpileplosetbackblackforewhitemag0.8ploaddselgeoidonsca.04magf1e3ploaddselfapplgreenmagf1e3ploaddselpilemommylbluelredmagf1e3axeyelsetmovieavistep100filepile.avimoviestartsolage1moviefinish11-1newconfdyn;打开动力计算功能genzonebricksize10510modelasmodnullrangex=0,5z=5,10;删除部分网格fixzrangex=-.1.1z=.110.1;设置静力边界条件fixzrangex=9.9,10.1z=.110.1fixyrangey=-.1.1fixyrangey=4.95.1propbulk2e8shear1e8;设置土体参数propbulk4e9shear2e9rangex=5,6z=5,10;设置墙体参数（土体参数的20倍）inidens2000;设置密度defsetup;动荷载中的变量赋值freq=1.0omega=2.0*pi*freqold_time=clockendsetup;执行变量赋值defwave;定义动荷载函数wave=sin(omega*dytime);定义动荷载变量endapplyxvel=1histwaverangez=-.1.1;施加动荷载applyzvel=0rangez=-.1.1histgpxvel5,2,0histgpxvel5,2,10histgpzvel5,2,10histdytimedeftim;估算程序运行的时间tim=0.01*(clock-old_time)endsetdynmultion;设置动态多步solveage1.0printtim;输出计算时间printdyn;输出动力计算相关信息savemult1.sav11-2newconfigdyngenzonebricksize1,1,50modelelaspropshear1e7bulk2e7inidens1000defsetupomega=2.0*pi*freqpulse=1.0/freqendsetfreq=4.0setupdefwaveifdytime>pulsewave=0.0elsewave=0.5*(1.0-cos(omega*dytime))endifendrangenamebottomz=-.1.1fixzrangez=.555;将上部网格都施加数值向约束applydquietsquietrangebottomapplysxz-2e5histwavesyz0.0szz0.0rangebottom;-2e5的系数来源于的值applynvel0planenorm0,0,1rangebottomhistgpxvel0,0,0histgpxvel0,0,25histgpxvel0,0,50histdytimehistwaveplotcreatehhhplotaddhist123vs4plotshowsolveage211-3new;第一步：静力计算阶段configdynsetdynoffgenzonebricksize632genzonebricksize232p0002genzonebricksize232p0402genzonewedgesize132p0202genzonewedgesize132p0432p1332p2402p3434&p4302p5404modelelasticpropbulk66667shear40000inidens0.0025setgrav00-10fixxrangex-0.010.01fixxrangex5.996.01fixyrangey-0.010.01fixyrangey2.993.01fixzrangez-0.10.1histunbalsolvesave11-3_1.sav;第二步：动力计算阶段setdynondefiniwaveper=0.01endiniwavedefwavewave=0.5*(1.0-cos(2*pi*dytime/per))endfreexyzranz-0.10.1;去掉模型底部原有的静力条件applynquietsquietdquietranz-0.10.1;静态边界条件applydstress1.0histwaveranz-0.10.1;加动力荷载applyff;施加自由场边界条件groupff_cornergroupff_sideranx06groupff_siderany03groupmain_gridranx06y03setdyntime=0;设置动力计算从0s开始histreset;清空已有的历史信息histunbalhistdytime;主体网格histgpxvel210histgpxvel215.0;柱体网格histgpxvel-1-10histgpxvel-1-15;平行于y方向的二维自由场网格histgpxvel-100histgpxvel-105.0;平行于x方向的二维自由场网格histgpxvel2-10histgpxvel2-15.0solveage0.015save11-3_2.sav11-4confdygenzonebricksize3,3,3modelelaspropbulk1e8shear0.3e8inidens1000fixzrangez-.1,.1setdyn=on,grav00-10,hist_rep=1histgpzdisp3.0,1.5,3.0histdytimeplotcreatehhplotaddhis1vs2save11-4_1.sav;保存文件，在后续计算中需调用该文件cyc15011-5;（1）质量分量和刚度分量共同作用rest11-4_1.savsetdyndamprayleigh122.8solveage=0.2titleverticaldisplacementversustime(mass&stiffnessdamping)plotshowpause;（2）只有质量分量rest11-4_1.savsetdyndamprayleigh222.8masssolveage=0.08titleverticaldisplacementversustime(massdampingonly)plotshowpause;（3）只有刚度分量rest11-4_1.savsetdyndamprayleigh222.8stiffnesssolveage=0.08titleverticaldisplacementversustime(stiffnessdampingonly)plotshow11-6rest11-4_1.savsetdyndamprayleigh0.0522.8sethist_rep=5solveage=0.5titleverticaldisplacementversustime(5%Rayleighdamping)plotshowpauserest11-4_1.savsetdyndamplocal0.1571;=pi*0.05sethist_rep=5solveage=0.5titleverticaldisplacementversustime(5%Localdamping)plotshow11-7;振动台试验的例子newconfigdynamicfluiddefmodel_dimh_R=0h_R1=h_R+1.0endmodel_dimgenzonbrip000-10p1300-10p201-10p3000p4301-10p5010p6300h_Rp7301h_Rsize30110groupsandgenzonbrip0000p1300h_Rp2010p3001p4301h_Rp5011p6300h_R1p7301h_R1size3011grouptop;genzonbrip000-.5p130-.5p201-.5p3000p431-.5p5010p6300p7310size30110modelelasticpropbulk=3e7shear=1e7fric=35inidens2000modelfl_isopropporo0.5perm1e-8inifmod2e8fdens1000inipp0grad00-10e3ranz0-10.0fixzranz-9.9-10.1fixxranx-.1.1fixxranx29.930.1fixysetgrav10setfluidoffdynoffinifmod0setmechrat1e-6solvedefini_conf_k0=1.0pnt=zone_headloopwhilepnt#nullval=_k0*z_szz(pnt)+(_k0-1.)*z_pp(pnt)z_sxx(pnt)=valz_syy(pnt)=valpnt=z_next(pnt)end_loopendini_confsolvesave11-7.sav11-8rest11-7.savsetdynonfluidoninifmod2e8setfluidpcutonmodelfinnrangrosandpropbulk=3e7shear=1e7fric=35rangrosandinidens2000rangrosandpropff_latency=50ff_switch=0ff_c1=0.8ff_c2=0.79ff_c3=0.45ff_c4=0.73rangrosand;扭剪试验结果defsetupfreq=5.0ampl=2omega=2.0*pi*freqendsetupdefsine_wavevv=9.36e-2*sin(omega*dytime)ifdytime<2.0sine_wave=dytime/2.0*vvelseifdytime<20.0sine_wave=vvelseifdytime<=30.sine_wave=(30.0-dytime)/10.0*vvendifendifendififdytime>30.0sine_wave=0.0endifendfreexapplyxvel=1.0histsine_waveranz-9.9-10.1applyxvel=1.0histsine_waveranx-.1.1applyxvel=1.0histsine_waveranx29.930.1setdyndamplocal.314callppr.datsethist_rep100setlarge;setdyndt3e-4setmechrat1e-20defsolve_agesloopn(1,39)save_file='11-8_'+string(n)+'s.sav'commandsolagen;savesave_fileendcommandendloopendsolve_agessavesave_filehistwrite20212223243031323334vs2file10-8_Outfile_pp.txthistwrite4041424344vs2file10-8_Outfile_xdis.txt12-1genzonbrisize333inipp30e3grad00-10e3plotconppoutlineon;outline的作用是在云图中显示网格的轮廓12-2genzonbrisize333waterdensity1000setgrav10watertableori003norm00112-3genzonbrisize333inipp30e3grad00-10e3ranx-.1.1fixppranx-.1.1plotconppoutlineon：显示节点上的孔压云图plotbconppoutlineon：显示单元中心点处平均孔压块图12-4genzonbrisize333apppp30e3grad00-10e3ranx-.1.1plotconppoutlineon：显示节点上的孔压云图plotbconppoutlineon：显示单元中心点处平均孔压块图12-5configfluidgenzonebricksize20110modelmohrpropbulk8.33e6shear3.85e6fric15coh10e3tens1e10;fixxrangex-.1.1fixxrangex19.920.1fixxyzrangez-.1.1fixy;---applyloadslowly---deframpramp=min(1.0,float(step)/200.0)endapplynstress=-40e3histramprangex-.13.1z9.910.1;---fluidflowmodel---modelfl_isoinifmod2e9;---porepressurefixedatzeroatthesurface---fixpp0rangez9.910.1;---settings---setfloff;---histories---histgppp2,.5,9;---test---step750saveload.savplotsetplaneori00.50norm010plotconppplaneouonplotaddfapredplane12-6configfluidsetfluidoffgenzonbrickp000-10size20110genzonbrickp0500p11500p2510p3905p41510p5915p61105p71115&size1015groupsoilgroupdamranx57z-50groupdamranid201aid211aid221aid231aid241agroupdamranid202aid212aid222aid232aid242amepropbu3e7sh1e7inipp0grad00-10e3ranz0-10inidens2000modelfl_isoproppor0.5perm1e-10inifden1000ften-1e10inisat0.0ranz05modelfl_nullrangrodam;inipp0rangrodamfixzranz-10fixxranx0fixxranx20fixysetgrav10solvesaveelastic.sav12-7restelastic.savinixd0yd0zd0xv0yv0zv0appnstress-40e3grad0010e3ranz04x09solvesavepressure.sav12-8restpressure.savsetfluidonmechoffinifmod2e3ften0.0rangrosoilinixd0yd0zd0xv0yv0zv0apppp40e3grad00-10e3ranz04x09apppp0ranz0x1520histid=10zoneppid215solve12-9configfluidsetfluidbiotongenzonebrip0000p1500p2020p30020size3110groupFacerangez1820groupSoftrangez1018;softsoil8mgroupClayrangez010modelmohrpropbulk33e6shear7e6c0f30rangegroupFacepropbulk.83e6shear.17e6c7.0e3f5rangegroupSoftpropbulk6.66e6shear1.42e6c8.0e3f20rangegroupClayinidens1500modelfl_isoiniftens-5e5inifdens1000inibiot_mod4e9propperm10e-13biot_c1propperm10e-10rangrouFacefixxyzrangez-0.10.1fixyrangy-.010.01fixyrangy1.92.1fixxrangex-.01.01fixxrangex4.95.1;InitialStressK0=1setgrav00-10iniszz-400e3grad0020e3inisyy-400e3grad0020e3inisxx-400e3grad0020e3inipp200e3grad00-10e3setfluidoffstep5000inixd0yd0zd0xv0yv0zv0;forfixedpponthePVDpositiondeffixPPloopi(1,10)pv=-80e3+(i-1)*10e3zt=20.1-(i-1)zb=20.1-icommandfixpppvrangx-.1.1zzbztend_commandend_loopendfixppstep5000setfluidoninibiot_mod4.3e7ranggroupSoftnot;modifybiot_modinibiot_mod9.0e6ranggroupSoftsetmechforce1.5e3setmechsubs100auto;slavesetfluidsubs30;---histories---defdayday=fltime/24/3600enddayhistid=1dayhistid=2gpzd2.5,0,20histid=3gpzd2.5,0,10solage7.776e613-1configcppudm;必须设置cppudm的选项才可以进行模型载入modelloadDuncanChang.dllgenzonbrip1.600p20.60p300.6size111;三轴试验尺寸：0.6*0.6*0.6mmodelduncanpropcohesion110e3friction48.5fricDel0.0ratiofail0.79ke704ne0.38kb303mb0.18kur844.8fixzranz-.01.01;模型底部边界的竖直向速度约束defsigma3;定义一个sigma变量便于进行多工况计算sigma3=-600e3endsigma3;设置初始应力appnstresssigma3ranx-.01.01appnstresssigma3ranx.59.61appnstresssigma3rany-.01.01appnstresssigma3rany.59.61appnstresssigma3ranz.59.61;iniden2190;为了与理论值作比较，不考虑重力、密度的作用iniszzsigma3inisyysigma3inisxxsigma3solve;得到加载前的初始应力inixdis0ydis0zdis0xvel0yvel0zvel0;位移和速度清零tab1nameloads;建立一个空表，用来保存荷载-沉降曲线plotab1;第1次加卸载defload1p_gp=gp_near(0,0,0.6)loopn(1,50)zss_load=sigma3-float(n)*20e3;加载1000kPacommandappnstresszss_loadranz.59.61solveendcommandz_dis=-1.0*gp_zdisp(p_gp)/0.6;得到应变z_load=(sigma3-zss_load);得到主应力差commandtab1z_disz_load;分别保存应变和主应力差endcommandend_loopendload1save13-1.sav13-2rest13-1.savdefunload1p_gp=gp_near(0,0,0.6)loopm(1,25)zss_load=sigma3-1000e3+float(m)*2e3;卸载500kPacommandappnstresszss_loadranz.59.61solveendcommandz_dis=-1.0*gp_zdisp(p_gp)/0.6;得到应变z_load=(sigma3-zss_load);得到主应力差commandtab1z_disz_load;分别保存应变和主应力差endcommandendloopendunload1pausetlogonsetlogfile1.logpritab1setlogoff;第2次加卸载defload2p_gp=gp_near(0,0,0.6)loopn(1,75)zss_load=sigma3-1000e3+500e3-float(n)*20e3;加载至2000kPacommandappnstresszss_loadranz.59.61solveendcommandz_dis=-1.0*gp_zdisp(p_gp)/0.6;得到应变z_load=(sigma3-zss_load);得到主应力差commandtab1z_disz_load;分别保存应变和主应力差endcommandend_loopendload2save13-2.sav14-1n;===================================;建立网格模型genzonebrick&p0000p1200p200.50p3003size313genzonebrick&p0200p12000p220.50p3203&size1713ratio1.0311genzonebrick&p0203p12003p220.53p312013&p4200.53p5120.513p620013&p7200.513size17117ratio1.0311;====================================;设置边界条件fixxyzrangez-0.10.1fixxrangex19.920.1fixxrangex-0.10.1fixy;====================================;初始地应力的生成modelelaspropdensity2000bulk3e9shear1e9setgravity00-10.0solveinixdisp0ydisp0zdisp0inixvel0yvel0zvel0;===================================;安全系数求解modelmohrpropdensity2000bulk1.0e8shear3.0e7&coh12380.0tens1.0e6fric20dila20solvefosfileslope3dfos1.savassociated14-2n;=====================================;建立网格模型genzonebrick&p0000p1200p200.50p3003size313genzonebrick&p0200p12000p220.50p3203size1713&ratio1.0311genzonebrick&p0203p12003p220.53p312013&p4200.53p5120.513p620013&p7200.513size17117ratio1.0311;**********************************************;自定义强度折减法defSSR;=====================================;定义有关参数及循环终止条件ait1=0.02k11=1.0k12=2.0ks=(k11+k12)/2loopwhile(k12-k11)>ait1coh1=12380/ksfri1=(atan((tan(20*pi/180))/ks))*180/pidila1=20.0ten1=1e6grav0=-10dens1=2000K1=1e8G1=3e7;=====================================;折减的实现过程commandmodelnull;初始应力场的生成modelelasticprobulk1e10she3e9densdens1fixxyzrangez-0.10.1fixxrangex19.920.1fixxrangex-0.10.1fixysetgrav00grav0solveinixdisp0ydisp0zdisp0inixvel0yvel0zvel0;塑性阶段求解modelmohrprobulkK1sheG1densdens1cohcoh1&frictionfri1dildila1tensten1setmechratio9.8e-6solvestep10000endcommand;二分法的实现过程ifmech_ratio<1.0e-5k11=ksk12=k12elsek12=ksk11=k11endifks=(k11+k12)/2endloop;=====================================;计算结果的保存fosfile0='_fos'+'.sav'commandsavefosfile0endcommandend;**********************************************;程序执行及结果显示SSRprks15-1n;导入分组网格impgrid1groupsoilrangegrouprocknot;定义材料性质modelmohrproden2200bulk4e7she1.5e7co8e4&fric28ten9e3dil10rangegroupsoilproden2900bulk2.5e10she1.2e10co4e6&fric45ten2e6dil15rangegrouprock;设置边界条件fixxyzrangez-0.10.1fixxyzrangez7.98.1fixy;加载appnstress-5e5rangex-0.10.1appnstress-5e5rangex4.95.1;设置初始条件inizvel-2.5e-5rangez7.98.1;定义轴向应力求解函数defaxi_stressf_accum=0.0pnt=gp_headloopwhilepnt#nullifgp_zpos(pnt)<0.15thenf_accum=f_accum+gp_zfunbal(pnt)endifpnt=gp_next(pnt)endloopaxi_stress=f_accum/0.5end;定义侧向应力求解函数deflat_stressg_accum=0.0pnt=gp_headloopwhilepnt#nullifgp_xpos(pnt)<0.1theng_accum=g_accum+gp_xfunbal(pnt)endifpnt=gp_next(pnt)endlooplat_stress=g_accum/0.8end;跟踪记录有关重要变量histn1histaxi_stresshistlat_stresshistgpzdis008histgpxdis004histgpxdis054histunbal;绘制应力-位移曲线plothist-1vs-3;axialstressvsaxialdisplace;求解step15000;结果保存及数据输出save1.savhistwrite1fileaxi_stress.dathistwrite2filelat_stress.dathistwrite3fileaxi_dis.dathistwrite4filelat_dis1.dathistwrite5filelat_dis2.dathistwrite6fileunbalance.dat15-2n;===================================;导入网格数据并生成“辅助”界面单元impgridwwgrowaterinterface1faceranx-4.9604.9y0.1300z-4.9704.9setgrav0-9.810waterden1000;===================================;定义生成潜水面的函数defwater_tablep_i=i_headp_ie=i_elem_head(p_i)　loopwhilep_ie#null;返回三个相邻界面单元的三个顶点的地址p_gp1=ie_vert(p_ie,1)　p_gp2=ie_vert(p_ie,2)p_gp3=ie_vert(p_ie,3);将顶点坐标赋予给网格节点x1=in_pos(p_gp1,1)y1=in_pos(p_gp1,2)z1=in_pos(p_gp1,3)x2=in_pos(p_gp2,1)y2=in_pos(p_gp2,2)z2=in_pos(p_gp2,3)x3=in_pos(p_gp3,1)y3=in_pos(p_gp3,2)z3=in_pos(p_gp3,3);以这三个节点，生成潜水面commandwatertablefacex1,y1,z1x2,y2,z2x3,y3,z3endcommandp_ie=ie_next(p_ie)endloopend;===================================;进行边坡分析impgridaagroupsoilrangegroup2anygroup4anygrouprockrangegroup1anygroup3any;生成潜水面water_table;删除为“辅助”网格单元和界面单元int1deledelerangegroupwater;初始化材料密度:饱和的和非饱和的defini_densp_z=zone_headloopwhilep_z#nullifz_group(p_z)='soil'thenifz_pp(p_z)#0.0then　　　　z_density(p_z)=2120elsez_density(p_z)=2020endifendififz_group(p_z)='rock'thenifz_pp(p_z)#0.0thenz_density(p_z)=3460elsez_density(p_z)=3360endifendifp_z=z_next(p_z)endloopendini_dens;施加边界约束条件fixxyzrangy-0.10.1fixzrangex0600z-0.10.1fixxrangex-0.10.1z0700fixzrangex0600z699.9700.1fixxrangex599.9600.1z0700;通过弹性求解生成初始应力场modelelasprobulk1e10she1e10rangegroupsoilprobulk1e10she1e10rangegrouprocksolvefo10;位移和速度归零inixdisp0ydisp0zdisp0inixvel0yvel0zvel0;定义材料特性;冰碛土物理力学参数指标modelmohrprobulk3.0e8she1.8e8co2.10e5fric25ten1e4&dil10rangegroupsoil;强节理化辉长岩物理力学参数指标probulk8e8she6e8co4.60e5fric37ten4e3&dil15rangegrouprockhistn=5histunbalplohist1setmechratio1.0e-6solve