Open Submitted: 22 March 2021 Accepted: 11 June 2021 Published Online: 02 July 2021
AIP Advances 11, 075009 (2021); https://doi.org/10.1063/5.0051080
more...View Affiliations
View Contributors
• George Loubimov
• Michael Kinzel

Early aerodynamic mathematical methods rely on potential flow concepts that formally isolate aerodynamic drag into a profile and an induced drag component. The more recent evolution of Navier–Stokes-based Computational Fluid Dynamics (CFD) methods typically directly computes aerodynamic forces. It does so using surface integration of pressure and viscous forces, which does not readily enable conventional separation of profile and induced drag. Isolating induced drag from aerodynamic drag is not well developed using CFD, leading to the present effort that derives a mathematical framework to extract induced drag from CFD model results. The present approach builds on mass, momentum, and energy equation control volume analysis performed within the CFD results. We find that the energy equation provides the necessary means for the closure of quantifying induced drag within the context of minor assumptions. In addition, the results indicate an interesting character in the development of energy losses in the context of induced drag associated with flow reorganization into a tip vortex. The results from the approach indicate accuracy in the method as displayed with good correlation to predictions from analytical and potential flow methods for a variety of wing planforms. Results indicate a novel and efficient method to extract induced drag from CFD models.
It is of interest to engineers to identify individual contributions of drag to an arbitrary body because specific design changes can be applied to lower these individual aerodynamic losses. A prime example of this is the finite wing, which experiences a combination of viscous and inviscid drag. These pressure and viscous forces can readily be computed from Computational Fluid Dynamics (CFD) simulations at a discrete level, where surface integration provides the bulk aerodynamic forces. In classical aerodynamics, it is conventional to decompose the two-dimensional, airfoil, forces from the three-dimensional, induced, effects. However, in the context of CFD, such decoupling is non-trivial. This presents a gap for CFD simulations in the inability to isolate effects of cross-sectional shape and planform shape that are critical in the design of wings. This study aims at closing that gap through a semi-analytical method to directly compute induced and profile drag components on wings through viscous CFD simulations.
A plethora of methods exists to estimate individual drag contributions. In the context of numerical methods, potential flow methods serve as an excellent, computationally inexpensive option for estimating induced drag effects. These range from low order analytical models such as Prandtl’s Lifting Line Theory (LLT) and11. J. D. Anderson and J. D. Anderson, Jr., A History of Aerodynamics: And its Impact on Flying Machines (Cambridge University Press, 1998), Vol. 8. lifting surface methods to vortex lattice methods.2,32. R. M. James, “On the remarkable accuracy of the vortex lattice method,” Comput. Methods Appl. Mech. Eng. 1, 59–79 (1972). https://doi.org/10.1016/0045-7825(72)90021-73. J. A. Cole, M. D. Maughmer, M. Kinzel, and G. Bramesfeld, “Higher-order free-wake method for propeller–wing systems,” J. Aircr. 56, 150–165 (2019). https://doi.org/10.2514/1.c034720 These methods provide fidelity to address how shed vorticity interacts with the circulation about lifting surfaces and enables the calculation of induced effects. These methods are also coupled to integral-boundary-layer solvers4–64. H. Schlichting and K. Gersten, Boundary-Layer Theory (Springer, 2016).5. K. Stewartson, “Further solutions of the Falkner-Skan equation,” Math. Proc. Cambridge Philos. Soc. 50, 454–465 (1954). https://doi.org/10.1017/S030500410002956X6. B. Thwaites, “Approximate calculation of the laminar boundary layer,” Aeronaut. Q. 1, 245–280 (1949). https://doi.org/10.1017/s0001925900000184 to solve for the corresponding adjusted profile drag and provide an estimate to the total drag distribution. A key attribute to potential flow methods is its foundation upon elementary functions that can be superimposed, which enables the direct assessment of induced and profile drag. Other force-decomposition approaches include the works of Chang77. C.-C. Chang, “Potential flow and forces for incompressible viscous flow,” Proc. R. Soc. London, Ser. A 437, 517–525 (1992). https://doi.org/10.1098/rspa.1992.0077 and Niu.88. Y.-Y. Niu, T.-S. Hsu, C. T. Hsieh, C. C. Chang, and C. C. Chu, “How does a gurney flap enhance the aerodynamics forces?” AIAA J. 48, 2710–2714 (2010). https://doi.org/10.2514/1.j050437 Both authors explored identifying drag-elements in viscous flows using potential flow models.
Navier–Stokes-based CFD methods directly include induced-drag contributions, but identifying each component is not straightforward due to the underlying formulation approach. One estimation approach utilizes CFD under inviscid approximations where the net body force should be equal to the induced drag and, in the context of inviscid assumptions, provides an estimate of the inviscid drag behavior.99. W. F. Phillips, S. R. Fugal, and R. E. Spall, “Minimizing induced drag with wing twist, computational-fluid-dynamics validation,” J. Aircr. 43, 437–444 (2006). https://doi.org/10.2514/1.15089 This particular approach is computationally inefficient as two individual simulations need to be completed and compared. The present method outlined by Phillips99. W. F. Phillips, S. R. Fugal, and R. E. Spall, “Minimizing induced drag with wing twist, computational-fluid-dynamics validation,” J. Aircr. 43, 437–444 (2006). https://doi.org/10.2514/1.15089 requires two simulations: one simulation requires inviscid assumptions to estimate the induced drag while the second solves the problem with viscous assumptions and by default provides one force that is composed of viscous and inviscid drag. Second, this is only an approximation to the behavior of the induced drag as the effects of the viscous conditions are neglected. Finally, solving Euler equations also requires added dissipation terms1010. T. H. Pulliam, “Artificial dissipation models for the Euler equations,” AIAA J. 24, 1931–1940 (1986). https://doi.org/10.2514/3.9550 as Euler equations do not naturally possess dissipation terms like Navier–Stokes equations do. This modeled dissipation possesses the capability to incur numerical error into the inviscid analysis. Hence, it is of value to explore the possibility of extracting both the profile and induced drag components from a single viscous simulation to reduce computational cost and have an understanding of viscous effects on the induced drag.
Novel approaches on estimating individual drag contributions have stemmed from early experimental work by Betz,1111. A. Betz, “A method for the direct determination of wing-section drag” Report No. 337, November 1925, available at https://digital.library.unt.edu/ark:/67531/metadc65142/m1/3/. Jones,1212. B. M. Jones, The Measurement of Profile Drag by the Pitot-Traverse Method (HM Stationery Office, 1936). and Maskell.1313. E. Maskell, “Progress towards a method for the measurement of the components of the drag of a wing of finite span,” RAE Technical Report No. 72232, 1972. In their work, far-field wake survey methods were used to estimate the profile and induced drag coefficients in a control volume (CV) approach. Using high-fidelity CFD simulations with well refined wake regions, it is possible to replicate these approaches with accurate results.14,1514. C. Marongiu and R. Tognaccini, “Far-field analysis of the aerodynamic force by lamb vector integrals,” AIAA J. 48, 2543–2555 (2010). https://doi.org/10.2514/1.j05032615. M. Méheut and D. Bailly, “Drag-breakdown methods from wake measurements,” AIAA J. 46, 847–862 (2008). https://doi.org/10.2514/1.29051 Destarac and Van Der Vooren1616. D. Destarac and J. Van Der Vooren, “Drag/thrust analysis of jet-propelled transonic transport aircraft; definition of physical drag components,” Aerosp. Sci. Technol. 8, 545–556 (2004). https://doi.org/10.1016/j.ast.2004.03.004 were able to conduct a thorough drag decomposition by individually modeling each source of drag (i.e., viscous, wave, and induced drag) using boundary and volume integrals. Cummings et al.1717. R. Cummings, M. Giles, and G. Shrinivas, “Analysis of the elements of drag in three-dimensional viscous and inviscid flows,” in 14th Applied Aerodynamics Conference (AIAA, 1996), p. 2482. explored utilizing a momentum balance in CFD for both near- and far-field to decompose drag. However, a downside to using these approaches (in the context of CFD) is the requirement of a highly refined wake region. It is of interest to alleviate this computational load by using a method that can estimate the individual drag coefficients in a near-field approach.
Within CFD, conventional models employ the coupled solutions of the mass, momentum, and energy equations (and mass and momentum for incompressible flow). Forces are gathered either by force integration or through control volume analysis using the momentum equations. When using the momentum deficit, one determines the reactionary forces on an arbitrary body of interest. The energy equation, however, is critical to these assessments and can provide additional insight. Consider the example of the propeller; although forces are estimated by conducting a momentum balance, the increase in momentum in the wake is also captured by the energy increase at the propeller disk. In the context of induced drag, the energy equation can also be utilized, and it is hypothesized that the energy budget can directly elucidate induced drag and that the character aligns with the flow redistribution in the tip vortices. Hence, it is proposed that using a combination of the momentum and energy equations, additional insight into wing loading can be developed.
There is a long history of utilizing energy equation aerodynamics; however, it appears to be mostly limited to early theoretical development. For example, by processing through standard texts, we can see the application of energy equation in the context of propellers,1818. W. Barnes and W. McCormick, Aerodynamics Aeronautics and Flight Mechanics (Wiley, 1995). helicopters,19,2019. B. W. McCormick, Aerodynamics of V/STOL Flight (Dover Publishing, 1999).20. G. J. Leishman, Principles of Helicopter Aerodynamics with CD Extra (Cambridge University Press, 2006). and other propulsion systems. The evaluation of the energy equation to extract insight from CFD for aerodynamic losses is much less studied and is partially developed. Extensive work was carried out by Drela2121. M. Drela, “Power balance in aerodynamic flows,” AIAA J. 47, 1761–1771 (2009). https://doi.org/10.2514/1.42409 who explored decomposing tightly coupled boundary layer ingestion propulsion. Drela balanced sources and sinks of power to further understand interference drag. In contrast to using the fluid momentum equations, the energy equation can be utilized to provide a concise accounting of fluid energy and allows for the decomposition of forces in tightly integrated systems. In his work, various energy sources and sinks are estimated and balanced, which gives insight into the efficiency of individual systems operating in a tightly coupled setting. Arntz et al.22,2322. A. Arntz, O. Atinault, and A. Merlen, “Exergy-based formulation for aircraft aeropropulsive performance assessment: Theoretical development,” AIAA J. 53, 1627–1639 (2015). https://doi.org/10.2514/1.j05346723. A. Arntz and D. Hue, “Exergy-based performance assessment of the NASA common research model,” AIAA J. 54, 88–100 (2016). https://doi.org/10.2514/1.j054127 expanded on the work by Drela by considering the exergy involving mechanical and thermal-energy balances. In his work, exergy balances were conducted using CFD methods on the NASA common research model. Additional studies by Williams et al.2424. D. M. Williams, D. S. Kamenetskiy, and P. R. Spalart, “On stagnation pressure increases in calorically perfect, ideal gases,” Int. J. Heat Fluid Flow 58, 40–53 (2016). https://doi.org/10.1016/j.ijheatfluidflow.2015.12.005 evaluated increases in stagnation pressure above inlet conditions. Williams utilized the energy equation to observe how stagnation pressure can increase on the outer edges of shear layers. The works of Drela, Arntz, and Williams depict an energy tracking approach that allows for a refined analysis of fluid physics, which a momentum deficit cannot provide. Furthermore, in the context of CFD, studies of the energy budget have not been explicitly used to predict induced drag and, thus, are of interest for exploration.
In the case of induced drag, control volume analysis can provide a variety of flow insights. The induced drag is characterized by reorganization of fluid momentum, which results in an induced velocity that creates lift-induced drag. The reorganization of this fluid momentum implies that the energy equation has the potential to separate these effects from viscous losses. Through force integration, computing by integrating local surface pressure and viscous forces demands integration over a closed body and cannot isolate induced from profile drag. The proposed methodology is fundamentally different as it isolates profile drag using energy-equation-based control volume assessments. The profile drag can be subtracted from the total drag (computed from pressure/shear integration) to elucidate the induced drag. Through this approach, we alleviate far-field resolution requirements and isolate mesh resolution requirements to the near field. This paper approaches this goal as follows: We initiate with the model formulation. We then develop a simple CFD model to highlight the effectiveness of the derived methodology. This is followed by wing-shaped model variations that span the wing planform. The results are then compared to potential flow models.
A. Mathematical formulation
Computing the net force on an arbitrary body through surface integration in CFD is well defined and is relatively straight forward. However, it is not possible to decompose the components into profile and induced components. This can be recovered through control volume analysis,19,2519. B. W. McCormick, Aerodynamics of V/STOL Flight (Dover Publishing, 1999).25. I. G. Currie, Fundamental Mechanics of Fluids (CRC Press, 2016). as detailed below. A variable control volume is used to assess the momentum and energy equations within the context of a steady CFD solution without gravity. The variable aspect of the control volume fixes the inlet control surface while evaluating control surface outlets across various cross-sectional planes within the computational domain (see Fig. 1). The cross-sectional plane crosses through the wing (on body assessments) as well as the vortex wake (off-body). Discrete evaluations of the momentum and energy on these variable control volumes enable spatial analysis of the integral effects on fluid momentum and energy. They will be referred to as spatial momentum conservation (SMC) and spatial energy conservation (SEC). In the context of our variable control volume, it is important to note that no-penetration conditions are placed on far-field surfaces parallel to the free-stream velocity (much like a large wind tunnel). This control volume simplifies the assessment to only the outlet. When utilizing this method with true far-field boundary conditions, the approach must also include the energy flux in/out of the domain. In this study, we focus on the x-direction momentum (aligning with the free-stream velocity) as follows:
 $D(x)=∬CS(x)ρ(V∞−Vx)V⃗⋅n⃗dA+∬CS(x)(p∞−p)dA,$ (1)
