Thrust Characteristics of Ducted Propeller and Hydrodynamics of an Underwater Vehicle in Control Motions

: A numerical technique to simulate the hydrodynamic behavior of ducted propellers attached to an underwater vehicle traveling under the mutually interacting ﬂow ﬁelds of the vehicle and the propellers is proposed; the hydrodynamic performance of the propellers and the hydrodynamic loading on the main body of the vehicle when it is in different kinds of motion is investigated numerically. In the research, 3D geometric models of the duct, propeller, and main body of the vehicle are ﬁrst constructed according to their geometrical features. A computational ﬂuid dynamics (CFD) technique based on the hybrid algorithm of dynamic mesh-nested sliding mesh is applied to solve the Navier–Stokes equations that govern the ﬂuid motion around the duct, propeller, and main body of the vehicle when it is in motion. These equations are solved numerically with the CFD code Fluent. With the proposed numerical simulation technique, the hydrodynamic characteristics of the thrusts generated by the ducted propellers and the loading on the main body in the vehicle system under the mutually interacting ﬂow ﬁelds are observed. The results of our numerical simulation indicate that the hybrid algorithm of dynamic mesh-nested sliding mesh can simulate multiple degrees of freedom of motion in underwater vehicle systems. In different motion states, the main body exerts a signiﬁcant inﬂuence on the investigated ﬂow ﬁelds; during the vehicle motions, negative wakes are formed on both sides of the main body, which lead to a decrease in the thrusts generated by the propellers on both sides. The thrust of the middle propeller is greater than that of the normal single one because of the obstructing effect in the tunnel caused by the main body.


Introduction
The underwater robot, known as a remotely operated vehicle, is a working robot. It can do the work normally completed by humans by diving into the water to complete underwater operations; thus, it is also known as a diving instrument. An underwater robot usually comprises a main body and several ducted propellers. In some cases, the horizontal and vertical wings are set up to control the balance, which increases its stability [1]. Ducted propellers are the primary actuators controlling the motion of the underwater robot. At present, the research on underwater vehicles is mainly carried out through model tests and numerical simulation. Model tests can directly reflect the performance of the underwater vehicle, and many scholars have carried out research on the underwater robot through model tests [2]. However, model tests can be limited by a test site such as a towing tank, and making a test model is a complex job, costing time and money. At the same time, it is restricted by monitoring equipment and technology such as sensors. Sometimes we cannot get the data we want. The reliability of CFD simulation needs to be verified by physical tests. This is because the accuracy of the calculation results of CFD technology usually depends on whether the definitions of the boundary conditions and physical parameters in the pre-processing of complex flow field simulations are truly consistent with reality, and whether the calculation methods and post-processing are accurate. In a numerical simulation, a model can be quickly and easily established and modified, and CFD software has a substantial number of post-processing modules that acquire the data of specific areas at specific times easily and accurately. Scholars have implemented additional functions with UDFs and other secondary developments of CFD software [3,4].
In current numerical simulation methods for the underwater robot system, using the principle of relativity to simulate motion is the most common practice. The underwater robot is placed in the center of the computational domain, and the rotation motion of the propellers is simulated by the sliding mesh technique. The motion of the underwater robot is simulated by changing the velocity conditions. Scholars have studied various types of propellers using sliding mesh, and the results showed that sliding mesh technology can effectively simulate the working conditions of a propeller at different rotation speeds [5]. However, as a quasi-static numerical analysis method, sliding mesh technology is an equivalent motion realized by a relative motion principle; the object is stationary and the grid rotates around the object, so it cannot reflect the actual movement. Dynamic mesh technology can realize the multi-degree-of-freedom motion of the robot itself in water [4], nevertheless, as the propeller speed increases and the mesh deformation rate increases, negative volume occurs easily, which leads to calculation failure. Using the dynamic mesh technique in conjunction with sliding mesh is another method of realizing the linear motion of the underwater robot and the rotating motion of the propeller. It is a more practical and feasible numerical simulation method that overcomes some of the shortcomings of dynamic mesh or sliding mesh alone [6]. However, laminating technology is used to update the mesh in this method, which can only realize forward and backward movement [7].
In this paper, three-dimensional underwater robot motion with multiple degrees of freedom was simulated through the hybrid technique of dynamic mesh and sliding mesh. The sliding mesh was used in the propeller rotation area, which solved the problem of negative volume of the dynamic mesh encountered in past research, and through spring smoothing and local remeshing, the mesh could be updated to overcome the limitation of the robot to only linear motion and realize its multi-degree-of-freedom motion in the fluid domain. According to the method, the hydrodynamics of the underwater robot system under different motion states and the thrust characteristics of each position propeller were studied. Firstly, the 3D geometric model was constructed according to selected geometric elements of the underwater robotic system [8]. Secondly, a hybrid computational domain of structural and unstructured meshes was constructed, and then the N-S equations were solved with the hybrid mesh technique and computational fluid dynamics method in the whole computational domain [9]. The hydrodynamic characteristics of the underwater robot system under different motion states were simulated and observed.

Governing Equations and Turbulence Models
The fluid in the research was assumed to be incompressible viscous fluid. The equations governing the fluid motions around the duct, propeller, and underwater vehicle in an unsteady motion are given as follows: 1.
Momentum equation: ρ ∂u y ∂t + u x ∂u y ∂x + u y ∂u y ∂y + u z ∂u y ∂z where x, y, z are the three-dimensional space coordinates of the flow field, u x , u y , u z the three-dimensional velocity components of the fluid, p the pressure of the flow field, ρ the density of the fluid, and ν the viscosity coefficient of the fluid motion. A standard RNG k-ε turbulence model was applied to describe turbulence within the flow field; the turbulence equations are as follows: In these equations, G k represents the generation of turbulence kinetic energy due to the mean velocity gradients, calculated as described in Modeling Turbulent Production in the k-ε Models. G b is the generation of turbulence kinetic energy due to buoyancy, calculated as described in Effects of Buoyancy on Turbulence in the k-ε models. Y M represents the contribution of the fluctuating dilatation in compressible turbulence to the overall dissipation rate, calculated as described in Effects of Compressibility on Turbulence in the k-ε models. The quantities α k and α ε are the inverse effective Prandtl numbers for k and ε, respectively. S K and S ε are user-defined source terms.

