Numerical Simulation on the Resistance Performance of Ice-Going Container Ship Under Brash Ice Conditions

Ice resistance prediction is a critical issue in the preliminary design of ships navigating brash ice conditions, which is closely related to the safety of a ship to navigate encounter brash ice, and has significant effects on the kinds of propellers and motor power needed. In research on this topic, model tests and full-scale tests on ships have thus far been the primary approaches. In recent years, the application of the finite element method (FEM) has also attracted interest. Some researchers have conducted numerical simulations on ship–ice interactions using the fluid–structure interaction (FSI) method. This study used this method to predict and analyze the resistance of an ice-going ship, and compared the results with those of model ship tests conducted in a towing tank with synthetic ice to discuss the feasibility of the FEM. A numerical simulation and experimental methods were used to predict the brash ice resistance of an ice-going container ship model in a condition with three concentrations of brash ice (60%, 80%, and 90%). A comparison of the results yielded satisfactory agreement between the numerical simulation and the experiments in terms of both observed phenomena and resistance values, indicating that the proposed numerical simulation has significant potential for use in related studies in the future.


Introduction
As the basic equipment required for polar shipping, polar expedition, and polar resource exploitation, ice-going ships are associated with several fundamental issues in the area, such as ice resistance, ice load forecasting, and studies on the structural and overall performance of such vessels, that directly affect the manufacture and research capabilities of ice-going ships (Zhu et al., 2015). Of these issues, predicting the resistance performance of ice-going ships is a major subject in preliminary design, as it is related to the navigation performance of ships in ice and has significant effects in determining the kinds of propellers and motor power needed.
Since the first model test of icebreaker resistance in 1964, empirical formulae and model ship tests have gradually become the primary methodology of studies on the resistance performance of ice-going ships (Guo et al., 2016). Empirical formulae are mathematical formulae based on test data from model or full-scale ships. For most empirical formulae employed to predict the sailing resistance of ships in ice, resistance in ice has been divided into separate parts based on ship-ice interaction (Lindqvist, 1989;Riska, 1997;Spencer and Jones, 2001;Jeong and Choi, 2008). For example, Spencer divided the ice-breaking resistance of a ship into four parts: open-water resistance, ice buoyancy resistance, ice-clearing resistance, and ice-breaking resistance. These empirical formulae have different forms and influencing factors; for most of them, the primary factors influencing the resistance of a ship to ice are its main dimensions, the speed of the ship, and the thickness, density, and the bending strength of ice. Colbourne believed that the resistance to brash ice is related to the width of the ship as well as the density, thickness, speed, and distribution concentration of brash ice. In Spencer's empirical formula, the factors influencing the resistance are ice density, water density, ship width, ice thickness, ship speed, and water speed. In gener-al, an empirical formula can only be applied to a certain type or several similar types of ships. Model ship tests are the most important means of investigating the resistance of a ship to ice (Kim et al., 2006;Molyneux and Kim, 2007;Jeong and Choi, 2008;Jeong et al., 2010), as they allow for direct and effective observation of the ship as it navigates icy waters, and the measurement of the resistance performance of various ship types, ice conditions, and ship speeds. This yields sufficient data pertaining to the performance of the ship in ice. This method also has non-substitutable applications in the fields of hull optimization and novel ship design. Since the tests in model ice basins are restricted by such conditions as the field used, cost, and technical means, some researchers have conducted model ship tests in normal-temperature basins using synthetic model ice to study the resistance performance of ships in ice. Although such an approach cannot simulate all real ice-breaking scenarios, it remains useful for research in the field because of its high practical feasibility Lee et al., 2009;Song et al., 2007).
In recent years, with increasing computational power of computers and widespread application of numerical simulations to engineering, the fluid-structure interaction (FSI) algorithm based on the finite element method (FEM). Li et al. (2018) selected methods for the ship performance in level ice, ridged ice and channel ice are evaluated based on fullscale measurement data of two ships, the modelling of subprocesses and ice properties led to certain scatter in applying these methods for speed prediction. Liu et al. (2018) introduced the use of peridynamics to calculate ice loads and simulate the ship-ice interaction process without giving the gravity and buoyancy of ice, which obtained dynamic crack generation and propagation of ice rubble, and sliding, rotation and accumulation of broken ice. Aquelet et al. (2003) has been successfully applied to the analysis of the mechanisms of structures in ice and the ocean. Past applications include the sloshing of a liquefied natural gas (LNG) tank in an LNG carrier, the destructive effect of level ice on offshore platforms (Wang and Derradji-Aouat, 2011), and the collision interaction between ice-going ships and small icebergs (Gagnon and Wang, 2012). Wang et al. (2018) used the LS-DYNA finite element code to simulate the interaction between sloping marine structure and level ice, but the simulations only modeled the initial part of ice-structure interaction. Studies on ship-ice interaction using the FEM began much later. Wang and Derradji-Aouat (2010), using the FEM, studied the resistance performance of a model of the standard Canadian icebreaker Terry Fox in brash ice and correctly simulated the buoyancy effects of brash ice in water. Resistance in brash ice was measured under ship speeds of 0.4 m/s and 0.6 m/s, respectively, where the results at the former speed were closer to experimental data. The potential factors leading to this, according to them, were the rigid body mechanics, geometry of ice, boundary effects, and mesh dependency. Lee et al. (2013) calculated brash ice resistance of two ice cargo vessels with different stem styles in a pre-cut ice test. In comparison with the results of model tests, they discovered that the movement of irregular thin ice in the numerical simulation was rather unstable, and from the numerical simulations for both vessels were slightly higher in values than the test values. Kim et al. (2013) performed comprehensive comparative analysis on numerically simulated resistance values, and refrigerated and unrefrigerated ice test data at three waterline angles (30°, 40°, and 50°), three concentrations of brash ice (60%, 80%, and 90%), and three speeds (0.1 m/s, 0.3 m/s, and 0.6 m/s). The results showed that the numerical simulation and test results were in good agreement, both qualitatively and quantitatively, thus showing that numerical simulation can be applied to studies to predict the brash ice resistance of ice-going ships and those on hull form optimization.
The FE approach is still in an exploratory phase in the context of predicting the resistance of ships to brash ice. Some researchers have used both experiments and numerical simulations for comparative studies and analyses. However, most of these studies only considered a limited number of conditions in simulations, and the analyses were often rather simple. Kim et al. (2014) conducted simulations under three concentrations of brash ice and three speeds, which had already been considered through extensive numerical work in similar studies; their analysis of the simulation results was also fairly simple. A large number of researchers have employed a uniform shape and size of the models of ice in their tests and numerical simulations, without considering the change in the dimensions of brash ice (Wang and Derradji-Aouat, 2010;Kim et al., 2013;Kim et al., 2014). Jeong et al. (2017) described the model test and analysis methodology for a brash ice condition in a square-type ice model basin, and the ice operating condition by the installed engine was verified to evaluate the ship's performance in a brash ice channel. At present, numerical simulations cannot be used to study the resistance performance of ships in ice independently of experiments and empirical formulae. However, compared with empirical formulae, numerical methods can obtain not only the average resistance of brash ice, but also the details of the ship-ice interaction, the instantaneous force curve, and the local load on the ship in ice. Moreover, it is not confined to a particular area of the sea or ship type, and is useful in ship optimization and the design of novel ships. Compared with the experimental method, numerical simulations are more efficient in terms of cost and speed of implementation, and are free of field or staff limitations. They offer the ideal complement to experiments. As more numerical cases are accumulated and the relevant algorithms are optimized, numerical simulations will become ever more accurate. Hence, numerical methods have good prospects in analyzing the resistance performance of ships in ice.
By considering the requirements of ships in order to navigate in ice as well as brash ice size, brash ice concentration, and ship speed, this study simulates brash ice conditions and comprehensively analyzes the outcomes, with a focus on the average resistance and the instantaneous load. By investigating the resistance performance of ships in ice under various conditions, this work has more general and practical significance than the above specific aim. An experiment involving synthetic model ice is conducted in a towing tank at the Harbin Engineering University to measure the resistance of a model container ship navigating under brash ice condition. The test matrix consisted of three concentrations of brash ice (60%, 80%, and 90%), seven synthetic model ices with different sizes (5 cm×5 cm, 10 cm× 10 cm, 15 cm×15 cm, 20 cm×20 cm, 25 cm×25 cm, 30 cm× 30 cm, and 35 cm×35 cm), and six ship speeds (0.2, 0.4, 0.6, 0.8, 1.0, and 1.2 m/s). The numerical simulations are performed under the same conditions. The experimental phenomena, average resistance, and instantaneous load curve are then analyzed based on both the experiment and simulation, and the reliability and practical usefulness of the numerical simulation is confirmed from this comparison.

