the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
Thermo-hydrological and thermo-mechanical modeling of freezing soil and frost quake occurrence
Abstract. Frost quakes are seismic events originating in frozen ground, traditionally attributed to ice expansion when air temperature decreases rapidly, the soil is saturated and has little or no snow cover. However, novel observations presented here question the necessity of some of the previous assumption and meteorological conditions, driving us to consider alternative mechanisms. This study investigates frost quake formation through numerical modeling and seismological and hydrological observations in Tähtelä, Finland, during winter 2022–2023. We analyzed soil and atmospheric conditions during frost quake occurrences, noting a strong correlation with rapid air temperature decrease below -20 °C, and with varying snow cover. We modeled the thermo-hydrological (TH) processes, such as cryosuction-driven ice lens growth, using Amanzi-ATS, while the thermo-mechanical (TM) evolution was modeled with OpenGeoSys (OGS). The ATS TH simulation results suggest an important role of cryosuction for the appearance of frost quakes. Focusing on volumetric effects, the OGS TM simulation results reveal tensional and shear stress rates in the soil, both of which are able to cause fracturing leading to frost quakes. Future work should integrate fully-coupled thermo-hydro-mechanical simulations and laboratory experiments to refine predictive models and assess infrastructure risks in cold climates.
- Preprint
(6335 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
- CC1: 'Comment on egusphere-2026-2611', Ivo Baselt, 08 Jul 2026
-
RC1: 'Comment on egusphere-2026-2611', Kehua You, 18 Jul 2026
This study challenges the conventional understanding of frost quakes by demonstrating, through the integration of field observations, thermo-hydrological modeling, and thermo-mechanical modeling, that frost quakes can occur under conditions that have not been recognized in previous studies. The work highlights the importance of using a fully coupled thermo-hydro-mechanical model to investigate the mechanisms governing frost quakes. Overall, the manuscript is well written, and I only have a few minor comments.
- Line 178: What does (L_f) represent? Please define it when it is first introduced.
- Most seismic events were observed in wetlands and irrigated channels. A common characteristic of these environments is their relatively high water content. Does this observation suggest that water content is an important controlling factor for frost quake occurrence?
- Line 203: The simulated ice content decreases from 0.42 to 0.044 (approximately one order of magnitude) when the cell size increases from 2 cm to 5 cm. Is this result correct? If so, could the authors explain the strong sensitivity to spatial resolution?
- Line 211: The simulation domain extends from the ground surface to a depth of 25 m. Why is a boundary condition specified at a depth of 5 m rather than at the bottom of the model domain?
- Table 1: Both lambda_s and lambda_sR are mentioned. Are these two parameters different? If so, please clarify their definitions.
- Please include the equation used to calculate the bulk thermal conductivity. Because bulk thermal conductivity is a key parameter in the model and can be estimated using different mixing models.
- Line 245: The thermo-mechanical model uses the parallel model to calculate bulk thermal conductivity. Is the same approach also used in the thermo-hydrological model? If different formulations are adopted, please explain the rationale and discuss any potential impact on the simulation results.
- Figure 5a: I suggest adding a horizontal line indicating the -20 oC temperature to facilitate interpretation of the results.
- Line 313: What does "temperature rate" refer to? Would "rate of temperature change" or "temperature change rate" be a clearer description?
Citation: https://doi.org/10.5194/egusphere-2026-2611-RC1 -
RC2: 'Comment on egusphere-2026-2611', Anonymous Referee #2, 23 Jul 2026
The manuscript presents a numerical study using thermo-hydrological and thermo-mechanical modelling of frozen soil to investigate the mechanisms that induce frost quakes. Data collected at a site (Tähtelä, Finland) during extreme winter events in 2022-2023 include local seismological observations, air temperature, snow depth, soil water content and soil temperature at various depths.
Thermo-hydrological modelling was first conducted to test the hypothesis that ice lens growth driven by cryosuction is the dominant mechanism of frost heave. The soil column model was detailed, but the boundary conditions at the soil surface (soil-air interaction, the role of the snow layer, etc.) were not clarified. In addition, the heat transfer model (related to the ground thermal conductivity) was not explicitly stated. For instance, the soil’s thermal conductivity appears to be constant (Table 1). The results of this part are inconclusive.
Thermo-mechanical modelling was then performed to assess the role of volumetric thermal and phase-change strains in frost quake. As water drainage during freezing was ignored in this part, the frost heave estimated directly from phase-change strain would be overestimated. The role of thermal contraction, associated with the cooling of ice, is not clearly identified. As a consequence, significant compressive lateral stress was obtained (Figure 10e), whereas tensile stress would be expected to explain the cracks that induce frost quake. For this reason, the results of this part are equally inconclusive.
Further details comments
- Table 1: The parameters lambda_S and alpha^S_T were not explained.
- Table 1: It seems that c_pS is considered constant?
- Figure 9: It seems the Soil freezing characteristic curve was not considered in the model.
Citation: https://doi.org/10.5194/egusphere-2026-2611-RC2 -
RC3: 'Comment on egusphere-2026-2611', Anonymous Referee #3, 03 Aug 2026
This manuscript combines valuable seismological observations with thermo-hydrological and thermo-mechanical simulations to investigate frost-quake formation. However, the main mechanistic conclusions are not clearly supported by the models employed. In particular, cryosuction is defined and interpreted inconsistently, and ATS does not simulate discrete ice-lens formation (as far as I know); the OGS model imposes a fully saturated and laterally confined state that contradicts both the ATS results and the field conditions. Moreover, the apparent associations with frost quakes largely reflect their common dependence on the imposed temperature forcing rather than demonstrated causal mechanisms. These are fundamental issues requiring substantial new modeling and analysis.
- The definition of cryosuction and the sign of the Clapeyron relation appear incorrect. Eq. (3) is negative below freezing. Under the conventional compressive-positive pressure definition, the linear ice–water Clapeyron relation instead gives PI-PW∝ (Tm-T)>0. If the authors use a different stress convention, it must be explicitly stated and applied consistently. At present, Eqs. (3)–(5), the discussion of increasing cryosuction, and the positive values plotted in Fig. 8c are mutually inconsistent. The sign of the theta_f in Eq. (4), and the placement of the interfacial-tension ratio (beta), also require re-derivation.
- The quantity plotted as “cryosuction” in Fig. 8 is apparently not the quantity defined in (3). Equation (3) predicts approximately 1.22 MPa per degree of undercooling, and therefore about -3.7 MPa at -3°C. Figure 8 instead reaches approximately +8 MPa, which closely matches the ATS liquid–ice capillary-pressure evaluator after including a surface-tension factor. Thus, the plotted field appears to be an ATS constitutive diagnostic rather than PI-PW as defined. The exact ATS output variable, sign convention, equation units, and post-processing must be reported.
- ATS does not simulate ice-lens formation or growth. The Painter–Karra formulation partitions water among liquid, gas, and pore ice within a continuum cell and allows temperature-induced water redistribution through Richards flow. It does not represent segregated ice lenses, lens nucleation, a frozen fringe, particle rejection, lens spacing, opening of a flaw, soil displacement, or fracture propagation. Indeed, the authors explicitly state that the solid phase remains mechanically implicit and that the available pore volume caps combined water and ice saturation. Figure 6 therefore shows distributed pore-ice accumulation, not ice lenses. All claims of “cryosuction-driven ice-lens growth” must either be removed or supported using an actual frost-heave/ice-lens model.
- No quantitative evidence is presented that cryosuction caused meaningful water influx or excess-ice growth. The plotted “cryosuction” is essentially prescribed by temperature through the Clapeyron constitutive relation. Its increase during cold spells is therefore mathematically inevitable and cannot constitute independent mechanistic evidence. The manuscript does not show liquid pressure, hydraulic gradients, Darcy fluxes, cumulative water influx, water mass balance, or excess ice relative to freezing of in-situ pore water. A control simulation with cryosuction disabled is also absent. Without these diagnostics, the results cannot distinguish cryosuction-driven migration from simple local freezing.
- The Clapeyron equation provides an equilibrium thermodynamic relation between ice pressure, liquid pressure, and temperature; it does not independently determine PI. The manuscript nevertheless treats PI as a proxy for mechanical stress without showing PW, mechanical confinement, or the effective-stress relation. ATS does not solve an independent ice momentum equation in this configuration. Therefore, the plotted liquid–ice pressure difference cannot be interpreted directly as a simulated soil stress or as evidence that a fracture threshold was reached.
- ATS predicts a strongly partially saturated soil, with gas volume fractions of approximately 0.09–0.21. OGS instead assumes complete saturation and converts the entire pore volume to ice essentially by -2°C. Only the ATS surface temperature and initial temperature profile are transferred; water, ice, gas, liquid pressure, and water flux are not transferred. Calling this “one-way coupling” obscures that the two models represent different materials and initial states. The OGS phase-expansion results therefore cannot be interpreted as the mechanical response of the ATS-simulated soil.
- The sides and bottom boundaries are mechanically fixed, which force zero lateral strain and consequently produce large compressive lateral stress during freezing. Figure 10 reports stresses approaching 100MPa, despite a soil Young’s modulus of only 53MPa. Such stresses are implausible for unconsolidated sandy soil and would induce granular rearrangement, plasticity, damage, or failure long before the reported elastic state was reached. The statement that plasticity may be regarded as merely “post-critical” does not resolve this problem; once yielding occurs, the calculated elastic stresses and stored energy are no longer physically meaningful.
- A positive stress rate does not mean that the stress state is tensile; it may simply indicate unloading from a highly compressive state. Figure 10e appears predominantly compressive, while the red regions in Fig. 10f show changes in stress, not necessarily tensile stress. No tensile strength, shear strength, fracture toughness, damage law, failure criterion, energy-release rate, or pre-existing flaw geometry is included. Consequently, the claim that the calculated rates are “able to cause fracturing” is unsupported. Likewise, Eq. (10) gives the maximum shear stress on a rotated plane, not a modeled shear-stress component. Mode-III fracture cannot be inferred from this quasi-1D/2D calculation, and buckling is not represented at all.
- Frost-quake sources occur in spatially distinct wetland, channel, and observatory clusters. In contrast, soil temperature, water content, and snow depth are measured at one station and modeled as one homogeneous sandy column. Some events are hundreds of meters to more than a kilometer from that station. The reported 78 cm snow depth therefore cannot automatically be assigned to the quake source locations, particularly where snow redistribution, drainage channels, and soil properties likely differ. The conclusion of frost quakes occurring beneath thick snow requires source-location-specific snow and subsurface information.
- The CS615 is a dielectric reflectometer, not a direct measurement of total water, liquid water, or ice content. Freezing strongly alters the dielectric response, and no frozen-soil calibration, temperature correction, uncertainty range, or observation operator is described. The apparent decline in measured “water content” during freezing may largely reflect liquid-to-ice phase change rather than drainage or cryosuction. Agreement between this signal and simulated liquid water therefore cannot validate the proposed mechanism without an appropriate frozen-soil calibration.
- The “best-fit” parameter set is selected without an objective function, goodness-of-fit statistics, parameter ranges, or uncertainty analysis. No statistical test quantifies the claimed association between simulated variables and frost-quake occurrence, and no non-event cold periods, lead–lag relationships, or event-detection completeness are evaluated. Because temperature, ice content, and the plotted cryosuction are mathematically coupled, visual coincidence of their extrema is not evidence of independent causation.
- The reported mean ice content of 0.42 for the 2 cm grid exceeds the stated porosity of 0.3 and differs by an order of magnitude from the 5 cm result of 0.044; this is either a serious numerical sensitivity or a consequential error. No time-step convergence is provided even though the conclusions depend on rapid temperature and stress rates, and the authors themselves acknowledge sensitivity to forcing-sampling interval. The OGS mesh resolution is not reported. I found key parameter values are missing, including residual saturation, the phase-expansion coefficient, the sigmoid parameter (k), reference states, initial conditions, the spin-up procedure, and several constitutive options. The unexplained three-year OGS simulation, the internal groundwater boundary at 5 m within a 25 m domain, incorrect dates, errors in Table 1, and the mismatch between the stated ATS version and cited software archive further prevent reproduction.
- Can double-check what the most accurate definition of cryosuction is? Looks like the cryosuction referred to mostly occurs during the freezing process. Will it also exist during the thawing process?
Citation: https://doi.org/10.5194/egusphere-2026-2611-RC3
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 62 | 38 | 11 | 111 | 12 | 7 |
- HTML: 62
- PDF: 38
- XML: 11
- Total: 111
- BibTeX: 12
- EndNote: 7
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
I would like to congratulate the authors on this very interesting and timely contribution. I particularly appreciate the attempt to combine seismological observations with thermo-hydrological and thermo-mechanical modelling to reassess the mechanisms behind frost quake occurrence.
One aspect that I found especially stimulating is the potential link between frost-quake-related crack formation and preferential-flow concepts in frozen soils. The authors already highlight the broader relevance of their work for wintertime hydrological processes, especially where rain, snowmelt, freezing, and thawing interact in the Earth´s Critical Zone. However, the manuscript mainly focuses on the formation mechanism of cracks and does not further discuss the possible thermo-hydraulic consequences of such cracks once they have formed.
This connection seems highly relevant because there is a related branch of cryosphere and frozen-soil research that investigates how macropores and preferential pathways affect infiltration, heat transport, refreezing, drainage onset, and runoff generation in seasonally frozen soils. Many experimental studies necessarily rely on idealised representations of such preferential structures (e.g. Watanabe and Kugisaki, 2017; Mohammed et al., 2018; Pittman et al., 2020; Bauer et al., 2026). These experiments provide important benchmark data for dual-porosity, dual-permeability, and dual-domain modelling approaches, in which macropores are represented as a separate flow or continuum domain (e.g. Larsbo et al., 2019; Heinze and Blöcher, 2019; Mohammed et al., 2021; Heinze, 2021; Heinze, 2025; Khanahmadi et al., 2026).
This is where I see a particularly interesting contribution of the present manuscript. The crack formation discussed here could provide a mechanistic bridge between idealised macropore concepts and naturally generated preferential structures in frozen ground. In other words, frost-quake-related cracks may not only be interpreted as a mechanical consequence of freezing, but potentially also as transient macropore-like structures that influence subsequent water flow, heat transport, refreezing, and cryosuction during rainfall or snowmelt events. The causal direction may even be twofold: macropores may lose functionality when they refreeze, while newly formed cracks may subsequently create a new macropore-like system that becomes hydraulically active.
In this context, it would be very helpful if the manuscript could further discuss, or at least provide order-of-magnitude estimates for, the expected geometry of frost-quake-related cracks. Possible crack aperture, penetration depth, lateral extent, spacing, orientation, and connectivity would be highly relevant parameters. I fully understand that these quantities may not be directly observable from the present data set. Nevertheless, even a discussion of plausible ranges inferred from frost depth, ice-body thickness, seismic source locations, or the simulated stress and ice zones would substantially increase the value of the study for future modelling efforts.
Such information would allow two complementary modelling strategies. At the continuum scale, frost-quake-induced cracks could be represented in dual-porosity or dual-permeability frameworks, where the frozen matrix and the crack or macropore domain are treated as interacting domains with different hydraulic and thermal properties. For this, approximate geometric information would be required to estimate the corresponding macroporosity or crack-domain volume fraction. At a more explicit scale, crack aperture, length, spacing, and connectivity could be used to resolve individual cracks as discrete preferential pathways in numerical models. This would be a natural next step beyond idealised cylindrical macropores and would directly connect thermo-mechanical crack formation with thermo-hydrological function.
I therefore suggest extending the introduction and/or the discussion by explicitly addressing the possible hydrological relevance of frost-quake-related cracks. In particular, it would be valuable to clarify whether the authors view these cracks primarily as mechanical failure features, or whether they may also act as transient preferential-flow structures after their formation. This could open an interesting path for future work linking frost-quake mechanics, crack geometry, dual-domain modelling, and discrete crack-scale simulations in seasonally frozen soils.
Let me also add three minor comments which might help to improve the manuscript:
Line 125: From my point of view, it would be useful to cite the original work for the Kozeny-Carman method, Carman (1956), in addition to the more recent source.
Line 296: Shouldn´t the year in 17.11.2023 be 2022? Also, the text between 294 to 298 refer to dates between 2022 (I guess) until May 2024. However, the caption in Fig. 7 shows only the period until 2023. So, is Fig. 7 only one example?
Figure 1: According to the caption, subfigure c shows the soil station installation. However, I can only identify a borehole with some wires. Perhaps the subfigure could be improved by using a clearer image or by adding text arrows that indicate the relevant components of the installation.
References:
Bauer, Julian; Müller, Sebastian; Heinze, Thomas; Khanahmadi Bafghi, Homa; Baselt, Ivo (2026): Thermohydraulic experiments on water infiltration into frozen slopes: the role of macropores and initial water content. In: The Cryosphere 20 (6), S. 3483–3509. DOI: 10.5194/tc-20-3483-2026.
Carman, Philip Crosbie (1956): Flow of Gases Through Porous Media: Academic Press.
Heinze, Thomas; Blöcher, Johanna R. (2019): A model of local thermal non-equilibrium during infiltration. In: Advances in Water resources 132, S. 103394. DOI: 10.1016/j.advwatres.2019.103394.
Heinze, Thomas (2021): A Multi‐Phase Heat Transfer Model for Water Infiltration Into Frozen Soil. In: Water Resour. Res. 57 (10). DOI: 10.1029/2021WR030067.
Heinze, Thomas (2025): A local thermal non-equilibrium model for rain-on-snow events. In: Hydrol. Earth Syst. Sci. 29 (8), S. 2059–2080. DOI: 10.5194/hess-29-2059-2025.
Khanahmadi, Homa; Bauer, Julian; Baselt, Ivo; Heinze, Thomas (2026): The influence of macropores on the thermal state of soil during infiltration in the absence of thermal equilibrium. In: Journal of Hydrology 668, S. 134983. DOI: 10.1016/j.jhydrol.2026.134983.
Larsbo, Mats; Holten, Roger; Stenrød, Marianne; Eklo, Ole Martin; Jarvis, Nicholas (2019): A Dual‐Permeability Approach for Modeling Soil Water Flow and Heat Transport during Freezing and Thawing. In: Vadose zone j. 18 (1), S. 1–11. DOI: 10.2136/vzj2019.01.0012.
Mohammed, Aaron A.; Kurylyk, Barret L.; Cey, Edwin E.; Hayashi, Masaki (2018): Snowmelt Infiltration and Macropore Flow in Frozen Soils: Overview, Knowledge Gaps, and a Conceptual Framework. In: Vadose zone j. 17 (1), S. 1–15. DOI: 10.2136/vzj2018.04.0084.
Mohammed, Aaron A.; Cey, Edwin E.; Hayashi, Masaki; Callaghan, Michael V. (2021): Simulating preferential flow and snowmelt partitioning in seasonally frozen hillslopes. In: Hydrol. Process. 35 (8), Artikel e14277, e14277. DOI: 10.1002/hyp.14277.
Pittman, Freda; Mohammed, Aaron; Cey, Edwin (2020): Effects of antecedent moisture and macroporosity on infiltration and water flow in frozen soil. In: Hydrol. Process. 34 (3), S. 795–809. DOI: 10.1002/hyp.13629.
Watanabe, Kunio; Kugisaki, Yuki (2017): Effect of macropores on soil freezing and thawing with infiltration. In: Hydrol. Process. 31 (2), S. 270–278. DOI: 10.1002/hyp.10939.