Numerical investigation of wind influences on photovoltaic arrays mounted on roof

The wind-induced response of photovoltaic (PV) panel installed on building roof is influenced by the turbulence induced by the pattern of both panels and roofs. Different roof types cause different flow patterns around PV panels, thus change the flow mechanism exerted on PV panels. In this study, the effects of roof types, heights and the PV array layouts on the net wind loads of the PV panel is investigated. The software Fluent is adopted and the three-dimensional Reynolds-averaged Navier–Stokes (RANS) method is used for numerical analysis. The rationality and accuracy of the numerical results obtained from the current study are verified through comparison with the results of wind tunnel experiments. The maximum wind uplift on the PV panels increases with the panel tilt angle for two types of roofs, but decreases with the increase of the PV array edge setback. Moreover, the maximum wind uplift also increases with the roof clearance within a certain range. Both the row spacing of PV array and the building height have a great effect on the maximum aerodynamic uplift for gable roofs, whereas both of them have slight influence for flat roofs.


Introduction
With the shortage and the pollution caused by traditional energy sources, the clean energy has been vigorously used in the world in recent years. Solar energy is abundant in China, which is utilized through installing the PV array power systems on buildings roofs. The PV arrays installed on flat and gable roofs are presented in Figure 1.
The flow field above the roof may be disturbed significantly by the PV array which makes the value of wind loads of the PV panels completely different compared to the bare roof. The interaction between the turbulence induced near the roof edge and the turbulent flow created by the PV array edge may result in a large uplift on the PV array and cause damage to the PV panels.
At present, both ground-mounted and roof-mounted PV array have been investigated to estimate wind pressure on PV panels. The wind pressure on the groundmounted PV panel is mainly affected by PV array parameters, while the roof-mounted PV panel is also affected by the building dimensions and the roof types. This study focuses on the PV array mounted on roof.
There are many proprietary studies concerning the effect of PV array parameters on the aerodynamic loads of the PV panel, but there are few investigations considering the effect of roof types. The shading effect resulted from the first row of PV arrays was studied by Radu et al. (1986) through the wind tunnel test. The negative net pressure coefficients of the PV panel were lower than those on the roof without PV panels mounted through wind pressure tests by Wood et al. (2001). The wind loads of the PV array were influenced significantly by the PV panel tilt angle and the PV array setback from the roof leading edge. The wind flow mechanism related to the wind loads of the roof-mounted PV array was researched by Kopp et al. (2012) taking into consideration of two panel tilt angles. A wind tunnel experiment conducted by Cao et al. (2013) evaluates the wind loads on PV panels located on a flat roof. They have pointed out that the turbulence generated by the PV panel edge became predominant as the PV panel tilt angle increased, and the wind uplift on the PV panels became large. The wind uplift also increased with the distance between the adjacent PV arrays. A wind tunnel experiment on PV panels was implemented by Aly and Bitsuamlak (2014). It was found that the wind pressure on the PV panel depends on the location of panels. Generally, the PV panels close to the roof corners were subjected to larger wind uplifts. Kopp (2014) carried out wind tunnel experiments to find out the influences of PV panel tilt angle and row spacing on the aerodynamic pressure of PV panels fixed to a flat roof. It was found that there was an obvious increase in the pressure coefficient only for PV panel tilt angles ranging from 2°t o 10°. Mou et al. (2017) investigated the influence of building geometric size on the mean aerodynamic pressure of square-shaped tall building by steady Reynoldsaveraged Navier-Stokes method. The Realizable k-ε turbulence model (Shih, Liou, Shabbir, Yang, & Zhu, 1995) was adopted and their numerical results were validated by wind tunnel experimental data. Dong et al. (2018) investigated the mechanism of the voice origin on the landing gear employing the Delayed Detached Eddy Simulation method and the Unsteady Reynolds-Averaged Navier-Stokes method using the software ANSYS-CFX, the Large Eddy Simulation (LES) was also applied. For the mean wind condition, the URANS and LES gave reasonable prediction and the grid density effect was weak.
For the flow fluctuations involving separation and attachment, the grid resolution had a significant effect on the predicted results and DDES had a great advantage in predicting pressure fluctuations compared to URANS and LES. The influence of building height on the wind uplifts of PV arrays was investigated by Ginger et al. (2011) for the flat roof configuration. The PV panel tilt angle was set at 30°and two sizes (2.7, 10 m) of building roof height were considered. The influence of the building height was found unimportant for the cases studied. Ferreira et al. (2018) studied the influences of panels on the wind load of the roof throug both wind tunnel tests and numerical method. While, the wind loads exerted on the panels were not analyzed. The parameters concerned in previous studies are summarized in Table 1, N/M presents that item was not mentioned in their paper.
The uplift on the PV panels is resulted from the interaction between the building-generated turbulence and the PV panels. Different roof types cause different types of flow pattern surrounding the PV panels, thus change flow mechanism of the maximum wind uplift on PV panels. However, none of the foregoing studies investigates the effect of roof types on the wind loads exerted on the panel. This paper is intended to address the following issues: (1) investigated the influence of panel tile angle, (2) investigated the roof clearance, (3) investigated the array edge setback. Based on CFD approach, the flow mechanisms that produces the maximum wind uplift on PV panels is analyzed for two different types of gable and flat roofs. To minimize the total mesh and resolve the near wall flow accurately, the nesting style mesh is used in the flow domain with prismatic cells created around the PV panels and building surfaces. The number of computational cells is set at the magnitude of 1.2 × 10 7 . The RANS method is employed using the commercial software Fluent. Finally, the reliability of the CFD simulation is verified by means of the comparison with the others' wind tunnel experiments.

