Skip to content

Commit

Permalink
(*)Avoiding using RHO_0 in non-Boussinesq debugging
Browse files Browse the repository at this point in the history
  When in non-Boussinesq mode, write out mass-based tracer inventories from
calls to MOM_tracer_chkinv and mean thicknesses in mass per unit area from calls
to MOM_state_stats, rather than the approximate volumes that depend on the
Boussinesq reference density.  Similarly the thicknesses and transports output
by write_u_accel and write_v_accel are in mass per unit area or mass flux per
unit face length in non-Boussinesq mode.  This is done specifically by using
GV%H_to_MKS in place of GV%H_to_m; these are identical in Boussinesq mode, but
GV%H_to_MKS avoids rescaling by the Boussinesq reference density in
non-Boussinesq mode.  Several comments were also updated to reflect these
changes.  All solutions are identical, but there are changes in the units of
several diagnostic output fields that the model can write out when debugging
non-Boussinesq mode runs.
  • Loading branch information
Hallberg-NOAA committed Mar 21, 2024
1 parent 4fd0162 commit 8210921
Show file tree
Hide file tree
Showing 3 changed files with 29 additions and 22 deletions.
13 changes: 7 additions & 6 deletions src/core/MOM_checksum_packages.F90
Original file line number Diff line number Diff line change
Expand Up @@ -268,10 +268,11 @@ subroutine MOM_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, pe
! Local variables
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] (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).
tmp_V, & ! The column-integrated volume [m3] or mass [kg] (unscaled to permit reproducing sum),
! depending on whether the Boussinesq approximation is used
tmp_T, & ! The column-integrated temperature [degC m3] or [degC kg] (unscaled to permit reproducing sum)
tmp_S ! The column-integrated salinity [ppt m3] or [ppt kg] (unscaled to permit reproducing sum)
real :: Vol, dV ! The total ocean volume or mass and its change [m3] or [kg] (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]
real :: T_scale ! The scaling conversion factor for temperatures [degC C-1 ~> 1]
Expand All @@ -284,7 +285,7 @@ subroutine MOM_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, pe
! assumption we will not turn this on with threads
type(stats), save :: oldT, oldS
logical, save :: firstCall = .true.
real, save :: oldVol ! The previous total ocean volume [m3]
real, save :: oldVol ! The previous total ocean volume [m3] or mass [kg]