Numerical method
When simulating ships navigating in ice, two effects should be considered simultaneously: the ship-ice interaction, and fluid-structure interaction (FSI). The former can be effectively analyzed via the FE software ANSYS, and is a dominant factor. FSI is also essential, as the simulation of ship-ice interaction and the calculation of the resistance are adversely affected if the initial state and the action of ice in water are not properly simulated. Therefore, the arbitrary Lagrange-Euler algorithm (ALE) in the software LS-DYNA was used to analyze the FSI process, and correctly simulated the mechanical properties of ice in water.
The ALE algorithm integrates the advantages of Lagrange's and Euler's algorithms. Specifically, to process the motion of structural boundaries, it can effectively track the motion of the structural boundaries of substances; the internal meshing is independent of the physical substance, thus preventing grids from severe distortion. The use of the ALE algorithm requires no special modeling of fluid domains. During meshing, hexahedral grids are normally generated for fluid domains such as water and air. In the calculation files, these hexahedral grids are defined as ALE grids. During calculation, the ALE algorithm first executes one or more Lagrange time step operations when the grids deform with material flow. It then executes ALE time step operations to: (1) maintain the boundary conditions of the material following deformation and restructure the internal units of the grids without changing the topological relations among grids (called the smooth step); and (2) transport element variables (density, energy, stress tensor, etc.) from deformed grids and the nodal velocity vector into the restruc-tured grids (called the advection step). The completion of these two operations constitutes that of one step of ALE calculation.

