VIRTUAL SIMULATION OF PMM TESTS INDEPENDENT OF TEST PARAMETERS

Summary The dynamic planar motion mechanism (PMM) tests are simulated numerically using computational fluid dynamics for a prolate spheroid underwater vehicle (PSUV) to find the effects of test parameters: the amplitude, frequency and flow velocity, and make the simulation independent of them. An amplitude of the sinusoidal path of the pure sway and heave tests less than 0.01L and a frequency less than 0.03 Hz are necessary to find accurate results for the maneuvering hydrodynamic derivatives. A ratio of angular frequency to the flow velocity equal to one and an amplitude of less than 0.03L provide relatively accurate results for pure yaw and pitch tests. The calculated test parameters are validated using them for the simulation of the PMM tests for two SUBOFF submarine models to control and compare with the experimental results.


Introduction
Manoeuvrability of a marine vehicle is usually obtained by course-keeping, coursechange and speed-changing abilities. Various mathematical models have been developed to estimate the manoeuvrability of a surface or an underwater vehicle at different design stages, including those based on derivative models or hydrodynamic coefficients. These models have been developed based on the original Abkowitz [1] model and the use of the Taylor series extension. Determination of the hydrodynamic coefficients or derivatives is a prerequisite to predict the vessel manoeuvrability using these models. The hydrodynamic coefficients in these models are divided into three categories: restoring coefficients, damping coefficients and added mass coefficients. Restoring coefficients depend on the position of the body and are usually due to the interaction of the gravity and the buoyancy forces. Damping coefficients are the function of the body velocity and are created due to the fluid viscosity and formation of waves on the free surface of the water. Added mass coefficients have resulted from body motion acceleration. The fluid particles are moving with a variable speed due to the acceleration of the body and this causes the energy dissipation of the body.
There are three different approaches to find the damping and added mass coefficients for manoeuvrability of marine vehicles, analytical, experimental and numerical approaches. The analytical approach is mostly based on the potential flow theory and can be applied for simple geometrical models. The analytical results are obtained by solution of boundary value problems governing by the Laplace equation subjected to different sets of linearized boundary conditions. The experimental methods include execution of a set of tests on a down-scaling model of the body at a towing tank or a manoeuvring basin. Conventional experiments are known as "captive model tests" which include rotating arm (RA) tests and planar motion mechanism (PMM) tests. Rotating arm tests are carried out in a maneuvering basin and are used to find the damping coefficients. Planar motion mechanism tests are done in a towing tank and are used to obtain all damping and added mass coefficients. A huge marine laboratory facilities is required to obtained high-precision coefficients with minimum scale effects [2].
The PMM tests may be done in a straight path when the model has a drift angle with the fluid flow. Such a test is called as Oblique Towing Test (OTT) and provides the damping coefficients depending on the translational velocities. The oblique towing test is also called as the static test since the model is running with a constant velocity in a straight path and therefore, there is no acceleration involved. The PMM tests may also be done in a sinusoidal path with various orientation of the body. These types of tests are also called as dynamic tests since the body is acted by inertia forces. The dynamic tests can provide the acceleration and velocity dependent hydrodynamic coefficients. Two of the most important PMM tests are pure sway and pure yaw tests. These types of model tests need special equipment and are expensive, time consuming and their results include the scaling effects due to inconsistency of Reynolds number between the ship and the model, Hajivand and Mousavizadegan [2].
There are two approaches to find damping and added mass coefficients using numerical calculations. The first approach is based on the potential flow theory. It provides all added mass coefficients but cannot predict the damping coefficients completely. This approach may predict the part of the damping coefficients due to the formation of waves on the free surface of water by the motion of the body. The second approach is based on the real fluid flow theory and virtual simulation of RA or PMM tests using computational fluid dynamics (CFD). Nowadays, CFD is a crucial tool for various aspects of marine vehicle hydrodynamics such as resistance, propeller performance and manoeuvrability, not only for research but also as a design tool [3]. Virtual simulation of different types of PMM tests, including oblique towing tests, pure sway and pure yaw tests, provides all necessary coefficients for the maneuvering mathematical model equations to predict the marine vehicle manoeuvres.
Zhang et al [4] simulated PMM test to obtain all of hydrodynamic coefficients of an underwater vehicle using Fluent. They didn't present any validation for calculations. Lee et al [5] obtained added mass and inertia coefficient of an unmanned underwater vehicle (UUV) by simulation of PMM test using CFD analysis. They validated the CFD calculation with experimental data and concluded that virtual simulation of PMM test can provided the hydrodynamic coefficients with a reasonable accuracy. Yu-cun et al. [6]calculated the hydrodynamic coefficients of a SUBOFF submarine by numerical simulations. They used Reynolds Averaged Navier-Stokes (RANS) simulations to simulate the static and dynamic PMM tests of the SUBOFF submarine model. The computational results were verified by comparison with experimental data, which show that this method can be used to estimate the hydrodynamic derivatives of a fully appended submarine. It is indicated that the OTT tests is more effective to find damping coefficients. There are also some literatures address the computation of hydrodynamic coefficients for surface ships such, Stern et al. [7], Simonsen et al. [8] Hajivand and Mousavizadegan [9], [10].
The frequency and the amplitude of the pure sway (heave) and pure yaw (pitch) tests may call as the laboratory parameters. The values of these parameters have considerable effects on the accuracy of the results. A quite simple method is developed to make the simulation of these dynamic PMM tests, independent from amplitude and frequency of the test. The method is developed based on the simulation of dynamic PMM tests on a Prolate Spheroid Underwater Vehicle (PSUV). The prolate spheroid is widely used to estimate the hydrodynamic characteristics of symmetric underwater vehicles such as torpedoes, submarines [11]. The effects of amplitude and frequency on the results have been studied and analysed, then the relationships for the frequency and amplitude of the motion are derived. The aspect ratio of the prolate spheroid is selected equal to = 1 6 to have an elongated axisymmetric body to resemble an underwater vehicle. The notations and are the polar and equatorial radiuses of the prolate spheroid, respectively.
The derived relationships are applied for two SUBOFF models to validate them. The SUBOFF submarine model is designed and tested by David Taylor [12], [13]. In order to validate, the fully appended model has been simulated as well as the model with 4 hydroplanes, as shown in Fig. 1. The total length of the SUBOFF model is 4.356 , the distance between forward and the aft perpendicular is 4.261 and the maximum diameter of the hull is 0.508 . The center of forces calculation is measured on the center line of the body and at a distance of 2.013 from forward perpendicular. The positive direction of the coordinate axes, the x-axis is in the direction of the head, the y-axis in the left direction of the body and the z-axis is considered to be the downward direction. The second coordinate system is the inertial coordinates system, which is used to define the translation and rotation motions of the body-fixed coordinates system to the earth-fixed coordinates as shown in Fig. 1. In this system, the vehicle location is expressed in ξ, η, ζ coordinates after being defined in the body-fixed reference frame. The orientations of the body-fixed references frame are expressed by Euler angles (yaw), θ (pitch) and φ (roll).

