The Theoretical Research for the Rotor/Fuselage Unsteady Aerodynamic Interaction Problem

J. Aerosp. Technol. Manag., São José dos Campos, Vol.8, No 3, pp.281-288, Jul.-Sep., 2016 ABSTRACT: Based on coupled unsteady panel/free-wake method, a universal analysis model was established, which provides a good prediction for the rotor/fuselage unsteady aerodynamic interaction. Considering the deficiencies of the traditional time-marching rotor free-wake algorithms, notably on stability and efficiency, the CB3D algorithm with 3rd-order accuracy is proposed. For solving the problem that part of the wake vortices may penetrate the fuselage, a “material line” rectification method with 3rd-order accuracy is proposed. An analysis for the model accuracy was then conducted to validate the accuracy of the new model, and a comparison against the available experimental data is performed. The simulated results show a good agreement with these experimental data. With the new model, several simulations are conducted for the typical rotor/fuselage aerodynamic interaction, and the results are analyzed.


INTRODUCTION
The rotor unsteady aerodynamic interaction is one of the research hotspots of rotorcraft.It not only has a great influence on aerodynamic general layout designing of helicopter, but also plays a more and more important role in helicopter flight dynamic modeling and simulation.Many scholars from different countries have carried out researches on this issue.This paper focuses on the prediction model about the typical aerodynamic interaction of rotor-fuselage interaction.
Helicopter rotor flow-field simulation methods can be divided into 2 categories: Computational Fluid Dynamics (CFD) and free-wake.Serious non-physical dissipation happens when the CFD method is used to simulate the blade vortex interaction.In addition, a huge amount of computation is needed for CFD method to accurately predict the rotor flow-field.By contrast, the problems mentioned will not happen in free-wake model, so it is more suitable for the calculation of unsteady aerodynamics interaction in engineering application.
In the free-wake prediction approaches, comparing to relaxation iteration algorithm (Bagai 1995), the time-stepping free-wake method can not only simulate the steady flight, but also predict the flow-field of a transient flying rotor.However, at present, the widely used PC2B algorithm (Bhagwat 2003) is a kind of predict-correction method.It has a relatively lower efficiency when considering the influence of other boundary such as the fuselage.Recently a one-step CB2D algorithm was proposed (Li and Chen 2012) and it has a higher efficiency.But divergence of free-wake structure may happen when it is used to simulate the rotor interaction problems.So, considering the deficiencies of traditional CB2D algorithm on stability, an explicit CB3D algorithm with 3rd-order accuracy is proposed in this paper.CB3D eliminates the anti-dissipation effect caused by 2nd-order error terms in CB2D algorithm, which could improve the numerical stability, analysis accuracy and computational efficiency in the simulation of rotor aerodynamic interaction problems.
In the research of rotor-fuselage interaction, the literature (Bi and Leishman 1990;Leishman and Bi 1990;Smith and Betzina 1986) introduced several experiments on this issue, and the pressure distribution on fuselage surface was measured.
In the theoretical analysis field, some papers (Rand 1988;Lorber and Egolf 1990;Crouse et al. 1992) used a prescribed wake method to predict the interaction flow-field.In view of the particularity of rotor-fuselage aerodynamic interaction, many empirical and ideal components were added into the analytical model, which limited its universal property and accuracy.The literature (Mavris et al. 1989;Xu et al. 2002;Quackenbush et al. 1999) used the more advanced free-wake method to simulate the rotor flow-field.In the free-wake model, only a few used tip vortices (Katz 2001) and others used a full-span free-wake model.For solving the problem that some of the wake vortices discrete points may move not physically into the fuselage surface, Lorber and Egolf (1990) placed the vortices near the fuselage as half-circle vortex segments around it, and some scholars placed the wake points which move into the fuselage in a mirror-image location.These rectification approaches lacked theoretical basis and then had a low rectification accuracy, which may lead to a less accurate aerodynamic interaction analysis model.
A coupled unsteady panel/time-stepping free-wake method is established to predict the helicopter rotor interaction flowfield.The panel method is used to simulate the influence of the helicopter fuselage.A "material line" rectification method is proposed to solve the problem that part of wake vortices may penetrate the fuselage surface.To validate the accuracy of this model in rotor interaction, a comparison of free-wake structure, rotor aerodynamic performance and the pressure on fuselage against the available experimental data is carried out, and the simulated results agree well with the experimental data.

