Skip to content
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
23 changes: 12 additions & 11 deletions JOSS/paper.md
Original file line number Diff line number Diff line change
Expand Up @@ -29,13 +29,13 @@ Analyzing the aerodynamics of unsteady airfoils is essential for understanding t

We introduce PANKH **(Panel Analysis of uNsteady Kinematics of Hovering airfoils)**, an open-source C++ tool that employs the unsteady vortex panel method to evaluate aerodynamic forces on airfoils in arbitrary motion. Its flexible, modular design enables users to specify custom kinematic patterns (impulsive motion, pitching, plunging), ideal for exploring bio-inspired flapping flight and gust responses.

The solvers source code, validation examples, and comprehensive Doxygen-generated API documentation are hosted on GitHub, ensuring accessibility and reproducibility. Future enhancements may include expanded input capabilities, improved wake modeling techniques for more accurate unsteady flow prediction, and the integration of viscous effects along with two-way fluid–structure interaction (FSI) to simulate flexible airfoils and their coupling with the flow in a strongly coupled manner.
The solver's source code, validation examples, and comprehensive Doxygen-generated API documentation are hosted on GitHub, ensuring accessibility and reproducibility. Future enhancements may include expanded input capabilities, improved wake modeling techniques for more accurate unsteady flow prediction, and the integration of viscous effects along with two-way fluid–structure interaction (FSI) to simulate flexible airfoils and their coupling with the flow in a strongly coupled manner.

# Statement of Need

For low-Mach, incompressible applications—where viscous phenomena are largely confined to thin boundary layers and wakes—potential-flow solvers provide a computationally efficient alternative to high-fidelity CFD. Representative studies and benchmarks report computational cost reductions varying from tens to several orders of magnitude depending on solver formulation and problem size [@katz2001low;@dimitri;@persson2012numerical]. Narrowing down to the case of a pitching and plunging airfoil, PANKH demonstrated a runtime reduction on the order $10^3$ compared to the commercial package ANSYS, when executed on identical single-core hardware. High-fidelity CFD simulations, executed on high-performance computing (HPC) clusters comprising multiple nodes and cores, often demand days of runtime due to their intensive computational requirements and parallel processing across distributed architectures. Experimental aerodynamic investigations, on the other hand, involve considerable capital investment and extended timelines. Designing wind/water tunnel setups, procuring flow visualization tools, and conducting controlled tests may take several months to years. Moreover, accommodating different flow conditions often requires substantial modifications to the experimental apparatus. To mitigate these limitations, scaling laws and dimensional analysis are commonly used but come with their own assumptions and constraints. In contrast, PANKH is a modular and open-source aerodynamic solver written in C++. It is designed to offer a rapid and accessible alternative for potential flow analysis in low-speed aerodynamic applications. Although still in its early development stage, PANKH is capable of delivering accurate estimates of aerodynamic loads within minutes. It runs efficiently even on standard single-core desktop systems. Future enhancements of the software will include parallelization via OpenMP or MPI, significantly improving performance for large-scale unsteady simulations. Additionally, the planned development of a graphical user interface (GUI) will streamline usability across platforms, including Android-based systems, transforming PANKH into a portable and interactive educational tool for both undergraduate and graduate-level aerodynamics courses.

