diff --git a/src/core/MOM_checksum_packages.F90 b/src/core/MOM_checksum_packages.F90 index 917a4afdc3..d9855a98d3 100644 --- a/src/core/MOM_checksum_packages.F90 +++ b/src/core/MOM_checksum_packages.F90 @@ -253,8 +253,8 @@ subroutine MOM_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, pe real, dimension(G%isc:G%iec, G%jsc:G%jec) :: & tmp_A, & ! The area per cell [m2] (unscaled to permit reproducing sum). tmp_V, & ! The column-integrated volume [m3] (unscaled to permit reproducing sum) - tmp_T, & ! The column-integrated temperature [degC m3] - tmp_S ! The column-integrated salinity [ppt m3] + tmp_T, & ! The column-integrated temperature [degC m3] (unscaled to permit reproducing sum) + tmp_S ! The column-integrated salinity [ppt m3] (unscaled to permit reproducing sum) real :: Vol, dV ! The total ocean volume and its change [m3] (unscaled to permit reproducing sum). real :: Area ! The total ocean surface area [m2] (unscaled to permit reproducing sum). real :: h_minimum ! The minimum layer thicknesses [H ~> m or kg m-2] @@ -294,6 +294,8 @@ subroutine MOM_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, pe T%average = T%average + dV*Temp(i,j,k) S%minimum = min( S%minimum, Salt(i,j,k) ) ; S%maximum = max( S%maximum, Salt(i,j,k) ) S%average = S%average + dV*Salt(i,j,k) + tmp_T(i,j) = tmp_T(i,j) + dV*Temp(i,j,k) + tmp_S(i,j) = tmp_S(i,j) + dV*Salt(i,j,k) endif if (h_minimum > h(i,j,k)) h_minimum = h(i,j,k) endif diff --git a/src/diagnostics/MOM_spatial_means.F90 b/src/diagnostics/MOM_spatial_means.F90 index 7969ee11f8..7fc83f9b40 100644 --- a/src/diagnostics/MOM_spatial_means.F90 +++ b/src/diagnostics/MOM_spatial_means.F90 @@ -19,7 +19,7 @@ module MOM_spatial_means public :: global_i_mean, global_j_mean public :: global_area_mean, global_area_mean_u, global_area_mean_v, global_layer_mean public :: global_area_integral -public :: global_volume_mean, global_mass_integral +public :: global_volume_mean, global_mass_integral, global_mass_int_EFP public :: adjust_area_mean_to_zero contains @@ -234,6 +234,49 @@ function global_mass_integral(h, G, GV, var, on_PE_only, scale) end function global_mass_integral +!> Find the global mass-weighted order invariant integral of a variable in mks units, +!! returning the value as an EFP_type. This uses reproducing sums. +function global_mass_int_EFP(h, G, GV, var, on_PE_only, scale) + 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_(GV)), & + intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & + optional, intent(in) :: var !< The variable being integrated + logical, optional, intent(in) :: on_PE_only !< If present and true, the sum is only done + !! on the local PE, but it is still order invariant. + real, optional, intent(in) :: scale !< A rescaling factor for the variable + type(EFP_type) :: global_mass_int_EFP !< The mass-weighted integral of var (or 1) in + !! kg times the units of var + + ! Local variables + real, dimension(SZI_(G), SZJ_(G)) :: tmpForSum + real :: scalefac ! An overall scaling factor for the areas and variable. + integer :: i, j, k, is, ie, js, je, nz, isr, ier, jsr, jer + + is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke + isr = is - (G%isd-1) ; ier = ie - (G%isd-1) ; jsr = js - (G%jsd-1) ; jer = je - (G%jsd-1) + + scalefac = GV%H_to_kg_m2 * G%US%L_to_m**2 + if (present(scale)) scalefac = scale * scalefac + + tmpForSum(:,:) = 0.0 + if (present(var)) then + do k=1,nz ; do j=js,je ; do i=is,ie + tmpForSum(i,j) = tmpForSum(i,j) + var(i,j,k) * & + ((scalefac * h(i,j,k)) * (G%areaT(i,j) * G%mask2dT(i,j))) + enddo ; enddo ; enddo + else + do k=1,nz ; do j=js,je ; do i=is,ie + tmpForSum(i,j) = tmpForSum(i,j) + & + ((scalefac * h(i,j,k)) * (G%areaT(i,j) * G%mask2dT(i,j))) + enddo ; enddo ; enddo + endif + + global_mass_int_EFP = reproducing_sum_EFP(tmpForSum, isr, ier, jsr, jer, only_on_PE=on_PE_only) + +end function global_mass_int_EFP + !> Determine the global mean of a field along rows of constant i, returning it !! in a 1-d array using the local indexing. This uses reproducing sums. diff --git a/src/diagnostics/MOM_sum_output.F90 b/src/diagnostics/MOM_sum_output.F90 index 668c297658..a7cae98620 100644 --- a/src/diagnostics/MOM_sum_output.F90 +++ b/src/diagnostics/MOM_sum_output.F90 @@ -733,10 +733,6 @@ subroutine write_energy(u, v, h, tv, day, n, G, GV, US, CS, tracer_CSp, dt_forci enddo ; enddo ; enddo call sum_across_PEs(CS%ntrunc) - ! Sum the various quantities across all the processors. This sum is NOT - ! guaranteed to be bitwise reproducible, even on the same decomposition. - ! The sum of Tr_stocks should be reimplemented using the reproducing sums. - if (nTr_stocks > 0) call sum_across_PEs(Tr_stocks,nTr_stocks) call max_across_PEs(max_CFL, 2) diff --git a/src/tracer/MOM_CFC_cap.F90 b/src/tracer/MOM_CFC_cap.F90 index 7296f1d469..c174fe4c39 100644 --- a/src/tracer/MOM_CFC_cap.F90 +++ b/src/tracer/MOM_CFC_cap.F90 @@ -4,6 +4,7 @@ module MOM_CFC_cap ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_coms, only : EFP_type use MOM_diag_mediator, only : diag_ctrl, register_diag_field, post_data use MOM_error_handler, only : MOM_error, FATAL, WARNING use MOM_file_parser, only : get_param, log_param, log_version, param_file_type @@ -14,6 +15,7 @@ module MOM_CFC_cap use MOM_io, only : vardesc, var_desc, query_vardesc, stdout use MOM_open_boundary, only : ocean_OBC_type use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP use MOM_time_manager, only : time_type use time_interp_external_mod, only : init_external_field, time_interp_external use MOM_tracer_registry, only : register_tracer, tracer_registry_type @@ -341,14 +343,13 @@ end subroutine CFC_cap_column_physics !> Calculates the mass-weighted integral of all tracer stocks, !! returning the number of stocks it has calculated. If the stock_index !! is present, only the stock corresponding to that coded index is returned. -function CFC_cap_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function CFC_cap_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), & intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]. - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each - !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(CFC_cap_CS), pointer :: CS !< The control structure returned by a !! previous call to register_CFC_cap. character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated. @@ -357,11 +358,6 @@ function CFC_cap_stock(h, stocks, G, GV, US, CS, names, units, stock_index) !! stock being sought. integer :: CFC_cap_stock !< The number of stocks calculated here. - ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - real :: mass ! The cell volume or mass [H L2 ~> m3 or kg] - integer :: i, j, k, is, ie, js, je, nz - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke CFC_cap_stock = 0 if (.not.associated(CS)) return @@ -377,15 +373,8 @@ function CFC_cap_stock(h, stocks, G, GV, US, CS, names, units, stock_index) call query_vardesc(CS%CFC12_desc, name=names(2), units=units(2), caller="CFC_cap_stock") units(1) = trim(units(1))//" kg" ; units(2) = trim(units(2))//" kg" - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 - stocks(1) = 0.0 ; stocks(2) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - mass = G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k) - stocks(1) = stocks(1) + CS%CFC11(i,j,k) * mass - stocks(2) = stocks(2) + CS%CFC12(i,j,k) * mass - enddo ; enddo ; enddo - stocks(1) = stock_scale * stocks(1) - stocks(2) = stock_scale * stocks(2) + stocks(1) = global_mass_int_EFP(h, G, GV, CS%CFC11, on_PE_only=.true.) + stocks(2) = global_mass_int_EFP(h, G, GV, CS%CFC12, on_PE_only=.true.) CFC_cap_stock = 2 diff --git a/src/tracer/MOM_OCMIP2_CFC.F90 b/src/tracer/MOM_OCMIP2_CFC.F90 index 5fe55b896b..28a9501d51 100644 --- a/src/tracer/MOM_OCMIP2_CFC.F90 +++ b/src/tracer/MOM_OCMIP2_CFC.F90 @@ -3,25 +3,28 @@ module MOM_OCMIP2_CFC ! This file is part of MOM6. See LICENSE.md for the license. -use MOM_coupler_types, only : extract_coupler_type_data, set_coupler_type_data -use MOM_coupler_types, only : atmos_ocn_coupler_flux -use MOM_diag_mediator, only : diag_ctrl -use MOM_error_handler, only : MOM_error, FATAL, WARNING -use MOM_file_parser, only : get_param, log_param, log_version, param_file_type -use MOM_forcing_type, only : forcing -use MOM_hor_index, only : hor_index_type -use MOM_grid, only : ocean_grid_type -use MOM_io, only : file_exists, MOM_read_data, slasher, vardesc, var_desc, query_vardesc -use MOM_open_boundary, only : ocean_OBC_type -use MOM_restart, only : query_initialized, MOM_restart_CS -use MOM_sponge, only : set_up_sponge_field, sponge_CS -use MOM_time_manager, only : time_type +use MOM_coms, only : EFP_type +use MOM_coupler_types, only : extract_coupler_type_data, set_coupler_type_data +use MOM_coupler_types, only : atmos_ocn_coupler_flux +use MOM_diag_mediator, only : diag_ctrl +use MOM_error_handler, only : MOM_error, FATAL, WARNING +use MOM_file_parser, only : get_param, log_param, log_version, param_file_type +use MOM_forcing_type, only : forcing +use MOM_hor_index, only : hor_index_type +use MOM_grid, only : ocean_grid_type +use MOM_io, only : file_exists, MOM_read_data, slasher +use MOM_io, only : vardesc, var_desc, query_vardesc +use MOM_open_boundary, only : ocean_OBC_type +use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP +use MOM_sponge, only : set_up_sponge_field, sponge_CS +use MOM_time_manager, only : time_type use MOM_tracer_registry, only : register_tracer, tracer_registry_type use MOM_tracer_diabatic, only : tracer_vertdiff, applyTracerBoundaryFluxesInOut -use MOM_tracer_Z_init, only : tracer_Z_init -use MOM_unit_scaling, only : unit_scale_type -use MOM_variables, only : surface -use MOM_verticalGrid, only : verticalGrid_type +use MOM_tracer_Z_init, only : tracer_Z_init +use MOM_unit_scaling, only : unit_scale_type +use MOM_variables, only : surface +use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -478,14 +481,13 @@ end subroutine OCMIP2_CFC_column_physics !> This function calculates the mass-weighted integral of all tracer stocks, !! returning the number of stocks it has calculated. If the stock_index !! is present, only the stock corresponding to that coded index is returned. -function OCMIP2_CFC_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function OCMIP2_CFC_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), & intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]. - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each - !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(OCMIP2_CFC_CS), pointer :: CS !< The control structure returned by a !! previous call to register_OCMIP2_CFC. character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated. @@ -494,11 +496,6 @@ function OCMIP2_CFC_stock(h, stocks, G, GV, US, CS, names, units, stock_index) !! stock being sought. integer :: OCMIP2_CFC_stock !< The number of stocks calculated here. - ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - real :: mass ! The cell volume or mass [H L2 ~> m3 or kg] - integer :: i, j, k, is, ie, js, je, nz - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke OCMIP2_CFC_stock = 0 if (.not.associated(CS)) return @@ -514,15 +511,8 @@ function OCMIP2_CFC_stock(h, stocks, G, GV, US, CS, names, units, stock_index) call query_vardesc(CS%CFC12_desc, name=names(2), units=units(2), caller="OCMIP2_CFC_stock") units(1) = trim(units(1))//" kg" ; units(2) = trim(units(2))//" kg" - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 - stocks(1) = 0.0 ; stocks(2) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - mass = G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k) - stocks(1) = stocks(1) + CS%CFC11(i,j,k) * mass - stocks(2) = stocks(2) + CS%CFC12(i,j,k) * mass - enddo ; enddo ; enddo - stocks(1) = stock_scale * stocks(1) - stocks(2) = stock_scale * stocks(2) + stocks(1) = global_mass_int_EFP(h, G, GV, CS%CFC11, on_PE_only=.true.) + stocks(2) = global_mass_int_EFP(h, G, GV, CS%CFC12, on_PE_only=.true.) OCMIP2_CFC_stock = 2 diff --git a/src/tracer/MOM_generic_tracer.F90 b/src/tracer/MOM_generic_tracer.F90 index f8c0f6ac06..31acb51160 100644 --- a/src/tracer/MOM_generic_tracer.F90 +++ b/src/tracer/MOM_generic_tracer.F90 @@ -29,7 +29,7 @@ module MOM_generic_tracer use g_tracer_utils, only: g_tracer_get_pointer,g_tracer_get_alias,g_tracer_set_csdiag use MOM_ALE_sponge, only : set_up_ALE_sponge_field, ALE_sponge_CS - use MOM_coms, only : max_across_PEs, min_across_PEs, PE_here + use MOM_coms, only : EFP_type, max_across_PEs, min_across_PEs, PE_here use MOM_diag_mediator, only : post_data, register_diag_field, safe_alloc_ptr use MOM_diag_mediator, only : diag_ctrl, get_diag_time_end use MOM_error_handler, only : MOM_error, FATAL, WARNING, NOTE, is_root_pe @@ -40,7 +40,7 @@ module MOM_generic_tracer use MOM_io, only : file_exists, MOM_read_data, slasher use MOM_open_boundary, only : ocean_OBC_type use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS - use MOM_spatial_means, only : global_area_mean + use MOM_spatial_means, only : global_area_mean, global_mass_int_EFP use MOM_sponge, only : set_up_sponge_field, sponge_CS use MOM_time_manager, only : time_type, set_time use MOM_tracer_diabatic, only : tracer_vertdiff, applyTracerBoundaryFluxesInOut @@ -568,13 +568,12 @@ end subroutine MOM_generic_tracer_column_physics !! being requested specifically, returning the number of stocks it has !! calculated. If the stock_index is present, only the stock corresponding !! to that coded index is returned. - function MOM_generic_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_index) + function MOM_generic_tracer_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each - !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(MOM_generic_tracer_CS), pointer :: CS !< Pointer to the control structure for this module. character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated. character(len=*), dimension(:), intent(out) :: units !< The units of the stocks calculated. @@ -584,14 +583,12 @@ function MOM_generic_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_ !! number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] type(g_tracer_type), pointer :: g_tracer, g_tracer_next real, dimension(:,:,:,:), pointer :: tr_field real, dimension(:,:,:), pointer :: tr_ptr character(len=128), parameter :: sub_name = 'MOM_generic_tracer_stock' - integer :: i, j, k, is, ie, js, je, nz, m - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke + integer :: m MOM_generic_tracer_stock = 0 if (.not.associated(CS)) return @@ -605,7 +602,6 @@ function MOM_generic_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_ if (.NOT. associated(CS%g_tracer_list)) return ! No stocks. - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 m=1 ; g_tracer=>CS%g_tracer_list do call g_tracer_get_alias(g_tracer,names(m)) @@ -613,12 +609,8 @@ function MOM_generic_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_ units(m) = trim(units(m))//" kg" call g_tracer_get_pointer(g_tracer,names(m),'field',tr_field) - stocks(m) = 0.0 tr_ptr => tr_field(:,:,:,1) - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + tr_ptr(i,j,k) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, tr_ptr, on_PE_only=.true.) !traverse the linked list till hit NULL call g_tracer_get_next(g_tracer, g_tracer_next) diff --git a/src/tracer/MOM_lateral_boundary_diffusion.F90 b/src/tracer/MOM_lateral_boundary_diffusion.F90 index 4a98aa1934..227e3ffb06 100644 --- a/src/tracer/MOM_lateral_boundary_diffusion.F90 +++ b/src/tracer/MOM_lateral_boundary_diffusion.F90 @@ -8,16 +8,17 @@ module MOM_lateral_boundary_diffusion use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end use MOM_cpu_clock, only : CLOCK_MODULE use MOM_checksums, only : hchksum -use MOM_domains, only : pass_var, sum_across_PEs +use MOM_domains, only : pass_var use MOM_diag_mediator, only : diag_ctrl, time_type use MOM_diag_mediator, only : post_data, register_diag_field use MOM_diag_vkernels, only : reintegrate_column -use MOM_error_handler, only : MOM_error, FATAL, is_root_pe +use MOM_error_handler, only : MOM_error, MOM_mesg, FATAL, is_root_pe use MOM_file_parser, only : get_param, log_version, param_file_type use MOM_grid, only : ocean_grid_type use MOM_remapping, only : remapping_CS, initialize_remapping use MOM_remapping, only : extract_member_remapping_CS, remapping_core_h use MOM_remapping, only : remappingSchemesDoc, remappingDefaultScheme +use MOM_spatial_means, only : global_mass_integral use MOM_tracer_registry, only : tracer_registry_type, tracer_type use MOM_unit_scaling, only : unit_scale_type use MOM_verticalGrid, only : verticalGrid_type @@ -169,13 +170,11 @@ subroutine lateral_boundary_diffusion(G, GV, US, h, Coef_x, Coef_y, dt, Reg, CS) real, dimension(SZK_(GV)) :: tracer_1d !< 1d-array used to remap tracer change to native grid real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: tracer_old !< local copy of the initial tracer concentration, !! only used to compute tendencies. - real, dimension(SZI_(G),SZJ_(G)) :: tracer_int !< integrated tracer before LBD is applied - !! [conc H L2 ~> conc m3 or conc kg] - real, dimension(SZI_(G),SZJ_(G)) :: tracer_end !< integrated tracer after LBD is applied. - !! [conc H L2 ~> conc m3 or conc kg] - integer :: i, j, k, m !< indices to loop over + real :: tracer_int_prev !< Globally integrated tracer before LBD is applied, in mks units [conc kg] + real :: tracer_int_end !< Integrated tracer after LBD is applied, in mks units [conc kg] real :: Idt !< inverse of the time step [T-1 ~> s-1] - real :: tmp1, tmp2 !< temporary variables [conc H L2 ~> conc m3 or conc kg] + character(len=256) :: mesg !< Message for error messages. + integer :: i, j, k, m !< indices to loop over call cpu_clock_begin(id_clock_lbd) Idt = 1./dt @@ -236,22 +235,11 @@ subroutine lateral_boundary_diffusion(G, GV, US, h, Coef_x, Coef_y, dt, Reg, CS) if (CS%debug) then call hchksum(tracer%t, "after LBD "//tracer%name,G%HI) - tracer_int(:,:) = 0.0; tracer_end(:,:) = 0.0 - ! tracer (native grid) before and after LBD - do j=G%jsc,G%jec ; do i=G%isc,G%iec - do k=1,GV%ke - tracer_int(i,j) = tracer_int(i,j) + tracer_old(i,j,k) * & - (h(i,j,k)*(G%mask2dT(i,j)*G%areaT(i,j))) - tracer_end(i,j) = tracer_end(i,j) + tracer%t(i,j,k) * & - (h(i,j,k)*(G%mask2dT(i,j)*G%areaT(i,j))) - enddo - enddo; enddo - - tmp1 = SUM(tracer_int) - tmp2 = SUM(tracer_end) - call sum_across_PEs(tmp1) - call sum_across_PEs(tmp2) - if (is_root_pe()) write(*,*)'Total '//tracer%name//' before/after LBD:', tmp1, tmp2 + ! tracer (native grid) integrated tracer amounts before and after LBD + tracer_int_prev = global_mass_integral(h, G, GV, tracer_old) + tracer_int_end = global_mass_integral(h, G, GV, tracer%t) + write(mesg,*) 'Total '//tracer%name//' before/after LBD:', tracer_int_prev, tracer_int_end + call MOM_mesg(mesg) endif ! Post the tracer diagnostics diff --git a/src/tracer/MOM_tracer_flow_control.F90 b/src/tracer/MOM_tracer_flow_control.F90 index 2ae72a3270..ce747bba01 100644 --- a/src/tracer/MOM_tracer_flow_control.F90 +++ b/src/tracer/MOM_tracer_flow_control.F90 @@ -3,21 +3,22 @@ module MOM_tracer_flow_control ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_coms, only : EFP_type, assignment(=), EFP_to_real, real_to_EFP, EFP_sum_across_PEs use MOM_diag_mediator, only : time_type, diag_ctrl use MOM_error_handler, only : MOM_error, FATAL, WARNING -use MOM_file_parser, only : get_param, log_version, param_file_type, close_param_file -use MOM_forcing_type, only : forcing, optics_type -use MOM_get_input, only : Get_MOM_input -use MOM_grid, only : ocean_grid_type -use MOM_hor_index, only : hor_index_type +use MOM_file_parser, only : get_param, log_version, param_file_type, close_param_file +use MOM_forcing_type, only : forcing, optics_type +use MOM_get_input, only : Get_MOM_input +use MOM_grid, only : ocean_grid_type +use MOM_hor_index, only : hor_index_type use MOM_open_boundary, only : ocean_OBC_type -use MOM_restart, only : MOM_restart_CS -use MOM_sponge, only : sponge_CS -use MOM_ALE_sponge, only : ALE_sponge_CS +use MOM_restart, only : MOM_restart_CS +use MOM_sponge, only : sponge_CS +use MOM_ALE_sponge, only : ALE_sponge_CS use MOM_tracer_registry, only : tracer_registry_type -use MOM_unit_scaling, only : unit_scale_type -use MOM_variables, only : surface, thermo_var_ptrs -use MOM_verticalGrid, only : verticalGrid_type +use MOM_unit_scaling, only : unit_scale_type +use MOM_variables, only : surface, thermo_var_ptrs +use MOM_verticalGrid, only : verticalGrid_type #include ! Add references to other user-provide tracer modules here. @@ -582,8 +583,8 @@ subroutine call_tracer_stocks(h, stock_values, G, GV, US, CS, stock_names, stock type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stock_values !< The integrated amounts of a tracer - !! on the current PE, usually in kg x concentration [kg conc]. + real, dimension(:), intent(out) :: stock_values !< The globally mass-integrated + !! amount of a tracer [kg conc]. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(tracer_flow_control_CS), pointer :: CS !< The control structure returned by a !! previous call to @@ -611,8 +612,10 @@ subroutine call_tracer_stocks(h, stock_values, G, GV, US, CS, stock_names, stock ! Local variables character(len=200), dimension(MAX_FIELDS_) :: names, units character(len=200) :: set_pkg_name - real, dimension(MAX_FIELDS_) :: values - integer :: max_ns, ns_tot, ns, index, pkg, max_pkgs, nn + ! real, dimension(MAX_FIELDS_) :: values + type(EFP_type), dimension(MAX_FIELDS_) :: values_EFP + type(EFP_type), dimension(MAX_FIELDS_) :: stock_val_EFP + integer :: max_ns, ns_tot, ns, index, pkg, max_pkgs, nn, n if (.not. associated(CS)) call MOM_error(FATAL, "call_tracer_stocks: "// & "Module must be initialized via call_tracer_register before it is used.") @@ -625,59 +628,59 @@ subroutine call_tracer_stocks(h, stock_values, G, GV, US, CS, stock_names, stock ! Add other user-provided calls here. if (CS%use_USER_tracer_example) then - ns = USER_tracer_stock(h, values, G, GV, US, CS%USER_tracer_example_CSp, & + ns = USER_tracer_stock(h, values_EFP, G, GV, CS%USER_tracer_example_CSp, & names, units, stock_index) - call store_stocks("tracer_example", ns, names, units, values, index, stock_values, & + call store_stocks("tracer_example", ns, names, units, values_EFP, index, stock_val_EFP, & set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif ! if (CS%use_DOME_tracer) then ! ns = DOME_tracer_stock(h, values, G, GV, CS%DOME_tracer_CSp, & ! names, units, stock_index) -! call store_stocks("DOME_tracer", ns, names, units, values, index, stock_values, & +! do n=1,ns ; values_EFP(n) = real_to_EFP(values(n)) ; enddo +! call store_stocks("DOME_tracer", ns, names, units, values_EFP, index, stock_val_EFP, & ! set_pkg_name, max_ns, ns_tot, stock_names, stock_units) ! endif if (CS%use_ideal_age) then - ns = ideal_age_stock(h, values, G, GV, US, CS%ideal_age_tracer_CSp, & + ns = ideal_age_stock(h, values_EFP, G, GV, CS%ideal_age_tracer_CSp, & names, units, stock_index) - call store_stocks("ideal_age_example", ns, names, units, values, index, & - stock_values, set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + call store_stocks("ideal_age_example", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_regional_dyes) then - ns = dye_stock(h, values, G, GV, US, CS%dye_tracer_CSp, & - names, units, stock_index) - call store_stocks("regional_dyes", ns, names, units, values, index, & - stock_values, set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + ns = dye_stock(h, values_EFP, G, GV, CS%dye_tracer_CSp, names, units, stock_index) + call store_stocks("regional_dyes", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_oil) then - ns = oil_stock(h, values, G, GV, US, CS%oil_tracer_CSp, & - names, units, stock_index) - call store_stocks("oil_tracer", ns, names, units, values, index, & - stock_values, set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + ns = oil_stock(h, values_EFP, G, GV, CS%oil_tracer_CSp, names, units, stock_index) + call store_stocks("oil_tracer", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_OCMIP2_CFC) then - ns = OCMIP2_CFC_stock(h, values, G, GV, US, CS%OCMIP2_CFC_CSp, names, units, stock_index) - call store_stocks("MOM_OCMIP2_CFC", ns, names, units, values, index, stock_values, & - set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + ns = OCMIP2_CFC_stock(h, values_EFP, G, GV, CS%OCMIP2_CFC_CSp, names, units, stock_index) + call store_stocks("MOM_OCMIP2_CFC", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_CFC_cap) then - ns = CFC_cap_stock(h, values, G, GV, US, CS%CFC_cap_CSp, names, units, stock_index) - call store_stocks("MOM_CFC_cap", ns, names, units, values, index, stock_values, & - set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + ns = CFC_cap_stock(h, values_EFP, G, GV, CS%CFC_cap_CSp, names, units, stock_index) + call store_stocks("MOM_CFC_cap", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_advection_test_tracer) then - ns = advection_test_stock( h, values, G, GV, US, CS%advection_test_tracer_CSp, & + ns = advection_test_stock( h, values_EFP, G, GV, CS%advection_test_tracer_CSp, & names, units, stock_index ) - call store_stocks("advection_test_tracer", ns, names, units, values, index, & - stock_values, set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + ! do n=1,ns ; values_EFP(n) = real_to_EFP(values(n)) ; enddo + call store_stocks("advection_test_tracer", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_MOM_generic_tracer) then - ns = MOM_generic_tracer_stock(h, values, G, GV, US, CS%MOM_generic_tracer_CSp, & + ns = MOM_generic_tracer_stock(h, values_EFP, G, GV, CS%MOM_generic_tracer_CSp, & names, units, stock_index) - call store_stocks("MOM_generic_tracer", ns, names, units, values, index, stock_values, & - set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + call store_stocks("MOM_generic_tracer", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) nn=ns_tot-ns+1 nn=MOM_generic_tracer_min_max(nn, got_min_max, global_min, global_max, & xgmin, ygmin, zgmin, xgmax, ygmax, zgmax ,& @@ -685,20 +688,26 @@ subroutine call_tracer_stocks(h, stock_values, G, GV, US, CS, stock_names, stock endif if (CS%use_pseudo_salt_tracer) then - ns = pseudo_salt_stock(h, values, G, GV, US, CS%pseudo_salt_tracer_CSp, & + ns = pseudo_salt_stock(h, values_EFP, G, GV, CS%pseudo_salt_tracer_CSp, & names, units, stock_index) - call store_stocks("pseudo_salt_tracer", ns, names, units, values, index, & - stock_values, set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + call store_stocks("pseudo_salt_tracer", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif if (CS%use_boundary_impulse_tracer) then - ns = boundary_impulse_stock(h, values, G, GV, US, CS%boundary_impulse_tracer_CSp, & + ns = boundary_impulse_stock(h, values_EFP, G, GV, CS%boundary_impulse_tracer_CSp, & names, units, stock_index) - call store_stocks("boundary_impulse_tracer", ns, names, units, values, index, & - stock_values, set_pkg_name, max_ns, ns_tot, stock_names, stock_units) + call store_stocks("boundary_impulse_tracer", ns, names, units, values_EFP, index, stock_val_EFP, & + set_pkg_name, max_ns, ns_tot, stock_names, stock_units) endif - if (ns_tot == 0) stock_values(1) = 0.0 + ! Sum the various quantities across all the processors. + if (ns_tot > 0) then + call EFP_sum_across_PEs(stock_val_EFP, ns_tot) + do n=1,ns_tot ; stock_values(n) = EFP_to_real(stock_val_EFP(n)) ; enddo + else + stock_values(1) = 0.0 + endif if (present(num_stocks)) num_stocks = ns_tot @@ -713,11 +722,13 @@ subroutine store_stocks(pkg_name, ns, names, units, values, index, stock_values, intent(in) :: names !< Diagnostic names to use for each stock. character(len=*), dimension(:), & intent(in) :: units !< Units to use in the metadata for each stock. - real, dimension(:), intent(in) :: values !< The values of the tracer stocks + type(EFP_type), dimension(:), & + intent(in) :: values !< The values of the tracer stocks integer, intent(in) :: index !< The integer stock index from !! stocks_constants_mod of the stock to be returned. If this is !! present and greater than 0, only a single stock can be returned. - real, dimension(:), intent(inout) :: stock_values !< The master list of stock values + type(EFP_type), dimension(:), & + intent(inout) :: stock_values !< The master list of stock values character(len=*), intent(inout) :: set_pkg_name !< The name of the last tracer package whose !! stocks were stored for a specific index. This is !! used to trigger an error if there are redundant stocks. diff --git a/src/tracer/advection_test_tracer.F90 b/src/tracer/advection_test_tracer.F90 index 8fdb525b4a..b37822823a 100644 --- a/src/tracer/advection_test_tracer.F90 +++ b/src/tracer/advection_test_tracer.F90 @@ -3,16 +3,18 @@ module advection_test_tracer ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_coms, only : EFP_type use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux use MOM_diag_mediator, only : diag_ctrl use MOM_error_handler, only : MOM_error, FATAL, WARNING -use MOM_file_parser, only : get_param, log_param, log_version, param_file_type -use MOM_forcing_type, only : forcing +use MOM_file_parser, only : get_param, log_param, log_version, param_file_type +use MOM_forcing_type, only : forcing use MOM_grid, only : ocean_grid_type use MOM_hor_index, only : hor_index_type use MOM_io, only : slasher, vardesc, var_desc, query_vardesc use MOM_open_boundary, only : ocean_OBC_type use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP use MOM_sponge, only : set_up_sponge_field, sponge_CS use MOM_time_manager, only : time_type use MOM_tracer_registry, only : register_tracer, tracer_registry_type @@ -75,8 +77,8 @@ function register_advection_test_tracer(HI, GV, param_file, CS, tr_Reg, restart_ ! Local variables character(len=80) :: name, longname -! This include declares and sets the variable "version". -#include "version_variable.h" + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "advection_test_tracer" ! This module's name. character(len=200) :: inputdir character(len=48) :: flux_units ! The units for tracer fluxes, usually @@ -344,13 +346,12 @@ end subroutine advection_test_tracer_surface_state !> Calculate the mass-weighted integral of all tracer stocks, returning the number of stocks it has calculated. !! If the stock_index is present, only the stock corresponding to that coded index is returned. -function advection_test_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function advection_test_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each + type(EFP_type), dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(advection_test_tracer_CS), pointer :: CS !< The control structure returned by a previous !! call to register_advection_test_tracer. character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated. @@ -359,7 +360,6 @@ function advection_test_stock(h, stocks, G, GV, US, CS, names, units, stock_inde integer :: advection_test_stock !< the number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -374,14 +374,9 @@ function advection_test_stock(h, stocks, G, GV, US, CS, names, units, stock_inde return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 do m=1,CS%ntr call query_vardesc(CS%tr_desc(m), name=names(m), units=units(m), caller="advection_test_stock") - stocks(m) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + CS%tr(i,j,k,m) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, CS%tr(:,:,:,m), on_PE_only=.true.) enddo advection_test_stock = CS%ntr diff --git a/src/tracer/boundary_impulse_tracer.F90 b/src/tracer/boundary_impulse_tracer.F90 index ea60a09608..44423b5650 100644 --- a/src/tracer/boundary_impulse_tracer.F90 +++ b/src/tracer/boundary_impulse_tracer.F90 @@ -3,24 +3,26 @@ module boundary_impulse_tracer ! This file is part of MOM6. See LICENSE.md for the license. -use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux -use MOM_diag_mediator, only : diag_ctrl -use MOM_error_handler, only : MOM_error, FATAL, WARNING -use MOM_file_parser, only : get_param, log_param, log_version, param_file_type -use MOM_forcing_type, only : forcing -use MOM_grid, only : ocean_grid_type -use MOM_hor_index, only : hor_index_type -use MOM_io, only : vardesc, var_desc, query_vardesc -use MOM_open_boundary, only : ocean_OBC_type -use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS -use MOM_sponge, only : set_up_sponge_field, sponge_CS -use MOM_time_manager, only : time_type +use MOM_coms, only : EFP_type +use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux +use MOM_diag_mediator, only : diag_ctrl +use MOM_error_handler, only : MOM_error, FATAL, WARNING +use MOM_file_parser, only : get_param, log_param, log_version, param_file_type +use MOM_forcing_type, only : forcing +use MOM_grid, only : ocean_grid_type +use MOM_hor_index, only : hor_index_type +use MOM_io, only : vardesc, var_desc, query_vardesc +use MOM_open_boundary, only : ocean_OBC_type +use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP +use MOM_sponge, only : set_up_sponge_field, sponge_CS +use MOM_time_manager, only : time_type use MOM_tracer_registry, only : register_tracer, tracer_registry_type use MOM_tracer_diabatic, only : tracer_vertdiff, applyTracerBoundaryFluxesInOut -use MOM_tracer_Z_init, only : tracer_Z_init -use MOM_unit_scaling, only : unit_scale_type -use MOM_variables, only : surface, thermo_var_ptrs -use MOM_verticalGrid, only : verticalGrid_type +use MOM_tracer_Z_init, only : tracer_Z_init +use MOM_unit_scaling, only : unit_scale_type +use MOM_variables, only : surface, thermo_var_ptrs +use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -287,13 +289,12 @@ end subroutine boundary_impulse_tracer_column_physics !> This function calculates the mass-weighted integral of the boundary impulse, !! tracer stocks returning the number of stocks it has calculated. If the stock_index !! is present, only the stock corresponding to that coded index is returned. -function boundary_impulse_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function boundary_impulse_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), intent(in ) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent( out) :: stocks !< the mass-weighted integrated amount of each - !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in ) :: US !< A dimensional unit scaling type + type(EFP_type), dimension(:), intent( out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(boundary_impulse_tracer_CS), pointer :: CS !< The control structure returned by a previous !! call to register_boundary_impulse_tracer. character(len=*), dimension(:), intent( out) :: names !< The names of the stocks calculated. @@ -302,14 +303,8 @@ function boundary_impulse_stock(h, stocks, G, GV, US, CS, names, units, stock_in !! being sought. integer :: boundary_impulse_stock !< Return value: the number of stocks calculated here. -! This function calculates the mass-weighted integral of all tracer stocks, -! returning the number of stocks it has calculated. If the stock_index -! is present, only the stock corresponding to that coded index is returned. - ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - integer :: i, j, k, is, ie, js, je, nz, m - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke + integer :: m boundary_impulse_stock = 0 if (.not.associated(CS)) return @@ -322,15 +317,10 @@ function boundary_impulse_stock(h, stocks, G, GV, US, CS, names, units, stock_in return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 do m=1,1 call query_vardesc(CS%tr_desc(m), name=names(m), units=units(m), caller="boundary_impulse_stock") units(m) = trim(units(m))//" kg" - stocks(m) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + CS%tr(i,j,k,m) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, CS%tr(:,:,:,m), on_PE_only=.true.) enddo boundary_impulse_stock = CS%ntr diff --git a/src/tracer/dye_example.F90 b/src/tracer/dye_example.F90 index dca01e974a..d7c7a7bad3 100644 --- a/src/tracer/dye_example.F90 +++ b/src/tracer/dye_example.F90 @@ -3,6 +3,7 @@ module regional_dyes ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_coms, only : EFP_type use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux use MOM_diag_mediator, only : diag_ctrl use MOM_error_handler, only : MOM_error, FATAL, WARNING @@ -13,6 +14,7 @@ module regional_dyes use MOM_io, only : vardesc, var_desc, query_vardesc use MOM_open_boundary, only : ocean_OBC_type use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP use MOM_sponge, only : set_up_sponge_field, sponge_CS use MOM_time_manager, only : time_type use MOM_tracer_registry, only : register_tracer, tracer_registry_type @@ -74,13 +76,13 @@ function register_dye_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS) !! structure for the tracer advection and diffusion module. type(MOM_restart_CS), target, intent(inout) :: restart_CS !< MOM restart control struct -! Local variables -! This include declares and sets the variable "version". -#include "version_variable.h" + ! Local variables character(len=40) :: mdl = "regional_dyes" ! This module's name. character(len=200) :: inputdir ! The directory where the input files are. character(len=48) :: var_name ! The variable's name. character(len=48) :: desc_name ! The variable's descriptor. + ! This include declares and sets the variable "version". +# include "version_variable.h" real, pointer :: tr_ptr(:,:,:) => NULL() logical :: register_dye_tracer integer :: isd, ied, jsd, jed, nz, m @@ -325,13 +327,12 @@ end subroutine dye_tracer_column_physics !> This function calculates the mass-weighted integral of all tracer stocks, !! returning the number of stocks it has calculated. If the stock_index !! is present, only the stock corresponding to that coded index is returned. -function dye_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function dye_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of - !! each tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(dye_tracer_CS), pointer :: CS !< The control structure returned by a !! previous call to register_dye_tracer. character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated. @@ -342,9 +343,7 @@ function dye_stock(h, stocks, G, GV, US, CS, names, units, stock_index) !! calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - integer :: i, j, k, is, ie, js, je, nz, m - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke + integer :: m dye_stock = 0 if (.not.associated(CS)) return @@ -357,15 +356,10 @@ function dye_stock(h, stocks, G, GV, US, CS, names, units, stock_index) return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 do m=1,CS%ntr call query_vardesc(CS%tr_desc(m), name=names(m), units=units(m), caller="dye_stock") units(m) = trim(units(m))//" kg" - stocks(m) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + CS%tr(i,j,k,m) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, CS%tr(:,:,:,m), on_PE_only=.true.) enddo dye_stock = CS%ntr diff --git a/src/tracer/ideal_age_example.F90 b/src/tracer/ideal_age_example.F90 index d5c813b3d0..5913251b14 100644 --- a/src/tracer/ideal_age_example.F90 +++ b/src/tracer/ideal_age_example.F90 @@ -3,6 +3,7 @@ module ideal_age_example ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_coms, only : EFP_type use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux use MOM_diag_mediator, only : diag_ctrl use MOM_error_handler, only : MOM_error, FATAL, WARNING @@ -13,6 +14,7 @@ module ideal_age_example use MOM_io, only : file_exists, MOM_read_data, slasher, vardesc, var_desc, query_vardesc use MOM_open_boundary, only : ocean_OBC_type use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP use MOM_sponge, only : set_up_sponge_field, sponge_CS use MOM_time_manager, only : time_type, time_type_to_real use MOM_tracer_registry, only : register_tracer, tracer_registry_type @@ -78,8 +80,8 @@ function register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS) !! diffusion module type(MOM_restart_CS), target, intent(inout) :: restart_CS !< MOM restart control struct -! This include declares and sets the variable "version". -#include "version_variable.h" + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "ideal_age_example" ! This module's name. character(len=200) :: inputdir ! The directory where the input files are. character(len=48) :: var_name ! The variable's name. @@ -369,14 +371,13 @@ end subroutine ideal_age_tracer_column_physics !> Calculates the mass-weighted integral of all tracer stocks, returning the number of stocks it !! has calculated. If stock_index is present, only the stock corresponding to that coded index is found. -function ideal_age_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function ideal_age_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), & intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each + type(EFP_type), dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(ideal_age_tracer_CS), pointer :: CS !< The control structure returned by a previous !! call to register_ideal_age_tracer. character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated. @@ -386,7 +387,6 @@ function ideal_age_stock(h, stocks, G, GV, US, CS, names, units, stock_index) integer :: ideal_age_stock !< The number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -401,15 +401,10 @@ function ideal_age_stock(h, stocks, G, GV, US, CS, names, units, stock_index) return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 do m=1,CS%ntr call query_vardesc(CS%tr_desc(m), name=names(m), units=units(m), caller="ideal_age_stock") units(m) = trim(units(m))//" kg" - stocks(m) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + CS%tr(i,j,k,m) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, CS%tr(:,:,:,m), on_PE_only=.true.) enddo ideal_age_stock = CS%ntr diff --git a/src/tracer/oil_tracer.F90 b/src/tracer/oil_tracer.F90 index 6f690ab760..0c5a4e6e8d 100644 --- a/src/tracer/oil_tracer.F90 +++ b/src/tracer/oil_tracer.F90 @@ -3,24 +3,27 @@ module oil_tracer ! This file is part of MOM6. See LICENSE.md for the license. -use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux -use MOM_diag_mediator, only : diag_ctrl -use MOM_error_handler, only : MOM_error, FATAL, WARNING -use MOM_file_parser, only : get_param, log_param, log_version, param_file_type -use MOM_forcing_type, only : forcing -use MOM_grid, only : ocean_grid_type -use MOM_hor_index, only : hor_index_type -use MOM_io, only : file_exists, MOM_read_data, slasher, vardesc, var_desc, query_vardesc -use MOM_open_boundary, only : ocean_OBC_type -use MOM_restart, only : query_initialized, MOM_restart_CS -use MOM_sponge, only : set_up_sponge_field, sponge_CS -use MOM_time_manager, only : time_type, time_type_to_real +use MOM_coms, only : EFP_type +use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux +use MOM_diag_mediator, only : diag_ctrl +use MOM_error_handler, only : MOM_error, FATAL, WARNING +use MOM_file_parser, only : get_param, log_param, log_version, param_file_type +use MOM_forcing_type, only : forcing +use MOM_grid, only : ocean_grid_type +use MOM_hor_index, only : hor_index_type +use MOM_io, only : file_exists, MOM_read_data, slasher +use MOM_io, only : vardesc, var_desc, query_vardesc +use MOM_open_boundary, only : ocean_OBC_type +use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP +use MOM_sponge, only : set_up_sponge_field, sponge_CS +use MOM_time_manager, only : time_type, time_type_to_real use MOM_tracer_registry, only : register_tracer, tracer_registry_type use MOM_tracer_diabatic, only : tracer_vertdiff, applyTracerBoundaryFluxesInOut -use MOM_tracer_Z_init, only : tracer_Z_init -use MOM_unit_scaling, only : unit_scale_type -use MOM_variables, only : surface, thermo_var_ptrs -use MOM_verticalGrid, only : verticalGrid_type +use MOM_tracer_Z_init, only : tracer_Z_init +use MOM_unit_scaling, only : unit_scale_type +use MOM_variables, only : surface, thermo_var_ptrs +use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -81,7 +84,7 @@ function register_oil_tracer(HI, GV, US, param_file, CS, tr_Reg, restart_CS) ! Local variables character(len=40) :: mdl = "oil_tracer" ! This module's name. -! This include declares and sets the variable "version". + ! This include declares and sets the variable "version". # include "version_variable.h" real, dimension(NTR_MAX) :: oil_decay_days !< Decay time scale of oil [days] character(len=200) :: inputdir ! The directory where the input files are. @@ -402,13 +405,12 @@ end subroutine oil_tracer_column_physics !> Calculate the mass-weighted integral of the oil tracer stocks, returning the number of stocks it !! has calculated. If the stock_index is present, only the stock corresponding to that coded index is returned. -function oil_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function oil_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each - !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(oil_tracer_CS), pointer :: CS !< The control structure returned by a previous !! call to register_oil_tracer. character(len=*), dimension(:), intent(out) :: names !< the names of the stocks calculated. @@ -418,9 +420,7 @@ function oil_stock(h, stocks, G, GV, US, CS, names, units, stock_index) integer :: oil_stock !< The number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - integer :: i, j, k, is, ie, js, je, nz, m - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke + integer :: m oil_stock = 0 if (.not.associated(CS)) return @@ -433,15 +433,10 @@ function oil_stock(h, stocks, G, GV, US, CS, names, units, stock_index) return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 do m=1,CS%ntr call query_vardesc(CS%tr_desc(m), name=names(m), units=units(m), caller="oil_stock") units(m) = trim(units(m))//" kg" - stocks(m) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + CS%tr(i,j,k,m) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, CS%tr(:,:,:,m), on_PE_only=.true.) enddo oil_stock = CS%ntr diff --git a/src/tracer/pseudo_salt_tracer.F90 b/src/tracer/pseudo_salt_tracer.F90 index c441e519be..6c22daa150 100644 --- a/src/tracer/pseudo_salt_tracer.F90 +++ b/src/tracer/pseudo_salt_tracer.F90 @@ -3,6 +3,7 @@ module pseudo_salt_tracer ! This file is part of MOM6. See LICENSE.md for the license. +use MOM_coms, only : EFP_type use MOM_debugging, only : hchksum use MOM_diag_mediator, only : post_data, register_diag_field, safe_alloc_ptr use MOM_diag_mediator, only : diag_ctrl @@ -14,6 +15,7 @@ module pseudo_salt_tracer use MOM_io, only : vardesc, var_desc, query_vardesc use MOM_open_boundary, only : ocean_OBC_type use MOM_restart, only : query_initialized, MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP use MOM_sponge, only : set_up_sponge_field, sponge_CS use MOM_time_manager, only : time_type use MOM_tracer_registry, only : register_tracer, tracer_registry_type @@ -253,13 +255,12 @@ end subroutine pseudo_salt_tracer_column_physics !> Calculates the mass-weighted integral of all tracer stocks, returning the number of stocks it has !! calculated. If the stock_index is present, only the stock corresponding to that coded index is returned. -function pseudo_salt_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function pseudo_salt_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each !! tracer, in kg times concentration units [kg conc] - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(pseudo_salt_tracer_CS), pointer :: CS !< The control structure returned by a previous !! call to register_pseudo_salt_tracer character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated @@ -269,10 +270,6 @@ function pseudo_salt_stock(h, stocks, G, GV, US, CS, names, units, stock_index) integer :: pseudo_salt_stock !< Return value: the number of !! stocks calculated here - ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - integer :: i, j, k, is, ie, js, je, nz - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke pseudo_salt_stock = 0 if (.not.associated(CS)) return @@ -285,14 +282,9 @@ function pseudo_salt_stock(h, stocks, G, GV, US, CS, names, units, stock_index) return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 call query_vardesc(CS%tr_desc, name=names(1), units=units(1), caller="pseudo_salt_stock") units(1) = trim(units(1))//" kg" - stocks(1) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(1) = stocks(1) + CS%diff(i,j,k) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(1) = stock_scale * stocks(1) + stocks(1) = global_mass_int_EFP(h, G, GV, CS%diff, on_PE_only=.true.) pseudo_salt_stock = 1 diff --git a/src/tracer/tracer_example.F90 b/src/tracer/tracer_example.F90 index a41f0ab76d..3848b84eff 100644 --- a/src/tracer/tracer_example.F90 +++ b/src/tracer/tracer_example.F90 @@ -3,22 +3,25 @@ module USER_tracer_example ! This file is part of MOM6. See LICENSE.md for the license. -use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux -use MOM_diag_mediator, only : diag_ctrl -use MOM_error_handler, only : MOM_error, FATAL, WARNING -use MOM_file_parser, only : get_param, log_param, log_version, param_file_type -use MOM_forcing_type, only : forcing -use MOM_grid, only : ocean_grid_type -use MOM_hor_index, only : hor_index_type -use MOM_io, only : file_exists, MOM_read_data, slasher, vardesc, var_desc, query_vardesc -use MOM_open_boundary, only : ocean_OBC_type -use MOM_restart, only : MOM_restart_CS -use MOM_sponge, only : set_up_sponge_field, sponge_CS -use MOM_time_manager, only : time_type +use MOM_coms, only : EFP_type +use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux +use MOM_diag_mediator, only : diag_ctrl +use MOM_error_handler, only : MOM_error, FATAL, WARNING +use MOM_file_parser, only : get_param, log_param, log_version, param_file_type +use MOM_forcing_type, only : forcing +use MOM_grid, only : ocean_grid_type +use MOM_hor_index, only : hor_index_type +use MOM_io, only : file_exists, MOM_read_data, slasher +use MOM_io, only : vardesc, var_desc, query_vardesc +use MOM_open_boundary, only : ocean_OBC_type +use MOM_restart, only : MOM_restart_CS +use MOM_spatial_means, only : global_mass_int_EFP +use MOM_sponge, only : set_up_sponge_field, sponge_CS +use MOM_time_manager, only : time_type use MOM_tracer_registry, only : register_tracer, tracer_registry_type -use MOM_unit_scaling, only : unit_scale_type -use MOM_variables, only : surface -use MOM_verticalGrid, only : verticalGrid_type +use MOM_unit_scaling, only : unit_scale_type +use MOM_variables, only : surface +use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -64,8 +67,8 @@ function USER_register_tracer_example(HI, GV, param_file, CS, tr_Reg, restart_CS ! Local variables character(len=80) :: name, longname -! This include declares and sets the variable "version". -#include "version_variable.h" + ! This include declares and sets the variable "version". +# include "version_variable.h" character(len=40) :: mdl = "tracer_example" ! This module's name. character(len=200) :: inputdir character(len=48) :: flux_units ! The units for tracer fluxes, usually @@ -358,14 +361,13 @@ end subroutine tracer_column_physics !> This function calculates the mass-weighted integral of all tracer stocks, !! returning the number of stocks it has calculated. If the stock_index !! is present, only the stock corresponding to that coded index is returned. -function USER_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_index) +function USER_tracer_stock(h, stocks, G, GV, CS, names, units, stock_index) 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_(GV)), & intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2] - real, dimension(:), intent(out) :: stocks !< the mass-weighted integrated amount of each - !! tracer, in kg times concentration units [kg conc]. - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(EFP_type), dimension(:), intent(out) :: stocks !< The mass-weighted integrated amount of each + !! tracer, in kg times concentration units [kg conc] type(USER_tracer_example_CS), pointer :: CS !< The control structure returned by a !! previous call to register_USER_tracer. character(len=*), dimension(:), intent(out) :: names !< The names of the stocks calculated. @@ -376,9 +378,7 @@ function USER_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_index) !! stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] - integer :: i, j, k, is, ie, js, je, nz, m - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke + integer :: m USER_tracer_stock = 0 if (.not.associated(CS)) return @@ -390,15 +390,10 @@ function USER_tracer_stock(h, stocks, G, GV, US, CS, names, units, stock_index) return endif ; endif - stock_scale = US%L_to_m**2 * GV%H_to_kg_m2 do m=1,NTR call query_vardesc(CS%tr_desc(m), name=names(m), units=units(m), caller="USER_tracer_stock") units(m) = trim(units(m))//" kg" - stocks(m) = 0.0 - do k=1,nz ; do j=js,je ; do i=is,ie - stocks(m) = stocks(m) + CS%tr(i,j,k,m) * (G%mask2dT(i,j) * G%areaT(i,j) * h(i,j,k)) - enddo ; enddo ; enddo - stocks(m) = stock_scale * stocks(m) + stocks(m) = global_mass_int_EFP(h, G, GV, CS%tr(:,:,:,m), on_PE_only=.true.) enddo USER_tracer_stock = NTR