Response to review of Daniele Pedretti

DP: The manuscript ("Preferential Pathways for Fluid and Solutes in Heterogeneous Groundwater Systems: Self-Organization, Entropy, Work") provides a new framework that combines the use of free energy and entropy to characterize and quantify the emergence of preferential flow and channelled transport in heterogeneous media. Although the specific components of this framework are not novel themselves ,as correctly acknowledged by the authors, their combined use makes it a novel way that help disentangling open questions regarding the mechanisms of transport in heterogenous porous media. The paper is well written. The objective, methodologies and conclusions are clear. I have summarized a line-by-line set of comments for the authors, which makes me recommending me accepting this manuscript after major revisions.


RESPONSE:
We agree that a transfer of the proposed assessment to a three 3D stochastic media would yield interesting insights, which will for sure differ somewhat from what we found in 2D. However, we respectfully disagree that this might imply a different qualitatively behavior, as long as we work in a confined system (no flow boundary conditions for the upper, lower, inlet and outlet boundaries). The local changes in power arise from the local feedback on the pressure head gradient in front of the low conductivity bottlenecks. Gradients steepen ahead of the bottlenecks, which implies a higher power, locally. This feedback will also occur in a 3D confined system, as it is a direct result of the boundary conditions. It would not likely occur in an aquifer with a free surface neither in 2D or 3D. In the revised manuscript we will state that an analysis in 3 D is the next step of a forthcoming study, and explain why think that we expect qualitatively similar results due to the above stated reasons.
We furthermore like to stress, that percolation considerations are not relevant here, as the domains are well connected and well above a percolation threshold. Moreover, the original study by Edery et al. (2014), which we cite, presented already a critical path analysis in reference to the percolation threshold, based on the common assumption that preferential flows are a manifestation of percolation, controlled by the lower cutoff for the hydraulic conductivity from which a path is possible. The limitations of percolation theory in evaluating the preferential flow are presented therein, and as such, the equivalence or connection between percolation and entropy is not straightforward.

