Data mining

Read

Local data

All data of the simulation are accessible from python; when you open the Inspector, blue labels of various data can be clicked – left button for getting to the documentation, middle click to copy the name of the object (use Ctrl-V or middle-click to paste elsewhere). The interesting objects are among others (see Omega for a full list):

  1. O.engines

    Engines are accessed by their index (position) in the simulation loop:

    O.engines[0]      # first engine
    O.engines[-1]     # last engine
    

    Note

    The index can change if O.engines is modified. Labeling introduced in the section below is a better solution for reliable access to a particular engine.

  2. O.bodies

    Bodies are identified by their id, which is guaranteed to not change during the whole simulation:

    O.bodies[0]                                                         # first body
    [b.shape.radius for b in O.bodies if isinstance(b.shape,Sphere)]    # list of radii of all spherical bodies
    sum([b.state.mass for b in O.bodies])                               # sum of masses of all bodies
    numpy.average([b.state.vel[0] for b in O.bodies])                   # average velocity in x direction
    

    Note

    Uniqueness of Body.id is not guaranteed, since newly created bodies might recycle ids of deleted ones.

  3. O.forces

    Generalized forces (forces, torques) acting on each particle. They are (usually) reset at the beginning of each step with ForceResetter, subsequently forces from individual interactions are accumulated in InteractionLoop. To access the data, use:

    O.forces.f(0)     # force on #0
    O.forces.t(1)     # torque on #1
    
  4. O.interactions

    Interactions are identified by ids of the respective interacting particles (they are created and deleted automatically during the simulation):

    O.interactions[0,1]   # interactions of #0 with #1
    O.interactions[1,0]   # the same object
    O.bodies[0].intrs()   # all interactions of body #0
    for i in O.bodies[12].intrs(): print (i.isReal,i.id1,i.id2)    # get some info about interactions of body #12
    [(i.isReal,i.id1,i.id2) for i in O.bodies[12].intrs()]         # same thing, but make a list
    

Labels

Engines and functors can be labeled, which means that python variable of that name is automatically created.

Yade [1]: O.engines=[
   ...:    NewtonIntegrator(damping=.2,label='newtonCustomLabel')
   ...: ]
   ...: 

Yade [2]: newtonCustomLabel.damping=.4

Yade [3]: O.engines[0].damping              # O.engines[0] and newtonCustomLabel are the same objects
Out[3]: 0.4

Yade [4]: newtonCustomLabel==O.engines[0]   # O.engines[0] and newtonCustomLabel are the same objects
Out[4]: True

Exercises

  1. Find meaning of this expression:

    max([b.state.vel.norm() for b in O.bodies])
    
  2. Run the Gravity deposition script, pause after a few seconds of simulation. Write expressions that compute

    1. kinetic energy \(\sum \frac{1}{2} m_i |v_i| ^2\)
    2. average mass (hint: use numpy.average)
    3. maximum \(z\)-coordinate of all particles
    4. number of interactions of body #1

Global data

Useful measures of what happens in the simulation globally:

unbalanced force
ratio of maximum contact force and maximum per-body force; measure of staticity, computed with unbalancedForce.
porosity
ratio of void volume and total volume; computed with porosity.
coordination number
average number of interactions per particle, avgNumInteractions
stress tensor (periodic boundary conditions)
averaged force in interactions, computed with normalShearStressTensors
fabric tensor
distribution of contacts in space (not yet implemented); can be visualized with plotDirections

Energies

Evaluating energy data for all components in the simulation (such as gravity work, kinetic energy, plastic dissipation, damping dissipation) can be enabled with

O.trackEnergy=True

Subsequently, energy values are accessible in the O.energy; it is a dictionary where its entries can be retrived with keys() and their values with O.energy[key].

Save

PyRunner

To save data that we just learned to access, we need to call Python from within the simulation loop. PyRunner is created just for that; it inherits periodicy control from PeriodicEngine and takes the code to run as text (must be quoted, i.e. inside '...') attributed called command. For instance, adding this to O.engines will print the current step number every one second wall clock time:

O.engines=O.engines+[ PyRunner(command='print(O.iter)',realPeriod=1) ]

Writing complicated code inside command is awkward; in such case, we define a function that will be called:

def myFunction():
        '''Print step number, and pause the simulation is unbalanced force is smaller than 0.05.'''
        print(O.iter)
        if utils.unbalancedForce()<0.05:
                print('Unbalanced force is smaller than 0.05, pausing.')
                O.pause()

Now this function can be added to O.engines:

O.engines+=[PyRunner(command='myFunction()',iterPeriod=100)]

or, in general, like that:

O.engines=[
        # ...
        PyRunner(command='myFunction()',iterPeriod=100) # call myFunction every 100 steps
]

Warning

If a function was declared inside a live yade session (ipython) then an error NameError: name 'myFunction' is not defined will occur unless python globals() are updated with command

globals().update(locals())

Exercises

  1. Run the Gravity deposition simulation, but change it such that:
    1. utils.unbalancedForce is printed every 2 seconds.
    2. check every 1000 steps the value of unbalanced force
      • if smaller than 0.2, set damping to 0.8 (hint: use labels)
      • if smaller than 0.1, pause the simulation