Geometry details
The geometries of buildings and PV array were set up according to the tests dimension of Kopp et al. (2012) using a size ratio of 1:30. Figure 2 shows the geometric size details. Two roof types common for residential building were considered. The size of the panel in the model scale were 670 (Length) × 33 (Width) × 1.7 (Thickness) mm to be the same with those used in the tests of Kopp et al. (2012). The roof width B was set as 22.5 m in full scale. Four PV panel tilt angles (2°, 10°, 20°, and 30°) common for rooftop installations were constructed in ICEM. For the purpose of academic research, the PV panels on the gable roof are symmetrically arranged on both sides of the ridge in this study. The H is the building roof height. Table 2 provides the geometric details. The row spacing was set according to practical application to minimize the shading effect. The minimum value of PV panel clearance was chosen based on typical values used in practice. Both of these parameters were considered from characteristic values used in practice up to possible, for example, worst case or upper bounds according to Kopp (2014). A value of 1.2 m is typical for PV array edge setback in practice and two other values were considered ranging from the worst case to possible upper bounds. The lowest value of roof height was set to be same as the wind tunnel experiment by Kopp et al. (2012). The other two values were set as multiples of 7.3 m. The flow domain dimensions are presented in Figure 2(a) with an distance equals to 5H from inlet to building based on the suggestions (Franke, Hellsten, Schlunzen, & Carissimo, 2007). The considered blockage ratio was less than approximately 5% (Holmes, 2007).

Mesh details
The computational domain and mesh for the CFD simulation were created based on the recommendations proposed by the Working Group of the Architectural Institute of Japan (AIJ, 2004). The grids were set up in software of ICEM. The flow field was meshed by nesting style ( Figure 3). There were prismatic cells generated near the panels and building adopting a refinement ratio of 1.2. For the inner domain, there were tetrahedral cells generated surrounding the prismatic cells. For the outer domain, it was filled by hexahedral cells (Figure 3(a)). The interface (Figure 3(b)) between outer and inner domain was defined to send messages. The inner domain and the close-up grid near the PV panel are presented in Figure 3(c-e).

