JFCMV  Vol.8 No.3 , July 2020
Numerical Simulation of Wake Deflection Control around NACA0012 Airfoil Using Active Morphing Flaps
Abstract: This study demonstrates an active flow control for deflecting a direction of wake vortex structures behind a NACA0012 airfoil using an active morphing flap. Two-dimensional direct numerical simulations are performed for flows at the chord Reynolds number of 10,000, and the vortex pattern in the controlled and noncontrolled wakes as well as the effect of an actuation frequency on the control ability are rigorously investigated. It is found that there is an optimum actuation-frequency regime at around F + = 2.00 which is normalized by the chord length and freestream velocity. The wake vortex pattern of the well-controlled case is classified as the 2P wake pattern according to the Williamson’s categorization [1] [2], where the forced oscillation frequency corresponds to the natural vortex shedding frequency without control. The present classification of wake vortex patterns and finding of the optimum frequency regime in the wake deflection control can lead to a more robust design suitable for vortex-induced-vibration (VIV) related engineering systems.

1. Introduction

When a body is placed in the vortex-dominated flows such as wakes of cylinders or airfoils at a nonzero angle of attack (AoA), the body experiences a fluctuating lift force and vibrates due to the unsteady vortex motion interacting with the body itself, which is known as a vortex-induced vibration (VIV) and regarded as an important problem ranging from a design of civil engineering structures to a flutter constraint of aircraft wing designs [3] [4] [5] [6]. The VIV typically occurs when the vortex shedding frequency locks onto the structural natural frequency, which results in the large amplitude structural vibration [7]. The prediction of VIV and lock-in phenomena has been of practical interests and thus extensively studied, e.g., so called the Griffin plot [8] [9] for predicting maximum structural displacement and the other semi-empirical wake oscillator models [4] [5] [10]. These models have been kept improved [11] [12] and effectively used in the engineering design to avoid destruction of the system due to excessive VIV phenomena. However, the use of these prediction models is often imperfect as the vortex motion is approximated by the potential flow (i.e., point vortex model) along with the linearized boundary conditions and also small velocity perturbation, which cannot account for highly unsteady and nonlinear vortex behaviors such as rolled-up wake patterns in the practical engineering problems [13]. Therefore, it is expected that another new technology can be incorporated with those VIV prediction models so that a more robust engineering design is established for avoiding a VIV-related destruction.

In the present study, we focus on the active flow control for deflecting a vortex shedding direction in the wake such that the wake-body interaction does not occur. There were many studies for classifying the wake vortex pattern behind oscillatory cylinders or pitching airfoils in the forced motion. Koochesfahani has experimentally investigated the wake vortex structures behind a pitching airfoil of the chord Reynolds number 12,000 by changing the amplitude, frequency, and shape of the oscillation waveform [13]. It reported that the wake vortex pattern as well as the mean velocity deficit is drastically changed according to the oscillation parameters and waveforms. Williamson and Roshko [1] [2] have experimentally conducted a parametric study of the cylinder oscillating in the direction perpendicular to a flow, and provided a map of wake vortex synchronization patterns in the amplitude-wavelength plane (Figure 3(a) in [1]). They categorized the wake vortex structure into different patterns, e.g., “S”, “P”, and “2P”, which stand for the single vortex, paired vortex, and formation of vortex pairs, respectively (Figure 3(b) in [1]). The vortex types mentioned above have been later widely refined in terms of the wake vortex classification [2], and also the theoretical background of the 2P vortex has been recently represented based on the point vortex theory [14]. However, the aforementioned studies focused on the wake vortex structures symmetric in the freestream direction (i.e., the wake vortices shedding in the freestream direction), and thus the control of the vortex shedding direction by forced oscillation (wake deflection control) has not been adequately investigated. Indeed, there are many other experimental and numerical studies of active flow control with respect to suppressing the flow separation around airfoils at stalled conditions or cylinders in the freestream [15] [16] [17]; nevertheless, they have mostly discussed the flow field in the vicinity of the airfoils and cylinders, whereas the wake structures have not been the main focus and thus inadequately discussed.

