Skip to content

Commit

Permalink
Use reproducing_sum with unscale in 5 places
Browse files Browse the repository at this point in the history 
  Use reproducing_sum with unscale to calculate the global mass diagnostic,
globally integrated ocean surface area, offline tracer transport residuals and
in calculating the OBC inflow area in tidal_bay_set_OBC_data().  This change allows
for the elimination or replacement of 6 rescaling factors and one added instance
with a consistent conversion factor and diagnostic units occur on the same line.
All answers are bitwise identical.
  • Loading branch information
@Hallberg-NOAA
Hallberg-NOAA committed Dec 6, 2024
1 parent 1c0d3d1 commit d6f594e
Show file tree
Hide file tree
Showing 4 changed files with 21 additions and 31 deletions.
9 changes: 5 additions & 4 deletions src/diagnostics/MOM_diagnostics.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -332,9 +332,9 @@ subroutine calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, &
if (CS%id_masso > 0) then
mass_cell(:,:) = 0.0
do k=1,nz ; do j=js,je ; do i=is,ie
mass_cell(i,j) = mass_cell(i,j) + (GV%H_to_kg_m2*h(i,j,k)) * US%L_to_m**2*G%areaT(i,j)
mass_cell(i,j) = mass_cell(i,j) + (GV%H_to_RZ*h(i,j,k)) * G%areaT(i,j)
enddo ; enddo ; enddo
masso = reproducing_sum(mass_cell)
masso = reproducing_sum(mass_cell, unscale=US%RZL2_to_kg)
call post_data(CS%id_masso, masso, CS%diag)
endif

Expand Down Expand Up @@ -1644,8 +1644,9 @@ subroutine MOM_diagnostics_init(MIS, ADp, CDp, Time, G, GV, US, param_file, diag
Time, 'Mass per unit area of liquid ocean grid cell', 'kg m-2', conversion=GV%H_to_kg_m2, &
standard_name='sea_water_mass_per_unit_area', v_extensive=.true.)

CS%id_masso = register_scalar_field('ocean_model', 'masso', Time, &
diag, 'Mass of liquid ocean', 'kg', standard_name='sea_water_mass')
CS%id_masso = register_scalar_field('ocean_model', 'masso', Time, diag, &
'Mass of liquid ocean', units='kg', conversion=US%RZL2_to_kg, &
standard_name='sea_water_mass')

