CC BY-NC-ND
0+ MODEL
|jNAO^ Available online at www.sciencedirect.com
ScienceDirect
Iggftr
® International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7
http://www.journals.elsevier.com/international-journal-of-naval-architecture-and-ocean-engineering/
Time domain simulation for icebreaking and turning capability of bow-first
icebreaking models in level ice
Donghyeong Ko, Kyung-Duk Park*, Kyoungsoo Ahn
Hyundai Heavy Industries Co, Ltd., Ulsan, South Korea
Received 15 January 2016; revised 16 February 2016; accepted 24 February 2016 Available online ■ ■ ■
Abstract
Recent icebreaking ships need to be designed to enhance not only icebreaking capability but also turning ability. For the evaluation of ice resistance induced by an icebreaking hull form, HHI (Hyundai Heavy Industries) has developed the hybrid empirical formulas (Park et al., 2015) by considering the geometrical hull shape features, such as waterline and underwater sections. However, the empirical formulas have inherent limits to the precise estimation of the icebreaking and turning ability because the breaking process and the resulting pattern are ignored. For this reason, numerical calculation in time domain is performed to predict the icebreaking process and pattern. In the simulation, varying crushing stress according to velocity vectors and contact areas between hull and ice is newly introduced. Moreover, the simulation results were verified by comparing them with the model test results for three different bow-first icebreaking models.
Copyright © 2016 Society of Naval Architects of Korea. Production and hosting by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
Keywords: Ice resistance; Icebreaking and turning capability; Empirical formulas; Level ice
1. Introduction
Recently, Hyundai Heavy Industries (HHI) has been involved in many projects relating to the transportation of natural resources by icebreaking vessels from or through the arctic region, starting from the hull form development for a 190,000 dwt class arctic ore carrier (Park et al., 2011). For icebreaking ship design, it is essential to know how to evaluate and minimize the ice resistance as well as to confirm whether the ship satisfies the specific requirements, such as minimum icebreaking speed and turning ability given the ice conditions and the engine power.
For the evaluation of ice resistance, several empirical formulas have been introduced by many researchers, but the typical method proposed by Lindqvist (1989) has been widely used among ice model basins and ship yards. This method uses
* Corresponding author.
E-mail address: parkkd@hhi.co.kr (K.-D. Park).
Peer review under responsibility of Society of Naval Architects of Korea.
the basic hull information such as length, breadth, depth, hull angles, etc. Recently HHI had developed the hybrid methods (Park et al., 2015) to get a relatively high level of accuracy by considering geometrical hull shape, such as waterline at Design Load Water Line (DLWL) and underwater hull sections. In detail, ice resistance was assumed to comprise of three components: breaking, clearing and buoyancy resistances. The Shimanskii icebreaking resistance (1938), the Enkvist buoyancy resistance (1972) and the newly modified clearing resistance from Ionov (Poznyak and Ionov, 1981) were combined.
However, the empirical formulas have inherent limits to precise estimation of the icebreaking and turning ability because the icebreaking process and its pattern are ignored. That is, although the empirical formulas give a certain level of accuracy in ice resistance estimation, Shimanskii's formula does not consider the icebreaking process and its pattern around the ship's shoulder in sufficient detail, so it is not appropriate to use when judging the quality of the ship's turning ability.
http://dx.doi.org/10.1016/j.ijnaoe.2016.02.004
2092-6782/Copyright © 2016 Society of Naval Architects of Korea. Production and hosting by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).
+ MODEL
2 D. Ko et al. / International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7
On the other hand, Su (2011) suggested a time-domain numerical simulation method based on an icebreaking procedure with constant crushing stress. This method well represents the icebreaking patterns and estimates the turning ability of icebreaking ships. In addition, it shows a good agreement with full scale measurements for a small size icebreaker shown by Riska et al. (2001). Through this study, the authors tried to apply this technique to commercial icebreaking ships and checked whether it is suitable for the estimation of icebreaking patterns and turning capability or not.
However, it was found that the commercial icebreaking ships typically have long parallel middle body, so the ice resistance in time domain simulation increases dramatically due to the constant crushing stress around this area. Therefore, a new method with varying crushing stress is proposed to avoid the problem. In this study, the crushing stress is assumed to vary with not only the contact area between hull and ice (Ashby et al., 1986) but also the linear velocity vector and the angular velocity vector of ship. Moreover, its applicability to ice resistance estimation for three different bow-first ice-breaking models is tested comparing with the results of model tests. For the simulation of ship motions, the modular type mathematical model with multiple POD propulsors (Kim et al., 2006) is used.
2. Numerical simulation
2.1. Varying crushing stress
Ice resistance is comprised of three parts: breaking, buoyancy and clearing parts. This study is mainly focused on the breaking part in the numerical simulation, and the others are simply treated by the empirical formulas of Lindqvist (1989). The overall process of the numerical simulation is denoted by the flow chart in Fig. 1.
—e——&— -— —- -J,---'
—e—©—e— —e>-
6800 7000 7200 7400 7600 7800 X(m)
Fig. 2. Description of ice nodes generation.
Initial inputs are hull and ice information such as waterline shape, propeller open water characteristics and ice properties. The waterline of the hull and ice are marked in discretized nodes as shown in Fig. 2. As a ship moves forward, the contact area between the hull and ice is calculated, and whether its vertical load exceeds the bending failure load or not is checked according to the same logical process proposed by Su (2011). If bending failure occurs, the ice geometry is regenerated and then the forces and moment acting on the ship are calculated. During the process, the nodes which are not in contact with hull are excluded from the calculation to enhance the computational efficiency and the number of ice nodes in the far field from hull is reduced to avoid a waste of memory.
Su (2011) proposed using the constant crushing stress to calculate ice load in level ice, but it was found that this constant value makes the ice load increase sharply around the long parallel middle body of ship. To solve this problem, the concept of varying crushing stress defined by the contact area between hull and ice is newly introduced in this study.
Fig. 1. Flow chart of numerical simulation.
Detection contact area around the hull
Vertical load > Bending failure load
Ice breaking and regeneration
Calculation of forces and moment working on the ship
+ MODEL
D. Ko et al. / International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7
According to Ashby et al. (1986), the crushing stress sare be defined as a function of contact area A with Eq. (1).
where, Li is an idealized size of cubical independent cells, DL is a moved distance before the cell fails and PL is the average force during the period of contact.
The relation between crushing stress and contact areas is given as a theoretical curve in Fig. 3, where Li = 1, PL = 15 MN and DL = 0.02 m are used. In addition, the peak stresses measured in the field are assumed to be approximately 1—2 MPa for contact areas of Az 100 m2. By slightly moving the theoretical curve to the left using the assumption that the peak stress is approximately 1 MPa for contact areas of A = 100 m2, the resulting theoretical curve is defined by Eq. (2) which is valid for A > > L2 and higher than the minimum level of 0.33 MPa for very large contact areas.
: 0.33(l + 20.3/VA)
On the other hand, for very small contact areas, the crushing stress curve is constrained to be bounded by the observed data for the maximum stress, 10.5 MPa, to avoid the unrealistic situation in which the stress soars infinitely as contact areas decreases to zero (Sanderson, 2014).
When calculating the resulting contact area between hull and ice, the velocity vector of the ship (U) needs to be considered, which is comprised of two components, equivalent linear velocity ( Ui) and angular velocity ( Ur) as in Eq. (3). Then, as shown in Fig. 4, all the nodes along the waterline of hull have equivalent linear velocity in the case of going straight but different angular velocity at every position in the case of making turns.
U = U,+ Ur
Once the contact areas are determined, the resulting crushing stress is calculated according to the modified theoretical curve of Ashby et al. (1986) and divided into two
10' 105
Fig. 3. Theoretical curve of Ashby et al. (1986) bounded by the observed data.
Fig. 4. Schematic diagram of coordinate systems.
components: normal and vertical to the ice sheet, considering the slope angles between hull and ice (a, b in Fig. 5). As shown in Fig. 5, vertical crushing stress and normal crushing stress sn act only as the ship contacts the ice, but they are zero if the ship moves far from the ice or if such a contact between hull and ice happens. Therefore, the resulting crushing stresses are given as in Eqs. (4) and (5), where the normal stress is the cosine a component and the vertical stress is the cosine b component of sarea in the side view of Fig. 5. They are also computed by the inner product between velocity vector (U) and inward normal vector of ice sheet (sn) as in the top view of Fig. 5. However, the two crushing stresses are not to exceed sarea and its minimum is to be zero with t\ and t2, as shown in Eq. (5).
Sn = S area ' tl ' COS a S = S area ' ^2 ' COS b
( U 'On
if > 0
U 'O„
if , n < 0
U' sn 1 if , n > 0
U' sn 0 if —-— < 0
To see how the varying crushing stress and the inner product mentioned above affect the numerical results, time series ice resistances for the same model are compared in Fig. 6. The black solid line shows the results with constant crushing stress (2.5 MPa proposed by Su, 2011) and red dotted line shows results with the varying crushing stress and the inner product proposed in this study. In the case of constant crushing stress, very many peak values higher than 10,000 kN are observed during the simulation and thus the resulting ship speed gradually drops to zero. However, the results of the varying crushing stress in the study show reasonable peak values lower than black solid line and thus are nearing the constant speed quickly and smoothly.
+ MODEL
4 D. Ko et al. / International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7
Side view
Fig. 5. Description of contact area and crushing stress.
For more detailed comparison, time series ice resistances with and without the effect of inner product U • sn are shown in Fig. 7. Both simulations are carried out with the varying crushing stress suggested in this study for the same model. The black solid line is the result with the inner product and the red dotted line without it. From the comparison of numerical results, it is known that the inner product also plays a role in reducing the peak values of ice resistance to some extent.
Fig. 7. Comparison of time series ice resistances between with and without the effect of inner product U • an:.
2.2. Simulation of ship performance based on modular model
For icebreaking ships, the assessment of turning maneuver is very important. For the precise consideration of icebreaking procedure and the resulting pattern, the step-by-step time domain simulation is necessary. This simulation is consistent with the maneuvering simulation model. The only difference is the small time step for icebreaking procedure.
In this study, the modular type mathematical model is used. This model is based on the twin POD propulsor model (Kim et al., 2006), and then slightly modified for three POD pro-pulsors as shown in Fig. 8. The coordinate system is already shown in Fig. 4. XYZ is a right-handed coordinate system with
Time (s)
Fig. 6. Comparison of time series ice resistances between constant and varying crushing stress.
Z vertically downward and its origin is located at the center of gravity. By the Newton's second law, the ship's equations of motion of (u, v, r) can be described in Eq. (6).
M$ (u - vr) = XH + Xp + Xr + Xice
M$(v + ur) = Yh + Xp + Yr + Yce (6)
Izz$r = Nh + Np + Nr + Nice
where, M is the mass of ship, Mx is the added mass in surge motion. The terms (X, Y and N) with subscript (H, p and R) represent the hydrodynamic forces and moments acting on
Fig. 8. A example of stern hull form with three POD propulsors.
+ MODEL
D. Ko et al. / International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7
hull and those induced by POD propellers and POD struts respectively. Xice and Yice are the forces induced by ice in surge and sway motion, and Nice is the induced moment by ice in yaw motion respectively, which are described in the previous section in detail.
Hydrodynamic force and moment acting on hull can be written as follows.
Xh = -MxU + (My + Xvr) ■ vr + X„ $ v2 + Xrr ■ r2 + X(u) Yh = -MyV + Yv ■ v + Yvvv ■ V3 + Yr ■ r + Yrrr ■ r3 + Yvrr ■ vr2 + Yn
Nh = -Jzz ■r+Nv^v + Nvvv ■v + Nr ■ r+Nrrr ■ r + Nvr
■ vr
■ rv
■ rv
where, the hydrodynamic coefficients acting on hull are followed by the equation of Lee et al. (2003). Mx and My are the added mass in surge motion and sway motion respectively, Jzz is the added moment of inertia in yaw motion and X(u) represents ship resistance which varies as a function of u. The partial derivatives of forces and moment with respect to motion parameters are denoted by the subscript notation as below:
dY 9v'
d3Y dvdv2
Hydrodynamic force and moment induced by POD propellers is modeled as follows.
XP = (1 - t)(TP + TC + TS) cos d YP = (1 -t) (TP + TC + TS) sin d
NP = xp • (1 -1) (TP + TC + T^ sin d - yp • (1 -1) (TP - TS) cos d
where Tp C S are the thrust forces induced by port, center and starboard side POD propellers which are positioned at xp and yp respectively, t is the average thrust deduction factor measured in model testing, and d is rudder angle.
Hydrodynamic force and moment induced by POD struts can be written as follows.
Xr = -(1-tR)(FN+Fcn)sind Yr = (1+ aH) (FP+FC+FN)cosd
Nr = (xr+aHXH) (FN+FC+Fsn)cosd-yR • (1 -r) (fp-F^) sind
where Fpp'C'S are the forces induced by port, center and starboard side POD struts respectively, tR is the rudder deduction fraction, aH is the coefficient of hydrodynamic force acting on hull induced by rudder action, xH is the x coordinate on which aH acts, xR and yR are the x and y coordinates of each POD strut.
3. Numerical results and verifications
3.1. Icebreaking capability
For three different ship models shown in Fig. 9, the numerical simulations for bow-first icebreaking were carried out
and their results are compared with model test results. Ship A and Ship B are both characterized by 20 degrees of the bow angle (stem profile angle) with an ice forefoot but L/B values are slightly different. Ship C is characterized by 25 degrees of the bow angle without ice forefoot but its L/B is the same to ship A. Ice model tests for Ship A and Ship B are conducted at the ice model basin of AARC in Finland, and the one for Ship C at the ice model basin of HSVA in Germany. The principal dimensions of the ships and ice model test conditions are summarized in Table 1.
The numerical simulations have been carried out with 10,000 s of computation time and 0.05 s of time step for all ships. The target thickness of level ice is 1.5 m for all ships, but the flexural strengths are 500 kPa for Ship A and Ship B, but 600 kPa for Ship C according to ice model test condition. The comparison of icebreaking patterns between Ship A and Ship C is plotted in Fig. 10. Because Ship C has larger bow angle and harder flexural strength than Ship A, the ice cusp of Ship C along the center line is bigger than the one of Ship A and this brings about the worst icebreaking ability (highest ice resistance and lowest resulting ship speed) among three ships.
Three vessel speeds at 100% MCR power are compared in Fig. 11. Although the net thrust (total propeller thrust minus the open water resistance at the same speed) and the propeller revolution speed were different between two model basins because different stock propellers were applied, the same propeller characteristics of AARC ice model basin were all applied in this numerical simulation for direct comparison. The resulting average speed at 100% MCR power in the simulation is 2.4 m/s for Ship A, 2.6 m/s for Ship B and 1.1m/ s for Ship C respectively and they show a good agreement with the model test results qualitatively. For reference, model test results were 2.5 m/s for Ship A, 2.7 m/s for Ship B and 1.3 m/s for Ship C.
To verify the relation between ship speed and ice resistance, the calculated ice resistances for three ships are compared with the model test results and presented in Fig. 12. The ice resistance values in the numerical simulation agree well with the model test results in qualitative terms, but the gradient according to ship speed is relatively lower than model test results.
3.2. Turning capability
Following the numerical simulation for going straight ahead, the turning simulation was carried out on the assumption that all POD propulsors turn clockwise with 35° simultaneously. The target ice thickness and flexural strength during turning test were 1.7 m and 500 kPa for Ship A and Ship B in AARC ice model basin, but 1.5 m and 600 kPa for Ship C in HSVA ice model basin. The diameters of the turning circle were obtained from 360° turning, but those of the model tests were extrapolated from the measured initial turning path which was much less than a quarter of a circle because of the limitation of tank size.
Approximate turning diameters measured in model test results were 9300 m for ship A, 15,000 m for ship B and
ARTICLE IN PRESS
+ MODEL
6 D. Ko et al. / International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7
Fig. 9. Photos of ship models (Ship A, B and C from left).
Table 1
Principal dimension of ships and ice model test conditions.
Ship A Ship B Ship C
L/B 5.5 6.0 5.5
B/T 4.3 3.9 4.3
Bow angle (degree) 20 20 25
Ice model basin AARC AARC HSVA
Ice thickness (m) 1.5 1.5 1.5
Flexural strength (kPa) 500 500 600
6000 m for ship C respectively. During the simulation of initial turning with 100% MCR power applied, all ships were stuck in the ice due to the excessive ice loads caused by buoyancy and clearing resistances estimated by empirical formulas. Thus, 150% MCR power was applied to all turning simulations to match up turning speed with model test results. As a result, the estimated turning diameters from the numerical simulations were 9500 m for ship A, 15,600 m for ship B and 7800 m for ship C respectively as shown in Fig. 13.
4. Conclusions and future works
Time domain simulation for the estimation of icebreaking patterns and turning ability in level ice is presented. In this study, varying crushing stress is newly suggested to solve the problem of excessive crushing stress around long parallel middle body in commercial icebreaking ships. The values of crushing stress are calculated by a modified version of the theoretical curve given by Ashby et al. (1986) but bounded by the maximum value of 10.5 MPa referring to the experimental data. In addition, the inner product of the ship's velocity vector and the normal vector of contact area is superimposed on the crushing stress calculation to regulate inappropriate interaction between hull nodes and ice nodes. For its verification and
Fig. 11. Time history of resulting icebreaking speeds at 100% MCR power.
validation, the numerical results are compared with model test results for bow-first icebreaking models. Although there is a quantitative difference in gradient, the simulations well reflect the change of hull form in the estimation of icebreaking and turning ability in level ice.
For future works, this time-domain simulation will be modified by considering the Enkvist buoyancy resistance
7000 6000 ¡Z 5000 | 4000 '1 3000 о 2000 1000
°0 0.5 1 1.5 2 2.5 3 3.5 4
Speed(m/s)
Fig. 12. Comparison of ice resistances between simulation and model test.
+ MODEL
D. Ko et al. / International Journal of Naval Architecture and Ocean Engineering xx (2016) 1—7 7
---Simulation (ship A) -Simulation (ship C)
X (m) x ю"
Fig. 13. Turning simulations of three model ships.
(1972) and the modified clearing resistance (Park et al., 2015) to enhance its accuracy instead of using Lindqvist's empirical formula (1989). Additionally the developed numerical code needs to be expanded to be utilized for stern-first icebreaking ships as well as bow-first, and then various commercial icebreaking ship models and its model test results will be compared in order to confirm the availability of the numerical simulation.
References
Ashby, M.F., Palmer, A.C., Thouless, M., Goodman, D.J., Howard, M.W., Hallam, S.D., Murrell, S.A.F., Jones, N., Sanderson, T.J.O., Ponter, A.R.S., 1986. Nonsimultaneous failure and ice loads on Arctic structures. In: Offshore Technology Conference 1986, NO. OTC 5127, Houston, USA, 399-404.
Enkvist, E., 1972. On the Ice Resistance Encountered by Ships Operating in the Continuous Mode of Icebreaking. Report No.24. The Swedish Academy of Engineering Science in Finland, Helsinki, Finland.
Kim, Y.G., Kim, S.Y., Park, Y.H., Yu, B.S., Lee, S.W., 2006. Prediction of maneuverability of a ship with POD propulsion system. J. Soc. Nav. Archit. Korea 43 (2), 164-170.
Lee, T.I., Ahn, K.S., Lee, H.S., Yum, D.J., 2003. On an empirical prediction of hydrodynamic coefficients for modern ship hulls. In: Proceedings of MARSIM'03, vol. III. Kanazawa, Japan, 25-28th August, RC-1-1-RC-1-8.
Lindqvist, G., 1989. A straightforward method for calclation of ice resistance of ships. In: Proceedings of 10th International Conference on Port and Ocean Engineering under Arctic Conditions (POAC), Lulea, Sweden.
Park, K.D., Chung, Y.K., Jang, Y.S., Kim, H.S., Molyneux, D., 2011. Development of hull forms for a 190,000 dwt icebreaking ore carrier. In: Preceedings of the ASME 2011 30th International Conferences on Ocean, Offshore and Arctic Engineering (OMAE2011), Rotterdam, The Netherlands.
Park, K.D., Kim, M.C., Kim, H.S., 2015. Calculation of ice clearing resistance using normal vector of hull form and direct calculation of buoyancy force under the hull. Int. J. Nav. Archit. Ocean Eng. 7, 699-707.
Poznyak, I.I., Ionov, B.P., 1981. The division of icebreaking resistance into components. In: Preceedings 6th STAR Symposium. SNAME, New York, 249-252.
Riska, K., Leiviska, T., Nyman, T., Fransson, L., Lehtonen, J., Eronen, H., Backman, A., 2001. Ice performance of the Swedish multi-purpose icebreaker Tor Viking II. Proc. PAC 849-865.
Sanderson, T.J.O., 2014. Ice Mechanics. Springer, Lexington, USA.
Shimanskii, Y.A., 1938. Conditional Standards of Ice Qualities of a Ship. Trans. Arctic Research Institute, Northern Sea Route Administration Publishing House, vol. 130. Translation T-381-01 by Engineering Consulting and Translation Center (ECTC), P.O. Box 1377, Jackson Heights, New York, NY 11372. Leningrad.
Su, B., 2011. Numerical Predictions of Global and Local Ice Loads on Ships. (Ph.D. thesis). Department of Marine Technology, Norwegian University of Science and Technology, Norway.