As such, the forced oscillation can drastically change the wake vortex structures, and thus it is expected that the vortex shedding direction can be appropriately deflected by controlling the vortex-vortex interaction in the wake, which leads to more robust engineering design for avoiding VIV phenomena. This study numerically investigates the wake deflection control for a two-dimensional NACA0012 airfoil at AoA = 5 deg. in the chord Reynolds number 10,000. In specific, the aft-portion of the airfoil is deformed in the direction perpendicular to the flow in a periodic motion, which we call the active morphing flap hereinafter. The primary objective of this study is to clarify whether the active morphing flap can effectively control and deflect the vortex shedding direction in the wake. Furthermore, the effect of actuation frequency on the control ability and associated wake-vortex structures is discussed. To the best of our knowledge, the classification of vortex patterns in such wake deflection control had never been focused with respect to the Williamson’s map [7]. For these objectives, a parametric study was performed for various actuation frequencies while keeping the same input momentum coefficient. Note that although the targeted Reynolds number is relatively low compared to the large industrial flows (such as aircraft wings of O(106 - 107)), the main focus of the present study is to investigate the modification of two-dimensional vortex structures which could be involved as large-scale flow motions in practical engineering problems. The rest of this paper is organized as follows: Section 2 specifies the problem settings; Section 3 describes the numerical methodology including a grid and validation studies; Section 4 discusses the results, and Section 5 concludes the paper.

2. Problem Specification

2.1. Flow and Airfoil Settings

In this study, the flow over a NACA0012 airfoil is considered for the wake deflection control, where the chord Reynolds number and AoA are set to be 10,000 and 5 deg., respectively. The aft-portion of the airfoil is assumed to be smoothly deformed as an active morphing flap with a length of l and amplitude of A (A is defined at the trailing edge). Figure 1 shows the schematic of the airfoil with the active morphing flap. The origin of the coordinate (x, z) is located at the leading edge. The trailing edge oscillates as z = A sin(ωt) where z = 0 represents the airfoil chord line without control and also corresponds to the original trailing edge position. ω is normalized by a uniform flow velocity U and the chord length h as F + = ωh/(2πU). F + is so-called the actuation frequency which is widely used in the relevant studies [15] [16] [17]. Let (xs(t), zs(t)) represent the coordinates of a grid point on the airfoil surface at time t. The deformation of the airfoil surface is defined as follows.

Figure 1. Schematic of the active morphing flap.

x s ( t ) = x s ( 0 ) , z s ( t ) = w 1 z s ( 0 ) + { z s ( 0 ) + w 2 A sin ( ω t ) } ( 1 w 1 ) , w 1 = 1 2 { 1 + k w tanh ( s / R 1 / 2 ) } , w 2 = 1 cos ( π 2 x s ( 0 ) x 1 x t e x 1 ) , s = ( x s ( 0 ) x t e ) 2 + ( z s ( 0 ) z t e ) 2

The parameters kw and R control a deformation gradient and set to be kw = 5.0 and R = 0.9h. (xte, zte) represents the coordinate of the original trailing edge and is defined as (xte, zte) = (h, 0). x1 indicates the location where the deformation starts, which is set to be x1 = h − l in this study. Equations above indicate that the airfoil surface deforms only vertically in time, and the weight functions of w1 and w2 provide a smooth deformation of the airfoil surface with keeping the amplitude of A in the trailing-edge movement. Note that the equations for the surface deformation above can be also applied to determine the deformation of the entire grid. The present study mainly focuses on the basic capability of the wake deflection control under forced oscillation and associated wake vortex modification. Therefore, the forced motion affecting the wake vortex structure is hoped to be as simple as possible, and thus we do not adopt a conventional pitching airfoil that may strongly oscillate leading-edge separated shear layers or a rigidly-oscillating trailing-edge flap that may cause a sharp perturbation in the flow passing over the kinked airfoil surface in the vicinity of the flap.

2.2. Examined Cases

To clarify the effect of F + on the wake deflection control, an induced flow momentum is kept constant throughout this study. Cμ is so called the momentum coefficient [15] [16] [17], which represents the ratio of the freestream and the induced momentum per unit time as follows:

C μ = ρ s a c t 2 l a c t ρ U 2 h ,

where sact is a root mean square (RMS) value of the trailing edge velocity; l a c t is the representative length scale of the morphing flap; ρ is the freestream density. sact and l a c t are defined as follows:

s a c t = A ω / 2 , l a c t = D a c t , D a c t = A l / 2 .

In this study, the momentum coefficient is fixed as 0.01, which indicates the momentum induced by the flap can be regarded as 1% of that of the freestream and is relatively small compared to the existing study [15]. Table 1 summarizes all the cases examined in this paper. According to the definition of the momentum coefficient, the oscillation amplitude is varied along with the change of the actuation frequency.

Table 1. Cases examined in this study.

3. Methodology and Validation Study

