# Collapse process of pipe-framed greenhouses under snow loading

Kazuya Takahashi¹, Yasushi Uematsu¹,ᵃ

¹Department of Architecture and Building Science, Tohoku University, Sendai 980-8579, Japan

**Abstract.** The present paper investigates the collapse process of pipe-framed greenhouses under snow loading based on a non-linear finite element analysis, in which not only geometrical but also material non-linearity are taken into account. Both even and uneven distributions of snow load are considered. The non-linear behavior of the frame up to the collapse is numerically analysed. Buckling analysis of the frame is also carried out. The objective of the present paper is to establish more reasonable snow resistant design of pipe-framed greenhouses commonly used in the agricultural and horticultural industries in Japan. First, two-dimensional analysis is carried out with a basic model without reinforcement in order to understand the fundamental behavior of the frame up to the collapse, when a part of the frame touches the ground. Then, some measures of reinforcement to improve the snow resistance of the frame are discussed. The effect of reinforcement on the maximum load capacity is made clear. It is found that some types of reinforcement are not effective to the uneven load distribution, although they are effective to the even distribution. Finally, a three-dimensional analysis is carried out with a three dimensional frame model. The structure consists of arch pipes arranged in parallel and connected with each other by horizontal tie-beams (circular pipes). Discussion is made of the effect of three-dimensionality of the frame on the maximum load capacity.

## 1 Introduction

Pipe-framed greenhouses (Figure 1(a)), widely used in the agricultural and horticultural industries in Japan, are generally designed to a lower level of structural safety than conventional buildings, because of the need to minimize capital costs, the demand for a higher level of light transmission and so on. Being inherently light and flexible, they are vulnerable to snow loading and often experience damage due to heavy snowfall (see Figure 1(b), provided by Dr. Moriyama, National Institute for Rural Engineering, Japan). Indeed, many pipe-framed greenhouses were collapsed by heavy snowfall in Kanto District, Japan, during February 14 to 16 of 2014. Depression of arch pipes is the most popular type of collapse [1, 2].

<figure><img src="image_274.png" /><figcaption>(a) Example</figcaption></figure>

<figure><img src="image_275.png" /><figcaption>(b) Collapse due to heavy snow</figcaption></figure>

<caption>**Figure 1.** Pipe-framed greenhouse.</caption>

In Japan, pipe-framed greenhouses are usually designed based on allowable stress or deformation limit following the structural guidelines [3, 4]. In the allowable deformation design, it is assumed that the structure will not lose the resistance because of the restraint effect given by the surrounding members, even if the stresses of some members exceed the elastic range. However, such an assumption has not been verified yet. Furthermore, it is not clear how much safety margin they have against the collapse, when they are designed based on the allowable stress.

The present paper investigates the collapse process of pipe-framed greenhouses under snow loading based on a non-linear finite element analysis. In the analysis, both the geometrical and the material non-linearity are considered because of large deformation and plastic range of stress. Uneven snow accumulation on the roof may be induced by wind (snow drift). Furthermore, snow sliding may occur on the south roof due to solar radiation. Therefore, both even and uneven distributions of snow load are considered.

First, two-dimensional analysis is carried out in order to understand the fundamental characteristics of the collapse process of the basic frame. Then, the effect of reinforcement on the maximum load capacity of the frame is discussed; six types of reinforcements are investigated.

It should be mentioned that the present paper is an extended version of our previous paper [6], focusing on the snow resistance. Furthermore, discussion is made of the effect of reinforcement on the strength of the frame.

## 2 Finite element analysis

### 2.1. Analytical model

The analytical model is schematically illustrated in Figure 2. The span $S$, height $H$, eaves height $h$, and length $L$ are 5.4 m, 3.0 m, 1.55 m and 14.4 m, respectively. The arch pipes are carbon steel tubes for general structural purpose (JIS G3445 STKM11A) with Young's modulus of $E =$

<figure><img src="image_276.png" /><figcaption>**Figure 2.** Analytical Model (Perspective) : unit mm.</figcaption></figure>

<figure><img src="image_277.png" /></figure>

<footnote>ᵃ Corresponding author: yu@archi.tohoku.ac.jp</footnote>
2.05 × 10⁵ N/mm², Poisson's ratio of ν = 0.33, and yield stress of σ_y = 175 N/mm². The external diameter φ and thickness t of the arch pipes are 22.2 mm and 1.2 mm, respectively. The arch pipes are arranged in parallel with a spacing of d = 0.45 m and connected with each other by horizontal tie-beams (circular pipes with the same dimension as that of the arch pipes). The connection of the arch pipes and tie-beams is assumed 'rigid' or 'hinged' for the purpose of simplicity; the actual condition is complicated depending on the construction method. The column base is assumed to be fixed to the ground.

The resistance of pipe-framed greenhouses against
snow load can be improved by adding members to the
frame. Figure 3 shows a schematic illustration of the
reinforcements made on the frame. The members for
reinforcement are the same as that of the arch pipes. The
boundary condition of these added members are assumed
'rigid' or 'hinged'; the practical condition is thought to be
in-between. In Case 5, temporary vertical supports are
installed at intervals of 3.15 m (= 7d) when heavy snowfall
is forecast. Diagonal arches are installed in Case 6, as
shown in Figure 4. Note that this figure shows the
development view of the frame. In Cases 1 to 4, both two-
dimensional and three-dimensional analyses are made, and
the effect of three-dimensionality on the maximum load
capacity is discussed. In Cases 5 and 6, on the other hand,
three-dimensional analysis is made, because the structures
are basically three dimensional.

<figure><img src="image_278.png" /><figcaption>**Figure 3.** Reinforcement methods (Cases 1 to 5).</figcaption></figure>

<figure><img src="image_279.png" /><figcaption>Figure 4. Reinforcement method (Case 6).</figcaption></figure>

<figure><img src="image_280.png" /><figcaption>Figure 5. Finite element model of the arch pipe.</figcaption></figure>

**2.2. Method of analysis**

A nonlinear finite element analysis is made by using a
commercial software ‘ABAQUS 6.13’. The arch pipe is
divided into 58 beam elements (see Figure 5). The
relationship between stress (σ) and strain (ε) of the
material is assumed to be represented by a bi-linear model
as shown in Figure 6, according to Ogawa et al [5]. The
slope of the second polygonal line is approximately 1/420
of that of the first polygonal line. In the analysis, both the
geometric and the material non-linearity are considered.
The arc-length method is used, in which the arc length is
controlled so that the equilibrium path of the structure up
to the collapse can be obtained.

Buckling may occur before reaching the maximum
load. Therefore, buckling analysis is also carried out by
solving an eigenvalue problem. Buckling load and mode
are given by the eigenvalue and eigenvector, respectively.

<figure><img src="image_281.png" /><figcaption>**Figure 6.** Stress-strain relationship of the material.</figcaption></figure>

**2.3. Load distribution**

As mentioned above, uneven snow load distribution may
occur due to the effects of wind and snow sliding.
Therefore, not only even but also uneven distributions of
snow load are considered. The snow loads per unit area, S<sub>W</sub>
and S<sub>L</sub>, on the windward and leeward roofs are assumed
constant over each area, but the ratio of S<sub>W</sub> and S<sub>L</sub> is
changed so that the total load is kept constant. The snow
load per unit area is represented by the product of average
snow density ρ<sub>s</sub> (assumed 9.8 N/m<sup>2</sup>·cm according to Ref.
[3]) and snow depth d<sub>s</sub>. Figure 7 shows three typical
distributions of snow loads; i.e. S<sub>W</sub> : S<sub>L</sub> = 1 : 1 (Type 1),
0.5 : 1.5 (Type 2) and 0 : 2 (Type 3). The total snow load
on the structure is the same for all distribution types. The
snow loads are applied to the nodes as concentrated loads
evaluated by considering the tributary area of each node.
For the purpose of simplicity, the loads are assumed to act
in the vertical direction despite the deformation of the
frame. This assumption is not practical when the
deformation becomes large, resulting in the redistribution
of snow accumulation due to sliding and others.

<figure><img src="image_282.png" /><figcaption>**Figure 7.** Typical load distributions.</figcaption></figure>
# 3 Results on 2-D frames

## 3.1. Collapse process of the basic frame

First, the fundamental characteristics of collapse process of the basic two-dimensional frame with no reinforcement are investigated. Figures 8 shows the transition of the equilibrium configuration of the frame under uniform snow load (Type 1) when the snow load is increased up to

<figure><img src="image_283.png" /><figcaption>**Figure 8.** Collapse process for even load distribution (Type 1, $S_W : S_L = 1 : 1$).</figcaption></figure>

<figure><img src="image_284.png" /><figcaption>**Figure 9.** Collapse process for even load distribution with consideration of initial imperfection (Type 1, $S_W : S_L = 1 : 1$).</figcaption></figure>

<figure><img src="image_285.png" /><figcaption>**Figure 10.** Collapse process for uneven load distribution (Type 2, $S_W : S_L = 0.5 : 1.5$).</figcaption></figure>

<figure><img src="image_286.png" /><figcaption>**Figure 11.** Collapse process for uneven load distribution (Type 3, $S_W : S_L = 0 : 2$).</figcaption></figure>

the maximum value, beyond which the structure collapses. The figure with small open circles represents the deformation of the frame when the maximum stress induced in the frame attains the allowable stress (yield stress). The deformation is relatively small. The dashed line shows the trace of the ridge. In this case, the deformation is always symmetric with respect to the centerline. However, when a small geometric imperfection similar to the buckling mode is introduced, the deflection becomes asymmetric with respect to the centerline, as shown in Figure 9. Note that the maximum imperfection is as small as $0.01\phi$, with $\phi$ being the external diameter of the arch pipes. Similar collapse process is observed for uneven snow load, as shown in Figures 10 and 11. These results imply that the collapse pattern may become asymmetric with respect to the centerline in most cases. However, almost symmetric patterns are often observed in the damage investigations. This discrepancy may be related to the restraint effect of horizontal members connecting the arch pipes and the gable walls (a kind of 3D effect).

## 3.2. Maximum load capacity of the basic frame

Because the snow density $\rho_s$ is assumed constant despite the snow depth $d_s$, the snow load $S$ per unit area is proportional to $d_s$. Therefore, the snow load is represented by a reference snow depth $d_{ref}$ in the present paper. That is, the frame is subjected to uniform snow load of $\rho_s \times d_{ref}$ in Type 1. On the other hand, the windward and leeward halves are respectively loaded by $0.5\rho_s \times d_{ref}$ and $1.5\rho_s \times d_{ref}$ in Type 2. Only the leeward half is loaded by $2\rho_s \times d_{ref}$ in Type 3 (see Fig. 7).

Table 1 shows the snow depths, $d_{ref,a}$ and $d_{ref,u}$, providing the allowable stress (yield stress) and the ultimate strength, respectively. The flowing features can be observed from the results.

1) The $d_{ref,u}/d_{ref,a}$ ratio for the frames subjected to uniform snow load (Type 1) is relatively small compared with those for uneven distributions (Types 2 and 3). This feature implies that the safety margin of this frame up to the collapse is not so large when the frame is designed based on the yield stress as the allowable stress. It may be more reasonable to design the pipe-framed greenhouses based on the ultimate strength.

2) The value of $d_{ref,a}$ for uniform snow load is decreased by a small initial imperfection. This is due to the effect of additional bending moment induced by the initial imperfection. Because some initial imperfection is unavoidable in practice, the snow resistance of the structures should be evaluated by considering such a feature. On the other hand, the $d_{ref,u}$ value seems less sensitive to the initial imperfection. This feature implies again that the design based on the ultimate strength is more reasonable.

3) The unevenness in load distribution affects the $d_{ref,a}$ value significantly. Such an uneven distribution of snow load induces larger bending moments, resulting in a decrease in $d_{ref,a}$. However, the effect on $d_{ref,u}$ is relatively small. It is interesting to note that the $d_{ref,u}/d_{ref,a}$ ratio is almost the same, approximately 1.5, both for Types 2 and 3.
<table><caption>**Table 1.** Snow depth providing the allowable stress and the ultimate strength.</caption><thead><tr><th rowspan="2">Load distribution</th><th colspan="2">Snow depth (cm)</th><th rowspan="2">d<sub>ref,u</sub><br>d<sub>ref,a</sub></th></tr><tr><th>d<sub>ref,a</sub></th><th>d<sub>ref,u</sub></th></tr></thead><tbody><tr><td>Type 1<br>(No initial imperfection)</td><td>14.71</td><td>19.82</td><td>1.35</td></tr><tr><td>Type 1<br>(Initial imperfection)</td><td>14.66</td><td>19.74</td><td>1.33</td></tr><tr><td>Type 2<br>(No initial imperfection)</td><td>11.91</td><td>18.11</td><td>1.52</td></tr><tr><td>Type 3<br>(No initial imperfection)</td><td>10.13</td><td>14.97</td><td>1.48</td></tr></tbody></table>

## 3.2. Effect of reinforcement

As mentioned above, the boundary condition of the added members for reinforcement is assumed ‘rigid’ or ‘hinged’. It was found that the boundary condition does not affect the results significantly. This is because the members are relatively thin and subjected to tension in most cases. Therefore, the results for ‘rigid’ will be shown hereafter.

Table 2 summarizes the results of $d_{ref,a}$ and $d_{ref,u}$ for all reinforcement methods and load distributions obtained from the two-dimensional analysis. It is clear that the snow resistance of pipe-framed greenhouses can be improved by adding members to the frame. The improvement is fairly significant in Cases 1 and 2, particularly when the snow load is uniform, due to the effect of the horizontal member that is subjected to tension. The effect becomes smaller with an increase in the unevenness of load distribution. The horizontal member may interference with the work inside the greenhouse. Then, the horizontal member is moved upward in Case 3. It is found that the improvement

<table><caption>**Table 2.** Effect of reinforcement on the ultimate snow depth.</caption><thead><tr><th rowspan="2">Reinforcement</th><th rowspan="2">Load type<br>(S<sub>w</sub> : S<sub>L</sub>)</th><th colspan="2">Snow depth<br>(cm)</th><th rowspan="2">d<sub>ref,u</sub> / d<sub>ref,a</sub></th></tr><tr><th>d<sub>ref,a</sub></th><th>d<sub>ref,u</sub></th></tr></thead><tbody><tr><td rowspan="3">Case 1</td><td>1<br>(1.0:1.0)</td><td>47.10</td><td>58.73</td><td>1.25</td></tr><tr><td>2<br>(0.5:1.5)</td><td>29.53</td><td>35.40</td><td>1.20</td></tr><tr><td>3<br>(0.0:2.0)</td><td>19.41</td><td>24.07</td><td>1.24</td></tr><tr><td rowspan="3">Case 2</td><td>1<br>(1.0:1.0)</td><td>56.46</td><td>65.10</td><td>1.16</td></tr><tr><td>2<br>(0.5:1.5)</td><td>34.33</td><td>41.18</td><td>1.20</td></tr><tr><td>3<br>(0.0:2.0)</td><td>21.41</td><td>29.05</td><td>1.36</td></tr><tr><td rowspan="3">Case 3</td><td>1<br>(1.0:1.0)</td><td>28.52</td><td>33.67</td><td>1.18</td></tr><tr><td>2<br>(0.5:1.5)</td><td>20.95</td><td>24.03</td><td>1.15</td></tr><tr><td>3<br>(0.0:2.0)</td><td>16.08</td><td>18.93</td><td>1.18</td></tr><tr><td rowspan="3">Case 4</td><td>1<br>(1/0:1.0)</td><td>28.12</td><td>35.77</td><td>1.27</td></tr><tr><td>2<br>(0.5:1.5)</td><td>28.67</td><td>35.06</td><td>1.22</td></tr><tr><td>3<br>(0.0:2.0)</td><td>27.77</td><td>34.21</td><td>1.23</td></tr></tbody></table>

on the strength is limited compared with that in Cases 1 and 2. Case 4 seems to be the best among the four cases tested from the viewpoints of both strength and workability. This structure provides similar maximum load capacity for all load distributions. Even in this case, the $d_{ref,u}d_{ref,a}$ ratio is not so large, approximately 1.2.

## 4 Results on 3-D frames

In the three-dimensional model, the arch pipes are arranged in parallel with a spacing of $d = 0.45$ m and connected with each other by horizontal tie-beams at the ridge and eaves (see Figure 2). As mentioned above, two kinds of joint conditions are considered; i.e. ‘rigid’ and ‘hinged’. Buckling analysis is also carried out.

The results for the two joint conditions are summarized in Tables 3 and 4. The results for the basic frame with no reinforcement are almost the same as those obtained from the two-dimensional analysis despite the boundary condition. When the joint is ‘rigid’, the critical loads of buckling, not shown in the table, are generally much larger than the maximum loads obtained from the non-linear analysis. When the joint is ‘hinged’, on the other hand, the buckling load becomes smaller than the maximum load in some cases. The snow depth corresponding to the buckling is represented by $d_{ref,b}$ in the table. Figure 12 shows the buckling mode in Case 4. The buckling mode is dominated by the displacements in the longitudinal direction with almost no in-plane deformation of the frame. This is due to low horizontal stiffness of the structure. Similar failure modes are often observed in the snow-induced damage to greenhouses, as shown in Figure 13. This greenhouse is not pipe-framed greenhouse but a reinforced pipe-framed greenhouse, in which the frame is constructed of square steel tubes (about 50 mm × 500 mm) and the roof is constructed of arch pipes. Because the horizontal stiffness of this kind of structure is rather low, the structure tends to be collapsed by snow load, as shown in Figure 13.

<figure><img src="image_287.png" /><figcaption>**Figure 12.** Buckling mode (Case 4).</figcaption></figure>

<figure><img src="image_288.png" /><figcaption>**Figure 13.** Collapse of a reinforced pipe-framed greenhouse due to snow load.</figcaption></figure>
<table><caption>**Table 3.** Snow resistance of 3D frames: the connection of arch pipes and tie-beams is assumed 'rigid'.</caption><thead><tr><th rowspan="2">Frame</th><th rowspan="2">Load type<br>(Sw: SL)</th><th colspan="2">Snow depth<br>(cm)</th><th rowspan="2">d<sub>ref,u</sub> / d<sub>ref,a</sub></th></tr><tr><th>d<sub>ref,a</sub></th><th>d<sub>ref,u</sub></th></tr></thead><tbody><tr><td rowspan="3">Basic</td><td>1 (1.0:1.0)</td><td>14.56</td><td>20.44</td><td>1.40</td></tr><tr><td>2 (0.5:1.5)</td><td>11.81</td><td>18.85</td><td>1.60</td></tr><tr><td>3 (0.0:2.0)</td><td>10.03</td><td>15.40</td><td>1.54</td></tr><tr><td rowspan="3">Case 1</td><td>1 (1.0:1.0)</td><td>45.32</td><td>58.89</td><td>1.30</td></tr><tr><td>2 (0.5:1.5)</td><td>29.52</td><td>35.73</td><td>1.21</td></tr><tr><td>3 (0.0:2.0)</td><td>19.36</td><td>24.24</td><td>1.25</td></tr><tr><td rowspan="3">Case 2</td><td>1 (1.0:1.0)</td><td>55.15</td><td>65.45</td><td>1.19</td></tr><tr><td>2 (0.5:1.5)</td><td>33.30</td><td>41.46</td><td>1.25</td></tr><tr><td>3 (0.0:2.0)</td><td>22.53</td><td>29.27</td><td>1.30</td></tr><tr><td rowspan="3">Case 3</td><td>1 (1.0:1.0)</td><td>28.14</td><td>33.79</td><td>1.20</td></tr><tr><td>2 (0.5:1.5)</td><td>20.88</td><td>24.47</td><td>1.17</td></tr><tr><td>3 (0.0:2.0)</td><td>14.84</td><td>19.39</td><td>1.31</td></tr><tr><td rowspan="3">Case 4</td><td>1 (1.0:1.0)</td><td>14.94</td><td>36.06</td><td>2.41</td></tr><tr><td>2 (0.5:1.5)</td><td>12.20</td><td>35.55</td><td>2.91</td></tr><tr><td>3 (0.0:2.0)</td><td>10.43</td><td>34.92</td><td>3.35</td></tr><tr><td rowspan="3">Case 5</td><td>1 (1.0:1.0)</td><td>12.52</td><td>39.65</td><td>3.17</td></tr><tr><td>2 (0.5:1.5)</td><td>12.25</td><td>28.44</td><td>2.32</td></tr><tr><td>3 (0.0:2.0)</td><td>12.12</td><td>19.81</td><td>1.63</td></tr><tr><td rowspan="3">Case 6</td><td>1 (1.0:1.0)</td><td>15.19</td><td>22.08</td><td>1.45</td></tr><tr><td>2 (0.5:1.5)</td><td>13.24</td><td>20.05</td><td>1.51</td></tr><tr><td>3 (0.0:2.0)</td><td>11.35</td><td>16.33</td><td>1.44</td></tr></tbody></table>

