Consistently use frac_sno_eff rather than frac_sno in qg calculations

The change here only changes answers for (1) urban pervious road, and
(2) more broadly with use_subgrid_fluxes = .false.

Previously, code in CalculateSurfaceHumidity was inconsistent regarding
the use of frac_sno vs. frac_sno_eff. I wanted to make this consistent
to facilitate upcoming cleanup of CalculateSurfaceHumidity. I went back
and forth as to which one we should use here. (See also
#822).

Argument for frac_sno_eff: This is consistent with the use of
frac_sno_eff elsewhere in the calculation of surface fluxes.

Argument for frac_sno: Even in CLM4.0 (prior to the introduction of
frac_sno_eff and frac_h2osfc), frac_sno was used to calculate surface
relative humidity. That makes me think that the parameterizations were
improved by using frac_sno here, even though frac_sno at that point was
mainly just used for albedo / radiation calculations. So sticking with
frac_sno would be more consistent with what has been done for a long
time, and would presumably be a bit better scientifically.

However, when I tried changing this to use frac_sno rather than
frac_sno_eff:

    qg(c) = frac_sno_eff(c)*qg_snow(c) + (1._r8 - frac_sno_eff(c) - frac_h2osfc(c))*qg_soil(c) &
          + frac_h2osfc(c) * qg_h2osfc(c)

I got a crash with use_subgrid_fluxes false (test
SMS_D_Ld1_P4x1.f10_f10_musgs.I2000Clm45BgcCropQianRsGs.bishorn_gnu.clm-no_subgrid_fluxes),
due to a check I put in relatively recently to ensure that the top-layer
soil liquid never goes significantly negative:

WJS 218: frac_sno, frac_sno_eff, frac_h2osfc, snl =    6.6447358408778726E-005   1.0000000000000000        0.0000000000000000               -1
 qg_snow, qg_soil, qg_h2osfc, qg =    3.8414454029812475E-003   3.6634834595404343E-003   3.6635176705822736E-003   3.6634952846414733E-003
 ERROR: In UpdateState_TopLayerFluxes, h2osoi_liq has gone significantly negative
 Bulk/tracer name = bulk
 c, lev_top(c) =          218           0
 h2osoi_liq_top_orig  =    6.6448689668604415E-004
 h2osoi_liq           =   -2.6754040213948407E-004
 frac_sno_eff         =    1.0000000000000000
 qflx_liq_grnd*dtime  =    0.0000000000000000
 qflx_dew_grnd*dtime  =    0.0000000000000000
 qflx_evap_grnd*dtime =    9.3202729882552822E-004
 ENDRUN:
 ERROR: In UpdateState_TopLayerFluxes, h2osoi_liq has gone significantly negative

What seems to be going on is: With a little bit of snow present, when I
use frac_sno rather than frac_sno_eff, we end up with a qg value very
close to that of soil. But I think the later evaporation happens from
the snow, at which point I guess the dryer air above leads to there
being too large evaporation from the snow.

So in the end, I decided to consistently use frac_sno_eff, even though
this might be a small step backwards in terms of the accuracy of surface
humidity for urban pervious road and with use_subgrid_fluxes = .false.

src/biogeophys/SurfaceHumidityMod.F90

@@ -93,7 +93,6 @@ subroutine CalculateSurfaceHumidity(bounds, &
 
          frac_h2osfc      =>    waterdiagnosticbulk_inst%frac_h2osfc_col    , & ! Input:  [real(r8) (:)   ] fraction of ground covered by surface water (0 to 1)
          frac_sno_eff     =>    waterdiagnosticbulk_inst%frac_sno_eff_col   , & ! Input:  [real(r8) (:)   ] eff. fraction of ground covered by snow (0 to 1)
-         frac_sno         =>    waterdiagnosticbulk_inst%frac_sno_col       , & ! Input:  [real(r8) (:)   ] fraction of ground covered by snow (0 to 1)
          h2osoi_ice       =>    waterstatebulk_inst%h2osoi_ice_col          , & ! Input:  [real(r8) (:,:) ] ice lens (kg/m2)
          h2osoi_liq       =>    waterstatebulk_inst%h2osoi_liq_col          , & ! Input:  [real(r8) (:,:) ] liquid water (kg/m2)
          qg_snow          =>    waterdiagnosticbulk_inst%qg_snow_col        , & ! Output: [real(r8) (:)   ] specific humidity at snow surface [kg/kg]
@@ -139,8 +138,8 @@ subroutine CalculateSurfaceHumidity(bounds, &
                psit = max(smpmin(c), psit)
                ! modify qred to account for h2osfc
                hr   = exp(psit/roverg/t_soisno(c,1))
-               qred = (1._r8 - frac_sno(c) - frac_h2osfc(c))*hr &
-                    + frac_sno(c) + frac_h2osfc(c)
+               qred = (1._r8 - frac_sno_eff(c) - frac_h2osfc(c))*hr &
+                    + frac_sno_eff(c) + frac_h2osfc(c)
                soilalpha(c) = qred
 
             else if (col%itype(c) == icol_road_perv) then
@@ -158,7 +157,7 @@ subroutine CalculateSurfaceHumidity(bounds, &
                   hr_road_perv = hr_road_perv + rootr_road_perv(c,j)
                end do
                ! Allows for sublimation of snow or dew on snow
-               qred = (1.-frac_sno(c))*hr_road_perv + frac_sno(c)
+               qred = (1.-frac_sno_eff(c))*hr_road_perv + frac_sno_eff(c)
 
                ! Normalize root resistances to get layer contribution to total ET
                if (hr_road_perv > 0._r8) then
@@ -195,8 +194,8 @@ subroutine CalculateSurfaceHumidity(bounds, &
             if (snl(c) < 0) then
                call QSat(t_soisno(c,snl(c)+1), forc_pbot(c), eg, degdT, qsatg, qsatgdT_snow)
                qg_snow(c) = qsatg
-               dqgdT(c) = frac_sno(c)*qsatgdT_snow + &
-                    (1._r8 - frac_sno(c) - frac_h2osfc(c))*hr*qsatgdT_soil
+               dqgdT(c) = frac_sno_eff(c)*qsatgdT_snow + &
+                    (1._r8 - frac_sno_eff(c) - frac_h2osfc(c))*hr*qsatgdT_soil
             else
                ! To be consistent with hs_top values in SoilTemp, set qg_snow to qg_soil
                ! for snl = 0 case. This ensures hs_top_snow will equal hs_top_soil.
@@ -212,7 +211,6 @@ subroutine CalculateSurfaceHumidity(bounds, &
                qg_h2osfc(c) = qg_soil(c)
             end if
 
-            !          qg(c) = frac_sno(c)*qg_snow(c) + (1._r8 - frac_sno(c) - frac_h2osfc(c))*qg_soil(c) &
             qg(c) = frac_sno_eff(c)*qg_snow(c) + (1._r8 - frac_sno_eff(c) - frac_h2osfc(c))*qg_soil(c) &
                  + frac_h2osfc(c) * qg_h2osfc(c)