Equation of continuity
The equation of continuity in fluid mechanics is In material, this can be expressed in the following form: Applying the Jacobin matrix yields:

Equation of momentum
The momentum balance on the system boundary is where ρ, t, and b denote density, surface tension, and body force per unit mass, respectively.

Equation of energy
The total energy, E, is the sum of kinetic energy and internal energy, which can be written in the following form: Thus, the energy balance equation is V p q where , , and denote the relative volume, hydrostatic pressure, and bulk viscosity, respectively. According to the finite element theory, the energy balance can be rewritten as: The central difference method is used for the time stepping: ∆t where is the time step size.
In ALE calculations, the material domain and reference domain are defined. The speed of the material, the speed of the reference system, and the speed difference are denoted by , , and , respectively. Therefore, the ALE governing equation is written as: J where is the Jacobian matrix, which represents the relat-ive volume derivative of the two domains.
The material time derivatives in ALE formulation are written as: b r b where means that is expressed as a function of the reference domain.
Hence, the mass, momentum, and energy conversation equations of the ALE method can be expressed as: 2.4 Penalty method The penalty method first requires defining the two contacting surfaces as the master surface and slave surface. When penetration occurs, a series of normal springs are set between the master and slave surface. The spring is intended to confine the penetration. The force of the spring is called the interface force. The interface force of the penalty method can be calculated using the following formula: where k is the stiffness factor and δ, the penetration distance.
The stiffness factor k can be determined by the formulation where is the scale factor for the interface stiffness, is the bulk modulus of the contact unit, is the contact area, and V is the volume of the master segment.

