A Study on the Ice Resistance Characteristics of Ships in Rafted Ice Based on the Circumferential Crack Method

Abstract

In previous studies of ship–ice interactions, most studies focused on ship–level ice interactions, overlooking potential rafted ice conditions in extreme ice conditions. The purpose of this study is to develop a numerical model for predicting ship resistance in rafted ice regions. Numerical modeling of rafted ice was carried out using preset grid cells. By comparing the model test results, the accuracy and reliability of the numerical model are verified. On this basis, we undertook the analysis of the impacts of different ice thicknesses, ship speeds, bending strengths, and crushing strengths on the ice resistance of ships under level and rafted ice conditions. The results show that the ice resistance of ships is significantly higher than that of rafted ice under the condition of level ice; however, level ice and rafted ice have different effects on ship ice resistance. Compared with level ice, the ice resistance of ships navigating in rafted ice is more concentrated. The findings of the present research can serve as a technical reference for studies focused on predicting ship resistance in rafted ice regions.

1. Introduction

Typical features of sea ice in the polar regions include brash ice, floating ice, layered ice, and rafted ice. Rafted ice is one of the specific ice formations in the polar regions, especially during the initial and final sea ice periods. The dynamic effects of the fractures, extrusion, and accumulation of sea ice cause an increase in sea ice thickness. During navigation, a ship is subjected to non-linear solid ice resistance, which significantly challenges a ship’s safe navigation.

For the navigation safety design of polar ships, researchers have proposed various ship performance prediction methods under sea ice conditions, which can be divided into experimental [1,2,3,4,5], analytical [6,7,8,9,10], and numerical methods [11,12,13,14,15]. Experimental methods include full-scale measurements and model tests. Full-scale ship trials are challenging to replicate and involve high costs, while model tests impose strict requirements on the experimental equipment and methodology. Empirical formula methods involve theoretical analyses of ship–ice interaction processes but often simplify the ship and sea ice models, which has particular limitations for complex sea ice and ship models. In recent years, with the rapid improvement of computer performance, numerical methods have been effectively applied. Numerical methods have significant development potential compared to experimental and empirical analytical methods. Numerical methods applied to ship–ice interactions mainly include the finite element method (FEM), the discrete element method (DEM), and the circumferential crack method (CCM). The FEM has been widely applied in the context of ship–ice interaction problems over the years.

Yu et al. [16] employed the finite element method to numerically simulate periodic ice loads in the interaction between sea ice and conical structures, and the calculated sea ice bending damage process was similar to that of the results of full-scale measurements. Feng et al. [17] used the cohesive element method to simulate the interaction between ice and structures and conducted with an analysis of parameter sensitivity. It was found that the structural response was very sensitive to changes in the fracture energy, and the stress–strain curve of the body unit had a significant effect on the simulation. This method was also used by Wang et al. [18] to simulate the continuous icebreaking process of ships at different heeling angles, and they analyzed the continuous icebreaking process of different ships with different transverse inclinations, with the results showing that the ice resistance capacity of the ship and the extension length of the sea ice crevasse increased with the increase in the ship’s transverse inclination angle. Lee et al. [19] proposed a method to analyze the ice load in the frequency domain, and the trend of the overall power spectral density with the bow angle was analyzed using different regression methods (linear interpolation, support vector machine, random forest, and deep neural network), and it was found that the deep neural network method performed the best. Shi et al. [20] proposed an elastic-plastic iceberg material model with temperature gradient effect to study the dynamic collision process between a floating production storage and offloading vessel (FPSO) and an iceberg. The simulation results are compared with the design specification to verify the validity of the iceberg model, and the effects of different iceberg shapes and temperatures on the collision process are analyzed. The results show that the structural damage of a floating production storage and offloading vessel (FPSO) is affected by the structural strength, the iceberg strength, and the localized shape of the iceberg.Based on the interaction process between ships and ice as well as the theory of sea ice fracturing, Lu et al. [21] proposed an edge-crack theory model. Using the extended finite element method, the mechanism of long crack propagation between parallel ice-breaking channels was studied. The maximum distance between parallel channels without sea ice fracturing was investigated and validated against experimental results.