<table><caption>**Table 4.** Snow resistance of 3D frames: the connection of arch pipes and tie-beams is assumed 'hinged'.</caption><thead><tr><th rowspan="2">Frame</th><th rowspan="2">Load type<br>(Sw: SL)</th><th colspan="3">Snow depth<br>(cm)</th><th rowspan="2">d<sub>ref,u</sub> / d<sub>ref,a</sub></th></tr><tr><th>d<sub>ref,a</sub></th><th>d<sub>ref,u</sub></th><th>d<sub>ref,b</sub></th></tr></thead><tbody><tr><td rowspan="3">Basic</td><td>1 (1.0:1.0)</td><td>13.99</td><td>19.61</td><td>28.41</td><td>1.40</td></tr><tr><td>2 (0.5:1.5)</td><td>11.17</td><td>17.82</td><td>27.50</td><td>1.60</td></tr><tr><td>3 (0.0:2.0)</td><td>9.40</td><td>14.82</td><td>25.39</td><td>1.58</td></tr><tr><td rowspan="3">Case 1</td><td>1 (1.0:1.0)</td><td>41.83</td><td>51.45</td><td>26.24</td><td>*</td></tr><tr><td>2 (0.5:1.5)</td><td>21.70</td><td>25.07</td><td>25.57</td><td>1.16</td></tr><tr><td>3 (0.0:2.0)</td><td>13.96</td><td>17.07</td><td>23.93</td><td>1.22</td></tr><tr><td rowspan="3">Case 2</td><td>1 (1.0:1.0)</td><td>42.19</td><td>50.72</td><td>26.23</td><td>*</td></tr><tr><td>2 (0.5:1.5)</td><td>21.66</td><td>25.07</td><td>25.56</td><td>1.16</td></tr><tr><td>3 (0.0:2.0)</td><td>13.95</td><td>17.08</td><td>23.92</td><td>1.22</td></tr><tr><td rowspan="3">Case 3</td><td>1 (1.0:1.0)</td><td>23.56</td><td>30.82</td><td>27.00</td><td>**</td></tr><tr><td>2 (0.5:1.5)</td><td>17.60</td><td>22.91</td><td>25.99</td><td>1.30</td></tr><tr><td>3 (0.0:2.0)</td><td>13.68</td><td>17.83</td><td>23.74</td><td>1.30</td></tr><tr><td rowspan="3">Case 4</td><td>1 (1.0:1.0)</td><td>14.34</td><td>35.26</td><td>26.65</td><td>**</td></tr><tr><td>2 (0.5:1.5)</td><td>11.57</td><td>34.81</td><td>25.25</td><td>**</td></tr><tr><td>3 (0.0:2.0)</td><td>9.79</td><td>34.08</td><td>22.65</td><td>**</td></tr><tr><td rowspan="3">Case 5</td><td>1 (1.0:1.0)</td><td>10.44</td><td>26.56</td><td>29.06</td><td>2.54</td></tr><tr><td>2 (0.5:1.5)</td><td>10.25</td><td>22.90</td><td>28.36</td><td>2.23</td></tr><tr><td>3 (0.0:2.0)</td><td>10.16</td><td>17.59</td><td>26.55</td><td>1.73</td></tr><tr><td rowspan="3">Case 6</td><td>1 (1.0:1.0)</td><td>13.71</td><td>19.41</td><td>43.21</td><td>1.42</td></tr><tr><td>2 (0.5:1.5)</td><td>11.34</td><td>17.75</td><td>43.20</td><td>1.57</td></tr><tr><td>3 (0.0:2.0)</td><td>9.42</td><td>14.67</td><td>43.08</td><td>1.56</td></tr></tbody></table>

Because the diagonal arches introduced in Case 6 increases the horizontal stiffness of the frame with 'hinged' joints, the value of *d*<sub>ref,b</sub> increases significantly.

The following features can be observed from Tables 1 to 3:

1) The values of *d*<sub>ref,a</sub> and *d*<sub>ref,u</sub> for 'rigid' connection are generally larger than those for 'hinged' connection. This is related to the stiffening effect of tie-beams on the arch pipes in the case of 'rigid' connection.

2) In the case of 'hinged' connection, the value of *d*<sub>ref,b</sub> is generally small except for Case 6 where the horizontal stiffness of the frame is increased by the diagonal braces. The value of *d*<sub>ref,b</sub> is comparative to or smaller than *d*<sub>ref,u</sub> in many cases.

3) The value of *d*<sub>ref,u</sub> for the 'rigid' connection is generally larger than that for the 2D frame. However, the difference is only a little, which indicates that the 3D effect on the maximum load capacity is relatively small.

4) The value of *d*<sub>ref,a</sub> for the 3D frame is not necessarily larger than that for the corresponding 2D frame. This feature is contrary to our expectation that the constraint of tie-beams would decrease the stress induced in the members and increase the load capacity. Careful investigation of the results on the stresses induced in the members of the 3D model indicates that the location of *d*<sub>ref,a</sub> is different from that for the 2D model. For example, in Case 4, *d*<sub>ref,a</sub> occurs at the top of the arch pipe in the 2D model, while at the column base of an arch pipe in the 3D model. In Case 5, *d*<sub>ref,a</sub> occurs in the tie-beam at the top of the arch pipes in the 3D model. The 3D effect induces additional stresses in the members and the maximum stress is induced at some point different from that in the 2D frame. In such a case, 2D analysis is not appropriate, because it will underestimate the snow resistance of the structure, when designing it based on the allowable stress.

5) The value of *d*<sub>ref,u</sub> is less sensitive to the analytical model. Therefore, it is more desirable to design the structure based on the maximum load capacity not on the allowable stress.

## 5 Concluding remarks

First, the two dimensional analysis shows the fundamental characteristics of the collapse process of pipe-framed greenhouses under snow loading. Both initial imperfection and uneven load distribution decrease the snow depth providing the allowable stress (yield stress) in the member, significantly. Therefore, the snow resistance of the structure should be evaluated by considering these factors, when the structural design is carried out based on the allowable stress. On the other hand, the maximum load capacity of the frame seems less sensitive to these factors.

Then, the effects of various reinforcements on the snow loads providing the allowable stress and the maximum load were investigated. The boundary condition of the added members minutely affects the load capacity of the frame. Generally speaking, the reinforcements tested in the present paper increases the snow resistance of

<footnote>* *d*<sub>ref,b</sub> is the lowest.</footnote>

<footnote>** *d*<sub>ref,b</sub> is larger than *d*<sub>ref,a</sub> but smaller than *d*<sub>ref,u</sub>.</footnote>
the structures significantly. However, some of them are
not so effective for uneven load distribution.

Finally, three-dimensional frames with and without reinforcement were analysed. The structure consists of arch pipes arranged in parallel and connected with each other by horizontal tie-beams. When the horizontal stiffness of the frame is low, the buckling provides the critical condition. However, in the practical condition, where the connection between the arch pipe and tie beam is regarded as ‘rigid’, the buckling load is much larger than the maximum load obtained from the non-linear analysis. The maximum load is not so sensitive to the analytical model that the snow resistance of the frames can be evaluated based on the maximum load not on the allowable stress.

**6 Acknowledgement**

The present study is financially supported by The Japan
Iron and Steel Federation Grant (2016). The authors thanks
Dr. Moriyama at National Institute for Rural Engineering
of Japan, for providing us useful information on the snow
resistance of greenhouses.

**References**

1. Y. Uematsu, Damage to structures due to heavy snow during February 14 to 16, 2014 – Mainly on the damage to carports and greenhouses in Gunma Prefecture -, Tohoku Journal of Disaster Science, Vol.51, pp.165-170 (2015). (in Japanese)

2. M. Hideki, S. Inoue, I. Kamiishi, Damage to Greenhouses by Heavy Snow in Feb. 2014, JSSI & JSSE Joint Conference – 2014/Hachinohe (2014). (in Japanese)

3. Japan Greenhouse Horticulture Association, Horticultural facility safety structural criteria, (1997). (in Japanese)

4. Japan Greenhouse Horticulture Association, Underground push-fit pipe-framed greenhouse safety structure guideline, (1988). (in Japanese)

5. H. Ogawa, I. Tsuge, Y. Sato, S. Hoshiba, S. Yamashita, Experimental analysis on strength of pipe-houses with ground anchoring (1) – Actual size experiment -, The Journal of the Society of Agricultural Structures, Japan, Vol. 9, No.3, pp.29-38 (1989). (in Japanese)

6. K. Takahashi, Y. Uematsu, Collapse processes of pipe-framed greenhouses under wind and snow loads, The Journal of the Society of Agricultural Structures, Japan, Vol.47, No.1, pp. 108 (2016). (in Japanese)
# Study on evaluation of roof snow load considering rain-on-snow surcharge: Statistical evaluation of snow cover and precipitation in winter in Japan

Masaya Otsuki<fnref n="1" />,ª, Toru Takahashi<fnref n="2" />, Yoshihiko Saito<fnref n="1" />, Takuya Tsutsumi<fnref n="3" /> and Kikitsu Hitomitsu<fnref n="4" />

**Abstract.** To examine regions for which the roof snow load, influenced by rain after snow, has to be estimated, a statistical evaluation of snow cover and rainfall during the winter period in Japan was conducted to examine the calculation method. In a comparison of daily and two-day rainfall, the latter was slightly greater without substantial difference in the distribution. The two-day rainfall is considered to be an appropriate data for considering the rain load. Investigation of areal distribution of rainfall revealed that the rainfall for an entire winter was greater in regions that received large amounts of rainfall in summer. The extreme values of rainfall on snow cover were found to be higher in regions with a long snow-cover period, ranging from Hokuriku to Tohoku and southern Hokkaido as well as the Kanto-Koshin region. The load conversion value of the rainfall is often relatively higher in general areas with snow cover depth of less than 1 m, such as the Pacific Ocean side of Tohoku and the Kanto-Koshin region, compared with their roof snow load. In such regions, there is a concern that the rain-on-snow surcharge, when added to the roof snow load, will exceed the design load.

## 1 Introduction

From February 14–16, 2014, the Kanto region of Japan was hit by heavy snow caused by a rapidly developing low-pressure system. A large number of buildings were damaged. The roof snow-load increase owing to a rain-on-snow surcharge was indicated as the main cause of the damage. Figure 1 shows the weather time series observed at Tokyo District Meteorological Observatory in one of the regions that sustained damage. The snow started on February 14, accumulated 27 cm, and then turned to rain after 2:00 a.m. on February 15. The snowfall was approximately 42 mm, and the rainfall was 54 mm, for a total precipitation of 96 mm. Since the design load of the region is 0.6 kN/m<fnref n="2" />, a roof snow load of approximately 0.42 kN/m<fnref n="2" /> was within the permissible range. Subsequently, as shown in the time-series graph, a rainfall of approximately 20 mm (0.2 kN/m<fnref n="2" />) led to the collapse of buildings. It is suspected that the snow on the roofs absorbed the rain and that a load of 0.62 kN/m<fnref n="2" />, exceeding the design load, was applied to the roof and led to the collapse. (Since this precipitation does not consider the capture rate owing to the wind velocity, the actual precipitation may have been greater).

<figure><img src="image_289.png" /><figcaption>**Figure 1.** Weather time series of Tokyo at time of damage owing to rainfall increasing the load (February 15, 2014).</figcaption></figure>

In the United States, the rain-on-snow surcharge of a roof snow load is specified in ASCE7 [1]. In regions where the roof snow load is 0.96 kN/m<fnref n="2" /> or less, all roofs with a slope angle of less than W/15.2° are required to allow for a surcharge of 0.24 kN/m<fnref n="2" />. [W refers to the roof span (m).] In Canada, NBC sets the roof snow load S<sub>S</sub> and rain load S<sub>R</sub> for each city, and the sum of these is to be used [2]. Values for both countries are based on a 50-year return period, but they can be considered as precautionary values if snowfall and rainfall are perceived not as individual events but as interdependent events. On the other hand, there has not been any clear regulation in Japan. This study examines the effect of the rain-on-snow surcharge on the roof snow load, and analyzes the characteristics of snow cover and rainfall in winter using the daily snowfall data. The goal is to establish a calculation method for the surcharge load.

## 2 Data used and methods of analysis

Precipitation, air temperature, and snow depth data measured by meteorological stations of the Japan Meteorological Agency were used for analysis. The meteorological stations include manned meteorological offices (142 locations) and 1500 unmanned AMeDAS (195 of which measures snow depth) (Figure 2). All data from the beginning of the observation period until 2014 (from December to March only) were used. Okinawa was excluded because snow does not fall in this region.

Rainfall and snowfall were distinguished according to the daily average temperature based on the precipitation data to extract the daily precipitation data of the days when rain was determined to have fallen. Winter daily rainfall was then calculated using the assumption that snow changes to rain when the daily average air temperature rises above 2°C [3] among various proposed methods of rain/snow determination. Similarly, the

<figure><img src="image_290.png" /></figure>

<footnote n="1">¹Yukiken Snow Eaters Co., Ltd., 5-6, Nishi-7, Minami-2, Chuo-ku, Sapporo 060-0062, Japan</footnote>

<footnote n="2">²Department of Architecture, Chiba University, 1-33 Yayoicho, Inage-ku, Chiba 263-8522, Japan</footnote>

<footnote n="3">³Northern Regional Building Research Institute, Higashi-1-3-1-20, Midorigaoka, Asahikawa 078-8801, Japan</footnote>

<footnote n="4">⁴National Institute for Land and Infrastructure Management, Tachihara-1, Tsukuba 305-0802, Japan</footnote>

<footnote>ª Corresponding author: otsuki@snow-eaters.com</footnote>
maximum rainfall for each winter was obtained by extracting the rainfall precipitation data of only the stations where snow cover was observed (snow depth of more than 1 cm). To take into account the nighttime rainfall events that may continue to the next day, two-day rainfall data were also compiled. Further, the extreme value of the snow depth during rain was statistically analyzed to evaluate the roof snow load, as it is necessary to understand the snow depth characteristics during rain in winter.

<figure><img src="image_291.png" /><figcaption>**Figure 2.** Meteorological stations whose data were analyzed.</figcaption></figure>

# 3 Statistical evaluation of daily and two-day rainfall

## 3.1 Occurrence frequency of daily and two-day rainfall

Figures 3 and 4 show the ranked averaged winter maximum daily and two-day rainfall by number of stations. In the graph, the white bars represent the values for the entire winter, and the black bars represent the collected statistics of the rainfall on days when there was snow cover. The number of observation stations with snow cover is less than the total number for the entire winter because the former is the statistic for only the points where snow cover was observed. According to the graph, the greatest number of stations observed 10–20 mm as the multiyear average of the maximum winter daily rainfall in the case of snow cover, with the highest values at 40–50 mm. In case of the entire winter period, the highest number of stations observed 30–40 mm. The data for which the number of observation stations peak for the two-day rainfall are 10–20 mm and 40–50 mm, respectively. Compared with the daily data, two-day rainfall precipitation tends to be higher, but the difference is not significant. This means that a rainfall event usually

takes no more than two days, and a substantial amount of rain does not fall on both days.

Figures 5 and 6 show the ranked daily and two-day extreme rainfall by number of observation stations. In the case of snow cover, the extreme value distribution of the daily rainfall peaks at 50–60 mm, while the two-day data is distributed throughout a range of 40–90 mm without a distinguishing peak. For the entire winter period, the extreme values of the daily rainfall are in the range of 60–70 mm. For two-day rainfall, the extreme values are in the range of 80–90 mm. This shows a tendency that the peak values are higher by approximately 30–40 mm over

<figure><img src="image_292.png" /><figcaption>**Figure 3.** Ranked averaged winter maximum daily rainfall by number of observation stations.</figcaption></figure>

<figure><img src="image_293.png" /><figcaption>**Figure 4.** Ranked averaged winter maximum two-day rainfall by number of observation stations.</figcaption></figure>

<figure><img src="image_294.png" /><figcaption>**Figure 5.** Ranked extreme winter maximum daily rainfall by number of observation stations.</figcaption></figure>

<figure><img src="image_295.png" /><figcaption>**Figure 6.** Ranked extreme winter maximum two-day rainfall by number of observation stations.</figcaption></figure>
those of the multiyear averaged values. In terms of the peak occurrence of observation stations, the rainfall for the entire winter period tended to be greater by a range of 10–20 mm for daily rainfall and approximately 30 mm (although the distribution is unclear) for two-day rainfall than that with snow cover.

## 3.2 Areal distribution of daily and two-day rainfall

### 3.2.1 Precipitation distribution in Japan

Regional distribution of rainfall in winter is examined in this section. First, Figure 7 shows the mesh climatic data of precipitation from the Japan Meteorological Agency. Values tend to be higher on the Japan Sea side from

<figure><img src="image_296.png" /></figure>

<figure><img src="image_297.png" /><figcaption>**Figure 7.** Mesh climatic data of monthly precipitation (Japan Meteorological Agency).</figcaption></figure>

Tohoku to the Hokuriku-Chugoku regions in winter (January), and extremely low along the coast of the Pacific Ocean from Hokkaido to Kyushu. This is a result of the typical winter pressure pattern (east: high pressure, west: low pressure), where snow clouds bring heavy snow to the coast of the Japan Sea. In summer (August), on the other hand, the coast of the Pacific Ocean from Kanto to Kyushu regions receives considerable rain owing to substantial rainfall from typhoons and stationary fronts.

### 3.2.2 Rainfall distribution for the entire winter period

With the above facts in mind, the multiyear average and extreme distribution of winter rainfall (two-day data, the entire winter period), shown in Figures 8 and 9 are examined. The distribution of average values overall is smaller in northern Japan and greater in southern Japan. Values tend to be higher on the Pacific Ocean side from the Tohoku to Shikoku regions, Kyushu, and the Japan Sea side from the Hokuriku to Chugoku regions (in the range of 70–150 mm). On the other hand, the values are smaller in Hokkaido and the inland part of Tohoku in the range of 20–30 mm, unlike the distribution of the winter precipitation according to the mesh climatic data. It is thought that most of the precipitation in winter is snowfall, that most of it falls on the Japan Sea side, and that two-day rainfall is often observed on the Pacific Ocean side as well, similar to the precipitation characteristics in summer.

Extreme values appear to indicate that the overall rainfall is increasing compared with the average, and that the difference in the area distribution of rainfall is narrowing. The values for inland Hokkaido are similarly low at approximately 20 mm as is the multiyear average, but those on the Pacific Ocean side from Tohoku to Kyushu exceed 100 mm. The rainfall on the Pacific Ocean side of Hokkaido and inland Tohoku where the

<figure><img src="image_298.png" /><figcaption>**Figure 8.** Average distribution of winter maximum rainfall (two-day data, entire winter period).</figcaption></figure>
<figure><img src="image_299.png" /><figcaption>**Figure 9.** Extreme value distribution of winter maximum rainfall (two-day data, entire winter period).</figcaption></figure>

average values were relatively small just as high as other regions with regard to the extreme values. The reason for this is probably that substantially large rainfall events are recorded in these regions depending on the year. Such records become apparent by plotting the extreme values, but they do not stand out when averaged over multiple years.

### 3.2.3 Distribution of rainfall in snow cover period

In this section, the distribution of the multiyear average and extreme values of two-day rainfall on snow cover are examined (Figures 10 and 11). The number of plots is smaller than that for the entire winter because only the stations that measure snow depth are included, and the distribution trend is completely different. Regarding the multiyear average, the values are high at 30–60 mm on the Japan Sea side from the Tohoku to Chugoku regions, but small at 30 mm or less in other regions. The regions where the average values are high roughly overlap the regions of high rainfall in the winter period, and the rainfall during this period is more likely to be reflected in the statistics because the snow cover exists for a longer time in these regions. In regions with less precipitation and almost no snow cover period in winter, on the other hand, the two-day rainfall is expected to be extremely small when there is snow cover.

According to the extreme value results, the majority of stations recorded 50 mm or more, and values for the Pacific Ocean side of the Hokkaido, Tohoku, and Kanto-Koshin regions tended to be high, in addition to the stations on the Japan Sea side where the average values were high. While the values for the average were high only on the Japan Sea side, plotting the extreme values reveals large values in the range of 100 mm even in the Tohoku and Kanto-Koshin regions. This suggests that heavy rainfall events occasionally occur once every few years in such regions. From the extreme value distribution of two-day rainfall on snow cover, heavy

rainfall events are likely to be extracted on the Japan Sea coast because the snow cover exists for a longer time between December and March. On the other hand, conditions for snow cover forming are rare in Kanto-Koshin and the Pacific Ocean side of the Tohoku region, but it is thought that meteorological disturbances are likely to occur in these areas, causing heavy rainfalls during the snow cover periods.

<figure><img src="image_300.png" /><figcaption>**Figure 10.** Average distribution of winter maximum rainfall (two-day data, snow cover days only).</figcaption></figure>

<figure><img src="image_301.png" /><figcaption>**Figure 11.** Extreme value distribution of winter maximum rainfall (two-day data, snow cover days only).</figcaption></figure>

## 3.3 Summary

To summarize up to this point, the regions with more rainfall in summer tend to have higher daily and two-day rainfall during the entire winter period, while the
multiyear average is high in the Hokuriku region for snow cover. Extreme values are recorded in a relatively wide area, including the Pacific Ocean side of the Tohoku and Kanto-Koshin regions, where it snows little. In other words, snowfall and rainfall events occur occasionally in the regions with little snow cover and a small roof snow load, and dangerous conditions easily exceeding the roof snow load could occur depending on the amount of snow cover and rainfall. Thus, for the regions that need to allow for a rain-on-snow surcharge, the rainfall ratio for the original roof snow load is important.

## 4 Extreme value analysis of snow depth and rain load ratio

According to the 50-year probability [3] distribution of the annual maximum snow depth, snow cover is as deep as 150–200 cm in the Hokuriku region and on the Japan Sea side of Hokkaido (Figure 12). The roof snow load for these regions is known to be very high. Therefore, even if a rainfall event occurs, the relative load would not increase much. By contrast, the roof snow load is small in the coastal part of Tohoku and in the areas west of Kanto-Koshin. There, the relative load would increase substantially during rainfall. Considering this, the two-day rainfall load conversion value ratio (multiyear average and 50-year return period) in relation to the roof snow load is shown in Figures 13 and 14. The roof snow load here is the product of the snow cover depth shown in Figure 12 and the snow cover density of 0.3 kN/m³ in heavy snow areas that recorded a snow depth of 1 m or more, and of 0.2 kN/m³ in other general areas. The unit load of snow cover can be specified as 20 N/m²/cm in general areas, and more for the heavy snow areas [4]. The ratio is 1 or more at locations with a two-day rain load exceeding the roof snow load. These regions are consistent with those having less snow cover such as Kanto on the Pacific Ocean side of Kinki, and the Shikoku and Kyushu regions.

