Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Contents lists available at ScienceDirect

Journal of Wind Engineering & Industrial Aerodynamics

journal homepage: www.elsevier.com/locate/jweia

Computational evaluation of wind loads on sun-tracking ground-mounted

photovoltaic panel arrays
Giovanni Paolo Reina, Giuliano De Stefano *
Dipartimento di Ingegneria Industriale e dell’Informazione, Universit
a della Campania, Aversa, Italy

A R T I C L E I N F O A B S T R A C T

Keywords: Computational fluid dynamics is employed to evaluate the mean wind loads on sun-tracking ground-mounted
Computational fluid dynamics (CFD) solar photovoltaic panel arrays. Reynolds-averaged Navier-Stokes simulations are performed using a finite
Solar panel array volume-based numerical method. The mean turbulent flow around the panel is simulated by following two
Wind loading
different approaches, by considering either the full three-dimensional model or the reduced model with periodic
Unsteady RANS
boundary conditions in the spanwise homogeneous direction. The periodic model is demonstrated to well
reproduce the aerodynamics of the panel with a considerable reduction of the computational cost. The wind
loading history due to the continuous rotation of the system is directly simulated by means of the dynamic
meshing technique, which allows for a further savings of computational resources with respect to static
calculations.

1. Introduction demonstrated that it is possible to increase the overall daily output power
gain by more than 20% compared to a fixed mounted system, e.g.
In recent years, solar power has been strongly emerging as the first (Al-Mohamad, 2004).
source of alternative energy for industrial applications. Besides the One of the most important goals of the engineering research on PV
relative simplicity of the electricity production process and the possibility systems, either fixed or mobile, is the analysis of the aerodynamic loads
of reducing the effect on global warming and air pollution, the main acting on the solar panels, and indirectly on the support structures. A
reason for such a trend lies in the continuously decreasing price of solar large number of works on roof-mounted panels exist in the literature
photovoltaic (PV) panels, whereas the cost of traditional energy sources (Wood et al., 2001; Chung et al., 2008; 2011; Kopp and Banks, 2013;
has been increasing. However, the productivity of PV cells still represents Pratt and Kopp, 2013; Banks, 2013; Browne et al., 2013; Cao et al., 2013;
a challenging issue, since the maximum efficiency reached by the large- Aly and Bitsuamlak, 2014; Stathopoulos et al., 2014; Warsido et al.,
area commercial cells is about 24% (Blakers et al., 2013). Furthermore, 2014), whereas studies on ground-mounted panels have been appearing
the highest performance of the PV system is achieved when the solar rays only recently.
are approximately perpendicular to the given panel surface, and thus From the experimental point of view, the study of ground-mounted
only for few hours per day. systems is very complex since the characteristic flow spatial scales that
In the industrial field, a widespread solution to the problem of low are involved require the realization of very small models for boundary
efficiency is the use of solar farms, which are large-scale PV power sta- layer wind tunnel tests. Nevertheless, among the others, important
tions that are capable of generating large quantities of electricity. In these experimental findings were obtained by Stathopoulos et al. (2014), who
plants, even thousands of ground-mounted solar panels are connected to studied the aerodynamic loads on a stand-alone PV panel for different
form systems of length equal to tens of meters and beyond. As to the configurations and inclinations. In particular, the maximum peak of the
panel surface orientation, recent studies have produced some smart so- pressure coefficient for a ground-mounted panel was demonstrated for an


lutions to this problem, such as the application of solar trackers, e.g. inclination angle of 135 with respect to the wind direction. The pressure
(Ghosh et al., 2010; Rohr et al., 2015). The latter are devices that allow distribution on the upper and the lower surfaces of a single panel for four
the solar panels to rotate around one or two axes, while following the different wind directions was investigated by Abiola-Ogedengbe et al.
different inclination of the solar rays during the day, so as to extend the (2015), while the effect of lateral and longitudinal spacing between
temporal period of maximum efficiency of the PV cells. This way, it was

* Corresponding author.
E-mail addresses: giovannipaolo.reina@unicampania.it (G.P. Reina), giuliano.destefano@unicampania.it (G. De Stefano).

http://dx.doi.org/10.1016/j.jweia.2017.09.002
Received 2 May 2017; Received in revised form 2 August 2017; Accepted 2 September 2017
Available online 15 September 2017
0167-6105/© 2017 Elsevier Ltd. All rights reserved.
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

