Merge branch 'dev/master' of github.com:ESMG/MOM6 into dev/master

src/core/MOM_continuity_PPM.F90

@@ -1017,10 +1017,12 @@ subroutine set_zonal_BT_cont(u, h_in, hL, hR, BT_cont, uh_tot_0, duhdu_tot_0, &
 
   do k=1,nz ; do I=ish-1,ieh ; if (do_i(I)) then
     visc_rem_lim = max(visc_rem(I,k), min_visc_rem*visc_rem_max(I))
-    if (u(I,j,k) + duR(I)*visc_rem_lim > -du_CFL(I)*visc_rem(I,k)) &
-      duR(I) = -(u(I,j,k) + du_CFL(I)*visc_rem(I,k)) / visc_rem_lim
-    if (u(I,j,k) + duL(I)*visc_rem_lim < du_CFL(I)*visc_rem(I,k)) &
-      duL(I) = -(u(I,j,k) - du_CFL(I)*visc_rem(I,k)) / visc_rem_lim
+    if (visc_rem_lim > 0.0) then ! This is almost always true for ocean points.
+      if (u(I,j,k) + duR(I)*visc_rem_lim > -du_CFL(I)*visc_rem(I,k)) &
+        duR(I) = -(u(I,j,k) + du_CFL(I)*visc_rem(I,k)) / visc_rem_lim
+      if (u(I,j,k) + duL(I)*visc_rem_lim < du_CFL(I)*visc_rem(I,k)) &
+        duL(I) = -(u(I,j,k) - du_CFL(I)*visc_rem(I,k)) / visc_rem_lim
+    endif
   endif ; enddo ; enddo
 
   do k=1,nz
@@ -1770,10 +1772,12 @@ subroutine set_merid_BT_cont(v, h_in, hL, hR, BT_cont, vh_tot_0, dvhdv_tot_0, &
 
   do k=1,nz ; do i=ish,ieh ; if (do_i(i)) then
     visc_rem_lim = max(visc_rem(i,k), min_visc_rem*visc_rem_max(i))
-    if (v(i,J,k) + dvR(i)*visc_rem_lim > -dv_CFL(i)*visc_rem(i,k)) &
-      dvR(i) = -(v(i,J,k) + dv_CFL(i)*visc_rem(i,k)) / visc_rem_lim
-    if (v(i,J,k) + dvL(i)*visc_rem_lim < dv_CFL(i)*visc_rem(i,k)) &
-      dvL(i) = -(v(i,J,k) - dv_CFL(i)*visc_rem(i,k)) / visc_rem_lim
+    if (visc_rem_lim > 0.0) then ! This is almost always true for ocean points.
+      if (v(i,J,k) + dvR(i)*visc_rem_lim > -dv_CFL(i)*visc_rem(i,k)) &
+        dvR(i) = -(v(i,J,k) + dv_CFL(i)*visc_rem(i,k)) / visc_rem_lim
+      if (v(i,J,k) + dvL(i)*visc_rem_lim < dv_CFL(i)*visc_rem(i,k)) &
+        dvL(i) = -(v(i,J,k) - dv_CFL(i)*visc_rem(i,k)) / visc_rem_lim
+    endif
   endif ; enddo ; enddo
   do k=1,nz
     do i=ish,ieh ; if (do_i(i)) then

src/diagnostics/MOM_diag_to_Z.F90

@@ -850,6 +850,8 @@ subroutine register_Z_tracer(tr_ptr, name, long_name, units, Time, G, CS, standa
   character(len=256) :: posted_cmor_standard_name
   character(len=256) :: posted_cmor_long_name
 
+  if (CS%nk_zspace<1) return
+
   if (present(standard_name)) then
     posted_standard_name = standard_name
   else
@@ -858,7 +860,7 @@ subroutine register_Z_tracer(tr_ptr, name, long_name, units, Time, G, CS, standa
 
   call register_Z_tracer_low(tr_ptr, name, long_name, units, trim(posted_standard_name), Time, G, CS)
 
- if (present(cmor_field_name)) then
+  if (present(cmor_field_name)) then
     ! Fallback values for strings set to "NULL"
     posted_cmor_units         = "not provided"   !
     posted_cmor_standard_name = "not provided"   ! values might be replaced with a CS%missing field?
@@ -966,9 +968,10 @@ subroutine MOM_diag_to_Z_init(Time, G, GV, param_file, diag, CS)
 
   character(len=40)  :: mod = "MOM_diag_to_Z" ! module name
   character(len=200) :: in_dir, zgrid_file    ! strings for directory/file
-  character(len=48)  :: flux_units
+  character(len=48)  :: flux_units, string
   integer :: z_axis, zint_axis
   integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, nk, id_test
+  logical :: diag_mediator_is_using_z
   isd  = G%isd   ; ied = G%ied  ; jsd  = G%jsd  ; jed  = G%jed ; nk = G%ke
   IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB
 
@@ -991,7 +994,19 @@ subroutine MOM_diag_to_Z_init(Time, G, GV, param_file, diag, CS)
                  "The file that specifies the vertical grid for \n"//&
                  "depth-space diagnostics, or blank to disable \n"//&
                  "depth-space output.", default="")
+
+  ! Check that the diag_mediator z-sapce remapping is not using the same module name
+  string = ''
+  call get_param(param_file, mod, "DIAG_REMAP_Z_MODULE_SUFFIX", string, &
+                 default='_z_new', do_not_log=.true.)
+  diag_mediator_is_using_z = .false.
+  if (trim(string) == '_z') diag_mediator_is_using_z = .true.
+
   if (len_trim(zgrid_file) > 0) then
+    if (diag_mediator_is_using_z) call MOM_error(FATAL, "MOM_diag_to_Z_init:"// &
+           "Z_OUTPUT_GRID_FILE can not be used when DIAG_REMAP_Z_MODULE_SUFFIX='_z'." // &
+           " Z_OUTPUT_GRID_FILE='"//trim(zgrid_file)//"'")
+
     call get_param(param_file, mod, "INPUTDIR", in_dir, &
                  "The directory in which input files are found.", default=".")
     in_dir = slasher(in_dir)  
@@ -1004,7 +1019,7 @@ subroutine MOM_diag_to_Z_init(Time, G, GV, param_file, diag, CS)
                  "from the size of the variable zw in the output grid file.", &
                  units="nondim")
   else
-    in_dir = "" ; CS%nk_zspace = -1
+    CS%nk_zspace = -1
   endif
 
   if (CS%nk_zspace > 0) then
@@ -1041,7 +1056,7 @@ subroutine MOM_diag_to_Z_init(Time, G, GV, param_file, diag, CS)
         missing_value=CS%missing_trans)
     if (CS%id_vh_z>0) call safe_alloc_ptr(CS%vh_z,isd,ied,JsdB,JedB,CS%nk_zspace)
 
-  else
+  elseif (.not. diag_mediator_is_using_z) then
 
     ! Check whether diag-table is requesting any z-space files; issue a warning if it is.
 
@@ -1229,6 +1244,7 @@ function ocean_register_diag_with_z(tr_ptr, vardesc_tr, G, Time, CS)
   isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed ; nk = G%ke
   if (.not.associated(CS)) call MOM_error(FATAL, &
          "register_Z_tracer: Module must be initialized before it is used.")
+  if (CS%nk_zspace<1) return
 
   if (CS%num_tr_used >= MAX_FIELDS_) then
     call MOM_error(WARNING,"ocean_register_diag_with_z:  Attempted to register and use "//&

src/framework/MOM_diag_mediator.F90

@@ -28,6 +28,8 @@ module MOM_diag_mediator
 !********+*********+*********+*********+*********+*********+*********+**
 
 use MOM_coms,             only : PE_here
+use MOM_cpu_clock,        only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
+use MOM_cpu_clock,        only : CLOCK_MODULE, CLOCK_ROUTINE
 use MOM_error_handler,    only : MOM_error, FATAL, is_root_pe
 use MOM_file_parser,      only : get_param, log_param, log_version, param_file_type
 use MOM_grid,             only : ocean_grid_type
@@ -162,13 +164,16 @@ module MOM_diag_mediator
   ! Number of z levels used for remapping
   integer :: nz_remap
 
-  ! Define z star on u, v, T grids, these are the interface positions
+  ! Output grid thicknesses
   real, dimension(:,:,:), allocatable :: h_zoutput
 
   ! Keep track of which remapping is needed for diagnostic output
   logical :: do_z_remapping_on_u, do_z_remapping_on_v, do_z_remapping_on_T
   logical :: remapping_initialized
 
+  !> String appended to module name for z*-remapped diagnostics
+  character(len=6) :: z_remap_suffix = '_z_new'
+
   ! Pointer to H and G for remapping
   real, dimension(:,:,:), pointer :: h => null()
   type(ocean_grid_type), pointer :: G => null()
@@ -181,6 +186,9 @@ module MOM_diag_mediator
 
 end type diag_ctrl
 
+! CPU clocks
+integer :: id_clock_diag_mediator, id_clock_diag_z_remap, id_clock_diag_grid_updates
+
 contains
 
 !> Sets up diagnostics axes
@@ -312,6 +320,17 @@ subroutine set_axes_info(G, GV, param_file, diag_cs, set_vertical)
                             "Depth of cell center", direction=-1)
     id_zzi = diag_axis_init('zi_remap', diag_cs%zi_remap, "meters", "z", &
                             'Depth of interfaces', direction=-1)
+    call get_param(param_file, mod, "DIAG_REMAP_Z_MODULE_SUFFIX", diag_cs%z_remap_suffix, &
+                 'This is the string attached to the end of "ocean_model"\n'// &
+                 'for use in the model column of the diag_table to indicate\n'// &
+                 'a diagnostic should be remapped to z*-coordinates.', &
+                 default='_z_new')
+    if (trim(diag_cs%z_remap_suffix) == '_z') then
+      ! This will conflict with the older MOM_diag_to_Z module for z-output
+      call get_param(param_file, mod, "Z_OUTPUT_GRID_FILE", string, default="", do_not_log=.true.)
+      if (len(trim(string))>0) call MOM_error(FATAL,"MOM_diag_mediator, set_axes_info: "// &
+           "Z_OUTPUT_GRID_FILE must be blank to use DIAG_REMAP_Z_MODULE_SUFFIX='_z'")
+    endif
   else
     ! In this case the axes associated with these will never be used, however
     ! they need to be positive otherwise FMS complains.
@@ -471,6 +490,7 @@ subroutine post_data_0d(diag_field_id, field, diag_cs, is_static, mask)
   logical :: used, is_stat
   type(diag_type), pointer :: diag => null()
 
+  if (id_clock_diag_mediator>0) call cpu_clock_begin(id_clock_diag_mediator)
   is_stat = .false. ; if (present(is_static)) is_stat = is_static
 
   ! Iterate over list of diag 'variants', e.g. CMOR aliases, call send_data
@@ -487,6 +507,7 @@ subroutine post_data_0d(diag_field_id, field, diag_cs, is_static, mask)
     diag => diag%next
   enddo
 
+  if (id_clock_diag_mediator>0) call cpu_clock_end(id_clock_diag_mediator)
 end subroutine post_data_0d
 
 subroutine post_data_1d_k(diag_field_id, field, diag_cs, is_static)
@@ -507,6 +528,7 @@ subroutine post_data_1d_k(diag_field_id, field, diag_cs, is_static)
   integer :: isv, iev, jsv, jev
   type(diag_type), pointer :: diag => null()
 
+  if (id_clock_diag_mediator>0) call cpu_clock_begin(id_clock_diag_mediator)
   is_stat = .false. ; if (present(is_static)) is_stat = is_static
 
   ! Iterate over list of diag 'variants', e.g. CMOR aliases.
@@ -522,6 +544,7 @@ subroutine post_data_1d_k(diag_field_id, field, diag_cs, is_static)
     diag => diag%next
   enddo
 
+  if (id_clock_diag_mediator>0) call cpu_clock_end(id_clock_diag_mediator)
 end subroutine post_data_1d_k
 
 subroutine post_data_2d(diag_field_id, field, diag_cs, is_static, mask)
@@ -541,6 +564,7 @@ subroutine post_data_2d(diag_field_id, field, diag_cs, is_static, mask)
 
   type(diag_type), pointer :: diag => null()
 
+  if (id_clock_diag_mediator>0) call cpu_clock_begin(id_clock_diag_mediator)
   ! Iterate over list of diag 'variants' (e.g. CMOR aliases) and post each.
   call assert(diag_field_id < diag_cs%next_free_diag_id, &
               'post_data_2d: Unregistered diagnostic id')
@@ -550,6 +574,7 @@ subroutine post_data_2d(diag_field_id, field, diag_cs, is_static, mask)
     diag => diag%next
   enddo
 
+  if (id_clock_diag_mediator>0) call cpu_clock_end(id_clock_diag_mediator)
 end subroutine post_data_2d
 
 subroutine post_data_2d_low(diag, field, diag_cs, is_static, mask)
@@ -657,6 +682,7 @@ subroutine post_data_3d(diag_field_id, field, diag_cs, is_static, mask)
   type(diag_type), pointer :: diag => null()
   real, allocatable :: remapped_field(:,:,:)
 
+  if (id_clock_diag_mediator>0) call cpu_clock_begin(id_clock_diag_mediator)
   ! Iterate over list of diag 'variants', e.g. CMOR aliases, different vertical
   ! grids, and post each.
   call assert(diag_field_id < diag_cs%next_free_diag_id, &
@@ -671,8 +697,10 @@ subroutine post_data_3d(diag_field_id, field, diag_cs, is_static, mask)
         call MOM_error(FATAL,"post_data_3d: no mask for regridded field.")
       endif
 
+      if (id_clock_diag_z_remap>0) call cpu_clock_begin(id_clock_diag_z_remap)
       allocate(remapped_field(DIM_I(field),DIM_J(field), diag_cs%nz_remap))
       call remap_diag_to_z(field, diag, diag_cs, remapped_field)
+      if (id_clock_diag_z_remap>0) call cpu_clock_end(id_clock_diag_z_remap)
       if (associated(diag%mask3d)) then
         ! Since 3d masks do not vary in the vertical, just use as much as is
         ! needed.
@@ -681,12 +709,15 @@ subroutine post_data_3d(diag_field_id, field, diag_cs, is_static, mask)
       else
         call post_data_3d_low(diag, remapped_field, diag_cs, is_static)
       endif
+      if (id_clock_diag_z_remap>0) call cpu_clock_begin(id_clock_diag_z_remap)
       deallocate(remapped_field)
+      if (id_clock_diag_z_remap>0) call cpu_clock_end(id_clock_diag_z_remap)
     else
       call post_data_3d_low(diag, field, diag_cs, is_static, mask)
     endif
     diag => diag%next
   enddo
+  if (id_clock_diag_mediator>0) call cpu_clock_end(id_clock_diag_mediator)
 
 end subroutine post_data_3d
 
@@ -817,6 +848,7 @@ subroutine diag_update_target_grids(diag_cs)
   if (.not. allocated(diag_cs%zi_remap)) then
     return
   endif
+  if (id_clock_diag_grid_updates>0) call cpu_clock_begin(id_clock_diag_grid_updates)
 
   if (.not. diag_cs%remapping_initialized) then
     call assert(allocated(diag_cs%zi_remap), &
@@ -854,6 +886,7 @@ subroutine diag_update_target_grids(diag_cs)
   ! when doing remapping.
   diag_cs%h_old(:,:,:) = diag_cs%h(:,:,:)
 #endif
+  if (id_clock_diag_grid_updates>0) call cpu_clock_end(id_clock_diag_grid_updates)
 
 end subroutine diag_update_target_grids
 
@@ -1127,7 +1160,7 @@ function register_diag_field(module_name, field_name, axes, init_time,         &
   ! Remap to z vertical coordinate, note that only diagnostics on layers
   ! (not interfaces) are supported, also B axes are not supported yet
   if (is_layer_axes(axes, diag_cs) .and. (.not. is_B_axes(axes, diag_cs)) .and. axes%rank == 3) then
-    if (get_diag_field_id_fms(trim(module_name)//'_z_new', field_name) /= DIAG_FIELD_NOT_FOUND) then
+    if (get_diag_field_id_fms(trim(module_name)//trim(diag_cs%z_remap_suffix), field_name) /= DIAG_FIELD_NOT_FOUND) then
       if (.not. allocated(diag_cs%zi_remap)) then
         call MOM_error(FATAL, 'register_diag_field: Request to regrid but no '// &
                        'destination grid spec provided, see param DIAG_REMAP_Z_GRID_DEF')
@@ -1139,7 +1172,7 @@ function register_diag_field(module_name, field_name, axes, init_time,         &
       call set_diag_mask(z_remap_diag, diag_cs, axes)
       call set_diag_remap_axes(z_remap_diag, diag_cs, axes)
       call assert(associated(z_remap_diag%remap_axes), 'register_diag_field: remap axes not set')
-      fms_id = register_diag_field_fms(module_name//'_z_new', field_name, &
+      fms_id = register_diag_field_fms(module_name//trim(diag_cs%z_remap_suffix), field_name, &
            z_remap_diag%remap_axes%handles, &
            init_time, long_name=long_name, units=units, missing_value=MOM_missing_value, &
            range=range, mask_variant=mask_variant, standard_name=standard_name, &
@@ -1160,13 +1193,13 @@ function register_diag_field(module_name, field_name, axes, init_time,         &
     if (is_root_pe() .and. diag_CS%doc_unit > 0) then
       msg = ''
       if (present(cmor_field_name)) msg = 'CMOR equivalent is "'//trim(cmor_field_name)//'"'
-      call log_available_diag(associated(z_remap_diag), module_name//'_z_new', field_name, &
+      call log_available_diag(associated(z_remap_diag), module_name//trim(diag_cs%z_remap_suffix), field_name, &
                               cm_string, msg, diag_CS, long_name, units, standard_name)
     endif
 
     ! Remap to z vertical coordinate with CMOR names and attributes
     if (present(cmor_field_name)) then
-      if (get_diag_field_id_fms(module_name//'_z_new', cmor_field_name) /= DIAG_FIELD_NOT_FOUND) then
+      if (get_diag_field_id_fms(module_name//trim(diag_cs%z_remap_suffix), cmor_field_name) /= DIAG_FIELD_NOT_FOUND) then
         if (.not. allocated(diag_cs%zi_remap)) then
           call MOM_error(FATAL, 'register_diag_field: Request to regrid but no '// &
                          'destination grid spec provided, see param DIAG_REMAP_Z_GRID_DEF')
@@ -1178,7 +1211,7 @@ function register_diag_field(module_name, field_name, axes, init_time,         &
         call set_diag_mask(cmor_z_remap_diag, diag_cs, axes)
         call set_diag_remap_axes(cmor_z_remap_diag, diag_cs, axes)
         call assert(associated(cmor_z_remap_diag%remap_axes), 'register_diag_field: remap axes not set')
-        fms_id = register_diag_field_fms(module_name//'_z_new', cmor_field_name, &
+        fms_id = register_diag_field_fms(module_name//trim(diag_cs%z_remap_suffix), cmor_field_name, &
              cmor_z_remap_diag%remap_axes%handles, &
              init_time, long_name=trim(posted_cmor_long_name), units=trim(posted_cmor_units), missing_value=MOM_missing_value, &
              range=range, mask_variant=mask_variant, standard_name=trim(posted_cmor_standard_name), &
@@ -1198,7 +1231,7 @@ function register_diag_field(module_name, field_name, axes, init_time,         &
       endif
       if (is_root_pe() .and. diag_CS%doc_unit > 0) then
         msg = 'native name is "'//trim(field_name)//'"'
-        call log_available_diag(associated(cmor_z_remap_diag), module_name//'_z_new', cmor_field_name, &
+        call log_available_diag(associated(cmor_z_remap_diag), module_name//trim(diag_cs%z_remap_suffix), cmor_field_name, &
                                 cm_string, msg, diag_CS, posted_cmor_long_name, posted_cmor_units, &
                                 posted_cmor_standard_name)
       endif
@@ -1629,6 +1662,10 @@ subroutine diag_mediator_init(G, param_file, diag_cs, err_msg, doc_file_dir)
   character(len=240) :: doc_file, doc_file_dflt, doc_path
   character(len=40)  :: mod  = "MOM_diag_mediator" ! This module's name.
 
+  id_clock_diag_mediator = cpu_clock_id('(Ocean diagnostics framework)', grain=CLOCK_MODULE)
+  id_clock_diag_z_remap = cpu_clock_id('(Ocean diagnostics remapping)', grain=CLOCK_ROUTINE)
+  id_clock_diag_grid_updates = cpu_clock_id('(Ocean diagnostics grid updates)', grain=CLOCK_ROUTINE)
+
   call diag_manager_init(err_msg=err_msg)
 
   ! Allocate and initialise list of all diagnostics (and variants)