Mesh Update Mode
In our research, simulations of the underwater vehicle motion were accomplished by applying dynamic mesh-nested sliding mesh in the computational domain; that is, the area where the propellers are located was set as sliding mesh, while the area outside the propellers was set as dynamic mesh. The dynamic mesh technique mainly utilizes mesh stretching, compressing, increasing, reducing, and remeshing to change the computational domain [10,11]. In updating dynamic mesh, three techniques are usually used, i.e., spring smoothing, layering, and remeshing [12].
Selection of mesh update methods was based on the physical nature of our simulation regarding the number and type of mesh, accuracy of calculation, time spent on calculation, and sub-regional division in the computational domain [13]. After several trial calculations, we decided to combine the spring smoothing and remeshing methods, forming a composite dynamic updating mesh system to achieve accurate and fast computation.
In the simulation, the CFD code Fluent was used for determining the thrusts generated by the ducted propellers and the loading on the main body at every time step. There are two main ways to specify the motion of the boundary of dynamic mesh in Fluent: using transient Profile files or UDFs (User-Defined Functions). Some simple motion forms or constant motion can be specified directly using Profile files, while UDFs are needed to describe more complex functional motions [14]. For consistency, all motion assignments in the simulation were described by UDFs.

Comparison between Numerical and Experimental Results of a Single Ducted Propeller Hydrodynamics
Due to the high speed rotation of the propeller, the deformation and update speed of the grid around the propeller far exceeds that of other parts of the underwater vehicle. Thus, the rotation of the propeller is the most significant and critical influence factor on the calculation of the underwater vehicle system based on a lot of previous research. Limited by test conditions, funding, time, etc., we could not conduct a model test of the underwater robot. Although the validation of the underwater vehicle calculation cannot be validated from the model test, it is a good option to validate the most critical and significant part of the validation. Therefore, in this paper, we compared the calculation of a single ducted propeller with the model test, and carried out the calculation of the underwater vehicle based on it.
To validate the effectiveness of the hybrid algorithm of dynamic mesh-nested sliding mesh that was used in the analysis of the underwater vehicle system, existing experimental data on the representative ducted propeller of type Ka 4-70/19A [15] were taken as a reference. Comparisons between the numerical results of the hydrodynamic characteristics of the ducted propeller and those of the laboratory experimental data in [15] under the same condition were made. In computation, the parameters of the ducted propeller and the operation condition were made consistent with those of model test, and the primary parameters of the test ducted propeller are given in Table 1. Catheter inlet diameter (mm) 61

Computational Domain Construction and Boundary Conditions
To perform the comparison computation, four kinds of computational domains were constructed: Domains I to IV. The four domains were all constructed in cylindrical form sharing the same axis. Domains I to III were nested together by volume in descending order from Domain III to Domain I, constituting an integrated motion domain with the centers located at the midpoint of their common axis. In numerical computations, the integrated motion domain consisting of Domains I to III traveled along the axis of Domain IV to simulate the linear motion of the ducted propeller, and Domain IV acted as the background domain.
The four computational domains possessed different structures: Domain I comprised a cylindrical rotating flow field induced by the rotating propeller in the duct with its domain axis and domain length consistent with those of the duct, its outer boundary equal to 105% diameter of the propeller, and the fluid rotating speed in the domain the same as that of the propeller; Domain II comprised a thin flow field surrounding Domain I describing the detail of the flow field; Domain III was an outer flow field of Domains I and II; and Domain IV acted as the background domain. The relationships of the four domains are shown in To ensure continuity of the flow fields across the boundaries of Domains Ⅰ to Ⅳ, the interface technique was adopted along the boundaries to combine the domains, forming an integral computational one. In this boundary condition technique, the original outer boundary of Domain Ⅳ became the boundary of the integral domain while the original boundary conditions of Domains Ⅰ and Ⅱ and those between Domains Ⅱ and Ⅲ were determined automatically with the interface technique. In the four domains, the sliding mesh technique was adopted in Domain Ⅰ, while the dynamic mesh technique was adopted in Domains Ⅱ to Ⅳ. The fundamental parameters and techniques applied in the four domains are given in Table 2, and the boundary conditions in the computational domains applied in this research are given in Table 3.

Results of Comparison
To examine the independence of the mesh number in the calculation results, different numbers and types of meshes in computation were used to compare the numerical results with the laboratory experimental data to find a suitable number of meshes to achieve the goals of accuracy and low computation time. The mesh information of the three computational cases is provided in Table 4. Table 4. Mesh information in the three cases.

Domain
Number of Meshes (10 5 ) To ensure continuity of the flow fields across the boundaries of Domains I to IV, the interface technique was adopted along the boundaries to combine the domains, forming an integral computational one. In this boundary condition technique, the original outer boundary of Domain IV became the boundary of the integral domain while the original boundary conditions of Domains I and II and those between Domains II and III were determined automatically with the interface technique. In the four domains, the sliding mesh technique was adopted in Domain I, while the dynamic mesh technique was adopted in Domains II to IV. The fundamental parameters and techniques applied in the four domains are given in Table 2, and the boundary conditions in the computational domains applied in this research are given in Table 3.

