Electromagnetic Scattering of Near-Field Turbulent Wake Generated by Accelerated Propeller

: The electromagnetic scattering study of the turbulent wake of a moving ship has important application value in target recognition and tracking. However, to date, there has been insufﬁcient research into the electromagnetic characteristics of near-ﬁeld propeller turbulence. This study presents a new procedure for evaluating the electromagnetic scattering coefﬁcient and imaging characteristics of turbulent wakes in the near ﬁeld. By controlling the different values of the net momenta, a turbulent wake was generated using the large-eddy simulation method. The results show that the net momentum transferred to the background ﬂow ﬁeld determines the development of the turbulent wake, which explains the formation mechanism of the turbulence. Combined with the turbulent energy attenuation spectrum, the electromagnetic scattering characteristics of the turbulent wake were calculated using the two-scale facet mode. Using this method, the impact of different parameters on the scattering coefﬁcient and the electromagnetic image of the turbulence wake were investigated, to explain the modulation mechanism and electromagnetic imaging characteristics of the near-ﬁeld turbulent wake. Moreover, an application for estimating a ship’s heading is proposed based on the electromagnetic imaging characteristics of the turbulent wake.


Introduction
The movement of a ship changes the texture of the sea surface at a given time, making it possible to detect a moving ship [1][2][3][4]. At present, many studies have focused on wakes produced by towed bodies and self-propelled bodies in a steady state. However, a ship usually breaks the steady state by frequently accelerating and decelerating, transmitting the additional momentum to the background flow field, which leads to complex physical processes on the ocean surface, such as the formation of a jet, free-surface turbulence, and a large-scale vortex street structure [5]. The strong turbulence and nonlinear structure generated by the extra momentum will affect the electromagnetic scattering of the sea surface. Therefore, it is of great practical significance to study the electromagnetic scattering of the near-field turbulent wake generated by a propeller, as shown in Figure 1.
The imaging characteristics and detection of a far-field turbulent wake have been studied extensively [6][7][8], but little research has been conducted on the near-field turbulent wake. In a previous study, a numerical model of turbulence with small momentum was built, and the relationship between the excess momentum and flow characteristics of a turbulent wake was studied [9,10]. Wang et al. [11] researched the turbulence generated by the propellers of underwater targets. The interaction between underwater circular jets and free surfaces was investigated in the works of Madnia and Bernal [12]. These studies demonstrated that strong turbulence has an obvious influence on the free surface. In terms of the electromagnetic scattering of turbulence, Wang et al. [13] used the forward-backward iterative physical optics (FBIPO) method to calculate the coupled scattering characteristics of ships and waves in the near field. Luo [14] used the FLUENT commercial software to study the electromagnetic scattering in a turbulent structure caused by a propeller. However, this research was limited to two-dimensional wave models at low grazing angles. By using the vector radiation equation and the volume scattering theory for foam, Sun et al. [15] modified the traditional two-scale electromagnetic scattering method to create a turbulent scattering model of the sea surface. However, the calculation process used in this study is based on empirical constants and basic data obtained by solving the equation of motion. It is difficult to deal with nonlinear effects of turbulence. In summary, there is still no suitable method for modeling turbulent wakes, and the explanation of the electromagnetic scattering characteristics and imaging mechanism of turbulent wakes is insufficient. This paper demonstrates a simulation method to capture the turbulence wake of a ship with a propeller, using the STAR_CCM+ fluid dynamics computing software. The simulation results show the structure of the propeller turbulence, and the wake width is verified. A two-scale facet model (TSFM) which can calculate arbitrarily inclined small facet elements, is presented. Together with the turbulent energy attenuation spectrum, the electromagnetic scattering characteristics of the turbulent wake were fully investigated under different radar parameters and other variables. These results show the electromagnetic characteristics of a propeller turbulent wake and imaging mechanisms that are not possessed by other wakes. Finally, an application for estimating a ship's heading is proposed based on the characteristics of the turbulent wake.

