Calculation Methods of Icebreaking Capability for a Double-Acting Polar Ship

: As a key parameter, icebreaking capability is often used to judge whether a polar ship could navigate in level ice at a certain speed. This paper presents two methods to calculate icebreaking capability. The ﬁrst one is a static method based on the estimation of ice resistance under di ﬀ erent ice thicknesses and ship speeds. The second is a dynamic method that involves solving the equation of motion. A series of model tests with a double-acting icebreaking tanker were also carried out in the ice basin of the Krylov State Research Center to measure ice resistances. The simulated ice resistances were compared with model tests results for both ahead and astern running operations. The calculated icebreaking capability based on static and dynamic methods was validated with the model test result. A good agreement was achieved between measurement and simulation. The discrepancy between the model test result and the result simulated by the static or dynamic method was minor.


Introduction
For safe operation, ship performance in ice should be taken into account seriously.The ability of a polar ship to break ice is measured in uniform ice conditions by the speed that the ship can attain in ice of a certain thickness [1].The speed that the ship can attain in ice is determined by the propeller thrust available to overcome the ice resistance.It is thus clear that the performance in ice is influenced by the resisting forces and the propulsive forces, which can be improved (resisting forces minimized and propulsive forces maximized) by the hull shape and propulsion design, respectively.
It is important to make a reasonable prediction of ice loads on a polar ship.Many researchers have calculated ice loads numerically with the ice resistance formula [2][3][4] or numerical simulations [5][6][7][8][9][10][11][12].However, the icebreaking capability for a polar ship has been seldom studied.Su [13] used a circumferential crack method to calculate the icebreaking capability of icebreaker Tor Viking II and compared the simulated result with full-scale data.
Zhou et al. [14] proposed a squared ice grid method to calculate dynamic ice forces.The ice grid could be taken as a simple static finite element.It would fail against the structure when a certain failure criterion was met.Similarly to the method developed by Su et al. [9], it is difficult to choose the speed-dependent icebreaking length coefficient and the bending capacity coefficient, which should be taken from model tests or full-scale trials.In the present paper, we proposed a formula to calculate the two coefficients, which are functions of ice thickness and relative velocity of ice and structure, based on which the ice resistance is calculated.The simulated ice resistances were then compared with the model test data.

Mathematical Model
There are two methods which could be applied to estimate icebreaking capability for a polar ship.The first method, the static method, considers the static balance between the ice resistance and available thrust from prolusion system of the polar ship.The second method, the dynamic method, obtains the steady speed of the polar ship under the available thrust in level ice by solving the equation of motion.The principles of two methods are described in the following section.

Static Method
Using the static method, we need to evaluate the ice resistance of polar ship and available thrust under certain ice conditions.A method was developed to calculate ice resistance by Zhou et al. [14], which was used in the present study.The ice force consists of the icebreaking force and the ice submersion force.The icebreaking is mainly due to the local ice crushing force during individual bending failure of an ice wedge.To calculate the icebreaking force, the hull was discretized into many nodes at the waterline (as shown in Figure 1).J. Mar.Sci.Eng.2020, 8, 179 2 of 24

Mathematical Model
There are two methods which could be applied to estimate icebreaking capability for a polar ship.The first method, the static method, considers the static balance between the ice resistance and available thrust from prolusion system of the polar ship.The second method, the dynamic method, obtains the steady speed of the polar ship under the available thrust in level ice by solving the equation of motion.The principles of two methods are described in the following section.

Static Method
Using the static method, we need to evaluate the ice resistance of polar ship and available thrust under certain ice conditions.A method was developed to calculate ice resistance by Zhou et al. [14], which was used in the present study.The ice force consists of the icebreaking force and the ice submersion force.The icebreaking is mainly due to the local ice crushing force during individual bending failure of an ice wedge.To calculate the icebreaking force, the hull was discretized into many nodes at the waterline (as shown in Figure 1).The icebreaking force occurs when ice collides with the hull.The local icebreaking force on each hull node was estimated at each time step.The detailed information on the calculation of local ice force Fbi was previously described by Zhou et al. [14].The total icebreaking force were calculated as the sum of all local forces Fbi on each panel or hull node on the waterline at a certain time instance: The loop of calculating the icebreaking force in the time domain is presented in Figure 2. The ice submersion force occurs when the broken ice floe fragmented from intact ice pile around the hull underwater.When a ship moves forward against the ice, the ice breaks at the bow area and the broken ice floe becomes submerged and slides along the hull.There were almost no ice floes sliding at the bottom of the ship model, and most of the ice floes moved laterally beneath the bordering ice sheet at the low drift speed range [15].Based on the ice resistance formula described by Lindqvist [2], the speed dependence of ice submersion resistance due to the loss of potential energy of the submerged ice floes and friction between the hull and ice floes was written as In early design, the net thrust was estimated by Juva and Riska [16] as the following ice rule: Figure 3 shows two typical curves of ice resistance and net thrust as a factor of ship speed with a certain ice thickness.The intersection of the two curves was taken as the expected steady speed ve of the ship.By changing the ice thickness, the h-v curve could be calculated.The icebreaking force occurs when ice collides with the hull.The local icebreaking force on each hull node was estimated at each time step.The detailed information on the calculation of local ice force F bi was previously described by Zhou et al. [14].The total icebreaking force were calculated as the sum of all local forces F bi on each panel or hull node on the waterline at a certain time instance: The loop of calculating the icebreaking force in the time domain is presented in Figure 2. The ice submersion force occurs when the broken ice floe fragmented from intact ice pile around the hull underwater.When a ship moves forward against the ice, the ice breaks at the bow area and the broken ice floe becomes submerged and slides along the hull.There were almost no ice floes sliding at the bottom of the ship model, and most of the ice floes moved laterally beneath the bordering ice sheet at the low drift speed range [15].Based on the ice resistance formula described by Lindqvist [2], the speed dependence of ice submersion resistance due to the loss of potential energy of the submerged ice floes and friction between the hull and ice floes was written as In early design, the net thrust was estimated by Juva and Riska [16] as the following ice rule: Figure 3 shows two typical curves of ice resistance and net thrust as a factor of ship speed with a certain ice thickness.The intersection of the two curves was taken as the expected steady speed v e of the ship.By changing the ice thickness, the h-v curve could be calculated.