In the discrete element method (DEM) realm, Hanse et al. [22] employed a two-dimensional discrete circular-disc viscoelastic model to simulate broken ice and adjusted numerical model calculation parameters according to ice tank experiments. Lau et al. [23] conducted a series of numerical simulations on the interaction between ice offshore structures and ice ships using the three-dimensional block discrete element model. Liu et al. [24] calculated the impact of factors such as ship speed, ice thickness, and ship width on the ice resistance of ships using the DEM. Dong et al. [25] established an ice channel model based on the discrete features of broken ice. Using image segmentation methods to extract ice channel regions and introducing intelligent corner regression networks to accurately delineate ice channel boundaries, this method has shown good accuracy in real ice channel recognition. Xie et al. [26] simulated the ship–water interaction using a coupled CFD-DEM method and established a discretized propeller model (DPM) and a body force model (BFM). The results indicate that the BFM method can be used effectively for the assessment of the main engine power and hull profile optimization during the ship development and design stages. Regarding the circumferential crack method, Zhou et al. [27] proposed a method based on the circumferential crack approach to distinguish the forms of sea ice damage according to the ship’s heel angle, and they compared the numerical simulation results with the model test results, which achieved a good consistency. Moreover, Gu et al. [28] predicted the slewing motion of a polar ship in horizontal ice, considered the effect of hull camber on different damage modes of sea ice, and analyzed the results in comparison with the results of real ruler measurements, and the two results are in good agreement.

Several scholars have also worked on rafted ice material modeling. Hopkins et al. [29] utilized the discrete element approach using circular-disk and block models to validate that the relative motion of two flat ice blocks can result in either overlapping or crushing and breaking. The former leads to rafted ice, and the latter is the initial process of ice ridge formation. After observing the natural appearance of rafted ice through experiments, Leppäranta et al. [30] found that the ice crystals at the contact point between the two layers of level ice in the rafted ice formed a granular structure, and the shear strength of rafted ice was thus lower than that of level ice. Bailey et al. [31] found that the shear force at the adhesive interface of artificially created rafted ice was approximately 30% lower than that of level ice through experiments. Parmerter et al. [32] established a numerical sea ice rafting model capable of calculating the bending stress during the ice rafting process. The results showed that the increase in the bending stress of sea ice is proportional to the square of the ice thickness.

Although these methods have been applied to study ship interactions with level ice, ice floes, broken ice, ice ridges, etc., most of the research on rafted ice has focused on its mechanical properties and physical models. There has been relatively less exploration on ship collisions with rafted ice. This paper combines a preset grid method with the circumferential crack icebreaking assumption to establish a numerical model for rafted ice. The model will be used to predict the resistance of ships in the rafted ice region and compare the numerical simulation results with the model test results. On this premise, the effects of different ice thicknesses, ship speeds, and sea ice characteristics on the level of ship ice resistance and rafted ice are studied. This study supports subsequent ship resistance predictions in rafted ice regions more effectively and holds a specific engineering application value.

2. Numerical Model

2.1. Circumferential Crack Method

When ships navigate in polar regions, the interaction between the ship and the sea ice leads to localized compression and fragmentation of the free edge of the ice when the ship’s bow comes into contact with the ice. With the increase in the contact area between the ship and the sea ice, there is a corresponding increase in crushing force, resulting in circumferential cracks parallel to the contact area or radial cracks perpendicular to the contact area. Based on the physical phenomena observed in full-scale measurements and model tests, the hypothesis of circumferential crack occurrence is adopted in this study. The geometric shape of the fractured floating ice is assumed to be wedge-shaped, with the ice wedge angle denoted as $\theta$. The icebreaking radius of the shape of the ice wedge is R, as expressed in the literature [33].

$$R=C_l\times l\left(1.0+C_v\times v_{n,2}\right)$$

(1)

where $C_l$ and $C_v$ represent empirical parameters, $v_{n,2}$ denotes the relative normal velocity between the ship and the sea ice, $l$ refers to the characteristic length of the ice, which can be expressed as follows:

$$l={\left[\frac{E_i{h}_i^3}{12\left(1-v^2\right){\rho }_{w}g}\right]}^{1/4}$$