Experiment setup
The experiments are conducted in a towing tank at Harbin Engineering University, where the model of an ice-going container ship is used as a test target and paraffin is selected as the material for synthetic model ice. The experiments are aimed at measuring the resistance performance of the ice-going container ship in brash ice, exploring the methodology of testing the resistance of the model ship in synthetic model ice, accumulating data for subsequent serial tests, and establishing the database for numerical simulations to analyze the resistance performance of the ship in ice. This is expected to provide a comprehensive approach to study the resistance performance of ships in brash ice.
Paraffin is used as the material for synthetic model ice, as mentioned above. The ice-breaking effect in actual navigation is ignored in our test. The sea ice simulated in the ex-periment is new ice or one-year ice with a thickness of approximately 1.05 m. According to the similarity criterion for the ice resistance model test, the thickness of synthetic model ice is 2 cm. There are seven sizes of brash ice, with the concentrations of 60%, 80%, and 90%. Six speed points are set in the experiment for the model ship and a total of 18 navigation scenarios are considered. An open-water test is first performed at each speed point, followed by repeated resistance tests under the brash ice condition. The resistance gap between the two types of tests is defined as brash ice resistance at the given speed point.

Brash ice conditions
The towing tank at Harbin Engineering University is 108 m long, 7 m wide, and 3.5 m deep. It costs extra labor and funds to lay synthetic model ice in the entire tank. Thus, with the aim of reducing the effects of the walls of the tank, a 3-m wide and 20-m long fence consisting of buoys is setup in the tank as a brash ice condition for the test. The model ship navigation is therefore confined to this condition, as shown in Fig. 1.

Model ship
By using an ice-going transport ship as the test model, the experiment focuses only on the resistance performance of the ship in ice, and ignores the effects of ice breaking. By considering commercial transport requirements of an ice-going ship, whereby the ship should demonstrate good resistance performance in both ice and open water, a standard container ship with a bulbous bow is chosen as the test ship  GUO Chun-yu et al. China Ocean Eng., 2018, Vol. 32, No. 5, P. 546-556 model, which, in the experiment, is made of fiberglass.

Properties of ice
The synthetic ice employed in the experiment is composed of semi-refined paraffin wax. A certain mass of granular paraffin wax materials is melted and poured into the mold before being divided into thin square blocks after cooling and solidifying. The model ice has an average thickness of approximately 2 cm and an approximate relative density of 0.9001 g/cm 3 , close to the density of real ice 0.917 g/cm 3 . The friction coefficient of the model ice is 0.63, close to that of real ice 0.3-0.45. There are seven sizes of the 2-cm-thick synthetic ice: 5 cm×5 cm, 10 cm×10 cm, 15 cm×15 cm, 20 cm×20 cm, 25 cm×25 cm, 30 cm×30 cm, and 35 cm×35 cm, as shown in Fig. 2. The quantities of synthetic ice of various sizes are determined according to the statistics of Arctic ice sizes. The dimensional distribution of brash ice in Arctic brash ice conditions follows a log-normal distribution (Tuovinen, 1979):

Simulation setup
A numerical simulation is conducted to examine the resistance of an ice-going ship under a brash ice condition by using the FE-based FSI method. A comparative analysis is then performed using the results of both the numerical simulations and the experiments on the synthetic model ice. The numerical simulations share the ship model with the experiments as well as a highly similar test setup.
Numerical modeling primarily involves four aspects, i.e., ship, air, water, and ice, as shown in Table 2. The material of the ship is deemed a rigid body; the ice is of the same shape and size as that used in the experiment, and its material is set as an isotropic elastic body. Ice breaking is ignored in the simulation, with only the ship-ice interaction taken into account in order for it to resemble the experimental scenario. Models of the ship body and the ice are shown in Fig. 3.
Throughout the numerical simulation process, accurate simulation of the pressure in the fluid domain is critical because it determines whether the floating ice experiences the correct buoyant force and, hence, whether precise kinetic characteristics can be captured during the ship-ice interactions. Fig. 4 shows the simulation results for the hydrostatic pressure. The pressure is stable and is the same as the true hydrostatic pressure at one atmospheric pressure, so that the model ice can float in water under the action of the gravity and buoyancy.
The preprocessing is conducted on ANSYS, where the fluid domains and ice models are generated, and is followed by meshing. The element types and properties of the materials as well as the load, the constraint, and the collision and boundary conditions are set. A velocity load along the X direction is added to the ship, and varied as 0.2, 0.4, 0.6, 0.8, 1.0 and 1.2 m/s, respectively, which forces the ship to move along only the X-axis, without the displacement along any other degree of freedom. A gravitational acceleration load is exerted on the model ice to restrict the ice    boundary and keep the ice within the fluid domain. Finally, ship-ice and ice-ice contacts are set as automatic surfaceto-surface contact, and the water boundary as non-reflecting. By completing all the above, the numerical conditions of the model are successfully setup.

