Numerical modeling of the stress-strain state of the ice beam by specified constitutive model

3.1. The demonstration of cantilever beam failure experiments

The cantilever beam failure test performed in the open air, with outdoor temperature at –2.0 °C in Novikbay in February, 2016 (Vladivostok), therefore, during the experiment process, the continuous increasing loads functioned on the ‘P’ point of cantilever beam (Fig. 1) till the beam failed around the fixed point (A-A section), the relevant parameters and beam dimensions recorded in Table 1.

Fig. 1Cantilever beam failure test scheme

Due to the experiments performed in the static state (quasi-static), the implicit algorithm will be applied during the numerical simulation process.

3.2. The constitutive model and relevant theories demonstrations:

General descriptions about the ice microstructure constitutive properties and characteristics lead researchers regard this material like the geotechnical material categories – ‘Coulomb-Mohr’, this expression is also explained by the CSA [9], Referring the ANSYS material library (ANSYS 2020) [10], for researchers modelling, it is possible that not only adopt ‘Coulomb-Mohr’ material model, but also choose ‘User Defined Model’ (temporary name ‘Exponent’) to be deployed on the ice beams, their basic mechanic characteristics and key parameters are expressed below respectively:

1) About the ‘Coulomb-Mohr’ model [10]. when the combination of pressure and shear stress reaches the cohesion of the material particles. Yielding occurs when the shear stress on any plane in the material reaches this criterion (on yield surface), the equation expressed below Eq. (1):

where $\tau$: shear strength, $C$: apparent cohesion, ${\sigma }_m$: mean stress, $\varphi$: angle of internal friction, when material failed ($\tau +{\sigma }_m\mathrm{t}\mathrm{a}\mathrm{n}\varphi \ge C$), considering the failure of beam caused by maximum shear stress on the beam, also means, when the beam failure, the maximum shear stress must be on the break surface, ${\tau }_{max}\ge C-{\sigma }_m\mathrm{t}\mathrm{a}\mathrm{n}\varphi$.

2) About the ‘User Defined Model’ [10]. the power law equation has a user-defined initial yield stress ${\sigma }_0$ and exponent $N$. the current yield stress ${\sigma }_y$ is given by solving the following equation Eq. (2):

where $G$ is the shear modulus determined from the user-defined elastic constants and ${\widehat{\epsilon }}^{pl}$ is the accumulated equivalent plastic strain.

3) About the flexural strength, Westergaard [11] proposed the flexural strength (${\sigma }_f$) formulations from the perspective of geometric dimensions on the circle area ice plate, Timco & Brien [12] also expressed ${\sigma }_f$ equations from the perspective of the brine volume fraction. In this study, we import the experiment break loads ($P$) to confirm ${\sigma }_f$, which is exactly conducted by the material mechanical principle Eq. (3):

where $I$ is the inertia moment, $b$ is the breadth of beam, $h$ is thickness of beam, $P$, exert loads, $L$ is the length of beam.

4) About the deflection, the floating ice plate can be thought of as an elastic plate resting on an elastic foundation, and its deflection is governed by the differential equation which highly depended on the ‘flexural rigidity of the plate’ and ‘biharmonic operator ${\nabla }^4$’, here we simplify the equation by classic material mechanic description as below Eq. (4):

Overall, the full set of ice cantilever beam mechanical and physical parameters extracted in Table 2.

4.1. Stress state analysis

Functioned numerical model with the ‘Coulomb-Mohr’ parameters from Table 1, when loads (P) reached the measured failure loads, the simulation max failure shear stress and mean stress occurred on Section A-A (Fig. 3(a), 3(b)), substituting the simulation values (read from Fig. 3(a) and Fig. 3(b)) into Eq. (1), results come that ${\tau }_{max}+{\sigma }_m\mathrm{t}\mathrm{a}\mathrm{n}\varphi =$ 395.56 KPa, compared with initial condition $C_1=$ 410 KPa, means ${\tau }_{max}+{\sigma }_m\mathrm{t}\mathrm{a}\mathrm{n}\varphi \approx C_1$, beam is at the critical point of failure at section A-A.