<figure><img src="image_302.png" /><figcaption>**Figure 12.** 50-year probability of snow depth distribution (according to Japan Meteorological Agency) [3].</figcaption></figure>

<figure><img src="image_303.png" /><figcaption>**Figure 13.** Ratio between two-day rainfall-equivalent load(average) and roof snow load.</figcaption></figure>

<figure><img src="image_304.png" /><figcaption>**Figure 14.** Ratio between two-day rainfall-equivalent load(50-year return period) and roof snow load.</figcaption></figure>

By contrast, the ratio is relatively small in regions with heavy snowfall such as the Hokkaido, Tohoku, and Hokuriku regions. Figure 15 shows the number of occurrences ranked by ratio. In heavy snow areas, the ratio for many stations is 0.3 or less in the 50-year probability values for two-day rainfall precipitation, and that for most stations is 2.0 or less for multiyear average values. On the other hand, the ratio for 70% and 50%, respectively, of the observation stations exceeds 1.0 (the rain load exceeding the roof snow load) in the general areas with less snow cover. Further, although the occurrence numbers are small, there are occurrences that appear evenly in the 0.3 to 0.9 range, and it is possible that the risk of rain-on-snow events applying a load far greater than the bearing strength of the structure is quite high. An inference can be made from these data that rain-on-snow surcharge must be considered in general areas
with a snow cover depth of less than 1 m in Japan. However, because snow is not thought to absorb much rain unless there is deep snow cover, general areas with less snow cover, many of which are expected to be stations with the above ratio of 1 or higher, should be excluded from consideration.

<figure><img src="image_305.png" /><figcaption>**Figure 15.** Occurrence percentage ranked by ratio between two-day rainfall-equivalent load and roof snow load.</figcaption></figure>

## 5 Conclusion

In this study, the snowfall and rainfall characteristics of winter in Japan were investigated.

- Statistics of maximum annual daily rainfall and two-day rainfall in the winter period based on precipitation data from the Japan Meteorological Agency revealed that the occurrence count peaked at 30-40 mm for the average daily rainfall during the entire winter period, and at 10-20 mm when there was snow cover. With regard to extreme values, peaks were observed at 60-70 mm and 50-60 mm, respectively.

- A comparison of daily and two-day rainfall showed that two-day rainfall was slightly greater without a substantial difference in the distribution. Thus the two-day rainfall is considered to be appropriate for examining the rain load.

- The number of data for the entire winter period naturally tended to be greater than those with snow cover in the processing of maximum values for statistics. The statistical values of rainfall during the entire winter period as a precautionary measure may be applied to a surcharge load. If there is a strong dependency and connection between snowfall and rainfall, however, the use of rainfall on snow cover is a possibility, and further investigation into the relationship between rainfall and snowfall is necessary.

- An investigation of the areal distribution of rainfall revealed that the rainfall for the entire winter was greater in regions where summer rainfall was high. The extreme values of rainfall on snow cover were

found to be higher in regions with a long snow cover period, ranging from Hokuriku to Tohoku and southern Hokkaido as well as the Kanto-Koshin regions.

- The load conversion value of the two-day rainfall is often relatively higher in the general areas on the Pacific Ocean side with snow cover depth of less than 1m, such as the Kanto-Koshin region, in relation to their snow load. In such regions, there is a concern that the roof snow load exceeds the design load owing to the rain-on-snow surcharge.

## 6 Future Challenges

The actual degree of load increase during continuous rainfall on snow on sloped roofs is still not known. According to Otsuki (2015) [5], the rain-on-snow load peaks when rainfall becomes heavy, the water starts to drain, and the speed of the rainfall becomes roughly equal to that of the drainage. Figure 16 shows how the load increases when a constant strength of rain (10 mm/h) falls on snow on four different roof models. The dashed line in the figure indicates the load conversion of the cumulative rainfall, and the load is plotted on this line when there is no drainage.

According this figure, the greater the span and the smaller the slope of the roof, the less rainfall it takes to reach the peak. The load on all roofs shifts slightly to the right of the dashed line by a few tens of millimeters before reaching the peak, indicating that the water started to drain during the rain. The drainage further accelerated after the peak load was reached on all roofs, and the load continued to decrease (the 50-m roof reached its peak at approximately 70 mm). In other words, reflecting all of the rainfall in the surcharge is too much. By continuing this experiment, the roof shape, snow depth, peak load, and their relationship to rainfall precipitation at that time will be clarified. Further, the amount of this rainfall will be compared with the rainfall distribution during the winter period to examine the regions where it is

<figure><img src="image_306.png" /><figcaption>**Figure 16.** Load change over time owing to rainfall (example of snow depth at 50 cm) [5].</figcaption></figure>
necessary to consider the effect of rain on the roof snow load. In addition, the roof-snow-load calculation method considering rain on snow will be examined.

**Acknowledgments**

This study was conducted as a part of 2014 Building Standard Provision Promotion Projects (S17: Evaluation of Roof Snow Load Considering Rain-on-Snow Surcharge) by the Japanese Ministry of Land, Infrastructure, Transport and Tourism, and the authors would like to express their appreciation for the funding.

**References**

1. ASCE, *Minimum Design Loads for Buildings and Other Structures* (ASCE/SEI 7-10) (2010, 2013)

2. *National Building Code of Canada* (NRCC53301) (2010)

3. Architectural Institute of Japan, *Recommendations for Loads on Buildings* (2015)

4. *The Building Standard Law of Japan* (2013)

5. M. Otsuki et al., *Evaluation of Roof Snow Load Considering Rain-on-Snow Surcharge; Summaries of JSSI & JSSE Joint Conference on Snow and Ice Research*, 283 (2015)
# Structural damage caused by rain-on-snow load in Japan

Toru Takahashi<fnref n="1" />,ᵃ, Takahiro Chiba<fnref n="2" /> and Kazuki Nakamura<fnref n="3" />

**Abstract.** In February 2014, heavy snow and additional rain led many large span structures and carports to collapse in Kanto region, Japan. We had almost the same depth of heavy snow one week before the event, but no structural damage was reported. Therefore, the authors focused on precipitation data and estimated roof snow load for each structure that collapsed time was clarified. As a result, it was clarified that designers would better to use precipitation data for design snow load.

## 1 Foreword

In February 2014, we got heavy snow twice with interval of one week in Japan. Both of snowfall was caused by low-pressure system that approached over the south coast of Japan, but only the second one led many large span structures to collapse.

In this paper, the authors checked time history data of weather stations and estimated the snow load for each structure that collapsed time was clarified. Then, problem of ordinary design method for snow load in Japan was discussed.

## 2 Outline of damage

Table 1 shows outline of damaged buildings during the disaster and Figure 1 shows distribution of collapsed large span structures in Kanto region. This was the first time in Kanto region to have such a number of structures' collapse. Besides, large numbers of carports with lightweight and simply supported as cantilever were collapsed but the authors do not discuss the issue.

<table><caption>Table 1. Outline of snow disaster in February 2014 in Japan<fnref n="1" />)</caption><thead><tr><th rowspan="2">Prefecture</th><th colspan="3">residential house</th><th colspan="2">non-residential buildings</th></tr><tr><th>fully damaged</th><th>half damaged</th><th>partially damaged</th><th>public buildings</th><th>other</th></tr></thead><tbody><tr><td>Hokkaido</td><td></td><td></td><td>2</td><td></td><td>1</td></tr><tr><td>Aomori</td><td></td><td></td><td>2</td><td></td><td></td></tr><tr><td>Iwate</td><td></td><td></td><td>3</td><td>1</td><td>7</td></tr><tr><td>Miyagi</td><td></td><td></td><td>25</td><td>6</td><td>60</td></tr><tr><td>Akita</td><td></td><td></td><td></td><td></td><td>1</td></tr><tr><td>Fukushima</td><td>1</td><td></td><td>3</td><td>7</td><td>34</td></tr><tr><td>Tochigi</td><td></td><td>2</td><td>24</td><td></td><td></td></tr><tr><td>Gumma</td><td>2</td><td></td><td>27</td><td></td><td>16</td></tr><tr><td>Saitama</td><td></td><td>1</td><td>19</td><td>3</td><td>45</td></tr><tr><td>Tokyo</td><td></td><td>1</td><td>76</td><td>20</td><td>72</td></tr><tr><td>Kanagawa</td><td></td><td></td><td>21</td><td>2</td><td>41</td></tr><tr><td>Yamanashi</td><td>13</td><td>42</td><td>357</td><td></td><td></td></tr><tr><td>Nagano</td><td></td><td></td><td>14</td><td></td><td>24</td></tr><tr><td>Gifu</td><td></td><td></td><td></td><td></td><td>1</td></tr><tr><td>Shizuoka</td><td></td><td></td><td>4</td><td></td><td>39</td></tr><tr><td>Aichi</td><td></td><td></td><td>3</td><td></td><td>1</td></tr><tr><td>Kyoto</td><td></td><td></td><td>3</td><td></td><td>1</td></tr><tr><td>Oita</td><td></td><td></td><td>2</td><td>1</td><td>4</td></tr><tr><td>Miyazaki</td><td></td><td></td><td></td><td></td><td>1</td></tr></tbody></table>

As shown in Figure 1, it looks like there is a zone of collapse. There seems to have had heavier roof snow load than design value in the zone.

<figure><img src="image_307.png" /><figcaption>Figure 1. Distribution of collapsed structures in Kanto region.</figcaption></figure>

<table><caption>Table 2. List of structure whose collapsed time was clarified.</caption><thead><tr><th>No.</th><th>Date</th><th>Time</th><th>Place</th><th>Damaged system</th></tr></thead><tbody><tr><td>1</td><td>14th Feb.</td><td>23:00</td><td>Hachioji</td><td>Arcade</td></tr><tr><td>2</td><td>15th Feb.</td><td>05:10</td><td>Yamato</td><td>Diamond truss</td></tr><tr><td>3</td><td>15th Feb.</td><td>05:30</td><td>Kumagaya</td><td>Membrane</td></tr><tr><td>4</td><td>15th Feb.</td><td>05:30</td><td>Yokohama</td><td>Prefab. truss</td></tr><tr><td>5</td><td>15th Feb.</td><td>06:00</td><td>Mitaka</td><td>Parallel truss</td></tr><tr><td>6</td><td>15th Feb.</td><td>07:00</td><td>Yokohama</td><td>Cantilever roof</td></tr><tr><td>7</td><td>15th Feb.</td><td>08:00</td><td>Fujimi</td><td>Parallel truss</td></tr><tr><td>8</td><td>15th Feb.</td><td>09:00</td><td>Yamanashi</td><td>Parallel truss (old)</td></tr><tr><td>9</td><td>15th Feb.</td><td>09:15</td><td>Ome</td><td>Full web frame</td></tr><tr><td>10</td><td>15th Feb.</td><td>09:30</td><td>Nikko</td><td>Parallel truss</td></tr><tr><td>11</td><td>15th Feb.</td><td>10:20</td><td>Fukaya</td><td>Unknown (factory)</td></tr></tbody></table>

<figure><img src="image_308.png" /></figure>

<footnote n="1">¹Department of Architecture, Chiba University, 1-33 Yayoicho, Inage-ku, Chiba, Japan</footnote>

<footnote n="2">²Department of Architecture, Hokkaido University of Science, 7-15-4-1, Maeda, Teine-ku, Sapporo, Japan</footnote>

<footnote n="3">³National Research Institute for Earth Science and Disaster Prevention, 3-1 Tennodai, Tsukuba, Japan</footnote>

<footnote>ª Corresponding author: takahashi.toru@faculty.chiba-u.jp</footnote>
Table 2 shows the list of collapsed large span structure whose collapsed time was clarified²). The snowfall began in Kanto region in the morning of February 14th, 2014 and changed to rain in the next morning. Structure No.1 collapsed midnight of February 14th and the others collapsed in the morning of February 15th. Photos 1, 2, 3 and 4 show examples of collapsed buildings. Investigation reports of the accidents say these structures were safe for design snow load level (600 N/m²)³). In the figure, collapsed time shown in table 2 and expected accumulated precipitation i.e. snow load on the roof are plotted with triangle mark. Because 1.0 mm of precipitation means 9.8 N/m² of roof snow load if its shape factor is 1.0, expected snow load on the roof might be heavier than design snow load level and it might be 800 to 900 N/m² as shown in Figure 2, except of structure No.1. In Figure 2, catch ratio of precipitation gauge was corrected after Yokoyama et al⁴). In the figure, it can be seen that snow load on the roof might exceed at 4 o'clock in the morning at Tokyo, and after that, additional rain-on-snow load might be affected.

This time series investigation results of roof snow loads are consistent with results of structural analyses in the accidents investigation reports of collapsed structures.

<figure><img src="image_309.png" /><figcaption>**Photo 1.** Example of collapsed building. (No.7 in Table 1.)</figcaption></figure>

<figure><img src="image_310.png" /><figcaption>**Photo 2.** Internal view of Photo 1.</figcaption></figure>

<figure><img src="image_311.png" /><figcaption>**Photo 3.** Example of collapsed building. (No.9 in Table 1.)</figcaption></figure>

<figure><img src="image_312.png" /><figcaption>**Photo 4.** Internal view of Photo 3.</figcaption></figure>

## 3 Statistical data study

Conventionally, in Japan, roof snow load of the general area, which is defined by building code, is calculated by multiplying the annual maximum snow depth of return period 50 years in the region and the unit snow weight (2.0 kN/m³). By the damage of this heavy snow accident, the influence of rainfall has become clear. Therefore, in the following clause, the authors verify the influence of rainfall using a method of snow layer model, which has been proposed by the authors to estimate the ground snow depth and the ground snow weight by using precipitation and temperature⁶), and extreme value statistics.

### 3.1 Longitudinal data

Figure 3 shows longitudinal data of annual maximum snow depth (left hand) and annual maximum snow weight on the ground (right hand), which were estimated by the snow layer model, for Tokyo, Yokohama, Kumagaya and Kofu from 1962 to 2014 seasons. In the figure, design levels of snow depth and snow weight are also shown for the convenience.

When we focus attention on snow depth, the latest value was not the largest one, for Tokyo and Yokohama. In the case of Tokyo, snow depth in 2014 (0.27 m) was second deepest value and was less than design snow
depth (0.3 m). On the other hand, at Kumagaya and Kofu, snow depth in 2014 was enormously high. The value reaches almost twice of design value. If we consider only this value, it is difficult to describe why no structure (except of the deteriorated old one) was collapsed around there.

Therefore, next, the authors evaluated annual maximum ground snow weight by snow layer model and compared them with snow depths. Even in Tokyo and Yokohama, the value of 2014 has become the largest. In Kofu, there appeared large value in 1998. There was a snow disaster in 1998 around Kofu region⁷). The values are consistent with structural damage.

## 3.2 Extreme value statistics

Figure 4 shows plots of annual maximum ground snow depth (left hand) and annual maximum ground snow weight (right hand) for Tokyo, Yokohama, Kumagaya and Kofu on the Gumbel probability paper. In the figure, design values of snow depth and snow weight are also shown for the convenience. Triangle plots are data based on 1962 to 1999, which was used in the AIJ recommendation 2004⁸) version, and circle plots are data based on 1962 to 2014⁹).

Comparing with these plots, the effect of 2014 data was clarified. For snow depth of Tokyo and Yokohama, there is no difference between both regression formula,

can also be understood that the snow depth 0.3 m that has been used in the design is a value that corresponds to return period 50 years. In Kumagaya, although 0.3 m was a value exactly corresponding in recurrence interval 50 years in the data by 1999, by regression including the data in 2014, it would be the evaluation 17.8 years. On the other hand, if observed value 0.62 m in 2014 was applied to this regression, it would become the evaluation of recurrence interval 272 years. In Kofu, although 0.5 m was a value exactly corresponding in recurrence interval 74 years in the data by 1999, by regression including the data in 2014, it would be the evaluation 18.2 years. On the other hand, if observed value 1.14 m in 2014 was applied to this regression, it would become the evaluation of recurrence interval 353 years.

When attached to the ground snow weight, as long as it has been based on the data by 1999, it has turned out that using 0.6 kN/m² for a design value of recurrence interval 50 years in Tokyo, Yokohama, and Kumagaya has been almost comparable. However, if based on the data by 2014, the recurrence interval that 0.6 kN/m² corresponds would become 24 years in Tokyo, 26 years in Yokohama and 16 years in Kumagaya. On the other hand, in Kofu, if based on the data by 1999, 1.0 kN/m² corresponds in 138 years, and even if it uses the data by 2014, it would become 67 years. This result was also considered as the circumstantial evidence that the damage of large span structure was little in Yamanashi.

<figure><img src="image_313.png" /><figcaption>**Figure 2.** Example of weather data during the disaster⁵).</figcaption></figure>
<figure><img src="image_314.png" /><figcaption>**Figure 3.** Longitudinal data of annual maximum snow depth and annual maximum snow weight on the ground.</figcaption></figure>
<figure><img src="image_315.png" /></figure>

<figure><img src="image_316.png" /></figure>

<figure><img src="image_317.png" /></figure>

<figure><img src="image_318.png" /></figure>

<figure><img src="image_319.png" /></figure>

<figure><img src="image_320.png" /></figure>

<figure><img src="image_321.png" /></figure>

<figure><img src="image_322.png" /><figcaption>Figure 4. Plots on Gumbel probability paper of annual maximum snow depth and annual maximum snow weight on the ground.</figcaption></figure>
## 4 Equivalent snow density

Figure 5 shows the relation between annual maximum snow depth and equivalent unit snow weight (snow density). Ordinary Japanese design sequence (building code) needs snow depth first, then requests equivalent unit snow weight. 2.0 kN/m³ is normally used in general region and 3.0 kN/m³ is used for heavy snow region. AIJ (Architectural Institute of Japan) recommends equation 1 for the equivalent unit snow weight $p_s$ [kN/m³] based on observations in heavy snow region⁸⁹⁾.

<figure><img src="image_323.png" /><figcaption>Figure 5. Relation between annual maximum snow depth and equivalent unit snow weight.</figcaption></figure>

$$p_s = 0.72\sqrt{d/d_{ref}} + 2.32 \quad (1)$$

where, $d$ [m] is snow depth on the ground and $d_{ref}$ [=1 m] is standard snow depth for make it non-dimensional number. Figure 5 shows that snow density for small snow depth could be very large based on precipitation data. In such a case, even if we use equation 1, it could not be expounded. It seems better to use precipitation data directly. For 2014 value, it could be said that the ground snow weight might be around 1.0 kN/m². It includes rain-on-snow surcharge. The value exceeded design load in Tokyo, Yokohama and Kumagaya, and led to structural damage. On the other hand, in Kofu, it was almost even with design load, and little damage for large span structure was reported.

## 5 Results

Almost of structural damage in 2014 heavy snow might occur because snow load on the roof exceeded their design snow load. It could be said that total precipitation during the event largely 1.0 kN/m² in Kanto region. It could be said that we would be better to use precipitation data directly to decide design snow load value that including rain-on-snow surcharge.

In Kofu, it was happy that design snow load had sufficient margin (long return period) against actual load, and it led little damage for large span structures.

## References

1. Fire and Disaster Management Agency: Heavy Snow Disaster Report during February 14 to 16, (2014) http://www.fdma.go.jp/ (in Japanese)

2. Asahi-shimbun newspaper, Feb. 17th, 2014, Nikkei-shimbun newspaper, Feb. 17th, 2014, etc. (2014)

3. For example, Investigation Board on Collapse of Fujimi Arena: Investigation Report on Collapse of Fujimi Arena (2014) (in Japanese)

4. Yokoyama K., Ohno H., Kominami Y., Inoue S. and Kawakata T.: Performance of Japanese Precipitation Gauges in Winter, Seppyo (Journal of Snow and Ice), 65, 303-316 (2003)

5. Japan Meteorological Agency: Observations at each weather station (as of February 2014) http://www.data.jma.go.jp/obd/stats/etrn/index.php

6. Takahashi, T., Kawamura, T., Kuramoto K.: Estimation of Ground Snow Load Using Snow Layer Model, Journal of Structural and Construction Engineering (Transactions of AIJ), 66, No.545, 35-41 (2001)

7. Takahashi, T., Kawaguchi, K., Ohi, K.: Snow Disaster in Yamanashi Prefecture in January, 1998, Journal of Snow Engineering of Japan, 14, 141-144 (1998)

8. Architectural Institute of Japan: Recommendations for Loads on Buildings (2004), AIJ, Tokyo (2004)

9. Architectural Institute of Japan: Recommendations for Loads on Buildings (2015), AIJ, Tokyo (2015)
# A new method of predicting slide snow load for sloped roofs

Xuanyi Zhou¹, a, Jialiang Li¹, Peng Huang¹, Ming Gu¹ and Lulu Sun¹

¹State Key Laboratory of Disaster Reduction in Civil Engineering, Tongji University, Shanghai, China

**Abstract.** This study applies a new method to simulate sliding snow load on sloped roofs based on the snowmelt model for two regions in China. The positive energy absorbed by snowpacks on sloped roofs is regarded as a key indicator of snow sliding based on a previous snowmelt model for building roofs and based on field observation results from the previous literature. The impacts of roof slope, shielding effect of neighbouring buildings and heat gained from within buildings are analysed. And, this study presents a simplified formula of sliding snow load as a function of roof slope, shielding effect of neighbouring buildings, and heat transfer coefficients for convenient application of structural engineers.