Dynamic Method
The dynamic method is based on solving the equation of motion including both the net thrust and ice resistance at each time step.The steady speed could be achieved when the net thrust and ice resistance came to a balance.According to Newton's second law, the linear coupled differential equations of motion in the body-fixed coordinate could be written in the following form: The Runge-Kutta method was used to solve the resulting equations of motion.

Model Tests in the Ice Basin
The model tests were carried out in the ice basin of the Krylov State Research Center (KSRC).The length of the tank is 45 m including the trimming area (as shown in Figure 4).The width and depth are 6 and 1.75 m respectively.The refrigeration system ensures uniform cooling of air over the tank down to -29 °C .The towing carriage could move by the rails running along both sides of the tank at a speed of 0.0005~1 m/s.This carriage possesses a number of openings serving to spray sea water uniformly over the test tank to generate the model ice.Usually, it takes two days to generate

Dynamic Method
The dynamic method is based on solving the equation of motion including both the net thrust and ice resistance at each time step.The steady speed could be achieved when the net thrust and ice resistance came to a balance.According to Newton's second law, the linear coupled differential equations of motion in the body-fixed coordinate could be written in the following form: The Runge-Kutta method was used to solve the resulting equations of motion.

Model Tests in the Ice Basin
The model tests were carried out in the ice basin of the Krylov State Research Center (KSRC).The length of the tank is 45 m including the trimming area (as shown in Figure 4).The width and depth are 6 and 1.75 m respectively.The refrigeration system ensures uniform cooling of air over the tank down to -29 °C .The towing carriage could move by the rails running along both sides of the tank at a speed of 0.0005~1 m/s.This carriage possesses a number of openings serving to spray sea water uniformly over the test tank to generate the model ice.Usually, it takes two days to generate

Dynamic Method
The dynamic method is based on solving the equation of motion including both the net thrust and ice resistance at each time step.The steady speed could be achieved when the net thrust and ice resistance came to a balance.According to Newton's second law, the linear coupled differential equations of motion in the body-fixed coordinate could be written in the following form: The Runge-Kutta method was used to solve the resulting equations of motion.

Model Tests in the Ice Basin
The model tests were carried out in the ice basin of the Krylov State Research Center (KSRC).The length of the tank is 45 m including the trimming area (as shown in Figure 4).The width and depth are 6 and 1.75 m respectively.The refrigeration system ensures uniform cooling of air over the tank down to −29 • C. The towing carriage could move by the rails running along both sides of the tank at a speed of 0.0005~1 m/s.This carriage possesses a number of openings serving to spray sea water uniformly over the test tank to generate the model ice.Usually, it takes two days to generate an ice sheet.The ice thickness ranges from 10 to 80 mm with an error less than ±1.5 mm.The bending strength of the model ice is 10~70 kPa with relative error of no more than ±20%.
an ice sheet.The ice thickness ranges from 10 to 80 mm with an error less than ±1.5 mm.The bending strength of the model ice is 10~70 kPa with relative error of no more than ±20%.
An artic tanker model with PC4 ice class was used in the model test.The dimensions of the tanker were reduced to the model scale with a geometric scale factor of λ = 42.5.The particulars of the full-scale vessel in the design draft are given in Table 1. Figure 5 presents the bow and stern of the ship model used for the tests.An artic tanker model with PC4 ice class was used in the model test.The dimensions of the tanker were reduced to the model scale with a geometric scale factor of λ = 42.5.The particulars of the full-scale vessel in the design draft are given in Table 1. Figure 5 presents the bow and stern of the ship model used for the tests.an ice sheet.The ice thickness ranges from 10 to 80 mm with an error less than ±1.5 mm.The bending strength of the model ice is 10~70 kPa with relative error of no more than ±20%.An artic tanker model with PC4 ice class was used in the model test.The dimensions of the tanker were reduced to the model scale with a geometric scale factor of λ = 42.5.The particulars of the full-scale vessel in the design draft are given in Table 1. Figure 5 presents the bow and stern of the ship model used for the tests.The ship model was constrained to the carriage, then the force could be measured.A dynamometer was installed and applied to measure the ice force acting on the model ship.The general setups of the model tests are shown in Figure 6.The ship model was constrained to the carriage, then the force could be measured.A dynamometer was installed and applied to measure the ice force acting on the model ship.The general setups of the model tests are shown in Figure 6.The ice sheet was made by cooling, spraying, freezing and tempering.Two ice sheets with thicknesses of 1.3 and 1.6 m in full scale were generated for these tests.Two operation modes were tried in the ice.The first mode was ahead running as normal.The other one was astern running with operating propellers.The constant speeds used to tow the ship model were 1.27, 1.90, 2.54 and 3.80 kn.In total, there were 16 tests performed with the specific information shown in Table 2.A picture of the model test is given in Figure 7.The ice sheet was made by cooling, spraying, freezing and tempering.Two ice sheets with thicknesses of 1.3 and 1.6 m in full scale were generated for these tests.Two operation modes were tried in the ice.The first mode was ahead running as normal.The other one was astern running with operating propellers.The constant speeds used to tow the ship model were 1.27, 1.90, 2.54 and 3.80 kn.In total, there were 16 tests performed with the specific information shown in Table 2.A picture of the model test is given in Figure 7.