Dynamic maneuvering equations of a marine vehicle and PMM test
The dynamic equations of a marine vehicle with six degrees of freedom may be written as follows [14], [15]: Where ν = {u, v, w, p, q, r} T denotes the linear and angular velocity vector with coordinates in the body-fixed frames and denotes the position vector with coordinates in the earth-fixed frame. The notation M RB is the inertia matrix of the body and its elements depend on the mass and mass distribution of the body and the position of reference frame system. The notation M A is the added mass matrix and its elements depend on the inertia of the surrounding fluid. Based on the Newton second law, it is expected that the inertia terms from the fluid to be linearly dependent on accelerations. Using Taylor's series expansion, the elements of M A , such as X u̇, Y v̇, N ṙ and etc., may be expressed as the first derivative of hydrodynamic forces and moments with respect to the acceleration in initial conditions that the body is moving with a constant speed along a straight course in longitudinal direction. The notation ( ) is the Coriolis and centripetal matrix of the body and depend on the mass and mass distribution of the body. The notation C A is the hydrodynamic Coriolis and centripetal matrix and its elements is expressed using the added mass coefficients of the body.
The symbol D(ν) is the hydrodynamic damping matrix. The hydrodynamic damping for a body is caused by wave formation, skin friction, separation of fluid and vortex shedding. It may be expressed as D = D l + D nl , where D l and D nl are linear and non-linear damping matrices. Using Taylor's series expansion, the elements of D l and D nl , such as , | | , , | | , , | | and etc., may be expressed as the first and higher order derivatives of hydrodynamic forces and moments with respect to the velocities in initial conditions that the body is moving with a constant speed in along a straight course in longitudinal direction. The notation G(η) is the restoring forces and moments vector, E and are the forces and moments of the environment and control and propulsion systems, respectively. A complete explanation of equation (1) and the relationships for the notations may be found in Fossen [14], [15]and Khac Duc Do and Jie Pan [16].
The elements of matrices M A and D are called as the hydrodynamic coefficients or derivatives. Determination of hydrodynamic coefficients is a prerequisite to solve the equation (1) and obtain the manoeuvrability characteristics of a vessel. These coefficients may be obtained by dynamic PMM tests which are pure sway, pure yaw, pure heave and pure pitch tests. Descriptions of these tests are summarized in Table 1 by providing schematic illustrations, assumptions, the relationships and explanation of the symbols of the tests.

