Merge branch 'Hallberg-NOAA-revise_Kd_extra' into dev/gfdl

.testing/Makefile

@@ -445,13 +445,16 @@ $(foreach d,$(DIMS),$(eval $(call CMP_RULE,dim.$(d),symmetric dim.$(d))))
 
 # Regression testing only checks for changes in existing diagnostics
 %.regression.diag: $(foreach b,symmetric target,work/%/$(b)/chksum_diag)
-	@! diff $^ | grep "^[<>]" | grep "^>" \
+	@! diff $^ | grep "^[<>]" | grep "^>" > /dev/null \
 	  || ! (\
 	    mkdir -p results/$*; \
 	    (diff $^ | tee results/$*/chksum_diag.regression.diff | head -n 20) ; \
 	    echo -e "$(FAIL): Diagnostics $*.regression.diag have changed." \
 	  )
-	@diff $^ || echo -e "$(WARN): New diagnostics in $<"
+	@cmp $^ || ( \
+	  diff $^ | head -n 20; \
+	  echo -e "$(WARN): New diagnostics in $<" \
+	)
 	@echo -e "$(PASS): Diagnostics $*.regression.diag agree."
 
 
@@ -480,7 +483,7 @@ work/%/$(1)/ocean.stats work/%/$(1)/chksum_diag: build/$(2)/MOM6 $(VENV_PATH)
 	fi
 	mkdir -p $$(@D)/RESTART
 	echo -e "$(4)" > $$(@D)/MOM_override