The wake deflection control is analyzed by a direct numerical simulation (DNS) of the two-dimensional (2D) compressible Navier-Stokes equations in this study. The governing equations are expressed in the body-fitted coordinate system, and the spatial derivatives of the convective and viscous terms, metrics, and Jacobian are evaluated by a sixth-order compact scheme with a tenth-order spatial filtering with a filtering coefficient 0.495 [18]. Note that the metrics and Jacobian for the coordinate transformation are evaluated by their symmetric conservative forms which can avoid the freestream preservation error on moving and deforming grids [19] [20]. The time integration is approximated by the second-order backward difference scheme, and the lower-upper symmetric alternating direction implicit and symmetric Gauss-Seidel (ADI-SGS) method is employed with five sub-iterations in each time step. The time step size is 0.0005 normalized by the sound speed of the ambient condition and the chord length, and the maximum Courant number is 1.21. The Mach number is set to be 0.2 where the compressible effect is negligible [16] [17]. At the outflow boundary, all variables are extrapolated from one point in front of the boundary. A no-slip adiabatic condition is adopted on the surface of the airfoil. A C-type grid is used around the airfoil with the outer boundary located at 25 times the chord length away from the leading edge, and the zonal grid approach [21] is adopted to increase the spatial resolution of the wake region. In Figure 2, the blue and red grids show the C-type airfoil and refined grids in the wake, respectively. The wall-normal grid spacing at the wall surface is 0.03/(Re)1/2 (or 0.12% of the chord length), which is of DNS resolution. The number of grid points is 362,103 (=1203 × 301) for the airfoil grid and 324,742 (=638 × 509) for the refined wake grid. The simulation starts from the uniform flow. When the controlled simulation is performed, the morphing flap starts oscillation after 1 flow pass over the chord length, and the simulation continues for 10 flow passes to reach quasi-steady state. Then, a time average is taken for 10 flow passes.

Next, the validation study is performed. Figure 3 shows the lift curve at the AoA between 0 and 6 deg. The red circle and blue cross represent two- and three-dimensional computations with the present solver, and the green cross represents the reference experimental data [22]. Here, the three-dimensional (3D) simulation is performed by extruding the 2D grid in the spanwise direction by 20% of the chord length. The results show that the targeted low-AoA (5 deg) flow condition can be accurately computed; furthermore, the 3D simulation also corresponds to the 2D simulation results, which indicates that the present flow is mostly dominated by 2D flow motion and thus the 2D simulation is sufficient to capture the flow physics. Note that we have also confirmed that the 3D simulation of the controlled case (F + = 2.00) results in approximately the same mean flow as the 2D simulation. Figure 4 shows the grid convergence of the surface pressure coefficient (Cp) at AoA = 5 deg. Four different grids with an ultra-coarse to fine density (9,605 to 686,485 points) were examined, which shows an excellent convergence. The 3D simulation with the extruded fine grid also corresponds to the 2D results. The fine 2D grid is hereafter used in this study.

Figure 2. Computational grids.

Figure 3. Validation of the lift coefficient: red circle and blue cross represent two- and three-dimensional computations with the present solver; green cross represents the reference experimental data [22].

Figure 4. Grid convergence of the time-averaged surface pressure coefficient (Cp).

4. Results and Discussion

First, the noncontrolled case (case 0) is discussed. In Figure 5(a), the contour color shows an instantaneous y-vorticity field, where red and blue regions represent clockwise and counter-clockwise motions, respectively. The flow is separated in the aft-portion of the suction side, and a pair of counter-rotated vortices is periodically released from the trailing edge. The present vortex pair is so-called the “P” mode (Figure 3(b) in [1]), which is fed into the downstream direction and forms the wake structure aligned in the freestream direction (AoA = 5 deg.). In Figure 5(a), one-dimensional (1D) profiles colored by black in the wake region show the RMS of the vertical (z-direction) velocity fluctuation w’ along grid lines at x/h = 1.5, 1.75, and 2.0. The 1D profiles of RMS of w’ show local peaks in the vicinity of the instantaneous wake vortices, and thus it is confirmed that the location of wake vortex structures corresponds to the local peak of RMS of w’ which moves upward by 5 deg. as going downstream. Figure 5(b) shows the power spectrum density of Cp (pressure coefficient) at (x/h, z/h) = (2.0, 0.0) with respect to a non-dimensional frequency (Strouhal number: St) based on the chord length and freestream velocity. There are clear peaks at St = 1.0 and its harmonic and higher modes. In the present Reynolds number, a single vortex generated from the separated shear layer in the suction surface (at around x/h = 0.75) directly becomes a clockwise vortex in the P mode wake. Therefore, it is expected that the shear layer instability on the suction surface strongly affects the formation of the P mode wake in addition to the wake instability. For this reason, the active morphing flap that directly oscillates the separated shear layer in the aft-portion of the airfoil is expected to have a strong effect on the wake modification.

