the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
A Minimal Bonded Discrete Element Model for Sea Ice Breakup
Abstract. Current climate models simulate sea ice as a continuous medium despite the inherent discontinuities at the scale of individual floes. Discrete element models (DEMs) can directly represent these discontinuities, making them valuable tools for understanding subgrid-scale sea ice dynamics. However, DEMs are typically limited by their high computational cost and large parameter spaces, which hinder their validation and interpretability. In this work, we developed, calibrated, and evaluated the performance of a low-complexity, two-dimensional bonded DEM for process-based studies of ice flow and fracture. We modified the linear bonded particle model implemented in the molecular dynamics software, LAMMPS (Large-Scale Molecular Massively Parallel Simulator), to prevent failure in compression because undamaged sea ice is significantly stronger under compression than tension. Through an ensemble of simple shear experiments, we related the effective elastic modulus and shear strength to the bond stiffness, critical strain, and particle size. We further observed that the simulated ice deformation was localized in space and intermittent in time, in agreement with observations. Using a canonical geometry idealizing sea ice flow through the Nares Strait, we demonstrated that the model could reproduce ice arch formation and collapse previously observed in this region and reproduced in higher-complexity models. We found that our model could represent four possible regimes: no fracturing, stable ice arch formation, the formation of an ice arch that subsequently collapses, and the absence of an ice arch buttressing flow through the channel. These regimes collapsed onto a single control parameter given by the product of bond stiffness and critical strain. In summary, this low-complexity DEM offers a computationally efficient framework for evaluating short-timescale, discrete sea ice processes, which may be used in conjunction with more comprehensive models to better understand key ice deformation processes.
- Preprint
(46393 KB) - Metadata XML
- BibTeX
- EndNote
Status: open (until 13 Oct 2026)
- RC1: 'Comment on egusphere-2026-3166', Anonymous Referee #1, 04 Oct 2026 reply
-
RC2: 'Comment on egusphere-2026-3166', Anonymous Referee #2, 07 Oct 2026
reply
The authors present a discrete element model of sea ice without tangential forces and show numerically that at least some of the macroscale sea ice dynamics can be described using granular mechanics. While the latter idea is not new, such numerical demonstrations are welcome, and the paper is interesting. The paper would, however, benefit from engaging with literature beyond sea ice, particularly on granular materials and discrete element modelling, where many of the modelling choices made here have been studied extensively. Further rigour in terms of solid and continuum mechanics is also needed, as spring stiffnesses and critical strains are treated as material properties, which they are not. Comparisons to the same model with tangential interactions would demonstrate the claimed advantages of the approach and strengthen the conclusions, although this reviewer understands that this may require extensive additional work. This reviewer considers that the manuscript requires major revisions but is likely to become a nice contribution to the sea ice community. The main comments are:
Computational efficiency: The authors claim several times that their model requires only a fraction of the computational cost of a model with particle orientations. This claim remains unsubstantiated and authors can simply profile their code and present data on efficiency. As all particles are circular and all interactions are central, no torques arise, and integrating rotations would not change the results; omitting them therefore saves only a negligible amount of computation. A model accounting for orientation would also include tangential interactions such as friction, so any savings stem from omitting tangential force calculations, which changes the physics of the model rather than merely its cost. Moreover, in this reviewer's experience, contact- and bond-related distance and overlap calculations dominate the computational time in DE models, while time integration takes only a very small fraction. The authors should substantiate the claim, at minimum by profiling their code, or otherwise revise or soften it.
Absence of tangential interactions: Interparticle friction and resistance to particle rotation are known to be key controls of the bulk shear strength of granular materials, and rolling resistance is commonly introduced to represent particle shape effects while retaining the low computational cost of circular or spherical particles. Frictionless circular particles, in contrast, offer little shear resistance, so the shear strength of the model must arise mainly from the bonds. Once bonds fail, the broken ice behaves as a frictionless granular material, whereas sea ice is known to be a frictional material. Tuning other parameters may yield a strength representative of sea ice, but these parameters then partly compensate for the missing tangential interactions and lose their physical interpretation, making them difficult to constrain. The authors should place their approach in the context of the granular materials literature and discuss its implications for parameter calibration.
Particle size dependence: The macroscopic shear modulus appears independent of particle size (Fig. 3), and the same largely applies to strength. As the number of bonds crossing a given section scales inversely with particle size, this independence depends on how the bond stiffness and failure criterion are defined with respect to particle size, for example whether the bond cross-sectional area is implicitly scaled with particle size. The authors should clarify this explicitly and relate k_b to an effective modulus (see the comment on the failure criterion). Demonstrating the behaviour of samples with different particle sizes under uniaxial tension would also help the reader understand the mechanical behaviour of the model.
Failure criterion: The authors use a critical strain to define bond failure and recommend its use. This choice should be discussed further. Critical strain is not a fundamental material parameter and cannot be constrained independently of stiffness from material experiments. The same applies to the spring stiffness k_b, for which an effective (Young's) modulus would be the appropriate material parameter. The relevant failure parameter is the bond strength, here k_b·ε_c, which in a 2D model can be related to the macroscopic (vertically integrated) strength. The authors' own results support this: as stated in the abstract, the results fall on a master curve defined by k_b·ε_c, indicating that bond strength governs the behaviour of the modelled ice. This is to be expected, as it defines the load-carrying capacity of the bonds. The paper currently alternates between discussing k_b and ε_c separately and discussing strength. This reviewer recommends that the authors parameterise failure in terms of bond strength (critical force or stress), or at least revise their analysis, discussion and recommendations accordingly.
Arching: The arching observed in the simulations is a well-known phenomenon in granular materials, and relating the results to this literature would strengthen the paper. A comparable geometry is granular flow through a silo outlet, where non-cohesive materials are known to form stable arches in which contact forces are considerably larger than in the surrounding force chains. Illustrating the contact forces in the arching cases would be beneficial, as these are readily available from the simulations. Arch stability is known to depend on interparticle friction and particle shape, so arching with frictionless circular particles warrants explanation: is it sustained by intact bonds, by boundary conditions, or by geometry alone? The results suggest that bond strength dominates. A comparison to a model with tangential interactions would clarify this, and if the arching behaviour differed substantially, this would be a strong demonstration of the need to account for friction and floe shape in sea ice models.
More detailed comments (as the revisions requested above will likely lead to considerable changes in the manuscript, the comments below are brief and focus on other topics):
- It would add clarity to state clearly from the beginning that the model contains no tangential forces or torques that would cause particle rotations. Orientations are then not merely left unintegrated; integrating them is unnecessary.
- L35-43: The authors omit the extensive ice engineering literature in which bonded DEM has been used for decades. Contrary to the statement here, small-scale DEM simulations of ice can be parameterised using measured material properties, as done in recent ice engineering studies. Large-scale DEM does not share this advantage.
- L68-69: It is not clear how only half of the computations are needed. As no torques arise, omitting rotations saves only a negligible amount of computation, and the most demanding contact and bond calculations are still required. Any further savings arise from omitting tangential forces, which changes the model itself. The reader should also be reminded that the absence of tangential forces is due to the model not including friction.
- L83: "Sensitive to initial particle assembly": this holds for the details of the response, but the overall behaviour and statistics should be insensitive to the initial assembly.
- L97-99: Please provide an illustration or close-up of the bonded specimens. The description suggests that bonds may be created between particles that are not immediate neighbours.
- L104: Suggest "magnitude" rather than "strength", as strength refers to a material property.
- L105-106: Why is the critical damping coefficient not defined for each particle pair using the effective mass, which would be the more standard approach? Currently the damping ratio varies throughout the simulated field, although this is likely not an issue in quasi-static simulations.
- L108-110: Is shear failure under compression allowed? This is a very common failure mode for sea ice. The use of critical strain could also be addressed here.
- Eq. 3: The model is described as minimalistic, yet the contact model includes several parameters (various stiffnesses and a smoothing factor) not commonly used in ice modelling. Can these be measured, or are they tuning parameters? Also, the last term is a scalar, not a vector.
- L156-159: Could the authors elaborate on why they want to avoid dilatancy?
- L220-222: The authors apply shear deformation, not shear forcing. If this reviewer understands correctly, the model does not allow shear failure. Is the local failure therefore due to tension?
- L225-227: It is not immediately clear why G (and in practice also τ) is constant for all Λ. This should be explained to the reader (see the main comment on particle size dependence).
- Fig. 4: This figure illustrates the challenge of using critical strain as the failure criterion, as the slope of the fitted line depends on both k_b and ε_c. It would be informative to vary k_b and ε_c while keeping their product constant, to show whether τ_max depends on the local deformation response or only on the local strength.
- L246-247: Based on the discussion above, this reviewer does not consider it justified to recommend critical strain for defining failure.
- L338-339: The sentence "The symmetry about the line…" in fact states that the regimes depend on the bond strength.
- Sec. 3.2: This section largely focuses on arching, discussed above. Relating it to the granular materials literature may affect how the results are presented. As strength again appears to be the dominant parameter, the separate discussion of k_b and ε_c could be condensed.
- L361-365: This result in particular would benefit from discussion in light of earlier work on granular materials.
- L394-395: This result is to be expected, as stiffness governs the elastic response and the failure criterion governs failure. The authors could consider shortening or rephrasing this statement.
- L407-409: "Smaller parameter space makes calibration tractable": this holds only if the parameters can be measured independently, and it is not immediately clear how this could be done with the chosen parameters.
- L421-426 and L434-436: This reviewer believes the local response in the model accounts for ridging, although not for ice thickening. As the contact stiffness applies at a larger scale than in some of the cited models, it is natural for it to be lower, as it reflects not only the material response but also, for example, the stiffness of ridges at particle interfaces.
Citation: https://doi.org/10.5194/egusphere-2026-3166-RC2
Data sets
Data for: A Minimal Bonded Discrete Element Model for Sea Ice Breakup Alexis Arlen https://doi.org/10.5281/zenodo.20435014
Model code and software
lexiarlen/P1-model-dev: v1.0.2-submission: Sea Ice DEM analysis code Alexis Arlen https://doi.org/10.5281/zenodo.20434932
Custom LAMMPS code for bonded discrete element sea ice simulations Lexi Arlen https://doi.org/10.5281/zenodo.20434107
Video supplement
Video Supplement for: A Minimal Bonded Discrete Element Model for Sea Ice Breakup Alexis Arlen https://doi.org/10.5281/zenodo.20435748
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 340 | 119 | 53 | 512 | 51 | 46 |
- HTML: 340
- PDF: 119
- XML: 53
- Total: 512
- BibTeX: 51
- EndNote: 46
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1
See attached pdf with comments.