Simulation Setup
The ice field was composed of squared ice grids with the same size.The length of each ice grid

Simulation Setup
The ice field was composed of squared ice grids with the same size.The length of each ice grid was set to be the breaking length.The initial positions of the ship and ice field are shown in Figure 8.

Simulation Setup
The ice field was composed of squared ice grids with the same size.The length of each ice grid was set to be the breaking length.The initial positions of the ship and ice field are shown in Figure 8.
The general ice data with respects to bending strength, crushing strength, elastic modulus, water and ice density, ice-hull friction coefficient, Poisson ratio and gravity acceleration are given in Table 3.The bending strength, crushing strength, elastic modulus and ice-hull friction coefficient were measured directly from the model test.The unknown parameters were Poisson ratio, and water and ice density.The Poisson ratio was selected based on previous study [17,18].The water contains many chemical components including salts.Herein, the density was taken as sea water density with 1025 kg.m −3 .The ice density varies from 720 to 940 kg.m −3 with an average of 910 kg.m −3 [19].The ice density was taken as 900 kg.m −3 in the present study.All data shown in the present paper are in full scale if not specified otherwise.The general ice data with respects to bending strength, crushing strength, elastic modulus, water and ice density, ice-hull friction coefficient, Poisson ratio and gravity acceleration are given in Table 3.The bending strength, crushing strength, elastic modulus and ice-hull friction coefficient were measured directly from the model test.The unknown parameters were Poisson ratio, and water and ice density.The Poisson ratio was selected based on previous study [17,18].The water contains many chemical components including salts.Herein, the density was taken as sea water density with 1025 kg.m −3 .The ice density varies from 720 to 940 kg.m −3 with an average of 910 kg.m −3 [19].The ice density was taken as 900 kg.m −3 in the present study.All data shown in the present paper are in full scale if not specified otherwise.

Simulation Examples
Test 11 was simulated using an ice drift velocity of 1.27 kn and ice thickness of 1.3 m.The speed-dependent icebreaking length coefficient C l was selected to be 0.24 while the bending capacity coefficient C f was set to be 2.5.The longitudinal ice force exposed to the ship was of interest for analysis of both operation modes herein.Therefore, only the longitudinal ice force was considered and taken as the ice force for short.The time series of the ice force components due to breaking and submersion are given in Figure 9.It can be seen that the ice submersion force increases steadily to a constant value.The breaking component is nonlinear and changing rapidly, significantly depending on the failure of the ice grid around the waterline of the hull.
analysis of both operation modes herein.Therefore, only the longitudinal ice force was considered and taken as the ice force for short.The time series of the ice force components due to breaking and submersion are given in Figure 9.It can be seen that the ice submersion force increases steadily to a constant value.The breaking component is nonlinear and changing rapidly, significantly depending on the failure of the ice grid around the waterline of the hull.The process would be discussed briefly at several time intervals.The time intervals t1, t2 and t3 were selected for study, which corresponded to the first three peaks at the early stage of the ice-hull interaction.The time intervals are labelled in Figure 9.The icebreaking snapshots at each time instant are presented in Figure 10.At time t1, the central grid is breaking against the hull, then the ice grid fails in bending and the icebreaking force turns to zero.When it comes to time t2, two ice grids neighboring the failed ice grid at time instant t1, are colliding with the hull.Then the three ice grids failed after t2 and the fourth ice grid starts to intrude the hull and the three peak occurs at t3.
Figure 11 shows a part of the ice force based on Figure 9 at the steady state when the whole bow part is travelled into ice field.The time intervals t4~t7 where the peaks occur are included and the ice-ship interaction snapshots are shown in Figure 12.There are seven ice grids interacting with the hull at time instant t4.Then the foremost ice grid which are colliding with the ship fails, but the total ice force continues to augment as the overall interaction area increases at t5.When it evolves to t6, the peak decreases though the number of icebreaking grid is the same as that at t5.Only four ice grids are breaking against the hull as t7 and the force tends to decrease afterwards.The process would be discussed briefly at several time intervals.The time intervals t1, t2 and t3 were selected for study, which corresponded to the first three peaks at the early stage of the ice-hull interaction.The time intervals are labelled in Figure 9.The icebreaking snapshots at each time instant are presented in Figure 10.At time t1, the central grid is breaking against the hull, then the ice grid fails in bending and the icebreaking force turns to zero.When it comes to time t2, two ice grids neighboring the failed ice grid at time instant t1, are colliding with the hull.Then the three ice grids failed after t2 and the fourth ice grid starts to intrude the hull and the three peak occurs at t3. Figure 11 shows a part of the ice force based on Figure 9 at the steady state when the whole bow part is travelled into ice field.The time intervals t4~t7 where the peaks occur are included and the ice-ship interaction snapshots are shown in Figure 12.There are seven ice grids interacting with the hull at time instant t4.Then the foremost ice grid which are colliding with the ship fails, but the total ice force continues to augment as the overall interaction area increases at t5.When it evolves to t6, the peak decreases though the number of icebreaking grid is the same as that at t5.Only four ice grids are breaking against the hull as t7 and the force tends to decrease afterwards.