Analysis of the numerical phenomena
The numerical simulation is repeated for 18 navigation scenarios constituted by the six speed points under three concentrations of brash ice. The results of the simulations are processed on the post-processing software LS-PRE-POST to obtain animations of the model ship navigating the brash ice conditions. Fig. 5a shows the beginning of the simulation, Fig. 5b shows its end, and the right side of Fig.  6 shows snapshots of the interaction between the model ship and ice during navigation. Similarly, ship navigation in the experiment is video-recorded to obtain snapshots corresponding to the ship-ice interaction before, during, and after the experiment, as shown in Figs. 5c and 5d, and the left side of Fig. 6. Based on the experimental phenomena, the snapshots acquired through the simulations are analyzed to observe whether the numerical method has correctly simulated the various scenarios. Fig. 5 illustrates the comparison between brash ice distribution before and after navigation under a brash ice concentration of 90%. According to Figs. 5c and 5d, the model  GUO Chun-yu et al. China Ocean Eng., 2018, Vol. 32, No. 5, P. 546-556 551 ship drove the brash ice to move and constantly pushed it aside during the experiment; thus, a large amount of brash ice is accumulated in the advance direction and on both sides under the brash ice conditions, forming a sparsely distributed trail behind the ship. Figs. 5a and 5b show that the numerical simulation achieves similar patterns. Hence, in terms of simulating the overall phenomena, the numerical simulation and experiment are in good agreement. In addition to simulating the overall phenomenon of navigation, the details of the ship-ice interaction are important for simulation. A comparison in terms of the details of the ship-ice interaction between the numerical simulation and the experiment is shown in Fig. 6. When the model ship interacts with the brash ice, the ice rolling, accumulating, slipping, and sinking, all of which are typical phenomena associated with ship-ice interaction, and directly affect the resistance of the ship model to brash ice. The experiments reveal that the ship model navigating at a lower speed pushed aside the brash ice at the bow in a steady and slow fashion, and drove some of it forward. At this stage, the brash ices mostly accumulate up and slip; in case of high brash ice concentration, accumulating is particularly prominent. With a slight impact on the instantaneous brash ice load, these two phenomena significantly influence the average brash ice resistance. A massive amount of brash ice is accumulated at the bow in particular once the ship has navigated for some distance, thus significantly increase the average resistance. At a higher ship speed, ship-ice collisions intensified. As a result, the brash ice turns over more often upon impact with the ship. In the experiments, the heaving motion of the ship is significant, and causes the brash ice to frequently sink. The rolling and sinking of the brash ice not only affect the average resistance, but also significantly increase the instantaneous load during the ship-ice interaction. Observation has shown that neither of these phenomena occurs alone, in either the experiment or numerical simulation, in a manner that it is separable. In the numerical simulations, these typical physical phenomena of brash ice, i.e., rolling, slipping, accumulating, and sinking, are distinctly observed as well.
In summary, the FE simulation based on the FSI method successfully simulates numerically the entire process of an ice-going ship navigating in brash ice, and delivers the results that are consistent with experimental results in terms of the overall phenomenon of ship model navigation and the detailed traits of ship-ice interaction.

