**Development of an Advanced Fluid-Structure-Acoustics Framework for Predicting and Controlling the Noise Emission from a Wind Turbine under Wind Shear and Yaw**

#### **Mingyue Zhou 1, Matias Sessarego 1, Hua Yang <sup>2</sup> and Wen Zhong Shen 1,\***


Received: 28 September 2020; Accepted: 26 October 2020; Published: 28 October 2020

**Abstract:** Noise generated from wind turbines is a big challenge for the wind energy industry to develop further onshore wind energy. The traditional way of reducing noise is to design low noise wind turbine airfoils and blades. A wind turbine operating under wind shear and in yaw produces periodic changes of blade loading, which intensifies the amplitude modulation (AM) of the generated noise, and thus can give more annoyance to the people living nearby. In this paper, the noise emission from a wind turbine under wind shear and yaw is modelled with an advanced fluid-structure-acoustics framework, and then controlled with a pitch control strategy. The numerical tool used in this study is the coupled Navier–Stokes/Actuator Line model EllipSys3D/AL, structure model FLEX5, and noise prediction model (Brooks, Pope and Marcolini: BPM) framework. All simulations and tests were made on the NM80 wind turbine equipped with three blades made by LM Wind Power. The coupled code was first validated against field load measurements under wind shear and yaw, and a fairly good agreement was obtained. The coupled code was then used to study the noise source control of the turbine under wind shear and yaw. Results show that in the case of a moderate wind shear with a shear exponent of 0.3, the pitch control strategy can reduce the mean noise emission about 0.4 dB and reduce slightly the modulation depth that mainly occurs in the low-frequency region.

**Keywords:** aeroacoustics; wind turbine; noise modelling; noise control

#### **1. Introduction**

As a clean and renewable energy, wind power has the advantages of low cost, low environmental pollution and wide availability. These unique advantages of wind power make it an important part of the sustainable energy mix in many countries. However, wind energy also has some drawbacks that hinder its global use. The noise caused by wind turbines has become one of the primary sources polluting the urban environment today. Wind turbine noise caused by the movement of the blades through the air is often seen as an essential aspect causing great annoyance as compared with other noise sources [1]. In Danish regulation [2], a modern wind turbine should be set up at least four times its tip height away from residential areas. Even segregated by such a long distance, a turbine still produces a sound pressure level (SPL) of more than 40 dB (A), which almost equals some common home appliances [3]. A wind turbine operating under wind shear and in yaw produces periodic changes of blade loading, which intensifies the amplitude modulation (AM) of the generated noise, and thus can give more annoyance to the people living nearby. If the aerodynamic noise of a wind turbine can be reduced without limiting its rotational speed, wind energy can be utilized with a high efficiency. Therefore, in this paper, aerodynamic characteristics and noise features of a wind turbine in wind shear and yaw is studied.

Noise from a wind turbine or wind farm often has a feature of amplitude modulation with a varying sound pressure amplitude during its operation. There are two types of AM noise: the one is created by changing the blade noise directivity to receiver [4] and the other is caused by non-uniform inflow and operation, for example the operation in wind shear and yaw. There are three methods (Japanese F-S method [5], UK method [6] and min-max method [7]) to characterize the modulation depth. Recently, Barlas et al. [8] studied the amplitude modulation caused by a wind turbine wake using an advanced Computational Fluid Dynamics (CFD) (EllipSys3D/AL) and a sound propagation model based on solving parabolic wave equation and found that the modulation depth can reach 4–5 dB even with the same strength of noise source.

The prediction formulae of trailing edge noise were summarized first using data from an experiment by Schlinker and Amiet in 1981 [9]. In 1985, Grosveld [10] used the helicopter noise formulae synthesized by Schlinker to calculate the emission noise of several two-blade, low-power downwind horizontal axis wind turbines, and the results were in good agreement with the measured data. In 1986, De Wolf continued the research based on a similar model [11]. Back to 1981, Viterna [12] utilized an approach for the low-frequency noise of a wind turbine. A rather complex airfoil self-noise prediction model (known as the Brooks, Pope and Marcolini: BPM model) using wind tunnel acoustic measurements of National Advisory Committee for Aeronautics (NACA) 0012 airfoils was published by Brooks et al. [13] in 1989. In 1996, Fuglsang and Madsen [14] used the BPM model and the inflow noise model by Lowson [15] modified from Amiet's model [16] to calculate the noise emission of a Vestas V27 wind turbine. In 2003, the National Renewable Energy Laboratory (NREL) [17] developed its NREL AirFoil Noise (NAFNoise) based on the BPM model. The BPM was improved by including the real blade geometry for predicting wind turbine noise generation by Zhu et al. from Technical University of Denmark (DTU) in 2005 [18]. In 2009, Bowdler [19] studied wind shear effects on noise at different locations and at different time instants of the day and of the year using field data collected at wind farm sites. Oerlemans [20] investigated the influence of wind shear on the amplitude modulation of wind turbine noise in 2015. Wind shear and atmospheric turbulence effects on wind turbine noise using the Monin–Obukhov similarity theory was studied in 2016 by Yuan and Cotté [21]. To predict the trailing-edge noise, a more advanced trailing-edge noise model (TNO) was developed by Parchen [22] using a relation between the sound pressure level at far field and the surface pressure spectrum at trailing edge and related the far field sound spectrum as a function of turbulent boundary layer quantities. A refined TNO model using CFD was made by Kamruzzaman et al. [23] to consider the non-isotropic issue in the trailing-edge boundary layer. Tian [24] extended Amiet's model for simulating noise from a wind turbine.

This paper deals with the development of an advanced fluid-structure-noise technique for simulating noise generated from a wind turbine and controlling the noise generation of the wind turbine in wind shear and yaw by adjusting the pitch angle of its blades. A coupled flow and acoustics code that simultaneously predicts the noise and aerodynamic outputs of a wind turbine is applied in this study.

The paper is organized as follows. Section 2 presents the related theory and methodology. Validation and control results are presented in Section 3. Finally, conclusions are drawn in Section 4.

#### **2. Theory and Methodology**

The flow past a wind turbine is modelled with the actuator line (AL) method introduced by Sørensen and Shen [25], which was implemented in DTU's in-house finite volume code EllipSys3D. The actuator line method was coupled successfully with a high-order spectral method by Kleusberg et al. [26,27]. EllipSys3D was developed by the co-operation of the Department of Mechanical Engineering at DTU [28] and the Department of Wind Energy at Risø National Laboratory [29] that is now the Department of Wind Energy at DTU. The AL method [25,30] represents the aerodynamic loads of

wind turbine blades by using a body force distributed along rotating lines. The body force imposed in the Navier–Stokes equations is computed by the aero-elastic code FLEX5 [31] also developed at DTU. The coupling of FLEX5 and EllipSys3D/AL methods is performed at every time-step. The flow data at the blades from EllipSys3D/AL are interpolated from the flow mesh and then fed into FLEX5. Using these information, the angle of attack and relative velocity including blade deformation/motions are calculated in FLEX5 and the aerodynamic loads on the blades are obtained further by using tabulated airfoil lift, drag and moment vs. angle of attack (AoA) data including dynamic stall. The new blade positions and loads from FLEX5, which includes the effect of blade bending and motion of the rotor (up/down/side-to-side/rotations) are then fed into EllipSys3D/AL as a body force on the moving actuator lines. With these new blade positions, the new loads are redistributed on the CFD mesh and a new solution is obtained in EllipSys3D/AL. In a standard aero-elastic code, these effects and motions are calculated by using the blade momentum theory (BEM) as the aerodynamic model. Here, the high-fidelity Navier–Stokes AL method was used to solve the fluid-structure interaction problem instead. To perform with large eddy simulation, the mixed scale model of Ta Phuoc [32,33] was used. To calculate the noise emission from the blades, the BPM code was used. The variables of angle of attack and relative velocity were calculated in FLEX5 and fed to BPM. A flow chart of the combined EllipSys3D/AL, FLEX5 and BPM framework is shown in Figure 1.

**Figure 1.** Flow chart of the combined EllipSys3D/AL, FLEX5 and BPM program.

In EllipSys3D, the solution to the incompressible Navier–Stokes equations is advanced in time using an iterative time-stepping method. The velocity-pressure coupling equations at each time-step are solved iteratively by sub-iterations with the usage of under-relaxation. First, the code solves the momentum equations as a predictor so that the solution can be advanced in time. The pressure correction equation, i.e., the rewritten continuity equation in satisfying the local mass flux conservation in the discretized form, is solved as a corrector making the final flow field satisfy the continuity constraint. This is a two-step procedure corresponding to a single sub-iteration, and the process is repeated until the solution becomes convergent within sub-iterations. The variables are then updated after the convergent solution is updated, and then followed by the next time step.

The aero-elastic code FLEX5 was designed to simulate a wind turbine's dynamic behavior regarding different wind conditions. FLEX5 operates in the time domain, and the output is the time-series of simulated loads and deflections. The detailed information of nine modules in FLEX5 can be found in [34].

The acoustic prediction BPM model used in this work takes the detailed geometrical and flow information into account, e.g., tip shape, blade geometry (chord and twist distributions) as well as instantaneous wind speed and direction. The velocity at each blade segment is computed by accounting

for the induction effect and the vibration velocity of the blades. Furthermore, the geometry contours of airfoil sections of blade segments are also included in the prediction as important inputs, such as the boundary thickness and displacement thickness. The boundary layer parameters for a NACA 0012 was benchmarked with the viscous-inviscid interaction code XFOIL [35] and the strong viscous-inviscid interactive coupling code Q3UIC [36], and compared with the measured values [13]. The comparisons (not shown) indicated that XFOIL under-predicts the boundary layer thicknesses about 30% and Q3UIC gave a good agreement with the measurement. In this study, we compute the boundary quantities of airfoils using Q3UIC. For more advanced calculations, the structure deformation in the blade torsion direction can also be considered, which will lead to changes in angle of attack. More details about the noise prediction calculation can be found in [18,37].

