Aerodynamic analysis on unsteady characteristics of a ducted fan hovering in ceiling effect

Ducted fans have been widely used in unmanned aerial vehicles (UAVs) for various operations due to the high efficiency, low noise and high safety. Proximity effects caused by the confined environment, including the ceiling effect (CE), bring significant interference to the aerodynamic performance of ducted fans, which serves as a main challenge to stability. In this study, the computational fluid dynamic simulation with the sliding mesh technique and the Unsteady Reynolds Averaged Navier-Stokes (URANS) method is conducted to evaluate the ceiling effect. Time-averaged simulation results show that the ceiling effect results in the increase of both rotor thrust and duct thrust. Transient simulation results show that there is a critical distance of 0.4 rotor radius. The ceiling effect leads to little fluctuation in thrust for distances greater than 0.4 rotor radius, but substantial unsteadiness in the flow field for distances less than 0.4 rotor radius. The unsteady movement of vortex structures results in thrust fluctuations with frequencies typically lower than the blade passing frequency of 200 Hz. The fluctuation accounts for up to 40% of the total thrust. The results are essential to the flight controller design and UAV path planning for enhancement of flight stability in ceiling effect.


Introduction
Due to high maneuverability and strong adaptability, unmanned aerial vehicles (UAV) with vertical take-off and landing functions have been extensively used for civil and military applications (Ai et al., 2020;Bingöl & Güzey, 2022;Wei et al., 2022).Ducted fans as potential thrustors of UAVs have received widespread attention (Aboelezz et al., 2022;Hu et al., 2022;Luo et al., 2021).Compared with propeller propulsion, ducted fans can generate additional thrust while reducing runtime noise with the help of the duct.Above characteristics ensure that the ducted fan propulsion system has larger thrust at the same size or a more compact structure at the same thrust, which is suitable for the UAV application scenarios, especially narrow terrain operations.When the ducted-fan drones carry out emergency tasks in narrow zones, such as the Search and Rescue in the event of a disaster (Husain et al., 2022), the complex and changeable environment will cause aerodynamic interference to the ducted fan.In case that the disturbance is strong enough and develops over time, the stability of the ducted fan decreases sharply, and thrust fluctuation or even loss of control may occur.Therefore, it is of great significance to investigate the aerodynamic performance of ducted fan in restricted environment to improve its stability and broaden the application scenarios of ducted-fan UAVs.
Researches on ducted fans began in the 1930s (Stipa, 1932).Early researchers mainly investigated the aerodynamic performance of ducted fans in hover, crosswinds and forward flight by experimental methods, and revealed the impact of structural parameters on the performance.This period of research was reported by Sack et al. (Sacks & Burnell, 1962) and Pereira (Pereira, 2008).In the twentieth century, with the development of computational fluid dynamics, researchers adopted experimental and numerical methods to study the aerodynamic performance of ducted fans, and proposed new optimization configurations and flow control methods.Relevant development in this period was comprehensively described by Zhang et al. (Zhang & Barakos, 2020).In recent years, with the development of electrification technology, distributed electric propulsion technology with high efficiency, low noise and high safety has become an inevitable trend of future aviation development, and the research of ducted fans has sprung up once again.The review article by Qian et al. (2022) pointed out that the aerodynamic characteristics of ducted fans in proximity effects were paid highly attention by researchers in the recent.The proximity effect, brought by the restricted environment, includes ground effect (GE), ceiling effect (CE) and wall effect (WE).Different wall surfaces with proximity effect have different impacts on the aerodynamic performance of the ducted fan.For example, the ceiling or ground mainly induce flow distortion at the inlet or outlet of the ducted fan, and the wall will cause the ducted fan to produce lateral force.
Research focus on the proximity effect has lain on the ground effect.Mi (2020) investigated the effects of rigid ground, water surface and dynamic waves on the ducted fan by CFD simulation.The results showed that the rigid ground had the greatest influence on thrust, followed by the water surface, and dynamic waves would cause unpredictable and periodic changes.Zhao et al. (2022) studied the aerodynamic characteristics of a ducted fan hovering and transition in ground effect.They found that the ground effect was significant when the height above the ground was less than 1.5D (D is the diameter of the fan).The increase in effective angle of attack and ambient pressure led to the increase of rotor thrust, while the decrease of effective angle of attack led to the decrease of duct thrust.The results are consistent with the studies of Ai (Ai et al., 2020), Han (Han et al., 2021), and Deng (Deng et al., 2020).
Compared with the ground effect, the ceiling effect is equally important for UAVs operating in restricted environments, because it will bring energy savings due to the increase of thrust, and in the meanwhile, it may also lead to collision between the UAV and the ceiling.Related research on propellers under ceiling effect has been rather mature (Hsiao & Chirarattananon, 2018;Jimenez-Cano et al., 2019;Powers et al., 2013;Sanchez-Cuevas et al., 2017).Few relevant studies investigated the aerodynamic characteristics of ducted fans in ceiling effect.Han et al. (2019) carried out steady numerical simulation to investigate the performance of the ducted fan in the constrained environment.Results showed that the ceiling significantly improved the aerodynamic performance of the ducted fan at distance below 1.5 fan rotor radius with enhanced figure of merit.To be specific, in ceiling effect, the rotor thrust increased because of the increase in the pressure difference between the pressure side and suction side, and the duct thrust increased due to the decrease in the static pressure at the duct lip.The aerodynamic efficiency of the ducted fan was therefore improved as a result of the increase in total thrust and the decrease in power consumption.Ai et al. (2021) established a thrust mechanism model of coaxial ducted fan under ceiling effect based on the momentum theory and the blade element theory.The model took the vertical climb rate, induced velocity of coaxial rotors, and the velocity interaction between the coaxial rotors into consideration, which was in good agreement with the experimental results.In addition, they studied the ceiling effect under different vertical velocities and vertical heights through numerical calculation methods.They found that as the ceiling height or the climbing velocity decreased, the total thrust increased due to the expansion of the low-pressure area at the duct lip.
It has been reached a consensus that the range of ceiling effects on the aerodynamic performance of ducted fans is within 4R (R is the radius of the fan rotor), and significant effect is observed below 1.5R approximately.As the distance from the ceiling decreases, the duct thrust, rotor thrust and total thrust are all improved due to the ceiling effect.However, current researches on the ceiling effect were conducted by steady simulations without exception.The unsteady characteristics of flow field caused by ceiling has not been captured and analysed, especially at an extremely close distance from the ceiling.Ceiling effect acts as a threat to UAV operations, not only because it will lead to the shifts in thrust, but also bring unsteady changes, which has not been addressed in any prior research.
In this paper, the unsteady CFD method with sliding mesh technique is applied to investigate the timeaveraged and transient aerodynamic characteristics of ducted fan hovering in ceiling effect, and the physical principles behind are compared and analysed.The main contributions of this paper include two parts.On the one hand, we verify the trend of ducted fan thrust under ceiling effect, and thoroughly analyse the causes accounting for the change.On the other hand, we reveal for the first time that there is a critical distance for the ducted fan hovering in ceiling effect.When the distance is greater than the critical value, the ducted fan is in stable operation, and the thrust changes regularly with the decrease of distance without significant fluctuation.When the distance is less than the critical value, the ducted fan enters an unstable zone, and the thrust fluctuates greatly out of accordance with the existing model predictions.We hope this study provides a reference boundary for stable operations in terms of the path planning of ducted fan UAV in restricted environments.
The remaining parts are organized as follows: Section 2 introduces the geometric model and numerical simulation methods.Section 3 presents the experimental validation of the CFD method.Section 4 compares and analyses the time-averaged and transient characteristics of the ducted fan in ceiling effect.Section 5 draws the conclusions and provides future perspectives.

