Merge branch 'dev/gfdl' into wave_speeds_cleanup

src/core/MOM.F90

@@ -1456,7 +1456,7 @@ subroutine step_MOM_thermo(CS, G, GV, US, u, v, h, tv, fluxes, dtdia, &
     if (CS%debug) then
       call MOM_thermo_chksum("Pre-oda ", tv, G, US, haloshift=0)
     endif
-    call apply_oda_tracer_increments(US%T_to_s*dtdia, Time_end_thermo, G, GV, tv, h, CS%odaCS)
+    call apply_oda_tracer_increments(dtdia, Time_end_thermo, G, GV, tv, h, CS%odaCS)
     if (CS%debug) then
       call MOM_thermo_chksum("Post-oda ", tv, G, US, haloshift=0)
     endif

src/core/MOM_density_integrals.F90

@@ -1550,10 +1550,10 @@ subroutine find_depth_of_pressure_in_cell(T_t, T_b, S_t, S_b, z_t, z_b, P_t, P_t
                                             !! are anomalous to [R ~> kg m-3]
   real,                  intent(in)  :: G_e !< Gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
   type(EOS_type),        intent(in)  :: EOS !< Equation of state structure
-  type(unit_scale_type), intent(in)  :: US !< A dimensional unit scaling type
+  type(unit_scale_type), intent(in)  :: US  !< A dimensional unit scaling type
   real,                  intent(out) :: P_b !< Pressure at the bottom of the cell [R L2 T-2 ~> Pa]
   real,                  intent(out) :: z_out !< Absolute depth at which anomalous pressure = p_tgt [Z ~> m]
-  real, optional,        intent(in)  :: z_tol !< The tolerance in finding z_out [Z ~> m]
+  real,                  intent(in)  :: z_tol !< The tolerance in finding z_out [Z ~> m]
 
   ! Local variables
   real :: dp    ! Pressure thickness of the layer [R L2 T-2 ~> Pa]
@@ -1583,8 +1583,7 @@ subroutine find_depth_of_pressure_in_cell(T_t, T_b, S_t, S_b, z_t, z_b, P_t, P_t
   Pa_left = P_t - P_tgt ! Pa_left < 0
   F_r = 1.
   Pa_right = P_b - P_tgt ! Pa_right > 0
-  Pa_tol = GxRho * 1.0e-5*US%m_to_Z
-  if (present(z_tol)) Pa_tol = GxRho * z_tol
+  Pa_tol = GxRho * z_tol
 
   F_guess = F_l - Pa_left / (Pa_right - Pa_left) * (F_r - F_l)
   Pa = Pa_right - Pa_left ! To get into iterative loop

src/core/MOM_dynamics_split_RK2.F90

@@ -1545,10 +1545,10 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param
   CS%id_PFv = register_diag_field('ocean_model', 'PFv', diag%axesCvL, Time, &
       'Meridional Pressure Force Acceleration', 'm s-2', conversion=US%L_T2_to_m_s2)
   CS%id_ueffA = register_diag_field('ocean_model', 'ueffA', diag%axesCuL, Time, &
-       'Effective U-Face Area', 'm^2', conversion = GV%H_to_m*US%L_to_m, &
+       'Effective U-Face Area', 'm^2', conversion=GV%H_to_m*US%L_to_m, &
        y_cell_method='sum', v_extensive=.true.)
   CS%id_veffA = register_diag_field('ocean_model', 'veffA', diag%axesCvL, Time, &
-       'Effective V-Face Area', 'm^2', conversion = GV%H_to_m*US%L_to_m, &
+       'Effective V-Face Area', 'm^2', conversion=GV%H_to_m*US%L_to_m, &
        x_cell_method='sum', v_extensive=.true.)
   if (GV%Boussinesq) then
     CS%id_deta_dt = register_diag_field('ocean_model', 'deta_dt', diag%axesT1, Time, &

src/core/MOM_dynamics_unsplit.F90

@@ -692,10 +692,10 @@ subroutine initialize_dyn_unsplit(u, v, h, Time, G, GV, US, param_file, diag, CS
   CS%id_PFv = register_diag_field('ocean_model', 'PFv', diag%axesCvL, Time, &
       'Meridional Pressure Force Acceleration', 'm s-2', conversion=US%L_T2_to_m_s2)
   CS%id_ueffA = register_diag_field('ocean_model', 'ueffA', diag%axesCuL, Time, &
-       'Effective U Face Area', 'm^2', conversion = GV%H_to_m*US%L_to_m, &
+       'Effective U Face Area', 'm^2', conversion=GV%H_to_m*US%L_to_m, &
        y_cell_method='sum', v_extensive=.true.)
   CS%id_veffA = register_diag_field('ocean_model', 'veffA', diag%axesCvL, Time, &
-       'Effective V Face Area', 'm^2', conversion = GV%H_to_m*US%L_to_m, &
+       'Effective V Face Area', 'm^2', conversion=GV%H_to_m*US%L_to_m, &
        x_cell_method='sum', v_extensive=.true.)
 
   id_clock_Cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=CLOCK_MODULE)

src/core/MOM_dynamics_unsplit_RK2.F90

@@ -655,10 +655,10 @@ subroutine initialize_dyn_unsplit_RK2(u, v, h, Time, G, GV, US, param_file, diag
   CS%id_PFv = register_diag_field('ocean_model', 'PFv', diag%axesCvL, Time, &
       'Meridional Pressure Force Acceleration', 'meter second-2', conversion=US%L_T2_to_m_s2)
   CS%id_ueffA = register_diag_field('ocean_model', 'ueffA', diag%axesCuL, Time, &
-       'Effective U-Face Area', 'm^2', conversion = GV%H_to_m*US%L_to_m, &
+       'Effective U-Face Area', 'm^2', conversion=GV%H_to_m*US%L_to_m, &
        y_cell_method='sum', v_extensive=.true.)
   CS%id_veffA = register_diag_field('ocean_model', 'veffA', diag%axesCvL, Time, &
-       'Effective V-Face Area', 'm^2', conversion = GV%H_to_m*US%L_to_m, &
+       'Effective V-Face Area', 'm^2', conversion=GV%H_to_m*US%L_to_m, &
        x_cell_method='sum', v_extensive=.true.)
 
   id_clock_Cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=CLOCK_MODULE)

src/core/MOM_variables.F90

@@ -168,6 +168,10 @@ module MOM_variables
     PFv => NULL(), &       !< Meridional acceleration due to pressure forces [L T-2 ~> m s-2]
     du_dt_visc => NULL(), &!< Zonal acceleration due to vertical viscosity [L T-2 ~> m s-2]
     dv_dt_visc => NULL(), &!< Meridional acceleration due to vertical viscosity [L T-2 ~> m s-2]
+    du_dt_visc_gl90 => NULL(), &!< Zonal acceleration due to GL90 vertical viscosity
+                           ! (is included in du_dt_visc) [L T-2 ~> m s-2]
+    dv_dt_visc_gl90 => NULL(), &!< Meridional acceleration due to GL90 vertical viscosity
+                           ! (is included in dv_dt_visc) [L T-2 ~> m s-2]
     du_dt_str => NULL(), & !< Zonal acceleration due to the surface stress (included
                            !! in du_dt_visc) [L T-2 ~> m s-2]
     dv_dt_str => NULL(), & !< Meridional acceleration due to the surface stress (included

src/core/MOM_verticalGrid.F90

@@ -105,7 +105,7 @@ subroutine verticalGridInit( param_file, GV, US )
                    log_to_all=.true., debugging=.true.)
   call get_param(param_file, mdl, "G_EARTH", GV%g_Earth, &
                  "The gravitational acceleration of the Earth.", &
-                 units="m s-2", default = 9.80, scale=US%Z_to_m*US%m_s_to_L_T**2)
+                 units="m s-2", default=9.80, scale=US%Z_to_m*US%m_s_to_L_T**2)
   call get_param(param_file, mdl, "RHO_0", GV%Rho0, &
                  "The mean ocean density used with BOUSSINESQ true to "//&
                  "calculate accelerations and the mass for conservation "//&

src/diagnostics/MOM_obsolete_params.F90

@@ -95,8 +95,9 @@ subroutine find_obsolete_params(param_file)
   call obsolete_real(param_file, "MIN_Z_DIAG_INTERVAL")
   call obsolete_char(param_file, "Z_OUTPUT_GRID_FILE")
 
-  ! This parameter is on the to-do list to be obsoleted.
-  ! call obsolete_logical(param_file, "NEW_SPONGES", hint="Use INTERPOLATE_SPONGE_TIME_SPACE instead.")
+  call read_param(param_file, "INTERPOLATE_SPONGE_TIME_SPACE", test_logic)
+  call obsolete_logical(param_file, "NEW_SPONGES", warning_val=test_logic, &
+                        hint="Use INTERPOLATE_SPONGE_TIME_SPACE instead.")
 
   ! Write the file version number to the model log.
   call log_version(param_file, mdl, version)

src/framework/MOM_domains.F90

@@ -180,7 +180,7 @@ subroutine MOM_domains_init(MOM_dom, param_file, symmetric, static_memory, &
   !$ if (.not.MOM_thread_affinity_set()) then
   !$   call get_param(param_file, mdl, "OCEAN_OMP_THREADS", ocean_nthreads, &
   !$              "The number of OpenMP threads that MOM6 will use.", &
-  !$              default = 1, layoutParam=.true.)
+  !$              default=1, layoutParam=.true.)
   !$   call get_param(param_file, mdl, "OCEAN_OMP_HYPER_THREAD", ocean_omp_hyper_thread, &
   !$              "If True, use hyper-threading.", default=.false., layoutParam=.true.)
   !$   call set_MOM_thread_affinity(ocean_nthreads, ocean_omp_hyper_thread)

src/ice_shelf/MOM_ice_shelf.F90

@@ -1468,7 +1468,7 @@ subroutine initialize_ice_shelf(param_file, ocn_grid, Time, CS, diag, forces_in,
 
   call get_param(param_file, mdl, "G_EARTH", CS%g_Earth, &
                  "The gravitational acceleration of the Earth.", &
-                 units="m s-2", default = 9.80, scale=US%m_s_to_L_T**2*US%Z_to_m)
+                 units="m s-2", default=9.80, scale=US%m_s_to_L_T**2*US%Z_to_m)
   call get_param(param_file, mdl, "C_P", CS%Cp, &
                  "The heat capacity of sea water, approximated as a constant. "//&
                  "The default value is from the TEOS-10 definition of conservative temperature.", &

src/ice_shelf/MOM_ice_shelf_dynamics.F90

@@ -407,7 +407,7 @@ subroutine initialize_ice_shelf_dyn(param_file, Time, ISS, CS, G, US, diag, new_
                  fail_if_missing=.true.)
     call get_param(param_file, mdl, "G_EARTH", CS%g_Earth, &
                  "The gravitational acceleration of the Earth.", &
-                 units="m s-2", default = 9.80, scale=US%m_s_to_L_T**2*US%Z_to_m)
+                 units="m s-2", default=9.80, scale=US%m_s_to_L_T**2*US%Z_to_m)
 
     call get_param(param_file, mdl, "GLEN_EXPONENT", CS%n_glen, &
                  "nonlinearity exponent in Glen's Law", &

src/initialization/MOM_fixed_initialization.F90

@@ -174,14 +174,13 @@ subroutine MOM_initialize_fixed(G, US, OBC, PF, write_geom, output_dir)
 
 end subroutine MOM_initialize_fixed
 
-!> MOM_initialize_topography makes the appropriate call to set up the bathymetry.  At this
-!! point the topography is in units of [Z ~> m] or [m], depending on the presence of US.
+!> MOM_initialize_topography makes the appropriate call to set up the bathymetry in units of [Z ~> m].
 subroutine MOM_initialize_topography(D, max_depth, G, PF, US)
   type(dyn_horgrid_type),           intent(in)  :: G  !< The dynamic horizontal grid type
   real, dimension(G%isd:G%ied,G%jsd:G%jed), &
-                                    intent(out) :: D  !< Ocean bottom depth [Z ~> m] or [m]
+                                    intent(out) :: D  !< Ocean bottom depth [Z ~> m]
   type(param_file_type),            intent(in)  :: PF !< Parameter file structure
-  real,                             intent(out) :: max_depth !< Maximum depth of model [Z ~> m] or [m]
+  real,                             intent(out) :: max_depth !< Maximum depth of model [Z ~> m]
   type(unit_scale_type),            intent(in)  :: US !< A dimensional unit scaling type
 
   ! This subroutine makes the appropriate call to set up the bottom depth.