Results of Comparison
To examine the independence of the mesh number in the calculation results, different numbers and types of meshes in computation were used to compare the numerical results with the laboratory experimental data to find a suitable number of meshes to achieve the goals of accuracy and low computation time. The mesh information of the three computational cases is provided in Table 4.  Table 5 shows the numbers and types of meshes. The mesh in case 1 was sparse and the mesh in case 2 was fine. The finer mesh was used in this section. Figures 2-4 show the results of the numerical simulation and model test of the ducted propeller when the rotational speed was set at 1500 rpm and the incoming velocity was zero for the mesh numbers 1.76 × 10 5 , 2.78 × 10 5 , and 4.12 × 10 5 , respectively [15]. Figure 5 shows the variation curve of the ducted propeller torque coefficient and efficiency when the ducted propeller had different Va speeds with zero velocity incoming flow, and the propeller speed was 1500 rpm. In Figures 2-4 we see that the numerical simulation calculation error decreased gradually with the increase in the number of meshes: the average error was 30% when the total number of meshes was 1.76 × 10 5 ; the average error was 25% when the total number of meshes was 2.78 × 10 5 ; and the average error was 18% when the total number of meshes was 4.12 × 10 5 . Therefore, the accuracy of the numerical simulation increased with the increase in mesh fineness.  Table 5 shows the numbers and types of meshes. The mesh in case 1 was sparse and the mesh in case 2 was fine. The finer mesh was used in this section. Figures 2-4 show the results of the numerical simulation and model test of the ducted propeller when the rotational speed was set at 1500 rpm and the incoming velocity was zero for the mesh numbers , respectively [15]. Figure 5 shows the variation curve of the ducted propeller torque coefficient and efficiency when the ducted propeller had different Va speeds with zero velocity incoming flow, and the propeller speed was 1500 rpm. In Figures 2-4, we see that the numerical simulation calculation error decreased gradually with the increase in the number of meshes: the average error was 30% when the total number of meshes was 5  . Therefore, the accuracy of the numerical simulation increased with the increase in mesh fineness. As shown in Figure 5, the propeller's torque coefficient gradually decreased with the increase in the advance speed, which means that the propeller generated more thrust by absorbing the same torque when the advance speed increased. The efficiency curve is consistent, and the efficiency increases with the increase in inlet velocity.                 As shown in Figure 5, the propeller's torque coefficient gradually decreased with the increase in the advance speed, which means that the propeller generated more thrust by absorbing the same torque when the advance speed increased. The efficiency curve is consistent, and the efficiency increases with the increase in inlet velocity. Figures 6 and 7 show the numerical simulation results of the thrust coefficient of the ducted propeller and model test results when the rotation speed was 500 rpm and 1000 rpm, respectively, and the incoming flow velocity was zero while the inlet velocity (Va) differed [15].
By comparing the results of Figures 4, 6 and 7, we found that the numerical calculation results from the ducted propeller were slightly smaller than the experimental results. The average error was 18%, and with the inlet velocity increased, the error decreased in a general trend, and when the inlet velocity was 1 kn, the error was less than 8%. Furthermore, the thrust measured in the model test included some components and connecting components attached to the duct and the propeller. Unfortunately, we could not obtain the detailed dimensional data of these components, so these components were not included in the simulation in this paper; therefore, the actual error should have been smaller. Nevertheless, the curve trends of numerical simulation and the model test results were consistent. Thus, the results reflected the thrust variation and trend of the propeller, which could be used as a supplement to the model test. It also indicated that the hybrid technique of dynamic mesh and sliding mesh used in this paper was suitable for the numerical simulation with spring smoothing and remeshing.  Figures 6 and 7 show the numerical simulation results of the thrust coefficient of the ducted propeller and model test results when the rotation speed was 500 rpm and 1000 rpm, respectively, and the incoming flow velocity was zero while the inlet velocity (Va) differed [15].  By comparing the results of Figures 4, 6, and 7, we found that the numerical calculation results from the ducted propeller were slightly smaller than the experimental results. The average error was 18%, and with the inlet velocity increased, the error decreased in a general trend, and when the inlet velocity was 1 kn, the error was less than 8%. Furthermore, the thrust measured in the model test included some components and connecting components attached to the duct and the propeller. Unfortunately, we could not obtain the detailed dimensional data of these components, so these components were not included in the simulation in this paper; therefore, the actual error should have been smaller. Nevertheless, the curve trends of numerical simulation and the model test results were consistent. Thus, the results reflected the thrust variation and trend of the propeller, which could be used as a supplement to the model test. It also indicated that the hybrid technique of dynamic mesh and sliding mesh used in this paper was suitable for the numerical simulation with spring smoothing and remeshing.
The longitudinal comparison is made of the numerical simulation results in Figures 4, 6, and 7: with the increase in inlet velocity, the thrust value of the ducted propeller decreased gradually under the same rotational speed. For example, when the rotational speed was 1000 rpm, the thrust coefficient was 0.3257 at 0.2 kn and 0.0687 at 1.0 kn. It was reduced by 78.9%. This was mainly because the thrust of the ducted propeller is directly related to the velocity at the surface of the propeller and the axial-induced velocity of the propeller. With the increase in inlet velocity, the axial induction velocity decreases, resulting in less propeller thrust. Table 6 shows the corresponding data of axial-induced velocity and inlet velocity produced by the propeller when it was moving straight at different speeds.  The longitudinal comparison is made of the numerical simulation results in Figures 4, 6, and 7: with the increase in inlet velocity, the thrust value of the ducted propeller decreased gradually under the same rotational speed. For example, when the rotational speed was 1000 rpm, the thrust coefficient was 0.3257 at 0.2 kn and 0.0687 at 1.0 kn. It was reduced by 78.9%. This was mainly because the thrust of the ducted propeller is directly related to the velocity at the surface of the propeller and the axial-induced velocity of the propeller. With the increase in inlet velocity, the axial induction velocity decreases, resulting in less propeller thrust. Table 6 shows the corresponding data of axial-induced velocity and inlet velocity produced by the propeller when it was moving straight at different speeds.  Figure 8 shows the mesh update process when the inlet velocity was 0.4 kn and the rotational speed of the ducted propeller was 1000 rpm.  Figure 8 shows that the initial mesh of the computational domain was regular and high-quality. As the time steps progressed, the motion calculation of Domain II containing the ducted propeller moved in the calculation Domain III. The regularity of the mesh was destroyed and the quality began to worsen. The mesh distortion rate and skewness (the parameter characterizing the mesh quality and the range is between 0 and 1; the smaller the parameter, the better the mesh quality) was always controlled below 0.7, which indicated that the quality of the updated mesh met the calculation requirement as long as the mesh quality control parameters were set reasonably. For the numerical result error, the following factors are relevant: (1) In this paper, we essentially discretized the governing equation of the fluid and solved the discrete governing equation. The obtained numerical solution met the convergence condition, and there was an error between the numerical solution and the analytical solution.
(2) The numerical calculation was based on the fluid micelle-mesh, the minimum computing unit. The quality of the mesh directly affected the accuracy of the numerical results. (3) Considering the time cost of calculation, the number of meshes cannot be too large, so the mesh quality needed to be further improved. Moreover, due to the mesh updates in the form of spring smoothing and remeshing, the mesh quality was reduced compared with the initial mesh over time.

