{'fraction_non_alphanumeric': 0.036274673338893525, 'fraction_numerical': 0.017559077008618293, 'mean_word_length': 4.328316643953309, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Modelling iceberg capsize in the open ocean

At near-grounded glacier termini, calving can lead to the capsize of kilometer-scale (i.e. gigatons) unstable icebergs. The transient contact force applied by the capsizing iceberg on the glacier front generates seismic waves that propagate over teleseismic distances. The inversion of this seismic signal is of great interest to get insight into actual and past capsize dynamics. However, the iceberg size, which is of interest for geophysical and climatic studies, cannot be recovered from the seismic amplitude alone. This is because the capsize is a complex process involving interactions between the iceberg, the glacier and the surrounding water. For this reason, we propose an approach based on the computation of the source force with a mechanical model for various iceberg geometries and the estimation of the iceberg size from a comparison of the computed and seismically inverted time histories of the force. The capsize dynamics is captured by computational fluid dynamic (CFD) simulations for a simple scenario of capsize in the open ocean without glacier front. This approach allows assessing the complexity of the fluid motion around an iceberg. Expressing the results in terms of appropriate dimensionless variables, we show that laboratory scale and field scale capsizes can be directly compared. The capsize dynamics is found to be highly sensitive to the iceberg aspect ratio and to the fluid and ice densities. We present a semi-analytical simplified fluid-structure model (SAFIM) that can be implemented in solid mechanics models dealing with contact dynamics. This model accounts for hydrodynamic forces through calibrated drag and added-mass effects. We show that SAFIM significantly improves the accuracy of the modeled horizontal force and horizontal iceberg velocity compared with existing simplified models. SAFIM has been calibrated against the reference CFD simulations...


Introduction
A current concern in climate science is the estimation of the mass balances of glaciers and ice sheets. The Greenland ice sheet mass balance contributes significantly to sea level rise, accounting for about 15% of the annual global sea level rise between [Jay Zwally et al., 2011. However, it is difficult to draw conclusions on general trends given the high uncertainties [Lemke et al., 2007] in these estimations, notably due to difficulties in estimating and partitioning the ice-sheet mass losses [van den Broeke et al., 2009]. Ice mass balance can be determined by calculating the difference between the (i) surface mass balance, mainly determined by inland ice gains minus ice losses and (ii) dynamic ice discharge, mainly made up of submarine melting and iceberg calving [Van Den Broeke et al., 2016]. One third to one half of the ice mass losses of the Greenland ice sheet are due to dynamic ice discharge [Enderlin et al., 2014]. Note that dynamic ice discharge is a complex phenomenon, influenced by ocean and atmospheric forcing and by glacier geometry and dynamics [Benn et al., 2017].
When a glacier marine terminus approaches a near-grounded position, calving typically occurs through the capsize of glacier-thickness icebergs. An iceberg will capsize under buoyancy for a width-to-height ratio below a critical value [MacAyeal and Scambos, 2003]. In Greenland, buoyancy-driven iceberg capsize can cause a horizontal and vertical movement of the glacier tongue [Murray et al., 2015], potentially destabilizing or restabilizing the ice tongue depending on variations of the glacier's flotation state [Cassotto et al., 2019]. Iceberg capsize can also trigger ice-shelf break-off as observed in Antarctica [MacAyeal andScambos, 2003, Burton et al., 2013]. Capsizing icebergs can slowly push ice-mélange away from the glacier front [Amundson et al., 2010] or squeeze it and rigidify it [Peters et al., 2014], which can affect the stabilizing buttressing stresses produced by ice-mélange at the calving front [Todd and Christoffersen, 2014]. Note that during capsize, important volumes of fjord water are mixed [Burton et al., 2012], which can affect fish populations [Moon et al., 2016]. Moreover, icebergs deteriorate through various processes such as break-up, wave erosion and solar or submarine convection melting [Job, 1978, Savage, 2007, which releases freshwater that can potentially affect overturning ocean circulation [Vizcaino et al., 2015, Marsh et al., 2015. [Wagner et al., 2017] explain that icebergs mainly melt through wind-driven wave erosion that leads to lateral thinning and thus eventually a buoyancy driven capsize of the iceberg (width-to-height ratio below a critical value [MacAyeal and Scambos, 2003]). Moreover, iceberg drift model simulations have shown that capsizing icebergs live longer and drift further than non-capsizing icebergs [Wagner et al., 2017] .
When they capsize right after calving, icebergs exert an almost horizontal transient contact force on the glacier front. This force is responsible for the generation of magnitude M w = 5 earthquakes that are recorded globally [Ekström Göran, 2003] and can be recovered from seismic waveform inversion [Walter et al., 2012, Sergeant et al., 2016. Monitoring iceberg calving and capsizing with such "glacial earthquakes" is a promising alternative to satellite imagery or airborne optical and radar sensors that suffer from limited time and spatial coverage. There is no direct link between the size of an iceberg and the seismic signal [Tsai et al., 2008, Walter et al., 2012, Sergeant et al., 2016, Sergeant et al., 2018 and many different processes involving the interactions between the iceberg, glacier, bedrock, water and ice mélange contribute to the type of calving, earthquake magnitude and seismic waveform [Amundson et al., 2010, Tsai et al., 2008.
To investigate in detail the link between iceberg volume, contact force and the generated seismic signals, the use of a hydrodynamic model coupled with a dynamic solid mechanics model is mandatory. The iceberg-ocean interaction governs the iceberg capsize dynamics and thus controls the time-evolution of the contact force which is responsible for the seismic waveform and amplitudes. Full modelling of the glacier-ocean-bedrock-iceberg-ice-mélange system is beyond capabilities of most existing models because it requires complexe and costly coupling between solid mechanics, contact dynamics and fluid dynamics. The simplified models proposed in the literature to describe capsizing icebergs assume a 2D rectangular rigid iceberg under gravity and account for the buoyancy force, iceberg-glacier contact force and simplified hydrodynamic effects using either added-masses [Tsai et al., 2008] or pressure drag [Burton et al., 2012, Sergeant et al., 2018, Sergeant et al., 2019. These models have been proposed to describe a specific aspect of the capsize: its vertical and rotational motion [Burton et al., 2012] validated against laboratory experiments, or the horizontal force that icebergs exert on the glacier fronts [Tsai et al., 2008, Sergeant et al., 2018. To build a complete catalogue of seismogenic calving events that can be used for seismic inversion and iceberg characterization, the model must accurately describe the interaction between the iceberg and the ocean. At the same time, its formulation should either remain sufficiently simple or be based on the interpolation on the response surface constructed on numerous full-model simulations. In particular, the horizontal force and the torque exerted by the fluid on the iceberg should be modelled as accurately as possible, since it controls the horizontal contact force [Tsai et al., 2008, Sergeant et al., 2018.
The present paper aims to provide an understanding the complex interactions between a capsizing iceberg and the surrounding water using a reference fluid dynamics solver and to capture the main feature of this interaction in a simplified model formulation suitable for further computations of seismic wave sources. For this, we use a Computational Fluid Dynamics solver (ISIS-CFD Software for Numerical Simulations of Incompressible Turbulent Free-Surface Flows) to generate reference results for the capsizing motion. This model solves Reynolds Averaged Navier-Stokes equations (RANSE) and handles interactions between rigid solids and fluids with a free surface. This model is not yet able to compute dynamic contacts between solids. This state-of-art model has been largely validated for marine engineering cases [Visonneau, 2005, Queutey et al., 2014, Visonneau et al., 2016 but not yet applied to kilometresize objects like capsizing icebergs in the open ocean, which give rise to high vorticity. Before applying ISIS-CFD to field-scale iceberg capsize, we will evaluate how well it can reproduce smallscale laboratory experiments (typical dimension of 10 cm). We will compare ISIS-CFD simulations to laboratory experiments conducted by [Burton et al., 2012] using a polymer proxy for icebergs capsizing in a narrow fjord-like water tank.
To obtain a model that can be easily coupled with a solid mechanics model, we will then propose a simplified formulation (SAFIM, for Semi-Analytical Floating Iceberg Model) for the interactions between iceberg and water. In this model, introduced hydrodynamic forces account for water drag and added-masses, these two effects being considered uncoupled and complementary. Such a complementary description of hydrodynamic forces was initially proposed for modelling the effect of waves on vertical piles [Morison et al., 1950] and has been widely used for modelling the effect of waves and currents on bulk structures [Venugopal et al., 2009, Tsukrov et al., 2002. The SAFIM hydrodynamic forces involve some coefficients that need to be calibrated to represent the effects of the hydrodynamic flow on the capsize motion. To do so, SAFIM parameters are fitted against reference results provided by ISIS-CFD.
When an iceberg capsizes in contact with the glacier front, it exerts a force on the glacier until it starts to drift away which happens before the tilt angle reaches 90 • . Therefore, in the goal of future investigation of seismic waves generated during the iceberg-toglacier contact, we are mainly interested in the dynamics of a capsizing iceberg in the open ocean only until this value of the tilt angle -i.e. θ < 90 • .
The paper is organized as follows. Section 2 presents the reference ISIS-CFD fluid dynamics model and its results are compared with those of laboratory experiments from [Burton et al., 2012]. The complexity of the fluid motion surrounding the iceberg and the pressure on the iceberg are then discussed. The similarities between the laboratory-scale and field-scale simulations are then presented. In section 3, we present the SAFIM simplified model and we discuss the assumptions with respect to those of other models from the literature. In Section 4, ISIS-CFD and SAFIM are compared for different case-studies and error quantification and fitting of parameters are discussed. Section 5 is an overall discussion: comparison of previous models, the new SAFIM model and the reference ISIS-CFD model, followed by a comparison of the results of drag and added-masses models compared with values computed by ISIS-CFD and discussion of the sensitivity of the model to geophysical variables such as the iceberg geometry and the density of water and ice.
2 Insights into iceberg capsize physical processes from computational fluid dynamics simulations 2.1 Considered geometry In this paper, we investigate the capsize of unstable icebergs in the open ocean, i.e. away from the glacier front. Water and ice densities are noted ρ w and ρ i . For the field scale, we take the typical values of ρ w = 1025 kg/m 3 and ρ i = 917 kg/m 3 . For the laboratory experiments, the densities are set to ρ w = 997 kg/m 3 and ρ i = 920 kg/m 3 for direct comparison with [Burton et al., 2012]. The considered geometry is shown in Fig. 1. The upward vertical direction is denoted e z and the horizontal direction e x . The out-of-plane dimension L of the iceberg (i.e. along e y ) is assumed to be sufficiently large compared to the height H and the width W such that the problem can be considered as two-dimensional. This two-dimensional assumption is discussed in Section 5.6. : Iceberg geometrical characteristics with the iceberg centre of gravity G and the centre of buoyancy B. The forces acting on the moving iceberg are: the gravity force F g , the hydrostatic fluid force F s and the hydrodynamic pressure indicated by the pink shaded area.
In this two-dimensional set-up, icebergs are assumed to be rectangular with in-plane dimensions H and W and an aspect ratio defined as ε = W/H. Rectangular icebergs in a vertical position are unstable, i.e. will capsize spontaneously, for aspect ratios smaller than a critical value [MacAyeal and Scambos, 2003]: ε < ε c = 6 ρ i ρw−ρi ρ 2 w . This critical aspect ratio is ε c ≈ 0.75 for densities in the field and ε c ≈ 0.65 for the densities in laboratory experiments. For ε > ε c , icebergs are vertically stable and will not capsize unless initially tilted sufficiently [Burton et al., 2012].
The iceberg is assumed to be rigid, i.e. it does not deform elastically. The mass of an iceberg per unit of thickness along e y is given by m = ρ i H 2 ε. Points G and B are the centre of gravity of the iceberg and the centre of buoyancy, respectively. The iceberg position is described by the horizontal and vertical positions of G, noted x G and z G respectively, and by the inclination θ of the iceberg's longer side with respect to the vertical. H w is the water depth.

Presentation of the ISIS-CFD solver
The ISIS-CFD solver, developed by LHEEA in Nantes (France), is a state-of-the-art solver for the dynamics of mixed turbulent flows [Guilmineau et al., 2017, Guilmineau et al., 2018 interacting with solid and flexible bodies [Hay et al., 2006, Leroyer and Visonneau, 2005, Durand et al., 2014 and with a free surface Visonneau, 2007, Leroyer et al., 2011]. Today, it is one of the very few available software products capable of solving problems as complex as interactions between solids and a fluid with a free surface. The target applications of ISIS-CFD are in the field of marine engineering, e.g. modelling the dynamic interactions between a ship and surface waves [Visonneau, 2005, Queutey et al., 2014, Visonneau et al., 2016 or the complex flows and interactions involved in the global hull-oars-rower system in Olympic rowing [Leroyer et al., 2012]. ISIS-CFD solves the Reynolds Averaged Navier-Stokes Equations (RANSE) [Robert et al., 2019] and uses a few different turbulence models. For the specific application of iceberg capsize (Fig. 1), two different turbulence equations were tested and were found to give similar results: k-w [Menter, 1993] and Spalart-Allmaras [Spalart and Allmaras, 1992]. The code uses an adaptive grid refinement [Wackers et al., 2012] or an overset meshing method (mandatory to deal with large amplitude body motion close to a wall for example) to connect two non-conforming meshes. The mesh used here is a converged mesh with n = 43000 elements. An example of a typical mesh in illustrated The iceberg (in white) is in the centre of a squared domain (fine mesh in red) moving with the iceberg inside a background domain (coarser mesh in blue). The mesh is automatically refined around the air-water interface and also around the interface between the two domains to obtain an equivalent cell size around the overlapping region so as to minimize interpolation errors between domains) in Fig. 2 for an iceberg with ε = 0.5. The coupling between the solid and the fluid is stabilized with a relaxation method based on the estimation of the periodically updated added-mass matrix [Yvin et al., 2018].
In the field of application of the ISIS-CFD flow solver, the typical range of Reynolds numbers (Re, a dimensionless ratio of inertial forces to viscous forces) extends from 10 6 for model-scale ship flow to 10 9 for full-scale ships [Visonneau et al., 2006]. In this field of ship hydrodynamics, the viscous contribution to the force is as high as 50% of the total force and therefore the scale effect focuses strongly on the prediction of the turbulent boundary layer (relatively thinner in full-scale flows than in model test conditions with delayed flow separation) whereas the scale effect on the free-surface deformation is generally smaller except in the aft-body region due to modification of the pressure field [Mizine et al., 2009].
Here, the ISIS-CFD solver is applied to a significantly different geometry, type of motion and dimensions i.e. to the modelling of the capsize of rectangular structures with a typical dimension of one kilometre, such as icebergs. This application provides new challenges for the ISIS-CFD solver: high Reynolds numbers Re ≈ 10 9 , greater lengths and velocities together with massive separations, high vorticity and free surface motion. [Burton et al., 2012] Since ISIS-CFD simulations will be compared with laboratory data (Sections 2.4 and 2.5), we summarize here the technical details of the corresponding experiments conducted by [Burton et al., 2012]. The experimental results presented below were kindly provided by J.C. Burton.

Laboratory experiments conducted by
The laboratory experiments consist in the capsize of parallelepipedic plastic iceberg proxies of density ρ i = 920±1 kg/m 3 in a long and narrow fjord-like tank 244 cm long, 28 cm wide and 30 cm tall, filled with water at room temperature (ρ w = 997 kg/m 3 ). To assess the effect of water depth H w on iceberg-capsize dynamics, two types of experiments were conducted in which the water height was varied from 11.4 cm to 24.3 cm. The iceberg height was H = 10.3 cm, (H < H w ) and the width varied between W = 2.5 cm and W = 10.2 cm, corresponding to aspect ratios ε = W/H ranging between 0.25 and 1 (see Fig. 1). The length of the iceberg was L = 26.6 cm, which is slightly smaller than the tank width to reduce edge effects so that the flow can be considered as two-dimensional. The plastic icebergs were initially placed slightly tilted with respect to the vertical position and were held by hand close to hydrostatic equilibrium. When the surface of the water became still, several drops of dye were introduced around the plastic iceberg to visualize the water flow. Then the icebergs were released to capsize freely. The capsizes were recorded by a camera located outside the tank. Snapshots are shown in Fig. 3 (top row). Further experimental details can be found in [Burton et al., 2012]. A selection of four experiments are presented here, corresponding to aspect ratios ε = 0.246, ε = 0.374, ε = 0.496 and ε = 0.639.  [Burton et al., 2012] in the top row (a, b, c); colour map of the dimensionless velocity u = u/ √ gH, where u is the velocity field of the fluid, in the water with ISIS-CFD simulations in the middle row (d, e, f); dimensionless velocity along the surface of the iceberg with SAFIM simulations in the bottom row (g, h, i), for three moments of the capsize (indicated on top of each column). The time scale is calibrated such that t = 0 s corresponds to the first time when the iceberg reaches θ = 90 • , as in [Burton et al., 2012]. The iceberg aspect ratio is ε = 0.496. The floor and walls are not shown in the snapshots. An animated figure computed by ISIS-CFD in available in supplementary material [S1].

Description of fluid motion around capsizing icebergs with ISIS-CFD and laboratory experiments
We assess the patterns of the fluid motion around a capsizing iceberg. Laboratory experiments show, to some extent, the fluid motion using dye at some specific locations. The ISIS-CFD computational fluid dynamics model makes a valuable contribution to the understanding of the complex motion of the fluid surrounding a capsizing iceberg since it computes the whole velocity field in the fluid. Fluid velocity colour maps computed with ISIS-CFD are qualitatively compared with the images of the laboratory experiments in Fig. 3. We also show maps of the iceberg velocity computed along its surface using the calibrated SAFIM model (see Section 3). Results for the capsize of an iceberg of aspect ratio ε = 0.496 are shown for three different times: in Figs 3(a), 3(d) and 3(g) during capsize; in Figs 3(b), 3(e) and 3(h) just at the end of the capsize when the iceberg reaches the horizontal position for the first time (θ = 90 • ); and in Figs 3(c), 3(f) and 3(i) some time after the capsize. The arrows represent the dimensionless velocity |u | = |u|/ √ gH, where u is the velocity field of the fluid. We observe a good qualitative agreement between the position and inclination of the iceberg obtained by ISIS-CFD and the laboratory experiments. Note that the iceberg is submarine when it reaches θ ≈ 90 • for the first time (Figs 3b and 3e). The motion of the fluid -initially at rest -is visible all around the capsizing iceberg. At a distance of H/2, the fluid velocity is typically 1.5 cm/s (Fig. 3d). Large vortices, associated with the iceberg motion, are clearly visible throughout the capsize in Fig 3 (top and middle row). The intense fluid motion represents an important source of energy dissipation and thus slows down the motion of the iceberg. Note also that the iceberg moves towards the left during capsize. Quantitative observations are given in the next section.

Comparison of ISIS-CFD with laboratory experiments
We will now quantitatively compare ISIS-CFD with the laboratory experiments by [Burton et al., 2012]. In Fig. 4, results are provided for three different unstable icebergs: a thin  [Burton et al., 2012] ε = 0.246 iceberg, a medium ε = 0.374 iceberg and a thick ε = 0.639 iceberg. The horizontal position x G , vertical position z G and tilt angle θ are plotted against time. As the plastic icebergs were initially positioned by hand, [Burton et al., 2012] observed some variability in the experimental initial conditions and the overall results. To provide an estimate of the variability of the protocol, three experiments with identical plastic icebergs and the same (nominal) initial conditions were conducted for each aspect ratio. We selected these three aspect ratios because of the consistency of the experimental results. For example, for the ε = 0.496 iceberg (studied in Section 2.4), one of the three laboratory experiments showed a change in direction of the horizontal motion during the capsize, probably due to an initial non-zero horizontal velocity in a direction opposite to that induced by the rotational motion. For this reason, the set of experiments for the iceberg with ε = 0.496 was not considered here.
The initial conditions in ISIS-CFD were chosen to fit within the variability of the laboratory experiments. In ISIS-CFD simulations, the icebergs were tilted by a small angle of 0.5 • for the thin and medium iceberg (black curves in Figs 4d, 4e) and a larger angle of 15 • for the thicker iceberg (black curves in Fig. 4i). The icebergs were initially placed in hydrostatic equilibrium. The water level in the tank H w = 11.4 cm or H w = 24.3 cm was found to have a negligible effect on the iceberg motion: results are within the data variability shown in Fig. 4. Therefore, the experiments with a constant water depth H w = 24.3 cm are presented and this water level was set as a boundary condition in ISIS-CFD simulations.
We first analyse the motion of the iceberg during capsize. The iceberg is initially slightly tilted with respect to the vertical position and is in hydrostatic equilibrium. Once released, it tilts to reach the horizontal position with associated upward and sideward motion. It then rises out of the water in a rocking motion superimposed with a continuous horizontal displacement (Fig. 4). The thinner the iceberg, the quicker it moves in the horizontal direction with a quasi constant velocity at least for the first 1.5 s. Indeed, between time t = −0.5 s (beginning of the horizontal motion) and time t = 1.0 s, the icebergs move horizontally by a distance of about d = 5 cm for a ε = 0.246 iceberg Fig. 4(a), d = 3 cm for a ε = 0.374 iceberg Fig. 4(b) and d = 1.5 cm for a ε = 0.639 iceberg Fig. 4(c). This horizontal motion is an important aspect of the iceberg capsize that we would like to focus on here. Note that, besides gravity and buoyancy, the only horizontal external force acting on the capsizing iceberg is a force coming from the relative motion of water around the iceberg (air has a negligible effect here). These hydrodynamic forces are responsible for the horizontal iceberg motion. They need to be captured accurately by the model as they contribute considerably to the contact force generated between the iceberg and the glacier front when a just-calved iceberg capsizes [Tsai et al., 2008, Sergeant et al., 2018. Fig. 4 shows that ISIS-CFD results are in very good agreement with the laboratory data. For the studied icebergs with three different aspect ratios, all the values of the horizontal movement of the centre of gravity of the icebergs Figs. 4(a), 4(b) and 4(c) computed by ISIS-CFD fit within the range of the laboratory data. Slightly less accurate results are obtained for the vertical movement. However, results are still very close to the laboratory data. We quantify the discrepancy by calculating the relative error on ISIS-CFD data with respect to the closest laboratory data, each time the ISIS-CFD data is outside the laboratory data range. The maximum error is 2.7% and is reached for the vertical motion of the thinnest iceberg (Fig. 4d). For the vertical motion of the two thicker icebergs (e) and (f) and for the inclination of all three icebergs (Figs 4g, 4h and 4i), the error is lower than 1%. This error quantification is conservative since we do not take into account uncertainties on each laboratory measurement.
This slight discrepancy on the vertical and rotational motion computed by ISIS-CFD could be due to differences between the laboratory and simulation set-ups with regards to the 2-dimensional approximation (see section 5.6) and the initial conditions as discussed above. Another reason for this slight discrepancy could be related to the turbulence model treated by a RANS approach. Even if the generation of large vortices and separations are not initially induced by turbulent phenomena (Euler approach captures similar flow topologies), their evolution can be affected by turbulent effects for which the RANS approach is not specifically designed for.

From laboratory scale to field scale
In the previous section, ISIS-CFD simulations were shown to fit laboratory experiments very well. However, our aim is to reproduce the dynamics of the capsize of field-scale icebergs with dimensions of several hundred metres, i.e. four orders of magnitude larger than for the laboratory scale. Also, as pointed out by [Sergeant et al., 2018], the laboratory-scale Reynolds number Re = LU/ν ≈ 10 3 , is 6 orders of magnitude smaller than the characteristic Reynolds number Re ≈ 10 9 for the field scale (with L the typical length, U the typical speed and ν the dynamic viscosity of the fluid). Global viscous effects are expected to be more pronounced for laboratory-scale than for field-scale capsize. Therefore, the question is whether laboratory-scale experiments can be used to understand the kinematics of the field-scale iceberg calving.
We compare the horizontal force applied to the iceberg during its capsize computed by ISIS-CFD for the two cases: (1) a field-scale iceberg of height H = 800 m and (2) a laboratory-scale iceberg of height H = 0.1 m, all other parameters being the same: aspect ratio ε = 0.25, infinite water pool, same densities of the water, the ice and the air, taken to be equal to the field densities (Section 2.1). Results are given in Fig. 5 using dimensionless variables, i.e. a dimensionless horizontal force F x = F x /(mg) acting on the capsizing iceberg and a dimensionless time t = t/ H/g. We observed that the two curves corresponding to the two scales are very similar from the beginning of the movement until t ≈ 15.3, which corresponds to θ ≈ 90 • . This similarity between the forces at the laboratory and field scales can be explained using the Vaschy-Buckingham π theorem, assuming that the effect of viscosity is negligible, as detailed in Appendix A. After the buoyancy driven capsize, i.e. for times larger than t ≈ 15.3, the discrepancy between laboratory and field scales is larger and dimensions start to play a more important role. This discrepancy probably originates from the fact that after the buoyancy driven capsize, the iceberg motion is driven by the evolution of complex vortices and different gravity-waves dynamics. The other variables of the system (vertical force and torque acting on the iceberg and horizontal and vertical displacement and inclination of the iceberg) are also very similar for the laboratory and field scales.
Since it was demonstrated that the laboratory and field scales produce the same horizontal dimensionless force, in the remaining simulations we will present only the results for the laboratoryscale iceberg with H = 0.103 m and for densities in the field and in absence of sea floor. The laboratory scale was chosen because numerical convergence is easier to achieve in ISIS-CFD for the laboratory scale than for the field scale. The sensitivity of the capsize to the densities will t' 3=15.6 t' 4=17.5 Figure 5: The horizontal dimensionless force F = F x /(mg) that the water exerts on the iceberg is plotted against a dimensionless time t = t/ H/g for both laboratory (H = 0.1 m, in red) and field (H = 800 m, in blue) scales. The iceberg inclination θ (grey curve with scale on the right) is plotted against dimensionless time t . The vorticity fields around the iceberg at four different times are also shown: (a) corresponds to the initial phase of iceberg's motion t/ H/g ≈ 7.3, (b) corresponds to the maximal force in the left direction t/ H/g ≈ 12.2, (c) corresponds to the maximal force in the right direction t/ H/g ≈ 15.6 and (d) at time t/ H/g ≈ 17.1, corresponds to the oscillations of the iceberg after complete capsize. The pink vectors represents the local hydrodynamic pressure. The colour maps are discussed in Section 5.3 Variable name Dimensionless variable Table 1: Table of dimensionless variables, with N denoting the dimensionless variable related to N and with the iceberg linear mass density m = ρ i H 2 ε, G the center of gravity of the iceberg and the characteristic time T H = H/g. Note that the dimensionless forces and moments can also be written through a normalization by the characteristic mass m, length H and time T H with the following formulas: . See Section. 5.5 for a discussion on a non-dimensionalization with a dimensionless time T ρ,H that depends on the densities. be discussed in Section 5.5. Also, the depth of the sea floor was observed to have no effect on the capsize dynamics. Results for icebergs with different heights can be deduced with a factor of proportionality given by the non-dimensionalization in Table. 1.
In ISIS-CFD simulations and laboratory experiments, we observe five stages during the capsize: 1. In the initial phase (0 < t < 6), the horizontal force F x oscillates around zero with a negligible amplitude: about 1% of its extremum amplitude. This stage is the initiation of the capsize with buoyancy and gravity forces making the iceberg rotate and rise.
2. Then the absolute value of F x increases, first slowly and then faster until the first extremum at t ≈ 12.2. This is explained by the fact that the induced vertical and rotational velocities and accelerations of the iceberg produce a resisting hydrodynamic fluid force that has a nonzero horizontal component. The horizontal component of this hydrodynamic force is the only horizontal force acting on the iceberg and it induces a horizontal motion of the iceberg.
3. The absolute value of F x decreases to F x = 0 before t ≈ 14.4, which corresponds to θ ≈ 70 • (Fig. 5). This can be explained by the fact that the horizontal motion of the iceberg (explained in phase ii) triggers a horizontal resisting fluid force that compensates part of the hydrodynamic force (mentioned in phase ii).
4. The force F x then becomes positive and increases to an extremum at t ≈ 15.6, where the iceberg is horizontal θ = 90 • .
5. In the post-capsize phase (after t ≈ 15.6), F x oscillates around zero and is slowly damped to vanish at t 15. The iceberg rocks around θ = 90 • and drifts to the left.
The highest water velocities in the surrounding ocean are reached when the iceberg is close to θ = 90 • . The values of the dimensionless velocities around the iceberg are shown in Fig 6. We observe that for an iceberg of height H (here 800 m) : • the water flows at an average speed of ≈ 0.01 √ g H (here ≈ 88 cm/s) at a distance of about H (here 800 m) from the iceberg.
• the water flows at the average speed of ≈ 0.0005 √ g H (here ≈ 4.4 cm/s) at a distance of about 3.5 H (here 2800 m) of the iceberg.
• the highest velocity in the fluid during capsize is ≈ 0.5 √ g H (here ≈ 42 m/s) and is reached by the fluid close to the iceberg.
Moreover, we observe that the iso-lines for the velocities are roughly semi-circles centred on the iceberg, with a radius slightly higher in horizontal direction. The reference ISIS-CFD model has the advantage of being very accurate for fluid-structure interactions but it cannot readily model contacts between deformable solids. As explained in the introduction, we aim to construct a simpler fluid-structure interaction model that can be more easily coupled with dynamic solid mechanics models. Thus we propose a simple empirical model that can be used to estimate the horizontal force applied to a capsizing iceberg. This model was initiated by the work of [Sergeant et al., 2018, Sergeant et al., 2019 and is extended and validated in this study. One possible way to construct such a simple model of a capsizing iceberg, used in [Tsai et al., 2008], [Burton et al., 2012] and [Sergeant et al., 2018], consists in solving equations of a rigid iceberg motion subject to relevant forces and moments while discarding water motion. The problem is then greatly simplified at the expense of this strong approximation. The general equations of iceberg motion for such simplified models can be written in two dimensions as: (1) where I = ρ i H 4 ε(1 + ε 2 ) is the moment of inertia of the iceberg with respect to its centre of mass G and around an axis parallel to e y (multiplying by the iceberg thickness along e y gives the inertia for the three-dimensional case). Such a formulation accounts for the hydrostatic force F s and the corresponding torque M s computed at the centre of gravity G, the gravity force F g , overall hydrodynamic (or drag) effects expressed by the force F d and the associated torque M d and so-called added-masses m xx , m zz , I θθ , m xz , J zθ and J xθ that account for the mass of water that must be displaced during the iceberg motion (the added-masses will be discussed below). We also define the total forces and torque (including the action of the added-mass) on each coordinate: the horizontal force is F x = mẍ G , the vertical force is F z = mz G and the torque is M y = Iθ.
Hydrodynamically induced forces that oppose the motion of solids in a fluid are commonly called "drag forces" F d and the corresponding drag moment M d accounts for a particular distribution of the drag pressure along the solid. The drag force is usually normalized by the squared relative velocity between the fluid and the solid with a factor of fluid density and it acts in the opposite Added-mass (AM):  [Tsai and Ekström, 2007] and [Burton et al., 2012] for iceberg capsize modelling direction of this velocity. A comparison between the pressure drag and the negligible friction drag is discussed in Appendix 5.3.
When the fluid motion is not computed, the added-mass (AM) should also be included in the model. Added-masses add some additional inertia to the moving solid. This effect is known to be significant when the density of the solid and the fluid are comparable, such as for ice and water. The matrix of added-masses, which is symmetrical [Yvin et al., 2018, Molin, 2002, has the following form: Added-mass effects are of two types. First, a force associated with an added-mass can arise in a given direction due to an acceleration in that direction, which corresponds to the diagonal terms m xx , m zz and I θθ in Eq. (1-3). Second, an added-mass force can arise in a given direction due to an acceleration in another direction, which is accounted for by the coupled terms m xz , J zθ , J xθ .
Within this framework, models proposed by [Burton et al., 2012] and [Tsai et al., 2008], summarized in Table 2, differ in the way they account for the drag and the added-mass. In the formulation proposed by [Tsai et al., 2008], pressure drag and the associated moment are not considered and only some of diagonal terms are taken into account in the added-mass matrix. In the formulation by [Burton et al., 2012], added-mass effects are neglected. As for the drag effects, they are assumed to depend only on individual components of the velocity of the centre of gravity G and on the angular velocity, for example the drag along e x only depends on the velocityẋ G . As a consequence, both models predict that an iceberg initially at rest (ẋ G = 0,ẍ G = 0) will not experience any horizontal movement along e x during its capsize. As discussed above, this result contradicts experimental and ISIS-CFD results.

SAFIM: a new capsize model
To improve the above mentioned simplified models, we have developed a new model that continues along the lines of previous propositions by our group [Sergeant et al., 2018, Sergeant et al., 2019. In addition to the previously used drag formulation, this model uses a drag coefficient varying with the aspect ratio and integrates added-mass effects. As will be shown, SAFIM reproduces the main results of ISIS-CFD. For example, it predicts the horizontal movement of capsizing icebergs initially at rest.
A particular feature of the drag model in SAFIM is that it is based on the local velocity of the iceberg's surface. Therefore, the total drag force and the drag moment, evaluated as integrals of local pressure and moment over the submerged part of the iceberg Γ s , is described as: where n is the local normal to the surface of the iceberg and v n = |v · n| is the normal velocity, r is the local position vector and r G is the position vector of the iceberg's centre of mass; the wedge sign ∧ denotes a vector product. We consider a quadratic dependence of the local drag force on the normal velocity in the following. See Section 5.3 for a discussion of various drag models. The integral expressions for the drag force and moment are given in Appendix B. The factor α is the only adjustable parameter of the drag model. It can be used to adjust the model with respect to the reference ISIS-CFD simulations. We recall that in the original papers of [Sergeant et al., 2018, Sergeant et al., 2019, this factor α was set to α = 1 by default. However, due to the complexity of the fluid flow, the optimal value of α may change with the geometry of the iceberg.
This formulation, rather than attempting to describe the local pressure accurately, which is difficult based on geometrical considerations only (see [Sergeant et al., 2018] and Section 5.3), aims at providing a good approximation of the global forces and moments acting on the rotating iceberg. As opposed to the simplified drag model of [Burton et al., 2012], in which the drag force and moment depend only on the velocity of the centre of gravity, in this formulation the hydrodynamic forces F d and the torque M d depend on the iceberg's current configuration x G , z G , θ and on the three velocitiesẋ G ,ż G ,θ, since the local velocity v n is a function of all six state variables. This makes it possible to produce a horizontal force acting on the iceberg during capsize even for icebergs initially at rest. Another advantage is that a unique fit-parameter α is required to represent the drag effect, contrary to three independent fit parameters used in [Burton et al., 2012]. This makes it possible to easily generalize the model to more complex iceberg geometries. The closed form expressions of the drag forces and the moment are given in Appendix B.
As for the added-masses in SAFIM, we will consider two different possibilities: simplified addedmasses and computed added-masses.
The simplified added-masses option is based on analytical formulas for the diagonal terms of the added-mass matrix. The coupled terms of added-mass are taken to be equal to zero: m xz = 0, J zθ = 0, J xθ = 0. The formulas used are taken from [Wendel, 1956] for fully or partly submerged solids and were adapted to a capsizing body. The horizontal and vertical added-masses take the following forms: where H eff and W eff are the effective height and width defined as the projection of the submerged part of the iceberg along the vertical and horizontal axes, respectively (see Fig. 1). The integral equations of the simplified added-masses m xx and m zz are given in Appendix B. The effective height and width depend on the vertical and angular positions of the iceberg that evolve during the capsize. Therefore m xx and m zz also evolve during the capsize. On the other hand, the added moment of inertia I θθ is assumed to depend only on the height of the iceberg, so it remains constant during the capsize: In order to adjust the added-mass effect used in SAFIM to reproduce the reference ISIS-CFD results, we introduce three factors in the above equations: C x , C z and C θ . In the next section, we compute the optimal values of these coefficients for various aspect ratios by minimizing the mismatch between SAFIM and ISIS-CFD.
Computed added-masses are calculated for each position of the iceberg. This is done by applying a unit acceleration on the iceberg for the considered degree of freedom and solving eqs. (18) and (16) of [Yvin et al., 2018] on the fluid domain using a numerical method such as the finiteelement, finite-volume or boundary element method. The integration of the induced pressure field on the body gives the complete column of the symmetrical added-mass matrix corresponding to the considered degree of freedom, with both diagonal and coupled entries. Similarly to the simplified added-mass method, the values of the computed added-mass also depend on the iceberg position and they therefore evolve during the capsize. For the computed added-masses, the coupled terms are non-zero, giving rise to a coupling between horizontal, vertical and rotational accelerations.
To solve the motion equations (1-3) with SAFIM, the Störmer-Verlet integration scheme is used. Since SAFIM has only three degrees of freedom, the integration over time is very fast, only a few seconds compared to few hours for ISIS-CFD on a single CPU. The time step in SAFIM that ensures a sufficiently accurate results is ∆t = 0.1 s in the field scale and ∆t = 0.001 s in the laboratory scale. In both cases, this corresponds to a dimensionless time step of ∆t = ∆t/ H/g ≈ 0.01.

Studied cases and methodology for SAFIM coefficient calibrations
The validation of the proposed SAFIM model should be suited to the final objectives: (1) accurate reproduction of the forces exerted on the iceberg during capsize, especially the horizontal component and the moment since they determine the resulting contact force when the iceberg is in contact with the glacier; (2) ease of implementation in a finite element solver for simulation of the whole iceberg-glacier-bedrock-ocean system, (3) suitability of the model for the entire range of possible geometries of icebergs encountered in the field.
In this context, we consider 2D icebergs with rectangular cross-sections (Fig. 1). 3D effects are discussed in Section 5.6. We use typical densities observed in the field (Section 2.1). As will be discussed in Section 5, the considered density has a non-negligible effect on the calving kinetics.
For this study, we select four iceberg geometries with a wide range of aspect ratios: exactly the same four geometries as described in Section 2.3.
A small initial tilt angle is set to 0.5 • (except for the thickest iceberg). This does not change the overall motion of the iceberg compared to a zero tilt -i.e. initially vertical iceberg-but it affects the "destabilization" period during which the iceberg slowly rotates or oscillates around its vertical position [Sergeant et al., 2018]. For the thickest iceberg with ε = 0.639, the initial tilt angle is 15 • . This angle produced ISIS-CFD simulations that best matched the laboratory experiments. Moreover, for the same angle of inclination and height, the initial torque is smaller for larger icebergs [Wagner et al., 2017]. Also, ε = 0.639 is close to the critical value ε c = 0.65 corresponding to the limit for iceberg stability (Section. 2.1).
To compare SAFIM and ISIS-CFD results, we compute the mismatch in the time-evolution of the horizontal forces F x (L 2 norm) during the capsize. The phases of the capsize that we focus on are phases (ii) and (iii) (defined in Section 2.6). The reason we do not seek to perfectly model phase (i) with SAFIM is discussed in Section 5.2. Also, SAFIM is designed to model the capsize phase but cannot model the very end of the capsize (θ > 80 • ), i.e. phase (iv), nor the post-capsize phase (v). In these phases, forces induced by complex fluid motion, which are difficult to parametrize, are expected to dominate gravity and buoyancy forces.
For SAFIM with a drag force and no added-masses, the mismatch can be written: with F x the total horizontal force acting on the iceberg, t 1 such that F x (t 1 ) = 1/6F min and t 1 < t Fmin and t 2 such that F x (t 2 ) = 1/6F max and t 2 > t Fmin , with F min the first extremum of the force and t Fmin the time at which it occurs. The SAFIM response is shifted artificially in time so that the first extremum of the curve fits that of ISIS-CFD. This time shift ∆t is discussed in Section 5.2. For the SAFIM model with a drag force and added-masses, the mismatch is written: where t = 0 is the beginning of the simulation and t 3 the first time for which the horizontal force crosses zero F x (t 3 ) = 0 after t Fmin . The curve is not shifted in time as for E 1 .
Equations (10) and (11) are computed for a parametric space of drag and added-mass coefficient detailed in Appendix C. The optimized values of the coefficients α, C x , C z and C θ are chosen such that the errors E 1 and E 2 are minimum.

Effect of drag and added-mass on the iceberg horizontal force
We analyse the horizontal force produced by SAFIM in order to understand the effects of the drag force and added-masses on the dynamics of the iceberg capsize. To do so, the dynamics of an iceberg with an aspect ratio ε = 0.246 were simulated by SAFIM with or without drag and with or without added-masses (Fig. 7) and compared to the results of ISIS-CFD. We also assessed the sensitivity of the total horizontal force F x to the drag coefficient and added-masses in Fig. 7.
First, without the drag and without non-diagonal terms of the added-masses, the horizontal force predicted by SAFIM was equal to zero F x = 0, as expected. Next, we analysed the results of SAFIM with a drag force and without added-masses. Fig. 7 shows the horizontal force computed by SAFIM without added-mass for two values of the drag coefficient: α = 1 (purple curve) and α = 0.85 (blue curve). The value α = 0.85 is the optimized value of α obtained by minimizing the error E 1 (Eq. 10). The force has a slightly higher amplitude (around t ≈ 8.5) and duration for SAFIM with α = 1 than with the optimized drag coefficient α = 0.85. Whatever the value of the drag coefficient, the SAFIM horizontal force curve does not fit the ISIS-CFD curve over time: the amplitude and shape are very similar to the ISIS-CFD results, but the initiation process occurs faster with SAFIM. When the SAFIM curves are shifted in time by ∆t = 3.0 (cyan curve in Fig. 7), the previous SAFIM curve (with drag coefficient α = 0.85 and without added-masses) fits the ISIS-CFD curve well. The shape is the same and the error on the waveform is E 1 = 5.2%, with a relative error on the first force extremum of 4%. Then, we analysed the effects of simplified added-masses. The SAFIM simulation with drag and simplified added-mass coefficients all equal to 1 gives horizontal force results that are very different from the reference ISIS-CFD results: the duration of the capsize is largely overestimated and the amplitude is strongly underestimated (orange curve in Fig. 7). The optimized drag and added-mass coefficients that give a minimum error E 2 are α = 1.1 for the drag, C θ = 0.75 for the added moment of inertia and zero added-masses (C x , C z ) = (0, 0). The corresponding results (yellow curve in Fig. 7) are in a very good agreement with the reference results (E 2 = 10%) both for the shape and for the time corresponding to the force extremum t ≈ 11.45 (same as for ISIS-CFD). The added moment of inertia (associated with coefficient C θ ) slows down the initial rotation of the iceberg. However, the amplitude of the extremum is slightly underestimated by 8% compared to that predicted by ISIS-CFD. The accuracy of the formula of the simplified added-masses with coefficients equal to 1 is discussed in section 5.4.
The curve for SAFIM with computed added-masses and no drag force (dark green curve in Fig. 7) gives a curve that fits the reference results quite well in amplitude but not in time and predicts a huge second extremum with the wrong sign (negative instead of positive). However, when correcting the drag coefficient to α = 3.0, which minimizes the error E 2 , the curve fits better in time, reproducing the initial slow change of the force, but the amplitude and the shape still do not fit (black curve in Fig. 7). Similarly to the simplified added-masses, the computed added-masses slow down the initial rotation of the iceberg.
This analysis suggests that the drag force has mainly an effect on the amplitude of the first force extremum and that the added-masses have an effect on the duration of the initiation of the capsize. Also, optimized coefficients of drag and added-masses improve the model significantly compared with the case with all coefficients set to 1. Further discussions on the pros and cons of the SAFIM models are given in Section 5.1.

Effect of the iceberg's aspect ratio on SAFIM overall performance
We will now analyse the forces and the torque acting on the four selected geometries of icebergs computed by ISIS-CFD and SAFIM. The evolution of the dimensionless horizontal force F x , vertical force F z , torque M θ , horizontal displacement x G , vertical displacement z G and inclination θ obtained by ISIS-CFD and SAFIM are plotted in Fig. 8 for SAFIM results obtained with drag and without added-masses and in Fig. 9 for SAFIM results with drag and simplified added-masses. SAFIM models use optimized parameters indicated in Table 3 for each aspect ratio.
We will first discuss the sensitivity of the computed forces to the aspect ratios ε. We observe that the amplitude of the first extremum of both the horizontal force F x and the vertical force F z decreases with increasing aspect ratio. Also, the amplitude of the horizontal acceleration is equal toẍ G = F x /m = gF x , so the horizontal acceleration of the iceberg decreases with increasing aspect ratio. This is consistent with the observed slower horizontal displacement of icebergs with larger aspect ratios as reported in Section 2.5. Also, the durations of the first negative part of the horizontal force and first positive part of the vertical forces are similar for the four aspect ratios: about 5 units of dimensionless time.
The minimal error E 1 of SAFIM with the drag and without added-masses increases with the aspect ratio (from 5% for ε = 0.246 and up to 24% for ε = 0.639). We also found that the optimal drag coefficient α increases in an approximately affine way with the aspect ratio (Fig. 10a) as: α opt ( ) ≈ −1.6 + 8.8 with the coefficient of determination equal to R 2 = 0.98. This linear regression is valid within the range of studied aspect ratios 0.246 ≤ ≤ 0.639. Note that this formula should not be used for < 0.18 for which the drag coefficient would be negative, which is physically meaningless.
The minimal errors E 2 for SAFIM with the drag and the simplified added-masses are greater than the errors E 1 for SAFIM with drag without added-masses and with the artificial time shift as plotted in Fig. 10(b) for each aspect ratio. As for E 1 , the error E 2 increases with the aspect ratio (from 10% for ε = 0.246 up to 26% for ε = 0.639). Note that for all the four studied cases, optimization of the error requires keeping only one non-zero added-mass coefficient: the added moment of inertia coefficient. Indeed, the simplified added-masses are added to reproduce the duration of the force's growth as accurately as possible, without changing the shape of the rest of the force curve that is already well represented by the drag force model. Therefore, the simplified added-masses should allow the slow initiation of the rotation, which can be explained by an added moment of inertia of the surrounding fluid. Moreover, the torque M θ and the vertical force F z are better represented by SAFIM with simplified added-masses (Fig. 9a) than without added-mass and with an artificial time shift (Fig. 8a) for the four studied geometries.
We also compared SAFIM including drag and the computed added-masses with the reference ISIS-CFD model. This SAFIM simulation predicts the time and amplitude of the extremum of the force and the torque less accurately then the two previous SAFIM simulations: the error E 2 for SAFIM with drag and computed added-masses is E 2 > 34% for all the four studied cases (Fig. 10b). The curves for this less accurate model are not shown.

Discussion
First we will discuss the performance of SAFIM in Sections 5.1 and 5.2, then SAFIM modelling choices in Sections 5.3 and 5.4 and finally the sensitivity of the model to the main geophysical constraints in Sections 5.5, 5.6 and 5.7.

SAFIM performance and comparison with previous models from the literature
The advantages of the formulated and validated SAFIM model with drag and without addedmasses is that (1) it can be readily implemented in a finite element model like the one in [Sergeant et al., 2018], (2) it requires only one parameter, the drag coefficient α opt ( ) ≈ −1.6+8.8ε, it quite accurately reproduces the shape and amplitude of the horizontal force with an error that increases with the iceberg aspect ratio from 5% for ε = 0.246 up to 25% for ε = 0.639. The    drawback of this model is that it does not correctly simulate the kinetics of the iceberg capsize, especially the time needed to reach the peak force (see discussion in Section 5.2). In addition, the evolution of the torque and vertical force is not well reproduced. The advantage of SAFIM with drag and simplified added-masses is that it correctly reproduces the time of the force extremum (no shift in time is needed) and it reproduces the torque and the vertical force better than SAFIM with drag and no added-masses. The drawbacks of SAFIM with drag and simplified added-masses is that it underestimates the amplitude of the first force extremum by 10%. SAFIM with drag and computed added-masses gives less accurate results than the two other SAFIM versions. Considering that the computed added-masses are physical and accurate [Yvin et al., 2018], the drag force given in SAFIM is physically lacking since it does not make it possible to reproduce the dynamics of the iceberg.
Here, we discuss the possible implementation of added-mass in a finite element solver with the aim of reproducing glacier calving behaviours as discussed in [Sergeant et al., 2018, Sergeant et al., 2019. While adding an arbitrary added-mass matrix is feasible in the framework of dynamics of rigid solids, its use in a finite element framework is not straightforward. As the optimization analysis showed, the only term that needs to be added in the simplified added-masses is the added moment of inertia (Table 3), keeping the total mass and the outer shape of the iceberg unchanged. It is a challenging, if even possible, task to compute such a mass distribution in order to change the inertia of the iceberg accordingly. If such a redistribution of the iceberg density is not possible, the use of simplified added-masses might therefore require adding an external moment on the iceberg at each time step proportional to the acceleration weighted by a factor of the added moment of inertia. Using computed added-masses might also require adding external forces and moments on the iceberg at each time step proportional to the coupled terms of the added-mass matrix. One advantage of computed added-masses over simplified added-masses is that there is no adjustable parameter. However, at every time step, the added-mass matrix should be computed which considerably slows down and complexifies the simulation.
Finally, the proposed SAFIM model well predicts the first part of the horizontal force applied by the fluid on the iceberg, either when using a drag force only (i.e. no added-masses) and shifting the curve in time or when using a drag force and simplified added-masses (and no shift in time). However, the evolution of the force after the capsize (θ > 80 • ) is not well modelled by SAFIM. This is probably due to the fact that the evolution of the local pressure force is governed by a complex established fluid motion around the iceberg (see Section 5.3) which is hard to parametrize without full fluid dynamics computations. Given the global aim of modelling seismic waves generated by capsizing icebergs, this part should not play any role as it happens after the loss of the iceberg contact with the glacier front. An advantage of SAFIM over previous models [Tsai et al., 2008 is that, thanks to a special form of the drag force, it can describe the horizontal movement of a capsizing iceberg triggered by its rotation. As shown in [Sergeant et al., 2018], SAFIM can distinguish between a top-out and bottom-out capsize, when used to simulate the contact force between a capsizing iceberg and a rigid wall. A qualitative explanation of the emerging non-zero horizontal force given by the drag force is discussed in Appendix B.
We calculate the error in SAFIM with a drag coefficient α = 1 for all aspect ratios and without added-mass (Fig. 10). Exactly the same model was used in [Sergeant et al., 2018, Sergeant et al., 2019, but for modelling an iceberg capsizing in contact with a glacier. The error E 1 is about two times greater than when taking the optimum value of the α coefficient for each aspect ratio, and the amplitude and duration of the first negative part of the force is underestimated -except for the thinnest iceberg with ε = 0.246, for which the opposite is true, as mentioned in Section 4.2. Note that this error is relevant for a freely capsizing iceberg. In future work, the errors for an iceberg capsizing in contact with a glacier should be estimated but this will require a reference model for fluid-structure interactions that can track the contact between solids, which is a highly challenging problem.

Initiation phase of the capsize
In the previous sections, the drag parameter α for SAFIM with a drag force and without addedmasses was optimised by implementing an artificial time shift of the SAFIM force curve with respect to the ISIS-CFD force curve. This was done because, as already mentioned, SAFIM without addedmasses is not able to predict the accurate duration of the initiation of the capsize, where the motion is slow and the horizontal force is close to zero.
Various reasons suggest that this initial phase may not be relevant in the global context of the ultimate objective of the project, i.e. estimation of the volume loss on marineterminating icebergs. To achieve this objective, we need to compare the modelled contact force with the inverted seismic source force. The very beginning of the seismic force has too small a signal-to-noise ratio, therefore it is the first peak of the force that is used as a reference to compare the seismic force and the modelled force. Also, because this force evolves very slowly at the beginning, it will not be responsible for the generation of seismic waves with a period of 50 s that is predominantly observed on glacial earthquake seismograms [Ekström Göran, 2003, Tsai and Ekström, 2007, Tsai et al., 2008, Sergeant et al., 2018. Another reason for ignoring the beginning of the capsize is that the duration of the initial slow rotation (phase i) of the iceberg is strongly dependent on the initial angle of inclination of the iceberg which is hard to constrain in the field data and has little effect on the capsize (phases ii, iii, iv) when it is sufficiently small.
Nevertheless, if we consider a complete bedrock-glacier-iceberg-ocean system, the initial detachment of the iceberg can result in various other effects such as basal sliding or vertical oscillations of the glacier tongue, which can produce a seismic signal. Therefore the superposition of these phenomena can be erroneous if the timing is not well reproduced.
To solve this issue, simulations of the complete bedrock-glacier-iceberg-ocean system with a full fluid dynamics model coupled with a model for dynamics of deformable solids would seem to be unavoidable, however, as already discussed, it lies beyond actual computational possibilities.

Drag force and local pressure field
Following [Burton et al., 2012], a linear drag model with a local pressure proportional to the normal velocity |v n | was also tested in SAFIM. It results in the following modification to equations (5) and (6): Such a drag model yields poorer results than the original model with quadratic dependency when compared with the reference ISIS-CFD model. In addition, other drag models were tested with linear and quadratic pressure dependency on the velocity, with a non-uniform parameter α on the surface of the iceberg and with drag depending on the sign of the local normal velocity v n . Of all drag models tested, the most accurate was the model with quadratic dependency on the normal velocity and with a constant α-factor over the whole surface of the iceberg. However, to better fit the reference results, the α-factor was made dependent on the iceberg's aspect ratio ε, which is an important difference with the original model presented in [Sergeant et al., 2018].
To further our understanding of the forces generated by the fluid, we analyse the hydrodynamic pressure distribution on the sides of the iceberg, computed by ISIS-CFD and defined as p dyn = p tot − p sta , with p tot the total fluid pressure and p sta the hydrostatic pressure computed for the reference still water level (z = 0). In particular, we attempted to establish a link between the spatial distribution of the hydrodynamic pressure on the iceberg and the local features of the fluid flow, notably with the normalized vorticity (see Fig. 5) which is defined as: ω = − H/g e y · (∇ ∧ u), with a negative value (blue) accounting for a vortex rotating clockwise and a positive (red) value for a counter-clockwise vortex. On the four snapshots presented in Fig. 5, we also plot the dimensionless hydrodynamic pressure p dyn = p dyn /(ρ i H g). The hydrodynamic pressure vectors are plotted in shaded pink and are shown with a dense array of vectors with the origin on the iceberg's surface and pointing outwards from the iceberg for a negative pressure and inwards for a positive pressure. Note that these values are about two orders of magnitude lower than the average hydrostatic pressure. The dynamic pressure is higher at locations where there is a vortex close to the surface of the iceberg such as on the corner furthest right in Figs. 5(b) and 5(c), on the bottom part of the left side and in the middle on the right side of the iceberg. This correlation suggests that the dynamic pressure field is highly dependent on the vortices in the fluid. Such an evolution of complex vortex motion cannot be reproduced within SAFIM and requires the resolution of the equations of fluid motion as in ISIS-CFD. Note that the high values of the pressure on the top side of the iceberg in Figs 5(c) and 5(d) are due to an additional hydrostatic pressure produced by the wave that is above the reference sea level.
We also attempted to correlate, in the full CFD simulations, the local hydrodynamic pressure p dyn with the normal velocity v n via a power law as is the case in SAFIM: Note that in SAFIM, the coefficients are taken to be equal to a = 2 and b = −α sign(v n ) ρ w /2. The coefficients a and b that minimize the error between p dyn and b|v n | a (with a least square method on the log of equation 14) were calculated on portions of length W/2 along the sides of the iceberg and at various times during the simulated capsize (θ < 90 • ). We observed that the values of a and b vary along the sides of the iceberg and with time. Also, the accuracy of the fit varies along the sides: the accuracy is good on the bottom part of the long sides of the iceberg (least square coefficient R 2 > 0.9) and poor on the top part of the long sides of the iceberg and close to the corners of the iceberg (least square coefficient R 2 < 0.5). Also, we tried to correlate the dynamic pressure for a = 2, as in SAFIM, but there is no local correlation between the squared normal velocity and the dynamic pressure. Nevertheless, this choice ensures rather accurate overall drag forces and moments acting on the iceberg due to dynamic pressure.

Accuracy of added-mass values
The simplified added-masses, describing the diagonal terms of the added-mass matrix, will now be compared with the reference computed added-masses. Both added-mass matrices depend on the current configuration, more precisely on the z-coordinate of the centre of gravity G and on the tilt angle θ, and therefore they should be updated at every time step. These matrices are calculated for the iceberg's kinetics computed by ISIS-CFD. We show the time evolution of the diagonal terms: horizontal added-mass m xx in Fig. 11(a), the vertical added-mass m zz in Fig. 11(b) and the added moment of inertia I θθ in Fig. 11(c), for the capsize of a laboratory-scale iceberg with ε = 0.246 and H = 0.103m. The simplified horizontal and vertical added-masses are in very good agreement with the corresponding computed added-masses: relative error (with the L 2 norm) of 21% on m xx for ε = 0.246 and 23% for ε = 0.496; relative error of 11% on m zz for ε = 0.246 and 13% for ε = 0.496. The simplified added moment of inertia I θθ is assumed to be constant in our model whereas the computed added moment of inertia varies in time and has a smaller value: relative error (with the L 2 norm) of 48% for ε = 0.246 and 52% for ε = 0.496.
For the aspect ratio ε = 0.246, the horizontal added-mass m xx decreases from 2.5 m at the beginning down to ≈ 0.23 m, where m is the iceberg mass, whereas, symmetrically, the vertical added-mass m zz increases from 0.23 m at the beginning up to ≈ 2.5 m at the end of the capsize. The horizontal added-mass m xx measures the resistance of the fluid to a horizontal acceleration x G of the iceberg. The iceberg has a longer submerged vertical extension (of the order of H) at the beginning than at the end (of the order of W ) of the capsize, thus it needs to displace a greater volume of the fluid in a horizontal motion at the beginning than at the end of the capsize (θ > 90 • ). Therefore, a greater added-mass m xx is expected at the beginning of the capsize. The vertical added-masses m zz , sensitive to the horizontal extension of the iceberg, experience the opposite variations in time. Also, the values of m xx at the beginning of the capsize are similar to those of m zz at the end of the capsize because they is closely linked to the iceberg length in the x and z directions as mentioned above. Similarly, the values of m xx at the end of the capsize are similar to those of m zz at the beginning of the capsize. For the aspect ratio ε = 0.246, the computed added moment of inertia is equal to the moment of inertia of the iceberg I at the beginning. Then it increases to 1.25 I and decreases below the iceberg's moment of inertia to 0.9 I (see Fig. 11).
For the aspect ratio ε = 0.496, the variations of the dimensionless added-masses are different: m xx varies from 1.3 m to ≈ 0.4 m and m zz varies from 0.4 m to ≈ 1.4 m. As for the ε = 0.246 case, the value of m xx at the beginning (resp. end) of the capsize is close to that of m zz at the end (resp. beginning) of the capsize.
The added-masses (resp. inertia) are also of the same order of magnitude as the masses (resp. inertia) of the iceberg for other geometries. For the case of a two-dimensional thin ellipse with an aspect ratio of b/a, with a the along-flow dimension and b the cross-flow dimension, [Newman, 1999] gives the added-masses and added moment of inertia. For an ellipse with b/a = 0.2, the transverse added-mass is equal to 0.9 times the mass of the displaced volume of fluid (i.e. the submerged volume of the solid times the density of the fluid) and the added moment of inertia is equal to 0.7 times the inertia of the displaced volume of fluid. For similar densities for the fluid and the solid, the added-masses are close to 0.9 times the masses of the solid and the added moment of inertia is close to 0.7 times the inertia for the solid.
The optimized values of the coefficients for the added-masses in SAFIM (in terms of minimization of the error E 2 , see Section 4.1) are: C x = 0, C z = 0 for all ε, C θ = 0.75 for ε = 0.246 and C θ = 0.5 for ε = 0.496 (see Table 3). The optimized value C z = 0, is consistent with the choice of m zz = 0 in [Tsai et al., 2008], however it is not equal to the computed vertical added-mass. The optimized coefficient C x = 0 gives m xx = 0. The horizontal added-mass m xx from [Tsai et al., 2008] varies similarly to the computed added-mass but is less accurate than the simplified added-mass with C xx = 1 and totally different from the optimized value m xx = 0. The added moment of inertia I θθ with the formula in [Tsai et al., 2008] is constant throughout the capsize and different from the optimized added moment of inertia. However, the formula for added-masses and added moment of inertia from [Tsai et al., 2008] were given for the simulation of an iceberg capsizing in contact with a wall, which may significantly affect the values of the added-masses. For ε = 0.246, note that the optimized simplified added moment of inertia (C θ = 0.75) is close to the computed added moment of inertia. However, the simplified added moment of inertia is not in agreement with the computed added moment of inertia for higher ε. In other words, the values of the simplified added-mass matrix in SAFIM with the optimized parameters C x , C z and C θ are not in agreement with the reference computed added-masses. However, as discussed in section 5.1, the SAFIM model with the simplified added-mass matrix gives better results than SAFIM with the computed added-mass matrix.

Effect of water/ice densities
The laboratory experiments were conducted with water and ice densities slightly different from those in the field (see Section. 2.1). As shown in Fig. 12, the dynamics of the iceberg computed by ISIS-CFD with field densities is significantly different from those obtained with laboratory densities: the amplitude, duration and initiation of the capsize are very sensitive to changes in densities. The results for both field and laboratory densities were obtained for the laboratory scale (H = 0.103 cm and ε = 0.246). This sensitivity is also well demonstrated by SAFIM with a drag model and without added-masses. Note that no change in the drag coefficient α is needed to accurately reproduce this effect with SAFIM.
In Section 2.6, we pointed out the similarity between laboratory scale and field scale simulations if the same water and ice densities were used in both. To obtain the dimensionless variables, we used the time scale T H = H/g, length scale H and mass scale m (see Table. 1). However, as shown in Fig. 12, using different densities yields great differences in the horizontal force. Here, we explain how a simulation of a laboratory-scale iceberg with laboratory densities can be related to a simulation of a field-scale iceberg with field densities. We use the same approach as in Section 2.6, with length scale H and mass scale m; but we introduce a time scale depending on the densities  field SAFIM drag ( = 0.85), no AM, field ISIS-CFD, lab SAFIM drag ( = 0.85), no AM, lab Figure 12: Horizontal force acting on a capsizing iceberg computed by ISIS-CFD and SAFIM for two different sets of densities: laboratory densities (red and yellow curves, "ρ lab" in the legend) ρ w = 997 kg/m 3 and ρ i = 920 kg/m 3 and field densities (blue and cyan curves, "ρ field" in the legend) ρ w = 1025 kg/m 3 and ρ i = 917 kg/m 3 . No change in the optimized parameter α is needed to reproduce ISIS-CFD results by SAFIM.
as proposed by [Tsai et al., 2008] (ignoring the factor 2π): In Fig. 13, we plot the dimensionless horizontal force F x = F x T 2 /(mH) with respect to the dimensionless time t/T for time scale T = T H (Figs 13 a, b, c) and for time scale T = T ρ,H (Figs 13 d, e, f) and for three aspect ratios ε = 0.25 (Figs 13 a, d), ε = 0.375 (Figs 13 b, e) and ε = 0.5 (Figs 13 c, f). For a time scale T H , which does not involve densities, the dimensionless curves differ considerably whereas for T ρ,H , which takes the densities into account, the agreement is largely improved, especially for small aspect ratios.
Using a shift in time, the fit can be improved even further. Therefore, to upscale the laboratoryscale experiments to the field scale, a dimensionless time scale T ρ,H involving the difference in densities should be used rather than a simple scaling T H . To attempt to reduce the remaining mismatch, a more elaborate time scale, length scale and mass scale -involving the initial angle of inclination or another formula dependent on the dimensions and densities -should be introduced.
The effect of the densities is taken into account in the time scale T ρ,H through the square root of the ratio of the ice density divided by the water/ice density difference: for laboratory densities : for field densities : As densities seem to have a relatively large impact on capsize dynamics, more realistic water and ice densities, including their spatial heterogeneity, should probably be considered in future capsize models. Water density depends on salinity and temperature. For example, in the fjord of the Bowdoin glacier (northwest Greenland), water density may range between 1015 kg/m 3 and 1028 kg/m 3 [Ohashi et al., 2019, Sejr et al., 2017, Middelbo et al., 2018, Holding et al., 2019. Ice density is more difficult to evaluate as in situ measurements are rare. It depends on the volume fraction of air bubbles, which is for example around 20 − 30% for firn at ≈ 40 m in depth [Herron and Langway, 1980]. Ice density can also be increased if the iceberg contains melt pockets. The density of the iceberg may then be heterogeneous and can probably range between ≈ 600 kg/m 3 and ≈ 930 kg/m 3 (the density of pure ice at 10 • C being 918.7 kg/m 3 ).
With these ranges of ice and water field densities, the factor ρ i /(ρ w − ρ i ) varies between the extreme values ≈ 1.18 and ≈ 3.31, which corresponds to an even greater spread than in our lab/field comparison (equations 15 and 16). Therefore consideration of the effect of density and its variability has to be integrated in the inverse problem for iceberg volume estimation based on the seismic signal inversion.

3D effects
Capsizing icebergs have the following typical dimensions: full-glacier-height 500 m H 1000 m, width in the glacier's flow direction W 0.75H [MacAyeal and Scambos, 2003], width along the glacier's coast line, which is generally greater than the iceberg's height H L, with an upper limit equal to the glacial fjord width. However, as discussed above, in our modelling we neglect the effect of the third dimension on the dynamics of the capsizing iceberg. The first argument to support this simplification is that iceberg capsizing in a narrow fjord-like tank (laboratory experiments of [Burton et al., 2012]) is very well reproduced with the two-dimensional ISIS-CFD model (Section 2). In the field, icebergs capsize in fjords with much more complex geometries. For example, the fjord may be much wider (in the e y direction) than the iceberg which would yield a truly three-dimensional motion of the fluid. This motion of the fluid in the third direction can induce motion such as vortices on each side of the iceberg which may have an effect on the motion of the iceberg that has not yet been evaluated.

Effect of the iceberg geometry
This study was conducted with the assumption that the icebergs have a perfectly rectangular (parallelepipedic) shape and smooth surface. However, icebergs in the field have much more complex shapes. The freeboard of an iceberg has irregularities that can range from a scale larger than 100 m down to a scale less than 0.1 m [Landy et al., 2015]. The roughness of the submerged part of icebergs is poorly documented because of the difficulty in conducting suitable measurements. In future work, we could estimate the roughness of some well documented icebergs, such as the PII-B-1 tabular iceberg in Northwest Greenland scanned with a Reson 8125 multibeam sonar by [Wagner et al., 2014]. In fluid mechanics modelling, surface features have a great impact on the boundary layer close to the surface and in some cases also on the whole flow [Krogstad and Antonia, 1999]. A sensitivity analysis would be needed to assess the influence of the surface features and surface roughness on the dynamics of capsizing icebergs.
Furthermore, in our simulations, icebergs were initially in hydrostatic equilibrium. In [Sergeant et al., 2018], the effect of hydrostatic imbalance of the iceberg at the initiation of the capsize was assessed by varying the vertical position of the iceberg with respect to the water level. Hydrostatic imbalance results in a different evolution of the contact force and different dominant frequencies of generated seismic waves. This is supported by seismic observations of calving events.

Conclusion
In this study, we have improved the understanding of free iceberg capsize in open water through fluid-dynamics simulations (ISIS-CFD solver) validated against laboratory experiments [Burton et al., 2012]. In particular, we have shown the complexity of the fluid motion and the dynamics of the iceberg during capsize: vortices around the iceberg during and after capsize, motion of the fluid around the iceberg (velocity of ≈ 88 cm/s for a H = 800 m high iceberg at a distance H from the iceberg) , wave generation, iceberg submergence when reaching the horizontal position and a significant horizontal displacement of the iceberg during capsize. Moreover, we have shown that the non-dimensionalized horizontal force F x = F x /(ρ i H 2 εg) is invariant with the height H of the iceberg.
We have presented here a Semi-Analytical Floating Iceberg Model (SAFIM) and demonstrated its accuracy for various geometries and water/ice densities by comparing the results with the direct numerical CFD simulations made by ISIS-CFD. Our simple model is slightly more complex but more accurate than the one used in our previous study [Sergeant et al., 2018]: the new feature is that the drag parameter depends on the iceberg aspect ratio (affine function) to minimize the error with the reference CFD simulations. SAFIM has an error of 5% to 20% on the horizontal force F x (without added-masses) during the capsize phase for different geometries. An extension of this model to more complex iceberg shapes and to three dimensions is relatively straightforward. Different options are offered by SAFIM. For accurate modelling of the amplitude of the fluid forces, SAFIM should be used with drag and without added-masses. For accurate modelling of the time of the peak force, it should be used with a drag force and an added moment of inertia. In the global context of estimations of iceberg volume by analysis of seismic signals generated during iceberg capsize in contact with a glacier front, based on the discussion on the time-shift Section 5.2, SAFIM should be used with an optimized drag coefficient α and no added-masses. However, this model has been validated only for the case of capsize of an iceberg in the open ocean. Further validation should be conducted for the simulation of the capsize of an iceberg in contact with a glacier. In the geophysical context of modelling seismogenic iceberg capsize, further studies would help improve the model accuracy. Examples of such studies include (i) modelling of the full iceberg-ocean-glacier-bedrock system, which is computationally very challenging and (ii) sensitivity analysis of the iceberg dynamics to the iceberg shape, surface roughness and fjord geometry, which may be very complex. The characteristic variables of the system are the dimensions, H and W , the densities ρ i and ρ w , the water viscosity µ w and gravity g, so n = 6.
The n − p = 3 dimensionless ratios are chosen here to be: The calculation of the horizontal force F x (t) from the n = 6 independent characteristic variables of the problem can be written as: The Vashy-Buckingham -π theorem states that the problem can be written: To estimate the effect of viscosity, we compare the pressure and viscous forces. The fluid force on the surface of the iceberg calculated by ISIS-CFD is the sum of a friction-induced force (locally tangent to the fluid/solid interface) and a pressure-induced force (normal to this interface). In the case of an iceberg with aspect ratio = 0.25, the friction force is found to be ≈ 300 times smaller than the pressure force for the field-scale case (H = 800 m) and ≈ 10 times smaller for the laboratory case (H = 0.103 m) as illustrated in Fig. 14. Therefore, viscous effects can be reasonably neglected in both scales. This leads to the following approximation for equations (17): i.e. for similar initial conditions and boundary conditions, the evolution with time of the dimensionless force F x = F x mg only depends on the aspect ratio ε and the density ratio ρ w ρ i . However, the function f remains unknown and is investigated in Section 2.6 and Section 5.5.
B Integrated expressions of SAFIM hydrostatic force, drag force, torque and simplified added-masses B.1 Integrated expressions of SAFIM hydrostatic and drag forces and torque The integrated expressions for hydrostatic and drag forces and the associated moments are given below for SAFIM and for a rectangular iceberg as in Fig. 1. All these expressions are implemented in the Python code available online at [Yastrebov and Bonnet, 2020]. The effect of the hydrostatic pressure is given by the following integral: F s = −ρ w g Γs zn dΓ, where n is the outward surface normal and Γ s is the submerged part of the iceberg. The torque induced by this pressure distribution with respect to the centre of gravity G at position r G is given by: The drag force is: where v n is the local normal velocity component. The drag moment with respect to G is: The calculation of the integral of the pressure drag is split into integration over all submerged or partly submerged sides of the rectangular iceberg. Consider a partly submerged side S = AB and let us assume that corner A is a submerged corner and B is a corner outside the water. Then the velocity v M of a point M ∈ [AB] is: where r M = r A + ξ(r B − r A ) with ξ ∈ [0, ξ i ] and ξ i defining the position of the point on the given side located at the sea level. So for the side AB, the contribution of the drag force is given by where v n =ṙ G · n +θe y ∧ (r A + ξ(r B − r A ) − r G ) · n.
For the case of a totally submerged side then ξ i = 1. For the case of a side totally outside the water, the contribution to the drag force is zero.

B.2 Qualitative explanation of horizontal motion of the capsizing iceberg
With the formulation of the drag force given above, we can reproduce the horizontal motion of a freely capsizing iceberg, which is observed experimentally and reproduced with the accurate ISIS-CFD simulations. Obtaining a closed form solution of SAFIM equations Eqs. (1), (2) and (3) is out of reach. We wish to give here some intuitive explanation of the horizontal motion of the iceberg. The resultant of the buoyancy and gravity forces moves the iceberg upwards and makes it rotate: these two effects initiate the vertical and rotational motion of the iceberg. The induced velocity produces a force with a non-zero horizontal component. We now explain why these two initial motions -upwards and rotation-, together generate a horizontal drag force, in the framework of SAFIM. We draw the velocity v (triple red arrow) of several points on the surface of the iceberg and its normal component v n n (dashed red arrow). In SAFIM, the elementary drag force dF d (solid black arrow) is collinear with n and opposes the normal velocity v n n. The projection of the elementary drag forces on the horizontal axis dF dx = (dF d · e x )e x is shown by a dashed green arrow if it is leftward and dashed blue arrow if it is rightward. The integral of these horizontal elemental forces results in the global horizontal force F d · e x . For the case of upward motion Fig. 15(a) of the iceberg (Fig. 15(a)), the vertical local velocity is constant along the iceberg surface. First, we show that the contribution to the horizontal drag force of all the points on the two long sides of the iceberg is a leftward elementary force dF dx . As shown in Fig. 15, on both long sides on the iceberg, the projection of the velocity on the normal vectors v n n is rightward and the opposite vector, which gives the direction of dF d , has a negative horizontal component, i.e. dF dx is a leftward force. Second, on the small submerged side CD, the horizontal elementary force dF dx is rightward but with a smaller amplitude than on the long side because of the different inclination of that side compared to the long sides. This is valid as long as the iceberg has an aspect ratio smaller than 1 and is tilted less than θ = 45 • . This is the case here since we assess the initialization of the horizontal motion for a small initial inclination of the iceberg (θ i ≤ 15 in the simulations) and for unstable icebergs, i.e. icebergs with aspect ratios ε < ε c < 1 (see Section. 2.1). Additionally, the surface of the iceberg contributing to a leftward drag force (in green) is longer and has larger elementary horizontal forces than the surface of the iceberg contributing to a rightward force (in dark blue). Therefore, the horizontal force resulting from the upward motion is leftward. For the case of rotational motion ( Fig. 15(b)), the velocity increases with the distance to the centre of rotation. Because we assume no initial horizontal motion of the iceberg, G is initially immobile and the iceberg rotates around G. The further away the point is, the more it contributes to the drag force. Two points located at the same distance from G, but with opposite normal velocity vectors v n n (i.e. one point on the blue line and one point on the green line) have the same absolute contribution to the drag force but in opposite directions. Thus the drag force on the part of the surface coloured solid green compensates the drag force on the part coloured blue. The remaining part of the surface illustrated by dashed green lines, induces a leftward total horizontal force. Therefore, by superposing vertical and rotation motion, we obtain a net drag force in the direction of the initial tilt of the iceberg's top, here to the left.

B.3 Integrated expressions of simplified added-masses
The calculation of the simplified added-masses given in equation (7), requires the calculation of the effective height H eff and effective width W eff of the submarine part of the iceberg (see Fig. 1). To calculate them, we use the positions of the four corners: C, D, E, F (Fig. 15). The coordinates of a corner P ∈ {C, D, E, F } have the general expression : x P = z G + δ P 1 H 2 cos(θ) + δ P 2 L 2 sin(θ) with (δ P 1 ,δ P 2 ,δ P 3 ,δ P 4 ) defined as follows for the four corners The effective height can be calculated with the following expressions : H eff = max ((z w − z C ) ; (z w − z D ) ; (z w − z E ) ; (z w − z F )) where z w is the water level.
The effective width, defined as the distance between the leftmost and the rightmost points of the submerged part of the iceberg, is calculated similarly, but after checking which are the submerged corners and the geometrical intersection between the water surface and the iceberg sides.  Table 3: First two columns : geometrical characteristics and initial conditions of the studied icebergs. Laboratory-scale iceberg simulations have height H = 0.1 m and field-scale iceberg simulations have height H = 800 m. The density of the water is ρ w = 1025 kg m -3 and the density of the ice is ρ i = 917 kg m -3 . Next columns : parameters minimizing the error on F x and the corresponding error for SAFIM without added-masses, SAFIM with computed added-masses and SAFIM with simplified added-masses.