Prediction of Vertical Tail Buffet Using CFD/CSD Coupling Method

.e vertical tail buffet induced by the vortex breakdown flow is numerically investigated. .e unsteady flow is calculated by solving the RANS equations..e structural dynamic equations are decoupled in the modal coordinates. .e radial basis functions (RBFs) are employed to generate the deformation mesh. .e buffet response of the flexible tail is predicted by coupling the three sets of equations. .e results show that the presence of asymmetry flow on the inner and outer surface of the tail forced the structural deflection offsetting the outboard..e frequency of the 2nd bendingmode of the tail structure meets the peak frequency of the pressure fluctuation upon the tail surface, and the resonance phenomenon was observed. .erefore, the 2nd bending responses govern the flow field surrounding the vertical tail. Finally, the displacement of the vertical tail is small, while the acceleration with a large quantitation forces the vertical tail undergoing severe addition inertial loads.


Introduction
LARGE swept delta wing combined with vertical tail configuration is widely used by modern flight vehicles. is kind of aerodynamic configuration is designed to reduce drag and improve the manipulated performance at a large attacking angle. However, the vortices emanating from the leading edge of the delta wing might breakdown before reaching the vertical tail at a certain attacking angle. In such cases, the vertical tails suffer from the wideband frequency excitation of the highly turbulent flow, which would produce a severe buffet on the tails. Wentz [1] conducted a wind tunnel experiment to investigate F/A-18 tail buffet problem, and the results revealed that the buffet is a phenomenon of resonance induced by the excitation of vortex breakdown flow, while the frequency bandwidth of the excitation covers the natural frequency of tails. e characteristics of buffeting are different from the traditional flutter or limit cycle oscillation (LCO). e buffet displacement with slight amplitude produces large acceleration responses that make the tail suffer from severe inertia force loads and fatigue failure.
Numerical studies for the aeroelastic responses in the vortical flow conducted before are mainly focused on the aerodynamic load, flutter, or LCO responses for lightweight materials. However, studies on the time/frequency characteristic of buffet responses are relatively inadequate. Findlay [2] simulated leading-edge vortex breakdown flow nearby vertical tails by solving thin-layer Navier-Stokes equations. It was found that the computed pressure distributions on the surface of the tails agree well with the experimental measurements. Kang [3] calculated the vertical tail buffeting loads at large attacking angle and obtained the buffeting characteristics. Kandil [4][5][6][7][8] investigated the bending-torsion interaction effect on the tail's buffet, as well as the effect of freestream Mach number on the tail's buffet. e results revealed that the deflections and loads of the coupled bending-torsion response case are substantially lower than those of the uncoupled bending-torsion response case, and the buffet loads in the transonic flow are lower than those in the subsonic flow. Guillaume [9] simulated the F/A-18 tail buffet in the transonic flow and obtained the displacement responses of the tails. Ahmed Elmekawy [10] conducted numerical simulations for F/A-18 tail buffet with four turbulence models and compared the computational results with the experiment data. It was found that the buffet problem could be simulated as a two-way fluid-structure interaction, and by using the nonlinear eddy viscosity model, since these give better results when compared to the other considered models. Attar [11] conducted numerical simulations for the full-span delta wing buffet in the vortex breakdown flow. Attar believes the lift enhancement phenomena are due to the reorganization of the flow and reformation of a leading-edge vortex structure. e present study conducts a numerical simulation to investigate the tail buffet characteristics. Roe [12] proposed a finite different scheme for obtaining numerical solutions to Riemann problems. Results show that this fast scheme is very accurate and might be used to reduce the number of iterations. is method was employed in this paper to discretize the Euler flux terms of N-S equations. Yoon [13] developed a vectorizable and unconditionally stable scheme to solving the Euler and N-S equations. Results show that this method is efficient and robust. is method was used for the time integration of the N-S equations in this paper. Rodden [14] demonstrated a fast interpolation scheme to fitting the payload data of an infinite plane. is method was used here to perform the transformation of structural deformation from the CSD grid points to the CFD grid points. Jakobsson [15] proposed an interpolation scheme to propagate the deformations from the boundaries to the interior of the CFD mesh. Both structured and unstructured meshes can be treated by this method with low computational costs. is method was employed here to generate the deformation mesh.
e unsteady Reynolds-averaged Navier-Stokes (RANS) equations with k-w SST turbulent model were used to model the vortical flow. e unsteady RANS solver was validated using the data from the Hummel's 76-degree delta wing [16]. A structural model of the vertical tail was established by the finite element analysis. e flexible tail was positioned after the leading-edge vortex breakdown point to investigate the buffet characteristics. e deflection and acceleration responses are calculated by solving RANS equations coupled with structural dynamic equations.

Numerical Method
e numerical method for simulating the tail buffet was used to solve the aeroelastic equations, which consisted of the Navier-Stokes equations and structural dynamic equations. ese equations coupled and advanced in time. Figure 1 shows the time stepping method.
Details of the methods are described as follows.