Next, the controlled cases are discussed. Figure 6 summarizes instantaneous y-vorticity fields with control. In the wake region of each case, the black 1D lines show the RMS of w’ in the noncontrolled case (the same in Figure 5(a)), and the colored 1D lines show the controlled case. The local peaks of RMS of w’ approximately correspond to the location of wake vortex structures as discussed in the noncontrolled case, which differs depending on the actuation frequency. In the case of F + = 8.00 in Figure 6(a), small vortices are generated in the wake due to the high frequency of the flap oscillation, which does not result in a significant modification of the RMS of w’ (red line). In the moderate frequency cases (F + = 4.00 and 2.00 in Figure 6(b) and Figure 6(c)), the local peaks of RMS of w’ are moved downward in comparison with the noncontrolled case, which indicates that the direction of the shedding vortices is effectively deflected using F + = 4.00 and 2.00. Interestingly, the peak value of the RMS of w’ becomes slightly stronger in the F + = 4.00 case than F + = 2.00, and thus the F + = 2.00 case should be preferable in terms of decreasing the VIV if a body is located in the wake. It is noteworthy that vortex patterns seem to be different between these two cases although the peak locations of velocity fluctuation are similarly modified. In the low frequency cases (F + = 1.59, 1.00, 0.50), the wake deflection is not achieved. In specific, the w’-RMS of F + = 1.00 and 0.50 broadly distribute in the vertical direction, which corresponds to larger vertical motions of the wake vortices. As such, the capability of the wake deflection control strongly depends on the actuation frequency. The frequency needs to be around F + = 2.00 for the effective control, which corresponds to the shedding frequency of each single vortex in the P-mode wake without control.

(a) instantaneous y-vorticity field (b) PSD of Cp at (x/h, z/h) = (2.0, 0.0)

Figure 5. Noncontrolled case: (a) shows the instantaneous y-vorticity field with 1D profiles (black lines) of RMS of vertical (z-direction) velocity fluctuation w’; (b) shows the PSD of Cp with respect to St based on the chord length and the freestream velocity.

(a) F+ = 8.00 (case 1) (b) F+ = 4.00 (case 2) (c) F+ = 2.00 (case 3) (d) F+ = 1.59 (case 4) (e) F+ = 1.00 (case 5) (f) F+ = 0.50 (case 6)

Figure 6. Controlled cases: contour color shows the instantaneous y-vorticity field. Black and colored 1D profiles represent RMS of vertical (z-direction) velocity fluctuation w’ in the noncontrolled and controlled cases, respectively.

The previous paragraph discussed the effect of actuation frequency on the capability of the wake deflection control. Next, we focus on the vortex pattern in the noncontrolled and controlled cases of F + = 8.00, 2.00, and 1.59. Figure 7 shows the broad view of the wake vortex structures with and without control. The noncontrolled case (Figure 7(a)) shows the P mode wake as was already discussed, where the pair of counter rotating vortices is shedding in the freestream direction. In the controlled cases, the highest frequency case (F + = 8.00) does not show a clear synchronized pattern but small vortices coalesce into a large-scale structure, which is called “C” in [1]. In the most effective controlled case with F + = 2.00, the wake vortex is split into two directions. One of them comprises a series of strong vortex pairs and is horizontally shedding downstream; the other line consists of a pair of small clockwise and weak counter-clockwise vortices. This is approximately categorized as “2P” (or “P + S”), which typically appears in a fundamental lock-in flow field [1] [2] and thus leads to the significant modification of the wake direction. The low frequency case (F + = 1.59) shows the wake structure consisting of single vortex and a pair of counter-rotating vortices (P + S). Such a low frequency regime is not controllable although the clear synchronized P + S mode is observed. To summarize, in the present flow condition, the wake vortex pattern is originally the P mode, and then it is transient to P + S, 2P, and C modes as F + increases. Figure 8 schematically illustrates vortex patterns based on our observation and the sketch in the Williamson and Roshko’s map [1] [2]. Note that in Figure 8(e), the “2S” (two single vortices) mode which emerges in the typical Karman vortex street is included for reference although it is not observed in this study. Furthermore, the C mode consists of small vortices and does not show a clear synchronized pattern, and thus the sketch in Figure 8(c) is only based on the case 1 of this study and does not specify generalized structures. The more detailed classification of the C mode has been recently proposed in [23].

