THE APPLICATION OF MESHLESS CYLINDER CONTROL SURFACE IN RANKINE-KELVIN HYBRID METHOD UDC 629

The solution of the three-dimensional seakeeping problem with forward speed in the frequency domain still has some well-known problems. In this paper, a Rankine-Kelvin hybrid method benefiting the merits of both the Rankine source and Kelvin source is presented and demonstrated by studying the wave diffraction/radiation problem with zero forward speed as an example. A meshless cylindrical surface is selected to be the control surface dividing the fluid domain into two regions, and the velocity potential and its normal derivative on the control surface are represented by series expansions. The present RankineKelvin hybrid method is validated by the added mass and damping associated with the linear radiation force, and comparison is made with the documented analytical solution. All the work in the paper sheds light on solving the forward speed hydrodynamic problems.


Introduction
As it is known, when solving the hydrodynamic problem with forward speed, the implementation of the Kelvin source method using the Green function associated with a translating and oscillating source bring along several well-known problems concerning its wave component.In particular, it behaves with complex singularities and high oscillations when both the field and source points approach to the free surface [1].In the earlier work, numerical results were in very good agreement for immersed bodies like sphere and ellipsoid, or surface piecing body of very slender shape like Wigley hull [2].Unfortunately, for the realistic container ships, the results present large discrepancy from those of model tests.Ten and Chen (2010) [3] came up with a new method to solve the forward speed problem, and Hui Li, Lizhu Hao The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Huilong Ren, Bo Tian Hybrid Method 88 they try to benefit the merits of the Rankine source method and the Kelvin source method by using a hemisphere as a control surface.However, the method made the calculation of the integration over the control surface complex and may lead to some troubles about singularities.
In this work, we use a cylinder as a control surface instead of a semi-sphere in the same spirit as Liang and Chen's work [4].The Rankine source method is used in the internal domain, and the fundamental solution is in the form of 1 r , with r standing for the distance between the field point and the source point.Panels are distributed over the body surface and the free surface in the internal domain.Kelvin source method is used in the external domain and the Green function which satisfies the free-surface boundary condition and radiation condition is adopted.On the meshless cylindrical control surface, the velocity potential and its normal derivative are expanded into Fourier-Laguerre series and Lv et al. [5] have proved the effectiveness of approximation by Laguerre series for diffraction waves on the cylinder.
The objective of this work is to solve the zero speed seakeeping problem as an example to show the effectiveness of Rankine-Kelvin hybrid method with a meshless cylindrical control surface.In the numerical implementation, we consider a hemisphere floating at the free surface, and the forward speed is out ruled at the moment.The added mass and damping coefficients are calculated with the present Rankine-Kelvin hybrid method, and the convergence test associated with the terms of Laguerre functions and Fourier series is made.In addition, the influence of different control surface radius on the numerical results is studied through comparing with the analytical solution, and eventually a satisfactory accuracy is obtained.Due to the fact that there are many new defined coefficients concerning multi-fold integrals in this method, Chebyshev expansions are utilized to approximate the resultant integral to improve computational efficiency.The successful application of the Rankine-Kelvin hybrid method in zero speed problem has established a foundation for solving the forward speed problem.

The series expansion method of velocity potential and its normal derivative
In the Rankine-Kelvin hybrid method, as a meshless circular cylinder surface is selected as control surface, the velocity potential and its normal derivative on the control surface need to be expressed analytically.In this paper, they are expanded into Fourier-Laguerre series.
The velocity potential on the control surface is a function associated with the polar angle  and vertical coordinate z positively downward, which can be expressed as   , z  .
As the velocity potential   , z  decreases exponentially with the increase of z, and when z tends to infinity, the velocity potential tends to zero.The velocity potential is expanded into series of Laguerre as follows The orthogonal property of Laguerre polynomials is considered in Eq. ( 2) and Eq.(3).
The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Hui Li, Lizhu Hao Hybrid Method Huilong Ren, Bo Tian 89 where nk  is Kronecker delta function.
Thus,   n z is multiplied on the two sides of the Eq. ( 1) simultaneously and we can integrate the both sides from zero to infinity with respect to z , and the coefficient of the series expansion can be expressed as follows As the velocity potential on the control surface is a periodic function about  of which the period is The substitution Eq. ( 5) into Eq.(1) yields The substitution Eq. (4) into Eq.( 6) yields In the same way, the normal derivative of velocity potential   , z  of an arbitrary point on the control surface can be also expanded into Laguerre-Fourier series as follows In (7) and (9), nm  and nm  are coefficients of the series expansion and defined by (8) and (10), respectively.