Sensitivity Study
To give insight into the present method, sensitivity study has been performed with respect to influence of simulation time interval, speed-dependent icebreaking length coefficient Cl and bending capacity coefficient Cf on the ice loads.
The specified ice thickness was 1.3 m and the speed was 1.27 kn (as same as that used in Test 11).Different values of time step lengths (0.03, 0.05, 0.1, 0.2, 0.5, 1 s) were applied.The ice resistance encountered by the ship as a function of time interval is plotted in Figure 13.As shown in Figure 13, the simulated ice load converges quite well when the time step length decreases.It comes to be steady at around 2400 kN when the time step drops to 0.2s.The computation time versus time step length is also given in Figure 14.From Figure 14, it is clear that the computation time increases exponentially with decreasing time interval.Considering the balance between computation time and the accuracy of the numerical simulation, the time step length of 0.1 s was appropriate and used in the following

Sensitivity Study
To give insight into the present method, sensitivity study has been performed with respect to influence of simulation time interval, speed-dependent icebreaking length coefficient C l and bending capacity coefficient C f on the ice loads.
The specified ice thickness was 1.3 m and the speed was 1.27 kn (as same as that used in Test 11).Different values of time step lengths (0.03, 0.05, 0.1, 0.2, 0.5, 1 s) were applied.The ice resistance encountered by the ship as a function of time interval is plotted in Figure 13.As shown in Figure 13, the simulated ice load converges quite well when the time step length decreases.It comes to be steady at around 2400 kN when the time step drops to 0.2s.The computation time versus time step length is also given in Figure 14.From Figure 14, it is clear that the computation time increases exponentially with decreasing time interval.Considering the balance between computation time and the accuracy of the numerical simulation, the time step length of 0.1 s was appropriate and used in the following simulations.

Sensitivity Study
To give insight into the present method, sensitivity study has been performed with respect to influence of simulation time interval, speed-dependent icebreaking length coefficient Cl and bending capacity coefficient Cf on the ice loads.
The specified ice thickness was 1.3 m and the speed was 1.27 kn (as same as that used in Test 11).Different values of time step lengths (0.03, 0.05, 0.1, 0.2, 0.5, 1 s) were applied.The ice resistance encountered by the ship as a function of time interval is plotted in Figure 13.As shown in Figure 13, the simulated ice load converges quite well when the time step length decreases.It comes to be steady at around 2400 kN when the time step drops to 0.2s.The computation time versus time step length is also given in Figure 14.From Figure 14, it is clear that the computation time increases exponentially with decreasing time interval.Considering the balance between computation time and the accuracy of the numerical simulation, the time step length of 0.1 s was appropriate and used in the following simulations.In addition, effects of speed-dependent icebreaking length coefficient Cl and bending capacity coefficient Cf were studied.In total, 18 icebreaking length coefficients ranging from 0.18 to 0.35 were considered in the simulation.The corresponding results are presented in Figure 15.It is shown in Figure 15 that as the icebreaking parameter Cl increases, the ice resistance decreases monotonously.Increasing Cl means large ice floe broken from intact ice sheet and long distance of travelling is needed before consequent collision.The impact of bending capacity coefficient Cf on ice resistance is given in Figure 16.It shows that the ice resistance increases as Cf increases in general.The relationship is almost linear.The ice resistance rises from 2730 kN to 3200 kN when the bending capacity coefficient increases from 0.20 to 0.35.In addition, effects of speed-dependent icebreaking length coefficient C l and bending capacity coefficient C f were studied.In total, 18 icebreaking length coefficients ranging from 0.18 to 0.35 were considered in the simulation.The corresponding results are presented in Figure 15.It is shown in Figure 15 that as the icebreaking parameter C l increases, the ice resistance decreases monotonously.Increasing C l means large ice floe broken from intact ice sheet and long distance of travelling is needed before consequent collision.The impact of bending capacity coefficient C f on ice resistance is given in Figure 16.It shows that the ice resistance increases as C f increases in general.The relationship is almost linear.The ice resistance rises from 2730 kN to 3200 kN when the bending capacity coefficient increases from 0.20 to 0.35.Increasing Cl means large ice floe broken from intact ice sheet and long distance of travelling is needed before consequent collision.The impact of bending capacity coefficient Cf on ice resistance is given in Figure 16.It shows that the ice resistance increases as Cf increases in general.The relationship is almost linear.The ice resistance rises from 2730 kN to 3200 kN when the bending capacity coefficient increases from 0.20 to 0.35.As shown in Figs. 15 and 16, ice resistance is very sensitive to the parameters Cf and Cl.This remains a question.In previous research, the parameters Cf and Cl were often selected based on measurements or ice resistance by some empirical formulas [9,14,17,18].In fact, measurements are sometimes unavailable due to limited condition.Empirical formula brings different uncertainties.Therefore, we proposed a new formula which could be suitable for the simulation.The length of ice grid is equal to the icebreaking radius.By including two parameters Cvs and Cls, the icebreaking radius was expressed as  =      (5) (5) where = 1 − (   3√ℎ  ) 0.7 The ice bending capacity coefficient could be written as a function of relative velocity between ice and hull: As shown in Figures 15 and 16, ice resistance is very sensitive to the parameters C f and C l .This remains a question.In previous research, the parameters C f and C l were often selected based on measurements or ice resistance by some empirical formulas [9,14,17,18].In fact, measurements are sometimes unavailable due to limited condition.Empirical formula brings different uncertainties.Therefore, we proposed a new formula which could be suitable for the simulation.The length of ice grid is equal to the icebreaking radius.By including two parameters C vs and C ls , the icebreaking radius was expressed as where The ice bending capacity coefficient could be written as a function of relative velocity between ice and hull: The following simulation was carried out based on the proposed method to determine the key factors C f and C l .