Finally, the PSD of Cp is discussed in Figure 9. The location of pressure sampling is set to be in the downstream of the trailing edge at (x/h, z/h) = (2.0, 0.0), where the grid is not deformed in time. In the case of F + = 8.00, the PSD of controlled case (red line) becomes broadly higher than that of the noncontrolled case (black line). The higher frequency modes (St > 10) are enhanced broadly associated with a few peaks at harmonic modes of the input frequency such as St = 8, 16, and 32. This indicates that the control using a high-frequency input results in enhancement of nonharmonic modes as well as harmonics of the input frequency, and thus no clear synchronized structure appears. On the other hand, the effectively controlled case (F + = 2.00, Figure 8(b)) shows clear harmonics of an input frequency, i.e., St = 2.0, 4.0, and 8.0 and so on. This is considered as a typical lock-in phenomenon where the natural vortex shedding frequency without control (St = 2.0 as discussed in Figure 4(b)) corresponds to the oscillation frequency of the morphing flap. The F + = 1.59 case shows peaks of PSD at the input frequency and its harmonics; however the peak at St = 1.59 does not originally exist in the noncontrolled case so that the strong lock-in phenomena unlikely occurs. As such, the well-controlled case (F + = 2.00) introduces vortex perturbation close to the natural shedding frequency, which results in the clearly synchronized wake-vortex structures, and the wake deflection control is achieved.

(a) noncontrolled (b) F + = 8.00 (case 1) (c) F + = 2.00 (case 3) (d) F + = 1.59 (case 4)

Figure 7. Vortex pattern of noncontrolled and controlled cases: contour color shows the instantaneous y-vorticity field.

(a) P (b) C (c) P + S (d) 2P (e) 2S

Figure 8. Schematics of wake vortex patterns: (a) P mode (paired vortices); (b) C mode (coalesced vortices) mode; (c) P + S mode (pair and single vortices); (d) Two P mode; (e) Two S mode.

(a) F + = 8.00 (case 1) (b) F + = 2.00 (case 3) (c) F + = 1.59 (case 4)

Figure 9. PSD of Cp at (x/h, z/h) = (2.0, 0.0): black line shows the noncontrolled case.

5. Conclusion

This paper has demonstrated the wake deflection control around a NACA0012 airfoil at the chord Reynolds number of 10,000 with an active morphing flap. A series of two-dimensional DNS was performed, and the vortex pattern in the controlled wake as well as the effect of the actuation frequency on the control ability has been rigorously investigated. It is found that there is an optimum actuation frequency regime around nondimensional frequency of F + = 2.00, which produces the 2P wake pattern based on Williamson’s categorization [1] [2]. In this case, the input frequency corresponds to the natural vortex shedding frequency, and thus the typical lock-in flow field can be obtained. The other controlled flows also exhibit synchronized vortex structures such as P + S, which, however, does not result in significant modification of the wake-vortex shedding direction. Although further investigation is required to address the mechanism of wake-deflection control based on such strong vortex-vortex interaction phenomena more in detail, the present classification of wake vortex patterns and finding of the optimum frequency regime in the wake deflection control can lead to a more robust design suitable for VIV-related engineering systems.


The present study was partly supported by JSPS KAKENHI (Grant Number K19K152050). This work was also in part supported by Council for Science, Technology and Innovation (CSTI), Cross-ministerial Strategic Innovation Promotion Program (SIP), “Materials Integration” for revolutionary design system of structural materials (Funding agency: JST). Numerical simulations of this study were performed on the Supercomputer system “AFI-NITY” at the Advanced Fluid Information Research Center, Institute of Fluid Science, Tohoku University.

Cite this paper: Abe, Y. , Konishi, T. and Okabe, T. (2020) Numerical Simulation of Wake Deflection Control around NACA0012 Airfoil Using Active Morphing Flaps. Journal of Flow Control, Measurement & Visualization, 8, 121-133. doi: 10.4236/jfcmv.2020.83007.

[1]   Williamson, C.H.K. and Roshko, A. (1988) Vortex Formation in the Wake of an Oscillating Cylinder. Journal of Fluids and Structures, 2, 355-381.

[2]   Morse, T.L. and Williamson, C.H.K. (2009) Fluid Forcing, Wake Modes, and Transitions for a Cylinder Undergoing Controlled Oscillations. Journal of Fluids and Structures, 25, 697-712.