Geometric model
Overall, the ducted fan with short duct, which is effective in reducing the air resistance and additional weight and is easy to be integrated with the fuselage, is a widely used configuration with application prospects.The investigated ducted fan, as is illustrated in Figure 1(a), has three blades with uniform spacing.The radii of propeller and duct are R = 103 mm and 146 mm respectively, and the duct length is h = 50 mm.The tip clearance measures δ = 0.91 mm. Figure 1(b) displays the specific parameters of the fan.The pitch angle varies almost linearly from 40.63°at the hub to 20.35°at the tip.The maximum chord length is observed at around 0.63R.

Domain division and mesh generation
Figure 2 shows the scheme of the entire computational domain.The entire computational domain is divided into a static domain (including b1 and b2) and a rotating domain (c).These two domains are connected by two pairs of interfaces, namely Interface_up & Inter-face_top and Interface_bottom & Interface_down.The static domain is a cylinder with a diameter of 20D and a height of 13D, which is used to simulate the flow field outside the ducted fan.The top surface of the cylinder is used to model the ceiling.The rotating domain contains the inner field of the ducted fan and some extended parts above and below the inlet and outlet to investigate the internal flow, of which the rotational speed is set to be identical to that of the ducted fan.For different cases, the distance of the ducted fan from the ceiling changes accordingly, leaving the rest of the settings unchanged.
To cover all important flow regions including surface boundary regions, blade tip regions and regions near the ceiling in particular, the structural and hexahedron grids were delicately constructed.The construction of the mesh in static domain was based on ANSYS ICEM CFD.Since the flow near the ducted fan and the ceiling is complex and the flow in the domain away from them is the opposite, we did not adopt a uniform structural grid distribution but a scheme that the grid size changed exponentially, which implied that grids close to the ducted fan and ceiling were denser, while grids away from the ducted fan were sparser.This was achieved by setting the distribution laws of grid nodes to Geometric 1 and Geometric 2 in ICEM.Such grid arrangement could accurately and effectively simulate the complex aerodynamic effects of the interference between the ducted fan and the ceiling.The construction of the mesh in rotating domain was based on NUMECA AutoGrid5, which is a powerful software for processing the mesh of turbomachines such as axial compressors.Figure 3 shows the resulting mesh on the surface of the ducted fan, especially mesh details for the blade tip and root regions.The minimum skewness of three-dimensional mesh was over 19, maximum aspect ratio was close to 600 and maximum expansion ratio was 3.85.The grid quality could well meet the calculation requirements.The overall grid number was approximately 12 million, including 6 million for each of the rotating and static domains.