panels on the wind loading of a ground-mounted solar array was studied introduced. The geometry of the solar PV panel system under investi-
by Warsido et al. (2014). gation and the numerical settings are presented, together with the free-
Also due to the difficulties encountered in wind tunnel testing, the stream turbulence conditions that are assumed.
Computational Fluid Dynamics (CFD) analysis has recently become an
important predictive tool for ground-mounted PV systems. For instance,
Aly and Bitsuamlak (2013) investigated the sensitivity of wind loads on 2.1. Computational settings
the geometric scale of a stand-alone panel, by demonstrating that the
mean loads are not significantly affected by the model size, while the The geometric model examined in this work corresponds to a ground-
peak loads are indeed varying. Shademan et al. (2014) performed mounted solar panel array, on which a single horizontal axis tracker is
three-dimensional Reynolds-Averaged Navier-Stokes (RANS) simulations installed. The sub-panel geometry is similar to that one described in
to determine the effect of lateral spacing between sub-panels on wind Jubayer and Hangan, 2014, which represents the reference study that is
loads in a PV system, as well as the influence of longitudinal spacing used for validating the present approach, as discussed in the following
between solar panels in the arrayed configuration. Unsteady RANS section. In fact, the present PV system, which is typical of solar farms,
models were also employed by Jubayer and Hangan (2014, 2016) to consists of 36 panels arranged in portrait orientation to form a single row.
investigate the wind effects on a stand-alone panel as well as a solar panel The numerical simulation of more complex systems with multiple rows of
array immersed in the Atmospheric Boundary Layer (ABL) flow, with a panels will be the subject of our future work. Each panel has the chord
variable wind direction. length C ¼ 2 m, the width of 1.2 m and the thickness of 0.007 m. It is
The main goal of this work is the computational evaluation of the worth noting that the low thickness of the present model does not
mean wind loading on sun-tracking ground-mounted PV panel arrays, represent the dimension of the frame of real panels, which is here not
which are commonly used in solar farms. To fulfill the aim, two different simulated. In fact, including the frame would cause a significant
numerical approaches are followed. First, the entire row is directly complication for the mesh generation, which is out of the scope of the
simulated in a fully three-dimensional computational domain. The present analysis. The transverse length of the whole system is L ¼ 43.2 m.
alternative approach consists in the simulation of a small portion of the The existing gap between the single panels is actually neglected in the
solar panel row by imposing periodic boundary conditions in the span- realization of the CFD model, because its effect upon the aerodynamics of
wise homogeneous direction. The present RANS solutions are examined the PV system is considered to be negligible, as suggested by Wu et al.,
by correlating the wind flow field past the obstacle with the panel surface 2010. Therefore, the row of panels is modeled as a single plate with an
pressure distributions. In addition, the dynamic mesh technique is aspect ratio of L/C ¼ 21.6. The PV structure is supported by six vertical
applied to determine the history of the wind loading on the panel columns that are equally spaced by 7.2 m. Each column is modeled as a
structures due to the action of a single axis sun-tracking device. The bar with the height of 1 m and the square cross-section of side length
numerical simulations are performed using the commercial solver ANSYS 0.1 m.
Fluent, which is commonly and successfully employed in the industrial The solar panel array is placed within a computational domain whose
aerodynamics research, e.g. Liu et al., 2011; Moonen and Carmeliet, spatial dimensions are the following: 21.6 C (longitudinal length) along
2012; Tan et al., 2016. the x-axis, which is aligned with the wind direction, 6.3 C (vertical
The remainder of this paper is organized as follows. In x2, the overall length) along the y-axis, and 32 C (spanwise length) along the z-axis. The
computational methodology is described, with a particular focus on the upstream boundary is placed at a distance of 5 C from the obstacle, the
inlet turbulent flow boundary conditions. After the preliminary valida- downstream boundary at 16 C, the top boundary at 5 C and the lateral
tion of the proposed approach for a stand-alone panel, the results of the boundaries at 10.5 C. The shortness of the gap between the upstream
present numerical experiments for a sun-tracking ground-mounted panel boundary and the panel array reflects the fact that suitable profiles for the
array are presented and critically discussed in x3. Concluding remarks are velocity and the turbulence intensity of the upcoming wind are explicitly
provided in x4. imposed at the inlet boundary, as discussed in the following. This way, a
substantial savings of computational resources can be yielded. In fact, the
2. Computational modeling approach computational domain is created also following the practical guidelines
provided by Franke et al. (2007). In Fig. 1, the domain is depicted for a


In this section, the overall computational modeling approach is given angle of inclination of the panel, which is θ ¼ 60 , with respect to
the oncoming wind direction. However, a number of three-dimensional

Fig. 1. Computational domain for the PV panel system: full model.

284
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

RANS calculations are conducted for various panel inclinations, in order

