# Fluid-structure modeling of iceberg capsize

Nicolas de Pinho Dias, Alban Leroyer, Anne Mangeney, Olivier Castelnau

► To cite this version:

Nicolas de Pinho Dias, Alban Leroyer, Anne Mangeney, Olivier Castelnau. Fluid-structure modeling of iceberg capsize. Ocean Engineering, 2025, 336, pp.121765. 10.1016/j.oceaneng.2025.121765. hal-05161111

HAL Id: hal-05161111

https://hal.science/hal-05161111v1

Submitted on 16 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

Public Domain
Fluid-structure modeling of iceberg capsize

Nicolas De Pinho Dias ID a,b,*, Alban Leroyer ID b, Anne Mangeney ID a, Olivier Castelnau ID c

ARTICLE INFO

Keywords:

Computational fluid dynamics

Fluid-structure interactions

Glacial earthquake

Glacier ice-mass-loss

Iceberg capsize

A B S T R A C T

Ice sheet evolution models require consideration of ocean-glacier processes such as iceberg calving. The capsize of large icebergs from marine-terminating glaciers as a result of calving generates teleseismic waves that contain information about the physical processes involved. Deriving the calved iceberg volume from the seismic signal requires the coupling of seismic inversion with numerical simulation of the capsize, which currently lacks high-fidelity hydrodynamic effects. Therefore, a Computational Fluid Dynamics model combining a Reynolds Averaged Navier-Stokes Equations solver with a volume penalization method and a spring-damper contact force formulation is used to simulate the capsize of an iceberg colliding with a glacier front. The modeled capsizes are in good agreement with laboratory-scale experimental results from the literature, for many configurations (capsize in the open ocean or in contact with a floating or grounded glacier terminus, for bottom-out or top-out configurations and several iceberg aspect ratios). In grounded glacier configurations, we show that the lateral flow confinement effects are significant. The initial iceberg tilt angle and the iceberg/glacier friction coefficient have little influence. Turbulence mainly affects the post-capsize iceberg drift velocity. In addition, a field-scale 200 m high iceberg capsize simulation shows near-field water velocities in excess of 10 m/s.

# 1. Introduction

Global warming is the cause of increasing melting of polar ice sheets. Ice sheets contain about 70% of Earth fresh water and are a major contributor to sea level rise. More specifically, the contribution of Greenland is estimated to be between 20 and 37% of sea level rise, corresponding to 0.47 to 0.77 mm/y (van den Broeke et al., 2016; Bamber et al., 2018; Zemp et al., 2019; Shepherd et al., 2020; IPCC, 2022). Several processes are involved in the mass loss of the Greenland ice sheet such as surface melting (Zwally et al., 2002; Mote, 2007; Zheng et al., 2022). Subglacial runoffs accelerate the glacier sliding by reducing friction between the ice and the bedrock (Nick et al., 2012; Seddik et al., 2019) and provide an important subglacial heat source which increases ice melting (Young et al., 2022). At marine-terminating glaciers, the ice mass loss is mainly due to submarine melting (Rignot et al., 2010; Haseloff and Sergienko, 2022; Slater and Straneo, 2022) and iceberg calving, i.e. the process through which pieces of ice break off marine-terminating glaciers (van der Veen, 1996; Amundson et al., 2008, 2010). The contributions of sub-marine melting and iceberg calving to Greenland ice mass loss are not easily determined and these processes are related to each other (Wagner et al., 2019). A recent study estimated that ice mass flux into the

ocean, which includes iceberg calving and submarine melting, accounts for 66% of Greenland ice mass loss (Mouginot et al., 2019).

Multiple factors are responsible for the stresses leading to iceberg calving. On top of seasonality (Nettles and Ekström, 2010; Nick et al., 2012; Bartholomaus et al., 2015; Gräff and Walter, 2021), studies highlight heat transfers linked to oceanic (Walter et al., 2012; Sugiyama et al., 2014; Truffer and Motyka, 2016; Todd et al., 2019) and atmospheric circulations (Hanna et al., 2011, 2013; Lea et al., 2014; Mattingly et al., 2018; Wang et al., 2021; Slater and Straneo, 2022). Mechanical sources of stress such as the hydrostatic pressure (Hanson and Hooke, 2000; James et al., 2014; Wagner et al., 2016), enhanced by tidal fluctuations leads to crevasse opening and eventually to ice breaking (Seddik et al., 2019; van Dongen et al., 2019; Christmann et al., 2021; Kneib-Walter et al., 2022). Crevasses which opened up-glacier due to bedrock topography are advected by glacier flow and produce large calvings (Jouvet et al., 2017; van Dongen et al., 2020; Berg and Bassis, 2022). Subglacial channels delivering freshwater to the ocean (plume) can form an undercutting at the base of the glacier front (van Dongen et al., 2020). Due to this undercutting, the glacier is thinner and crevasses are more likely to propagate over the whole glacier reduced thickness, leading to a calving (Benn et al., 2023). Overall, it seems that

<footnote>a Institut de Physique du Globe (IPGP), Université Paris Cité, 1 Rue Jussieu, Paris, 75005, France</footnote>

<footnote>b Nantes Université, École Centrale de Nantes, CNRS, LHEEA, UMR 6598, 1 Rue de la Noë, Nantes, 44000, France</footnote>

<footnote>c PIMM, CNRS, Arts et Métiers Science and Technology, CNAM, 155 Bd de l'Hôpital, Paris, 75013, France</footnote>
icebergs tend to calve near the grounding line of their glacier (van der Veen, 1996; van der Veen, 2002; James et al., 2014; Benn et al., 2023).

One can find different types of calvings depending on the size of the
iceberg and whether it is submarine (i.e. a piece of ice detaching from
a submerged part of the glacier front and rising to the water surface) or
aerial (i.e. a piece of ice detaching from an emerged part of the glacier
and falling into the water) (van der Veen, 2002; Glowacki et al., 2015;
Sugiyama et al., 2019; Glowacki, 2022). The type of calving is influenced
by the shape of the iceberg, its position when the fracture occurs, the
angle of the iceberg with respect to the vertical axis and the competition
between buoyancy and gravity. The presence of ice-melange can also
inhibit top-out calvings by strengthening the upper part of the glacier
(Walter et al., 2012).

These events can be followed by a capsize of the iceberg under the
torque created by the competition between gravity and buoyancy forces.
Capsize criteria have been established mathematically for parallelepiped
iceberg (MacAyeal et al., 2003; Burton et al., 2012; Wagner et al., 2017).
Depending on which edge remains in contact with the glacier during
the process, the capsize is called bottom-out (the bottom part of the
iceberg moves away from the calving front) or top-out (the top part
of the iceberg moves away from the calving front) Tsai et al. (2008);
Amundson et al. (2010).

Different measurement techniques have been developed to observe
iceberg capsizes and associated glacier dynamic response. Some rely
on imagery to count events and estimate iceberg volume: satellite, un-
manned aerial vehicle (UAV), time-lapse camera (Jouvet et al., 2017;
Minowa et al., 2018; Bézu and Bartholomaus, 2024) but depending on
the equipment used, the snapshot frequency can be as long as a few
days. Laser interferometry and scanning can also be used but like im-
agery techniques, they depend on the weather (Cook et al., 2021b). As
iceberg calvings make loud noises, Glacial EarthQuakes, denoted GEQ
(Podolskiy and Walter, 2016; Podolskiy et al., 2022) and pressure waves
that can be recorded under water using ocean bottom seismometers and
hydrophones (Podolskiy et al., 2021, 2022).

An alternative strategy is to use the data from the global seismic
network data recorded by stations in Greenland (GLISN, IRIS) to de-
tect glacial earthquakes. Podolskiy and Walter (2016) reviewed a few
cryo-seismic events occurring on glaciers and sorted them by their main
frequency content. The authors show that GEQ contain frequencies rang-
ing from 0.001 to 0.1 Hz while icequakes (i.e. vibrations generated by
fracture in the ice) contains frequencies above 0.1 Hz. The difference
in frequency content allows to identify and distinguish the processes
at stake. Sergeant et al. (2016) also shows, on the case of a capsize
at Jakobshavn glacier, that the ground displacement in the frequency
band of ice fracture is about one order of magnitude smaller than the
ground displacement in the GEQ frequency band. Therefore, the avail-
able 30 years of data can give us information about iceberg calving and
its evolution in time (Tsai and Ekström, 2007; Veitch and Nettles, 2012;
Olsen and Nettles, 2017). Seismic inversion is a way to decipher the
seismic signal of GEQ that has been previously explored by Tsai and Ek-
ström (2007); Veitch and Nettles (2012); Sergeant et al. (2016); Veitch
and Nettles (2017); Sergeant et al. (2018, 2019); Olsen et al. (2021).
The assumption made is that the iceberg exerts a time-varying contact
force on the glacier which generates seismic waves during the capsize.
The pressure field and water waves generated by the iceberg motion
also produce a fluctuating force on the glacier and on the fjord bed.
The recorded seismic signal is a mix of all phenomena occurring dur-
ing the capsize. Each of the time-varying forces, convolved with Earth
Green's functions (i.e. Earth impulse response), results in the recorded
seismic signal. For each individual glacier, the comparison between
forces inverted from 30 years of seismic signal and simulated forces al-
low to infer the calved-iceberg volume. As a result, it will be possible
to know the spatio-temporal evolution of the Greenland ice-mass-loss
through iceberg calving Sergeant et al. (2019). It will then be possible
to link this evolution to the glaciers dynamics and to external forcings
described above (oceanic and atmospheric circulation, heat transfers,

mechanical sources of stress). However, Sergeant et al. (2019) and Olsen
et al. (2021) observe discrepancies in shape and intensity between syn-
thetic forces (from numerical modeling and laboratory experiments, re-
spectively) and forces inverted from the seismic signal. The model used
by these authors only includes added masses and a simplistic drag force
to model the effect of water on the iceberg dynamics. Adding hydrody-
namic pressure forces within the whole water domain could change the
dynamical response of both the iceberg and the glacier.

In the case of large icebergs of cubic-kilometer (and gigatons) scale,
calving-induced stress variations have been suggested to create fractures
in the ice, weakening the glacier and potentially leading to more calv-
ings (Murray et al., 2015c; Amundson et al., 2022). Another feedback in-
volves ice melange, which is a mixture of ice debris floating on the ocean
and that acts like a cohesive material. Thus, the ice-melange supports
the glacier and helps in preventing iceberg calving (Walter et al., 2012;
Cassotto et al., 2015, 2021). However, a recent study reports that ice-
berg calving weakens melange, which, in turn, cannot fortify the glacier
and leads to more calvings (Wehrlé et al., 2023).

Current models of glacier front evolution (from a few hours to up
to hundreds of years) including calving focus on the stresses acting on
the glacier and include only the buoyancy component or heat exchanges
due to plumes or ocean circulation (Gagliardini et al., 2013; Cook et al.,
2021a; Jones et al., 2021). Thus they do not take into account the Fluid-
Structure Interaction (FSI) which involves a complex and turbulent flow
during the few minutes of a capsize.

Previous capsize models are either experimentally based (Amundson
et al., 2012b; Burton et al., 2012; Murray et al., 2015c) or simplified ana-
lytical model including added mass (Sergeant et al., 2018, 2019; Bonnet
et al., 2020). Their results reasonably match high-fidelity Computational
Fluid Dynamics (CFD) modeling after parametrization for some simple
geometrical configuration, like in an open ocean (Bonnet et al., 2020)
but do not reproduce the end of the capsize dynamics and can hardly be
generalized to any glacier and bedrock geometries. Therefore, the use
of CFD is necessary. Bonnet et al. (2020) used CFD to compute iceberg
capsize in open ocean and obtained results consistent with laboratory
experiments of Burton et al. (2012). The next step is to add a contact
between the glacier front and the iceberg.

Several studies were conducted on the tsunami generated by the
calving of icebergs from simple analytical capsize models (Levermann,
2011; MacAyeal et al., 2011) to experiments and numerical simulations
of falling/capsizing icebergs (Heller et al., 2019; Lüthi and Vieli, 2016;
Chen et al., 2020; Wolper et al., 2021). However, these studies focus on
the far-field hydrodynamics at the water surface and disregard the near-
field and deep-water effects. Tsunamis can also provide information on
the calving rate (Minowa et al., 2019).

Depending on the formulation used, solid-solid (iceberg-glacier) con-
tacts can be implemented more or less naturally in numerical models,
even with different laws for the type of contact. Using the finite ele-
ment method, Fluid-Structure Interactions (FSI) studies have proposed
a contact model with different boundary conditions in the arteries (Bur-
man et al., 2020, 2022). Other methods are based on particles inter-
acting through the forces they exert on each other. This is the case of
Smooth Particle Hydrodynamics (Heller et al., 2016), Discrete Element
Method (Crawford et al., 2021) and Material Point Method (Wolper
et al., 2021). Immersed boundary methods (Peskin, 2002) are quite com-
mon in FSI and have been used to simulate tsunami generation by ice-
berg fall (Chen et al., 2020). The model of Chen et al. (2020) includes
water viscosity, hydrodynamic pressure and solves Reynolds-averaged
Navier-Stokes equations (RANSE) but do not focus on the near-field ice-
berg area where all the glacier-iceberg-ocean interactions take place, see
Murray et al. (2015c).

In this study, we have implemented a solid-solid/iceberg-glacier con-
tact model (see Fig. 1) in the CFD finite-volume solver used by Bonnet
et al. (2020). Main novelties of this work include:

* a glacier/iceberg contact model in a finite-volume solver
<figure><img src="image_1.png" /><figcaption>**Fig. 1.** 3D sketch of a large iceberg capsizing bottom-out against a glacier terminus. The system is symmetrical and the symmetry plane (in red) corresponds to $y = 0$, therefore cutting both the glacier and the iceberg in half. $L$, $W_i$ and $H$ are the iceberg length, width and height, respectively and we define the aspect ratio $\epsilon = L/H$. $W_g$ is the glacier width and $\Delta W = (W_g - W_i)/2$ is the distance between the fjord/tank lateral side and the iceberg lateral face. $\theta(t)$ is the iceberg rotation angle. When $\theta = 0$, the iceberg is parallel to the glacier front and covers the area delimited by the dashed lines on the glacier wall. $\Delta z_{bot}$ is the distance between the iceberg bottom and the bedrock at $t = 0$. As a first step we assume the bedrock to be flat and horizontal. $d_w$ stands for water depth (at rest). The water flow paths is represented in purple in the case of a 2D computation and in orange in the case of a 3D computation.</figcaption></figure>

* numerical simulations of lab-scale experiments of iceberg capsize against a glacier terminus

* an estimate of near-field and deep-water flow resulting from the capsize

The method is presented in Section 2. Our model is validated against laboratory experiments from the literature (Amundson et al., 2012b; Burton et al., 2012; Murray et al., 2015c) presented in Section 3. Effects of numerical parameters are discussed in details in Section 4. The influence of physical parameters is quantified in Section 5. Section 6 provides a detailed comparison to the experimental results. Then, an application to a field-scale situation is presented in Section 7 (Table 1).

## 2. Numerical methods

### 2.1. Fluid/structure interaction model

We use the numerical solver ISIS-CFD developed at Ecole Centrale de Nantes (France) which is integrated and distributed by Cadense Design Systems as the software suite Fine/Marine. Initially designed for naval hydrodynamic purposes, it contains features to model FSI. ISIS-CFD can solve unsteady incompressible Reynolds-Averaged Navier-Stokes Equations with numerous turbulence models including the widely used $k-\omega$ model (Menter, 1994). It is based on a finite-volume method with an unstructured mesh to handle complex geometries.