Computational method
In this paper, the flow can be regarded as incompressible because the blade tip Mach number is close to 0.1 when   the fan rotating speed is = 4000 r/min.From this, we adopted the pressure-based solver in a mature commercial software ANSYS Fluent for numerical investigations (Halfi et al., 2020;Menni et al., 2021).The classical shear stress transport (SST) k-ω turbulence model was chosen because it has high computational accuracy for predictions of free shear turbulence, boundary layer turbulence and moderate separation turbulence (Shadab et al., 2022), which satisfies the simulation of ducted fans.The robust Coupled algorithm is used to solve the continuity equation and momentum equation simultaneously, and the second-order upwind scheme is used to improve the calculation accuracy.The steady initial field is calculated by multiple reference frame (MRF) method, and the transient flow field is calculated by sliding mesh technique.
Boundary conditions for simulation of a ducted fan in ceiling effect are implemented according to Figure 2(a).In the static domain, the ceiling is set as the non-slip wall, and the rest of the outer boundary 'Far' is set as the velocity inlet with velocity amplitude zero.There exist four pairs of interfaces in the static domain, among which the two pairs mentioned above are interfaces between the static and rotating domains.In the rotating domain, the inner duct wall is set as the non-rotating wall in the absolute coordinate system.

Validation of turbulence model
For verification purposes, we conducted bench experiments in Chongqing Innovation Center of Beijing Institute of Technology.Figure 4 shows the composition of the entire experimental bench, including the aluminium bracket for mechanical support, the tested ducted fan, and associated controllers and sensors.Four types of sensors are installed on the bench to collect thrust, torque, motor speed, and the current and voltage of the driving motor in real time.The measuring range and accuracy of each sensor are listed in Table 1.
Figure 5 shows a schematic diagram of the data acquisition system.The force sensor, torque sensor, and current and voltage sensor are connected to the data collector ADS1256, which is mounted on the main control board stm32f103 with the rotary encoder.The serial communication between the main control board on the test bench and the personal computer is implemented by 485 transmission modules.We conducted three repeated hover experiments to reduce random errors, and the uncertainty is shown in Table 2.Among all the operating conditions, type B uncertainty caused by instrumental error is dominant especially at low rotating speeds.The relative uncertainty in cases above 3000 r/min is below 20%, indicating that the experimental results have good reliability.The simulation results based on the URANS method are compared with experimental data in Figure 6.It can be seen that the computed rotor thrust overpredicts about 20% more than the measurement, while the computed duct thrust underpredicts about 10%.The error in total thrust predictions is below 5%.The torque agreement is relatively good, with an error of less than 5%.Given the discrepancies between tests and numerical modelling, the agreement between numerical and experimental results is acceptable.The overall error may be caused by the differences between the numerical and experimental setups, for example, the aerodynamic interference of the test bench support and the driven module installed at the fan bottom.
We considered the influence of different turbulence models on the results as well.Generally speaking, the large eddy simulation (LES) and direct numerical simulation (DNS) can improve the calculation accuracy, but also bring huge computational cost.As an alternative, URANS method can provide fast aerodynamic prediction though the calculation accuracy is doubtful.For illustration purposes, we compared our results with a hybrid method, Delayed Detached Eddy Simulation (DDES), which combines the advantages of the   high precision of LES and the low computational cost of RANS.SST turbulence model was adopted for both methods.Table 3 shows that the URANS method can well predict the performance of the investigated ducted fan because the differences between two methods are below 1.5%.Figure 7 demonstrates that compared with the DDES method, the URANS method is accurate in calculation of the velocity and pressure distributions in the internal flow field despite little difference in calculating the transient evolution process of rotor wake.Based on the above mentioned, we can conclude that it is reliable to use the URANS method to study the performance of the ducted fan and reveal flow characteristics.