where D(x) is the drag evaluation on the cross section CS(x) as indicated in Fig. 1. The normal vector perpendicular to the surface of CS(x) is denoted as $n⃗$, and the area of the cross section is given by A. The fluid density is given by ρ, $V⃗$ is the velocity vector (Vx being the velocity component in the x-axis), and p is the pressure. Equation (1) describes the change in forces along the x-axis based on the change in x-axis fluid momentum. This can be expressed in terms of a drag coefficient (by dividing by $12ρV∞2S$), with an induced and profile drag component, as follows:
 $CD=2DρV∞2S.$ (2)
We next consider the drag coefficient as the decomposition of profile ($CD0$) and induced drag ($CDi$) as follows:
 $CD=CD0+CDi.$ (3)
The underlying issue is that these provide two equations [Eqs. (2) and (3)] and three unknowns (CD, $CD0$, and $CDi$) when developed to solve for the induced and profile drag components. Therefore, there is a gap in the ability to decompose induced and profile drag components.
The present application of the energy equation starts here; we hypothesize that the energy equation can serve to decompose the profile and induced drag components. In order to provide closure, we introduce the energy equation. The energy equation is, similarly, evaluated on various control surfaces,
 $QCV−WCV=∬CS(x)e+ρ(V∞2−V2)2V⃗⋅n⃗dA+∬CS(x)(p∞−p)V⃗⋅n⃗dA.$ (4)