Icebreaking Ahead
The ahead running tests were first considered with input parameters shown in Tables 2 and 3.By using different towing speeds of the hull, four tests in ice sheet with a thickness of 1.3 m were simulated.The transition stage of ice force was ignored.The steady state was studied after the bow was completely inside the ice sheet.The time series of the ice force are shown in Figure 17.The time series of the ice force present a similar trend.The ice force is cyclic with several typical peaks.As the speed increases, the number of peaks increases because of the higher number of ice grid failures within the same time length.The maximum ice forces are around 5000 kN for Test 11~13 and 6000 kN for Test 14. Test 21~24 were also simulated by increasing ice thickness to 1.6m and the numerical results are displayed in Figure 18.
was completely inside the ice sheet.The time series of the ice force are shown in Figure 17.The time series of the ice force present a similar trend.The ice force is cyclic with several typical peaks.As the speed increases, the number of peaks increases because of the higher number of ice grid failures within the same time length.The maximum ice forces are around 5000 kN for Test 11~13 and 6000 kN for Test 14. Test 21~24 were also simulated by increasing ice thickness to 1.6m and the numerical results are displayed in Figure 18.The ice force was averaged as ice resistance.The ice resistances derived from numerical simulation as well as model test were compared herein and are presented in Figure 19.It is indicated that the general trend of two curves agree well with each other.The ice resistance curves go up as the towing speed increases.The relative errors between the calculated and measured ice resistance values are discussed in a subsequent section.The ice force was averaged as ice resistance.The ice resistances derived from numerical simulation as well as model test were compared herein and are presented in Figure 19.It is indicated that the general trend of two curves agree well with each other.The ice resistance curves go up as the towing speed increases.The relative errors between the calculated and measured ice resistance values are discussed in a subsequent section.

Icebreaking Astern
The astern running tests were included with the input parameters shown in Table 2 and Table 3  as

Icebreaking Astern
The astern running tests were included with the input parameters shown in Tables 2 and 3 as well.Eight tests were simulated in ice sheets with thicknesses of 1.3 and 1.6 m with different towing speeds of the hull.The time series of the ice force at the steady state are shown in Figure 20 for Tests 15~18 and Figure 21 for Tests 25~28.The simulated and measured ice resistances are compared in Figure 22.

Icebreaking Astern
The astern running tests were included with the input parameters shown in Table 2 and Table 3 as well.Eight tests were simulated in ice sheets with thicknesses of 1.3 and 1.6 m with different towing speeds of the hull.The time series of the ice force at the steady state are shown in Figure 20 for Tests 15~18 and Figure 21 for Tests 25~28.The simulated and measured ice resistances are compared in Figure 22.The discrepancy between the model test and simulation results are summarized in Table 4 for all the tests.The minimum absolute error is 0.7% while the maximum is up to 13.7%.The average absolute error is 6.7%.With the assistance of aft propellers, the astern running operation causes a lower ice resistance than the ahead running operation in a 1.3-m-thick ice sheet.The discrepancy between the model test and simulation results are summarized in Table 4 for all the tests.The minimum absolute error is 0.7% while the maximum is up to 13.7%.The average absolute error is 6.7%.With the assistance of aft propellers, the astern running operation causes a lower ice resistance than the ahead running operation in a 1.3-m-thick ice sheet.As described above, the static method is based on the static balance of the ice resistance and thrust delivered by propellers.The bollard for the heading forward operation is 3305 kN and the open water speed of the ship could reach up to 16 kn.The thrust curve as a function of ship speed is calculated according to Equation (3) and presented in Figure 23.The simulated ice resistances are also plotted in Figure 23 for all the ahead running tests.For the curve composed of tests with an ice thickness of 1.3 m, the ice resistance counteracts 3261 kN thrust.The ship speed vsh1 at balance is 2.39 kn.For the curve of the tests with an ice thickness of 1.6m, there is no intersection with the thrust curve.The ice resistance has to be extrapolated as a dashed line in Figure 23.The corresponding speed vsh2 is 0.59 kn at the intersection.also plotted in Figure 23 for all the ahead running tests.For the curve composed of tests with an ice thickness of 1.3m, the ice resistance counteracts 3261 kN thrust.The ship speed vsh1 at balance is 2.39 kn.For the curve of the tests with an ice thickness of 1.6m, there is no intersection with the thrust curve.The ice resistance has to be extrapolated as a dashed line in Figure 23.The corresponding speed vsh2 is 0.59 kn at the intersection.
The measured ice resistances and evaluated thrust force are shown in Figure 24 for all the ahead running tests.For the curve composed of the tests with an ice thickness of 1.3m, the ice resistance counteracts thrust with 3179 kN.The ship speed vmh1 at balance is 1.54 kn.For the curve with an ice thickness of 1.6m, there is no intersection with the thrust curve.The ice resistance has to be extrapolated as a dashed line in Figure 24.The corresponding speed vmh2 is 0.81 kn at the intersection.
The ice thickness-speed curve (h-v curve) is derived based on the obtained ship speeds.The results of the h-v curve are shown in Figure 25.It is indicated that at the astern running speed of 1 knot, the ice-breaking capability is 1.52 m based on the measurement and 1.53 m from the simulation.The simulated result is close to the measured value.The measured ice resistances and evaluated thrust force are shown in Figure 24 for all the ahead running tests.For the curve composed of the tests with an ice thickness of 1.3 m, the ice resistance counteracts thrust with 3179 kN.The ship speed vmh1 at balance is 1.54 kn.For the curve with an ice thickness of 1.6m, there is no intersection with the thrust curve.The ice resistance has to be extrapolated as a dashed line in Figure 24.The corresponding speed vmh2 is 0.81 kn at the intersection.The simulated ice resistances are given in Figure 26 for all astern running tests.For the curve composed of tests with an ice thickness of 1.3m, the ice resistance counteracts the thrust with 2795 kN.The ship speed vss1 at balance is 2.01 kn.For the curve of tests with an ice thickness of 1.6m, there is no intersection with the thrust curve.The ice resistance has to be extrapolated as a dashed line in Figure 26.The corresponding speed vss2 is 0.64 kn at the intersection.
The measured ice resistances and evaluated thrust force are shown in Figure 27 for all ahead running tests.For both ice resistance curves with ice thicknesses of 1.3 and 1.6 m, there is no intersection with the thrust curve.The ice resistances have to be extrapolated.At the intersections, the corresponding speeds vms1 and vms2 are 1.16 kn, 0.87 kn respectively.The h-v curves for both The simulated ice resistances are given in Figure 26 for all astern running tests.For the curve composed of tests with an ice thickness of 1.3 m, the ice resistance counteracts the thrust with 2795 kN.The ship speed vss1 at balance is 2.01 kn.For the curve of tests with an ice thickness of 1.6m, there is no intersection with the thrust curve.The ice resistance has to be extrapolated as a dashed line in Figure 26.The corresponding speed vss2 is 0.64 kn at the intersection.
intersection with the thrust curve.The ice resistances have to be extrapolated.At the intersections, the corresponding speeds vms1 and vms2 are 1.16 kn, 0.87 kn respectively.The h-v curves for both the simulation and the measurement are expressed in Figure 28.The results show that at the astern running speed of 1 knot, the ice-breaking capability is 1.47 m based on the measurement and 1.52 m from the simulation.The discrepancy was minor and could be accepted.The measured ice resistances and evaluated thrust force are shown in Figure 27 for all ahead running tests.For both ice resistance curves with ice thicknesses of 1.3 and 1.6 m, there is no intersection with the thrust curve.The ice resistances have to be extrapolated.At the intersections, the corresponding speeds vms1 and vms2 are 1.16 kn, 0.87 kn respectively.The h-v curves for both the simulation and the measurement are expressed in Figure 28.The results show that at the astern running speed of 1 knot, the ice-breaking capability is 1.47 m based on the measurement and 1.52 m from the simulation.The discrepancy was minor and could be accepted.

