A Body-Nonlinear Green ’ s Function Method with Viscous Dissipation Effects for Large-Amplitude Roll of Floating Bodies

A novel time-domain body-nonlinear Green’s function method is developed for evaluating large-amplitude roll damping of two-dimensional floating bodies with consideration of viscous dissipation effects. In the method, the instantaneous wetted surface of floating bodies is accurately considered, and the viscous dissipation effects are taken into account based on the “fairly perfect fluid” model. As compared to the method based on the existing inviscid body-nonlinear Green’s function, the newly proposed method can give a more accurate damping coefficient of floating bodies rolling on the free surface with large amplitudes according to the numerical tests and comparison with experimental data for a few cases related to ship hull sections with bilge keels.


Introduction
The seakeeping performance of floating bodies such as ships in harsh sea conditions is one of most important hydrodynamic characteristics that need to be carefully investigated by engineers.Generally, ship motions can be acceptably predicted using inviscid linear potential theories except for roll response.The difficulty of predicting roll response originates from the nonlinearity of roll motion and viscous effects of water.As compared to other motions, the roll motion is inherently highly nonlinear when the roll amplitude is medium to large.The viscous effects always accompany with the roll motion, and they become increasingly significant when ships are fitted with bilge keels.To accurately predict the roll motion, both nonlinear and viscous effects should be considered in the analysis.
As far back as a century ago, William Froude discovered that the roll motion of ships in regular beam waves is nonlinear in certain conditions and should be described using nonlinear equations [1].Systematic studies on the nonlinearity of ship roll motion only started from 1950s, when nonlinear dynamics theories began to be applied to investigate the nonlinear roll characteristics.Typically, the roll response of ships in regular beam waves can be modeled using a single degree of freedom differential equation, which includes linear terms and nonlinear damping and stiffness terms.The nonlinear stiffness terms mainly characterize the shape of the static stability diagram, while the nonlinear damping terms reflect the viscous effects as well as the nonlinearity.The effect of these nonlinear terms is enhanced as the wave amplitude increases.However, it is hard to obtain the exact analytical solution for the nonlinear differential equation with respect to roll response.Therefore, various approximate analytical methods have been developed for solving the nonlinear roll equation to this day, such as the multi-scale method [2], perturbation-based method [3], Krylov-Bogoliubov asymptotic method [4], and so on.
To model and solve the nonlinear differential equation for roll response, the nonlinear damping and stiffness terms must be determined beforehand.There are two approaches to achieve this: one is the parametric identification approach [5], and the other is the non-parametric identification approach [6].In the parametric identification approach, one needs to presuppose the expressions of damping and stiffness, which are generally assumed to be odd-ordered polynomials, e.g., linear plus cubic damping ), and linear plus cubic, quintic stiffness (c 1 θ + c 3 θ 3 + c 3 θ 3 ) [7], where θ, . θ are roll angle and roll angular velocity, respectively, and b i , c i are damping and stiffness coefficients, respectively.However, as in some works, the even-ordered terms might also be adopted, such as linear plus quadric damping [1] or linear plus quadric, cubic, quintic stiffness [8].The nonlinear stiffness could be evaluated using hydrostatic restoring curve fitting method, while nonlinear damping coefficients are difficult to obtain using the classical potential theories.Thereby, most studies that made parametric identification for nonlinear roll models concentrated on investigating the nonlinear damping coefficients.Generally, roll damping can be identified using either the extinction curve method [5] or the energy method [9] from the free roll decay curves.To enhance the numerical stability, some ingenious data processing tools such as support vector regression [10], kalman filter [11], and wavelets [12] are employed in the damping identification.In contrast to the parametric identification approach, which depends on the given nonlinear damping and stiffness models, the non-parametric identification approach does not make any priori assumption of the damping form.A pioneer work using this approach was made by Jang et al. [6], which identified the nonlinear damping by solving the nonlinear Volterra integral equation of the first kind.
All abovementioned studies and methods were based on the differential equation of roll with nonlinear damping models.In fact, in the ship engineering field, sometimes the nonlinear damping is approximated to an equivalent linear one for simplicity.Normally, the equivalent linear roll damping coefficient depends upon not only the roll frequency and forward speed, but also on the roll amplitude.Therefore, it is vital to correctly determine the equivalent linear roll damping coefficient for the large-amplitude roll analysis.
Currently, there are two main ways to obtain the equivalent roll damping coefficient.One is from carrying out ship model tank tests, which can provide reliable resulting data.The other is computational fluid dynamics (CFD) methods, which are preferred due to their lower cost and shorter time consumption.Yeung et al. [13] experimentally and numerically investigated hydrodynamic coefficients of various forced rolling sections and observed that viscous effects were significantly present for sharply edged sections.Korpus and Falzarano [14] compared the roll moment from the potential theory with that from using Reynolds-averaged Navier-Stokes equations (RANS) and found that the shear roll moment was negligible, while the vortex effect was identified to have significant influence on the phase and magnitude.Yıldız et al. [15] demonstrated that the unsteady RANS (URANS) calculations are in a good agreement with experimental results under large amplitudes and shallow draft conditions.Lavrov et al. [16] clearly captured the vortex shedding around a forced rolling hull section with sharp keel by using the RANS solver in OpenFOAM, which suggests that the viscous effects mainly originate from the vortex shedding.Similar results were achieved by using the mesh-free numerical model [17].Actually, the vortex shedding around sharp edges is the major source of viscous effects in the roll damping.
However, CFD methods still require significant computational resources, which might not meet the engineering demand.An alternative approach is to exploit semi-empirical methods.One of the most common semi-empirical methods is Ikeda's estimation method [18], which divides the total equivalent linear roll damping coefficient into skin friction damping, eddy damping, wave damping, lift damping, and bilge keel damping.Amongst these, only wave damping and lift damping can be calculated using potential theories, while the rest are associated with fluid viscous effects and can be obtained using empirical approaches.Beyond Ikeda's estimation method, Brown and Patel [19] incorporated the Direct Vortex Method (DVM) into a potential theory to model viscous effects within an inviscid model, in which the roll center of bilge keels needs to be estimated.Based on the work of Brown and Patel [19], Downie et al. [20] developed a method using DVM to assess the roll damping of ships.
More recently, Guo et al. [21] first proposed a linear time-domain Green's function method (TGF_V) for assessing the interaction between floating bodies and water waves with viscous dissipation effects, which can accurately evaluate the hydrodynamic coefficients of sharply edged floating bodies with small roll amplitude.The TGF_V [21] can consider the full viscous dissipation effects of water, while other similar Green's functions [22][23][24] can only take the partial viscous dissipation effects into account.However, the linear TGF_V or other Green's functions cannot be employed to evaluate the equivalent damping coefficient of large-amplitude roll motion, which is known to be highly nonlinear.
In this paper, a body-nonlinear time-domain Green's function method ( TGF_N _V) is first developed for evaluating the equivalent damping of large-amplitude roll motion for two-dimensional floating bodies with consideration of viscous dissipation effects of water.The TGF_N _V is a significant development of the TGF_N [25] with both nonlinear and viscous dissipation effects of the large-amplitude roll motion into account.
The paper is organized as follows.The mathematical model of the TGF_N _V is developed in Section 2.Then, in Section 3, the TGF_N _V is employed to calculate the equivalent roll damping of hull sections with bilge keels under different roll amplitudes, and the numerical results are compared with experimental ones and those using TGF_N.Finally, the conclusions are drawn in Section 4.

