In snowy winter regions, the combined effect of wind and snow loads is a critical factor in the structural design of solar panel support systems. However, current design codes lack specific combination factors for wind and snow loads on solar panels. This paper investigates the combination factor for a single-axis tracking solar panel in Harbin, China, based on joint wind-snow hazard characteristics. We employ a multi-layer snowmelt model to compute historical ground snow pressure values and then fit probability distributions to wind speed and snow pressure data pairs using four different sampling methods. Linear regression analysis reveals negligible correlation between the two variables. Joint hazard contours are constructed, and finite element analysis of a typical solar panel support structure is performed to determine the load effects. The results indicate that when both wind and snow act as downward pressures on the solar panel and the design considers column axial force and main beam bending moment, the recommended wind-snow load combination factor is 0.7.
1. Introduction
Photovoltaic systems, consisting of lightweight panels mounted on support structures, are vulnerable to combined wind and snow actions. In many cold regions, the simultaneous occurrence of strong winds and heavy snowfall can lead to structural failures. Although several studies have investigated wind loads on solar panels – including wind tunnel tests for pressure coefficients and interference effects – and snow loads on roofs with solar arrays, the combined wind-snow hazard has received much less attention. Existing building codes (e.g., ASCE 7, NBCC, GB 50009) and the Chinese technical specification for photovoltaic support structures (NB/T 10115) prescribe a load combination factor for wind and snow that is identical to that used for general building structures. This generic factor may be overly conservative or unconservative for solar panel systems. Therefore, there is a need to derive a dedicated combination factor based on the probabilistic joint occurrence of wind and snow.
This study focuses on a single-axis tracking solar panel installed in Harbin, a representative city in northeastern China with a snowy climate. The methodology comprises: (i) calculating daily ground snow pressure using a multi-layer energy-balance snowmelt model driven by meteorological data from 1951 to 2016; (ii) extracting wind speed and snow pressure data pairs using four different sampling strategies; (iii) fitting probability models and selecting the optimal ones via K-S tests and AIC; (iv) constructing joint wind-snow hazard contours for different return periods; (v) performing finite element analysis of a solar panel support structure to compute load effects (column axial force, main beam bending moment, and torsion); and (vi) deriving the combination factor that equates the combined load effect to the maximum effect from the joint hazard contour.
2. Methodology
2.1 Multi-Layer Snowmelt Model
To obtain reliable ground snow pressure data, we employ a physically based multi-layer snowmelt model that simulates the accumulation, compaction, and melting of snow on an open field. The model solves the energy balance equation for each snow layer:
\[
\frac{dU_i(t)}{dt} = L_n(t) + S_n(t) + H_s(t) + E_l(t) + Q_c(t) + Q_p(t) + Q_g(t)
\]
where \(U_i\) is the internal energy of layer \(i\); \(L_n\) is net longwave radiation; \(S_n\) is net shortwave radiation; \(H_s\) is sensible heat flux; \(E_l\) is latent heat flux; \(Q_c\) is heat conduction between adjacent layers; \(Q_p\) is heat advected by precipitation; and \(Q_g\) is heat conduction from the ground (assumed to be 0°C). The mass balance for the top layer and internal layers are given by:
\[
\frac{dW_m(t)}{dt} = P_{rain}(t) + P_{snow}(t) – M_{out}(t) – W_e(t)
\]
\[
\frac{dW_i(t)}{dt} = M_{out,i}(t) – M_{out,i+1}(t)
\]
where \(W\) is snow water equivalent; \(P_{rain}\) and \(P_{snow}\) are rainfall and snowfall rates; \(M_{out}\) is meltwater outflow; and \(W_e\) is sublimation/evaporation. Hourly meteorological inputs (temperature, precipitation, wind speed, relative humidity) are derived from daily observed data via linear interpolation. The model has been validated in previous studies. For Harbin, we run the model for each winter from 1951 to 2016, producing a time series of ground snow pressure at hourly intervals.
2.2 Sampling of Wind–Snow Data Pairs
Wind speed data (daily maximum wind speed) are obtained from the same meteorological station. Since the snow accumulation and melting processes extend over several days, we define a “snow event” as a continuous period starting from the first snowfall until the last snow on the ground completely melts. To simplify, we create an envelope of ground snow pressure by taking the maximum value between two consecutive snowfall events. Four methods are then used to pair wind speed with ground snow pressure for each winter:
- Method 1: For each snow event, pair the maximum ground snow pressure (between two snowfall events) with the maximum wind speed occurring in the same interval.
- Method 2: Pair the maximum ground snow pressure (between two snowfall events) with the maximum wind speed observed within the first three days after the snowfall event (since snow particles consolidate after three days).
- Method 3: Pair the maximum ground snow pressure over the entire snow event with the maximum wind speed during that snow event.
- Method 4: Pair the winter’s maximum ground snow pressure with the winter’s maximum wind speed.
These methods yield different sample sizes and capture different physical interpretations of concurrent wind and snow.
2.3 Probability Models and Goodness-of-Fit
For each data pair sample, we fit three probability distributions to the marginal variables (wind speed \(V\) and ground snow pressure \(S\)): Gumbel (Type I), Lognormal, and Generalized Extreme Value (GEV). Parameter estimation uses Maximum Likelihood Estimation (MLE), Method of Moments (MOM), or Least Squares Method (LSM) as appropriate. The Kolmogorov–Smirnov (K-S) test at significance level 0.05 is applied to assess goodness-of-fit. When multiple models pass the K-S test, the Akaike Information Criterion (AIC) is used to select the best model:
\[
AIC = N \ln(RSS) + 2k
\]
\[
RSS = \sum_{i=1}^N (x_i – \hat{x}_i)^2
\]
where \(k\) is the number of parameters (2 for Gumbel and Lognormal, 3 for GEV). The model with the smallest AIC is preferred.
Table 1: K-S test statistic \(D_n\) for wind speed and ground snow pressure (Harbin, all methods).
| Model | Method 1 Wind | Method 1 Snow | Method 2 Wind | Method 2 Snow | Method 3 Wind | Method 3 Snow | Method 4 Wind | Method 4 Snow |
|---|---|---|---|---|---|---|---|---|
| Lognormal (MOM) | 0.034 | 0.089 | 0.048 | 0.093 | 0.045 | 0.170 | 0.144 | 0.078 |
| Lognormal (MLE) | 0.033 | 0.046 | 0.043 | 0.048 | 0.039 | 0.057 | 0.142 | 0.061 |
| GEV (MLE) | 0.036 | 0.044 | 0.040 | 0.044 | 0.038 | 0.086 | 0.140 | 0.067 |
| Gumbel (LSM) | 0.038 | 0.098 | 0.062 | 0.091 | 0.047 | 0.232 | 0.119 | 0.067 |
| Gumbel (MOM) | 0.039 | 0.086 | 0.046 | 0.081 | 0.052 | 0.200 | 0.151 | 0.065 |
| Gumbel (MLE) | 0.036 | 0.044 | 0.047 | 0.044 | 0.038 | 0.175 | 0.140 | 0.075 |
| Critical \(K_\theta\) | 0.056 | 0.056 | 0.099 | 0.167 | ||||
Table 2: AIC values for each probability model (Harbin).
| Model | Method 1 Wind | Method 1 Snow | Method 2 Wind | Method 2 Snow | Method 3 Wind | Method 3 Snow | Method 4 Wind | Method 4 Snow |
|---|---|---|---|---|---|---|---|---|
| Lognormal (MOM) | 1264.6 | -1234.5 | 1645.2 | -1132.4 | 356.9 | -443.2 | 232.2 | -221.3 |
| Lognormal (MLE) | 1243.0 | 171.4 | 1856.5 | 240.3 | 341.4 | -226.9 | 231.3 | -215.7 |
| GEV (MLE) | 1236.1 | -64.2 | 1448.4 | -49.9 | 363.3 | 662.1 | 237.8 | -208.5 |
| Gumbel (LSM) | 1251.5 | -1374.6 | 1623.4 | -1508.6 | 406.2 | -438.2 | 241.0 | -225.9 |
| Gumbel (MOM) | 1248.7 | -1342.6 | 1648.9 | -1480.1 | 409.7 | -417.5 | 235.9 | -222.2 |
| Gumbel (MLE) | 1236.7 | -830.7 | 1931.8 | -969.1 | 413.9 | -244.4 | 235.1 | -228.5 |
The selected optimal models are summarized in Table 3.
Table 3: Optimal probability models for wind speed and ground snow pressure (Harbin).
| Variable | Method 1 | Method 2 | Method 3 | Method 4 |
|---|---|---|---|---|
| Wind speed | GEV (MLE) | GEV (MLE) | Lognormal (MLE) | Lognormal (MLE) |
| Ground snow pressure | GEV (MLE) | GEV (MLE) | Lognormal (MLE) | GEV (MLE) |
2.4 Joint Hazard Contours
We first examine the correlation between wind speed and ground snow pressure using linear regression. The coefficient of determination \(R^2\) (square of Pearson correlation) is computed for each data pair method (Table 4). All values are well below 0.24, indicating negligible correlation. Therefore, wind speed and ground snow pressure can be treated as independent variables.
Table 4: Coefficient of determination \(R^2\) for different data pair methods.
| Method | 1 | 2 | 3 | 4 |
|---|---|---|---|---|
| \(R^2\) | 0.0022 | 0.0134 | 0.0351 | 0.0140 |
For independent variables, the joint exceedance probability for wind speed \(v\) and snow pressure \(s\) is given by:
\[
P(V > v, \; S > s) = \big[1 – F_V(v)\big] \cdot \big[1 – F_S(s)\big]
\]
where \(F_V\) and \(F_S\) are the marginal cumulative distribution functions. The joint return period \(T\) is:
\[
T = \frac{1}{[1 – F_V(v)]\,[1 – F_S(s)]}
\]
For a given return period \(T\), the hazard contour consists of all pairs \((v, s)\) satisfying the above equation. We generate contours for \(T = 10, 25, 50, 100\) years. The 25-year contour is of particular interest because the design service life of a solar panel support structure is typically 25 years. The hazard contours show that as the return period increases, the magnitudes of both wind speed and snow pressure increase.
3. Case Study: Single-Axis Tracking Solar Panel
3.1 Finite Element Model and Load Application
We consider a single-axis tracking solar panel system consisting of a row of panels 100 m long, 1.5 m high, inclined at 30°, supported by a main beam (torque tube) and 13 steel columns spaced 8 m apart. A finite element model of the structure is created. The wind load on the solar panel is calculated using the standard formula:
\[
w_k = \beta_z \mu_z \mu_s w_0
\]
where \(w_0 = 0.5 \rho v^2\) (with \(\rho = 1.25\ \text{kg/m}^3\)), \(\beta_z = 1.2\), \(\mu_z = 1.0\) (height < 10 m, terrain B), and \(\mu_s\) is the shape coefficient. For a 30° inclined panel, the pressure distribution is trapezoidal: at the windward edge \(\mu_{s1} = 11.4\), at the leeward edge \(\mu_{s2} = 10.6\) (pressure), and for suction \(\mu_{s3} = -1.4\), \(\mu_{s4} = -0.6\). The distributed wind load is decomposed into a uniform component plus a torque about the main beam axis. The snow load on the panel is uniformly distributed with a shape coefficient \(\mu_r = 0.85\) (from code for 30° roof). Thus, the snow pressure is \(s_k = \mu_r s_0\).
3.2 Load Combinations and Structural Responses
Two load cases are considered:
- Case 1: Wind pressure + snow pressure (both acting downward).
- Case 2: Wind suction + snow pressure (opposite directions).
For Case 2, the effects partially cancel, so the combined effect is less critical. Therefore, we focus on Case 1. Using the 25-year hazard contours from all four sampling methods, we extract a series of (wind speed, snow pressure) pairs along the contour. For each pair, we compute the structural responses: maximum column axial force \(P_{\max}\), maximum main beam bending moment \(M_{\max}\), and maximum main beam torque \(T_{\max}\). We then identify the maximum combined response \(S_m\) (e.g., \(P_{m}\) or \(M_{m}\)) over the contour. Additionally, we compute the responses due to the 25-year wind load alone (\(S_{25,w}\)) and the 25-year snow load alone (\(S_{25,s}\)) (these correspond to the endpoints of the contour where the other variable is zero).
Table 5: Column axial force responses (kN) for Case 1.
| Method | \(S_{25,w}\) (kN) | \(S_{25,s}\) (kN) | \(S_{P,m}\) (kN) |
|---|---|---|---|
| 1 | -3.60 | -6.22 | -6.39 |
| 2 | -2.87 | -6.16 | -6.28 |
| 3 | -5.16 | -6.17 | -6.47 |
| 4 | -9.85 | -7.00 | -10.81 |
Table 6: Main beam bending moment responses (kN·m) for Case 1.
| Method | \(S_{25,w}\) (kN·m) | \(S_{25,s}\) (kN·m) | \(S_{M,m}\) (kN·m) |
|---|---|---|---|
| 1 | -2.40 | -3.83 | -3.98 |
| 2 | -1.92 | -3.80 | -3.90 |
| 3 | -3.45 | -3.80 | -4.07 |
| 4 | -6.58 | -4.31 | -7.21 |
Note that the torque response depends only on wind speed because snow load is uniform and does not produce net torque. Therefore, the maximum torque is always at the highest wind speed on the contour; no combination effect exists for torque.
3.3 Determination of Combination Factor
The wind-snow load combination factor \(\psi\) is defined as the ratio of the maximum combined load effect to the sum of the individual 25-year load effects:
\[
\psi = \frac{S_m}{S_{25,w} + S_{25,s}}
\]
Since the relationship between loads and effects is linear (static analysis), \(\psi\) equally applies to the loads themselves. We compute \(\psi\) separately for column axial force (denoted \(\psi_P\)) and for main beam bending moment (\(\psi_M\)).
Table 7: Combination factors \(\psi_P\) and \(\psi_M\) from four sampling methods.
| Factor | Method 1 | Method 2 | Method 3 | Method 4 | Mean |
|---|---|---|---|---|---|
| \(\psi_P\) | 0.7 | 0.7 | 0.6 | 0.6 | 0.7 |
| \(\psi_M\) | 0.6 | 0.7 | 0.6 | 0.7 | 0.7 |
The mean values for both load effects round to 0.7. Therefore, we recommend a combination factor of 0.7 for wind and snow loads acting simultaneously as downward pressures on a solar panel. This result is consistent across the different sampling methods, indicating robustness.
For design, the factored load combination can be written as:
\[
\text{Design effect} = 0.7 \times (\text{wind load effect}_{25} + \text{snow load effect}_{25})
\]
or equivalently, the characteristic values of wind and snow loads can be multiplied by 0.7 before summing. The torque response does not require a combination factor because it is generated only by wind load; snow load contributes negligible torque.
4. Conclusion
This study provides a probabilistic framework to determine the wind-snow load combination factor for solar panel support structures. Using Harbin as a case study, we obtained the following key findings:
- The optimal probability models for wind speed and ground snow pressure vary with the sampling method: GEV (MLE) for Methods 1 and 2, and Lognormal (MLE) or GEV (MLE) for Methods 3 and 4. The variables are statistically independent.
- When wind and snow loads both act as downward pressures (Case 1), the combined effect on column axial force and main beam bending moment is critical. When wind acts as suction (Case 2), the effects partially cancel and no combination is needed.
- Based on 25-year joint hazard contours and finite element analysis, the recommended combination factor is 0.7 for both axial force and bending moment. This factor can be applied to the characteristic wind and snow loads in the limit state design of single-axis tracking solar panels in snowy regions.
- The torque response is dominated by wind load alone; hence, no wind-snow combination is required for torsion design.
The combination factor derived here may serve as a reference for future revisions of solar panel design codes. Further research should consider wind-induced snow drifting and variable panel geometries, which may alter the snow load distribution and thus the combination factor.