Overview
EDEMpy is a Python library for post-processing and analyzing EDEM simulation data. It leverages the open-source hdf5 file format, which has been used to store EDEM data since EDEM 2017. With EDEMpy, extracting specific data from a simulation deck and processing it in a customizable and reusable manner becomes effortless. For instance, users can extract the force acting on a specific geometry over time, compare results across multiple simulation decks, track particles' residence time throughout a simulation, or calculate post-processing properties like the segregation index using EDEM data.
Since the release of version 1.1.0, users can now write data to the simulation without relying on the EDEM graphical user interface (GUI). This allows them to automate the process of creating an EDEM simulation and makes EDEMpy more powerful than ever before.
In this tutorial, you will learn how to utilize the new setting functions in EDEMpy to construct a complete simulation from scratch. We will cover topics such as creating materials, setting material properties and interactions, adding particles, defining geometries, creating factories, and configuring model parameters.
1. Importing the necessary libraries and opening the deck
The first step is to import the libraries that we'll be using throughout this tutorial: edempy and numpy. The edempy library enables us to work with EDEM simulation data, while numpy provides useful numerical computing functionalities. Once we have imported the required libraries, we need to open our simulation deck in write mode, denoted by mode='w', which allows us to modify and write data to the deck.
from edempy import Deck
import numpy as np
deck = Deck("Setting_Functions_Example.dem", mode='w')
2. Creating the bulk material
2.1 Setting the bulk material properties
The first step in our simulation is to create the bulk material. We can achieve this by utilizing the createMaterial function. To designate it as a bulk material, we need to set the last argument of the createMaterial function to True. In this case, we will name our bulk material as BulkMaterial, and set its properties as per the following table:
| Property | Value |
|---|---|
| Material Name | BulkMaterial |
| Poisson Ratio | 0.25 |
| Shear Modulus | 10 MPa |
| Density | 4100 kg/m³ |
| Work function | 0 J |
# To create a bulk material, the last argument must be set to True
deck.createMaterial('BulkMaterial', 0.25, 10_000_000, 4100, 0, True)
2.2 Setting, particle size and shape
With the bulk material created, we can now add a particle to it. This can be achieved using the createMultisphereParticleType function. To utilize this function, we need to specify the bulk material and define the particle name. In this case, we will name it as MyParticle. By default, the createMultisphereParticleType function creates a triple-sphere shape.
Once the particle is created, we can adjust the radius and position of the spheres that constitute our particle using the setMultisphereSpherePosition and setMultispherePhysicalRadius functions. Then, we use the setAutoCalculateParticleProperties function to ensure the accurate calculation of the particle's properties.
# Add a particle to the bulk material
deck.createMultisphereParticleType('BulkMaterial', particle_name = 'MyParticle')
# Set the position and radius of each sphere in MyParticle
deck.setMultisphereSpherePosition('MyParticle', 0, [0.0, 0.004, 0.0])
deck.setMultispherePhysicalRadius('MyParticle', 0, 0.005)
deck.setMultisphereSpherePosition('MyParticle', 1, [-0.00346, -0.002, 0.0])
deck.setMultispherePhysicalRadius('MyParticle', 1, 0.005)
deck.setMultisphereSpherePosition('MyParticle', 2, [0.00346, -0.002, 0.0])
deck.setMultispherePhysicalRadius('MyParticle', 2, 0.005)
# Enable auto-calculation of particle properties
deck.setAutoCalculateParticleProperties('MyParticle', True)
3. Defining the equipment material
After creating the bulk material and particles, the next step is to create the equipment material. We will name the equipment material as EquipMaterial. To create the equipment material, we will use the createMaterial function, but this time we need to set the last argument to False. The properties of EquipMaterial are shown in the table below:
| Property | Value |
|---|---|
| Material Name | EquipMaterial |
| Poisson Ratio | 0.3 |
| Shear Modulus | 70 GPa |
| Density | 7800 kg/m³ |
| Work function | 0 J |
# To create the equipment material, the last argument must be set to False
deck.createMaterial('EquipMaterial', 0.30, 70_000_000_000, 7800, 0, False)
3.1. Defining the interaction parameters
Now that all the materials in our simulation have been created we can set the interaction parameters between them. This can be achieved using the createInteraction function. The interaction parameters are specified in the table below:
| Contact parameter | Particle-Particle | Particle-Geometry |
|---|---|---|
| Coefficient of restitution | 0.36 | 0.40 |
| Coefficient of static friction | 0.64 | 0.45 |
| Coefficient of rolling friction | 0.13 | 0.15 |
# Define particle-particle interactions
deck.createInteraction('BulkMaterial', 'BulkMaterial', 0.36, 0.64, 0.13)
# Define particle-geometry interactions
deck.createInteraction('BulkMaterial', 'EquipMaterial', 0.40, 0.45, 0.15)
4. Defining the geometries
To create the geometries, we need to specify the coordinates of each node as well as the connectivity matrix. The connectivity matrix represents the relationships between the vertices (points) and the triangles in the mesh. Each row of the connectivity matrix corresponds to a triangle in the mesh.
In this example, we create two geometries: a virtual square named MySquare and a physical cube named MyCube. After specifying the coordinates and connectivity matrix for each geometry, we can use the createGeometry function to create them. Additionally, if you have an existing geometry, you can import it into your simulation using Python libraries like trimesh. These libraries provide access to the geometry's nodes and connectivity list, allowing you to import the geometry into your simulation seamlessly.
# Define the 4 corners of the square
square_coords = np.array([[-0.5, -0.5, 0.5],[0.5, -0.5, 0.5],[0.5, 0.5, 0.5],[-0.5, 0.5, 0.5]])
# Define the 2 triangles that make up the square
square_triangles = np.array([[0,3,1], [1,3,2]])
# Define the 8 corners of the cube
cube_coords= np.array([[-1, -1, -1], [1, -1, -1], [1, 1, -1], [-1, 1, -1], [-1, -1, 1], [1, -1, 1], [1, 1, 1], [-1, 1, 1]])
# Define the 12 triangles that make up the cube
cube_triangles= np.array([[0,3,1],[1,3,2],[0,4,7],[0,7,3],[4,5,6],[4,6,7],[5,1,2],[5,2,6],[2,3,6],[3,7,6],[0,1,5],[0,5,4]])
# Create the square
deck.createGeometry('MySquare', 'EquipMaterial', square_coords, square_triangles, geometry_type='Virtual')
# Create the cube
deck.createGeometry('MyCube', 'EquipMaterial', cube_coords, cube_triangles, geometry_type='Physical')
4.1. Defining geometry motion, particle factory and particle size distribution
Now that the geometries are created, we can add a linear rotation to MyCube and a particle factory to MySquare. To add a linear rotation to MyCube, we use the createKinematic function and specify the parameters listed in the table below:
| Name | Box_rotation |
| Type | Linear Rotation |
| Start Time | 2 s |
| Initial Velocity | 10 rpm |
| Axis Start Point | (0, 0, -1) |
| Axis End Point | (0, 0, 1) |
To add a particle factory to MySquare, we use the createFactory function and specify the parameters listed in the table below:
| Name | BulkMaterial_Factory |
| Material Name | BulkMaterial |
| Velocity | fixed |
| Fixed Velocity | (0, 0, -1) m/s |
| Start Time | 0 s |
| Total to Create | 30000 particles |
| Creation Rate | 15000 particles/s |
Finally, we can define the particle size distribution using the setParticleUserDefinedSizeDistribution function. Here, we specify the scale and the percentage mass in each size class. The last argument of the function is set to True so that the scales are applied to the particle radius.
# Set geometry motion - linear rotation
rotation_in_rpm = 10
rotation_in_rad_s = rotation_in_rpm * (2*np.pi/60)
deck.createKinematic('MyCube', kinematic_name= 'Cube_rotation', kinematic_type= 'Linear Rotation', start_time= 2.0,
rotation_init_angular_velocity_value= rotation_in_rad_s, rotation_start_point= [0.0, 0.0, -1.0],
rotation_end_point= [0.0, 0.0, 1.0])
# Create particle factory
deck.createFactory('MySquare', 'BulkMaterial_Factory', 'BulkMaterial', velocity= 'fixed', fixed_velocity=[0.0, 0.0, -1.0],
start_time= 0, is_mass_creation= False, total_to_create= 30_000, is_mass_generation= False,
creation_rate= 15_000)
# Define the particle size distribution
scales = np.array([1.0, 1.2, 1.5, 1.7])
mass_percentages = np.array([45, 30, 5, 20])
deck.setParticleUserDefinedSizeDistribution('MyParticle', scales, mass_percentages, is_scaled_by_radius= True)
5. Defining the physics model
Finally, we need to define the physics models that we'll use in our simulation. For particle-particle contacts, we'll use the Hertz-Mindlin model with JKR cohesion. Since the default Hertz-Mindlin model is already present, we first need to delete it using the deletePhysicsModelSurfSurf function. Then, we specify the surface energy value for the particle-particle contacts and add the model to our simulation using the createPhysicsModelSurfSurf function. For the particle-geometry contacts, we'll keep the default Hertz-Mindlin model and add the Relative Wear model. Since the Relative Wear model doesn't require any additional parameters, we can simply add it to our simulation using the createPhysicsModelSurfGeom function.
# Delete the default Hertz-Mindling model
deck.deletePhysicsModelSurfSurf('Hertz-Mindlin (no slip)')
# Define the JKR surface energy
model_data=np.array([('BulkMaterial', 'BulkMaterial', 2)], dtype=([('name1', 'S256'), ('name2', 'S256'), ('value', '<f8')]))
# Add the particle-particle contact model
deck.createPhysicsModelSurfSurf('Hertz-Mindlin with JKR',{"data": model_data})
# Add the particle-geometry contact model
deck.createPhysicsModelSurfGeom('Relative Wear',{})
6. Setting the environment and saving the deck
In the final step of our simulation setup, we need to define the minimum and maximum values of the simulation domain. This can be done using the setDomainMin and setDomainMax functions. The domain represents the bounding box within which the simulation will take place. We also need to save our work using the save function to store the simulation deck for future use.
# Set the minimum values of the simulation domain
deck.setDomainMin([-1.01, -1.01, -1.01])
# Set the maximum values of the simulation domain
deck.setDomainMax([1.01, 1.01, 1.01])
# Save the simulation deck
deck.save()Pre-Requisite
EDEM 2023.1 or later installed
Usage/Installation Instructions
Download the attached zip file and extract all files. Open EDEM and go to the Analyst tab. Click on File, Run EDEMpy script and select the script.
Note: the deck and the script should be in the same folder.