DP:
Thanking the authors for considering our 2017 WRR publication, I also suggest having a look at our follow-up manuscript (Bianchi and Pedretti 2018 WRR https://doi.org/10.1029/2018WR022827) where we extend our previous theory by computing the geological entrogram on evolving sampling scales. I think that most conclusions we got in those studies there are very much in line with those obtained through this study. Indeed, in the 2018 paper we also address the question of 2D vs 3D models, and described at page 4444 how solute particles tend to sample specific K clusters when travelling in the heterogeneous media.

RESPONSE:
We will be pleased to read your 2018 paper and refer to the study and conclusions in our revised manuscript, if/where appropriate. Best Regards

Daniele Pedretti
University of Milan -UNIMI (Italy). Line-specific comments DP: L72-76 I wouldn't be so strict. Someone has succeeded in this task (Zhang et al 2013 JH for instance). What is really complicated is finding a universal way to predict solute transport based on the aquifers geological structures. In Bianchi and Pedretti's works on geological entropy we found an explanation for that: the lower the structure's Shannon entropy, the more organized the flow and transport patterns in the field. In that set of works our aim was to start from the geology and not from self-customed flow fields (e.g. power-law distributed seepage velocities).

RESPONSE:
We thank DP for pointing this out, and will revise this passage to be less strict. DP: L97 "the probability of solutes to pass through HIGH (not low) conductivity regions". Please , fix it.

RESPONSE:
We will fix this; thanks for noting the typo. DP: L98 Please consider also our follow-up study (Bianchi and Pedretti 2018 WRR), where we study evolving scales of geological entropy rather than studying fixed-size blocks (as we did in the 2017 paper). In the 2018 paper we developed the concept of entrogram scale, which is also nicely correlated with the emergence of preferential flow and solute channelling.

RESPONSE:
We will be pleased to read your 2018 paper and refer to the study and conclusions in our revised manuscript, if/where appropriate. DP: L101 enigmatic OK, emergent not really I would say.

RESPONSE:
Agreed. We will delete the term emergent. DP: L104 Again, I wouldn't be so strict ("virtually impossible"). I'd rather just say that such predictions remain challenging. For instance, Bianchi and Pedretti works or Zhang et al 2013 showed that it can be done. There is also a set of works by Rizzo and de Barros showing that predictions can be made starting from the aquifer structures.

RESPONSE:
Agreed. We will revise this passage accordingly. We note, too, that we can also achieve useful predictions of preferential transport of solutes in the partially saturated zone, in cases where detailed information about the pdf of macropores in soil and their connectivity are available. DP: L158 please see comment at L98.

RESPONSE:
As for the L98 comment, we will read your 2018 paper and refer to the study and conclusions in our revised manuscript, if/where appropriate. DP: L165-on. Rather than Objectives, these are Results and Conclusions.

RESPONSE:
We agree that we provide here some foreshadowing on the results. The idea is to interest and motivate the reader to study, in particular, the theory section. In the revised manuscript, from the end of line 164 ("Specifically"), we will define the text as a new paragraph.
DP: L185 I wonder if all these nice concepts can be exported directly to 3D models, considering the different percolation thresholds between 2D and 3D models. Could the authors discuss on this?

RESPONSE:
We think that the above mentioned local feedback on the head gradient will also occur in a 3D confined system. As such we expect qualitatively similar behavior in three 3D, as noted in our response to the second comment above. In the revised manuscript, we will include a brief consideration of these comments/responses in the conclusion and outlook section. DP: L216 why no local dispersion? This is a physical mechanism, which can substantially modify the solute pathways by increasing mixing and coalescence among the so-called "lamellas". Why neglecting it? I think this should have been investigated, from low to high Peclet numbers.

RESPONSE:
The modeling approach used here is well accepted and accurate/sufficient for our purposes. Hydrodynamic dispersion is a macroscale fingerprint of diffusive transversal mixing of solute between flow lines of different fluid velocities. Adding a local dispersion term would incorporate implicit assumptions about subscale heterogeneity in velocities. DP: L256 I totally agree, but again, I think that transversal dispersion could have a big impact here.

RESPONSE:
We agree that transversal dispersion would of course work against steepening of lateral gradients. But lateral dispersion is the result of variability in the transversal flow field. We do not see much option for that, because the system is confined and flow is at steady state.

DP: L388 see comment at L98
RESPONSE: We will address this point, as noted in our response to the L98 comment. DP: L395 how this new theory could be connected to previously developed ones, such as the "lamella" description of solute transport in heterogenous media (e.g. https://doi.org/10.1017/jfm.2015.117) which also strongly depends on concentration gradients transversal to the main flow directions? or the concepts of "least-resistance paths" (e.g. https://doi.org/10.1002/2017WR020418)?

RESPONSE:
We think these concepts are well connected. At the end of the day, an increase in spatial organization implies to a steepen gradient. There are many ways to express this, but there is only one mechanism to explain why this implies ordering and this is the second law of thermodynamics. This is because production of physical entropy implies essentially to deplete gradients. DP: L403-405 this looks like a conclusion of this work, rather than a result. Consider moving it to the appropriate sections RESPONSE: Indeed we do a little bit of foreshadowing here, but we think this is helpful to keep the reader on track. We state clearly that "In the following, we demonstrate". DP: L442 this is very similar to a conclusion by Bianchi and Pedretti 2018 (page 4444), which reads " These observations are further confirmed by the CDFs of the subsamples of K values, which show a significantly higher variability for 2-D fields compared to 3-D (Figures 6c and 6d). These results strengthen the knowledge that solutes tend to travel in the upper 15-20% of the K distribution (K classes 32-40 in Figure 6), which tend to fully percolate in 3-D correlated random fields regardless of their overall structure (Fogg et al., 2000;Fogg & Zhang, 2016;Harter, 2005). " Consider commenting on that RESPONSE: We will happily include this point in our discussion here. DP: L457 I find it a bit complicated to relate "watts" to something related to groundwater. I mean, the entire derivation of power is clearly described in the previous sections, but as a hydrogeologist I have some problem to understand, for instance, if this is a high/low power. For instance would 2 watts be a high or lower power in this context?

RESPONSE:
We absolutely agree that power and W/m is not very common groundwater and vadose zone hydrology. Zehe et al. (2013) analyzed energy conversion associated with infiltration and found that macropores increase power in the infiltrating water flux. The maximum values during a rainfall event were of order 2 W/m2. The use of W/m2 is much more common when dealing with the landsurface energy balance, as evaporation as water flux can be expressed as energy flux as well. In this context it is interesting to recall that the climatological land-surface energy balance is of order 100 W/m2. By comparison, a difference of 2 W/m2 is hence quite significant. In the revised manuscript, we will add text to include this background information, and comment on the significance of a 2 W per unit width increase as noted on L457. DP: L455-478: these are great findings. It is particularly interesting that at some point the system behaves effectively as a 1D system. I'm wondering if they also hold for a 3D system, which is in general more percolated than a 2D one.

RESPONSE:
Thanks. You raise a very good point. We expect that the 3D system will deviate even more from the 1D approximation, but the local feedback on the head gradient will remain in a confined aquifer. As mentioned above, we will address this point in the conclusions and outlook section of the revised manuscript. DP: L487-488 for the same token, then solute injection mode (i.e resident vs fluxaveraged) could be also important to control the "effective" system power, right? Even if the flow field is the same, you change the way particles are already injected into certain