Mathematical Model
A right-hand Cartesian coordinate system o − xy is defined with ox axis placed on the undisturbed water surface and the oy axis orienting positively upward.Figure 1   The differential equation of roll with equivalent linear damping can be written as where a 44 is the added mass, b 44 is the equivalent linear damping coefficient, C(θ) is the stiffness, and is the external roll moment.In this paper, the attention is mainly focused on the evaluation of the equivalent linear damping coefficient b 44 .

Definite Problem of the TGF_N_V
Let ψ(x, y, t) be the disturbed velocity potential of water due to the roll motion of a floating body.Let Ω, S F , S B (t), S ∞ be the water domain, undisturbed water surface, instantaneous wetted surface of the floating body, and flow boundary at infinity, respectively.The definite conditions for the disturbed velocity potential ψ(x, y, t) can be formulated as where is the non-dimensional viscosity coefficient, ω 0 is the reference frequency, and V n (t) is the normal velocity of the instantaneous wetted surface.
From Equation ( 2), one can observe that in the definite conditions of the TGF_N _V, the normal velocity V n (t) is satisfied on the instantaneous wetted surface of the floating body, while the normal velocity in TGF_V [21] is given on the mean wetted surface S B (t) of the floating body.This difference makes the definite problem of the TGF_N _V to be body-nonlinear.