Low-to-medium fidelity solvers like PANKH solve Laplace’s equation to model inviscid, incompressible, and irrotational flows, minimizing computational complexity while delivering precise lift estimations for aerodynamic analysis. By solving Laplace’s equation with appropriate boundary conditions, PANKH enables rapid aerodynamic analysis, making it particularly useful for preliminary design, scaling studies, and parametric investigations. Despite their advantages, open-source tools for unsteady potential flow analysis remain scarce. Established tools such as [XFOIL](https://web.mit.edu/drela/Public/web/xfoil/) [@drela1989xfoil], NASA's FoilSIM III [@nasa_foilsim], JavaFoil [@hepperle_javafoil] are tailored for steady-state aerodynamics and lack the capability to address unsteady flow phenomena. Existing research codes for unsteady vortex panel methods are often proprietary, poorly documented, or no longer maintained. In contrast, PANKH is purpose-built for the unsteady aerodynamics of hovering airfoils, offering advanced features tailored for applications such as fixed-wing aircraft, flapping-wing micro-air vehicles (MAVs), ornithopters, and hovering rotorcraft. The software repository includes validation cases comparing PANKH's results with experimental studies by [@anderson1998oscillating] and [@floryan2017scaling], as well as with numerical simulations reported by [@dimitri]. For three-dimensional wing design, understanding two-dimensional cross-sectional properties through airfoil analysis is critical for selecting optimal shapes.
Low-to-medium fidelity solvers like PANKH solve Laplace’s equation to model inviscid, incompressible, and irrotational flows, minimizing computational complexity while delivering precise lift estimations for aerodynamic analysis. By solving Laplace’s equation with appropriate boundary conditions, PANKH enables rapid aerodynamic analysis, making it particularly useful for preliminary design, scaling studies, and parametric investigations. Despite their advantages, open-source tools for unsteady potential flow analysis remain scarce. Established tools such as [XFOIL](https://web.mit.edu/drela/Public/web/xfoil/) [@drela1989xfoil], NASA's FoilSIM III [@nasa_foilsim], JavaFoil [@hepperle_javafoil] are tailored for steady-state aerodynamics and lack the capability to address unsteady flow phenomena. Existing research codes for unsteady vortex panel methods are often proprietary, poorly documented, or no longer maintained. In contrast, PANKH is purpose-built for the unsteady aerodynamics of hovering airfoils, offering advanced features tailored for applications such as fixed-wing aircraft, flapping-wing micro-air vehicles (MAVs), ornithopters, and hovering rotorcraft. The software repository includes validation cases comparing PANKH's results with experimental studies by @anderson1998oscillating and @floryan2017scaling, as well as with numerical simulations reported by @dimitri. For three-dimensional wing design, understanding two-dimensional cross-sectional properties through airfoil analysis is critical for selecting optimal shapes.

PANKH shows strong potential for real-world applications across research, education, and hobbyist domains. In the aviation industry, it can be integrated with structural solvers to perform coupled fluid–structure interaction analyses, aiding the design of flexible lifting surfaces. Its modular and procedural backend makes it suitable for classroom use, allowing students to simulate and visualize various unsteady flow scenarios with ease. PANKH can also support hobbyists and independent experimentalists working on flapping-wing vehicles, to estimate aerodynamic forces on key sections—supporting early design decisions and material selection.

Expand All @@ -61,11 +61,11 @@ The wake vorticity sheet is discretized as follows:

In the next time step, the current wake panel is convected downstream with the flow and becomes part of the known wake vorticity field, modeled as a point vortex of equivalent strength.
Consequently, the system has a total of $n+1$ unknowns, which are as follows:
- **(a)** $n$ bound vortex strengths $\gamma_j\,(1 \leq j \leq n)$ on the airfoil, and
- **(b)** $(n+1)^{\text{th}}$ unknown, $\gamma_{wp}$, representing the strength of the latest shed wake panel.
- $n$ bound vortex strengths $\gamma_j\,(1 \leq j \leq n)$ on the airfoil, and
- $(n+1)^{\text{th}}$ unknown, $\gamma_{wp}$, representing the strength of the latest shed wake panel.

Accordingly, a system of $n+1$ equations is required to determine the $n+1$ unknowns. The following section describes the formulation of these $n+1$ equations in detail.
For additional methodological details, refer to [@vezza1985method] and [@basu_1978].
For additional methodological details, refer to @vezza1985method and @basu_1978.

Each airfoil panel has a control point (or collocation point) at its midpoint, where the no-penetration boundary condition is enforced, yielding $n-1$ equations:
$$
Expand All @@ -87,14 +87,14 @@ V_{y_i}
\gamma_{j+1}
\end{bmatrix}
$$
where $[\mathbf{P}]_{ji}$ is the $2 \times 2$ *panel coefficient matrix*, describing the coefficients of the velocity induced by $j^{\text{th}}$ panel on $i^{\text{th}}$ control point. The expression for the $[\mathbf{P}]_{ji}$ is derived from [@phillips2004mechanics].
where $[\mathbf{P}]_{ji}$ is the $2 \times 2$ *panel coefficient matrix*, describing the coefficients of the velocity induced by $j^{\text{th}}$ panel on $i^{\text{th}}$ control point. The expression for the $[\mathbf{P}]_{ji}$ is derived from @phillips2004mechanics.

The velocity contribution due to the wake, $V_{wake}$, comes from the most recently shed wake panel and all previously shed vortices. This is more complex for the current wake panel, since its strength $\gamma_{wp}$, length $l_{wp}$, and orientation $\theta_{wp}$ are initially unknown. To compute its velocity contribution, a guess for length and orientation is made. A $2\times2$ Jacobian matrix is formed along with a residual from two equations (one for length, one for orientation), which is solved using Newton's method. The corrected length and orientation are then used to compute $V_{wake}$ for the current time step.
The velocity contribution due to the wake, $V_{\text{wake}}$, comes from the most recently shed wake panel and all previously shed vortices. This is more complex for the current wake panel, since its strength $\gamma_{wp}$, length $l_{wp}$, and orientation $\theta_{wp}$ are initially unknown. To compute its velocity contribution, a guess for length and orientation is made. A $2\times2$ Jacobian matrix is formed along with a residual from two equations (one for length, one for orientation), which is solved using Newton's method. The corrected length and orientation are then used to compute $V_{wake}$ for the current time step.

### Other Physical Considerations

#### Trailing Edge Condition (Kutta Condition)
Satisfying the boundary conditions alone does not yield a unique solution for the bound vortex strengths, $\gamma_j\,(1 \leq j \leq n)$. To obtain a unique solution, the flow must leave the airfoil's sharp trailing edge smoothly along the bisector line. This requirement is known as the Kutta Condition [@eldredge2019mathematical]. Mathematically, this condition is expressed as:
Satisfying the boundary conditions alone does not yield a unique solution for the bound vortex strengths, $\gamma_j\,(1 \leq j \leq n)$. To obtain a unique solution, the flow must leave the airfoil's sharp trailing edge smoothly along the bisector line. This requirement is known as the Kutta condition [@eldredge2019mathematical]. Mathematically, this condition is expressed as:
$$
\gamma_{\text{TE}}(t_k) = 0, \quad \text{where} \quad
\gamma_{\text{TE}}(t_k) = \gamma_1(t_k) + \gamma_n(t_k) + \gamma_{wp}(t_k)
Expand All @@ -116,19 +116,20 @@ Here, $\Gamma$ represents the total circulation in the flow domain enclosed by t

This equation provides the $(n+1)^{\text{th}}$ constraint, required to solve for the n+1 unknowns in the system.

The resulting system, formulated as shown in **Figure 2**, is expressed as a single matrix equation and solved using the PartialPivLU decomposition solver available in the *Dense Linear Problems and Decompositions* module of the Eigen Library [@eigenweb]. In-depth insights into post-processing steps, encompassing the solution of the unsteady Bernoulli equation and computations of aerodynamic loads, are detailed in [@chowdhury2025fmfp].
The resulting system, formulated as shown in **Figure 2**, is expressed as a single matrix equation and solved using the PartialPivLU decomposition solver available in the *Dense Linear Problems and Decompositions* module of the Eigen Library [@eigenweb]. In-depth insights into post-processing steps, encompassing the solution of the unsteady Bernoulli equation and computations of aerodynamic loads, are detailed by @chowdhury2025fmfp.

![Equations formulated to determine the set of unknowns at a time instant $t_k$.](ax=b.jpg)

### Wake Modeling

PANKH utilizes a robust time-stepping wake model to capture the intricate dynamics of vortex shedding, which is a critical aspect of unsteady aerodynamics. Shed vortices convect, interact, and induce velocities on both the airfoil and neighboring wake elements, significantly influencing the pressure distribution and resulting aerodynamic loads. While some existing methods, such as the one used in ULVPC [@prosser2011applicability], model the most recently shed vortex as a point vortex for simplicity, this approximation may limit the physical fidelity of wake dynamics. In contrast, PANKH adopts the vortex panel-based approach of [@basu_1978]. In this approach, the newly shed wake vortex is modeled as a constant-strength vortex panel whose strength, orientation, and length are not known a priori. This leads to a nonlinear system of equations, which is solved iteratively using Newtons method.
PANKH utilizes a robust time-stepping wake model to capture the intricate dynamics of vortex shedding, which is a critical aspect of unsteady aerodynamics. Shed vortices convect, interact, and induce velocities on both the airfoil and neighboring wake elements, significantly influencing the pressure distribution and resulting aerodynamic loads. While some existing methods, such as the one used in ULVPC [@prosser2011applicability], model the most recently shed vortex as a point vortex for simplicity, this approximation may limit the physical fidelity of wake dynamics. In contrast, PANKH adopts the vortex panel-based approach of @basu_1978. In this approach, the newly shed wake vortex is modeled as a constant-strength vortex panel whose strength, orientation, and length are not known a priori. This leads to a nonlinear system of equations, which is solved iteratively using Newton's method.

The entire numerical framework discussed above is implemented in the code. The implementation is structured into separate .cpp files, each containing functions dedicated to specific numerical tasks. The function and variable names are chosen to align with their respective roles in the numerical framework, ensuring clarity and ease of use for the user. A detailed guide on input parameters, motion specification, and output interpretation is provided in the project's [README](https://github.com/coding4Acause/2d_UnsteadyVortexPanel/blob/main/README.md). An extensive testing and validation suite has been used to ensure the accuracy of the results obtained from the source code.
The entire numerical framework discussed above is implemented in the code. The implementation is structured into separate `.cpp` files, each containing functions dedicated to specific numerical tasks. The function and variable names are chosen to align with their respective roles in the numerical framework, ensuring clarity and ease of use for the user. A detailed guide on input parameters, motion specification, and output interpretation is provided in the project's [README](https://github.com/coding4Acause/2d_UnsteadyVortexPanel/blob/main/README.md). An extensive testing and validation suite has been used to ensure the accuracy of the results obtained from the source code.

# Acknowledgements
This work was supported by the Defence Research and Development Organization (DRDO) under Grant No. DFTM/14/3958/IITJ/P/04/FOMS-02/1034/D(R&D)/2024.


# References