In wind turbine control systems, the proportional–integral–derivative (PID) controller is widely applied to control the generator speed and blade pitch angle. Obviously, the ideal controller should be developed by using the PID controller theory. In the present study, given the practical complexity, two control strategies, of (1) full step pitch input control and (2) interpolated pitch input control, are developed, which are P-type controls.

For the full step pitch input control, the blade pitch β is controlled at each time-step with an objective of AoA at one or more selected cross-sections, α*i*, equal to its mean value in a cycle at one or more selected cross-sections, α*i,mean*, as follows

$$\beta = \frac{1}{N} \sum\_{i=\overline{1,\ldots,N}} \alpha\_i - \alpha\_{i,\text{mean}}$$

where *N* is the number of selected cross-sections. After approximately 100 s, the solution becomes stable. Alternatively, the interpolation pitch input control can be used and Figure 2 illustrates the concepts. The idea is to control the blade pitch in the way that (i) at the peaks or troughs of one or more cross-sections, the blade pitch is controlled as the one in the full step pitch input control; (ii) at the azimuth positions in between the blade pitch is controlled by a pitch angle interpolated from the pitch angles at the peaks and troughs in (i) according to its azimuth position. Concretely, the peaks and troughs are extracted first with their corresponding azimuth angles and the mean AoAs are obtained from the maximums and minimums (red circles in Figure 2):

**Figure 2.** Pitch input setup based on angle of attack (AoA).

#### **3. Numerical Results**

In this section, the accuracy of the used numerical code is first assessed by performing a validation against field measurements performed in the DanAero project [38]. Results from the control study under the same conditions are presented afterwards.

#### *3.1. Validation of the Simulation Framework*

Before going further to discuss the performance of a wind turbine in different operation setups, the accuracy of the simulation code should be verified first. The 2.3 MW NM80 wind turbine was used here.

The design geometry of the LM 38.8 blade and aerodynamic airfoil data for the generic 2.3 MW variable speed and pitch controlled NM80 wind turbine were provided through a confidential agreement of using the DANAERO database. The airfoil data were used directly in AL/FLEX5, and the structural data originally prepared for HAWC2 were used to produce an equivalent FLEX5 input file. We assume the wind turbine operates in an ideal situation, so the tilt angle and tower shadow effects are neglected. Some key parameters of the wind turbine are listed in Table 1.


**Table 1.** Key wind turbine parameters.

Specific parameters for the two verification cases are listed below in Table 2.

**Table 2.** Description of Case I and Case II.


Figure 3 describes the configuration of the wind turbine and some parameters defined in Table 2.

**Figure 3.** Configuration of a wind turbine under wind shear and yaw.

In the computations, a relatively fine mesh of about 19 million cells was used with 192 × 192 × 512 cells in the transversal, vertical and streamwise directions in a domain of [−16R, 16R] × [0R, 19R] × [−16R, 34R] (R is the rotor radius), respectively. A time-step based on rotor radius and inflow wind speed of 0.002 was used. A resolution of 30 cells per rotor radius was used in the rotor plane and 18 cells per rotor radius were uniformly distributed from 1D (rotor diameter) in front of the turbine to 11D behind the turbine in order to have a good resolution in the wake region. The airfoil data used in the computations were 2-dimensional airfoil data provided from the International Energy Agency (IEA) Task 29 consortium. To take into account the dynamic stall effects, the Øye dynamic stall model [34] was used as this is a part of FLEX5.

The boundary condition used in the computations is that at the inlet (min z), spanwise boundary (x direction) and up-boundary (max y), a sheared wind velocity profile is used, at the ground (y = 0) no-slip condition is used, and at the outlet (max z) the convective boundary condition is used. In this study, no synthetic inflow turbulence is used.

Figures 4 and 5 illustrate the comparisons between the measured data and simulation results. The simulation results were extracted from a cycle after the flow was stabilized. The forces normal and tangential to the local airfoil chord at 33%, 48%, 76% and 92% rotor radius noted as Fn33, Ft33, Fn48, Ft48, Fn76, Ft76, Fn92, Ft92, respectively, are shown below:

**Figure 4.** *Cont*.

**Figure 4.** Forces normal and tangential to the local chord at 33%, 48%, 76% and 92% rotor radius for Case I (**a**) Fn33; (**b**) Ft33; (**c**) Fn48; (**d**) Ft48; (**e**) Fn76; (**f**) Ft76; (**g**) Fn92; (**h**) Ft92.

**Figure 5.** *Cont*.

**Figure 5.** Forces normal and tangential to the local chord at 33%, 48%, 76% and 92% rotor radius for Case II (**a**) Fn33; (**b**) Ft33; (**c**) Fn48; (**d**) Ft48; (**e**) Fn76; (**f**) Ft76; (**g**) Fn92; (**h**) Ft92.

As observed, the magnitudes between the measurements and simulations are similar for most of the quantities. In Case I, the computed normal force does not follow with the measured one and this is probably due to the Øye dynamic stall model used in the computations, which does not perform well at small angles of attack. On the other hand, the tangential force agrees well with the measurements. In Case II, the normal force is seen to be very well captured, while the tangential force is slightly over-predicted. The trends in function of azimuth angle (with 0◦ when the blade is down-pointing) followed the numerical results well. Since the tangential force is sensitive to the local angle of attack and airfoil data, calibrating airfoil data and its 3-dimensional corrections to rotational effects is a challenge for accurately predicting the performance of a wind turbine under wind shear and yaw. From the comparisons, it is also seen that the agreement is fairly good for both small (6◦) and large (38◦) yaw angles under a moderate wind shear. Therefore, it can be concluded from the validation that the simulation technique used in the present paper can predict the wind turbine performance relatively well.

#### *3.2. Power Performance and Load Simulation*

Since a wind turbine is rarely working with a large yaw angle, we define scenarios for the control study according to Case I in the validation. The specifics of the four different scenarios are described below:


Figure 6 illustrates the streamwise velocity and pressure field for Case 1 and Case 4. Horizontal axis x and vertical axis y are non-dimensional with the rotor radius of R = 40 m, which means the real values of them should equal to the axis values times R. The rotor center is located at x = 1.5 and y = 1.5. In Figure 6a, a velocity slowdown (induction effect) created by the rotor when exacting wind energy can be observed in the areas both in front of and behind the rotor, where the impact of induction effect is larger in the wake. When no shear is included, no wind speed stratification above the ground can be observed. A significant wind speed stratification above the ground can be observed in Figure 6b, while the wake is also mitigated by the lower longitudinal wind speed due to the shear. Similarly, Figure 6c,d illustrate the pressure field for Case 1 and Case 4, respectively. In general, the pressure in front of the rotor is larger than the pressure behind the rotor. By introducing the wind shear and yaw, the pressure both in front of and behind the rotor is changed, especially for the pressure in front of the rotor, which becomes more homogeneous.

**Figure 6.** Streamwise velocity field and pressure field for Case 1 and Case 4: (**a**) streamwise velocity for Case 1; (**b**) streamwise velocity for Case 2; (**c**) pressure for Case 1; (**d**) pressure for Case 2.

The angle of attack (AoA) at 90% length of the blade for the four scenarios with no control is investigated. The comparison of the AoA can be found in Figure 7. In general, it can be observed that the fluctuations of AoA become much larger by introducing yaw or shear. In this study, the 0.3 wind shear power law exponent leads to a more significant amplitude change than the 10-degree yaw angle.

**Figure 7.** Angle of attack in a cross-section of 35.24 m (90%R) without control strategy.

Figure 8 illustrates the comparisons of AoA under pitch control with respect to different scenarios. It is noted that the pitch angle control in this study is performed using the AoA at 90%R where the wind turbine noise source center is located [37]. In the beginning, due to inexact interpolation values, the interpolated pitch control strategy results in more throbbing of AoA than the full step control strategy. However, when it is converged, both control techniques lead to a significant reduction of AoA fluctuations in all the four scenarios, while the mean of AoA remains the same.

**Figure 8.** Angle of attack in section 35.24 m (90%R) with control strategy: (**a**) Case 1; (**b**) Case 2; (**c**) Case 3; (**d**) Case 4.

Figure 9 shows the pitch input values based on the results of AoA mentioned above. Pitch angles for two control strategies are calculated by using the angle of attack at 90%R of the no control case. It can be observed that the pitch is constant with 0.15◦ pitch offset when no control is performed. In terms of pitch in all the four cases, no significant difference can be found between control strategies I and II.

**Figure 9.** *Cont.*

**Figure 9.** Pitch in a cross-section of 35.24 m (90%R) with and without control strategy: (**a**) Case 1; (**b**) Case 2; (**c**) Case 3; (**d**) Case 4.

In Figure 10, power and thrust coefficients, noted as *Cp* and *CT*, respectively, are plotted for the four scenarios. It can be seen that *Cp* is similar in all the four scenarios whereas introducing the shear and yaw slightly reduces *Cp* especially when *Cp* becoming more stationary, i.e., after 110 s. The reduction of the mean *Cp* is mainly due to the yaw misalignment, where the axial wind speed is smaller. It also should be noted that the scenarios with a yaw angle, *Cp* has more periodical fluctuations due to the yaw angle. A similar trend like *Cp* can be found in the thrust coefficient *CT*, shown in Figure 10b, where the reduction of *CT* is realized by the shear and yaw. The most significant difference happened in the scenario with both yaw and wind shear. Moreover, it can be observed that the difference between the yaw versus the yaw plus shear scenarios is relatively small. The reduction of *CT* in Case 2, Case 3 and Case 4 is due to the reduction of the axial wind speed caused by the yaw.

**Figure 10.** *Cp* and *CT* without control strategy: (**a**) power coefficient; (**b**) thrust coefficient.