Green's Function for the TGF_N_V
Analogous to the TGF_V [21], the Green's function method is employed to solve the Equation ( 2).Let G(p, t; q, τ; ) be the Green's function, which can be expressed as the combination of instantaneous term G(p, q; ) and free-surface memory term G(p, t; q, τ; ).
where δ(•) is the Dirac function, H(•) is the Heaviside function, p(x, y) is the field point, q(ξ, η) is the source point, ξ, η is the coordinate of q at ox and oy axis, respectively, τ is the time when q pulses, and is the non-dimensional viscosity coefficient.
According to the definition of Green's function, the free-surface memory term should satisfy the following definite conditions: where the subscript q from ∇ q , ∇ 2 q means that the operation is taken with respect to variable q(ξ, η).In the second equation of ( 4), the sign of the viscous term changes to minus due to the relation Obviously, the definite conditions given in Equation ( 4) are the same as those in Guo et al. [21], i.e., they are not related to the body condition.Thereby, according to Guo et al. [21], the expression of the free-surface memory term G(p, t; q, τ; where k 0 = ω 2 0 /g is the reference wave number in deep water conditions.On the other hand, the definite conditions for the instantaneous term G are associated with the free-surface memory term G.The derivation of definite conditions for G is detailed in Appendix A, and the final results are shown in Equation (A8), which are the same as those in Guo et al. [21].According to Appendix B of Guo et al. [21], the unknown constant is β = 1.Therefore, the final definite conditions for G are as follows: where r pq is the distance between field point p and source point q.
Since the definite conditions for the instantaneous term G and the expression of free-surface memory term G are the same as those in Guo et al. [21], the expression of the instantaneous term G must also be the same as that in Guo et al. [21], i.e., where E 1 (z) is the complex exponential integral function defined as

Boundary Integral Equation for the TGF_N_V
According to (A7) in Appendix A, one can obtain the boundary integral equation 2πψ(t, p) + S B (t) G ∂ψ(t,q) ∂n q − ψ(t, q) ∂G ∂n q ds q = t 0 dτ S B (t) G ∂ψ(τ,q) where A and B are water surface-wetted surface of the hull intersections (see Figure 1).The boundary integral Equation ( 10) is different from that in Guo et al. [21].There is an addition term at the right end of Equation ( 10) that reveals the body-nonlinear effect, while there is no similar term in boundary integral equation from Guo et al. [21].
The boundary integral Equation ( 10) is a mixed source and dipole distribution model.One might prefer the pure source distribution model in numerical practice.Let ψ i be the velocity potential of the flow in the interior of the floating body.Through the same procedure, one can obtain the boundary integral equation with respect to ψ i as follows: It is assumed that the following conditions are satisfied on the wetted surface of the hull: σ(τ, q) = ∂ψ(τ,q) ∂n q − ∂ψ i (τ,q) ∂n q , on S B (t) (12) where σ is the source density.
Subtracting Equation (11) from Equation (10) and taking Equation ( 12) into consideration yields 2πψ(t, p) + S B (t) σ(t, q)Gds q = t 0 dτ S B (t) σ(τ, q) Gds q − 1 g t 0 ∂ψ(τ,q) ∂τ According to Duan [25], where n x is the horizontal component of the normal vector on the intersection of the wetted surface.Substituting Equations ( 14) into (13), one gets a boundary integral equation only with source points distributing on the instantaneous wetted body surface.
Taking the derivative of ( 15) with respect to normal vector n p on point p, the source density boundary integral equation is obtained.
where V n (t, p) is the velocity of point p in the outward normal direction.

