Development of numerical method to simulate flows around a ship in regular waves including the effect of ship propulsion plant model

Development of Numerical Method to Simulate Flows around a Ship in Regular Waves including the Effect of Ship Propulsion Plant Model VII International Conference on Computational Methods in Marine Engineering MARINE 2017 M. Visonneau, P. Queutey and D. Le Touz´e (Eds) DEVELOPMENT OF NUMERICAL METHOD TO SIMULATE FLOWS AROUND A SHIP IN REGULAR WAVES INCLUDING THE EFFECT OF SHIP PROPULSION PLANT MODEL KUNIHIDE OHASHI∗ ∗National Maritime Research Institute 6-38-1 Shinkawa, Mitaka, Tokyo, Japan e-mail: [email protected] Key words: Ship Propulsion Plant Model, RANS, Overset-Grid Method Abstract. A numerical method to simulate flows with propeller effects including the responce of a ship propulsion plant has been developed. The dynamics of a ship propulsion plant is modeled by the function of the diesel engine cntrol system. Propeller torque which is computed by the propeller model with the interaction of flow fields is put to the ship propusion plant model, then the propeller rotational speed is obtained by solving the equation of the rotational motion of the propeller shaft line. Present method can reproduce the time history of propeller rotational speed and torque in the condition with the regular waves. The amplitude of fluctuations shows agreement with the measured data. The detail analysis of flow fields and ship motions which is difficult to be obtained at the experiment is also carried out. 1 INTRODUCTION Recently, Reynolds Averaged Navier-Stokes(RANS) simulations are utilized at the design stage of ship performance, and the numerical simulations are gradually applied to the more complex problems. URANS solver which can cope with the overset-grid method is coupled with a ship propulsion plant model through a propeller model which compute the body forces which express the propulsion effect in the present research. The present method can treat the response of the ship propulsion plant and the fluctuation of the propeller torque with ship motions in waves. The ship propulsion plant model[1] based on the mathematical equation of diesel engine components is employed. Computational results are validated with the experimental data[2]. Additionally, the detail analysis for the flows around the hull with motions are carried out. The relation between the response of the ship propulsion plant and flow fields are examined. 1 975 First A. Author, Second B. Author and Third C. Author 2 COMPUTATIONAL METHOD 2.1 Base solver An in-house structured CFD solver[3] is employed. The governing equation is 3D RANS equation for incompressible flows. Artificial compressibility approach is used for the velocitypressure coupling. Spatial discretization is based on a finite-volume method. A cell centered layout is adopted in which flow variables are defined at the centroid of each cell and a control volume is a cell itself. Inviscid fluxes are evaluated by the third-order upwind scheme based on the flux-difference splitting of Roe. The evaluation of viscous fluxes is second-order accurate. For unsteady flow simulations, a dual time stepping approach is used in order to recover incompressibility at each time step. It is consisted from the second order two-step backward scheme for the physical time stepping and the first order Euler implicit scheme for the pseudo time. The linear equation system is solved by the symmetric Gauss-Seidel (SGS) method. For free surface treatment, an interface capturing method with a single phase level set approach is employed. Incoming Regular headsea waves are generated at the region inside of the computational domain[4]. Body motions are obtained by solving the equations of motion, and motions are taken into account by a moving grid technique with the grid deforming methodology. Grid velocities are contained in the inviscid terms to satisfy the geometrical conservation law. The grid velocities are derived from the volume where an each cell face sweeps. The boundary condition on a body is given as the velocities of the body motion. The regions where the overset relations are composed deform with the body motions to maintain the overset information, and the amount of deformations gradually decreases with the distance from the body surfaces. Such the way is adapted to be able to avoid the computational load with using the dynamic overset-grid method. 2.2 Overset-grid method The weight values for the overset-grid interpolation are determined by an in-house system[5]. The detail of the system can be found on [5], the summary is described. 1. The priority of the computational grid is set. 2. The cells of a lower priority grid and inside a body is identified (called as in-wall cell in here). 3. Receptors cells which the flow variables have to be interpolated from donor cells are defined. Two cells on a higher priority grid and facing to the outer boundary are set as receptor cells to satisfy the third order discretization of NS solver. Additionally, two cells neighborhood of in-wall cells, the cells of a lower priority grid and inside the domain of a higher priority grid are also set as the receptor cell. 4. The weight values for the overset interpolation are determined by solving the inverse problem based on Ferguson spline interpolation. Flow variables of the receptor cell are updated when the boundary condition is set. The forces and moments are integrated on the higher priority grid to eliminate the lapped region on body 2 976 First A. Author, Second B. Author and Third C. Author surfaces. At first, the cell face of the lower priority grid is divided into small pieces. Secondly, the small piece is projected to the cell face of the higher priority grid by using the normal vector of the higher priority face. Then the 2D solid angle is computed and the small piece is decided in or out of the higher priority face. Once the small piece is in the higher priority face, the area ratio of the piece is set to zero. Finally, the area ratio is integrated on the lower priority face, then we have the ratio to integrate the forces and moments on lower priority face. 2.3 Propeller model The propeller model based on the potential theory[6] is applied to archive the propulsion condition. The propeller effect is taken account by the body forces which are computed by the following equations. fx(r, θ)=Γ(r, θ)Vθ r− 1 2CDNc(r)1+(h(r)/r)2VoxVoθ 2πr (1) ft(r, θ)=Γ(r, θ)Vx r+ 1 2CDNc(r)1+(h(r)/r)2V2 oθ 2πr (2) fy=ftcosθ, fz=ftsinθ (3) where rand θare the cylindrical coordinate at the propeller plane. The propeller circle is divided into the fan-shaped segments at (r, θ). Γ(r, θ) is the vortex strength, and Γ(r, θ) is determined by the equation which is based on the boundary condition at the segment (r, θ). Vx,Vθare the inflow velocities to the segment, Vox,Voθ are the circumferential averaged velocities. CDis a drag coefficient which is given by the empirical formula, Nis the blade number of a propeller, c(r) is a chord length at each radial direction, h(r) is the pitch of a free stream vortex. The coupling between the flow field and the propeller model is made by the velocities Vx,Vθ. At first, the velocity at the segment center is interpolated from the computational grid. Then, the Γ(r, θ) is determined, and body forces are computed. Finally, the cross section between the computational grid and the propeller segment is searched, and the forces at the cell is derived by the multiplication of body forces and the cross section area. The coordinate transfer which the body fixed coordinate to the earth fixed coordinate or its inverse are made by using Euler angles. 2.4 Ship Propulsion Plant Model The ship propulsion plant model[1] based on the mathematical equation of diesel engine components is employed. The model can simulate the behavior of the diesel engine in the transient conditions interacting with the major components of the ship propulsion plant. At first, the speed sensor gain is determined by the sensor model. ¯ Zs=Ks(¯nsp −¯ne) (4) where Ksis speed sensor gain coefficient, ¯nsp is target engine rotational speed, ¯neis actual speed of engine. The next equation is the summation unit, and ¯ Zip is isodromic feedback amount, ¯ Zfb is proportional feed back amount. 3 977 First A. Author, Second B. Author and Third C. Author ¯ ψ=¯ Zs−¯ Zip −¯ Zfb (5) Dead band is applied by following equation. ¯ Z0=   ¯ ψ−ǫ/2if ¯ ψ > ǫ/2 0if |¯ ψ|≤ǫ/2 ¯ ψ+ǫ/2if ¯ ψ<−ǫ/2 (6) where ǫis dead band. The fuel flow rate is determined by the power piston model. Tpp d¯ hp dt =¯ Z0(7) where Tpp is the time constant of the power piston. Eq.(8) and eq.(9) are the feedback mechanism model to determine the isodromic feedback gain ¯ Zip and the proportional feedback gain ¯ Zfb. Ti d¯ Zip dt +¯ Zip =KiTi d¯ hp dt (8) ¯ Zfb =Kfb¯ hp(9) where Tiis the time constant of the isodromic feedback, Kiis the gain coefficient of isodromic feedback, Kfb is gain coefficient of the proportional feedback. The engine torque ¯ Qecan be computed by the equation based on the relationship between the fuel flow rate ¯ hpand engine revolution number ¯ne. ¯ Qe=0.5¯ h2/3 p+1.5¯ h1/3 p¯ne−¯n2 e(10) The variables ¯ne,¯ hp,¯ Qeare non-dimensionalized by the engine revolution number nMCR, the fuel flow hMCR and the torque QeMCR at the Maximum Continuous Rating(MCR) ¯ne=ne/nMCR,¯ hp=hp/hMCR,¯ Qe=Qe/QeMCR (11) Finally, the shaft rotational speed is obtained by solving the equation of the rotational motion of a shaft line. 2πIp dn dt =Qe−Qp(12) where Ipis moment of inertia of whole shaftline, Qpis propeller torque. The equations of the ship propulsion plant model are discretized by the first order. The propeller torque Qpis obtained by the simplified propeller model. The actual speed of engine neequals to the propeller rotational speed nin the similar way of the measurement condition. 4 978 First A. Author, Second B. Author and Third C. Author 3 COMPUTED RESULTS The container hull form which is utilized on the experiment[2] is selected. The model ship length is about 4.0m and Reynolds number is R=2.8×106. Froude number is Fn =0.128. The particular of the propeller is same with the experiment. The wave length ratio of incoming regular headsea waves is λ/L =0.9,1.0,1.1, and the wave height ratio hw/L is 0.025. The motions are free to pitch and heave. The non-dimensional physical time step size is set to ∆t=0.002 Table 1 shows the division number of computational grids arranged with the priority of the overset-grid method. Figure 1 shows the grids near the hull body. The minimum spacing on the wall is set as y+= 1 to apply k−ωSST turbulence model. The hull and refined rectangular grid are overlapped to the rectangular grid to generate the incoming regular headsea waves. Table 1: Division number of computational grid Grid IM×JM×KM Refined Rect. 65×65×65 Hull 193×217×65 Rect. 217×105×57 Figure 1: Computational grids Table 2 shows coefficients and time constants of ship propulsion plant model, target engine rotational speed. All the coefficients and time constants are used at the experiment and nondimensionalized. Figure 2 shows the comparison of time history of propeller rotational speed, Figure 3 shows the comparison of time history of propeller torque, and Figure 4 depicts the comparison of the non-dimensionalized fuel consumption. All the figures are at the condition λ/L =1.1, the propeller rotational speed and propeller torque are dimensionalized for the direct comparison with the experimental data. As similar with the experiment data, the time phase is corrected 5 979 First A. Author, Second B. Author and Third C. Author Table 2: Coefficients and time constants of ship propulsion plant model Ks12.5 Tpp 2.771 ×10−3 Ki1.4 Ti5.773 ×10−2 Kfb 0.75 nMCR 75.05 nsp 50.0 QMCR 4.573 ×10−6 Ip2.655 ×10−10 by the propeller rotational speed to start from the trough of fluctuations. The computational results can reproduce the fluctuations of the variables. The differences of the time averaged value between the measurement and computation are mainly caused by the difference of the time averaged value of the inflow to the propeller. The computational results are able to simulate the characteristics which the propeller rotational speed is decreasing to the contrast of the propeller torque, additionally the fuel consumption is increasing to maintain the target engine rotational speed. TIme[sec.] Prop. Revolution Number[rps] 0123456 9.4 9.5 9.6 9.7 9.8 9.9 Comp. Meas. Figure 2: Comparison of time history of propeller rotational speed TIme[sec.] Prop. Torque[Nm] 0123456 0.08 0.09 0.1 0.11 0.12 0.13 0.14 Comp. Meas. Figure 3: Comparison of time history of propeller torque From Figure 5 to Figure 7 show the comparisons of the double amplitude with changing λ/L. The computational results show agreement with the measured results. The double amplitude when the propeller rotational speed is fixed at the time averaged value with the ship propulsion 6 980 First A. Author, Second B. Author and Third C. Author TIme[sec.] Non-dimensional Fuel Consumption 0123456 0.25 0.3 0.35 0.4 0.45 0.5 0.55 Comp. Meas. Figure 4: Comparison of time history of fuel consumption plant model is also shown in Figure 7. The present computational result can simulate the decreasing of the amplitude with the ship propulsion plant model. Wave length Prop. Rotational Speed[rps] 0.7 0.8 0.9 1 1.1 1.2 1.3 0.1 0.2 0.3 0.4 Meas. Comp. Figure 5: Comparison of double amplitude of propeller rotational speed Wave length Prop. Toque[Nm] 0.7 0.8 0.9 11.1 1.2 1.3 0.01 0.02 0.03 0.04 Meas. Meas.(Prop. Rot. Speed Fixed) Comp. Comp.(Prop. Rot. Speed Fixed) Figure 6: Comparison of double amplitude of propeller torque Wave length Non-dimensional Fuel Consumption 0.7 0.8 0.9 1 1.1 1.2 1.3 0 0.05 0.1 0.15 0.2 Meas. Comp. Figure 7: Comparison of double amplitude of fuel consumption Figure 8 shows the time histories of the propeller rotational speed, ship motions and wave 7 981 First A. Author, Second B. Author and Third C. Author height within the one encounter period at λ/L =1.1 condition. The propeller rotational speed closely fluctuates with the heave motion of the ship. Figure 9 shows the free surface near the hull from T1 to T4 in Figure 8, and Figure 10 shows the axial velocity contour and cross flow vectors at the propeller plane in ship fixed coordinate. T1 The heave motion reaches maximum value, and the pitch motion takes bow up position. The ship nose appears from the free surface, and inflow velocity wto the propeller takes negative value. The propeller rotates in clockwise direction seeing from the stern, thus the axial velocity uof the starboard side is faster than the velocity of the port side at a calm water condition. On the contrary, the axial velocity of the port side is faster than the velocity uof the starboard side at the time T1. T2 The heave motion becomes negative value, and pitch motion takes bow down position. The inflow velocity wto the propeller still takes negative value. The magnitude of the velocity wis smaller than the value of the time T1. The axial velocity uis accelerated in both sides. T3 The heave motion reaches minimum value, and the pitch motion takes from bow down position to the bow up position. The inflow velocity wto the propeller takes positive value. The axial velocity uis accelerated in the starboard side. T4 The heave motion takes almost zero position, and the pitch motion takes bow up position again. The ship nose also appears from the free surface. The inflow velocity wto the propeller still takes positive value. The magnitude of the velocity wis smaller than the value of the time T3. The axial velocity uis accelerated as similar with the calm water condition. The inflow velocities are changed by the incoming waves and ship motions, then the propeller torque which is the major input value to the ship propulsion plant model varies. The propeller revolution number is changed with the behavior of the ship propulsion plant model. The interactions between the change of the propeller characteristics in waves and the vary of the propeller revolution number which is caused by the behavior of the ship propulsion plant model is major reason for the time fluctuation of the variables. 4 CONCLUSIONS - The numerical method including the behavior of ship propulsion plant model by using the URANS solver and the propeller model has developed. - From the comparisons with the experimental data, the present method can reproduce the time fluctuations of the fuel flow rate, the propeller torque and rotational speed. - The detail flow analysis including the effect of the ship propulsion plant model is become possible by the present method. 5 ACKNOWLEDGEMENT This work has been supported by JSPS KAKENHI Grant Number JP16K06919. 8 982 First A. Author, Second B. Author and Third C. Author Time[sec.] Heave[mm], Pitch[deg.] Prop. Revolution Number[rps] 0 0.2 0.4 0.6 0.8 1 1.2 -25 -20 -15 -10 -5 0 5 10 15 20 25 9.4 9.5 9.6 9.7 9.8 9.9 Prop. Revolution Number Heave Pitch T4 T1 T2 T3 Figure 8: Time history of propeller rotational speed and motions T1T2 T3T4 Figure 9: Instantaneous view of free surface and ship with ship propulsion plant model REFERENCES [1] O., Bondarenko, M., Kashiwagi, M., Naito, Dynamics of Diesel Engine in the Framework of Ship Propulsion Plant, Proc. of JASNAOE spring meeting, Vol.8, 2009. [2] Y., Kitagawa, K., Tanizawa, Y., Tsukada, M., Ueno, Development of a Methodology on Tank Model Test Measuring Speed Drop of Actual Ship under Waves, Journal of the Japan Society of Naval Architects and Ocean Engineers, Vol.22, pp.21-34, 2015. [3] K., Ohashi, T., Hino, N., Hirata, H., Kobayashi, Development of NS solver with a structured overset grid method The 28th Computational Fluid Dynamics Symposium, Japan, 2014. 9 983