Virtual Simulation Of PMM Tests Independent Of Test Parameters
Ardeshiri, Mousavizadegan, Kheradmand.  • Subscripts in and out are relation to the force and moment that are in and out of translation phase, respectively.

•
The parameter 0 is the amplitude of the movement, is the oscillation frequency of the arms.

•
The notations and are the forces, and are the velocities in and directions, respectively.
• The parameters and are the moments around and axes, respectively.

Governing equations of fluid flow modelling
The flow of water around a moving underwater vehicle may be modelled by the Navier Stokes equations which can be used for both laminar and turbulent flow. It is necessary to use very fine meshing to capture all the turbulence effects for a turbulent flow regime using the Navier Stokes equations. Alternatively, the Reynolds-averaged Navier Stokes (RANS) equations may also be applied to model the turbulent flow. The RANS equations are obtained based on statistical tools known as Reynolds decomposition where the flow parameters are decomposed into time-averaged and fluctuation parts, i.e. u = u ̅ + u ′ and p = p ̅ + p ′ where u ̅ and p ̅ are the time-averaged and u′ and p′ are the fluctuation velocity and pressure, respectively.
The Reynolds-averaged Navier-Stokes (RANS) equations are applied to model the water flow around the bodies to simulate the PMM tests. The RANS equations may be given as follows for an incompressible flow (Ferziger and Peric, 2002).
Where μ is the fluid dynamic viscosity of the fluid, ρ and u i define the density and velocity component of i (i = 1, 2, 3), respectively. The notation g i is the components of gravitational acceleration, −ρu i ′ u j ′ = τ ij is the Reynolds stress tensor and shows the turbulence effects. For modeling the turbulence effects in RANS equations, various models have been developed as turbulence models. The k-ε Realizable turbulence model is used to model the turbulence effects and find the Reynolds stress tensor, [17] and [18].

Numerical implementation
The RANS equations are numerically solved by several methods such as Finite Difference, Finite Volume and Finite element method [19]. Finite volume method (FVM) is the most widely used method. A finite number of control volumes that are called cells in FVM generates in computational fluid domain. The conservation laws that are expressed as RANS equation should be satisfied for the cell. The Star CCM+ software which is developed based on FVM, is used to do the computations. The collocated point, Considered at the center of each cell, all field variables except velocity are calculated at this point, which makes the fields of velocity and pressure independent from each other. In order to solve this problem, the Rhie and Chow interpolation algorithm [20] is used. The SIMPLE (semi implicit method for pressure-linked equations) algorithm is used to couple the pressure and velocity. The advection terms are discretized by a high-resolution numerical scheme [21]. A linear interpolation scheme is used for interpolating the pressure, while the velocity is interpolated using a trilinear numerical scheme. The Star CCM+ manual [22] is followed and the problem is solved in the time domain until a steady state is reached.