Equivalent Roll Damping Using the TGF_N_V
The Bernoulli's equation with consideration of viscous dissipation effects of water can be written as ∂ψ ∂t where P is the overall pressure, which can be decomposed into the sum of hydrostatic pressure p 0 and hydrodynamic pressure p P = p 0 + p (18) The hydrostatic pressure p 0 relates only to the immersion depth of the wetted surface, while the roll damping associates only with roll angular velocity, according to Equation (1).Therefore p 0 does not make contribution to the roll damping.
Subtracting Equation (19) from Equation ( 17) yields ∂ψ ∂t Then the roll moment of the floating body due to hydrodynamic pressure is where n c is the generalized normal vector on the wetted surface defined as To enhance the numerical stability, the first term on the right-hand side of Equation ( 21) can be transformed using Reynold's transport theorem [25] S B (t) The calculation of r × ∇ψ and |∇ψ| 2 was not detailed in Duan [25].In this paper we would like to provide a possible approach to calculate them.As shown in Figure 1, V τ is the tangential velocity of water on the wetted surface, which can be obtained through the numerical difference where ∆s is the segment length, i is the segment index for the wetted surface.Now, one can obtain the following expressions: where n d is the dot product of r and n .Substituting Equations ( 23) into ( 21), and then substituting Equations ( 25) and ( 26) into the resulting equation, one gets roll moment Once the time series of the roll moment is obtained, the equivalent roll damping can be calculated using the cosine Fourier transform [16] where T is the period of roll, Im{•} takes the imaginary part of the complex number.The integral in Equation (28) starts from the time nT, when the time series has reached a steady state.The imaginary part of f (t) is selected for calculation due to that the imaginary part of roll (Im{ .θ}) reaches to its maximum angular velocity at time nT and (n + 1)T.
Generally, the quantities related to roll damping need to be made dimensionless for the convenience of use or comparison.In this paper the roll frequency, roll damping, and roll moment are nondimensionalized as follows: where b, d are beam and mean draft of the hull section.

Application of the TGF_N_V for Solving Roll Damping of Hull Sections with Bilge Keels
A novel body-nonlinear Green's function method with viscous dissipation effects ( TGF_N _V) was developed in Section 2. In this section, TGF_N _V is employed to solve the roll damping of hull sections with bilge keels under small and large roll amplitude.For comparison purpose, numerical results from linear inviscid time-domain Green's function method (TGF) [25], body-nonlinear inviscid time-domain Green's function method (TGF_N) [25], and experimental fluid dynamic (EFD) results [13] are also provided in the cases.