By implementing the pitch adjustment, more fluctuations in Figure 11 are observed in the *Cp* and *CT* signals while the mean *Cp* and *CT* almost remain the same. By introducing the control methods, the fluctuation magnitudes of both control strategies increase, especially in control technique I. Some jumps of *Cp* and *CT* values occur in the transient period as well. Moreover, the transient periods when performing both control techniques are longer, as compared to no control. For control technique I, the unexpected large fluctuations in *CT* and *Cp* might be explained by the large pitch changes at each time step. For control technique II, similar results are observed when the fluctuation magnitudes are much smaller. The possible explanation is that there are fewer input pitch values and the corresponding interpolation enables better convergence of the elastic code.

**Figure 11.** Power and thrust coefficients with a pitch adjustment: (**a**) *Cp* for Case 1; (**b**) *CT* for Case 1; (**c**) *Cp* for Case 2; (**d**) *CT* for Case 2; (**e**) *Cp* for Case 3; (**f**) *CT* for Case 3; (**g**) *Cp* for Case 4; (**h**) *CT* for Case 4.

#### *3.3. Acoustical Analysis*

Effects of Wind Shear and Yaw on Sound Emission

Figure 12 illustrates the polar plot of the sound pressure level in the function of azimuth angle of blade 1 for the four cases, which is called the amplitude modulation (AM) noise. AM noise is often a big issue because it gives more perception to the people living nearby. However, it should be noted that there are two different phenomena that make AM noise: one is caused by the blade noise directivity change during its operation [4] and another is caused by the non-uniform flow and operation, such as wind shear and yaw. In the present study, the latter is focused. The observer is located at 1.5 m height and a distance of 500 m downstream from the turbine (500 m away from the turbine is often a place where people are living although the sound propagation effects due to ground and atmospheric conditions are not considered here). The choice of the 500 m distance downstream location is to minimize the AM noise created by the blade noise directivity change. The speed of sound used in the calculations is 340 m/s. The turbulence level and turbulence length scale are set as 3% and 100 m, respectively. This choice was made in accordance with the 0.3 wind shear exponent defined in the scenarios. A more realistic choice can be made when the detailed flow information is measured. The 0◦ azimuth angle is defined when blade 1 is vertically down and coincides with the tower. All these four scenarios are performed with no control. As can be observed in Cases 1 and 3, the sound pressure level is slightly larger when the shear is included, while it significantly decreases with a 10◦ yaw, i.e., for both Cases 2 and 4. Normally, the SPL from one blade is highest when it is pointing upwards. The reason can be seen from the BPM sound prediction equations that the SPL is related to the fifth power of the oncoming velocity. For the case of three blades, at zero azimuth angle (blade 1 is pointing downwards), the other two blades are pointing at 120◦ and 240◦ when the noise is almost highest during its operation for Cases 3 and 4.

**Figure 12.** Sound pressure level in function of blade azimuth angle without control at 1.5 m height and 500 m distance downstream from the turbine at a wind speed of 6 m/s.

Figure 13 illustrates the comparison of SPL regarding different control scenarios. It can be observed that compared with the non-controlled case, the controlled SPL results, i.e., with control techniques I and II, are similar for Case 1 and Case 2. For Case 1, the peaks for both control strategies are in step, whereas the magnitude of the no control results is slightly smaller. For Case 2, the fluctuations are also generally in step, and both the time and magnitude of the peaks are consistent. Therefore, the implementation of the control techniques may not lead to noise reductions for the no-shear scenarios. The main reason for the negligible difference in Cases 1 and 2 is that the modifications of the pitch and the corresponding AoA are relatively small where the fluctuations of noise level are mainly due to blade deflections. A larger difference is expected for a larger pitch. It should be noted that Case 2 corresponds to a slightly lower SPL as compared to Case 1, due to the yaw-caused lower axial wind speed. For Cases 3 and 4, a larger decrease of SPL, i.e., about 0.4 dB, can be observed by introducing either control strategy I or II, whereas the trend remains the same. No significant difference can be found between the results of the two control techniques. It can be concluded that by introducing the control techniques, the reduction of SPL can be realized for the shear scenarios. In comparison to the angle of attack plots in Figure 8 with those plots of SPL in Figure 13, a non-linear effect between these two quantities is seen due to the blade structural responses. The AoA variation of the interpolation control in Figure 8d is reduced about 1◦. According to the rule of a 1◦ AoA reduction resulting in a 1 dB noise reduction [39], the noise change in the present study is within this range. For comparison, it can be observed that the differences between the controlled and uncontrolled results of A-weighted SPL regarding all the cases are negligible. The reason is that the

A-weighted SPL accounts for more of the high-frequency noise to which humans are more sensitive and de-emphasizes the low-frequency contribution. In our case, the reduction of the noise mainly occurs in the low-frequency part and thus the A-weighting SPL is similar in all the cases. It should be noted that the low-frequency noise is important in far-field as it can propagate in a long distance with very small air absorption. To check the AM noise, the modulation depth for Case 4 is calculated according to the method proposed in [7] and shown in Table 3. From the table, it is seen that the modulation depth is slightly reduced.

**Figure 13.** Sound pressure level and A-weighted sound pressure level of the turbine at 6 m/s wind speed: (**a**) sound pressure level (SPL) for Case 1; (**b**) A-weighted SPL for Case 1; (**c**) SPL for Case 2; (**d**) A-weighted SPL for Case 2; (**e**) SPL for Case 3; (**f**) A-weighted SPL for Case 3; (**g**) SPL for Case 4; (**h**) A-weighted SPL for Case 4.


**Table 3.** Modulation depth for Case 4.

Moreover, this phenomenon can be verified in Figure 14, where the contribution of SPL in different frequencies at the time instant (120 s) is shown. From the figure, it is seen that at the low frequency (which is less than 315 Hz), SPL is reduced obviously by the pitch control strategy.

**Figure 14.** SPL spectra of case 4 of the turbine at 6 m/s wind speed at a time instant of 120 s after control.

Figure 15 illustrates the SPL and A-weighted SPL in Cases 1, 2 and 4, but with a lager inflow wind speed of 10 m/s. With the lager wind speed, the SPL increases about 8.5 dB in each case. The interesting thing deserved to be noticed is that the pitch control methods at 10 m/s display a similar effect in Case 4 as at 6 m/s but the variations in amplitude are more important. Since the A-weighted SPL was designed according to human hearing and emphasizes the contribution from the high frequency sound centered at 1000 Hz, it can give an impression that the control methods give a small effect on A-weighted SPL. The modulation depth is listed in Table 3 and a small reduction is seen. It is worth noting that the yaw and shear considered in this study are relatively small and these effects can be much bigger in large shear and yaw cases. Moreover, the noise propagation effects due to changing atmospheres and wakes [8] are not included here. In the future, the present noise source-modelling framework will be further coupled with the propagation model [8] in order to control wind turbine noise at far field.

**Figure 15.** Sound pressure level and A-weighted Sound pressure level from the wind turbine at 10 m/s wind speed: (**a**) SPL for Case 1; (**b**) A-weighted SPL for Case 1; (**c**) SPL for Case 2; (**d**) A-weighted SPL for Case 2; (**e**) SPL for Case 4; (**f**) A-weighted SPL for Case 4.

#### **4. Conclusions**

In this paper, an advanced flow-structure-acoustics framework has been developed and validated against load field measurements for the case of the 2.3 MW NM80 wind turbine under a moderate wind shear and a small and a large yaw angle. The code has been further used for studying and controlling noise generation from a wind turbine under a moderate wind shear and a small yaw, which is often the case when a wind turbine operates. Through analysis of the results, several conclusions are summarized here:

1. Simulations based on the coupled Ellipsys3D/AL/Flex5 framework agreed well with field load measurements on the 2.3 MW NM80 turbine.


**Author Contributions:** Conceptualization, W.Z.S.; Formal analysis, M.Z.; Investigation, M.Z.; Methodology, W.Z.S.; Project administration, W.Z.S.; Resources, W.Z.S.; Software, M.S. and W.Z.S.; Supervision, M.S., H.Y. and W.Z.S.; Visualization, M.Z.; Writing—original draft, M.Z.; Writing—review & editing, M.S., H.Y. and W.Z.S. All authors have read and agreed to the published version of the manuscript.

**Funding:** This paper was partially funded by the Danish Energy Agency in the project (EUDP2018, Journal nr. 64018-0084) and the Innovation Fund Denmark in the project entitled "DecoWind—Development of low-noise and cost-effective wind farm control technology" under project number (8055-00041B). The DanAero projects were funded partly by the Danish Energy Agency, (EFP2007. Journal nr. 33033-0074 and EUDP 2009-II. Journal nr: 64009-0258) and partly by self-funding from the project partners Vestas, Siemens Gamesa Renewable Energy, LM Wind Power, Dong Energy and DTU Wind Energy.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Abbreviations**


#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Performance Characteristics of an Orthopter-Type Vertical Axis Wind Turbine in Shear Flows**

**Rudi Purwo Wijayanto 1,2,\*, Takaaki Kono 3,\* and Takahiro Kiwata <sup>3</sup>**


Received: 9 January 2020; Accepted: 2 March 2020; Published: 4 March 2020

