2.1 Mathematical formulation

In the first step, the potential solution of the linear vortex lattice is calculated, by solving the linear system that imposes the non-penetration boundary condition on the vortex lattice panel as follows:

(1) \begin{align}{\hat{\textbf n}} \cdot {{\textbf{u}}_\phi } = {\hat{\textbf n}} \cdot \left( {{{\textbf{u}}_b} - {{\textbf{U}}_\infty } - {{\textbf{u}}_\psi }} \right),\end{align}

where ${{\textbf{u}}_\phi }$ is the potential velocity, ${{\textbf{u}}_b}$ is the body velocity, ${{\textbf{U}}_\infty }$ is the free-stream velocity ${{\textbf{u}}_\psi }$ is the rotational perturbation velocity. In the vortex lattice method, a surface is discretised as a sheet of vortex wing with intensity ${{\rm{\Gamma }}_{{i_v}}}$ , equivalent to a piecewise-uniform surface doublet distribution with the same intensity ${\mu _{{i_v}}}$ . The boundary condition is written for each panel collocation point and the problem can be rewritten in a linear system form as

(2) \begin{align}{\textbf{A}}{\rm{\Gamma }} = {\rm{RHS}},\end{align}

where each term ${a_{ij}}$ of the matrix ${\textbf{A}}$ is defined as the velocity component induced by unit strength ${{\rm{\Gamma }}_i}$ , normal to the surface at collocation point $i$ . Then, a fixed-point iterative problem is solved by imposing a convergence between the lift coefficient derived from the aerodynamic 2D look-up table and the one calculated from the stripe intensity. The total velocity $\boldsymbol{u}_{{v_{nl}}}^i$ acting on the control point (centre) at time $t$ and iteration $i$ of each stripe can be calculated as follows:

(3) \begin{align}\boldsymbol{u}_{{v_{nl}}}^i = {\boldsymbol{U}_\infty } - {\boldsymbol{u}_{{v_{nl}},b}} + {\boldsymbol{u}_{tot\backslash {v_{nl}}}} + \boldsymbol{u}_{{v_{nl}}}^{i - 1}\end{align}

where ${\boldsymbol{U}_\infty }$ is the free stream velocity, ${\boldsymbol{u}_{{v_{nl}},b}}$ is the body velocity, ${\boldsymbol{u}_{tot\backslash {v_{nl}}}}$ is the velocity induced by all elements (potential wake, vortical particles, and lifting body except the non-linear vortex lattice elements), $\boldsymbol{u}_{{v_{nl}}}^{i - 1}$ is the induced velocity of the stripe on the other stripes which is calculated as the velocity induced by a surface doublet whose vertices are the leading and trailing edge nodes and the intensity is obtained as follows:

(4) \begin{align}{\rm{\Gamma }}_{{v_{nl}}}^{i - 1} = d{\rm{\Gamma }}_{{v_{nl}}}^1 + \mathop {\mathop \sum \limits_{k = 2} }\limits^N d{\rm{\Gamma }}_{{v_{nl}}}^k - d{\rm{\Gamma }}_{{v_{nl}}}^{k - 1}\end{align}

where $d{\rm{\Gamma }}_{{v_{nl}}}^k$ is the intensity panel $k$ -th of the stripe calculated at iteration $i - 1$ , and $N$ is the number of chordwise elements. Then, the velocity is projected on the plane formed by the stripe normal $\hat{\boldsymbol{n}}$ and chord $\hat{\boldsymbol{t}}$ :

(5) \begin{align}\boldsymbol{u}_{2D}^i = \left( {\boldsymbol{u}_{{v_{nl}}}^i \cdot \hat {\boldsymbol{n}}} \right)\hat{\boldsymbol{n}} + \left( {\boldsymbol{u}_{{v_{nl}}}^i \cdot \hat{\boldsymbol{t}}} \right)\hat{\boldsymbol{t}}\end{align}

Angle-of-attack $\alpha $ is calculated as:

(6) \begin{align}{\alpha ^i} = {\tilde \alpha ^i} - \alpha _{2D}^i\end{align}

where ${\tilde \alpha ^i}$ is calculated as:

(7) \begin{align}{\tilde \alpha ^i} = {\rm{atan}}\!\left( {\frac{{\boldsymbol{u}_{2D}^i \cdot \hat{\boldsymbol{n}}}}{{\boldsymbol{u}_{2D}^i \cdot \hat{\boldsymbol{t}}}}} \right)\end{align}

This angle must be corrected by taking out the velocity induced by the stripe itself. According to Ref. [Reference Piszkin and Levinsky20], this correction is equal to:

(8) \begin{align}\alpha _{2D}^i = \frac{{2 \cdot {\rm{\Gamma }}_l^{i - 1}}}{{{C_{L/\alpha }}{\rm{\;}}c\left| {\boldsymbol{u}_{2D}^i} \right|}}\end{align}

where ${C_{L/\alpha }}$ is the slope of the lift curve in the linear range obtained from the tabulated data, and $c$ is the stripe chord. For each iteration, the vortex lattice linear system is solved by modifying only the right-hand side and taking as a result the updated circulations. To numerically stabilise the method two kinds of relaxation are available:

• • Constant relaxation, where the updated right-hand side of the system is

  (9) \begin{align}{\rm{RH}}{{\rm{S}}_i} = {\rm{RH}}{{\rm{S}}_{i - 1}} + \alpha \cdot {\textbf{r}}\end{align}
  with $\alpha $ being the constant relaxation factor and ${\textbf{r}}$ is the residual vector defined as the difference between the lift coefficient solved from the linear system $C_L^{{\rm{inv}}}$ and the one obtained for the aerodynamic 2D look-up table at iteration $i$ $C_L^{{\rm{visc}}}$ :
  (10) \begin{align}{r_s} = \left( {C_L^{{\rm{inv}}} - C_L^{{\rm{visc}}}} \right) \cdot {\rm{sign}}\!\left( {C_L^{{\rm{inv}}} - C_L^{{\rm{visc}}}} \right)\end{align}

• • Aitken acceleration [Reference Küttler and Wall21], where the updated right-hand side of the system is

  (11) \begin{align}{\alpha _i} = - {\alpha _{i - 1}}\frac{{{{\textbf{r}}_{i - 1}} \cdot {\rm{\Delta }}{\textbf{r}_i}}}{{{\rm{\Delta }}{{\textbf{r}}_i} \cdot {\rm{\Delta }}{{\textbf{r}}_i}}}\end{align}
  with ${\alpha _i}$ being the Aitken relaxation factor for iteration $i$ , ${\textbf{r}}$ is the residual vector and ${\rm{\Delta }}{{\textbf{r}}_i}$ is the difference between the residuals at iteration $i$ and $i - 1$ . Then the ${\rm{RHS}}$ is updated as described in Equation (9).

For non-linear vortex lattice elements, aerodynamic loads are calculated using the unsteady formulation of the Kutta-Joukowsky theorem, as for the linear case plus a viscous correction on the drag term:

(12) \begin{align}\boldsymbol{d}{\boldsymbol{F}_{i \in s}} = \boldsymbol{dF}_{i \in s}^{sys} + \frac{1}{2}{{ho} _\infty }({\boldsymbol{U}^2}\!\left( {\boldsymbol{r},t} \right)\!d{A_{i \in s}}c_d^s\!\left( {{\rm{cos}}\!\left( \alpha \right)\!\hat{\boldsymbol{t}} + {\rm{sin}}\!\left( \alpha \right)\!\hat{\boldsymbol{n}}} \right)\end{align}

where $i \in s$ stands for a panel $i$ that belongs to the stripe $s$ and $\boldsymbol{dF}_{i \in s}^{sys}$ is the force calculated from the system.

2.1.1 Coupled governing equations

The potential contribution of the velocity field results from the free-stream velocity’s superposition and the induced velocity from the doublet and source singularities, modeling the bound vorticity on the surface of the investigated bodies. Given the free stream velocity, the rotational component of the velocity field, the state of the wake, and the motion of the solid components of the model, the intensity of the surface doublets are the actual unknowns of the potential problem obtained from the surface panel discretisation, lifting line discretisation and the subdividing vortex lattice in the linear $\left( {{v_l}} \right)$ and non-linear case $\left( {{v_{nl}}} \right)$ :

(13) \begin{align}\left\{\begin{array}{l} \left[\begin{array}{c@{\quad}c@{\quad}c@{\quad}c} \textbf{A}_{s,s} & \textbf{A}_{s,v} & \textbf{A}_{s,v_{nl}} & \textbf{A}_{s,l} \\[5pt] \textbf{C}_{v,s} & \textbf{C}_{v,v} & \textbf{C}_{v,v_{nl}} & \textbf{C}_{v,l} \\[5pt] \textbf{C}_{v_{nl},s} & \textbf{C}_{v_{nl}, v} & \textbf{C}_{v_{nl}, v_{nl}} & \textbf{C}_{v_{nl}, l} \end{array} \right] \left\{ \begin{array}{c} \boldsymbol{\Gamma}_s \\[5pt] \boldsymbol{\Gamma}_{vl} \\[5pt] \boldsymbol{\Gamma}_{nl} \end{array}\right\} = -\left[\begin{array}{c@{\quad}c@{\quad}c} \textbf{B}_{s, s} & \textbf{0}_{s, v_{l}} & \textbf{0}_{s, v_{nl}} \\[5pt] \textbf{D}_{v, s} & -\textbf{I}_{v, v} & \textbf{0}_{v, v_{nl}} \\[5pt] \textbf{D}_{v_{nl}, s} & \textbf{0}_{v_{nl}, v} & -\textbf{I}_{v_{nl}, v_{nl}}\\[5pt] \end{array}\right] \left\{ \begin{array}{c} \sigma_s \\[5pt] \sigma_{vl} \\[5pt] \sigma_{nl} \end{array}\right\} -\left[\begin{array}{c} \textbf{A}_{sw} \\[5pt] \textbf{C}_{v_lw} \\[5pt] \textbf{C}_{v_{nl}w} \end{array}\right] \boldsymbol{\Gamma}_w \\[5pt] \boldsymbol{\Gamma}_l = \boldsymbol{f}_l ( \boldsymbol{\Gamma}, \boldsymbol{\sigma}, \boldsymbol{\Gamma}_w )\end{array}\right.\end{align}

where the vector ${{\boldsymbol{\Gamma }}_w}$ collects the intensity of the wake doublets, ${\boldsymbol{\Gamma }} = \left( {{{\boldsymbol{\Gamma }}_s},{{\boldsymbol{\Gamma }}_{{v_l}}},{{\boldsymbol{\Gamma }}_{{v_{nl}}}},{{\boldsymbol{\Gamma }}_l}} \right)$ collects the intensity of the body doublets on the surface panel, the vortex lattice and the lifting line elements, the vector $\boldsymbol\sigma = \left( {{\boldsymbol\sigma _s},{\boldsymbol\sigma _{{v_l}}},{\boldsymbol\sigma _{{v_{nl}}}}} \right)$ contains the unperturbed relative normal velocity at the collocation points of the surface panels and the vortex lattice elements, ${\textbf{0}_{sv}}$ is a null matrix ${N_s} \times {N_v}$ , ${{\textbf{I}}_{vv}}$ is the identity matrix ${N_v} \times {N_v}$ . The non-linear vortex lattice is considered implicit in the sense that its solution is obtained by iteratively solving the linear system that contains all the elements except the lifting line, whereas the solution of the lifting line problem is decoupled from the other element type.

2.2 Validation of the method over fixed-wing test cases

In order to validate the implementation of the numerical element, this approach was tested on fixed-wing benchmark cases available from literature related to aerodynamic and aeroelastic problems.

Firstly, a simple finite wing benchmark case was considered consisting of a rectangular wing with an aspect ratio equal to 12, with a constant chord length equal to 1 m and a constant NACA 4415 aerofoil section. Figure 1 shows the comparison between the wing lift coefficients computed by DUST and high-fidelity CFD simulations results [Reference Petrilli, Paul, Gopalarathnam and Frink22]. In particular, DUST simulation results are related to the use of both the classical LL and VL elements and the non-linear vortex lattice method (NL-VL).

Figure 1. Comparison of the NACA 4415 wing lift coefficient ( ${C_L}$ ) as a function of angle-of-attack: high-fidelity CFD simulation results from Ref. [Reference Petrilli, Paul, Gopalarathnam and Frink22] (CFD), (DUST VL) LL and NL-VL simulations results. Mach number = 0.2.

The effectiveness of the non-linear vortex-lattice method is evident in the linear range as well as in the post-stall region, as shown by the remarkable agreement with high-fidelity CFD simulations. In particular, the NL-VL results are quite similar to the ones obtained using the DUST LL approach considering viscous contributions for airloads evaluation. This result confirms the correctness of the implementation based on the integration of the viscous aerodynamic 2D data and shows the capabilities of the method to capture stalled and post-stalled conditions with respect to the purely potential VL model. In addition, thanks to the possibility of the NL-VL elements to evaluate the chordwise pressure distribution, Fig. 2 shows the ${\rm{\Delta }}{C_p}$ distribution evaluated between upper and lower surface of the NACA4415 aerofoil by DUST simulations compared to high-fidelity CFD simulations results from [Reference Petrilli, Paul, Gopalarathnam and Frink22] for a nearly stalled and post-stalled condition.

Figure 2. Comparison of ${\rm{\Delta }}{C_p}$ distribution evaluated between upper and lower surface of the NACA4415 aerofoil for high-fidelity CFD simulations from Ref. [Reference Petrilli, Paul, Gopalarathnam and Frink22] (CFD) and DUST NL-VL simulations. Mach number = 0.2.

This comparison shows an accuracy similar to CFD of DUST NL-VL method in terms of chordwise pressure difference over the aerofoil for the nearly stalled condition characterised by a low degree of flow separation. On the other hand, due to the limitation of the VL approach for highly separated flow conditions, pressure distribution comparison at very high angle-of-attack shows higher discrepancy with respect to CFD.

Secondly, the Lovell [Reference Lovell23] wing test case was considered for a further validation of the DUST methodology. This test case consists in a 30.5 $^{\circ }$ leading edge swept wing (see Fig. 3) used for wind tunnel investigation of the effect of flap deflections. The geometry is simple, with a constant aerofoil geometry along the span and experimental results are available for the isolated wing with different flap deflections. DUST numerical models were built using 50 spanwise elements for both LL and NL-VL while 5 elements in chord were used for the NL-VL model. DUST simulations results are compared to VLM coupled with RANS solutions [Reference Parenteau, Laurendeau and Carrier16] and experimental data for lift, drag and pitching moment coefficients, see Figs. 4, 5, and 6. Results comparison showed that DUST models reproduce the behaviour of the lift and polar curves obtained by Parenteau et al. [Reference Parenteau, Laurendeau and Carrier16] with a similar numerical model, thus further validating the implementation of NL-VL model. In particular, the comparison of the ${C_L} - {C_M}$ curve shows a better agreement of the solution obtained with DUST with respect to experimental data.

Figure 3. DUST simulation results of the Lovell wing test case [Reference Lovell23]: wake visualisation and sectional load distribution.

Figure 4. Lovell clean configuration: CL. Comparison between experimental results, Parenteau et al. VLM (2D RANS), DUST LL, and DUST NLVL.

Figure 5. Lovell clean configuration: CD. Comparison between experimental results, Parenteau et al. VLM (2D RANS), DUST LL, and DUST NLVL.

Figure 6. Lovell clean configuration: CM. Comparison between experimental results, Parenteau et al. VLM (2D RANS), DUST LL, and DUST NLVL.

Thirdly, an aeroelastic benchmark test was considered, i.e. the Goland’s wing [Reference Goland24]. This test case was widely used in literature to validate aeroelastic codes. The study of this low-aspect-ratio wedged wing with aspect ratio $ \approx 3.33$ is also of interest because it highlights the effect of using a 2D or a 3D aerodynamic model on flutter computation. This aeroelastic model was investigated in a previous work by Savino et al. [Reference Savino, Cocco, Zanotti, Tugnoli, Masarati and Muscarello12] using DUST aerodynamic models built with classical VL and SP elements coupled to the multibody solver MBDyn code. Indeed, as a matter of fact, the use of LL approach for this kind of coupled simulation failed due to the implicit numerical scheme that characterises this element. The present activity aimed to extend the DUST-MBDyn capabilities by including the use of the non-linear vortex lattice method for the evaluation of wing flutter. In particular, to study the flutter instability, a non-zero angle-of-attack of 0.05° was introduced as a perturbation, as done by Murua et al. [Reference Murua, Palacios and Graham25]. The frequency and damping of the model response were identified from the time history of the wing-tip deflection using the matrix pencil estimation (MPE) method [Reference Hua and Sarkar26]. Table 1 reports a comparison of flutter speed and frequency computed by several authors for this test case. Results from the coupled DUST-MBDyn code are in quite good agreement with those obtained by similar codes using 3D aerodynamic models [Reference Murua, Palacios and Graham25, Reference Patil, Hodges and Cesnik27, Reference Wang, Chen, Liu, Mook and Patil28], while a certain difference of flutter speed occurs with respect to calculations obtained with 2D aerofoil models. In particular, results obtained using the NL-VL validated the suitability of this approach to investigate aeroelastic problems. In detail, the discrepancy with results obtained with the same multibody structural model MBDyn, but using its built-in aerodynamic model based on two-dimensional unsteady strip theory, indicates the higher capability of the coupled code based on surface elements as SP and VL for the investigation of aeroelastic problems over low aspect ratio wings.

Table 1. Comparison of flutter speed and frequency computed for Goland’s wing

Further validation of the coupled code DUST-MBDyn approach based on the use of different aerodynamic models including the NL-VL approach is provided in Figs. 7 and 8 showing the computed frequency and damping of the first beam torsional mode of the wing as functions of the free-stream speed.

Figure 7. Frequency $f$ of the first torsional modes vs speed ${\boldsymbol{U}_\infty }$ for Goland’s wing. Coupled DUST-MBDyn simulation results (VL, NL-VL and SP mesh) and MBDyn results with 2D strip theory aerodynamic model.

Figure 8. Damping ( $\xi $ ) of the first torsional modes vs speed ${\boldsymbol{U}_\infty }$ for Goland’s wing. Coupled DUST-MBDyn simulation results (VL, NL-VL and SP mesh) and MBDyn results with 2D strip theory aerodynamic model.

Numerical results of the coupled simulations obtained using a panel mesh (SP) show slightly higher aerodynamic damping than those obtained using both VL approaches. An increase in the predicted flutter speed of approximately 3.7% is also observed. Given these minor differences in the results obtained with the two models, a vortex lattice mesh appears to be more convenient than a surface panel one, as the computational cost is reduced significantly with no significant loss in accuracy. In particular, the use of the non-linear VL elements could provide better performance for test cases characterised by more consistent viscous effects. In the considered case, as angles of attack are small and within the linear range of the aerofoil, the viscous correction is not significant. It is instead important to consider a correction for compressibility since flutter occurs at a Mach number around 0.5. In particular, for surface panels and classical vortex lattice models, Prandtl-Glauert correction was used, while for NL-VL model compressibility effect is evaluated considering 2D aerodynamic tables.

3.1 XV-15 proprotor

The proprotor of the XV-15 tiltrotor with metal blades was simulated using NL-VL elements in different flight conditions, i.e. hover condition, forward flight in helicopter mode, and aircraft mode. For the present test case, the three propeller blades were modeled as rigid blades using non-linear vortex lattice elements only. The tabulated 2D aerodynamic data of the propeller aerofoils were taken from Felker et al. [Reference Signor, Felker and Betzina29]. DUST simulations were performed considering 10 propeller revolutions with a time discretisation of 5 ${^\circ}$ of blade azimuthal angle. In the following, DUST simulations results obtained using the NL-VL approach are compared to the recent high-fidelity numerical simulations results obtained by Jia et al. [Reference Jia, Moore and Wang30] using a Detached Eddy Simulation (DES) approach and to the experimental data collected in the test campaigns described by Felker et al. [Reference Signor, Felker and Betzina29] and Betzina [Reference Betzina31].

3.1.1 Hover flight condition

Figures 9 and 10 show the comparison of the rotor thrust coefficient ( ${C_T}$ ) as a function, respectively, of blade collective angle and of rotor torque coefficient ( ${C_Q}$ ). Moreover, Fig. 11 shows the comparison of the rotor figure of merit ( $FM$ ).

Figure 9. Comparison of the ${C_T}/\sigma $ vs collective angle for the XV-15 proprotor in hover. Experimental data from Ref. [Reference Signor, Felker and Betzina29] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

Figure 10. Comparison of the ${C_Q}/\sigma $ vs ${C_T}/\sigma $ for the XV-15 proprotor in hover. Experimental data from Ref. [Reference Signor, Felker and Betzina29] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

Figure 11. Comparison of the ${F_M}$ vs ${C_T}/\sigma $ for the XV-15 proprotor in hover. Experimental data from Ref. [Reference Signor, Felker and Betzina29] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

The performance curves obtained with DUST describe quite well the behaviour of the experimental data in the whole range of blade collective angles tested. Moreover, the comparison of the curves shows that the DUST approach provides similar capabilities to a DES approach in terms of aerodynamic performance evaluation for such a case, but requires a quite lower amount of computational effort. The figure of merit comparison depicted in Fig. 11 shows a very good agreement with experiments of the overall rotor performance computed by DUST using NL-VL, particularly at low ${C_T}/\sigma $ . Moreover, the lower computational cost required by the DUST approach allows us to perform a decidedly higher amount of simulations with respect to high-fidelity CFD, thus covering with a finer step of collective angle the whole operational range described by the experimental curve. Indeed, the computational time required to complete the simulation of the rotor configuration was about 8 minutes using a workstation with a Dual Intel Xeon Gold 6230R @2.10Ghz with 104 cores processor.

A flow field representation of the present test condition is presented in Fig. 12, showing the helical vortical structure of the proprotor wake in hover computed by DUST highlighted by the iso-surfaces of Q-criterion.

Figure 12. Wake visualisation of the XV-15 proprotor in hover at ${\theta _{075}} = {15^ \circ }$ by means of iso-surfaces of Q-criterion computed by DUST coloured by Mach number.

3.1.2 Forward flight condition

In forward flight conditions, helicopter mode configurations of the XV-15 proprotor were investigated with DUST considering three shaft angle attitudes, i.e. $\alpha = - {5^ \circ }$ , $\alpha = {0^ \circ }$ and $\alpha = {5^ \circ }$ , at advance ratio 0.17. This choice enabled to investigate DUST capabilities in both propulsive and descending forward flight conditions, analogously to what was done in the work by Jia et al. [Reference Jia, Moore and Wang30] that was here considered together with experimental data for the assessment of NL-VL performance capabilities. Figures 13 and 14 show the comparison of the rotor torque coefficient ( ${C_Q}$ ) as a function, respectively, of rotor thrust coefficient ( ${C_T}$ ) and of lift coefficient ( ${C_L}$ ).

Figure 13. Comparison of the ${C_Q}/\sigma $ vs ${C_T}/\sigma $ for the XV-15 proprotor in forward flight. Experimental data taken from Ref. [Reference Betzina31] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

Figure 14. Comparison of the ${C_Q}/\sigma $ vs ${C_L}/\sigma $ for the XV-15 proprotor in forward flight. Experimental data from Ref. [Reference Betzina31] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

A quite good agreement of DUST NL-VL simulations results with experiments is found also for this flight condition. In particular, DUST results slightly underestimate the experimental performance curves in the whole range of collective blade angles considered. Nevertheless, the discrepancies of DUST results from experimental data are relatively lower with respect to CFD simulation results, particularly in the lower range of collective angles tested. On the other hand, at higher thrust levels or lift coefficient ranges, larger discrepancies with respect to experimental data are found for DUST representation of the proprotor aerodynamic performance. This behaviour could be related to the limitation of the DUST approach to reproduce accurately high blade loading conditions characterised by significantly separated flow regions. Moreover, these discrepancies could be also due to high three-dimensional flow effects occurring in this configuration.

Indeed, these flight conditions are also characterised by relevant aerodynamic interactional effects, as shown by the flow visualisation presented in Fig. 15 highlighting the interaction of the tip vortices with the downstream blades.

Figure 15. Wake visualisation of the XV-15 proprotor in advanced flight at 5 ${^\circ }$ collective, and $\alpha = {5^ \circ }$ by means of iso-surfaces of Q-criterion computed by DUST coloured by Mach number.

3.1.3 Aircraft mode flight condition

Aircraft mode configurations of the XV-15 proprotor were simulated using DUST at an advance ratio of 0.337 for different collective blade pitch angles. The capabilities of the DUST NL-VL approach to reproduce the aerodynamic performance of the proprotor in this flight condition are evaluated by comparing the torque coefficient ( ${C_Q}$ ) and propulsive efficiency ( $\eta $ ) as a function of the thrust coefficient ( ${C_T}$ ), as can be seen respectively in Figs. 16 and 17.

Figure 16. Comparison of the ${C_Q}/\sigma $ vs ${C_T}/\sigma $ for the XV-15 proprotor in aircraft mode flight. Experimental data from Ref. [Reference Betzina31] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

Figure 17. Comparison of the propulsive efficiency eta vs ${C_T}/\sigma $ for the XV-15 proprotor in aircraft mode flight. Experimental data from Ref. [Reference Betzina31] (Exp), numerical data from Ref. [Reference Jia, Moore and Wang30] (CFD) and DUST (NL-VL).

A very good agreement between DUST simulation results and experimental data is found also for this flight conditions. In particular, the discrepancies with respect to experimental curves exposed by DUST simulations results are quite lower with respect to the ones obtained by high-fidelity CFD simulations, thus confirming the suitability of DUST NL-VL approach for an accurate evaluation of propellers’ aerodynamic performance. Indeed, a slight overestimation of the propulsive efficiency evaluated by experiments can be observed from DUST simulation results over almost the whole range of rotor thrust conditions tested, while a higher discrepancy in the order of 20% is found only for the lowest blade loading condition tested.

The flow field representation of an aircraft mode flight condition computed by DUST is shown in Fig. 18 by means of iso-surfaces of Q-criterion, highlighting a quite coherent helical vortical structure of the proprotor wake without interactions due to the free-stream velocity convection characterising this flight condition.

Figure 18. Wake visualisation of the XV-15 proprotor in aircraft mode by means of iso-surfaces of Q-criterion computed by DUST coloured by Mach number.

3.2 Tandem propellers

The second rotorcraft application concerns the aerodynamic interaction between multiple propeller configurations typical of eVTOL vehicles. In particular, two propellers in tandem with different overlapping of the rotor disk were simulated by DUST using LL, SP and the NL-VL elements. DUST simulation results were compared with experimental data available from a wind tunnel campaign performed at the S. De Ponte wind tunnel of Politecnico di Milano. The experimental setup of the tandem propeller models inside the wind tunnel test section is shown in Fig. 19. Details about the propeller models and the experimental campaign can be found in the following works [Reference Zanotti32, Reference Zanotti and Algarotti33].

Figure 19. Layout of the experimental setup of the tandem propellers at S. De Ponte wind tunnel of POLIMI [Reference Zanotti and Algarotti33].

The axial distance ( ${L_x}$ ) between the two propeller disks was set equal to five rotor radii, while several lateral separation distances ( ${L_y}$ ) between propeller rotation axes were considered. A scheme of the tandem propellers layout is depicted in Fig. 20. DUST simulations reproduced wind tunnel test conditions that consisted of runs performed with tandem co-rotating clockwise propellers both with a rotational speed controlled to 7050 RPM. This RPM target value was considered to reproduce a typical tip Mach number, i.e. ${M_t} = 0.325$ , of full-scale eVTOL aircraft propellers in cruise flight conditions.

Figure 20. Layout of the tandem propellers numerical model built for DUST simulations.

The blade geometry used to build the numerical model of the propeller was digitally created by means of a 3D scanning of the blade model. In particular, CAD software was used to generate the blade geometry from the surfaces provided by the scanning system. The maximum difference between the reconstructed blade CAD geometry and the 3D scanned surfaces was below 0.1 mm. A total number of 12 sections were extracted along the span direction. For each section, the aerofoil geometry was extracted and the distribution of twist, chord, sweep and dihedral was derived along the blade radial coordinate (r), as reported in Table 2. The blade numerical model was built considering aerofoils of the GOE, NACA and MH series reproducing the sectional geometries derived from the scan. The selected aerofoils as well as their spanwise positions along the blade are reported in Table 2. The 2D aerodynamic coefficients of the selected aerofoils were calculated by XFOIL simulations [Reference Drela34] in the range of angles of attack before stall. The method presented in Ref. [Reference Battisti, Zanne, Castelli, Bianchini and Brighenti35] was used to calculate the post-stall behaviour of the sectional aerodynamic load coefficients in the angle-of-attack range between $ \pm 180^\circ $ .

As previously stated three different numerical models were used for DUST simulations, i.e. using LL elements, surface panels (SP), and NL-VL. For each of the three blades, both LL, NL-VL and SP models have a total of 50 elements in the spanwise direction. In particular, for the NL-VL model, a chordwise discretisation of 5 elements was used, while 30 elements were used for the SP case. The spinner-nacelle surface was modeled with 1212 surface panel elements. In particular, Fig. 21 shows the layout of the propeller mesh built for DUST simulations with blades discretised with vortex lattice elements.

All DUST simulations were performed considering a length of 10 propeller revolutions with a time discretisation of 4 ${^\circ }$ of blade azimuthal angle. The computational time required to complete the simulation of a single propellers configuration was about 10 minutes for the LL mesh and about 20 minutes for the NL-VL using a workstation with a Dual Intel^® Xeon Gold 6230R @2.10GHz processor with $52$ physical cores and $2$ threads for each core. A tandem propellers configuration simulation required about 17 minutes for the LL mesh and about 40 minutes for the NL-VL. Thus, the use of the NL-VL method does not introduce a consistent increase of computational effort, particularly if compared to the time required for the simulation of the same propeller configurations by high-fidelity CFD approach [Reference Caccia, Abergo, Savino, Morelli, Zhou, Gori, Zanotti, Gibertini, Vigevano and Guardone36].

Before analysing the tandem propeller configuration, DUST simulation results obtained over simulations of the single propeller are presented. Figure 22 shows the comparison of the thrust coefficient ${C_T}$ and power coefficient ${C_P}$ as a function of advance ratio $J$ obtained with the different DUST elements compared to experimental measurements.

Table 2. Aerofoils sections, chord and twist distributions along the propeller blade span

Figure 21. Layout of the propeller model mesh with blades discretised with vortex lattice elements.

Figure 22. Single propeller performance evaluated by DUST compared to experimental data, ${M_t} = 0.325$ .

The performance curves of the single propeller obtained with LL and NL-VL elements are in quite good agreement with experimental results, while a net overestimation of both thrust and power is obtained with the blade SP elements model. This underlines the importance of viscosity effects for this propeller flight condition that are well captured by LL and NL-VL methods. In particular, the NL-VL model performs slightly better than the LL model in terms of thrust coefficient capture at high advance ratios, probably due to a better representation of the pressure distribution on blade surface thanks to chord-wise discretisation, while a very similar behaviour of the models can be found considering the evaluation of the power coefficient in the whole range of advance ratio tested.

Then, the tandem propellers test case was investigated, particularly by focusing the attention on an advance ratio $J$ equal to 0.8 representing a typical eVTOL aircraft cruise flight condition. The rear propeller performance computed by DUST for different overlapping degrees of the propellers disks is presented in Fig. 23 as a function of the lateral distance ${L_y}$ and compared with experimental data [Reference Zanotti32]. In particular, the average values of the rear propeller thrust and power coefficients were normalised with respect to the values obtained from the simulations of the single propeller configuration. This representation enables to highlight the interactional effect provided by the front propeller slipstream on the rear propeller performance.

Figure 23. Rear propeller performance evaluated by DUST compared to experimental data as a function of the lateral distance ${L_y}$ and ${L_X} = 5$ R, $J = 0.8$ , ${M_t} = 0.325$ .

The experimental data indicate that when the lateral separation distance between propellers is set to ${L_y}$ = 2R, the normalised performance coefficients tend to converge towards unity, indicating that aerodynamic interactions have a negligible effect on rear propeller performance in this configuration. As the degree of overlapping between the propeller disks is increased to ${L_y} = 1R$ , the experimental results show only minor losses in performance with respect to single propeller performance, i.e. below 5%. When the lateral distance between propellers decreases below one radius, the curves’ slope increases significantly due to the negative impact of aerodynamic interaction effects on the rear propeller’s performance. In particular, for the co-axial propeller configuration, the performance of the rear propeller decreases by nearly 30% and 20% for thrust and power coefficients, respectively. The physical explanation of the occurring phenomena causing this loss of rear propeller performance can be provided by the analysis of the local velocities computed on rear propeller blade as described in detail in Ref. [Reference Zanotti and Algarotti33]. In particular, this analysis showed that the ingestion of front propeller slipstream provides a remarkable increase of the axial velocity component experienced by the rear propeller blade as well as a slight variation of tangential velocity with respect to single propeller configuration. The combination of these effects are responsible of a decrease of the local effective angle-of-attack occurring on rear propeller blade and responsible of the large performance loss observed for this co-axial configuration. This effect on rear propeller performance becomes larger for co-axial configuration, as rear propeller disk is fully invested by front propeller disk, while increasing the vertical distance the effect on rear propeller performance decreases due to the lower degree or overlapping between propeller disks.

The comparison of DUST simulation results shows a quite good agreement of the LL and NL-VL computations with experimental data along the entire range of lateral separation distances tested. On the other hand, the SP model presents greater discrepancies with respect to experimental values, particularly for the configurations characterised by the higher degree of disks overlapping. This confirms the suitability of the models including viscous corrections for the performance evaluation of rotorcraft test cases with massive aerodynamic interactions between wakes. In particular, the NL-VL model performs better than the LL model for the evaluation of the rear propeller thrust loss in co-axial configuration, where the higher degree of interactional effects due to the front propeller slipstream occurs. On the other hand, LL and NL-VL models behave almost similarly concerning the evaluation of rear propeller power loss.

Further assessment of DUST capabilities to investigate the present problem can be provided by the comparison of the flow fields computed using the different numerical models built for the tandem propellers. In particular, Fig. 24 shows the comparison of the average freestream ( $u$ ) velocity components computed for the co-axial propellers configuration, representing the most demanding test case in terms of the interactional effect on performance. The numerical flow fields are here compared also to experimental data obtained by stereo PIV surveys [Reference Zanotti and Algarotti33].

Figure 24. Comparison of the averaged axial velocity components for tandem propeller configurations ${L_x} = 5R$ and ${L_y} = 0R$ at $J = 0.8$ and ${M_t} = 0.325$ . PIV (left) and DUST (right).

PIV results obtained for the co-axial tandem propeller configuration ( ${L_y} = 0$ ) reveal that the slipstream generated by the front propeller causes the acceleration of the rear propeller outer wake region (see Fig. 24). The acceleration provided by both propellers is accurately captured by DUST simulations with LL and NL-VL models. On the other hand, the SP numerical model provides a quite higher value of flow acceleration related to propeller blade rotation. Thus, flow field comparison further confirms the need for a numerical model including viscous corrections as the NL-VL element for an accurate investigation of complex interactional aerodynamics phenomena typical of novel rotorcraft vehicle configurations. An additional aspect to be considered is the capability of the NL-VL element to provide accurate pressure distributions along both spanwise and chordwise coordinates of the blades. This represents a promising feature that could improve the performance of aeroacoustics computation made with FWH approach with respect to the use of 1-D numerical elements as lifting lines or potential 2D elements as surface panels [Reference Caccia, Abergo, Savino, Morelli, Zhou, Gori, Zanotti, Gibertini, Vigevano and Guardone36].