Dynamic Method
By solving the equation of motion (Equation (4)), the ship speed could be calculated in time domain.The Dof of the equation was reduced to surge.This is enough for an icebreaking capability analysis.The added mass of the ship in surge is taken as 3.5% of the initial mass.The potential damping and restoring terms are set to be zero.The thrust due to propellers is included by application of Equation.(3).The properties of the ice sheet in Table 3 are used.
For ahead running, there are eight ice thicknesses ranging from 0.8 to 1.55 m considered in the simulation.The simulation stops when the ship speed becomes relatively stable.The simulated time series of ship speeds under different ice thicknesses are shown in Figure 29.The ship starts to accelerate under maximum thrust and zero ice resistance at the beginning.Then, the thrust decreases vms1 vms2 Figure 28.h-v curves from the measurement and the simulation using the static method for astern running.

Dynamic Method
By solving the equation of motion (Equation (4)), the ship speed could be calculated in time domain.The Dof of the equation was reduced to surge.This is enough for an icebreaking capability analysis.The added mass of the ship in surge is taken as 3.5% of the initial mass.The potential damping and restoring terms are set to be zero.The thrust due to propellers is included by application of Equation (3).The properties of the ice sheet in Table 3 are used.
For ahead running, there are eight ice thicknesses ranging from 0.8 to 1.55 m considered in the simulation.simulation stops when the ship speed becomes relatively stable.The simulated time series of ship speeds under different ice thicknesses are shown in Figure 29.The ship starts to accelerate under maximum thrust and zero ice resistance at the beginning.Then, the thrust decreases and the ice resistance increases eventually as the ship speed increases.If the resistance is lower than the ice resistance before reaching the maximum speed, the ship speed will increase to a steady value in ice thicknesses of 0.8 m and 1.0, as shown in Figure 29.If the ice resistance is larger than the thrust when the speed reaches the peak, the ship speed will drop afterwards and come to a relatively steady value under ice thicknesses of 1.1, 1.2, 1.3, 1.4, 1.45, and 1.55 m.For the simulation in 1.55 m ice, the ship comes to stop completely and gets stuck in ice because the maximum thrust is not enough to overcome the ice resistance.As the ice thickness increases, the ship tends to vibrate severely.The amplitude of the speed vibration could be up to 0.17 m/s for the test in 1.5 m ice.
The attainable speeds of the ship were collected based on Figure 29.Then, h-v curves from the dynamic simulation and the model test are compared and given in Figure 30.The results show that the icebreaking capability of the ship is 1.45 m at a speed of 1 kn by using the dynamic analysis method.The simulated result is slightly lower than the model test result at 1.52.
thrust when the speed reaches the peak, the ship speed will drop afterwards and come to a relatively steady value under ice thicknesses of 1.1, 1.2, 1.3, 1.4, 1.45, and 1.55 m.For the simulation in 1.55 m ice, the ship comes to stop completely and gets stuck in ice because the maximum thrust is not enough to overcome the ice resistance.As the ice thickness increases, the ship tends to vibrate severely.The amplitude of the speed vibration could be up to 0.17 m/s for the test in 1.5 m ice.The attainable speeds of the ship were collected based on Figure 29.Then, h-v curves from the dynamic simulation and the model test are compared and given in Figure 30.The results show that the icebreaking capability of the ship is 1.45 m at a speed of 1 kn by using the dynamic analysis method.The simulated result is slightly lower than the model test result at 1.52.For astern running, there are also eight ice thicknesses ranging from 0.8 to 1.5 m considered in the simulation.The simulation stops when the ship speed becomes relatively stable.The simulated time series of the ship speeds under different ice thicknesses are shown in Figure 31.The phenomenon is quite similar to that that occurred to the ahead running tests in Figure 29.We are not going into details herein.For the simulation in 1.5 m ice, the ship comes to halt completely because of the inefficient thrust.The attainable ship speeds are collected as a factor of ice thickness.Then, hv curves obtained from dynamic simulation and model test are c given in Figure 32.The results show that the icebreaking capability of the ship is 1.42 m at a speed of 1 kn for astern running.The simulated result is slightly lower than the model test result at 1.47.For astern running, there are also eight ice thicknesses ranging from 0.8 to 1.5 m considered in the simulation.The simulation stops when the ship speed becomes relatively stable.The simulated time series of the ship speeds under different ice thicknesses are shown in Figure 31.The phenomenon is quite similar to that that occurred to the ahead running tests in Figure 29.We are not going into details herein.For the simulation in 1.5 m ice, the ship comes to halt completely because of the inefficient thrust.The attainable ship speeds are collected as a factor of ice thickness.Then, h-v curves obtained from dynamic simulation and model test are c given in Figure 32.The results show that the icebreaking capability of the ship is 1.42 m at a speed of 1 kn for astern running.The simulated result is slightly lower than the model test result at 1.47.
phenomenon is quite similar to that that occurred to the ahead running tests in Figure 29.We are not going into details herein.For the simulation in 1.5 m ice, the ship comes to halt completely because of the inefficient thrust.The attainable ship speeds are collected as a factor of ice thickness.Then, hv curves obtained from dynamic simulation and model test are c given in Figure 32.The results show that the icebreaking capability of the ship is 1.42 m at a speed of 1 kn for astern running.The simulated result is slightly lower than the model test result at 1.47.The calculated and icebreaking capabilities of the polar ship at the speed of 1 kn are summarized in Table 5.