Hui Li, Lizhu Hao
The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Huilong Ren, Bo Tian Hybrid Method 90 3.The application of rankine-kelvin hybrid method to zero speed problem

Definition of parameters
We consider the ship floating in the sea of infinite depth with zero speed and use a cylinder to separate the whole fluid domain into two domains.In the external domain, the normal vectors of the boundary surfaces positively inward to the domain.In the internal domain, the normal vectors point positively outwards of the domain.The normal vectors and the coordinate system are shown in Figure 1.It is convenient to use the cylindrical coordinate system as a cylindrical surface is chosen to be the control surface.Here is the transformation between the two coordinates.
where h is the radius of the cylinder, The boundary conditions in the internal domain are given as followings where  denotes the oscillation frequency.

External problem
In the external domain, the Green function used here is given in the form of [6] is the wave number, 0 () is the zeroth order Bessel function of the first kind, .. PV is principle-value integral.
For the velocity potential at an arbitrary field point P on the control surface in the external domain, application of the Green's second identity provides Application of Eq. ( 8), we have We define new coefficients So Eq. ( 17) can be rewritten as ,, 00 As the subscripts n and m change, we can get a system of equations which can be expressed as follow

 
, Finally we can get the relationship between the series expansion coefficients of the velocity potential and its normal derivative on the control surface.
And we define the matrix DN as As the application of the free-surface Green function in the external domain problem will result in the occurrence of irregular frequencies, extended boundary integral equation method presented in [7] is adopted to remove irregular frequencies, and the interior free surface is divided into panels shown in Figure 2. The boundary integral equations are given as follows The dipole distribution  on each panel of the interior free surface is assumed constant.

The validation of the external domain
To verify the relationship between the velocity potential and its normal derivative obtained from the external domain, we have calculated the diffraction potential of an infinitely long vertical circular cylinder.
We consider a regular wave as the incident wave, which propagates along the positive axis of x , and the corresponding first order incident wave potential is as follows Among them, A is the amplitude of the incident wave, and   ,, Rz  is the cylindrical coordinates of an arbitrary point.According to the body surface condition, we can get Applying Eqs.( 9) and (10), the series expansion coefficient of the normal derivative of diffraction potential can be obtained Appendix 1 can be referred for the calculation details of Eq. (37) Substitute Eq. (37) into Eq.( 22) and (33), we can get the series expansion coefficients of the diffraction potential, then numerical solution can be obtained based on Eq. (7).
At the same time, the diffraction force of an infinitely long vertical circular cylinder d wj F can be calculated based on the following Eq.0 Analytical solution of the diffraction potential [6] is as follow Substitute Eq. (42) into Eq.( 39), analytical solution to the diffraction force can be expressed as follow Comparison of the numerical solutions for the diffraction potential and diffraction force with the analytical ones are presented in Section 4.

Internal problem
In the internal domain, Rankine source is adopted which is given as follow The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Hui Li, Lizhu Hao Hybrid Method Huilong Ren, Bo Tian 95 For the velocity potential at an arbitrary field point P on the boundary surface of internal domain consisting of control surface C S free surface F S and body surface B S , application of the Green's second identity yields The free surface and body surface are divided by B N and F N panels, respectively.Thus, expression (45) becomes Substitute the boundary conditions Eqs. ( 12) and (13) into Eq.( 46), where n  is the component of the normal vector on B S .In Case 1, where the field point P is on the control surface