to simulate the panel rotation around the horizontal axis by means of the
sun tracker, as discussed in the next section.
Based on the above geometric model, a reduced computational
domain with a spanwise extent of 3.6 C is obtained by cutting the original
domain with two vertical planes at equal distances from a supporting
column. The reduced model is shown in Fig. 2, where the geometric scale
is enlarged with respect to the previous Fig. 1, for the sake of clarity. For
the reduced model, due to the high aspect ratio of the PV system under
study, periodic boundary conditions are imposed to the flow in the ho-
mogeneous spanwise direction.
Concerning the discretization of the full computational domain, an
unstructured mesh is generated throughout the fluid region with a
minimum cell size of 105 m, except in the wall regions that are close to
the panel, the columns and the ground. In fact, the wall boundaries are
Fig. 4. Detail of the inflation on the surface of the panel.
inflated with 20 layers of clustered hexahedral cells in order to ensure
that the yþ value for the first grid point close to the walls is everywhere of
unitary order, as empirically verified. That allows for the proper nu-
merical resolution of the boundary layer regions. The overall numerical in order to impose the periodic boundary flow conditions, the two lateral
grid involves about 12  106 computational cells and the maximum sides of the computational domain are directly matched. Specifically, the
skewness factor is 0.84. For example, in Fig. 3, the global view of the two surface meshes on the periodic faces shown in Fig. 2 are conforming


computational mesh generated for the full model at θ ¼ 60 is shown, and the nodes on one side are coincident with those of the other side.
while the inflation on the PV panel is highlighted in Fig. 4. The same The present numerical meshes are obtained as a result of a grid
mesh settings are used for the reduced periodic model, which results in sensitivity study that is performed involving different grids with three
utilizing a hybrid computational grid of about 2  106 cells. In this case, different minimum cell sizes that are 104, 105 and 106 m, respec-
tively. By controlling the aerodynamic force coefficients, the second
resolution is chosen as the fair compromise between the accuracy of the
solution and the computational cost.

2.2. Freestream turbulence modeling

The turbulent incompressible flow around the PV system described in

the previous section is simulated by using the three-dimensional RANS
modeling approach, supplied with the Shear Stress Transport kω tur-
bulence model (Menter, 1994), which is commonly exploited in
computational wind engineering, e.g. Jubayer and Hangan, 2016.
The upcoming wind flow, which is given and known, corresponds to
the atmospheric turbulence model for open terrain provided by the En-
gineering Science Data Unit (ESDU) 1982, 1983. By assuming the wind
speed at the height of 10 m from the ground as the reference velocity,
which is Vref ¼ 26 m/s, and the panel chord C as the reference length, the
flow Reynolds number is Re ¼ VrefC/ν ¼ 3.56  106, where ν stands for
the kinematic viscosity of the air. The oncoming wind is simulated by
Fig. 2. Computational domain for the PV panel system: reduced periodic model. imposing the following logarithmic law boundary layer velocity profile at

Fig. 3. Computational grid for the PV panel system: full model.

285
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

the inlet of the computational domain


u* y þ y0
u¼ ln (1)
κ y0

where u is the mean streamwise velocity and y the distance form the
ground, u* ¼ 1.822 m/s is the frictional velocity, κ ¼ 0.41 is the von
Karman constant, and y0 is the aerodynamic roughness length, which is
0.03 m for open terrain. The present mean wind velocity profile at the
inlet boundary is shown in Fig. 5, compared with the theoretical profile
derived from ESDU 1982.
As suggested by Yang et al. (2009), in order to satisfy the requirement
of modeling an equilibrium ABL flow at the inlet boundary, the wind
velocity profile must be coupled with suitable profiles of resolved tur-
bulence variables, which are the turbulent kinetic energy and the specific
turbulence dissipation rate for the present SST k–ω model. Therefore, the
following inlet profiles of these two variables are considered
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi

 ffi
u2* y þ y0
k ¼ pffiffiffiffiffi
ffi C1 ln þ C2 (2)
Cμ y0 Fig. 6. Wind turbulence intensity profile at the inlet boundary.

u* 1 3. Results
ω ¼ pffiffiffiffiffiffi ; (3)
κ Cμ y þ y0
As common practice, the wind loadings on the PV panel are expressed
where C1 ¼ 0.025 and C2 ¼ 0.41 are two constants (experimentally in terms of the drag and the lift force coefficients, namely, CD ¼ D/(qref
determined) that describe the inflow turbulence level, and Cμ ¼ 0.04. The Aref) and CL ¼ L/(qref Aref), where D and L are the mean aerodynamic
outlet boundary of the computational domain is set as pressure-outlet, force components along the x-axis and the y-axis directions, respectively,
while the lateral sides and the upper boundary are imposed to obey a while qref ¼ 12 ρV 2ref represents the dynamic pressure of the air stream and
symmetry condition. As regards the solid walls, the panel and the col- Aref corresponds to the surface area of the panel. In addition, the mean
umns surfaces are all set as no-slip smooth walls, while the bottom of the pitching moment M about the z-axis, evaluated with respect to the center
domain, which represents the terrain surface, is modeled as no-slip rough of the panel, is considered by defining the moment coefficient, namely,
wall. In order to correctly simulate the ground by maintaining the desired CM ¼ M/(qrefArefC), where C stands for the panel chord length. The
equilibrium ABL flow condition, the roughness height is set as ks ¼ 4.9  surface pressure distribution on the panel is expressed in terms of the
106 m. This value corresponds to a roughness constant Cs ¼ 6  104, in mean pressure coefficient, that is Cp ¼ (p  pref)/qref, where pref repre-
accordance with the relationship ks ¼ Ey0/Cs (Blocken et al., 2007), sents the ambient atmospheric pressure.
where E ¼ 9.793 is an integration constant and y0 is the aerodynamic In the following, before showing the results of the current study for a
roughness length already used. The present wind turbulence intensity sun-tracking PV system, the preliminary validation of the overall
(TI) profile at the inlet boundary is shown in Fig. 6, compared with the computational modeling approach is presented for a fixed ground-
theoretical profile derived from ESDU 1983. mounted stand-alone panel. Specifically, the case recently studied by
For the reduced periodic case, the boundary conditions and the nu- Jubayer and Hangan (2014) as well as by Abiola-Ogedengbe et al.
merical settings are practically the same as for the full model. The only (2015), which we try to replicate, is considered.
very important difference is in the lateral sides of the domain, where
periodic boundary conditions are imposed. All the numerical simulations
3.1. Preliminary validation
are performed by using the commercial CFD code ANSYS Fluent R16.2.