Conclusions
This paper presents a static method and a dynamic method to calculate icebreaking capability numerically.The simulated ice results were compared with the model tests results.The following conclusions could be derived based on the present study.
(1) The ice resistance is very sensitive to the speed-dependent icebreaking length coefficient and the bending capacity coefficient.
(2) The scatter between the measured and simulated ice resistance is small and less than 13.7%.This means that the proposed formula to calculate the speed-dependent icebreaking length coefficient and bending capacity coefficient are reasonable.
(3) Considering the balance between ice resistance and thrust available, the icebreaking The calculated and measured icebreaking capabilities of the polar ship at the speed of 1 kn are summarized in Table 5.

Conclusions
This paper presents a static method and a dynamic method to calculate icebreaking capability numerically.The simulated ice results were compared with the model tests results.The following conclusions could be derived based on the present study.
(1) The ice resistance is very sensitive to the speed-dependent icebreaking length coefficient and the bending capacity coefficient.
(2) The scatter between the measured and simulated ice resistance is small and less than 13.7%.This means that the proposed formula to calculate the speed-dependent icebreaking length coefficient and bending capacity coefficient are reasonable.
(3) Considering the balance between ice resistance and thrust available, the icebreaking capability for ahead running is relatively high compared to astern running based on both the model test and the simulation results.
(4) The static method gives a higher estimation on the icebreaking capability than the model test, but the difference is within 3.4%.
(5) The static method underestimates the icebreaking capability by 3.4%~4.6%.This could be accepted in early design of polar ship.
It should be noted that the effect of the propeller washing on icebreaking capability is not considered in both present numerical simulation and model test.The high-speed wake may assist the polar ship to break ice and thus, the capability.All conclusions are made based on the present polar tanker herein.More studies are needed to validate the present method further.

Figure 3 .
Figure 3. Illustration of the ice resistance and net thrust.

Figure 4 .
Figure 4. Ice tank of the Krylov State Research Center (KSRC).

Figure 4 .
Figure 4. Ice tank of the Krylov State Research Center (KSRC).

Figure 5 .
Figure 5. Photographs of the ship model.(a) Bow view (b) Stern view.

Figure 6 .
Figure 6.Side view illustration of the test system.

Figure 7 .
Figure 7. Model test of ahead running.

Figure 8 .
Figure 8.Initial positions of the ship and ice field.

Figure 8 .
Figure 8.Initial positions of the ship and ice field.

Figure 9 .
Figure 9.Time series of the total ice force and ice force components.

Figure 9 .
Figure 9.Time series of the total ice force and ice force components.

Figure 11 .Figure 11 .Figure 12 .
Figure 11.Part of the total ice force and ice force components.

Figure 13 .
Figure 13.Simulated ice resistance as a function of the time step length

Figure 14 .
Figure 14.Computation time as a function of the time step length.