Motion Control of Propeller
A ship is assumed to be subjected to two types of forces as it passes through the water: the resistance of the incoming flow and the thrust of the propeller.
When the ship sails in a horizontal direction at speed U, the resistance F R generated by the flow field and the thrust F T generated by the propeller correspond to the forward resistance momentum J R and backward propulsion momentum J T , respectively. According to [16], the resistance momentum J R of a ship is defined as: where R S is the drag coefficient of the ship, U is the towing speed magnitude in the x-direction, S D is the cross-sectional area, and ρ is the water density. In this study, to obtain the propulsion momentum J T , an indirect measurement method is adopted. The magnitude of the net force was applied as a criterion to estimate the difference between the thrust and drag. When the towing speed U is constant, the propulsion momentum J T is equal to the resistance momentum J R provided the net force is close to 0. Otherwise, the self-propelled speed U P was adjusted by increasing the speed of the propeller until the net force was neglected.
The magnitude of the net force when the speed of the ship is 2.2 ms −1 is shown in Figure 2a. The oscillation of the net force is due to the wave-building resistance of the ship and which gradually stabilises. The corresponding relationship between different propeller speeds n and the ship's self-propelled speed U P is shown in Figure 2b. When the net force is 0, the resistance momentum and propulsion momentum are equal in magnitude and opposite in direction. Therefore, the propulsion momentum J T at equilibrium can be expressed as: where U P is the magnitude of the self-propelled speed in the x-direction. In the simulation, the level of resistance momentum J R can be controlled by adjusting the speed of the ship, and the level of the propulsion momentum J T can be controlled by adjusting the speed of the propeller.
To evaluate the effect of the propulsion momentum J T , using the resistance momentum J R , the dimensionless J D can be defined as: If J D = 1, the propulsion momentum is equal to the resistance momentum, and the propeller has no additional net momentum transferred to the background flow field. If J D < 1, the propulsion momentum is greater than the drag momentum, and the propeller has additional net momentum transferred to the background flow field.

Geometric Modeling
Considering the computational efficiency, the entire simulation region was limited to a cuboid measuring 54 × 36 × 24 m. The distance from the stern to the velocity inlet was 18 m, the distance from the stern to the pressure outlet was 36 m, the width of both sides was 18 m, and the water depth was 18 m.
The hull model used in the simulation is a zoomed KRISO container ship (KCS) developed by the Korea Maritime and Marine Engineering Research Institute [17], as shown in Figure 3. The propeller type is KP505. The main parameters of the ship and propeller are listed in Table 1.  The KCS was subdivided by the cutting body mesh generator, to make the boundary layer more consistent with the flow field in a single direction. The KP505 propeller is subdivided by the polyhedral mesh generator, which makes the boundary layer more suitable for turbulence development. Finally, the number of ship grid points was 680,000, the number of propeller grid points was 96,000, and the total number of grid points in the entire simulation region was 6.82 million. Most of these are used to capture the propeller turbulence and the free surface. The grid discretization diagram of the KCS with the KP505 propeller is shown in Figure 4.

Governing Equation and Watershed Division
Owing to the existence of turbulence caused by the motion of the propeller, the large eddy simulation (LES) method was selected to calculate the flow field, given that it can handle the relationship between the different scale vortices and accurately capture the turbulent structure.
The LES method employed in this study is based on a numerical model for transient hydrodynamic flows. Its equations were solved by spatial filtering instead of by obtaining an average value [18]. The transport equations of the filtering mass, momentum, and energy can be written as: ∂ρ ∂t ∂ ∂t where ρ is the density, ∼ v is the filtering velocity, ∼ p is the filtering pressure, I is the unit tensor,T is the filtering stress tensor, q is the filtering heat flux, and f b is the resultant force. T SGS (subgrid scale (SGS)) is implemented by the Boussinesq approximation and is written as: where µ t is the turbulent viscosity and S is the strain rate tensor.
To obtain a closed-form solution of the governing equations, the turbulent viscosity µ t must be obtained at the sub-grid scale. The WALE sub-grid model [19] was used to solve µ t which can be expressed as: where ∆ is the width of the grid filter and S W is the deformation parameter. The width of the grid filter ∆ is defined by the cell volume V: where C W = 0.544 is the coefficient of the model, κ = 0.41 is von Karman's constant, and d is the distance to the ship surface.
In the simulation, the characteristics of the free surface were captured using the volume of fluid (VOF) method [20]. The VOF transport equation between different phase fractions α can be expressed as: According to the phase fraction, the density ρ and dynamic viscosity µ of the different positions in the stratified water are defined as: where ρ i is the density of the phases and µ i is the dynamic viscosity of the different phases.
To accurately capture the free surface characteristics, the adaptive mesh refinement (AMR) method was added to the simulation. This method, which dynamically regenerates the mesh and executes linear interpolation with the deformation of the free surface, efficiently takes advantage of the grid resources and reduces the sharpening of the free surface.
For the solution process, the second-order implicit scheme was selected for time discretization, the second-order upwind scheme was selected for space discretization, and in the turbulence solver, the AMG linear solver equation was used to solve the flow field.

