MOIST: a MATLAB-based fully coupled one-dimensional isotope and soil water transport model

Fu, Han; Neil, Eric J.; Li, Huijie; Si, Bingcheng

doi:https://doi.org/10.5194/hess-2022-422

Preprints

https://doi.org/10.5194/hess-2022-422

Preprints

16 Jan 2023

| 16 Jan 2023

Status: this discussion paper is a preprint. It has been under review for the journal Hydrology and Earth System Sciences (HESS). The manuscript was not accepted for further review after discussion.

MOIST: a MATLAB-based fully coupled one-dimensional isotope and soil water transport model

Han Fu, Eric J. Neil, Huijie Li, and Bingcheng Si

Abstract. Modeling water stable isotope transport in soil is crucial to sharpen our understanding of water cycles in terrestrial ecosystems. However, isotope and soil water transport are not fully coupled in current models. In this study, we developed MOIST: a MATLAB-based one-dimensional isotope and soil water transport model, a program that solves one-dimensional water, heat, and isotope transport equations simultaneously. Results showed that the MOIST model has good agreements to the theoretical tests and semi-analytical solutions of isotope transport under fixed boundary conditions. Furthermore, we validated the program with short- and long-term measurements from lysimeters studies. The overall Nash-Sutcliff efficiency coefficient (NSE) of soil water and deuterium (²H) transport for the short-term measurements are 0.66 and 0.69, respectively, with respective determine coefficient (R²) of 0.82 and 0.70, mean absolute error (MAE) of 0.02 m³ m^-3 and 11.84 ‰. For the long-term lysimeter study, the overall NSE, R², and MAE of simulated δ¹⁸O are 0.47, 0.49, and 0.92 ‰, respectively. These indices indicated the excellent performance of the MOIST model in simulating water flow and isotope transport in simplified ecosystems, suggesting a great potential of our program in promoting understandings of ecohydrological processes in terrestrial ecosystems.

Received: 16 Dec 2022 – Discussion started: 16 Jan 2023

Download & links

Han Fu et al.

Status: closed

RC1:
'Comment on hess-2022-422', Anonymous Referee #1, 03 Feb 2023

I suggest rejection for the following reasons:
Equation (1) is the water content-based Richards equation. It is well known that this kind of formulation cannot handle saturated problems and is not well posed at the interface between two layers, because water content is discontinuous. The mixed form of Richards equation should be used.
The heat transfer equation (2) is not correct. Csoil depends on the water content and should be embed in the time derivative.
More information should be provided concerning the solvers (ode113, ode23tb).
The tests 1 to 6 are very qualitatively discussed. Only different types of processes are checked. Physical processes can be verified but it does not mean that the computed variables and the process kinetics are correct. Moreover, these tests are development tests. They do not provide any new information on processes and therefore, should not be part of the manuscript. It is expected that models overcome these kinds of tests before publication.
L479-480: the reason for poor MAE value is unclear to me.
L521-540: The analysis of the difference between fully coupled or sequential approach (segregation) is convincing but it applies for an explicit time scheme discretization whereas HYDRUS and SiSPat use an implicit scheme. Moreover, the flow equation is written in terms of water content for MOIST, the other codes are using pressure based or a mixed form of Richards equation.
L608-610: The discussion about boundary conditions and intermodal conductivity is very popular. There are key papers not cited in the manuscript that review some of the techniques (see for example Belfort et al.,. On equivalent hydraulic conductivity for oscillation–free solutions of Richards equation. Journal of Hydrology, 2013, 505, pp.202-217).
MOIST was used to simulate two types of experiments and the authors concluded that MOIST is more accurate and reliable. This is not supported by the provided results. These results only show that MOIST might be better calibrated not that the numerical scheme – fully coupled- is better than other schemes. Parameters used by MOIST and the other models should be given.
The comparisons do not provide any information on the code accuracy and efficiency. To demonstrate the ‘excellent performance of the MOIST,’ the authors should compare their code with other existing codes (for example looking at breakthrough curves at different locations) and check detailed mass balances, time and space discretization sensitivity and computer time.