Numerical Simulation Analysis of the Underwater Robotic System
The thrust characteristics of the ducted propeller were analyzed as above, and the influence of the main body of the underwater robot on the propeller was not considered. In fact, this effect was significant for the thrust of the propeller [4]. The whole underwater robotic system was numerically simulated to obtain the coupling relationship between the main body of the underwater robot and the propeller of the control mechanism.

Construction of the Geometric Model of the Underwater Robotic System
In this paper, the main body of the underwater robot was a water droplet streamline cylinder structure. The overall length of the main body was L = 510 mm. The maximum width was 320 mm, about 3/5 L from the stern, and the total height was 240 mm. The port and starboard were symmetrically arranged with a ducted propeller, which mainly controlled the direction of motion. The stern was arranged with a ducted propeller, mainly to produce forward or backward thrust. The ducted propellers adopted the Ka 4-70/19A standard for ducted propellers, as shown in Figure 9, and the specific parameters are  Figure 8 shows that the initial mesh of the computational domain was regular and high-quality. As the time steps progressed, the motion calculation of Domain II containing the ducted propeller moved in the calculation Domain III. The regularity of the mesh was destroyed and the quality began to worsen. The mesh distortion rate and skewness (the parameter characterizing the mesh quality and the range is between 0 and 1; the smaller the parameter, the better the mesh quality) was always controlled below 0.7, which indicated that the quality of the updated mesh met the calculation requirement as long as the mesh quality control parameters were set reasonably. For the numerical result error, the following factors are relevant: (1) In this paper, we essentially discretized the governing equation of the fluid and solved the discrete governing equation. The obtained numerical solution met the convergence condition, and there was an error between the numerical solution and the analytical solution.
(2) The numerical calculation was based on the fluid micelle-mesh, the minimum computing unit. The quality of the mesh directly affected the accuracy of the numerical results. (3) Considering the time cost of calculation, the number of meshes cannot be too large, so the mesh quality needed to be further improved. Moreover, due to the mesh updates in the form of spring smoothing and remeshing, the mesh quality was reduced compared with the initial mesh over time.

Numerical Simulation Analysis of the Underwater Robotic System
The thrust characteristics of the ducted propeller were analyzed as above, and the influence of the main body of the underwater robot on the propeller was not considered. In fact, this effect was significant for the thrust of the propeller [4]. The whole underwater robotic system was numerically simulated to obtain the coupling relationship between the main body of the underwater robot and the propeller of the control mechanism.

Construction of the Geometric Model of the Underwater Robotic System
In this paper, the main body of the underwater robot was a water droplet streamline cylinder structure. The overall length of the main body was L = 510 mm. The maximum width was 320 mm, about 3/5 L from the stern, and the total height was 240 mm. The port and starboard were symmetrically arranged with a ducted propeller, which mainly controlled the direction of motion. The stern was arranged with a ducted propeller, mainly to produce forward or backward thrust. The ducted propellers adopted the Ka 4-70/19A standard for ducted propellers, as shown in Figure 9, and the specific parameters are shown in Table 7. The three-dimensional model of the ducted propeller was constructed using the method of converting the shape value to the three-dimensional space coordinate point [16]. The 3D geometric model of the underwater robotic system is shown in Figure 10. shown in Table 7. The three-dimensional model of the ducted propeller was constructed using the method of converting the shape value to the three-dimensional space coordinate point [16]. The 3D geometric model of the underwater robotic system is shown in Figure  10. The coordinate system was defined uniformly here, as shown in Figure 10b: the direction of the z-axis was defined as the forward and backward direction of the underwater robot, where the positive value was forward and the negative value backward. The Y direction was the left and right direction of the underwater robot, where the positive value was right and the negative value was left. The X direction was the heave direction of the underwater robot, where the positive value was rising and the negative value was sinking.

Computational Domain and Boundary Conditions
The computational domain, the position of ducted propellers, and the main body of the underwater robot are shown in Figure 11.  shown in Table 7. The three-dimensional model of the ducted propeller was constructed using the method of converting the shape value to the three-dimensional space coordinate point [16]. The 3D geometric model of the underwater robotic system is shown in Figure  10. The coordinate system was defined uniformly here, as shown in Figure 10b: the direction of the z-axis was defined as the forward and backward direction of the underwater robot, where the positive value was forward and the negative value backward. The Y direction was the left and right direction of the underwater robot, where the positive value was right and the negative value was left. The X direction was the heave direction of the underwater robot, where the positive value was rising and the negative value was sinking.