Model resolution investigation
To make sure the independent meshing was obtained, different grids number for 0°wind angle with M1 (8.75M cells), M2 (12.13M cells), and M3 (17.20M cells) were created. The PV array configuration is consistent with that used in the tests of Kopp et al. (2012). The surface pressure coefficients (C p = (P − P ∞ )/(0.5ρU 2 ∞ )) were calculated. The reference velocity U ∞ was measured at the height of roof in the undisturbed upstream flow. The net pressure coefficient C p was calculated by the difference of wind pressures between the PV panels upper and lower sides. The wind load predictions of three meshes were contrasted to experimental data (Kopp et al., 2012) (Figure 4). The largest negative net C p value of the PV array was listed to compare with the test for the design consideration. The results obtained from Grids M2 and M3 agree with each other well. The panels of row 7-12 were experiencing positive aerodynamic loads, and there exists deviations between the simulations and the experiments. This is similar with previous study (Tominaga, Mochida, Murakami, & Sawaki, 2008) because the flow unsteadiness (reattachment, etc.) caused by vortex shedding may not be predicted well by RANS method. Integrating the computational devices and the precision of simulation, the grid number M2 was selected for subsequent simulation. The minimum cell size was 0.15 mm. The number of y + = y √ ρτ w /μ (τ w denotes the wall shear stress, μ is the molecular viscosity, ρ refers to the air density, and y is the height of the cell nearest to the wall) is lower than 5. The software FLUENT was used to simulate the flow structure surrounding the roof-mounted PV arrays. The second order upwind scheme was used for the spatial discretisation of the momentum, turbulence kinetic energy, and specific dissipation. Convergence criteria of 10 −5 were employed for all variables. During the simulation, the static pressure monitored at a point upstream of the building reached steady state.

