Conditional thermodynamic diags; Flux diag bugfix

Several thermodynamic diabatic diagnostics in the legacy_diabatic update
could be registered and computed when the thermodynamics was disabled,
so these are now only conditionally checked and registered when
`ENABLE_THERMODYNAMICS` is true.

We also fixed a bug in the flux diagnostic registration by splitting the
logical association and memory access into two separate stages, since it
was raising errors in builds with stricter memory access.

src/core/MOM_forcing_type.F90

@@ -2239,14 +2239,19 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles)
     if (handles%id_net_massout > 0 .or. handles%id_total_net_massout > 0) then
       do j=js,je ; do i=is,ie
         res(i,j) = 0.0
-        if (associated(fluxes%lprec) .and. fluxes%lprec(i,j) < 0.0) &
-          res(i,j) = res(i,j) + fluxes%lprec(i,j)
-        if (associated(fluxes%vprec) .and. fluxes%vprec(i,j) < 0.0) &
-          res(i,j) = res(i,j) + fluxes%vprec(i,j)
-        if (associated(fluxes%evap) .and. fluxes%evap(i,j) < 0.0) &
-          res(i,j) = res(i,j) + fluxes%evap(i,j)
-        if (associated(fluxes%seaice_melt) .and. fluxes%seaice_melt(i,j) < 0.0) &
-          res(i,j) = res(i,j) + fluxes%seaice_melt(i,j)
+        if (associated(fluxes%lprec)) then
+          if (fluxes%lprec(i,j) < 0.0) res(i,j) = res(i,j) + fluxes%lprec(i,j)
+        endif
+        if (associated(fluxes%vprec)) then
+          if (fluxes%vprec(i,j) < 0.0) res(i,j) = res(i,j) + fluxes%vprec(i,j)
+        endif
+        if (associated(fluxes%evap)) then
+          if (fluxes%evap(i,j) < 0.0) res(i,j) = res(i,j) + fluxes%evap(i,j)
+        endif
+        if (associated(fluxes%seaice_melt)) then
+          if (fluxes%seaice_melt(i,j) < 0.0) &
+            res(i,j) = res(i,j) + fluxes%seaice_melt(i,j)
+        endif
       enddo ; enddo
       if (handles%id_net_massout > 0) call post_data(handles%id_net_massout, res, diag)
       if (handles%id_total_net_massout > 0) then
@@ -2269,15 +2274,20 @@ subroutine forcing_diagnostics(fluxes, sfc_state, dt, G, diag, handles)
         if (associated(fluxes%frunoff)) &
           res(i,j) = res(i,j) + fluxes%frunoff(i,j)
 
-        if (associated(fluxes%lprec) .and. fluxes%lprec(i,j) > 0.0) &
-          res(i,j) = res(i,j) + fluxes%lprec(i,j)
-        if (associated(fluxes%vprec) .and. fluxes%vprec(i,j) > 0.0) &
-          res(i,j) = res(i,j) + fluxes%vprec(i,j)
+        if (associated(fluxes%lprec)) then
+          if (fluxes%lprec(i,j) > 0.0) res(i,j) = res(i,j) + fluxes%lprec(i,j)
+        endif
+        if (associated(fluxes%vprec)) then
+          if (fluxes%vprec(i,j) > 0.0) res(i,j) = res(i,j) + fluxes%vprec(i,j)
+        endif
         ! fluxes%cond is not needed because it is derived from %evap > 0
-        if (associated(fluxes%evap) .and. fluxes%evap(i,j) > 0.0) &
-          res(i,j) = res(i,j) + fluxes%evap(i,j)
-        if (associated(fluxes%seaice_melt) .and. fluxes%seaice_melt(i,j) > 0.0) &
-          res(i,j) = res(i,j) + fluxes%seaice_melt(i,j)
+        if (associated(fluxes%evap)) then
+          if (fluxes%evap(i,j) > 0.0) res(i,j) = res(i,j) + fluxes%evap(i,j)
+        endif
+        if (associated(fluxes%seaice_melt)) then
+          if (fluxes%seaice_melt(i,j) > 0.0) &
+            res(i,j) = res(i,j) + fluxes%seaice_melt(i,j)
+        endif
       enddo ; enddo
       if (handles%id_net_massin > 0) call post_data(handles%id_net_massin, res, diag)
       if (handles%id_total_net_massin > 0) then

src/parameterizations/vertical/MOM_diabatic_driver.F90

@@ -2938,32 +2938,34 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di
     enddo
   endif
 