Domain and boundary conditions
The overset meshing is used to simulate the sinusoidal motion of the bodies in dynamic PMM tests. The grid components can have any size, scale, and shape, but should be in such a way to interact with other components and to create complete overlap on the computing domain. In addition, the domain of computation should be selected large enough to satisfy the immersed body conditions and considering the effect of the simulation runtime. In other word, the computation domain dimension should be selected so that the next period of motion of body can also be verified. Accordingly, the computational domain dimensions depend on the fluid velocity 0 , the body motion amplitude a, the body motion frequency and the required time, , for simulation. However, the boundary condition should also be satisfied. According to Gao et al [23], the computational domain dimensions are selected as follows. The domain of computation and the notations are depicted in Fig. 2. The boundary conditions are also shown in Fig.2 and chosen such that the opposite wall of the body moving and the top and bottom walls, have the velocity inlet boundary condition. The pressure outlet boundary condition is for upstream boundary, and non-slip boundary condition is applied to the immersed body in the calculations.

Mesh generation
The structured trimmer meshing model is used to discretize the computational domain. This method is an optimal and powerful way to produce structured mesh for a grid with minimal symmetry. Fig. 3 shows the mesh around the PSUV and the value of y+ on the body to simulate the pure yaw test. In the generation of mesh elements, it is considered that the minimum number of elements satisfying the conditions of residuals less than 10 −5 and + = 30. The requirement of overset mesh generation to prevent negative finite volume errors is also controlled during the computations. The Roach method [ [24]- [25]- [26]] was used with the GCI (Grid Convergence Index) criterion to control the independence of the results from the number of grids. In each simulation, the three key parameters (f ) are the force along the x − axis, , the force along the or z − axis ( in the horizontal plane or in the vertical plane) and the moment around the local coordinates reference ( in the horizontal plane or on the vertical plane) at 4° angle for analysing the sensitivity of the results to the gridding in the simulations. The results of PSUV pure sway test simulation are presented in Table 2. The following relationships provide the steps to control the criterion for each parameter. Where fi is the key parameter, n is the number of cells, r is the ratio of the cells number, D = 3 is for 3D geometry, ε is the relative error of the key parameter between two levels of grid, F s = 1.25 is the safety factor. If the parameter α is close to one then will be independent of grid. The notation is the apparent order of the method, which is calculated for the three levels of the mesh as follows: Which in the above equations is φ 12 = f 1 − f 2 and φ 23 = f 2 − f 3 .

Development of test parameters
The manoeuvring coefficients for sway, heave, pitch and yaw motions may be obtained by simulation of the dynamic PMM tests which are pure sway, pure heave, pure pitch and pure yaw tests. These tests with the relevant parameters and formula are illustrated in Table  1. The accuracy of the result depends on the some test parameters which are the amplitude and the frequency of the motion in each test. The effects of variation of these parameters are studied on a prolate spheroid underwater vehicle (PSUV) by variation of these parameters for the pure sway and yaw motion. Since the PSUV have a symmetrical shape with its longitudinal axis, the results for the sway and heave, and pitch and yaw are identical. The PSUV is selected due to the fact that their analytical solutions of the added mass coefficients are available and can be used for validation and development of the test results, Korotkin However, it cannot be zero due to the fact that the motion will be steady and there is no added mass involved. The value of 0.012 for is good enough and the errors is less than 10%. This value is equal to 0.01 where is the length of the prolate spheroid.  Fig. 5. The analytical solution is ̇′ (Theory) = -0.02668 and the numerical solution for ̇′ is increased with increasing the frequency of the body motion. The velocity of body has also a detrimental effect on the accuracy of the solution. It may be due to the influence of the boundary layer and phenomena such as flow separation and their contribution to calculation, and thus the probability of bias from the theoretical value has been increased. The numerical computations show that the fluid flow and the frequency should be as low as possible to find an accurate result. If it is assumed that an error less than 10% is acceptable, the body motion frequency should be less than 0.2 to guarantee the desired accuracy regardless of the velocity of body in longitudinal direction. However, taking into account the computational runtime and cost within the range of acceptable error less than 10%, the value of 0.2 for the body motion frequency is more acceptable. The achievement of the defined residual values (less than 10 −5 ) in the time of flow solving cycles with the increasing velocity of the body was achieved over a longer time, which, in turn, increased the runtime. The body motion in the pure sway test is shown in Fig. 6. The body is moved with a speed of 1 / , the amplitude of 0.12 and the frequency0.2 . The vorticity field around the body is shown in the simulation of the motion of a pure sway at a speed of 1 m/s. We can see the time history of the wake for pure sway test. The wake flow behind the body is also shows a sinusoidal motion due to the motion of the bodies. The wake after the body fades after about two cycles behind the body.