Simulation Results of Turbulent Wakes
According to references (Walker et al. [21]), the half-width of the turbulence caused by the propeller jet should satisfy |y/x| = 0.1. The width at x = 0 is caused by the propeller radius. In Figure 5, the development trend of the turbulent wake in the simulation is the same as the theoretical turbulence width, which demonstrates the reliability of the simulation. The effect of propeller turbulence on the free surface at different net momentums was studied when F r = U/ gL = 0.29, and L is the length of the ship in Table 1. Considering that the ship cannot be ignored in the actual situation, the simulation results also include the Kelvin wake generated by the volume effect of the ship. The results are presented for different values of J D at 1, 0.173, and 0.097, as shown in Figure 6. When J D = 1, there is no net momentum generated by the propeller, as shown in Figure 6a, so the turbulent wake is not obvious at the free surface. As J D decreases, the propulsion momentum increases. Excess energy is transferred to the flow field, and turbulent wakes begin to dominate, as shown in Figure 6b,c. The periodicity of the transverse wave train disappeared, and randomness and asymmetry of the turbulence wake appeared. At the junction of the turbulent wake and the transverse wave, there are small breaking water ripples, which is a manifestation of the dissipation of the turbulent wake.
The transformation process from the Kelvin wake to the turbulent wake is more intuitively depicted in Figure 7. As the net momentum increases, the wave pattern changes irregularly, and the nonlinear characteristics are enhanced. The wave height decreased rapidly. These behaviours coincide with the characteristics of the turbulence.

Construction of Composite Sea Surface
To analyse the electromagnetic scattering of the turbulent wake, first, a composite sea surface with turbulent wake had to be constructed.
The sea surface can be formed by the linear superposition of monochromatic waves with different frequencies and phases. By integrating a wave of random phase in the frequency and direction [22], the sea surface height can be expressed as: where ψ nm is the initial phase, A nm = 2S(ω n , φ m )∆ω∆φ, and S(ω n , φ m ) is the power density function of the wave spectrum. The ELH (Elfouhaily-Longuet-Higgins) spectrum proposed by Elfouhaily et al. [23,24] was chosen for this study. ω n is the angle frequency, k n is the wavenumber, and φ m is the azimuth angle. Although the wake and sea surface were obtained, nonlinear wave and wave-current coupling should be considered for the wake and sea surface. Bretherton and Garrett [25] proposed a scheme to deal with the relationship between waves and streams. The sea spectrum S was replaced by a wave-action spectrum H based on the wave-action balance equation. The sea spectrum affected by wave modulation can be described by the following equation [25,26]: where ω 0 is the intrinsic frequency of the wave component, . x = dω 0 /dk + U and x = (x, y) represent the two-dimensional spatial position, and U is the mean current of the wake. . k = dk/dt denotes the propagation speed in wavenumber space, . θ = dθ/dt is the propagation speed in the direction of motion, and θ is the direction angle. Q tot is the source function.
To close the equation, only the effect of wind on the deep sea was considered in this study, and the source term was simplified to three terms [26]: where Q in is the wind input term, Q ds is the dissipative term, and Q nl is the wave-wave nonlinear interaction term. The parameterisation of the wind input term Q in , as proposed by Donelan et al. [27], can be expressed as: where ρ a and ρ w are the densities of air and water, respectively; U w is the wind speed; c is the phase speed; B n (k) and is the spectral saturation which is local in the wavenumber space.
The dissipative term Q ds can be calculated using the following equation [28]: where T 1 (k, θ) is the inherent breaking term and T 2 (k, θ) represents the cumulative effect of short-wave breaking. The wave-wave nonlinear interaction term can be expressed by the Boltzmann integral [29,30]: where G k represents the coefficient between wave components, and subscripts 1 to 4 indicate the four wave components satisfying the resonance conditions. After determining the source function, the new action spectrum H modulated by the wake is obtained from the integral solution, and the new sea spectrum S calculated with Equation (15) is obtained. The wave height of the modulated sea surface can be expressed as: (24) A comparison of the wave heights between the modulated and original sea surfaces is shown in Figure 8. At x = 7 m, the sea surface begins to be modulated by the velocity fields generated by the ship and turbulence. The amplitude of the modulated sea surface was significantly lower in the crest and trough, indicating that the horizontal velocity of the turbulence would inhibit the sea surface wave height. With a decrease in the turbulence, the modulation effect also decreases. Through this process, the problem of excessive wave amplitude caused by the linear superposition of the wake and sea surface can be solved. Finally, the composite sea surface is expressed as: where z wake (x, y, t) is the wave height field containing the Kelvin wake and turbulent wake in the simulation.