## 1 Introduction

Snow can accumulate easily on the roofs of buildings in cold regions. When the weather warms up, snowpack often becomes unstable sliding drops. When two buildings with different heights are adjacent to each other, the snowpack sliding from the upper roof can aggravate the snow load on the lower one. The prediction of slide snow/ice loads on roofs is vital for structure design in heavy snowfall regions. Some studies have been conducted on the snow-slide loads on roofs. Sack et al. (1987) [1] developed a method to calculate the accumulation of snow on unobstructed, slippery roofs and pointed out that snow sliding on the roof—snow interface is initiated by either a reduction in resistance or shear failure. Lepage and Schuyler (1988) [2] presented a numerical simulation method that is combined with a scale model test to predict the quantity of snow that could eventually slide off sloped roofs. The simulation method considered several years of historical observations of wind, temperature, precipitation, and solar radiation. The scale model tests included a water flume and a wind tunnel experiment. Takakura et al. (2000) [3] reported that snowslide is affected not only by roof slope, but also by the temperature characteristics in each region. Based on the available historical weather data, Isyumov and Mikitiuk (2008) [4] developed a numerical method to predict the frequency of snow sliding from the surfaces of sloped buildings. In this method, 0 °C was the prerequisite air temperature in snowpack sliding. To explore the relationship between temperature characteristics and snow-slide, Chiba et al. (2012) [5] presented the snow-slide probabilities of different roof slopes as a function of the accumulated temperature above freezing. However, most of them do not provide a practical method for structural engineers. Current load codes and Zhou et al. (2013) [6] established calculation methods that disregard some crucial effects of key factors. For example, these load codes generally regard the slide snow/ice load coefficient as fixed. Zhou et al. (2013) considered the effects of heat obtained from within buildings and the shielding effect of neighbouring buildings, but they ignored the effects of roof slope. Further, Zhou et al. (2015) [7] developed a new method to simulate sliding snow load on sloped roofs. The

positive energy absorbed by snowpacks on sloped roofs was regarded as a key indicator of snow sliding based on a previous snowmelt model for building roofs and based on field observation results from the previous literature. In this study we adopted the new method developed by Zhou et al. (2015) to analyse the characteristics of slide snow in some regions in China

## 2 Snowmelt model on building roofs

The snowmelt model of Zhou et al. (2015) includes the energy and mass balance equations.

The energy balance equation can be described as follows,

$$ \frac{dU}{dt} = S_n + L_a - L_t + H + E_l + Q_p + Q_r \quad (1) $$

where, $S_n$ denotes the solar radiation absorbed by snow on roofs; $L_a$ the atmospheric longwave radiation, $L_t$ the terrestrial longwave radiation; $H$ sensible heat; $E_l$ latent heat; $Q_p$ the heat advected by precipitation (rainfall or snowfall); and $Q_r$ the heat obtained from within the building. The unit of the energies above is $kJm^{-2}hr^{-1}$.

The mass balance equation in the snow layer (Zhou et al., 2015) can be expressed as follows,

$$ \frac{dW}{dt} = P_r + P_s - M_{out} - W_E \quad (2) $$

where $P_r$ denotes the precipitation as rainfall; $P_s$ the precipitation as snowfall; $M_{out}$ the rate of snowmelt; $W_E$ the sublimation or condensation from the snowpack. The unit of the variables on the right side of Equation (2) is $mhr^{-1}$.

## 3 Method of predicting slide snow load for sloped roofs

The concept of accumulated temperature above freezing was proposed by Chiba et al. (2012) [5]. Chiba et al.

<figure><img src="image_324.png" /></figure>

<footnote>a Corresponding author: zhouxytj@tongji.edu.cn</footnote>
believed that the data of the accumulated temperature above freezing follows log-normal distribution and summarized the relationship between snow-slide probability and the accumulated temperature above freezing. In Zhou et al.'s (2015) study, snow sliding is closely related to the internal physical states of snowpack, which is determined by the absorption of energy by snowpack on roofs. The energy absorbed by the snowpack reflects the essence of snow slide more effectively than air temperature does. In especial, the positive energy absorbed from the surrounding environment can alter snow phases and induces further snow slide. The concept of the positive energy that is shadowed is shown in Fig. 1.

<figure><img src="image_325.png" /><figcaption>Fig. 1. Accumulated positive energy</figcaption></figure>

<table><caption>**Table 1** Comparison of sufficient condition between field observation and simulated result (Chiba et al., 2012; Zhou et al., 2015)</caption><thead><tr><th>Roof slope α (1)</th><th>Critical accumulated temperature above freezing (Chiba et al., 2012) (2)</th><th>Critical accumulated positive energy (simulated results) (3)</th><th>Critical accumulated temperature above freezing (simulated results) (4)</th></tr></thead><tbody><tr><td>1.4</td><td>33°C·hr</td><td>1.0 MJ/m²</td><td>33°C·hr</td></tr><tr><td>1.0</td><td>46°C·hr</td><td>1.3 MJ/m²</td><td>43°C·hr</td></tr><tr><td>0.7</td><td>72°C·hr</td><td>2.0 MJ/m²</td><td>70°C·hr</td></tr><tr><td>0.5</td><td>110°C·hr</td><td>2.8 MJ/m²</td><td>115°C·hr</td></tr><tr><td>0.4</td><td>161°C·hr</td><td>4.1 MJ/m²</td><td>165°C·hr</td></tr><tr><td>0.3</td><td>314°C·hr</td><td>8.8 MJ/m²</td><td>310°C·hr</td></tr><tr><td>0.2</td><td>480°C·hr</td><td>15.0 MJ/m²</td><td>482°C·hr</td></tr></tbody></table>

Chiba et al. (2012) built several roof models with different slopes in Sapporo, Japan. Field observations were performed from 2001 to 2008. For every snow-siding event, the accumulated temperature above freezing was analyzed based on eight years' data of field observation. Chiba believed that the data of the accumulated temperature above freezing follows log-normal distribution and summarized the relationship between snow-slide probability and the accumulated temperature above freezing. Column 2 of Table 1 lists the accumulated temperatures above freezing under a 95% sliding probability for different roof slopes. The accumulated temperature above freezing can be regarded as a sufficient condition for snow sliding on sloped roofs because the values correspond to a high sliding

probability. Thus, they are labeled in Table 1 as critical accumulated temperature above freezing. The assumed critical values in Column 3 of Table 1 are adjusted through trial and error until the differences between Column 2 and Column 4 are minimized. After continuous trials, the preset critical accumulated positive energy depicted in Column 3 reduces the differences between Columns 4 and 2, which suggests that Column 3 of Table 1 can reflect the results of the field observation in Sapporo. The critical accumulated positive energy increases significantly as roof slope decreases; thus, snowpacks on low-sloped roofs can be expected to slide down with more difficulty. In the following simulations, the values of critical accumulated positive energy as listed in Column 3 of Table 1 are applied as the thresholds of sufficient condition for snow sliding on differently sloped roofs.

# 4 Simulation of sliding snow loads for two regions in China

The following part would apply the method described above to predict the slide snow loads for sloped roofs in China.

## 4.1 Effect of geographic differences on the energy exchange

This study analysed the meteorological data of several representative regions with snowfall in Jiamusi and Nanchang, which represent the northern region and southern region of China, respectively. The meteorological data of the two regions were collected from the China Meteorological Data Sharing Service System and included the daily average wind speed, average temperature, daily maximum temperature, daily minimum temperature, relative humidity, and amount of snowfall or rainfall during 62 winters from November 1, 1951 to March 31, 2013. Figure 2 depicts the meteorological data on the representative regions in a certain winter. All meteorological data include four measured values at some fixed time every day and the data were linearly interpolated to obtain the hourly data. As shown in Figure 2, Jiamusi experiences evenly distributed snowfall and little rainfall in winter, whereas Nanchang encounters concentrated snowfall and heavy rainfall. Moreover, the winter temperature in Jiamusi mainly ranges between -30 °C and -10 °C, whereas that in Nanchang usually exceeds 0 °C. Even during snowfall, the temperature is greater than -5 °C in Nanchang. In addition, wind speed in Jiamusi is slightly higher and fluctuates more than that in Nanchang. The humidities of both regions are similar, although that in Nanchang fluctuates more.

Figure 3 and Figure 4 exhibit the energy absorbed/released by snowpacks on roofs in Jiamusi and Nanchang during a certain winter. These figures were drawn based on the hypothesis that the snow on the roof does not slide down. Thus, they can reflect the entire process of energy absorption/release in winter. In these figures, short wave radiation Sn, sensible heat, and heat
transfer are positive. Net longwave radiation is generally negative, but it may also be positive on occasion. The latent heat and energy from precipitation are positive at times and negative at others. Their absolute values are smaller than those of other energies. Radiation energy including short wave and longwave radiation dominates the energy, followed by that obtained from sensible heat and the heat transfer from within buildings. The remaining energies only contribute a small proportion. The proportion of Sn in total energy is lower in Jiamusi than in Nanchang mainly because of the low latitude of Nanchang, which enhances solar radiation in Nanchang. The proportion of the net longwave radiation is lower in Nanchang compared with that in Jiamusi. This is because southern China experiences rainy days in winter more frequently than northern China, resulting in an increase in cloudy sky emissivity and decrease in the absolute value of net longwave radiation.

<figure><img src="image_326.png" /></figure>

<figure><img src="image_327.png" /></figure>

<figure><img src="image_328.png" /></figure>

<figure><img src="image_329.png" /><figcaption>Figure 2. Meteorological data in a certain winter</figcaption></figure>

<figure><img src="image_330.png" /><figcaption>Figure 3. Energy absorbed/released in a certain winter of Jiamusi</figcaption></figure>

<figure><img src="image_331.png" /><figcaption>Figure 4. Energy absorbed/released in a certain winter of Nanchang</figcaption></figure>

Figure 5 depicts the air and snow temperatures in Jiamusi and Nanchang. It can be seen that snow temperatures variation differs from that of external air temperature to a certain extent. This result can be attributed to the influence of radiation, wind speed, and humidity on snow energy content, although air temperature does affect it significantly. Throughout winter, the air temperature in Jiamusi is significantly lower than 0 °C according to Figure 5, which facilitates the accumulation of snow on roofs. The mean temperature of snow fluctuates with air temperature. Snow is a poor heat conductor. When external temperature is lower than 0 °C, snow temperature is relative higher for most of time. From the figure, Nanchang has a short snow accumulation period for its winter temperature is high.
<figure><img src="image_332.png" /></figure>

<figure><img src="image_333.png" /><figcaption>**Figure 5.** Air temperature and snow temperature in a certain winter</figcaption></figure>

## 4.2 Characteristics of snow sliding

The developed method is used to simulate the slide snow loads on several sloped roofs. The differences in snowslide mechanisms between Jiamusi and Nanchang were analysed. Fig. 6 and Fig. 7 present the simulated results for accumulated positive energy and events of snow sliding in two representative regions. From Figure 6 and Figure 7, we can compare the snow-sliding events in Jiamusi and Nanchang under the same roof slope. The sliding events in Jiamusi were found more evenly distributed than those in Nanchang. And these sliding events of Jiamusi were uncorrelated with snowfall. However, the sliding events in Nanchang were relatively concentrated and mainly occurred several days after snowfall. After snowfall, the air temperature in Nanchang generally increases, which induces the snow-slide. In addition, high-sloped roofs easily facilitated snow-sliding events so that the mass of sliding decreases significantly. The phenomenon corresponds to the proposed sufficient condition on sloped roofs, in which roofs with high slopes facilitate sliding events for small accumulated positive energy is needed.

As roof slope increases, sliding events remain still concentrated in Nanchang compared with those in Jiamusi. This difference is not only related to the concentrated snowfall in Nanchang, but also to their varied positive energy absorptions. According to Figure 6 and Figure 7, the positive energy absorption of snow after snowfall is much more in Nanchang than that in Jiamusi. Thus, snow accumulates for only a short period in Nanchang, whereas Jiamusi enjoys a long accumulation period.

<figure><img src="image_334.png" /><figcaption>**Figure 6.** Accumulated positive energy and events of snow sliding in Jiamusi</figcaption></figure>

<figure><img src="image_335.png" /><figcaption>**Figure 7.** Accumulated positive energy and events of snow sliding in Nanchang</figcaption></figure>

## 5 Simplified formula of sliding snow load coefficient

The sliding snow load coefficient is defined as the ratio of the snow load that slides from the upper roof to the total snow load on the upper roof. Figure 8 depicts the effect of roof slope on the sliding snow load coefficient. The heat transfer coefficient of the roof is 0.4 W/m²·k, and PSA is 0. PSA is defined as the ratio of the blocked area on a roof to the entire roof area. The coefficients obtained from the Canadian code (National Research Council of Canada, 2005), European code (British Standard Institution, 2003), and ASCE standard (American Society of Civil Engineers, 2010) are constant at 0.5, 0.5 and 0.4, respectively. The sliding snow load coefficients decrease quickly with the increase in roof slope when roof slope is lower than 0.5 (26.6°); when roof slope is larger than 0.5 (26.6°), these values decrease gradually. In general, the sliding snow load coefficients of both the northern and southern regions vary significantly with the change of roof slope.

<figure><img src="image_336.png" /><figcaption>**Figure 8.** Sliding load coefficient as a function of slope</figcaption></figure>
Given the differences in the sliding snow load coefficients between the northern and southern regions of China, the following simplified formula was obtained [3]. For northern regions,

$$ \begin{cases} (0.45-0.35 \times \alpha) \times (1.35+0.15 \times \text{PSA}) \times (1.7-0.3 \times K) & 0 \le \alpha \le 0.5 \\ (2.88-0.25 \times \alpha) \times (0.43+0.02 \times \text{PSA}) \times (0.5-0.1 \times K) & \alpha > 0.5 \end{cases} \quad (3) $$

For southern regions,

$$ \begin{cases} (0.4-0.34 \times \alpha) \times (1.8+0.4 \times \text{PSA}) \times (1.15-0.2 \times K) & 0 \le \alpha \le 0.5 \\ 0.23 \times (1.8+0.4 \times \text{PSA}) \times (1.15-0.2 \times K) & \alpha > 0.5 \end{cases} \quad (4) $$

where $\alpha$ is the roof slope, which is defined as the change in the vertical direction divided by the corresponding change in the horizontal direction; PSA is the shield rate at which adjacent buildings shield the upper roof; K is the heat transfer coefficient of the roof (W/m²k).

## 6 Conclusions

This study applies a new method to simulate sliding snow load on sloped roofs based on the snowmelt model for building roofs presented by Zhou et al. (2015). The simulation of sliding snow loads for two regions in China was performed and their characteristics of representative regions were analyzed. The sliding events in the northern regions, which were uncorrelated with snowfall, were more evenly distributed than those in the southern regions. However, the sliding events in the southern regions were relatively concentrated and mainly occurred several days after snowfall. The distribution of these events in the southern regions remained concentrated even with an increase in roof slope.

## Acknowledgements

This project is jointly supported by the National Natural Science Foundation (51278368, 51478359), the Ministry of Science and Technology of China (SLDRCE14-B-10) and the Fundamental Research Funds for the Central Universities, which are gratefully acknowledged.

## References

1. Sack RL, Arnholtz D, Haldeman JS, Sloped roof snow loads using simulation. Journal of Structural Engineering, ASCE, 113 (8): 1829-1833 (1987)

2. Lepage MF, Schuyler GD, A Simulation to Predict Snow Sliding and Lift-off on Buildings, Proceedings of the Engineering Foundation Conference on a Multidisciplinary Approach to Snow Engineering. Santa Barbara, California (1988)

3. Takakura M, Chiba T, Ito T, Tomabechi T, On the Term Controlled by Snow Load: Practical Use of Snow Sliding on Pitched Roof. Journal of Structural and Construction Engineering, Transactions of AIJ 528: 53-57. (2000)

4. Isyumov N, Mikitiuk M, Sliding Snow and Ice from Sloped Building Surfaces: Its Prediction, Potential Hazards and Mitigation. The 6th Snow Engineering Conference, Whistler B.C., Canada (2008)

5. T. Chiba, T. Tomabechi, T. Takahashi, Study on evaluation of snow load considering roof snow-slide on gable roofs. *Snow Engineering VII*, Fukui, Japan, 231-241, (2012).

6. X.Y. Zhou, Y.Q. Zhang, M. Gu, J.L. Li, Simulation Method of sliding Snow Load on Roofs and its application in Some Representative Regions of China. Nat. Hazards 67, 2: 295-320, (2013).

7. X.Y. Zhou, J.L. Li, M. Gu, L.L.Sun, A new simulation method on sliding snow load on sloped roofs, Nat. Hazards, 77, 1: 39-65, (2015).
# Experimental study of the distribution of snow deposits on the surface of structures with complex three-dimensional shape of the roof

Poddaeva Olga¹, Pavel Churin¹ª

¹Moscow State University of Civil Engineering (Mgsu) Yaroslavskoe shosse, 26, Moscow, 129337, Russian Federation

**Abstract.** Currently, due to the modern trends of architecture in our country unique buildings such as trade, exhibition and business centres, sports arenas, airport complexes, and others are being actively constructed. One of the key features of such structures is a complex spatial geometric shape of the roof, which together with the peculiarities of the Russia climate make the need for accurate accounting snow load at design significantly increase. This article describes an experimental study of the snow deposits distribution on the structures surface with complex three-dimensional roof shape by the example of the airport complex in the Rostov-on-Don city, which is being built for the World Cup.

## 1 Investigated object and problem analysis

A site is located at Russian Federation, Rostov region. A new terminal construction is associated with the air network expansion and provides international passenger transport services.

As it seen from Figure 1, the most difficult and interesting airport terminal part in terms of aerodynamics and the snow load studying is the roof as it has a complex three-dimensional shape with a bending in all the directions, the exact reproduction of such forms is possible only by using the advantage of 3D printing.

<figure><img src="image_337.png" /><figcaption>**Figure 1.** Airport complex in the city of Rostov-on-Don</figcaption></figure>

Snow cover consists of solid precipitation in winter and ice crust and the water formed during the thaw. Depending on climatic conditions, snow can have a

different thickness and density. The product of thickness and density gives the snow mass on area unit, i.e. snow load on the ground that defines a basis for snow load on the roof.

The normative document that regulates the accounting impact of snow on structures in the Russian Federation is the SP 20.13330.2011 "Loads and effects"[1]

Statistical methods allow establishing design values of snow load on the ground by the results of meteorological observations. Snow load on the roof is formed by a number of additional factors, them interconnection is needed to be taken inter account. Full snow load acts on buildings and structures, where snow accumulates throughout the winter. Let's consider the four main factors that determine the value of the snow load on the buildings surfaces:

- amount of falling in the winter solid precipitation (details are provided in Annex G SP 20.13330.2011);

- snow falling from the inclined surfaces;

- melting snow on the surfaces of heated buildings

- transport of snow as a result of wind action, leading to non-uniform deposits on the coating surface and the demolition of some fallen snow from a coating.

Moving snow under the wind influence is a major factor, due to which the level of the snow load is not the same for the cover area. As a result of snow transport the places where snow imposed and places with snow accumulation (the so-called snow bags) arise. The location of these sites is due primarily to the configuration of the surface coating (his profile) and with the wind direction on the other hand. These phenomena are studied by simulation in the wind tunnel, as well as through natural observations of snow deposits on real roofs. The research results are displayed in the norms given load factor of transition from the ground to the load to the load on the cover.

Schemes of snow deposits distribution which listed in Annex G to SP 20.13330.2011 (as, indeed, in other standards, such as Eurocodes [2]) cover only a limited set of the most common forms of coatings, while often problems arise on account of the snow load on other forms of roofing or clarify the data provided in norms due

<figure><img src="image_338.png" /></figure>

<footnote>ª Corresponding author: pashok @inbox.ru</footnote>
to the need of specification of bearing ability. It is necessary to take into account two factors that have been used in the preparation contained in the norms unified design schemes of loading the snow:

- it assumed approximately uniform distribution of the winds from the characteristic directions;

- it used the balance ratio, which consists in the fact that the amount of snow in a snow bag equal to the volume of snow, demolished with the overlying section.

Based on these assumptions and the provisions of the regulations set out below in difficult cases we can recommend an approach to the estimation of the maximum possible value of the snow bag, based on the completely filled it with snow.

The obvious unpredictability of the snow cover value has led to an understanding of its random nature and the need to use statistical methods for the design values of snow loads determination on the basis of natural observations and modeling in a wind tunnel.

Current regulations SP 20.13330.2011 "Loads and effects" were developed on the methodological basis, where as the raw statistical data used yearly highs weight of snow cover according to the data snow surveys at sites protected from the effects of wind. Full standard value of snow load on the horizontal projection of the cover is equal S = S₀μ, where S₀ is a standard value of snow cover weight per 1 m² of horizontal earth surface; μ is a conversion factor from the weight of snow cover on the ground surface to the snow load on the cover.

In the formation of snow load even on the covers single-span buildings there is a nonuniform deposition of snow on the bi-coatings, as well as deflation and falling from the snow surface coatings. These phenomena are taken into account by means of the coefficients μ listed in the regulations of most countries, depending mainly on the slopes of coatings and present differences of heights.

Determination of snow loads on structures nontraditional forms and/or unusual scales may lead to a number of new issues. A wind tunnel simulation is used for the snow load determination. Based on the results of physical experiment a zoning of complex roof for primitives that are presented in the regulations [1] is performed.

According to SP 20.13330.2011 - "Loads and effects." P. 10 "Snow load", p. 10.1, Note 2: for spatial coatings of complex geometric shapes and coatings having the largest characteristic dimension in the plane of more than 100 m, their snow load distribution scheme is necessary to determine with test data based on a specially developed recommendations.