The three-dimensional RANS simulation of the flow around a ground-

mounted solar panel is conducted for the geometric model illustrated in
Fig. 7, which is directly taken from Jubayer and Hangan, 2014. The solar
PV panel, which has the chord C ¼ 2.48 m and the width of 7.29 m, is


inclined by θ ¼ 25 with respect to the wind direction, with a distance
of 1.65 m between the trailing edge and the ground. The resulting dis-
tributions of the mean pressure coefficient along the midline of the upper
and the lower surfaces of the PV panel are depicted in Fig. 8. Here, we

Fig. 5. Mean wind velocity profile at the inlet boundary. Fig. 7. Preliminary validation: geometric model (Jubayer and Hangan, 2014).

286
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Fig. 8. Mean pressure coefficient for a stand-alone panel: present results compared to
reference experiments (square symbols) and numerical data (dashed line).

compare the mean surface pressures from our simulation against the
numerical data in Jubayer and Hangan, 2014 and the wind tunnel
measurements in Abiola-Ogedengbe et al., 2015. The present results are
in good agreement with the experiments for both the surfaces of the
panel, differently from the numerical study of Jubayer and Hangan,
2014, which provides some discrepancies for the lower surface when
compared to the experiments. In particular, as expected, the maximum
wind load takes place at the leading edge region.
When making a comparison in terms of the overall mean aerodynamic
coefficients, which are presented in Table 1, the difference with reference
numerical data appears fully acceptable so that the proposed computa-
tional modeling approach can be considered fully validated.

3.2. Full versus periodic model

In this section, the results of the two different computational

modeling procedures for evaluating the wind loads on the sun-tracking Fig. 9. Mean pressure coefficient on the panel section for the two different computational



approaches at the inclination angles (a) θ ¼ 25 and (b) 25 .
ground-mounted PV system described in x2.1 are presented. As already
mentioned, the ground-mounted solar panel array under examination is

equipped with a sun-tracking device, which allows a complete rotation of discrepancies exist for θ ¼ ∓60 . This is due to the fact that the periodic


120 about a fixed horizontal axis. Practically, the panel rotates from the model unavoidably alters the flow evolution in the z-axis direction. In


initial inclination of θ ¼ 60 , with respect to the oncoming flow di- particular, the existence of the vortices developing at the ends of a real

panel, which are more important for high tilt angles, causes a reduction
rection, up to θ ¼ 60 . Five different numerical calculations are con-



of the surface pressures that cannot be reproduced by the periodic model.
ducted at different angles of inclination, namely, θ ¼ ∓60 , ∓25 and 0 .
This difference between the two solutions also affects the mean aero-
Regarding the effect of the supporting columns on the panel aero-
dynamic coefficients, which are slightly overestimated in the periodic
dynamics, as empirically verified, the computations performed without
case, as illustrated in the following.
the supporting bars substantially give the same results of those ones with
The mean aerodynamic coefficients CD, CL and CM are reported
the bars included, which are actually considered in the following dis-
against the panel inclination θ in Figs. 11, 12 and 13, respectively, for the
cussion. For the full modeling approach, each CFD calculation takes
two different approaches that are proposed. As expected, the angle
about 12 h by employing a parallel workstation supplied with two Intel
variation causes a change of sign for CL and CM, while CD shows a min-
Xeon 2.2 GHz processors, whereas the periodic modeling approach is
imum when the panel is aligned with the wind direction. The force co-
about eight times faster.
efficients appear substantially symmetric with respect to zero angle of
The mean pressure distributions for the panel section corresponding
inclination, which demonstrates the negligibility of the ground effect for
to the centerline of a supporting bar are drawn in Figs. 9 and 10, for the