Analysis of the numerical data
Data from the experiments and numerical simulations are analyzed. At a data sampling frequency of 50 Hz in the experiment, the ship navigates several times for each condition to yield the improved accuracy. The average brash ice resistance for single navigation is calculated based on the recorded data. Data from all instances of navigation under the same conditions are then further averaged to obtain the average brash ice resistance for each condition. At a sampling frequency of 200 Hz in the numerical simulations, the curves of the average brash ice resistance and instantaneous load over time in the numerical simulation are obtained directly from LS-PREPOST. By calculating the gap between the corresponding numerical and experimental values, the relative error in the numerical simulation is analyzed. The average brash ice resistance values obtained through experiments and numerical simulation are listed in Table 3.
According to the table, a relative error of less than 50% commonly occurred between the numerical simulation and experiment, which is slightly better than results obtained by Kim. Even so, compared with the accuracy of CFD prediction in the context of the resistance performance in open waters, the FE-based numerical simulation method is still not sufficiently accurate to independently predict the resistance of the ship in brash ice, but can serve as a supplement or reference for other methods.
Ship-ice collisions are the major factor affecting the resistance of ice-going ships navigating under brash ice conditions, but the ship speed and brash ice concentration also influence the net resistance of ice-going ships. However, for the same speed and brash ice concentration, the results show a difference due to the random brash ice distribution and collision processes. Therefore, while the average brash ice resistance does display certain patterns as a whole, these are still substantially different from the pattern in which the ship resistance varies with the velocity in open water. That is, the variation in the brash ice resistance with speed under the same brash ice concentration cannot be mistakenly described as a smooth and steady curve. In order to better compare the results of the numerical simulations and the experiments, brash ice resistance data from the experiment and numerical simulation under various sailing situations are fitted linearly, as shown in Fig. 7. Fig. 7 shows a small gap between the linear fitting trends of the experimental results and numerical simulation data. The slope of the linear fit from the numerical simulations is higher than that from the experiments, which is consistent with the fact that the numerical simulation error at a concentration of 60% is relatively higher than that at other concentrations. Under 80% and 90% intensity, the results of the experiments and numerical simulations are satisfactorily consistent with regard to the linear fitting trend. As we only calculated for one ship type and three intensities, this is insufficient to prove, using the numerical simulation prediction, that the resistance in brash ice is better for high brash ice intensity than that for low brash ice intensity. Table 3 shows that there is no error distribution law for the resistance in brash ice. On the whole, the numerical simulation perfectly reflects the trend of linear variation of the resistance in brash ice, and has important reference value.
When sailing in ice, in addition to the average resistance of brash ice, the instantaneous load generated by the ship-ice interaction is an important physical parameter that can be simulated and calculated through the FEM. In this example, the changes in intensity of the brash ice and the speed of the ship influenced the ship-ice interaction, leading to corresponding variations in the instantaneous load curve. To observe whether the numerical simulations correctly predicted such effects, this study extracted the unsteadiness curves of the resistance in brash ice from the numerical simulation and compared them under two scenarios-the same brash ice intensity and different speeds (Fig. 8), and different brash ice intensities and the same speed (Fig. 9).
To obtain more significant comparison results, the unsteadiness curves from the numerical simulations for 0.4 m/s and 1.0 m/s, two speed points with a relatively large gap, are compared under brash ice intensities of 60%, 80%, and 90% as shown in Figs. 8a, 8b, and 8c, respectively. As shown in these figures, regardless of the brash ice intensity, the mean value of the unsteadiness curve at 1 m/s is larger than that at 0.4 m/s, with more intensive peaks and more severe oscillations. This is because for the same brash ice intensity, as the speed increases, the ice and the ship interacted more often, resulting in a greater acting force and stronger effects. This is seen as the brash ice rolling and sinking more frequently under ship-ice interactions, increases the average resistance of the brash ice, and an elevated curve, with more peaks and more severe oscillations in the unsteadiness curve. All these show that the results of the numerical simulations reflected how the increase in speed  GUO Chun-yu et al. China Ocean Eng., 2018, Vol. 32, No. 5, P. 546-556 affected the ship-ice interaction.
Since three curves plotted in one figure would cause confusion, only the unsteadiness curves for two brash ice intensities at the same speed from the numerical simulations are selected for comparison. In the case of a ship speed of 0.6 m/s, the unsteadiness curves for the brash ice intensities of 60% and 90% are shown in Fig. 9a for comparison; in the case of a speed of 0.8 m/s, the unsteadiness curves for the intensities of 80% and 90% are shown in Fig. 9b. From the figures, it is clear that under the same ship speed, as brash ice became more intensive, the mean value increased while the peak value, oscillation amplitude, and peak intensity slightly increased. This is because at the same ship speed, as the brash ice became more intensive, the ship and the ice interacted more often, leading to a more significant ship-ice effect. In the experiments, when the amount of brash ice increases, there is slightly more accumulating, rolling, slipping, and sinking of the brash ice and larger average resistance. For the unsteadiness curves, this is represented by the elevated curves, a few more peaks, and slightly larger oscil-lation amplitudes. It should be noted that the ship-ice interactions mainly occur in the brash ice around the body of the ship, where an increase in the intensity significantly affects the accumulating of the brash ice, leading to more obvious changes in the mean values. Owing to the random distribution of brash ice, the peak value, oscillation frequency, and peak intensity of the curve only show the trends of overall increase rather than a significant change.
To analyze the characteristics of the unsteadiness curves from the numerical simulations in greater detail, we compare the curves obtained in the experiments with those obtained in the numerical simulations under the same conditions. It should be noted that due to the random distribution of brash ice, and the random and complex interaction between ship and ice, an analysis of the similarity in values is not required for the unsteadiness curves generated by the two methods. It can be seen that in the absence of the experiments, numerical simulations cannot comprehensively predict the resistance in brash ice independently. By comparing the unsteadiness curves, we have hoped to observe the features of the curves in both methods and demonstrate how the numerical simulations reflect the patterns of variations in the instantaneous load generated by the ship-ice interaction.
In order to generate a clearer comparison, the unsteadiness curve from the numerical simulations is filtered, and a data sampling frequency similar to that used in the experiment is used to add the mean values of both curves. By analyzing the four conditions at two speeds of 0.4 and 1.0 m/s, the unsteadiness curves for the brash ice intensities of 60% and 90% are generated as shown in Figs. 10a-10d, respectively. Fig. 10 shows that when the sampling frequency decreased, the curves of the numerical simulation and experiment show a similar pattern of variation that reflected the randomness and oscillations of the ship-ice interaction. By analyzing the unsteadiness curves, we observe that the results of the numerical simulations can adequately reflect the instantaneous load generated by the changes in the ship-ice interaction, and hence can be used as a supplemental meth-od to experimentation to gauge the resistance of brash ice. Numerical simulations are important for guidance prior to experiments, or when some experimental conditions are lacking.

Conclusions
Using the FEM-based FSI algorithm, a numerical simulation is conducted in this study to predict the brash ice resistance of an ice-going ship navigating a broken ice condition. The results are compared with those of the model ship experiments using synthetic ice in a towing tank. The experiments and numerical simulations are carried out under eighteen working conditions: the product of three brash ice concentrations and six different navigation speeds. Various ship-ice interaction phenomena, the average resistance values, the linear trend of resistance relative to speed, and instantaneous load curves are obtained both experimentally and numerically. The effects of the navigation speed of the model ship and brash ice concentration on instantaneous load are analyzed. The results indicate that the numerical simulations exhibited good agreement with the experiments in terms of both the phenomena and results. Although it cannot be directly applied to the prediction of the average resistance value, a numerical simulation is significant as a reference for the prediction of the linear trend of the average resistance relative to speed and the instantaneous load curve. It complements the model ship tests, and hence is promising for application. Many deficiencies still exist in the FEM for application in this area, which suggests the need for more numerical simulations under various conditions as well as comparisons with the experimental results, including the data from full-scale ship tests.