**Abstract:** To properly conduct a micro-siting of an orthopter-type vertical axis wind turbine (O-VAWT) in the built environment, this study investigated the effects of horizontal shear flow on the power performance characteristics of an O-VAWT by performing wind tunnel experiments and computational fluid dynamics (CFD) simulations. A uniform flow and two types of shear flow (advancing side faster shear flow (ASF-SF) and retreating side faster shear flow (RSF-SF)) were employed as the approaching flow to the O-VAWT. The ASF-SF had a higher velocity on the advancing side of the rotor. The RSF-SF had a higher velocity on the retreating side of the rotor. For each type of shear flow, three shear strengths (Γ = 0.28, 0.40 and 0.51) were set. In the ASF-SF cases, the power coefficients (*Cp*) were significantly higher than the uniform flow case at all tip speed ratios (λ) and increased with Γ. In the RSF-SF cases, *CP* increased with Γ. However, when Γ = 0.28, the *CP* was lower than the uniform flow case at all λ. When Γ = 0.51, the *CP* was higher than the uniform flow case except at low λ; however, it was lower than the ASF-SF case with Γ = 0.28. The causes of the features of *CP* were discussed through the analysis of the variation of blade torque coefficient, its rotor-revolution component and its blade-rotation component with azimuthal angle by using the CFD results for flow fields (i.e., horizontal velocity vectors, pressure and vorticity). These results indicate that a location where ASF-SFs with high Γ values dominantly occur is ideal for installing the O-VAWT.

**Keywords:** orthopter; vertical axis wind turbine; power coefficient; torque coefficient; shear flow; wind tunnel; CFD; delayed detached-eddy simulation

#### **1. Introduction**

Since the 2000s, interest in installing small wind turbines (SWTs) in the built environment has been growing [1–9]. Wind conditions in the built environment are complex in nature and are characterized by lower wind speeds and higher turbulence because of the presence of obstructions [8,9]. For SWTs to be able to make up their costs within their lifetimes, they should have high efficiency and be placed at sites with high wind speeds, such as coastal sites or high-elevation inland sites. However, in the built environment, keeping the rotational speed of an SWT's rotor as low as possible is preferable from the viewpoint of aerodynamic noise [10,11]. Therefore, the optimal tip-speed ratio of an SWT in the built environment should be as low as possible, while the maximum power coefficient of the SWT should be as high as possible.

Wind turbines are classified into horizontal-axis wind turbines (HAWTs) and vertical-axis wind turbines (VAWTs), based on the orientation of their rotation axes. Generally, in the built environment, VAWTs are preferable to HAWTs because VAWTs do not suffer, as much as HAWTs, from reduced energy outputs from frequent wind direction changes [12]. Wind turbines are further classified into lift-type wind turbines and drag-type wind turbines, based on the aerodynamic force component that

acts on a blade and dominantly contributes to the rotor rotation. With regard to lift-type VAWTs, a lot of research on Darrieus-type VAWTs, including straight-bladed and helical-bladed ones, has been conducted [13–16]. With regard to drag-type VAWTs, a lot of research on Savonius-type VAWTs has been conducted [17–21]. In general, the optimal tip-speed ratio of a drag-type VAWT is less than 1.0 [11], which is much smaller than that of a lift-type VAWT. Moreover, although the maximum power coefficient of a drag-type VAWT is generally much smaller than that of a lift-type VAWT, the power coefficient of a drag-type VAWT is generally greater than that of a lift-type VAWT at a low tip-speed ratio, of less than 1.0. Therefore, a drag-type VAWT is favorable in the built environment and was researched by our research group.

Our group [22,23] researched a drag-type VAWT called the orthopter-type VAWT (O-VAWT). The O-VAWT is a variable-pitch VAWT; each of the flat-plate blades not only revolves around the main shaft but also rotates around its own blade axis, which is rotationally supported by a pair of connecting arms. We investigated the effects of the number and aspect ratio of the flat-plate blades on the power performance of the O-VAWT in a uniform flow by conducting wind tunnel experiments with an open test section and three-dimensional computational fluid dynamics (CFD) simulations. When the number of the blades was three and the aspect ratio of the blades was 1:1, the maximum power coefficient was 0.25 and the optimal tip-speed ratio was 0.4 [22]. Here, the optimal tip-speed ratios of Savonius-type VAWTs are in the range of 0.45 to 1.0 [21]. Except for several studies that were conducted using a wind tunnel with a very high blocking ratio of the closed test section, the maximum power coefficient of the Savonius-type VAWT was, at most, 0.25 [21]. That is, the O-VAWT has a lower optimal tip-speed ratio than Savonius-type VAWTs, although the maximum power coefficient is relatively high. Therefore, from the viewpoint of aerodynamic noise, the O-VAWT can be more favorable in the built environment as compared to Savonius-type VAWTs. Except for our studies, studies on the power performance of O-VAWTs are very limited. Shimizu et al. [24] investigated the effects of the aspect ratio of the blade on the power performance of an O-VAWT with two blades whose cross-sectional shape was an ellipse by conducting wind tunnel experiments. Bayeul-Line et al. [25] examined the effects of the blade's cross-sectional shape (elliptical and straight) and the initial blade stagger angle on the performance of an O-VAWT by conducting 2-dimensional CFD simulations. Cooper and Kennedy [26] examined the power performance of an O-VAWT with three blades whose cross-sectional shape was the upstream half of a NACA0010-65 section reflected about the mid-chord by conducting theoretical analysis with a multiple-stream tube model and field measurements. Our group [23] conducted wind tunnel experiments to compare the performance of an O-VAWT with elliptic blades and one with flat-plate blades. By considering the mechanical loss torque, we obtained the maximum power coefficient of 0.246 at a tip-speed ratio of 0.4 for the O-VAWT with elliptic blades and of 0.288 at a tip-speed ratio of 0.4 for the O-VAWT with flat-plate blades. It should be noted that these studies on O-VAWTs were conducted in conditions where the approaching flows had uniform distribution.

To properly conduct a micro-sitting of an O-VAWT in the built environment, it is important to understand the effects of the strong shear approach flow with on the performance of the O-VAWT. Figure 1 illustrates the approaching wind flow to a building. As the wind flow approaches the building, the wind speed decreases, and the pressure increases. Then, the wind flow proceeds along the upwind face of the building and separates at the corners on the roof and the side walls. As the separated wind flow is not obstructed by the building, the pressure decreases, and the wind speed increases. Near the upwind corners, the wind speed increases more than that of the approaching wind. In addition, reverse flow regions are formed between the separated shear layer and the building's walls. As a result, strong shear flows are formed vertically over the roof surface and horizontally over the side walls. Due to the mixing of momentum, the shear becomes weaker as the flow proceeds downstream. To utilize the increased wind speed over the roof of a building, the effects of building shapes and wind directions on the wind conditions have been investigated (e.g., [27]). Furthermore, the effects of wind conditions, such as wind speed, turbulent intensity, and skew angles, on the potential energy yield and the power performance of a wind turbine have been studied (e.g., [28,29]). In this study, we

investigated the effects of horizontal shear flow on the performance of the O-VAWT by conducting wind tunnel experiments and three-dimensional CFD simulations. A uniform flow and two types of shear flow were employed as the approaching flow to the O-VAWT. One type had a higher velocity on the advancing side of the rotor. The other type had a higher velocity on the retreating side of the rotor. For each type of shear flow, we set three different shear strengths.

**Figure 1.** Illustration of the approaching wind flow to the wind turbine near upwind corners of a building. (**a**) Bird view; (**b**) enlarged top view.

#### **2. Experimental Approach**

#### *2.1. Wind Turbine Model*

In this paper, we employed a right-handed Cartesian coordinate system (*x*1, *x*2, *x*3) = (*x*, *y*, *z*), in which the *z*-direction was aligned in the vertical direction. The wind turbine used in this study was an O-VAWT with three flat-plate blades as shown in Figure 2. The blade had a height of *<sup>h</sup>* <sup>=</sup> 4.00 <sup>×</sup> <sup>10</sup>−<sup>1</sup> m, chord length of *<sup>c</sup>* <sup>=</sup> 4.00 <sup>×</sup> <sup>10</sup>−<sup>1</sup> m and thickness of 4.0 <sup>×</sup> <sup>10</sup>−<sup>3</sup> m. Each of the blades not only revolved around the main shaft but also rotated around its own blade axis, which was rotationally supported by a pair of connecting arms. The distance between the main shaft and one of the blade axes was *R* = 2.55 <sup>×</sup> 10−<sup>1</sup> m. In addition, each of the blade axes was connected with the main shaft by a chain via sprockets. Since the ratio of the number of teeth on the sprocket of the main shaft to that of the blade axis was 1:2, each of the blades rotated around the own blade axis a half time while the rotor revolved around the main shaft one time. When seen from the top, the rotor revolved around the main shaft counterclockwise and each of the blades rotated around the own blade axis clockwise as shown in Figure 2c. Therefore, by using the angular velocity of the rotor revolution (ω), the angular velocity of the blade rotation is expressed as −ω/2. Furthermore, as shown in Figure 2c, according to

the azimuthal angle of a blade (ϕ), we call the range of 90◦ < ϕ < 270◦ the "upwind region" of the rotor, 0◦ < ϕ < 90◦ and 270◦ < ϕ < 360◦ the "downwind region" of the rotor, 180◦ < ϕ < 360◦ the "advancing side of the rotor" and 0◦ < ϕ < 180◦ the "retreating side" of the rotor. The O-VAWT is designed so that the drag force on a blade is large on the advancing side of the rotor while being small on the retreating side of the rotor.

**Figure 2.** Orthopter-type vertical-axis wind turbine (O-VAWT). (**a**) A photograph of O-VAWT with three flat blades, main shaft, arm, the chains and sprockets; (**b**) an isometric view of O-VAWT with three flat blades; (**c**) motion of rotor and blades viewed from the top; (**d**) a projected swept area of the rotor viewed from the upwind side.

Figure 2d shows a projected swept area of the rotor (*A*), which is defined as:

$$A = \left(2R + 0.5\,\text{c}\right)h.\tag{1}$$

We define the diameter of the rotor as:

$$D = 2R.\tag{2}$$

The O-VAWT had a rotor's diameter of *D* = 5.1 <sup>×</sup> 10−<sup>1</sup> m and a projected rotor s swept area of *A* = 2.84 <sup>×</sup> 10−<sup>1</sup> m<sup>2</sup> can be considered as a micro wind turbine. A small scale wind turbine that has a diameter up to 1.25 m and the swept area up to 1.2 m<sup>2</sup> is categorized as a micro wind turbine [4].

