砂土地基碾压加密的模拟
碾压法和夯实法是常用的地基处理方式，碾压法是采用碾压、振动压实机械，如压路机、振动压路机、推土机、羊足碾等机械，来回反复碾压、震动使地基土密实、强度提高、压缩性降低，从而使地基得到处理的一种密实方法。
本文主要针对于砂土碾压法进行建模，建模步骤主要分为两步，首先生成一个比较松的砂土地基。之后在砂土地基上生成一个碾轮，并让其来回滚动，使得砂土加密。
先看一下地基的生成方式，为了使研究结果比较好，这里的砂土没有经历过成样、预压然后自重，而是直接生成带有摩擦系数的颗粒，在自重下沉降。
newdef parrdmin=6e-3rdmax=9e-3width=4height=width/3.0end@pardomain extent [-width] [width]wall generate box [-width*0.5] [width*0.5] [-height*0.5] [height*0.5] expand 1.5ball distribute porosity 0.5 radius [rdmin] [rdmax] box [-width*0.5] [width*0.5] [-height*0.5] [height*0.5]cmat default model linearcmat add 1 model linear method deformability emod 100e6 kratio 1.5 property fric 0.5 range contact type ball-ballcmat add 2 model linear method deformability emod 100e6 kratio 1.5 range contact type ball-facetball attribute density 2.7e3 damp 0.7cycle 1000 calm 50cycle 2000set gravity 9.8solvesave sample
注意，需要使用calm命令清除一下颗粒的速度，沉降后的试样如下图所示：
可以看到图中的颗粒因为自重沉降，导致地面不是很平衡，所以我们对其进行削面处理，删除表面的颗粒，然后平衡。
restore sample[diding=height*0.5*0.22]ball delete range y [diding] [height*0.5]cycle 1solvesave diji
这步运行完后如下图所示：
之后就是碾轮的生成了，我们使用wall来模拟碾轮，对碾轮施加竖向的伺服力来模拟碾轮对地基的重力。可以根据密度以及体积来计算重力的大小。这时候也可以删除上部的墙了。这里使用wall sevro命令来施加伺服力，这其实是对之前伺服程序的一个包装，里面的参数都可以在伺服程序中找到。
restore diji[nianlun_rad=height*0.1][force=3.8e3*9.8*math.pi*nianlun_rad^2]wall generate id 10 name nianlun circle position [-width*0.5 nianlun_rad] [diding nianlun_rad*1.01] radius [nianlun_rad]wall delete range id 3wall servo activate on gainupdate 100 vmax [height] yforce [-force] range id 10ball property fric 0.5wall property fric 0.5 range id 5 68cmat modify 2 model linear method deformability emod 100e6 kratio 1.5 range contact type ball-facetcycle 1000solve aratio 1e-5save create_nianlun
注意这里的碾轮和砂土地基间是有摩擦的，不然轮子会整个陷进地基中。
最后一步就是进行碾轮的反复碾压了。
时间和位移场清零肯定是必须的。
关键的一个点是根据移动速度去计算滚动速度，这样保证轮子和颗粒间发生的是静摩擦。还需要注意时刻更新碾轮的滚动中心，不然轮子就飞了。
这里设置了测量圆去监测地基孔隙率的变化，并用history进行记录。
还需要注意这里的循环方式，速度有了，碾轮的居留我们也是有的，那么碾压一次需要的时间是可以计算的。那么根据时间来实现循环碾压也就是可行的了，每次到达碾压时间后，反向移速和滚动速度即可，也可以顺便保存个save文件进行后处理。
restore create_nianlunset mech age 0ball attribute displacement multiply 0[yisu=0.2][spin=yisu/nianlun_rad][cishu=5][mea_rad=width*0.05]wall attribute xvel [yisu] spin [-spin] range id 10[wp_nianlun=wall.find(10)]def updata_centerwhilesteppingwall.rotation.center(wp_nianlun)=wall.pos(wp_nianlun)end[onetime=(width-2*nianlun_rad)/yisu][timerecord=0][file_count=1]measure deletemeasure create id 1 position [-width*0.5*0.2] [-height*0.5*0.3] radius [mea_rad]measure create id 2 position [-width*0.5*0.1] [-height*0.5*0.3] radius [mea_rad]measure create id 3 position [width*0.5*0.1] [-height*0.5*0.3] radius [mea_rad]measure create id 4 position [width*0.5*0.2] [-height*0.5*0.3] radius [mea_rad]def mp_getmp1 = measure.find(1)mp2 = measure.find(2)mp3 = measure.find(3)mp4 = measure.find(4)end@mp_getdef jiancewhilesteppingmeasure_poros1=measure.porosity(mp1)measure_poros2=measure.porosity(mp2)measure_poros3=measure.porosity(mp3)measure_poros4=measure.porosity(mp4)measure_av=(measure_poros1 measure_poros2 measure_poros3 measure_poros4)*0.25forceLunY = wall.force.contact.y(wp_nianlun)forceLunX = wall.force.contact.x(wp_nianlun)time_now = mech.ageweiyiX = yisu * mech.agenianlun_pos = wall.pos.x(wp_nianlun)end@jiancedef circle_nianyawhilesteppingif mech.age-timerecord>onetime thenwall.vel.x(wp_nianlun)=-wall.vel.x(wp_nianlun)wall.spin(wp_nianlun)=-wall.spin(wp_nianlun)timerecord=mech.agefilename=string.build("jieguo%1",file_count)commandsave @filenameendcommandfile_count =1endifendhistory deletehistory id 1 @time_nowhistory id 2 @weiyiXhistory id 3 @forceLunYhistory id 4 @forceLunXhistory id 5 @measure_poros1history id 6 @measure_poros2history id 7 @measure_poros3history id 8 @measure_poros4history id 9 @measure_avsolve time [cishu*onetime]
碾压开始前状态为：
图中红色圈圈为测量圆的位置。
第一次碾压：
第二次碾压：
第三次碾压：
第四次碾压
可以看出随着碾压次数的提高，位移区域是在增加的。从位移场中也可以看出这个颗粒位移的趋势。
显示力链图，可以看到碾轮向地基的力链延伸。
这里显示记录的孔隙率的变化：
四个孔隙率随着碾压的进行都在减小，并且碾压的越多，后面减小的量越小，也是比较常见的指数曲线。
这里受限于计算效率，只算了五次，读者可以算十次或者二十次看看效果。
有的碾轮需要加激振力，只需要不停的改变碾轮伺服力的大小即可。
