Better document parameter units in 4 user modules

  Add units arguments to 17 unlogged get_param calls in 4 user modules
(DOME2d_initialization, ISOMIP_initialization, Kelvin_initialization and
seamount_initialization) to help detect inconsistent units and scaling factors.
Also added comments describing many internal real variables and their units in
the DOME2d_initialization module and seamount_initialize_temperature_salinity.
All answers and output are bitwise identical.

src/user/DOME2d_initialization.F90

@@ -45,9 +45,13 @@ subroutine DOME2d_initialize_topography( D, G, param_file, max_depth )
   real,                    intent(in)  :: max_depth !< Maximum ocean depth in arbitrary units
 
   ! Local variables
+  real    :: bay_depth           ! Depth of shelf, as fraction of basin depth [nondim]
+  real    :: l1, l2              ! Fractional horizontal positions where the slope changes [nondim]
+  real    :: x                   ! Fractional horizontal positions [nondim]
+  real    :: dome2d_width_bay    ! Width of shelf, as fraction of domain [nondim]
+  real    :: dome2d_width_bottom ! Width of deep ocean basin, as fraction of domain [nondim]
+  real    :: dome2d_depth_bay    ! Depth of shelf, as fraction of basin depth [nondim]
   integer :: i, j
-  real    :: x, bay_depth, l1, l2
-  real    :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
   ! This include declares and sets the variable "version".
 # include "version_variable.h"
 
@@ -106,28 +110,30 @@ subroutine DOME2d_initialize_thickness ( h, depth_tot, G, GV, US, param_file, ju
   real :: e0(SZK_(GV))     ! The resting interface heights, in depth units [Z ~> m], usually
                            ! negative because it is positive upward.
   real :: eta1D(SZK_(GV)+1)! Interface height relative to the sea surface
-                           ! positive upward, in depth units [Z ~> m].
-  integer :: i, j, k, is, ie, js, je, nz
-  real    :: x
-  real    :: min_thickness
-  real    :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
+                           ! positive upward, in depth units [Z ~> m]
+  real    :: x                      ! Fractional horizontal positions [nondim]
+  real    :: min_thickness          ! Minimum layer thicknesses [Z ~> m]
+  real    :: dome2d_width_bay       ! Width of shelf, as fraction of domain [nondim]
+  real    :: dome2d_width_bottom    ! Width of deep ocean basin, as fraction of domain [nondim]
+  real    :: dome2d_depth_bay       ! Depth of shelf, as fraction of basin depth [nondim]
   character(len=40) :: verticalCoordinate
+  integer :: i, j, k, is, ie, js, je, nz
 
   is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke
 
   if (.not.just_read) &
     call MOM_mesg("MOM_initialization.F90, DOME2d_initialize_thickness: setting thickness")
 
-  call get_param(param_file, mdl,"MIN_THICKNESS",min_thickness, &
+  call get_param(param_file, mdl,"MIN_THICKNESS", min_thickness, &
                  default=1.e-3, units="m", do_not_log=.true., scale=US%m_to_Z)
   call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE", verticalCoordinate, &
                  default=DEFAULT_COORDINATE_MODE, do_not_log=.true.)
   call get_param(param_file, mdl, "DOME2D_SHELF_WIDTH", dome2d_width_bay, &
-                 default=0.1, do_not_log=.true.)
+                 units="nondim", default=0.1, do_not_log=.true.)
   call get_param(param_file, mdl, "DOME2D_BASIN_WIDTH", dome2d_width_bottom, &
-                 default=0.3, do_not_log=.true.)
+                 units="nondim", default=0.3, do_not_log=.true.)
   call get_param(param_file, mdl, "DOME2D_SHELF_DEPTH", dome2d_depth_bay, &
-                 default=0.2, do_not_log=.true.)
+                 units="nondim", default=0.2, do_not_log=.true.)
 
   if (just_read) return ! All run-time parameters have been read, so return.
 