(2)

where $E_i$ signifies the elastic modulus, $h_i$ represents the ice thickness, $v$ represents the poisson ratio, ${\rho }_w$ signifies the density of water, and $g$ denotes the acceleration in gravity.

This paper converts the fan-shaped ice wedge to a square through area-equivalent treatment, as shown in Figure 1. Assuming that the icebreaking radius of the ice wedge is equal to that of the side length of the square grid cell and that the areas are equal, the icebreaking angle $\theta$ of the two satisfy the following:

$$R^2=\frac{\theta }{2}{R}^2\to \theta =2rad$$

(3)

2.2. Icebreaking Force

During polar ship icebreaking navigation, the compressive force gradually increases as the contact area between sea ice and the ship’s hull increases. Before bending failure, the icebreaking force $F_{cr}$ generated by the compression between the ship’s hull and sea ice, with the force being perpendicular to the contact area, can be expressed as follows [34]:

$$F_{cr}={\sigma }_c\cdot A_c$$

(4)

where $A_c$ is the contact area and ${\sigma }_c$ is the crushing strength of the ice.

2.3. Contact Area

When $L_d\cdot \tan \left(\phi \right)\le h_i$, the contact area between the ship and the sea ice has not reached the bottom of the sea ice at this point, resulting in a triangular contact area as follows:

$$A_c=\frac{1}{2}{L}_h\frac{L_d}{\cos \left(\phi \right)}$$

(5)

When $L_d\cdot \tan \left(\phi \right)\ge h_i$, the contact area between the ship and the sea ice reaches the bottom, resulting in a quadrilateral contact as follows:

$$A_c=\frac{1}{2}\left(L_h+L_h\frac{L_d-h_i/\tan \left(\phi \right)}{L_d}\right)\frac{h_i}{\sin \left(\phi \right)}$$

(6)

where $A_c$ represents the contact area, $L_h$ denotes the contact length, $L_d$ refers to the contact length, $h_i$ represents the ice thickness, and $\phi$ signifies the outward tilt angle at different ship nodes.

2.4. Ice Failure Model

During the interaction between ships and sea ice, the failure mode of sea ice is influenced by various factors, including ship angle, ice thickness, and the relative velocity between the ice and the ship. The ice failure model includes both bending failure and crushing failure in the present study. According to the research by Zhou et al. [27], different sea ice damage modes were used to distinguish the relationship between bending and crushing damage, and they found that the hull camber angle produces different sea ice damage modes and that the ship–ice friction coefficient affects the ultimate hull camber angle of the sea ice failure mode. On this basis, Gu et al. [28] assumed that the friction coefficient between the ship and the sea ice was 0.1 and calculated that the limiting angle of the ship–ice failure mode was 84.2894°, which means that, during the icebreaking voyage, if the angle between the ship and the ice is more than 84.2894°, this will lead to crushing damage, whereas if it is less than 84.2894°, this will lead to bending damage.

According to Kerr [35], the expression for the ultimate load of ice bending failure is given as follows:

$$P_f=C_f{\left(\frac{\theta }{\pi }\right)}^2{\sigma }_f{h}_i^2$$

(7)

where $C_f$ signifies the empirical parameter, $\theta$ signifies the idealized ice fracturing angle, ${\sigma }_f$ denotes the bending strength of ice, and $h_i$ represents the ice thickness.

When ice experiences crushing failure, the localized icebreaking force acting on the hull, according to ISO/FDIS 19906-2019 [36], can be expressed as follows:

$$F_{cr}=P_G·A_c$$

(8)

$$P_G=C_R\left[{\left(\frac{h}{h_1}\right)}^n{\left(\frac{L}{h}\right)}^m+f_{AR}\right]$$

(9)

where $C_R$ represents the ice strength coefficient, $h$ signifies the ice thickness, $h_1$ signifies the reference ice thickness of 1 m, $m$ and $n$ are the empirical coefficients, and $f_{AR}$ represents the ice strength coefficient.

2.5. Rafted Ice Model