-	rm -rf results/$$*/std.$(1).{out,err}
+	rm -f results/$$*/std.$(1).{out,err}
 	cd $$(@D) \
 	  && $(TIME) $(5) $(MPIRUN) -n $(6) ../../../$$< 2> std.err > std.out \
 	  || !( \
@@ -543,7 +546,7 @@ work/%/restart/ocean.stats: build/symmetric/MOM6 $(VENV_PATH)
 	  && halfperiod=$$(printf "%.f" $$(bc <<< "scale=10; 0.5 * $${daymax} * $${timeunit_int}")) \
 	  && printf "\n&ocean_solo_nml\n    seconds = $${halfperiod}\n/\n" >> input.nml
 	# Remove any previous archived output
-	rm -rf results/$*/std.restart{1,2}.{out,err}
+	rm -f results/$*/std.restart{1,2}.{out,err}
 	# Run the first half-period
 	cd $(@D) && $(TIME) $(MPIRUN) -n 1 ../../../$< 2> std1.err > std1.out \
 	  || !( \

src/core/MOM_variables.F90

@@ -241,15 +241,6 @@ module MOM_variables
   real, pointer, dimension(:,:,:) :: &
     Ray_u => NULL(), & !< The Rayleigh drag velocity to be applied to each layer at u-points [Z T-1 ~> m s-1].
     Ray_v => NULL()    !< The Rayleigh drag velocity to be applied to each layer at v-points [Z T-1 ~> m s-1].
-  real, pointer, dimension(:,:,:) :: Kd_extra_T => NULL()
-                !< The extra diffusivity of temperature due to double diffusion relative to the
-                !! diffusivity of density [Z2 T-1 ~> m2 s-1].
-  real, pointer, dimension(:,:,:) :: Kd_extra_S => NULL()
-                !< The extra diffusivity of salinity due to double diffusion relative to the
-                !! diffusivity of density [Z2 T-1 ~> m2 s-1].
-  ! One of Kd_extra_T and Kd_extra_S is always 0. Kd_extra_S is positive for salt fingering;
-  ! Kd_extra_T is positive for double diffusive convection.  They are only allocated if
-  ! DOUBLE_DIFFUSION is true.
   real, pointer, dimension(:,:,:) :: Kd_shear => NULL()
                 !< The shear-driven turbulent diapycnal diffusivity at the interfaces between layers
                 !! in tracer columns [Z2 T-1 ~> m2 s-1].

src/parameterizations/vertical/MOM_bkgnd_mixing.F90

@@ -237,16 +237,13 @@ subroutine bkgnd_mixing_init(Time, G, GV, US, param_file, diag, CS)
                  units="nondim", default=1.0)
 
   if (CS%Bryan_Lewis_diffusivity .or. CS%horiz_varying_background) then
-
     prandtl_bkgnd_comp = CS%prandtl_bkgnd
-    if (CS%Kd /= 0.0) prandtl_bkgnd_comp = Kv/CS%Kd
+    if (CS%Kd /= 0.0) prandtl_bkgnd_comp = Kv / CS%Kd
 
     if ( abs(CS%prandtl_bkgnd - prandtl_bkgnd_comp)>1.e-14) then
-      call MOM_error(FATAL,"set_diffusivity_init: The provided KD, KV,"//&
-                           "and PRANDTL_BKGND values are incompatible. The following "//&
-                           "must hold: KD*PRANDTL_BKGND==KV")
+      call MOM_error(FATAL, "bkgnd_mixing_init: The provided KD, KV and PRANDTL_BKGND values "//&
+                            "are incompatible. The following must hold: KD*PRANDTL_BKGND==KV")
     endif
-
   endif
 
   call get_param(param_file, mdl, "HENYEY_IGW_BACKGROUND", CS%Henyey_IGW_background, &
@@ -264,7 +261,7 @@ subroutine bkgnd_mixing_init(Time, G, GV, US, param_file, diag, CS)
 
   if (CS%Kd>0.0 .and. (trim(CS%bkgnd_scheme_str)=="BRYAN_LEWIS_DIFFUSIVITY" .or.&
                           trim(CS%bkgnd_scheme_str)=="HORIZ_VARYING_BACKGROUND" )) then
-    call MOM_error(WARNING, "set_diffusivity_init: a nonzero constant background "//&
+    call MOM_error(WARNING, "bkgnd_mixing_init: a nonzero constant background "//&
          "diffusivity (KD) is specified along with "//trim(CS%bkgnd_scheme_str))
   endif
 
@@ -540,7 +537,7 @@ subroutine check_bkgnd_scheme(CS, str)
   if (trim(CS%bkgnd_scheme_str)=="none") then
     CS%bkgnd_scheme_str = str
   else
-    call MOM_error(FATAL, "set_diffusivity_init: Cannot activate both "//trim(str)//" and "//&
+    call MOM_error(FATAL, "bkgnd_mixing_init: Cannot activate both "//trim(str)//" and "//&
                    trim(CS%bkgnd_scheme_str)//".")
   endif
 

src/parameterizations/vertical/MOM_diabatic_aux.F90

@@ -20,7 +20,7 @@ module MOM_diabatic_aux
 use MOM_opacity,       only : optics_type, optics_nbands, absorbRemainingSW, sumSWoverBands
 use MOM_tracer_flow_control, only : get_chl_from_model, tracer_flow_control_CS
 use MOM_unit_scaling,  only : unit_scale_type
-use MOM_variables,     only : thermo_var_ptrs, vertvisc_type! , accel_diag_ptrs
+use MOM_variables,     only : thermo_var_ptrs ! , vertvisc_type, accel_diag_ptrs
 use MOM_verticalGrid,  only : verticalGrid_type
 use time_interp_external_mod, only : init_external_field, time_interp_external
 use time_interp_external_mod, only : time_interp_external_init
@@ -223,15 +223,23 @@ end subroutine make_frazil
 
 !> This subroutine applies double diffusion to T & S, assuming no diapycal mass
 !! fluxes, using a simple triadiagonal solver.
-subroutine differential_diffuse_T_S(h, tv, visc, dt, G, GV)
+subroutine differential_diffuse_T_S(h, T, S, Kd_T, Kd_S, dt, G, GV)
   type(ocean_grid_type),   intent(in)    :: G    !< The ocean's grid structure
   type(verticalGrid_type), intent(in)    :: GV   !< The ocean's vertical grid structure
   real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
                            intent(in)    :: h    !< Layer thicknesses [H ~> m or kg m-2]
-  type(thermo_var_ptrs),   intent(inout) :: tv   !< Structure containing pointers to any
-                                                 !! available thermodynamic fields.
-  type(vertvisc_type),     intent(in)    :: visc !< Structure containing vertical viscosities, bottom
-                                                 !! boundary layer properies, and related fields.
+  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
+                           intent(inout) :: T    !< Potential temperature [degC].
+  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
+                           intent(inout) :: S    !< Salinity [PSU] or [gSalt/kg], generically [ppt].
+  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), &
+                           intent(inout)    :: Kd_T !< The extra diffusivity of temperature due to
+                                                 !! double diffusion relative to the diffusivity of
+                                                 !! diffusivity of density [Z2 T-1 ~> m2 s-1].
+  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), &
+                           intent(in)    :: Kd_S !< The extra diffusivity of salinity due to
+                                                 !! double diffusion relative to the diffusivity of
+                                                 !! diffusivity of density [Z2 T-1 ~> m2 s-1].
   real,                    intent(in)    :: dt   !<  Time increment [T ~> s].
 
   ! local variables
@@ -250,27 +258,12 @@ subroutine differential_diffuse_T_S(h, tv, visc, dt, G, GV)
                        ! interface [H-1 ~> m-1 or m2 kg-1].
   real :: b_denom_T    ! The first term in the denominators for the expressions
   real :: b_denom_S    ! for b1_T and b1_S, both [H ~> m or kg m-2].
-  real, dimension(:,:,:), pointer :: T=>NULL(), S=>NULL()
-  real, dimension(:,:,:), pointer :: Kd_T=>NULL(), Kd_S=>NULL() ! Diffusivities [Z2 T-1 ~> m2 s-1].
   integer :: i, j, k, is, ie, js, je, nz
 
   is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = G%ke
   h_neglect = GV%H_subroundoff
 
-  if (.not.associated(tv%T)) call MOM_error(FATAL, &
-      "differential_diffuse_T_S: Called with an unassociated tv%T")
-  if (.not.associated(tv%S)) call MOM_error(FATAL, &
-      "differential_diffuse_T_S: Called with an unassociated tv%S")
-  if (.not.associated(visc%Kd_extra_T)) call MOM_error(FATAL, &
-      "differential_diffuse_T_S: Called with an unassociated visc%Kd_extra_T")
-  if (.not.associated(visc%Kd_extra_S)) call MOM_error(FATAL, &
-      "differential_diffuse_T_S: Called with an unassociated visc%Kd_extra_S")
-
-  T => tv%T ; S => tv%S
-  Kd_T => visc%Kd_extra_T ; Kd_S => visc%Kd_extra_S
-!$OMP parallel do default(none) shared(is,ie,js,je,h,h_neglect,dt,Kd_T,Kd_S,G,GV,T,S,nz) &
-!$OMP                          private(I_h_int,mix_T,mix_S,h_tr,b1_T,b1_S, &
-!$OMP                                  d1_T,d1_S,c1_T,c1_S,b_denom_T,b_denom_S)
+  !$OMP parallel do default(private) shared(is,ie,js,je,h,h_neglect,dt,Kd_T,Kd_S,G,GV,T,S,nz)
   do j=js,je
     do i=is,ie
       I_h_int = 1.0 / (0.5 * (h(i,j,1) + h(i,j,2)) + h_neglect)