Pressure-velocity coupling is performed using a Rhie and Chow formulation (Rhie and Chow, 1982) and a SIMPLE-type algorithm (Patankar, 1980). An interface capturing method (Queutey and Visonneau, 2007) allows to account for air/water interactions. An Arbitrary Lagrangian–Eulerian (Leroyer et al., 2008) formulation allows for deformation of the fluid domain during the iceberg motion. In unsteady

configuration, the time discretization is performed with a Backward Differentiation Formula of order 2 (BDF2). Each time step contains an inner loop of about 20 iterations (called nonlinear loop) used to treat nonlinearities of the system with a Picard linearization method. The fluid-structure coupling is implemented in the nonlinear loop where the solid dynamics is solved and coupled to the fluid through a physically-based relaxation using an evaluation of the added-mass coefficients updated every five time steps (Yvin et al., 2018). ISIS-CFD is parallel and uses the Message Passing Interface (MPI) protocol. Computations were performed on the GLiCID supercomputer (Ecole Centrale de Nantes, France).

The considered geometry of the capsizing iceberg is illustrated in Fig. 1. The tank bottom is a rigid solid with a slip condition (no shear stress $\tau_{wall} = 0$ Pa). The top of the tank is subjected to an imposed hydrostatic pressure such that, at initialization, $p = 0$ at the water free surface. The iceberg is rigid with a boundary condition of type “wall function”. This boundary condition are used to model the near-wall region without fully resolving the viscous sub-layer, which would require a very fine mesh. It reduces computational cost while maintaining acceptable accuracy by approximating velocity profiles and other quantities near walls based on the known log region properties in turbulent boundary layer. The tank is large enough to avoid reflections of water waves generated by the capsize. In addition, far-field boundary conditions and elongated cells at $x \to \infty$ provide a good numerical wave damping (see Appendix A, Fig. A.18). The glacier wall ($x = 0$) is modeled as a porous medium described in Section 2.2. Meanwhile, the tank boundary is located, depending on the mesh, in $x < 0$ ($x \approx -0.02$ m on Fig. 2) and slip boundary conditions are applied on this surface. Computations are performed for 2D and 3D geometries. In 3D, the plane $y = 0$ is a symmetry plane and only half of the space is considered ($y > 0$),
<table><caption>**Table 1**
Symbols and abbreviations used in this paper.</caption><thead><tr><th>Notation</th><th>Quantity</th><th>Unit</th></tr></thead><tbody><tr><td>a</td><td>"stiffness" of the spring-damper contact force model</td><td>—</td></tr><tr><td>A</td><td>Area of the iceberg perpendicular to the flow during the post-capsize drift</td><td>m<sup>2</sup></td></tr><tr><td>b</td><td>"damping" of the spring-damper contact force model</td><td>—</td></tr><tr><td>d<sub>w</sub></td><td>Water depth</td><td>m</td></tr><tr><td>F<sub>glac→iceb</sub></td><td>Contact force exerted by the glacier on the iceberg</td><td>N</td></tr><tr><td>F<sub>r</sub></td><td>Froude number</td><td>—</td></tr><tr><td>F<sub>W</sub></td><td>Iceberg weight</td><td>N</td></tr><tr><td>g</td><td>Acceleration of the gravity</td><td>m/s<sup>2</sup></td></tr><tr><td>h</td><td>Cell size (see Appendix B)</td><td>m</td></tr><tr><td>H</td><td>Iceberg height</td><td>m</td></tr><tr><td>k</td><td>Parameter used in the VPM source term</td><td>—</td></tr><tr><td>L</td><td>Iceberg length</td><td>m</td></tr><tr><td>M3, M4, M5</td><td>Levels of mesh refinement</td><td>—</td></tr><tr><td>n</td><td>Exponent used in the VPM source term</td><td>—</td></tr><tr><td>N<sub>x</sub>, N<sub>y</sub>, N<sub>z</sub></td><td>Positive integers used to build the mesh before refinement (see Appendix B)</td><td>—</td></tr><tr><td>P<sub>dyn</sub></td><td>Hydrodynamic pressure</td><td>Pa</td></tr><tr><td>r<sub>d</sub></td><td>Refinement diffusion (see Appendix B)</td><td>—</td></tr><tr><td>Re</td><td>Reynolds number</td><td>—</td></tr><tr><td>S<sub>VPM</sub></td><td>Volume Penalization Method source term</td><td>N/m<sup>3</sup></td></tr><tr><td>t</td><td>Time coordinate</td><td>s</td></tr><tr><td>t'</td><td>Non-dimensional time coordinate</td><td>—</td></tr><tr><td>T<sub>x,y,z</sub></td><td>Iceberg translation along axis x, y, z</td><td>m</td></tr><tr><td>T<sub>x,y,z</sub></td><td>Iceberg translation along axis x, y, z of the most left corner (mlc)</td><td>m</td></tr><tr><td>T<sub>x,y,z</sub><sup>mlc,slc</sup></td><td>or top left corner (tlc) of the iceberg</td><td>m</td></tr><tr><td>u<sub>x,y,z</sub></td><td>Water velocity along axis x, y, z</td><td>m/s</td></tr><tr><td>u<sub>x,y,z</sub><sup>mlc,slc</sup></td><td>Velocity along axis x, y, z of the most left corner (mlc)</td><td>m/s</td></tr><tr><td>V<sub>S</sub></td><td>Submerged volume of the iceberg</td><td>m<sup>3</sup></td></tr><tr><td>V<sub>iol</sub></td><td>Total Volume of the iceberg</td><td>m<sup>3</sup></td></tr><tr><td>x<sub>1</sub>, y<sub>1</sub>, z</td><td>Usual cartesian coordinates</td><td>m</td></tr><tr><td>x<sub>ivall</sub></td><td>Position of the glacier wall</td><td>m</td></tr><tr><td>x<sub>ivall</sub><sup>mlc,slc</sup></td><td>Position along axis x of the most left corner (mlc)</td><td>m</td></tr><tr><td>W<sub>z</sub></td><td>Glacier width</td><td>m</td></tr><tr><td>W<sub>t</sub></td><td>Iceberg width</td><td>m</td></tr><tr><td>δ</td><td>name given to the simulation described in Appendix D</td><td>—</td></tr><tr><td>Δt</td><td>Numerical time step</td><td>s</td></tr><tr><td>ΔW</td><td>Distance between the fjord/tank lateral side and the iceberg lateral face</td><td>m</td></tr><tr><td>Δz<sub>bot</sub></td><td>Vertical distance between the iceberg bottom and the bedrock at t = 0</td><td>m</td></tr><tr><td>c</td><td>Aspect ratio of the iceberg</td><td>—</td></tr><tr><td>η</td><td>Water surface elevation</td><td>m</td></tr><tr><td>θ</td><td>Iceberg rotation angle with respect to the vertical</td><td>°</td></tr><tr><td>θ<sub>0</sub></td><td>Iceberg initial rotation angle with respect to the vertical</td><td>°</td></tr><tr><td>μ<sub>C</sub></td><td>Coulomb friction parameter</td><td>—</td></tr><tr><td>ρ<sub>l</sub></td><td>Iceberg density (for plastic or ice depending on the case)</td><td>kg/m<sup>3</sup></td></tr><tr><td>ρ<sub>w</sub></td><td>Water density</td><td>kg/m<sup>3</sup></td></tr><tr><td>τ</td><td>Time scale τ = √H/g</td><td>s</td></tr><tr><td>τ<sub>0</sub></td><td>Time scale used in the VPM source term</td><td>s</td></tr><tr><td colspan="3"><strong>Abbreviations</strong></td></tr><tr><td>BO</td><td>Bottom-Out</td><td></td></tr><tr><td>CFD</td><td>Computational Fluid Dynamics</td><td></td></tr><tr><td>CoG</td><td>Center of Gravity</td><td></td></tr><tr><td>CPU</td><td>Central Processing Unit</td><td></td></tr><tr><td>FG</td><td>Floating Glacier</td><td></td></tr><tr><td>FSI</td><td>Fluid-Structure Interaction</td><td></td></tr><tr><td>GEQ</td><td>Glacial EarthQuake</td><td></td></tr><tr><td>GG</td><td>Grounded Glacier</td><td></td></tr><tr><td>OO</td><td>Open Ocean</td><td></td></tr><tr><td>TO</td><td>Top-Out</td><td></td></tr><tr><td>VPM</td><td>Volume Penalization Method</td><td></td></tr></tbody></table>

see Fig. 1. Thus, mirror boundary conditions are applied on the plane $y = 0$ and “wall function” boundary conditions are applied on the other lateral side.

Since an iceberg capsize is a situation including two solids (glacier and iceberg) of which the relative rotation can exceed 90°, we use overlapping grids (overset) as shown in Fig. 2. Two domains are defined: one for the water tank/fjord which is fixed and one attached to the iceberg (respectively in white and blue on Fig. 2). Information crossing these two domains is interpolated and numerical errors due to interpolation are reduced thanks to adaptive mesh refinement (Wackers et al., 2012).

This technique is applied at the interface between the two numerical domains but also in the free surface vicinity and allows to perform accurate computations at a reasonable cost.

Despite adapting the time step to one of the two Courant numbers (one for the free surface and one for the overset), a particularly complex and unstable location is where the air-water interface meets the overset one. We noticed that numerical instabilities could appear in this region but they are quickly dissipated and neither influence the iceberg dynamics nor the water flow.
<figure><img src="image_2.png" /><figcaption>Fig. 2. Snapshot of mesh domains (iceberg at 45° here). The rigid wall is represented in red ($x < 0$). The iceberg solid body is in white. The tank domain is in white with the mesh. The iceberg-related domain is shown in blue. The lighter blue corresponds to the overlapping layer between iceberg and tank domains. The free surface (initially at $z = 0$) is shaded in dark. Refined mesh areas are the wall, the ground ($z = -0.2$ m here), the iceberg contour and the free surface. Some streamlines are also represented in black lines, they highlight the flow deviation due to the wall presence (volume penalization method).</figcaption></figure>

Refined areas include the iceberg contour, the glacier wall ($x = 0$ plane) and the tank bottom in the vicinity of the iceberg, where the generated flow structures can have a large impact. In order to capture the fluid dynamics and permit a good interpolation between the over-set grids in narrow spaces (under the iceberg), a sufficient number of cells must be placed between the iceberg and the tank bottom boundary. For this reason, decreasing $\Delta z_{bol}$ while keeping $H$ (the reference length of our simulations) constant triggers the need to refine more the region under the iceberg. Thus it leads to more computationally expensive simulations.

Overall, the number of cells ranges from about $2 \times 10^4$ (open ocean in 2D) to about $3 \times 10^6$ for the most refined simulations presented in this paper (Grounded Glacier in 3D, see Section 7). Cases of icebergs with large aspect ratios require more cells to cover the whole iceberg perimeter/surface. In the cases involving narrow lateral or vertical spaces $\Delta W$ and $\Delta z_{bol}$, a large number small cells is necessary to capture the fine details. These features increase the computational cost.