[3]   Williamson, C.H.K. and Govardhan, R. (2004) Vortex-Induced Vibrations. Annual Review of Fluid Mechanics, 36, 413-455.

[4]   Sarpkaya, T. (1979) Vortex-Induced Oscillations—A Selective Review. Journal of Applied Mechanics. Transactions ASME, 46, 241-258.

[5]   Parkinson, G. (1989) Phenomena and Modelling of Flow-Induced Vibrations of Bluff Bodies. Progress in Aerospace Science, 26, 169-224.

[6]   Theodorsen, T. (1935) General Theory for Aerodynamic Instability and the Mechanism of Flutter. NACA TR 496.

[7]   Blevins, R.D. (1990) Flow-Induced Vibrations. Van Nostrand Reinhold, New York.

[8]   Griffin, O.M. (1980) Vortex-Excited Cross Flow Vibrations of a Single Cylindrical Tube. ASME Journal of Pressure Vessel Technology, 102, 158-166.

[9]   Khalak, A. and Williamson, C.H.K. (1999) Motions, Forces and Mode Transitions in Vortex-Induced Vibrations at Low Mass-Damping. Journal of Fluids and Structures, 13, 813-851.

[10]   Bishop, R.E.D. and Hassan, A.Y. (1963) The Lift and Drag Forces on a Circular Cylinder in a Flowing Fluid. Proceedings of the Royal Society Series A, 277, 32-50.

[11]   Klamo, J.T. (2007) Effects of Damping and Reynolds Number on Vortex-Induced Vibrations. PhD Thesis, California Institute of Technology, Pasadena.

[12]   Farshidianfar, A. and Zanganeh, H. (2010) A Modified Wake Oscillator Model for Vortex-Induced Vibration of Circular Cylinders for a Wide Range of Mass-Damping Ratio. Journal of Fluids and Structures, 26, 430-441.

[13]   Koochesfahani, M.M. (1988) Vortical Patterns in the Wake of an Oscillating Airfoil. AIAA Journal, 27, 1200-1205.

[14]   Basu, S. and Stremler, M. (2015) On the Motion of Two Point Vortex Pairs with Glide-Reflective Symmetry in a Periodic Strip. Physics of Fluids, 27, Article ID: 103603.

[15]   Corke, T.C. and Post, M.L. (2004) Separation Control Using Plasma Actuators—Stationary and Oscillating Airfoils. AIAA 004-0841.

[16]   Abe, Y., Okada, K., Nonomura, T. and Fujii, K. (2015) The Effects of Actuation Frequency on the Separation Control over an Airfoil Using a Synthetic Jet. Progress in Flight Physics, 7, 147-168.

[17]   Sato, M., Nonomura, T., Okada, K., Asada, K., Aono, H., Yakeno, A., Abe, Y. and Fujii, K. (2015) Mechanisms for Laminar Separated-Flow Control Using Dielectric-Barrier-Discharge Plasma Actuator at Low Reynolds Number. Physics of Fluids, 27, Article ID: 117101.

[18]   Lele, S.K. (1992) Compact Finite Difference Schemes with Spectral-Like Resolution. Journal of Computational Physics, 103, 16-42.

[19]   Abe, Y., Iizuka, N., Nonomura, T. and Fujii, K. (2013) Conservative Metric Evaluation for High-Order Finite Difference Schemes with the GCL Identities on Moving and Deforming Grids. Journal of Computational Physics, 232, 14-21.

[20]   Abe, Y., Iizuka, N., Nonomura, T. and Fujii, K. (2014) Geometric Interpretations and Spatial Symmetry Property of Metrics in the Conservative Form for High-Order Finite-Difference Schemes on Moving and Deforming Grids. Journal of Computational Physics, 260, 163-203.

[21]   Fujii, K. (1995) Unified Zonal Method Based on the Fortified Solution Algorithm. Journal of Computational Physics, 118, 92-108.

[22]   Ohtake, T., Nakae, Y. and Motohashi, T. (2007) Nonlinearity of the Aerodynamic Characteristics of NACA0012 Airfoil at Low Reynolds Numbers. Journal of the Japan Society for Aeronautical and Space Sciences, 55, 439-445. (In Japanese)

[23]   Hu, Z., Liu, J., Gan, L. and Xu, S. (2019) Wake Modes behind a Streamwisely Oscillating Cylinder at Constant and Ramping Frequencies. Journal of Visualization, 22, 505-527.