Currently, there are two main types of rafted ice. One type is that of finger-rafted ice, where the ice body does not move as a rigid body, accumulating internal stresses that are then released a shear force, thereby forming finger-rafted ice. Another type is that of layered rafted ice, where one ice layer fractures and climbs onto another under external dynamic forces, becoming rafted and forming a layered structure [37]. Figure 2 illustrates the phenomenon of vertical layering in the material thickness direction of consolidated ice, including the smooth level ice layer, the consolidation layer, and the submerged layer. The submerged layer is formed by immersion in water, interacts with the level of ice under the action of buoyancy, and ultimately forms the intermediate consolidation layer. Shafrova et al. [38] conducted experiments on the freezing process of first-year ice ridges and noted that various factors, such as seawater infiltration, temperature, salinity, and pressure, affect the strength of the frozen bond between ice bodies during the formation process. While studying the material properties of consolidated ice, Chen et al. [39] obtained a fragment function relationship between its compressive strength and the strain rate. The calculation formula is as follows:

$${\sigma }_c^{\prime }=\left\{\begin{array}{l}0.37{\dot{\epsilon }}^{0.2} \dot{\epsilon }>4.6\times {10}^{-4}\\ 53 {\dot{\epsilon }}^{0.2} \dot{\epsilon }\le 4.6\times {10}^{-4}\end{array}\right\}$$

(10)

where ${\sigma }_c^{\prime }$ is the revised compressive strength.

This study focuses on layered rafted ice, which belongs to composite ice formations. The stacking of two level ice layers primarily develops it. Currently, numerous scholars have proposed corresponding numerical models based on the characteristics exhibited by ice. Ni et al. [40] introduce cohesive elements to numerically model the intra- and interlayer structures of the rafted ice layers, respectively. By randomly deleting the cohesive elements within the model, the porosity of the natural rafted ice was successfully simulated, and numerical calculations of the collision between a ship and rafted ice were carried out, and it was found that the method was excellent in simulating the crack extension of rafted ice. In this study, based on the significant vertical layered structure of rafted ice, the solidified and submerged layers are collectively referred to as the second layer of rafted ice in the numerical model. And a correction factor is introduced to define the constitutive parameters of the ice layer, which can include the mechanical properties of the low sea ice in the condensation layer, reasonably expressing the differences between the rafted ice layers, with the formula provided below.

$$\sigma =\left\{\begin{array}{l}C{\sigma }_1 n=1 \\ C{\sigma }_2 n=2\end{array}\right\}$$

(11)

where $\sigma$ is the ice strength, $C$ is the correction factor, and $n$ is the number of layers of rafted ice.

In the process of ship–ice interactions, the action of icebreaking forces causes the formation of ice cracks. As these cracks spread, the ice gradually breaks and destroys. The formation of rafted ice crevasses is simulated using grid cells, and the sea ice failure model is introduced. The rafted ice is separated into isolated grid cells. The side length of the grid cell is related to the icebreaking radius $R$, as shown in Figure 3. In the numerical model, when the icebreaking force reaches the load-bearing limit of the grid cells, the rafted ice is destroyed. Introducing this failure model allows for a more detailed consideration of the rafted ice layer’s destruction process and simulates the ice layer’s fracture behavior in numerical simulations. Each grid cell represents a discrete unit of the ice layer, and by monitoring the impact of icebreaking forces on these units, it is possible to track the real-time generation and propagation of cracks, ultimately simulating the complex failure process of the rafted ice under the action of icebreaking forces. Figure 4 illustrates the computational flow of the ship icebreaking simulation, which mainly includes the numerical model and numerical process.

3. Model Test in Rafted Ice

3.1. Experimental Description

A relevant model test of a ship model sailing in rafted ice was performed in an outdoor ice tank of Harbin Engineering University [41]. The ice water tank is 20 m in length, 2 m in width, and 1.5 m deep, and it can naturally make different ice features in winter. The purpose of this is to conduct a towing test on a ship operating in rafted ice, and the test follows Froude similarity and Cauchy similarity using a polar ship with a model scale of 1:60. The main particulars of the full-scale and model-scale ships are listed in

Table 1. Main particulars of the ship.