Establishment of Electromagnetic Model
To study the contribution of turbulent wake to sea surface electromagnetic scattering and its imaging characteristics, the electromagnetic model should be able to calculate local scattering. The scattering coefficients calculated according to the first-order perturbation method proposed by Furks and Voronovich can be expressed as follows [31,32]: where P, Q = h, v represents the polarization mode, k i and k s represent the unit vectors of the incident and scattering directions, respectively, ε is the complex permittivity of seawater, q l is the projection vector of q = (k s − k i ) on the mean surface z = 0, S(q l ) is the sea spectrum calculated above, and F PQ represents the polarization factor which can be calculated by: where R VV and R HH are the Fresnel reflection coefficients of different polarities, and θ i , θ s , and φ s represent the incident angle, scattering angle, and scattering azimuth angle, respectively.
Considering the tilt modulation of large-scale waves to small-scale waves, a facet scattering model suitable for an arbitrarily tilted surface should be established [33,34]. As shown in Figure 9, (x,ŷ,ẑ) is the global coordinate system, and the local coordinate system is built on an arbitrary inclined micro-rough mean plane O l . The corresponding unit polarization vectors are represented by (ĥ i ,v i ,ĥ s ,v s ) and, respectively, (ĥ l i ,v l i ,ĥ l s ,v l s ). Their relations can be expressed as follows: According to Equations (31), the polarization factorF PQ in the global coordinate system can be expressed as: Given that the capillary spectrum is used to describe the microscale structure of the sea surface, the ELH spectrum is divided into large-scale gravity waves S grav E and small-scale capillary waves S capi E according to the number of truncated waves k cut : where k cut = k/6. According to Equation (26), the sea surface scattering coefficients of arbitrarily microrough tilted facets can be expressed as: For the turbulent wake scattering coefficients of an arbitrarily micro-rough tilted facet, the small-scale capillary spectrum S capi E is replaced by the semi-empirical turbulent energy attenuation spectrum S Tur : where the turbulent energy attenuation spectrum S Tur can be expressed as [35]: where U is the amplitude of the ship speed, L is the length of the ship, l = 1/k 0 is the integral scale, and k 0 is the wavenumber corresponding to the maximum turbulent energy. W(x) is the width of the propeller turbulence. E(k) is the turbulent energy density spectrum [36], which can be expressed as: where γ = 0.25 is the empirical parameter.

Influence of the Turbulent Wake under Different Parameters
To better understand the effect of the turbulent wake on electromagnetic scattering, ∆σ was defined as the difference between the scattering coefficients, and its expression is: where σ Tur PQ (k i , k s ) is the scattering coefficient of the sea surface with the turbulent wake, and σ Sea PQ (k i , k s ) is the scattering coefficient of the sea surface without a turbulent wake. First, ∆σ, that is, the backscattering coefficients on the y = 0 cross-section are shown in Figures 10 and 11. The curve is not smooth because we have extracted the scattering coefficient at each point rather than the total normalised scattering coefficient of the sea surface.  The regularities for the changes in ∆σ under different radar parameters are shown in Figure 10. Only HH polarization is given because the VV polarization is similar. The basic parameters are as follows: radar frequency of 1 GHz, radar incidence angle of 30 • , incidence azimuth angle of 0 • , sea surface wind speed of 3 ms −1 , and J D = 0.097. The variable parameters are shown in each figure. Figure 10a compares the results of ∆σ with different incident angles. The incidence angle has an effect on the Bragg scattering and the tilting effect produced by the turbulent wake. Compared with 30 • , the scattering coefficient ∆σ at 60 • has a larger amplitude and better reflects the characteristics of the surface slope generated by the turbulent wake.
The results for ∆σ with different radar frequencies are shown in Figure 10b. The radar frequency only affects the Bragg scattering of the turbulent wake but does not change the tilting modulation of the turbulence on the sea surface. The amplitude of the scattering coefficient is related to the value of the wave number k 0 corresponding to the maximum turbulent energy. The wavenumber k 0 = 11.75 was selected for this study. As the radar frequency increases, the wavenumber k satisfying the Bragg resonance shifts far away from the wavenumber k 0 corresponding to the maximum peak value of turbulent energy; therefore, the Bragg resonance weakens and the scattering coefficient decreases.
The regularities for the change in ∆σ with different wind speeds and values of J D are shown in Figure 11. The radar incidence angle changed to 60 • and the other parameters are consistent with those shown in Figure 10.
The results for ∆σ at different wind speeds are shown in Figure 11a. The Bragg scattering of the turbulent wake based on the turbulence spectrum is not affected by the wind speed, but the higher sea waves make the turbulent wakes on the composite sea surface less visible, weakening the tilting modulation of the turbulence on the surface and resulting in a significant reduction in ∆σ. Figure 11b compares the results for ∆σ with different magnitudes of J D . The intensity of the turbulence increases with a decrease in the magnitude of J D . The turbulent flow further exacerbates the slope change on the water surface, which causes ∆σ to increase in the HH polarization image.