Citation: https://doi.org/10.5194/hess-2022-422-RC1
- AC1: 'Reply on RC1', Han Fu, 10 Feb 2023
  
  For the responses to reviewer 1, please see the attatchment.
  
  Citation: https://doi.org/10.5194/hess-2022-422-AC1
RC2:
'Comment on hess-2022-422', Anonymous Referee #2, 27 Feb 2023

The paper has some interesting aspects. A fully coupled isotope transport model in the soil-plant-atmosphere continuum needs to be improved in the recent literature. Some approximations have been made to simulate isotope transport in soil using HYDRUS for example, but the results could be better. This is an intricate problem that must consider water content and movement influences water temperature, and both influence isotope transport and fractionation, and temperature may also affect water movement.

The paper claims to solve the transport equations simultaneously. I would like to see the numeric scheme that shows this back-forward process. And yet more information about it is needed. Unfortunately, numerical implementation has only one paragraph.
Another question is, how do equations 41 and 42 take isotope fractionation from temperature variation into account?
It has the potential for a good paper, but I recommend rejection and resubmission at this stage because I need to determine whether the paper answers the proposed objectives.

Citation: https://doi.org/10.5194/hess-2022-422-RC2
- AC2: 'Reply on RC2', Han Fu, 13 Mar 2023
  
  For the responses to reviewer 2, please see the attatchment.
  
  Citation: https://doi.org/10.5194/hess-2022-422-AC2
RC3:
'Comment on hess-2022-422', Anonymous Referee #3, 03 Mar 2023
General comments
The paper presents a model simulating the coupled 1D transport of water, energy and stable isotope in soils. It is centered on the presentation of the model formulation and on the evaluation of the accuracy of the numerical solution using theoretical tests, the comparison with analytical solutions and finally, the comparison with published field data. The main claim of the authors concerning the originality of their approach is that the transport of stable isotopes is fully coupled with the water and energy transport equations, as compared to other existing models where the water stable isotope transport equation is solved after resolution of the water and energy transport equations. The authors argue that, with their solution, coarser vertical spatial discretization can be used, without losing accuracy; non-coupled models requiring a very fine vertical discretization to ensure the same accuracy. The new model has therefore a serious advantage in terms of computing times.
However, the presentation made by the authors of this full coupling is not detailed and clear enough to be fully convincing (see detailed comments below). The comparison with the analytical simulations lacks quantitative measures. The comparison with the field measurement is interesting but not sufficient to claim that the MOIST model performs better than the other models. The moisture conditions are not very dry so the model is only tested in a limited range of conditions and not in the most sensitive ones when water vapor becomes dominant at the surface. The discussion about the performance of the model and the comparison with the other models is not supported by robust results, but appears more as assumptions made by the authors.
Two options for this paper could be considered:
The main objective is to show that the numerical methods is robust and accurate. If this is the case, the paper fails to demonstrate it because there is no robust comparison between other models and the accuracy with a coarser resolution is not really shown. Furthermore, I miss a sensitivity analysis on how the upper boundary condition is prescribed and on the choice of the kinetic fractionation factor formulation. Several options exist in the literature and the impact of this choice is not discussed in the paper. Such sensitivity analysis could complement the analysis of the model performance

The purpose could be to show that the model provides good simulations of existing experiments and provides new insights on these experiments. However, it is not the case, as the comparison is restricted to the statistical measures on the comparison between simulated and observed values.

