sensible+latent heatfluxes using linear bulk formula (#371)

* 'heatflux_linear' flag: sensible+latent heatfluxes using traditional linear bulk formula

* ERROR correction: heatflux_linear is logical

* ERROR/syntax correction.

* Add option atmbndy='mixed' boundary layer condition

* For 'atmbndy': change 'default' (obsolete) to 'similarity'. Same physics

* For 'atmbndy': change 'default' (obsolete) to 'similarity'. Same physics

* New options for 'atmbndy'

* "atmbndy" options: 'similarity', 'constant', 'mixed'

* Abort if "atmbndy" option is unknown + set default=similarity for backward compability

* Update doc/source/science_guide/sg_boundary_forcing.rst

Co-authored-by: Philippe Blain <levraiphilippeblain@gmail.com>

* atmbndy option 'ccsm' not an option. Previous given as 'default'

* namelist parameter plus string inside double quotes

* Re-interduce file, after accidentally git rm command

* Remove "ccsm" as an option to atmbndy

* atmbndy: 'constant' or 'mixed'. Anything else end in 'similarity' behaviour

* Error correction. Revert atmbndy to constant as before

* Put constants intp icepack_parameters

* Introduce p0012, p0015 from icepack_parameters

* Changed p0012/p0015 to specific constants senscoef/latncoef

* senscoef/latncoef introduced

Co-authored-by: Philippe Blain <levraiphilippeblain@gmail.com>

columnphysics/icepack_atmo.F90

@@ -15,6 +15,7 @@ module icepack_atmo
       use icepack_kinds
       use icepack_parameters,  only: c0, c1, c2, c4, c5, c8, c10
       use icepack_parameters,  only: c16, c20, p001, p01, p2, p4, p5, p75, puny
+      use icepack_parameters,  only: senscoef, latncoef
       use icepack_parameters,  only: cp_wv, cp_air, iceruf, zref, qqqice, TTTice, qqqocn, TTTocn
       use icepack_parameters,  only: Lsub, Lvap, vonkar, Tffresh, zvir, gravit
       use icepack_parameters,  only: pih, dragio, rhoi, rhos, rhow
@@ -359,8 +360,16 @@ subroutine atmo_boundary_layer (sfctype,            &
       ! as in Jordan et al (JGR, 1999)
       !------------------------------------------------------------
 
-      shcoef = rhoa * ustar * cp * rh + c1
-      lhcoef = rhoa * ustar * Lheat  * re
+      if (trim(atmbndy) == 'mixed') then
+         !- Use constant coefficients for sensible and latent heat fluxes
+         !    similar to atmo_boundary_const but using vmag instead of wind
+         shcoef = senscoef*cp_air*rhoa*vmag
+         lhcoef = latncoef*Lheat *rhoa*vmag
+      else ! 'similarity'
+         !- Monin-Obukhov similarity theory for boundary layer
+         shcoef = rhoa * ustar * cp * rh + c1
+         lhcoef = rhoa * ustar * Lheat  * re
+      endif
 
       !------------------------------------------------------------
       ! Compute diagnostics: 2m ref T, Q, U
@@ -503,8 +512,8 @@ subroutine atmo_boundary_const (sfctype,  calc_strair, &
       ! coefficients for turbulent flux calculation
       !------------------------------------------------------------
 
-      shcoef = (1.20e-3_dbl_kind)*cp_air*rhoa*wind
-      lhcoef = (1.50e-3_dbl_kind)*Lheat *rhoa*wind
+      shcoef = senscoef*cp_air*rhoa*wind
+      lhcoef = latncoef*Lheat *rhoa*wind
 
       end subroutine atmo_boundary_const
 
@@ -952,7 +961,7 @@ subroutine icepack_atm_boundary(sfctype,                   &
                                    delt,     delq,     &
                                    lhcoef,   shcoef    )
          if (icepack_warnings_aborted(subname)) return
-      else ! default
+      else
          call atmo_boundary_layer (sfctype,                 &
                                    calc_strair, formdrag,   &
                                    Tsf,      potT,          &

columnphysics/icepack_parameters.F90

@@ -247,10 +247,12 @@ module icepack_parameters
          qqqice = 11637800._dbl_kind ,&! for qsat over ice
          TTTice = 5897.8_dbl_kind    ,&! for qsat over ice
          qqqocn = 627572.4_dbl_kind  ,&! for qsat over ocn
-         TTTocn = 5107.4_dbl_kind      ! for qsat over ocn
+         TTTocn = 5107.4_dbl_kind    ,&! for qsat over ocn
+         senscoef= 0.0012_dbl_kind   ,&! Sensible heat flux coefficient for constant-based boundary layer
+         latncoef= 0.0015_dbl_kind     ! Latent heat flux coefficient for constant-based boundary layer
 
       character (len=char_len), public :: &
-         atmbndy = 'default' ! atmo boundary method, 'default' ('ccsm3') or 'constant'
+         atmbndy = 'similarity'        ! atmo boundary method, 'similarity', 'constant' or 'mixed'
 
       logical (kind=log_kind), public :: &
          calc_strair     = .true.  , & ! if true, calculate wind stress
@@ -646,12 +648,12 @@ subroutine icepack_init_parameters(   &
          TTTocn_in        ! for qsat over ocn
 
       character (len=*), intent(in), optional :: &
-         atmbndy_in ! atmo boundary method, 'default' ('ccsm3') or 'constant'
+         atmbndy_in   ! atmo boundary method, 'similarity', 'constant' or 'mixed'
 
       logical (kind=log_kind), intent(in), optional :: &
-         calc_strair_in, &  ! if true, calculate wind stress components
-         formdrag_in,    &  ! if true, calculate form drag
-         highfreq_in        ! if true, use high frequency coupling
+         calc_strair_in,     & ! if true, calculate wind stress components
+         formdrag_in,        & ! if true, calculate form drag
+         highfreq_in           ! if true, use high frequency coupling
 
       integer (kind=int_kind), intent(in), optional :: &
          natmiter_in        ! number of iterations for boundary layer calculations
@@ -1326,12 +1328,12 @@ subroutine icepack_query_parameters(   &
          TTTocn_out        ! for qsat over ocn
 
       character (len=*), intent(out), optional :: &
-         atmbndy_out ! atmo boundary method, 'default' ('ccsm3') or 'constant'
+         atmbndy_out   ! atmo boundary method, 'similarity', 'constant' or 'mixed'
 
       logical (kind=log_kind), intent(out), optional :: &
-         calc_strair_out, &  ! if true, calculate wind stress components
-         formdrag_out,    &  ! if true, calculate form drag
-         highfreq_out        ! if true, use high frequency coupling
+         calc_strair_out,     & ! if true, calculate wind stress components
+         formdrag_out,        & ! if true, calculate form drag
+         highfreq_out           ! if true, use high frequency coupling
 
       integer (kind=int_kind), intent(out), optional :: &
          natmiter_out        ! number of iterations for boundary layer calculations

configuration/driver/icedrv_init.F90

@@ -563,8 +563,8 @@ subroutine input_data
       if (formdrag) then
       if (trim(atmbndy) == 'constant') then
          write (nu_diag,*) 'WARNING: atmbndy = constant not allowed with formdrag'
-         write (nu_diag,*) 'WARNING: Setting atmbndy = default'
-         atmbndy = 'default'
+         write (nu_diag,*) 'WARNING: Setting atmbndy = similarity'
+         atmbndy = 'similarity'
       endif
 
       if (.not. calc_strair) then

configuration/scripts/icepack_in

@@ -97,7 +97,7 @@
 
 &forcing_nml
     formdrag        = .false.
-    atmbndy         = 'default'
+    atmbndy         = 'similarity'
     calc_strair     = .true.
     calc_Tsfc       = .true.
     highfreq        = .false.

doc/source/icepack_index.rst

@@ -52,7 +52,7 @@ either Celsius or Kelvin units).
    "atmiter_conv", ":math:`\bullet` convergence criteria for ustar", "0.0"
    "atm_data_format", ":math:`\bullet` format of atmospheric forcing files", ""
    "atm_data_type", ":math:`\bullet` type of atmospheric forcing", ""
-   "atmbndy", ":math:`\bullet` atmo boundary layer parameterization (‘default’ or ‘constant’)", ""
+   "atmbndy", ":math:`\bullet` atmo boundary layer parameterization (‘similarity’,‘constant’ or 'mixed')", ""
    "awtidf", "weighting factor for near-ir, diffuse albedo", "0.36218"
    "awtidr", "weighting factor for near-ir, direct albedo", "0.00182"
    "awtvdf", "weighting factor for visible, diffuse albedo", "0.63282"

doc/source/science_guide/sg_boundary_forcing.rst

@@ -120,7 +120,12 @@ parameterizations. Rain and all melted snow end up in the ocean.
 Wind stress and transfer coefficients for the
 turbulent heat fluxes are computed in subroutine
 *atmo\_boundary\_layer* following :cite:`Kauffman02`, with additions and changes as detailed in Appendix A of :cite:`Roberts15` for high frequency coupling (namelist variable ``highfreq``).
-The resulting equations are provided here.
+The resulting equations are provided here for the default boundary layer 
+scheme, which is based on Monin-Obukhov theory (``atmbndy = ‘stability’``). 
+Alternatively, ``atmbndy = ‘constant’`` provides constant coefficients for
+wind stress, sensible heat and latent heat calculations (computed in subroutine
+*atmo\_boundary\_const*); ``atmbndy = ‘mixed’`` uses the stability based 
+calculation for wind stress and constant coefficients for sensible and latent heat fluxes.
 
 The wind stress and turbulent heat flux calculation accounts for both
 stable and unstable atmosphere–ice boundary layers. We first define the

doc/source/user_guide/interfaces.include

@@ -944,7 +944,7 @@ icepack_init_parameters
            TTTocn_in        ! for qsat over ocn
 
         character (len=char_len), intent(in), optional :: &
-           atmbndy_in ! atmo boundary method, 'default' ('ccsm3') or 'constant'
+           atmbndy_in ! atmo boundary method, 'similarity', 'constant' or 'mixed'
 
         logical (kind=log_kind), intent(in), optional :: &
            calc_strair_in, &  ! if true, calculate wind stress components
@@ -1306,7 +1306,7 @@ icepack_query_parameters
            TTTocn_out        ! for qsat over ocn
 
         character (len=char_len), intent(out), optional :: &
-           atmbndy_out ! atmo boundary method, 'default' ('ccsm3') or 'constant'
+           atmbndy_out ! atmo boundary method, 'similarity', 'constant' or 'mixed'
 
         logical (kind=log_kind), intent(out), optional :: &
            calc_strair_out, &  ! if true, calculate wind stress components

doc/source/user_guide/ug_case_settings.rst

@@ -314,8 +314,10 @@ forcing_nml
    :widths: 15, 15, 30, 15 
 
    "", "", "", ""
-   "``atmbndy``", "``constant``", "bulk transfer coefficients", "``default``"
-   "", "``default``", "stability-based boundary layer", ""
+   "``atmbndy``", "string", "bulk transfer coefficients", "``similarity``"
+   "", "``similarity``", "stability-based boundary layer", ""
+   "", "``constant``", "constant-based boundary layer", ""
+   "", "``mixed``", "stability-based, but constant for sensible+latent heatfluxes", ""
    "``atmiter_conv``", "real", "convergence criteria for ustar", "0.0"
    "``atm_data_file``", "string", "file containing atmospheric data", "' '"
    "``atm_data_format``", "``bin``", "read direct access binary forcing files", "``bin``"