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.

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

Momentum budget terms multiplied by fractional layer-thicknesses #1163

Merged
merged 27 commits into from
Jul 29, 2020
Merged
Show file tree
Hide file tree
Changes from 23 commits
Commits
Show all changes
27 commits
Select commit Hold shift + click to select a range
14531c2
New diagnostics for barotropic momentum budget calculations are creat…
hmkhatri May 28, 2020
3124b4d
More thickness weighted diagnostics
hmkhatri May 30, 2020
d9e1ec2
Implemented Andrew Shao's suggestions
hmkhatri Jun 1, 2020
2a0ee28
More fractional thickness multiplied diagnostics
hmkhatri Jun 2, 2020
a18cd08
Some barotropic diagnostics obtained from fractional thickness multip…
hmkhatri Jun 11, 2020
bd29c37
Removed trailing spaces
hmkhatri Jun 11, 2020
192fb9b
All fractional-thickness multiplied diagnostics implemented
hmkhatri Jun 19, 2020
fb119b0
define diagnostic variables as pointers
hmkhatri Jun 22, 2020
a7173e0
Merge branch 'dev/gfdl' into diag_TW
hmkhatri Jun 22, 2020
6de102b
Shorter initlization of 2D arrays
hmkhatri Jun 22, 2020
8f30470
Modifications in vertical friction diagnostics
hmkhatri Jun 23, 2020
3a9a944
Diagnotics initilization as pointers to save memory allocation
hmkhatri Jun 24, 2020
53b5c77
Vertical friction module
hmkhatri Jun 28, 2020
ff54401
Removed commented lines and trailing spaces
hmkhatri Jun 30, 2020
d212cc6
Diagnostic description change
hmkhatri Jul 2, 2020
3d819c8
Diagnostic description change
hmkhatri Jul 2, 2020
7d2b052
Fractional thickness-weighted diagnostics for accelaration due to rel…
hmkhatri Jul 3, 2020
4dc88c3
Modifications in vertical friction diagnostics (fractional thickness-…
hmkhatri Jul 4, 2020
20fe0e5
Merge branch 'dev/gfdl' into diag_TW
hmkhatri Jul 10, 2020
ba0aa2f
Modifications in diag_hfrac_u and diag_hfrac_v calls
hmkhatri Jul 10, 2020
c815434
Made allocation of hfrac at u/v points optional. These arrays are onl…
hmkhatri Jul 14, 2020
306e18b
Merge branch 'dev/gfdl' into diag_TW
hmkhatri Jul 18, 2020
428e61f
Fractional-thickness weighted 3D diagnostics are now commented out as…
hmkhatri Jul 22, 2020
e2ede33
Initilization of 2D diagnostic arrays changed from pointer stype to a…
hmkhatri Jul 27, 2020
082f3b5
Trailing spaces removed and line lengths reduced to be under 120 char…
hmkhatri Jul 27, 2020
8067c5b
Comment modified
hmkhatri Jul 27, 2020
b5b3030
Merge branch 'dev/gfdl' into diag_TW
adcroft Jul 27, 2020
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
150 changes: 150 additions & 0 deletions src/core/MOM_CoriolisAdv.F90
Original file line number Diff line number Diff line change
Expand Up @@ -74,6 +74,10 @@ module MOM_CoriolisAdv
!>@{ Diagnostic IDs
integer :: id_rv = -1, id_PV = -1, id_gKEu = -1, id_gKEv = -1
integer :: id_rvxu = -1, id_rvxv = -1
! integer :: id_hf_gKEu = -1, id_hf_gKEv = -1
integer :: id_hf_gKEu_2d = -1, id_hf_gKEv_2d = -1
! integer :: id_hf_rvxu = -1, id_hf_rvxv = -1
integer :: id_hf_rvxu_2d = -1, id_hf_rvxv_2d = -1
!>@}
end type CoriolisAdv_CS

Expand Down Expand Up @@ -211,6 +215,16 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS)
real :: QUHeff,QVHeff ! More temporary variables [H L2 T-1 s-1 ~> m3 s-2 or kg s-2].
integer :: i, j, k, n, is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz

! Diagnostics for fractional thickness-weighted terms
real, pointer, dimension(:,:) :: &
hf_gKEu_2d => NULL(), hf_gKEv_2d => NULL(), & ! Depth sum of hf_gKEu, hf_gKEv [L T-2 ~> m s-2].
hf_rvxu_2d => NULL(), hf_rvxv_2d => NULL() ! Depth sum of hf_rvxu, hf_rvxv [L T-2 ~> m s-2].
!real, pointer, dimension(:,:,:) :: &
! hf_gKEu => NULL(), hf_gKEv => NULL(), & ! accel. due to KE gradient x fract. thickness [L T-2 ~> m s-2].
! hf_rvxu => NULL(), hf_rvxv => NULL() ! accel. due to RV x fract. thickness [L T-2 ~> m s-2].
! 3D diagnostics hf_gKEu etc. are commented because there is no clarity on proper remapping grid option.
! The code is retained for degugging purposes in the future.

! To work, the following fields must be set outside of the usual
! is to ie range before this subroutine is called:
! v(is-1:ie+2,js-1:je+1), u(is-1:ie+1,js-1:je+2), h(is-1:ie+2,js-1:je+2),
Expand Down Expand Up @@ -828,6 +842,78 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS)
if (CS%id_gKEv>0) call post_data(CS%id_gKEv, AD%gradKEv, CS%diag)
if (CS%id_rvxu > 0) call post_data(CS%id_rvxu, AD%rv_x_u, CS%diag)
if (CS%id_rvxv > 0) call post_data(CS%id_rvxv, AD%rv_x_v, CS%diag)

! Diagnostics for terms multiplied by fractional thicknesses

! 3D diagnostics hf_gKEu etc. are commented because there is no clarity on proper remapping grid option.
! The code is retained for degugging purposes in the future.
!if (CS%id_hf_gKEu > 0) then
! call safe_alloc_ptr(hf_gKEu,G%IsdB,G%IedB,G%jsd,G%jed,G%ke)
! do k=1,nz ; do j=js,je ; do I=Isq,Ieq
! hf_gKEu(I,j,k) = AD%gradKEu(I,j,k) * AD%diag_hfrac_u(I,j,k)
! enddo ; enddo ; enddo
! call post_data(CS%id_hf_gKEu, hf_gKEu, CS%diag)
!endif

!if (CS%id_hf_gKEv > 0) then
! call safe_alloc_ptr(hf_gKEv,G%isd,G%ied,G%JsdB,G%JedB,G%ke)
! do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
! hf_gKEv(i,J,k) = AD%gradKEv(i,J,k) * AD%diag_hfrac_v(i,J,k)
! enddo ; enddo ; enddo
! call post_data(CS%id_hf_gKEv, hf_gKEv, CS%diag)
!endif

if (CS%id_hf_gKEu_2d > 0) then
call safe_alloc_ptr(hf_gKEu_2d,G%IsdB,G%IedB,G%jsd,G%jed)
hf_gKEu_2d(:,:) = 0.0
do k=1,nz ; do j=js,je ; do I=Isq,Ieq
hf_gKEu_2d(I,j) = hf_gKEu_2d(I,j) + AD%gradKEu(I,j,k) * AD%diag_hfrac_u(I,j,k)
enddo ; enddo ; enddo
call post_data(CS%id_hf_gKEu_2d, hf_gKEu_2d, CS%diag)
endif

if (CS%id_hf_gKEv_2d > 0) then
call safe_alloc_ptr(hf_gKEv_2d,G%isd,G%ied,G%JsdB,G%JedB)
hf_gKEv_2d(:,:) = 0.0
do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
hf_gKEv_2d(i,J) = hf_gKEv_2d(i,J) + AD%gradKEv(i,J,k) * AD%diag_hfrac_v(i,J,k)
enddo ; enddo ; enddo
call post_data(CS%id_hf_gKEv_2d, hf_gKEv_2d, CS%diag)
endif

!if (CS%id_hf_rvxv > 0) then
! call safe_alloc_ptr(hf_rvxv,G%IsdB,G%IedB,G%jsd,G%jed,G%ke)
! do k=1,nz ; do j=js,je ; do I=Isq,Ieq
! hf_rvxv(I,j,k) = AD%rv_x_v(I,j,k) * AD%diag_hfrac_u(I,j,k)
! enddo ; enddo ; enddo
! call post_data(CS%id_hf_rvxv, hf_rvxv, CS%diag)
!endif

!if (CS%id_hf_rvxu > 0) then
! call safe_alloc_ptr(hf_rvxu,G%isd,G%ied,G%JsdB,G%JedB,G%ke)
! do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
! hf_rvxu(i,J,k) = AD%rv_x_u(i,J,k) * AD%diag_hfrac_v(i,J,k)
! enddo ; enddo ; enddo
! call post_data(CS%id_hf_rvxu, hf_rvxu, CS%diag)
!endif

