# Snow Load Shape Coefficients and Snow Prevention Method for Stepped Flat Roofs
## 3 Numerical Simulation  
### Simulation Scheme  
In this study, the Euler-Euler method was adopted in which both the air phase and snow phase are regarded as a continuum. The proportions of the two phases in the mixture multiphase model are defined by the volume fraction. The continuity equation and momentum equation of the mixed phase are solved by the mixture model.  
The continuity equation of the mixed phase is Equation (3):  
\[\frac{\partial}{\partial t}(\rho_{m})+\
abla\cdot\left(\rho_{m}\overset{ \rightarrow}{v}_{m}\right)=0. \tag{3}\]  
The mixed phase density (\(\rho_{m}\)) and the average mass velocity (\(\overset{\rightarrow}{v}_{m}\)) are determined using Equations (4) and (5):  
\[\rho_{m}=\sum_{k=1}^{2}\alpha_{k}\rho_{k}, \tag{4}\]  
\[\overset{\rightarrow}{v}_{m}=\frac{\sum\limits_{k=1}^{2}\alpha_{k}\rho_{k} \overset{\rightarrow}{v}_{k}}{\rho_{m}}, \tag{5}\]  
where \(k=1\), \(2\) represents the principal term (air phase) and secondary term (snow phase), \(\alpha_{k}\) is the volume fraction of the \(k\)th term, \(\rho_{k}\) is the density of the \(k\)th term, and \(\overset{\rightarrow}{v}_{k}\) is the velocity of the \(k\)th term.  
The momentum equation for the mixture model is expressed as Equation (6):  
\[\frac{\partial}{\partial t}\left(\rho_{m}\overset{\rightarrow}{v}_{m} \right)+\
abla\cdot\left(\rho_{m}\overset{\rightarrow}{v}_{m}\overset{ \rightarrow}{v}_{m}\right)=-\
abla p+\
abla\cdot\left[\mu_{m}\bigg{(}\
abla \overset{\rightarrow}{v}_{m}+\
abla\overset{\rightarrow}{v}_{m}^{T}\bigg{)} \right]+\rho_{m}\overset{\rightarrow}{g}+\overset{\rightarrow}{F}-\
abla \cdot\left(\sum_{k=1}^{2}\alpha_{k}\rho_{k}\overset{\rightarrow}{v}_{dr,k} \overset{\rightarrow}{v}_{dr,k}\right), \tag{6}\]  
where \(F\) is the volume force, and the mixed viscosity and the drift velocity of the \(k\)th term are calculated using Equations (7) and (8):  
\[\mu_{m}=\sum_{k=1}^{2}\alpha_{k}u_{k}, \tag{7}\]  
\[\overset{\rightarrow}{v}_{dr,k}=\overset{\rightarrow}{v}_{k}-\overset{ \rightarrow}{v}_{m}. \tag{8}\]  
The air phase is solved using the Reynolds equation [19]. The turbulence model adopts the realizable \(k\)-\(\varepsilon\) model [20,21], which improves the dissipation rate transport equation compared with the standard \(k\)-\(\varepsilon\) model. [PERSON] et al. (2016) [22] successfully simulated the wind field and snow distribution on stepped flat roofs with the help of the realizable \(k\)-\(\varepsilon\) model but did not consider the reverse impact of the snowpack on the wind field.  
The management of the snow transport and movement in the computational domain is achieved using the snow phase volume fraction (Equation (9)) as well as the erosion and deposition models (Equations (10) and (11)) [23, 24, 25]. The governing equations are as follows:  
\[\frac{\partial}{\partial t}(\alpha_{2}\rho_{2})+\
abla\cdot\left(\alpha_{2} \rho_{2}\overset{\rightarrow}{v}_{m}\right)=-\
abla\cdot\left(\alpha_{2} \rho_{2}\overset{\rightarrow}{v}_{dr,2}\right), \tag{9}\]  
\[q_{ero}=A_{ero}\left(u_{*}^{2}-u_{st}^{2}\right)\quad u_{*}\geq u_{*t}, \tag{10}\]  
\[q_{ap}=\phi w_{f}\frac{u_{*t}^{2}-u_{*}^{2}}{u_{*t}^{2}}\quad u_{*}<u_{st}, \tag{11}\]  
where \(A_{ero}\) is the erosion constant, equal to 7.0 \(\times\) 10\({}^{-4}\)[23]; \(u_{*}\) is the wall friction velocity; \(u_{*t}\) is the threshold friction velocity of snow drifting; and \(\Phi\) is the snow mass concentration.