Validation of mesh and time step independence
To verify grid convergence, we chose the hovering condition supported by experimental data and altered the mesh size of stationary and rotating domains.Three mesh grids with different qualities, namely coarse (6 million cells), medium (9 million cells), and fine (12 million cells), were generated for comparison.As is shown in Table 4, the calculated thrust coefficient using URANS in three cases were 5.4695 N, 5.5946 N, and 5.6363 N, with corresponding error compared with fine mesh being 2.96%, 0.74%, and 0%.This indicates that the medium and fine quality mesh had enough computational accuracy, and we  selected the finest mesh grids with 12 million cells in the follow-up study.Time step size is an important parameter in transient simulation, which has a significant effect on accuracy and performance.It must be not only small enough in theory to solve the time-related features, but also balanced with the computational cost.We carried out the time step independence verification on the basis of the finest grid.According to our experience, we chose four different time steps as dt = 0.00005 s, 0.0001 s, 0.00025 s, and 0.0005 s, with corresponding rotating angles of the ducted fan per step being 1.2°, 2.4°, 6°and 12°.The test case of the ducted fan hovering at a distance of z = 0.5R from the ceiling was selected, and the results are shown in Figure 8.It is apparent that that dt = 0.00005 s is accurate and effective in capturing the transient change of the total thrust, so the time step of dt = 0.00005 s is used in the subsequent research.

Results and discussions
Ducted fans are efficient in generating thrust in the hovering state due to the additional thrust produced by the duct.When the ducted fan operates near the ceiling, however, the ceiling turns the direction of incoming flow at the inlet and leads to changes in the performance of both the duct and the rotor.The impact mainly lies on the thrust and torque in the vertical direction rather than in other directions due to the symmetric configuration.In this section, the influence of the ceiling on the aerodynamic characteristics of a ducted fan is numerically investigated based on the sliding mesh technique.Related operating parameters and physical parameters in this study are listed in Table 5.
In the following parts, Section 4.1 aims to reveal the relationship between the thrust and torque characteristics of the ducted fan in ceiling effect and the distance from the ceiling.Section 4.2 provides an explanation for the time-averaged and transient changes in thrust from the perspective of three-dimensional flow fields.

Aerodynamic performance in ceiling effect
Previous studies have shown that, the aerodynamic performance of ducted fans is more sensitive to wall constraints compared with that of propellers.Taking the ground effect as an example, the ground effect has a significant impact on the propeller in the range of 0 to 2R, while the range changes to 0 to 4R for the ducted fan (Han et al., 2021;Matus-Vargas et al., 2021;Mi, 2020).
The influence of ceiling effect on thrust characteristics of the ducted fan is illustrated in Figure 9, where T ICE denotes the thrust in ceiling effect and T OCE denotes the thrust out of ceiling effect.Experimental data of cases with z/R = 0.4, 0.5, 1, and 2 well supports the simulation results.The figure reveals two major features of the ducted fan performance.
(a) Time-averaged thrust characteristics The ceiling effect is significant for a distance up to 4R, and the aerodynamic interaction can be negligible when the distance exceeds 4R.With the decrease of distance, the duct thrust increases significantly, the fan rotor thrust increases slightly, and the total thrust shows a tendency towards growth.The maximum total thrust ratio 1.74 is achieved at z/R = 0.4.In order to quantitatively characterize the influence of the ceiling on the total thrust of the ducted fan, a classic prediction model given by the reference (Jardin et al., 2017) is adopted, which writes  This model shows that the total thrust changes with distance according to the power law in ceiling effect, which is effective in the flight controller design for better stability of the ducted fan aircraft.When the distance is less than 0.4R, the ducted fan enters an unstable operating state and the total thrust starts to change irregularly.Large thrust variations indicate a serious deterioration in the ducted fan's performance in this situation.