#### *2.2. Experimental Setup for Uniform Flow Case*

Figure 3a shows the experimental setup for the uniform flow case. The experiments were conducted using a closed circuit wind tunnel with an open test section. The size of the cross section of the wind tunnel outlet was 1.25 m × 1.25 m. The blockage ratio which is defined as the ratio of the projected rotor s swept area to the wind tunnel outlet area was approximately 18%. The O-VAWT was set in the test section so that the rotor center was at the center of the cross-section of the wind tunnel outlet and 0.850 m downwind of the wind tunnel outlet. Here, we defined the rotor center as the point on the rotational axis of the rotor and at the mid-height of the blades. In addition, we set the origin of the coordinate system at the rotor center, as shown in Figure 2c,d. The rotor was driven by a motor (Mitsubishi Electric, GM-S) and its rotational speed (ω) was monitored by using a digital tachometer (Ono Sokki, HT-5500) and controlled by using an inverter (Hitachi, SJ200). The rotor torque was measured by using a torque meter (TEAC, TQ-AR), which was connected to the motor and shaft via couplings. The output signal of the torque meter was converted by a 16-bit analog-to-digital converter with a sampling interval of 0.5◦, and 36,000 items (50 revolutions) of data were stored. To measure the reference wind speed *U*∞, an ultrasonic anemometer (Kaijo Sonic, DA-650-3TH and TR-90 AH) was set approximately 2 m upwind of the wind tunnel outlet. The value of *U*<sup>∞</sup> was kept at 8 m/s. The value of tip speed ratio λ, which is defined as:

$$
\lambda = \text{RayOL}\_{\text{os}\_{\lambda}} \tag{3}
$$

was varied from 0.1 to 0.8 with an increment of 0.1.

#### *2.3. Experimental Setup for Shear Flow Cases*

To generate a horizontal shear flow, a perforated panel and a splitter plate were installed at the outlet of the wind tunnel as shown in Figure 3b,c and Figure 4. The perforated panel had a width of 1.0 m, a height of 1.5 m and a thickness of 2 <sup>×</sup> 10−<sup>3</sup> m, and was set so that it covered half of the wind tunnel outlet in the horizontal wind direction. Due to the existence of the perforated panel, the pressure upwind of the panel increased and the wind flow rate through the wind-tunnel-outlet area covered by the panel decreased while that through the uncovered wind-tunnel-outlet area increased. When the perforated panel covered the wind-tunnel-outlet area upwind of the retreating side of the rotor, the wind speed of the generated shear flow was higher on the advancing side of the rotor. Hereafter, this type of shear flow is referred to as "advancing side faster shear flow" (ASF-SF). On the other hand, when the perforated panel covered the wind-tunnel-outlet area upwind of the advancing side of the rotor, the wind speed of the generated shear flow was higher on the retreating side of the rotor. Hereafter, this type of shear flow is referred to as "retreating side faster shear flow" (RSF-SF). The splitter plate was set vertically, parallel to the wind tunnel wall and at the center of the wind tunnel outlet to avoid the horizontal component of the wind velocity in the generated shear flow becoming significant. The splitter plate had a width of 0.88 m, a height of 1.25 m and a thickness of 4 <sup>×</sup> <sup>10</sup>−<sup>3</sup> m.

(**a**)

(**c**)

**Figure 3.** The experimental apparatus and measurement devices; (**a**) in uniform flow, (**b**) in shear flows and (**c**) the porous plate position at the nozzle exit of the wind tunnel in case of shear flows.

**Figure 4.** The perforated panel and splitter plate set at the outlet of the wind tunnel.

To investigate the effects of the strength of the shear flow on the performance of the O-VAWT, we generated three kinds of shear flows by using three perforated panels shown in Table 1. With regard to a staggered round-hole perforated panel, the shielding ratio Φ, which is the ratio of the area that shields the airflow to the whole area of the perforated panel, can be computed by:

$$\Phi = 1 - \frac{\pi \, d^2}{2 \sqrt{3} \, L^{2'}} \tag{4}$$

where *d* is a diameter of a hole and *L* is the distance between the centers of adjacent holes.


**Table 1.** Perforated panels. Here, *d* is the diameter of a hole, *L* is the distance between the centers of adjacent holes and Φ is the shielding ratio.

Except for the installation of the perforated panels and the splitter plate, the experimental setup for the measurement of the performance of the O-VAWT was the same as the uniform flow case. Prior to the measurement of the performance of the O-VAWT, we measured the horizontal profiles of the generated shear flows at 0.10 m downwind of and at the center height of the wind tunnel outlet by using an x-type hot-wire probe (Kanomax, 0252R-T5).

#### *2.4. Torque and Power Coe*ffi*cients*

Due to the difficulty of evaluating the mechanical losses of the bearings, the sprockets and the chains, this study considers only the aerodynamic torque generated by the blades as the rotor torque of the O-VAWT. The aerodynamic torque generated by the blades was computed by:

$$T\_B = T\_{\text{uvB}} - T\_{\text{uvB}\_{\text{V}}} \tag{5}$$

where *TwB* was the measured aerodynamic torque generated by the rotor when the blades were not removed; and *TwoB* was the measured aerodynamic torque generated by the rotor when the blades were removed. It is worth noting that, generally, when the blades were not removed from the rotor, the rotor generated positive torque while the motor acted as a load to keep the value of ω constant. Conversely, at all tip speed ratios, when the blades were removed from the rotor, the rotor generated negative torque while the motor acted as the driving force of the rotor revolution. Therefore, at all tip speed ratios, *TB* was higher than *TwB*.

The power coefficient describes that fraction of the power in the wind that may be converted by the turbine into mechanical work [30] and is defined in this study as:

$$\mathcal{C}p = \frac{T\_{B}\omega}{0.5\rho A l l\_{0}^{3}}\tag{6}$$

and the torque coefficient is defined as:

$$C\_T = \frac{T\_B}{0.5 \,\rho A \,\mathrm{ul}\_0^2 \mathrm{R}}.\tag{7}$$

Here, the *U*<sup>0</sup> is the time-mean stream-wise velocity, *u*(*x*, *y*, *z*), averaged over the projected rotor's swept area at *x* = 0 and is computed by:

$$
\Omega I\_0 = \frac{\int\_{-\left(R+0.5c\right)}^{R} \overline{u}(0, y, 0) dy}{2R + 0.5c}. \tag{8}
$$

#### **3. Numerical Approach**

The CFD software utilized to simulate the wind flow field was ANSYS Fluent 17.2 [31,32]. The numerical approach was based on our previous paper [22,33], in which the CFD simulations with the delayed detached eddy simulation (DDES) turbulence model of flow around the O-VAWT were conducted and the validities of the grid resolution and the time-step size were confirmed.

#### *3.1. Governing Equations and Discretization Method*

The flow field around the wind turbine was assumed to be incompressible and isothermal. The DDES turbulence mode treats near-wall region in a manner like a Reynolds-averaged Navier–Stokes (RANS) turbulence model and treats the rest of the flow field in a manner like a large-eddy simulation (LES) turbulence model [34]. This model has the potential to achieve higher accuracy than RANS models and save a large number of computing resources compared with pure LES models. The governing equations for the CFD simulation with the DDES turbulence model based on the Spalart–Allmaras (SA) model are the continuity equation:

*Appl. Sci.* **2020**, *10*, 1778

$$\frac{\partial \mathbf{u}\_i}{\partial \mathbf{x}\_i} = \mathbf{0},\tag{9}$$

the Navier–Stokes equation:

$$\frac{\partial u\_i}{\partial t} + \frac{\partial u\_i u\_j}{\partial \mathbf{x}\_j} = -\frac{1}{\rho} \frac{\partial p}{\partial \mathbf{x}\_i} + \nu \frac{\partial}{\partial \mathbf{x}\_j} \left( \frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i} - \frac{2}{3} \delta\_{ij} \frac{\partial u\_i}{\partial \mathbf{x}\_j} \right) - \frac{\partial}{\partial \mathbf{x}\_j} \left[ \hat{\mathbf{v}} \left( \frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i} \right) - \frac{2}{3} k \delta\_{ij} \right],\tag{10}$$

and the transport equation for the kinematic eddy viscosity <sup>∼</sup> ν:

$$\frac{\partial \check{\boldsymbol{\nu}}}{\partial t} + \frac{\partial \check{\boldsymbol{\nu}} \boldsymbol{u}\_{l}}{\partial \boldsymbol{x}\_{l}} = \mathbb{C}\_{\mathbb{H}1} \boldsymbol{\hat{\nu}} \left( \mathbb{S} + \frac{\check{\boldsymbol{\nu}}}{\kappa^{2} d\_{\mathrm{DDES}}} \Big( \mathbbm{1} - \frac{\boldsymbol{\chi}}{1 + \boldsymbol{\chi} f\_{l1}} \Big) \right) + \frac{1}{\check{\boldsymbol{\nu}} \boldsymbol{\hat{\nu}}} \Big[ \frac{\partial}{\partial \boldsymbol{x}\_{l}} \Big\{ (\boldsymbol{\nu} + \boldsymbol{\hat{\nu}}) \frac{\partial \check{\boldsymbol{\nu}}}{\partial \boldsymbol{x}\_{l}} \Big\} + \mathbb{C}\_{\mathbb{H}2} \Big( \frac{\partial \check{\boldsymbol{\nu}}}{\partial \boldsymbol{x}\_{l}} \Big)^{2} \right] - \mathbb{C}\_{\mathbb{u}1} f\_{\mathbb{U}} \left( \frac{\check{\boldsymbol{\nu}}}{d\_{\mathrm{DRES}}} \right)^{2} \tag{11}$$