5.2
The effect of the test parameters for the pure yaw test of the PSUV The same procedure as used for pure yaw test as applied for pure sway test to find the effects of the body motion amplitude and frequency on the accuracy of the numerical results. The PSUV is moved in longitudinal direction ( −direction) with a speed between 1 to 4 / . At first the effect of the variation of frequency is studied by setting a constant amplitude of = 0.01 with different frequencies. Then, the effect of variation of amplitude is investigated by setting = 1 with four different amplitudes of 0.012 , 0.024 , 0.036 and 0.12 .
The results of simulations of pure yaw test for N r′ at different frequencies and velocities are given in Table 3. The ratio of and the error are also shown in the table 3. The error (E) is calculated as follows.
For all cases, the error is minimum when the value of 0 = 1. Therefore, the effects of amplitude variations are investigated when 0 = 1.   Fig 7 when the frequency in each run is set to follow the relationships = 1. The error with respect to the analytical value that is ̇′ (Theory) = −0.0011. The error is between 3.5% to 6% when the frequency of the motion is set using the relationships = 1. The error remains constant when the amplitude of the body motion is more than 0.02. The velocity field around the body during the pure yaw test simulation is shown in Fig. 8. The PSUV axis is tangent to the path in pure yaw simulation. The numerical results are sensitive to the values of amplitude and frequency. However, the pure sway numerical simulation is more sensitive to the amplitude variation but the pure yaw numerical simulation is more sensitive to the frequency variation.
The added mass force and coefficient are due to the acceleration of imparted to the fluid by the body and therefore it has inertial characterise. They can be calculated using the potential theory as used to obtain the analytical solution. The fluid viscosity may affect the result if the real fluid flow simulation is used. The effect may be increased by increasing the amplitude of the motion due the separation phenomena. This effect is controlled by changing the wall boundary condition of PSUV from wall-no slip to wall-slip condition. The numerical solutions for both cases are shown in Table 4 for minimum and maximum amplitude. The error is significantly decreased by eliminating the viscosity and separation phenomena with wall-slip condition.  The transverse velocity should be zero in the pure yaw test that is v = ẏ0 cos ψ i − U 0 sin ψ i = 0. The value of is varied during the test between 0 and − 0 , Fig. 9. If the value of 0 is small then cos 0 ≈ 1 and 0 ≈ 0 and it may be written that = 0 cos − 0 sin =̇0 − 0 . The transverse motion of the body is 0 = 0 and then ẏ0 = −a 0 ω cos ωt. The transverse velocity may be written as: The value of 0 = 0 if the ratio 0 = 1.