RANS Equations.
e unsteady compressible Navier-Stokes (N-S) equations in the conservative form can be written as where t V represents the computational domain and V represents the computational boundary. e vector of conservative variables Q � [ρ, ρu, ρv, ρw, ρe]/J, J represents the Jacobian determinant, ρ is density, u, v, w represents the velocity magnitude in three coordinate directions, e represents the total energy per unit mass, and t is time. All the variables mentioned above are dimensionless. A semidiscrete finite-volume formulation is used to discretize the N-S equations. e viscous flux terms Γ i are discretized by the central difference scheme, and the inviscid flux vector F i is discretized by FDS method of Roe [12], where the subscript variable i � 1, 2, 3, represents three coordinate directions. e LU-SGS implicit scheme [13] is used for time integration. e turbulence model was the k-ω SST model, which solves a two-equation turbulence model for the kinetic turbulent energy k and specific dissipation rate ω: where G k represents the generation of turbulence kinetic energy due to mean velocity gradient. G ω represents the generation of ω. Y k and Y ω represent the dissipation of k and ω due to turbulence. D ω represents the cross-diffusion term. S k and S ω are user-defined source terms. σ k and σ ω are the Prandtl numbers for k and ω. μ t is the turbulence viscosity: where S is the strain rate magnitude. e modal constant a 1 � 0.31, Re t � ρ·k/(μ·ω), and the blending function F is given by where y is the distance to the closest surface. Note that the constant coefficient α in the k − ω SST model is defined as α � α 1 ·F + α 2 ·(1 − F); the model constant α 1 and α 2 are 5/9 and 0.44, respectively.