Computational Domain and Boundary Conditions
The computational domain, the position of ducted propellers, and the main body of the underwater robot are shown in Figure 11. The coordinate system was defined uniformly here, as shown in Figure 10b: the direction of the z-axis was defined as the forward and backward direction of the underwater robot, where the positive value was forward and the negative value backward. The Y direction was the left and right direction of the underwater robot, where the positive value was right and the negative value was left. The X direction was the heave direction of the underwater robot, where the positive value was rising and the negative value was sinking.

Computational Domain and Boundary Conditions
The computational domain, the position of ducted propellers, and the main body of the underwater robot are shown in Figure 11.
According to the characteristics of the problems studied in this paper, the computational domain was divided into four computational domains. Domain I was the rotating flow field with the propeller shaft as the rotating axis, its length Lc and its diameter up to 105% of the propeller diameter. In this paper, I-I, I-II and I-III represent the computational domains of the three ducted propellers, where I-I and I-II are the computational domains of the left-and right-ducted propellers, and I-III is the computational domain of the middleducted propeller. The motion in computational Domain II includes the robot body and the ducted propellers; Domain III is made up of unstructured mesh and includes Domain II; Domain IV is composed of the structural mesh including Domain III. Table 8 Table 8 shows the types of boundary conditions between each computational domain. Transition boundary between regions I and IV Interface Figure 12 shows the overall computational domain body mesh and the surface mesh of the underwater robotic system. Table 9 shows the types and numbers of meshes of the computational domains.    Figure 12 shows the overall computational domain body mesh and the surface mesh of the underwater robotic system. Table 9 shows the types and numbers of meshes of the computational domains.  Table 8 shows the types of boundary conditions between each computational domain. Transition boundary between regions I and IV Interface Figure 12 shows the overall computational domain body mesh and the surface mesh of the underwater robotic system. Table 9 shows the types and numbers of meshes of the computational domains.     Figure 13 shows the thrust of the three ducted propellers and the fluid resistance curve of the main body with the coupling of the main body and ducted propeller when the underwater robotic system advanced at different speeds; the rotational speed of each propeller was 1200 rpm.   Figure 13 shows the thrust of the three ducted propellers and the fluid resistanc curve of the main body with the coupling of the main body and ducted propeller when the underwater robotic system advanced at different speeds; the rotational speed of each propeller was 1200 rpm. As seen in Figure 13, no matter how the forward speed changed, the thrust of th propeller located behind the underwater robotic system was greater than the thrust of single propeller. Moreover, the proportion of increase grew sharply, as the advance speed increased. It increased from 9.9% at 0.2 kn to 349% at 2.0 kn. Meanwhile, as the forward speed increased, the thrust generated by the rear propeller decreased from 7.462 N at 0. kn to 5.426 N at 2.0 kn, a decrease of 27.3%.

Linear Motion
The reason is that under the influence of the main body, the water flow in the are behind the main body is slow. It leads to the decrease in the axial average advance speed at the rear propeller surface and the increase in the axial-induced velocity produced b the propeller. Thus, the thrust of the rear propeller is greater than that of a single propeller At the same time, when the underwater robotic system advances in a straight line at dif ferent speeds, the flow field behind the underwater robot is extremely uneven because o the influence of the main body. The uneven variation of the flow field causes the thrust o the rear propeller to fluctuate greatly, more than that of a single propeller. When the speed is high, the flow behind the main body will produce a reflux phenomenon, which inten sifies the decrease in the axial flow velocity in this area and reduces the velocity of wate flow perceived at the surface of propeller. Finally, the rear propeller sends out greate thrust than a single propeller.
The propellers on the left and right sides of the underwater robotic system have th same thrust curve height, but they do not completely coincide. Moreover, each thrust wa As seen in Figure 13, no matter how the forward speed changed, the thrust of the propeller located behind the underwater robotic system was greater than the thrust of a single propeller. Moreover, the proportion of increase grew sharply, as the advance speed increased. It increased from 9.9% at 0.2 kn to 349% at 2.0 kn. Meanwhile, as the forward speed increased, the thrust generated by the rear propeller decreased from 7.462 N at 0.2 kn to 5.426 N at 2.0 kn, a decrease of 27.3%.
The reason is that under the influence of the main body, the water flow in the area behind the main body is slow. It leads to the decrease in the axial average advance speed at the rear propeller surface and the increase in the axial-induced velocity produced by the propeller. Thus, the thrust of the rear propeller is greater than that of a single propeller. At the same time, when the underwater robotic system advances in a straight line at different speeds, the flow field behind the underwater robot is extremely uneven because of the influence of the main body. The uneven variation of the flow field causes the thrust of the rear propeller to fluctuate greatly, more than that of a single propeller. When the speed is high, the flow behind the main body will produce a reflux phenomenon, which intensifies the decrease in the axial flow velocity in this area and reduces the velocity of water flow perceived at the surface of propeller. Finally, the rear propeller sends out greater thrust than a single propeller.
The propellers on the left and right sides of the underwater robotic system have the same thrust curve height, but they do not completely coincide. Moreover, each thrust was smaller than that of a single propeller. The decrease proportion rose with the increase in forward speed, from 2.6% at 0.2 kn to 135% at 2.0 kn.
The reason for this result is that the underwater robotic system is a symmetrical arrangement. Although the underwater robotic system moves forward at different speeds, the flow field between the left and right propellers is almost identical. Moreover, because of the complexity of the flow field and the error in the numerical calculation itself, the calculation results cannot be completely consistent with experimental results. As the main body of the underwater robot system moves forward, the flow near the two sides of the main body has backward velocity, forming a negative wake. Therefore, the axial velocities of the left and right propellers increase, and the thrusts of the propellers decrease. Moreover, with the rise in forward speed, the velocity of negative wake also increases, which makes the proportional reduction of propeller thrust gradually increase.
From the resistance curve of the main body in Figure 13, we see that with the increase in forward velocity, the fluid resistance of the main body increased sharply. This was consistent with the geometric relationship where the fluid resistance was proportional to the square of the velocity of the object in the fluid. Figure 14 shows the thrust curve of the ducted propellers and the ducted propellers under the influence of the main body, when the underwater robotic system heaved at different speeds; the rotational speed of each propeller was 1200 rpm.