As initial conditions, we impose a fluid domain at rest. The iceberg of height $H$ is tilted to a small angle $\theta(t = 0s) = \theta_0$ with respect to the vertical axis and placed in hydrostatic position with one of its corners in contact with the glacier depending on the configuration: bottom-out (see Figs. 1 and 2) or top-out (the iceberg's bottom left corner is in contact with the wall while its top left corner is moving away from the wall). Values of $\theta_0$ depend on the study case and are reported in the corresponding results section.

The CPU time ranges from 6 h in 2D to more than 96 h in 3D but is very dependent on the case studied, the average computation time is about 20 h on the GLiCID supercomputer.

## 2.2. Solid-solid contact with large motion

The iceberg is considered to be a rigid body and follows Newton's second law of dynamics. Forces acting on the iceberg are gravity, buoyancy and the water/air pressure forces, computed at each time step.

Modeling the glacier wall, however, raises difficulties because solid-solid contact is still an open question in the finite-volume framework. By combining two numerical techniques, we model the contact between a fixed wall and a moving iceberg (see Fig. 2). The first one is a Volume Penalization Method (VPM), which is quite common in fluid-structure interactions (Hou et al., 2012; Engels et al., 2015; Hester et al., 2021). It mimics the behavior of a porous medium and is based on the Darcy-Brickmann equations. The wall domain ($x < 0$) is filled with fluid that has velocities close to zero. This corresponds to the red region in the Fig. 2. Note that the streamlines represented on this figure are indeed deflected by the wall. This is simply done here by adding, to the momentum equations, a source term of the form

$$S_{VPM,i} = -\rho_w \frac{u_i}{\tau_0} \text{ with } \tau_0 = t_0 k \left( \frac{\Delta t}{t_0} \right)^n, \quad (1)$$

with $u_i$ the flow velocity component in direction $i = x, y, z, k$ a constant parameter ($k = 1$ in our case), $\Delta t$ the current time step, $t_0 = 1$ s and $n$ a parameter. This technique works on the condition that $S_{VPM}$ is large compared to other momentum terms. Since $\Delta t$ is small and $k$ is fixed, the condition leads to “n is large enough”. In practice $n = 1$ is the limit between a porous (the fluid can still move inside the medium) and solid behavior (fluid being immobile). However, if this term is too large, the computation can become unstable, even if the additional term is treated implicitly.

Note that the time step is variable, unknown a priori and depends on the Courant number. Therefore, there is no rule for choosing the $n$ value. We used $n = 1.5$ for the 2D cases and $n = 1.2$ for the 3D cases which was found to be a good compromise between efficiency and convergence, see Section 4.3 for a sensitivity study.

In our model, there is no way to prevent the iceberg from falling into the wall. However, it is possible to impose a force on the iceberg that will keep its corner close to the wall edge. We therefore impose a spring-damper contact force on any iceberg mesh node located in the negative-$x$ domain, i.e. the wall domain, given by Eq. (1)

$$\frac{\vec{F}_{glac \to iceb} \cdot \vec{x}}{F_W} = \begin{cases} a \frac{|x^{mlc} - x^{wall}|^2}{H^2} - b \frac{|u_x^{mlc}|^2}{gH} \operatorname{sign}(u_x^{mlc}), & \text{if } x^{mlc} - x^{wall} < 0 \\ 0, & \text{otherwise} \end{cases} \quad (2)$$

with the iceberg weight $F_W = \rho_l g \epsilon H^2 W_l$ and $\epsilon = L/H$. Subscripts “glac”, “iceb” and “mlc” stand for “glacier”, “iceberg” and “most left corner” respectively. The iceberg most left corner is in contact with the glacier. It corresponds to the top left corner in the bottom-out case (see Fig. 2) and to the bottom left corner in the top-out case (see Fig. 12). In expression (2), $a$ and $b$ are dimensionless parameters tuned to limit the iceberg penetration inside the wall (see Section 4), $x^{mlc} - x^{wall}$ is the distance between the iceberg most left corner and the wall ($x^{wall} = 0$) and $u_x^{mlc}$ is the iceberg most-left-corner horizontal velocity. This model has been described and was initially calibrated by Bonnet (2021) in a dimensional form. The formulation presented in Eq. (2) is adapted to any length scale $H$.

As for the friction law between solids, we add a simple Coulomb-type friction force to the wall $F_C$ (along tangential vector $\vec{z}$) given by Eq. (3).

$$F_C = \pm \vec{F}_{glac \to iceb} \cdot \vec{z} = \pm \mu_C \vec{F}_{glac \to iceb} \cdot \vec{x}, \quad (3)$$

with the + sign in bottom-out configuration and - sign in top-out. By varying $\mu_C$, it is possible to compute a contact between smooth or rough surfaces which can be different in the laboratory or in the field.

## 3. Experimental database

Three sets of iceberg capsize experiments were conducted in a water tank made of plexiglass with plastic icebergs and published in Amundson et al. (2012b); Burton et al. (2012) and Murray et al. (2015c), see
<figure><img src="image_3.png" /><figcaption>Fig. 3. Sketches of the three laboratory experiments with (a) the Open Ocean (OO) case (Burton et al., 2012) and (b, c) the contact cases: (b) Floating Glacier (FG) (Amundson et al., 2012b) and (c) Grounded Glacier (GG) (Murray et al., 2015c). Hatched areas correspond to the glacier wall and the ground.</figcaption></figure>

Fig. 3. Each of them corresponds to a different configuration but the experimental process was similar. The authors immersed a plastic iceberg of density $\rho_i = 920 \text{ kg/m}^{-3}$ into a tank filled with fresh water of density $\rho_w = 997 \text{ kg/m}^{-3}$. According to Burton et al. (2012), tank outer dimensions are length × width × height = 244 cm × 30 cm × 30 cm. The lateral walls of the tank have a thickness of 1.3 cm.

Initially held by hand in a hydrostatic, near vertical position ($\theta_0 \approx 0$), the iceberg is released and capsizes. This feature makes repeatability difficult and initial conditions are uncertain. Depending on the configuration, the iceberg height $H$ is equal to 10.3 or 20.3 cm, the aspect ratio $\epsilon = L/H$ varies between 0.22 and 0.63. The transverse iceberg length is $W_i = 26.7$ cm. More specifically:

* In (Burton et al., 2012) the iceberg is allowed to capsize in “Open Ocean” (OO) and the water depth $d_w$ varies from 11.4 cm to 24.3 cm which is about 1.1H to 2.4H (see Fig. 3 (a)). Therefore $\Delta z_{bot}$ varies between 0.18H and 1.44H.

* In the case described by Amundson et al. (2012b), there is a fixed “Floating Glacier” (FG) made of plastic with a height of 10.3 cm (equal to the height of the iceberg). The water depth $d_w$ was approximately equal to 13 cm $\approx 1.3H$ (see Fig. 3 (b)). Thus $\Delta z_{bot} \approx 0.34H$. The plastic glacier had a length of 15 cm.

* Finally, in Murray et al. (2015c), the 20.3-cm-high iceberg is placed flush against a wall that extends to the tank bottom. A narrow channel (flow section of about $2 \times 260 \text{ mm}^2$) allows water to flow under the plastic glacier but is not expected to have a significant effect on the capsize dynamics. Therefore, we will ignore this channel and refer to this case as “Grounded Glacier” (GG) in the following sections. The water depth is $d_w = 20$ cm which leaves about $\Delta z_{bot} = 2$ cm = 0.1H (see Figs. 1 and 3(c)). In this last configuration, both bottom-out and top-out iceberg capsizes were performed.

These laboratory experiments were scaled based on Froude number ($Fr$) conservation which is more practical than Reynolds number ($Re$) scaling for large-scale hydrodynamic applications. Usually,

$$Re = \frac{\rho_w u_{ref} L_{ref}}{\mu}, \quad Fr = \frac{u_{ref}}{\sqrt{g L_{ref}}}$$

with $u_{ref}$, a reference flow velocity and $L_{ref}$, a reference length. Here $u_{ref}$ is unknown a priori so, in Burton et al. (2012), the Froude number is computed as

$$Fr = \sqrt{8\pi \left(1 - \frac{\rho_i}{\rho_w}\right) (1 - \epsilon)} \approx 1$$

with a difference of about 10% between the field and the laboratory due to different densities of the water (sea water versus fresh water) and iceberg (ice versus plastic). The Reynolds numbers

$$Re = Fr \frac{\rho_w}{\mu} \sqrt{\frac{gH^3}{8\pi}}$$

are above $Re_{field} = 10^7$ for field conditions and about $Re_{lab} = 10^4$ for lab conditions (Burton et al., 2012). Therefore, turbulence might not develop as much in the experiment as it does at the field scale. Still, some eddies appear in the water as shown in Burton et al. (2012).

In (Burton et al., 2012) paper, a gap $\Delta W$ of 3.5 mm was maintained between the iceberg and each lateral wall (see Fig. 1). This gap is meant to prevent contact of the iceberg with the tank side during motion. However, it is possible that some rounding errors may have been introduced in Burton et al. (2012) and having $\Delta W = 6$ mm on each side is more realistic (J.C. Burton, personal communication). We checked this last statement by computing relative distances from the video presentation made by Amundson et al. (2012a).

Among the various quantities measured, translations and rotations were captured by tracking the black spots attached to the plastic (see Fig. 3). This will serve as the basis for validating our numerical model.

## 4. Numerical accuracy

In this section section we investigate the numerical convergence as numerical parameters are varied. This includes the mesh size, time step and coefficients involved in Eqs. (1) and (2).

We define the dimensionless time $t' = t/\tau$, where $\tau = \sqrt{H/g}$ is the time scale of gravity (with $H$ the iceberg height and $g$ the gravity). The variable $t'$ will be employed until Section 6.

### 4.1. Mesh refinement effect

The objective of this study is to present the results of a mesh convergence analysis of the 3D model in grounded glacier configuration. Our results are compared to experimental results from Murray et al. (2015c) for aspect ratios $\epsilon$ of 0.22, 0.28 and 0.43. Three levels of refinement were selected and designated M3, M4 and M5 from the coarser to the finer following the nested mesh refinement characteristics described in Appendix B. A fine mesh contains more and smaller cells than a coarse one. The number of cells used during a computation is reported in Table 2.

Since the time step depends on the cell size through the Courant number, a convergence in time is concomitant with a convergence in space.

<table><caption>**Table 2**
Number of million cells in a 3D simulation for each mesh level and each iceberg aspect ratio $\epsilon$.</caption><thead><tr><th rowspan="2">Level of<br>refinement</th><th colspan="3">Number of millions of cells</th></tr><tr><th>ε = 0.22</th><th>ε = 0.28</th><th>ε = 0.43</th></tr></thead><tbody><tr><td>M3</td><td>0.55</td><td>0.56</td><td>0.61</td></tr><tr><td>M4</td><td>1.57</td><td>1.70</td><td>1.97</td></tr><tr><td>M5</td><td>3.87</td><td>3.72</td><td>5.1</td></tr></tbody></table>
<table><caption>**Table 3**</caption><caption>Overview of cases run in the Section 4.2. The wall penetration coefficients (a, b) come from Eq. (2). Simulations are performed with an iceberg of aspect ratio ε = 0.22. Bold characters correspond to default values used in Section 6.</caption><thead><tr><th>Wall pen. coeff. (a,b)</th><th>(14.10<sup>3</sup>, 7.10<sup>1</sup>)</th><th>(14.10<sup>4</sup>, 7.10<sup>2</sup>)</th><th>(14.10<sup>5</sup>, 7.10<sup>3</sup>)</th></tr></thead></table>

Results are displayed on Fig. 4a and b. Mesh levels are plotted from the coarser to the finer in dotted, dashed and solid lines respectively. Blue, orange and green correspond to $\epsilon = 0.22$, 0.28 and 0.43 respectively. Circles show experimental results. One can see on Fig. 4a that finer-mesh solutions converge as expected. Indeed, M3 (dotted) curves exhibit a distinctive shape from other simulations and experiments while M4 and M5 are almost overlapping. Fig. 4b is a zoom around $\theta = 90^{\circ}$. The slight differences in arrival time are more clearly discernible. Still, amplitude and frequency of oscillations are very well resolved and M4

solutions seem to be sufficiently converged as it is almost superimposed to the M5 curve. In order to avoid excessive CPU cost, a M4 mesh level is chosen for all other simulations.

## 4.2. Wall penetration

Tuning the couple of parameters (a, b) from Eq. (2) enables reducing the wall penetration. The selected values of (a, b) are summarized in Table 3. We can track the iceberg top-left-corner position, which is in contact with the wall in the bottom-out configuration. Ideally, this corner must not enter into the wall and thus remain on the plane $x = 0$. Fig. 5 shows, in blue, the top left corner position along $x$ denoted $T_x^{tlc}$, expressed as a fraction of the typical wall cell size equal to 1.2 mm for mesh level M4. Three couples of parameters (a, b) are considered and the blue curves follow a similar pattern: for $t' \in [0, 2]$, $T_x^{tlc}$ decreases abruptly to negative values, meaning that the iceberg corner is penetrating the wall. Since we are using a spring-damper contact

<figure><img src="image_4.png" /><figcaption>**Fig. 4.** Iceberg rotation angle as a function of the dimensionless time $t' = t/\tau = t/\sqrt{H/g}$ (for mesh levels Mi, $i = 3, 4, 5$ (see Table 2)) compared to experiments from Murray et al. (2015c) for three aspect ratios $\epsilon = 0.22$, 0.28 and 0.43. Upper row: complete time series, lower row: zoom around $\theta(t) = 90^{\circ}$ (end of the capsize).</figcaption></figure>
force, the wall/iceberg model behaves as a mass-spring-damper model and exhibits a transient regime for $t' \in [0, 2]$ with damped oscillations which frequency depends on a (stiffness) and b (damping). It is followed by a prolonged increase of the wall penetration from $t' \approx 2$ to $t' \approx 65$. $T_x^{tlc}$ then increases quickly around $t' \approx 65$ when the iceberg starts drifting away from the wall. As the values of a and b increase, there is a clear convergence toward a no-penetration condition. For the couple $(a, b) = (14.10^5, 7.10^3)$, the iceberg corner penetrates a distance of only 3% of a cell size, i.e. $\approx 36 \mu m$.

On Fig. 5, the contact force value divided by the iceberg weight $F_W$ is represented in red. A transient regime is also visible for $t' \in [0, 2]$. One noteworthy aspect of this model is that changing $(a, b)$ does not affect the contact force exerted by the glacier on the iceberg which reaches 3% of the iceberg weight at $t' \approx 60$. The variation of $(a, b)$ only affects the top left corner penetration into the wall. One can check that according to Eq. (2), increasing a and b by 100 leads to a wall penetration divided by $\sqrt{100} = 10$. Notes that $F_{glac→iceb}$ reaches about 3% of the iceberg weight.

As the three tested couples exhibit reasonable wall penetration (less than 30% of the typical most refined cell size, i.e. less than 0.36 mm) and very similar contact forces, any of these could be retained for subsequent analysis. However, the couple $(a, b) = (14.10^5, 7.10^3)$, which gives the least wall penetration, makes the solver stiffer and produces artifacts visible on Fig. 5 at $t' = 23$. To circumvent this issue and maintain the wall penetration as small as possible, we selected $(a, b) = (14.10^4, 7.10^2)$ for all subsequent simulations.

### 4.3. Sensitivity of the volume penalization method (VPM)

We checked the efficiency of the porous medium implementation by varying the exponent $n$ in Eq. (1) (see Table 4 for the values). It is also possible to vary $k$, however, given that $\Delta t/t_0$ is of the order of $10^{-3}$, the effect of $n$ is dominant. In these simulations, the contact force model was kept unchanged $(a, b) = (14.10^4, 7.10^2)$. The first case is $n = 0$ which means water can freely flow inside the wall volume even though the contact force pushing the iceberg away is still active. Extending the numerical domain in the negative $x$ can have a significant effect in this case as it provides more space for water to flow. In these simulations, the wall domain had a thickness of 2 cm. Fig. 6 displays the iceberg absolute

<caption>**Table 4**</caption>

Overview of cases run in the Section 4.3. The VPM exponent $n$ comes from Eq. (1). Simulations are performed with an iceberg of aspect ratio $\epsilon = 0.22$. Bold characters correspond to default values used in Section 6.

<table><thead><tr><th>VPM exponent n</th><th>0</th><th>0.5</th><th>1</th><th>1.2</th><th>1.5</th></tr></thead></table>

rotation angle in simulations of 3D grounded glacier with $\theta_0 = 2^\circ$ and $\epsilon = 0.22$ for the different values of $n$.

As can be observed, there is a discernible contrast in the anticipated trajectory between the cases where $n < 1$ and those where $n \ge 1$. For $n = 0$, the wall is fully porous and for $n = 0.5$, partially porous. Hence, water can flow in and out of the wall volume. Curves $\theta(n = 0)$ and $\theta(n = 0.5)$ reach 90° before cases with $n \ge 1$. All the curve slopes are very similar which means that the porosity only plays a role in the triggering of the capsize. A high porosity wall leads to a faster capsize because the water initially located inside such wall can participate in pushing the iceberg.

The VPM source term in Eq. (1) becomes significant when $n \ge 1$. In this case, the fluid located in the wall is “locked” ($||\vec{u}|| \approx 0$ m/s). This is verified in the simulations as $n = 1, 1.2, 1.5$ curves are almost identical. The small differences between these three cases might come from the first mesh layer at $x \lesssim 0$ which exhibits a slight leakage due to the porous formulation. Indeed, within this mesh layer, velocities are not exactly zero ($||\vec{u}|| \approx 10^{-2}$ m/s). They quickly converge towards smaller values as $x$ decreases, for example, in the second layer of cells, we observe ($||\vec{u}|| \approx 10^{-4}$ m/s). Further increase of $n$ does not make a strong difference in the dynamics even though it slightly reduces the “leaks”. The drawback of high $n$ values is the more difficult solver convergence which increases the CPU time, see Section 2.2.

In the following, we present results obtained for $n = 1.2$.

### 5. Geometrical & physical parametric study

In this section, we present the effect of geometrical and physical parameters on the capsize dynamics. Three different sets of simulations are performed. The selected parameters are the initial tilt angle $\theta_0$, the lateral space $\Delta W$ between the plastic iceberg and the tank lateral wall

<figure><img src="image_5.png" /><figcaption>**Fig. 5.** Effect of the parameters of the contact force model, Eq. (2). (Blue left axis) the dimensionless horizontal position of the iceberg top left corner $T_x^{tlc}$ versus dimensionless time $t' = t/\tau$ for three couples of contact force parameters $(a, b)$. Negative values of $T_x^{tlc}$ mean that the iceberg corner is inside the wall. (Red right axis) Dimensionless contact forces are indicated. Grounded glacier configuration with $\epsilon = 0.22$ and M4 mesh is selected. In this case, the typical wall-cell size along $x$ is 1.2 mm.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 6.** Effect of the exponent *n* (see Eq. (1)) on the iceberg rotation angle *θ* as a function of the dimensionless time *t'* = *t*/τ. Increasing *n* means the wall is less porous.</figcaption></figure>

(see Fig. 1) and the Coulomb friction coefficient $\mu_C$. Studied cases are reported in Table 5.

## 5.1. Influence of iceberg initial inclination

In the comparison between experiments and simulation, the initial iceberg inclination is the main source of uncertainty. We investigate the impact of varying $\theta_0$ on the overall dynamics. We selected $\theta_0 = [0.5, 1, 2, 4, 8]°$ and a grounded glacier case with aspect ratio 0.22. On Fig. 7, we display the iceberg rotation as a function of the time. It was decided to take $t' = 0$ when the angle $\theta$ reaches 20°. On the figure, it is obvious that for smaller $\theta_0$, the capsize takes more time to complete. The limit would be, if $\theta_0 = 0°$ and the fluid is at rest everywhere, corresponding to an (instable) equilibrium and no capsize.

Interestingly, for $t' > 0$ s, all the curves are almost superimposed. Only small differences of about 0.1° appear in the oscillations at $\theta \approx 90°$. It means that there is no significant influence of the initial tilt angle on the overall dynamics. Indeed, even for $\theta_0 = 8°$, the iceberg angular velocity quickly reaches the same value as the other cases. This phenomenon can be explained by the fact that the capsize process starts at a very slow rate. In the initial phase of the motion, up to an angle of 20°, the build-up of inertia is insufficient to exert a notable influence on

<caption>**Table 5**</caption>

Overview of cases run in the Sections 5.1–5.3. All simulations are performed with an iceberg of aspect ratio $\epsilon = 0.22$ and each line corresponds to a different study. Bold characters correspond to default values used in Section 6. $\theta_0$ is the initial tilt imposed on the iceberg. $\Delta W$ is the lateral gap between iceberg and tank. Coulomb friction with $\mu_C = 0.2$ was added to the case (*).

<table><thead><tr><th></th><th>0.5</th><th>1</th><th>2</th><th>4</th><th>8</th></tr></thead><tbody><tr><td>Initial angle θ<sub>0</sub></td><td>0.5</td><td>1</td><td>2</td><td>4</td><td>8</td></tr><tr><td>ΔW (mm)</td><td>0 (2D)</td><td>3.5</td><td>3.5*</td><td>6</td><td>13</td></tr><tr><td>Friction coefficient μ<sub>C</sub></td><td>0</td><td>0.05</td><td>0.2</td><td>1</td><td>1.5</td></tr></tbody></table>

<figure><img src="image_7.png" /><figcaption>**Fig. 7.** Iceberg rotation during a capsize for different initial angles as a function of the dimensionless time $t' = t/\tau$. The Grounded Glacier configuration with $\epsilon = 0.22$ is selected. Curves are shifted in time to collapse at $\theta(t' = 0 s) = 20°$. The inset shows a zoom on the curves around $\theta \approx 90°$.</figcaption></figure>