## 2 Layout

For experimental studies a model of the passenger terminal of the airport complex was designed and manufactured. Considering the wind tunnel working part size, the best possible scale of layout was chosen that is 1: 250.

<table><caption>**Table 1.** Comparison of the main characteristic of the model and the size of the full-scale object.</caption><thead><tr><td>Overall dimension</td><td>Model</td><td>Full-scale object</td></tr></thead><tbody><tr><td>length, mm</td><td>1322</td><td>331200</td></tr><tr><td>width, mm</td><td>906,5</td><td>226700</td></tr><tr><td>height, mm</td><td>132</td><td>33213</td></tr></tbody></table>

The main requirement to the model for aerodynamic testing is a strict observance of the geometric similarity of the object under study. To meet this requirement at the stage of layout design based on the raw data provided by the customer the album design drawings was created.

Walls of model are made of plywood. Columns and roof support are made of a steel wire and aluminum plates (Figure 2).

<figure><img src="image_339.png" /><figcaption>**Figure 2.** Layout framework</figcaption></figure>

Elements of the roof due to the complex spatial geometric forms were made on 3D printer Cube Pro (Figure 3).

<figure><img src="image_340.png" /><figcaption>**Figure 3.** 3D Printer Cube Pro and roofing element.</figcaption></figure>

After the assembly of all the elements of the layout, the layout was painted matt black to enhance the results of photographic phenomena snow transport and snow deposits (Figure 4).

<figure><img src="image_341.png" /><figcaption>**Figure 4.** Painting layout.</figcaption></figure>
## 3 Test procedure

Research was carried out in the Boundary Layer Wind Tunnel of Moscow State University of Civil Engineering (Figure 5). This plant has a closed circulation loop and modular fan unit of the nine units, with the length of the working zone of 18.9 meters.

<figure><img src="image_342.png" /><figcaption>**Figure 5.** Boundary Layer Wind Tunnel of Moscow State University Of Civil Engineering</figcaption></figure>

The research method of the main regularities experimental determination of snow transport and snow deposits formation on roofing construction objects coverings is based on the following actions:

1. The model is rigidly mounted in the working area of the wind tunnel on the turntable;

2. Rotation angle of the relative to the direction of the oncoming flow is set to mark 0°;

3. The surface model is cleaned of extraneous elements such as dust, sawdust, wood flour residues and so forth.

4. The uniform layer of model material (brand 200 wood flour, bulk density – 200 kg/m³, humidity – 3,5-4%) with 1 to 1,5 mm thick is applied on model;

<figure><img src="image_343.png" /><figcaption>**Figure 6.** Model of the airport complex coated with a layer of modeling material</figcaption></figure>

5. Engines start of a wind tunnel is made. The flow velocity increases uniformly from 0.1 to 7.6 m/s After the beginning of the phenomenon emergence the transport model material the speed of a stream is established by a constant and the purge proceeds before formation of characteristic snow deposits. The control of the beginning transport moment and the moment of snow deposits characteristic forms formation is carried out visually;

6. The work of engine wind tunnel stops, the experimental group captures results using a camera;

7. Position model changes according to the experimental plan. Repeat steps 3-6.

Characteristic speeds of a wind at which investigated snow transport on model in the range from 3,4 to 8 m/s. Modeling of snow is carried out by wood flour with a size of particles from 50 micron to 250 microns at humidity of 3,5-4% which is displaced from a smooth colored surface at a speed of wind of 3,4 m/s [3]. As a result of a long exposition of the model covered with a thin layer of wood flour in an air stream with a speed of 6-7 m/s the picture of snow deposits is formed. An analysis of the snow bags forms was performed for the zoning of the roof at the elementary entities for which the standard calculation of snow load [4].

Tests were carried out on a model scale of 1: 250 (Figure 3). When the model was under blowing from different directions the wind effect on creating of snow deposits was simulated in accordance with the wind rose for the region. In the experiment complete circular exposure to the wind direction with a step 45° was performed. Pictures at different wind directions were proposed (Figure 7-8), that can be used for data obtaining about the snow deposits phenomenon, directions and volumes of transported snow.

<figure><img src="image_344.png" /><figcaption>**Figure 7.** Model of airport complex</figcaption></figure>

<figure><img src="image_345.png" /><figcaption>**Figure 8.** Typical areas of snow deposits, the angle of the incoming flow is 45°</figcaption></figure>
## 4 Analysis of the results

Analysis of snow deposits allows to select specific zones as a result of increased snow deposits. On the surface of the roof there are always areas with elevated snow deposits determined by wind direction. Any portion of the roof may be in the area of high snow deposits depending on wind direction.

On the basis of the data provided by the Customer and data of the carried-out climatic analysis, the characteristic directions of wind in the Rostov region during the winter period are eastern (the angle of statement of model is 315°) and western (the angle of statement of model is 135°). For those directions and also for the direction 0°, schemes of snow deposits at the airport complex roof surface are made.

<figure><img src="image_346.png" /><figcaption>**Figure 9.** The scheme of distribution of snow deposits, flow direction - south-east - 0°</figcaption></figure>

Analysis of the results shows that for all elements of the roofing there is a flow direction in which this element will carry the maximum of snow load. Based on this, for the calculation of the snow load is necessary to consider the full picture of the snow covered roof. Figure 6 shows the zones indicating the wind direction at which the element receives the maximum snow load.

<figure><img src="image_347.png" /><figcaption>**Figure 10.** Zones indicating the wind direction at which the element receives the maximum snow load</figcaption></figure>

## Acknowledgement

This work was financially supported by the Ministry of Education and Science of the Russian Federation within the framework of the state order to MGSU №2014/107 project "Fundamental studies of wind effects (including extreme) to the unique building and bridge structures"

## References

1. SP 20.13330.2011 - "Loads and effects."

2. EN 1991 Eurocode 1: Actions on structures

3. Xuanyi Zhou, Jinhai Hu, Ming Gu. Wind tunnel test of snow loads on a stepped flat roof using different granular materials // Natural Hazards. 2014. Volume 74, Issue 3, pp 1629-1648

4. E. Hjorth-Hansen, I. Holand, S. Loset, H. Norem Snow Engineering 2000: Recent Advances and Developments // CRC Press, 2000
# Theoretical and experimental study of ice accretion due to freezing rain on an inclined cylinder

Krzysztof Szilder

Aerospace, National Research Council, Ottawa, ON, K1A 0R6, Canada

**Abstract.** This paper reports on an investigation of ice formation on a cylinder exposed to freezing rain. Freezing rain icing is characterized by supercooled raindrop impingement, with vertical and horizontal velocity components, water flow over the surface driven by aerodynamic and gravity forces, and freezing or shedding of the surface water. In the first part of the paper, a novel, three-dimensional analytical model of ice growth on a cylinder is derived and the results are discussed. This is followed by a discussion of experimental tests in the NRC Altitude Icing Wind Tunnel. The experiments were constrained to horizontal drop impingement, and we encountered difficulties delivering large drops at low speed to the test section. Hence, the experiments were used primarily to validate the NRC morphogenetic ice accretion model. The morphogenetic model is a discrete element, random walk model that emulates the motion and freezing of individual fluid elements on the accretion surface. After validating the numerical methodology against the wind tunnel experiments, the model could be used to simulate ice formation under more realistic and general geometrical configurations and environmental conditions.

## 1 Introduction

Freezing rain consists of raindrops that become supercooled while passing through a sub-freezing layer of air. When these drops impact a surface at a sub-zero temperature, they partly freeze. The resulting ice can accumulate to a thickness of many centimetres during severe freezing rain events called ice storms. Although ice storms are not violent, they often cause power outages and tree damage, because of the weight of the ice accretions and wind-on-ice loads. In addition, the shape of ice accretions may alter the aerodynamics of transmission lines and bridge cables, causing large amplitude vibrations called galloping. Finally, pieces of ice falling from structures such as overhead transmission lines, wind turbines and telecommunication masts can cause serious injuries and property damage.

In order to better understand the physics of freezing rain and its consequences for engineered structures, more observations and controlled experiments are needed. In aeronautical icing wind tunnels, tiny drops are delivered to the test section by the horizontal flow and the drop trajectories are almost horizontal. In nature, however, freezing precipitation drops are larger and have a substantial vertical velocity component. Thus, there are limits to how accurately one can simulate ice formation due to freezing rain in a horizontal icing wind tunnel, designed for simulating in-flight icing. This discrepancy leads to different ice shapes in a wind tunnel and nature. In addition, in aeronautical icing wind tunnels, one encounters great difficulty delivering large drops characteristic of freezing rain to the test section. These differences lead to a different drop impingement domain and distribution than in the natural environment. Consequently, wind tunnel experiments can simulate only the “boundary conditions” of the environmental parameter space to which engineered structures are exposed.

Corresponding author: Krzysztof.Szilder@nrc-cnrc.gc.ca

Existing freezing rain ice accretion models are quite simplistic, often consisting of a couple of algebraic equations [1] and assuming simple two-dimensional geometry [2]. These models do not reflect the complexity of the ice accretion process. Instead of computing it, the ice accretion shape is assumed to have a simple geometrical form. Recognizing the scarcity of ice accretion experimental data [3-4] and the oversimplification of existing predictive models, the objective of this work is to develop a validated numerical model that can predict realistic ice accretion shapes under diverse freezing rain conditions. The focus of our current analysis is ice formation on a simple, but common, three-dimensional cylindrical geometry. Such a geometry can be found in a variety of engineered structures, such as transmission lines, telecommunication masts and bridge cables. In the future, the model can be readily extended to include more complex substrate shapes.

This contribution consists of three sections: development and analysis of an analytical ice accretion shape model that yields back-of-the-envelope estimates, experimental data collection in the National Research Council (NRC) Altitude Icing Wind Tunnel (AIWT), and validation of NRC's 3D morphogenetic icing model for experimental conditions.

<figure><img src="image_348.png" /></figure>
## 2 Analytical model of ice accretion on a cylinder in 3D space

In this section, we examine the ice mass, its shape and its location on a cylinder under freezing rain conditions, in a three-dimensional configuration. We assume that the freezing rain drops are large enough to follow straight-line trajectories. In addition, the impinging precipitation is assumed to freeze promptly upon impact. There is no surface flow of liquid. These assumptions give a zero-order model of impingement and the ice accretion processes. This model can help identify where the freezing rain impinges and estimate extreme ice shapes, sizes and loads, under general geometrical configurations. Comparison between the predictions of this simple model and experimental or numerical results from more complex models will allow identification of the influence of physical processes, such as surface water flow, gradual water freezing, water shedding, and the formation of rough ice surfaces or voids in the ice structure. Such processes might occur in nature but they are not included in this zero-order approximation.

### 2.1 Freezing precipitation with vertical and horizontal velocity components

The drop trajectory direction is given by the droplet velocity, $\vec{V}$, which is the vector sum of the drop terminal velocity, $V_z$, and the horizontal wind velocity with components $V_x$ and $V_y$, Fig. 1. The effect of turbulence and slower changes in wind speed and direction are neglected here but they could be incorporated in the future. The present values could be considered as mean values, and a fluctuating or time dependent term could be added to each velocity component.

In order to represent the most general configuration of the cylinder relative to the precipitation impingement direction, we establish the coordinate system shown in Fig. 1. Without loss of generality, the axis of the cylinder lies in the y-z plane, but it is inclined to the x-y plane at the inclination angle $\omega$. It varies from 0° to 90°. The angle between the drop velocity vector and the vertical is called the precipitation angle, $\epsilon$.

The impinging freezing rain mass per unit length of cylinder, m (kg/m), is given by:

$$m = P \sin \delta 2R \rho_w \quad (1)$$

where: P is the equivalent depth of the approaching precipitation, (m) (depth if it accumulated uniformly on a plane perpendicular to its velocity vector); $\delta$ is the approach angle between the drop impingement direction and cylinder axis; R is cylinder radius, (m); and $\rho_w$ is water density, (kg/m³).

<figure><img src="image_349.png" /><figcaption>Figure 1. A sketch of the general geometry for freezing rain impinging on a cylinder.</figcaption></figure>

The equivalent depth of precipitation P can be considered to be a vector with direction along the precipitation trajectory. It has vertical, $P_v$, and horizontal, $P_H$, components. The horizontal and vertical components are related by the precipitation angle:

$$P = \sqrt{P_v^2 + P_H^2} = \frac{P_v}{\cos \epsilon} \quad (2)$$
$$\text{where} \quad \tan \epsilon = \frac{\sqrt{V_x^2 + V_y^2}}{V_z}$$

It can be shown that the negative unit drop velocity vector, $-\vec{V}$, and unit cylinder axis vector, $\vec{C}_A$, are given by:

$$-\vec{V} = [\sin \epsilon \sin \gamma, \sin \epsilon \cos \gamma, \cos \epsilon] \quad (3)$$
$$\text{and} \quad \vec{C}_A = [0, \cos \omega, \sin \omega]$$

where the yaw angle, $\gamma$, is defined as the angle between the projection of the precipitation velocity vector and the projection of the cylinder axis onto a horizontal plane. It varies from 0° to 180°.

The dot product of these two unit vectors gives the approach angle:

$$\cos \delta = \sin \epsilon \cos \gamma \cos \omega + \cos \epsilon \sin \omega \quad (4)$$

Combining Equations (1) and (2) gives the expression for ice accretion mass as a function of precipitation and cylinder characteristics:

$$m = \lambda P_v 2R \rho_w \quad \text{where} \quad \lambda = \frac{\sin \delta}{\cos \epsilon} \quad (5)$$

where $\lambda$ is the overall impingement coefficient and the approach angle is given by Eq. (4). It should be noted that the impingement coefficient has the same value as the local collision efficiency at the stagnation point.
The cylinder area where ice forms can also be identified. Since the drop trajectories are assumed to be straight lines, impingement occurs on the half of the cylinder circumference that is exposed to drops. We will call the angle between the uppermost location on the cylinder cross-section, perpendicular to its axis, and the middle of the ice semicircle the maximum impingement angle, α, see Fig. 1. This angle is in a plane perpendicular to the cylinder axis and it should be also noted that ice accretion is homogenous along the length of the cylinder.

In order to find the maximum impingement angle, we examine the components of the negative drop velocity vector along the cylinder axis, $\vec{C}_A$, and perpendicular to it, $\vec{C}_R$:

$$-\vec{V} = -[V_X, V_Y, V_Z] = \vec{C}_A + \vec{C}_R \quad (6)$$

Using second term of Eq. (3) and eliminating the magnitude of the cylinder axis vector gives components of the vector perpendicular to the cylinder axis:

$$\vec{C}_R = [-V_X, (V_Z \cos \omega - V_Y \sin \omega) \sin \omega, -(V_Z \cos \omega - V_Y \sin \omega) \cos \omega] \quad (7)$$

The dot product of $\vec{C}_R$ and a vector perpendicular to the cylinder axis in the y-z plane $\vec{C}_U = [0, -\sin\omega, \cos\omega]$ gives, with the help of first term of Eq. (3), the maximum impingement angle:

$$\cos \alpha = \operatorname{sign} \left[ \left( \frac{\sin \gamma}{\cos \gamma \sin \omega - \cot \varepsilon \cos \omega} \right)^2 + 1 \right]^{-0.5} \quad (8)$$
$$\operatorname{sign} = \begin{cases} -1 & \text{when } \tan \varepsilon \cos \gamma \tan \omega > 1 \\ 1 & \text{when } \tan \varepsilon \cos \gamma \tan \omega < 1 \end{cases}$$

From this location, the mass impingement amount decreases, following a cosine function of an angle in a plane perpendicular to cylinder axis, to vanish at the ±90° locations.

<figure><img src="image_350.png" /><figcaption>**Figure 2.** Impingement coefficient, λ, (left column) and maximum impingement angle, α, (right column) as a function of yaw, γ, and inclination, ω, angle.</figcaption></figure>

(a) vertical precipitation, ε = 0°.

(b) horizontal speed equal to terminal speed, ε = 45°.

(c) horizontal speed twice terminal speed, ε = 63.4°.
Equations (5) and (8) have been used to examine ice accretion mass, shape and its location as a function of the precipitation angle, ε, yaw angle, γ, and inclination angle, ω. Three cases have been considered: vertical precipitation, ε = 0°, Fig. 2a; horizontal velocity equals terminal velocity, ε = 45°, Fig. 2b; and horizontal velocity twice the terminal velocity, ε = 63.4°, Fig. 2c. Since the typical fall speed for drizzle is around 3.5 m/s and for rain is 6.5 m/s, in the extreme case, the wind speed would be approximately 50 km/h. For vertical precipitation, the precipitation yaw angle is meaningless and λ = cos ω, Fig. 2a. In addition, the maximum impingement angle is always 0°. When the cylinder is horizontal, ω = 0°, the precipitation impinges perpendicular to the cylinder axis, and consequently, the impingement coefficient is unity. An increasing cylinder inclination angle leads to a decrease in the impingement mass per unit length of the cylinder. When the cylinder is vertical, ω = 90°, drops do not impinge on the cylindrical surface (impingement on the cylinder ends is not included in this analysis), resulting in an impingement coefficient of zero.

Horizontal airflow enhances the impinging precipitation amount, P, for a given value of raingauge-measured precipitation, Pv, as given by Eq. (2). The case when the precipitation horizontal speed (which is assumed to equal the wind speed) equals the vertical terminal drop speed has been examined in Figure 2b. In addition, the ice shapes have been plotted for four inclination angles, Fig. 3, for a vertical precipitation amount equal to the cylinder radius, $P_v = R$. In order to obtain the ice shape it has been assumed that the initial growth rate does not change with time and that the ice grows radially away from the cylinder surface. For given inclination angle, the ice shape is displayed for varying yaw angle with five enhanced cross-sections at 0°, 45°, 90°, 135° and 180° that are depicted with a different circumference colour. The ice surface colour, with corresponding contour lines, depicts the local ice thickness. It should be noted that the accretion shape and thickness for given conditions is uniform along the length of the cylinder and the apparent longitudinal variations on Fig. 3 are the result of yaw angle variations.

<figure><img src="image_351.png" /><figcaption>**Figure 3.** Schematic illustration of the influence of yaw angle on ice accretion shape for four inclination angles, assuming ε = 45°. H/R is the relative ice accretion thickness.</figcaption></figure>

(a) ω = 0°, (b) ω = 30°, (c) ω = 60°, (d) ω = 90°
<figure><img src="image_352.png" /><figcaption>Figure 4. The impingement coefficient, λ, as a function of yaw angle, γ, and precipitation angle, ε. (a) horizontal cylinder, ω = 0°, (b) vertical cylinder, ω = 90°.</figcaption></figure>

For a precipitation angle of 45°, Fig. 2b, and for a horizontal cylinder, Fig. 3a, the solution is symmetrical about a 90° yaw angle. When the yaw angle is 0° or 180°, the impingement coefficient is unity since the impinging mass flux is enhanced by $\sqrt{2}$ but diminished by $\sqrt{2}$ because of the approach angle. For these two limiting cases, maximum impingement occurs at the uppermost point of the cylinder, Fig. 3a. When the yaw angle is 90°, the impingement coefficient is $\sqrt{2}$, because the approach angle is 90°, but since $\epsilon = 45°$, the vertical precipitation is enhanced by a factor of $\sqrt{2}$. For this case, the maximum impingement occurs 45° from the vertical, as can be seen in the right-hand image of Fig. 2b. The ice shapes and thickness for all yaw angles are depicted in Fig. 3a.

For a cylinder inclination angle of 30°, Fig. 3b, changes of yaw angle have a more pronounced effect on accretion mass and shape than in the horizontal cylinder case, Fig. 3a. When the yaw angle is 0°, the impinging water mass is small, since the incoming drop direction lies only 15° from the cylinder axis. Increasing the yaw angle leads to an ice mass increase; it reaches a maximum for a yaw angle of approximately 125°, Fig. 2b, when the maximum impingement angle is approximately 35°. A further increase of yaw angle diminishes the impinging mass; for a yaw angle of 180°, the ice accretion is symmetrical about a vertical plane. Associated ice shape changes are depicted in Fig. 3b. For a cylinder inclination angle of 60°, variations in the impinging mass are somewhat similar to the $\omega = 30°$ case, but the ice location is markedly different for yaw angles smaller than 90°. For low yaw angles, ice forms underneath the cylinder. Finally, for a vertical cylinder, a change of yaw angle leads only to a change in the ice location, but the ice mass is independent of yaw angle, Fig. 3d.

When the horizontal velocity is increased twofold, $\epsilon = 63.4°$, the impingement coefficient increases, Fig. 2c. The precipitation amount is enhanced by a factor of $\sqrt{5}$, see Eq. (2). When the cylinder is horizontal and the yaw angle is 90°, the impingement coefficient has a maximum value of $\sqrt{5}$. For a vertical cylinder, the impingement coefficient decreases to 2. For a yaw angle of 180°, the maximum precipitation occurs when $\omega = \epsilon$. Conversely, for a yaw angle of 0°, zero precipitation occurs when $\omega = 90°-\epsilon$. These last two statements are valid for any precipitation angle.

Horizontal and vertical cylinder configurations have also been analysed as a function of yaw and precipitation angle, Fig. 4. The horizontal case is symmetrical about a yaw angle of 90°, for which $\lambda = 1/\cos \epsilon$. A high value of horizontal velocity leads to a precipitation angle increase and rapidly increasing impingement coefficient. However, when the wind blows along the cylinder, corresponding to yaw angle of 0° or 180°, the impingement coefficient remains unity. For a vertical configuration, the impingement coefficient is not a function of yaw and $\lambda = \tan(\epsilon)$, Fig. 4b. Consequently, an increase of precipitation angle leads to an initially gradual but then a more dramatic increase in impinging precipitation mass.

