Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

+Add and use post_product_[uv] #22

Merged
merged 3 commits into from
Dec 7, 2021
Merged

+Add and use post_product_[uv] #22

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
3aad8ca
+Add post_product_[uv] and post_product_sum_[uv]
Hallberg-NOAA Dec 3, 2021
2d053e6
Use post_product_u for momentum budget diagnostics
Hallberg-NOAA Dec 3, 2021
ef795df
Merge branch 'dev/gfdl' into post_product
Hallberg-NOAA Dec 6, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
390 changes: 87 additions & 303 deletions src/core/MOM_dynamics_split_RK2.F90
View file
Open in desktop

Large diffs are not rendered by default.

96 changes: 21 additions & 75 deletions src/diagnostics/MOM_diagnostics.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,8 @@ module MOM_diagnostics
use MOM_coupler_types, only : coupler_type_send_data
use MOM_density_integrals, only : int_density_dz
use MOM_diag_mediator, only : post_data, get_diag_time_end
use MOM_diag_mediator, only : post_product_u, post_product_sum_u
use MOM_diag_mediator, only : post_product_v, post_product_sum_v
use MOM_diag_mediator, only : register_diag_field, register_scalar_field
use MOM_diag_mediator, only : register_static_field, diag_register_area_ids
use MOM_diag_mediator, only : diag_ctrl, time_type, safe_alloc_ptr
Expand Down Expand Up @@ -226,8 +228,6 @@ subroutine calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, &
!! previous call to diagnostics_init.

! Local variables
real, dimension(SZIB_(G),SZJ_(G),SZK_(G)) :: usq ! squared eastward velocity [L2 T-2 ~> m2 s-2]
real, dimension(SZI_(G),SZJB_(G),SZK_(G)) :: vsq ! squared northward velocity [L2 T-2 ~> m2 s-2]
real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: uv ! u x v at h-points [L2 T-2 ~> m2 s-2]

integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state
Expand All @@ -238,12 +238,6 @@ subroutine calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, &
real :: work_2d(SZI_(G),SZJ_(G)) ! A 2-d temporary work array.
real :: rho_in_situ(SZI_(G)) ! In situ density [R ~> kg m-3]

real, allocatable, dimension(:,:) :: &
hf_du_dt_2d, hf_dv_dt_2d ! z integeral of hf_du_dt, hf_dv_dt [L T-2 ~> m s-2].

real, allocatable, dimension(:,:,:) :: h_du_dt ! h x dudt [H L T-2 ~> m2 s-2]
real, allocatable, dimension(:,:,:) :: h_dv_dt ! h x dvdt [H L T-2 ~> m2 s-2]

! tmp array for surface properties
real :: surface_field(SZI_(G),SZJ_(G))
real :: pressure_1d(SZI_(G)) ! Temporary array for pressure when calling EOS [R L2 T-2 ~> Pa]
Expand Down Expand Up @@ -278,70 +272,32 @@ subroutine calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, &
call diag_copy_storage_to_diag(CS%diag, diag_pre_sync)

if (CS%id_h_pre_sync > 0) &
call post_data(CS%id_h_pre_sync, diag_pre_sync%h_state, CS%diag, alt_h = diag_pre_sync%h_state)
call post_data(CS%id_h_pre_sync, diag_pre_sync%h_state, CS%diag, alt_h=diag_pre_sync%h_state)

if (CS%id_du_dt>0) call post_data(CS%id_du_dt, CS%du_dt, CS%diag, alt_h = diag_pre_sync%h_state)
if (CS%id_du_dt>0) call post_data(CS%id_du_dt, CS%du_dt, CS%diag, alt_h=diag_pre_sync%h_state)

if (CS%id_dv_dt>0) call post_data(CS%id_dv_dt, CS%dv_dt, CS%diag, alt_h = diag_pre_sync%h_state)
if (CS%id_dv_dt>0) call post_data(CS%id_dv_dt, CS%dv_dt, CS%diag, alt_h=diag_pre_sync%h_state)

if (CS%id_dh_dt>0) call post_data(CS%id_dh_dt, CS%dh_dt, CS%diag, alt_h = diag_pre_sync%h_state)
if (CS%id_dh_dt>0) call post_data(CS%id_dh_dt, CS%dh_dt, CS%diag, alt_h=diag_pre_sync%h_state)

