Abstract
Over recent years, considerable attention has been devoted to the optimization of energy production in wind farms, where yaw angles can play a significant role. In order to quantify and maximize such potential power, the simulation of wakes is vital. In the present study, an actuator line model code was implemented in the OpenFOAM flow solver. A tip treatment was applied to involve the tip effect induced by the pressure equalization from the suction and pressure sides. The Leishman–Beddoes dynamic stall (LB-DS) model modified by Sheng et al. was employed to consider the dynamic stall phenomenon. The developed ALM-CFD solver was validated for the NREL Phase VI wind turbine reference case. The solver was then used in simulating the yawed wind turbine, and power variation was compared with UBEM and CFD. Overall, according to the obtained data, the coupled solver compared well with CFD. There was an improvement in terms of prediction of the phase delay that is due to the dynamic stall. However, there was still negligible overestimation in deep stall conditions. Based on the obtained results, it is suggested that the reduction of power output follows a cosine to the power of X function of the yaw angle. In terms of visualizing wake, the results demonstrated that the current ALM code was satisfying enough to simulate skewed wake and vortices trajectory. The effect of advancing and retreating blade was captured. It was found that yaw led to the concentration of the induced velocity downstream, resulting in a lower velocity deficit on a broader area, which is essential for wind farm optimization.
1. Introduction
The efforts towards reducing the environmental impact caused by fossil fuel-based electricity production led to the definition of more sustainable and green power generation systems, among which wind turbines play one of the most relevant roles. During the last decades, wind has gone from an alternative to a mainstream energy source, which is currently driving the renewables evolution. The amount of electricity production covered by wind points out relevant shares in many countries, providing 15% of EU demand in 2019 [1] and showing a world average of 4.8% in 2018 [2]. With an average growth of installed power at 13.4% in the last ten years and the levelized cost of energy (LCOE) dropping for both onshore and offshore applications [1, 3, 4], wind turbines gather strong research interest. However, unlike the vertical axis wind turbines (VAWT) that are wind direction independent [5, 6], the effect of yaw angle in horizontal axis wind turbines (HAWT) is still an important issue.
In this context, the definition of numerical models to simulate wind turbines in different operating conditions is of utmost importance. Focusing on the aerodynamic analysis of the turbine wake, the actuator line method (ALM) is widely adopted, representing a valid compromise between low-order models and 3D calculations of the actual geometry. The method introduces a body-force scheme within a 2D or 3D flow solver. In ALM, the blades are replaced by rotating lines, along which a force distribution is imposed. This is determined through tabulated airfoil data. Since the turbine geometry is removed, no boundary layer resolution is required, yielding a lean grid and a much lower computational effort. Furthermore, if a wind farm installation is considered, the prediction of the wake interaction through an ALM implementation greatly simplifies the simulation complexity.
In recent years, several authors performed numerical analyses employing ALM. Mikkelsen et al. [7] employed an ALM combined with FLEX5 aeroelastic code and the LES solver EllipSys3D to study the wake interaction between three HAWTs in a row based on the Tjæreborg turbine geometry. Optimal pitch settings were investigated. The same flow solver was coupled to an ALM implementation by Shen et al. [8], who performed a validation considering Mexico and NREL Phase VI wind turbines under yaw conditions. Churchfield et al. [9] coupled an ALM code to an OpenFOAM LES solver. Improvements in velocity and body-force distribution function were introduced, and the code was validated on NREL Phase VI and NREL 5-MW wind turbines. Positive results were presented for the prediction of near wake behavior and tip vortices. Weihing et al. [10] developed an ALM running with the LES solver FLOWer. Two commercial wind turbines were considered under offshore and onshore conditions, and good agreement was achieved in comparison with 3D URANS simulations. Matiz-Chicacausa and Lopez [11] investigated the NREL Phase VI in downwind configuration, focusing on the tower shadow effect. An OpenFOAM URANS solver was employed, and the turbine was simulated using both ALM modeling and the full 3D geometry, pointing out matching with measurements. Wang et al. [12] studied the wake of two offshore NREL 5-MW wind turbines through ALM and OpenFOAM LES solver. Different yaw and tilt angles were tested, and a control strategy was proposed. Ravensbergen et al. [13] implemented an ALM within a variational multiscale flow solver, and NREL 5-MW, NTNU, and NREL Phase VI wind turbines were analyzed. The simulations were carried out in standalone, back-to-back, and complex terrain conditions.
Despite several benefits, the application of ALM in the literature has been limited to validation in yaw conditions [8]. In fact, further studies are required to reach more comprehensive and detailed understanding and clarify the weaknesses and powers of this method in yaw conditions.
The present paper aims to investigate the application of ALM in the prediction of wind turbine wake in yaw operation, which is vital to achieve reliable wind farm design. The method is implemented within the OpenFOAM flow solver. A relatively new dynamic stall model that is modified by Sheng for wind turbine conditions [14] is employed. The model is validated on the NREL Phase VI test case. The results are then compared with a lower-order method and a higher-fidelity approach, namely, 1D BEM and 3D URANS.
The outline of the remaining paper is as follows. In the next section, the computational model that includes the employed actuator line model, dynamic stall model, URANS governing equations, and geometric model are described. In the results and discussion section, after verification and validation of the current model, the result of yawed case is investigated. Finally, a summary and conclusions section is provided.
2. Computational Model
2.1. Actuator Line Model (ALM)
OpenFOAM package is a set of C++ libraries meant to solve partial differential equations. In the current paper, a C++ library was developed for the OpenFOAM toolbox to implement the actuator line model. The ALM employs a series of spanwise element points, representing sections with a constant airfoil, twist, chord, and oncoming wind. The volumetric forces caused by sections are projected into the flow field. This procedure is accomplished through source terms in the momentum equations and applies the same forces on the fluid as the turbine blades. The blade aerodynamic forces are calculated by coefficients from tabulated airfoil data. These forces are formulated as follows:where is the local angle of attack measured as the angle between the chord and the local relative flow, is the lift coefficient, and is the drag coefficient obtained from tabulated airfoil data. The lift and drag coefficients are linearly interpolated for the local Reynolds number. Here, denotes the relative velocity, resulting from inflow velocity and the element velocity . A sampling approach was employed to interpolate from the CFD resolved velocity field. A total of 16 sampling points were considered at a distance of from the location of elements. Then, a Gaussian function is utilized for the smooth projection of the calculated forces onto the flow field. It is defined as follows [15]:where is the volumetric force field in cells’ position, is the tip correction function, and characterizes the forces at actuator element points that are calculated by lift and drag forces. In this equation, is the body index (i.e., blades, hub, and tower), and describes the actuator element point index. A tip treatment was applied to consider the tip effect induced by the pressure equalization from the suction and pressure sides [16].where is the number of blades, is the pitch angle (including twist), is the blade radius, describes the location of the element point in the spanwise direction, and is the tip speed ratio.
2.2. Dynamic Stall Model
In this study, the modified Leishman–Beddoes dynamic stall (LB-DS) model by Sheng et al. [14] is utilized. Since the LB-DS model was first proposed for helicopter applications [17], Sheng has modified the model to coincide with the use of HAWTs. In this situation, flow is typically incompressible, and airfoils are relatively thick. The principles of this model are presented in this section. For more details, readers are referred to references [14, 18]. The LB-based models include three parts: the unsteady attached flow, stall onset, and separated flow. The air loads for the attached flow consist of impulsive and circulatory modules. The total normal force coefficient in the unsteady attached flow condition is equal to the sum of impulsive and circulatory normal force coefficient [14]:where is a deficiency function that is defined based on the LB 3G model [19], is the normal force curve slope, and is an equivalent of incidence, which is calculated by the geometrical angle of attack and three deficiency functions.
The deficiency functions of , , and are empirically driven by using the flow velocity and pitch rate. Additionally, there is a lag due to leading-edge pressure and boundary layer development, which is implemented by applying a first-order lag to :where is the delayed angle of attack, and the deficiency function is introduced as follows [14]:where is a time constant associated with pressure response, and it is principally independent of the airfoil shape. The nondimensional time step is
Sheng et al. defined a stall-onset criterion for the angle that the dynamic stall may begin. This criterion considers a reduced pitch rate [14].
The reduced pitch rate is expressed by the following equation:
The separated flow solution includes two indicial responses: (I) trailing-edge separation and (II) leading-edge vortex shedding. The trailing-edge separation module is associated with nonlinear effects of trailing-edge flow separation through a dimensionless separation point parameter. This parameter is attained by Kirchhoff flow approximation in unsteady conditions. It is described as follows [14]:where , , and are constants that change with Reynolds number. In addition to the pressure response delay, boundary layer (viscous) lag is formulated as follows [14]:where is the deficiency function, and it is defined as follows [14]:where is a time constant, which is dependent on the Mach number.
The normal force coefficient for dynamic separation point is superimposed with impulsive loading as follows [14]:
There is an additional overshoot in the normal force during vortex shedding. The vortex-induced normal force coefficient is expressed as follows [14]:where and are the dynamic and static separation points, respectively. is a modulation parameter, which is given by a periodic function that incorporates the track of the concentrated vorticity location. The above equation increases the vortex lift only when the leading-edge vortex (LEV) starts convection over the suction side.
By superposition of the loading, the total normal force coefficient is obtained as follows [14]:
The chordwise force coefficient, pitching moment coefficient, and constants are considered according to Sheng guidelines [14, 18].
2.3. Definitions
The nondimensional parameters of the tip speed ratio (), power coefficient (), and thrust coefficient () are expressed by equations (17)–(19), respectively.where is the rotor swept area, P describes the produced power from the product of the radius and the tangential force , and is the axial force acting on the rotor.
2.4. Governing Equations
The turbulent flow field around the wind turbine was computed by the finite volume method (FVM). The three-dimensional Reynolds-averaged Navier–Stokes (RANS) equations are formulated as follows [20]:
The Reynolds stress term was modeled based on the Boussinesq hypothesis according to equation (21), whereby the Reynolds stresses are concerned with the mean velocity gradients.
The realizable k-epsilon turbulence model (equations (22) and (23)) was used for its relatively low computational cost associated with the computation of the turbulent viscosity [21].
Regarding the temporal discretization, a second-order Crank–Nicolson method was applied for time integration. The PIMPLE algorithm was implemented for pressure-velocity coupling. This algorithm is a combination of SIMPLE and PISO algorithms. The matrices are solved by choosing the generalized geometric-algebraic multigrid (GAMG) for pressure and a preconditioned biconjugated gradient (PBiCG) for velocity [20]. The details of BEM and CFD data that are used in this paper for comparison are described in references [22, 23].
2.5. Geometric Model
In this study, the National Renewable Energy Laboratory (NREL) Phase VI reference wind turbine was modeled based on the sequence S settings represented in reference [24]. In this test sequence, an upwind, rigid turbine with a 0° cone angle was used, and the test was implemented for a wind speed, ranging from 5 m/s to 25 m/s. This case is comprised of a constant speed, constant pitch, and stall-regulated wind turbine. It has a two-bladed rotor with twisted, tapered blades and shaped using the S809 airfoil [25, 26]. The principal characteristics of the turbine are summarized in Table 1. The chord and twist distributions are illustrated in Figure 1.
(a)
(b)
(c)
3. Results and Discussion
3.1. Verification and Validation
Figure 2 shows the cuboid computational domain, which consists of nonconformal hexahedral (cubic) grids, which allows local refinement in the areas of interest that significantly reduce the computational resources required for the simulation. Figure 2(a) expresses the boundary conditions for the wind turbine simulation. For the inlet condition, a uniform freestream velocity boundary was employed, whose magnitude was set to the examined velocity at the NREL experiment [24], and its direction was set parallel to the X-direction. The outlet condition was set as the pressure outlet boundary. A standard atmospheric pressure was considered for the reference pressure and imposed on the boundary. Therefore, the dynamic pressure field is calculated considering this condition. No-slip boundary conditions were applied to the two lateral walls, the base, and the upper side of the domain.
(a)
(b)
(c)
All the studied cases were simulated by considering a domain that is geometrically similar to the wind tunnel used in reference [24]. For this reason, the width and height of the domain section were set to 24.4 m and 36.6 m. The effect of upwind and downwind distances (, ) on the of the turbine was investigated as described in Table 2. This examination was performed at a wind speed equal to 7 m/s. While the mesh size kept the same, the domain length was extended according to Table 2. A comparison was made between the data and the longest domain (case 7). The results revealed that if the cases were longer than case 1, the difference was less than 0.1%. A compromise was reached between computational cost and accuracy, and thus case 2 (bolded) was chosen for this study.
The effect of time step size was evaluated on the produced power coefficient (). The results of the evaluation are given in Table 3. The required data are obtained for cases with 0° and 30° of yaw angles. It was performed to ensure involving phenomenon associated with yaw condition. For these tests, the wind speed was set to 7 m/s. The shortest time step with a value of = 0.216° was used as the reference case. In order to make data comparable, the ratio of is expressed in Figure 3, where the subscript “ref” refers to case 1. Figure 3 depicts for time steps smaller than = 0.432°; there was a difference of less than 2% in compared to the shortest time step size.
As shown in Figure 2, a structured hexahedral mesh was generated for the ALM-CFD calculation considering the mesh block structure. The computational domain was divided into four zones (see Figure 2(c)) to provide coarser grids in the far region. The cell size increases two times between each zone. The investigated cell size at the rotor disc is reported in Figure 4. As shown in this figure, a case with mesh less than 0.125 m at the rotor disk can be chosen as the optimum mesh size with a balance between cost and accuracy.
The independence of the following simulation results from the mentioned numerical factors was verified; the results were compared to the available data in the literature. In Figure 5, the numerical results for NREL Phase VI are reported together with the NREL experimental data [24]. The produced torque is reported as a function of wind speed, and as it can be observed, the computed aerodynamic torque agrees with the standard deviations for the tested wind speeds [27, 28].
3.2. Yawed Case
Figure 6 shows the employed notation for yaw angle. Here, a positive yaw angle corresponds to a reduced inflow angle for the blade at 0° azimuth angle (12 o’clock). The right-handed Cartesian coordinate system is also depicted in this figure, whereby is in the -direction and tower is in the -direction.
In Figure 7(a), the produced aerodynamic torque of the turbine is compared for different yaw angles in = 10 m/s. It is evident that an increase in the yaw results in the oscillation of produced power. It can be observed that by increasing the yaw angle to 10°, 20°, and 30°, the oscillation amplitude of produced torque increases 5%, 18%, and 20% of their average, respectively. This oscillation is addressed by the angle of attack variation that is due to the turbine yaw misalignment. Although the wind speed during rotation is almost constant, the yaw misalignment leads to oscillation in the angle between wind chord and wind flow oscillate during a rotor revolution. This results in oscillation in both relative velocity and angle of attack.
(a)
(b)
The magnitude of the geometrical relative velocity and the geometrical angle of attack are calculated using equations (24) and (25), respectively, as follows [29]:where and denote the yaw angle and the azimuth angle and stands for the pitch angle that includes the twist. The effect of turbine on the flow field is absent in the geometrical relative velocity and the geometrical angle of attack. Figure 7(b) illustrates how the oscillation of geometrical angle of attack varies with yaw angle during a rotor revolution. It should also be mentioned, while blade 1 is at its maximum angle of attack, the second blade is at its minimum.
An increment in yaw misalignment raises the oscillation amplitude while the average decreases. This is a crucial phenomenon in wind turbine design. Nonetheless, simulation of it using a method like BEM is challenging at higher wind speed.
The geometric angle of attack at the blade tip for the yaw angle of 30° is plotted in Figure 8(a) for different wind speeds. As shown in the figure, both the amplitude and the average of grow with an increase in wind speed. The minimum angle of attack occurs at the azimuth angle of = 0°, after which it increases until its maximum at = 180°. The opposite is true when the yaw angle is in the negative direction. The geometrical relative velocity, , is also plotted in Figure 8(b). This figure shows that the variation of geometrical relative velocity varies inversely with the geometrical angle of attack, whereby the geometrical relative velocity reaches its maximum and minimum at azimuth angles of 0° and 180°, respectively.
(a)
(b)
As shown in Figure 8, turbines under yaw misalignment lead to a variation in the angle of attack. Hence, before employing the Sheng DS model in the study of yawed cases, an assessment was established between the Sheng DS model results and the experimental measurements reported in [14]. For this part, the simulation was performed for an isolated airfoil, whereas in the rest of this paper, an entire turbine has been simulated using the ALM. The results of the ramp-up test at a reduced pitch rate of 0.005 are plotted in Figure 9(a). The results for a single S809 produced by the employed DS model were very close to the experimental data in terms of the normal force. It predicted the stall onset just 2° earlier, and the effect of the vortex emitted at 27° was less than that in the experimental data; however, the results were still in good agreement. Furthermore, the results of a single oscillatory airfoil were compared with the experimental data, as reported in Figure 9(b). The hysteresis loop is obtained for an oscillatory airfoil with α = 15 ± 10° and a reduced frequency of . Again, the predicted stall was slightly earlier. It should be noted that the employed DS model was initially developed for helicopter applications, and even the modified version by Sheng inherited the main characteristics of the original model, e.g., employing the Kirchhoff equation, which is more suitable for thinner airfoils. Despite this, employing the Sheng DS model established significant capability in modeling delays and downstroke movement in the dynamic stall.
(a)
(b)
In Figure 10, the aerodynamic torque calculated using the current ALM code is compared with UBEM and CFD results, the details of which are available in [22, 23]. Data are presented for different wind speeds of 10 m/s, 13 m/s, and 15 m/s corresponding to TSR equal to 3.8, 2.9, and 2.5, respectively. All three figures were obtained for the NREL Phase VI wind turbine under 30° of yaw angle. In general, the accuracy and reliability of the presented ALM method were confirmed by CFD computation as the percentage error in the average aerodynamic torque between CFD, and the ALM calculation tools were very limited. The results can be compared in terms of amplitude and phase of the oscillating torque for further elaboration.
(a)
(b)
(c)
As can be observed, UBEM underestimated the torque for = 10 m/s and overestimated it for = 15 m/s; nonetheless, there was an outstanding improvement in the amplitude by ALM. It should be mentioned that UBEM uses yaw models to consider the nonorthogonality of flow to the rotor plane [30], whereas the ALM model obtains the inflow velocity from the resolved flow field; thus, all flow angles are inherently concerned.
Regarding the phase of oscillating torque, a phase shift error of about 30° and 15° appeared between UBEM and CFD results at 10 m/s and 13 m/s, respectively. Nevertheless, as the wind speed increased, the phase error decreased towards zero. Considering the angle of attack, the blades are subjected to dynamic stall at a wind speed of 10 m/s in some portions of a revolution. As the wind speed increases to 13 m/s, the blade is subjected to dynamic stall almost during the entire revolution. An outstanding improvement in decreasing the phase shift error appeared in the ALM results compared to UBEM. The phase was predicted accurately for = 10 m/s with a partial dynamic stall. However, the accuracy decreased slightly for = 13 m/s as there was an incremental oscillation after = 90°, caused by a delay in predicting the dynamic stall of the second blade. As the wind speed increased to 15 m/s, the high angle of attack during the revolution led to a complete flow separation and the formation of TEV. While the ALM uses the tabulated airfoil data, this three-dimensional phenomenon affects the inflow, resulting in a slightly inaccurate prediction by ALM.
The results obtained for the NREL Phase VI reference wind turbine with different static yaw misalignments are presented in Figure 11. Cases with yaw angles ranging from 0° to 30° were simulated in both negative and positive yaw directions. Square symbols show the average power, and the minimum and maximum values over the revolution are indicated by error bars. As shown, the results for negative and positive yaw angles are relatively symmetrical. The results indicate an average power reduction of 25% due to 30° yaw misalignment.
The variation of power with yaw angle was suggested to be approximated by [31]. The exponent was suggested to be equal to 3; however, it is only valid if the axial induction distribution is small compared to [31]. It was stated that measurements including wind tunnel models and full‐scale turbines have indicated the exponent could vary between 1.88 and 5.14 [32]. In the current case, = 2 shows to be in good agreement for the studied yaw misalignments by the ALM results, although it is not entirely fit.
Taking the amplitude of oscillations into consideration, it shows to grow continuously until 30°, after which there is no clear trend. This is related to the angles of attack as for cases with yaw angles greater than 30°, the blade is under stall condition in a portion of its revolution. The blades being under stall conditions lead to instabilities for further investigations. In other words, the amplitudes for these conditions are case-dependent, as they is a result of two stalled blades.
Figure 12 illustrates the wake structure as viewed from the side (x–z plane). The isosurface of Q visualizes the instantaneous vorticity, where the Q criterion is defined as follows [33]:where and represent the symmetric and antisymmetric parts of the velocity gradient tensor , respectively:
(a)
(b)
The flow field around a flow-aligned turbine and a 30° yawed case at a wind speed of 10 m/s is shown in Figures 12 and 13, respectively. For the flow-aligned turbine, the wake was initially developed behind the rotor blades in the form of a well-defined helical geometry. The wake structure was conserved until four rotor diameters downstream and then evolved into the turbulent wake. The strength of the helical wake vorticities was azimuthally symmetric; hence, the contour in Figure 12(b) was illustrated only for one time step.
(a)
(b)
For the yawed case, the contour at the blade section is illustrated for each 30° of revolution as indicated in Figure 13(b)I-XII. The comparison of yawed cases with the flow-aligned case indicates that the yaw misalignment yielded significantly different wake evolution. A periodic variation in the wake strength with respect to the azimuthal angle was manifested. The periodical deformation of the wake vorticities strength is addressed by the variation of AOA and relative velocity during revolutions. The corresponding AOA and relative velocity can be observed in Figure 8. As shown in Figure 13(b)I, a more robust wake was emitted when the blades were at the azimuth position of 0°. This position corresponds to the highest relative velocity during a revolution.
From a downstream wind turbine perspective, axial velocity distribution directly influences the inflow condition of the downstream wind turbines. The inflow condition has a crucial role in determining the aerodynamic and fatigue loads of the downstream wind turbine [34]. In Figures 14 and 15 , the streamwise velocity contours behind the aligned turbine are compared with the yawed wind turbine to ascertain the influence of emitted wake structure on the downstream velocity field. The data were reported at as indicated in Figures 14(a) and 15(a), where is the downstream distance and is the rotor diameter. Simulation indicated that the propagation of the wake lasted longer in the case of the aligned rotor. As it can be observed in contours at , the portion of contour in dark blue is broader for the aligned case.
(a)
(b)
(a)
(b)
The wake evolved in the form of a helical geometry with similar strength for the aligned case, as shown in Figure 14(a). Thus, the distribution of velocity deficit is axisymmetric concerning the wake age. As the wind turbine underwent 30° of yaw angle, a highly skewed deformed wake structure was developed (see Figure 15). According to Figure 15(b), the axial velocity distribution was no longer axisymmetric concerning the wake age because of the periodically oscillating wake vortices due to oscillating AOA and relative velocity. As strong wake vortices pass the downstream region, a high induced velocity region appeared at the azimuth angle of 240°. This agrees with the observation of Bastankhah and Porté-Agel [35] where the wake centre, defined as the point where the velocity deficit is maximum at each downwind location, was changed according to yaw angle. Thus, the downstream wind turbine receives significantly unsteady and asymmetric inflow conditions and blades suffer from a time-dependent aerodynamic load.
4. Conclusions
An actuator line model was implemented in OpenFOAM flow solver. The Shen’s tip treatment to consider the tip effect and Sheng’s LB dynamic stall model to consider the dynamic stall effects were applied. The ALM simulation results for the aligned case were in excellent agreement with available experiment, and the overall trend was reproduced for the output torque with respect to wind speed.
An increment in yaw angle leads to blade advancing and retreating effect whereby the angle of attack oscillates. The variation of produced power was investigated in the following three conditions:(1)The turbine at a wind speed of 10 m/s and under 30° of yaw represents an understall condition. In this condition, ALM results were in general matching in terms of both phase and amplitude.(2)The turbine at a wind speed of 13 m/s and under 30° of yaw is a condition that the blade is subjected to dynamic stall almost during the entire revolution. Employing the Sheng DS model established a significant capability in modeling delays and downstroke movement in the dynamic stall.(3)The turbine at a wind speed of 15 m/s and under 30° of yaw represents a deep stall-dominated flow. The phase was predicted accurately with a minor overestimation. Considering the fact that ALM reads airfoil coefficients from tabulated data, further modification in the model is required if this stall-dominated condition is in the scope of a study.
In general, despite the small number of cells, the ALM code is capable of predicting produced power with satisfying accuracy before deep stall-dominated conditions. Employing the DS model improved the prediction of phase delay in power variation.
Providing knowledge about turbine wakes is an advantage of ALM over BEM, which is a priority in wind farm design. According to the results, even using 400k cells, the ALM code could fairly simulate the skewed helical geometry of wake and capture the effect of advancing and retreating blade on the wake. It was found that the dominated flow by the helical wake lasts longer in the aligned case. On the other side, the deficit velocity is more concentrated for the yawed case. Considering these points can be exploited for wind farm optimization, ALM potentially promises to be utilized with less computational cost than CFD.
Nomenclature:
| : | Nondimensional constant for the vortex shedding (-) |
| : | Chord length (m) |
| ,: | Lift and drag coefficients (-) |
| : | Total normal force coefficient (-) |
| : | Slope of the normal force coefficient at the attached flow () |
| : | Circulatory normal force coefficient (-) |
| : | Impulsive normal force coefficient (-) |
| : | Normal force coefficient for the unsteady attached flow (-) |
| : | Normal force coefficient for the trailing-edge separation (-) |
| : | Normal force coefficient for the vortex shedding (-) |
| : | Power coefficient (-) |
| C1: | Realizable k-epsilon turbulence model coefficient (-) |
| C2: | Realizable k-epsilon turbulence model coefficient (-) |
| : | Effective diffusivity (-) |
| : | Deficiency function for separation point (-) |
| : | Deficiency function for geometrical angle of attack () |
| , : | Delayed flow separation points (-) |
| : | The volumetric force field in cell position () |
| , : | Lift and drag forces () |
| : | Deficiency function for the impulsive force coefficient (-) |
| : | Gravity acceleration () |
| : | Turbulent kinetic energy production rate due to the anisotropic part of the Reynolds stress tensor () |
| : | Mach number (-) |
| : | Pressure () |
| : | Distance from actuator point to the point where the force is applied (m) |
| : | Reduced pitch rate [-] |
| : | Nondimensional time [-] |
| : | Source terms |
| , : | Coefficients for the flow separation point () |
| : | Time (s) |
| : | Time constant, associated with the separation point delay (-) |
| : | Time constant, associated with the stall onset (-) |
| : | Flow velocity |
| : | Relative velocity |
| : | Freestream velocity |
| : | Modulation parameter for the vortex shedding |
| : | Volume of cell |
| : | Section width (). |
| : | Local angle of attack |
| : | Time derivative of angle of attack |
| : | Delayed angle of attack (rad) |
| : | Equivalent angle of attack (rad) |
| : | Critical stall-onset angle of attack (rad) |
| : | Constant critical stall-onset angle of attack (rad) |
| : | Static stall-onset angle (rad) |
| : | Angle of attack for the breakpoint of the flow separation (rad) |
| : | Pitch angle |
| : | Yaw angle |
| : | Step change in angle of attack or in time |
| : | Projection width () |
| : | Turbulent kinetic energy dissipation rate () |
| : | Angular displacement |
| : | Reduced frequency [-] |
| : | Tip speed ratio (-) |
| : | Dynamic viscosity (m2·s−1) |
| : | Dynamic eddy viscosity (m2·s−1) |
| : | The density of the fluid (kg·m−3) |
| : | Azimuth angle. |
| ALM: | Actuator line model |
| BEM: | Blade element momentum |
| UBEM: | Unsteady blade element momentum. |
Data Availability
Data are available on request. Request can be sent to Sahar Jannesarahmadi, [email protected], or Alireza Arabgolarcheh, [email protected].
Disclosure
This research received no specific grant from any funding agency in the public, commercial, or not-for-profit sectors.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.