ds n G n ds G ds n G ds G ds gn
,, 0 Eq. ( 49) can then be rewritten as follows In a similar way with Eq. ( 21), we can get a matrix form as follows ,, In Case 2, when the field point P is on the body surface B S or on the free surface F S , the panel including the field point P is numbered by  , then Eq. (47) becomes Then Eq. ( 54) can be rewritten as follows , , In a similar way with Eq. ( 21), we can get a matrix form as follows Here, E is a unit matrix with dimensions of   , Substitution of Eq. (33) into Eq.(57) yields Hui Li, Lizhu Hao The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Huilong Ren, Bo Tian Hybrid Method

98
Then the equations for the entire domain can be obtained by combining Eq. ( 53) and ( Once Eq. ( 59) is solved, one can obtain the velocity potential on the panels and the series expansion coefficients of the velocity potential on the control surface.The method for evaluating the coefficients we have defined for the multi-fold integrals is given in the appendices.From the expression of the multi-fold integrals, the Chebyshev expansion is used to approximate the integrations and improve the efficiency.

The solution of added mass and damping coefficient
The added mass and damping coefficient can be solved as follows where, ij a is added mass, ij b is the damping coefficient, as the body surface is discretized into panels, Eq. ( 60) can be rewritten as Then added mass and damping coefficients can be obtained by substituting the velocity potential on the body surface into Eq.(61).

Results for the diffraction potential
Numerical solution has been computed with the order of Laguerre function from 0 to 10 and the order of Fourier series from -10 to 10, the radius of the cylinder is 3.0 m, a fixed point P (3.0 1.5 1.0) is selected expressed in cylindrical coordinates on the circular cylinder, the diffraction potentials at the fixed point varying with wavenumber are shown in Figure 3, among them, method 1 refers to the original solution method of the exterior domain in subsection 3.2, Method 2 refers to the extended boundary integral equation method in subsection 3.3.
Fig. 3 The comparison results of diffraction potentials In Figure 3, the numerical solution agrees well with analytical solution for wavenumber varying from 0.1 to 1.6 except for 0  is 0.8.
As we can see, numerical solution has a fluctuation when wavenumber is 0.8, and this wave number is regarded as irregular frequency.Different points on the circular cylinder are chosen with h =3.0 and z =1.0.The diffraction potentials at the irregular frequency varying with different circumferential locations of the points are shown in Figure 4 - The comparison results of diffraction potentials at irregular frequency From Figure 4, we can see that the numerical solution obtained from extended boundary integral equations in method 2 is in good agreement with analytical solution, which illustrates the relationship between the velocity potential and its normal derivative on control surface is right and extended boundary integral equation method is capable of removing the irregular frequency.Diffraction force in the direction of x axis is given in Figure 5, and the numerical solution shows a satisfied accuracy.

The results of added mass and damping coefficients of a hemisphere
A hemisphere is chosen as the example for calculation, the numerical solutions about the added mass and damping coefficient are compared with the analytical solutions given by Hulme [8].For the convenience, we have defined the following conditions shown in table 1.The convergence test associated with the order of Fourier-Laguerre series has been made through condition 1, 2 and 3.The results are shown in the following figures, among them, L is the characteristic length which is equal to r the radius of hemisphere.From the results, we can see that the numerical results are in good agreement with the analytical solutions, in Figure 6 and Figure 8.The results of surge motion seem sensitive to the order of Fourier-Laguerre series at a high wave frequency, and the discrepancy between analytical solutions and numerical results in condition 1 are larger than that in condition 2 and condition 3. From Figure 7 and Figure 9, the results of heave added mass and damping coefficients have shown a good precision in condition 1.In conclusion, the numerical results show a good convergence with the increase of the order of Fourier-Laguerre series.
At the same time, the influence of different control surface radius to numerical results has also been studied through condition 3, 4 and 5.The results are shown in the following figures.From the results we can see, added mass and damping coefficients of surge motion are more sensitive to control surface radius than that of heave motion.In Figure 11 and Figure 13, the numerical results show a good agreement with the analytical solution in different conditions.In Figure 10 and Figure 12, the influence of the control surface radius to the numerical result is small at a low wave frequency, and with the increase of the wavenumber, the numerical results in condition 3 are more accurate than the other conditions overall.The main reason is more panels are needed on the free surface with the increase of the control surface radius.Consequently, considering the chosen of control surface radius, we prefer to a smaller one.

Conclusion
The Rankine-Kelvin hybrid method has been applied to solve zero speed seakeeping problems successfully.From this paper, we can draw the following conclusions 1.The application of Kelvin source in the external domain will result in irregular frequency, and an extended boundary integral equation method has been used to eliminate the irregular frequencies.2. The result is convergent with the increase of the order of Fourier-Laguerre series, and a better accuracy can be achieved when a smaller radius of the control surface is chosen with the same order of Fourier-Laguerre series.3.All the work in this paper has laid the foundation for solving the forward speed hydrodynamic problem.In this part, some integrals about Laguerre function will be calculated following his notes [9].
1 The calculation about the integration of Laguerre function 1.1 Eq. ( 1.1) can be transformed into recursion Eq. (1.2) Change the order of Laguerre polynomial from 1 k  into k and reuse Eq. ( 1.2), we can get the following result Finally, the integration about Laguerre function can be solved referring to Eq. (1.3) Eq. (1.5) can be transformed into recursion Eq. (1.6) Finally, the integration about Laguerre function can be solved referring to Eq. (1.6) When the integration is computed near Then we have


We use the result of Eq. (1.4), then Referring to Eq. ( 2), the analytical solution of Eq. (1.11) can be obtained as follow This part introduces the expression of Green function in the cylindrical coordinate system.
Kelvin source can be written as Eq.(2.1) r G is the Rankine source which is associated with the distance between the field point ( , , ) P x y z and the source point Q( , , ) 2) can be rewritten as In the cylindrical coordinate system, the field point P and the source point Q can be written as ( , , ) We can get The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Hui Li, Lizhu Hao Hybrid Method Huilong Ren, Bo Tian 107 In a similar way