In the context of an iso-thermal flow, the heat-generation term (QCV) can be discarded, leaving only the work term, WCV, on the left-hand side, which is intended to represent losses associated with drag. Under isothermal conditions, it is assumed that the internal energy per unit mass, e, is constant in this analysis; hence, its net flux is zero. By assuming an isothermal flow, the pressure changes will impact the enthalpy of the fluid; however, this purely becomes a function of pressure and density (where density is assumed to be constant). Even in the context of the isothermal model, the energy equation must still be satisfied in the context of any control volume. Although under isothermal conditions there is a decoupling of energy and state equations from momentum and mass, the energy equation must still be satisfied, and the present analysis exploits this for deeper understanding. Hence, the energy equation is still applicable and is applied in conventional aerodynamics, for example, a propeller in an incompressible flow (Lishman2020. G. J. Leishman, Principles of Helicopter Aerodynamics with CD Extra (Cambridge University Press, 2006). and McCcormick1919. B. W. McCormick, Aerodynamics of V/STOL Flight (Dover Publishing, 1999).). The energy equation can now be simplified to a control volume form assessed at each cross-sectional plane, CS(x), as follows:
 $−WCV=∬CS(x)ρ(V∞2−V2)2+(p∞−p)V⃗⋅n⃗dA.$ (5)
