This contribution presents Pysammos, a software that can be used to post-process the results of any Discrete Element Modelizations (DEM). Indeed, the discrete particle data obtained with the DEM can be transformed into continuous fields, more accessible for visualization purposes and deeper analysis. The architecture, the abilities, the validity, and the performances of Pysammos are depicted in this paper. This method appears user-friendly, computationally efficient, and flexible to the DEM software used.

The article is well written and well structured.

Even if I am positive about this article, I have a few queries to consider during this minor revision.

L. 55: I would change “force displacement law” to “contact law”.

L. 116: There is no verb in “Similarly, to enhance DEM analysis by extracting continuum fields without the need to handle inner-level source code”. Could you reformulate this sentence to ease the reading?

L. 149: I would not refer to Eq. 9 solely, as you use multiple functions for the coarse-graining function \Psi. Is a reference to Section 3.4.3 more accurate?

Eq. 5: Could you check this Equation? In the current version, an individual contact is computed twice. I guess that the second sum should start at j=i+1, in agreement with the formulation depicted in Weinhart et al. (2012).

L. 170: It appears in L. 385 that you made a different choice than other D2C software (MercuryDPM, for instance) for the particle search (particle vs contact point). A justification of this assumption can be relevant here.

L. 264: I would start with the Section “Data reading and handling” as the format of the input data is described. Please consider reorganizing your Sections.

Figure 6: Could you check this plot? Indeed, the Heavyside function does not verify the formulation detailed in Eq. 11. The function should be =0 for w<r, while the plot shows a non-null value for w<r<2w. This verification should be conducted for the Lucy and Gaussian functions.

Figure 6: I guess that W(r) should be Psi(r). If not, could you specify the definition of W(r)?

Figure B2: Same remark as Figure 6: I guess that W(r) should be Psi(r). If not, could you specify the definition of W(r)?

L. 529: This depth-independence of the vertical pressure is controversial. I agree with the Authors that Fig. 13 does not show a clear increase in the pressure. However, the relation P=rho.g.h is well established in the literature. I guess this relation is retrieved for the configuration depicted if a larger height h or a larger density rho is considered. Could you add limitations to the conclusion raised concerning the homogeneous vertical pressure with depth?

Figure 12: Could you consider modifying the scale used for the right map? The color assigned to the value -6.6 does not seem to be employed in the map. The range used for the scale can be shorter.

L. 574: As discussed with Section 4.4.3, the default value of w (=0.75*d_43) may not be the accurate choice to make. Maybe a default value of w that depends on the mean distance of the particles is a more resilient selection. In particular, if this code is used by a standard user (see the definition of the profile in Section 3.1.3). This remark is more a hint to the Authors, the text of the article does not have to be changed. I have not tried; maybe it is a bad idea.

L. 683: The consideration of the particle-wall interactions seems in Equation 5 to be already conducted (with the particles between N+1 and N+Nb). If it is finally not the case, could you delete this part of the Equation?

This paper describe a coarse-graining software, Pysammos, dedicated to the spatial averaging of discrete element modelling data, in particular the one produces by the software MFiX-DEM. The principle of the method and software are explained, as well as its specificities, and its performance are compared with other coarse-graining codes and with different geometries and simulation conditions.

The paper is overall well written and detailed, and can be useful as an entry point to coarse-graining, while proposing additional features not often present in already existing solutions. Some of the specific aspects, such as the phase identification and the validation with non-spherical particles, can be useful in a broad range of contexts. I only have some minor comments and questions below.

1/ In the simulations of particles suspended in a magma, it is not clear how the flow is created, and if it is a one-way or two-ways interaction between the fluid and particles. It would be interesting as well to provide the volume fraction field, which may have some spatial variations. Also, the simulations are indicated as 2D, but the particles seems to be three-dimensional?

2/ The results on the stress under the pile (4.4.4) are quite surprising. Experimental observations show essentially hydrostatic pressure, sometimes with a dip of pressure below the center, but both appear quite different from the results reported. Even if frictional redistribution of the stress occurs (similar to the Janssen effect), the weight of the grains still needs to be supported somehow. It would be good to also look at the shear stresses to assess whether this frictional redistribution does occur. Also, Fig 13a may be better plotted with a linear scale, since again the expectation would be a linear stress increase with depth on average.

3/ The calculation of the contact stress involves integrating the windowing function over the branch vector, as shown in equation 5. Is this approximated through averaging the values at the center, as is commonly done, or is the integral numerically solved? How? This may explain some of the differences pointed out in σ_c later on in the paper.

4/ There are many different clustering algorithm that exist, how was the particular one used chosen?

5/ line 144: the mathematical detail of the Dirac delta function is inaccurate. Technically, the Dirac function is infinity at the particle centre, and 0 otherwise, such as its integral over the whole space is 1.

6/ Equations 5 and C4: the tensor product, or Kronecker product that is used for the vector in those equations, should be specified. It is usually denoted with ⊗.

7/ line 193: the use of only equilibrium states to justify the replacement of ensemble averaging with time averaging appears quite restrictive. Most granular systems are out of equilibrium, and time averaging can still be meaningfully performed is the systems are at steady state.

8/ Fig. 11b: the scale shows velocities going up to 10^30 m/s, which is much higher than the speed of light. There is probably an issue with the scale.

9/ Fig 15 caption: the explanation of the "greyish" colour is incorrect. T_y is not zero at these locations, but undefined, since the volume fraction is zero (no particles). As indicated by equation C13, this would correspond to dividing a zero kinetic stress tensor trace by the zero density at that location, which numerically should result in a Not a Number.

10/ line 885: it would be good to define D right after eq. C13, rather than on line 895.

11/ eq C19: the equation uses d̄, which is not defined. Is it the same as D43 (eq. C6)?