yade-users team mailing list archive

[Anton <question706664@xxxxxxxxxxxxxxxxxxxxx>]

Question #706664 on Yade changed:
https://answers.launchpad.net/yade/+question/706664

    Status: Needs information => Open

Anton gave more information on the question:
Hello Jan and thanks for coming back to me.

Here is the code I am using right now:

# gravity deposition in box, showing how to plot and save history of data,
# and how to control the simulation while it is running by calling
# python functions from within the simulation loop

# import yade modules that we will use below
from yade import pack, plot
from yade import utils, export
from minieigen import Vector3

# create rectangular box from facets
O.bodies.append(geom.facetBox((.5, .5, .5), (.5, .5, .5), wallMask=31))

# create empty sphere packing
# sphere packing is not equivalent to particles in simulation, it contains only the pure geometry
sp = pack.SpherePack()
# generate randomly spheres with uniform radius distribution
sp.makeCloud((0, 0, 0), (1, 1, 1), rMean=.05, rRelFuzz=.5)
# add the sphere pack to the simulation
sp.toSimulation()

O.engines = [
        ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb(), Bo1_Facet_Aabb()]),
        InteractionLoop(
                # handle sphere+sphere and facet+sphere collisions
                [Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom()],
                [Ip2_FrictMat_FrictMat_FrictPhys()],
                [Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
        NewtonIntegrator(gravity=(0, 0, -9.81), damping=0.4),
        # call the checkUnbalanced function (defined below) every 2 seconds
        PyRunner(command='checkUnbalanced()', realPeriod=2),
        # call the addPlotData function every 200 steps
        PyRunner(command='addPlotData()', iterPeriod=100),
        VTKRecorder(iterPeriod=2000, recorders=['all'], fileName='/home/antonsie/DEM/Test/M01/-')
]
O.dt = .5 * PWaveTimeStep()

# enable energy tracking; any simulation parts supporting it
# can create and update arbitrary energy types, which can be
# accessed as O.energy['energyName'] subsequently
O.trackEnergy = True

# if the unbalanced forces goes below .05, the packing
# is considered stabilized, therefore we stop collected
# data history and stop
def checkUnbalanced():
	if unbalancedForce() < .05:
		O.pause()
		plot.saveDataTxt('bbb.txt.bz2')
		# plot.saveGnuplot('bbb') is also possible


# collect history of data which will be plotted
def addPlotData():
	# each item is given a names, by which it can be the unsed in plot.plots
	# the **O.energy converts dictionary-like O.energy to plot.addData arguments
	plot.addData(i=O.iter, unbalanced=unbalancedForce(), **O.energy)


# define how to plot data: 'i' (step number) on the x-axis, unbalanced force
# on the left y-axis, all energies on the right y-axis
# (O.energy.keys is function which will be called to get all defined energies)
# None separates left and right y-axis
plot.plots = {'i': ('unbalanced', None, O.energy.keys)}

# show the plot on the screen, and update while the simulation runs
plot.plot()

O.saveTmp()


# overall porosity
pa = utils.voxelPorosity(200,(0, 0, 0),(1, 1, 1))
pb = utils.porosity()
print('porosity of whole sample computed with voxel porosity', pa)
print('porosity of whole sample computed with porosity', pb)

# porosity from sub volumes (8 total)
p1 = utils.voxelPorosity(200,(0, 0, 0),(0.5, 0.5, 0.5))
p2 = utils.voxelPorosity(200,(0, 0.5, 0),(0.5, 1, 0.5))
p3 = utils.voxelPorosity(200,(0.5, 0, 0),(1, 0.5, 0.5))
p4 = utils.voxelPorosity(200,(0.5, 0.5, 0),(1, 1, 0.5))
p5 = utils.voxelPorosity(200,(0, 0, 0.5),(0.5, 0.5, 1))
p6 = utils.voxelPorosity(200,(0.5, 0, 0.5),(1, 0.5, 1))
p7 = utils.voxelPorosity(200,(0, 0.5, 0.5),(0.5, 1, 1))
p8 = utils.voxelPorosity(200,(0.5, 0.5, 0.5),(1, 1, 1))

# print porosity obtained by subvolumina
print(p1, p2, p3, p4, p5, p6, p7, p8)

#vtk = yade.export.VTKExporter('/home/antonsie/DEM/Test/M01/pa-')
#vtk.export('pa')

The code is working fine except for the last line. It returns the
porosity calculated with utils.voxelPorosity() for the whole sample.
>From the documentation I understood, that it is also calculating the
porosity for various subvolumina (in this case 200), but it is not clear
to me, how I can access this data. Can you help me out here?

So I turned to defining subvolumina in the code which I find to be a
tedious process. My resort would be to save them as a txt document and
convert it to vtk.

The VTKRecorder can't return the values I would need to calculate the porosity, which are the volumes of the subvolumina defined through utils.voxelPorosity() so I tried to use the VTKExporter.
The yade documentation says, with VTKExporter, the user is able to export any user defined variable. However, it is not clear to me, how to do that. I searched in the documentation and in examples on GitLab. Could you provide more infomation?

Thanks in advance!

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.