!! Diagnostics for terms multiplied by fractional thicknesses

! 3D diagnostics hf_du(dv)_dt are commented because there is no clarity on proper remapping grid option.
! The code is retained for degugging purposes in the future.
! The code is retained for debugging purposes in the future.
!if (CS%id_hf_du_dt > 0) then
! do k=1,nz ; do j=js,je ; do I=Isq,Ieq
! CS%hf_du_dt(I,j,k) = CS%du_dt(I,j,k) * ADp%diag_hfrac_u(I,j,k)
! enddo ; enddo ; enddo
! call post_data(CS%id_hf_du_dt, CS%hf_du_dt, CS%diag, alt_h = diag_pre_sync%h_state)
!endif

!if (CS%id_hf_dv_dt > 0) then
! do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
! CS%hf_dv_dt(i,J,k) = CS%dv_dt(i,J,k) * ADp%diag_hfrac_v(i,J,k)
! enddo ; enddo ; enddo
! call post_data(CS%id_hf_dv_dt, CS%hf_dv_dt, CS%diag, alt_h = diag_pre_sync%h_state)
!endif

if (CS%id_hf_du_dt_2d > 0) then
allocate(hf_du_dt_2d(G%IsdB:G%IedB,G%jsd:G%jed))
hf_du_dt_2d(:,:) = 0.0
do k=1,nz ; do j=js,je ; do I=Isq,Ieq
hf_du_dt_2d(I,j) = hf_du_dt_2d(I,j) + CS%du_dt(I,j,k) * ADp%diag_hfrac_u(I,j,k)
enddo ; enddo ; enddo
call post_data(CS%id_hf_du_dt_2d, hf_du_dt_2d, CS%diag)
deallocate(hf_du_dt_2d)
endif
! call post_product_u(CS%id_hf_du_dt, CS%du_dt, ADp%diag_hfrac_u, G, nz, CS%diag, alt_h=diag_pre_sync%h_state)
!if (CS%id_hf_dv_dt > 0) &
! call post_product_v(CS%id_hf_dv_dt, CS%dv_dt, ADp%diag_hfrac_v, G, nz, CS%diag, alt_h=diag_pre_sync%h_state)

if (CS%id_hf_dv_dt_2d > 0) then
allocate(hf_dv_dt_2d(G%isd:G%ied,G%JsdB:G%JedB))
hf_dv_dt_2d(:,:) = 0.0
do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
hf_dv_dt_2d(i,J) = hf_dv_dt_2d(i,J) + CS%dv_dt(i,J,k) * ADp%diag_hfrac_v(i,J,k)
enddo ; enddo ; enddo
call post_data(CS%id_hf_dv_dt_2d, hf_dv_dt_2d, CS%diag)
deallocate(hf_dv_dt_2d)
endif
if (CS%id_hf_du_dt_2d > 0) &
call post_product_sum_u(CS%id_hf_du_dt_2d, CS%du_dt, ADp%diag_hfrac_u, G, nz, CS%diag)
if (CS%id_hf_dv_dt_2d > 0) &
call post_product_sum_v(CS%id_hf_dv_dt_2d, CS%dv_dt, ADp%diag_hfrac_v, G, nz, CS%diag)

if (CS%id_h_du_dt > 0) then
allocate(h_du_dt(G%IsdB:G%IedB,G%jsd:G%jed,GV%ke))
h_du_dt(:,:,:) = 0.0
do k=1,nz ; do j=js,je ; do I=Isq,Ieq
h_du_dt(I,j,k) = CS%du_dt(I,j,k) * ADp%diag_hu(I,j,k)
enddo ; enddo ; enddo
call post_data(CS%id_h_du_dt, h_du_dt, CS%diag)
deallocate(h_du_dt)
endif
if (CS%id_h_dv_dt > 0) then
allocate(h_dv_dt(G%isd:G%ied,G%JsdB:G%JedB,GV%ke))
h_dv_dt(:,:,:) = 0.0
do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
h_dv_dt(i,J,k) = CS%dv_dt(i,J,k) * ADp%diag_hv(i,J,k)
enddo ; enddo ; enddo
call post_data(CS%id_h_dv_dt, h_dv_dt, CS%diag)
deallocate(h_dv_dt)
endif
if (CS%id_h_du_dt > 0) &
call post_product_u(CS%id_h_du_dt, CS%du_dt, ADp%diag_hu, G, nz, CS%diag)
if (CS%id_h_dv_dt > 0) &
call post_product_v(CS%id_h_dv_dt, CS%dv_dt, ADp%diag_hv, G, nz, CS%diag)