Principal Hull Data	Full-Scale	Model-Scale
Length between perpendiculars/m	122.5	2.04
Beam/m	23.32	0.38
Draught/m	7.8	0.13
Stem angle/°	20	20
Waterline angle/°	34	34

.

Xu et al. [41] selected three different speeds, 0.17, 0.27, and 0.37 m/s, for the ship model towing tests in the rafted ice region, and two tests were conducted for each speed, and the settings of these six ship model tests are listed in

Table 2. Parameters of rafted ice.

Case	Towing Speed/m/s	Bending Strength/MPa	Crushing Strength/MPa
1	0.17	0.85	1.73
2		0.91	1.40
3	0.27	0.84	1.76
4		0.97	1.86
5	0.37	0.72	1.14
6		0.74	1.06

.

In the numerical simulation, as shown in Figure 5, the rafted ice is divided into upper and lower layers, with each being composed of numerous square ice grids. Each ice grid’s length can be taken as the ship’s icebreaking radius during navigation. The model is divided into upper and lower layers that accurately reproduce the construction of rafted ice, with each layer consisting of several square grids. Fine discretization of the waterline at the ship’s draft ensures a reasonable simulation of the rafted ice damage and a precise description of the location of the ship–ice grid contact points. Figure 6 provides the initial top view of the ice–ship interaction.

3.2. Comparison of Model Tests

In the numerical simulation, each grid cell size represents the icebreaking radius. White grids indicate no contact between the ship and these grid cells. Red grids indicate interaction between the ship and the grid cells, while blue grids indicate grid cells that have failed. Figure 7 illustrates the interaction process between the ship and the rafted ice in the numerical simulation. With the continuous progress of the ship, the contact area between the hull and the grid cells gradually increases. When the icebreaking force exceeds the load-bearing limit of the grid cells, the grid cells fail, indicating the occurrence of fractures and a failure in the rafted ice.

In this paper, the test conditions of Case 5 are selected to analyze the numerically simulated time history curves of total resistance, as shown in Figure 8. The simulation results show obvious periodic characteristics with relatively stable peaks. This is due to the fact that the grid cell parameters of the rafted ice are fixed in the numerical simulation, and the ice resistance value oscillates and changes within a certain amplitude when the ship reaches the icebreaking stabilization stage. The comparison between the numerical simulation results and the experimental results is presented in

Table 3. Comparison between model test results and numerical results.

Case	Average in Experiment/N	Average in Simulation/N	Error/%
1	25.45	24.69	2.9
2	28.13	30.91	9.8
3	34.61	31.83	8.0
4	29.14	27.39	5.9
5	42.52	43.10	1.3
6	35.78	37.80	5.6

. It can be observed that with increasing ship speed, the resistance in the rafted ice also increases. By comparing six different experimental conditions, the error between the two ice resistance values is within 10%, indicating a good consistency and verifying the accuracy of the numerical model to some extent.

4. Sensitivity Analysis of Parameters

The numerical model simulating the ship’s motion in rafted ice established in this paper has been well validated through comparisons with previous experiments and numerical results. Additionally, this study focuses on the ice resistance characteristics of ships in the level and rafted ice regions. The model test study from the ice tank of Aalto University in Finland was chosen [42]. The model ship of MT Uikku had a scale of 1:31.6. The key parameters of both the full-scale and model-scale platforms are presented in

Table 4. Main particulars of MT Uikku.

Principal Hull Data	Full-Scale	Model-Scale
Length between perpendiculars/m	150	4.75
Beam/m	21.3	0.67
Draught/m	9.5	0.3
Stem angle/°	30	30
Waterline angle/°	21	21

. The ship model was towed by a trailer at a constant speed on level ice, and various experimental data were obtained by changing the ice thickness and velocity. Three test conditions were selected from this model test for comparisons, which were Case1, 2, and 3, with the parameters of the level ice being listed in

Table 5. Settings of level ice in model tests.

Case	Towing Speed/m/s	Towing Speed/kn	Ice Thickness/m	Bending Strength/MPa	Crushing Strength/MPa
1	0.09	0.97	0.76	0.844	2.192
2	0.09	0.97	1.03	0.669	2.485
3	0.09	0.97	0.63	1.029	5.389