(b) Transient thrust characteristics
In addition to the time-averaged characteristics, it's indicated in Figure 9 that the ceiling causes unstable thrust fluctuations of the ducted fan, especially at short distances.However, previous studies have made no effort in the investigation and interpretation of this transient phenomenon.
According to Figure 9, the error bar in the left panel indicates the extreme value of the thrust change, and the standard deviation (SD) in the right panel indicates the degree of dispersion of the thrust from the mean value over time.Thrust can be regarded as a constant when the distance from the ceiling is more than 0.4R.A threshold is reached at z/R = 0.4.The ceiling significantly interferes with the ducted fan's aerodynamics and causes a major fluctuation in the total thrust for distances below 0.4R.The maximum fluctuation occurs at z/R = 0.2, where the thrust ranges of the rotor and the duct are 1.3001 N and 1.7296 N, and corresponding SDs are 0.2827 and 0.4029 respectively.The range in total thrust is 2.7794 N, accounting for 40% of its average value, and the standard deviation is 0.6459.Moreover, it's clear by comparison that the duct is more susceptible to the ceiling effect and contributes more to fluctuations of total thrust.
For illustration purposes, Figure 10 compares the frequency spectrums of thrust at a different distance from the ceiling by Fast Fourier Transform (FFT).The x-coordinate represents the frequency of the thrust fluctuation, in which the blade passing frequency is 200 Hz at the speed of 4000 r/min.The y-coordinate refers to the amplitude of rotor thrust and duct thrust in the unit of Newtons.Low amplitude in the frequency spectrums in Figure 10(c,d) indicates that there are few fluctuations in thrust for cases above 0.4R.As the distance falls below 0.4R, thrust fluctuations are seen with a variety of frequencies that are typically lower than the blade passing   frequency.The predominant frequency fluctuations are at 6.66 Hz and 13.32 Hz.
The torque coefficient and figure of merit of the ducted fan are evaluated in Figure 11 in relation to the ceiling effect.As can be seen, the simulation results are in good agreement with experimental data, where the error in torque predictions is below 6% and that in FM predictions is below 10%.The torque is substantially unaffected by ceiling effect and remains close to the torque coefficient 0.01450 in hover except for slight variations at 0.4R and 0.5R.Besides, the FM increases from 0.85 at z/R = 5 to 1.79 at z/R = 0.4, which follows the power law with distance.This is largely attributed to the increase in total thrust induced by ceiling effect.

Flow characteristics
Changes in thrust are essentially due to the variations in the ducted fan flow field induced by the ceiling.For illustration purposes, this section investigates and analyses the time-averaged and transient flow characteristics in ceiling effect.