Keeping history

Yade provides the plot module used for storing and plotting variables (plotting itself will be discussed later). Let us start by importing this module and declare variable names that will be plotted:

from yade import plot
plot.plots={'t':('coordNum','unForce',None,'Ek')}                # kinetic energy will have legend on the right as indicated by None separator.

Periodic storing of data is done with PyRunner and the plot.addData function. Also let’s enable energy tracking:

O.trackEnergy=True
def addPlotData():
        # this function adds current values to the history of data, under the names specified
        plot.addData(t=O.time,Ek=utils.kineticEnergy(),coordNum=utils.avgNumInteractions(),unForce=utils.unbalancedForce())

Now this function can be added to O.engines:

O.engines+=[PyRunner(command='addPlotData()',iterPeriod=20)]

or, in general, like that:

O.engines=[  # ...,
        PyRunner(command='addPlotData()',iterPeriod=20)         # call the addPlotData function every 20 iterations
]

History is stored in plot.data, and can be accessed using the variable name, e.g. plot.data['Ek'], and saved to text file (for post-processing outside yade) with plot.saveDataTxt.

Plot

plot provides facilities for plotting history saved with plot.addData as 2d plots. Data to be plotted are specified using dictionary plot.plots

plot.plots={'t':('coordNum','unForce',None,'Ek')}

History of all values is given as the name used for plot.addData; keys of the dictionary are \(x\)-axis values, and values are sequence of data on the \(y\) axis; the None separates data on the left and right axes (they are scaled independently). The plot itself is created with

plot.plot()         # on the command line, F8 can be used as shorthand

While the plot is open, it will be updated periodically, so that simulation evolution can be seen in real-time.

Energy plots

Plotting all energy contributions would be difficult, since names of all energies might not be known in advance. Fortunately, there is a way to handle that in Yade. It consists in two parts:

  1. plot.addData is given all the energies that are currently defined:

    plot.addData(i=O.iter,total=O.energy.total(),**O.energy)
    

    The O.energy.total functions, which sums all energies together. The **O.energy is special python syntax for converting dictionary (remember that O.energy is a dictionary) to named functions arguments, so that the following two commands are identical:

    function(a=3,b=34)              # give arguments as arguments
    function(**{'a':3,'b':34})      # create arguments from dictionary
    
  2. Data to plot are specified using a function that gives names of data to plot, rather than providing the data names directly:

    plot.plots={'i':['total']+O.energy.keys()}
    

    where total is the name we gave to O.energy.total() above, while O.energy.keys() will always return list of currently defined energies.

Energy plot example

Plotting energies inside a live yade session, for example by launching examples/test/triax-basic-without-plots.py would look following:

from yade import plot
O.trackEnergy=True
O.step()                          # performing a single simulation step is necessary to populate O.energy.keys()
plot.plots={'t':O.energy.keys()+['total']}

def addPlotData():
        # this function adds current values to the history of data, under the names specified
        plot.addData( t=O.time , total=O.energy.total() , **O.energy )

O.engines+=[PyRunner(command='addPlotData()',iterPeriod=20)]

globals().update(locals())        # do this only because this is an example of a live yade session

Press F8 to show plot window and F11 to show 3D view, then press ▶ to start simulation.

Using multiple plots

It is also possible to make several separate plots, for example like this:

plot.plots={ 't':('total','kinetic') , 't ':['elastPotential','gravWork'] , 't  ':('nonviscDamp') }

Warning

There cannot be duplicate names declared in separate plots. This is why spaces were used above to indicate the same variable t.

With the caveat above, a following example inside a live yade session launched on examples/test/triax-basic-without-plots.py would look following:

from yade import plot
O.trackEnergy=True
plot.plots={ 't':('total','kinetic') , 't ':['elastPotential','gravWork'] , 't  ':('nonviscDamp') }

def addPlotData():
        # assign value to all three: 't', 't ' and 't  ' with single t=... assignment
        plot.addData( t=O.time , total=O.energy.total() , **O.energy )

O.engines+=[PyRunner(command='addPlotData()',iterPeriod=20)]

globals().update(locals())        # do this only because this is an example of a live yade session

plot.plot(subPlots=False)         # show plots in separate windows

plot.plot(subPlots=True)          # same as pressing F8: close current plot windows and reopen a single new one

Press F8 to show plot window and F11 to show 3D view, then press ▶ to start simulation, see video below:

Exercises

  1. Calculate average momentum in y direction.
  2. Run the Gravity deposition script, plotting unbalanced force and kinetic energy.
  3. While the script is running, try changing the NewtonIntegrator.damping parameter (do it from both Inspector and from the command-line). What influence does it have on the evolution of unbalanced force and kinetic energy?
  4. Think about and write down all energy sources (input); write down also all energy sinks (dissipation).
  5. Simulate Gravity deposition and plot all energies as they evolve during the simulation.

See also

Most Examples with tutorial use plotting facilities of Yade, some of them also track energy of the simulation.