.

Sequential ice resistance results through Case1, 2, and 3 were analyzed, as shown in Figure 9. The squares and triangles are the average and maximum ice resistance observed for Case1, 2, and 3, respectively, and the dots and pentagrams are the average and maximum ice resistance for Case1, 2, and 3 in the numerical simulation. In Case1, the numerical simulation’s average value is 598.59 kN, which closely matches the measurement of 560 kN from the ice tank experiment, with an error of approximately 6.8%. In Case2, the numerical simulation’s average value is 759.9 kN, while the measured value in the ice tank experiment is 830 kN, with an error of approximately 8%. In Case3, the average error between the two is 10%. By comparing the maximum ice resistance values under the three test conditions, it can be observed that only in Case2 is there a significant error between the numerical simulation and the ice tank test in the peak ice resistance. The main reason for this error is the inherent variability in ice parameters in different regions of level ice during the model test preparation process. However, the numerical simulation discretizes the ice field using grids coupled with predetermined ice parameters, which, to some extent, affects the peak ice resistance. Based on the above analysis, the numerical results are qualitatively and quantitatively consistent with the experimental data. The numerical method can also predict ice resistance for ships in level ice.

Numerical simulations were conducted under Case1, 2, and 3 to compare the average ice resistance of level and rafted ice of the same ice thicknesses for the three conditions, as shown in Figure 10. The ice resistance of rafted ice is significantly lower than that of level ice at the same speed, and the difference in ice resistance gradually increases as the ice thickness increases from 0.63 m to 1.03 m. Since the structure of rafted ice is composed of two thin layers of level ice that undergo secondary freezing, its overall strength is lower than that of level ice; therefore, rafted ice is more prone to damage during ship–ice interactions.

Due to the different formation mechanisms and internal structures of level ice and rafted ice, there are specific differences in their resistance characteristics. As shown in Figure 11, in Case1, for example, the time history curves of ice resistance for both the level ice and the rafted ice were analyzed. It can be observed that within the first 25 s, the two trends are similar, but the resistance is lower than that of the level ice during the same period. As the icebreaker keeps moving forward, both show an increasing trend in resistance. The ice resistance exhibits periodic fluctuations once the ship enters the stable icebreaking stage. The ice resistance fluctuates with more prominent peaks when interacting with level ice. However, due to the differences in the mechanical properties of the rafted ice layers in the numerical simulation, the ultimate loads on the ice grids are inconsistent. In Figure 12a, the first rafted ice layer where the ship’s bow is first contacted by the two grid cells has already failed, while in Figure 12b, the two grid cells of the second rafted ice layer at the same position are in action, and the forces generated by the different layers of ice cause the ship’s ice resistance to fluctuate more significantly. The peaks of the fluctuations are more significant.

4.1. Influence of Ice Thickness

The ice thickness is a crucial component influencing crushed ice. Diverse ice thicknesses of 0.6 m, 0.8 m, 1.0 m, and 1.2 m were selected to simulate the icebreaking of ships in the level and rafted ice regions, and the sailing speed was 1 kn, and the sea ice parameters were referred to in Table 5 and

Table 6. The mechanical parameters for numerical simulation of rafted ice.

Case	Velocity/kn	Bending Strength/MPa	Crushing Strength/MPa	Number of Layers	Correction Factor
1	0.97	0.844	2.192	1	1.0
		0.759	1.972	2	0.9
2	0.97	0.669	2.485	1	1.0
		0.602	2.236	2	0.9
3	0.97	1.029	5.389	1	1.0
		0.926	4.850	2	0.9

. Figure 13 compares the mean ice resistance of ships in level and rafted ice regions under varying ice thicknesses. As the ice thickness rises from 0.6 m to 1.0 m, the ship ice resistance of level ice increases from 558.02 kN to 1386.8 kN, while the ship ice resistance of the rafted ice increases from 335.20 kN to 819.14 kN. The ice resistance of ships in both the level and rafted ice regions increases with the growing ice thickness. Under the same ice thickness circumstances, the resistance of ships in the rafted ice region is relatively close to that of the level ice region at 0.6 m ice thickness, but at 1.2 m ice thickness, there is a significant difference between the two. It can be found that the ship ice resistance in level ice is more sensitive to the change in ice thickness compared to the ship ice resistance in the rafted ice area.