Time-averaged flow features
Figure 12 illustrates the velocity fields in the form of streamlines on a colored-contoured background and the tip flow field features.Four cases of z/R = 0.5, 1, 2, and hover are chosen for detailed comparison at a constant rotating speed = 4000 r/min, where velocity magnitudes are normalized by the blade tip velocity.For the case in hovering state shown in Figure 12(d), the air is sucked in from the open atmosphere, accelerated by the fan and creates a high-speed vertical wake below the outlet.Flow separation regions are observed at the top of the duct and the diffuser section due to the turning of radial airflow and adverse pressure gradient.In contrast, Figure 12(a-c) compares the influence of the ceiling on the aerodynamic performance of the ducted fan.As can be seen, the ceiling limits the direction of incoming flow, resulting in its acceleration in the entrance section.The flow separations at the top of the duct and in the diffuser section have somewhat expanded due to ceiling effect.Besides, obvious ceiling vortex gradually forms between the fan and the ceiling as the distance decreases to 0.4R in Figure 12(a).The stagnant area reduces the effective flow area of the fan and causes potential aerodynamic instability.
To quantitatively study the ceiling effect on flow, we collect the total mass flow rate through the ducted fan and total pressure rise with different distances in ceiling effect in Figure 13(a).The mass flow rate is calculated at the inlet and the total pressure rise is calculated to be the area-weighted average value between the outlet and the inlet.It seems that the ceiling and the duct together form an 'inlet throttle valve', and the proximity of the two acts as if the valve were closed.The total pressure rise and mass flow rate of the ducted fan are maintained at a constant value of 90 Pa and 0.46 kg/s due to the constant rotating speed when the ducted fan stays far from the ceiling.When the distance is reduced to below 0.4R, the ceiling seriously inhibits the inlet airflow due to significant ceiling effect, resulting in a sharp decrease in the mass flow rate and a decrease in the total pressure rise.The mass flow rate drops to less than 15% of its original value when the distance reaches 0.2R.In this case, the throttling effect becomes remarkable, and the ducted-fan performance is strictly inhibited.
Further, Figure 13(b) shows the mass flow distributions in ceiling effect, where the mainstream area is equally divided into five regions in the radial direction, and the axial velocity is normalized by the blade tip velocity.As can be seen, the major flow travels through the tip region, and the average axial velocity there is about 30% larger than that in the hub region at z/R = 2.With the decrease of distance, the ceiling effect results in the decrease of axial velocity in all regions, especially in the hub and tip regions, which are primarily responsible for the reduction in mass flow rate.The reverse flow begins to appear in the hub and tip regions when the  distance is below 0.3R.Only a small portion of airflow passes through the middle of each blade passage in this situation.
The change of axial velocity distributions induced by the ceiling has an impact on the blade attack angle, resulting in a change of its aerodynamic performance.Figure 14 shows the chordwise distributions of static pressure coefficients at 20%, 50% and 90% blade heights in ceiling effect, where the airfoil profiles at different blade heights are displayed in Figure 14(a).Through the comparison between Figure 14(b-d), it can be found that the low-pressure peak at the blade leading edge expands, indicating that the blade angles of attack at all blade heights increase with the enhancement of ceiling effect.The increase of angle of attack at the blade root is larger than that at the blade tip which is due in part to the significant reverse flow shown in Figure 13.However, compared with the blade root, the blade tip contributes more to the increased rotor thrust as the distance from the ceiling decreases.The improvements in thrust from z/R = 2 to z/R = 0.4 at 50% and 90% blade height are 30% and 22% respectively, while the thrust at 20% blade height even decreases by 15% as a result of large flow separation at the suction side induced by the excessive angle of attack.
In terms of the duct, Figure 15 reveals the impact of ceiling on the flow characteristics along the duct surface, where the separation points are established by the zero points of wall shear stress (WSS).According to Figure 15(a), there are two main flow separation regions on the surface, one at the top of the duct and the other at the inner duct.The flow separation region at the duct top is beneficial to the duct thrust because it contributes to a low-pressure zone.The radial displacement of point S in Figure 15(b) indicates that the ceiling effect expands the flow separation region, which enlarges the low-pressure zone and helps to increase the duct thrust.
Additionally, the flow profile near the inner duct surface is complicated due to the aerodynamic interference of the blade wake.The blade wake makes the airflow reattach to the duct surface below the flow separation at the duct lip and separate again at the duct diffuser, dividing the separation region into two zones.The upward movement of points X, Y and Z demonstrates that the ceiling effect expands the separation region at the inner duct as well, and the moving distance of point X is greater than that of point Z.The interaction mechanism between blade wake and flow separation near the duct surface requires further study.
The flow characteristics on the duct surface directly result in changes of the pressure distributions.For illustration purposes, Figure 16 shows the time-averaged static pressure distributions around the duct.The pressure of section A-B is slightly lower than the atmospheric pressure due to the flow separation caused by the airflow turning at the inlet of the ducted fan, and a maximum value point appears near point B. The significant pressure drop in the B-C section forms a low-pressure zone as a result of increased velocity at the duct lip according to the Bernoulli equation.The pressure rises briefly in the transition section (near point C) owing to the blade wake and then drops to form a second low-pressure zone in the C-D section due to flow separation in the diffusion section.In the D-A section, the pressure remains basically the same as the atmospheric pressure.Based on the geometric configuration of the duct, it's clear that the duct thrust is generated by the pressure difference between the D-E section with atmospheric pressure and the A-B, B-C sections with low pressure.When the ducted fan hovers near the ceiling, the ceiling accelerates the airflow passing by and enhances the flow separation according to Figures 12 and 15, resulting in the decrease of the static pressure in the A-B section.The pressure difference between D-E section and A-C section increases, and the duct thrust increases.When the distance is reduced to less than 0.4R, the flow rate and velocity is greatly limited in the internal flow field as shown in Figure 13, and as a consequence, the low-pressure zone in the A-B and B-C sections basically disappears.In this case, the average thrust generated by the duct is only around 50% of the thrust in hovering state.

