UCS Test with Voronoi Blocks
Problem Description
Note
The project file for this example may be viewed/run in 3DEC.[1] The data files used are shown at the end of this example.
This example describes how to perform an unconfined compression test on a cylindrical sample made up of bonded Voronoi blocks.
3DEC Model
The model is created by first importing a simple DXF of a parallelepiped. The shape is 6x6x12.5 cm. A voronoi assembly of blocks is then generated using the block generate voronoi
command with the fromgeometry keyword. A maximum edge length of 6 mm is specified, resulting in 3361 Voronoi blocks.
The model is then cut into a cylinder by importing a dxf of a cylinder and using the command block cut geometry
. This step is necessary because the Voronoi blocks along the edges of the parallelepiped exhibit unrealistic alignment, leading to increased strength resulting from the “crystalline” packing. By cutting away the edges, this alignment is disrupted.
The blocks are then zoned and the resulting 3DEC model is shown in Figure 1.
The zones are assumed to be elastic with a Young’s modulus of 76 GPa and a Poisson’s ratio of 0.25. The properties were chosen to mimic the numerical experiments described in Lan et al (2013). The contacts between the blocks were also assigned properties from Lan et al (2013) as shown in the following properties:


Normal stiffness 
7.6e4 GPa/m 
Shear stiffness 
7.6e4 GPa/m 
Cohesion 
40 MPa 
Friction 
10 degrees 
Tension 
14.4 MPa 
The block is loaded by applying a constant velocity of 0.05 m/s to the gridpoints on the top and bottom of the block. Axial stress is calculated by summing the vertical (z) reaction forces of all the top and bottom gridpoints and dividing by twice the crosssectional area. The axial strain is simply the sum of displacement of the top and bottom gridpoints divided by the original length.
The model is then cycled for 35,000 steps. The model solve
command is not used because the model will never reach equilibrium.
Results and Discussion
The axial stress versus axial strain curve from the 3DEC simulation is shown in Figure 2. A peak stress of ~752 MPa is observed, which is close to the value obtained in Lan et al (2013).
A crosssection of the “cracks” that form are shown in Figure 3. This figure is generated by plotting the joints, contouring by normal displacement and then checking the box to specify an aperture threshold of 1e5 m, so that only the open joints are plotted.
References
Lan, H., Martin, C.D., Qi, S. and Li, T., 2013. A 3D grain based model for characterizing the geometric heterogeneity of brittle rock, in Proceedings of the 47th U.S. Rock Mechanics / Geomechanics Symposium, San Francisco, CA, USA, 2326 June, American Rock Mechanics Association.
Data Files
sample.dat
model new
model random 10000
block tolerance 1e6
; make parallelepiped to fill with voronois
geometry import 'parallelepiped.dxf'
block generate voronoi maximumedge 0.006
; cut out cylinder to prevent alignment of voronois along the edges
model domain extent 1 1
geometry import 'cylinder.dxf'
block cut geometry 'cylinder'
block delete range cylinder end1 0 0 0 end2 0 0 0.125 radius 0.025 not
model save 'voro_rigid.sav'
;model restore 'voro_rigid.sav'
block zone generatenew maximumedge 0.002
; delete blocks that didn't zone
block delete range deformable not
model save 'voro_defor.sav'
ucs_test.dat
model restore 'voro_defor'
model largestrain on
; properties in MPa
block zone cmodel assign elastic
block zone property young 76e3 pois 0.25 density 0.002
block contact jmodel mohr
block contact property stiffnessnormal 7.6e7 stiffnessshear 7.6e7 ...
cohesion 40 tension 14.4 friction 10
block contact materialtable default property stiffnessnormal 7.6e7 ...
stiffnessshear 7.6e7 friction 10
block gridpoint group 'top' range positionz 0.125
block gridpoint group 'bottom' range positionz 0.
block gridpoint apply velocityz 0.05 range group 'top'
block gridpoint apply velocityz 0.05 range group 'bottom'
model cycle 1
; make lists of gridpoints forigin monitoring reaction forces
[global area = math.pi*0.025^2]
[global length = 0.125]
fish define findList
gp_list_top = list
gp_list_bottom = list
loop foreach gp block.gp.list
if block.gp.isgroup(gp,'top')
gp_list_top("end") = gp
endif
if block.gp.isgroup(gp,'bottom')
gp_list_bottom("end") = gp
endif
endloop
end
[findList]
[global gp_top = gp_list_top(1)]
fish define zzstress
; compression positive,unit MPa
local top_reaction = list.sum(block.gp.force.reaction.z(::gp_list_top))
local bottom_reaction = list.sum(block.gp.force.reaction.z(::gp_list_bottom))
zzstress = 0.5*(top_reaction + bottom_reaction)/area
zzstrain = 2.0*block.gp.disp.z(gp_top)/length
end
fish history zzstress
fish history zzstrain
model cycle 30000
model save 'ucs.sav'
Endnote
Was this helpful? ...  Itasca Software © 2024, Itasca  Updated: Sep 26, 2024 