This is a model description paper, describing innovations in sediment mixing and the treatment of CaCO3, which might impact the preservation of paleo proxies encoded in the CaCO3, in ways which the model would be intended to disentangle.  The results of the model presented here are focused on the effect of the sedimentary mixing dynamics, and how they might be coupled with the dynamics of CaCO3 chemical erosion. 

The paper could be much more interesting with the addition of some additional sensitivity experiments.  The authors cite the literature on the 14-C distribution in CaCO3 in the equatorial Pacific, which gets older as the %CaCO3 goes down, the opposite of what you would expect if the CaCO3 dissolved homogeneously – a low %CaCO3 would imply a short residence time and a low age.  Somehow shells become “armored” from dissolution if they survive an early stage.  I have attempted to replicate this by using multiple phases of CaCO3 with varying dissolution kinetics or solubilities, but I never managed to reproduce the trend in the observations.  I think this observation is mirrored in the 14-C age distribution of mollusk shells; it seems to be a general thing.  Perhaps attacking this problem is for a future study, while this is just a model development paper, but it doesn’t seem like another sensitivity plot or two would add too much baggage to the paper. 

I had some questions as I was reading, points of confusion.  For the Fickian diffusion, does the rate taper off exponentially with depth or is it an abrupt cutoff?  What is the difference between Ficking diffusion and homogeneous mixing? 

line 45: clarify what you mean; I would have thought that Fickian diffusion is random mixing.  In that section it might also be worth mentioning that some models use uniform mixing down to an abrupt cutoff, while others use an exponential dropoff in mixing rate.  

line 110.  Kudos for coding the model up in multiple languages!  But why python, when Julia seems just as elegant and flexible and also lots faster?  

line 225.  How can the initial condition have vanishinly small concentrations of all of the solid phases? Don’t they have to sum to fill the solid volume implied by the time-independent porosity?  (On further reading I understood that this is an initial state for an iteration, which by the time it converges will have solved the problem.  However, maybe a sentence here would help clarify.) 

line 260.  “time implicit method”. It took a bit of digging to figure out if the model is time-dependent or steady-state?  This was a clue in the text but it didn’t specify whether it applied to solid and dissolved species or what.  I figured it out from the figures, but it would have been useful to state it more explicitly earlier on.  (And on that, why bother with time dependence for the solutes?  It must slow things down a lot.)

Equations 23-24.  These are succinct descriptions of the matrices, but they are not very transparent as far as explaining what the mixing models do.  Why does homogeneous mixing use P rather than D?   (On subsequent rereading there is an extensive discussion on the formulation of homogeneous mixing, but a bit of summary here would be helpful.)

Would it be possible to make some kind of visualization of the transportation matrix, a heat map of some sort that would show how the mixing mechanisms differ? 

line 275.  I’m not familiar with LABS and would have appreciated a bit more. Is there reason to believe that the model organism is doing about the right amount of mixing?  (In further reading I found this in a comparison of results, but a statement here of where things are going would be useful, or the methodology of how the model was tuned.)  Is there also a good reason to simplify the spectrum of real critters into the single size?  How does the mixing differ from Fickian? 

line 500.  It would be interesting to integrate how much excess CaCO3 dissolution occurred due to the change in solubility (water depth) – how the buffering strength of the sea floor depends on the mixing model.

line 515.  Do the smaller particles have higher surface to volume, and also less mass, so they dissolve more quickly for those reasons also? It would be useful to add differences in kinetics or solubility here, and separate out the different effects. 

General Comments:

This paper presents a new model (IMP) that builds on previous reactive-transport models, aiming to better decipher proxy signals while addressing complications due to bioturbation and chemical erosion. The model includes complexity previously not addressed in similar models, as the authors break down CaCO3 particles into “classes” based on recorded proxy values and biologic and physiochemical characteristics. Further, they investigate the effects of sediment mixing and dissolution in their model using stable isotope records of a hypothetical environmental change event. Generally, the paper is very well written and undoubtedly of interest to the audience of Geoscientific Model Development. The authors do a good job laying out the motivation for their model and experiments and clearly outline the methods and assumptions. This work will help to improve understanding and interpretation of marine proxy records, providing a strong theoretical basis for disentangling geospatial variability in records from similar archives. I applaud the authors for providing such a well-revised, organized and thought-out manuscript. My critique is minimal and provided as a suggestion for two general areas of improvement: (1) expansion of explanation of transition matrices, namely automation-based and (2) expansion of the discussion.

Specific Comments:

The manuscript would benefit from an expansion on the LABS simulation approach; automation-based transition matrices described in the methods could be more thorough. As written, the paper requires unfamiliar readers to investigate this approach outside of the paper. This could also be achieved in the introduction.

While this model will be applied to interpreting archives of geologic events and such events are cited as motivation, there is little to no discussion later in the paper of the significance of their experiments with regard to these events. For example, what is the significance of the model result in which coarse fraction species become more dominant, in terms of records of past abrupt events of environmental change? Perhaps the authors decided to stay away from interpreting their model development results in terms of geologic applications, but some model-data comparison may be warranted in the discussion. This may or may not include a more representative simulation of an early Eocene hyperthermal event.

Following this more generally, the discussion section of the manuscript is slightly limited and could be expanded. For example, how does this new model and the results of experiments in this study inform understanding of examples of processes outlined in the introduction? How might the findings here bias proxies in specific geologic archives (e.g., size fraction differences previously unaccounted for in proxy records)? This type of expansion would not necessarily require re-interpretation of records of e.g. the PETM, but rather clearly lay out the implications of their experiments which may be significant to a proxy-based researcher in the field.