Fig. 3Stress simulation values with different constitutive models: a) max shear strength, b) mean stress, c) max normal stress

a)

b)

c)

For the user defined model (temporary name ‘Exponent’), it is easy to read the max normal stress on potential failure section A-A (Fig. 3(c)), ${\sigma }_{max}=$ 408 KPa, which is approximately equal to initial yield stress ${\sigma }_i=$ 410 KPa, the beam almost reached harden phrase, showed similar failure phenomenon like the previous model.

To be more precisely, focused on these models above, the marked section A-A is divided by several sub-sections showing in Fig. 4, extracting the normal stress distributions, then the 3D stress distributions charts can be established below (Fig. 5).

It is obvious that the 3D normal stress distribution charts reflected compressive and tensile values, meanwhile exhibited potential failure surface (section A-A), for these two models, most simulation stress values located in the range of (–230 KPa - 230 KPa). Apparently, due to the ‘Stress concentration effect’, a few exceeded or deviated simulation stress values appeared near the top and bottom areas. Overall, the main outline of stress distributions fit the material mechanic theory, albeit the peak values and outlines in the top (tensile) area showed different shapes.

Fig. 4Marked section A-A

Fig. 5The stress distribution (3D) from beam section A-A

4.2. Strain state analysis

Compared with swath of float ice plate, the beam scale is much smaller than that. Masterson, et al. have demonstrated the strain measurements [15], but not suitable for the small scale beams, in this article, the beam strain states can be expressed by equivalent and plastic state respectively from the ANSYS workbench simulation results. It is easily to observe that the maximum strain occurred on section A-A (Fig. 6) from general mechanical principle and simulation results. the dependency of Max strain with loads and flexural stress depicted in Fig. 7 respectively.

Fig. 6Max strain occurred on section A-A

In terms of the ‘Coulomb-Mohr’ model, equivalent strains (refer Fig. 7(a)), have steady increment till loads reached 1200 (N), around ${\sigma }_f$= 140 KPa, the precise values swelled from 0 to 2.8×10^-4, which indicates that strains performed elastic states, this state also can be verified in Fig. 7(b), where showed the plastic strain appeared when loads reached 1200 (N). during the loads continue growing up to the failure loads (1985 N), the irregular strain growth appeared, increased from 2.8×10^-4 up to 1.2×10^-3 till the critical point of failure came. All this procedure reflected the strain was non-linear and in irreversible plastic phrase. As for the ‘User-Defined Model’, from Fig. 7(b), it is clear that plastic strain started when loads equal to 1200 (N) as well, when loads hit the number of 1600 (N), ${\sigma }_f$= 170 KPa, the plastic strain increased from ‘0’ to 0.15×10^-4, which is much smaller than the ‘Coulomb-Mohr’ model, the later plastic strain behavior showed linear form with the consequent loads. Overall, these numerical models can reflect the basic strain state of ice material.

Fig. 7Dependency of loads and flexural stress with strain

a) Equivalent strain curves

b) Plastic strain curves

Fig. 8Dependance of loads with deformations

4.3. Deflection and strain-stress state analysis

As per Fig. 8, under the test loads, ‘Coulomb-Mohr’ model and ‘User Defined Model’ have the same deflection value ${∆}_s=$ 8.78 mm ,approximately equal to experimental value ${∆}_e=$ 9.00 mm. the simulation values fit the experimental data well, however, the theoretical calculation results ${∆}_t$= 5.03 mm have significant disparities with them, tracing back to the applied deflection equation (Eq. (3)), the possible main reason is conducted by ignoring the hydrostatic pressure.

Despite of the initial assumption that ice beam can be considered as homogenous elastic components, according to the depict from Fig. 9, when exert loads functioned on beam, the strain and stress state can reflect partial non-linear mechanic pattern, also indicate the variable of Young’s Modulus. It is easy to calculate the elastic modulus from chart below at the range from 0.59 GPa to 1.5 GPa (‘Coulomb-Mohr’), and 0.24 Gpa to 3.5 GPa (‘User defined’), fit the initial parameters from Table 3.

Fig. 9Strain and stress curves of beam on section A-A