the present geometry. The mean load coefficients have been verified to
inclination angles θ ¼ ∓25 and θ ¼ ∓60 , respectively. Apparently, for
be negligible at zero tilt angles. The corresponding data points are not
larger angles, the pressure loading on the panel is higher and more
reported in the above Figure.
uniform for both the ventral and the dorsal surface. When making a
In order to examine the local distribution of the surface pressure on
comparison between the two different computational models, the pres-

the panel, in Fig. 14, the contours of the mean pressure coefficient on the
sure distribution appears practically the same for θ ¼ ∓25 , while some
dorsal and the ventral surfaces of the plate are shown for a given tilt


angle, which is θ ¼ 60 . For the present case of straight wind, the nu-
Table 1
merical solution is symmetric about the spanwise midplane of the PV
Mean aerodynamic loads for a stand-alone panel.
system. From inspection of these pictures, it qualitatively appears that the
CD CL CM
periodic model locally matches quite well with the full model, for which
Present 0.56 1.20 0.11 only half the geometry is actually reported, at the midplane zone. The
Jubayer and Hangan (Jubayer and Hangan, 2014) 0.57 1.24 – highest pressure values are achieved on the dorsal (upper) surface, close
Error 1.7% 3.2% –
to the lower edge (Y/C ¼ 0.5) that corresponds to the leading edge of

287
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Fig. 12. Mean lift force coefficient versus the panel inclination for the two
different models.

Fig. 10. Mean pressure coefficient on the panel section for the two different computa-



tional approaches at the inclination angles (a) θ ¼ 60 and (b) 60 .
Fig. 13. Mean pitching moment coefficient versus the panel inclination for the two
different models.

1995). The method has been extensively and successfully utilized for the
numerical visualization of separated flows around obstacles, e.g. (De
Stefano and Vasilyev, 2014). By looking at this figure, where the
iso-surface of Q ¼ 30 s2 is reported, the development of an intense
vortex at the end of the panel is apparently well captured by the full
model, while this is impossible for the reduced one.
The RANS solution obtained with the periodic model is further
examined by considering the mean flow in the plane Z/C ¼ 0, which
corresponds to the middle section of a supporting column. The spanwise
vorticity contours are reported together with the streamlines in Figs. 16
and 17, for negative and positive inclinations, respectively. For high


negative tilt angles, as happens for θ ¼ 60 , the wind flow hits the
upper surface of the panel in the proximity of the leading edge, where it
separates. The downward separated fluid flow arrives on the downwind
surface by creating a very intense vortex in the area between the leading
edge and the support column, as recognized by looking at the vorticity
Fig. 11. Mean drag force coefficient versus the panel inclination for the two
different models. contours. On the other side of the panel, the upward flow reaches the
trailing edge, where it separates, while creating an extended recircula-
tion zone on the lower surface of the panel, downstream of the support


the panel at this tilt angle. For the full model, the surface pressure tends column. For low negative tilt angles, as happens for θ ¼ 25 , a similar
to decrease at the lateral zones, owing to the three-dimensionality of the flow configuration occurs, with the increased extension of the vortex
flow. The different behavior of the two models is evident when consid- close to the leading edge, while the opposite holds for the recirculation
ering the flow visualization reported in Fig. 15, where the main vortical region generated at the trailing edge of the panel. By reversing the
structure in the near wake of the obstacle is shown, for the same panel inclination of the plate, the flow configuration appears quite different, as