This relation provides the relation between fluid energy flux evaluated on CS(x) and the total work extracted within the control volume (CV). The work removed from the fluid energy can be evaluated as being proportional to the fluid power. The energy reduction is a result of work associated with the wing aerodynamics—this work term can be related to viscous dissipation or the work required to form the tip vortex. In order to interpret the work term as a force, we use the time-rate-of-change of work or the power (input or output). In this case, since a steady simulation is used, the time term is neglected, yielding the following equation:
 $WCV=PΔt∼−F⃗⋅V⃗,$ (6)
which can be used to extract the fluid force. Here, P is the power required. Using the inlet velocity, the accumulated force and drag coefficient can be solved for as
 $WCV=PΔt∼F⃗⋅V⃗≅DV∞=∬CS(x)ρ(V∞2−V2)2+(p∞−p)V⃗⋅n⃗dA.$ (7)
Note that this relation includes drag and is distinct from the momentum equation in which it utilizes the energy and not momentum flux. Rearranging this equation, in terms of drag coefficients, yields
 $CD=2ρV∞3S∬CS(x)ρ(V∞2−V2)2+(p∞−p)V⃗⋅n⃗dA.$ (8)
B. Example application
The numerical solutions are generated utilizing a commercial code, Star-CCM+.2626. STAR Siemens, CCM+ Users Manual, 2019. The Spalart–Allmaras one-equation turbulence model2727. P. Spalart and S. Allmaras, “A one-equation turbulence model for aerodynamic flows,” in 30th Aerospace Sciences Meeting and Exhibit (AIAA, 1992), p. 439. was used to model the turbulent flow. CFD results were bench-marked against lifting-line theory, and mesh-refinement studies were conducted. Mesh refinement studies suggested that the CFD solutions were minimized by numerical uncertainty. Direct comparisons against experimental data were not considered; however, the effort bench marks the CFD and theory against theoretical, text-book equations associated with lifting line theory that have been validated against a multitude of experiments over time. This validation approach is accepted in the context of verification and validation practice.2828. W. L. Oberkampf, T. G. Trucano, and C. Hirsch, “Verification, validation, and predictive capability in computational engineering and physics,” Appl. Mech. Rev. 57, 345–384 (2004). https://doi.org/10.1115/1.1767847
In order to verify that the application of the energy equation yields insight into induced drag, we first start by conducting an SMC and SEC assessment throughout a computational domain containing a wing of infinite span (or a 2D airfoil) in a control volume approach. The 2D airfoil used in the simulation is an NACA 0012, which operates at a Reynolds number of 675 000 (based on the chord length, c) and a ten-degree angle of attack (α = 10°). In Fig. 1, the cross-sectional planes represent control surfaces in which the momentum and energy equations are integrated using surface integrals. This approach allows for the observation of spatially varying momentum and energy and can be seen in Fig. 2. The plot indicates differences in the momentum and energy equation along the free-stream axis.
As shown in Fig. 2, the SEC and SMC assessments are evaluated and plotted in the context of a viscous CFD solution for a wing with infinite span. The plot in Fig. 2 is a simplified version of Eqs. (1) and (5); due to the two-dimensional simulation, an integral may be dropped in the control volume approach. In the comparison between the SEC and SMC assessments in Fig. 2, it is apparent that different methods provide different insights. The drag predicted from both methods are, as expected, virtually equal to each other as x/c (note that Cd represents the two-dimensional drag coefficient). Such an observation verifies that the approach of utilizing the energy equation to extract the drag is consistent with the momentum methods. Second, we highlight the distinct characteristics of the two curves. The momentum equation captures reactionary forces on the surface of the body. Hence, it provides a leading-edge suction force and the resulting negative drag coefficient on the front sections of the airfoil, followed by positive drag components that develop on the aft parts of the airfoil. The energy equation, conversely, displays the work carried out by the body onto the fluid. Thus, its evaluation displays the overall gradual losses along the airfoil from shear and pressure on the body. Hence, the SEC assessments eliminate the leading-edge suction force, thus enabling a spatial evaluation of the drag coefficient. In conclusion from this plot, although drag is predicted by both SMC and SEC, the resulting profiles from the evaluations differ, providing the potential for deeper insight into the wing forces.
In order to evaluate induced drag using the energy approach, a finite wing is evaluated and compared to results obtained for the infinite wing (note: CD is used to represent the three-dimensional drag coefficient). The finite wing has a rectangular planform with an aspect ratio of 10 while the free stream conditions remain identical, where Re = 675, 000 and α = 10°. The direct comparison uses the SEC approach in Fig. 3 while there is an SEC comparison to the SMC approach for only the finite wing in Fig. 4.
First, consider a comparison between the infinite (2D) and finite wing (3D) cases in Fig. 3. In this comparison, one can observe a difference that is driven by the 3D effects associated with the well understood vortical wake that drives induced drag. There is also a distinct trait for the finite wing associated with the gradual increase in the drag coefficient downstream of the trailing edge. In contrast, the infinite wing reached its loss associated with drag at the trailing edge. This additional, off-body loss mechanism arising in the SEC assessment is associated with the work required to transition from shed vortex filaments to a tip vortex.
An additional observation from Fig. 3 is that the energy associated with drag for the finite wing is less than that for the 2D configuration. Such an observation may be expected as suggested by 3D wing theory, implying that the 3D wing experiences an induced downwash and a resulting reduction in the effective angle of attack. Such reductions imply reduced sectional profile drag. In order to estimate this, we build on the previous observation from Fig. 2, i.e., that the energy loss from the infinite wing at the airfoil trailing edge is equivalent to the total drag coefficient. Continuing along this line of logic, the infinite wing is pitched down to match the effective angle of attack of the wing, under the assumption of a uniform downwash (consistent with an elliptical wing, which is not the case studied). Such a test is used to evaluate if the profile drag of a finite wing can be estimated using the energy loss from the leading edge (LE) to trailing edge (TE). The induced angle of attack is calculated using
 $αi=CLπAR$ (9)