-  CS%id_Tdif = register_diag_field('ocean_model',"Tflx_dia_diff",diag%axesTi, &
-      Time, "Diffusive diapycnal temperature flux across interfaces", &
-      "degC m s-1")
-  CS%id_Tadv = register_diag_field('ocean_model',"Tflx_dia_adv",diag%axesTi, &
-      Time, "Advective diapycnal temperature flux across interfaces", &
-      "degC m s-1")
-  CS%id_Sdif = register_diag_field('ocean_model',"Sflx_dia_diff",diag%axesTi, &
-      Time, "Diffusive diapycnal salnity flux across interfaces", &
-      "psu m s-1")
-  CS%id_Sadv = register_diag_field('ocean_model',"Sflx_dia_adv",diag%axesTi, &
-      Time, "Advective diapycnal salnity flux across interfaces", &
-      "psu m s-1")
-  CS%id_MLD_003 = register_diag_field('ocean_model', 'MLD_003', diag%axesT1, Time, &
-      'Mixed layer depth (delta rho = 0.03)', 'm', conversion=US%Z_to_m, &
-      cmor_field_name='mlotst', cmor_long_name='Ocean Mixed Layer Thickness Defined by Sigma T', &
-      cmor_standard_name='ocean_mixed_layer_thickness_defined_by_sigma_t')
-  CS%id_mlotstsq = register_diag_field('ocean_model','mlotstsq',diag%axesT1, Time, &
-      long_name='Square of Ocean Mixed Layer Thickness Defined by Sigma T', &
-      standard_name='square_of_ocean_mixed_layer_thickness_defined_by_sigma_t', &
-      units='m2', conversion=US%Z_to_m**2)
-  CS%id_MLD_0125 = register_diag_field('ocean_model','MLD_0125',diag%axesT1,Time, &
-      'Mixed layer depth (delta rho = 0.125)', 'm', conversion=US%Z_to_m)
-  CS%id_subMLN2  = register_diag_field('ocean_model','subML_N2',diag%axesT1,Time, &
-      'Squared buoyancy frequency below mixed layer', 's-2')
-  CS%id_MLD_user = register_diag_field('ocean_model','MLD_user',diag%axesT1,Time, &
-      'Mixed layer depth (used defined)', 'm', conversion=US%Z_to_m)
+  if (use_temperature) then
+    CS%id_Tdif = register_diag_field('ocean_model',"Tflx_dia_diff",diag%axesTi, &
+        Time, "Diffusive diapycnal temperature flux across interfaces", &
+        "degC m s-1")
+    CS%id_Tadv = register_diag_field('ocean_model',"Tflx_dia_adv",diag%axesTi, &
+        Time, "Advective diapycnal temperature flux across interfaces", &
+        "degC m s-1")
+    CS%id_Sdif = register_diag_field('ocean_model',"Sflx_dia_diff",diag%axesTi, &
+        Time, "Diffusive diapycnal salnity flux across interfaces", &
+        "psu m s-1")
+    CS%id_Sadv = register_diag_field('ocean_model',"Sflx_dia_adv",diag%axesTi, &
+        Time, "Advective diapycnal salnity flux across interfaces", &
+        "psu m s-1")
+    CS%id_MLD_003 = register_diag_field('ocean_model', 'MLD_003', diag%axesT1, Time, &
+        'Mixed layer depth (delta rho = 0.03)', 'm', conversion=US%Z_to_m, &
+        cmor_field_name='mlotst', cmor_long_name='Ocean Mixed Layer Thickness Defined by Sigma T', &
+        cmor_standard_name='ocean_mixed_layer_thickness_defined_by_sigma_t')
+    CS%id_mlotstsq = register_diag_field('ocean_model','mlotstsq',diag%axesT1, Time, &
+        long_name='Square of Ocean Mixed Layer Thickness Defined by Sigma T', &
+        standard_name='square_of_ocean_mixed_layer_thickness_defined_by_sigma_t', &
+        units='m2', conversion=US%Z_to_m**2)
+    CS%id_MLD_0125 = register_diag_field('ocean_model','MLD_0125',diag%axesT1,Time, &
+        'Mixed layer depth (delta rho = 0.125)', 'm', conversion=US%Z_to_m)
+    CS%id_subMLN2  = register_diag_field('ocean_model','subML_N2',diag%axesT1,Time, &
+        'Squared buoyancy frequency below mixed layer', 's-2')
+    CS%id_MLD_user = register_diag_field('ocean_model','MLD_user',diag%axesT1,Time, &
+        'Mixed layer depth (used defined)', 'm', conversion=US%Z_to_m)
+  endif
   call get_param(param_file, mdl, "DIAG_MLD_DENSITY_DIFF", CS%MLDdensityDifference, &
                  "The density difference used to determine a diagnostic mixed "//&
                  "layer depth, MLD_user, following the definition of Levitus 1982. "//&