Figure 15
Figure15that as the icebreaking parameter Cl increases, the ice resistance decreases monotonously.Increasing Cl means large ice floe broken from intact ice sheet and long distance of travelling is needed before consequent collision.The impact of bending capacity coefficient Cf on ice resistance is given in Figure16.It shows that the ice resistance increases as Cf increases in general.The relationship is almost linear.The ice resistance rises from 2730 kN to 3200 kN when the bending capacity coefficient increases from 0.20 to 0.35.

Figure 15 .
Figure 15.Simulated global ice load as a function of Cl.

Figure 16 .
Figure 16.Simulated global ice load as a function of C f.

Figure 17 .
Figure 17.Time series of the ice forces during the ahead running operation in ice sheet with a thickness of 1.3m.(a) Test 11 (b) Test 12. (c) Test 13 (d) Test 14.

Figure 17 .
Figure 17.Time series of the ice forces during the ahead running operation in ice sheet with a thickness of 1.3 (a) Test 11 (b) Test 12. (c) Test 13 (d) Test 14.

Figure 18 .
Figure 18.Time series of the ice forces during the ahead running operation in ice sheet with a thickness of 1.6 m.(a) Test 21 (b) Test 22. (c) Test 23 (d) Test 24.

Figure 18 .
Figure 18.Time series of the ice forces during the ahead running operation in ice sheet with a thickness of 1.6 m.(a) Test 21 (b) Test 22. (c) Test 23 (d) Test 24.

Figure 19 .
Figure 19.Comparison of the ice resistance for the ahead running mode.
well.Eight tests were simulated in ice sheets with thicknesses of 1.3 and 1.6 m with different towing speeds of the The time series of the ice force at the steady state are shown in Figure 20 for Tests 15~18 and Figure 21 for Tests 25~28.The simulated and measured ice resistances are compared in

Figure 19 .
Figure 19.Comparison of the ice resistance for the ahead running mode.

Figure 19 .
Figure 19.Comparison of the ice resistance for the ahead running mode.

Figure 20 .Figure 20 .Figure 20 .Figure 21 .
Figure 20.Time series of the ice forces during the astern running operation in ice sheet with a thickness of 1.3 m.(a) Test 15 (b) Test 16.(c) Test 17 (d) Test 18.

Figure 23 .
Figure 23.Simulated ice resistance and evaluated thrust as a factor of ship speed for ahead running.

Figure 23 .
Figure 23.Simulated ice resistance and evaluated thrust as a factor of ship speed for ahead running.

Figure 23 .
Figure 23.Simulated ice resistance and evaluated thrust as a factor of ship speed for ahead running.

Figure 24 .
Figure 24.Measured ice resistance and evaluated thrust as a factor of ship speed for ahead running.

Figure 24 .
Figure 24.Measured ice resistance and evaluated thrust as a factor of ship speed for ahead running.The ice thickness-speed curve (h-v curve) is derived based on the obtained ship speeds.The results of the h-v curve are shown in Figure 25.It is indicated that at the astern running speed of 1 knot, the ice-breaking capability is 1.52 m based on the measurement and 1.53 m from the simulation.The simulated result is close to the measured value.

Figure 25 .
Figure 25.h-v curves from the measurement and the simulation using the static method for ahead running.

Figure 25 .
Figure 25.h-v curves from the measurement and the simulation using the static method for ahead running.

Figure 26 .
Figure 26.Simulated ice resistance and evaluated thrust as a factor of ship speed for astern running.

Figure 26 .
Figure 26.Simulated ice resistance and evaluated thrust as a factor of ship speed for astern running.

J 24 Figure 27 .
Figure 27.Measured ice resistance and evaluated thrust as a factor of ship speed for astern running.

Figure 27 .
Figure 27.Measured ice resistance and evaluated thrust as a factor of ship speed for astern running.

Figure 27 .
Figure 27.Measured ice resistance and evaluated thrust as a factor of ship speed for astern running.

Figure 28 .
Figure 28.h-v curves from the measurement and the simulation using the static method for astern running.

Figure 29 .
Figure 29.Simulated time series of the dynamic ship speed in the ice sheets with different thicknesses for ahead running.

Figure 29 .
Figure 29.Simulated time of the dynamic ship speed in the ice sheets with different thicknesses for ahead running.J. Mar.Sci.Eng.2020, 8, 179 21 of 24

Figure 30 .
Figure30.curves from the measurement and the simulation using the dynamic method for ahead running.

Figure 31 .
Figure 31.Simulated time series of the dynamic ship speed in ice the sheets with different thicknesses for astern running.

Figure 31 .
Figure 31.Simulated time series of the dynamic ship speed in ice the sheets with different thicknesses for astern running.

Figure 32 .
Figure 32.h-v curves from measurement and simulation using the dynamic method for astern running.

Figure 32 .
Figure 32.h-v curves from measurement and simulation using the dynamic method for astern running.

Table 1 .
Main parameters of the arctic tanker.

(a) Figure 4. Ice
tank of the Krylov State Research Center (KSRC).

Table 1 .
Main parameters of the arctic tanker.

Table 1 .
Main parameters of the arctic tanker.
Figure 6.Side view illustration of the test system.

Table 3 .
Ice parameters used in the simulation.

Table 3 .
Ice parameters used in the simulation.

Table 4 .
Comparison of the resistance for all tests from the model test and simulation.
Figure 22.Comparison of the ice resistance for the astern running mode.

Table 4 .
Comparison of the ice resistance for all tests from the model test and simulation.

Table 5 .
Calculated and measured icebreaking capability.

Table 5 .
Calculated and measured icebreaking capability.