一个构造数值模拟试验的完整流程包括两个步骤:

  1. 生成初始模型。生成颗粒集合体,定义颗粒的材料参数,让其在重力作用下沉积,形成初始模型。
  2. 挤压。给定墙相应的速度,开始挤压。

这里是一个最简单的挤压计算实例,仅需30行命令。学会了该命令脚本,基本掌握VBOX的使用方法。

微观参数表(李长圣,2019)

颗粒直径         
d (m)
颗粒密度         
ρ (kg∙m-3)
颗粒摩擦系数         
缩短速率         
υ (m∙s-1)
时间步         
s
60,80 2500 0.3 2.0 0.05

李长圣 (2019) 基于离散元的褶皱冲断带构造变形定量分析与模拟. 博士论文. 南京大学. 推荐下载 最新修订版 提取码 zdem

登陆计算平台,运行 zdem push.py ,其中 push.py 中内容如下:

######################################
# title: 一个实例学会ZDEM
# date: 2021-04-26
# authors: 李长圣
# E-mail: sheng0619@163.com
# note:
# 括号内参数可根据模型大小及个人需要修改
# 脚本命令不区分大小写
# 用16个核心,实际用时<1小时
# 计算费用约<2元
# more info, see www.geovbox.com
#######################################
#程序初始化
START
#颗粒设为球,计算颗粒体积用4/3*pi*r^3计算
SET disk off
#设置研究范围 
BOX left 0.0 right 41000.0 bottom 0.0 height 11000.0 kn=0e10 ks=0e10 fric 0.00 
#设置挡板墙,这里模型采用hertz接触模型,挡板墙的kn ks无效,计算时取颗粒的参数
WALL ID 0, NODES (      0.0 ,     10.0 ) (  40000.0 ,     10.0 ), kn=0e10 ks=0e10 fric 0.0 COLOR black
WALL ID 1, NODES (     10.0 ,  10000.0 ) (     10.0 ,     10.0 ), kn=0e10 ks=0e10 fric 0.0 COLOR blue
WALL ID 2, NODES (  40000.0 ,     10.0 ) (  40000.0 ,  10000.0 ), kn=0e10 ks=0e10 fric 0.0 COLOR red
#在矩形范围内生成颗粒
GEN NUM 100000 rad discrete 60.0 80.0,  x ( 10.0, 40000.0), y ( 10.0, 10000.0), COLOR black GROUP ball_rand
#设置颗粒的微观参数 density 密度,firc 摩擦系数,shear 剪切模量,poiss 泊松比,damp 局部阻尼常数,heart Hertz-Mindlin接触模型
PROP density 2.5e3, fric 0.0, shear 2.9e9, poiss 0.2, damp 0.4, hertz
#设置时间步及重力加速度
SET  DT 5e-2,  GRAVITY  0.0,  -9.8 
#设置每1000步保存一次vtk格式的计算结果
SET  vtk 1000
#设置每1000步保存一次ps格式的计算结果
SET  ps 1000
#设置每1000步保存一次dat格式的计算结果
SET  print 1000
#沉积,计算5000步
CYC 5000
#删除4000米以上的颗粒
DEL RANGE y 4000.0 999000.0
#平衡,计算1000步
CYC 1000
#输出包含颗粒的[x y r]信息的初始模型 init_xyr.dat
#EXP init_xyr.dat

#设置bond粘结,使颗粒具有粘聚力,ebmod 杨氏模量, gbmod 剪切模量,tstrength 抗拉强度,sstrength 聚合强度, firc 摩擦系数
PROP ebmod 2e8 gbmod 2e8  tstrength 2e7 sstrength 4e7 fric 0.3 
#给地层赋上颜色
PROP COLOR lg          range y    0.0   500.0
PROP COLOR green       range y  500.0  1000.0
PROP COLOR yellow      range y 1000.0  1500.0
PROP COLOR red         range y 1500.0  2000.0
PROP COLOR black       range y 2000.0  2500.0
PROP COLOR mg          range y 2500.0  3000.0
PROP COLOR blue        range y 3000.0  3500.0
PROP COLOR gb          range y 3500.0  4000.0
PROP COLOR violet      range y 4000.0  4500.0

#设置挡板墙摩擦系数
WALL id 0 fric 0.3
WALL id 1 fric 0.3
WALL id 2 fric 0.3
#设置墙的挤压速度 x方向速度为2.0
WALL id 1 xv 2.0
#设置墙的挤压量x方向推进10000.0,每挤压2000.0保存一次计算结果
IMPLE wall id 1 xmove 10000.0 save 2000.0 print 1000.0 ps 1000.0 vtk 1000.0
#计算停止
STOP

单位 (km)

构造变形及应力应变的演化过程

挤压量10 km时,构造变形及应力应变