Control and validation of the results
The test parameters for pure sway and yaw tests are obtained by doing a huge CFD computations on a PSUV. The test parameters are amplitude and frequency in each velocity. The values are applied for simulation of PMM tests for two SUOFF submarine models as shown in Fig. 1. The requirements of the y + ≅ 30 and the achievement of the residual 10 −5 in the mesh generation are applied to simulate the pure sway and heave tests for the SUBOFF models.
According to the calculations on PSUV for pure sway test, the amplitude is set to be equal to 0.01L and the frequency is set 0.2 Hz at a speed of 4 m/s to obtain an accurate result. The results of the simulation for both SUBOFF models are shown in Table 5. It can be seen that the use of the PSUV simulation results provide a relatively accurate result for the coefficients ′ , ̇′ , ′ , ′ and ′ but very poor results for ̇′ and ′ .
The numerical results for the coefficients ′ and ̇′ have considerable error as shown in the table 5. These coefficients are obtained using following relationships in the model test in towing tank.  is the mass of the body, is the longitudinal centre of gravity, l is the distance between strut and centre of gravity, F s and F b are the measured force act on struts in stern and bow respectively, γ s and γ b are the phase angle in stern and bow, respectively. The experimental value of non-dimensional added moment about −axis due to the transverse acceleration is ̇′ = 0.000202. If it is considered that the SUBOFF appended model weight is m=756 kg and x G = − 0.003 , the relationship for the coefficient N′ v̇ may be given in the following form.
The coefficient ̇′ for different values of x G is given in Table 6. It shows that for small variation of x G (lower than 0.01L) the error is significantly increase then we can deduce that the error is may be due to the small variation of . This is also may happen for the coefficient ̇. In order to find the effect of ̇′ on prediction of maneuvering characteristics of the SUBOFF submarine, the straight ahead, turning circle, stopping manoeuvres for full appended SUBOFF model are simulated using experimental and calculated values for According to the calculations done on the PSUV for the pure yaw test, the amplitude is set to be equal to 0.01L and the value of frequency is set according to the ratio of / 0 = 1 or = 0 /2 . Therefore, the frequency is set to 0.637 Hz at a speed of 4 m/s to obtain the results. The results of simulation of pure yaw and pitch tests for both SUBOFF models are listed in Table 7. The experimental results are also given to provide a comparison and error estimation. It can be seen that the use of the data obtained from the PSUV simulation provide a relatively accurate result for almost all the coefficients. However, the results of computation for the horizontal plane, i.e. the pure yaw test, are more accurate than the result for the vertical plane, i.e. the pure pitch test. The numerical solutions for almost all coefficients are relatively accurate in compare with experimental data except for two added mass coefficients, ̇′ and ′ , and two damping coefficients, ′ and ′ . It indicate that the test parameters obtained using the PSUV can be applied and find the relatively accurate results for majority of the coefficients.
The velocity fields around the full appended SUBOFF model for pure heave and pitch tests are shown in Fig. 11 with respect to the body reference frame. The fluid velocity is identical to the velocity of the body just near the body surface due to the no flux and no slip condition. The wake flow behind the body has a sinusoidal pattern as the body motion.

Discussions and conclusions
We focused on the important question for virtual PMM tests simulation that "how the virtual simulation of PMM tests can be made independent to the test parameters, i.e. the amplitude, frequency and the fluid flow?". These tests are simulated for a PSUV for a wide range of amplitudes, frequencies and velocities and the proper values and relationships are obtained to find the accurate results for the simulation of tests. The test parameters are applied for simulation of two SUBOFF models to control and validate them. The results of computations and validation indicate that: − The smaller the amplitude of motions the more accurate results for the hydrodynamic derivatives; − At small values of motion frequency for pure sway and heave motion, i.e. 0.1, 0.2 and 0.3 Hz, the results are accurate an error less than 8.5%, and the results are almost independent to the frequency of motion. By increasing the frequency more than 0.3 HZ, the error increase significantly. to the amplitude, frequency and flow velocity. − The discrepancy of ̇′ in CFD computation is due to the effect of the fluid viscosity and the separation phenomena. It is controlled by changing the wall boundary condition of PSUV from wall-no slip to wall-slip condition. The error is significantly decreased using wall-slip boundary condition. − The error for the coefficients ̇′ and ′ may be due the effect of the longitudinal position of center of gravity. The computations show the different values of x G have a significant effect on the value of these two coefficients. However, these two coefficients has a little influence on the estimation of the maneuvering charateristics of the SUBOFF submarine model. − According to simulations and their results we can purpose the following value, table 8, for elimination or reduction the effect of amplitude and frequency in each velocity: The unsteady simulations are carried out on desktop PC with Intel Core i7 Processor (2x 3.4 GHz) and 16GB RAM. The average computational time for all bodies are given in table 9.