the subsequent phase, beyond 20°. This provides an explanation as to why it is possible to recover experimental results without knowing the exact initial iceberg inclination.

In terms of numerical computation, despite the adaptive and relatively large time step at the beginning of simulations, it is thus possible to save CPU time by taking a non-zero initial angle. Here, the CPU time needed for $\theta_0 = 0.5°$ was of 31 h and 22 h for $\theta_0 = 2°$. In the following we chose $\theta_0 = 2°$.

## 5.2. Importance of 3D effects

This section presents a discussion on the role of the lateral flow confinement on the capsize dynamics in a grounded glacier configuration (Murray et al., 2015c). During the capsize, water is expected to flow along the sides of the iceberg and under it. Fig. 1 shows, with purple and orange arrows, the water flow paths in a 2D or 3D model, respectively. If one of these paths is narrow, then we can say that the flow is confined. These flow paths are conditioned by $\Delta W$ and $\Delta z_{bot}$. As the values of these parameters increase, there is more space for the water to move and the capsize goes faster. Here, following the work of Murray et al. (2015c), $\Delta z_{bot}$ is fixed and imposed by the water depth $d_w = 20$ cm and iceberg height $H = 20.3$ cm. Thus, in this study, the flow confinement is varied using only $\Delta W$.

A 2D case is equivalent to having $\Delta W \to 0$. It is not possible to compute the 3D case with $\Delta W = 0$ due to numerical limitations. $\Delta W = 3.5$ mm corresponds to the setup described in Burton et al. (2012) and used in Murray et al. (2015c) experiments as well. $\Delta W = 6$ mm is the gap used in all other calculations we have conducted and which may be closer to the actual experimental configuration (J.C. Burton, personal communication). Lastly, $\Delta W = 16$ mm was selected arbitrarily as a larger value than the previous ones. The Coulomb-type friction with coefficient $\mu_C = 0.2$ (for plastic-plastic contact) was also added to a case with $\Delta W = 3.5$ mm to assess the relative importance of the iceberg-glacier friction and lateral flow section.
<figure><img src="image_8.png" /><figcaption>**Fig. 8.** Iceberg rotation angle during bottom-out capsizes as a function of the dimensionless time $t' = t/\tau$ and for different gaps between the iceberg and the tank sides $\Delta W$. Blue dots: Experiments from Murray et al. (2015c), Blue dash: 2D simulation, Pink: 3D simulation with $\Delta W = 3.5$ mm, Green: same with Coulomb friction $\mu_C = 0.2$ (plastic on plastic), Blue solid: 3D simulation with $\Delta W = 6$ mm, Brown: same with $\Delta W = 16$ mm. The curves are shifted such that $\theta(t' = 0) = 80^\circ$ to highlight dynamic differences.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Fig. 9.** Iceberg rotation angle as a function of the dimensionless time $t' = t/\tau$ for different Coulomb friction coefficients $\mu_C$ (see Eq. (3)). These 3D simulations were run with $\epsilon = 0.22$.</figcaption></figure>