Structural Dynamic Equations.
e structural dynamic equations in the form of finite element analysis can be written as [11] [ where N is the total number of node points of the structural model; u i and f i are vectors with three components in x, y, z directions: In a modal approach, the modal decomposition of the structure motion can be expressed as the following: where ω j is the natural frequency of jth mode, and the mode shape matrix [ϕ] � [ϕ 1 , ϕ 2 , . . ., ϕ 3N ]. Equation (8) can be solved by using a finite element solver (e.g., MSC.NASTRAN) to obtain [ϕ] and [Λ].
us, the displacement of the structure can be calculated as follows: e energy of high order modes is relatively small. In this paper, only the first four modes are used in equation (9). us, the vector [a] � [a 1 , a 2 , a 3 , a 4 ] T . Substituting equation (9) to equation (5) yields Multiplying equation (10) (c 1 , c 2 , . . . , c 3N ), and (k 1 , k 2 , . . . , k 3N ). (11) can be transformed in a state form as follows: where I represents the unit matrix. Equation (12) can be solved by the four-stage Runge-Kutta method. e infinite plane spline (IPS) method [14] was adopted as CFD-CSD interface data transformation method.

Mesh Deformation
Formulations. e radial basis function (RBF) [15] was used to transform the CFD mesh.
e interpolation function is defined as where F (r) denotes the interpolation function; N refers to the number of RBFs. In this paper, N is equal to the total number of the grid points on the surface. φ(·) denotes the RBF, w i refers to the weighting factor of the function φ(·), r refers to the position vector of grid point, r i refers to the central point position vector of the RBF, and d refers to the influencing radius of the interpolation function. e Wendland's C 2 function was employed as φ(·): e weighting factor w i can be obtained by solving the following equations: where ΔS i is a vector of the displacement magnitude of each grid point on the surface, the subscript i � 1, 2, 3 denotes the x, y, z direction, respectively, the weighting factor vector W � [w 1 , w 2 ,. . ., w N ] T , and the coefficient matrix Φ,

Results and Discussion
e computational configuration is a half-delta wing with a vertical tail, as shown in Figure 2. e x, y, z coordinate, the wing chord c, and spanwise length of the wing L are described in Figure 2. More details of the model geometry and the structural parameters are described in Section 3.1. e wind tunnel experiments of the rigid delta wing are employed to verify the computation code.
e buffet responses of the flexible tail as well as the instantaneous flowfields are demonstrated in Section 3.2.

Model, Grid, and Verification Test.
e geometry of the delta wing is shown in Figure 3. e aspect ratio of the wing A � L 2 /S � 1, where S represents the wing area. e leadingedge swept angle of the wing is 75.96°. e root chord length of the wing is 1 m. e maximum thickness of the wing is 0.021 m. e vertical tail is located at the 1.06 chord length behind the apex of the wing, with a leading-edge sweep of 62.5°. e structure model of the tail was a single aluminum spar hidden inside a balsa-wood covering, as shown in Figure 4. e aluminum spar has a constant thickness of 0.0017 m. e root chord of the aluminum spar is 0.04 m, and the tip chord is 0.012 m. e aluminum spar was constructed from 6061-T6 alloy with a density of ρ � 2693 kg/m 3 , a moduli of elasticity of E � 6.896 × 10 10 N/m, and a moduli of rigidity of G � 2.5925 × 10 10 N/m. e root chord of the balsa-wood covering was 0.232 m, the tip chord was 0.054 m, and the spanwise length was 0.204 m. e thickness of the balsa-wood gradually decreases from 0.02 m at the tail root to 0.01 m at the tail midspan. From the midspan to the tail tip, the thickness was constant. e tail cross-section has a semidiamond shape with a bevel angle of 20°. e balsa-wood density was ρ � 179.7 kg/m 3 , the moduli of elasticity was E � 6.896 × 10 10 N/m, and the moduli of rigidity was G � 2.5925 × 10 10 N/m. e fundamental frequencies of the first and second bending mode of the tail are 33.941 Hz and 152 Hz, respectively, while the torsional frequency was 241.09 Hz, as shown in Figure 5. Figure 6 shows the CFD-CSD interface interpolation effect for the mode shapes of the vertical tail.
An O-H grid of 87 × 107 × 51 grids that points in the wrap-around, streamwise, and normal directions, respectively, was used for the solution. Figure 7 compares the original mesh nearby the delta wing/tail surface with the mesh deformation. e flow conditions were Re c � 1.25 × 10 6 , the attacking angle was α � 40°, and the freestream Mach number was Ma � 0.3. e unsteady time step was chosen as 10 −5 s to capture the high frequency vortex breakdown flow. e experiment [16] for the vortex flow on the delta wing was employed to validate the computational capability of the CFD procedure. e freestream Mach number was 0.3, Re c � 2 × 10 6 , and the attacking angle was α � 20.5°. e wing surface pressure measurement positions were x/c � 0.3, 0.5, 0.7, and 0.9. Figure 8 compares the computed surface pressure with the experimental measurements. It was shown that the computed solutions agree well with the experimental measurements. Figure 9 compares the calculated surface flow patterns with the experiments. It was shown that the main and the secondary separation lines were both captured by the numerical procedure.
e results above show that the computational method presented in the present study is reliable.

Vertical Tail Buffet Results.
e computational instantaneous streamlines at 0.2 s and 0.25 s are shown in Figure 10, and the total pressure slices and tip deflection of the vertical tail at the position of x/c � 1.3 were also extracted. It is found that the leading-edge vortices proceed to the outward side of the vertical tail, resulting in the tail deflection towards the outward side. e buffet load coefficients ΔC P nearby the tip of the tail (at the position of 90% spanwise length of the tail, 50% local root length of the tail) and the root of the tail (at the position of 30% spanwise length of the tail, 50% local root length of the tail) were both extracted in Figure 11. It was found that the mean value and the fluctuation magnitude of the buffet loads at the tip of the tail were both larger than those at the root of the tail, and the frequency bandwidth of the fluctuation loads was from 100 Hz to 200 Hz. In addition, the peak frequency of the buffet loads at both the tip and root of the tail appears at 153 Hz, which is almost the same as the fundamental frequency of the second bending mode of the tail, demonstrating the tail structure resonance under the fluctuation load excitation of the vortex breakdown flow. Figure 12 compares the displacement responses with the acceleration responses on the midpoint of the leadingedge of the tail. It was observed that the amplitude of the displacement responses was small (magnitude: 10 −4 m). However, the amplitude of the acceleration responses was 5∼6 orders of magnitude greater than the displacement responses (magnitude: 10 3 m/s 2 ). Furthermore, the acceleration amplitude clearly diverges with time. e PSD curves in Figure 12 show that the peak frequency of both displacement and acceleration responses was 151.23 Hz, which meets the fundamental frequency of the second bending mode. e results above conclusively show that the buffet deflections with slight amplitude produce acceleration responses with a large amplitude, causing severe additional internal loads on the structure, and the resonance between the vortex breakdown flow and the second bending mode dominates the characteristics of the acceleration responses.

Conclusion
A computation method for solving RANS equations coupled with structural dynamic equations was presented in the present study to simulate the vertical tail buffet responses in the vortex breakdown flow. A validation computation of the vortex flow on Hummel's 76°delta wing was conducted to verify the numerical method. A delta wing/vertical tail configuration was modeled and simulated at a large angle of attack to obtain the buffet response of the tail. It was observed that the vortex breakdown flow on the inner and outer surface of the tail was asymmetry, and the vortices were located at the outer side of the tail, which makes the tail deflect to the outer side. e peak frequency of the buffet load was almost the same as the fundamental frequency of the second bending mode of the tail. e    Mathematical Problems in Engineering response amplitude of the second modal displacement and acceleration diverge with time. Meanwhile, the amplitudes of the first and third modal responses were irregular oscillations with time. e frequency of the modal displacement and acceleration responses was the fundamental frequency of each structural mode. In addition, the deflection amplitude of the tail is small (magnitude: 10 −4 m). However, the acceleration amplitude was much greater (magnitude: 10 3 m/s 2 ). e results above conclusively show that the resonance between the vortex breakdown flow and the second bending mode dominates the characteristics of the acceleration responses, and the buffet deflections with slight amplitude produce acceleration responses with large amplitude, causing severe additional internal loads on the structure.
Data Availability e data that support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.