call diag_restore_grids(CS%diag)

Expand All @@ -362,24 +318,14 @@ subroutine calculate_diagnostic_fields(u, v, h, uh, vh, tv, ADp, CDp, p_surf, &

if (CS%id_h > 0) call post_data(CS%id_h, h, CS%diag)

if (CS%id_usq > 0) then
do k=1,nz ; do j=js,je ; do I=Isq,Ieq
usq(I,j,k) = u(I,j,k) * u(I,j,k)
enddo ; enddo ; enddo
call post_data(CS%id_usq, usq, CS%diag)
endif
if (CS%id_usq > 0) call post_product_u(CS%id_usq, u, u, G, nz, CS%diag)

if (CS%id_vsq > 0) then
do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
vsq(i,J,k) = v(i,J,k) * v(i,J,k)
enddo ; enddo ; enddo
call post_data(CS%id_vsq, vsq, CS%diag)
endif
if (CS%id_vsq > 0) call post_product_v(CS%id_vsq, v, v, G, nz, CS%diag)

if (CS%id_uv > 0) then
do k=1,nz ; do j=js,je ; do i=is,ie
uv(i,j,k) = (0.5*(u(I-1,j,k) + u(I,j,k))) * &
(0.5*(v(i,J-1,k) + v(i,J,k)))
(0.5*(v(i,J-1,k) + v(i,J,k)))
enddo ; enddo ; enddo
call post_data(CS%id_uv, uv, CS%diag)
endif
Expand Down
103 changes: 103 additions & 0 deletions src/framework/MOM_diag_mediator.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,7 @@ module MOM_diag_mediator
#define MAX_DSAMP_LEV 2

public set_axes_info, post_data, register_diag_field, time_type
public post_product_u, post_product_sum_u, post_product_v, post_product_sum_v
public set_masks_for_axes
public post_data_1d_k
public safe_alloc_ptr, safe_alloc_alloc
Expand Down Expand Up @@ -1802,6 +1803,108 @@ subroutine post_data_3d_low(diag, field, diag_cs, is_static, mask)

end subroutine post_data_3d_low

!> Calculate and write out diagnostics that are the product of two 3-d arrays at u-points
subroutine post_product_u(id, u_a, u_b, G, nz, diag, mask, alt_h)
integer, intent(in) :: id !< The ID for this diagnostic
type(ocean_grid_type), intent(in) :: G !< ocean grid structure
integer, intent(in) :: nz !< The size of the arrays in the vertical
real, dimension(G%IsdB:G%IedB, G%jsd:G%jed, nz), &
intent(in) :: u_a !< The first u-point array in arbitrary units [A]
real, dimension(G%IsdB:G%IedB, G%jsd:G%jed, nz), &
intent(in) :: u_b !< The second u-point array in arbitrary units [B]
type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
real, optional, intent(in) :: mask(:,:,:) !< If present, use this real array as the data mask [nondim]
real, target, optional, intent(in) :: alt_h(:,:,:) !< An alternate thickness to use for vertically
!! remapping this diagnostic [H ~> m or kg m-2]

! Local variables
real, dimension(G%IsdB:G%IedB, G%jsd:G%jed, nz) :: u_prod ! The product of u_a and u_b [A B]
integer :: i, j, k

if (id <= 0) return

do k=1,nz ; do j=G%jsc,G%jec ; do I=G%IscB,G%IecB
u_prod(I,j,k) = u_a(I,j,k) * u_b(I,j,k)
enddo ; enddo ; enddo
call post_data(id, u_prod, diag, mask=mask, alt_h=alt_h)

end subroutine post_product_u

!> Calculate and write out diagnostics that are the vertical sum of the product of two 3-d arrays at u-points
subroutine post_product_sum_u(id, u_a, u_b, G, nz, diag)
integer, intent(in) :: id !< The ID for this diagnostic
type(ocean_grid_type), intent(in) :: G !< ocean grid structure
integer, intent(in) :: nz !< The size of the arrays in the vertical
real, dimension(G%IsdB:G%IedB, G%jsd:G%jed, nz), &
intent(in) :: u_a !< The first u-point array in arbitrary units [A]
real, dimension(G%IsdB:G%IedB, G%jsd:G%jed, nz), &
intent(in) :: u_b !< The second u-point array in arbitrary units [B]
type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output

real, dimension(G%IsdB:G%IedB, G%jsd:G%jed) :: u_sum ! The vertical sum of the product of u_a and u_b [A B]
integer :: i, j, k

if (id <= 0) return

u_sum(:,:) = 0.0
do k=1,nz ; do j=G%jsc,G%jec ; do I=G%IscB,G%IecB
u_sum(I,j) = u_sum(I,j) + u_a(I,j,k) * u_b(I,j,k)
enddo ; enddo ; enddo
call post_data(id, u_sum, diag)

end subroutine post_product_sum_u

!> Calculate and write out diagnostics that are the product of two 3-d arrays at v-points
subroutine post_product_v(id, v_a, v_b, G, nz, diag, mask, alt_h)
integer, intent(in) :: id !< The ID for this diagnostic
type(ocean_grid_type), intent(in) :: G !< ocean grid structure
integer, intent(in) :: nz !< The size of the arrays in the vertical
real, dimension(G%isd:G%ied, G%JsdB:G%JedB, nz), &
intent(in) :: v_a !< The first v-point array in arbitrary units [A]
real, dimension(G%isd:G%ied, G%JsdB:G%JedB, nz), &
intent(in) :: v_b !< The second v-point array in arbitrary units [B]
type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output
real, optional, intent(in) :: mask(:,:,:) !< If present, use this real array as the data mask [nondim]
real, target, optional, intent(in) :: alt_h(:,:,:) !< An alternate thickness to use for vertically
!! remapping this diagnostic [H ~> m or kg m-2]

! Local variables
real, dimension(G%isd:G%ied, G%JsdB:G%JedB, nz) :: v_prod ! The product of v_a and v_b [A B]
integer :: i, j, k

if (id <= 0) return

do k=1,nz ; do J=G%JscB,G%JecB ; do i=G%isc,G%iec
v_prod(i,J,k) = v_a(i,J,k) * v_b(i,J,k)
enddo ; enddo ; enddo
call post_data(id, v_prod, diag, mask=mask, alt_h=alt_h)

end subroutine post_product_v

!> Calculate and write out diagnostics that are the vertical sum of the product of two 3-d arrays at v-points
subroutine post_product_sum_v(id, v_a, v_b, G, nz, diag)
integer, intent(in) :: id !< The ID for this diagnostic
type(ocean_grid_type), intent(in) :: G !< ocean grid structure
integer, intent(in) :: nz !< The size of the arrays in the vertical
real, dimension(G%isd:G%ied, G%JsdB:G%JedB, nz), &
intent(in) :: v_a !< The first v-point array in arbitrary units [A]
real, dimension(G%isd:G%ied, G%JsdB:G%JedB, nz), &
intent(in) :: v_b !< The second v-point array in arbitrary units [B]
type(diag_ctrl), intent(in) :: diag !< regulates diagnostic output

real, dimension(G%isd:G%ied, G%JsdB:G%JedB) :: v_sum ! The vertical sum of the product of v_a and v_b [A B]
integer :: i, j, k

if (id <= 0) return

v_sum(:,:) = 0.0
do k=1,nz ; do J=G%JscB,G%JecB ; do i=G%isc,G%iec
v_sum(i,J) = v_sum(i,J) + v_a(i,J,k) * v_b(i,J,k)
enddo ; enddo ; enddo
call post_data(id, v_sum, diag)

end subroutine post_product_sum_v

!> Post the horizontally area-averaged diagnostic
subroutine post_xy_average(diag_cs, diag, field)
type(diag_type), intent(in) :: diag !< This diagnostic
Expand Down
Loading