APPENDIX 3
The computational methods of the integrals defined in the external domain.
, nm kl G can be divided into the following four parts Then we can calculate the four parts separately with the formulas mentioned in appendixes 1 and 2. Substitute Eq. (2.5) into Eq.(3.1), use the orthogonally of Fourier series, we can get Eq.(3.5) To make the calculation of numerical method accurate and efficient, we can deal with Eq. (3.9) as follows The normal derivative of Kelvin source on the control surface can be written as The result can be written as Hui Li, Lizhu Hao The Application of Meshless Cylinder Control Surface in Rankine-Kelvin Huilong Ren, Bo Tian Hybrid Method 110 In the extended boundary integral equation method, when the field point is on the control surface and the source point on the interior free surface, the normal derivative of Kelvin source on the free surface can be written as H defined in Eq. ( 27) can be divided into four parts , , When the field point is on the free surface and the source point on the control surface, , F kl G defined in Eq. ( 30) is divided into four parts , , where S   represents the area of panel on the body surface.Substitute the result of 1.1 and 1.2 in Appendix 1 into Eq.(4.1), then ,, B nm G  can be solved for.The normal derivative of Kelvin source on the body surface can be written as  When the field point is on the free surface and the source point is the control surface

Fig. 1
Fig. 1 Definition of the coordinate system and normal vectors

Fig. 2
Fig. 2 Schematic diagram of extended boundary integral equation method

Fig. 5
Fig. 5 The comparison results of diffraction force

Fig. 6 7 Fig. 8 9
Fig.6 Surge added mass coefficients Fig.7Heave added mass coefficients with different order of Fourier-Laguerre series with different order of Fourier-Laguerre series

Fig. 10 Fig. 12
Fig. 10 Surge added mass coefficients Fig. 11 Heave added mass coefficients for different control surface radius for different control surface radius Eq. (2.8) into (3.4),we can get

5 )
When the field point is on the body surface and the source point on the control surface

Table 1
When both the field point and source point are on the interior free surface, similar to the way of dealing with The calculation of the integrals defined in Eqs.(50) and (55)., when the source point Q is on the body surface or on the free surface, the source on each panel is constant.