Electromagnetic Imaging of Turbulent Wake
To represent the electromagnetic characteristics of turbulent wakes more intuitively, the backscattering coefficient distributions of the turbulence are presented in Figures 12 and 13, and the images are restricted to the same grey level scale between −50 and −17.5 dB. The radar frequency was 1 GHz, the radar incidence angle was 60 • , and the incident azimuth angle was 0 • .  The different scattering structures of the wake are clearly shown in Figure 12. The scattering coefficient was strongest at the stern, such that bright scattering points appeared. With an increase in the downstream distance, the scattering coefficient gradually decreased. At the same time, the turbulent scattering coefficient increases with an increase in J D , as discussed previously. Figure 13 shows the scattering coefficient distribution of the turbulent wake at different azimuth angles. For different incidence azimuth angles, the Kelvin wake can only preserve some characteristics. These conditions are frequently observed in SAR images [2]. However, the turbulent wake can maintain a relatively complete structure at different azimuth angles, which is one of the advantages of the turbulent wake compared with the Kelvin wake.
In the ocean, wind and ocean currents may affect the wake characteristics of ships and interfere with ship detection [37]. Referring to the reference (Fujimura et al. [38]), we added a flow field with a velocity of 0.4 ms −1 along the negative direction of the Y-axis to the simulation scene. At the same time, the simulation scene was further expanded, and the ship speed was 5 ms −1 . The radar frequency was 1 GHz, the radar incidence angle was 60 • , the incident azimuth angle was 90 • , and the wind speed was 5 ms −1 .
The distribution of the backscattering coefficients with an extra flow field is shown in Figure 14. And these images are restricted to the same grey level scale between −40 and −17.5 dB. The Kelvin wake and turbulent wake in the simulation results show different responses to the extra flow field. In Figure 14b,c, the red line denotes the outline of the ship wake, the blue line denotes the estimated ship heading, and the white line denotes the actual ship heading. The Kelvin wake is significantly affected by the flow field and the wake drifts. If we estimate the ship heading based on the Kelvin wake extension angle theory, there is bound to be a certain error, as shown in Figure 14c. These small errors in a large target tracking scene can have significant disadvantages.
Compared with the Kelvin wake, the near-field turbulent wake, owing to its fast speed and high momentum, has minimal drift. The estimate of the ship heading based on the Kelvin wake angle has an error of 6 • . However, the estimate of ship heading based on a strong turbulent wake has an error of less than 1 • , as shown in Figure 14b. Therefore, compared with the traditional estimate of ship heading based on the Kelvin wake, the estimate of ship heading based on the turbulent wake has better accuracy.

Conclusions
An accurate simulation method for the near-field turbulent wake was proposed based on a large eddy simulation method. The net momentum transferred to the background flow field determines the development of turbulence. Compared with the simulation based on the semi-empirical energy attenuation spectrum and width formula of turbulence, our method can better deal with the morphological changes of the model in the simulation, and the influence of nonlinear problems on the results. Therefore, the final results obtained abundant structures, including jet flow, turbulent wake, and broken water ripple, which exist in a real marine environment.
Combined with the turbulent energy attenuation spectrum, the electromagnetic scattering characteristics of the turbulent wake were calculated using the two-scale facet mode. First, the effects of different parameters on the scattering coefficient of the turbulent wake were analysed from the Bragg scattering and tilting modulation. Second, the backscattering coefficient distribution of the turbulence wake was calculated. The near-field turbulent wake has a strong scattering field and a certain angle. The structure of the turbulent wake is complete at different incidence azimuth angles. These are more conducive to the observation and recognition of turbulent wakes. Third, under the influence of the extra flow field, the turbulent wake has better stability than the Kelvin wake owing to its fast speed and high momentum. Because of these characteristics, compared with traditional estimates of ship heading based on the Kelvin wake, the estimate of ship heading based on the turbulent wake has better accuracy.

Conflicts of Interest:
The authors declare no conflict of interest.