The above analysis is valid for an outdoor environment where it can be assumed that vertical precipitation is the input parameter. However, in order to examine ice accretion in a horizontal wind tunnel configuration, the above analysis has been modified.

## 2.2 Freezing precipitation with only horizontal impingement

Since the experimental validation was performed in a horizontal icing wind tunnel, we will now consider the analytical model for horizontally incoming precipitation. We assume now that horizontal precipitation, $P_H$, is the input parameter and $\epsilon = 90°$. An analysis similar to that described in the previous section, which led to Eq. (5), now gives:

$$m = \lambda P_H 2R \rho_w \qquad (9)$$

where $\lambda = \sin \delta$ and $\cos \delta = \cos \gamma \cos \omega$

In addition, the maximum impingement angle function assumes the following simpler form:

$$\cos \alpha = \operatorname{sign} \left[ \left( \frac{\tan \gamma}{\sin \omega} \right)^2 + 1 \right]^{-0.5} \qquad (10)$$

where $\operatorname{sign} = \begin{cases} -1 & \text{when } 0° \le \gamma \le 90° \\ 1 & \text{when } 90° < \gamma \le 180° \end{cases}$

Graphics of the solutions to Eqs. (9) and (10) are given in Figures 5 and 6. When the yaw angle is 90°, the cylinder inclination angle does not influence the impingement coefficient. It is equal to its maximum value of unity,
since the impingement direction is always perpendicular to the cylinder axis, Fig. 5. The maximum impingement angle also remains constant at 90°. The resulting ice shapes can be seen in Fig. 6, for four different cylinder inclination angles, for the case when $P_H = R$. It may be noted that the ice shape is independent of inclination angle for a 90° yaw angle. Figure 5 shows that the solution is symmetrical about a yaw angle of 90°.

However, when the yaw angle is less than 90°, impingement occurs on the lower part of the cylinder. Consequently, the maximum impingement angle varies between 90° and 180°. For a yaw angle larger than 90°, impingement occurs on the upper surface. For a horizontal cylinder with a precipitation yaw angle of 0° or 180°, drops do not impinge on the cylinder and consequently, the impingement coefficient is zero.

<figure><img src="image_353.png" /><figcaption>**Figure 5.** The impingement coefficient, $\lambda$, (left) and maximum impingement angle, $\alpha$, (right), as a function of yaw, $\gamma$, and inclination, $\omega$, angle for horizontal precipitation.</figcaption></figure>

<figure><img src="image_354.png" /><figcaption>**Figure 6.** Schematic illustration of the influence of yaw angle on ice accretion shape, for four inclination angles, assuming horizontal precipitation impingement. H/R is the relative ice accretion thickness.</figcaption></figure>

(a) $\omega = 0°$, (b) $\omega = 30°$, (c) $\omega = 60°$, (d) $\omega = 90°$.
Figure 6 depicts the ice shape profiles/cross-sections for four cylinder inclination angles and a continuously varying yaw angle.

## 3 Ice accretion laboratory experiments

The ice accretion experiments were performed in the NRC Altitude Icing Wind Tunnel (AIWT), which is a specialized, closed-loop wind tunnel designed to simulate an in-flight atmospheric icing environment [5]. In the tunnel test section of 57 by 57 cm, a cylinder of diameter 5.08 cm was placed vertically and at an inclination angle of 30° with varying yaw angle with respect to practically horizontally incoming precipitation. All horizontal icing wind tunnels, including AIWT, have difficulties to deliver large drops that are characteristic of freezing rain to the test section. Consequently, smaller drops of volume median diameter of 120 µm were produced and delivered to the cylindrical surface in an airflow of 20 m/s. The air relative humidity was 80%, the airstream liquid water content was 0.26 g/m³, and the tests were run at atmospheric pressure. The ice was grown for 30 min and the ice accretion shapes were observed and documented. The aim of these tests was to collect a sample of data that could be used to validate our analytical and numerical models. A more extensive experimental campaign with varying inclination angle, yaw angle and larger air temperature range is planned in the future.

In this section, the focus is on the experimental rate of ice growth at the stagnation line and its comparison with a simple thermodynamic model outlined below. In the subsequent validation section, we will describe other characteristics of the experimental ice shapes. Since the analytical model described in the previous section assumes instant drop freezing upon impingement, it is applicable only for very cold conditions. However, when water flows on the surface during warmer conditions, a model that includes thermodynamics should be applied. Figure 7 shows the ice growth rate at the stagnation line for all of the experimental cases. For each experimental test, the ice thickness was measured at three locations spaced 5 cm apart along the cylinder axis. The experimental results for a vertical cylinder are depicted as points in Fig. 7a. Decreasing air temperature leads to a gradual increase of the ice growth rate, which is consistent with the observation that, at lower air temperatures, less water is shed from the cylinder. When the air temperature was -15°C, no water shedding was observed and the ice thickness reached a maximum value. The inclined cylinder at 30° was exposed to incoming precipitation at various yaw angles and tests were run for air temperatures of -2°C and -5°C, Fig. 7b. The experimental results show limited variability as a function of yaw angle. However, the strong influence of air temperature could be seen, lower temperatures being associated with greater ice growth rates.

These experimental results were compared with a simple thermodynamic ice accretion model outlined below. It has been assumed that freezing is determined by convective and evaporative heat fluxes. Hence the rate of

ice growth, $\dot{H}$ (m/s), is determined by the following heat balance equation:

$$ \dot{H} \rho_I L_F = h(T_S - T_A) + h \frac{\varepsilon L_V}{\rho c_P} (e_S - RH e_A) \quad (11) $$

where: $\rho_I$ is the glaze ice density, $L_F$ is the specific latent heat of freezing, $L_V$ is the specific latent heat of vaporization, $\varepsilon$ is the ratio of the molecular weights of water vapour and dry air, $p$ is the static pressure in the freestream, $c_P$ is the specific heat capacity of dry air at constant pressure, $T_S$ and $T_A$ are respectively the ice accretion surface temperature and the air temperature, $e_S$ and $e_A$ are the saturation vapour pressures at $T_S$ and $T_A$, and RH is the air relative humidity. The convective heat transfer coefficient, $h$, is computed using the Frossling number, $Fr$, defined as the ratio of the Nusselt number to the square root of the Reynolds number, $Fr = Nu/\sqrt{Re}$. Based on the experimental results [6], we assume that for comparatively low Reynolds number and a realistic range of surface roughness, the Frossling number is equal to unity. The sloped line in Fig. 7a and the horizontal lines in Fig. 7b are graphical representations of Eq. (11) assuming a glaze density of 900 kg/m³ and a surface temperature of 0°C. The agreement between the experimental points and this simple thermodynamic model is satisfactory.

<figure><img src="image_355.png" /><figcaption>Figure 7. Experimental ice growth rate (points) and simple thermodynamic model prediction (lines) for the stagnation point.</figcaption></figure>

(a) vertical cylinder, $\omega = 90°$

(b) inclined cylinder, $\omega = 30°$

Equation (11) gives the growth rate of ice thickness when there is sufficient water available to freeze. The maximum possible ice growth rate, $\dot{H}_M$, at the stagnation point is determined by the impinging liquid flux and given by:

$$ \dot{H}_M \rho_{IR} = \beta_M V LWC \quad (12) $$
where: $\rho_{IR}$ is the rime ice density at the stagnation point, $\beta_M$ is the collision efficiency at the stagnation point, and LWC is the liquid water content of the airstream. We have estimated a rime density of 750 kg/m³ [7] and calculated $\beta_M$ = 0.92 based on drop trajectory calculations for the experimental volume median drop diameter, air speed and cylinder radius. Equation (12) is depicted as the horizontal line on Fig. 7a. It, too, gives a reasonable agreement with experiments in the rime ice growth range.

We have already examined two models: the analytical model capable of ice mass and shape prediction when drops freeze on impact, and the thermodynamic model that predicts ice thickness at the stagnation line for rime and glaze conditions. Now we will discuss the NRC morphogenetic ice accretion model, which has a prediction domain that includes both mentioned models. In addition, the morphogenetic model can analyse cases with complex water flow and freezing on the surface, and it also predicts three-dimensional details of the ice accretion shape, mass and structure.

## 4 The morphogenetic model and its validation for horizontal precipitation

<figure><img src="image_356.png" /><figcaption>**Figure 8.** Comparison between morphogenetic model prediction and experimental ice accretion on a vertical cylinder. The first column shows three cross-sections, spaced 5 cm apart, produced by the morphogenetic model and depicted by the blueish lines. Three experimental ice cross-sections with the same spacing are represented by greenish lines. The red dot represents the stagnation line. The black line in (c) is the analytical model prediction. Side views of the numerical prediction and experimental ice are shown in the second and third columns, respectively.
(a) air temperature -2°C, (b) air temperature -5°C, (c) air temperature -15°C.</figcaption></figure>
The morphogenetic model was originally developed for aerospace in-flight icing applications (2D and 3D versions are described in [8] and [9], respectively). It has recently been modified to predict ice accretion shapes due to freezing rain precipitation. The main characteristics of our novel approach are described below. The morphogenetic model is a discrete element, random walk model that emulates the motion and freezing of individual fluid elements arriving at the accretion surface. Individual model fluid elements may be imagined to consist of an ensemble of drops, all of which undergo identical histories. However, each fluid element has a different history on the surface because of the stochastic nature of the model. A 3D rectangular lattice defines the accretion domain. Each fluid element begins a stochastic motion downstream along the surface from its initial impact location, which is also determined stochastically. At each step in the process, a random number is generated, and, according to its value, the element either freezes or moves along the surface. The model is sequential, so that as soon as a particular fluid element freezes, the behaviour of the next fluid element is considered.

Most of the results presented in this paper were obtained for a fluid element size of 0.4 × 0.4 × 0.4 mm and a morphogenetic model domain of 200 × 600 × 600 elements, giving 8 × 24 × 24 cm in the horizontal x and y, and vertical directions, respectively. The fluid element size was chosen as a compromise between model resolution and computational effort. An extensive sensitivity analysis to fluid element size has not yet been performed for the three-dimensional version of the model. However, a sensitivity analysis for the two-dimensional version of the model [8] suggests that the principal, large-scale features of the ice accretion are not sensitive to a reduction in fluid element size. However, with smaller elements, additional, fine-scale details of the ice structure are revealed.

In this section, the NRC morphogenetic model is validated for ice-on-cable applications, using the experimental results described above for horizontal airflow and drop impingement. In this version of the morphogenetic model, we assume that the impinging fluid elements move along straight lines. The precipitation flux is given by the right hand side of Eq. (12). Using our experimental values, the flux is 4.78 g/m²-s. The duration of the simulation for all examined cases is 30 min. We also assume a Frossling number of unity, which determines the convective heat transfer coefficient.

Results of the numerical prediction and laboratory tests of ice formation on a vertical cylinder for three air temperatures of -2°C, -5°C and -15°C are depicted in Fig. 8. The morphogenetic model shows that, when the air temperature is -2°C, the heat transfer is sufficient to freeze approximately 35% of the impinging liquid at the stagnation point. The unfrozen water flows downstream and freezes gradually. Because the aerodynamic drag acting on the surface liquid greatly exceeds gravity, the water flows nearly horizontally on the surface. Once the drops reach locations beyond 90° from the stagnation point, they are gradually shed. For the -2°C case, 72% of the impinging water mass freezes on the cylinder and the

remaining water is shed. Our comparison between the numerical ice accretion cross-sections and the experimental cross-sections shows good agreement. However, the experimental results exhibit higher variability than numerical prediction. The comparison between the side views of numerical and experimental ice accretions is also good. When the air temperature is -5°C, freezing is enhanced, leading to greater ice thickness on the upwind side of the cylinder. In the model, approximately 70% of the impinging water at the stagnation line freezes instantly, and only 0.2% of the impinging water is shed from the cylinder. Once again, the agreement between numerical and experimental ice cross-sections and side views is good. For an air temperature of -15°C, all impinging drops freeze instantly in the model, and there is no flow of liquid on the surface. In this case, the cross-section comparison is confirmed by the analytical model, where all impinging water freezes instantly. The agreement between the morphogenetic, analytical and experimental results is good.

Figure 9 shows a comparison between the model and experiments for ice forming on a 30° inclined cylinder, at a yaw angle of 0°. Since the collection efficiency at the stagnation point for this case is 0.5, see Fig. 5, the total impinging liquid mass is half that for the vertical case. When the air temperature is -2°C, 70% of the impinging water at the stagnation line freezes instantly. Drops that do not freeze immediately upon impingement freeze downstream and there is no water shedding. In addition, the rapid freezing and the low approach angle of 30° lead to the formation of a rough ice surface. The numerical and experimental ice cross-sections exhibit diverse ice surface irregularities, but the overall agreement between them is good. These ice shape anomalies can also be seen in the three-dimensional numerical and laboratory ice accretion views, the second and third columns in Fig. 9. Since the model ice features are best seen in profile, we show profile images, despite the different angle of view for the experimental ice accretions. The model accretions show the three-dimensionality of the ice accretion process, with ice rivulets aligned with the external airflow. They also show the variability of the ice accretion shape on different geometrical scales. When the air temperature is -5°C, all impinging drops freeze almost instantly upon impact. However, due to the shadowing effect of the ice irregularities, the resulting ice shapes are rather complex. In the aircraft icing literature, this type of accretion is sometimes called “lobster tail” [9]. It is encouraging to note that the predicted and observed 3D features are similar, inasmuch as the cross-sections are rather “jagged” in both the model and experiment. Both the morphogenetic model and experimental ice structures are full of cavities, which increase the ice structure thickness. However, the analytical model, depicted by the black line in the first (b) column, predicts a smaller ice thickness, since it assumes solid ice (high density) without voids.
<figure><img src="image_357.png" /><figcaption>**Figure 9.** Comparison between morphogenetic model prediction and experimental ice accretion on an inclined cylinder (30° inclination angle) and with a yaw angle of 0°. The first column is described in the caption of Figure 8. The second column depicts side views of a 10 cm segment of the morphogenetic model simulation, and the third column shows the corresponding experimental ice accretions, seen from below.</figcaption></figure>

(a) air temperature -2°C, (b) air temperature -5°C.

## 5 Conclusions

We have developed a morphogenetic numerical model that can simulate ice accretions on arbitrarily oriented cylinders exposed to freezing rain precipitation. This predictive model has been validated experimentally for a limited number of conditions in an aeronautical icing wind tunnel under horizontal droplet impingement. We have also developed a zero-order analytical model that identifies the location and size of freezing rain ice formation on a cylinder in a fully 3D geometrical space. This comparatively simple model assumes droplet freezing on impact, allowing swift identification of limiting ice shapes in a cold environment. The agreement between predictions of the morphogenetic and the analytical model for considered limiting conditions is satisfactory.

The validated morphogenetic model can be used to predict complex three-dimensional ice accretion shapes on overhead transmission lines, bridge cables, telecommunication masts or any other cylindrical objects exposed to freezing rain. This morphogenetic numerical model is not subject to the limitations of aircraft icing wind tunnels in trying to simulate freezing rain conditions. Although the model presented here was developed and validated for a cylindrical geometry, it could be easily extended to more complex geometries of engineering interest.

## Acknowledgement

The author wishes to thank Mr. G. Chevrette and Dr. D. Orchard for their help with experimental part of the project. Many thanks go to Prof. E. P. Lozowski for his help in improving the manuscript through his clarity of thinking and editorial skills.

## References

1. K. Jones, Ice Storms in the St. Lawrence Valley Region, ERDC/CRREL TR-03-1 (2003)

2. P. Fu, M. Farzaneh, G. Bouchard, Cold Reg. Sci. Technol. **46**, 132-146 (2006)

3. C. Demartino, H.H. Koss, C.T. Georgakis, F. Ricciardelli, J. Wind Eng. Ind. Aerodyn. **138**, 98 (2015)

4. L. E. Kollar, M. Farzaneh, Int J. Heat Mass Transfer, **53**, 849-861 (2010)

5. M. Oleskiw, F. H. Hyde, P.J. Penna. AIAA-2001-094 (2001)

6. E. Achenbach, Int. J. Heat Mass Transfer, **20**, 359-369 (1977)

7. W.C. Macklin, Q. J. R. Meteorol. Soc, **88**, 30-50 (1962)

8. K. Szilder, E.P. Lozowski. J. Aircraft, **41**, 834 (2004)

9. K. Szilder, E.P. Lozowski, SAE 2015 Int. Conf. Icing Aircraft, Engines, and Structures, 2015-01-2162 (2015)
# Parametric approach for assessing risks due to falling ice and snow

Jan Dale¹,ᵃ, Scott Gamble¹, Albert Brooks¹ and Jill Bond¹

¹Rowan Williams Davies & Irwin Inc., 650 Woodlawn Road West, Guelph, Ontario, Canada N1K 1B8

**Abstract.** Falling ice and snow from structures presents a significant challenge for building design in cold climate regions. The variables involved in the accumulation and subsequent shedding of snow and ice from a structure include building and feature geometry, thermal and surface characteristics, melt-water drainage paths, ambient air temperature fluctuations, solar radiation and wind forces. The transient and combined effect of these variables does not facilitate a discrete probabilistic prediction of the associated risks. As a result, there are no known standards or codes available to designers to provide guidance on addressing these concerns in design. This paper outlines a parametric approach for assessing the relative risks associated with falling ice and snow hazards for a range of façade details and climatic variables. The approach has been developed through a combination of scale model testing within a boundary layer wind tunnel, full scale testing in a climatic wind tunnel and numerical snow and ice accumulation models. The aim is to provide designers with a reliable approach to understand the risks associated with a particular design and the level of on-going maintenance efforts that will be required over the life of the structure to manage these risks and maintain a safe environment.

## 1 Introduction

Snow and ice accumulations which are shed from buildings and structures pose a significant hazard to pedestrians and property due to the associated large impact forces. However, snow and ice falling from structures into pedestrian areas remains a commonplace occurrence in cold climates. The key factors regarding whether the accumulations that are shed from a structure pose a serious risk are related to how massive the falling accumulations are and how often the events happen.

Ice and snow collection and the subsequent shedding can be influenced by, but is not limited to, air temperature, solar radiation and heat loss, insulation properties of accumulated snow, the cyclical freeze and thaw process created by daily temperature cycles and the orientation and surface friction properties of the materials. The transient and combined effects of these variables do not easily facilitate discrete quantification of all potential falling ice and snow scenarios. As a result, to the knowledge of the authors, there are no known criteria or guidelines for designers to assess acceptable performance considering falling ice and snow hazards. Rather, it is necessary to balance effective design detailing with on-going maintenance efforts over the life of the structure to effectively mitigate falling ice and snow concerns.

An approach for assessing the performance of a structure has been developed whereby building features that are identified through experience as sensitive to snow and ice accretion and shedding are evaluated. The approach relates the severity of a potential shedding scenario relative to that of a reference façade detail. The methodology involves an assessment of the magnitude of potential falling accumulations for a range of mean recurrence intervals. This is then combined with an estimate of the probability that a falling accumulation will impact either pedestrians or property before coming to rest to estimate the probability of a problematic shedding event. An estimate of the potential magnitude of this event relative to that of the reference façade detail is

used to determine a severity index associated with each façade detail. An example assessment result is shown in Figure 1. The following outlines the methodology proposed to provide design teams with adequate information to assess the performance of a structure.

<figure><img src="image_358.png" /><figcaption>**Figure 1.** Example assessment result.</figcaption></figure>

A combination of scale model testing within a boundary layer with tunnel, full scale testing in a climatic wind tunnel and numerical snow and ice accumulation models are used to provide insight into the relative risks.

## 2 Assessment of shedding severity

The random nature of many of the variables involved in the accumulation and subsequent shedding of snow from a building façade makes a deterministic solution to the mass and frequency of shedding almost impossible to achieve. Further, there needs to be the presence of pedestrians or property in the location of impact in order for an accumulation to pose a significant risk when shed from a structure. As such, the relative risks associated with falling accumulations are assessed through the introduction of a severity index. A probabilistic approach has been used to estimate the severity index at a level of

<figure><img src="image_359.png" /></figure>

<footnote>ᵃ Corresponding author: Jan.Dale@rwdi.com</footnote>
risk consistent with the design snow loading provisions for structural design.

## 2.1 Relative volume of accumulations shed from a façade surface

### 2.1.1 Accumulation magnitude

The potential risks associated with a falling snow event are directly related to the mass of the accumulation being shed. A relatively small accumulation falling from a building surface poses significantly less risk than a large dense accumulation. As a result, an estimate of the relative volume of potential falling accumulations has been used as the basis for the estimation of the severity of a shedding event.

As outlined in further detail in Sections 3 and 4, the effective depth of an accumulation is estimated based on the geometry and location of a particular façade detail on the building. Estimates of the effective depth of accumulation are derived based on the anticipated snow load combined with an average snow density of 3 kN/m³. This density is representative of an aged or wind driven snow accumulation. The resulting effective cross section is then used to estimate the potential volume of an accumulation that may be shed from the façade.

The key uncertainty associated with the volume of an accumulation that may be shed from a structure it the length of piece that falls. Given the range of conditions that may lead to the shedding scenario, accumulations may be shed in relatively small pieces, or as a consistent mass equal to the length of the façade detail. For the purposes of comparison, it has been assumed that the potential length of accumulation that is shed from a façade detail is proportional to the effective cross sectional area moment of inertia of the accumulation.

## 2.2 Reference severity index

