SecondBEST is a simplified implementation of the Bare Essentials for Surface Transfer (BEST) land surface model. The main simplifications are that SecondBEST does not incorporate the effects of vegetation, which reduces a lot of calculations within BEST. We also assume each grid cell is either fully bare soil or snow/ice, i.e, there are no fractional land types. The implementation here follows the description in Pitman et al.
This notebook is a ready reference to the equations used in SecondBEST so that interested users don't have to read the above reference if they don't wish to. However, if you plan on making heavy use of SecondBEST, you are encouraged to do so!
Surface types:
Spectral regions, phases of water
Water content
Certain soil properties are calculated initially based on the soil and area types.
Currently, SecondBEST recognises two soil types: clay and sand.
claysandclaysandland_iceland_ice.The albedo calculation is described in section 5.0 of Pitman et al. SecondBEST calculates albedo in the near-infraread and shortwave separately. This is also quite convenient since RRTMG, the radiative transfer code also uses two separate albedos.
The albedo for bare soil is calculated as (Eqns. 5.5, 5.6): $$\alpha_{SW,u} = 0.10 + 0.1clr + 0.06(1-W_{L,u})\\ \alpha_{LW,u} = 2\alpha_{SW,u}$$
where $clr$ is the soil colour defined above.
$W_{L,u}$ is the soil wetness factor, and is defined as the ratio of the liquid soil water content in the upper soil layer to its porosity.
$$W_{L,u} = X_{w,u}/X_v$$The albedo of Land Ice is calculated as (Eqns. 5.7, 5.8): $$\alpha_{SW,u} = 0.60 + 0.06(1-W_{L,u}) \\ \alpha_{LW,u} = \frac{1}{3}\alpha_{SW,u}$$
SecondBEST calculates the drag parameters as in section 6 of Pitman et al.
where $z_m$ is the hydrostatic mid-level height of the lowest atmosphere level and $\kappa$ is the von Karman constant.
Where $U_m$ is the magnitude of the wind speed at the lowest model level and the other notation is as before. $U_m$ has a minimum value that can be set during initialisation of SecondBEST.
A negative Richardson number indicates unstable conditions.
where $\epsilon$ is 0.01 for land and 1.0 for oceans.
The equations are from Section 8. in Pitman et al.
is given by the relatively well known bulk aerodynamic formula.
$$SHF = \rho_s C_{pd} U_m C_{Dh} (T_u - T_s)$$where $C_{pd}$ is the heat capacity at constant pressure of dry air, and $\rho_s$ is the density of air at the surface.
is also given by the bulk formula
$$LHF = L_v\rho_s U_m C_{Dh} \beta_u(q^*_u - q_s) $$where $q^*_u$ is the saturation specific humidity at the soil surface. A "wetness factor", $\beta_u$, is also required, defined as
$$\beta_u = \left \{ \begin{array}{1} 1,\ snow\\ \frac{W_{F,u}L_v}{L_i} + \frac{E_{usmax}}{E_{usP}},\ soil\end{array} \right .$$$L_v$ is the latent heat of vapourisation, $L_i$ is the latent heat of sublimation, $W_{F,u}$ is the frozen soil moisture in the upper soil layer. The first term containing these terms represents the rate at which frozen soil can sublimate.
$E_{usP}$ is the potential evaporation rate at the soil surface,
$$E_{usP} = \rho_s c_u (q^*_u - q_s)$$where $c_u$ is the soil conductance, defined as (Eq. 7.22)
$$c_u = C_{Dh}U_m.$$$E_{usmax}$ is the exfiltration rate given by
$$E_{usmax} = K_{HD}\Theta_u^{0.5B + 2} - K_{H0}\Theta_u^{2B + 3}$$Where
$$B = \left \{ \begin{array} 10,\ clay\\ 4,\ sand\end{array} \right .$$$$K_{H0} = \left \{ \begin{array} 0.001,\ clay\\ 0.1,\ sand \end{array} \right .$$and $$K_{HD} = \left ( \frac{-4K_{H0}B\Psi_0\rho_wX_v(1 - W_{F,u})}{\pi\Delta t}\right )$$
where
$\Theta_u$ is defined as
$$\Theta_u = \frac{W_{L,u} - 0.01}{1 - W_{F,u}}$$The wetness factors $W_{L,x}$ in the soil in the upper and lower layers is NOT allowed to fall below 0.01 due to physical and numerical issues.
are directly calculated from the momentum bulk coefficient $$U_s = -\rho_s C_{Dm}U_mu$$
The transfer of heat and moisture through the snow-ice-soil pack requires solving the diffusion equation along with local sources and sinks of heat and moisture which correspond to freezing/melting of water.
The conservation laws are:
$$c_vT_t = G_z + L_f\Gamma + c_{pw}RT_z - T_t(fict)$$$$\rho_w (X_w)_t = (R - (1-f_i)E)_z - \Gamma - \Pi$$$$\rho_i(X_i)_t = -(f_i E)_z + \Gamma$$where
$f_i = X_i/(X_w + X_i)$, $c_v$ is a volumetric heat capacity for soil.
$\Gamma = \frac{c_v}{L_f}\frac{(T_f - T)}{\Delta t}$ is the potential ice production rate, which is zero if there is no water for freezing or ice for melting.
$\Pi = A_vE_c/\Delta z$ is the water removed by transpiration, which is zero in this component (no vegetation effects)
And the flux densities are
$G = K_T T_z$
$R = \rho_w K_H(1 + \psi_z) = \rho_w K_H(z + \psi)_z$ for water
$E = \left \{ \begin{array}\\ 0,\ z < 0\\ E_{us},\ z=0 \end{array} \right .$ for water vapour
$T_t(fict)$ is a fictional source that is not necessary in our case since we are solving the full heat/mass transfer equations.
We use a Neumann condition for all variables $()_z=0$, representing a sink of energy and water to the deep soil.
We use an implicit scheme to step forward the diffusive (density) terms and an Euler scheme to step forward the other sources and sinks. The term $c_{pw}RT_z$ represents advection, and this will set the CFL condition $\Delta z/(c_{pw}R) > \Delta t/2$ where $\Delta t$ is the timestep to be used.