Vertical Deep Motion
smaller than that of a single propeller. The decrease proportion rose with the inc forward speed, from 2.6% at 0.2 kn to 135% at 2.0 kn.
The reason for this result is that the underwater robotic system is a symme rangement. Although the underwater robotic system moves forward at different the flow field between the left and right propellers is almost identical. Moreover, of the complexity of the flow field and the error in the numerical calculation it calculation results cannot be completely consistent with experimental results. As t body of the underwater robot system moves forward, the flow near the two sid main body has backward velocity, forming a negative wake. Therefore, the axial v of the left and right propellers increase, and the thrusts of the propellers decreas over, with the rise in forward speed, the velocity of negative wake also increase makes the proportional reduction of propeller thrust gradually increase.
From the resistance curve of the main body in Figure 13, we see that with the in forward velocity, the fluid resistance of the main body increased sharply. This w sistent with the geometric relationship where the fluid resistance was proportion square of the velocity of the object in the fluid. Figure 14 shows the thrust curve of the ducted propellers and the ducted pr under the influence of the main body, when the underwater robotic system he different speeds; the rotational speed of each propeller was 1200 rpm.  Figure 14 shows that the thrust curves of each propeller are very different, regularity is not apparent when the underwater robotic system heaves. The thrus of the left and right propellers on both sides of the main body also show great diff The fundamental reason for this is that the shape of the underwater robot is a flat line body, and the upper surface of the main body is an elliptical plane. When th water robot system heaves, the relative downward flow meets the resistance of t body and disperses randomly to both sides, which makes the flow field environm both sides of the main body clearly different.

Two-Dimensional Linear Motion
A numerical simulation of the two-dimensional linear motion of the underw botic system is discussed in this section. The UDF control program was written to the specified motion of the underwater robot. By changing the value of Va, the d two-dimensional linear motion velocity of the underwater robot system was alte speed of motion is shown in Formula (7) Figure 14 shows that the thrust curves of each propeller are very different, and the regularity is not apparent when the underwater robotic system heaves. The thrust curves of the left and right propellers on both sides of the main body also show great differences. The fundamental reason for this is that the shape of the underwater robot is a flat streamline body, and the upper surface of the main body is an elliptical plane. When the underwater robot system heaves, the relative downward flow meets the resistance of the main body and disperses randomly to both sides, which makes the flow field environments on both sides of the main body clearly different.

Two-Dimensional Linear Motion
A numerical simulation of the two-dimensional linear motion of the underwater robotic system is discussed in this section. The UDF control program was written to control the specified motion of the underwater robot. By changing the value of Va, the different two-dimensional linear motion velocity of the underwater robot system was altered. The speed of motion is shown in Formula (7) below: Figure 15 shows the thrust curve of the ducted propellers and the thrust curve of the ducted propellers in each position under the influence of the main body, when the underwater robot system made the two-dimensional linear motion at different speeds and the rotational speed of each propeller was 1200 rpm.
thrust of the single propeller. Additionally, the thrust of the rear propeller also d with the increase in the absolute speed. The reduction ratio had a decrease of 14. 3 7.131 N at an absolute speed of 0.282 kn to 6.112 N at 2.828 kn. The thrust curv left and right propellers coincide, which was reduced overall by 112.4%. With the in absolute speed, their thrust decreased from 6.612 N at an absolute speed of 0.2 -0.822 N at 2.828 kn. In addition, their combined thrust was less than the thrust of propeller. The reason for these outcomes is basically the same as those mentione as they are the result of the change in the movement speed and the influence of t body. Thus, the negative wake formed around the main body, which increased th ity of flow at the propeller surface and reduced the thrust from the propeller. As the rear propeller is behind the main body of the underwater robot, the locity of the flow field in this area was less affected by the change in the veloci underwater robotic system. It decreased the change in the velocity of flow at the r peller surface. Furthermore, it reduced the proportional change in the thrust of propeller beyond that of a single propeller.

Two-Dimensional Circular Motion
The circular motion velocity of the underwater robot system is changed by t of Va. The underwater robot was controlled by a UDF program to move along th fied circle. The speed of motion is as follows (8) Figure 16 shows the thrust curve of the ducted propellers in each position u influence of the main body and the fluid resistance curve of the main body, w underwater robotic system made a two-dimensional circular motion at differen and the rotational speed of propeller was 1200 rpm. It can be seen from Figure 15 that when the underwater robot system performed two-dimensional linear motion, the thrust of the rear propeller was always greater than the thrust of the single propeller. Additionally, the thrust of the rear propeller also decreased with the increase in the absolute speed. The reduction ratio had a decrease of 14.3%, from 7.131 N at an absolute speed of 0.282 kn to 6.112 N at 2.828 kn. The thrust curves of the left and right propellers coincide, which was reduced overall by 112.4%. With the increase in absolute speed, their thrust decreased from 6.612 N at an absolute speed of 0.282 kn to −0.822 N at 2.828 kn. In addition, their combined thrust was less than the thrust of a single propeller. The reason for these outcomes is basically the same as those mentioned above, as they are the result of the change in the movement speed and the influence of the main body. Thus, the negative wake formed around the main body, which increased the velocity of flow at the propeller surface and reduced the thrust from the propeller.
As the rear propeller is behind the main body of the underwater robot, the axial velocity of the flow field in this area was less affected by the change in the velocity of the underwater robotic system. It decreased the change in the velocity of flow at the rear propeller surface. Furthermore, it reduced the proportional change in the thrust of the rear propeller beyond that of a single propeller.