CS%id_thkcello = register_diag_field('ocean_model', 'thkcello', diag%axesTL, Time, &
long_name='Cell Thickness', standard_name='cell_thickness', &
Expand Down
9 changes: 3 additions & 6 deletions src/initialization/MOM_shared_initialization.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1314,17 +1314,14 @@ subroutine compute_global_grid_integrals(G, US)

! Local variables
real, dimension(G%isc:G%iec, G%jsc:G%jec) :: tmpForSumming ! Masked and unscaled cell areas [m2]
real :: area_scale ! A scaling factor for area into MKS units [m2 L-2 ~> 1]
integer :: i,j

area_scale = US%L_to_m**2
integer :: i, j

tmpForSumming(:,:) = 0.
G%areaT_global = 0.0 ; G%IareaT_global = 0.0
do j=G%jsc,G%jec ; do i=G%isc,G%iec
tmpForSumming(i,j) = area_scale*G%areaT(i,j) * G%mask2dT(i,j)
tmpForSumming(i,j) = G%areaT(i,j) * G%mask2dT(i,j)
enddo ; enddo
G%areaT_global = reproducing_sum(tmpForSumming)
G%areaT_global = US%L_to_m**2 * reproducing_sum(tmpForSumming, unscale=US%L_to_m**2)

if (G%areaT_global == 0.0) &
call MOM_error(FATAL, "compute_global_grid_integrals: zero ocean area (check topography?)")
Expand Down
24 changes: 10 additions & 14 deletions src/tracer/MOM_offline_main.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -625,27 +625,24 @@ real function remaining_transport_sum(G, GV, US, uhtr, vhtr, h_new)

! Local variables
real, dimension(SZI_(G),SZJ_(G)) :: trans_rem_col !< The vertical sum of the absolute value of
!! transports through the faces of a column, in MKS units [kg].
!! transports through the faces of a column [R Z L2 ~> kg].
real :: trans_cell !< The sum of the absolute value of the remaining transports through the faces
!! of a tracer cell [H L2 ~> m3 or kg]
real :: HL2_to_kg_scale !< Unit conversion factor to cell mass [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

HL2_to_kg_scale = GV%H_to_kg_m2 * US%L_to_m**2

trans_rem_col(:,:) = 0.0
do k=1,nz ; do j=js,je ; do i=is,ie
trans_cell = (ABS(uhtr(I-1,j,k)) + ABS(uhtr(I,j,k))) + &
(ABS(vhtr(i,J-1,k)) + ABS(vhtr(i,J,k)))
if (trans_cell > max(1.0e-16*h_new(i,j,k), GV%H_subroundoff) * G%areaT(i,j)) &
trans_rem_col(i,j) = trans_rem_col(i,j) + HL2_to_kg_scale * trans_cell
trans_rem_col(i,j) = trans_rem_col(i,j) + GV%H_to_RZ * trans_cell
enddo ; enddo ; enddo

! The factor of 0.5 here is to avoid double-counting because two cells share a face.
remaining_transport_sum = 0.5 * GV%kg_m2_to_H*US%m_to_L**2 * &
reproducing_sum(trans_rem_col, is+(1-G%isd), ie+(1-G%isd), js+(1-G%jsd), je+(1-G%jsd))
remaining_transport_sum = 0.5 * GV%RZ_to_H * reproducing_sum(trans_rem_col, &
is+(1-G%isd), ie+(1-G%isd), js+(1-G%jsd), je+(1-G%jsd), unscale=US%RZL2_to_kg)

end function remaining_transport_sum

Expand Down Expand Up @@ -876,8 +873,8 @@ subroutine offline_advection_layer(fluxes, Time_start, time_interval, G, GV, US,
real, dimension(SZI_(G),SZJB_(G),SZK_(GV)) :: vhtr_sub ! Remaining meridional mass transports [H L2 ~> m3 or kg]

real, dimension(SZI_(G),SZJB_(G)) :: rem_col_flux ! The summed absolute value of the remaining
! fluxes through the faces of a column or within a column, in mks units [kg]
real :: sum_flux ! Globally summed absolute value of fluxes in mks units [kg], which is
! mass fluxes through the faces of a column or within a column [R Z L2 ~> kg]
real :: sum_flux ! Globally summed absolute value of fluxes [R Z L2 ~> kg], which is
! used to keep track of how close to convergence we are.

real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: &
Expand All @@ -890,7 +887,6 @@ subroutine offline_advection_layer(fluxes, Time_start, time_interval, G, GV, US,
! Work arrays for temperature and salinity
integer :: iter
real :: dt_iter ! The timestep of each iteration [T ~> s]
real :: HL2_to_kg_scale ! Unit conversion factors to cell mass [kg H-1 L-2 ~> kg m-3 or 1]
character(len=160) :: mesg ! The text of an error message
integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
integer :: IsdB, IedB, JsdB, JedB
Expand Down Expand Up @@ -993,22 +989,22 @@ subroutine offline_advection_layer(fluxes, Time_start, time_interval, G, GV, US,
call pass_vector(uhtr,vhtr,G%Domain)

! Calculate how close we are to converging by summing the remaining fluxes at each point
HL2_to_kg_scale = US%L_to_m**2*GV%H_to_kg_m2
rem_col_flux(:,:) = 0.0
do k=1,nz ; do j=js,je ; do i=is,ie
rem_col_flux(i,j) = rem_col_flux(i,j) + HL2_to_kg_scale * &
rem_col_flux(i,j) = rem_col_flux(i,j) + GV%H_to_RZ * &
( (abs(eatr(i,j,k)) + abs(ebtr(i,j,k))) + &
((abs(uhtr(I-1,j,k)) + abs(uhtr(I,j,k))) + &
(abs(vhtr(i,J-1,k)) + abs(vhtr(i,J,k))) ) )
enddo ; enddo ; enddo
sum_flux = reproducing_sum(rem_col_flux, is+(1-G%isd), ie+(1-G%isd), js+(1-G%jsd), je+(1-G%jsd))
sum_flux = reproducing_sum(rem_col_flux, is+(1-G%isd), ie+(1-G%isd), js+(1-G%jsd), je+(1-G%jsd), &
unscale=US%RZL2_to_kg)

if (sum_flux==0) then
write(mesg,*) 'offline_advection_layer: Converged after iteration', iter
call MOM_mesg(mesg)
exit
else
write(mesg,*) "offline_advection_layer: Iteration ", iter, " remaining total fluxes: ", sum_flux
write(mesg,*) "offline_advection_layer: Iteration ", iter, " remaining total fluxes: ", sum_flux*US%RZL2_to_kg
call MOM_mesg(mesg)
endif

Expand Down
10 changes: 3 additions & 7 deletions src/user/tidal_bay_initialization.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -78,8 +78,7 @@ subroutine tidal_bay_set_OBC_data(OBC, CS, G, GV, US, h, Time)
real :: my_flux ! The vlume flux through the face [L2 Z T-1 ~> m3 s-1]
real :: total_area ! The total face area of the OBCs [L Z ~> m2]
real :: PI ! The ratio of the circumference of a circle to its diameter [nondim]
real :: flux_scale ! A scaling factor for the areas [m2 H-1 L-1 ~> nondim or m3 kg-1]
real, allocatable :: my_area(:,:) ! The total OBC inflow area [m2]
real, allocatable :: my_area(:,:) ! The total OBC inflow area [L Z ~> m2]
integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, n
integer :: IsdB, IedB, JsdB, JedB
type(OBC_segment_type), pointer :: segment => NULL()
Expand All @@ -94,8 +93,6 @@ subroutine tidal_bay_set_OBC_data(OBC, CS, G, GV, US, h, Time)

allocate(my_area(1:1,js:je))

flux_scale = GV%H_to_m*US%L_to_m

time_sec = US%s_to_T*time_type_to_real(Time)
cff_eta = CS%tide_ssh_amp * sin(2.0*PI*time_sec / CS%tide_period)
my_area = 0.0
Expand All @@ -105,12 +102,11 @@ subroutine tidal_bay_set_OBC_data(OBC, CS, G, GV, US, h, Time)
do j=segment%HI%jsc,segment%HI%jec ; do I=segment%HI%IscB,segment%HI%IecB
if (OBC%segnum_u(I,j) /= OBC_NONE) then
do k=1,nz
! This area has to be in MKS units to work with reproducing_sum.
my_area(1,j) = my_area(1,j) + h(I,j,k)*flux_scale*G%dyCu(I,j)
my_area(1,j) = my_area(1,j) + h(I,j,k)*(GV%H_to_m*US%m_to_Z)*G%dyCu(I,j)
enddo
endif
enddo ; enddo
total_area = US%m_to_Z*US%m_to_L * reproducing_sum(my_area)
total_area = reproducing_sum(my_area, unscale=US%Z_to_m*US%L_to_m)
my_flux = - CS%tide_flow * SIN(2.0*PI*time_sec / CS%tide_period)

do n = 1, OBC%number_of_segments
Expand Down

0 comments on commit d6f594e

Please sign in to comment.