The shedding performance model relates the severity of a potential shedding scenario relative to that of a reference façade detail. For the purposes of this research, a 10cm wide horizontal ledge surface has been used as the reference as accumulations shedding from ledges of this size have historically been considered acceptable by building owners. While surfaces of this nature will collect and potentially shed snow and ice accumulations, the size and consistency of these are not known to have been problematic. Larger accumulations that fall from larger ledges have increased tendency to be considered inappropriate, and mitigation investigations are commonly sought as a result.

As the potential risk associated with a falling accumulation is dependent on the presence of pedestrians or property in the vicinity of the impact area, probability of hazardous conditions is the combined probability of having an accumulation and an estimate of the probability of impacting something of significance. As a reasonable estimate of an average condition for reference, it has been assumed that the probability of occupancy in the impact area is 0.25. Given this, the resulting reference return

period for snow accumulation to achieve an effective 50-year return period falling snow hazard event is 12.5 years. This is determined as follows:

$$R_A = R_{FH} * 1 / P_P \quad (1)$$

Where:
* $R_A$ is the mean recurrence interval of the snow accumulation;

* $R_{FH}$ is the mean recurrence interval of a falling snow hazard; and,

* $P_P$ is the probability of persons or property in the impact area during a falling event.

The accumulation cross section for a 10cm wide ledge representative of a 12.5-year return period was determined based on historical meteorological information. The product of this cross section and associated moment of inertia was used as the reference for a severity index of unity. The ratio of the product of the effective cross sectional area and moment of inertia associated with other return periods and façade details to this reference magnitude is then used as the basis for comparison. This ratio is then scaled to account for the slope of the surface. A linear scale factor ranging from 1.0, for horizontal surfaces, to 10 for vertical surfaces has been introduced to account for the anticipated increase in shedding frequency for steeply sloped surfaces. The resultant is the severity index associated with a particular façade element.

Figure 2 shows a plot of the relative change in severity index for a 10cm ledge as a function of the return period of snow accumulation. As outlined, the effective 50-year return period falling snow event for the 10cm ledge, given a 12.5 year snow accumulation and 25% chance of having something of significance at the point of impact, has a severity index of 1.0. It should be noted that, the severity index of the same ledge detail would be increased if it were assumed that the ledge was located above an area where the probability of impacting something of significance was higher. For example, if it is assumed that an accumulation falling from the ledge would always impact something of significance, the severity index would be increased to approximately 3.5.

<figure><img src="image_360.png" /><figcaption>Figure 2. Severity index plot for reference surface.</figcaption></figure>

The severity index may be associated with approximate operational impacts for the purpose of comparison of the performance of various façade details. In general, severity indices of 1.0 or less may be
considered to be indicative of small falling accumulations that are likely to be considered a nuisance, but not likely to cause harm. Façade details with severity indices ranging from 1.0 to 5.0 are considered to be small enough to be handled through effective operational and maintenance procedures such as clearing operations and reduced access to the impact areas. Larger chunks that could be expected to cause some safety issues or property damage are expected from façade details with severity indices ranging from 5.0 through 10. Façade details with severity indices above 10 are expected to experience large chunks that may be considered a significant risk to life safety or extensive property damage. Clearly, the better performing building would shed smaller chunks and/or shed less frequently.

# 3 Model simulation techniques

## 3.1 Climatic wind tunnel simulations

Climatic wind tunnel simulations are conducted on building façade mock-up sections to observe the accumulation and shedding characteristics of snow on various details. A range of variables including temperature, snow and ice consistency, effects of wind and façade geometry can be varied to gain an understanding of the overall performance of a particular design. Figure 3 shows an example of wet snow accumulations on a horizontal sunshade element during a climatic wind tunnel simulation.

<figure><img src="image_361.png" /><figcaption>**Figure 3.** Example of snow accumulation on a sunshade during a climatic chamber simulation.</figcaption></figure>

It should be noted that although there is a difference in simulated to natural snow, ice particles from the snow gun were assumed to attach to surfaces in a similar manner to the ice at the end of the fingers of a snow flake. The ambient temperature, snow gun air pressure and water discharge were varied to give a range of conditions to simulate the various types of frozen accumulations. Observations of snow accumulations on various full scale mock-up façade details were used throughout this research.

### 3.1.1 Effect of snow consistency

As it is not practical to simulate all of the combinations of frozen precipitation that may occur, a range of

conditions that bound the general range of snow characteristics are modelled. Simulations considering relatively dry snow and wet snow are typically conducted to observe any variation in accumulation and shedding characteristics.

As shown in Figure 4, dry snow simulations do not show significant accumulations. While accumulations are often observed, they tend to be limited in depth and not well adhered to the façade surface. As a result, these accumulations are expected to be either melted or shed from the surface in a manner that does not pose a significant concern.

<figure><img src="image_362.png" /><figcaption>**Figure 4.** Example of dry snow accumulation on a mullion in a climatic chamber simulation.</figcaption></figure>

In contrast, the observations in wet snow simulations, as shown in Figure 5, are significantly different. Windward façade details readily accrete snow. Further, the accumulations tend to be relatively dense. When shed from the façade sections, the accumulations that form as a result of wet snow applications tend to remain relatively well intact. As a result, accumulations under these conditions are likely to be relatively more severe than those for dry snow events.

<figure><img src="image_363.png" /><figcaption>**Figure 5.** Example of wet snow accumulation (after partial shedding) on a mullion in a climatic chamber simulation.</figcaption></figure>

It should be noted that, the simulations shown in Figures 4 and 5 were conducted using similar snow application rates and durations. The only significant change to the simulation was the consistency of the applied snow.
### 3.1.2 Capture efficiency

While it is common knowledge that snow will accumulate on horizontal surfaces such as ledges and parapet caps, accumulations on steeply sloped or even vertical surfaces are often overlooked. This is likely due to the fact that most building codes and standards assume snow loading to be negligible for steeply sloped surfaces. However, the aerodynamics of many façade details, particularly those that are separated from the heated interior of a building by a thermal break, can accrete significant quantities of snow. The key difference is the duration that the accumulation will be present on the structure. For most cold climate regions, the structural loading due to snow is generally the combination of accumulations that occur over an entire winter. In the case of vertical façade details, the accumulations typically would only remain on a structure for a period of days or potentially weeks. As a result, the impact on structural loading may be considered to be negligible. However, falling accumulations can be of sufficient size to cause significant damage or injury and are shed from the structure potentially many times in a given winter.

The observations from climatic chamber simulations have been used as the basis for assessing the potential magnitude of accumulations on steeply sloped and vertical façade details. Figure 6 shows accumulations on a flat vertical façade during a climatic chamber simulation. It can be seen that, while accumulations do form on the façade surface, they are not densely packed and uniformly distributed. This is due to the redirection of wind flows along the surface of the façade. As a result, the shedding of these accumulations is generally not problematic.

<figure><img src="image_364.png" /><figcaption>**Figure 6.** Example of wet snow accumulation on a flat vertical façade during a climate chamber simulation.</figcaption></figure>

Alternatively, as seen previously in Figure 5, vertical projections from a façade surface do have the potential for significant, dense accumulations. These accumulations form on projections from the façade or other relatively narrow surfaces such as signs or other decorative façade elements as winds are not redirected to the same extent as with the main façade. As a result, there is the need to determine the potential magnitude of accumulation on these features to facilitate an estimate of the associated risks. The simulation parameters and observations from multiple climatic chamber simulations

of various façade projections were used to estimate the effective capture efficiency for this type of detail. This capture efficiency is analogous to the basic reduction factors used by most building codes and standards for deriving roof snow loading distributions from ground snow load estimates. The capture efficiency was defined as the ratio of the average depth of accumulation on the façade detail to the effective depth of snow falling on an unobstructed horizontal surface. Based on the observations, a capture efficiency of 0.8 was found to be a reasonable assumption.

## 3.2 Finite area element (FAE) modelling

The Finite Area Element (FAE) methodology [1] was employed to investigate the contribution to snow accumulations and potential shedding from independent variables. The FAE method includes a detailed computer simulation of the hour-by-hour deposition, drifting and melting of snow on building roofs and façade details. The absorption of rain and melt water into the snowpack as well as the thermal characteristics of the façade are also considered throughout each of the simulations. Entire winters are simulated, including the cumulative effects of successive storms and melting periods. Available meteorological records such as daily snowfall, hourly snowfall intensity, hourly wind speeds and directions, rainfall, temperature, solar radiation and thermal performance of the façade may be altered to observe their overall effect on snow accumulation magnitude.

### 3.2.1 Effects of wind directionality

Some façade details can be sensitive to wind exposure when considering the potential accumulations that may form. In many instances, such as with vertical mullions on a particular façade of a building, significant accumulations will only form for precipitation events associated with winds from a particular sector of directions. To account for this variability, the meteorological data used in the FAE assessment of snow accumulation were filtered to include only storm events from the sector of interest. In this way, the reduced frequency of a problematic shedding event was assessed.

The plot in Figure 7 shows the snowfall for the winter of 1995-96 based on meteorological data collected at LaGuardia International Airport in New York, NY. This figure is considering all snowfall events regardless of wind directionality. In comparison, the plot in Figure 8 shows the snow fall events for the same winter considering only events that occur when winds would impact the west façade of a subject building. It is evident in comparing these two plots that considering all snowfall events equally will over-estimate the total amount of snow when assessing the potential for wind driven accumulations.

The FAE methodology was used with the various filtered meteorological data sets to estimate the magnitude of potential accumulations. Figure 9 outlines the variation in snow load on a building surface if all winds are considered, while Figure 10 presents the
loading magnitude if only westerly winds are considered. It can be seen in comparing these two figures that there is a significant difference in potential loading for building features prone to accumulations from all wind directions compared to features limited to wind-driven accumulations from particular directional sectors.

<figure><img src="image_365.png" /><figcaption>**Figure 7.** Snowfall at LaGuardia International Airport during the winter of 1995-96 considering all winds.</figcaption></figure>

<figure><img src="image_366.png" /><figcaption>**Figure 8.** Snowfall at LaGuardia International Airport during the winter of 1995-96 considering only westerly winds.</figcaption></figure>

<figure><img src="image_367.png" /><figcaption>**Figure 9.** Snow load variation at LaGuardia International Airport during the winter of 1995-96 considering all winds.</figcaption></figure>

<figure><img src="image_368.png" /><figcaption>**Figure 10.** Snow load variation at LaGuardia International Airport during the winter of 1995-96 considering only westerly winds.</figcaption></figure>

### 3.2.2 Snow consistency consideration

Similar to the effect of wind directionality, some accumulations are sensitive to the consistency of the snow. This is particularly true for wind driven accumulations. As discussed based on observations in the climatic wind tunnel simulations, wind driven dry snow accumulations do not tend to build to problematic magnitudes and are generally shed from the façade in small amounts that would not pose a significant risk to pedestrians or property. However, wet snow accumulations tend to be more cohesive and can experience significant hardening through repetitive freeze thaw cycles before being shed from a building. As with the wind directionality assessments, the meteorological data were filtered to assess accumulations within a range of snow consistency. For the purposes of this research, snow occurring at temperatures less than -2.5° Celsius was considered to be dry, while that occurring at or above -2.5° Celsius was considered to be wet.

Similar to the consideration for wind directionality, the meteorological data were filtered to consider only wet snow precipitation events. Figure 11 shows an example of the wet snow precipitation events in New York. The corresponding snow accumulation magnitudes are presented in Figure 12.

<figure><img src="image_369.png" /><figcaption>**Figure 11.** Wet snow precipitation events at LaGuardia International Airport during the winter of 1995-96.</figcaption></figure>

<figure><img src="image_370.png" /><figcaption>**Figure 12.** Snow load variation at LaGuardia International Airport during the winter of 1995-96 considering only wet snow events.</figcaption></figure>
### 3.2.3 Effects of snow sliding

In addition to the ability to parametrically assess the impact of the climate variables associated with snow loading, the FAE methodology was used to assess the potential loading magnitude for sloped surfaces that will regularly shed accumulations. For simplification, the specific coefficient of friction and the resulting slipperiness of the roof were not included in these analyses. Rather, weather conditions such as rain and melting that encourage sliding from generally smooth surfaces were focussed on. This provided a starting point for the resulting snow load reduction due to sliding since the basic assumption was that low friction between the building surface and snow pack would be produced by liquid water at this interface, not the material itself. Sliding was assumed to occur whenever there was sufficient melt water produced in the snowpack or rain water draining through the snow pack. Note that this could happen multiple times a winter in some areas.

# 4 Parametric assessment of shedding magnitude

## 4.1 Façade geometry

To facilitate the parametric assessment, building details have been categorized as being horizontal, sloped or vertical. Each of these classifications has a fundamental difference in the relative amount of snow or ice that will accumulate during a precipitation event. The present work outlines the key concerns with respect to falling ice and snow for these generalized classifications and how the parametric analysis has been developed to estimate the risks associated with each.

### 4.1.1 Horizontal surfaces

The FAE methodology was employed directly to determine the magnitude of accumulations on horizontal building surfaces accounting for wind scouring and potential drift accumulations.

### 4.1.2 Sloped façade surfaces

Further, the effects of snow sliding from a building surface have been estimated. The formation of melt water is monitored throughout the FAE analysis process. An accumulation was considered to be able to slide given the condition that melt water is at the interface between the building façade and the accumulated snow. A probability of slide was applied in the analysis to determine the effects of sliding on sloped building façade details. An appropriate probability of slide was determined such that the 50-year accumulation of snow equalled that of a horizontal building surface, reduced by an appropriate slope reduction factor as obtained from building codes and standards.

### 4.1.3 Vertical façade surfaces

For vertical accumulations, the probability of slide given the condition of melt water at the façade-accumulation interface was set to be unity. In addition, the rate of accumulation was adjusted based on observations of vertical façade elements in climatic wind tunnel simulations.

# 5 Example assessment

A generic building incorporating various façade details and located in New York has been used as the basis for an example use of the proposed methodology. The building incorporates 0.3m wide sunshade elements on all facades; 0.3m wide parapet caps around the perimeter of the roof; 0.3m deep vertical mullions on the north façade; and 0.3m wide sills sloped at 45° on the north façade. The following outlines the use of the proposed assessment methodology.

Each of the building details highlighted in this example are sensitive to a different set of loading assumptions. The parapet cap would be susceptible to all snow consistencies and exposed to all wind directions. Similarly, the sunshade elements will collect snow during all snowfall events, but significant accumulations would be limited to events that occur from wind directions approaching the building façade. The sloped sill accumulations are similar to the sunshade elements, but are further reduced due to the relative increase in the frequency of shedding throughout the winter months. It has been assumed that the vertical mullion details would only have accumulations during wet snow events that are associated with winds blowing into the façade. These accumulations were assumed to shed when there was water present at the interface of the mullion and the snow. The resulting severity indices of these elements based on these assumptions is presented as a function of snow accumulation return period in Figure 13. It should be noted that, although the accumulation magnitude associated with the sloped surfaces are expected to be less than the horizontal surfaces, the relative severity indices could exceed that of other façade features as a result of the increased frequency of shedding.

With the severity index relationships determined for each of the façade details, one can readily assess the relative risks associated with falling snow. For example, if there are minimal pedestrian areas along the north, it may be assumed that the probability of persons or property in the impact area during a falling event is relatively low. If it is assumed to be 10%, the corresponding snow accumulation corresponding to the 5-year mean recurrence interval would yield the 50-year falling snow hazard. As the parapet cap would potentially be adjacent to a major urban pedestrian thoroughfare, the probability of persons in the impact area could be considered to be much higher. If it is assumed to be 100%, the severity index associated with the 50-year snow accumulation would be assessed. The assessment would yield severity indices of 1.2, 2.2, 6.2 and 32 for the vertical mullion, horizontal sunshade, sloped sill and parapet cap, respectively.
<figure><img src="image_371.png" /><figcaption>**Figure 13.** Effect of geometry on the severity index of a façade detail.</figcaption></figure>

In addition to the relative performance of different façade details, the proposed methodology may be used to optimize the design of the façade features and focus mitigation expenses on the areas that pose the most significant risk. As outlined in Figure 14, the same façade detail mounted on different building faces can yield significantly different severity indices. The plot is showing the variation in severity indices for a 0.3m wide sunshade mounted on the north, east, south, and west facades. If the probability of persons or property in the impact zone is assumed to be 50%, then the resulting severity indices would be 0.01, 0.3, 7.5 and 15.6 for the south, west, east and north facades, respectively. Effectively, the sunshades mounted on the south and west facades would be expected to shed accumulations that could be considered a nuisance, but not posing a significant risk to safety. However, the sunshades mounted on the east and north facades would be expected to pose a concern for public safety. As a result, it may be appropriate to focus mitigation efforts on the north and east facades, or eliminate these features from those facades all together. Ultimately, this methodology provides the design and ownership teams with the information required to make an educated decision on the risks associated with a particular design.

<figure><img src="image_372.png" /><figcaption>**Figure 14.** Effect of wind directionality on the severity index of a façade detail.</figcaption></figure>

# 6 Conclusions

Even though snow accumulations that fall from structures present a significant risk to pedestrians and property, there are no known standards or codes available to assist designers in addressing these concerns in design. The wide range of variables involved in the accumulation and subsequent shedding of snow from a structure does not facilitate a discrete probabilistic prediction of the associated risks. The presented methodology outlines a parametric approach for assessing the relative risks associated with falling snow hazards for a range of façade details and climatic variables. The approach was developed through a combination of boundary layer wind tunnel and climatic wind tunnel simulations in addition to numerical snow accumulation models. The result is a reliable approach for designers to understand the risks associated with a particular facade detail. Variations in façade details and pedestrian usage of the areas around a particular building are accounted for through the proposed methodology. This will allow building designers and owners to focus mitigation and maintenance efforts on the building details that pose the most significant risk of a falling snow hazard.

# References

1. Gamble, S.L., Kochanski, W.W. and Irwin, P.A. *Finite Area Element Snow Loading Prediction – Applications and Advancements*. Journal of Wind Engineering and Industrial Aerodynamics, 41-44: 1537-1548 (1992)

2. ASCE 7-10, 2010. *Minimum Design Loads for Buildings and Other Structures*. American Society of Civil Engineers, Reston, Virginia.

3. National Research Council of Canada. *National Building Code. 2010. Vol.2 Section 4.1.6*. Ottawa, ON.
Session 6

# Building / Simulation
# Prediction of snow loads: past, present and future

Peter A. Irwin¹a

¹Rowan Williams Davies and Irwin Inc., 650 Woodlawn Road West, Guelph, Ontario, N1K 1B8, Canada.

**Abstract.** : The evolution of building code provisions is reviewed, with a main focus on Canada, and the development of physical and computational modeling methods for studying snow accumulation, melting and sliding are described. Physical models in water flumes and wind tunnels using various particles to represent snow can provide useful guidance on the locations and forms of snow drifts, particularly when snow loads are dominated by a single event. Computational Fluid Dynamics with particles can provide similar guidance. To go beyond a single event it is necessary to numerically track the multiple processes of snowfall, rainfall, melting, and sliding over extended periods. Therefore numerical methods are needed that make full use of the meteorological records and the building's geometry and thermal characteristics. Some of these methods are described in the paper. With the trends towards better energy efficiency and more varied building shapes, buildings are being designed that increasing are in need of more sophisticated methods for predicting snow and ice behavior.

## 1 Introduction

In snowy countries there are many impacts of snow loads on buildings and on the infrastructure that we rely on. An important one is the extra gravity load that snow causes on a roof. The load can be relatively uniform or it may be redistributed by the wind and by ponding of melt water. In drifts it may result in loads that are several times the snow load on the ground. Therefore snow loading has always been of concern to designers, and building codes have long specified minimum design snow loads. In Canada, up until 1960, the specified roof loads were equal to the 30 year ground load applied uniformly over the roof (Allen [1]; Shriever et al [2]), but they evolved subsequently so as to make allowance for reduced loads due scouring of snow off roofs by wind, and increased loads that occur in areas of aerodynamic shade such as in roof steps, roof valleys, the downwind sides of gable, arch and dome roofs and around rooftop projections such as mechanical units (Taylor [3]). They also make allowance for sliding of snow off sloped roofs, both in reducing the snow on an upper roof and increasing loads on a lower roof.

Another effect of snow or ice is that it can alter the shape of a structure seen by the wind. The electrical power industry have learned how the accumulation of ice on power transmission lines can lead to large amplitude wind-induced galloping oscillations and have developed measures to deal with this problem. Another example is snow, ice or freezing rain accumulating on and blocking up an otherwise very porous fence at the edge of a bridge deck. Wind tunnel tests show this increases aerodynamic drag and can amplify the problems of vortex excitation, flutter and galloping. The reliability of the structure needs to be assessed with such possibilities in mind.

There are other practical issues that must be addressed in building design and operations. The sliding of snow and ice off a building can cause a hazard to pedestrians below. Lawsuits have been launched in some cases due to injury or death that resulted. The measures that have been employed to avoid the adverse effects of sliding are varied, ranging from melting the snow, to holding it on the roof, to treating it as an operational problem and

simply fencing off hazardous areas until the threat has diminished (Williams et al [4], Isyumov [5]). The same type of problem has occurred on cable supported bridges where snow and ice falling off the cables has impacted vehicle or pedestrian traffic below.

In this paper the focus will be on snow loading and how methods of prediction have evolved and will probably continue to evolve in future.

## 2 Snow loading viewed as a chain

In wind engineering the problem of wind loading is often represented as a chain. In fact the chain has been formally named the Alan G. Davenport wind loading chain. A somewhat similar concept can be applied to snow loading and is depicted in Figure 1.

<figure><img src="image_373.png" /><figcaption>**Figure 1.** Snow loading viewed as a chain</figcaption></figure>