The paper lies in between these two options but is not satisfactory for none of them. My feeling is that the paper presentation is centered on the model description and numerical performance and not on new scientific results. Therefore, the paper may be better suited to the Geoscientific Model Development journal (https://www.geoscientific-model-development.net/home.html ). Before being published either in GMD or HESS, the paper should strengthen the discussion about the model performance to get fully convincing results.

Specific comments
0/ It would be useful to provide a list of notations with the units.
1/L39: do you really address this problem with the MOIST model?
2/ L57: I do not consider that a Matlab program is something fully accessible as paying for an expensive Matlab licence is necessary to use the developments
3/ Eq. (4) D_v is not defined
4/ Eq. (5) looks strange
5/ Eq (10) in Braud et al. (2005) (their Equation (9) has an additional term). Why do you neglect it (as in Haverd and Cuntz, 2010)?
6/ Eq. (12) there are several options in the literature for the specification of n_D. Why do you chose this formulation?
7/ Eq. (24) and (25) (even if taken from Haverd and Cuntz, 2010) are strange as they are already provided in a discretized form, contrarily to the other equations. Same comment for Eq. (27)
8/ section 2.1.4: the formulation of the boundary conditions for the various equations is very important for getting correct results. The authors put much effort in fully coupling the equations in the soil but put much less attention in the formulation of the boundary conditions. The choices made would require more justifications even if they seem to be similar to Haverd and Cuntz (2010).
More generally, it seems that the authors have recoded the Haverd and Cuntz (2010) model without considering the litter. It would be relevant to say it if it the truth.
8/ L162. Should be Eq. (24) and (25)?

9/ L191: should be Eq. (19)?
10/ the description of the numerical implementation is too short. What are the variables that you are computing? I would be curious to see the three discretized coupled equations in order to see where is the benefit of having them fully coupled. In other words, on which system of equations do you implement the Matlab solver (in the equations presently given in the paper, there are more than three unknown sets of variables).
11/ L194: when you use the soil water pressure as variable, it is continuous at the interface between layers with different hydraulic properties, which limits the problem you mention here.
11/ L220. The reference with the DOI of the dataset should appear in the reference list, not only in the text. This is a reference like a standard paper and it should be cited as such.
12/ Eq. (34) and (35): did you checked the MOIST model behavior using the Van Genuchten (1980) model for the retention curve and the Brooks and Corey (1964) model for the hydraulic conductivity as done in Braud et al. (2005)?
13/Figure 5, 6, 7: You could zoom on the 0-0.5 cm layer and avoid the 0 y-axis at the very top of the figure to get more legible figures. Furthermore, the Test case 3 is almost invisible (at least on a printed version).
14/Section 3.2: graphical comparison is fine, but I would like to see a comparison between the simulated peak concentration and slope of the oxygen 18 – deuterium relationships (see Braud et al., 2005, Tables 3 and 5). This would provide a more robust evaluation of the model performance.
15/ the demonstration from Eq. (43) to (47) would require more details to be fully understandable.
16/ Figure 12 is not cited in the paper. Furthermore, if the purpose of the figure were to demonstrate that the coarser vertical spatial resolution provides results as accurate as in Figure 5, it would be more informative to show the difference between both simulations. Visually, it seems that the peak is simulated deeper in Figure 12 than in Figure 5, which would not be very satisfactory as the peak is located at the evaporation front.
17/ L588-593: this paragraph is not supported by the results. Other reasons than solving or not the fully coupled equations could explain discrepancies between the models, one of them being the specification of the boundary conditions or the specification of the kinetic fractionation factor, that may be different in the different models.
18/ L595-610: the question of the numerical method (cell-centered versus vertex-centered) is strange here as the method used in the MOIST model has not been presented before.

References

Brooks, R. H. and Corey, A. T.: Hydraulic properties of porous media, Colorado State University, Fort Collins, 27 pp.1964.
Braud, I., Bariac, T., Gaudet, J. P., and Vauclin, M.: SiSPAT-Isotope , a coupled heat, water and stable isotope (HDO and H218O) transport model for bare soil. Model description and first verification, Journal of Hydrology, 309, 277-300, 2005.
Haverd, V. and Cuntz, M.: Soil-Litter-Iso: A one-dimensional model for coupled transport of heat, water and stable isotopes in soil with a litter layer and root extraction, Journal of Hydrology, 388, 438-455, 2010..
Van Genuchten, M. T.: A closed-form equation for predicting the hydraulic conductivity of unsaturated soils, Soil Science Society of America Journal, 44, 892-898, 1980.
Citation: https://doi.org/10.5194/hess-2022-422-RC3
- AC3: 'Reply on RC3', Han Fu, 13 Mar 2023
  
  For the responses to reviewer 3, please see the attatchment.
  
  Citation: https://doi.org/10.5194/hess-2022-422-AC3

Status: closed

RC1:
'Comment on hess-2022-422', Anonymous Referee #1, 03 Feb 2023

I suggest rejection for the following reasons:
Equation (1) is the water content-based Richards equation. It is well known that this kind of formulation cannot handle saturated problems and is not well posed at the interface between two layers, because water content is discontinuous. The mixed form of Richards equation should be used.
The heat transfer equation (2) is not correct. Csoil depends on the water content and should be embed in the time derivative.
More information should be provided concerning the solvers (ode113, ode23tb).
The tests 1 to 6 are very qualitatively discussed. Only different types of processes are checked. Physical processes can be verified but it does not mean that the computed variables and the process kinetics are correct. Moreover, these tests are development tests. They do not provide any new information on processes and therefore, should not be part of the manuscript. It is expected that models overcome these kinds of tests before publication.
L479-480: the reason for poor MAE value is unclear to me.
L521-540: The analysis of the difference between fully coupled or sequential approach (segregation) is convincing but it applies for an explicit time scheme discretization whereas HYDRUS and SiSPat use an implicit scheme. Moreover, the flow equation is written in terms of water content for MOIST, the other codes are using pressure based or a mixed form of Richards equation.
L608-610: The discussion about boundary conditions and intermodal conductivity is very popular. There are key papers not cited in the manuscript that review some of the techniques (see for example Belfort et al.,. On equivalent hydraulic conductivity for oscillation–free solutions of Richards equation. Journal of Hydrology, 2013, 505, pp.202-217).
MOIST was used to simulate two types of experiments and the authors concluded that MOIST is more accurate and reliable. This is not supported by the provided results. These results only show that MOIST might be better calibrated not that the numerical scheme – fully coupled- is better than other schemes. Parameters used by MOIST and the other models should be given.
The comparisons do not provide any information on the code accuracy and efficiency. To demonstrate the ‘excellent performance of the MOIST,’ the authors should compare their code with other existing codes (for example looking at breakthrough curves at different locations) and check detailed mass balances, time and space discretization sensitivity and computer time.

Citation: https://doi.org/10.5194/hess-2022-422-RC1
- AC1: 'Reply on RC1', Han Fu, 10 Feb 2023
  
  For the responses to reviewer 1, please see the attatchment.
  
  Citation: https://doi.org/10.5194/hess-2022-422-AC1
RC2:
'Comment on hess-2022-422', Anonymous Referee #2, 27 Feb 2023

The paper has some interesting aspects. A fully coupled isotope transport model in the soil-plant-atmosphere continuum needs to be improved in the recent literature. Some approximations have been made to simulate isotope transport in soil using HYDRUS for example, but the results could be better. This is an intricate problem that must consider water content and movement influences water temperature, and both influence isotope transport and fractionation, and temperature may also affect water movement.

The paper claims to solve the transport equations simultaneously. I would like to see the numeric scheme that shows this back-forward process. And yet more information about it is needed. Unfortunately, numerical implementation has only one paragraph.
Another question is, how do equations 41 and 42 take isotope fractionation from temperature variation into account?
It has the potential for a good paper, but I recommend rejection and resubmission at this stage because I need to determine whether the paper answers the proposed objectives.

Citation: https://doi.org/10.5194/hess-2022-422-RC2
- AC2: 'Reply on RC2', Han Fu, 13 Mar 2023
  
  For the responses to reviewer 2, please see the attatchment.
  
  Citation: https://doi.org/10.5194/hess-2022-422-AC2
RC3:
'Comment on hess-2022-422', Anonymous Referee #3, 03 Mar 2023
General comments
The paper presents a model simulating the coupled 1D transport of water, energy and stable isotope in soils. It is centered on the presentation of the model formulation and on the evaluation of the accuracy of the numerical solution using theoretical tests, the comparison with analytical solutions and finally, the comparison with published field data. The main claim of the authors concerning the originality of their approach is that the transport of stable isotopes is fully coupled with the water and energy transport equations, as compared to other existing models where the water stable isotope transport equation is solved after resolution of the water and energy transport equations. The authors argue that, with their solution, coarser vertical spatial discretization can be used, without losing accuracy; non-coupled models requiring a very fine vertical discretization to ensure the same accuracy. The new model has therefore a serious advantage in terms of computing times.
However, the presentation made by the authors of this full coupling is not detailed and clear enough to be fully convincing (see detailed comments below). The comparison with the analytical simulations lacks quantitative measures. The comparison with the field measurement is interesting but not sufficient to claim that the MOIST model performs better than the other models. The moisture conditions are not very dry so the model is only tested in a limited range of conditions and not in the most sensitive ones when water vapor becomes dominant at the surface. The discussion about the performance of the model and the comparison with the other models is not supported by robust results, but appears more as assumptions made by the authors.
Two options for this paper could be considered:
The main objective is to show that the numerical methods is robust and accurate. If this is the case, the paper fails to demonstrate it because there is no robust comparison between other models and the accuracy with a coarser resolution is not really shown. Furthermore, I miss a sensitivity analysis on how the upper boundary condition is prescribed and on the choice of the kinetic fractionation factor formulation. Several options exist in the literature and the impact of this choice is not discussed in the paper. Such sensitivity analysis could complement the analysis of the model performance

The purpose could be to show that the model provides good simulations of existing experiments and provides new insights on these experiments. However, it is not the case, as the comparison is restricted to the statistical measures on the comparison between simulated and observed values.

The paper lies in between these two options but is not satisfactory for none of them. My feeling is that the paper presentation is centered on the model description and numerical performance and not on new scientific results. Therefore, the paper may be better suited to the Geoscientific Model Development journal (https://www.geoscientific-model-development.net/home.html ). Before being published either in GMD or HESS, the paper should strengthen the discussion about the model performance to get fully convincing results.

Specific comments
0/ It would be useful to provide a list of notations with the units.
1/L39: do you really address this problem with the MOIST model?
2/ L57: I do not consider that a Matlab program is something fully accessible as paying for an expensive Matlab licence is necessary to use the developments
3/ Eq. (4) D_v is not defined
4/ Eq. (5) looks strange
5/ Eq (10) in Braud et al. (2005) (their Equation (9) has an additional term). Why do you neglect it (as in Haverd and Cuntz, 2010)?
6/ Eq. (12) there are several options in the literature for the specification of n_D. Why do you chose this formulation?
7/ Eq. (24) and (25) (even if taken from Haverd and Cuntz, 2010) are strange as they are already provided in a discretized form, contrarily to the other equations. Same comment for Eq. (27)
8/ section 2.1.4: the formulation of the boundary conditions for the various equations is very important for getting correct results. The authors put much effort in fully coupling the equations in the soil but put much less attention in the formulation of the boundary conditions. The choices made would require more justifications even if they seem to be similar to Haverd and Cuntz (2010).
More generally, it seems that the authors have recoded the Haverd and Cuntz (2010) model without considering the litter. It would be relevant to say it if it the truth.
8/ L162. Should be Eq. (24) and (25)?

9/ L191: should be Eq. (19)?
10/ the description of the numerical implementation is too short. What are the variables that you are computing? I would be curious to see the three discretized coupled equations in order to see where is the benefit of having them fully coupled. In other words, on which system of equations do you implement the Matlab solver (in the equations presently given in the paper, there are more than three unknown sets of variables).
11/ L194: when you use the soil water pressure as variable, it is continuous at the interface between layers with different hydraulic properties, which limits the problem you mention here.
11/ L220. The reference with the DOI of the dataset should appear in the reference list, not only in the text. This is a reference like a standard paper and it should be cited as such.
12/ Eq. (34) and (35): did you checked the MOIST model behavior using the Van Genuchten (1980) model for the retention curve and the Brooks and Corey (1964) model for the hydraulic conductivity as done in Braud et al. (2005)?
13/Figure 5, 6, 7: You could zoom on the 0-0.5 cm layer and avoid the 0 y-axis at the very top of the figure to get more legible figures. Furthermore, the Test case 3 is almost invisible (at least on a printed version).
14/Section 3.2: graphical comparison is fine, but I would like to see a comparison between the simulated peak concentration and slope of the oxygen 18 – deuterium relationships (see Braud et al., 2005, Tables 3 and 5). This would provide a more robust evaluation of the model performance.
15/ the demonstration from Eq. (43) to (47) would require more details to be fully understandable.
16/ Figure 12 is not cited in the paper. Furthermore, if the purpose of the figure were to demonstrate that the coarser vertical spatial resolution provides results as accurate as in Figure 5, it would be more informative to show the difference between both simulations. Visually, it seems that the peak is simulated deeper in Figure 12 than in Figure 5, which would not be very satisfactory as the peak is located at the evaporation front.
17/ L588-593: this paragraph is not supported by the results. Other reasons than solving or not the fully coupled equations could explain discrepancies between the models, one of them being the specification of the boundary conditions or the specification of the kinetic fractionation factor, that may be different in the different models.
18/ L595-610: the question of the numerical method (cell-centered versus vertex-centered) is strange here as the method used in the MOIST model has not been presented before.

References

Brooks, R. H. and Corey, A. T.: Hydraulic properties of porous media, Colorado State University, Fort Collins, 27 pp.1964.
Braud, I., Bariac, T., Gaudet, J. P., and Vauclin, M.: SiSPAT-Isotope , a coupled heat, water and stable isotope (HDO and H218O) transport model for bare soil. Model description and first verification, Journal of Hydrology, 309, 277-300, 2005.
Haverd, V. and Cuntz, M.: Soil-Litter-Iso: A one-dimensional model for coupled transport of heat, water and stable isotopes in soil with a litter layer and root extraction, Journal of Hydrology, 388, 438-455, 2010..
Van Genuchten, M. T.: A closed-form equation for predicting the hydraulic conductivity of unsaturated soils, Soil Science Society of America Journal, 44, 892-898, 1980.
Citation: https://doi.org/10.5194/hess-2022-422-RC3
- AC3: 'Reply on RC3', Han Fu, 13 Mar 2023
  
  For the responses to reviewer 3, please see the attatchment.
  
  Citation: https://doi.org/10.5194/hess-2022-422-AC3

Han Fu et al.

Model code and software

MOIST Han Fu, Bingcheng Si https://github.com/HAN-2/MOIST

Han Fu et al.

Viewed

Total article views: 848 (including HTML, PDF, and XML)

HTML	PDF	XML	Total	BibTeX	EndNote
548	273	27	848	11	9

HTML: 548
PDF: 273
XML: 27
Total: 848
BibTeX: 11
EndNote: 9

Views and downloads (calculated since 16 Jan 2023)

Month	HTML	PDF	XML	Total
Jan 2023	232	71	4	307
Feb 2023	117	53	5	175
Mar 2023	84	45	8	137
Apr 2023	37	26	0	63
May 2023	24	12	0	36
Jun 2023	16	12	2	30
Jul 2023	21	31	6	58
Aug 2023	15	20	2	37
Sep 2023	2	3	0	5

Cumulative views and downloads (calculated since 16 Jan 2023)

Month	HTML	PDF	XML	Total
Jan 2023	232	71	4	307
Feb 2023	117	53	5	175
Mar 2023	84	45	8	137
Apr 2023	37	26	0	63
May 2023	24	12	0	36
Jun 2023	16	12	2	30
Jul 2023	21	31	6	58
Aug 2023	15	20	2	37
Sep 2023	2	3	0	5

Viewed (geographical distribution)

Total article views: 794 (including HTML, PDF, and XML) Thereof 794 with geography defined and 0 with unknown origin.

Country	#	Views	%

Latest update: 04 Sep 2023

Short summary

Current hydrological models segregate water and isotope transport within soil. Thus, MOIST, a MATLAB-based one-dimensional isotope and soil water transport model, was developed. Results indicated that MOIST had good performances on simulating transport of isotope and water within soil under theoretical and realistic conditions, even outperformed than HYDRUS-1D. Suggesting a great potential of MOIST in promoting understandings of ecohydrological processes in terrestrial ecosystems.


Total:	0
HTML:	0
PDF:	0
XML:	0