and is subtracted from the geometric angle of attack to determine the effective angle of attack. As observed in the 2D adjusted viscous theory case (in Fig. 3), the energy loss is nearly identical to that of the finite wing at the trailing edge, or
 $CD0=2ρV∞3S∬CS(x=TE)×ρ(V∞2−V2)2+(p∞−p)V⃗⋅n⃗dA.$ (10)
This result implies that the profile drag from the finite wing is well predicted using the evaluation of the energy equation at the trailing edge.
In order to extract the induced drag coefficient, the finite wing is evaluated using SMC and SEC assessments. The comparison is provided in Fig. 4, where the various drag components are highlighted. From the SMC assessment, the $xc→∞$ result yields the overall drag component (CD). The SEC assessment, from the wing LE to TE, yields the profile drag coefficient of the wing ($CD0$). The underlying deficit, i.e., $CD−CD0$, yields the induced drag on the wing. Hence, we are able to directly assess the induced drag, provided we apply the SEC assessment on a control volume spanning the LE and TE of the wing.
In the context of the aforementioned observation, it can be observed that the x/c result from the SEC assessment is not equal to CD from the SMC assessment (Fig. 4). This disparity demands additional attention and is associated with the fact that the energy balance [Eq. (7)] does not account for the work associated with lift. As presented, Eq. (7) accounts for both the viscous and pressure work on the airfoil through drag. However, work also appears in terms of an energy addition from lift being generated by the wing, which increases the energy such that the blue line (Fig. 4) is shifted biased downward relative to the red line. The overall energy budget can be considered as an energy conversion process, where the lift is being generated (energy increase) at the expense of drag-induced energy losses. The discussion of energy increases/losses refers to discrete changes in the total energy in the fluid, as calculated in the SEC drag assessment. Drag (both profile and induced) is an energy loss mechanism and is independent from an energy increase associated with lift. Thus, the overall energy budget in the SEC must also consider these mechanisms. The result implies the observed delta occurring between the overall predictions of CD from the SMC and SEC drag assessments.
In working with SEC, it is critical to quantify all forms of work. In order to consider the lift-induced work, an assumption extracted from wing theory is applied. It is assumed that the lifting power is estimated as the product of lift and downwash (w), i.e.,
 $Py,lift=12ρV∞2S∫−b/2b/2Cl(z)w(z)dz,$ (11)
and is computed using results from the VLM solution. This lift-added power can be described as an energy increase in the SEC assessment and is an error in the previous evaluation. In fact, this term must be added to Eq. (7) to complete the energy balance, which, in the context of drag coefficient, is given as
 $CDy,lift=Plift12ρV3S.$ (12)
The addition of this term to the energy budget is provided in Fig. 5, which corrects the deficit between the SEC and SMC assessments. Note that the $CDy,lift$ term is, in reality, distributed along the chord; however, here, it is added at the trailing edge leading to an immediate increase. This is preferred as we want to retain energy associated with lift while simultaneously retaining our closed-form solution of $CD0$ in Eq. (10). It is assumed that the downwash initiates at the trailing edge and is added there. Such an assumption is valid in the context of energy addition due to the fact that the wing acts as a no-flux boundary. While the aerodynamic lift is being distributed over the wing surface, the downwash caused by the vortex system is effectively zero up until the trailing edge.
Note how the final energy loss is now equal to the momentum yielded drag coefficient. This also implies that $CDi$ can now be explicitly solved for as follows:
 $CDi=2ρV∞3S∫∫CS(x=TE)∞ρ(V∞2−V2)2+(p∞−p)×V⃗⋅n⃗dA+12ρV∞2S∫−b/2b/2Clwdz,$ (13)
or in a simplified form,
 $CDi=2ρV∞3SCDi,wake−CDi,lift.$ (14)