where *ui* is the wind-velocity component in the *xi* direction; *p* is the pressure; ν is the kinematic viscosity; *t* is the time; ρ is the air density; *k* is the turbulence kinetic energy; δ*ij* is the Kronecker delta; *dDDES* is the DDES length scale; χ is (<sup>∼</sup> ν/ν); *S* is a scalar measure of the deformation tensor; *f*ν<sup>1</sup> and *fw* are damping functions; and *Cb*1, *Cb*2, *Cw*1, σ<sup>∼</sup> <sup>ν</sup> and κ are constants.

The DDES length scale is computed by:

$$d\_{\rm DDES} = d - f\_{\rm d} \max\left(0, d - c\_{\rm des} \Delta\_{\rm max}\right) \tag{12}$$

where *d* is the distance to the closest wall; *cdes* is the empirical constant; and Δ*max* is the maximum edge length of the local computational cell, i.e., Δ*max* = *max (*Δ*x,* Δ*y,* Δ*z)*. The switching between the RANS and the LES mode depends on the following shielding function;

$$f\_d = 1 - \tanh\left(\left(8r\_d\right)^3\right) \tag{13}$$

and

$$r\_d = \frac{\tilde{\nu}}{\sqrt{\mu\_{i,j}\mu\_{i,j}\kappa^2 d^2}}\tag{14}$$

where *ui*,*<sup>j</sup>* is the velocity gradient. The damping functions and closure coefficients are as follows:

$$f\_{\nu1} = \frac{\chi^3}{c\_{\nu1}^3 + \chi^{3'}}, \; f\_w = g \left[ \frac{1 + c\_{w3}^6}{g^6 + c\_{w3}^6} \right]^{1/b}, \; g = r + c\_{w2} (r^6 - r), \; c\_{w1} = \frac{c\_{b1}}{\chi^2} + \frac{(1 + c\_{b2})}{c\_{\nu}^2}, \tag{15}$$

$$\mathbf{c}\_{b1} = 0.1355, \; \mathbf{c}\_{b2} = 0.622, \; \mathbf{c}\_{\nu1} = 7.1, \; \mathbf{c}\_{\mathbf{u}2} = 0.3, \; \mathbf{c}\_{\mathbf{u}3} = 2.0, \; \sigma\_{\mathbf{v}} = \frac{2}{3}, \; \mathbf{c}\_{\text{des}} = 0.65, \; \mathbf{x} = 0.4187. \tag{16}$$

The governing equations are discretized by the finite-volume method. The advection terms of the Navier–Stokes equations are discretized by the bounded central-difference scheme. The advection term of the transportation equations for <sup>∼</sup> ν is discretized by a second-order upwind scheme. Other spatial derivatives are discretized by the central-difference scheme. The time integration is performed using the second-order implicit method.

#### *3.2. Numerical Setup*

The computational domain, the computational meshes and the boundary conditions are shown in Figure 5. As the same as the experimental setup, the origin of the coordinate system is defined at the center of the O-VAWT. The modeled O-VAWT was comprised of three blades, one main shaft and two sets of connecting arms. To reduce the computational cost, other components, such as the chains and sprockets are omitted. The sizes of these components are the same as those used in the experiment. The computational domain consists of three blade domains, one rotor domain and one far-field domain. The blade domain included one of the blades and rotates around each blade axis. The rotor domain

included these three blade domains, the main shaft and the connecting arms and rotates around the main shaft. The far-field domain was a stationary domain and its size was 23.5*D* × 17*D* × 5*D*. Except for uniform flow cases, a splitter plate with the same thickness as the experiment was set at *y* = 0 and *x* = −11.37*D* to −1.67*D*. In all domains, only unstructured meshes were used. Based on our previous mesh-resolution dependency tests [33], the number of computational cells in each domain were set as shown in Table 2. The total number of computational cells was approximately 10 million. All surfaces of the solid components were covered with boundary-layer meshes. The first grid nodes over the surface of the blades were *y*+ < 1 in all run cases.

**Figure 5.** Computational domain, computational mesh and boundary conditions. (**a**) Top view of the computational domain; (**b**) bird's eye view of the computational domain with computational mesh and boundary conditions; (**c**) modeled rotor; (**d**) blade domains and (**e**) rotor domain.



At the inlet boundary, the distributions of the stream-wise wind velocity shown in Table 3 were implemented. These distributions were set so that when the O-VAWT is absent, the distributions of the time-mean values of *u* at *x* ≈ −0.147 *D*, which corresponds to 0.1 m downwind of the wind tunnel outlet, matched well with those of the wind tunnel experiment, which were generated by using the three kinds of perforated panels with different shielding ratios Φ, as shown in Figure 6a,b. Here, *UH* and *UL* are the maximum and minimum streamwise velocities in the profile of a shear flow measured in the wind tunnel experiment; Γ is the velocity ratio defined as:

$$
\Gamma = \frac{\mathcal{U}\_H - \mathcal{U}\_L}{\mathcal{U}\_H + \mathcal{U}\_L} \tag{17}
$$

**Table 3.** Distribution of *u* at the inlet boundary. Here, Γ is the velocity ratio; Φ is the shielding ratio of a perforated panel; *UH* and *UL* are the maximum and minimum velocities in a shear flow.


**Figure 6.** Horizontal distribution of the time-mean values of streamwise velocity at the mid-height of the wind tunnel outlet at *x* ≈ −0.147*D*, which corresponds to 0.1 m downwind of the wind tunnel outlet, and *x* = 0, which corresponds to the position of the rotational axis of the O-VAWT, when the O-VAWT is absent. (**a**) ASF-SF cases at *x* ≈ −0.147*D*, (**b**) RSF-SF cases at *x* ≈ −0.147*D*, (**c**) ASF-SF cases at *x* = 0 and (**d**) RSF-SF cases *x* = 0.

It is worth noting that amplification factors of Γ = 0.28 and Γ = 0.51 are considered as 1.35 and 1.65, respectively, by computing *UH*/*U*∞. The amplification factor is defined as the ratio of wind speed in the case where there are buildings to wind speed in the case where these buildings are removed. The separated shear flow from a building with an aspect ratio of 1:1:2 reaches an amplification factor of 1.2 [35]. In addition, the separated shear flow from a building with an aspect ratio of 1:1:6 reaches an amplification factor of 1.7 [36]. At *x* = 0, the values of Γ do no change, as shown in Figure 6c,d. However, due to the momentum diffusion, the horizontal gradients of the time-mean streamwise velocity become weaker. One of the reasons for the relatively large discrepancies between the profiles obtained by the experiment and the CFDs can be that the turbulence intensities of the CFDs are significantly small compared to those of the experiments, as shown in Figure 7. Even though very high values of turbulence intensity were set at the inlet boundary, the turbulence intensities dissipated rapidly and became very small at *x* = 0 as compared to those of the experiments. As the setting of high turbulence intensity at the inlet boundary leads to computational instability, we set no perturbation condition at the inlet boundary. Table 4 shows the values of *U*<sup>0</sup> computed by Equation (8) and the Reynolds number which is defined as:

$$Re = \frac{\mathcal{U}\_0 D}{\nu} \tag{18}$$

**Figure 7.** Horizontal distribution of the turbulence intensity of ASF-SF cases at the mid-height of the wind tunnel outlet at *x* = 0, which corresponds to the position of the rotational axis of the O-VAWT, when the O-VAWT was absent. (**a**) Linear scale on the horizontal axis and (**b**) logarithmic scale in the horizontal axis.


**Table 4.** Experimental and computational fluid dynamics (CFD) results for the values of *U*<sup>0</sup> and *Re* for the uniform flow, ASF-SF and RSF-SF cases.

At the outlet boundary, the pressure outlet condition with *p* = 0 was imposed. On the surface of the O-VAWT and the splitter plate, the no-slip boundary conditions were set. The sliding mesh technique was used to couple the rotational domains and the stationary domain. The direction of the rotor and the blade rotations are counterclockwise and clockwise, respectively, when viewed from the top (Figure 5a). By changing the rotational speed of the rotor ω, the tip speed ratio λ was set at 0.2, 0.4, 0.5, 0.6, or 0.8. The time step sizes were set as *dt* = 0.5◦/ω.

#### *3.3. Torque and Power Coe*ffi*cients*

As mentioned in sub-Section 2.4, this study considers only the aerodynamic torque generated by the blades as the rotor torque of the O-VAWT. Since each of the blade axes was connected with the main shaft by a chain via sprockets, the aerodynamic torque on each of the blades about each of the blade axes was transmitted through the chain and contributed to the torque about the main shaft. Therefore, the rotor torque generated by a blade at an azimuthal angle ϕ is expressed as:

$$T\_B(\varphi) = T\_{B\\_rev}(\varphi) + T\_{B\\_rot}(\varphi),\tag{19}$$

where *TB\_rev*(ϕ) is the conventional blade torque that is calculated by multiplying the rotor radius and the component of the aerodynamic force on the blade at ϕ in the rotor-revolution direction; and *TB\_rot*(ϕ) is the torque generated by the component of the aerodynamic force on the blade at ϕ in the blade-rotation direction about the blade axis. Hereafter, we call *TB* (ϕ) the "blade torque," *TB\_rev*(ϕ) the "rotor-revolution torque" and *TB\_rot*(ϕ) the "blade-rotation torque." The blade torque coefficient (*CTB*), the rotor-revolution component (*CTB\_rev*) and the blade-rotation component (*CTB\_rot*) are defined as:

$$\mathcal{C}\_{TB}(\varphi) = \frac{T\_B \ (\varphi)}{0.5 \rho A l I\_0^2 R} \,\mathrm{'}\tag{20}$$

$$\mathcal{C}\_{TB\\_rev}(\boldsymbol{\varphi}) = \frac{T\_{B\\_rev}\left(\boldsymbol{\varphi}\right)}{0.5\rho A \mathcal{U}\_0^2 \mathcal{R}}\,\mathrm{}^\prime\tag{21}$$

and

$$C\_{TB\\_rot}(\varphi) = \frac{T\_{B\\_rot}\left(\varphi\right)}{0.5\rho A \,\mathrm{l}\,\mathrm{l}\_0^2 R}.\tag{22}$$

It should be noted that *CTB*, *CTB\_rev* and *CTB\_rot* are coefficients of one blade.

The CFD simulations were conducted for eight revolutions of the rotor. Using the data of the last two rotor revolutions, the torque coefficient *CT* was computed by the following formula:

$$\mathbf{C}\_{T} = \frac{n}{N\_{\varepsilon} - N\_{s}} \sum\_{N=N\_{s}+1}^{N\_{\varepsilon}} \int\_{0}^{2\pi} \mathbf{C}\_{TB} \, d\boldsymbol{q} \,. \tag{23}$$

Here, *n* (= 3) is the number of the blades; *Ns* (= 6) is the number of the rotor revolutions before starting the computation of *CT*; *Ne* (= 8) is the number of the rotor revolutions before finishing the computation of *CT*. The power coefficient of *CP* was computed by Equation (6).

#### **4. Results and Discussion**

In this section, the results of the wind tunnel experiments and the CFD simulations for the power performance of the O-VAWT, such as the dependency of the power and torque coefficients on the tip speed ratio, the variations of the torque coefficients with azimuthal angle, are presented for the uniform flow case and the shear flow cases. Subsequently, the causes of the features of the power performance of the O-VAWT are discussed based on the CFD results of the flow fields.

#### *4.1. Performance in Uniform Flow*

Figure 8a shows the power and torque coefficients (*CP* and *CT*) of the O-VAWT in the uniform flow. The CFD results are in good agreement with the experimental ones. The optimal tip speed ratio at which *CP* becomes the maximum is less than unity; 0.4 in the experiments and 0.5 in the CFD simulations. As mentioned in the introduction, this low optimal tip speed ratio is a favorable feature for the built environment from the viewpoint of aerodynamic noise. With increasing λ, *CT* decreases monotonically from a small tip speed ratio (λ = 0.2). These tendencies are commonly found among drag-type wind turbines. The wind-tunnel experimental results by Shimizu et al. [24] for *CP* and *CT* of an O-VAWT with two elliptical cross-sectional blades show the same tendencies as our results. The value of the maximum *CP* is 0.32 in our experiments and CFD simulations, while the value is 0.176 in Shimizu et al.'s experiments. The main factors for the better performance of our O-VAWT as compared to Shimizu et al.'s O-VAWT can be the number of blades and the cross-sectional shape of the blades. In our previous studies, the maximum Cp improved from 0.189 to 0.244 by changing the number of blades from two to three [33] and improved from 0.246 to 0.288 by changing the cross-sectional shape of blades from ellipse to rectangle [23].

**Figure 8.** Performance of the O-VAWT in the uniform flow; (**a**) the variation of power coefficient (*CP*) and torque coefficient (*CT*) with tip speed ratios (λ) by the experiments and CFD simulations, (**b**) rotor-revolution (*CT\_rev*) and blade-rotation (*CT\_rot*) components of torque coefficients (*CT*) computed by the CFD simulations. The wind-tunnel experimental results by Shimizu et al. [24] for *CP* and *CT* of an O-VAWT with two elliptical cross-sectional blades are added for reference.

Figure 8b shows *CT\_rev* and *CT\_rot* computed by the CFD results. The sum of *CT\_rev* and *CT\_rot* is *CT*. As well as *CT*, the value of *CT\_rev* decreases monotonically with an increase in λ*.* Conversely, the value of *CT\_rot* increases with an increase in λ. The values of *CT\_rev* are positive and larger than those of *CT\_rot* except for λ = 0.8. At λ = 0.8, *CT\_rev* is negative; however, *CT\_rot* is positive and its absolute value is larger than that of *CT\_rev*. As a result, *CT* is positive at λ = 0.8.

Figure 9 shows the variations of blade torque coefficient (*CTB*), its rotor-revolution component (*CTB\_rev*) and blade-rotation component (*CTB\_rot*) with respect to azimuthal angle (ϕ) at λ = 0.4 and 0.6. The value of *CTB* is significantly large in the upwind region of the advancing side (ϕ ≈ 180◦ to 270◦) of the rotor, being the maximum at ϕ ≈ 210◦. In the range of ϕ where *CTB* is significantly large, the contribution of *CTB\_rev* is dominant. Except for this range, *CTB\_rev* does not always positively contribute to *CTB*. At ϕ where the value of *CTB\_rev* is negative, the value of *CTB\_rot* is generally positive and the rotation of the blade positively contributes to *CTB*. Due to this positive contribution of *CTB\_rot* to *CTB*, the value of *CTB* is positive at almost all ϕ and the variation of *CTB* of the O-VAWT with respect to ϕ is smaller as compared to that of a Savonius-type VAWT. The variation of *CTB* of a Savonius-type VAWT with respect to ϕ in Figure 9 is a result of CFD simulation by Tian et al. [37]. The maximum *CP* and the optimal λ of the Savonius-type VAWT were 0.258 and 1.0, respectively.

**Figure 9.** The variation of blade torque coefficients (*CTB*), its rotor-revolution component (*CTB\_rev*) and its blade-rotation component (*CTB\_rot*) with azimuthal angle (ϕ) at; (**a**) λ = 0.4 and (**b**) λ = 0.6. Note that *CTB*, *CTB\_rev* and *CTB\_rot* are coefficients of one blade. The CFD simulation result by Tian et al. [37] for the variation of *CTB* of a Savonius-type VAWT with ϕ is added for comparison.

#### *4.2. Performance in Shear Flows*

Figure 10 compares the experimental results for *CP* of the O-VAWT in the cases of the advancing side faster shear flow (ASF-SF), the retreating side faster shear flow (RSF-SF) and the uniform flow. In the cases of the ASF-SF, *CP* is higher than in the case of the uniform flow at all λ and increases with an increase in Γ. In the cases of the RSF-SF, similar to the cases of the ASF-SF, *CP* increases with an increase in Γ. However, when Γ = 0.28, the values of *CP* are lower than those of the uniform flow case. Concerning the optimal λ, there is a trend that it shifts to higher λ in both cases of shear flows as compared to the uniform flow case, except for the RSF-SF with Γ = 0.28.

**Figure 10.** Power coefficient (*CP*) of O-VAWT in shear flows by the experiments.

Both in the cases of the ASF-SF (Figure 11a) and the RSF-SF (Figure 11b), the CFD results for *Cp*–λ curves are in good agreement with the experimental ones. In the following discussion, we use the CFD results for the torque of the O-VAWT and the flow field to explain the effects of shear flows on the characteristics of the power performance.

**Figure 11.** Power coefficient (*CP*) of O-VAWT in shear flows by the experiments and the CFD simulations in the cases of; (**a**) the ASF-SF and (**b**) the RSF-SF.

Figure 12 shows the variations of blade torque coefficient (*CTB*), its rotor-revolution component (*CTB\_rev*) and its blade-rotation component (*CTB\_rot*) with azimuthal angle (ϕ) in the cases of the ASF-SF. It is confirmed that these profiles are qualitatively the same between the cases of λ = 0.4 and λ = 0.6. Therefore, it is considered that the effects of the shear flow on the characteristics of the blade torque variations with ϕ do not significantly change around the optimal λ. As compared to the case of the uniform flow, *CTB* is higher on most of the advancing side of the rotor (ϕ ≈ 210◦ to 330◦) and lower in the upwind region of the retreating side of the rotor, specifically at ϕ ≈ 120◦ to 150◦. In particular, with an increase in Γ, *CTB* increases on most of the advancing side of the rotor (ϕ ≈ 210◦ to 330◦). The optimal ϕ at which *CTB* is the maximum shifts to the downwind direction on the advancing side of the rotor (ϕ ≈ 240◦) as compared to the case of the uniform flow. The effects of the shear flow on the variations of *CTB\_rev* with ϕ is almost the same as *CTB*. As compared to the case of the uniform flow, *CTB\_rev* is higher on most of the advancing side of the rotor (ϕ ≈ 210◦ to 330◦). In contrast, *CTB\_rot* is lower in most of the upwind region of the advancing side (ϕ ≈ 180◦ to 240◦) and slightly higher in the upwind region of the retreating side of the rotor, specifically at 120◦ to 150◦ as compared to the case of the uniform flow. With an increase in Γ, *CTB\_rot* decreases in the upwind region of the advancing side of the rotor (ϕ ≈ 180◦ to 240◦). Due to its negative values of *CTB\_rot*, the optimal ϕ at which *CTB* is the maximum slightly shifts to the downwind direction as compared to *CTB\_rev*.

Figure 13 shows the variation of blade torque coefficient (*CTB*), its rotor-revolution component (*CTB\_rev*) and its blade-rotation component (*CTB\_rot*) with azimuthal angle (ϕ) in the cases of the RSF-SF. Since these profiles are qualitatively the same between the cases of λ = 0.4 and λ = 0.6, it is considered that the effects of the shear flow on the characteristics of the blade torque variations with ϕ do not significantly change around the optimal λ. As compared to the case of the uniform flow, *CTB* is significantly higher in most of the upwind region of the retreating side (ϕ ≈ 120◦ to 180◦) and lower on most of the advancing side of the rotor (ϕ ≈ 210◦ to 330◦). In particular, with an increase in Γ, *CTB* increases in the upwind region of the retreating side of the rotor. The optimal ϕ at which *CTB* is the maximum shifts to the upwind region of the retreating side of the rotor (ϕ ≈ 150◦) as compared to the case of the uniform flow. The effects of the shear flow on the variation of *CTB\_rev* with ϕ is almost the same as *CTB*. As compared to the case of the uniform flow, *CTB\_rev* is higher in most of the upwind region of the retreating side (ϕ ≈ 100◦ to 160◦) and lower on most of the advancing side of the rotor (ϕ ≈ 210◦ to 330◦). In contrast, *CTB\_rot* is slightly lower in most of the upwind region of the retreating side (ϕ ≈ 90◦ to 150◦) and slightly higher in most of the upwind region of the advancing side of the

rotor (ϕ ≈ 180◦ to 240◦) as compared to the case of the uniform flow. Furthermore, with an increase in Γ, *CTB\_rot* decreases in most of the upwind region of the retreating side (ϕ ≈ 100◦ to 160◦) and increases in most of the upwind region of the advancing side of the rotor (ϕ ≈ 180◦ to 240◦).

**Figure 12.** The variation of blade torque coefficient (*CTB*), its rotor-revolution component (*CTB\_rev*) and its blade-rotation component (*CTB\_rot*) in the cases of the ASF-SF computed by the CFD simulations; (**a**) for λ = 0.4 and (**b**) for λ = 0.6. Note that *CTB*, *CTB\_rev* and *CTB\_rot* are coefficients of one blade.

**Figure 13.** The variation of blade torque coefficient (*CTB*), its rotor-revolution component (*CTB\_rev*) and its blade-rotation component (*CTB\_rot*) in the cases of the RSF-SF computed by the CFD simulations; (**a**) for λ = 0.4 and (**b**) for λ = 0.6. Noted that *CTB*, *CTB\_rev* and *CTB\_rot* are coefficients of one blade.

#### *4.3. Flow Characteristics*

Figures 14a, 15a and 16a show the temporal sequence of the horizontal distributions of normalized horizontal velocity vectors, normalized pressure and normalized vorticity, respectively, at the mid-height of the O-VAWT at λ = 0.4 in the case of the uniform flow. At ϕ = 180◦ to 270◦, the approaching flow to the blade has a large velocity component perpendicular to the blade, and the pressure on the upwind side of the blade is high. Due to this high pressure, *CTB\_rev* is significantly high in the range of ϕ = 180◦ to 270◦ in Figure 9. In addition, at ϕ = 210◦ to 270◦, due to the strong large vortex formed near the outer edge of the downwind side of the blade, the pressure is low near the vortex. This low pressure positively contributes to *CTB\_rev* (see Figure 9 at ϕ = 210◦ to 270◦). By contrast, this low pressure negatively contributes to *CTB\_rot* (see Figure 9 at ϕ = 210◦ to 270◦). At ϕ = 300◦ and 330◦, the approaching flow to the blade has a small velocity component perpendicular to the blade and the pressure on the upwind side of the blade is not high. Therefore, *CTB\_rev* is lower as compared to the upwind region of the advancing side of the rotor (see Figure 9 at ϕ = 180◦ to 270◦). At ϕ = 0◦ to 120◦, the attack angle of the blade is positive (here, the counterclockwise direction is defined as positive) and the flow separates over the outer side of the blade. Therefore, on the upwind edge of the blade and on the inner side of the blade near its upwind edge, the pressure is relatively high. In contrast, on the outer side of the blade near its upwind edge, the pressure is relatively lower. This pressure distribution contributes negatively to *CTB\_rev* and positively to *CTB\_rot*.

**Figure 14.** The temporal sequence of the horizontal distributions of normalized horizontal velocity vectors of the O-VAWT for λ = 0.4 in the case of; (**a**) the uniform flow, (**b**) the ASF-SF with Γ = 0.51 and (**c**) the RSF-SF with Γ = 0.51.

**Figure 15.** The temporal sequence of normalized pressure at the mid-height of the O-VAWT for λ = 0.4 in the case of; (**a**) the uniform flow, (**b**) the ASF-SF with Γ = 0.51 and (**c**) the RSF-SF with Γ = 0.51.

Figures 14b and 15b and Figure 16b show the temporal sequence of the horizontal distributions of normalized horizontal velocity vectors, normalized pressure and normalized vorticity, respectively, at the mid-height of the O-VAWT at λ = 0.4 in the case of the ASF-SF with Γ = 0.51. At ϕ = 210◦ to 330◦, the approaching flow to the blade has a larger velocity component perpendicular to the blade, and the pressure on the upwind side of the blade is higher as compared to the uniform flow case. Due to this higher pressure, *CTB\_rev* is higher than the uniform flow case (see Figure 12 at ϕ = 210◦ to 330◦). In addition, at ϕ = 210◦ to 240◦, wind speed is significantly increased at the outer edge of the blade and a significantly stronger and larger vortex is formed near the outer edge on the downwind side of the blade. Near the vortex, the pressure is lower as compared to the uniform flow case. This low pressure contributes to higher *CTB\_rev* and lower *CTB\_rot* as compared to the uniform flow case (see Figure 12 at ϕ = 210◦ to 240◦). At ϕ = 120◦ and 150◦, the pressure on the outer side of the blade is lower due to the lower speed approaching flow to the blade. Furthermore, the pressure on the inner side of the blade near its upwind edge is higher due to the existence of the blade at ϕ = 240◦ or 270◦, respectively. This pressure distribution contributes to lower *CTB\_rev* and higher *CTB\_rot* as compared to the uniform flow case (see Figure 12 at ϕ = 120◦ to 150◦).

**Figure 16.** The temporal sequence of normalized vorticity at the mid-height of the O-VAWT for λ = 0.4 in the case of; (**a**) the uniform flow, (**b**) the ASF-SF with Γ = 0.51 and (**c**) the RSF-SF with Γ = 0.51.

Figures 14c, 15c and 16c show the temporal sequence of normalized horizontal distributions of horizontal velocity vectors, pressure and vorticity, respectively, at the mid-height of the O-VAWT at λ = 0.4 in the case of the RSF-SF with Γ = 0.51. At ϕ = 120◦ and 150◦, the approaching flow to the blade has a larger velocity component perpendicular to the blade and the pressure on the outer side of the blade is higher as compared to the uniform flow case. Due to this higher pressure, *CTB\_rev* is higher than the uniform flow case (see Figure 13 at ϕ = 120◦ and 150◦). Furthermore, wind speed is increased at the upwind edge of the blade and a stronger vortex is formed near the upwind edge of the blade on its inner side. Near the vortex, the pressure is lower as compared to the uniform flow case. This low pressure contributes to higher *CTB\_rev* and lower *CTB\_rot* as compared to the uniform flow case (see Figure 13 at ϕ = 120◦ and 150◦). At ϕ = 210◦ to 300◦, except in the vicinity of the upwind side of the inner edge of the blade at ϕ = 210◦, the approaching flow to the blade has a smaller velocity component perpendicular to the blade and the pressure on the upwind side of the blade is lower as compared to the uniform flow case. Due to this lower pressure, *CTB\_rev* is lower than the uniform flow case (see Figure 13 at ϕ = 210◦ and 300◦). In addition, at ϕ = 210◦ and 240◦, the vortex formed near the outer edge of the downwind side of the blade is weaker and the pressure drop becomes smaller as

compared to the uniform flow case. This smaller pressure drop contributes to lower *CTB\_rev* and higher *CTB\_rot* as compared to the uniform flow (see Figure 13 at ϕ = 210◦ and 240◦).

#### **5. Conclusions**

We investigated the effects of horizontal shear flow on the performance characteristics of an orthopter-type vertical axis wind turbine (O-VAWT) by conducting wind tunnel experiments and computational fluid dynamics (CFD) simulations. In addition to a uniform flow, two types of shear flow were used as the approaching flow to the O-VAWT. One type was an advancing side faster shear flow (ASF-SF), which had a higher velocity on the advancing side of the rotor. The other type was a retreating side faster shear flow (RSF-SF), which had a higher velocity on the retreating side of the rotor. For each type of shear flow, we set three different velocity ratios (Γ = 0.28, 0.40 and 0.51), which were the ratios of the difference between the highest velocity and the lowest velocity in a shear flow to the sum of the highest and lowest velocities. The main findings are summarized as follows:


These findings are useful for micro-siting of an O-VAWT in the area where shear flows occur. A location where ASF-SFs with high Γ values dominantly occur is ideal for installing the O-VAWT. At a location where not only ASF-SFs but also RSF-SFs occur at high frequencies, higher Γ values are preferable. However, the shear flows utilized in this study are limited in their profiles and the relative positions to the rotor. To properly conduct the micro-siting of an O-VAWT in the area where various kinds of shear flows occur, such as the vicinity of a building, it is essential to understand the performance characteristics of the O-VAWT in the various kinds of shear flows. Therefore, in future research, we plan to investigate the effects of the broadness of the shear layer and the relative position

of the shear flow to the rotor on the O-VAWT's performance characteristics. In addition, we plan to investigate the effects of the turbulence intensity of the approaching flow on the CFD results for the O-VAWT's performance characteristics by setting obstacles, which emit eddies, upwind of the O-VAWT to avoid the rapid dissipation of high turbulence intensity.

**Author Contributions:** R.P.W. performed the numerical simulations and prepared this manuscript being supervised by T.K. (Takaaki Kono). All authors contributed to the analyses of the data. T.K. (Takaaki Kono) and T.K. (Takahiro Kiwata) supervised the entire work. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by World Bank Loan No. 8245-ID and a grant from Takahashi Industrial and Economic Research Foundation.

**Acknowledgments:** This research was supported by the Program for Research and Innovation in Science and Technology (RISET-Pro) Kemenristekdikti. The authors are thankful to technician Kuratani and student Taiki Sugawara for their help with the experiment.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Applied Sciences* Editorial Office E-mail: applsci@mdpi.com www.mdpi.com/journal/applsci

MDPI St. Alban-Anlage 66 4052 Basel Switzerland

Tel: +41 61 683 77 34 Fax: +41 61 302 89 18

www.mdpi.com