Two-Dimensional Circular Motion
The circular motion velocity of the underwater robot system is changed by the value of Va. The underwater robot was controlled by a UDF program to move along the specified circle. The speed of motion is as follows (8): vel [1] = Va * cos(M _PI * time) vel [2] = Va * sin(M _PI * time) (8) Figure 16 shows the thrust curve of the ducted propellers in each position under the influence of the main body and the fluid resistance curve of the main body, when the underwater robotic system made a two-dimensional circular motion at different speeds and the rotational speed of propeller was 1200 rpm. In Figure 16, we see that the thrust curve of each propeller evidently changes at different velocities.
(1) The thrust curve of the rear propeller.
At low speed, the thrust curve of the rear propeller shows obvious regularity. With the increase in linear velocity, the curve begins to fluctuate. Although the thrust curve of the rear propeller fluctuates sharply when the linear velocity exceeds 0.5 kn, the shape and trend of the curve does not change fundamentally. The reason is that the flow field of the underwater robot system changes fundamentally when the linear velocity of circular motion is less than 0.5 kn. With the increase in the linear velocity of circular motion, the influence on thrust regularity of propeller passes the inflection point. This means that the increase in linear velocity will not fundamentally change the thrust regularity of the propeller.
(2) The thrust curves of single, left, and right propellers.
At low speed, the thrust curves of the single propeller, left propeller, and right propeller are close in shape and trend. Nevertheless, according to the numerical simulation results in this paper, this consistency is destroyed when the linear velocity increases to a certain value (greater than 1.0 kn). When the angle of motion θ is greater than 1.5, each curve shows its own unique characteristics, which are obviously different from the other curves. This is mainly because when the linear velocity is high, the flow field becomes extremely irregular and the flow field environment of the area where the propeller is located in a different position results in different thrust regularities of the propellers. In Figure 16, we see that the thrust curve of each propeller evidently changes at different velocities.
(1) The thrust curve of the rear propeller.
At low speed, the thrust curve of the rear propeller shows obvious regularity. With the increase in linear velocity, the curve begins to fluctuate. Although the thrust curve of the rear propeller fluctuates sharply when the linear velocity exceeds 0.5 kn, the shape and trend of the curve does not change fundamentally. The reason is that the flow field of the underwater robot system changes fundamentally when the linear velocity of circular motion is less than 0.5 kn. With the increase in the linear velocity of circular motion, the influence on thrust regularity of propeller passes the inflection point. This means that the increase in linear velocity will not fundamentally change the thrust regularity of the propeller.
(2) The thrust curves of single, left, and right propellers.
At low speed, the thrust curves of the single propeller, left propeller, and right propeller are close in shape and trend. Nevertheless, according to the numerical simulation results in this paper, this consistency is destroyed when the linear velocity increases to a certain value (greater than 1.0 kn). When the angle of motion θ is greater than 1.5, each curve shows its own unique characteristics, which are obviously different from the other curves. This is mainly because when the linear velocity is high, the flow field becomes extremely irregular and the flow field environment of the area where the propeller is located in a different position results in different thrust regularities of the propellers. Figure 17 shows the mesh update image of the underwater robotic system when the linear velocity (Va) was 2.0 kn. It can be seen from Figure 17 that with the increase in the circular motion angle, the influence of the underwater robotic system on the flow field in its region gradually diffused from around the robot to the whole flow field area. As the angle θ of circular motion increased, the meshes became more and more irregular, which showed the deterioration of mesh quality. The deterioration of the mesh quality led directly to fluctuations in the numerical calculation results. Figure 17 shows the mesh update image of the underwater robotic system when the linear velocity (Va) was 2.0 kn. It can be seen from Figure 17 that with the increase in the circular motion angle, the influence of the underwater robotic system on the flow field in its region gradually diffused from around the robot to the whole flow field area. As the angle θ of circular motion increased, the meshes became more and more irregular, which showed the deterioration of mesh quality. The deterioration of the mesh quality led directly to fluctuations in the numerical calculation results.    linear velocity (Va) was 2.0 kn. It can be seen from Figure 17 that with the increase in the circular motion angle, the influence of the underwater robotic system on the flow field in its region gradually diffused from around the robot to the whole flow field area. As the angle θ of circular motion increased, the meshes became more and more irregular, which showed the deterioration of mesh quality. The deterioration of the mesh quality led directly to fluctuations in the numerical calculation results.    linear velocity (Va) was 2.0 kn. It can be seen from Figure 17 that with the increase in the circular motion angle, the influence of the underwater robotic system on the flow field in its region gradually diffused from around the robot to the whole flow field area. As the angle θ of circular motion increased, the meshes became more and more irregular, which showed the deterioration of mesh quality. The deterioration of the mesh quality led directly to fluctuations in the numerical calculation results.    The figures illustrate that the fluid resistance in both the Z direction and the Y direction increased with the increase in the linear velocity of the circular motion. Moreover, the range of increase widens. At low linear velocity, the resistance curve fluctuates slightly. In contrast, at high linear velocity, the fluctuation is noticeable. It is consistent with the physical law that states fluid resistance is proportional to the square of velocity. At the same time, we also see that the shape and trend of the fluid resistance curve at different linear velocities are basically the same. The difference is that the upper and lower amplitude of the curve are much larger at high speeds than that at low speeds.