Hull Section Under Small Amplitude Roll
A rectangular hull section with bilge keels is shown in Figure 2  As abovementioned, the novel Green's function method TGF_N _V and two contrast methods (TGF and TGF_N [25]) are employed for the roll damping calculation.In each Green's function method, the hull section is discretized into 160 segments (30 for each wall, 20 for each bilge keel, 60 for the flat bottom).In the linear method TGF, all segments are located on the mean wetted surface of the hull section and kept constant.On the other hand, in the nonlinear methods TGF_N and TGF_N _V, all segments move with roll motion, and the segments on two sidewalls need to be remeshed due to the change of draught.
Two viscosity distribution approaches are adopted for TGF_N _V in this case.As shown in Figure 3a, in the first approach a non-zero viscosity is set only on the segments in the vicinity of bilge keels (16 ≤ i ≤ 60, 101 ≤ i ≤ 145), while on the rest segments the viscosity is set to zero.In the second approach (see Figure 3b) a non-zero viscosity is set on all segments (1 ≤ i ≤ 160).The non-dimensional viscosity coefficient should be empirically selected or obtained through the CFD or EFD methods.It is found that desirable results can be achieved if taking the non-dimensional viscosity coefficient as = 0.25 for the first approach and = 0.17 for the second approach.The roll frequency ω is taken as the reference frequency ω 0 in TGF_N _V (see Equation ( 27)). Figure 4 portrays the equivalent roll damping using TGF, TGF_N, and TGF_N _V as compared with EFD.The lines labeled with "TGF_N_V 1 " and "TGF_N_V 2 " denote results from TGF_N _V by using the first and second viscosity distribution approach, respectively.Comparing the TGF_N_V 1 result with the TGF_N_V 2 one, one can find that both of two viscosity distribution approaches are effective for evaluating the roll damping of a hull section with bilge keels.Nonetheless, it can be observed that the second approach can yield smoother results.Thereby, in the next case, only the second viscosity distribution approach will be employed for study.From Figure 4, one can observe that, as compared with EFD data, the linear inviscid method TGF and body-nonlinear inviscid method TGF_N significantly underestimate the roll damping, though the discrepancy between TGF_N result and EFD data is smaller than that between the TGF result and EFD data, which suggests that the numerical results could be improved when the nonlinear effects are taken into account.Moreover, the TGF and TGF_N results decrease when ω > 0.8, which are in the opposite trend as compared to EFD data.On the other hand, the TGF_N _V results agree well with EFD data under the given non-dimensional viscosity coefficients, even when ω > 0.8.
The numerical results in this case suggest that, for the hull section with bilge keels under roll motion of small amplitude, the viscous effects in the roll damping are more important than the nonlinear effects.The study also confirms the possibility of employing TGF_N _V to estimate the viscous effects as well as nonlinear effects of a hull section with bilge keels rolling on the free surface.In contrast, the inviscid methods TGF and TGF_N are not appropriate for evaluating the roll damping of the hull section with bilge keels.As in the previous case, the S60 midship section is discretized into 160 segments (30 for each wall, 20 for each bilge keel, 60 for the flat bottom).It is worth mention that the bilge keels in this case are simplified to lines without thickness, and each side of the bilge keel is discretized into 10 segments.The numerical error may occur when the source point and field point located on the two sides of the bilge keel are overlapped.To avoid this, the value of the Green's function is set as zero once the source point and field point are overlapped.Note that for the body-nonlinear methods, the segments on walls could be stretched or compressed with the roll motion of the S60 midship section.The second viscosity distribution approach, i.e., a non-zero viscosity for all segments on the wetted surface, is employed for TGF_N _V in this case.It is found that for TGF_N _V the non-dimensional viscosity coefficient = 0.2, 0.25, 0.3 for the roll motion with amplitude θ = 0.15, 0.20, 0.25 (rad), respectively, could lead to desirable results.It can be noted that for larger roll amplitudes, it requires larger non-dimensional viscosity coefficient for TGF_N _V.The roll frequency ω is taken as the reference frequency ω 0 in TGF_N _V.

Hull Section Under Large Amplitude Roll
Figure 6 portrays the equivalent roll damping using TGF, TGF_N and TGF_N _V as compared with CFD [15] and EFD [15].One can find that, under all roll amplitudes, the linear inviscid method TGF and body-nonlinear inviscid method TGF_N dramatically underestimate the roll damping as compared with EFD data.Nonetheless, as compared to TGF, the TGF_N results are closer to EFD, which suggests that there exist nonlinear effects in the roll damping.On the other hand, the TGF_N _V results agree well with EFD data, though at high roll frequency TGF_N _V slightly underestimates the equivalent roll damping.In contrast, the CFD results also agree well with EFD data, though it overestimates the damping at high roll frequency.The numerical results in this case suggest that for the hull section with bilge keels under large amplitude roll, the viscous effects in the roll damping are much more important than the nonlinear effects, and one could obtain desirable numerical results only when both of nonlinear and viscous effects are considered.The study also confirms the effectiveness of TGF_N _V in evaluating the equivalent roll damping of a hull section with bilge keels rolling on the free surface with large amplitude.

