## 3 Computational Approach
### Computational Method  
We employ the flow solver LANS3D, which has been developed and verified to achieve high-fidelity simulations [38, 39, 40]. Generalized curvilinear coordinates \(\{\xi,\eta,\zeta\}\) were adopted to solve the governing equations. The spatial derivatives of the convective and viscous terms, metrics, and Jacobians were evaluated using a sixth-order compact difference scheme [41]. At the first and second points off the wall boundary, a second-order explicit difference scheme was adopted. Tenth-order low-pass filtering [41,42] was used with a filtering coefficient of \(0.495\). A backward second-order difference formula was used for time integration, and five sub-iterations [43] were adopted to ensure time accuracy. For time integration, the lower-upper symmetric alternating direction implicit and symmetric Gauss-Seidel (ADI-SGS) [44] methods were used. The time step size nondimensionalized by the free-stream velocity and the chord length was \(4\times 10^{-5}\). The maximum Courant number was approximately \(2.0\). The nondimensional time step based on the wall units was lower  
Figure 3: Distribution of body force in the \(x\) direction based on the Suzen–Huang model [34] and body force vectors near PA. Gray and orange areas represent the airfoil surface and exposed electrode, respectively.  
than 0.025 at the attached turbulent boundary layer. [PERSON] and [PERSON] [45] indicated that a time step of less than 0.4 is sufficient for the LES on a turbulent boundary layer. In the present study, wall-resolved LESs were conducted. The LES using the compact difference scheme with the high-order low-pass filter is well-validated and shows comparable results with standard LESs with explicit sub-grid-scale models when the computational grid resolution is fine enough [46; 47; 48]. The high-order low-pass filter adds numerical viscosity to computations at the only coarse grid region and implicitly acts as sub-grid scale models. Thus, the explicit sub-grid scale models were not used in the present study to avoid unpreferable turbulence decay. Utilized grid resolutions in the present study are explained in Section 3.4.  
All variables were extrapolated from the point inside the outflow boundary into the point at the boundary. At the outflow boundary, all variables were extrapolated from the grid points next to the boundary into the grid points at the boundary, and the static pressure was fixed as the free stream value. At the wall boundary, adiabatic and no-slip conditions were applied. For the boundaries in the spanwise direction, a periodic boundary condition was adopted. At the inflow boundary, a uniform free stream condition without disturbance was employed.  
All LESs were performed using the JAXA Supercomputer System Generation 2 and 3 (JSS2 and JSS3) at the Japan Aerospace Exploration Agency. Approximately 80 nodes (2560 cores) were used for each case.  
The computations were conducted until the aerodynamic coefficients became quasi-steady before obtaining the time-averaged flow field and aerodynamic coefficients. The computational time in nondimensional time before obtaining the data was longer than 30. The data obtaining duration was five in nondimensional time.