Three-Dimensional Circular Motion
The specified three-dimensional circular motion speed of the underwater robot is shown in the following Equation (9): vel[0] = 0.257 vel [1] = Va * cos(M _PI * time) vel [2] = Va * sin(M _PI * time) (9) By changing the value of Va, the different three-dimensional circular motion velocities of the underwater robot system can be altered. When the underwater robot system moved in a three-dimensional circular motion at different linear velocities, thrust curves of each propeller ensured that are shown in Figure 20.
tion increased with the increase in the linear velocity of the circular motion. Moreover, the range of increase widens. At low linear velocity, the resistance curve fluctuates slightly. In contrast, at high linear velocity, the fluctuation is noticeable. It is consistent with the physical law that states fluid resistance is proportional to the square of velocity. At the same time, we also see that the shape and trend of the fluid resistance curve at different linear velocities are basically the same. The difference is that the upper and lower amplitude of the curve are much larger at high speeds than that at low speeds.

Three-Dimensional Circular Motion
The specified three-dimensional circular motion speed of the underwater robot is shown in the following Equation (9) By changing the value of Va, the different three-dimensional circular motion velocities of the underwater robot system can be altered. When the underwater robot system moved in a three-dimensional circular motion at different linear velocities, thrust curves of each propeller ensured that are shown in Figure 20. It can be seen in Figure 20 that the linear velocity of circular motion had a major influence on the thrusts of the propellers at different positions. With the increase in the linear velocity of circular motion, the fluctuation amplitude of the thrust curve of each propeller begins to increase. When the linear velocity is high, the flow field is very disordered, which leads to the obvious fluctuation of the thrust curve. It also reflects the decisive influence of the velocity on the complexity of the flow field. Moreover, the thrust It can be seen in Figure 20 that the linear velocity of circular motion had a major influence on the thrusts of the propellers at different positions. With the increase in the linear velocity of circular motion, the fluctuation amplitude of the thrust curve of each propeller begins to increase. When the linear velocity is high, the flow field is very disordered, which leads to the obvious fluctuation of the thrust curve. It also reflects the decisive influence of the velocity on the complexity of the flow field. Moreover, the thrust curves of the single propeller, the left propeller, and right propeller cross many times, which means that the negative wake around the main body changed because of the three-dimensional circular motion. Additionally, the negative wake no longer reduced the thrust of the propeller at certain periods, and it, conversely, had a positive effect.
The thrust curve of the rear propeller is almost above the thrust curve of the single propeller, left propeller, and right propeller. This is mainly because the rear propeller is located behind the main body, where the blocking effect of the main body on the water flow still exists. Therefore, the flow in the rear propeller area becomes slower and the thrust from the propeller becomes larger.
Comparing the hydrodynamic characteristics of two-dimensional circular motion ( Figure 16) with those of three-dimensional circular motion (Figure 20), we found that the motion parameters of multiple dimensions had a significant effect on the thrust characteristics of the propeller control mechanism of the underwater robot system. They led to a drastic fluctuation of the curve. The main reason for this is that when the underwater robotic system makes a deep vertical motion, the upper surface of the main body hinders the flow of water, which eliminates the flow along both sides of the main body. Furthermore, it intensifies the complexity of the flow field in the propeller area, obviously causing the thrust to fluctuate. In addition, the disturbance of water flow is more obvious when the linear velocity increases.

Discussion and Conclusions
In this paper, three-dimensional motions with multiple degrees of freedom of the underwater robot system were numerically simulated and the thrust characteristics of the ducted propeller under different motion forms were studied and analyzed. Through the previous numerical simulation research and analysis, the following conclusions can be drawn: (1) The dynamic mesh was used in the motion of the underwater robot in the computational Domain III and the sliding mesh was used in the rotation motion of the ducted propeller. The technique of combined mesh allowed the underwater robot system to perform complex three-dimensional motions with multiple degrees of freedom. (2) When comparing the numerical simulation results with the experimental results, we determined that the hybrid technique of dynamic mesh and sliding mesh, using the mesh update method of spring smoothing and remeshing, met the engineering requirements. (3) Numerical simulation results showed that the flow field environment around the underwater robot system was obviously different when different motion modes were carried out at different motion speeds. Affected by this flow field environment, there were obvious differences in the hydrodynamic characteristics and the thrust characteristics of the control mechanism ducted propeller. At low speed, the fluid resistance of the main body increased sharply with the increase in motion speed. Additionally, the left and right propellers were influenced by the negative wake of the main body, and the thrust value was reduced; the thrust value of the rear propeller increased from the blocking effect of the main body on the flow. This obvious regularity changed fundamentally at high speed, and the thrust of the left, right, and rear propellers all fluctuated significantly. (4) Similar to the circumstance of thrust behaviors issued from a ducted propeller without the effects of the vehicle, for a given heading angle and a specific yawing motion, there is a roughly corresponding relationship between the thrust generated from the ducted propeller and the mean axial advance velocity on the propeller disk in the integrated system of the underwater vehicle and ducted propellers: the smaller the advance velocity, the greater the thrust. For the same rotating speed, the thrust generated by the propeller in the integrated system of underwater vehicle and ducted propellers is greater than that of the propeller without the effects of the vehicle, since the mean axial advance velocity met by the propeller disk in the integrated system under the influence in the vehicle's flow field is smaller than that without the effects of vehicle's flow field.
In this paper, because of the limitations of computer performance and calculation time cost, the quality of the mesh caused a large error in the results. In future research, we will refine the mesh size and improve the mesh type, for example, by using a polyhedral mesh instead of the tetrahedral mesh for more accurate results. In addition, the pressure distribution at the ducted propeller surface and the underwater influence on the flow near the propeller location are also worthy of study.