if (CS%id_hf_rvxv_2d > 0) then
call safe_alloc_ptr(hf_rvxv_2d,G%IsdB,G%IedB,G%jsd,G%jed)
hf_rvxv_2d(:,:) = 0.0
do k=1,nz ; do j=js,je ; do I=Isq,Ieq
hf_rvxv_2d(I,j) = hf_rvxv_2d(I,j) + AD%rv_x_v(I,j,k) * AD%diag_hfrac_u(I,j,k)
enddo ; enddo ; enddo
call post_data(CS%id_hf_rvxv_2d, hf_rvxv_2d, CS%diag)
endif

if (CS%id_hf_rvxu_2d > 0) then
call safe_alloc_ptr(hf_rvxu_2d,G%isd,G%ied,G%JsdB,G%JedB)
hf_rvxu_2d(:,:) = 0.0
do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie
hf_rvxu_2d(i,J) = hf_rvxu_2d(i,J) + AD%rv_x_u(i,J,k) * AD%diag_hfrac_v(i,J,k)
enddo ; enddo ; enddo
call post_data(CS%id_hf_rvxu_2d, hf_rvxu_2d, CS%diag)
endif
endif

end subroutine CorAdCalc
Expand Down Expand Up @@ -1087,6 +1173,70 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS)
'Zonal Acceleration from Relative Vorticity', 'm-1 s-2', conversion=US%L_T2_to_m_s2)
if (CS%id_rvxv > 0) call safe_alloc_ptr(AD%rv_x_v,IsdB,IedB,jsd,jed,nz)

!CS%id_hf_gKEu = register_diag_field('ocean_model', 'hf_gKEu', diag%axesCuL, Time, &
! 'Fractional Thickness-weighted Zonal Acceleration from Grad. Kinetic Energy', &
! 'm-1 s-2', v_extensive=.true., conversion=US%L_T2_to_m_s2)
!if (CS%id_hf_gKEu > 0) then
! call safe_alloc_ptr(AD%gradKEu,IsdB,IedB,jsd,jed,nz)
! call safe_alloc_ptr(AD%diag_hfrac_u,IsdB,IedB,jsd,jed,nz)
!endif

!CS%id_hf_gKEv = register_diag_field('ocean_model', 'hf_gKEv', diag%axesCvL, Time, &
! 'Fractional Thickness-weighted Meridional Acceleration from Grad. Kinetic Energy', &
! 'm-1 s-2', v_extensive=.true., conversion=US%L_T2_to_m_s2)
!if (CS%id_hf_gKEv > 0) then
! call safe_alloc_ptr(AD%gradKEv,isd,ied,JsdB,JedB,nz)
! call safe_alloc_ptr(AD%diag_hfrac_v,isd,ied,Jsd,JedB,nz)
!endif

CS%id_hf_gKEu_2d = register_diag_field('ocean_model', 'hf_gKEu_2d', diag%axesCuL, Time, &
'Depth-sum Fractional Thickness-weighted Zonal Acceleration from Grad. Kinetic Energy', &
'm-1 s-2', conversion=US%L_T2_to_m_s2)
if (CS%id_hf_gKEu_2d > 0) then
call safe_alloc_ptr(AD%gradKEu,IsdB,IedB,jsd,jed,nz)
call safe_alloc_ptr(AD%diag_hfrac_u,IsdB,IedB,jsd,jed,nz)
endif

CS%id_hf_gKEv_2d = register_diag_field('ocean_model', 'hf_gKEv_2d', diag%axesCvL, Time, &
'Depth-sum Fractional Thickness-weighted Meridional Acceleration from Grad. Kinetic Energy', &
'm-1 s-2', conversion=US%L_T2_to_m_s2)
if (CS%id_hf_gKEv_2d > 0) then
call safe_alloc_ptr(AD%gradKEv,isd,ied,JsdB,JedB,nz)
call safe_alloc_ptr(AD%diag_hfrac_v,isd,ied,Jsd,JedB,nz)
endif

!CS%id_hf_rvxu = register_diag_field('ocean_model', 'hf_rvxu', diag%axesCvL, Time, &
! 'Fractional Thickness-weighted Meridional Acceleration from Relative Vorticity', &
! 'm-1 s-2', v_extensive=.true., conversion=US%L_T2_to_m_s2)
!if (CS%id_hf_rvxu > 0) then
! call safe_alloc_ptr(AD%rv_x_u,isd,ied,JsdB,JedB,nz)
! call safe_alloc_ptr(AD%diag_hfrac_v,isd,ied,Jsd,JedB,nz)
!endif