Figure 14 illustrates the distribution of ship ice resistance in level and rafted ice conditions under four different ice thicknesses at a certain ship speed. It can be observed that the center of distribution of ship ice resistance in both level and rafted ice increases with the growth of ice thickness. In level ice, the distribution of ship ice resistance gradually transitions from a left-skewed distribution to a right-skewed distribution, indicating that increasing ice thickness increases the probability of encountering peak values in ship ice resistance. In addition, the distribution of ship ice resistance for rafted ice increases with the increase in ice thickness, especially in the conditions of 1 m and 2 m ice thickness, but the distribution of ship ice resistance is lower compared with that of level ice under the same ice thickness conditions. Notably, under different thickness conditions, the distribution of ship ice resistance in rafted ice is more concentrated compared to that of level ice. This suggests that the internal structure of level and rafted ice has distinct influences on the distribution of ice resistance. The structure of rafted ice, being more intricate, results in a more concentrated distribution of ice resistance, potentially leading to increased fatigue effects on the ship’s structure.

4.2. Influence of Ship Speed

This article selects speeds of 2 kn, 3 kn, 4 kn, and 5 kn to simulate ship icebreaking in the level ice and rafted ice regions. The sea ice parameters are shown in Table 5 and Table 6, with an ice thickness of 0.76 m. Figure 15 illustrates the variation in ship ice resistance trends for level and rafted ice at different ship speeds. With an almost linearly increasing relationship, ship speed highly influences the ship’s ice resistance in different ice conditions. A comparison between ship speeds of 2 kn and 5 kn shows that the ship ice resistance of level ice rises by 39.2%, while that of rafted ice increases by 38.4%. It can be observed that the resistance of level ice is more sensitive than that of rafted ice under the influence of ship speed.

The results in Figure 16 show that the median line of ship ice resistance increases with speed in all cases except for that of the case of a speed of 4 kn in level ice. The variation in ship speed directly influences the size of the grid cells. Figure 17 illustrates the icebreaking state of vessels in level ice simultaneously under four different ship speeds. It is observed that, under the condition of a 4 kn speed in level ice, continuous crushing occurs between the ship’s side and the grid cells. This leads to a more drastic variation in ship ice resistance, significantly increasing the probability of peak values and causing a more dispersed overall distribution of ice resistance. Therefore, speed not only impacts the peak magnitude of ice resistance but also significantly influences the distribution of ship ice resistance. This indicates that the change in speed may induce alterations in the interaction state between the vessel and ice, consequently affecting the overall characteristics of ice resistance.

4.3. Influence of Bending Strength

In the numerical simulations, bending strength is a crucial parameter directly influencing the load-bearing capacity of each ice grid. To analyze the impact of bending strength on ship ice resistance, especially considering the environmental differences in the growth of rafted ice, which are primarily reflected in the parameter variations of the lower sea ice layer, this paper selected the correction factors of 0.6, 0.7, 0.8, and 0.9 to set the bending strength for both the overall level ice and the lower layer of the rafted ice. Figure 18 illustrates the linearly increasing trend of ship ice resistance for the level and rafted ice as the correction factor for the bending strength rise from 0.6 to 0.9. The ship ice resistance for level ice increases from 539.15 kN to 701.42 kN, while for rafted ice, it increases from 373.14 kN to 431.87 kN. It is worth noting that the ice resistance generated by the ship in level ice remains higher than that in rafted ice under the influence of bending strength.

As shown in Figure 19, it can be observed that the median line of the ship ice resistance in level and rafted ice increases with the increasing bending strength, which implies that a higher peak in the ship ice resistance has an effect on the central tendency. Bending strength directly influences the load-bearing capacity of the ice layer, and with an increase in bending strength, the magnitude of ice resistance experienced by the vessel markedly rises. From the analysis of the overall distribution of ship ice resistance, it can be seen that the larger distance between the minima of the median line makes the ship ice resistance in level ice have a right-skewed distribution, while in the rafted ice, the median line and the minima are equally distant from each other, thereby showing a normal distribution. Despite an increase in bending strength, the overall trend changes relatively insignificantly. This suggests that bending strength does not significantly impact the distribution of ice resistance.