In Eq. (14), the first term in the bracket represents the work required to form the tip vortex throughout the wake, and the second term is the power generated by the wing in the form of lift. At this point, the utilization of the energy equation can be used independently of the momentum equations to decompose drag to viscous and inviscid components by conducting this energy balance.
Overall, the method applies classical fluid dynamics principles to modern CFD to interpret wing design. Conventionally, the net drag on a body, CD, can be computed through surface integration of the pressure and shear forces. The first approach shows how the profile drag, $CD0$, can be separated from CD. Second, it can be reworked to provide a closed-form solution to calculate $CDi$. In the context of mesh density, the application of Eqs. (3) and (10) may be more computationally efficient than Eq. (13) due to the fact that the wake must be well refined to calculate the final energy loss and $CDi$. Equation (10) purely depends on a well-refined grid near the wing to calculate the pressure and viscous work associated with the profile drag and allows for an implicit calculation of $CDi$. In this study, the induced drag coefficient is calculated by solving the system of equations utilizing Eq. (10) to close out the system.
The use of Eq. (13) depends on utilizing the downwash velocity, w, which is not trivial to determine using CFD exclusively. However, potential flow models serve as a viable and computationally affordable method to be coupled with CFD to close out this issue. With respect to the present effort, it is important to note that the purpose of Eq. (13) is to evaluate the accuracy of this concept and is used as supportive evidence.
An area of interest worth exploring is the application of an SEC assessment on an inviscid test case. Consider the SEC inviscid to viscous comparison in Fig. 6(a). The spatial variation in the energy loss (in the terms of CD) based on SEC assessments highlights insight into viscous assessments. First, note that the energy loss between the viscous and inviscid solutions until near the leading edge of the wing (x/c = 0.25) is nearly identical. This represents the wing applying a work term to the fluid in terms of pressure. Moving downstream from the leading edge, the inviscid CFD results divert from the viscous CFD due to the lack of viscous losses on the wing surface. In the viscous case, the pressure and shear forces dominate the fluid losses in the region between the leading and trailing edge. Compared to the viscous case, the energy loss in the inviscid case between x/c = 0.25 and x/c = 1 seem to be non-existent. This is expected as the inviscid simulation should be free of pressure and shear forces associated with profile drag. However, a caveat appears: why does the fluid energy decrease at the leading edge for the inviscid case? This discrepancy suggests that the inviscid case experiences energy losses along the foil boundary and poses a contradiction. It is suspected that this discrepancy is related to numerical dissipation, which is causing the inviscid model to incorrectly record losses, which represent profiled drag losses.
As one moves downstream from the trailing edge, the energy losses develop and are associated with the formation of the wing-tip vortex. One interesting observation is that the viscous and inviscid CFD profiles have an offset due to the discrepancy in viscous losses. Regardless, the energy losses in the wake mirror each other. The observation implies that the physical loss mechanisms in the wake are unaffected by losses upstream, which is consistent with potential flow theory and common decoupling of viscous mechanisms. This leads to the hypothesis that the dominant energy loss mechanism in the wakes is the formation of the wing-tip vortex.
To analyze the hypothesis that wake losses are associated with the reorganization of the flow into a tip vortex, the results in Fig. 6(a) are adjusted to highlight the energy losses aft of the trailing edge. This is carried out by evaluating the viscous and inviscid energy losses concerning the trailing edge. This approach isolates the energy losses in the wake and enables a comparison of the viscous and inviscid wakes. Results of the wake comparison are presented in Fig. 6(b).
The wake comparisons in Fig. 6(b) indicate that the presence of viscous shear in the wake results in only a subtle increase in energy loss. In fact, the slopes of both curves are nearly identical, highlighting that the work associated with the tip vortex formation is only occurring in the near wake. Due to the fact that the difference between the viscous and inviscid energy loss profiles is nearly constant, this also suggests that the viscous wake does not have a continual viscous shear loss term decreasing the fluid energy. If this were the case, the viscous energy loss profile would continue to diverge from the inviscid one. The difference in the adjusted profiles is likely due to numerical errors that increase with energy loss at the trailing edge. Ultimately, the inviscid problem results in re-organization of kinetic energy due to the induced flow. When viscosity is present, a work term associated with friction is introduced, which further affects the energy loss in the wake (as well as in the viscous boundary layer). This can be seen in Fig. 6(b); here, the difference in energy loss within the wake is caused by fluid shearing, which is not felt in the inviscid case. This difference is relatively small [based on Fig. 6(b)] and suggests that the primary energy loss is still an inviscid one (re-organization of kinetic energy). The added shearing mechanism results in a slight increase in induced-drag; however, the fact remains that there is a viscous component associated with this drag in real flows though the driver for this is still an inviscid mechanism.
Figure 7 shows the formation of the tip vortex. Here, the Q-criterion is used to capture the vortex structure and is colored by the energy equation used in this study. The color map helps visualize the increase in fluid energy; it is important to note that the energy equation is not integrated, which is why the fluid energy is increasing. The Q-criterion highlights flow fields that are likely being accelerated due to the interaction with the wing and not showing the viscous losses happening in the boundary layer or the surrounding flow. An increase in energy at the vortex is actually representative of the wing applying work to the fluid, which then formulates a tip-vortex and results in a net energy loss. Due to the plotted energy not being the integrated energy, observed increases in fluid energy can be balanced by another loss mechanism yet are not seen in this Q-criterion representation. This isosurface indicates the sharp energy loss associated with reorganizing the flow to form the tip vortices. The effect is specifically noticeable in the color contours of total energy on the tip vortices, where it can be observed that there is a quick decrease in total energy as the vortex forms followed by stabilization after it is formed.
A. Variations in the aspect ratio
In this section, the energy approach is used to extract the induced drag coefficients for rectangular wings with variations in the aspect ratio ($AR=b2S$ with variations of 5, 10, and 20). The changes in the AR are achieved while maintaining a constant chord so that the chordal Reynolds number is constant. Hence, span length (b) is varied. For completeness, the evaluations occur over a range of angles of attack (α) from zero to ten degrees. Results from the present method are presented and compared to predictions using a a vortex lattice method (VLM) (XFLR52929. See http://www.xflr5.tech/xflr5.htm for XFLR5 Panel Method; accessed 1 May 2020.) in Figs. 8 and 9. In general, it is clear that there is good agreement between the present predictions of $CDi$ for all ARs. There appears to be an increasing discrepancy in the prediction of CD0 with an increasing angle of attack that amplifies with lower AR wings. Due to this incorrect prediction given by the energy equation, the induced drag coefficient is over-predicted at larger angles of attack. Here, it is important to note that the global trend is well predicted for the numerical fidelity used in this study. Another worthy discussion point is the relationship between $CDi$ and CL; as shown in Fig. 10, the trend is captured well over multiple angles of attack, yet there are multiple outliers. These outliers suggest that either multiple points are required to fit a trend to estimate the span-wise efficiency or that further mesh refinement is required to capture a span-wise efficiency using one angle of attack. The latter point on mesh refinement will be explored in further detail in Sec. IV B.
B. Variations in the taper ratio
In this section, the variation in the taper ratio for a fixed aspect ratio of ten is examined and compared to theoretical calculations performed by McCormick1818. W. Barnes and W. McCormick, Aerodynamics Aeronautics and Flight Mechanics (Wiley, 1995). based on
 $CDi=CL2ARπ(1+δ).$ (15)
The tapered wings have a straight trailing edge to simplify the energy-based estimate of $CD0$; it was discovered that estimating $CD0$ for non-straight trailing edges was more convoluted due to the ambiguous location of the “mean trailing edge.” The SEC prediction of $CD0$ assumes that the energy post-trailing edge is a function of induced effects; this proved difficult to locate for non-straight trailing edges, which is why all tapered wings in this study have a swept leading edge and straight trailing edge.
Based on an earlier discussion, it was found that a range of angles of attack would need to be simulated to reach a better estimate of a span-wise efficiency. This was carried out for two wings of a constant aspect ratio with taper ratios of 1 and 0.5. The angle of attack was varied from zero to ten degrees. CFD predictions were again compared to VLM predictions. The results can be observed in Fig. 10, with $CDi$ plotted on the y-axis and CL plotted on the x-axis. Note that even with a slight delta between CFD and VLM predictions, the change in $CDi$ and $CDL$ due to taper modifications is almost identical. This suggests that at a bare minimum, the energy method can be used to capture span-wise efficiency by simulating multiple angles of attack.
The next step to refining this method is to utilize a single-point estimate of the span-wise efficiency. In Sec. IV A, it was observed that some outliers pose a threat to an accurate estimate of the span-wise efficiency. To further understand this behavior, a mesh refinement study is conducted. The results of this study are shown in Fig. 11, where the cell base size is plotted on the x-axis and the span-wise efficiency, δ, is plotted on the y-axis. This figure describes the predicted aerodynamic as a function of computational mesh resolution. The aerodynamics are plotted as a reference base size or a representative cell size for each computational mesh. It is clearly seen that there is a strong driver associated with span-wise efficiency and mesh refinement. This is most likely due to inaccuracy in CL as that term is squared in Eq. (15) and can cause large errors in predicting span-wise efficiency.
With a better understanding of how the mesh quality affects the accuracy of span-wise efficiency, single angles of attack will be utilized to predict span-wise efficiency for tapered wings. Two wings of a constant aspect ratio and with taper ratios of 0.5 and 0.25 are simulated at an attack of five degrees. The results are compared to McCormick’s analytical lifting line theory results in Fig. 12. The results show that CFD is capable of capturing changes in span-wise efficiency in a manner similar to that of LLT; however, the exact magnitudes vary in error. For example, at a taper ratio of 0.25, CFD significantly over-predicts the span-wise efficiency. Another, notable, point is the early prediction of diminishing returns with respect to a decreasing taper ratio. LLT predicts that a taper ratio of 0.33 yields the smallest δ; while there are only three CFD predictions, the trend points to a larger taper ratio, which lends a minimum value of δ.
This study systematically evaluates drag from a CFD prediction using both momentum and energy equations. The method uses this combination to separate the total drag coefficient into profile and induced drag components. It was shown that the energy equation provides a method to close out a system of equations that directly identifies the induced drag. The application of this approach provides a method to evaluate drag associated with the wing planform shape. Evaluating this procedure suggests that there are discrepancies between the final energy-drag and momentum drag. This brings in a need to also include work added from the wing in the form of lift. An estimate for this work was provided, which balanced out the final energy. This estimate also allowed for the induced drag coefficient to be explicitly solved for using only the energy equation. Thus, two methods were shown to provide estimates of the induced drag coefficient.
The first form involves computing the profile drag by measuring the energy loss from the leading to the trailing edge of the wing. This profile drag is then subtracted from the total force, which is a function of both pressure and shear forces acting on the wing surface. Thus, the induced drag is calculated using Eq. (3). The second form involves addressing the discrepancy between the SMC and SEC assessments downstream in the wake by estimating the work generated by the lifting wing. This is carried out by integrating the product of the sectional lift and downwash and subtracting it from the SEC assessment at the trailing edge. This consideration allows for a closed form solution of the induced drag without depending on the momentum equation.
The first form was utilized to estimate induced drag coefficients for various geometries. The induced drag coefficients were then calculated for three varying aspect ratios and a range of angles of attack and were shown to have good agreement with predictions made by the vortex lattice method. The wing geometry was again modified to include variations in the taper ratio. Instead of only comparing $CDi$, the span-wise efficiency δ was used to tie in CL predictions and was compared to lifting line theory predictions made by McCormick. Results showed that CFD has the capability to capture trends in δ but some uncertainty remains, primarily associated with estimates in $CD0$ at lower taper ratios.
The implications of this study suggest that the use of the energy equation can be viable to separate various drag components, specifically, through estimating inviscid criteria using viscous CFD simulations. While CFD-based design can be costly, this method holds value in high-fidelity simulations on aircraft at the end of the design stage. Furthermore, this method could be expanded to identify viscous and inviscid forces in aircraft geometries that are non-trivial and difficult to model using potential flow methods.
The authors acknowledge the University of Central Florida Stokes Advanced Research Computing Center for providing computational resources and support that have contributed to results reported herein. G.L. would also like to thank the DoD SMART Fellowship Program for funding this study.
The data that support the findings of this study are available from the corresponding author upon reasonable request.
1. 1. J. D. Anderson and J. D. Anderson, Jr., A History of Aerodynamics: And its Impact on Flying Machines (Cambridge University Press, 1998), Vol. 8. Google Scholar
2. 2. R. M. James, “On the remarkable accuracy of the vortex lattice method,” Comput. Methods Appl. Mech. Eng. 1, 59–79 (1972). https://doi.org/10.1016/0045-7825(72)90021-7, Google ScholarCrossref
3. 3. J. A. Cole, M. D. Maughmer, M. Kinzel, and G. Bramesfeld, “Higher-order free-wake method for propeller–wing systems,” J. Aircr. 56, 150–165 (2019). https://doi.org/10.2514/1.c034720, Google ScholarCrossref
4. 4. H. Schlichting and K. Gersten, Boundary-Layer Theory (Springer, 2016). Google Scholar
5. 5. K. Stewartson, “Further solutions of the Falkner-Skan equation,” Math. Proc. Cambridge Philos. Soc. 50, 454–465 (1954). https://doi.org/10.1017/S030500410002956X, Google ScholarCrossref
6. 6. B. Thwaites, “Approximate calculation of the laminar boundary layer,” Aeronaut. Q. 1, 245–280 (1949). https://doi.org/10.1017/s0001925900000184, Google ScholarCrossref
7. 7. C.-C. Chang, “Potential flow and forces for incompressible viscous flow,” Proc. R. Soc. London, Ser. A 437, 517–525 (1992). https://doi.org/10.1098/rspa.1992.0077, Google ScholarCrossref
8. 8. Y.-Y. Niu, T.-S. Hsu, C. T. Hsieh, C. C. Chang, and C. C. Chu, “How does a gurney flap enhance the aerodynamics forces?” AIAA J. 48, 2710–2714 (2010). https://doi.org/10.2514/1.j050437, Google ScholarCrossref
9. 9. W. F. Phillips, S. R. Fugal, and R. E. Spall, “Minimizing induced drag with wing twist, computational-fluid-dynamics validation,” J. Aircr. 43, 437–444 (2006). https://doi.org/10.2514/1.15089, Google ScholarCrossref
10. 10. T. H. Pulliam, “Artificial dissipation models for the Euler equations,” AIAA J. 24, 1931–1940 (1986). https://doi.org/10.2514/3.9550, Google ScholarCrossref
11. 11. A. Betz, “A method for the direct determination of wing-section drag” Report No. 337, November 1925, available at https://digital.library.unt.edu/ark:/67531/metadc65142/m1/3/. Google Scholar
12. 12. B. M. Jones, The Measurement of Profile Drag by the Pitot-Traverse Method (HM Stationery Office, 1936). Google Scholar
13. 13. E. Maskell, “Progress towards a method for the measurement of the components of the drag of a wing of finite span,” RAE Technical Report No. 72232, 1972. Google Scholar
14. 14. C. Marongiu and R. Tognaccini, “Far-field analysis of the aerodynamic force by lamb vector integrals,” AIAA J. 48, 2543–2555 (2010). https://doi.org/10.2514/1.j050326, Google ScholarCrossref
15. 15. M. Méheut and D. Bailly, “Drag-breakdown methods from wake measurements,” AIAA J. 46, 847–862 (2008). https://doi.org/10.2514/1.29051, Google ScholarCrossref
16. 16. D. Destarac and J. Van Der Vooren, “Drag/thrust analysis of jet-propelled transonic transport aircraft; definition of physical drag components,” Aerosp. Sci. Technol. 8, 545–556 (2004). https://doi.org/10.1016/j.ast.2004.03.004, Google ScholarCrossref
17. 17. R. Cummings, M. Giles, and G. Shrinivas, “Analysis of the elements of drag in three-dimensional viscous and inviscid flows,” in 14th Applied Aerodynamics Conference (AIAA, 1996), p. 2482. Google ScholarCrossref
18. 18. W. Barnes and W. McCormick, Aerodynamics Aeronautics and Flight Mechanics (Wiley, 1995). Google Scholar
19. 19. B. W. McCormick, Aerodynamics of V/STOL Flight (Dover Publishing, 1999). Google Scholar
20. 20. G. J. Leishman, Principles of Helicopter Aerodynamics with CD Extra (Cambridge University Press, 2006). Google Scholar
21. 21. M. Drela, “Power balance in aerodynamic flows,” AIAA J. 47, 1761–1771 (2009). https://doi.org/10.2514/1.42409, Google ScholarCrossref
22. 22. A. Arntz, O. Atinault, and A. Merlen, “Exergy-based formulation for aircraft aeropropulsive performance assessment: Theoretical development,” AIAA J. 53, 1627–1639 (2015). https://doi.org/10.2514/1.j053467, Google ScholarCrossref
23. 23. A. Arntz and D. Hue, “Exergy-based performance assessment of the NASA common research model,” AIAA J. 54, 88–100 (2016). https://doi.org/10.2514/1.j054127, Google ScholarCrossref
24. 24. D. M. Williams, D. S. Kamenetskiy, and P. R. Spalart, “On stagnation pressure increases in calorically perfect, ideal gases,” Int. J. Heat Fluid Flow 58, 40–53 (2016). https://doi.org/10.1016/j.ijheatfluidflow.2015.12.005, Google ScholarCrossref
25. 25. I. G. Currie, Fundamental Mechanics of Fluids (CRC Press, 2016). Google ScholarCrossref
26. 26. STAR Siemens, CCM+ Users Manual, 2019. Google Scholar
27. 27. P. Spalart and S. Allmaras, “A one-equation turbulence model for aerodynamic flows,” in 30th Aerospace Sciences Meeting and Exhibit (AIAA, 1992), p. 439. Google ScholarCrossref
28. 28. W. L. Oberkampf, T. G. Trucano, and C. Hirsch, “Verification, validation, and predictive capability in computational engineering and physics,” Appl. Mech. Rev. 57, 345–384 (2004). https://doi.org/10.1115/1.1767847, Google ScholarCrossref
29. 29. See http://www.xflr5.tech/xflr5.htm for XFLR5 Panel Method; accessed 1 May 2020. Google Scholar