 |
|
|
|
Coarse graining, as described on the previous page, has the disadvantage of smoothing out sharp boundaries. A solution to this is to calculate macroscopic quantities through a precise calculation of the overlap volumes between spheres and mesh elements, as described below.
MercuryCG offers this capability through the use of Severin Strobl's overlap library (https://github.com/severinstrobl/overlap), which implements the mathematical details of the precise overlap calculations outlined by Strobl et al. [1]. The overlap library is included in Mercury as a submodule and has to be activated by setting the CMake option MercuryDPM_OVERLAP to ON. Note that this option requires the CMake option MercuryDPM_EIGEN to be set to ON aswell, as the overlap library itself depends on the Eigen library for its calculations.
Once compiled correctly, macroscopic quantities based on exact overlap calculations can be selected by choosing the function ExactOverlap and the fields StandardFieldsBinning in MercuryCG. An example command line call could look like this:
More information on the syntax of MercuryCG's command line interface can be found in MercuryCG's documentation.
The binning process utilizes a grid of cubic cells. Given the minimum and maximum coordinates \((\mathbf{r}_\text{min} = (x_\text{min},\, y_\text{min},\, z_\text{min})^\text{T})\) and \((\mathbf{r}_\text{max} = (x_\text{max},\, y_\text{max},\, z_\text{max})^\text{T})\) of the considered space, along with the number of grid cells \((\mathbf{n} = (n_x,\, n_y,\, n_z)^\text{T})\), the vector \((\mathbf{d})\) representing the side lengths of the cubic cells in each dimension is determined as
\[ \mathbf{d} = \left(\frac{x_\text{max} - x_\text{min}}{n_x}, \frac{y_\text{max} - y_\text{min}}{n_y}, \frac{z_\text{max} - z_\text{min}}{n_z}\right)^\text{T}. \]
The center of each grid cell is calculated by
\[ \mathbf{r}_{ijk} = \mathbf{r}_\text{min} + \begin{bmatrix} (i + 0.5)d_x \\ (j + 0.5)d_y \\ (k + 0.5)d_z \end{bmatrix}, \]
where (d_x), (d_y), and (d_z) are the components of (\mathbf{d}), ensuring that the cells do not overlap and the space is fully covered by the grid. Here, the indices (i), (j), and (k) range from (0) to (n_x - 1), (0) to (n_y - 1), and (0) to (n_z - 1), respectively.
Given a mesh element or grid cell \(c\) of volume \(V_c\), the exact density \(j_c^q\) of a macroscopic quantity \(q_c\) within the space of the selected grid cell is calculated as
\[ j_c^q = \sum_{i=1}^{N} \frac{V_{c\cap i}}{V_cV_i} q_i. \]
Here, \(N\) is the number of particles in the simulation, \(q_i\) is the quantity of interest of particle \(i\), \(V_i\) is its volume \(i\), and \(V_{c \cap i}\) is the volume of the intersection of particle \(i\) with cell \(c\).
Figure 1 illustrates the intersection of one particle with a grid cell. The cell is colored in dark grey (its area corresponds to \(V_c\) in 3D); the particle is colored in light grey (its area corresponds to \(V_\text{i}\) in 3D). Their overlapping area is indicated by the hatched pattern (corresponding to \(V_{c\cap i}\) in 3D).
Given this method, the momentum density \(j_c^p\) is calculated as
\[ \mathbf{j}_c^p = \sum_{i=1}^{N} \frac{V_{c\cap i}}{V_cV_i} m_i \mathbf{v}_i, \]
where \(m_i\) is the mass and \(\mathbf{v}_i\) is the velocity of particle \(i\).
Thus, the velocity \(\mathbf{v}_c\) of cell \(c\) is calculated as
\[ \mathbf{v}_c = \frac{\mathbf{j}_c^p}{\rho_c}, \]
where
\[ \rho_c = \sum_{i=1}^{N} \frac{V_{c\cap i}}{V_cV_i} m_i \]
is the mass density of cell \(c\).
Using the exact binning method, interaction properties are also evaluated in a particle-centered way. For example, the contact force density \(\mathbf{j}_c^f\) is calculated as
\[ \mathbf{j}_c^f = \sum_{i=1}^{N} \frac{V_{c\cap i}}{V_cV_i} \sum_{j} \mathbf{f}_{ij}, \]
where the second sum runs over all particles \(j\) interacting with particle \(i\) and \(\mathbf{f}_{ij}\) is the contact force acting on particle \(i\).
In a similar fashion, the contact stress density is calculated as
\[ \mathbf{\sigma}_{c} = \sum_{i=1}^{N} \frac{V_{c\cap i}}{V_cV_i} \sum_{j} \mathbf{f}_{ij} \otimes \mathbf{l}_{ij}, \]
where \(\mathbf{l}_{ij}\) is the branch vector pointing from the center of particle \(i\) to the contact point with particle \(j\).
1.9.1