!CS%id_hf_rvxv = register_diag_field('ocean_model', 'hf_rvxv', diag%axesCuL, Time, &
! 'Fractional Thickness-weighted Zonal Acceleration from Relative Vorticity', &
! 'm-1 s-2', v_extensive=.true., conversion=US%L_T2_to_m_s2)
!if (CS%id_hf_rvxv > 0) then
! call safe_alloc_ptr(AD%rv_x_v,IsdB,IedB,jsd,jed,nz)
! call safe_alloc_ptr(AD%diag_hfrac_u,IsdB,IedB,jsd,jed,nz)
!endif

CS%id_hf_rvxu_2d = register_diag_field('ocean_model', 'hf_rvxu_2d', diag%axesCvL, Time, &
'Fractional Thickness-weighted Meridional Acceleration from Relative Vorticity', &
'm-1 s-2', conversion=US%L_T2_to_m_s2)
if (CS%id_hf_rvxu_2d > 0) then
call safe_alloc_ptr(AD%rv_x_u,isd,ied,JsdB,JedB,nz)
call safe_alloc_ptr(AD%diag_hfrac_v,isd,ied,Jsd,JedB,nz)
endif

CS%id_hf_rvxv_2d = register_diag_field('ocean_model', 'hf_rvxv_2d', diag%axesCuL, Time, &
'Depth-sum Fractional Thickness-weighted Zonal Acceleration from Relative Vorticity', &
'm-1 s-2', conversion=US%L_T2_to_m_s2)
if (CS%id_hf_rvxv_2d > 0) then
call safe_alloc_ptr(AD%rv_x_v,IsdB,IedB,jsd,jed,nz)
call safe_alloc_ptr(AD%diag_hfrac_u,IsdB,IedB,jsd,jed,nz)
endif

end subroutine CoriolisAdv_init

!> Destructor for coriolisadv_cs
Expand Down
15 changes: 14 additions & 1 deletion src/core/MOM_barotropic.F90
Original file line number Diff line number Diff line change
Expand Up @@ -27,6 +27,7 @@ module MOM_barotropic
use MOM_unit_scaling, only : unit_scale_type
use MOM_variables, only : BT_cont_type, alloc_bt_cont_type
use MOM_verticalGrid, only : verticalGrid_type
use MOM_variables, only : accel_diag_ptrs

implicit none ; private

Expand Down Expand Up @@ -394,7 +395,7 @@ module MOM_barotropic
subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, &
eta_PF_in, U_Cor, V_Cor, accel_layer_u, accel_layer_v, &
eta_out, uhbtav, vhbtav, G, GV, US, CS, &
visc_rem_u, visc_rem_v, etaav, OBC, BT_cont, eta_PF_start, &
visc_rem_u, visc_rem_v, etaav, ADp, OBC, BT_cont, eta_PF_start, &
taux_bot, tauy_bot, uh0, vh0, u_uh0, v_vh0)
type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
Expand Down Expand Up @@ -447,6 +448,7 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce,
real, dimension(SZI_(G),SZJB_(G),SZK_(G)), intent(in) :: visc_rem_v !< Ditto for meridional direction [nondim].
real, dimension(SZI_(G),SZJ_(G)), optional, intent(out) :: etaav !< The free surface height or column mass
!! averaged over the barotropic integration [H ~> m or kg m-2].
type(accel_diag_ptrs), optional, pointer :: ADp !< Acceleration diagnostic pointers
type(ocean_OBC_type), optional, pointer :: OBC !< The open boundary condition structure.
type(BT_cont_type), optional, pointer :: BT_cont !< A structure with elements that describe
!! the effective open face areas as a function of barotropic
Expand Down Expand Up @@ -2376,6 +2378,17 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce,
if (CS%id_frhatv1 > 0) CS%frhatv1(:,:,:) = CS%frhatv(:,:,:)
endif

if ((present(ADp)) .and. (associated(ADp%diag_hfrac_u))) then
do k=1,nz ; do j=js,je ; do I=is-1,ie
ADp%diag_hfrac_u(I,j,k) = CS%frhatu(I,j,k)
enddo ; enddo ; enddo
endif
if ((present(ADp)) .and. (associated(ADp%diag_hfrac_v))) then
do k=1,nz ; do J=js-1,je ; do i=is,ie
ADp%diag_hfrac_v(i,J,k) = CS%frhatv(i,J,k)
enddo ; enddo ; enddo
endif

if (G%nonblocking_updates) then
if (find_etaav) call complete_group_pass(CS%pass_etaav, G%Domain)
call complete_group_pass(CS%pass_ubta_uhbta, G%Domain)
Expand Down
Loading