diff --git a/src/ALE/MOM_ALE.F90 b/src/ALE/MOM_ALE.F90 index 33b498a60a..def4148f25 100644 --- a/src/ALE/MOM_ALE.F90 +++ b/src/ALE/MOM_ALE.F90 @@ -246,9 +246,9 @@ subroutine ALE_register_diags(Time, G, GV, US, diag, CS) ! These diagnostics of the state variables before ALE are useful for ! debugging the ALE code. CS%id_u_preale = register_diag_field('ocean_model', 'u_preale', diag%axesCuL, Time, & - 'Zonal velocity before remapping', 'm s-1') + 'Zonal velocity before remapping', 'm s-1', conversion=US%L_T_to_m_s) CS%id_v_preale = register_diag_field('ocean_model', 'v_preale', diag%axesCvL, Time, & - 'Meridional velocity before remapping', 'm s-1') + 'Meridional velocity before remapping', 'm s-1', conversion=US%L_T_to_m_s) CS%id_h_preale = register_diag_field('ocean_model', 'h_preale', diag%axesTL, Time, & 'Layer Thickness before remapping', get_thickness_units(GV), v_extensive=.true.) CS%id_T_preale = register_diag_field('ocean_model', 'T_preale', diag%axesTL, Time, & diff --git a/src/core/MOM.F90 b/src/core/MOM.F90 index 23c11cc05b..8a5f9bfe02 100644 --- a/src/core/MOM.F90 +++ b/src/core/MOM.F90 @@ -2371,7 +2371,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, & call register_surface_diags(Time, G, CS%sfc_IDs, CS%diag, CS%tv) call register_diags(Time, G, GV, US, CS%IDs, CS%diag) call register_transport_diags(Time, G, GV, US, CS%transport_IDs, CS%diag) - call register_tracer_diagnostics(CS%tracer_Reg, CS%h, Time, diag, G, GV, & + call register_tracer_diagnostics(CS%tracer_Reg, CS%h, Time, diag, G, GV, US, & CS%use_ALE_algorithm) if (CS%use_ALE_algorithm) then call ALE_register_diags(Time, G, GV, US, diag, CS%ALE_CSp) diff --git a/src/core/MOM_barotropic.F90 b/src/core/MOM_barotropic.F90 index 8d48ebbb0b..7b2f367487 100644 --- a/src/core/MOM_barotropic.F90 +++ b/src/core/MOM_barotropic.F90 @@ -4205,13 +4205,13 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, 'Barotropic meridional acceleration from baroclinic terms', 'm s-2', conversion=US%L_T2_to_m_s2) CS%id_eta_bt = register_diag_field('ocean_model', 'eta_bt', diag%axesT1, Time, & - 'Barotropic end SSH', thickness_units) + 'Barotropic end SSH', thickness_units, conversion=GV%H_to_m) CS%id_ubt = register_diag_field('ocean_model', 'ubt', diag%axesCu1, Time, & 'Barotropic end zonal velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_vbt = register_diag_field('ocean_model', 'vbt', diag%axesCv1, Time, & 'Barotropic end meridional velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_eta_st = register_diag_field('ocean_model', 'eta_st', diag%axesT1, Time, & - 'Barotropic start SSH', thickness_units) + 'Barotropic start SSH', thickness_units, conversion=GV%H_to_m) CS%id_ubt_st = register_diag_field('ocean_model', 'ubt_st', diag%axesCu1, Time, & 'Barotropic start zonal velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_vbt_st = register_diag_field('ocean_model', 'vbt_st', diag%axesCv1, Time, & @@ -4221,31 +4221,34 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, CS%id_vbtav = register_diag_field('ocean_model', 'vbtav', diag%axesCv1, Time, & 'Barotropic time-average meridional velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_eta_cor = register_diag_field('ocean_model', 'eta_cor', diag%axesT1, Time, & - 'Corrective mass flux', 'm s-1') + 'Corrective mass flux', 'm s-1', conversion=GV%H_to_m) CS%id_visc_rem_u = register_diag_field('ocean_model', 'visc_rem_u', diag%axesCuL, Time, & 'Viscous remnant at u', 'nondim') CS%id_visc_rem_v = register_diag_field('ocean_model', 'visc_rem_v', diag%axesCvL, Time, & 'Viscous remnant at v', 'nondim') CS%id_gtotn = register_diag_field('ocean_model', 'gtot_n', diag%axesT1, Time, & - 'gtot to North', 'm s-2', conversion=US%L_T_to_m_s**2) + 'gtot to North', 'm s-2', conversion=GV%m_to_H*(US%L_T_to_m_s**2)) CS%id_gtots = register_diag_field('ocean_model', 'gtot_s', diag%axesT1, Time, & - 'gtot to South', 'm s-2', conversion=US%L_T_to_m_s**2) + 'gtot to South', 'm s-2', conversion=GV%m_to_H*(US%L_T_to_m_s**2)) CS%id_gtote = register_diag_field('ocean_model', 'gtot_e', diag%axesT1, Time, & - 'gtot to East', 'm s-2', conversion=US%L_T_to_m_s**2) + 'gtot to East', 'm s-2', conversion=GV%m_to_H*(US%L_T_to_m_s**2)) CS%id_gtotw = register_diag_field('ocean_model', 'gtot_w', diag%axesT1, Time, & - 'gtot to West', 'm s-2', conversion=US%L_T_to_m_s**2) + 'gtot to West', 'm s-2', conversion=GV%m_to_H*(US%L_T_to_m_s**2)) CS%id_eta_hifreq = register_diag_field('ocean_model', 'eta_hifreq', diag%axesT1, Time, & - 'High Frequency Barotropic SSH', thickness_units) + 'High Frequency Barotropic SSH', thickness_units, conversion=GV%H_to_m) CS%id_ubt_hifreq = register_diag_field('ocean_model', 'ubt_hifreq', diag%axesCu1, Time, & 'High Frequency Barotropic zonal velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_vbt_hifreq = register_diag_field('ocean_model', 'vbt_hifreq', diag%axesCv1, Time, & 'High Frequency Barotropic meridional velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_eta_pred_hifreq = register_diag_field('ocean_model', 'eta_pred_hifreq', diag%axesT1, Time, & - 'High Frequency Predictor Barotropic SSH', thickness_units) + 'High Frequency Predictor Barotropic SSH', thickness_units, & + conversion=GV%H_to_m) CS%id_uhbt_hifreq = register_diag_field('ocean_model', 'uhbt_hifreq', diag%axesCu1, Time, & - 'High Frequency Barotropic zonal transport', 'm3 s-1') + 'High Frequency Barotropic zonal transport', 'm3 s-1', & + conversion=GV%H_to_m*US%L_to_m*US%L_T_to_m_s) CS%id_vhbt_hifreq = register_diag_field('ocean_model', 'vhbt_hifreq', diag%axesCv1, Time, & - 'High Frequency Barotropic meridional transport', 'm3 s-1') + 'High Frequency Barotropic meridional transport', 'm3 s-1', & + conversion=GV%H_to_m*US%L_to_m*US%L_T_to_m_s) CS%id_frhatu = register_diag_field('ocean_model', 'frhatu', diag%axesCuL, Time, & 'Fractional thickness of layers in u-columns', 'nondim') CS%id_frhatv = register_diag_field('ocean_model', 'frhatv', diag%axesCvL, Time, & @@ -4255,9 +4258,11 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, CS%id_frhatv1 = register_diag_field('ocean_model', 'frhatv1', diag%axesCvL, Time, & 'Predictor Fractional thickness of layers in v-columns', 'nondim') CS%id_uhbt = register_diag_field('ocean_model', 'uhbt', diag%axesCu1, Time, & - 'Barotropic zonal transport averaged over a baroclinic step', 'm3 s-1') + 'Barotropic zonal transport averaged over a baroclinic step', 'm3 s-1', & + conversion=GV%H_to_m*US%L_to_m*US%L_T_to_m_s) CS%id_vhbt = register_diag_field('ocean_model', 'vhbt', diag%axesCv1, Time, & - 'Barotropic meridional transport averaged over a baroclinic step', 'm3 s-1') + 'Barotropic meridional transport averaged over a baroclinic step', 'm3 s-1', & + conversion=GV%H_to_m*US%L_to_m*US%L_T_to_m_s) if (use_BT_cont_type) then CS%id_BTC_FA_u_EE = register_diag_field('ocean_model', 'BTC_FA_u_EE', diag%axesCu1, Time, & @@ -4269,9 +4274,9 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, CS%id_BTC_FA_u_W0 = register_diag_field('ocean_model', 'BTC_FA_u_W0', diag%axesCu1, Time, & 'BTCont type near west face area', 'm2', conversion=US%L_to_m*GV%H_to_m) CS%id_BTC_ubt_EE = register_diag_field('ocean_model', 'BTC_ubt_EE', diag%axesCu1, Time, & - 'BTCont type far east velocity', 'm s-1') + 'BTCont type far east velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_BTC_ubt_WW = register_diag_field('ocean_model', 'BTC_ubt_WW', diag%axesCu1, Time, & - 'BTCont type far west velocity', 'm s-1') + 'BTCont type far west velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_BTC_FA_v_NN = register_diag_field('ocean_model', 'BTC_FA_v_NN', diag%axesCv1, Time, & 'BTCont type far north face area', 'm2', conversion=US%L_to_m*GV%H_to_m) CS%id_BTC_FA_v_N0 = register_diag_field('ocean_model', 'BTC_FA_v_N0', diag%axesCv1, Time, & diff --git a/src/core/MOM_dynamics_unsplit.F90 b/src/core/MOM_dynamics_unsplit.F90 index 108f4c8943..9e8be65d7a 100644 --- a/src/core/MOM_dynamics_unsplit.F90 +++ b/src/core/MOM_dynamics_unsplit.F90 @@ -678,13 +678,17 @@ subroutine initialize_dyn_unsplit(u, v, h, Time, G, GV, US, param_file, diag, CS 'Meridional Thickness Flux', flux_units, x_cell_method='sum', v_extensive=.true., & conversion=H_convert*US%L_to_m**2*US%s_to_T) CS%id_CAu = register_diag_field('ocean_model', 'CAu', diag%axesCuL, Time, & - 'Zonal Coriolis and Advective Acceleration', 'meter second-2, conversion=US%L_T2_to_m_s2') + 'Zonal Coriolis and Advective Acceleration', 'm s-2', & + conversion=US%L_T2_to_m_s2) CS%id_CAv = register_diag_field('ocean_model', 'CAv', diag%axesCvL, Time, & - 'Meridional Coriolis and Advective Acceleration', 'meter second-2', conversion=US%L_T2_to_m_s2) + 'Meridional Coriolis and Advective Acceleration', 'm s-2', & + conversion=US%L_T2_to_m_s2) CS%id_PFu = register_diag_field('ocean_model', 'PFu', diag%axesCuL, Time, & - 'Zonal Pressure Force Acceleration', 'meter second-2', conversion=US%L_T2_to_m_s2) + 'Zonal Pressure Force Acceleration', 'm s-2', & + conversion=US%L_T2_to_m_s2) 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) + 'Meridional Pressure Force Acceleration', 'm s-2', & + conversion=US%L_T2_to_m_s2) id_clock_Cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=CLOCK_MODULE) id_clock_continuity = cpu_clock_id('(Ocean continuity equation)', grain=CLOCK_MODULE) diff --git a/src/diagnostics/MOM_diagnostics.F90 b/src/diagnostics/MOM_diagnostics.F90 index 54025a0ac0..4aa6fad710 100644 --- a/src/diagnostics/MOM_diagnostics.F90 +++ b/src/diagnostics/MOM_diagnostics.F90 @@ -92,17 +92,18 @@ module MOM_diagnostics ! The following arrays hold diagnostics in the layer-integrated energy budget. real, pointer, dimension(:,:,:) :: & - KE => NULL(), & !< KE per unit mass [m2 s-2] - dKE_dt => NULL(), & !< time derivative of the layer KE [m3 s-3] + KE => NULL(), & !< KE per unit mass [L2 T-2 ~> m2 s-2] + dKE_dt => NULL(), & !< time derivative of the layer KE [H L2 T-3 ~> m3 s-3] PE_to_KE => NULL(), & !< potential energy to KE term [m3 s-3] - KE_CorAdv => NULL(), & !< KE source from the combined Coriolis and advection terms [m3 s-3]. + KE_CorAdv => NULL(), & !< KE source from the combined Coriolis and + !! advection terms [H L2 T-3 ~> m3 s-3]. !! The Coriolis source should be zero, but is not due to truncation !! errors. There should be near-cancellation of the global integral !! of this spurious Coriolis source. - KE_adv => NULL(), & !< KE source from along-layer advection [m3 s-3] - KE_visc => NULL(), & !< KE source from vertical viscosity [m3 s-3] - KE_horvisc => NULL(), & !< KE source from horizontal viscosity [m3 s-3] - KE_dia => NULL() !< KE source from diapycnal diffusion [m3 s-3] + KE_adv => NULL(), & !< KE source from along-layer advection [H L2 T-3 ~> m3 s-3] + KE_visc => NULL(), & !< KE source from vertical viscosity [H L2 T-3 ~> m3 s-3] + KE_horvisc => NULL(), & !< KE source from horizontal viscosity [H L2 T-3 ~> m3 s-3] + KE_dia => NULL() !< KE source from diapycnal diffusion [H L2 T-3 ~> m3 s-3] !>@{ Diagnostic IDs integer :: id_u = -1, id_v = -1, id_h = -1 @@ -254,7 +255,7 @@ subroutine calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, & if (loc(CS)==0) call MOM_error(FATAL, & "calculate_diagnostic_fields: Module must be initialized before used.") - call calculate_derivs(dt, G, CS) + call calculate_derivs(US%s_to_T*dt, G, CS) if (dt > 0.0) then call diag_save_grids(CS%diag) @@ -916,8 +917,8 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%KE)) then do k=1,nz ; do j=js,je ; do i=is,ie - CS%KE(i,j,k) = US%L_T_to_m_s**2*((u(I,j,k)*u(I,j,k) + u(I-1,j,k)*u(I-1,j,k)) + & - (v(i,J,k)*v(i,J,k) + v(i,J-1,k)*v(i,J-1,k)))*0.25 + CS%KE(i,j,k) = ((u(I,j,k) * u(I,j,k) + u(I-1,j,k) * u(I-1,j,k)) & + + (v(i,J,k) * v(i,J,k) + v(i,J-1,k) * v(i,J-1,k))) * 0.25 ! DELETE THE FOLLOWING... Make this 0 to test the momentum balance, ! or a huge number to test the continuity balance. ! CS%KE(i,j,k) *= 1e20 @@ -936,19 +937,19 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%dKE_dt)) then do k=1,nz do j=js,je ; do I=Isq,Ieq - KE_u(I,j) = US%L_T_to_m_s**2*US%s_to_T*uh(I,j,k)*G%dxCu(I,j)*CS%du_dt(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * CS%du_dt(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie - KE_v(i,J) = US%L_T_to_m_s**2*US%s_to_T*vh(i,J,k)*G%dyCv(i,J)*CS%dv_dt(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * CS%dv_dt(i,J,k) enddo ; enddo do j=js,je ; do i=is,ie - KE_h(i,j) = CS%KE(i,j,k)*US%s_to_T*CS%dh_dt(i,j,k) + KE_h(i,j) = CS%KE(i,j,k) * CS%dh_dt(i,j,k) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%dKE_dt(i,j,k) = GV%H_to_m * (KE_h(i,j) + 0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1))) + CS%dKE_dt(i,j,k) = KE_h(i,j) + 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_dKEdt > 0) call post_data(CS%id_dKEdt, CS%dKE_dt, CS%diag) @@ -957,16 +958,16 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%PE_to_KE)) then do k=1,nz do j=js,je ; do I=Isq,Ieq - KE_u(I,j) = US%L_T_to_m_s**2*US%s_to_T*uh(I,j,k)*G%dxCu(I,j)*ADp%PFu(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * ADp%PFu(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie - KE_v(i,J) = US%L_T_to_m_s**2*US%s_to_T*vh(i,J,k)*G%dyCv(i,J)*ADp%PFv(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * ADp%PFv(i,J,k) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%PE_to_KE(i,j,k) = GV%H_to_m * 0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) + CS%PE_to_KE(i,j,k) = 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_PE_to_KE > 0) call post_data(CS%id_PE_to_KE, CS%PE_to_KE, CS%diag) @@ -975,20 +976,20 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%KE_CorAdv)) then do k=1,nz do j=js,je ; do I=Isq,Ieq - KE_u(I,j) = US%L_T_to_m_s**2*US%s_to_T*uh(I,j,k)*G%dxCu(I,j)*ADp%CAu(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * ADp%CAu(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie - KE_v(i,J) = US%L_T_to_m_s**2*US%s_to_T*vh(i,J,k)*G%dyCv(i,J)*ADp%CAv(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * ADp%CAv(i,J,k) enddo ; enddo do j=js,je ; do i=is,ie - KE_h(i,j) = -CS%KE(i,j,k) * G%IareaT(i,j) * & - US%s_to_T*(uh(I,j,k) - uh(I-1,j,k) + vh(i,J,k) - vh(i,J-1,k)) + KE_h(i,j) = -CS%KE(i,j,k) * G%IareaT(i,j) & + * (uh(I,j,k) - uh(I-1,j,k) + vh(i,J,k) - vh(i,J-1,k)) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%KE_CorAdv(i,j,k) = GV%H_to_m * (KE_h(i,j) + 0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1))) + CS%KE_CorAdv(i,j,k) = KE_h(i,j) + 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_KE_Coradv > 0) call post_data(CS%id_KE_Coradv, CS%KE_Coradv, CS%diag) @@ -1002,21 +1003,21 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS do k=1,nz do j=js,je ; do I=Isq,Ieq if (G%mask2dCu(i,j) /= 0.) & - KE_u(I,j) = US%L_T_to_m_s**2*US%s_to_T*uh(I,j,k)*G%dxCu(I,j)*ADp%gradKEu(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * ADp%gradKEu(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie if (G%mask2dCv(i,j) /= 0.) & - KE_v(i,J) = US%L_T_to_m_s**2*US%s_to_T*vh(i,J,k)*G%dyCv(i,J)*ADp%gradKEv(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * ADp%gradKEv(i,J,k) enddo ; enddo do j=js,je ; do i=is,ie - KE_h(i,j) = -CS%KE(i,j,k) * G%IareaT(i,j) * & - US%s_to_T*(uh(I,j,k) - uh(I-1,j,k) + vh(i,J,k) - vh(i,J-1,k)) + KE_h(i,j) = -CS%KE(i,j,k) * G%IareaT(i,j) & + * (uh(I,j,k) - uh(I-1,j,k) + vh(i,J,k) - vh(i,J-1,k)) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%KE_adv(i,j,k) = GV%H_to_m * (KE_h(i,j) + 0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1))) + CS%KE_adv(i,j,k) = KE_h(i,j) + 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_KE_adv > 0) call post_data(CS%id_KE_adv, CS%KE_adv, CS%diag) @@ -1025,16 +1026,16 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%KE_visc)) then do k=1,nz do j=js,je ; do I=Isq,Ieq - KE_u(I,j) = US%L_T_to_m_s**2*US%s_to_T*uh(I,j,k)*G%dxCu(I,j)*ADp%du_dt_visc(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * ADp%du_dt_visc(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie - KE_v(i,J) = US%L_T_to_m_s**2*US%s_to_T*vh(i,J,k)*G%dyCv(i,J)*ADp%dv_dt_visc(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * ADp%dv_dt_visc(i,J,k) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%KE_visc(i,j,k) = GV%H_to_m * (0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1))) + CS%KE_visc(i,j,k) = 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_KE_visc > 0) call post_data(CS%id_KE_visc, CS%KE_visc, CS%diag) @@ -1043,16 +1044,16 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%KE_horvisc)) then do k=1,nz do j=js,je ; do I=Isq,Ieq - KE_u(I,j) = US%s_to_T*uh(I,j,k)*US%L_to_m*G%dxCu(I,j)*US%L_T2_to_m_s2*ADp%diffu(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * ADp%diffu(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie - KE_v(i,J) = US%s_to_T*vh(i,J,k)*US%L_to_m*G%dyCv(i,J)*US%L_T2_to_m_s2*ADp%diffv(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * ADp%diffv(i,J,k) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%KE_horvisc(i,j,k) = GV%H_to_m * 0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) + CS%KE_horvisc(i,j,k) = 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_KE_horvisc > 0) call post_data(CS%id_KE_horvisc, CS%KE_horvisc, CS%diag) @@ -1061,20 +1062,20 @@ subroutine calculate_energy_diagnostics(u, v, h, uh, vh, ADp, CDp, G, GV, US, CS if (associated(CS%KE_dia)) then do k=1,nz do j=js,je ; do I=Isq,Ieq - KE_u(I,j) = US%L_T_to_m_s**2*US%s_to_T*uh(I,j,k)*G%dxCu(I,j)*ADp%du_dt_dia(I,j,k) + KE_u(I,j) = uh(I,j,k) * G%dxCu(I,j) * ADp%du_dt_dia(I,j,k) enddo ; enddo do J=Jsq,Jeq ; do i=is,ie - KE_v(i,J) = US%L_T_to_m_s**2*US%s_to_T*vh(i,J,k)*G%dyCv(i,J)*ADp%dv_dt_dia(i,J,k) + KE_v(i,J) = vh(i,J,k) * G%dyCv(i,J) * ADp%dv_dt_dia(i,J,k) enddo ; enddo do j=js,je ; do i=is,ie - KE_h(i,j) = CS%KE(i,j,k) * & - (CDp%diapyc_vel(i,j,k) - CDp%diapyc_vel(i,j,k+1)) + KE_h(i,j) = CS%KE(i,j,k) & + * (US%T_to_s * (CDp%diapyc_vel(i,j,k) - CDp%diapyc_vel(i,j,k+1))) enddo ; enddo if (.not.G%symmetric) & call do_group_pass(CS%pass_KE_uv, G%domain) do j=js,je ; do i=is,ie - CS%KE_dia(i,j,k) = KE_h(i,j) + GV%H_to_m * 0.5 * G%IareaT(i,j) * & - (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) + CS%KE_dia(i,j,k) = KE_h(i,j) + 0.5 * G%IareaT(i,j) & + * (KE_u(I,j) + KE_u(I-1,j) + KE_v(i,J) + KE_v(i,J-1)) enddo ; enddo enddo if (CS%id_KE_dia > 0) call post_data(CS%id_KE_dia, CS%KE_dia, CS%diag) @@ -1495,7 +1496,8 @@ subroutine MOM_diagnostics_init(MIS, ADp, CDp, Time, G, GV, US, param_file, diag CS%id_thkcello = register_diag_field('ocean_model', 'thkcello', diag%axesTL, Time, & long_name = 'Cell Thickness', standard_name='cell_thickness', units='m', v_extensive=.true.) CS%id_h_pre_sync = register_diag_field('ocean_model', 'h_pre_sync', diag%axesTL, Time, & - long_name = 'Cell thickness from the previous timestep', units='m', v_extensive=.true.) + long_name = 'Cell thickness from the previous timestep', units='m', & + v_extensive=.true., conversion=GV%H_to_m) ! Note that CS%id_volcello would normally be registered here but because it is a "cell measure" and ! must be registered first. We earlier stored the handle of volcello but need it here for posting @@ -1625,36 +1627,44 @@ subroutine MOM_diagnostics_init(MIS, ADp, CDp, Time, G, GV, US, param_file, diag ! terms in the kinetic energy budget CS%id_KE = register_diag_field('ocean_model', 'KE', diag%axesTL, Time, & - 'Layer kinetic energy per unit mass', 'm2 s-2') + 'Layer kinetic energy per unit mass', 'm2 s-2', & + conversion=US%L_T_to_m_s**2) if (CS%id_KE>0) call safe_alloc_ptr(CS%KE,isd,ied,jsd,jed,nz) CS%id_dKEdt = register_diag_field('ocean_model', 'dKE_dt', diag%axesTL, Time, & - 'Kinetic Energy Tendency of Layer', 'm3 s-3') + 'Kinetic Energy Tendency of Layer', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_dKEdt>0) call safe_alloc_ptr(CS%dKE_dt,isd,ied,jsd,jed,nz) CS%id_PE_to_KE = register_diag_field('ocean_model', 'PE_to_KE', diag%axesTL, Time, & - 'Potential to Kinetic Energy Conversion of Layer', 'm3 s-3') + 'Potential to Kinetic Energy Conversion of Layer', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_PE_to_KE>0) call safe_alloc_ptr(CS%PE_to_KE,isd,ied,jsd,jed,nz) CS%id_KE_Coradv = register_diag_field('ocean_model', 'KE_Coradv', diag%axesTL, Time, & - 'Kinetic Energy Source from Coriolis and Advection', 'm3 s-3') + 'Kinetic Energy Source from Coriolis and Advection', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_KE_Coradv>0) call safe_alloc_ptr(CS%KE_Coradv,isd,ied,jsd,jed,nz) CS%id_KE_adv = register_diag_field('ocean_model', 'KE_adv', diag%axesTL, Time, & - 'Kinetic Energy Source from Advection', 'm3 s-3') + 'Kinetic Energy Source from Advection', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_KE_adv>0) call safe_alloc_ptr(CS%KE_adv,isd,ied,jsd,jed,nz) CS%id_KE_visc = register_diag_field('ocean_model', 'KE_visc', diag%axesTL, Time, & - 'Kinetic Energy Source from Vertical Viscosity and Stresses', 'm3 s-3') + 'Kinetic Energy Source from Vertical Viscosity and Stresses', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_KE_visc>0) call safe_alloc_ptr(CS%KE_visc,isd,ied,jsd,jed,nz) CS%id_KE_horvisc = register_diag_field('ocean_model', 'KE_horvisc', diag%axesTL, Time, & - 'Kinetic Energy Source from Horizontal Viscosity', 'm3 s-3') + 'Kinetic Energy Source from Horizontal Viscosity', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_KE_horvisc>0) call safe_alloc_ptr(CS%KE_horvisc,isd,ied,jsd,jed,nz) if (.not. adiabatic) then CS%id_KE_dia = register_diag_field('ocean_model', 'KE_dia', diag%axesTL, Time, & - 'Kinetic Energy Source from Diapycnal Diffusion', 'm3 s-3') + 'Kinetic Energy Source from Diapycnal Diffusion', 'm3 s-3', & + conversion=GV%H_to_m*(US%L_T_to_m_s**2)*US%s_to_T) if (CS%id_KE_dia>0) call safe_alloc_ptr(CS%KE_dia,isd,ied,jsd,jed,nz) endif @@ -1832,10 +1842,10 @@ subroutine register_transport_diags(Time, G, GV, US, IDs, diag) standard_name='ocean_mass_y_transport_vertical_sum', x_cell_method='sum') IDs%id_dynamics_h = register_diag_field('ocean_model','dynamics_h', & diag%axesTl, Time, 'Change in layer thicknesses due to horizontal dynamics', & - 'm s-1', v_extensive = .true.) + 'm s-1', v_extensive=.true., conversion=GV%H_to_m) IDs%id_dynamics_h_tendency = register_diag_field('ocean_model','dynamics_h_tendency', & diag%axesTl, Time, 'Change in layer thicknesses due to horizontal dynamics', & - 'm s-1', v_extensive = .true.) + 'm s-1', v_extensive=.true.) end subroutine register_transport_diags @@ -1882,34 +1892,34 @@ subroutine write_static_fields(G, GV, US, tv, diag) 'Longitude of zonal velocity (Cu) points', 'degrees_east', interp_method='none') if (id > 0) call post_data(id, G%geoLonCu, diag, .true.) - id = register_static_field('ocean_model', 'area_t', diag%axesT1, & - 'Surface area of tracer (T) cells', 'm2', conversion=US%m_to_L**2, & + id = register_static_field('ocean_model', 'area_t', diag%axesT1, & + 'Surface area of tracer (T) cells', 'm2', conversion=US%L_to_m**2, & cmor_field_name='areacello', cmor_standard_name='cell_area', & - cmor_long_name='Ocean Grid-Cell Area', & + cmor_long_name='Ocean Grid-Cell Area', & x_cell_method='sum', y_cell_method='sum', area_cell_method='sum') if (id > 0) then call post_data(id, G%areaT, diag, .true.) call diag_register_area_ids(diag, id_area_t=id) endif - id = register_static_field('ocean_model', 'area_u', diag%axesCu1, & - 'Surface area of x-direction flow (U) cells', 'm2', conversion=US%m_to_L**2, & + id = register_static_field('ocean_model', 'area_u', diag%axesCu1, & + 'Surface area of x-direction flow (U) cells', 'm2', conversion=US%L_to_m**2, & cmor_field_name='areacello_cu', cmor_standard_name='cell_area', & - cmor_long_name='Ocean Grid-Cell Area', & + cmor_long_name='Ocean Grid-Cell Area', & x_cell_method='sum', y_cell_method='sum', area_cell_method='sum') if (id > 0) call post_data(id, G%areaCu, diag, .true.) - id = register_static_field('ocean_model', 'area_v', diag%axesCv1, & - 'Surface area of y-direction flow (V) cells', 'm2', conversion=US%m_to_L**2, & + id = register_static_field('ocean_model', 'area_v', diag%axesCv1, & + 'Surface area of y-direction flow (V) cells', 'm2', conversion=US%L_to_m**2, & cmor_field_name='areacello_cv', cmor_standard_name='cell_area', & - cmor_long_name='Ocean Grid-Cell Area', & + cmor_long_name='Ocean Grid-Cell Area', & x_cell_method='sum', y_cell_method='sum', area_cell_method='sum') if (id > 0) call post_data(id, G%areaCv, diag, .true.) - id = register_static_field('ocean_model', 'area_q', diag%axesB1, & - 'Surface area of B-grid flow (Q) cells', 'm2', conversion=US%m_to_L**2, & + id = register_static_field('ocean_model', 'area_q', diag%axesB1, & + 'Surface area of B-grid flow (Q) cells', 'm2', conversion=US%L_to_m**2, & cmor_field_name='areacello_bu', cmor_standard_name='cell_area', & - cmor_long_name='Ocean Grid-Cell Area', & + cmor_long_name='Ocean Grid-Cell Area', & x_cell_method='sum', y_cell_method='sum', area_cell_method='sum') if (id > 0) call post_data(id, G%areaBu, diag, .true.) diff --git a/src/framework/MOM_spatial_means.F90 b/src/framework/MOM_spatial_means.F90 index f7084ee7ea..f6282faa52 100644 --- a/src/framework/MOM_spatial_means.F90 +++ b/src/framework/MOM_spatial_means.F90 @@ -35,8 +35,8 @@ function global_area_mean(var,G) is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec tmpForSumming(:,:) = 0. - do j=js,je ; do i=is, ie - tmpForSumming(i,j) = ( var(i,j) * (G%US%L_to_m**2*G%areaT(i,j) * G%mask2dT(i,j)) ) + do j=js,je ; do i=is,ie + tmpForSumming(i,j) = var(i,j) * (G%areaT(i,j) * G%mask2dT(i,j)) enddo ; enddo global_area_mean = reproducing_sum( tmpForSumming ) * G%IareaT_global diff --git a/src/parameterizations/lateral/MOM_MEKE.F90 b/src/parameterizations/lateral/MOM_MEKE.F90 index d204db1305..6c038a4eb5 100644 --- a/src/parameterizations/lateral/MOM_MEKE.F90 +++ b/src/parameterizations/lateral/MOM_MEKE.F90 @@ -1148,7 +1148,7 @@ logical function MEKE_init(Time, G, US, param_file, diag, CS, MEKE, restart_CS) 'MEKE derived barotropic eddy-velocity scale', 'm s-1', conversion=US%L_T_to_m_s) if (.not. associated(MEKE%MEKE)) CS%id_Ut = -1 CS%id_src = register_diag_field('ocean_model', 'MEKE_src', diag%axesT1, Time, & - 'MEKE energy source', 'm2 s-3') + 'MEKE energy source', 'm2 s-3', conversion=(US%L_T_to_m_s**2)*US%s_to_T) CS%id_decay = register_diag_field('ocean_model', 'MEKE_decay', diag%axesT1, Time, & 'MEKE decay rate', 's-1', conversion=US%s_to_T) CS%id_GM_src = register_diag_field('ocean_model', 'MEKE_GM_src', diag%axesT1, Time, & diff --git a/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 b/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 index c4a2d0c38f..6543ab7af6 100644 --- a/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 +++ b/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 @@ -884,8 +884,8 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS, CS%diag => diag - if (GV%Boussinesq) then ; flux_to_kg_per_s = GV%Rho0*US%L_to_m**2*US%s_to_T - else ; flux_to_kg_per_s = US%L_to_m**2*US%s_to_T ; endif + if (GV%Boussinesq) then ; flux_to_kg_per_s = GV%H_to_kg_m2 * US%L_to_m**2 * US%s_to_T + else ; flux_to_kg_per_s = GV%H_to_m * US%L_to_m**2 * US%s_to_T ; endif CS%id_uhml = register_diag_field('ocean_model', 'uhml', diag%axesCuL, Time, & 'Zonal Thickness Flux to Restratify Mixed Layer', 'kg s-1', conversion=flux_to_kg_per_s, & diff --git a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 index cbf42d2b8b..2a17bfbd6f 100644 --- a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 +++ b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 @@ -3582,33 +3582,34 @@ subroutine bulkmixedlayer_init(Time, G, GV, US, param_file, diag, CS) Time, 'Surface mixed layer depth', 'm') CS%id_TKE_wind = register_diag_field('ocean_model', 'TKE_wind', diag%axesT1, & Time, 'Wind-stirring source of mixed layer TKE', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_RiBulk = register_diag_field('ocean_model', 'TKE_RiBulk', diag%axesT1, & Time, 'Mean kinetic energy source of mixed layer TKE', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_conv = register_diag_field('ocean_model', 'TKE_conv', diag%axesT1, & - Time, 'Convective source of mixed layer TKE', 'm3 s-3', conversion=US%Z_to_m*US%T_to_s**3) + Time, 'Convective source of mixed layer TKE', 'm3 s-3', & + conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_pen_SW = register_diag_field('ocean_model', 'TKE_pen_SW', diag%axesT1, & Time, 'TKE consumed by mixing penetrative shortwave radation through the mixed layer', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_mixing = register_diag_field('ocean_model', 'TKE_mixing', diag%axesT1, & Time, 'TKE consumed by mixing that deepens the mixed layer', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_mech_decay = register_diag_field('ocean_model', 'TKE_mech_decay', diag%axesT1, & Time, 'Mechanical energy decay sink of mixed layer TKE', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_conv_decay = register_diag_field('ocean_model', 'TKE_conv_decay', diag%axesT1, & Time, 'Convective energy decay sink of mixed layer TKE', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_TKE_conv_s2 = register_diag_field('ocean_model', 'TKE_conv_s2', diag%axesT1, & Time, 'Spurious source of mixed layer TKE from sigma2', & - 'm3 s-3', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'm3 s-3', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_PE_detrain = register_diag_field('ocean_model', 'PE_detrain', diag%axesT1, & Time, 'Spurious source of potential energy from mixed layer detrainment', & - 'W m-2', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'W m-2', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_PE_detrain2 = register_diag_field('ocean_model', 'PE_detrain2', diag%axesT1, & Time, 'Spurious source of potential energy from mixed layer only detrainment', & - 'W m-2', conversion=US%Z_to_m*US%L_to_m**2*US%T_to_s**3) + 'W m-2', conversion=US%Z_to_m*(US%L_to_m**2)*(US%s_to_T**3)) CS%id_h_mismatch = register_diag_field('ocean_model', 'h_miss_ML', diag%axesT1, & Time, 'Summed absolute mismatch in entrainment terms', 'm', conversion=US%Z_to_m) CS%id_Hsfc_used = register_diag_field('ocean_model', 'Hs_used', diag%axesT1, & diff --git a/src/parameterizations/vertical/MOM_vert_friction.F90 b/src/parameterizations/vertical/MOM_vert_friction.F90 index 1bed36e75e..b282995d3f 100644 --- a/src/parameterizations/vertical/MOM_vert_friction.F90 +++ b/src/parameterizations/vertical/MOM_vert_friction.F90 @@ -1712,16 +1712,20 @@ subroutine vertvisc_init(MIS, Time, G, GV, US, param_file, diag, ADp, dirs, & 'Meridional Viscous Vertical Coupling Coefficient', 'm s-1', conversion=US%Z_to_m*US%s_to_T) CS%id_h_u = register_diag_field('ocean_model', 'Hu_visc', diag%axesCuL, Time, & - 'Thickness at Zonal Velocity Points for Viscosity', thickness_units) + 'Thickness at Zonal Velocity Points for Viscosity', thickness_units, & + conversion=GV%H_to_m) CS%id_h_v = register_diag_field('ocean_model', 'Hv_visc', diag%axesCvL, Time, & - 'Thickness at Meridional Velocity Points for Viscosity', thickness_units) + 'Thickness at Meridional Velocity Points for Viscosity', thickness_units, & + conversion=GV%H_to_m) CS%id_hML_u = register_diag_field('ocean_model', 'HMLu_visc', diag%axesCu1, Time, & - 'Mixed Layer Thickness at Zonal Velocity Points for Viscosity', thickness_units) + 'Mixed Layer Thickness at Zonal Velocity Points for Viscosity', thickness_units, & + conversion=GV%H_to_m) CS%id_hML_v = register_diag_field('ocean_model', 'HMLv_visc', diag%axesCv1, Time, & - 'Mixed Layer Thickness at Meridional Velocity Points for Viscosity', thickness_units) + 'Mixed Layer Thickness at Meridional Velocity Points for Viscosity', thickness_units, & + conversion=GV%H_to_m) CS%id_du_dt_visc = register_diag_field('ocean_model', 'du_dt_visc', diag%axesCuL, & Time, 'Zonal Acceleration from Vertical Viscosity', 'm s-2', conversion=US%L_T2_to_m_s2) diff --git a/src/tracer/MOM_tracer_hor_diff.F90 b/src/tracer/MOM_tracer_hor_diff.F90 index 098a647ec8..5577115a48 100644 --- a/src/tracer/MOM_tracer_hor_diff.F90 +++ b/src/tracer/MOM_tracer_hor_diff.F90 @@ -449,20 +449,20 @@ subroutine tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online Coef_y(i,J) * (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i,j+1,k)))) enddo ; enddo if (associated(Reg%Tr(m)%df_x)) then ; do j=js,je ; do I=G%IscB,G%IecB - Reg%Tr(m)%df_x(I,j,k) = Reg%Tr(m)%df_x(I,j,k) + US%L_to_m**2*Coef_x(I,j) * & - (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i+1,j,k))*Idt + Reg%Tr(m)%df_x(I,j,k) = Reg%Tr(m)%df_x(I,j,k) + Coef_x(I,j) & + * (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i+1,j,k)) * Idt enddo ; enddo ; endif if (associated(Reg%Tr(m)%df_y)) then ; do J=G%JscB,G%JecB ; do i=is,ie - Reg%Tr(m)%df_y(i,J,k) = Reg%Tr(m)%df_y(i,J,k) + US%L_to_m**2*Coef_y(i,J) * & - (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i,j+1,k))*Idt + Reg%Tr(m)%df_y(i,J,k) = Reg%Tr(m)%df_y(i,J,k) + Coef_y(i,J) & + * (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i,j+1,k)) * Idt enddo ; enddo ; endif if (associated(Reg%Tr(m)%df2d_x)) then ; do j=js,je ; do I=G%IscB,G%IecB - Reg%Tr(m)%df2d_x(I,j) = Reg%Tr(m)%df2d_x(I,j) + US%L_to_m**2*Coef_x(I,j) * & - (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i+1,j,k))*Idt + Reg%Tr(m)%df2d_x(I,j) = Reg%Tr(m)%df2d_x(I,j) + Coef_x(I,j) & + * (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i+1,j,k)) * Idt enddo ; enddo ; endif if (associated(Reg%Tr(m)%df2d_y)) then ; do J=G%JscB,G%JecB ; do i=is,ie - Reg%Tr(m)%df2d_y(i,J) = Reg%Tr(m)%df2d_y(i,J) + US%L_to_m**2*Coef_y(i,J) * & - (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i,j+1,k))*Idt + Reg%Tr(m)%df2d_y(i,J) = Reg%Tr(m)%df2d_y(i,J) + Coef_y(i,J) & + * (Reg%Tr(m)%t(i,j,k) - Reg%Tr(m)%t(i,j+1,k)) * Idt enddo ; enddo ; endif do j=js,je ; do i=is,ie Reg%Tr(m)%t(i,j,k) = Reg%Tr(m)%t(i,j,k) + dTr(i,j) diff --git a/src/tracer/MOM_tracer_registry.F90 b/src/tracer/MOM_tracer_registry.F90 index 4680c058b4..41299be3e8 100644 --- a/src/tracer/MOM_tracer_registry.F90 +++ b/src/tracer/MOM_tracer_registry.F90 @@ -20,6 +20,7 @@ module MOM_tracer_registry use MOM_restart, only : register_restart_field, MOM_restart_CS use MOM_string_functions, only : lowercase use MOM_time_manager, only : time_type +use MOM_unit_scaling, only : unit_scale_type use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -53,13 +54,13 @@ module MOM_tracer_registry !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1] real, dimension(:,:,:), pointer :: df_x => NULL() !< diagnostic array for x-diffusive tracer flux - !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1] + !! [conc H L2 s-1 ~> conc m3 s-1 or conc kg s-1] real, dimension(:,:,:), pointer :: df_y => NULL() !< diagnostic array for y-diffusive tracer flux - !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1] + !! [conc H L2 s-1 ~> conc m3 s-1 or conc kg s-1] real, dimension(:,:), pointer :: df2d_x => NULL() !< diagnostic vertical sum x-diffusive flux - !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1] + !! [conc H L2 s-1 ~> conc m3 s-1 or conc kg s-1] real, dimension(:,:), pointer :: df2d_y => NULL() !< diagnostic vertical sum y-diffusive flux - !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1] + !! [conc H L2 s-1 ~> conc m3 s-1 or conc kg s-1] real, dimension(:,:), pointer :: df2d_conc_x => NULL() !< diagnostic vertical sum x-diffusive content flux !! [conc H m2 s-1 ~> conc m3 s-1 or conc kg s-1] real, dimension(:,:), pointer :: df2d_conc_y => NULL() !< diagnostic vertical sum y-diffusive content flux @@ -320,9 +321,10 @@ end subroutine lock_tracer_registry !> register_tracer_diagnostics does a set of register_diag_field calls for any previously !! registered in a tracer registry with a value of registry_diags set to .true. -subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, use_ALE) +subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, US, use_ALE) type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(tracer_registry_type), pointer :: Reg !< pointer to the tracer registry real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & intent(in) :: h !< Layer thicknesses @@ -396,10 +398,12 @@ subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, use_ALE) trim(flux_units), v_extensive = .true., x_cell_method = 'sum') Tr%id_dfx = register_diag_field("ocean_model", trim(shortnm)//"_dfx", & diag%axesCuL, Time, trim(flux_longname)//" diffusive zonal flux" , & - trim(flux_units), v_extensive = .true., y_cell_method = 'sum') + trim(flux_units), v_extensive = .true., conversion=US%L_to_m**2, & + y_cell_method = 'sum') Tr%id_dfy = register_diag_field("ocean_model", trim(shortnm)//"_dfy", & diag%axesCvL, Time, trim(flux_longname)//" diffusive zonal flux" , & - trim(flux_units), v_extensive = .true., x_cell_method = 'sum') + trim(flux_units), v_extensive = .true., conversion=US%L_to_m**2, & + x_cell_method = 'sum') else Tr%id_adx = register_diag_field("ocean_model", trim(shortnm)//"_adx", & diag%axesCuL, Time, "Advective (by residual mean) Zonal Flux of "//trim(flux_longname), & @@ -409,10 +413,12 @@ subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, use_ALE) flux_units, v_extensive=.true., conversion=Tr%flux_scale, x_cell_method = 'sum') Tr%id_dfx = register_diag_field("ocean_model", trim(shortnm)//"_diffx", & diag%axesCuL, Time, "Diffusive Zonal Flux of "//trim(flux_longname), & - flux_units, v_extensive=.true., conversion=Tr%flux_scale, y_cell_method = 'sum') + flux_units, v_extensive=.true., conversion=(US%L_to_m**2)*Tr%flux_scale, & + y_cell_method='sum') Tr%id_dfy = register_diag_field("ocean_model", trim(shortnm)//"_diffy", & diag%axesCvL, Time, "Diffusive Meridional Flux of "//trim(flux_longname), & - flux_units, v_extensive=.true., conversion=Tr%flux_scale, x_cell_method = 'sum') + flux_units, v_extensive=.true., conversion=(US%L_to_m**2)*Tr%flux_scale, & + x_cell_method='sum') endif if (Tr%id_adx > 0) call safe_alloc_ptr(Tr%ad_x,IsdB,IedB,jsd,jed,nz) if (Tr%id_ady > 0) call safe_alloc_ptr(Tr%ad_y,isd,ied,JsdB,JedB,nz) @@ -430,11 +436,13 @@ subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, use_ALE) Tr%id_dfx_2d = register_diag_field("ocean_model", trim(shortnm)//"_diffx_2d", & diag%axesCu1, Time, & "Vertically Integrated Diffusive Zonal Flux of "//trim(flux_longname), & - flux_units, conversion=Tr%flux_scale, y_cell_method = 'sum') + flux_units, conversion=(US%L_to_m**2)*Tr%flux_scale, & + y_cell_method='sum') Tr%id_dfy_2d = register_diag_field("ocean_model", trim(shortnm)//"_diffy_2d", & diag%axesCv1, Time, & "Vertically Integrated Diffusive Meridional Flux of "//trim(flux_longname), & - flux_units, conversion=Tr%flux_scale, x_cell_method = 'sum') + flux_units, conversion=(US%L_to_m**2)*Tr%flux_scale, & + x_cell_method='sum') if (Tr%id_adx_2d > 0) call safe_alloc_ptr(Tr%ad2d_x,IsdB,IedB,jsd,jed) if (Tr%id_ady_2d > 0) call safe_alloc_ptr(Tr%ad2d_y,isd,ied,JsdB,JedB)