@@ -229,15 +235,17 @@ subroutine DOME2d_initialize_temperature_salinity ( T, S, h, G, GV, US, param_fi
   logical,                                   intent(in)  :: just_read !< If true, this call will
                                                       !! only read parameters without changing T & S.
 
-  integer   :: i, j, k, is, ie, js, je, nz
-  real      :: x
-  integer   :: index_bay_z
-  real      :: delta_S
-  real      :: S_ref, T_ref         ! Reference salinity [S ~> ppt] and temperature [C ~> degC] within surface layer
-  real      :: S_range, T_range     ! Range of salinities [S ~> ppt] and temperatures [C ~> degC] over the vertical
-  real      :: xi0, xi1
+  real    :: x                  ! Fractional horizontal positions [nondim]
+  real    :: delta_S            ! Change in salinity between layers [S ~> ppt]
+  real    :: S_ref, T_ref       ! Reference salinity [S ~> ppt] and temperature [C ~> degC] within surface layer
+  real    :: S_range, T_range   ! Range of salinities [S ~> ppt] and temperatures [C ~> degC] over the vertical
+  real    :: xi0, xi1           ! Fractional vertical positions [nondim]
+  real    :: dome2d_width_bay   ! Width of shelf, as fraction of domain [nondim]
+  real    :: dome2d_width_bottom ! Width of deep ocean basin, as fraction of domain [nondim]
+  real    :: dome2d_depth_bay   ! Depth of shelf, as fraction of basin depth [nondim]
   character(len=40) :: verticalCoordinate
-  real    :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
+  integer :: index_bay_z
+  integer :: i, j, k, is, ie, js, je, nz
 
   is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke
 
@@ -370,10 +378,16 @@ subroutine DOME2d_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use_A
   real :: eta1D(SZK_(GV)+1)         ! Interface height relative to the sea surface
                                     ! positive upward [Z ~> m].
   real :: d_eta(SZK_(GV))           ! The layer thickness in a column [Z ~> m].
-  real :: dome2d_width_bay, dome2d_width_bottom, dome2d_depth_bay
+  real :: dome2d_width_bay          ! Width of shelf, as fraction of domain [nondim]
+  real :: dome2d_width_bottom       ! Width of deep ocean basin, as fraction of domain [nondim]
+  real :: dome2d_depth_bay          ! Depth of shelf, as fraction of basin depth [nondim]
   real :: dome2d_west_sponge_time_scale, dome2d_east_sponge_time_scale ! Sponge timescales [T ~> s]
-  real :: dome2d_west_sponge_width, dome2d_east_sponge_width
-  real :: dummy1, x, z
+  real :: dome2d_west_sponge_width  ! The fraction of the domain in which the western sponge for
+                                    ! restoring T/S is active [nondim]
+  real :: dome2d_east_sponge_width  ! The fraction of the domain in which the eastern sponge for
+                                    ! restoring T/S is active [nondim]
+  real :: dummy1, x                 ! Nondimensional local variables indicating horizontal positions [nondim]
+  real :: z                         ! Vertical positions [Z ~> m]
   integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
 
   is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke
@@ -405,15 +419,15 @@ subroutine DOME2d_initialize_sponges(G, GV, US, tv, depth_tot, param_file, use_A
      "DOME2d_initialize_sponges called with an associated ALE-sponge control structure.")
 
   call get_param(param_file, mdl, "DOME2D_SHELF_WIDTH", dome2d_width_bay, &
-                 default=0.1, do_not_log=.true.)
+                 units="nondim", default=0.1, do_not_log=.true.)
   call get_param(param_file, mdl, "DOME2D_BASIN_WIDTH", dome2d_width_bottom, &
                  default=0.3, do_not_log=.true.)
   call get_param(param_file, mdl, "DOME2D_SHELF_DEPTH", dome2d_depth_bay, &
                  default=0.2, do_not_log=.true.)
-  call get_param(param_file, mdl, "S_REF", S_ref, default=35.0, scale=US%ppt_to_S)
-  call get_param(param_file, mdl, "T_REF", T_ref, scale=US%degC_to_C, fail_if_missing=.false.)
-  call get_param(param_file, mdl, "S_RANGE", S_range, default=2.0, scale=US%ppt_to_S)
-  call get_param(param_file, mdl, "T_RANGE", T_range, default=0.0, scale=US%degC_to_C)
+  call get_param(param_file, mdl, "S_REF", S_ref, units="ppt", default=35.0, scale=US%ppt_to_S)
+  call get_param(param_file, mdl, "T_REF", T_ref, units="degC", scale=US%degC_to_C, fail_if_missing=.false.)
+  call get_param(param_file, mdl, "S_RANGE", S_range, units="ppt", default=2.0, scale=US%ppt_to_S)
+  call get_param(param_file, mdl, "T_RANGE", T_range, units="degC", default=0.0, scale=US%degC_to_C)
 
 
   ! Set the sponge damping rate as a function of position

src/user/ISOMIP_initialization.F90

@@ -343,8 +343,8 @@ subroutine ISOMIP_initialize_temperature_salinity ( T, S, h, depth_tot, G, GV, U
                   "A reference temperature used in initialization.", &
                   units="degC", scale=US%degC_to_C, fail_if_missing=.not.just_read, do_not_log=just_read)
       call get_param(param_file, mdl, "S_REF", S_Ref, &
-                  "A reference salinity used in initialization.", units="PSU", &
-                  default=35.0, scale=US%ppt_to_S, do_not_log=just_read)
+                  "A reference salinity used in initialization.", &
+                  units="PSU", default=35.0, scale=US%ppt_to_S, do_not_log=just_read)
       if (just_read) return ! All run-time parameters have been read, so return.
 
       ! write(mesg,*) 'read drho_dS, drho_dT', drho_dS1, drho_dT1