ROTOR AND FREE-WAKE AERODYNAMIC MODEL
The circulation of the blade bound vortex and the aerodynamic performance are calculated using the lifting surface theory.A series of vortex lattices are divided into the spanwise and chordwise direction, respectively.The sketch of blade lifting surface is illustrated in Fig. 1.From the experimental observations, the trailing vortices will roll up to a single tip vortex filament within a few chords behind the rotor blades (about 30 -60 degrees behind blade trailing edge).The tip vortex has been documented to be the most dominant structure in the rotor flow-field.The inboard trailers have a negligible influence on the rotor performance.
But in rotor-fuselage interaction flow-field, since most part of the slim fuselage is just below the rotor blade root, the free-wake model only keeping tip vortex will be less accurate in the prediction of the pressure on fuselage surface.Thus more sophisticated model is required to determine the rotor/fuselage interaction.Provision has been made to use a full-span free-wake model including 5 -7 trailing vortices, as shown in Fig. 2. Circulation of the vortices is determined by the maximum circulations on both sides of the vortex release point.In the free-wake method, the vortices are force-free and convect with the local air velocity.The motion is governed by the following equation: where: r is the position vector of a wake collocation point; (1) The Theoretical Research for the Rotor/Fuselage Unsteady Aerodynamic Interaction Problem ψ is the azimuth angle; ζ is the vortex age; Ω is the rotor rotational speed; u is the wake velocity; u ∞ is the free stream velocity; u freewake-induce (ψ, ζ) and u fuselage-induce (ψ, ζ) are the velocities induced by free-wake and fuselage, respectively, on vortices points.

