1. Introduction
Solar photovoltaic (PV) panels are being widely deployed around the globe in multi-row arrays as an essential step towards net-zero economies [
1]. As the economic performance of solar PV farms relies on maximising the energy yield during daylight, they have recently being equipped with tracking systems that vary the tilt angle (
) of the panel according to the position of the Sun to optimise the incidence angle [
2]. These can be single-axis with a horizontal or vertical rotation, or a combination of both for dual-axis systems. Whilst the latter can provide increases in annual energy yield of approximately 40% [
2], it is limited to small rows of panels as the swept area is very large, reducing power density. Thus, single-axis systems are mostly adopted by the industry with a torsional motion in the horizontal axis, as shown in
Figure 1. This tracking system comprises a torque tube that rotates the panel commonly in a range of tilt angles of
, which increases the energy yield by 20% compared to fixed solar PV panels [
3] with most improvements during the early or late hours of sunlight [
4].
The extreme slenderness of the solar PV panels, however, has brought aerodynamic challenges that need to be overcome as they can undergo extreme wind-induced motions that trigger aeroelastic phenomena which can eventually collapse the structure [
5]. The aerodynamic response of the structure mainly depends on the incoming wind conditions, namely speed and direction, and the tilt angle adopted by the panel. Recent research has suggested that the mechanical properties of the solar tracker, such as its inertia, damping and stiffness, determines whether a given wind speed triggers the phenomenon of flutter [
6,
7]. This fluid–structure interaction is complex due the unique characteristics of solar panels such as: (1) they operate in a wide range of tilt angles, including positions parallel to the ground up to
; (2) they have extremely slender structures (i.e., the length-to-thickness ratio can be well beyond 20) which decouples the turbulent flow structures generated at the leading and trailing edges; (3) they have a single degree of freedom that corresponds to the torsional moment without vertical or horizontal displacement; (4) the natural frequencies of the system derive from the stiffness, only provided by the torsional axis, and inertia, which only depends on that of the panel; (5) the distance to the bottom ground
h relative to the panel’s length
L is relatively small, e.g.,
is normally between 1.2 and 1.8, which generates an aerodynamic ground effect that changes the wake dynamics; and, (6) in large-scale farms, the different rows of panels induce an internal boundary layer that modify the onset flow conditions for each row [
8,
9].
Compared to other classic fluid mechanics problems, e.g., the bridge deck section or airfoil aerodynamics, the flow around solar PV panels presents notable differences that restrains us from applying theoretical concepts from the former fields to their aerodynamics. For instance, the flow around bridge decks has been an active topic of research in civil engineering but differs from the flow around solar PV panels as the angles of attack normally considered do not exceed
, the transverse sections have a length-to-thickness ratios below 10, and have at least three degrees of freedom, e.g., torsional motion with vertical and horizontal displacement [
10,
11,
12]. The thickness of the bridge decks induces recirculation areas on both the upper and lower surfaces near the leading edge which do not occur in solar panels as they are much more slender. Furthermore, in the flow around airfoils, large angles of attack can be studied, e.g., to account for dynamic stall, but this is not directly applicable to solar panels as they are curved surfaces designed to optimise the lift-to-drag ratio, with the capability of delaying flow separation to occur closer to the trailing edge [
13,
14,
15]. Finally, most studies looking at the flow around thin flat plates consider inclinations that are nearly parallel or perpendicular to the flow direction and, more importantly, do not account for the effect of the proximity to the ground [
16,
17].
To date, most research looking at the flow around solar PV panels has focused on structures with a short spanwise length, typically accounting for standalone panels. In these cases with a small chord-to-width aspect ratio, the shear generated over the lateral boundaries plays a significant role in generating vortices with the flow being highly three-dimensional [
6,
18,
19]. In solar farms, the PV panels are deployed in long rows, e.g., with 50 panels or more in the transverse direction, so that the vortical structures generated by the lateral ends have negligible influence on the overall aerodynamics [
20,
21]. At these solar farms, extreme wind events that surpass the design speed and lead to an unsteady response of the structure are likely to cause its structural collapse, partial fracture or fatigue loading [
22]. Such a fluid–structure interaction depends on the natural frequency of the panel, however, more importantly, it changes with the operational tilt angle and row number at which the panel sits [
7,
23]. Experimental tests in wind tunnels reveal that this aeroelastic phenomenon is at a fixed frequency but the underlying unsteady mechanism triggering it is to be fully understood yet for every tilt angle, which is the scope of this paper.
Earlier numerical studies of flow around solar PV panels adopted the Reynolds averaged Navier–Stokes (RANS) turbulent closures, as often being the industry standard. The steady nature of RANS introduces excessive numerical viscosity in the flow developed on the lee side of the panel, depicting a main recirculating bubble in the suction side with no unsteady vortices being shed from the shear layers [
19,
24,
25,
26], which have been partly seen in experiments of thin flat plates [
27,
28]. Due to the transient nature of the turbulent flow developed behind the panel, the eddy-resolving large-eddy simulation (LES) turbulence closure needs to be used as it resolves the flow scales smaller then a filter size, commonly the grid size, whilst modelling the smaller ones through a sub-grid scale [
16].
This paper aimed to investigate the nature and coherence of the vortical scales generated behind a single solar PV panel under uniform flow conditions varying the tilt angle in the range between
. We performed a three-dimensional high-fidelity LES to assess how the vortex shedding pattern and frequency varies for each tilt angle without inlet turbulence, and this can induce a faster loss of coherence of the flow structures. The paper is structured as follows:
Section 2 presents the in-house LES code used for the computations and the numerical setup.
Section 3 presents the results concerning the mean velocity field and recirculation regions, instantaneous vortical structures, and aerodynamic coefficients. Finally, the main conclusions are drawn in
Section 4.
2. Numerical Framework and Setup
Large-eddy simulation (LES) was performed, adopting the well-validated in-house code Hydro3D [
29,
30,
31,
32] that adopts the spatially filtered Navier–Stokes equations to resolve the largest flow structures whilst modelling the small scales with the wall-adapting local eddy-viscosity (WALE) sub-grid scale model [
33]. These equations, namely mass and momentum conservation, read:
Here, = is the spatially filtered velocity vector (for convenience, the symbol is omitted hereafter), the coordinates vector is = , denotes the fluid density, p is the relative pressure, is the kinematic viscosity of the fluid, and is a source term that accounts for the immersed boundary method (IBM) adopted to represent the tracker.
Hydro3D adopts a fractional-step method to advance the simulation in time, including a low-storage three-step Runge–Kutta method for calculating the convective and viscous terms that provide the non-divergence free velocity, and a multigrid method to resolve the Poisson pressure equation used to obtain the divergence-free velocity and pressure fields. Velocities are stored in a staggered arrangement in a rectangular Cartesian grid, as a computationally efficient methodology to perform simulations on multiple processors using the message passing interface (MPI) protocol. A local-mesh refinement (LMR) method allows to refine the grid in the vicinity of the area of interest [
29], which in this study is the solar PV panel.
The direct–forcing Immersed Boundary Method (IBM) used to compute the fluid–structure interaction ensures that a no-slip boundary condition at the solid interface is accomplished using an interpolation function in the communication between fluid and solid frameworks, which in our case is a delta function
from Yang et al. [
34]. Whilst details about the direct–forcing IBM implementation in Hydro3D and its validation can be found in Kara et al. [
35] and Ouro et al. [
36] for cylinder flows and in Ouro and Stoesser [
37] for tidal stream turbines, a brief description is provided below. First, the predicted fluid velocities are obtained in the Cartesian grid of the fluid which are interpolated onto the Lagrangian particles comprising the solid mesh. Then, the reaction forces on the Lagrangian grid (
=
,
,
) are computed to enforce the no-slip condition. In this case, a multi-direct forcing method is adopted to improve the accuracy of this method [
38] with two additional loops. As shown in
Figure 2, the tilt moment is obtained by integrating the streamwise and vertical forces at every Lagrangian point (L) according to:
Here, the volume assigned to every Lagrangian particle (
) is equal to that of a fluid cell (
) as required for numerical stability [
37], and
indicates the time step. The overall forces in the Cartesian directions (
,
) can be divided into tangential (
) and normal (
) components depending on the tilt angle (
), as:
The corresponding aerodynamic coefficients for these forces (
and
) and tilt moment (
) are:
In the following, positive tilt angles refer to those when the leading edge is above the height of the rotational centre, and vice versa, as indicated in
Figure 2.
Computational Setup
A rectangular solar PV panel that is 2.3 m long (
L) and has a thickness (
h) of 0.08 m, yielding a chord-to-thickness aspect ratio of approximately 30, is adopted as these are similar dimensions to those adopted by the industry in solar farms. The rotational axis is at the centre of the panel, around which the tilt moment is computed (Equation (
3)), at a height (
H) of 1.53 m from the bottom ground, i.e., resulting in a gap ratio
= 1.5. A total of 21 tilt angles are simulated with
= 0
, ± 2.5
, ± 5
, ± 7.5
, ± 10
, ± 15
, ± 20
, ± 25
, ± 30
, ± 35
, ± 40
, ± 45
, ± 50
, ± 55
, and ± 60
.
A side view of the computational domain is presented in
Figure 3, whose dimensions are 42 m (18.3
L) by 12 m (5.2
L) by 1 m (0.44
L) in the streamwise, vertical, and spanwise directions, respectively. The width of the computational domain is reduced as periodic boundary conditions are applied in this direction to simulate an infinitely long solar panel without wall effects as those that can be found in wind tunnels. Three levels of local-mesh refinement [
29] are used to optimise the use of computational resources whilst keeping a high spatial resolution near the panel. To avoid any influence from the inlet velocity conditions, e.g., vertical shear from a logarithmic boundary layer, a uniform velocity of
= 10 m/s is applied at the inflow boundary condition with a Neumann condition set at the outflow. No turbulence is included at the inlet so as to ease the depiction of the turbulent structures shed by the solar panel, similarly to when aerodynamic derivatives from sectional analysis are determined in wind tunnel testing [
23]. At the bottom, a wall function for hydrodynamically smooth walls is applied, whilst a no-shear slip condition is used at the top of the domain. The solar panel’s gravity centre is located at 11.85 m (5.1
L) from the inlet.
A mesh resolution of = = 0.01 m is adopted, i.e., 230 grid cells are distributed along the longitudinal direction and eight cells across its thickness, whilst in the spanwise direction, the resolution is = 0.025 m, totalling 40 cells over the transverse direction. The time step is variable according to a CFL condition of 0.6 to guarantee numerical stability. Each simulation runs for a physical time of approximately 180 s in total, with mean velocities averaged after the first 20 s and second-order moments after 50 s, adopting 36 CPUs with an equivalent computational load of 5000 CPUhours.
3. Results
In this section, we first present the results of the time-averaged streamwise velocities and instantaneous vorticity field, followed by the mean structural forces acting on the solar PV panel. The spectral analysis of the force time history is then adopted to compare the response of the panel for positive and negative tilt angles, and the Strouhal number obtained from the peak-energy frequency is presented and compared to the experimental data from the referenced studies.
3.1. Time-Averaged Flow Field
The time-averaged flow field generated on the suction side of the panels is characterised by a recirculation region whose pattern changes with the tilt angle and depending on whether this is positive or negative.
Figure 4 presents the contours of normalised mean streamwise velocity (
) with streamlines of mean flow for positive tilt angles ranging from
until
. For angles smaller than
, a recirculation bubble is generated over the suction side of the panel with an extension that increases vertically and longitudinally for larger tilt angles, as with
, this bubble extends over the whole upper surface but its centre is almost at
, i.e., at a similar location to the gravity and rotational centre. At
, the bubble’s centre is displaced further downwind and gets closer to the trailing edge.
As the tilt angle increases further, the recirculation region stemming from the leading edge grows in size and, after
, a smaller bubble starts to form at the trailing edge as a result of vortices being generated and shed at this extreme. As the trailing edge gets closer to the bottom surface, it generates a nozzle effect characterised by a high-momentum region underneath the panel. This, in turn, affects the lower recirculation area which impacts the development of the trailing edge vortices if compared to a boundless scenario, e.g., Fage and Johansen [
27], Chen and Fang [
28]. For tilt angles larger than
, both recirculation regions have a similar size but they are not fully symmetric.
The time-averaged streamwise velocity field obtained for the negative tilt angle cases is presented in
Figure 5. A similar recirculation pattern over the suction side for
is observed compared to their positive counterparts (
Figure 4) but with some differences. For
, the recirculating bubble does not extend entirely over the panel’s length as it did for
, and this is also notably reduced for
. With negative tilt angles, the leading edge approaches the bottom surface to trigger a diffuser-like effect, i.e., there is a flow acceleration underneath the solar panel after surpassing the leading edge that modifies the flow developed over the suction side for negative tilt angles. This constrains the behaviour of the vortex generation (as explained later) and the recirculation region. Such an effect is exacerbated for larger tilt angles, as the bottom gap in the becomes smaller, again characterised by two low-velocity recirculation regions, but these feature a longer extension in the streamwise direction owed to a faster bypass flow coming underneath the leading edge.
3.2. Instantaneous Velocity Field
To provide insights into the flow dynamics, the instantaneous vortical structures generated for the simulated tilt angles are presented in
Figure 6 and
Figure 7 for positive and negative angles, respectively, up to
35
. For tilt angles lower than
, vortices shed from the leading edge are generated and these move along the suction side of the solar PV panel. Note that the slenderness of this structure makes its thickness almost irrelevant to generate an significant pressure increase at the leading edge (or velocity reduction) when the flow impinges the structure, as in the case of bridge deck sections or elongated squared cylinders [
10,
11,
12,
39].
After
, a small recirculation region is formed at the leading edge, constrained by the laminar shear layer above it. Further increasing the tilt angle leads to the generation of an enclosed low-velocity area over the suction side that is not present for the smallest angles and which increases the lift force on the structure. It is after
that large leading-edge structures are generated and shed with an almost constant regular frequency, which we denote as a von Kármán vortex that generates a von Kármán vortex street [
27]. Note that previous studies looking at solar PV panel aerodynamics were performed with RANS models which damped the shear-layer instabilities observed in the present LES and failed to capture the von Kármán vortices due to the excessive numerical viscosity introduced [
19,
26].
For tilt angles larger than , the size of the von Kármán vortex is approximately that of the projected vertical length of the panel ( sin ), i.e., it grows in size with an increasing tilt angle. This large-scale structure is shed once it grows enough in size proportionally to and becomes unstable. However, once a tilt angle of is surpassed, a trailing-edge vortex starts to develop recurrently, which induces the shedding of the dominant von Kármán vortices. After , the interaction between these flow structures generated at the leading and trailing edge is enhanced and starts to shed at the same frequency, i.e., the vortex shedding mechanism behind the panel behaves as a classic bluff body.
Comparing the aerodynamic response of the flow around the solar panel for positive and negative angles, it is seen that below a tilt angle range of , few variations are observed, as the ground effect remains almost negligible. The latter effect starts to play an important role when the angle increases after this threshold, as the adopted ratio equal to 1.5 is relatively small. This leads to an asymmetric vortex shedding depending on the orientation of the panel. At negative tilt angles, the vortices developed over the suction side hit the bottom ground, which constrains their development. For positive angles, there is no physical boundaries for the vortices over the suction side to form, which increases their strength and coherence.
3.3. Mean Forces on the Solar PV Panel
The forces acting on the panel vary depending on the tilt angle magnitude and whether this is positive or negative. The immersed boundary method provides the reaction forces of the panel in the Cartesian coordinates (
and
) which are adopted to compute the torsional moment (
, Equation (
3)) as well as the drag and normal force directions.
Figure 8 presents the drag, normal force and torsional moment coefficients together with the force eccentricity normalised by the panel’s length (
with
). The current LES results are compared with those of Taylor and Browne [
23] who used a section model to determine the aerodynamic characteristics of a panel under no turbulent flow conditions. The drag force coefficient
indicates that, with increasing tilt angle, the force on the panel in the wind direction increases, which is in agreement with the experimental data [
23]. For panel positions at
,
values have a similar magnitude between positive and negative angles. However, the effect of the proximity to the ground starts to be noticeable at
with an asymmetric force increment, i.e., forces are larger with positive angles (the panel acts as a diffuser) compared to the results for negative values (when the leading edge is pointing towards the bottom surface). Note that the current results are for a solar PV panel configuration with a gap ratio
equal to 1.5 which can vary these coefficients if this parameter changes as in the case of circular cylinders [
36].
The coefficient of the force normal to the panel’s surface () features a linear increase up to . With an increasing tilt angle in the negative value range, the steady increment in continues almost linearly from until reaching a maximum value of 2.3 at = –60. For positive tilt angles, after , the lift coefficient rapidly increases until = –3.0 for after which there is a slight increase at larger angles up to –3.3 attained with . This asymmetric behaviour is again a consequence of the ground effect having a different impact in function of whether the suction side of the panel faces the bottom ground or not.
The normalised eccentricity reaches its maximum value when the panel is flat with a value of –0.19. This rapidly decreases to –0.1 at and remains nearly constant until , which coincides with a linear increase in and . Thereafter, for a larger , an asymmetric increase is observed with positive angles featuring lower values of , especially at with –0.01, whilst at negative angles, the minimum eccentricity is –0.075 for = –60.
Based on the distribution of force coefficients and eccentricity, three regions can be distinguished in the distribution of the aerodynamic forces. First, for
, the solar PV panel behaves accordingly with the thin airfoil theory in which the rate-of-change of
is constant, similarly to that of the
. This linear behaviour is derived from the steady increase in the recirculating bubble shown in
Figure 4 and
Figure 5. The second region is found in the range between
, corresponding to when there is flow reattachment over the panel’s suction side leading to an increase in lift forces but at a lower rate. Here, the moment coefficient is at its largest at
. In the last region, with
, the panel behaves similarly to a bluff body with a quasi-periodic vortex shedding with the turbulent structures shed from the trailing edge becoming notably larger than for smaller tilt angles.
3.4. Spectral Analysis of Forces and Strouhal Number Curve
The frequency of the vortex shedding is investigated with the Fourier transform using the
pwelch technique applied to the time history of the structure’s torsional moment (
) and vertical velocity fluctuation at a point located at
,
downwind. Spectra obtained for cases with
and
are shown in
Figure 9 as they represent the scenarios for the highest torsional moment and stall conditions, respectively. The most energetic peak is observed at a frequency close to the unity that corresponds to the shedding of von Kármán structures as their signature is also captured in the velocity spectra in the panel’s wake characterised by a single peak at the same frequency (
Figure 9c,d). Note that, in the latter plots, the Kolmogorov’s turbulence slope of –5/3 is observed over a frequency range starting at the most energetic peak up to approximately 30 Hz after which the energy decay rate increases.
The frequency of this dominant von Kármán flow structure is lower for positive angles than for negative ones (
Figure 9a,b), which is attributed to the proximity of the panel to the ground. For negative angles, this ground-effect impacts the wake development by constraining the free expansion of the vortex shedding normal to the panel’s surface, as seen in
Figure 6 and
Figure 7. Whilst the energy peaks for the von Kármán structures extend over a relatively broad range of frequencies, other series of distinct peaks appear in the spectra of
at higher frequencies. These occur at a narrower frequency band and correspond to the rapid shedding of vortices from the leading and trailing edge. The harmonics of these shear-layer vortices appear at all tilt angles, and after the third harmonic, the spectral energy of
rapidly decays with a slope of –3. Note that only for positive angles is there also a harmonic at twice the frequency of the von Kármán vortex, which is suppressed for negative angles.
The shedding frequency of the von Kármán vortices observed in
Figure 9 is used to calculate the reduced frequency based on the Strouhal number (
, with
f being the frequency corresponding to the most energetic energy peak captured in the spectra). The values of
obtained from the current numerical simulations are presented in
Figure 10 together with the experimental data from Fage and Johansen [
27] and Chen and Fang [
28] obtained in wind tunnel tests without bottom-wall effects. Both LES and experimental results from Fage and Johansen [
27] agree that the maximum
for positive tilt angles occurs for
with values of
= 0.30–0.34 whilst those from Chen and Fang [
28] continue to increase owing to the trapezoidal shape of the plate adopted in the experiments. For larger positive tilt angles,
steadily decreases until the value of approximately 0.15–0.16 for the maximum incidence of
[
27,
28], approaching the Strouhal value of 0.159 [
40]. This asymptotic
that would correspond to the panel acting as a bluff body perpendicular to the flow direction also agrees with
= 0.16, as found in other experiments [
28,
41]. Note that in all these experimental studies, the body was away from the bottom surface.
For tilt angles below
, there is a nearly linear increase in
from the minimum at
with
= 0.12. Over the range of negative tilt angles, i.e., with the leading edge closer to the bottom surface, the overall pattern is similar to the positive tilt angles but with larger Strouhal values which again peak at
=
, as depicted in
Figure 9a,b. This curve asymmetry between the positive and negative tilt angles is due to the proximity of the panel’s edges to the bottom surface. For negative angles, the suction surface faces the bottom ground and generates a high-speed region that accelerates the vortex shedding (
Figure 7). Conversely, with positive angles, it is the trailing edge which is the one closer to the ground but with the suction side looking in the opposite direction from the bottom ground, allowing the structure to develop and shed vortices unaffected by the ground (
Figure 6).
The vertical projection of the solar PV panel provides the characteristic length
sin
that can be considered as the scaling factor in the increase in vortex shedding frequency [
27]. The normalised Strouhal number
is presented in
Figure 11 to take into account such variation of the vertically projected panel length as the tilt angle varies. The
curve shows two distinct regions. First, from
,
increases as a function of
until the maximum
is reached, i.e., at
, and
is equal to 0.165 corresponding to the asymptotic
. The second region of this curve is characterised by an almost constant value of
, especially for positive tilt angles, as for negative ones, the bottom wall effect impacts the wake dynamics. The differences between positive and negative tilt angles are expected to reduce when increasing the clearance between the PV panel and bottom surface, i.e., the
ratio. However, typical values adopted by industrial solar PV panel designs are in the range of
, in which the ground effect is likely to be always present. These results suggest that once the vortex shedding excitation is at its maximum for
, the solar PV panel develops a shedding pattern such as a bluff body independent of the projected vertical length.
3.5. Effect of Ratio on the Vortex Shedding
Finally, we investigate the role of the relative distance of the solar PV panel to the ground, i.e.,
ratio, in the vortex shedding dynamics.
Figure 12 shows the Strouhal number for cases with
= 1.3, 1.5 and 1.7, which evidences that this ratio only changes the vortex shedding frequency for negative angles, i.e., when the vortical structures are shed in the suction side constrained by the ground surface.
results for the tilt angles of –
, –
, and –
, which show that a larger
, i.e., with the PV panel closer to the ground, increases the vortex shedding frequency, especially at
= –
, when the Strouhal number is maximum.
Figure 12 outlines that for larger negative tilt angles, the bluff body wake generated by the PV panel neglects the ground effect.
4. Conclusions
This paper studies the aerodynamics developed behind a single solar photovoltaic (PV) panel for a wide range of tilt angles up to at a relative distance to the ground of = 1.5, with H being the distance of the gravity centre to the bottom ground and L being the panel’s chord length. The results computed from the high-fidelity large-eddy simulations show that the proximity of the panel to the bottom surface impacted the vortex shedding and aerodynamic forces compared to classic flat plate studies, and its extreme slenderness did not generate a significant disruption in the leading edge flow as in the case of bridge decks.
Three regimes in the unsteady wake were observed: at low tilt angles, the flow is mostly attached over the suction side; in the range of to , more energetic and coherent vortical structures started to form a recirculating bubble; and, at larger tilt angles, the flow is fully separated and the panel’s wake is similar to that of a bluff body. These wake patterns varied slightly depending on the tilt angle being positive or negative due to the proximity to the ground. At positive tilt angles, the leading edge points upwards and the suction side of the panel allowed an unconstrained development and shedding of the vortical structures. Conversely, at negative tilt angles, the pressure side faced the bottom surface directly impacting the vortex generation and shedding, and the low position of the leading edge generated a diffuser-like effect.
The analysis of the body forces revealed that, at positive tilt angles, both vertical and horizontal forces were larger compared to their negative counterparts, except for the moment coefficient in which the maximum values were similar between positive and negative angles. The Strouhal number based on the peak frequency, obtained from the power spectral density of the moment coefficient and vertical velocity fluctuation at a point in the wake, attained larger values for negative angles due to the diffuser effect that enhanced vortex shedding. The range of Strouhal numbers varied with a tilt angle with an asymmetric "M" shape, reaching maxima of 0.30 and 0.38 at and , respectively. Varying the ratio, i.e., the relative distance of the structure to the bottom surface, changes the Strouhal number values for negative tilt angles between and , whilst this mostly remains the same for positive inclination angles.