Conclusions
This paper presents a novel body-nonlinear time-domain Green's function method ( TGF_N _V) for evaluating the equivalent roll damping of floating bodies on the free surface with viscous dissipation effects.The TGF_N _V is derived through the definite problems based on "fairly perfect fluid."The TGF_N _V is different from the existing body-nonlinear time-domain Green's function (TGF_N) [25] or linear time-domain Green's function method with viscous dissipation effects (TGF_V) [21] as follows.As compared to TGF_N, the viscous effects of the fluid can be approximately taken into account.As compared to TGF_V, the nonlinear effects of the large amplitude roll can be precisely considered.
The TGF_N _V is employed to evaluate the equivalent roll damping of hull sections with bilge keels, in which there exist both of viscous and nonlinear effects.The numerical results suggest that the viscous effects are more important than the nonlinear effects in the roll motion, and TGF_N _V can give desirable numerical results as compared with CFD or EFD results.In contrast, those inviscid Green's function methods (TGF or TGF_N) are not appropriate for estimating the equivalent roll damping of hull sections with bilge keels.∂ G ∂n q − G ∂ψ(τ, q) ∂n q ds q = 0 (A2) Using the boundary conditions for G on S ∞ , it can be deduced that the integral (A2) on S ∞ equals to 0. Using the initial and free-surface conditions for ψ(τ, q) and G(p, t; q, τ; ), the integral (A2) on S F can be transformed to t 0 dτ S F ψ ∂ G ∂n q − G ∂ψ ∂n q ds q = − 1 ξ(τ) is the horizontal velocity of intersection between water surface and the floating body, and A and B are the intersections located on the left and right side of the floating body (see Figure 1), respectively.
Substituting Equations (A3) into (A2) yields On the other hand, substituting Equations (A2) into (A1) and taking Equation (3) into account, one obtains t 0 dτ S F +S B +S ∞ ψ(τ, q) ∂ δ(t − τ)G ∂n q − δ(t − τ)G ∂ψ(τ, q) ∂n q ds q = 2πβψ(t, q) (A5) i.e., S B ψ(t, q) ∂G ∂n q − G ∂ψ(t, q) ∂n q ds q − 2πβψ(t, q) = S F +S ∞ G ∂ψ(t, q) ∂n q − ψ(t, q) ∂G ∂n q ds q (A6) Subtracting Equation (A4) from Equation (A6) yields ψ ∂G ∂n q − G ∂ψ ∂n q ds q − t 0 dτ S B ψ ∂ G ∂n q − G ∂ψ ∂n q ds q − 2πβψ If the right hand side of (A7) equals to 0, i.e., G = ∂G ∂n q + 1 g ∂ G ∂τ τ = t = 0 on S F , and G ∼ O 1/r pq on S ∞ , the singular points can be distributed only on the wetted surface S B .Therefore, the definite conditions for the instantaneous term G(p, q; ) can be concluded as following: As compared with Guo et al. [21], one can find that there is an addition term related to the water surface-floating body intersection in the relationship between G and G in Equation (A7).Nonetheless, the definite conditions for the instantaneous term G are the same as those in Guo et al. [21].
portrays the coordinate system o − xy and a hull section with bilge keels.

Figure 1 .
Figure 1.The coordinate system o − xy and a hull section with bilge keels.The hull section harmonically rolls on the water surface.A and B are water surface-hull wetted surface intersections.The wetted surface of the hull section is discretized into s segments.
. The beam and mean draught of the hull section are b = 0.4 m and d = 0.2 m, respectively.Two bilge keels are set at the corners that connect the wall and bottom of the hull section, and the bilge radius and bilge keel depth are R b = 0.02b and K d = 0.04b, respectively.The hull section harmonically rolls on the water surface with small amplitude θ = 0.05 rad ≈ 2.87 • .The axis of rotation is located at the center of waterplane of the hull section.The detailed description of the experiments can be found in Yeung et al. [13].

Figure 2 .
Figure 2. A rectangular hull section with bilge keels harmonically rolls on the water surface, with beam b = 0.4 m, mean draught d = 0.2 m, bilge radius R b = 0.02b, bilge keel depth K d = 0.04b, and roll amplitude θ = 0.05 rad.The rotation axis is located at the center of waterplane.

Figure 5
Figure 5 depicts an S60 midship section with bilge keels.The beam and mean draught of the hull section are b = 0.237 m and d = 0.057 m, respectively.The bilge radius and bilge keel depth are R b = 0.035 m and K d = 0.01 m, respectively.The single degree of freedom model experiment was carried out by Yıldız et al. [15], in which the S60 midship section model is fixed by a forced rolling device, and the roll moment is measured through the torque transducer.The S60 midship section model harmonically rolls on the water surface with large amplitudes θ = 0.15 rad (8.59 • ), 0.20 rad (11.6 • ), 0.25 rad (14.4 • ).The axis of rotation is located at the center of waterplane of the hull section.More detailed description on the model experiments can be found in Yıldız et al.[15].

Figure 5 .
Figure 5.A semi-submerged S60 midship section with bilge keels harmonically rolls in the presence of a free surface; beam b = 0.237 m, mean draught d = 0.057 m, bilge radius R b = 0.035 m, bilge keel depth K d = 0.01 m, roll amplitudes θ = 0.15, 0.20, 0.25 (rad).The rotation axis is located at the center of the waterplane.