Transient flow features
In the previous section, Figure 12 shows that the ceiling effect contributes to the increase in the number and intensity of vortex structures.The unsteady vortex may be the dominating reason accounting for the thrust fluctuation.Therefore, this section aims to reveal the development process of vortex structures from the perspective of a three-dimensional transient flow field.
Figure 17 visualizes the vortex structures of the flow field in ceiling effect by extracting iso-surfaces of Q-criterion = 1,000,000.It can be seen that vortex cores mainly appear inside the duct, demonstrating that the internal flow is essential to the ducted fan performance in ceiling effect.Figure 17(c,d) show typical flow fields above 0.4R., where the flow separation at the top of the duct and the tip leakage vortex can be observed.As the distance decreases, the ceiling effect leads to slight enhancement of the above vortices and the appearance of vortices near the hub region.Figure 17(a,b) show the unstable flow patterns induced by intense ceiling effect below 0.4R.On the one hand, the flow separation at the duct top is gradually weakened and finally disappears due to the decrease in mass flow rate.The tip leakage vortex expands upward as a result of abundant reversed flow.On the other hand, a large number of newly emerging large-scale vortices are observed in the flow field.These vortices, which are generally distributed in the blade root region and tip region, cluster in the position close to the blade surface.They form dead water zones and leads to the blockage of the mainstream.When the distance decreases further from z/R = 0.3 to z/R = 0.2, there is significant separated vortex on the suction side of the blade.Overall, unsteady vortices around the rotor account for the unsteady behaviour of the ducted fan in ceiling effect, and lead to the increase in the amplitude of thrust fluctuation.
For a clearer view of the development of the vortex structures, Figure 18 extracts and compares the circumferential motions of vortex structures inside the duct.The z-plane in the middle of the blade passage is chosen as the interception position.Figure 18(d) demonstrates that there are a pair of tip leakage vortices located at the blade tip and a pair of separated vortices located at the hub, which rotate synchronously with the fan blade (Fang & Tachie, 2019;Graziani et al., 2018).Figure 18(b,c) manifest that the ceiling effect is responsible for the enhancement in the hub separated vortex and tip leakage vortex.It is noteworthy that the vortex originally attached to the hub gradually falls off to the blade passage, forming a pair of large-scale counter-rotating vortices in the centre as the distance approaches 0.4R. Figure 18(a) shows that the internal flow in the ducted fan becomes entirely turbulent flow, and there is obvious stall at the blade root at a distance of 0.3R.The circles in the figure mark the moving trajectory of the hub vortices in one revolution.The vortex structures shed from the hub, interact with the passage vortex, move to the blade tip along the pressure side as a consequence of centrifugal force, and finally merge with the tip leakage vortex.The unsteady behaviour of vortices causes thrust fluctuation, which is challenging to forecast due to its weak periodicity.Operations under such conditions should be avoided as possible.
Vortex motion will cause local fluctuations in the flow field, and the intensity and frequency information of vortex structures can be obtained by using spectrum analysis.In the numerical simulation, four pressure monitoring points arranged at blade heights of 10% and 95% are used to quantitatively characterize the unsteady features in ceiling effect.For each blade height, the placement of monitoring points is shown in Figure 19.The monitoring points are created in a rotating reference frame and rotate synchronously with the ducted fan.For investigation purposes, the static pressure at each point is collected and processed from the 5th revolution to the 15th revolution.
Since the disturbance can propagate in the whole subsonic flow field, monitoring points arranged in any position are effective in reflecting the frequency information of the fluctuation despite the low amplitude when far from the disturbance source.In order to find an ideal position for monitoring point arrangement, Figure 20 shows the power spectral density at different monitoring points in the hub region at z/R = 0.2.According to the figure, Point 2 is the most appropriate monitoring point because its highest power of pressure fluctuation indicates it is the closest one to the main vortex structures.Thus Point 2 is adopted for the following studies.hub region becomes greater than that in the tip region due to the strong interaction between the reverse flow induced by the ceiling effect and the mainstream.
In summary, when the distance from the ceiling is greater than 0.4R, small-scale vortices are mainly concentrated in the hub and blade tip region, and their intensified motion due to ceiling effect leads to an increase of thrust fluctuation.When the distance is lower than 0.4R, however, the small-scale vortices develop into large-scale vortices in the whole flow field, resulting in significant thrust fluctuations.