The chain has five main links. Starting at the left hand end, the first is the climate. Right from this starting point one realises that snow loading is more complex in many ways than wind. We see that whereas for wind we only need concern ourselves with the statistics of wind speed and direction, for snow we need to know about the statistics of snow fall, wind, rain, temperature and solar radiation. Therefore there are many more variables to keep track of compared with wind loading. The second link in the chain is the terrain surrounding the site. Like in wind engineering the ground roughness and topography affect the wind's mean velocity profile and

<figure><img src="image_374.png" /></figure>

<footnote>a Corresponding author: peter.irwin@rwdi.com</footnote>
turbulence at the site which in turn affect air flow patterns over the building and the resulting drifting. The next link in the chain is the local aerodynamics of the building, which is a function of the building's shape and the shapes of other buildings in the immediate neighbourhood. One additional factor for snow in this link is that the formation of drifts or ice accumulations can alter the shape of the structure and this in turn affects the air flow. Some drifts such as those at roof steps "fill in" areas of aerodynamic shade and are thereby limited as to how deep they can get.

The fourth link in the chain for wind is the influence of the properties of the building. For wind loading these would primarily be the dynamic properties, i.e. the frequencies of the natural modes of vibration, the mass distribution and the damping. For snow the dynamic properties of the building are not relevant but several other properties can be very important. Two of these are the internal temperature and thermal resistance of the roof since these affect melting rate. Two more are the slope and slipperiness of the roof in so far as they affect re-distribution of the load through sliding, which may also result in impact loads on lower level roofs. For some buildings snow accumulations are countered by operational measures. Operational measures can involve surface heating to melt the snow, or procedures to remove snow manually or mechanically before it can accumulate excessively. Another very important property of the building is drainage. The drainage system needs to be able to be of an appropriate capacity to rid the building of melt water without blocking or re-freezing.

The final link in the chain is to meet criteria for structural safety and serviceability. For snow this involves synthesizing the effects of all the parameters that affect snow loading so that quantities such as the 50 or 100 year return period loads can be determined and serviceability problems kept to an appropriate low level of probability.

## 3 Codes and standards

Because of the large number of factors that influence snow loads to estimate loads appropriate for use in design with any precision is very complex. However, codes and standards must attempt to reduce the complexities and provide relatively simple rules for designers. The author's experience with snow loads is primarily in Canada and the USA. The National Building Code of Canada (NBCC [7]) has evolved considerably since the 1940's. In the 1940's simple uniform loads were specified as a blanket provision (Allen [1]). By 1980 drift loads and shifting of loads due to sliding had been introduced (Taylor [3]) and by 1995 the effect of roof size on both uniform and drift loads had been introduced (Irwin, et al [8]). Similar developments were taking place in the USA, (O'Rourke et al [9]). Through involvement in the ISO Standard the author has become more familiar with practices in other snowy countries such as Norway, Japan, Germany, Poland and Russia. Interesting differences in approach can be seen in different countries. For example in Japan snow removal is recognized as a method of controlling

snow loads whereas in Canada, for example, it is not. The ISO Standard was updated in 2013 and, as with the national codes and standards contains provisions for uniform loads, drift loads and sliding. Building codes and standards use the recorded ground snow loads (or depths as substitute, along with assumptions concerning density) and then develop factors that are applied for exposure, drifting, and sliding. In some codes the effect of the internal temperature of the building is accounted for by a thermal factor.

It should be mentioned that many building codes and standards permit and encourage special, project specific studies to determine wind loads by wind tunnel testing as an alternative to using the standard analytical provisions. In 2012 the American Society for Civil Engineers published a standalone standard, ASCE/SEI 49-12 [6], covering the requirements for the wind tunnel testing methods. It is noteworthy that the ASCE 49-12 Standard also covers use of scale model testing and associated analysis methods for establishing snow loads on buildings.

## 4 Physical modeling and computer modeling of particle drifting

Because drifting has been the cause of many roof failures (O'Rourke and Wikoff [10]) a number of researchers have investigated methods to simulate the particle drifting process using physical models in wind tunnels and water flumes (Isyumov [11], Iversen [12], Irwin and Williams [13], Anno [14], Peterka and Petersen [15], O'Rourke et al [16]. The scaling principles were discussed by Kind [17]. It is quite difficult to match all the required scaling parameters but nonetheless useful information can be obtained on typical drift patterns and on potential methods for mitigating drifts.

<figure><img src="image_375.png" /><figcaption>Figure 2. Physical modeling in water flume</figcaption></figure>

Figure 2 shows an example of a drifting study in RWDI's water flume using crushed glass as the model particles. The physical model approach is best suited to the situation where drifts are formed by wind from a single direction in a single event. They are also useful in showing the limiting shapes for drifts as they fill in areas of aerodynamic shade. However, modeling the cumulative effects of different wind directions and speeds becomes cumbersome to implement on physical models.
Computational Fluid Dynamics (CFD) methods have also be applied to the particle drifting process: Beyers et al [18], Thiis and Ramberg [19]. As with physical models the CFD approach is easiest to implement for a single event with a predominant wind direction. However, it does in principle have some flexibility to examine the cumulative effects of varying directions and speeds.

## 5 Full time history modeling

To estimate extreme snow loads corresponding to a given return period requires that the cumulative effects of all the parameters listed in Figure 1 be evaluated with reasonable precision over a sufficient number of winters to obtain statistically meaningful results. The statistics of ground loads have themselves been analyzed by many researchers and meaningful relationships developed for ground loads as a function of return period. Building codes and standards rely on these relationships. However, notwithstanding the assumptions of building codes, the true relationship between 50 year ground load and roof load is not one to one and varies widely depending on circumstances.

In Canada, Isyumov and Davenport [20] made the first serious attempt to simulate the detailed processes of snow accumulation and melting on roofs over protracted periods using statistical methods. They examined the case of flat roofs in different regions of Canada using a Monte Carlo simulation method based on the statistics of snow fall, wind and temperature. Further work along these lines was reported by Isyumov and Mikiutiuk [21]. The numerical approach, although approximate, allowed the importance of various parameters to be identified. For example in less cold areas peak snow loading is more likely to be associated with a single event, whereas in cold areas the prolonged sub-freezing temperatures meant that peak snow accumulations were more likely to be due to cumulative effects over many weeks.

Subsequently Irwin and Gamble [22] developed what they named the Finite Area Element (FAE) method, which allowed the numerical approach to be extended to roofs of any shape. They developed the method during the design of the 208 m span retractable roof of the Toronto SkyDome (now named the Rogers Centre), using wind tunnel measurements of wind velocities near the roof as input into their numerical simulations of snow fall, drifting, and melting. The following description of the FAE method is taken largely from Irwin et al [8]. In the FAE method the roof is divided into a grid forming a large number of small area elements, Fig. 3, and the mean wind velocity vectors at the grid node points at height about 1 m (or other similar value) above the roof are determined by wind tunnel tests, Figure 4 (or alternatively by CFD). This is done for wind direction intervals such as 10, 15, or 22.5 degrees and empirical relationships based on field data are used to relate drift rate to wind velocity and direction. All velocities are determined as a fraction of a reference value, high above the building. When the results of the wind tunnel tests are combined with the hourly meteorological wind records, a

set of velocities over the roof corresponding to any particular date and time can be established. Alternatively, for approximate studies on simple shapes such as flat roofs the velocities can be estimated, without wind tunnel tests, directly from the meteorological records and general knowledge of wind flow over buildings.

<figure><img src="image_376.png" /><figcaption>Figure 3. Finite area element (FAE) method</figcaption></figure>

<figure><img src="image_377.png" /><figcaption>Figure 4. Measurement of surface wind velocities in an RWDI wind tunnel for FAE method</figcaption></figure>

In the FAE method the computations are completed for every area element and, by evaluating the snow fluxes into and out of each area, the change in snow mass contained within each element during a selected time interval, usually 1 hour, is determined. This procedure is repeated for a whole series of time steps, covering first one winter, and then a whole series of winters. If the meteorological records show that, in addition to drifting, snow was falling during a particular time interval, the change in snow mass in each element due to this is incorporated. Likewise, if the heat balance within the snowpack permits melting, then the mass of snow turning to the liquid state is computed. If the snow is saturated,
liquid water is allowed to drain out; otherwise it is held within the snowpack. Rain falling onto the snow is treated as additional liquid water in this calculation. Refreezing when the temperature drops is also computed. There is also an on-off switch that monitors whether the snow has been rained upon or has melted, in which case drifting is turned off. If new snow subsequently arrives, then drifting is turned back on, at least until all the new snow has drifted away and the old layer of snow below is reached again. The end result of an FAE simulation is a complete time history of the average snow load in each area element on the roof for every hour for about 30 to 50 winters typically. Then the output is analyzed statistically to determine extreme snow loads or extreme values of other meteorological quantities. This allows computation of quantities such as 50 year return period uniform and drift loads. It is worth noting that if wind pressure coefficients are measured on the wind tunnel model in addition to wind velocities then the FAE method lends itself to the evaluation of wind and snow load combinations.

<figure><img src="image_378.png" /><figcaption>**Figure 5.** Heat transfer parameters for a roof</figcaption></figure>

The evaluation of melting and water retention in the snow pack is complex. Figure 5 illustrates the various factors affecting melting. They include direct and indirect solar radiation, long-wave radiation, convective heat transfer, internal temperature, roof insulation value, and various snow pack properties such as thermal conductivity, latent heat, liquid water retention, albedo, crust formation etc. Typically the FAE simulation is run on the ground load initially for a specific project and is calibrated to the code ground load. The calibration is accomplished through adjustments of the assumed snow pack computer model until a good alignment with the code ground load is achieved. Then the FAE simulation is applied to the roof and snow load factors developed. Figure 6 shows a time history of simulated ground load in Montreal for the winter of 1981 – 1982 which was a fairly typical winter. It can be seen that the ground load builds up and reduces several times during the winter due to snow falls followed by melting. In contrast the winter of 1977-1978 shown in Figure 7 was much colder and the ground load continued to build up over almost the entire winter. There were a large number of roof collapses in the Northeastern part of North America in that year. At locations where snow does not stay for long such as Denver, Colorado, studies using the FAE method have

shown that heat transfer through the roof has much less influence on the snow load since the snow arrives in a large fall and heat transfer, although it may reduce the time the snow stays on the roof, does not have time to have much influence the maximum load. In an area like Montreal heat transfer through the roof is more likely to affect the maximum load.

<figure><img src="image_379.png" /><figcaption>**Figure 6.** Ground load simulation for Montreal for the winter of 1981-1982</figcaption></figure>

<figure><img src="image_380.png" /><figcaption>**Figure 7.** Ground snow load simulation for Montreal for the winter of 1977-1978</figcaption></figure>

Figure 8 shows as an example the application of the FAE method to the roof of Terminal 1 at Toronto's Pearson International Airport. The snow load in each area element is colour coded so that areas of heavy snow load (red in the diagram) can be seen at a glance. Figure 9 shows a photograph of actual snow accumulations on the roof seen during construction. The FAE method gives a good indication of the location and magnitude of snow loads that a large area roof such as this can experience.

Although the FAE method was developed initially for large one-of-a-kind projects, it was subsequently applied to the development of improved provisions in the 1995 edition of the NBCC, specifically to account for the effect of roof size on loads on flat roofs and in roof steps (Irwin et al, [8]). Later the method was also used to improve the NBCC load provisions on arch roofs (Hochstenbach et al [23]). More recently Dale et al [24] have used the parametric modeling capabilities of the FAE method to examine the effect of building orientation and height on roof step loads, as well as the effects of
internal temperature and roof thermal resistance, and the variation of loads averaged over different roof areas.

<figure><img src="image_381.png" /><figcaption>**Figure 8.** Output from and FAE study on Terminal 1 at Pearson International Airport, Toronto</figcaption></figure>

<figure><img src="image_382.png" /><figcaption>**Figure 9.** Observed snow on Terminal 1 at Pearson International Airport, Toronto</figcaption></figure>

Figures 8 and 9 are both snapshots in time of roof loading, as simulated and as seen on a real roof. For design extreme value analysis of the annual maxima is needed. Figure 10 shows and example of extreme value analysis of the load in one of the valleys on the fabric roof of Denver International Airport. The analysis uses the simulated highest annual load from about 30 winters.

<figure><img src="image_383.png" /><figcaption>**Figure 10.** Example of extreme value analysis of snow load on a portion of roof</figcaption></figure>

## 6 Future directions

Because of the issue of climate change the energy consumption of buildings has become an increasing preoccupation of building designers. At the turn of this century buildings were responsible for approximately 35% of greenhouse gas emissions (Battles and Burns, [25]). Since then the rapid trend towards more energy efficient buildings implies reduced heat loss. This has the potential to increase snow loads and designers are becoming more aware of the need for more detailed studies to address this issue.

Another trend is the use of external shading devices to reduce solar energy input into buildings. These features can attract snow and ice accumulations and increase the likelihood of falling snow and ice. Care is required to examine this problem during design and, if necessary, develop appropriate solutions, such as heat tracing to cause melting, snow retention devices, and devices to break ice into small pieces less likely to cause injury or damage. Tests in climatic simulation chambers (such as the Jules Verne facility here in Nantes) that can generate both snow and wind can be very instructive in identifying problems and developing solutions for architects, Fig. 10 shows an example of such testing in the test facilities at the University of Ontario Institute of Technology.

<figure><img src="image_384.png" /><figcaption>**Figure 11.** Example of tests in a climatic simulation chamber for wind, snow and ice</figcaption></figure>

The clean energy trend has also led to the installation of solar panels on many buildings, in some cases in large arrays. The snow load patterns induced by the collectors are different from a bare roof and can be usefully studied using scale models, Irwin et al [26], Brooks et al. [27]. An array of solar panels greatly alters the air flow patterns on roofs which leads to different drift formations.
<figure><img src="image_385.png" /><figcaption>**Figure 12.** Water flume simulation of snow drifting on a bare flat roof for quartering winds</figcaption></figure>

<figure><img src="image_386.png" /><figcaption>**Figure 13.** Modification of drifts due to a solar panel array for quartering winds from the northwest</figcaption></figure>

Figure 12 shows a water flume simulation of snow drifting and scouring for quartering winds. A dominant feature for this wind direction is the scouring that occurs under the vortices generated by the upwind corner. Figure 13 shows the change in drift patterns due to the solar panel array for winds from the northwest. The single vortex on the northeast side of the roof breaks into multiple vortices and the vortex on the southwest half of the roof tends to be spread out more. The gap beneath a solar panel can be blocked off by sliding snow which can cause further changes in snow load distribution. Further research into the interaction of snow with solar panels is needed.

As indicated earlier the combination of wind and snow loads can benefit from time history simulation methods since both wind and snow loads can be tracked simultaneously. Bond et al [28] describe studies along these lines and more can be expected in future. Although snow loads on roofs act downwards and in most cases the dominant wind loads are upwards there can be some situations, such as a roof next to a taller building, or a roof with unbalance snow, where the snow and wind combination might govern.

<figure><img src="image_387.png" /><figcaption>**Figure 14.** Variation of critical speed for flutter with radial snow and ice thickness</figcaption></figure>

Combined wind and snow loads are an example of a multi-hazard situation. Besides loading, aerodynamic stability can also be influenced by snow and ice accumulations. For example the aerodynamic stability of flexible structures such as long span bridges can be affected by snow or ice accumulations. The current trend in North America is for long span bridges to have tall “anti-climb” fences at the deck edges. These are designed to be very porous, having no more than about 10% to 15% solidity, but snow or ice accumulations can greatly increase their solidity which can then adversely change the aerodynamic behaviour of the deck, making it far less stable against vortex excitation and flutter. This is illustrated in Fig.14 where the critical speed, as determined by sectional model tests, for the onset of flutter is plotted versus radial thickness of ice on the fence elements. It can be seen that the ice thickness can have a very adverse effect in reducing the flutter speed substantially.

In these situations a combined risk approach based on historical records of wind, snow and freezing rain, as well as wind tunnel tests of aerodynamic stability with and without fences blocked to varying degrees, can help determine the true risk of problems occurring. Table 1 illustrates the risk associated with ice being left on the fence for different numbers of days. The target return period for flutter was 10,000 years in this case. Therefore if a maintenance protocol can ensure ice does not remain for more than four days then the risk of flutter is acceptable.

<table><caption>**Table 1.** Effect of number of days of exposure to icing on return period for flutter of a bridge</caption><thead><tr><td>Assumed Exposure (days)</td><td>Return Period of Flutter (years)</td></tr></thead><tbody><tr><td>1</td><td>30,000</td></tr><tr><td>4</td><td>10,000</td></tr><tr><td>7</td><td>6,300</td></tr><tr><td>10</td><td>4,500</td></tr><tr><td>15</td><td>3,000</td></tr></tbody></table>
The sliding of snow off roofs is still an area of
significant uncertainty. The NBCC simply classifies roofs
as slippery or non-slippery. Further research is needed
into the sliding process. What may be a conservative
assumption for an upper roof is non-conservative for a
lower roof. Lepage and Schuyler [29] described the
development of computer and scale model simulation
methods for predicting the sliding of snow off sloped
roofs or the lifting off of snow mass by wind forces.
These phenomena are of particular concern on high rise
buildings because of the hazard to pedestrians below. The
simulation techniques enable the effectiveness of various
mitigation measures to be assessed. The sliding snow and
ice problem has been more recently discussed by
Williams et al [4] and by Isyumov [5].

An important continuing task is the validation of
simulation methods. These include scale model
simulations, CFD simulations, and full time history
methods such as FAE. Here field data is invaluable as a
means of calibrating their accuracy. Irwin [30] describes
comparison of the FAE method with field observations
on a large area stepped roof in Canada. The observations
were taken shortly after a heavy snow loading event that
caused part of the roof to collapse. During the
investigations the FAE method was run for the winter
months leading up to the collapse and was used to
estimate the loads in the collapsed part. It was calibrated
by comparing it with observations of load in roof steps in
un-collapsed parts of the roof. The indications were that
the FAE method tended to overestimate the capture
efficiency of the step by about 20% and this allowed the
results in the collapsed portion to be refined. This
illustrates the importance of capitalizing on any full scale
data that becomes available.

7 Conclusion

Snow loads continue to cause roof failures and we still
have much to learn. Snow engineering is in many ways
similar to wind engineering, only more complex due to
the increased number of variables involved and the fact
that in many cases cumulative effects over time must be
accounted for. In wind engineering the forces are
generated over a short interval but in snow engineering
the load is affected by the history of snow fall, wind
speed and direction, temperature, heat transfer through
the building, sliding and retention of melt water and rain
by the snow pack. In some cases the history leading up to
the maximum load can be many weeks long. The concept
of the A.G. Davenport wind loading chain used in wind
engineering can nonetheless serve as a good guide as to
what needs to be taken into account. There has been
progress in developing both physical and CFD models to
simulate single events of drifting, and numerical
techniques to track the multiple climatic variables and
their effects on the snow. These modeling techniques,
although often approximate can nonetheless provide good
guidance on probable extreme loads and on ways of
avoiding problems of overloading and sliding hazards.
Much of what has been learned about snow loading has

come from case histories and field observation. It cannot
be emphasized enough how important it is that
information continue to be gathered from case histories
and field observation. Only with such data can prediction
methods continue to be improved in future.

References

1. Allen, D.E., *Snow loads on roofs*, National Research Council of Canada, DBR Internal Report 106, 1956.

2. Schriever, W.R., Lutes, D.A., and Peter, B.G.W., *Snow Loads*, National Building Code of Canada, Supplement No. 4., Commentary No.2, pp. 568-573, 1970.

3. Taylor, D.A., *Roof snow loads in Canada*, Canadian Journal of Civil Engineering, Vol. 7, No.1, pp.1-18, 1980.

4. Williams, C.J., M. Carter, M., Hochstenbach, F., and Lovlin, T., *Sliding snow and ice on Buildings: a balance of risk cost and aesthetics*, Proceedings of Snow Engineering V, 5th International Conference on Snow Engineering, Davos Switzerland, p. 59, 2004

5. Isyumov, N., *Sliding Snow from Sloped Roofs: Its Prediction, Potential Hazard and Mitigation*, Proceedings of Snow Engineering VI sixth International Conference on Snow Engineering, Whistler, British Columbia, Canada, 2008.

6. ASCE/SEI 49-12, Wind tunnel testing for buildings and other structures, American Society for Civil Engineering/ Structural Engineering Institute, Standard 49-12, 2012.

7. NBCC, National Building Code of Canada, Various Editions, published by the National Research Council, 1970, 1975, 1980, 1985, 1990, 1995, 2005, 2010, 2015.

8. Irwin, P.A., Gamble, S.L and Taylor, D.A. 1995. *Effects of Roof Size, Heat Transfer, and Climate on Snow Loads: Studies for the 1995 NBC*. Canadian Journal of Civil Engineering, Vol. 22: 770-784.

9. O'Rourke, M., Tobiasson, W., and Wood, E., *Proposed code provisions for drifted snow loads*, ASCE J. Structural of Engineering, Vol. 112, No. 9, 1986.

10. O'Rourke, M., and Wikoff, J., *Snow related collapse and implications for building codes*, Structure Magazine, pp. 18-21, January 2013.

11. Isyumov, N., *An approach to the prediction of snow loads*, PhD Thesis, University of Western Ontario, 1971.

12. Iversen, J.D., *Drifting snow similitude – transport rate and roughness modeling*, Scientific Symposium on Snow Motion, Fort Collins, Colorado, August 12-17, 1979.

13. Irwin, P.A., and Williams, C.J., Application of snow simulation model tests to planning and design, Proceedings of the Eastern Snow Conference, Vol. 28, 40th Annual meeting, pp. 118-130, 1983.

14. Anno, Y., *Snow drift wind tunnels in Japan*, Proceedings First International Conference on Snow