4.4. Influence of Crushing Strength

Crushing strength is a crucial parameter in the interaction between ships and ice. In numerical simulations, correction factors of 0.6, 0.7, 0.8, and 0.9 were selected to set the compression strength for level ice and the lower layer of rafted ice. The sailing speed is 0.97 kn, and the ice thickness is 0.76 m. Figure 20 shows that with the increase of crushing strength, the average and maximum ship ice resistance tend to decline. The icebreaking force between the ship and the ice is affected by the crushing strength of the sea ice. Analysis of individual interactions between the ship and ice grids reveals that enhancing crushing strength shortens the time required for the breaking force to reach the ice’s load-bearing limit. The force remains zero until colliding with the next ice grid, leading to a decreasing trend in the mean ice resistance.

As shown in Figure 21, it can be observed that the median line of ship ice resistance decreases gradually with the increase in breaking strength in both the level and rafted ice. In the level ice, the median line of ship ice resistance is further away from the end of the minima, which makes the overall left-skewed distribution. In the rafted ice, the median line of ship ice resistance is further away from the extreme value end, and the overall distribution is right-skewed. This is due to the crushing strength directly affecting the icebreaking force between the ship and the ice. The icebreaking force between the ship and the ice increases with the increase in crushing strength. This means that the time of grid cell failure is accelerated under the condition of the same carrying capacity of the grid cells in the level ice. In the process of making contact with the next grid cell, no ice resistance is generated, which leads to a decrease in the trend of the ice resistance of the ship in the level ice. However, the process of collision between the ship and the rafted ice results in the failure of the upper grid cell, which does not mean that the lower grid cell will also fail due to differences in the mechanical properties of the layers of rafted ice. This alternating action leads to a significant difference in the distribution of the overall ice resistance from that in the level ice.

5. Conclusions

Based on the assumption of circumferential crack icebreaking, this paper employed a predefined grid method to numerically simulate the icebreaking navigation of ships in the rafted ice region. The numerical simulation results were compared with experimental data, demonstrating a good consistency. On this basis, the paper further investigates the influence of some key parameters on ship ice resistance and analyzes the resistance distribution characteristics of the level ice and rafted ice areas using probability density functions. The following conclusions are made:

• According to the structural characteristics and mechanical properties of rafted ice, this paper adopts a new numerical method to establish a numerical model and to simulate the icebreaking process of ships in the rafted ice area, which are the keys to success.
• Moreover, this paper utilizes the established numerical model for ship icebreaking in rafted ice areas to conduct numerical simulations and validate it against six operational conditions in an ice tank model experiment. The results demonstrate that this method can accurately predict the ice resistance experienced by ships in rafted ice areas, with the error between the simulated resistance value and the experimental value being within 10%.
• Therefore, the model accurately predicts ship ice resistance in level ice regions through numerical case analysis. It can compare the effects of ice thickness, ship speed, bending strength, and crushing strength on ship ice resistance in both the level and rafted ice areas. Simulation results indicate that the ship ice resistance in the level and rafted ice regions linearly increases with ice thickness, ship speed, and bending strength while linearly decreasing with crushing force. Comparing ship ice resistance in the level and rafted ice regions, the influences of ice thickness, ship speed, bending strength, and crushing strength on ship ice resistance are more sensitive in level ice than in rafted ice.
• Numerical simulations of ships operating in level and rafted ice show that the ice resistance generated in level ice is more significant than that in rafted ice. However, this does not imply that the potential damage to the vessel caused by rafted ice can be easily overlooked. In reality, the ice resistance from the interaction between the ship and the rafted ice is more concentrated than that in the level ice. This concentration can make the ship’s structure more susceptible to fatigue, increasing the risks associated with polar navigation.

The present method should be further verified with more measured data in the future. Moreover, variations in ship speed can impact the size of grid cells in the ice field, consequently influencing the distribution of ice resistance experienced by the ship, which could be studied further.