On Fig. 8, iceberg rotation angle $\theta$ is represented as a function of time. For this case, curves are shifted horizontally such that $\theta_0(t' = 0) = 80^\circ$. This representation highlights the differences in the iceberg kinematics at the beginning of the capsize. Here, aspect ratio is set to $\epsilon = 0.22$. This figure illustrates the range of behaviors that can be

<caption>**Table 6**</caption>

Overview of cases run in the Section 6. Configurations are Open Ocean (OO), Floating Glacier (FG), Grounded Glacier (GG) Bottom-out (BO) and Top-Out (TO). $\epsilon$ is the iceberg aspect ratio. Data sources are experiments (Exp) and simulations (Sim) in 2D or 3D.

<table><thead><tr><td>Name</td><td>Reference</td><td>ε</td><td>Data source</td></tr></thead><tbody><tr><td>OO</td><td>Burton et al. (2012)</td><td>0.5</td><td>Exp, Sim 2D, 3D</td></tr><tr><td>FG BO</td><td>Amundson et al. (2012b)</td><td>0.5</td><td>Exp, Sim 2D, 3D</td></tr><tr><td>GG BO</td><td>Murray et al. (2015c)</td><td>0.22, 0.28, 0.43</td><td>Exp, Sim 2D, 3D</td></tr><tr><td>GG TO</td><td>Murray et al. (2015c)</td><td>0.22, 0.28, 0.43</td><td>Exp, Sim 3D</td></tr></tbody></table>

obtained starting from a 2D case which displays the slowest capsize ($t' > 100$).

The results of this study demonstrate that $\Delta W$, i.e. the flow confinement, has a significant influence on the iceberg kinematics while the friction with $\mu_C = 0.2$ has almost no effect. Our numerical results perfectly match the experimental data for $\Delta W = 6$ mm and therefore we will use this gap value in the rest of this study.

### 5.3. Coulomb friction

As presented in Section 2.2, the sole parameter that can be adjusted to alter the wall/iceberg contact is the coefficient $\mu_C$, which is involved in the Coulomb friction law. To assess the impact of this parameter, a series of simulations varying $\mu_C$ were run. The values relevant to this study are $\mu_{C,ice-ice} = 0.05$ (order of magnitude given by Oksanen and Keinonen (1982); Mills (2008); Makkonen and Tikanmäki (2014) and $\mu_{C,plastic-plastic} = 0.2$ (Wróbel and Malgorzata, 2008). Fig. 9 shows the iceberg rotation angle as a function of time for different $\mu_C$. The shifts between curves $\mu_C = 0$ to $\mu_C = 0.2$ is very small (see Figs. 8 and 9). Curves with $\mu_C = 1$ and $\mu_C = 1.5$ clearly exhibit, as expected, slower rotations. These cases with very large $\mu_C$ values could be relevant in future studies to account for a rough contact between the iceberg and the glacier. In the following, we choose $\mu_C = 0$ as default value.

## 6. Comparison to laboratory experiments

In this section, we start by comparing simulations in 2D and 3D to experiments in the three configurations: Open Ocean (OO), Floating Glacier (FG) and Grounded Glacier (GG). It should be noted that for the rotations, in the following, only the absolute value of the rotation angle $\theta(t)$ is indicated. With respect to the frame displayed on Fig. 1, a bottom-out case leads to a negative $\theta$ and a top-out case to a positive one. Table 6 provides a summary of the investigated cases. For Bottom-Out simulations, Open Ocean and Floating Glacier cases are run with aspect ratio $\epsilon = 0.5$ while grounded glacier cases are run with $\epsilon = 0.22$, 0.28 and 0.43.

### 6.1. Iceberg rotation for bottom-out capsizes

Simulations can give access to many important quantities such as the iceberg translations and rotation during the capsize but also the velocity $\vec{u}$ and pressure fields inside the water. Fig. 10 shows a snapshot of the water velocity field at $t = 4.4$ s or $t' = 30.6$ in a 3D GG simulation with $\epsilon = 0.4$, iceberg height $H = 0.2$ m and $\Delta W \gg W_i$. This simulation was made for illustration purpose and is not related to the laboratory experiments where $\Delta W \ll W_i$. The velocity magnitude $||\vec{u}||$ is represented in colors and the flow direction is indicated with vectors of constant norm. The maximum water velocity is located on the iceberg bottom right edge and is due to its sharp edges. Note that, for $x < 0$, inside the porous medium (VPM), the water velocity $||\vec{u}||$ is of order 1 mm/s which is very small. The free surface is not maintained at $z = 0$ inside this region. This is due to the small leaks described in Section 4.3. The full time series can be found in the supplementary material.

Fig. 11 shows the iceberg rotation time series in the three configurations. Experimental, 2D and 3D results are represented with dots,
<figure><img src="image_10.png" /><figcaption>**Fig. 10.** Snapshot of the velocity field (magnitude in colors and direction in vectors) in a 3D GG simulation with $\epsilon = 0.4$ and $\Delta W \gg W_t$. The corresponding animated file (.mp4) is available in the supplementary material.</figcaption></figure>

dashed lines and solid lines respectively. Each color corresponds to a different aspect ratio. The first experimental data corresponds to $t' = 0$ and we shifted the time of all simulation curves such that they reach $\theta = 20^\circ$ at the same time as in the experiments. As shown in Section 5.1, this representation allows to ignore uncertainties on the iceberg initial tilt angle. The differences in capsize duration and in rotational velocity between 2D and 3D GG configurations are also highlighted by this representation.

This figure shows very good agreement between experimental and numerical results in open ocean and floating glacier cases. Fortunately, initial conditions only play a role in the triggering of the capsize and not in its full dynamics, as explained in Section 5.1. It can also be seen that, for the same aspect ratio, the presence of a glacier front slows down the process. In Section 5.2, we have shown that the presence of the solid glacier, compared to the Open Ocean case, makes the water motion more difficult by constraining it, thus the capsize takes more time. The OO iceberg capsizes faster than the GG ones for a similar reason: the gap underneath the iceberg ($\Delta z_{bot}$) is larger in the OO case (see Section 3) thus the water takes more time to move around the iceberg and the capsize is slowed down. In Appendix D we show that the capsize occurs as fast in a floating or grounded glacier case as long as the water depth is the same. It is difficult to compare the FG with the others because its initial iceberg tilt angle is $\theta_0 \approx 7.7^\circ$ (OO and GG cases start at $\theta_0 = 2^\circ$) Note that after reaching $\theta = 90^\circ$, the 20-cm-high plastic iceberg oscillates with a pitch period of $\approx 0.56$ s. Simulations give a very close period of $\approx 0.59$ s.

The grounded glacier model in 2D does not match very well with experimental data because the rotation is too slow. Therefore, we used these cases to investigate 3D water flow effects. In 2D, the rotation is about 20% longer to complete than in 3D (Fig. 11). At $20^\circ$, one can also notice the slope difference between 2D and 3D cases, indicating that less confined icebergs can capsize more easily as in 3D calculations. 3D simulations are able to very accurately reproduce experimental results as associated curves are almost superimposed.

After reaching $90^\circ$, the amplitude of oscillations of $\theta$ (pitch oscillations) increases with $\epsilon$ and this feature is also very well reproduced by the model. This might be due to the fact that icebergs with larger aspect ratios may produce larger waves but this falls out of the scope of this study.

## 6.2. Iceberg rotation for top-out capsizes

We also investigated top-out iceberg capsizes against a glacier. There, the volume between the glacier and the iceberg, below the free surface ($z = 0$) is filled with water (12). This situation is interesting because it highlights the limitations of a 2D model: as the iceberg bottom is in contact with the glacier wall, water cannot flow and only a depression appears in the water trapped above the contact point, maintaining the iceberg immobile. This situation is showed on Fig. 12 which is a snapshot of the hydrodynamic pressure field ($P_{dyn}$) in a GG case with $\epsilon = 0.22$ and $H = 0.203$ m. This snapshot has been taken at $t' = 139$ and $\theta_0 = 5^\circ$.

Thus for Top-Out events, a 3D modeling is mandatory. With the 3D model, the water flow uniquely takes the lateral path and the capsize is achieved. Fig. 13 shows good agreement between experiments (Mur-ray et al., 2015c) and 3D simulations. It can be seen that TO capsizes is faster than BO events. Plotting the bottom-out case with $\epsilon = 0.28$ on Fig. 13 highlights the sharper velocity increase in top-out cases. The oscillations of $\theta$ values after the complete capsize also differs significantly. The largest $\theta$ value reaches $97^\circ$ for the TO case versus $91^\circ$ for the BO case. It is worth noting that this feature is remarkably well reproduced by the CFD model.

## 6.3. Iceberg translations

We define the translations $T_x$ (along the x-axis) and $T_z$ (along the z) axis of the iceberg center-of-gravity during the capsize such as $T_x = x_{CoG}(t') - x_{CoG}(t' = 0)$ and $T_z = z_{CoG}(t') - z_{CoG}(t' = 0)$, where subscript CoG stands for Center of Gravity. Removing the initial positions of the center-of-gravity is convenient because $x_{CoG}(t' = 0)$ and $z_{CoG}(t' = 0)$ depend on the iceberg aspect ratio. Dimensionless translations are indicated $T_x/H$ and $T_z/H$.

On Fig. 14, one can observe the iceberg translation along the x-axis in the three configurations presented before.

* First of all, it can be observed that results of 2D Grounded Glacier models depart significantly from the 3D ones. We found earlier that 3D effects were also significant for the iceberg rotation in GG. 3D effects are limited for the Open Ocean and Floating Glacier cases

* Fig. 14 (right) represents the top-out x-translations. It shows a very good consistency between simulations and experiments, especially at the beginning of the capsize. The iceberg velocity is increasing slowly until an abrupt acceleration at $t' = 45$ for $\epsilon = 0.43$ and $t' = 65$ for $\epsilon = 0.22$. $\epsilon = 0.28$, just before $\theta$ has reached $90^\circ$. The abrupt acceleration is due to the change from a mainly upward translation to an horizontal translation of the center of gravity. Indeed the time slope change in $T_x$ corresponds to $\theta \approx 35^\circ$, see Fig. 13. The simulations slightly overestimate the drift after a complete capsize. It can be noted that top-out capsizes result in drifts that are two to three times larger than bottom-out capsizes at the same time.

* Translations $T_x/H$ in 3D bottom-out cases also very well reproduce experimental data for the three cases (OO, FG and GG).

* In the Open Ocean case, the iceberg pushes water on the positive x-side. Thus momentum conservation explains why the iceberg goes toward negative x.

* In the BO cases where the iceberg is in contact with the wall, i.e. FG and GG, similar patterns of $T_x$ are found with a plateau on the curves when the rotation reaches $90^\circ$. For example, the dark blue curve in Fig. 14 (left) has its plateau at $t' \approx 80$. It can be observed that these plateaus are also captures by the models. During the rotation, ($\theta < 90^\circ$), the iceberg corner remains in contact with the wall. When $\theta = 90^\circ$, the rotation is completed and the translation stops. This is why we observe a plateau on the x-translation curves. At the same time, the iceberg top left side dives into the water, generating an overpressure in the fluid. This overpressure would be the
<figure><img src="image_11.png" /><figcaption>**Fig. 11.** Absolute value of the iceberg rotation angle as a function of the dimensionless time $t' = t/\tau$ during the bottom-out capsize. For Grounded Glaciers (purple, dark blue, light blue), Floating Glacier (orange), Open Ocean (red). Aspect ratios $\epsilon$ are indicated. Circles, dashed and solid lines stand for experimental (Exp), 2D and 3D simulation results respectively. Simulation curves are shifted in time such that they cross the experimental one when $\theta = 20^\circ$. The insets show zooms on the curves around $\theta \approx 90^\circ$.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 12.** Hydrodynamic pressure field ($P_{dyn}$) in a 2D top-out GG case at $t' = 139$ or $t = 30$ s. The depression in the water between the glacier wall ($x = 0$) and the iceberg is strong enough to maintain it in this equilibrium position. In this simulation, $\epsilon = 0.22$ and, for visualization purposes, $\theta_0 = 5^\circ$.</figcaption></figure>

driver of the iceberg post-capsize drift (increase in $T_x$) observed for $\theta \ge 90^\circ$. This mechanism can also be seen on videos from the supplementary material of Murray et al. (2015b). Plotting translation for different aspect ratio in purple ($\epsilon = 0.28$) and green ($\epsilon = 0.43$) shows that smaller aspect ratios result in larger drifts and larger transport velocities.

The cause of this difference in drift is unclear. We do not know how the force pushing the iceberg away from the glacier front is linked to the other parameters and this is left for future work. It

<figure><img src="image_13.png" /><figcaption>**Fig. 13.** Iceberg absolute rotation angle as a function of the dimensionless time $t' = t/\tau$ during the top-out capsize in Grounded Glaciers configuration (purple, dark blue, light blue). Aspect ratios are indicated. Circles and solid lines stand for experimental (Exp) and 3D simulation results respectively. (Green): Computed response for the Bottom-Out $\epsilon = 0.28$ case is added for comparison. The inset shows a zoom on the curves around $\theta \approx 90^\circ$.</figcaption></figure>

is unlikely that drag forces are responsible for the drift evolution because drag forces are proportional to the area of the iceberg facing the water flow (here $A = L \times W_i$). A large flow-facing area would result in a large drag force and the iceberg would drift less than
<figure><img src="image_14.png" /><figcaption>**Fig. 14.** Dimensionless horizontal translation of the iceberg center-of-gravity during the bottom-out (left) and top-out (right) capsizes as a function of the dimensionless time $t' = t/\tau$. The grounded glacier cases are shown in purple, dark blue and light blue. The floating glacier and open ocean cases are shown in orange and red respectively. Aspect ratios are indicated. Circles and solid lines stand for experimental (Exp) and 3D simulation results respectively.</figcaption></figure>

with a small flow-facing area. The FG iceberg has $A_{FG} = 133.5 \text{ cm}^2$ and the GG case with $\epsilon = 0.43$ has a flow-facing area $A_{GG}(\epsilon = 0.43) = 233.1 \text{ cm}^2$. Thus the iceberg from the GG case at $\epsilon = 0.43$ has a larger flow-facing area but also drifts more. This is not consistent with the expectation, suggesting there are other mechanisms responsible for the larger drifts in the grounded glacier case.

* Unlike simulations, Amundson et al. (2012b) experiments on floating glacier do not show this behavior. We suspect that, perhaps the plastic iceberg hits one of the tank sides before stopping completely. $T_x$ thus remains constant, contrarily to the simulated iceberg which has no obstacles and drifts like in the GG cases.

* It seems that there is less energy dissipation in the simulations since we slightly overestimate the drift at large $t'$ in almost all cases. This might come from the fact that the iceberg can only drift in one direction: $\vec{x}$ while there might be some transversal drift in the experiment (along $\vec{y}$) and/or z-axis rotation.

**Fig. 15** represents the z-translation of iceberg center of gravity. There again, one can observe the perfect match between 3D modeling and experimental data. Interestingly, it can be seen that after complete capsize, $T_z$ does not oscillate around a stable position, but keeps increasing slowly. For a plastic parallelepiped of volume $V_{tot}$ float-

<figure><img src="image_15.png" /><figcaption>**Fig. 15.** Dimensionless vertical translation of the iceberg center-of-gravity during the bottom-out (left) and top-out (right) capsizes (distances / H) as a function of the dimensionless time $t' = t/\tau$. For grounded glaciers (purple, dark blue, light blue), floating glacier (orange), open ocean (red). Aspect ratios $\epsilon$ are indicated. Circles and solid lines stand for experimental (Exp) and 3D simulation results respectively.</figcaption></figure>
<table><caption>**Table 7**
Vertical displacement of iceberg center of gravity $\Delta z$ divided by $H$ (dimensionless) at the end of simulations for several aspect ratios $\epsilon$. Values are computed with Eq. (4).</caption><thead><tr><th>ε</th><th>0.22</th><th>0.28</th><th>0.43</th><th>0.5</th></tr></thead><tbody><tr><td>T<sub>z</sub>(t' → ∞)/H(.)</td><td>0.327</td><td>0.302</td><td>0.239</td><td>0.211</td></tr></tbody></table>

ing on water, the submerged volume $V_s$ is function of the densities: $V_s = V_{tot} \times \rho_l/\rho_w \approx V_{tot} \times 0.92$ in lab conditions. The center of mass is then located 0.42 $H$ under the water surface. Computing the difference between center of mass vertical position at rest before and after capsize leads to

$$T_z/H = 0.42(1 - \epsilon) \quad (4)$$

These lead to values displayed on Table 7 which are consistent with both experimental and numerical results. Physically, the slowly increasing value of $T_z$ observed at the end of the capsize is due to water remaining on the top of the iceberg, adding weight to the structure and being drained slowly. It is remarkable that this feature is also well captured by the solver.

Overall, our 3D simulations are in very good agreement with experimental results in every configuration. We also find that the 2D model is robust enough only in ungrounded-glacier cases and open ocean.

## 7. Towards capsize modeling of field-scale icebergs

In order to extend our computations to a field-scale scenario, we selected a large calving event from the literature (Jouvet et al., 2017; van Dongen et al., 2019; Kneib-Walter et al., 2022) which happened at Bowdoin glacier, North-West Greenland. Time-lapse imagery showed that the iceberg calved bottom-out against the glacier (van Dongen et al., 2019). Moreover, Bowdoin is a grounded glacier so the configuration is very similar to the one of Murray et al. (2015c) experiments. It is not possible to recover exactly the iceberg geometry but assuming a full-thickness calving, the iceberg dimensions were determined using aerial vehicles/satellite imagery and the fjord bathymetry: $L \times W_l \times H = 80 \times 650 \times 200 \text{ m}^3$ (Sugiyama et al., 2015; Minowa et al., 2019). The mass of this iceberg is of about $10^7$ tons. The main difference with the cases investigated earlier, besides dimensions, is that there are no lateral walls as Bowdoin glacier width is $W_g \approx 3$ km and the calving happened at about 1.5 km away from the fjord lateral side (van Dongen et al., 2020). With this geometry, the water can flow more on the side than under the iceberg which is a case we have not reproduced yet in this study.

At this location, the water depth reaches about 200 m which means $\Delta z_{bar} \approx 20$ m. Similar undercut geometries have been used in previous studies (van Dongen et al., 2019, 2020) and are based on field observations on other glaciers with meltwater plumes (Fried et al., 2015; How et al., 2019).

For the following field-scale 3D simulations, field densities are taken as $\rho_l = 917 \text{ kg/m}^3$ and $\rho_w = 1025 \text{ kg/m}^3$ (Cuffey and Paterson, 2010). The numerical mesh is similar to those of previous sections. We refine the area around the iceberg, close to the wall and at the fjord bottom as represented on Fig. 2. For the iceberg-glacier contact, the Coulomb friction coefficient is taken as $\mu_C = 0.05$ (ice-ice). In the field, this contact is a complex process including evolving surface geometries, roughness and lubrication by water. However, we have shown in Section 5.3 that friction does not significantly affect the capsize, even for very high coefficients ($\mu_C \ge 1$). The exact initial water velocity field is unknown but we expect it to be negligible compared to the flow generated by the capsize. This hypothesis is validated later as the water flow velocities in a fjord “at rest” is two orders of magnitude smaller

than the water flow velocities due to a capsize. Thus the initial water velocity field is set to 0 everywhere. The iceberg initial angle with respect to the vertical axis is also unknown but, we have shown that it does not influence the capsize dynamics (Section 5.1) thus we set $\theta_0 = 2^\circ$ and the iceberg velocities (in translations and rotations) are all equal to 0. Boundary conditions are the same as in the laboratory cases (see Section 2):

* wall function on the iceberg surface

* symmetry condition on the plane $y = 0$

* slip condition on the fjord bottom and its lateral side

* imposed hydrostatic pressure on the top of the domain

The glacier front is treated as low porosity medium with the Volume Penalization Method, see Section 2.

### 7.1. Water flow

Fig. 16 shows two snapshots extracted from the simulation at times 133 s and 190 s. The velocity magnitude is displayed on the ground which has a slip boundary condition (see Section 2) and on the symmetry plane $y = 0$. The glacier front consists of the surface $x = 0$ as on Fig. 2. For an iceberg of 200 m, one can notice water velocity magnitude of $\approx 10 \text{ m/s}$ in the close-iceberg vicinity as well as turbulent structures in red/dark red indicated by the black arrows on both snapshots. The maximum velocities are generated in a vortex under the iceberg on Fig. 16a and b (behind the transparent iceberg). Another vortex is visible under the green dot ($x = 400 \text{ m}, y = 0 \text{ m}$) due to high flow velocity on Fig. 16b and is advected towards positive $x$. As we move further away from the iceberg, the water velocity decreases down to $10 \text{ cm/s}$ at $\approx 800 \text{ m}$ from the origin, which is the typical velocity of currents in fjords (Johnson et al., 2011; Mortensen et al., 2014; Sutherland et al., 2014). We can say that water is disturbed up to more than twice the iceberg width in a roughly-circular area.

### 7.2. Tsunami genesis

Bottom-out/top-out capsize of large icebergs generate tsunamis which are recorded by two gauges in our simulations. These are located at ($x = 400 \text{ m}, y = 0 \text{ m}$) for gauge 1 and ($x = 0 \text{ m}, y = 400 \text{ m}$) for gauge 2 (purple and green dots on Fig. 16). The gauge recordings correspond to the local water surface elevation. The bottom-out simulation case is the same as aforementioned in this section and the top-out case has the same characteristics except the initial tilt angle $\theta_{0,bottom-out} = -\theta_{0,top-out}$. The water surface elevation $\eta$ is displayed on Fig. 17. This figure shows a wave of amplitude 2.5 m in bottom-out and 3 m in the top-out case recorded by the gauge located on the x-axis.

Note that the spike at $t = 220 \text{ s}$ corresponds to the iceberg passing under gauge 1. The emerged iceberg height is then of about 8.5 m, adding a thin layer ($\approx 1 \text{ m}$) of water trapped above the iceberg. Once this water layer falls back into the ocean, the signal vanishes because the gauge is blanked by the iceberg.

Smaller waves (1 to 2 m) are recorded close to the glacier wall (gauge 2). The wave produced is of quite small amplitude in this case compared to what was found in other studies (Chen et al., 2020; Wolper et al., 2021) but is dependent on the calving type. For example, Chen et al. (2020) found a 40-m high water elevation for an iceberg falling/sliding into the water. In the bottom-out case, at location (400 m, 0 m) (in light green on Fig. 17), two main frequencies can be observed. For $t \in [150, 200 \text{ s}$, the signal of frequency $f_v = 0.1 \text{ Hz}$ could correspond to the vortex passing under the wave gauge around $t \approx 180 \text{ s}$ and which is visible on Fig. 16b. The second main frequency, $t \in [200, 320 \text{ s}$ is the iceberg heave frequency $f_{heave} = 0.046 \text{ Hz}$ (vertical oscillations).
<figure><img src="image_16.png" /><figcaption>Fig. 16. Water velocity magnitude displayed on the ground (slip condition) and on the symmetry plane y = 0 (thus cutting the iceberg in half) and on the wall x = 0. In the half-domain simulation, the iceberg dimensions are L × Wᵢ × H = 80 × 325 × 200 m³. This shows that the water is affected by the capsize up to x = 700 m and y = ±800 m. Note that the color-scale is exponential. The glacier front is the plane x = 0. Purple and green dots indicate the positions of water wave gauge 1 (x = 400 m, y = 0 m) and wave gauge 2 (x = 0 m, y = 400 m), see Section 7.2. Black arrows show where the highest values of the velocity magnitude are located. An animated file (.mp4) is available in the supplementary material.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. 17. Water surface elevation η as a function of the time for the field-scaled capsize case. (TO) stands for Top-Out and (BO) for Bottom-Out. Gauge positions (x, y) are in meters and the colors correspond to the dots on Fig. 16.</figcaption></figure>

Our model is well-suited for estimating the generated-wave amplitude in the near field. Its propagation in fjords can then be studied with codes based on other formulations such as Saint-Venant (Popinet, 2012; Macías et al., 2016; Poulain et al., 2022a) or Boussinesq (Grilli et al., 2007; Poulain et al., 2022a; Svennevig et al., 2024) approximations.

Finally, it is worth noting that it would have been possible to run simulations at laboratory scale and rescale the results based on the Froude number which we want to preserve instead of performing field-scale simulations as above. In that case, the Reynolds number would be different but with a limited impact regarding the first stage of calving. In Appendix C, we show that a lab-scaled Reynolds number simulation ($Re_{lab} = 10^4$) with turbulence gives results close to an Euler formulation ($Re \to \infty$) for $\theta < 90^\circ$ and slightly differ in the iceberg drift velocity ($\theta > 90^\circ$). The ratio between field and laboratory length scales $H_{field}/H_{model}$ can be used to derive other quantities such as velocities u and characteristic times t.

$$Fr_{model} = \frac{u_{model}}{\sqrt{gH_{model}}} = \frac{u_{field}}{\sqrt{gH_{field}}} = Fr_{field}$$

which leads to $u_{field} = u_{model} \sqrt{\frac{H_{field}}{H_{model}}}$. Now, using the relationship between typical velocity, distance and time, we have

$$t_{field} = \frac{H_{field}}{u_{field}} = t_{model} \sqrt{\frac{H_{field}}{H_{model}}}$$

## 8. Conclusion

We have implemented a numerical model for fluid-structure interaction with solid-solid contact between a capsizing iceberg and a rigid glacier. The glacier wall is composed of a porous medium with a Volume Penalization Method and a spring-damper repulsive force to mimic the glacier/iceberg contact. The convergence of our solution is checked in analyses involving mesh refinement and wall-model parameters. The model is able to recover experimental results in various configurations and iceberg aspect ratio even though the experimental initial iceberg tilt angle was not precisely known. We have shown that this angle does not have a significant effect on the capsize dynamics.

In terms of CPU time, the 2D computation takes about 6 h on 2 CPU while the 3D computation takes about 24 h on 48 CPU. 2D simulations are numerically cheaper but they do not seem to represent the reality because no field study has reported a case of full- (or quasi-full-) glacier-width iceberg calving ($W_i \approx W_g$). We have shown that 3D water flows are important even in quasi-2D cases (see Section 5.2) and their effect on the capsize dynamics is not negligible. 3D simulations are then more accurate even though they are more expensive.

Our fluid-structure simulations are able to accurately capture details of the iceberg and flow dynamics such as:

* The influence of the iceberg aspect ratio $\epsilon$ on the timing of capsize.

* The plateau in x-translations when the iceberg angle reaches $\theta = 90^\circ$ during a bottom-out capsize and the absence of plateau in the top-out cases.

* The water remaining on the top of the horizontal iceberg after the capsize. This water is adding weight to the iceberg and is draining slowly, resulting in a slow emergence of the iceberg.

Application to a field-scale case illustrates the extent of flow disturbance following an 200-meters-high iceberg calving as well as its
intensity. Velocities over 10 m/s are estimated and the water velocity field is affected up to 800 m along the glacier front which could trigger ocean layer mixing (Grosfeld et al., 2001; Meredith et al., 2022) and biological fluxes similar to plume (Lydersen et al., 2014; Sugiyama et al., 2014; Nishizawa et al., 2019).

The current model is limited to iceberg capsizes in which a whole iceberg edge remains in contact with the glacier. One limitation is the fact that the iceberg is modeled as a parallelepiped-like rigid body. In reality, the iceberg and the glacier would deform and possibly break, thus changing their shape and dissipating energy in the process. The ice-melange is not modeled here and could play a significant role in the dynamics and energy dissipation as well.

Further work will involve more accurate geometries for the iceberg, glacier and fjord bathymetry. Additionally, the force histories computed in simulations will be useful to:

* Compute the mechanical response of the glacier with the solver ElmerIce (Gagliardini et al., 2013) to study the changes in local forces and the potential change in the glacier flow.

* Compute the iceberg volume by comparing forces inverted from the glacial earthquake signal to the force output by the simulations (Sergeant et al., 2019). Pirot et al. (2023) and Wetter et al. (2024) are currently working on machine-learning algorithms to detect and classify these glacial earthquake signals over the last 30 years (Tsai and Ekström, 2007; Veitch and Nettles, 2012; Olsen and Nettles, 2017), it will be possible to give more accurate estimation of calved-iceberg volume (Sergeant et al., 2019; Cheng et al., 2021; Köhler et al., 2022).

## CRediT authorship contribution statement

**Nicolas De Pinho Dias:** Writing – review & editing, Writing – original draft, Visualization, Validation, Software, Investigation; **Alban Leroyer:** Writing – review & editing, Validation, Supervision, Software, Resources, Methodology, Investigation, Funding acquisition, Formal analysis; **Anne Mangeney:** Writing – review & editing, Visualization, Validation, Supervision, Resources, Project administration, Methodology, Investigation, Funding acquisition; **Olivier Castelnau:** Writing – review & editing, Validation, Supervision, Methodology, Funding acquisition, Formal analysis, Conceptualization.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgments

We would like to thank J.C. Burton who greatly helped by sharing experimental details as well as Evgeniy Podolskiy and Stefan Frei for insightful discussions. We thank two anonymous reviewers for their precious comments which helped improving the quality of the manuscript. This work benefits from fundings from doctoral school STEP'UP and Agence Innovation Defense, financial support of the INSIS CNRS Institute through the project SIM-ICE, part of the grant PEPS INSIS 2023 Engineering to meet the challenge of climate change. It is granted access to the HPC resources of Nautilus through GLiCID, Ecole Centrale de Nantes, France. We would also acknowledge financial support from Horizon Europe MSCA Doctoral Network EnvSeis under grant agreement no. 101073148.

## Appendix A. Numerical wave damping

Placing elongated cells on the ocean side of our simulations provides a good water wave damping and reduces the reflections. Fig. A.18 shows the mesh used for the tank.

<caption>**Table B.8**</caption>

Characteristics of the mesh used. Three levels of refinement are selected and called M3, M4, M5. Refinement diffusion $r_d$ is a parameter linked to the number of number of cells in the most refined region. The typical cell size in each region depends on the mesh refinement level.

<table><thead><tr><th colspan="4">Mesh parameters</th></tr><tr><th>Case</th><th>r<sub>d</sub></th><th>Number of cells in the most refined layer</th><th>Typical cells size</th></tr></thead><tbody><tr><td>M3</td><td>2</td><td>6</td><td>4h/3</td></tr><tr><td>M4</td><td>3</td><td>8</td><td>h</td></tr><tr><td>M5</td><td>4</td><td>10</td><td>4h/5</td></tr></tbody></table>

## Appendix B. Nested mesh refinement

This section gives a general insight on the meshing procedure we used and how the refinement works for mesh convergence. We select three meshes called $M_i$ for $i \in [3, 4, 5]$. Each of them is initialized with a different number of cells which are multiples of $i$ in each direction. For example the mesh M4 can be initialized with $4N_x$ cells in the $x$ direction, $4N_y$ cells in the $y$ direction and $4N_z$ cells in the $z$ direction where $N_x, N_y, N_z$ are positive integers. Mesh M3 would have $3N_x, 3N_y$ and $3N_z$ cells in the respective directions. If we call $h_4 = h$ the typical size of a M4-mesh cell, then the typical size of a M3-cell is $h_3 = 4h/3$ and for the M5-mesh, $h = 4h/5$. The refinement is then performed by cutting the cells in 2 until the target cell size is reached. The most refined region contains $2(r_d + 1)$ cells where $r_d$ is a parameter. Table B.8 summarizes the mesh refinement parameters for the three selected meshes. On Fig. B.19, it can be seen that, as $i$ increases, we have more and smaller cells while the length of each region remains the same. This permits to have a space convergence in the simulations.

## Appendix C. Influence of the turbulence model

In order to quantify the influence of the turbulence on the capsize dynamics, an Euler formulation (no viscosity) was compared to the $k-\omega$ turbulence model. With both models, simulations were run in 2D with an iceberg of aspect ratio $\epsilon = 0.22$, height $H = 0.203$ m and in a Grounded Glacier case (GG). Fig. C.20 shows the dimensionless translation of the iceberg center of gravity $T_x/H$ (left) and its angle of rotation $\theta$ (right) as functions of the dimensionless time $t' = t/\sqrt{(H/g)}$. The Euler model is represented with the blue curve and the $k-\omega$ model is represented in orange. Both exhibit a very similar capsize dynamics. The main difference is noticeable after $t' \approx 87$ ($\theta$ has reached 90°). For $t' > 87$, the blue curve (Euler model) has a steeper slope than the orange curve ($k-\omega$ model). As explained in Section 6.3, for $t' > 87$, the iceberg is drifting away from the glacier. Thus, The observed difference in slope steepness indicates that the iceberg is drifting slower with the $k-\omega$ model which is expected when we introduce an energy dissipating mechanism like viscosity or turbulence.

## Appendix D. Grounded glacier with larger water depth

It is interesting to know if there is a difference in the iceberg capsize due to the presence of water under the glacier (floating glacier) or not (grounded glaceir).

A simulation is run almost like in the FG configuration (see, Section 3) with a the water depth $d_w = 1.3H$, the iceberg height $H$ and aspect ratio $\epsilon$ are $H = 10.3$ cm and $\epsilon = 0.5$ respectively. As in the FG case, the initial iceberg tilt angle is set to $\theta_0 = 7^\circ$. The only difference is that the glacier is not floating anymore but is grounded. This simulation is called case $\delta$ in this section. Only a 2D case is run because it was
<figure><img src="image_18.png" /><figcaption>**Fig. A.18.** Tank mesh used in our simulations with elongated cells on the right side. These provide a good numerical water wave damping and avoid reflections. The refined region on the left is where the capsize takes place.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Fig. B.19.** Nested mesh refinement features. M3, M4 and M5 are the names given to the meshes used in the convergence study (Section 4.1). M3 is the coarsest and M5 is the most refined mesh of the three. On this figure, each mesh contains three regions of equal length but with coarse, medium and small cells (from the left to the right). The refinement diffusion $r_d$ is a parameter used in the mesh generation such as the most refined region contains $2(r_d + 1)$ cells along the direction of interest.</figcaption></figure>

<figure><img src="image_20.png" /><figcaption>**Fig. C.20.** Dimensionless horizontal translation of the iceberg center of gravity $T_x/H$ (left) and its angle of rotation $\theta$ (right) as functions of the dimensionless time $t' = t/\sqrt{H/g}$. The Euler model is represented with the blue curve and the $k-\omega$ model is represented in orange.</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Fig. D.21.** Dimensionless horizontal translation of the iceberg center of gravity $T_x/H$ (left) and its angle of rotation $\theta$ (right) as functions of the dimensionless time $t' = t/\sqrt{H/g}$. In orange shows the experimental (dots), 2D (dashed line) and 3D (solid line) simulation results for the FG case. In green shows the results of the case $\delta$.</figcaption></figure>
shown that if $d_w$ is large enough $d_w \ge 1.3H$, 3D computations are not necessary.

Fig. D.21 shows the dimensionless horizontal translation of the ice-berg center of gravity $T_x/H$ (left) and its angle of rotation $\theta$ (right) as functions of the dimensionless time $t' = t/\sqrt(H/g)$. Orange dots show the experimental (dots), 2D (dashed line) and 3D (solid line) simulation results for the FG case. The green dashed line shows the results of the case $\delta$. It is clear that the $\theta$ curves are almost superimposed for $t' < 30$ (before the end of the capsize) meaning that the capsize dynamics is the same whether the glacier is floating or grounded.

The difference appears on the $T_x$ curves (Fig. D.21b) for $t' \ge 30$. In the case $\delta$, the iceberg post-capsize drift is smaller than in the FG case as the green dashed curve is below the orange dashed and solid curves. Since the iceberg drift seems to be linked to waves and hydrodynamic pressure fluctuations, the difference between 2D cases FG and $\delta$ might indicate a different behavior of the pressure field after the iceberg cap-size.

## Appendix E. Supplementary material

### E.1. Animation 1: Lab_scale_capsize_2dview.mp4

We provide an animation of the water velocity field in the plane $(x, y = 0, z)$. It is a side view of 3D simulation in grounded glacier configuration with iceberg aspect ratio $\epsilon = 0.4$ run at lab scale ($H = 0.2$ m). Colors represent the velocity magnitude and vectors represent the velocity direction. The domain $x < 0$ corresponds to the glacier wall.

### E.2. Animation 2: Field_scale_capsize_3dview.mp4

We provide an animation of the water velocity field from a 3D perspective. The case is a 3D grounded glacier with iceberg aspect ratio $\epsilon = 0.4$ run at field scale ($H = 200$ m). Colors represent the velocity magnitude. Note that the color scale is exponential with a minimum set at $||u|| = 0.1$ m/s which is the typical velocity of water currents in fjords (Johnson et al., 2011; Mortensen et al., 2014; Sutherland et al., 2014). The domain $x < 0$ corresponds to the glacier wall.

## Supplementary material

Supplementary material associated with this article can be found, in the online version, at 10.1016/j.oceaneng.2025.121765.

## References

* Amundson, J.M., Burton, J.C., Cathles, M., MacAyeal, D., Glotter, M., Zhang, W., McCormik, M., O'Hanrahan, E., 2012a. Hydrodynamics of capsize. Video trailer for a poster presentation at the 2012 International Glaciological Society meeting on 'Glaciers and ice sheets in a warming climate' in Fairbanks, Alaska. Online; accessed 02 July 2024. https://www.youtube.com/watch?v=WrJ8vHAG_al&t=3s

* Amundson, J.M., Burton, J.C., Correa-Legisos, S., 2012b. Impact of hydrodynamics on seismic signals generated by iceberg collisions. Ann. Glaciol. 53 (60), 106-112.

* Amundson, J.M., Fahnestock, M., Truffer, M., Brown, J., Lüthi, M.P., Motyka, R.J., 2010. Ice mélange dynamics and implications for terminus stability, Jakobshavn Isbræ, Greenland. J. Geophys. Res. 115. https://doi.org/10.1029/2009JF001405

* Amundson, J.M., Truffer, M., Lüthi, M.P., Fahnestock, M., West, M., Motyka, R.J., 2008. Glacier, fjord, and seismic response to recent large calving events, Jakobshavn Isbræ, Greenland. Geophys. Res. Lett. 35

* Amundson, J.M., Truffer, M., Zwinger, T., 2022. Tidewater glacier response to individual calving events. J. Glaciol. 68, 1117-1126.

* Bamber, J.L., Westaway, R.M., Marzeion, B., Wouters, B., 2018. The land ice contribution to sea level during the satellite era. Environ. Res. Lett. 13 (6), 063008. https://doi.org/10.1088/1748-9326/aac2f0

* Bartholomaus, T.C., Larsen, C.F., West, M.E., O'Neel, S., Pettit, E.C., Truffer, M., 2015. Tidal and seasonal variations in calving flux observed with passive seismology. J. Geophys. Res. 120 (11), 2318-2337. https://doi.org/10.1002/2015jf003641

* Benn, D.I., Todd, J., Luckman, A., Bevan, S., Chudley, T.R., Åström, J., Zwinger, T., Cook, S., Christoffersen, P., 2023. Controls on calving at a large Greenland tidewater glacier: stress regime, self-organised criticality and the crevasse-depth calving law. J. Glaciol., 1-16. https://doi.org/10.1017/jog.2023.81

* Berg, B., Bassis, J., 2022. Crevasse advection increases glacier calving. J. Glaciol., 1-10. https://doi.org/10.1017/jog.2022.10

* Bonnet, P., 2021. Mechanical modelling of the source of glacial earthquakes in polar regions. Ph.D. thesis. IPGP.

* Bonnet, P., Yastrebov, V.A., Queutey, P., Leroyer, A., Mangeney, A., Castelnau, O., Sergeant, A., Stutzmann, E., Montagner, J.-P., 2020. Modelling capsizing icebergs in the open ocean. Geophys. J. Int. 223, 1265-1287.

* van den Broeke, M.R., Enderlin, E.M., Howat, I.M., Kuipers Munneke, P., Noël, B. P.Y., Berg, W.J.vande, van Meijgaard, E., Wouters, B., 2016. On the recent contribution of the Greenland ice sheet to sea level change. Cryosphere 10 (5), 1933-1946. https://doi.org/10.5194/tc-10-1933-2016

* Burman, E., Fernández, M.A., Frei, S., Gerosa, F.M., 2022. A mechanically consistent model for fluid-structure interactions with contact including seepage. Comput. Methods Appl. Mech. Eng. 392, 114637. https://doi.org/10.1016/j.cma.2022.114637

* Burman, E., Fernández, M.A., Frei, S., 2020. A nitsche-based formulation for fluid-structure interactions with contact. ESAIM 54, 531-564. https://doi.org/10.1051/m2an/2019072

* Burton, J.C., Amundson, J.M., Abbot, D.S., Boghosian, A., Cathles, L.M., Correa-Legisos, S., Darnell, K.N., Guttenberg, N., Holland, D.M., MacAyeal, D.R., 2012. Laboratory investigations of iceberg capsize dynamics, energy dissipation and tsunamigenesis. J. Geophys. Res. 117. https://doi.org/10.1029/2011JF002055

* Bézu, C., Bartholomaus, T.C., 2024. Greenland ice sheet's distinct calving styles are identified in terminus change timeseries. Geophys. Res. Lett. 51 (21). https://doi.org/10.1029/2024g110224

* Cassotto, R., Fahnestock, M., Amundson, J.M., Truffer, M., Joughin, I., 2015. Seasonal and interannual variations in ice mélange and its impact on terminus stability, Jakobshavn Isbræ, Greenland. J. Glaciol. 61 (225), 76-88. https://doi.org/10.3189/2015jog13j235

* Cassotto, R.K., Burton, J.C., Amundson, J.M., Fahnestock, M.A., Truffer, M., 2021. Granular decoherence precedes ice mélange failure and glacier calving at jakobshavn isbræ. Nat. Geosci. 14 (6), 417-422. https://doi.org/10.1038/s41561-021-00754-9

* Chen, F., Heller, V., Briganti, R., 2020. Numerical modelling of tsunamis generated by iceberg calving validated with large-scale laboratory experiments. Adv. Water Resour. 142, 103647. https://doi.org/10.1016/j.adwatrems.2020.103647

* Cheng, D., Hayes, W., Larour, E., Mohajerani, Y., Wood, M., Velicogna, I., Rignot, E., 2021. Calving front machine (CALFIN): glacial termini dataset and automated deep learning extraction method for Greenland, 1972-2019. Cryosphere 15 (3), 1663-1675. https://doi.org/10.5194/tc-15-1663-2021

* Christmann, J., Helm, V., Khan, S.A., Kleiner, T., Müller, R., Morlighem, M., Neckel, N., Rückamp, M., Steinhage, D., Zeising, O., et al., 2021. Elastic deformation plays a non-negligible role in Greenland's outlet glacier flow. Commun. Earth Environ. 2 (1), 1-12.

* Cook, S.J., Christoffersen, P., Todd, J., 2021a. A fully-coupled 3d model of a large greenlandic outlet glacier with evolving subglacial hydrology, frontal plume melting and calving. J. Glaciol. 68 (269), 486-502. https://doi.org/10.1017/jog.2021.109

* Cook, S.J., Christoffersen, P., Truffer, M., Chudley, T.R., Abellán, A., 2021b. Calving of a large greenlandic tidewater glacier has complex links to meltwater plumes and mélange. J. Geophys. Res. 126 (4). https://doi.org/10.1029/2020JF006051

* Crawford, A.J., Benn, D.I., Todd, J., Aström, J.A., Bassis, J.N., Zwinger, T., 2021. Marine ice-cliff instability modeling shows mixed-mode ice-cliff failure and yields calving rate parameterization. Nat. Commun. 12 (1), 1-9.

* Cuffey, K.M., Paterson, W. S.B., 2010. The Physics of Glaciers. Academic Press.

* van Dongen, E., Jouvet, G., Walter, A., Todd, J., Zwinger, T., Asaji, I., Sugiyama, S., Walter, F., Funk, M., 2019. Tides modulate crevasse opening prior to a major calving event at Bowdoin glacier, Northwest Greenland. J. Glaciol. 66 (255), 113-123. https://doi.org/10.1017/jog.2019.89

* van Dongen, E. C.H., Åström, J.A., Jouvet, G., Todd, J., Benn, D.I., Funk, M., 2020. Numerical modeling shows increased fracturing due to melt-undercutting prior to major calving at Bowdoin glacier. Front. Earth Sci. 8. https://doi.org/10.3389/feart.2020.00253

* Engels, T., Kolomenskiy, D., Schneider, K., Sesterhenn, J., 2015. Numerical simulation of fluid-structure interaction with the volume penalization method. J. Comput. Phys. 281, 96-115. https://doi.org/10.1016/j.jcp.2014.10.005

* Fried, M.J., Catania, G.A., Bartholomaus, T.C., Duncan, D., Davis, M., Stearns, L.A., Nash, J., Shroyer, E., Sutherland, D., 2015. Distributed subglacial discharge drives significant submarine melt at a Greenland tidewater glacier. Geophys. Res. Lett. 42 (21), 9328-9336. https://doi.org/10.1002/2015gl065806

* Gagliardini, O., Zwinger, T., Gillet-Chaullet, F., Durand, G., Favier, L., de Fleurian, B., Greve, R., Malinen, M., Martín, C., Råback, P., Ruokolainen, J., Sacchettini, M., Schäfer, M., Seddik, H., Thies, J., 2013. Capabilities and performance of elmer/ice, a new-generation ice sheet model. Geosci. Model Dev. 6 (4), 1299-1318. https://doi.org/10.5194/gmd-6-1299-2013

* Glowacki, O., 2022. Distinguishing subaerial and submarine calving with underwater noise. J. Glaciol. 68 (272), 1185-1196. https://doi.org/10.1017/jog.2022.32

* Glowacki, O., Deane, G.B., Moskalik, M., Blondel, P., Tegowski, J., Blaszczyk, M., 2015. Underwater acoustic signatures of glacier calving. Geophys. Res. Lett. 42 (3), 804-812.

* Grilli, S.T., Ioualalen, M., Asavanant, J., Shi, F., Kirby, J.T., Watts, P., 2007. Source constraints and model simulation of the december 26, 2004, indian ocean tsunami. J. Waterway, Port, Coastal, Ocean Eng. 133 (6), 414-428. https://doi.org/10.1061/(asc)0733-950x(2007)133:6<1414>

* Grosfeld, K., Schröder, M., Fahrbach, E., Gerdes, R., Mackensen, A., 2001. How iceberg calving and grounding change the circulation and hydrography in the filchner ice shelf-ocean system. J. Geophys. Res. 106 (C5), 9039-9055. https://doi.org/10.1029/2000jc000601

* Gräff, D., Walter, F., 2021. Changing friction at the base of an alpine glacier. Sci. Rep. 11, https://doi.org/10.1038/s41598-021-90176-9

* Hanna, E., Fettweis, X., Mernild, S.H., Cappelen, J., Ribergaard, M.H., Shuman, C.A., Stefken, K., Wood, L., Mote, T.L., 2013. Atmospheric and oceanic climate forcing of the exceptional Greenland ice sheet surface melt in summer 2012. Int. J. Climatol. 34 (4), 1022-1037. https://doi.org/10.1002/joc.3743
Hanna, E., Huybrechts, P., Cappelen, J., Steffen, K., Bales, R.C., Burgess, E., McConnell, J.R., Peder Steffensen, J., Van den Broeke, M., Wake, L., Bigg, G., Griffiths, M., Savas, D., 2011. Greenland ice sheet surface mass balance 1870 to 2010 based on twentieth century reanalysis, and links with global climate forcing: Greenland ice sheet mass balance. J. Geophys. Res. 116 (D24), n/a-n/a. https://doi.org/10.1029/2011jd016387

Hanson, B., Hooke, R.L., 2000. Glacier calving: a numerical model of forces in the calving-speed/water-depth relation. J. Glaciol. 46 (153), 188-196. https://doi.org/10.3189/172756500781832792

Haseloff, M., Sergienko, O.V., 2022. Effects of calving and submarine melting on steady states and stability of buttressed marine ice sheets. J. Glaciol. 68 (272), 1149-1166. https://doi.org/10.1017/jog.2022.29

Heller, V., Bruggemann, M., Spinneken, J., Rogers, B.D., 2016. Composite modelling of subaerial landslide-tsunamis in different water body geometries and novel insight into slide and wave kinematics. Coastal Eng. 109, 20-41.

Heller, V., Chen, F., Brühl, M., Gabl, R., Chen, X., Wolters, G., Fuchs, H., 2019. Large-scale experiments into the tsunamigenic potential of different iceberg calving mechanisms. Sci. Rep. 9, 1-10.

Hester, E.W., Vasil, G.M., Burns, K.J., 2021. Improving accuracy of volume penalised fluid-solid interactions. J. Comput. Phys. 430, 110043. https://doi.org/10.1016/j.jcp.2020.110043

Hou, G., Wang, J., Layton, A., 2012. Numerical methods for fluid-structure interaction—A review. Commun. Comput. Phys. 12 (2), 337-377. https://doi.org/10.4208/cicp.291210.290411s

How, P., Schild, K.M., Benn, D.I., Noormets, R., Kirchner, N., Luckman, A., Vallot, D., Hulton, N. R.J., Borstad, C., 2019. Calving controlled by melt-under-cutting: detailed calving styles revealed through time-lapse observations. Ann. Jglaciol. 60 (78), 20-31. https://doi.org/10.1017/aog.2018.28

IPCC, 2022. Sea Level Rise and Implications for Low-Lying Islands, Coasts and Communities. Cambridge University Press. p. 321-446. p 336 each part contribution to SLR until 2015.

James, T.D., Murray, T., Selmes, N., Scharrer, K., O'Leary, M., 2014. Buoyant flexure and basal crevassing in dynamic mass loss at helheim glacier. Nat. Geosci. 7 (8), 593-596. https://doi.org/10.1038/ngeo2204

Johnson, H.L., Münchow, A., Falkner, K.K., Melling, H., 2011. Ocean circulation and properties in petermann fjord, Greenland. J. Geophys. Res. 116 (C1). https://doi.org/10.1029/2010jc006519

Jones, R.S., Gudmundsson, G.H., Mackintosh, A.N., McCormack, F.S., Whitmore, R.J., 2021. Ocean-driven and topography-controlled nonlinear glacier retreat during the holocene: southwestern ross sea, Antarctica. Geophys. Res. Lett. 48 (5). https://doi.org/10.1029/2020gl091454

Jouvet, G., Weidmann, Y., Seguinot, J., Funk, M., Abe, T., Sakakibara, D., Seddik, H., Sugiyama, S., 2017. Initiation of a major calving event on the Bowdoin glacier captured by UAV photogrammetry. Cryosphere 11 (2), 911-921. https://doi.org/10.5194/tc-11-911-2017

Kneib-Walter, A., Lüthi, M.P., Funk, M., Jouvet, G., Vieli, A., 2022. Observational constraints on the sensitivity of two calving glaciers to external forcings. J. Glaciol., 1-16 https://doi.org/10.1017/jog.2022.74

Köhler, A., Myklebust, E.B., Mæland, S., 2022. Enhancing seismic calving event identification in svalbard through empirical matched field processing and machine learning. Geophys. J. Int. 230, 1305-1317. https://doi.org/10.1093/gji/gga117

Lea, J.M., Mair, D.W.F., Nick, F.M., Rea, B.R., van As, D., Morlighem, M., Nienow, P.W., Weidick, A., 2014. Fluctuations of a greenlandic tidewater glacier driven by changes in atmospheric forcing: observations and modelling of kangiaia nunataa sermia, 1859-present. Cryosphere 8 (6), 2031-2045. https://doi.org/10.5194/tc-8-2031-2014

Leroyer, A., Barré, S., Kobus, J.-M., Visonneau, M., 2008. Experimental and numerical investigations of the flow around an oar blade. J. Mar. Sci. Technol. 13 (1), 1-15. https://doi.org/10.1007/s00773-007-0256-7

Levermann, A., 2011. When glacial giants roll over. Nature 472 (7341), 43-44.

Lüthi, M.P., Vieli, A., 2016. Multi-method observation and analysis of a tsunami caused by glacier calving. Cryosphere 10 (3), 995-1002.

Lydersen, C., Assmy, P., Falk-Petersen, S., Kohler, J., Kovacs, K.M., Reigstad, M., Steen, H., Strøm, H., Sundfjord, A., Varpe, Ø., Walczowski, W., Weslawski, J.M., Zajaczkowski, M., 2014. The importance of tidewater glaciers for marine mammals and seabirds in svalbard, norway. J. Mar. Syst. 129, 452-471. https://doi.org/10.1016/j.jmarsys.2013.09.006

MacAyele, D.R., Abbot, D.S., Sergienko, O.V., 2011. Iceberg-capsize tsunamigenesis. Ann. Jglaciol. 52 (58), 51-56.

MacAyele, D.R., Scambos, T.A., Hulbe, C.L., Fahnestock, M.A., 2003. Catastrophic ice-shelf break-up by an ice-shelf-fragment-capsize mechanism. J. Glaciol. 49 (164), 22-36. https://doi.org/10.3189/172756503781830863

Macías, J., Mercado, A., González-Vida, J.M., Ortega, S., Castro, M.J., 2016. Comparison and Computational Performance of Tsunami-HySEA and MOST Models for LANTEX 2013 Scenario: Impact Assessment on Puerto Rico Coasts. Springer International Publishing. pp. 3973-3997. https://doi.org/10.1007/978-3-319-55480-8_16

Makkonen, L., Tikanmäki, M., 2014. Modeling the friction of ice. Cold Reg. Sci. Technol. 102, 84-93. https://doi.org/10.1016/j.coldregions.2014.03.002

Mattingly, K.S., Mote, T.L., Fettweis, X., 2018. Atmospheric river impacts on Greenland ice sheet surface mass balance. J. Geophys. Res. 123 (16), 8538-8560. https://doi.org/10.1029/2018jd028714

Menter, F.R., 1994. Two-equation eddy-viscosity turbulence models for engineering applications. AIAA J. 32 (8), 1598-1605. https://doi.org/10.2514/3.12149

Meredith, M.P., Inall, M.E., Brearley, J.A., Ehmen, T., Sheen, K., Munday, D., Cook, A., Retallick, K., Van Landeghem, K., Gerrish, L., Annett, A., Carvalho, F., Jones, R., Naveira Garabato, A.C., Bull, C. Y.S., Wallis, B.J., Hogg, A.E., Scourse, J., 2022.

Internal tsunamigenesis and ocean mixing driven by glacier calving in Antarctica. Sci. Adv. 8 (47). https://doi.org/10.1126/sciadv.add0720

Mills, A., 2008. The coefficient of friction, particularly of ice. Phys. Educ. 43 (4), 392-395. https://doi.org/10.1088/0031-9210/43/4/006

Minowa, M., Podolskiy, E.A., Jouvet, G., Weidmann, Y., Sakakibara, D., Tsutaki, S., Genco, R., Sugiyama, S., 2019. Calving flux estimation from tsunami waves. Earth Planet. Sci. Lett. 515, 283-290. https://doi.org/10.1016/j.epsl.2019.03.023

Minowa, M., Podolskiy, E.A., Sugiyama, S., Sakakibara, D., Skvarca, P., 2018. Glacier calving observed with time-lapse imagery and tsunami waves at glaciar perito Moreno, Patagonia. J. Glaciol. 64 (245), 362-376.

Mortensen, J., Bendtsen, J., Lennert, K., Rysgaard, S., 2014. Seasonal variability of the circulation system in a west Greenland tidewater outlet glacier fjord, godthåbsfjord (64°n): godthåbsfjord. J. Geophys. Res. 119 (12), 2591-2603. https://doi.org/10.1002/2014jf003267

Mote, T.L., 2007. Greenland surface melt trends 1973-2007: evidence of a large increase in 2007. Geophys. Res. Lett. 34 (22). https://doi.org/10.1029/2007gl031976

Mouginot, J., Rignot, E., Bjork, A.A., van den Broeke, M., Millan, R., Morlighem, M., Noël, B., Scheuchl, B., Wood, M., 2019. Forty-six years of Greenland ice sheet mass balance from 1972 to 2018. Proc. Natl. Acad. Sci. 116 (19), 9239-9244. https://doi.org/10.1073/pnas.1904242116

Murray, T., Nettles, M., Selmes, N., Cathles, L.M., James, T.D., Edwards, S., Martin, I., O'Farrell, T., Aspey, R., Rutt, I., Baugé, T., 2015a. Reverse glacier motion during iceberg calving and the cause of glacial earthquakes. Sciences 349, 305-308. https://doi.org/10.1126/science.aab0460

Murray, T., Nettles, M., Selmes, N., Cathles, L.M., Burton, J.C., James, T.D., Edwards, S., Martin, I., O'Farrell, T., Aspey, R., Rutt, I., Baugé, T., 2015b. Reverse glacier motion during iceberg calving and the cause of glacial earthquakes. Science 349 (6245), 305-308. https://doi.org/10.1126/science.aab0460

Murray, T., Selmes, N., James, T.D., Edwards, S., Martin, I., O'Farrell, T., Aspey, R., Rutt, I., Nettles, M., Baugé, T., 2015c. Dynamics of glacial calving at the ungrounded margin of helheim glacier, southeast Greenland. J. Geophys. Res. 120, 964-982. https://doi.org/10.1002/2015jf003531

Nettles, M., Ekström, G., 2010. Glacial earthquakes in Greenland and Antarc- itane. Annu. Rev. Earth Planet. Sci. 38 (1), 467-491. https://doi.org/10.1146/ annurev-earth-040809-152414

Nick, F.M., Luckman, A., Vieli, A., Veen, C.J. V.D., As, D.V., De Wal, R. S. W.V., Pattyn, F., Hubbard, A.L., Floricioiu, D., 2012. The response of petermann glacier, Greenland, to large calving events, and its future stability in the context of atmospheric and oceanic warming. J. Glaciol. 58 (208), 229-239. https://doi.org/10.3189/ 2012jog11j242

Nishizawa, B., Kanna, N., Abe, Y., Ohashi, Y., Sakakibara, D., Asaji, I., Sugiyama, S., Yamaguchi, A., Watanuki, Y., 2019. Contrasting assemblages of seabirds in the subglacial meltwater plume and oceanic water of Bowdoin fjord, northwestern Greenland. JES J. Mar. Sci. 77 (2), 711-720. https://doi.org/10.1093/jejsj/fzs213

Oksanen, P., Keinonen, J., 1982. The mechanism of friction of ice. Wear 78 (3), 315-324. https://doi.org/10.1016/0043-1648(82)90242-3

Olsen, K.G., Nettles, M., 2017. Patterns in glacial-earthquake activity around Greenland, 2011-13. J. Glaciol. 63 (242), 1077-1089.

Olsen, K.G., Nettles, M., Cathles, L.M., Burton, J.C., Murray, T., James, T.D., 2021. Improved estimation of glacial-earthquake size through new modeling of the seismic source. J. Geophys. Res. 126 (12), e2021JF006384.

Patankar, S., 1980. Numerical Heat Transfer and Fluid Flow. CRC Press. P129: pressure-velocity coupling (checkerboard decoupling), SIMPLE, SIMPLER,...

Peskin, C.S., 2002. The immersed boundary method. Acta Numer. 11, 479-517. https://doi.org/10.1017/s0962492902000077

Pirot, E., Hibert, C., Mangeney, A., 2023. Enhanced glacial earthquake catalogues with supervised machine learning for more comprehensive analysis. Geophys. J. Int. 236 (2), 849-871. https://doi.org/10.1093/gji/ggad402

Podolskiy, E.A., Murai, Y., Kanna, N., Sugiyama, S., 2021. Ocean-bottom seismology of glacial earthquakes: the concept, lessons learned, and mind the sediments. Seismol. Res. Lett. 92 (5), 2850-2865. https://doi.org/10.1785/0220200465

Podolskiy, E.A., Murai, Y., Kanna, N., Sugiyama, S., 2022. Glacial earthquake-generating iceberg calving in a narwhal summering ground: the loudest underwater sound in the arctic? J. Acoust. Soc. Am. 151, 6-16. https://doi.org/10.1121/10.0009166

Podolskiy, E.A., Walter, F., 2016. Cryoseismology. Rev. Geophys. 54 (4), 708-758. https://doi.org/10.1002/2016rg000526

Popinet, S., 2012. Adaptive modelling of long-distance wave propagation and fine-scale flooding during the tohoku tsunami. Nat. Hazards Earth Syst. Sci. 12 (4), 1213-1227. https://doi.org/10.5194/nhess-12-1213-2012

Poulain, P., Friant, A.L., Mangeney, A., Viroulet, S., Fernandez-Nieto, E., Castro Diaz, M., Peruzzetto, M., Grandjean, G., Bouchut, F., Pedreros, R., Komorowski, J.-C., 2022a. Performance and limits of a shallow-water model for landslide-generated tsunamis: from laboratory experiments to simulations of flank collapses at Montagne Pelée (Martinique). Geophys. J. Int. 233 (2), 796-825. https://doi.org/10.1093/gji/ggac482

Poulain, P., Le Friant, A., Pedreros, R., Mangeney, A., Filippini, A.G., Grandjean, G., Lemoine, A., Fernández-Nieto, E.D., Castro Díaz, M.J., Peruzzetto, M., 2022b. Numerical simulation of submarine landslides and generated tsunamis: application to the on-going mayotte seismo-volcanic crisis. C. R. Géosci. 354 (S2), 361-390. https://doi.org/10.5802/crgeos.138

Queutey, P., Visonneau, M., 2007. An interface capturing method for free-surface hydrodynamic flows. Comput. Fluids 36 (9), 1481-1510. https://doi.org/10.1016/j.compfluid. 2006.11.007

Rhie, C., Chow, W., 1982. A numerical study of the turbulent flow past an isolated airfoil with tailing edge separation. In: 3rd Joint Thermophysics, Fluids, Plasma and Heat Transfer Conference. American Institute of Aeronautics and Astronautics. https://doi.org/10.2514/6.1982-998
Rignot, E., Koppes, M., Velicogna, I., 2010. Rapid submarine melting of the calving faces of west Greenland glaciers. Nat. Geosci. 3 (3), 187-191. https://doi.org/10.1038/ ngeo765

Seddik, H., Greve, R., Sakakibara, D., Tsutaki, S., Minowa, M., Sugiyama, S., 2019. Re-
sponse of the flow dynamics of Bowdoin glacier, northwestern Greenland, to basal
lubrication and tidal forcing. J. Glaciol. 65 (250), 225-238. https://doi.org/10.1017/
jog.2018.106

Sergeant, A., Mangeney, A., Stutzmann, E., Montagner, J., Walter, F., Moretti, L., Castelnau, O., 2016. Complex force history of a calving-generated glacial earthquake derived from broadband seismic inversion. Geophys. Res. Lett. 43 (3), 1055-1065.
https://doi.org/10.1002/2015gl066785

Sergeant, A., Mangeney, A., Yastrebov, V.A., Walter, F., Montagner, J.-P., Castelnau, O., Stutzmann, E., Bonnet, P., Ralaiarisoa, V. J.-L., Bevan, S., et al., 2019. Monitoring Greenland ice sheet buoyancy-driven calving discharge using glacial earthquakes. Ann. Glaciol. 60 (79), 75-95.

Sergeant, A., Yastrebov, V.A., Mangeney, A., Castelnau, O., Montagner, J.-P., Stutzmann, E., 2018. Numerical modeling of iceberg capsizing responsible for glacial earthquakes. J. Geophys. Res. 123 (11), 3013-3033. https://doi.org/10.1029/ 2018jf004768

Shepherd, A., Ivins, E., Rignot, J., E., Smith, B., van den Broeke, M., 2020. Mass balance of the Greenland ice sheet from 1992 to 2018. Nature 579 (7798), 233-239. https://doi.org/10.1038/s41586-019-1855-2

Slater, D.A., Straneo, F., 2022. Submarine melting of glaciers in Greenland amplified by atmospheric warming. Nat. Geosci. 15 (10), 794-799. https://doi.org/10.1038/
s41561-022-01035-9

Sugiyama, S., Minowa, M., Schaefer, M., 2019. Underwater ice terrace observed at the front of glacier grey, a freshwater calving glacier in Patagonia. Geophys. Res. Lett. 46 (5), 2602-2609. https://doi.org/10.1029/2018gl081441

Sugiyama, S., Sakakibara, D., Tsutaki, S., Maruyama, M., Sawagaki, T., 2015. Glacier dynamics near the calving front of Bowdoin glacier, northwestern Greenland. J. Glaciol. 61 (226), 223-232. https://doi.org/10.3189/2015jog14j127

Sugiyama, S., Sawagaki, T., Fukuda, T., Aoki, S., 2014. Active water exchange and life near the grounding line of an antarctic outlet glacier. Earth Planet. Sci. Lett. 399, 52-60.
https://doi.org/10.1016/j.epsl.2014.05.001

Sutherland, D.A., Straneo, F., Pickart, R.S., 2014. Characteristics and dynamics of two major Greenland glacial fjords. J. Geophys. Res. 119 (6), 3767-3791. https://doi.org/10.1002/2013jc009786

Svennevig, K., Hicks, S.P., Forbriger, T., Lecocq, T., Widmer-Schnidrig, R., Mangeney, A., Hibert, C., Korsgaard, N.J., Lucas, A., Satriano, C., et al., 2024. A rockslide-generated tsunami in a Greenland fjord rang earth for 9 days. Science 385 (6714), 1196-1205.

Todd, J., Christoffersen, P., Zwinger, T., Råback, P., Benn, D.I., 2019. Sensitivity of a calving glacier to ice-ocean interactions under climate change: new insights from a 3-d full-stokes model. Cryosphere 13 (6), 1681-1694. https://doi.org/10.5194/ tc-13-1681-2019

Truffer, M., Motyka, R.J., 2016. Where glaciers meet water: subaqueous melt and its rel-
evance to glaciers in various settings. Rev. Geophys. 54, 220-239.

Tsai, V.C., Ekström, G., 2007. Analysis of glacial earthquakes. J. Geophys. Res. 112, https://doi.org/10.1029/2006jf000596

Tsai, V.C., Rice, J.R., Fahnestock, M., 2008. Possible mechanisms for glacial earthquakes. J. Geophys. Res. 113.

van der Veen, C.J., 2002. Calving glaciers. Prog. Phys. Geogr. 26 (1), 96-122. https://doi.
org/10.1191/030913302pp327ra

Veitch, S.A., Nettles, M., 2012. Spatial and temporal variations in Greenland glacial- earthquake activity, 1993-2010. J. Geophys. Res. 117 (F4). https://doi.org/10.1029/ 2012jf002412

Veen, C.J., 1996. Tidewater calving. J. Glaciol. 42, 375-385.

Veitch, S.A., Nettles, M., 2017. Assessment of glacial-earthquake source parameters. J. Glaciol. 63 (241), 867-876. https://doi.org/10.1017/jog.2017.52

Wackers, J., Deng, G., Leroyer, A., Queutey, P., Visonneau, M., 2012. Adaptive grid refinement for hydrodynamic flows. Comput. Fluids 55, 85-100. https://doi.org/10.1016/j.compfluid.2011.11.004

Wagner, T. J.W., James, T.D., Murray, T., Vella, D., 2016. On the role of buoyant flexure in glacier calving. Geophys. Res. Lett. 43, 232. https://doi.org/10.1002/2015gl067247

Wagner, T. J.W., Stern, A.A., Dell, R.W., Eisenman, I., 2017. On the representation of capsizing in iceberg models. Ocean Model. 117, 88-96.

Wagner, T. J.W., Straneo, F., Richards, C.G., Slater, D.A., Stevens, L.A., Das, S.B., Singh, H., 2019. Large spatial variations in the flux balance along the front of a Greenland tide-water glacier. Cryosphere 13 (3), 911-925. https://doi.org/10.5194/tc-13-911-2019

Walter, J.I., Box, J.E., Tulaczyk, S., Brodsky, E.E., Howat, I.M., Ahn, Y., Brown, A., 2012. Oceanic mechanical forcing of a marine-terminating Greenland glacier. Ann. Glaciol. 53 (60), 181-192. https://doi.org/10.3189/2012aog60a083

Wang, W., Zender, C.S., van As, D., Fausto, R.S., Laffin, M.K., 2021. Greenland surface melt dominated by solar and sensible heating. Geophys. Res. Lett. 48 (7). https://doi.org/10.1029/2020gl090653

Wehrlé, A., Lüthi, M.P., Vieli, A., 2023. The control of short-term ice mélange weakening episodes on calving activity at major Greenland outlet glaciers. Cryosphere 17 (1), 309-326. https://doi.org/10.5194/tc-17-309-2023

Wetter, S., Hibert, C., Mangeney, A., Stutzmann, E., 2024. Constructing a new catalogue of Greenland's iceberg calving events through seismic data analysis and machine learning. EGU General Assembly 2024, Vienna, Austria, 14-19 Apr 2024, EGU24-7576.
https://doi.org/https://doi.org/10.5194/egusphere-egu24-7576

Wolper, J., Gao, M., Lüthi, M.P., Heller, V., Vieli, A., Jiang, C., Gaume, J., 2021. A glacier-ocean interaction model for tsunami genesis due to iceberg calving. Commun. Earth Environ. 2, https://doi.org/10.1038/s43247-021-00179-7

Wróbel, G., Malgorzata, S., 2008. Influence of temperature on friction coefficient of low density polyethylene. J. Achiev. Mater. Manuf. Eng. 28.

Young, T.J., Christoffersen, P., Bougamont, M., Tulaczyk, S.M., Hubbard, B., Mankoff, K.D., Nicholls, K.W., Stewart, C.L., 2022. Rapid basal melting of the Greenland ice sheet from surface meltwater drainage. Proc. Natl. Acad. Sci. 119 (10). https://doi.org/10.1073/pnas.2116036119

Yvin, C., Leroyer, A., Visonneau, M., Queutey, P., 2018. Added mass evaluation with a finite-volume solver for applications in fluid-structure interaction problems solved with co-simulation. J. Fluids Struct. 81, 528-546. https://doi.org/10.1016/j.jfluidstructs.2018.05.008

Zemp, M., Huss, M., Thibert, E., Eckert, N., McNabb, R., Huber, J., Barandun, M., Machguth, H., Nussbaumer, S.U., Gärtpner-Roer, I., Thomson, L., Paul, F., Maussion, F., Kutuzov, S., Cogley, J.G., 2019. Global glacier mass changes and their contributions to sea-level rise from 1961 to 2016. Nature 568 (7752), 382-386. https://doi.org/10.1038/s41586-019-1071-0

Zheng, L., Cheng, X., Shang, X., Chen, Z., Liang, Q., Wang, K., 2022. Greenland ice sheet daily surface melt flux observed from space. Geophys. Res. Lett. 49 (6). https://doi.org/10.1029/2021gl096690

Zwally, H.J., Abdalati, W., Herring, T., Larson, K., Saba, J., Steffen, K., 2002. Surface melt-induced acceleration of Greenland ice-sheet flow. Science 297 (5579), 218-222.
https://doi.org/10.1126/science.1072708