position. Here, the vortex detection is carried out by applying the Q- illustrated in Fig. 17. For low positive tilt angles, as happens for θ ¼ 25 ,
criterion, which identifies vortices as connected regions with a positive a single recirculation region is created on the upper surface of the panel
second invariant of the velocity gradient tensor (Jeong and Hussain, that extends down to the trailing edge. On the contrary, for high positive

288
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Fig. 14. Contours of the mean pressure coefficient on the upper and the lower surfaces of the panel at θ ¼ 60 for the two different models.

tilt angles, as happens for θ ¼ 60 , the flow past the solar panel generates determined from the analysis of the fluid flow around a smaller part,
two very large counter rotating vortices. The highest vorticity zones in while assuming periodic boundary conditions in the homogeneous
the flow field correspond to the shear layers developing at the spanwise direction, regardless of the tilt angle. The main consequence is a
panel edges. considerable reduction of the computational cost of the CFD analysis, in
When looking at the mean loadings, the present results demonstrate terms of both calculation time and memory. However, it is worth noting
that the aerodynamic behavior of a single row of arrayed panels can be that the present periodic method can be only applied for straight winds.

289
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Fig. 15. Main vortical structure in the near wake of the panel at θ ¼ 60 , identified by
the Q-criterion, for (a) the full and (b) the periodic models.

Corner vortices at oblique wind directions, which usually produce the

maximum force on the panel, could only be estimated by modeling the
full panel array.

3.3. Dynamic approach

As discussed above, in order to increase its efficiency, the PV panel

array under study is equipped with a solar tracking system, which allows


the complete rotation of 120 around the z-axis. When evaluating the
history of the aerodynamic loads during the period of rotation, the dy-
namic meshing technique can be effectively exploited. This methodology
allows to obtain an accurate time-dependent numerical solution by
continuously varying the computational mesh during the CFD calcula- Fig. 16. Streamlines and spanwise vorticity contours in the midplane for two different
tion, consistently with the changing position of the moving obstacle. To


negative inclinations that are (a) θ ¼ 60 and (b) 25 .
prevent the deterioration of the mesh quality and/or the degeneration of
existing cells, due to the geometry continuous modification, two different
methods can be used, which are referred to as smoothing and remeshing.
The smoothing technique consists in moving the interior nodes of the different flow configurations by means of a single relatively short
computational grid without changing their number and connectivity. The calculation.
remeshing technique allows for the local update of the numerical mesh The time step size for the dynamic simulation is chosen equal to
by either adding or deleting cells, where the boundary displacement 103 s, in order to maintain a Courant-Friedrichs-Lewy number of unitary
would be otherwise locally too large with respect to the mesh size. In this order and, thus, achieve an accurate and time consistent numerical so-
study, after some preliminary testing of the above methods for the pre- lution. The whole calculation, which is made by employing the same
sent geometry, both remeshing and diffusion-based smoothing proced- parallel workstation utilized for the previous cases, takes about 24 h.
ures are actually utilized. Even considering the time needed for the mesh generation and the pre-
In order to assess the capabilities of the dynamic approach for the processing stage, a substantial saving of computational time with
present case, unsteady RANS simulations are performed for the two- respect to the equivalent number of static calculations is basi-
dimensional model corresponding to the cross-section of the panel, by cally obtained.


The solution provided by the dynamic approach is compared to those
imposing a complete rotation, from θ ¼ 60 to 60 , in 10 s. Usually, the

ones obtained with the reduced periodic model, without considering the
rotational speed of a solar tracker is lower, namely, about 15 per hour, to
supporting column, for four different tilt angles. When looking at the
transfer from a tilted to stow position during a high wind event. In some
continuous histories of the mean aerodynamic coefficients, which are
particular circumstances, the panel can be moved at much higher rota-


reported in Figs. 18, 19 and 20, the comparison with the different static
tional speeds, ranging between 5 and 18 per second. However, even if
solutions is quite successful. As expected, the mean drag force achieves
not fully representative of typical operations of a single axis solar array,
its minimum value when the panel gets the position parallel to the wind
the speed is also chosen to make it possible to simulate a number of

290
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Fig. 19. Mean lift force coefficient for the dynamic approach, compared with static
simulations at different inclinations.

Fig. 20. Mean pitching moment coefficient for the dynamic approach, compared with
static simulations at different inclinations.

direction so that the impact surface is limited. The mean lift force and the
mean pitching moment show a complex behavior, which reflects the
different flow configurations that occur for the various panel inclinations,
Fig. 17. Streamlines and spanwise vorticity contours in the midplane for two different
as discussed in the previous section. A further analysis of the dynamic



positive inclinations that are (a) θ ¼ 25 and (b) 60 . solution is presented in Figs. 21 and 22, where the mean pressure coef-
ficient distribution is compared to that one obtained with the different
static calculations on the panel mid section. Here, the same four different
angles of inclination contemplated above are considered, namely, θ ¼



∓60 and θ ¼ ∓25 . The comparison between the static and the dynamic
approaches is quite successful. However, a marked difference exists at the


tilt angle θ ¼ 25 for the lower surface of the panel, where the mean
flow around the obstacle seems not very well reproduced by the dynamic
solution. This is likely to be caused by the deterioration of the instanta-
neous numerical mesh, which could be avoided by more stringently
enforce the parameters that control the mesh quality, like the skewness
factor and the cell aspect ratio. This issue however requires further
detailed investigation.

3.4. Discussion

Since the supporting structures are indirectly affected by the pressure

forces acting on the panel array, also due to the large exposed surface
area, the knowledge of the wind loadings is of great importance for the
design of the PV system. In particular, the net pressure coefficient
Fig. 18. Mean drag force coefficient for the dynamic approach, compared with static resulting from the combination of the pressure coefficients at the upper
simulations at different inclinations. and the lower surfaces of the panel has to be considered and the panel

291
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

Fig. 22. Mean pressure coefficient on the panel section for the dynamic and the static
Fig. 21. Mean pressure coefficient on the panel section for the dynamic and the static




approaches at (a) θ ¼ 25 and (b) 25 .
approaches at (a) θ ¼ 60 and (b) 60 .

inclinations associated with the maximum aerodynamic loads have to be 4. Conclusions

identified. The values of the wind loads for the PV panel array are
however different from those for the support structures, due to the In this work, computational fluid dynamics is employed for the
different tributary areas. evaluation of wind loads on sun-tracking ground-mounted photovoltaic
It is worth emphasizing that the present work based on the RANS panels under atmospheric boundary layer flow. The main goal is the
approach, like similar studies on the subject, (e.g. Jubayer and Hangan, numerical analysis of the mean turbulent flow around panel arrays,
2016), focuses on the analysis of the mean aerodynamic loads, whose which are commonly used in solar farms. Two different modeling ap-
knowledge would not be yet sufficient for the design process. In fact, for proaches are considered, which consist in simulating either the whole
the practical design of a sun-tracking ground-mounted PV system, not three-dimensional system or a reduced portion, while assuming periodic
only the static peak loads should be determined but, in addition, the boundary conditions in the spanwise direction. The flow is simulated for
torsional instability issue should be addressed, e.g. (Rohr et al., 2015). various panel inclinations with respect to the wind direction by using the
The use of more reliable numerical techniques such as Large Eddy Reynolds-Averaged Navier-Stokes approach, supplied with the Shear
Simulation (LES) and Detached Eddy Simulation (DES) would provide Stress Transport k–ω turbulence model.
more detailed solutions, where important flow characteristics such as For straight wind, the reduced periodic model provides results that
flow separation and vortex shedding would be simulated, e.g. (Breuer are fully acceptable when compared to the fully three-dimensional so-
et al., 2003). This way, the numerical results would be improved and a lution, with a great advantage in terms of computational time and power.
better insight for the design process would be achieved. However, given Also, the dynamic meshing technique is applied to determine the history
the actual flow Reynolds-number, the application of the LES method to of the aerodynamic loads on the panel due to the action of a horizontal
the present PV system is practically unfeasible, because of the enormous axis sun-tracking device. The flow around the cross-section of the panel is
number of grid points that would be necessary to yield an accurate so- simulated for a period of time corresponding to its complete rotation. The
lution. On the contrary, the less expensive DES method seems promising resulting mean aerodynamic loads are in good agreement with the results
for computational wind engineering and its application to the present of static simulations at given tilt angles, with a further reduction of the
case is the subject of our ongoing research. overall computational costs.
Finally, the numerical results should be viewed in conjunction with The use of more sophisticated numerical methods such as detached-
wind tunnel testing in order to improve the understanding of fluid- eddy simulation is the subject of our ongoing research on computa-
structure interaction and reduce the costs of research and develop- tional wind engineering for the prediction of ground-mounted solar panel
ment, which is particularly sensible for the wind engineering applica- aerodynamics.
tions, e.g. (Meroney, 2016).

292
G.P. Reina, G. De Stefano Journal of Wind Engineering & Industrial Aerodynamics 170 (2017) 283–293

References Jeong, J., Hussain, F., 1995. On the identification of a vortex. J. Fluid Mech. 285, 69–94.
Jubayer, C.M., Hangan, H., 2014. Numerical simulation of wind effects on a stand-alone
ground mounted photovoltaic (PV) system. J. Wind Eng. Ind. Aerodyn. 134, 56–64.
Abiola-Ogedengbe, A., Hangan, H., Siddiqui, K., 2015. Experimental investigation of wind
Jubayer, C.M., Hangan, H., 2016. A numerical approach to the investigation of wind
effects on a standalone photovoltaic (pv) module. Renew. Energy 78, 657–665.
loading on an array of ground mounted solar photovoltaic (PV) panels. J. Wind Eng.
Al-Mohamad, A., 2004. Efficiency improvements of photo-voltaic panels using a sun-
Ind. Aerodyn. 153, 60–70.
tracking system. Appl. Energ 79, 345–354.
Kopp, G.A., Banks, D., 2013. Use of wind tunnel test method for obtaining design wind
Aly, A.M., Bitsuamlak, G., 2013. Aerodynamics of ground-mounted solar panels: test
loads on roof-mounted solar arrays. J. Struct. Eng. 139, 284–287.
model scale effects. J. Wind Eng. Ind. Aerodyn. 123, 250–260.
Liu, B., Qu, J., Zhang, W., Qian, G., 2011. Numerical simulation of wind flow over
Aly, A.M., Bitsuamlak, G., 2014. Wind induced pressures on solar panels mounted on
transverse and pyramid dunes. J. Wind Eng. Ind. Aerodyn. 99, 879–888.
residential homes. J. Archit. Eng. 20, 1–12.
Menter, F., 1994. Two-equation eddy-viscosity turbulence models for engineering
Banks, D., 2013. The role of corner vortices in dictating peak wind loads on tilted flat
applications. Am. Inst. Aeronaut. Astronaut. J. 32, 1598–1605.
solar panels mounted on large, flat roofs. J. Wind Eng. Ind. Aerodyn. 123, 192–201.
Meroney, R., 2016. Ten questions concerning hybrid computational/physical model
Blakers, A., Zin, N., McIntosh, K., Fong, K., 2013. High efficiency silicon solar cells.
simulation of wind flow in the built environment. Build. Environ. 96, 12–21.
Energy Procedia 33, 1–10.
Moonen, P., Dorer, V., J., Carmeliet, 2012. Effect of flow unsteadiness on the mean wind
Blocken, B., Stathopoulos, T., Carmeliet, J., 2007. CFD simulation of the atmospheric
flow pattern in an idealized urban environment. J. Wind Eng. Ind. Aerodyn. 104-106,
boundary layer: wall function problems. Atmos. Environ. 41, 238–252.
389–396.
Breuer, M., Jovicic, N., Mazaev, K., 2003. Comparison of DES, RANS and LES for
Pratt, R.N., Kopp, G.A., 2013. Velocity measurements around low-profile, tilted, solar
separated flow around a flat plate at high incidence. Int. J. Numer. Meth. Fluids 41,
arrays mounted on large flat-roofs, for wall normal wind directions. J. Wind Eng. Ind.
357–388.
Aerodyn. 123, 226–238.
Browne, M.T.L., Gibbons, M.P.M., Gamble, S., Galsworthy, J., 2013. Wind loading on
Rohr, C., Bourke, P., Banks, D., 2015. Torsional instability of single-axis solar tracking
tilted roof-top solar arrays: the parapet effect. J. Wind Eng. Ind. Aerodyn. 123,
systems. In: Proceedings of the 14th Int. Conf. Wind. Eng., June 21-26, pp. 1–7. Porto
202–213.
Alegre Brazil.
Cao, J., Yoshida, A., Saha, P., Tamura, Y., 2013. Wind loading characteristics of solar
Shademan, M., Barron, R.M., Balachandar, R., Hangan, H., 2014. Numerical simulation of
arrays mounted on flat roofs. J. Wind Eng. Ind. Aerodyn. 123, 214–225.
wind loading on ground-mounted solar panels at different flow configurations. Can.
Chung, K., Chang, K., Liu, Y., 2008. Reduction of wind uplift of a solar collector model.
J. Civ. Eng. 41, 728–738.
J. Wind Eng. Ind. Aerodyn. 96, 1294–1306.
Stathopoulos, T., Zisis, I., Xypnitou, E., 2014. Local and overall wind pressure and force
Chung, K., Chang, K., Chou, C., 2011. Wind loads on residential and large-scale solar
coefficients for solar panels. J. Wind Eng. Ind. Aerodyn. 125, 195–206.
collector models. J. Wind Eng. Ind. Aerodyn. 99, 59–64.
Tan, L., Zhang, W., Bian, K., An, Z., Zu, R., Qu, J., 2016. Numerical simulation of three-
De Stefano, G., Vasilyev, O.V., 2014. Wavelet-based adaptive simulations of three-
dimensional wind flow patterns over a star dune. J. Wind Eng. Ind. Aerodyn. 159,
dimensional flow past a square cylinder. J. Fluid Mech. 748, 433–456.
1–8.
ESDU, 1982. Strong winds in the atmospheric boundary layer. Part 1: mean hourly wind
Warsido, W.P., Bitsuamlak, G.T., Barata, J., 2014. Influence of spacing parameters on the
speeds. Eng. Sci. Data Unit Number 82026.
wind loading of solar array. J. Fluid Struct. 48, 295–315.
ESDU, 1983. Strong winds in the atmospheric boundary layer. Part 2: discrete gust
Wood, G., Denoon, R., Kwok, K., 2001. Wind loads on industrial solar panel arrays and
speeds. Eng. Sci. Data Unit Number 82045.
supporting roof structure. Wind. Struct. 4, 481–494.
Franke, J., Hellsten, A., Schlunzen, H., Carissimo, B., 2007. Best practice guideline for the
Wu, Z., Gong, B., Wang, Z., Li, Z., Zang, C., 2010. An experimental and numerical study of
CFD simulation of flows in the urban environment. Qual. Assur. Improv. Microscale
the gap effect on wind load on heliostat. Renew. Energy 35, 797–806.
Meteor. Models Cost Action 732.
Yang, Y., Gu, M., Chen, S.Q., Jin, X.Y., 2009. New inflow boundary conditions for
Ghosh, H., Bhowmik, N., Hussain, M., 2010. Determining seasonal optimum tilt angles,
modeling the neutral equilibrium atmospheric boundary layer in computational wind
solar radiations on variously oriented, single and double axis tracking surfaces at
engineering. J. Wind Eng. Ind. Aerodyn. 97, 88–95.
Dhaka. Renew. Energy 35, 1292–1297.

293