character(len=80) :: lMsg
integer :: is, ie, js, je, nz, i, j, k
Expand All @@ -308,7 +309,7 @@ subroutine MOM_state_stats(mesg, u, v, h, Temp, Salt, G, GV, US, allowChange, pe
h_minimum = 1.E34*GV%m_to_H
do k=1,nz ; do j=js,je ; do i=is,ie
if (G%mask2dT(i,j)>0.) then
dV = US%L_to_m**2*G%areaT(i,j)*GV%H_to_m*h(i,j,k)
dV = US%L_to_m**2*G%areaT(i,j)*GV%H_to_MKS*h(i,j,k)
tmp_V(i,j) = tmp_V(i,j) + dV
if (do_TS .and. h(i,j,k)>0.) then
T%minimum = min( T%minimum, T_scale*Temp(i,j,k) ) ; T%maximum = max( T%maximum, T_scale*Temp(i,j,k) )
Expand Down
18 changes: 10 additions & 8 deletions src/diagnostics/MOM_PointAccel.F90
Original file line number Diff line number Diff line change
Expand Up @@ -95,9 +95,10 @@ subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st
real :: du ! A velocity change [L T-1 ~> m s-1]
real :: Inorm(SZK_(GV)) ! The inverse of the normalized velocity change [L T-1 ~> m s-1]
real :: e(SZK_(GV)+1) ! Simple estimates of interface heights based on the sum of thicknesses [m]
real :: h_scale ! A scaling factor for thicknesses [m H-1 ~> 1 or m3 kg-1]
real :: h_scale ! A scaling factor for thicknesses [m H-1 ~> 1] or [kg m-2 H-1 ~> 1]
real :: vel_scale ! A scaling factor for velocities [m T s-1 L-1 ~> 1]
real :: uh_scale ! A scaling factor for transport per unit length [m2 T s-1 L-1 H-1 ~> 1 or m3 kg-1]
real :: uh_scale ! A scaling factor for transport per unit length [m2 T s-1 L-1 H-1 ~> 1]
! or [kg T m-1 s-1 L-1 H-1 ~> 1]
real :: temp_scale ! A scaling factor for temperatures [degC C-1 ~> 1]
real :: saln_scale ! A scaling factor for salinities [ppt S-1 ~> 1]
integer :: yr, mo, day, hr, minute, sec, yearday
Expand All @@ -108,7 +109,7 @@ subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st
integer :: file

Angstrom = GV%Angstrom_H + GV%H_subroundoff
h_scale = GV%H_to_m ; vel_scale = US%L_T_to_m_s ; uh_scale = GV%H_to_m*US%L_T_to_m_s
h_scale = GV%H_to_mks ; vel_scale = US%L_T_to_m_s ; uh_scale = h_scale*vel_scale
temp_scale = US%C_to_degC ; saln_scale = US%S_to_ppt

! if (.not.associated(CS)) return
Expand Down Expand Up @@ -232,7 +233,7 @@ subroutine write_u_accel(I, j, um, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st
do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hv(I,j,k) ; enddo
endif
if (present(str)) then
write(file,'(/,"Stress: ",ES10.3)', advance='no') vel_scale*US%Z_to_m * (str*dt / GV%Rho0)
write(file,'(/,"Stress: ",ES10.3)', advance='no') (uh_scale*GV%RZ_to_H) * (str*dt)
endif

if (associated(CS%u_accel_bt)) then
Expand Down Expand Up @@ -435,9 +436,10 @@ subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st
real :: dv ! A velocity change [L T-1 ~> m s-1]
real :: Inorm(SZK_(GV)) ! The inverse of the normalized velocity change [L T-1 ~> m s-1]
real :: e(SZK_(GV)+1) ! Simple estimates of interface heights based on the sum of thicknesses [m]
real :: h_scale ! A scaling factor for thicknesses [m H-1 ~> 1 or m3 kg-1]
real :: h_scale ! A scaling factor for thicknesses [m H-1 ~> 1] or [kg m-2 H-1 ~> 1]
real :: vel_scale ! A scaling factor for velocities [m T s-1 L-1 ~> 1]
real :: uh_scale ! A scaling factor for transport per unit length [m2 T s-1 L-1 H-1 ~> 1 or m3 kg-1]
real :: uh_scale ! A scaling factor for transport per unit length [m2 T s-1 L-1 H-1 ~> 1]
! or [kg T m-1 s-1 L-1 H-1 ~> 1]
real :: temp_scale ! A scaling factor for temperatures [degC C-1 ~> 1]
real :: saln_scale ! A scaling factor for salinities [ppt S-1 ~> 1]
integer :: yr, mo, day, hr, minute, sec, yearday
Expand All @@ -448,7 +450,7 @@ subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st
integer :: file

Angstrom = GV%Angstrom_H + GV%H_subroundoff
h_scale = GV%H_to_m ; vel_scale = US%L_T_to_m_s ; uh_scale = GV%H_to_m*US%L_T_to_m_s
h_scale = GV%H_to_mks ; vel_scale = US%L_T_to_m_s ; uh_scale = h_scale*vel_scale
temp_scale = US%C_to_degC ; saln_scale = US%S_to_ppt

! if (.not.associated(CS)) return
Expand Down Expand Up @@ -576,7 +578,7 @@ subroutine write_v_accel(i, J, vm, hin, ADp, CDp, dt, G, GV, US, CS, vel_rpt, st
do k=ks,ke ; if (do_k(k)) write(file,'(ES10.3," ")', advance='no') h_scale*hv(i,J,k) ; enddo
endif
if (present(str)) then
write(file,'(/,"Stress: ",ES10.3)', advance='no') vel_scale*US%Z_to_m * (str*dt / GV%Rho0)
write(file,'(/,"Stress: ",ES10.3)', advance='no') (uh_scale*GV%RZ_to_H) * (str*dt)
endif

if (associated(CS%v_accel_bt)) then
Expand Down
20 changes: 12 additions & 8 deletions src/tracer/MOM_tracer_registry.F90
Original file line number Diff line number Diff line change
Expand Up @@ -786,14 +786,16 @@ subroutine tracer_array_chkinv(mesg, G, GV, h, Tr, ntr)
integer, intent(in) :: ntr !< number of registered tracers

! Local variables
real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1 or m3 kg-1]
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: tr_inv ! Volumetric tracer inventory in each cell [conc m3]
real :: total_inv ! The total amount of tracer [conc m3]
real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1] or cell
! masses to kg [kg H-1 L-2 ~> 1], depending on whether the Boussinesq approximation is used
real :: tr_inv(SZI_(G),SZJ_(G),SZK_(GV)) ! Volumetric or mass-based tracer inventory in
! each cell [conc m3] or [conc kg]
real :: total_inv ! The total amount of tracer [conc m3] or [conc kg]
integer :: is, ie, js, je, nz
integer :: i, j, k, m

is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke
vol_scale = GV%H_to_m*G%US%L_to_m**2
vol_scale = GV%H_to_MKS*G%US%L_to_m**2
do m=1,ntr
do k=1,nz ; do j=js,je ; do i=is,ie
tr_inv(i,j,k) = Tr(m)%conc_scale*Tr(m)%t(i,j,k) * (vol_scale * h(i,j,k) * G%areaT(i,j)*G%mask2dT(i,j))
Expand All @@ -814,16 +816,18 @@ subroutine tracer_Reg_chkinv(mesg, G, GV, h, Reg)
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]

! Local variables
real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1 or m3 kg-1]
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: tr_inv ! Volumetric tracer inventory in each cell [conc m3]
real :: total_inv ! The total amount of tracer [conc m3]
real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1] or cell
! masses to kg [kg H-1 L-2 ~> 1], depending on whether the Boussinesq approximation is used
real :: tr_inv(SZI_(G),SZJ_(G),SZK_(GV)) ! Volumetric or mass-based tracer inventory in
! each cell [conc m3] or [conc kg]
real :: total_inv ! The total amount of tracer [conc m3] or [conc kg]
integer :: is, ie, js, je, nz
integer :: i, j, k, m

if (.not.associated(Reg)) return

is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke
vol_scale = GV%H_to_m*G%US%L_to_m**2
vol_scale = GV%H_to_MKS*G%US%L_to_m**2
do m=1,Reg%ntr
do k=1,nz ; do j=js,je ; do i=is,ie
tr_inv(i,j,k) = Reg%Tr(m)%conc_scale*Reg%Tr(m)%t(i,j,k) * (vol_scale * h(i,j,k) * G%areaT(i,j)*G%mask2dT(i,j))
Expand Down

0 comments on commit 8210921

Please sign in to comment.