CB3D NUMERICAL INTEGRATION SCHEME
The original CB2D algorithm is shown in Eq. 2, where r l,k is the position vector of the k th vortex discrete point on a vortex line shed by blade located at the l-th azimuth angle; u(r l,k ) is the velocity of the vortex point, and the term with a γ coefficient is the artificial dissipative one proposed by Li and Chen ( 2012): (2) Equation 2 is obtained based on the assumption that Δψ = Δζ, and its modified equation (expanded at the point Since the principal error terms on the right side of the equation are the 2nd-order ones, the CB2D is a 2nd-order accuracy algorithm.
According to Eq. 3 we can transform the error terms As analyzed in the literature (Li and Chen 2012; Warming and Hyett 1974), the terms containing the 2nd-order partial derivatives r ζζ and r ζ 2 , which mainly play an anti-dissipation role in the free-wake iteration, because of their coefficient variables u r , u rr , often become minus.So in the literature the artificial dissipative term γΔψ 2 r ζζ (γ ∈ (0, Δψ)) with a constant positive coefficient is added to strengthen the numerical dissipative effect, improving the stability of the numerical scheme.In the present study, an approach is proposed to further improve stability of the scheme as follows.
Note the following existing relations: The 2nd-order terms on the right side of the equations are just the error terms in Eq. 3. If the polynomials on the left side of Eq. 4 are added to the right side of Eq. 2, according to the coefficients of the corresponding terms in Eq. 3, then a new algorithm is obtained: The modified equation (expanded at the point (l + 1/2, k + 1/2)) of Eq. 5 is as follows: Given that γ ∈ (0, Δψ) (Li and Chen 2012), Eq. 5 is an algorithm of 3rd-order accuracy.This algorithm is referred to as CB3D scheme in this study.Since the 2nd-order antidissipative terms in Eq. 3 have been eliminated, the stability of the new scheme should be strengthened.
Next the CB2D and CB3D algorithms will be used to determine the flow-field of an isolated rotor (Li and Chen 2012), and the comparisons of iteration convergence speed and wake geometry are shown in Figs.3a and 3b.
From Fig. 3 we can see that the convergence speed of CB3D is faster than CB2D, and the wake geometry predicted by the 2 algorithms have got a little difference, which could be attributed to the different error terms in the modified equations of the 2 algorithms.

FUSELAGE MODEL
The interaction effects of fuselage on rotor flow-field are all added by the panel method.According to the low-speed aerodynamics (Katz 2001), the velocity potential ϕ at a point in space induced by the panels is determined by the distribution of doublets μ and sources σ, as follows:  According to the research made by Katz (2001), the control point is set at the centroid of every element.The flow-field velocity on every control point needs to satisfy the no-penetration boundary conditions.For solving the doublet and source strength of every panel, the source panel strength can be written as Eq.8: where: V → ∞ is the free stream velocity; V → freewake-induce is the induced velocity at control point induced by the far-wake vortex; V → fuselage is the motion velocity of fuselage at control point.Based on the relationship between velocity potential and velocity, combined with the no-penetration boundary condition on solid surface, the velocity at the i th control point should satisfy the equation: where: B i,k and C i,k are the induced velocity at the i th panel element induced by the unit strength source or doublet distributed at the k th panel element.
We can get N equations of N panel elements like Eq. 9.The source strengths σ k can be achieved from Eq. 8, so the doublet strengths can be obtained by solving the N-dimensional linear equations.Then the induced velocity and the induced velocity potential at any point in the flow-field can be calculated.
Considering that the flow on the fuselage is unsteady under the influence of rotor interaction, the pressure and load on the fuselage is obtained through solving the unsteady Bernoulli Equation in:

"MATERIAL LINE" RECTIFICATION METHOD
The description of flow-field by the free-wake mathematical model is based on the Lagrange method.When updating the free-wake geometry on discrete temporal and spatial step, part of the wake vortices may not physically penetrate the fuselage surface as shown in Fig. 5.As shown in Fig. 6, in the rectification, first, at the last discrete time step when the wake collocation point is still outside the solid surface, a panel control point near the collocation point should be selected as an auxiliary point.The auxiliary point and the collocation point are connected with a straight line, then extend the line segment to another space auxiliary point, which is demanded to not penetrate the solid surface after one time step and get as close to the collocation point as possible.
Considering the extended line segment as a fluid material line, since one discrete time step is very short, the material line segment can be thought as a straight line segment after one time step (Fig. 7).So, if the 2 auxiliary points do not penetrate the solid surface, all the points on the material line segment do not penetrate the solid surface, and the rectified location of collocation point could be predicted by interpolation between 2 auxiliary points on the line segment.
In view of this problem, a "material line" rectification method is proposed in this section.The method is based on the fact that the fluid micelles at the panel control points will

Surface auxiliary point
Auxiliary point

Collocation point
Fuselage surface If we write the collocation point coordinates as (x p0 , y p0 , z p0 ), the fuselage surface auxiliary point coordinates as (x 10 , y 10 , z 10 ), and the auxiliary point coordinates on the extended line as (x 20 , y 20 , z 20 ), after one time step, the coordinates of each of the abovementioned 3 points can be written as (x pt , y pt , z pt ), (x 1t , y 1t , z 1t ), and (x 2t , y 2t , z 2t ), and the collocation point coordinates (x pt , y pt , z pt ) after rectification can be predicted by the following functions: Next, an accuracy analysis for the "material line" rectification method will be given.Let us choose the rectification formula as a representative sample for the analysis.If we consider k = (z p0 -z 10 )/(z 20 -z 10 ), then we have: At the last time step, the velocity on z direction can be written as V z , then we have: So the interpolation for collocation point location is an interpolation for its velocity in essence.Next, this section will give an accuracy analysis for point velocity interpolation.
If we set the surface auxiliary point as original point (subscript 1) and set l axis along the "material line", the coordinates of space auxiliary point (subscript 2) in this coordinate system can be written as Δl.According to the proportional relation in Eq. 13, expanding V z2 and V z1 at location of the collocation point with a Taylor series, we can get: Adding Eq. 15 and Eq. 16, we obtain: (12) According to the former auxiliary points method chosen, the magnitude orders of the material line length Δl and (V 20 gΔt), that is (V 20 g(Δψ/Ω), are close.From Eq. 17, we can see that the interpolation accuracy of "material line" rectification method for velocity is of 2nd-order.So the interpolation accuracy of the method for location is O(Δl 2 )gΔt, whose 3rd-order is in temporal discrete scale.

SIMULATION AND ANALYSIS FOR ROTOR-FUSELAGE AERODYNAMIC INTERACTION
Leishman team conducted a series of experiments (Bi and Leishman 1990; Leishman and Bi 1990) about rotor-fuselage interaction phenomenon; the outlook of fuselage model they used is similar to the actual helicopter fuselage as shown in the left part of Fig. 4. In this section, the rotor/fuselage aerodynamic interaction of Maryland model in forward flight is simulated and compared with measured data.Finally, the influences of advanced ratio on unsteady aerodynamics interaction are analyzed.The simulated typical free-wake structure in interaction is shown in Fig. 8.
As can be seen in Fig. 8, in forward flight situation, the free-wake vortex approaches the rear of fuselage.This indicates that the vortex-panel interactions become more intensive as the free-wake moves towards the backward of fuselage.
The distribution of pressure coefficient (C p ) along the fuselage left centerline at various advanced ratios is shown in Fig. 9.The abscissa in Fig. 9 means the dimensionless lengthwise position on the left centerline -X f is the lengthwise position which is 0 at the airframe nose and L is the total length of the fuselage.As can be seen, the simulated results compare well with the experimental data, which further validates the accuracy of aerodynamic interaction model developed in this paper.From Fig. 9 we can also point out that, as the free stream velocity increases, the pressure coefficient on the fuselage nose increases.The negative peak of pressure coefficient emerges at the rear of the fuselage, which reflects the location of vortex/ panel interaction and moves backward on the fuselage as the advanced ratio increases.

CONCLUSIONS
In this paper, an unsteady panel/time-stepping free-wake method is established to predict the helicopter rotor flow-field in aerodynamic interaction conditions.The panel method is applied to simulate the influence of helicopter fuselage on flow-field.A "material line" rectification method is proposed to solve the problem that part of the wake vortices may penetrate the fuselage in numerical simulation.Considering the problem that the stability of traditional time-stepping free-wake algorithms is deficient for the simulation of unsteady aerodynamics in rotor interaction, an explicit CB3D algorithm with 3rd-order accuracy is proposed.According to this study, a few general conclusions can be drawn:

•
The coupled panel/time-stepping free-wake method is proved to be accurate and efficient for predicting the rotor flow-field of rotor-fuselage interaction.The simulated results are in good agreement with available experimental data.

•
The "material line" rectification method is flexible and has a 3rd-order precision.It is applicable to solve the nonphysical motion problem encountered by the interaction between rotor wake geometry and solid surface.
• Full-span free-wake model improves the simulation accuracy of rotor/fuselage interaction phenomenon.As the free stream velocity increases, the pressure coefficient on the head of fuselage increases, and the negative peak of pressure coefficient at the rear of fuselage moves backward.

Figure 1 .
Figure 1. Figure of blade lifting and tip vortices.

Figure 2 .
Figure 2. Full-span free-wake geometry trailing from one blade.

Figure 3 .
Figure 3. Comparisons of iterated convergence speed and wake geometry predicted by CB2D and CB3D.(a) Convergence speed of CB2D and CB3D; (b) Wake geometry predicted by CB2D and CB3D.
Sf represents the whole boundary surface; r is the distance from the doublet or source to the measured point; n → is the local outward normal vector of the boundary surface.The sources and doublets are uniformly distributed on every surface panel element.The discrete fuselage model is shown in Fig. 4. where: P, V → and ϕ are, respectively, the pressure, velocity and potential of the calculating point; ρ is the air density,The Theoretical Research for the Rotor/Fuselage Unsteady Aerodynamic Interaction Problem in kg•m -3 ; R is the blade radius, in m; P ∞ is far-field pressure; C p is defined as the non-dimensional pressure coefficient in this paper.The studies ofCrouse et al. (1992) andQuackenbush et al. (1999) demonstrated that, in fuselage surface pressure determination, the unsteady term cannot be ignored.So, the unsteady term ∂ϕ/∂t is the aporia in unsteady flow-field pressure prediction.If r is the distance from the potential panel to the predicted point, that is r = √ (x -x΄) 2 + (y -y΄) 2 + (z -z΄) 2 , considering, ∂ϕ/∂x = (∂ϕ/∂r) g (∂r/∂x), we have, ∂ϕ/∂x΄ = (∂ϕ/∂r) g (∂r/∂x΄) = 1g (∂ϕ/∂r) g (∂r/∂x), and the calculation formula of ∂ϕ/∂t can be written as: is the motion velocity of the potential panel.So, in the case of a stationary fuselage, we can get ∂ϕ/∂t = -1g

Figure 6 .
Figure 6.Normal located wake collocation point and chosen auxiliary point at the last time step.

Figure 7 .
Figure 7. Location of the rotor wake collocation point after rectification.

Figure 8 .
Figure 8. Space diagram and side elevation of rotor wake geometry at an advance ratio of 0.2.

Figure 9 .
Figure 9.Comparison of predicted and experimental time average pressure on the left centerline of the fuselage surface at advance ratios of 0.075, 0.10 and 0.15 (C T /σ = 0.085, α s = −6°, where C T is the rotor thrust coefficient and α s is the rotor shaft tilt angle).
The Theoretical Research for the Rotor/Fuselage Unsteady Aerodynamic Interaction Problem