diff --git a/src/core/MOM_CoriolisAdv.F90 b/src/core/MOM_CoriolisAdv.F90 index 8dab711d32..e13bd492b8 100644 --- a/src/core/MOM_CoriolisAdv.F90 +++ b/src/core/MOM_CoriolisAdv.F90 @@ -76,8 +76,10 @@ module MOM_CoriolisAdv 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_intz_gKEu_2d = -1, id_intz_gKEv_2d = -1 ! integer :: id_hf_rvxu = -1, id_hf_rvxv = -1 integer :: id_hf_rvxu_2d = -1, id_hf_rvxv_2d = -1 + integer :: id_intz_rvxu_2d = -1, id_intz_rvxv_2d = -1 !>@} end type CoriolisAdv_CS @@ -219,12 +221,17 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS) real, allocatable, dimension(:,:) :: & hf_gKEu_2d, hf_gKEv_2d, & ! Depth sum of hf_gKEu, hf_gKEv [L T-2 ~> m s-2]. hf_rvxu_2d, hf_rvxv_2d ! Depth sum of hf_rvxu, hf_rvxv [L T-2 ~> m s-2]. + !real, allocatable, dimension(:,:,:) :: & ! hf_gKEu, hf_gKEv, & ! accel. due to KE gradient x fract. thickness [L T-2 ~> m s-2]. ! hf_rvxu, hf_rvxv ! 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. +! Diagnostics for depth-integrated momentum budget terms + real, dimension(SZIB_(G),SZJ_(G)) :: intz_gKEu_2d, intz_rvxv_2d ! [L2 T-2 ~> m2 s-2]. + real, dimension(SZI_(G),SZJB_(G)) :: intz_gKEv_2d, intz_rvxu_2d ! [L2 T-2 ~> m2 s-2]. + ! 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), @@ -883,6 +890,22 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS) deallocate(hf_gKEv_2d) endif + if (CS%id_intz_gKEu_2d > 0) then + intz_gKEu_2d(:,:) = 0.0 + do k=1,nz ; do j=js,je ; do I=Isq,Ieq + intz_gKEu_2d(I,j) = intz_gKEu_2d(I,j) + AD%gradKEu(I,j,k) * AD%diag_hu(I,j,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_gKEu_2d, intz_gKEu_2d, CS%diag) + endif + + if (CS%id_intz_gKEv_2d > 0) then + intz_gKEv_2d(:,:) = 0.0 + do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie + intz_gKEv_2d(i,J) = intz_gKEv_2d(i,J) + AD%gradKEv(i,J,k) * AD%diag_hv(i,J,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_gKEv_2d, intz_gKEv_2d, CS%diag) + endif + !if (CS%id_hf_rvxv > 0) then ! allocate(hf_rvxv(G%IsdB:G%IedB,G%jsd:G%jed,GV%ke)) ! do k=1,nz ; do j=js,je ; do I=Isq,Ieq @@ -918,6 +941,22 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS) call post_data(CS%id_hf_rvxu_2d, hf_rvxu_2d, CS%diag) deallocate(hf_rvxu_2d) endif + + if (CS%id_intz_rvxv_2d > 0) then + intz_rvxv_2d(:,:) = 0.0 + do k=1,nz ; do j=js,je ; do I=Isq,Ieq + intz_rvxv_2d(I,j) = intz_rvxv_2d(I,j) + AD%rv_x_v(I,j,k) * AD%diag_hu(I,j,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_rvxv_2d, intz_rvxv_2d, CS%diag) + endif + + if (CS%id_intz_rvxu_2d > 0) then + intz_rvxu_2d(:,:) = 0.0 + do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie + intz_rvxu_2d(i,J) = intz_rvxu_2d(i,J) + AD%rv_x_u(i,J,k) * AD%diag_hv(i,J,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_rvxu_2d, intz_rvxu_2d, CS%diag) + endif endif end subroutine CorAdCalc @@ -1163,24 +1202,24 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) 'Potential Vorticity', 'm-1 s-1', conversion=GV%m_to_H*US%s_to_T) CS%id_gKEu = register_diag_field('ocean_model', 'gKEu', diag%axesCuL, Time, & - 'Zonal Acceleration from Grad. Kinetic Energy', 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'Zonal Acceleration from Grad. Kinetic Energy', 'm s-2', conversion=US%L_T2_to_m_s2) if (CS%id_gKEu > 0) call safe_alloc_ptr(AD%gradKEu,IsdB,IedB,jsd,jed,nz) CS%id_gKEv = register_diag_field('ocean_model', 'gKEv', diag%axesCvL, Time, & - 'Meridional Acceleration from Grad. Kinetic Energy', 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'Meridional Acceleration from Grad. Kinetic Energy', 'm s-2', conversion=US%L_T2_to_m_s2) if (CS%id_gKEv > 0) call safe_alloc_ptr(AD%gradKEv,isd,ied,JsdB,JedB,nz) CS%id_rvxu = register_diag_field('ocean_model', 'rvxu', diag%axesCvL, Time, & - 'Meridional Acceleration from Relative Vorticity', 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'Meridional Acceleration from Relative Vorticity', 'm s-2', conversion=US%L_T2_to_m_s2) if (CS%id_rvxu > 0) call safe_alloc_ptr(AD%rv_x_u,isd,ied,JsdB,JedB,nz) CS%id_rvxv = register_diag_field('ocean_model', 'rvxv', diag%axesCuL, Time, & - 'Zonal Acceleration from Relative Vorticity', 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'Zonal Acceleration from Relative Vorticity', 'm 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) + ! 'm 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) @@ -1188,7 +1227,7 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) !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) + ! 'm 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) @@ -1196,7 +1235,7 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) CS%id_hf_gKEu_2d = register_diag_field('ocean_model', 'hf_gKEu_2d', diag%axesCu1, Time, & 'Depth-sum Fractional Thickness-weighted Zonal Acceleration from Grad. Kinetic Energy', & - 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'm 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) @@ -1204,12 +1243,28 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) CS%id_hf_gKEv_2d = register_diag_field('ocean_model', 'hf_gKEv_2d', diag%axesCv1, Time, & 'Depth-sum Fractional Thickness-weighted Meridional Acceleration from Grad. Kinetic Energy', & - 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'm 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_intz_gKEu_2d = register_diag_field('ocean_model', 'intz_gKEu_2d', diag%axesCu1, Time, & + 'Depth-integral of Zonal Acceleration from Grad. Kinetic Energy', & + 'm2 s-2', conversion=GV%H_to_m*US%L_T2_to_m_s2) + if (CS%id_intz_gKEu_2d > 0) then + call safe_alloc_ptr(AD%gradKEu,IsdB,IedB,jsd,jed,nz) + call safe_alloc_ptr(AD%diag_hu,IsdB,IedB,jsd,jed,nz) + endif + + CS%id_intz_gKEv_2d = register_diag_field('ocean_model', 'intz_gKEv_2d', diag%axesCv1, Time, & + 'Depth-integral of Meridional Acceleration from Grad. Kinetic Energy', & + 'm2 s-2', conversion=GV%H_to_m*US%L_T2_to_m_s2) + if (CS%id_intz_gKEv_2d > 0) then + call safe_alloc_ptr(AD%gradKEv,isd,ied,JsdB,JedB,nz) + call safe_alloc_ptr(AD%diag_hv,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) @@ -1227,8 +1282,8 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) !endif CS%id_hf_rvxu_2d = register_diag_field('ocean_model', 'hf_rvxu_2d', diag%axesCv1, Time, & - 'Fractional Thickness-weighted Meridional Acceleration from Relative Vorticity', & - 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'Depth-sum Fractional Thickness-weighted Meridional Acceleration from Relative Vorticity', & + 'm 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) @@ -1236,12 +1291,28 @@ subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) CS%id_hf_rvxv_2d = register_diag_field('ocean_model', 'hf_rvxv_2d', diag%axesCu1, Time, & 'Depth-sum Fractional Thickness-weighted Zonal Acceleration from Relative Vorticity', & - 'm-1 s-2', conversion=US%L_T2_to_m_s2) + 'm 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 + CS%id_intz_rvxu_2d = register_diag_field('ocean_model', 'intz_rvxu_2d', diag%axesCv1, Time, & + 'Depth-integral of Meridional Acceleration from Relative Vorticity', & + 'm2 s-2', conversion=GV%H_to_m*US%L_T2_to_m_s2) + if (CS%id_intz_rvxu_2d > 0) then + call safe_alloc_ptr(AD%rv_x_u,isd,ied,JsdB,JedB,nz) + call safe_alloc_ptr(AD%diag_hv,isd,ied,Jsd,JedB,nz) + endif + + CS%id_intz_rvxv_2d = register_diag_field('ocean_model', 'intz_rvxv_2d', diag%axesCu1, Time, & + 'Depth-integral of Fractional Thickness-weighted Zonal Acceleration from Relative Vorticity', & + 'm2 s-2', conversion=GV%H_to_m*US%L_T2_to_m_s2) + if (CS%id_intz_rvxv_2d > 0) then + call safe_alloc_ptr(AD%rv_x_v,IsdB,IedB,jsd,jed,nz) + call safe_alloc_ptr(AD%diag_hu,IsdB,IedB,jsd,jed,nz) + endif + end subroutine CoriolisAdv_init !> Destructor for coriolisadv_cs diff --git a/src/core/MOM_barotropic.F90 b/src/core/MOM_barotropic.F90 index f30b86f57a..17e7ebee40 100644 --- a/src/core/MOM_barotropic.F90 +++ b/src/core/MOM_barotropic.F90 @@ -2675,6 +2675,17 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, enddo ; enddo ; enddo endif + if ((present(ADp)) .and. (present(BT_cont)) .and. (associated(ADp%diag_hu))) then + do k=1,nz ; do j=js,je ; do I=is-1,ie + ADp%diag_hu(I,j,k) = BT_cont%h_u(I,j,k) + enddo ; enddo ; enddo + endif + if ((present(ADp)) .and. (present(BT_cont)) .and. (associated(ADp%diag_hv))) then + do k=1,nz ; do J=js-1,je ; do i=is,ie + ADp%diag_hv(i,J,k) = BT_cont%h_v(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) diff --git a/src/core/MOM_dynamics_split_RK2.F90 b/src/core/MOM_dynamics_split_RK2.F90 index a89514cde4..66f142eb21 100644 --- a/src/core/MOM_dynamics_split_RK2.F90 +++ b/src/core/MOM_dynamics_split_RK2.F90 @@ -161,14 +161,17 @@ module MOM_dynamics_split_RK2 integer :: id_CAu = -1, id_CAv = -1 ! integer :: id_hf_PFu = -1, id_hf_PFv = -1 integer :: id_hf_PFu_2d = -1, id_hf_PFv_2d = -1 + integer :: id_intz_PFu_2d = -1, id_intz_PFv_2d = -1 ! integer :: id_hf_CAu = -1, id_hf_CAv = -1 integer :: id_hf_CAu_2d = -1, id_hf_CAv_2d = -1 + integer :: id_intz_CAu_2d = -1, id_intz_CAv_2d = -1 ! Split scheme only. integer :: id_uav = -1, id_vav = -1 integer :: id_u_BT_accel = -1, id_v_BT_accel = -1 ! integer :: id_hf_u_BT_accel = -1, id_hf_v_BT_accel = -1 integer :: id_hf_u_BT_accel_2d = -1, id_hf_v_BT_accel_2d = -1 + integer :: id_intz_u_BT_accel_2d = -1, id_intz_v_BT_accel_2d = -1 !>@} type(diag_ctrl), pointer :: diag !< A structure that is used to regulate the @@ -336,6 +339,12 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s hf_CAu_2d, hf_CAv_2d, & ! Depth integeral of hf_CAu, hf_CAv [L T-2 ~> m s-2]. hf_u_BT_accel_2d, hf_v_BT_accel_2d ! Depth integeral of hf_u_BT_accel, hf_v_BT_accel + real, dimension(SZIB_(G),SZJ_(G)) :: & + intz_PFu_2d, intz_CAu_2d, intz_u_BT_accel_2d ! [L2 T-2 ~> m2 s-2]. + + real, dimension(SZI_(G),SZJB_(G)) :: & + intz_PFv_2d, intz_CAv_2d, intz_v_BT_accel_2d ! [L2 T-2 ~> m2 s-2]. + real :: dt_pred ! The time step for the predictor part of the baroclinic time stepping [T ~> s]. logical :: dyn_p_surf @@ -897,6 +906,21 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s ! enddo ; enddo ; enddo ! call post_data(CS%id_hf_PFv, hf_PFv, CS%diag) !endif + if (CS%id_intz_PFu_2d > 0) then + intz_PFu_2d(:,:) = 0.0 + do k=1,nz ; do j=js,je ; do I=Isq,Ieq + intz_PFu_2d(I,j) = intz_PFu_2d(I,j) + CS%PFu(I,j,k) * CS%ADp%diag_hu(I,j,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_PFu_2d, intz_PFu_2d, CS%diag) + endif + if (CS%id_intz_PFv_2d > 0) then + intz_PFv_2d(:,:) = 0.0 + do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie + intz_PFv_2d(i,J) = intz_PFv_2d(i,J) + CS%PFv(i,J,k) * CS%ADp%diag_hv(i,J,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_PFv_2d, intz_PFv_2d, CS%diag) + endif + if (CS%id_hf_PFu_2d > 0) then allocate(hf_PFu_2d(G%IsdB:G%IedB,G%jsd:G%jed)) hf_PFu_2d(:,:) = 0.0 @@ -930,6 +954,21 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s ! enddo ; enddo ; enddo ! call post_data(CS%id_hf_CAv, hf_CAv, CS%diag) !endif + if (CS%id_intz_CAu_2d > 0) then + intz_CAu_2d(:,:) = 0.0 + do k=1,nz ; do j=js,je ; do I=Isq,Ieq + intz_CAu_2d(I,j) = intz_CAu_2d(I,j) + CS%CAu(I,j,k) * CS%ADp%diag_hu(I,j,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_CAu_2d, intz_CAu_2d, CS%diag) + endif + if (CS%id_intz_CAv_2d > 0) then + intz_CAv_2d(:,:) = 0.0 + do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie + intz_CAv_2d(i,J) = intz_CAv_2d(i,J) + CS%CAv(i,J,k) * CS%ADp%diag_hv(i,J,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_CAv_2d, intz_CAv_2d, CS%diag) + endif + if (CS%id_hf_CAu_2d > 0) then allocate(hf_CAu_2d(G%IsdB:G%IedB,G%jsd:G%jed)) hf_CAu_2d(:,:) = 0.0 @@ -963,6 +1002,21 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s ! enddo ; enddo ; enddo ! call post_data(CS%id_hf_v_BT_accel, hf_v_BT_accel, CS%diag) !endif + if (CS%id_intz_u_BT_accel_2d > 0) then + intz_u_BT_accel_2d(:,:) = 0.0 + do k=1,nz ; do j=js,je ; do I=Isq,Ieq + intz_u_BT_accel_2d(I,j) = intz_u_BT_accel_2d(I,j) + CS%u_accel_bt(I,j,k) * CS%ADp%diag_hu(I,j,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_u_BT_accel_2d, intz_u_BT_accel_2d, CS%diag) + endif + if (CS%id_intz_v_BT_accel_2d > 0) then + intz_v_BT_accel_2d(:,:) = 0.0 + do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie + intz_v_BT_accel_2d(i,J) = intz_v_BT_accel_2d(i,J) + CS%v_accel_bt(i,J,k) * CS%ADp%diag_hv(i,J,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_v_BT_accel_2d, intz_v_BT_accel_2d, CS%diag) + endif + if (CS%id_hf_u_BT_accel_2d > 0) then allocate(hf_u_BT_accel_2d(G%IsdB:G%IedB,G%jsd:G%jed)) hf_u_BT_accel_2d(:,:) = 0.0 @@ -1388,6 +1442,16 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param conversion=US%L_T2_to_m_s2) if(CS%id_hf_PFv_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hfrac_v,isd,ied,Jsd,JedB,nz) + CS%id_intz_PFu_2d = register_diag_field('ocean_model', 'intz_PFu_2d', diag%axesCu1, Time, & + 'Depth-integral of Zonal Pressure Force Acceleration', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if(CS%id_intz_PFu_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hu,IsdB,IedB,jsd,jed,nz) + + CS%id_intz_PFv_2d = register_diag_field('ocean_model', 'intz_PFv_2d', diag%axesCv1, Time, & + 'Depth-integral of Meridional Pressure Force Acceleration', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if(CS%id_intz_PFv_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hv,isd,ied,Jsd,JedB,nz) + CS%id_hf_CAu_2d = register_diag_field('ocean_model', 'hf_CAu_2d', diag%axesCu1, Time, & 'Depth-sum Fractional Thickness-weighted Zonal Coriolis and Advective Acceleration', 'm s-2', & conversion=US%L_T2_to_m_s2) @@ -1398,6 +1462,16 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param conversion=US%L_T2_to_m_s2) if(CS%id_hf_CAv_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hfrac_v,isd,ied,Jsd,JedB,nz) + CS%id_intz_CAu_2d = register_diag_field('ocean_model', 'intz_CAu_2d', diag%axesCu1, Time, & + 'Depth-integral of Zonal Coriolis and Advective Acceleration', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if(CS%id_intz_CAu_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hu,IsdB,IedB,jsd,jed,nz) + + CS%id_intz_CAv_2d = register_diag_field('ocean_model', 'intz_CAv_2d', diag%axesCv1, Time, & + 'Depth-integral of Meridional Coriolis and Advective Acceleration', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if(CS%id_intz_CAv_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hv,isd,ied,Jsd,JedB,nz) + CS%id_uav = register_diag_field('ocean_model', 'uav', diag%axesCuL, Time, & 'Barotropic-step Averaged Zonal Velocity', 'm s-1', conversion=US%L_T_to_m_s) CS%id_vav = register_diag_field('ocean_model', 'vav', diag%axesCvL, Time, & @@ -1428,6 +1502,16 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param conversion=US%L_T2_to_m_s2) if(CS%id_hf_v_BT_accel_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hfrac_v,isd,ied,Jsd,JedB,nz) + CS%id_intz_u_BT_accel_2d = register_diag_field('ocean_model', 'intz_u_BT_accel_2d', diag%axesCu1, Time, & + 'Depth-integral of Barotropic Anomaly Zonal Acceleration', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if(CS%id_intz_u_BT_accel_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hu,IsdB,IedB,jsd,jed,nz) + + CS%id_intz_v_BT_accel_2d = register_diag_field('ocean_model', 'intz_v_BT_accel_2d', diag%axesCv1, Time, & + 'Depth-integral of Barotropic Anomaly Meridional Acceleration', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if(CS%id_intz_v_BT_accel_2d > 0) call safe_alloc_ptr(CS%ADp%diag_hv,isd,ied,Jsd,JedB,nz) + id_clock_Cor = cpu_clock_id('(Ocean Coriolis & mom advection)', grain=CLOCK_MODULE) id_clock_continuity = cpu_clock_id('(Ocean continuity equation)', grain=CLOCK_MODULE) id_clock_pres = cpu_clock_id('(Ocean pressure force)', grain=CLOCK_MODULE) diff --git a/src/core/MOM_variables.F90 b/src/core/MOM_variables.F90 index d81cf28e17..886ee77510 100644 --- a/src/core/MOM_variables.F90 +++ b/src/core/MOM_variables.F90 @@ -188,6 +188,9 @@ module MOM_variables real, pointer :: diag_hfrac_u(:,:,:) => NULL() !< Fractional layer thickness at u points real, pointer :: diag_hfrac_v(:,:,:) => NULL() !< Fractional layer thickness at v points + real, pointer :: diag_hu(:,:,:) => NULL() !< layer thickness at u points + real, pointer :: diag_hv(:,:,:) => NULL() !< layer thickness at v points + end type accel_diag_ptrs !> Pointers to arrays with transports, which can later be used for derived diagnostics, like energy balances. diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index fc6d97040a..5bbc495e93 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -188,6 +188,7 @@ module MOM_hor_visc integer :: id_diffu = -1, id_diffv = -1 ! integer :: id_hf_diffu = -1, id_hf_diffv = -1 integer :: id_hf_diffu_2d = -1, id_hf_diffv_2d = -1 + integer :: id_intz_diffu_2d = -1, id_intz_diffv_2d = -1 integer :: id_Ah_h = -1, id_Ah_q = -1 integer :: id_Kh_h = -1, id_Kh_q = -1 integer :: id_GME_coeff_h = -1, id_GME_coeff_q = -1 @@ -276,8 +277,9 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, grad_d2vel_mag_h, & ! Magnitude of the Laplacian of the velocity vector, squared [L-2 T-2 ~> m-2 s-2] boundary_mask_h ! A mask that zeroes out cells with at least one land edge [nondim] - real, allocatable, dimension(:,:) :: hf_diffu_2d ! Depth sum of hf_diffu [L T-2 ~> m s-2] - real, allocatable, dimension(:,:) :: hf_diffv_2d ! Depth sum of hf_diffv [L T-2 ~> m s-2] + real, allocatable, dimension(:,:) :: hf_diffu_2d, hf_diffv_2d ! Depth sum of hf_diffu, hf_diffv [L T-2 ~> m s-2] + real, dimension(SZIB_(G),SZJ_(G)) :: intz_diffu_2d ! Depth-integral of diffu [L2 T-2 ~> m2 s-2] + real, dimension(SZI_(G),SZJB_(G)) :: intz_diffv_2d ! Depth-integral of diffv [L2 T-2 ~> m2 s-2] real, dimension(SZIB_(G),SZJB_(G)) :: & dvdx, dudy, & ! components in the shearing strain [T-1 ~> s-1] @@ -1662,6 +1664,21 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, deallocate(hf_diffv_2d) endif + if (present(ADp) .and. (CS%id_intz_diffu_2d > 0)) then + intz_diffu_2d(:,:) = 0.0 + do k=1,nz ; do j=js,je ; do I=Isq,Ieq + intz_diffu_2d(I,j) = intz_diffu_2d(I,j) + diffu(I,j,k) * ADp%diag_hu(I,j,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_diffu_2d, intz_diffu_2d, CS%diag) + endif + if (present(ADp) .and. (CS%id_intz_diffv_2d > 0)) then + intz_diffv_2d(:,:) = 0.0 + do k=1,nz ; do J=Jsq,Jeq ; do i=is,ie + intz_diffv_2d(i,J) = intz_diffv_2d(i,J) + diffv(i,J,k) * ADp%diag_hv(i,J,k) + enddo ; enddo ; enddo + call post_data(CS%id_intz_diffv_2d, intz_diffv_2d, CS%diag) + endif + end subroutine horizontal_viscosity !> Allocates space for and calculates static variables used by horizontal_viscosity(). @@ -2378,6 +2395,20 @@ subroutine hor_visc_init(Time, G, GV, US, param_file, diag, CS, MEKE, ADp) call safe_alloc_ptr(ADp%diag_hfrac_v,G%isd,G%ied,G%JsdB,G%JedB,GV%ke) endif + CS%id_intz_diffu_2d = register_diag_field('ocean_model', 'intz_diffu_2d', diag%axesCu1, Time, & + 'Depth-integral of Zonal Acceleration from Horizontal Viscosity', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if ((CS%id_intz_diffu_2d > 0) .and. (present(ADp))) then + call safe_alloc_ptr(ADp%diag_hu,G%IsdB,G%IedB,G%jsd,G%jed,GV%ke) + endif + + CS%id_intz_diffv_2d = register_diag_field('ocean_model', 'intz_diffv_2d', diag%axesCv1, Time, & + 'Depth-integral of Meridional Acceleration from Horizontal Viscosity', 'm2 s-2', & + conversion=GV%H_to_m*US%L_T2_to_m_s2) + if ((CS%id_intz_diffv_2d > 0) .and. (present(ADp))) then + call safe_alloc_ptr(ADp%diag_hv,G%isd,G%ied,G%JsdB,G%JedB,GV%ke) + endif + if (CS%biharmonic) then CS%id_Ah_h = register_diag_field('ocean_model', 'Ahh', diag%axesTL, Time, & 'Biharmonic Horizontal Viscosity at h Points', 'm4 s-1', conversion=US%L_to_m**4*US%s_to_T, &