Conclusions
Based on URANS and the sliding mesh technique, this paper investigates the unsteady aerodynamic performance of the ducted fan in ceiling effect, and compares and analyses the influence of ceiling on the time-averaged and transient characteristics.Conclusions are drawn as follows.
(1) Time-averaged characteristics When the ducted fan is hovering in ceiling effect, the ceiling forces the airflow to pass along the radial direction and turn axially to the rotor-disk plane.The ceiling limits the flow through the ducted fan as the distance from the ceiling decreases, resulting in a decrease in the axial velocity of the mainstream and an increase in the blade angle of attack.The rotor thrust therefore increases.In the meanwhile, the low-pressure region in the duct inlet section is expanded due to the acceleration of the airflow, resulting in the increase of the duct thrust.The total thrust shows an upward trend.
(2) Transient characteristics The stable flow patterns are maintained when the ducted fan hovers far away from the ceiling though flow separations at the hub, the duct top, and the diffuser are enhanced by ceiling effect.However, the ceiling leads to a significant reduction in mass flow rate and arouses substantial unsteadiness in the flow field inside the duct when the distance is less than the critical value 0.4R.Unsteady vortices are primarily located between the fan rotor and the ceiling, which shed from the rotor surface, merge with each other and travel between blade passages.The unsteady movement of vortex structures results in thrust fluctuations with a variety of frequencies that are typically lower than the blade passing frequency 200 Hz.
Current research comprehensively reveals the aerodynamic characteristics of the ducted fan in ceiling effect and the potential range (z/R ≥ 0.4) for stable operations, which provides the key theoretical basis and detailed data reference for the operation of UAV in ceiling effect.However, our research does not consider the dynamic characteristics of single ducted fans in proximity effect, nor does it consider the mutual aerodynamic interference between several ducted fans applied to multi-ducted fan UAV.The aerodynamic characteristics under these conditions will be much more complicated, which requires further exploration and investigation.

Figure 1 .
Figure 1.Description of the investigated ducted fan: (a) overall scheme, (b) blade chord and pitch angle distributions.

Figure 2 .
Figure 2. The boundary conditions and computational mesh adopted: (a) overall settings, (b) static zone settings, (c) rotating zone settings.

Figure 3 .
Figure 3. Detailed computational mesh on the surface of the ducted fan.

Figure 4 .
Figure 4. Schematics diagram of the test bench.

Figure 5 .
Figure 5. Data collection system of the test bench.

Figure 6 .
Figure 6.Comparisons between computational and experimental results of the ducted fan in open space at different rotating speeds: (a) rotor thrust, duct thrust and total thrust, (b) rotor torque.

Figure 7 .
Figure 7. Flow field comparisons between the DDES and URANS simulations in hover: (a) velocity magnitude profiles, (b) pressure profiles.

Figure 8 .
Figure 8.Time step independence test for the model in ceiling effect (z = 0.5R, = 4000 r/min).

Figure 9 .
Figure 9. Thrust characteristics of the ducted fan hovering in ceiling effect at 4000 r/min: (a) value and SD of rotor thrust, (b) value and SD of duct thrust, (c) value and SD of total thrust.

Figure 11 .
Figure 11.Torque coefficient and figure of merit of the ducted fan in ceiling effect.

Figure 13 .
Figure 13.Mainstream velocity characteristics in ceiling effect: (a) total mass flow rate and total pressure rise versus distance, (b) mass flow distributions in different regions.

Figure 15 .
Figure 15.Flow characteristics on the duct surface: (a) diagram of wall shear stress at z/R = 2, (b) positions of flow separation points at different distances from the ceiling.

Figure 16 .
Figure 16.Time-averaged pressure coefficient distribution around the duct section for different positions below the ceiling.

Figure 18 .
Figure 18.Diagram of normalized z-vorticity in the middle of the blade passage in one revolution for different positions below the ceiling: (a) z/R = 0.3, (b) z/R = 0.4, (c) z/R = 0.5, (d) z/R = 1.

Figure 19 .
Figure 19.Schematic diagram of distributed monitoring points for pressure signals.

Figure 21
Figure 21 shows the pressure power spectrum of Point 2 at different blade heights in different cases, where the peak signal is aroused by vortex motion.At z/R = 1 and z/R = 2, there is minimal fluctuation in pressure in the hub and blade tip regions, and their frequencies lie mainly within 100 Hz.The high-frequency signal is weak.With the enhancement of the ceiling effect, the movement of the hub and tip vortices is enhanced, and the

Figure 20 .
Figure 20.Power spectral density of monitoring points of 10% blade height at z/R = 0.2.

Table 1 .
Details of sensors used in the experiment.

Table 2 .
Uncertainty for different parameters in repeated hover tests.

Table 3 .
Comparisons of computational results based on DDES and URANS methods ( = 4000 r/min).

Table 4 .
Grid convergence tests for the free space model ( = 4000 r/min).

Table 5 .
Operating parameters and physical parameters used in the study.