@@ -479,7 +479,7 @@ subroutine ISOMIP_initialize_sponges(G, GV, US, tv, depth_tot, PF, use_ALE, CSp,
             default=DEFAULT_COORDINATE_MODE)
 
   call get_param(PF, mdl, "ISOMIP_TNUDG", TNUDG, "Nudging time scale for sponge layers (days)", &
-                 default=0.0, scale=86400.0*US%s_to_T)
+                 units="days", default=0.0, scale=86400.0*US%s_to_T)
 
   call get_param(PF, mdl, "T_REF", t_ref, "Reference temperature", default=10.0, &
                  do_not_log=.true.)

src/user/Kelvin_initialization.F90

@@ -92,11 +92,11 @@ function register_Kelvin_OBC(param_file, CS, US, OBC_Reg)
   endif
   if (CS%mode /= 0) then
     call get_param(param_file, mdl, "DENSITY_RANGE", CS%rho_range, &
-                   default=2.0, do_not_log=.true., scale=US%kg_m3_to_R)
+                   units="kg m-3", default=2.0, do_not_log=.true., scale=US%kg_m3_to_R)
     call get_param(param_file, mdl, "RHO_0", CS%rho_0, &
-                   default=1035.0, do_not_log=.true., scale=US%kg_m3_to_R)
+                   units="kg m-3", default=1035.0, do_not_log=.true., scale=US%kg_m3_to_R)
     call get_param(param_file, mdl, "MAXIMUM_DEPTH", CS%H0, &
-                   default=1000.0, do_not_log=.true., scale=US%m_to_Z)
+                   units="m", default=1000.0, do_not_log=.true., scale=US%m_to_Z)
   endif
 
   ! Register the Kelvin open boundary.

src/user/seamount_initialization.F90

@@ -205,9 +205,19 @@ subroutine seamount_initialize_temperature_salinity(T, S, h, G, GV, US, param_fi
                                                       !! only read parameters without changing T & S.
 
   ! Local variables
+  real :: xi0, xi1  ! Fractional positions within the depth range [nondim]
+  real :: r         ! A nondimensional sharpness parameter with an exponetial profile [nondim]
+  real :: S_Ref     ! Default salinity range parameters [ppt].
+  real :: T_Ref     ! Default temperature range parameters [degC].
+  real :: S_Light, S_Dense, S_surf, S_range ! Salinity range parameters [S ~> ppt].
+  real :: T_Light, T_Dense, T_surf, T_range ! Temperature range parameters [C ~> degC].
+  real :: res_rat   ! The ratio of density space resolution in the denser part
+                    ! of the range to that in the lighter part of the range.
+                    ! Setting this greater than 1 increases the resolution for
+                    ! the denser water [nondim].
+  real :: a1, frac_dense, k_frac  ! Nondimensional temporary variables [nondim]
   integer :: i, j, k, is, ie, js, je, nz, k_light
-  real    :: xi0, xi1, r, S_surf, T_surf, S_range, T_range
-  real    :: T_ref, T_Light, T_Dense, S_ref, S_Light, S_Dense, a1, frac_dense, k_frac, res_rat
+
   character(len=20) :: verticalCoordinate, density_profile
 
   is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke
@@ -235,15 +245,16 @@ subroutine seamount_initialize_temperature_salinity(T, S, h, G, GV, US, param_fi
       ! These parameters are used in MOM_fixed_initialization.F90 when CONFIG_COORD="ts_range"
       call get_param(param_file, mdl, "T_REF", T_ref, default=10.0, do_not_log=.true.)
       call get_param(param_file, mdl, "TS_RANGE_T_LIGHT", T_light, &
-                 default=T_Ref, scale=US%degC_to_C, do_not_log=.true.)
+                 units="degC", default=T_Ref, scale=US%degC_to_C, do_not_log=.true.)
       call get_param(param_file, mdl, "TS_RANGE_T_DENSE", T_dense, &
-                 default=T_Ref, scale=US%degC_to_C, do_not_log=.true.)
+                 units="degC", default=T_Ref, scale=US%degC_to_C, do_not_log=.true.)
       call get_param(param_file, mdl, "S_REF", S_ref, default=35.0, do_not_log=.true.)
       call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", S_light, &
-                 default = S_Ref, scale=US%ppt_to_S, do_not_log=.true.)
+                 units="1e-3", default=S_Ref, scale=US%ppt_to_S, do_not_log=.true.)
       call get_param(param_file, mdl, "TS_RANGE_S_DENSE", S_dense, &
-                 default = S_Ref, scale=US%ppt_to_S, do_not_log=.true.)
-      call get_param(param_file, mdl, "TS_RANGE_RESOLN_RATIO", res_rat, default=1.0, do_not_log=.true.)
+                 units="1e-3", default=S_Ref, scale=US%ppt_to_S, do_not_log=.true.)
+      call get_param(param_file, mdl, "TS_RANGE_RESOLN_RATIO", res_rat, &
+                 units="nondim", default=1.0, do_not_log=.true.)
       if (just_read) return ! All run-time parameters have been read, so return.
 
       ! Emulate the T,S used in the "ts_range" coordinate configuration code