Boundary conditions
The upcoming speed U(z) and the wind turbulence intensity I(z) are exerted at the flow inlet according to Engineering Science Data Unit (ESDU, 1982(ESDU, , 1983 speed and turbulence intensity profiles. The equation of wind velocity is based on Equation (1) considering the height of the PV array is less than 30 m, the wind intensity is exerted at inlet profile base on Equation (2): Where: z represents the height from the ground. A terrain roughness parameter z 0 of 0.03 m is used which represents open terrain. The friction velocity u * is calculated following the relationship: The reference wind velocity U 10 is 26 m/s at a 10 m height from the ground. σ (z) is standard deviation of U at height z and can be expressed as: where: f is Coriolis parameter which equals 0.000115. The turbulence kinetic energy k(z) and turbulence dissipation rate ε(z) are resolved as follow: Where: κ represents Von Karman constant which equals 0.41 in this study. The User Defined Function (UDF) was employed to define the inlet conditions. The outlet was defined as outflow with a zero gauge pressure. Slip walls were applied at the sides of the domain and no slip walls were employed at the bottom of the domain, the PV panels and building surfaces. During the CFD simulation, considering that the RNG k-ε is a High-Reynolds number model, the Enhanced Wall Treatment was used. Since the SST k-ω is a Low-Reynolds number model, the near-wall flow is resolved directly.
Modeling of wind flow inside an empty flow domain was conducted to check horizontal homogeneity of the flow (Blocken, Stathopoulos, & Carmeliet, 2007). The wind speed and turbulence kinetic energy conditions of

RANS model validation
A RANS simulation is conducted using a steady solver and the resulted mean aerodynamic loads of the PV array are calculated. The second order upwind scheme is applied for discretising the convective terms, and the standard scheme is employed for the pressure interpolation. The semi-implicit method for the pressure-linked equations consistent (SIMPLEC) scheme was applied for the pressure-velocity coupling. Both the RNG k-ε and SST k-ω models are used. The RNG k-ε model is proposed by Yakhot and Smith (1992) and the equations corresponding to turbulence kinetic energy k and turbulence dissipation rate ε are presented by Equations (8) and (9): Where: ρ represents the air density, t represents the time, u i denotes the space time-averaged velocity. X i and X j are spatial coordinate components, μ is the wind kinematic viscosity, μ t refers to the kinematic eddy viscosity. G k represents the turbulence kinetic energy produced by the laminar velocity gradient. α k and α ε are turbulence Prandtl constants for equation k and equation ε respectively. C 1ε and C 2ε are model constants and considered to be 1.42 and 1.68 in this study. The SST k-ω model is proposed by Menter (1993). The transport formulas corresponding to turbulence kinetic energy k and turbulence dissipation rate ω are presented by Equations (10) and (11): Where: Y k and Y ω are the dssipations of k and ω, respectively, σ k and σ ω are the turbulence Prandtl constants Figure 6. Contrast of aerodynamic coefficients on panels obtained using RNG and SST models for different wind directions.
for k and ω. G ω is the generation of ω. D ω is the crossdiffusion term. The distance between the panels is fixed to 56 mm with an PV array tilt angle of 20°which is same with the existing experiments conducted by Kopp et al. (2012). For 0°wind direction, the coefficients of wind pressure predicted by the SST k-ω model match with the test results (Kopp et al., 2012) (Figure 6) well for the PV panels upstream. The net C p values for upstream rows were overestimated by the RNG. For the downstream PV panels, the deviation between both models was prominent caused by the insufficient of steady RANS simulations (Tominaga et al., 2008). For 180°wind angle, the predictions of SST k-ω model match with the experimental results rather well except for the rows 11 and 12. The net mean C p values calculated from the RNG k-ε deviated deviating from the experimental data obviously.
The mesh counts and analysis expense corresponding to two RANS models are shown in Table 3. The precision of the RANS predictions has positive correlation with correlation coefficient. SST k-ω model has better accuracy compared with RNG k-ε model and the computing consumption was acceptable, so the SST k-ω model was applied for further parameter analysis.

Parameter analysis
In this section, the influences of various PV array parameters are investigated considering a flat roof and a gable roof. The wind direction is set as 0°for all the following numerical models. The PV arrays were modeled using a range of configurations and arrangements, and the resultant aerodynamic load distribution of the PV arrays are analyzed in correlation to the flow fields around the PV arrays. The geometric details of the models and the analyzed parameters are presented in Table 2.

Effect of panel tilt angle
The flow structures surrounding the array fixed to flat roof for various panel tilt angles are presented in Figure 7. The turbulence induced by the roof edging is obvious for the 2°tilted PV array. However, the vortices resulted from panel edging becomes predominant for the 30°tilt angle PV array configuration. Increasing the PV panel tilt angle from 2°to 20°results in a significant increase in the largest uplifts on the PV array. However, this increase is not apparent as the PV panel tilt angle increases from 20°t o 30° (Figure 9(a)). The flow fields for the gable roof cases are shown in Figure 8. The local turbulence induced by the PV panel edging become predominant as the panel tilt angle increases from 2°to 30°. For the arrays having a tilt angle of 2°, 10°, and 20°, the main eddy caused by the ridge dominates the flow surrounding the PV arrays on the leeside of the roof. For the 30°tilted PV array, the ridge vortex induced by ridge is separated into small eddies by the panel edges.
The net wind load coefficients on the PV array for various tilt angles are presented in Figure 9. The PV panels at the windward side of the roof are mainly experiencing positive wind loads. However, the PV panels put on the roof leeside are mainly suffered from wind uplift. Row 7, which is closest to the ridge at the roof leeside, is subjected to the largest wind uplift. It can be concluded that the largest uplift of the PV array generally increases with increasing PV array tilt angle for both the roof types.

Effect of roof clearance
The flow fields for various model configurations are modeled to evaluate the effect of the roof clearance on the aerodynamic loads of the panels for both the flat and gable roofs. The panel tilt angle equals to 20°. The pressure contour is presented with the velocity streamlines in Figures 10 and 11. The localized vortices between the PV panels can be observed clearly. In the case of the flat roof, increasing the roof clearance amplifies the vortices caused by flow separation at panel edge. The largest negative net mean pressure coefficient decreases from −0.14 to −0.24 as the clearance increases from 0.09 m to 0.41 m (Figure 14(a)). The pressure coefficients on the upper and lower sides of row 1 are presented in Figure 12. It can be observed that the pressure coefficients of panel lower surface increase as the clearance increases, while the pressure on the upper side remains constant. The local turbulence effect caused by the PV panel is enhanced with the increasing of the roof clearance. As the clearance increases from 0.41 m to 1.02 m, the increase of   Figure 12. This demonstrates that the enhancement of the local turbulence effect is partially offset by the pressure equalization effect. Generally, the interaction between the pressure equalization effect and the local turbulence effect exerts influence on the aerodynamic loads as the roof clearance increasing.
In the case of the gable roof, the flow structures around the PV arrays for three clearances are presented (Figures 8(c), 11(a), and 11(b)). The net pressure coefficients on the PV arrays corresponding to these three clearances are compared in Figure 14(b). It can be observed that row 7, which is close to the roof ridge, experiences the largest wind uplift, and the pressure coefficients on the upper and lower sides of row 7 are presented in Figure 13. The increase in the largest negative net pressure coefficient is negligible as the clearance increases from 0.09 m to 0.41 m. The PV panel at roof leeside is predominantly affected by the ridge induced   turbulence, and the flow field remains unchanged with the increasing of the clearance. As the clearance increases from 0.41 m to 1.02 m, the largest negative net pressure coefficient decreases from −0.19 to −0.26. It can be observed that the flow beneath the PV panels is enhanced with the clearance increasing (Figures 8(c), 11(b)). The lower surface of row 7 is attached by the flow across the ridge. This results to an enhancement of wind uplift in row 7. Overall, the largest aerodynamic suction of the PV array tends to increase as the clearance increases for the PV panel mounted on both flat and gable roofs.

Effect of array edge setback
The PV array edge setback from the roof leading edge may be the parameter that most significantly influences the wind loads on the PV array (Kopp, 2014). The influence of the PV array edge setback is investigated by considering three PV array setback values (0.5, 1.2, and 2.1 m), and three PV array tilt angles (10°, 20°, and 30°) for both the flat and gable roofs. Since the roof pressure is more dependent on roof height than building length (Lin & Surry, 1998), the size and strength of large vortices generated by building edge are determined mainly by roof height rather than building depth. The influence of changing building depth D on the pressure distribution is ignored.
Through calculations, the results of different PV panel tilt angles are consistent. Considering the length of the paper, only the results for the 20°tilted PV array are presented in this section. The flow fields around the PV array are presented in Figures 15 and 16. The uplift of the PV array tend to decrease with the increase in the edge setback (Figure 17(a)). For the flat roof, the largest negative net wind load coefficient of the PV array tends to decrease from −0.12 to −0.23 as the PV array edge setback decreases from 2.1 m to 1.2 m. The PV array can be affected by the vortex separated from the roof leading edge significantly as the PV array edge setback decreases. It should also be noted that the decrease in the largest negative wind load coefficient is   subtle as the array edge setback decreases from 1.2 m to 0.5 m. Figure 17(b) shows the wind load coefficients for the gable roof configurations. The largest negative net pressure coefficients of the PV array decrease significantly as the setback decreases from 2.1 m to 1.2 m, while this decrease is negligible when the setback value decreases from 1.2 m to 0.5 m. The turbulence induced by the roof edge has a strong effect on the PV array when the PV array edge setback decreases to a certain range. To take advantage of the roof space and reduce the wind uplift acting on the PV array, a value of 2.1 m is therefore recommended for the PV array edge setback in PV array installations.

Effect of array row spacing
Three values of row spacing are considered for the PV array tilt angle of 2°and 30°while considering both the flat and gable roofs. The flow fields around the PV array on the gable roof are presented in Figures 18 and 19. The resultant net pressure coefficients of the PV array corresponding to both the flat and gable roof cases are presented in Figures 22 and 23, respectively. The variation in the net wind pressure distribution on the array fixed to flat roof is subtle over the range of row spacing considered. For the 2°tilt angle array, the largest negative net pressure coefficient on the PV array decreases from −0.057 to −0.085 as the row spacing increases   from 0.135 m to 1.12 m. The pressure coefficients of row 1 with the 30°PV panel tilt angle for the flat roof are presented in Figure 20. It can be observed that the fluctuations on both the upper and lower sides of the PV panel remain consistent with the increase in the row spacing between the panels. Therefore, the change in the negative net pressure coefficients on the PV array is negligible with the variation in the row spacing for the flat roof configuration.
However, there are significant differences in the net pressure coefficients on the arrays mounted on the gable roof over the range of the row spacing considered ( Figure 23). The largest negative net mean pressure coefficient tends to decrease from −0.024 to −0.165 for the 2°tilted PV array as the row spacing increases from 0.135 m to 1.12 m. The change in the largest negative net aerodynamic load coefficient of the array tilted to 30°angle is also apparent with the variation in the row spacing. This is because the fluctuation corresponding to the upper side of the panel is apparent with the variation in the row spacing, while no significant difference is observed corresponding to the lower side of the panel (Figure 21). An increase in the wind uplift is therefore resulted. This demonstrates that the flow separation at roof ridge is enhanced with the increase in the row spacing for the gable roof configuration; therefore, wind uplift on top surface of row 7 is increased.

Effect of building height
The interaction between the vortex generated by buildings and the PV array can influence the uplift of PV panel. The size and strength of the roof edging induced eddy are controlled by the building height rather than by the building length as Section 4.3 stated. The building width was also reported to be an insignificant parameter about wind-induced pressure distributed on PV panel (Saha, Yoshida, & Tamura, 2011;Stathopoulos, Zisis, & Xypnitou, 2014). The effect of three values of building heights on wind pressure are examined for the flat and gable roof configurations.
The obtained results demonstrate that there are obvious distinctions in the aerodynamic coefficient of the array for the two roof types. For the flat roof, the variations in the net wind pressure coefficients are small for the four tilt angles simulated. Figure 26(a) shows the wind load coefficients on the 20°tilted array. A change in the building height has a negligible effect on the aerodynamic coefficients of the array, which is consistent with the results of the study conducted by Saha et al. (2011) (Figure 24).
In the case of the gable roof, the building height has a significant effect on the pressure coefficient distributions on the PV array. The largest negative net mean pressure coefficients tend to decrease with the increase in the building height, especially for the 2°tilt angle PV array.      ilt angle PV array. The largest negative net mean pressure coefficient on the PV array is significantly decreased with the increase in the building height from 7.3 to 14.6 m. The roof height effect on the largest negative aerodynamic coefficient is weakened as the panel angle increases and almost vanishes for the configuration corresponding to the 30°PV array tilt angle. This is because the localized vortices between the PV panels are predominant for the 30°tilted PV array (Figure 25), and the influence of the turbulence generated by the roof leading edge is insignificant.

Conclusion
A detailed investigation of the wind load characteristics for roof-mounted PV arrays is provided employing the RANS method. Combined with array parameters and roof height, the impact of changing roof types on wind pressure of the PV panel is thoroughly studied. Both flat and gable roofs are considered. The applicability of SST k-ω model to this type of flow has been verified. The main conclusions can be summarized as follows: (1) For both the flat and gable roof configurations, the local vortices become dominant as the PV panel tilt angle increases, and the pressure gradient between the upper and lower sides of the PV array becomes large. The negative net pressure coefficients on the PV array generally decrease with the increment of panel tilt angle.
(2) The local vortices dominate the flow around the PV array in the case of the 20°tilt angle when the roof clearance is small for both the flat and gable roof configurations. With the increase in the roof clearance, the local turbulence is offset by the pressure equalization partially for the flat roof case, while for the gable roof case, the flow across roof ridge attaches to the lower surface of row 7, results to an enhancement of the wind suction of row 7. (3) The wind uplift of the array has a trend of increasing with the decrease in the edge setback for both roof types. The PV array may be subjected to a strong turbulence generated by the roof edge in a certain roof zone. A PV array setback value of 2.1 m in full scale is recommended for PV array installations. (4) The effect of the row spacing is not apparent in the case of a flat roof, however, the negative net pressure coefficients on the PV array are affected significantly by changing the row spacing in the case of a gable roof. The flow separation by the roof ridge is enhanced with the row spacing increasing, which results in higher uplift of row 7. (5) The largest wind uplift at the 2°tilt angle of PV array increase significantly with an increase of the building height in the case of the gable roof, but decrease as the PV array tilt angle increases and almost vanish at 30°tilt angle of PV array. The effect of building height can be neglected in the case of a flat roof corresponding to the four PV panel tilt angles considered.

Future work
This study investigated the aerodynamic structure surrounding the roof-mounted PV array and the net mean C p on PV panels by means of the RANS approach, and mainly analyzing the mean wind loads of panels. The simulated results of downstream panels deviate from the wind tunnel tests apparently due to the limitation of RANS. Large Eddy Simulation (LES) with finer grids is ongoing for investigating dynamic wind pressure of PV panels and the effect of vortex shedding. The research results will be published in future.

Disclosure statement
No potential conflict of interest was reported by the authors.