From c3c1106377efc414dac392a29ec2719741f4b130 Mon Sep 17 00:00:00 2001 From: Alistair Adcroft Date: Sun, 19 Jun 2016 11:03:59 -0400 Subject: [PATCH] Split OBC_kind_u into OBC_kind_u and OBC_direction_u - OBC_kind_u/v used to indicate whether an open boundary was OBC_SIMPLE or OBC_FLATHER_* and if the latter would also indicate the outward direction. This was not extensible and also made rational reconciliation of masks awkward. - OBC_kind_u/v can now be only OBC_NONE, OBC_SIMPLE, OBC_FLATHER and OBC_direction_u/v has been introduced to indicate outward direction. Although direction is not needed for OBC_SIMPLE it does no harm AND it does help when generalizing the setup that needs to know which side of the open boundary bathymetry could be altered. - No answer changes. --- src/core/MOM_barotropic.F90 | 180 +++++++------- src/core/MOM_continuity_PPM.F90 | 20 +- src/core/MOM_legacy_barotropic.F90 | 232 ++++++++++-------- src/core/MOM_open_boundary.F90 | 116 ++++++--- .../lateral/MOM_hor_visc.F90 | 9 +- src/tracer/MOM_tracer_advect.F90 | 12 +- src/tracer/MOM_tracer_hor_diff.F90 | 2 - src/user/DOME_initialization.F90 | 2 + src/user/user_initialization.F90 | 4 +- 9 files changed, 322 insertions(+), 255 deletions(-) diff --git a/src/core/MOM_barotropic.F90 b/src/core/MOM_barotropic.F90 index a359a9f90f..3fb030c3ce 100644 --- a/src/core/MOM_barotropic.F90 +++ b/src/core/MOM_barotropic.F90 @@ -105,9 +105,9 @@ module MOM_barotropic use MOM_grid, only : ocean_grid_type use MOM_hor_index, only : hor_index_type use MOM_io, only : vardesc, var_desc -use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE, OBC_NONE -use MOM_open_boundary, only : OBC_FLATHER_E, OBC_FLATHER_W -use MOM_open_boundary, only : OBC_FLATHER_N, OBC_FLATHER_S +use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE, OBC_NONE, OBC_FLATHER +use MOM_open_boundary, only : OBC_DIRECTION_E, OBC_DIRECTION_W +use MOM_open_boundary, only : OBC_DIRECTION_N, OBC_DIRECTION_S use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS use MOM_tidal_forcing, only : tidal_forcing_sensitivity, tidal_forcing_CS use MOM_time_manager, only : time_type, set_time, operator(+), operator(-) @@ -362,6 +362,8 @@ module MOM_barotropic OBC_mask_u => NULL(), & OBC_mask_v => NULL() integer, dimension(:,:), pointer :: & + OBC_direction_u => NULL(), & + OBC_direction_v => NULL(), & OBC_kind_u => NULL(), & OBC_kind_v => NULL() real, dimension(:,:), pointer :: & @@ -2378,44 +2380,46 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans, uhbt(I,j) = BT_OBC%uhbt(I,j) ubt(I,j) = BT_OBC%ubt_outer(I,j) vel_trans = ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then - cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL - u_inlet = cfl*ubt_old(I-1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 - ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external - h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i-1,j)) ! internal - - H_u = BT_OBC%H_u(I,j) - vel_prev = ubt(I,j) - ubt(I,j) = 0.5*((u_inlet + BT_OBC%ubt_outer(I,j)) + & - (BT_OBC%Cg_u(I,j)/H_u) * (h_in-BT_OBC%eta_outer_u(I,j))) - - vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then - cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL - u_inlet = cfl*ubt_old(I+1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 - ! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external - h_in = eta(i+1,j) + (0.5-cfl)*(eta(i+1,j)-eta(i+2,j)) ! internal - - H_u = BT_OBC%H_u(I,j) - vel_prev = ubt(I,j) - ubt(I,j) = 0.5*((u_inlet+BT_OBC%ubt_outer(I,j)) + & - (BT_OBC%Cg_u(I,j)/H_u) * (BT_OBC%eta_outer_u(I,j)-h_in)) - - vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_N) then - if ((vbt(i,J-1)+vbt(i+1,J-1)) > 0.0) then - ubt(I,j) = 2.0*ubt(I,j-1)-ubt(I,j-2) - else - ubt(I,j) = BT_OBC%ubt_outer(I,j) - endif - vel_trans = ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_S) then - if ((vbt(i,J)+vbt(i+1,J)) > 0.0) then - ubt(I,j) = 2.0*ubt(I,j+1)-ubt(I,j+2) - else - ubt(I,j) = BT_OBC%ubt_outer(I,j) + elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then + cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL + u_inlet = cfl*ubt_old(I-1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 + ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external + h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i-1,j)) ! internal + + H_u = BT_OBC%H_u(I,j) + vel_prev = ubt(I,j) + ubt(I,j) = 0.5*((u_inlet + BT_OBC%ubt_outer(I,j)) + & + (BT_OBC%Cg_u(I,j)/H_u) * (h_in-BT_OBC%eta_outer_u(I,j))) + + vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then + cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL + u_inlet = cfl*ubt_old(I+1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 + ! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external + h_in = eta(i+1,j) + (0.5-cfl)*(eta(i+1,j)-eta(i+2,j)) ! internal + + H_u = BT_OBC%H_u(I,j) + vel_prev = ubt(I,j) + ubt(I,j) = 0.5*((u_inlet+BT_OBC%ubt_outer(I,j)) + & + (BT_OBC%Cg_u(I,j)/H_u) * (BT_OBC%eta_outer_u(I,j)-h_in)) + + vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_N) then + if ((vbt(i,J-1)+vbt(i+1,J-1)) > 0.0) then + ubt(I,j) = 2.0*ubt(I,j-1)-ubt(I,j-2) + else + ubt(I,j) = BT_OBC%ubt_outer(I,j) + endif + vel_trans = ubt(I,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_S) then + if ((vbt(i,J)+vbt(i+1,J)) > 0.0) then + ubt(I,j) = 2.0*ubt(I,j+1)-ubt(I,j+2) + else + ubt(I,j) = BT_OBC%ubt_outer(I,j) + endif + vel_trans = ubt(I,j) endif - vel_trans = ubt(I,j) endif if (BT_OBC%OBC_kind_u(I,j) /= OBC_SIMPLE) then @@ -2436,52 +2440,54 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans, vhbt(i,J) = BT_OBC%vhbt(i,J) vbt(i,J) = BT_OBC%vbt_outer(i,J) vel_trans = vbt(i,J) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then - cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL - v_inlet = cfl*vbt_old(i,J-1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl<1 - ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external - h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i,j-1)) ! internal - - H_v = BT_OBC%H_v(i,J) - vel_prev = vbt(i,J) - vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & - (BT_OBC%Cg_v(i,J)/H_v) * (h_in-BT_OBC%eta_outer_v(i,J))) - - vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then - cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL - v_inlet = cfl*vbt_old(i,J+1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl <1 - ! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external - h_in = eta(i,j+1) + (0.5-cfl)*(eta(i,j+1)-eta(i,j+2)) ! internal - - H_v = BT_OBC%H_v(i,J) - vel_prev = vbt(i,J) - vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & - (BT_OBC%Cg_v(i,J)/H_v) * (BT_OBC%eta_outer_v(i,J)-h_in)) - - vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_E) then - if ((ubt(I-1,j)+ubt(I-1,j+1)) > 0.0) then - vbt(i,J) = 2.0*vbt(i-1,J)-vbt(i-2,J) - else - vbt(i,J) = BT_OBC%vbt_outer(i,J) - endif - vel_trans = vbt(i,J) + elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then + cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL + v_inlet = cfl*vbt_old(i,J-1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl<1 + ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external + h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i,j-1)) ! internal + + H_v = BT_OBC%H_v(i,J) + vel_prev = vbt(i,J) + vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & + (BT_OBC%Cg_v(i,J)/H_v) * (h_in-BT_OBC%eta_outer_v(i,J))) + + vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then + cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL + v_inlet = cfl*vbt_old(i,J+1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl <1 + ! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external + h_in = eta(i,j+1) + (0.5-cfl)*(eta(i,j+1)-eta(i,j+2)) ! internal + + H_v = BT_OBC%H_v(i,J) + vel_prev = vbt(i,J) + vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & + (BT_OBC%Cg_v(i,J)/H_v) * (BT_OBC%eta_outer_v(i,J)-h_in)) + + vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_E) then + if ((ubt(I-1,j)+ubt(I-1,j+1)) > 0.0) then + vbt(i,J) = 2.0*vbt(i-1,J)-vbt(i-2,J) + else + vbt(i,J) = BT_OBC%vbt_outer(i,J) + endif + vel_trans = vbt(i,J) !!!!!!!!!!!!!!!!!!! CLAMPED !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(i,J) ! ! vbt(i,J) = (vbt(i-1,J) + CFL*vbt(i,J)) / (1.0 + CFL) ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_W) then - if ((ubt(I,j)+ubt(I,j+1)) < 0.0) then - vbt(i,J) = 2.0*vbt(i+1,J)-vbt(i+2,J) - else - vbt(i,J) = BT_OBC%vbt_outer(i,J) - endif - vel_trans = vbt(i,J) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_W) then + if ((ubt(I,j)+ubt(I,j+1)) < 0.0) then + vbt(i,J) = 2.0*vbt(i+1,J)-vbt(i+2,J) + else + vbt(i,J) = BT_OBC%vbt_outer(i,J) + endif + vel_trans = vbt(i,J) !!!!!!!!!!!!!!!!!! CLAMPED !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(i,J) ! ! vbt(i,J) = (vbt(i-1,J) + CFL*vbt(i,J)) / (1.0 + CFL) ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + endif endif if (BT_OBC%OBC_kind_v(i,J) /= OBC_SIMPLE) then @@ -2532,8 +2538,8 @@ subroutine apply_eta_OBCs(OBC, eta, ubt, vbt, BT_OBC, G, MS, halo, dtbt) if ((OBC%apply_OBC_u_flather_east .or. OBC%apply_OBC_u_flather_west) .and. & associated(BT_OBC%OBC_mask_u)) then - do j=js,je ; do I=is-1,ie ; if (BT_OBC%OBC_mask_u(I,j)) then - if (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then + do j=js,je ; do I=is-1,ie ; if (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL u_inlet = cfl*ubt(I-1,j) + (1.0-cfl)*ubt(I,j) ! Valid for cfl <1 ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external @@ -2542,7 +2548,7 @@ subroutine apply_eta_OBCs(OBC, eta, ubt, vbt, BT_OBC, G, MS, halo, dtbt) H_u = BT_OBC%H_u(I,j) eta(i+1,j) = 2.0 * 0.5*((BT_OBC%eta_outer_u(I,j)+h_in) + & (H_u/BT_OBC%Cg_u(I,j))*(u_inlet-BT_OBC%ubt_outer(I,j))) - eta(i,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then cfl = dtbt*BT_OBC%Cg_u(I,j)*G%IdxCu(I,j) ! CFL u_inlet = cfl*ubt(I+1,j) + (1.0-cfl)*ubt(I,j) ! Valid for cfl <1 ! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external @@ -2557,8 +2563,8 @@ subroutine apply_eta_OBCs(OBC, eta, ubt, vbt, BT_OBC, G, MS, halo, dtbt) if ((OBC%apply_OBC_v_flather_north .or. OBC%apply_OBC_v_flather_south) .and. & associated(BT_OBC%OBC_mask_v)) then - do J=js-1,je ; do i=is,ie ; if (BT_OBC%OBC_mask_v(i,J)) then - if (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then + do J=js-1,je ; do i=is,ie ; if (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then cfl = dtbt*BT_OBC%Cg_v(i,J)*G%IdyCv(i,J) ! CFL v_inlet = cfl*vbt(i,J-1) + (1.0-cfl)*vbt(i,J) ! Valid for cfl <1 ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external @@ -2567,7 +2573,7 @@ subroutine apply_eta_OBCs(OBC, eta, ubt, vbt, BT_OBC, G, MS, halo, dtbt) H_v = BT_OBC%H_v(i,J) eta(i,j+1) = 2.0 * 0.5*((BT_OBC%eta_outer_v(i,J)+h_in) + & (H_v/BT_OBC%Cg_v(i,J))*(v_inlet-BT_OBC%vbt_outer(i,J))) - eta(i,j) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then cfl = dtbt*BT_OBC%Cg_v(i,J)*G%IdyCv(i,J) ! CFL v_inlet = cfl*vbt(i,J+1) + (1.0-cfl)*vbt(i,J) ! Valid for cfl <1 ! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external @@ -2635,6 +2641,7 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D allocate(BT_OBC%eta_outer_u(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%eta_outer_u(:,:) = 0.0 allocate(BT_OBC%OBC_mask_u(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%OBC_mask_u(:,:)=.false. allocate(BT_OBC%OBC_kind_u(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%OBC_kind_u(:,:)=OBC_NONE + allocate(BT_OBC%OBC_direction_u(isdw-1:iedw,jsdw:jedw)); BT_OBC%OBC_direction_u(:,:)=OBC_NONE allocate(BT_OBC%Cg_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%Cg_v(:,:) = 0.0 allocate(BT_OBC%H_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%H_v(:,:) = 0.0 @@ -2643,11 +2650,13 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D allocate(BT_OBC%eta_outer_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%eta_outer_v(:,:)=0.0 allocate(BT_OBC%OBC_mask_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%OBC_mask_v(:,:)=.false. allocate(BT_OBC%OBC_kind_v(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%OBC_kind_v(:,:)=OBC_NONE + allocate(BT_OBC%OBC_direction_v(isdw-1:iedw,jsdw:jedw)); BT_OBC%OBC_direction_v(:,:)=OBC_NONE if (associated(OBC%OBC_mask_u)) then do j=js-1,je+1 ; do I=is-1,ie BT_OBC%OBC_mask_u(I,j) = OBC%OBC_mask_u(I,j) BT_OBC%OBC_kind_u(I,j) = OBC%OBC_kind_u(I,j) + BT_OBC%OBC_direction_u(I,j) = OBC%OBC_direction_u(I,j) enddo ; enddo if (OBC%apply_OBC_u) then do k=1,nz ; do j=js,je ; do I=is-1,ie @@ -2683,6 +2692,7 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D do J=js-1,je ; do i=is-1,ie+1 BT_OBC%OBC_mask_v(i,J) = OBC%OBC_mask_v(i,J) BT_OBC%OBC_kind_v(i,J) = OBC%OBC_kind_v(i,J) + BT_OBC%OBC_direction_v(i,J) = OBC%OBC_direction_v(i,J) enddo ; enddo if (OBC%apply_OBC_v) then do k=1,nz ; do J=js-1,je ; do i=is,ie @@ -2732,6 +2742,7 @@ subroutine destroy_BT_OBC(BT_OBC) if (associated(BT_OBC%OBC_mask_u)) deallocate(BT_OBC%OBC_mask_u) if (associated(BT_OBC%OBC_kind_u)) deallocate(BT_OBC%OBC_kind_u) + if (associated(BT_OBC%OBC_direction_u)) deallocate(BT_OBC%OBC_direction_u) deallocate(BT_OBC%Cg_u) deallocate(BT_OBC%H_u) deallocate(BT_OBC%uhbt) @@ -2740,6 +2751,7 @@ subroutine destroy_BT_OBC(BT_OBC) if (associated(BT_OBC%OBC_mask_v)) deallocate(BT_OBC%OBC_mask_v) if (associated(BT_OBC%OBC_kind_v)) deallocate(BT_OBC%OBC_kind_v) + if (associated(BT_OBC%OBC_direction_v)) deallocate(BT_OBC%OBC_direction_v) deallocate(BT_OBC%Cg_v) deallocate(BT_OBC%H_v) deallocate(BT_OBC%vhbt) diff --git a/src/core/MOM_continuity_PPM.F90 b/src/core/MOM_continuity_PPM.F90 index 8a31121f1c..2562ca995f 100644 --- a/src/core/MOM_continuity_PPM.F90 +++ b/src/core/MOM_continuity_PPM.F90 @@ -48,8 +48,8 @@ module MOM_continuity_PPM use MOM_error_handler, only : MOM_error, FATAL, WARNING, is_root_pe use MOM_file_parser, only : get_param, log_version, param_file_type use MOM_grid, only : ocean_grid_type -use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE -use MOM_open_boundary, only : OBC_FLATHER_E, OBC_FLATHER_W, OBC_FLATHER_N, OBC_FLATHER_S +use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE, OBC_FLATHER +use MOM_open_boundary, only : OBC_DIRECTION_E, OBC_DIRECTION_W, OBC_DIRECTION_N, OBC_DIRECTION_S use MOM_variables, only : BT_cont_type use MOM_verticalGrid, only : verticalGrid_type @@ -229,11 +229,11 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC, if (apply_OBC_u_flather_east .or. apply_OBC_u_flather_west) then do k=1,nz ; do j=LB%jsh,LB%jeh do I=LB%ish,LB%ieh+1 - if (OBC%OBC_mask_u(I-1,j) .and. (OBC%OBC_kind_u(I-1,j) == OBC_FLATHER_E)) & + if (OBC%OBC_kind_u(I-1,j) == OBC_FLATHER .and. (OBC%OBC_direction_u(I-1,j) == OBC_DIRECTION_E)) & h(i,j,k) = h_input(i-1,j,k) enddo do i=LB%ish-1,LB%ieh - if (OBC%OBC_mask_u(I,j) .and. (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W)) & + if (OBC%OBC_kind_u(I,j) == OBC_FLATHER .and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W)) & h(i,j,k) = h_input(i+1,j,k) enddo enddo ; enddo @@ -257,11 +257,11 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC, if (apply_OBC_v_flather_north .or. apply_OBC_v_flather_south) then do k=1,nz do J=LB%jsh,LB%jeh+1 ; do i=LB%ish-1,LB%ieh+1 - if (OBC%OBC_mask_v(i,J-1) .and. (OBC%OBC_kind_v(i,J-1) == OBC_FLATHER_N)) & + if (OBC%OBC_kind_v(i,J-1) == OBC_FLATHER .and. (OBC%OBC_direction_v(i,J-1) == OBC_DIRECTION_N)) & h(i,j,k) = h_input(i,j-1,k) enddo ; enddo do J=LB%jsh-1,LB%jeh ; do i=LB%ish-1,LB%ieh+1 - if (OBC%OBC_mask_v(i,J) .and. (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S)) & + if (OBC%OBC_kind_v(i,J) == OBC_FLATHER .and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S)) & h(i,j,k) = h_input(i,j+1,k) enddo ; enddo enddo @@ -284,11 +284,11 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC, if (apply_OBC_v_flather_north .or. apply_OBC_v_flather_south) then do k=1,nz do J=LB%jsh,LB%jeh+1 ; do i=LB%ish-1,LB%ieh+1 - if (OBC%OBC_mask_v(i,J-1) .and. (OBC%OBC_kind_v(i,J-1) == OBC_FLATHER_N)) & + if (OBC%OBC_kind_v(i,J-1) == OBC_FLATHER .and. (OBC%OBC_direction_v(i,J-1) == OBC_DIRECTION_N)) & h(i,j,k) = h_input(i,j-1,k) enddo ; enddo do J=LB%jsh-1,LB%jeh ; do i=LB%ish-1,LB%ieh+1 - if (OBC%OBC_mask_v(i,J) .and. (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S)) & + if (OBC%OBC_kind_v(i,J) == OBC_FLATHER .and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S)) & h(i,j,k) = h_input(i,j+1,k) enddo ; enddo enddo @@ -312,11 +312,11 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC, if (apply_OBC_u_flather_east .or. apply_OBC_u_flather_west) then do k=1,nz ; do j=LB%jsh,LB%jeh do I=LB%ish,LB%ieh+1 - if (OBC%OBC_mask_u(I-1,j) .and. (OBC%OBC_kind_u(I-1,j) == OBC_FLATHER_E)) & + if (OBC%OBC_kind_u(I-1,j) == OBC_FLATHER .and. (OBC%OBC_direction_u(I-1,j) == OBC_DIRECTION_E)) & h(i,j,k) = h_input(i-1,j,k) enddo do i=LB%ish-1,LB%ieh - if (OBC%OBC_mask_u(I,j) .and. (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W)) & + if (OBC%OBC_kind_u(I,j) == OBC_FLATHER .and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W)) & h(i,j,k) = h_input(i+1,j,k) enddo enddo ; enddo diff --git a/src/core/MOM_legacy_barotropic.F90 b/src/core/MOM_legacy_barotropic.F90 index c6f077ec6a..d8232d1b77 100644 --- a/src/core/MOM_legacy_barotropic.F90 +++ b/src/core/MOM_legacy_barotropic.F90 @@ -107,9 +107,9 @@ module MOM_legacy_barotropic use MOM_grid, only : ocean_grid_type use MOM_hor_index, only : hor_index_type use MOM_io, only : vardesc, var_desc -use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE, OBC_NONE -use MOM_open_boundary, only : OBC_FLATHER_E, OBC_FLATHER_W -use MOM_open_boundary, only : OBC_FLATHER_N, OBC_FLATHER_S +use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE, OBC_NONE, OBC_FLATHER +use MOM_open_boundary, only : OBC_DIRECTION_E, OBC_DIRECTION_W +use MOM_open_boundary, only : OBC_DIRECTION_N, OBC_DIRECTION_S use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS use MOM_tidal_forcing, only : tidal_forcing_sensitivity, tidal_forcing_CS use MOM_time_manager, only : time_type, set_time, operator(+), operator(-) @@ -352,6 +352,8 @@ module MOM_legacy_barotropic OBC_mask_u => NULL(), & OBC_mask_v => NULL() integer, dimension(:,:), pointer :: & + OBC_direction_u => NULL(), & + OBC_direction_v => NULL(), & OBC_kind_u => NULL(), & OBC_kind_v => NULL() real, dimension(:,:), pointer :: & @@ -2236,44 +2238,46 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans, uhbt(I,j) = BT_OBC%uhbt(I,j) ubt(I,j) = BT_OBC%ubt_outer(I,j) vel_trans = ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then - cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL - u_inlet = cfl*ubt_old(I-1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 - ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external - h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i-1,j)) ! internal + elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then + cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL + u_inlet = cfl*ubt_old(I-1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 + ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external + h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i-1,j)) ! internal + + H_u = BT_OBC%H_u(I,j) + vel_prev = ubt(I,j) + ubt(I,j) = 0.5*((u_inlet + BT_OBC%ubt_outer(I,j)) + & + (BT_OBC%Cg_u(I,j)/H_u) * (h_in-BT_OBC%eta_outer_u(I,j))) - H_u = BT_OBC%H_u(I,j) - vel_prev = ubt(I,j) - ubt(I,j) = 0.5*((u_inlet + BT_OBC%ubt_outer(I,j)) + & - (BT_OBC%Cg_u(I,j)/H_u) * (h_in-BT_OBC%eta_outer_u(I,j))) - - vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then - cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL - u_inlet = cfl*ubt_old(I+1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 - ! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external - h_in = eta(i+1,j) + (0.5-cfl)*(eta(i+1,j)-eta(i+2,j)) ! internal + vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then + cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL + u_inlet = cfl*ubt_old(I+1,j) + (1.0-cfl)*ubt_old(I,j) ! Valid for cfl<1 + ! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external + h_in = eta(i+1,j) + (0.5-cfl)*(eta(i+1,j)-eta(i+2,j)) ! internal - H_u = BT_OBC%H_u(I,j) - vel_prev = ubt(I,j) - ubt(I,j) = 0.5*((u_inlet+BT_OBC%ubt_outer(I,j)) + & - (BT_OBC%Cg_u(I,j)/H_u) * (BT_OBC%eta_outer_u(I,j)-h_in)) + H_u = BT_OBC%H_u(I,j) + vel_prev = ubt(I,j) + ubt(I,j) = 0.5*((u_inlet+BT_OBC%ubt_outer(I,j)) + & + (BT_OBC%Cg_u(I,j)/H_u) * (BT_OBC%eta_outer_u(I,j)-h_in)) - vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_N) then - if ((vbt(i,J-1)+vbt(i+1,J-1)) > 0.0) then - ubt(I,j) = 2.0*ubt(I,j-1)-ubt(I,j-2) - else - ubt(I,j) = BT_OBC%ubt_outer(I,j) - endif - vel_trans = ubt(I,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_S) then - if ((vbt(i,J)+vbt(i+1,J)) > 0.0) then - ubt(I,j) = 2.0*ubt(I,j+1)-ubt(I,j+2) - else - ubt(I,j) = BT_OBC%ubt_outer(I,j) + vel_trans = (1.0-bebt)*vel_prev + bebt*ubt(I,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_N) then + if ((vbt(i,J-1)+vbt(i+1,J-1)) > 0.0) then + ubt(I,j) = 2.0*ubt(I,j-1)-ubt(I,j-2) + else + ubt(I,j) = BT_OBC%ubt_outer(I,j) + endif + vel_trans = ubt(I,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_S) then + if ((vbt(i,J)+vbt(i+1,J)) > 0.0) then + ubt(I,j) = 2.0*ubt(I,j+1)-ubt(I,j+2) + else + ubt(I,j) = BT_OBC%ubt_outer(I,j) + endif + vel_trans = ubt(I,j) endif - vel_trans = ubt(I,j) endif if (BT_OBC%OBC_kind_u(I,j) /= OBC_SIMPLE) then @@ -2294,52 +2298,54 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans, vhbt(i,J) = BT_OBC%vhbt(i,J) vbt(i,J) = BT_OBC%vbt_outer(i,J) vel_trans = vbt(i,J) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then - cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL - v_inlet = cfl*vbt_old(i,J-1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl<1 - ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external - h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i,j-1)) ! internal - - H_v = BT_OBC%H_v(i,J) - vel_prev = vbt(i,J) - vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & - (BT_OBC%Cg_v(i,J)/H_v) * (h_in-BT_OBC%eta_outer_v(i,J))) + elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then + cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL + v_inlet = cfl*vbt_old(i,J-1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl<1 + ! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external + h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i,j-1)) ! internal + + H_v = BT_OBC%H_v(i,J) + vel_prev = vbt(i,J) + vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & + (BT_OBC%Cg_v(i,J)/H_v) * (h_in-BT_OBC%eta_outer_v(i,J))) - vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then - cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL - v_inlet = cfl*vbt_old(i,J+1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl <1 - ! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external - h_in = eta(i,j+1) + (0.5-cfl)*(eta(i,j+1)-eta(i,j+2)) ! internal + vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then + cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(I,j) ! CFL + v_inlet = cfl*vbt_old(i,J+1) + (1.0-cfl)*vbt_old(i,J) ! Valid for cfl <1 + ! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external + h_in = eta(i,j+1) + (0.5-cfl)*(eta(i,j+1)-eta(i,j+2)) ! internal - H_v = BT_OBC%H_v(i,J) - vel_prev = vbt(i,J) - vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & - (BT_OBC%Cg_v(i,J)/H_v) * (BT_OBC%eta_outer_v(i,J)-h_in)) + H_v = BT_OBC%H_v(i,J) + vel_prev = vbt(i,J) + vbt(i,J) = 0.5*((v_inlet+BT_OBC%vbt_outer(i,J)) + & + (BT_OBC%Cg_v(i,J)/H_v) * (BT_OBC%eta_outer_v(i,J)-h_in)) - vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_E) then - if ((ubt(I-1,j)+ubt(I-1,j+1)) > 0.0) then - vbt(i,J) = 2.0*vbt(i-1,J)-vbt(i-2,J) - else - vbt(i,J) = BT_OBC%vbt_outer(i,J) - endif - vel_trans = vbt(i,J) + vel_trans = (1.0-bebt)*vel_prev + bebt*vbt(i,J) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_E) then + if ((ubt(I-1,j)+ubt(I-1,j+1)) > 0.0) then + vbt(i,J) = 2.0*vbt(i-1,J)-vbt(i-2,J) + else + vbt(i,J) = BT_OBC%vbt_outer(i,J) + endif + vel_trans = vbt(i,J) !!!!!!!!!!!!!!!!!!! CLAMPED !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(i,J) ! ! vbt(i,J) = (vbt(i-1,J) + CFL*vbt(i,J)) / (1.0 + CFL) ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_W) then - if ((ubt(I,j)+ubt(I,j+1)) < 0.0) then - vbt(i,J) = 2.0*vbt(i+1,J)-vbt(i+2,J) - else - vbt(i,J) = BT_OBC%vbt_outer(i,J) - endif - vel_trans = vbt(i,J) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_W) then + if ((ubt(I,j)+ubt(I,j+1)) < 0.0) then + vbt(i,J) = 2.0*vbt(i+1,J)-vbt(i+2,J) + else + vbt(i,J) = BT_OBC%vbt_outer(i,J) + endif + vel_trans = vbt(i,J) !!!!!!!!!!!!!!!!!! CLAMPED !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! ! cfl = dtbt * BT_OBC%Cg_v(i,J) * G%IdyCv(i,J) ! ! vbt(i,J) = (vbt(i-1,J) + CFL*vbt(i,J)) / (1.0 + CFL) ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + endif endif if (BT_OBC%OBC_kind_v(i,J) /= OBC_SIMPLE) then @@ -2391,24 +2397,26 @@ subroutine apply_eta_OBCs(OBC, eta, ubt, vbt, BT_OBC, G, MS, halo, dtbt) if ((OBC%apply_OBC_u_flather_east .or. OBC%apply_OBC_u_flather_west) .and. & associated(BT_OBC%OBC_mask_u)) then do j=js,je ; do I=is-1,ie ; if (BT_OBC%OBC_mask_u(I,j)) then - if (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then - cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL - u_inlet = cfl*ubt(I-1,j) + (1.0-cfl)*ubt(I,j) ! Valid for cfl <1 -! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external - h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i-1,j)) ! internal - - H_u = BT_OBC%H_u(I,j) - eta(i+1,j) = 2.0 * 0.5*((BT_OBC%eta_outer_u(I,j)+h_in) + & - (H_u/BT_OBC%Cg_u(I,j))*(u_inlet-BT_OBC%ubt_outer(I,j))) - eta(i,j) - elseif (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then - cfl = dtbt*BT_OBC%Cg_u(I,j)*G%IdxCu(I,j) ! CFL - u_inlet = cfl*ubt(I+1,j) + (1.0-cfl)*ubt(I,j) ! Valid for cfl <1 -! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external - h_in = eta(i+1,j) + (0.5-cfl)*(eta(i+1,j)-eta(i+2,j)) ! internal - - H_u = BT_OBC%H_u(I,j) - eta(i,j) = 2.0 * 0.5*((BT_OBC%eta_outer_u(I,j)+h_in) + & - (H_u/BT_OBC%Cg_u(I,j))*(BT_OBC%ubt_outer(I,j)-u_inlet)) - eta(i+1,j) + if (BT_OBC%OBC_kind_u(I,j) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then + cfl = dtbt * BT_OBC%Cg_u(I,j) * G%IdxCu(I,j) ! CFL + u_inlet = cfl*ubt(I-1,j) + (1.0-cfl)*ubt(I,j) ! Valid for cfl <1 +! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i+1,j) ! external + h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i-1,j)) ! internal + + H_u = BT_OBC%H_u(I,j) + eta(i+1,j) = 2.0 * 0.5*((BT_OBC%eta_outer_u(I,j)+h_in) + & + (H_u/BT_OBC%Cg_u(I,j))*(u_inlet-BT_OBC%ubt_outer(I,j))) - eta(i,j) + elseif (BT_OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then + cfl = dtbt*BT_OBC%Cg_u(I,j)*G%IdxCu(I,j) ! CFL + u_inlet = cfl*ubt(I+1,j) + (1.0-cfl)*ubt(I,j) ! Valid for cfl <1 +! h_in = 2.0*cfl*eta(i+1,j) + (1.0-2.0*cfl)*eta(i,j) ! external + h_in = eta(i+1,j) + (0.5-cfl)*(eta(i+1,j)-eta(i+2,j)) ! internal + + H_u = BT_OBC%H_u(I,j) + eta(i,j) = 2.0 * 0.5*((BT_OBC%eta_outer_u(I,j)+h_in) + & + (H_u/BT_OBC%Cg_u(I,j))*(BT_OBC%ubt_outer(I,j)-u_inlet)) - eta(i+1,j) + endif endif endif ; enddo ; enddo endif @@ -2416,24 +2424,26 @@ subroutine apply_eta_OBCs(OBC, eta, ubt, vbt, BT_OBC, G, MS, halo, dtbt) if ((OBC%apply_OBC_v_flather_north .or. OBC%apply_OBC_v_flather_south) .and. & associated(BT_OBC%OBC_mask_v)) then do J=js-1,je ; do i=is,ie ; if (BT_OBC%OBC_mask_v(i,J)) then - if (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then - cfl = dtbt*BT_OBC%Cg_v(i,J)*G%IdyCv(i,J) ! CFL - v_inlet = cfl*vbt(i,J-1) + (1.0-cfl)*vbt(i,J) ! Valid for cfl <1 -! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external - h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i,j-1)) ! internal - - H_v = BT_OBC%H_v(i,J) - eta(i,j+1) = 2.0 * 0.5*((BT_OBC%eta_outer_v(i,J)+h_in) + & - (H_v/BT_OBC%Cg_v(i,J))*(v_inlet-BT_OBC%vbt_outer(i,J))) - eta(i,j) - elseif (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then - cfl = dtbt*BT_OBC%Cg_v(i,J)*G%IdyCv(i,J) ! CFL - v_inlet = cfl*vbt(i,J+1) + (1.0-cfl)*vbt(i,J) ! Valid for cfl <1 -! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external - h_in = eta(i,j+1) + (0.5-cfl)*(eta(i,j+1)-eta(i,j+2)) ! internal - - H_v = BT_OBC%H_v(i,J) - eta(i,j) = 2.0 * 0.5*((BT_OBC%eta_outer_v(i,J)+h_in) + & - (H_v/BT_OBC%Cg_v(i,J))*(BT_OBC%vbt_outer(i,J)-v_inlet)) - eta(i,j+1) + if (BT_OBC%OBC_kind_v(i,J) == OBC_FLATHER) then + if (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then + cfl = dtbt*BT_OBC%Cg_v(i,J)*G%IdyCv(i,J) ! CFL + v_inlet = cfl*vbt(i,J-1) + (1.0-cfl)*vbt(i,J) ! Valid for cfl <1 +! h_in = 2.0*cfl*eta(i,j) + (1.0-2.0*cfl)*eta(i,j+1) ! external + h_in = eta(i,j) + (0.5-cfl)*(eta(i,j)-eta(i,j-1)) ! internal + + H_v = BT_OBC%H_v(i,J) + eta(i,j+1) = 2.0 * 0.5*((BT_OBC%eta_outer_v(i,J)+h_in) + & + (H_v/BT_OBC%Cg_v(i,J))*(v_inlet-BT_OBC%vbt_outer(i,J))) - eta(i,j) + elseif (BT_OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then + cfl = dtbt*BT_OBC%Cg_v(i,J)*G%IdyCv(i,J) ! CFL + v_inlet = cfl*vbt(i,J+1) + (1.0-cfl)*vbt(i,J) ! Valid for cfl <1 +! h_in = 2.0*cfl*eta(i,j+1) + (1.0-2.0*cfl)*eta(i,j) ! external + h_in = eta(i,j+1) + (0.5-cfl)*(eta(i,j+1)-eta(i,j+2)) ! internal + + H_v = BT_OBC%H_v(i,J) + eta(i,j) = 2.0 * 0.5*((BT_OBC%eta_outer_v(i,J)+h_in) + & + (H_v/BT_OBC%Cg_v(i,J))*(BT_OBC%vbt_outer(i,J)-v_inlet)) - eta(i,j+1) + endif endif endif ; enddo ; enddo endif @@ -2493,6 +2503,7 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D allocate(BT_OBC%eta_outer_u(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%eta_outer_u(:,:) = 0.0 allocate(BT_OBC%OBC_mask_u(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%OBC_mask_u(:,:)=.false. allocate(BT_OBC%OBC_kind_u(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%OBC_kind_u(:,:)=OBC_NONE + allocate(BT_OBC%OBC_direction_u(isdw-1:iedw,jsdw:jedw)); BT_OBC%OBC_direction_u(:,:)=OBC_NONE allocate(BT_OBC%Cg_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%Cg_v(:,:) = 0.0 allocate(BT_OBC%H_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%H_v(:,:) = 0.0 @@ -2501,11 +2512,13 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D allocate(BT_OBC%eta_outer_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%eta_outer_v(:,:)=0.0 allocate(BT_OBC%OBC_mask_v(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%OBC_mask_v(:,:)=.false. allocate(BT_OBC%OBC_kind_v(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%OBC_kind_v(:,:)=OBC_NONE + allocate(BT_OBC%OBC_direction_v(isdw-1:iedw,jsdw:jedw)); BT_OBC%OBC_direction_v(:,:)=OBC_NONE if (associated(OBC%OBC_mask_u)) then do j=js-1,je+1 ; do I=is-1,ie BT_OBC%OBC_mask_u(I,j) = OBC%OBC_mask_u(I,j) BT_OBC%OBC_kind_u(I,j) = OBC%OBC_kind_u(I,j) + BT_OBC%OBC_direction_u(I,j) = OBC%OBC_direction_u(I,j) enddo ; enddo if (OBC%apply_OBC_u) then do k=1,nz ; do j=js,je ; do I=is-1,ie @@ -2541,6 +2554,7 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D do J=js-1,je ; do i=is-1,ie+1 BT_OBC%OBC_mask_v(i,J) = OBC%OBC_mask_v(i,J) BT_OBC%OBC_kind_v(i,J) = OBC%OBC_kind_v(i,J) + BT_OBC%OBC_direction_v(i,J) = OBC%OBC_direction_v(i,J) enddo ; enddo if (OBC%apply_OBC_v) then do k=1,nz ; do J=js-1,je ; do i=is,ie @@ -2590,6 +2604,7 @@ subroutine destroy_BT_OBC(BT_OBC) if (associated(BT_OBC%OBC_mask_u)) deallocate(BT_OBC%OBC_mask_u) if (associated(BT_OBC%OBC_kind_u)) deallocate(BT_OBC%OBC_kind_u) + if (associated(BT_OBC%OBC_direction_u)) deallocate(BT_OBC%OBC_direction_u) deallocate(BT_OBC%Cg_u) deallocate(BT_OBC%H_u) deallocate(BT_OBC%uhbt) @@ -2598,6 +2613,7 @@ subroutine destroy_BT_OBC(BT_OBC) if (associated(BT_OBC%OBC_mask_v)) deallocate(BT_OBC%OBC_mask_v) if (associated(BT_OBC%OBC_kind_v)) deallocate(BT_OBC%OBC_kind_v) + if (associated(BT_OBC%OBC_direction_v)) deallocate(BT_OBC%OBC_direction_v) deallocate(BT_OBC%Cg_v) deallocate(BT_OBC%H_v) deallocate(BT_OBC%vhbt) diff --git a/src/core/MOM_open_boundary.F90 b/src/core/MOM_open_boundary.F90 index c036da9935..2e740f09ad 100644 --- a/src/core/MOM_open_boundary.F90 +++ b/src/core/MOM_open_boundary.F90 @@ -31,8 +31,11 @@ module MOM_open_boundary public set_Flather_data integer, parameter, public :: OBC_NONE = 0, OBC_SIMPLE = 1, OBC_WALL = 2 -integer, parameter, public :: OBC_FLATHER_E = 4, OBC_FLATHER_W = 5 -integer, parameter, public :: OBC_FLATHER_N = 6, OBC_FLATHER_S = 7 +integer, parameter, public :: OBC_FLATHER = 3 +integer, parameter, public :: OBC_DIRECTION_N = 100 !< Indicates the boundary is an effective northern boundary +integer, parameter, public :: OBC_DIRECTION_S = 200 !< Indicates the boundary is an effective southern boundary +integer, parameter, public :: OBC_DIRECTION_E = 300 !< Indicates the boundary is an effective eastern boundary +integer, parameter, public :: OBC_DIRECTION_W = 400 !< Indicates the boundary is an effective western boundary !> Open-boundary data type, public :: ocean_OBC_type @@ -50,13 +53,18 @@ module MOM_open_boundary OBC_mask_u => NULL(), & !< True at zonal velocity points that have prescribed OBCs. OBC_mask_v => NULL() !< True at meridional velocity points that have prescribed OBCs. ! These arrays indicate the kind of open boundary conditions that are to be applied at the u and v - ! points, and can be OBC_NONE, OBC_SIMPLE, OBC_WALL, or one of OBC_FLATHER_[EWNS]. Generally these + ! points, and can be OBC_NONE, OBC_SIMPLE, OBC_WALL, or OBC_FLATHER. Generally these ! should be consistent with OBC_mask_[uv], with OBC_mask_[uv] .false. for OBC_kind_[uv] = NONE ! and true for all other values. integer, pointer, dimension(:,:) :: & OBC_kind_u => NULL(), & !< Type of OBC at u-points. OBC_kind_v => NULL() !< Type of OBC at v-points. - ! The following apply at points with OBC_kind_[uv] = OBC_FLATHER_x. + ! These arrays indicate the outward-pointing orientation of the open boundary and will be set to + ! one of OBC_DIRECTION_N, OBC_DIRECTION_S, OBC_DIRECTION_E or OBC_DIRECTION_W. + integer, pointer, dimension(:,:) :: & + OBC_direction_u => NULL(), & !< Orientation of OBC at u-points. + OBC_direction_v => NULL() !< Orientation of OBC at v-points. + ! The following apply at points with OBC_kind_[uv] = OBC_FLATHER. real, pointer, dimension(:,:,:) :: & rx_old_u => NULL(), & !< The rx_old_u value for radiation coeff for u-velocity in x-direction ry_old_v => NULL(), & !< The ry_old_v value for radiation coeff for v-velocity in y-direction @@ -220,17 +228,17 @@ subroutine open_boundary_impose_normal_slope(OBC, G, depth) if (.not.associated(OBC)) return - if (associated(OBC%OBC_kind_u)) then + if (associated(OBC%OBC_direction_u)) then do j=G%jsd,G%jed ; do I=G%isd,G%ied-1 - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) depth(i+1,j) = depth(i,j) - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) depth(i,j) = depth(i+1,j) + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) depth(i+1,j) = depth(i,j) + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) depth(i,j) = depth(i+1,j) enddo ; enddo endif if (associated(OBC%OBC_kind_v)) then do J=G%jsd,G%jed-1 ; do i=G%isd,G%ied - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) depth(i,j+1) = depth(i,j) - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) depth(i,j) = depth(i,j+1) + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) depth(i,j+1) = depth(i,j) + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) depth(i,j) = depth(i,j+1) enddo ; enddo endif @@ -294,7 +302,7 @@ subroutine Radiation_Open_Bdry_Conds(OBC, u_new, u_old, v_new, v_old, & if (OBC%apply_OBC_u_flather_east .or. OBC%apply_OBC_u_flather_west) then do k=1,nz ; do j=js,je ; do I=is-1,ie ; if (OBC%OBC_mask_u(I,j)) then - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then dhdt = u_old(I-1,j,k)-u_new(I-1,j,k) !old-new dhdx = u_new(I-1,j,k)-u_new(I-2,j,k) !in new time backward sasha for I-1 rx_new = 0.0 @@ -311,7 +319,7 @@ subroutine Radiation_Open_Bdry_Conds(OBC, u_new, u_old, v_new, v_old, & ! OBC%rx_old_h(I,j,k) = rx_avg ! h_new(I+1,j,k) = (h_old(I+1,j,k) + rx_avg*h_new(I,j,k)) / (1.0+rx_avg) !original endif - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then dhdt = u_old(I+1,j,k)-u_new(I+1,j,k) !old-new dhdx = u_new(I+1,j,k)-u_new(I+2,j,k) !in new time backward sasha for I+1 rx_new = 0.0 @@ -333,7 +341,7 @@ subroutine Radiation_Open_Bdry_Conds(OBC, u_new, u_old, v_new, v_old, & if (OBC%apply_OBC_v_flather_north .or. OBC%apply_OBC_v_flather_south) then do k=1,nz ; do J=js-1,je ; do i=is,ie ; if (OBC%OBC_mask_v(i,J)) then - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then dhdt = v_old(i,J-1,k)-v_new(i,J-1,k) !old-new dhdx = v_new(i,J-1,k)-v_new(i,J-2,k) !in new time backward sasha for J-1 rx_new = 0.0 @@ -351,7 +359,7 @@ subroutine Radiation_Open_Bdry_Conds(OBC, u_new, u_old, v_new, v_old, & ! h_new(i,J+1,k) = (h_old(i,J+1,k) + rx_avg*h_new(i,J,k)) / (1.0+rx_avg) !original endif - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then dhdt = v_old(i,J+1,k)-v_new(i,J+1,k) !old-new dhdx = v_new(i,J+1,k)-v_new(i,J+2,k) !in new time backward sasha for J+1 rx_new = 0.0 @@ -391,12 +399,18 @@ subroutine set_Flather_positions(G, OBC) if (.not.associated(OBC%OBC_mask_u)) then allocate(OBC%OBC_mask_u(G%IsdB:G%IedB,G%jsd:G%jed)) ; OBC%OBC_mask_u(:,:) = .false. endif + if (.not.associated(OBC%OBC_direction_u)) then + allocate(OBC%OBC_direction_u(G%IsdB:G%IedB,G%jsd:G%jed)) ; OBC%OBC_direction_u(:,:) = OBC_NONE + endif if (.not.associated(OBC%OBC_kind_u)) then allocate(OBC%OBC_kind_u(G%IsdB:G%IedB,G%jsd:G%jed)) ; OBC%OBC_kind_u(:,:) = OBC_NONE endif if (.not.associated(OBC%OBC_mask_v)) then allocate(OBC%OBC_mask_v(G%isd:G%ied,G%JsdB:G%JedB)) ; OBC%OBC_mask_v(:,:) = .false. endif + if (.not.associated(OBC%OBC_direction_v)) then + allocate(OBC%OBC_direction_v(G%isd:G%ied,G%JsdB:G%JedB)) ; OBC%OBC_direction_v(:,:) = OBC_NONE + endif if (.not.associated(OBC%OBC_kind_v)) then allocate(OBC%OBC_kind_v(G%isd:G%ied,G%JsdB:G%JedB)) ; OBC%OBC_kind_v(:,:) = OBC_NONE endif @@ -421,11 +435,18 @@ subroutine set_Flather_positions(G, OBC) do j=G%jsd,G%jed ; do I=G%IsdB,G%IedB if ((I+G%idg_offset) == east_boundary) then !eastern side OBC%OBC_mask_u(I,j) = .true. - OBC%OBC_kind_u(I,j) = OBC_FLATHER_E + OBC%OBC_direction_u(I,j) = OBC_DIRECTION_E + OBC%OBC_kind_u(I,j) = OBC_FLATHER OBC%OBC_mask_v(i+1,J) = .true. - if (OBC%OBC_kind_v(i+1,J) == OBC_NONE) OBC%OBC_kind_v(i+1,J) = OBC_FLATHER_E + if (OBC%OBC_direction_v(i+1,J) == OBC_NONE) then + OBC%OBC_direction_v(i+1,J) = OBC_DIRECTION_E + OBC%OBC_kind_v(i+1,J) = OBC_FLATHER + endif OBC%OBC_mask_v(i+1,J-1) = .true. - if (OBC%OBC_kind_v(i+1,J-1) == OBC_NONE) OBC%OBC_kind_v(i+1,J-1) = OBC_FLATHER_E + if (OBC%OBC_direction_v(i+1,J-1) == OBC_NONE) then + OBC%OBC_direction_v(i+1,J-1) = OBC_DIRECTION_E + OBC%OBC_kind_v(i+1,J-1) = OBC_FLATHER + endif endif enddo ; enddo endif @@ -435,11 +456,18 @@ subroutine set_Flather_positions(G, OBC) do j=G%jsd,G%jed ; do I=G%IsdB,G%IedB if ((I+G%idg_offset) == west_boundary) then !western side OBC%OBC_mask_u(I,j) = .true. - OBC%OBC_kind_u(I,j) = OBC_FLATHER_W + OBC%OBC_direction_u(I,j) = OBC_DIRECTION_W + OBC%OBC_kind_u(I,j) = OBC_FLATHER OBC%OBC_mask_v(i,J) = .true. - if (OBC%OBC_kind_v(i,J) == OBC_NONE) OBC%OBC_kind_v(i,J) = OBC_FLATHER_W + if (OBC%OBC_direction_v(i,J) == OBC_NONE) then + OBC%OBC_direction_v(i,J) = OBC_DIRECTION_W + OBC%OBC_kind_v(i,J) = OBC_FLATHER + endif OBC%OBC_mask_v(i,J-1) = .true. - if (OBC%OBC_kind_v(i,J-1) == OBC_NONE) OBC%OBC_kind_v(i,J-1) = OBC_FLATHER_W + if (OBC%OBC_direction_v(i,J-1) == OBC_NONE) then + OBC%OBC_direction_v(i,J-1) = OBC_DIRECTION_W + OBC%OBC_kind_v(i,J-1) = OBC_FLATHER + endif endif enddo ; enddo endif @@ -449,11 +477,18 @@ subroutine set_Flather_positions(G, OBC) do J=G%JsdB,G%JedB ; do i=G%isd,G%ied if ((J+G%jdg_offset) == north_boundary) then !northern side OBC%OBC_mask_v(i,J) = .true. - OBC%OBC_kind_v(i,J) = OBC_FLATHER_N + OBC%OBC_direction_v(i,J) = OBC_DIRECTION_N + OBC%OBC_kind_v(i,J) = OBC_FLATHER OBC%OBC_mask_u(I,j+1) = .true. - if (OBC%OBC_kind_u(I,j+1) == OBC_NONE) OBC%OBC_kind_u(I,j+1) = OBC_FLATHER_N + if (OBC%OBC_direction_u(I,j+1) == OBC_NONE) then + OBC%OBC_direction_u(I,j+1) = OBC_DIRECTION_N + OBC%OBC_kind_u(I,j+1) = OBC_FLATHER + endif OBC%OBC_mask_u(I-1,j+1) = .true. - if (OBC%OBC_kind_u(I-1,j+1) == OBC_NONE) OBC%OBC_kind_u(I-1,j+1) = OBC_FLATHER_N + if (OBC%OBC_direction_u(I-1,j+1) == OBC_NONE) then + OBC%OBC_direction_u(I-1,j+1) = OBC_DIRECTION_N + OBC%OBC_kind_u(I-1,j+1) = OBC_FLATHER + endif endif enddo ; enddo endif @@ -463,11 +498,18 @@ subroutine set_Flather_positions(G, OBC) do J=G%JsdB,G%JedB ; do i=G%isd,G%ied if ((J+G%jdg_offset) == south_boundary) then !southern side OBC%OBC_mask_v(i,J) = .true. - OBC%OBC_kind_v(i,J) = OBC_FLATHER_S + OBC%OBC_direction_v(i,J) = OBC_DIRECTION_S + OBC%OBC_kind_v(i,J) = OBC_FLATHER OBC%OBC_mask_u(I,j) = .true. - if (OBC%OBC_kind_u(I,j) == OBC_NONE) OBC%OBC_kind_u(I,j) = OBC_FLATHER_S + if (OBC%OBC_direction_u(I,j) == OBC_NONE) then + OBC%OBC_direction_u(I,j) = OBC_DIRECTION_S + OBC%OBC_kind_u(I,j) = OBC_FLATHER + endif OBC%OBC_mask_u(I-1,j) = .true. - if (OBC%OBC_kind_u(I-1,j) == OBC_NONE) OBC%OBC_kind_u(I-1,j) = OBC_FLATHER_S + if (OBC%OBC_direction_u(I-1,j) == OBC_NONE) then + OBC%OBC_direction_u(I-1,j) = OBC_DIRECTION_S + OBC%OBC_kind_u(I-1,j) = OBC_FLATHER + endif endif enddo ; enddo endif @@ -598,10 +640,10 @@ subroutine set_Flather_data(OBC, tv, h, G, PF, tracer_Reg) call pass_var(tv%S, G%Domain) do k=1,nz ; do j=js,je ; do I=is-1,ie if (OBC%OBC_mask_u(I,j)) then - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then OBC_T_u(I,j,k) = tv%T(i,j,k) OBC_S_u(I,j,k) = tv%S(i,j,k) - elseif (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then + elseif (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then OBC_T_u(I,j,k) = tv%T(i+1,j,k) OBC_S_u(I,j,k) = tv%S(i+1,j,k) elseif (G%mask2dT(i,j) + G%mask2dT(i+1,j) > 0) then @@ -621,10 +663,10 @@ subroutine set_Flather_data(OBC, tv, h, G, PF, tracer_Reg) do k=1,nz ; do J=js-1,je ; do i=is,ie if (OBC%OBC_mask_v(i,J)) then - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then OBC_T_v(i,J,k) = tv%T(i,j,k) OBC_S_v(i,J,k) = tv%S(i,j,k) - elseif (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then + elseif (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then OBC_T_v(i,J,k) = tv%T(i,j+1,k) OBC_S_v(i,J,k) = tv%S(i,j+1,k) elseif (G%mask2dT(i,j) + G%mask2dT(i,j+1) > 0) then @@ -651,28 +693,28 @@ subroutine set_Flather_data(OBC, tv, h, G, PF, tracer_Reg) call add_tracer_OBC_values("S", tracer_Reg, OBC_in_u=OBC_S_u, & OBC_in_v=OBC_S_v) do k=1,nz ; do j=jsd,jed ; do I=isd,ied-1 - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) then + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) then tv%T(i+1,j,k) = tv%T(i,j,k) ; tv%S(i+1,j,k) = tv%S(i,j,k) - elseif (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) then + elseif (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) then tv%T(i,j,k) = tv%T(i+1,j,k) ; tv%S(i,j,k) = tv%S(i+1,j,k) endif enddo ; enddo ; enddo do k=1,nz ; do J=jsd,jed-1 ; do i=isd,ied - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) then + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) then tv%T(i,j+1,k) = tv%T(i,j,k) ; tv%S(i,j+1,k) = tv%S(i,j,k) - elseif (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) then + elseif (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) then tv%T(i,j,k) = tv%T(i,j+1,k) ; tv%S(i,j,k) = tv%S(i,j+1,k) endif enddo ; enddo ; enddo endif do k=1,nz ; do j=jsd,jed ; do I=isd,ied-1 - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) h(i+1,j,k) = h(i,j,k) - if (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W) h(i,j,k) = h(i+1,j,k) + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E) h(i+1,j,k) = h(i,j,k) + if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) h(i,j,k) = h(i+1,j,k) enddo ; enddo ; enddo do k=1,nz ; do J=jsd,jed-1 ; do i=isd,ied - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) h(i,j+1,k) = h(i,j,k) - if (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S) h(i,j,k) = h(i,j+1,k) + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N) h(i,j+1,k) = h(i,j,k) + if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) h(i,j,k) = h(i,j+1,k) enddo ; enddo ; enddo end subroutine set_Flather_data diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index 191eb62005..897810c527 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -93,8 +93,7 @@ module MOM_hor_visc use MOM_grid, only : ocean_grid_type use MOM_lateral_mixing_coeffs, only : VarMix_CS use MOM_MEKE_types, only : MEKE_type -use MOM_open_boundary, only : ocean_OBC_type, OBC_FLATHER_E, OBC_FLATHER_W -use MOM_open_boundary, only : OBC_FLATHER_N, OBC_FLATHER_S +use MOM_open_boundary, only : ocean_OBC_type, OBC_FLATHER use MOM_verticalGrid, only : verticalGrid_type implicit none ; private @@ -587,8 +586,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, G%IareaCu(I,j)) / (0.5*(h(i+1,j,k) + h(i,j,k)) + h_neglect) if (apply_OBC) then ; if (OBC%OBC_mask_u(I,j)) then - if ((OBC%OBC_kind_u(I,j) == OBC_FLATHER_E) .or. & - (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W)) diffu(I,j,k) = 0.0 + if ((OBC%OBC_kind_u(I,j) == OBC_FLATHER)) diffu(I,j,k) = 0.0 endif ; endif enddo ; enddo @@ -600,8 +598,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, CS%DX2h(i,j+1)*str_xx(i,j+1))) * & G%IareaCv(i,J)) / (0.5*(h(i,j+1,k) + h(i,j,k)) + h_neglect) if (apply_OBC) then ; if (OBC%OBC_mask_v(i,J)) then - if ((OBC%OBC_kind_v(i,J) == OBC_FLATHER_N) .or. & - (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S)) diffv(I,j,k) = 0.0 + if ((OBC%OBC_kind_v(i,J) == OBC_FLATHER)) diffv(I,j,k) = 0.0 endif ; endif enddo ; enddo diff --git a/src/tracer/MOM_tracer_advect.F90 b/src/tracer/MOM_tracer_advect.F90 index 949039f7e5..674a5f1e41 100644 --- a/src/tracer/MOM_tracer_advect.F90 +++ b/src/tracer/MOM_tracer_advect.F90 @@ -14,8 +14,8 @@ module MOM_tracer_advect use MOM_error_handler, only : MOM_error, FATAL, WARNING, MOM_mesg, is_root_pe use MOM_file_parser, only : get_param, log_version, param_file_type use MOM_grid, only : ocean_grid_type -use MOM_open_boundary, only : ocean_OBC_type, OBC_FLATHER_E -use MOM_open_boundary, only : OBC_FLATHER_W, OBC_FLATHER_N, OBC_FLATHER_S +use MOM_open_boundary, only : ocean_OBC_type, OBC_DIRECTION_E +use MOM_open_boundary, only : OBC_DIRECTION_W, OBC_DIRECTION_N, OBC_DIRECTION_S use MOM_tracer_registry, only : tracer_registry_type, tracer_type, MOM_tracer_chksum use MOM_verticalGrid, only : verticalGrid_type @@ -478,9 +478,9 @@ subroutine advect_x(Tr, hprev, uhr, uh_neglect, OBC, domore_u, ntr, Idt, & ! Tracer fluxes are set to prescribed values only for inflows ! from masked areas. if (((uhr(I,j,k) > 0.0) .and. ((G%mask2dT(i,j) < 0.5) .or. & - (OBC%OBC_kind_u(I,j) == OBC_FLATHER_W))) .or. & + (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W))) .or. & ((uhr(I,j,k) < 0.0) .and. ((G%mask2dT(i+1,j) < 0.5) .or. & - (OBC%OBC_kind_u(I,j) == OBC_FLATHER_E))) ) then + (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_E))) ) then do_i(I) = .true. ; do_any_i = .true. uhh(I) = uhr(I,j,k) endif @@ -738,9 +738,9 @@ subroutine advect_y(Tr, hprev, vhr, vh_neglect, OBC, domore_v, ntr, Idt, & ! Tracer fluxes are set to prescribed values only for inflows ! from masked areas. if (((vhr(i,J,k) > 0.0) .and. ((G%mask2dT(i,j) < 0.5) .or. & - (OBC%OBC_kind_v(i,J) == OBC_FLATHER_S))) .or. & + (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S))) .or. & ((vhr(i,J,k) < 0.0) .and. ((G%mask2dT(i,j+1) < 0.5) .or. & - (OBC%OBC_kind_v(i,J) == OBC_FLATHER_N))) ) then + (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_N))) ) then do_i(i) = .true. ; do_any_i = .true. vhh(i,J) = vhr(i,J,k) endif diff --git a/src/tracer/MOM_tracer_hor_diff.F90 b/src/tracer/MOM_tracer_hor_diff.F90 index 73b581d299..ab45688a66 100644 --- a/src/tracer/MOM_tracer_hor_diff.F90 +++ b/src/tracer/MOM_tracer_hor_diff.F90 @@ -21,8 +21,6 @@ module MOM_tracer_hor_diff use MOM_neutral_diffusion, only : neutral_diffusion_init, neutral_diffusion_end use MOM_neutral_diffusion, only : neutral_diffusion_CS use MOM_neutral_diffusion, only : neutral_diffusion_calc_coeffs, neutral_diffusion -use MOM_open_boundary, only : ocean_OBC_type, OBC_FLATHER_E -use MOM_open_boundary, only : OBC_FLATHER_W, OBC_FLATHER_N, OBC_FLATHER_S use MOM_tracer_registry, only : tracer_registry_type, tracer_type, MOM_tracer_chksum use MOM_variables, only : thermo_var_ptrs use MOM_verticalGrid, only : verticalGrid_type diff --git a/src/user/DOME_initialization.F90 b/src/user/DOME_initialization.F90 index d0337bcc9d..6c1723be71 100644 --- a/src/user/DOME_initialization.F90 +++ b/src/user/DOME_initialization.F90 @@ -327,6 +327,7 @@ subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg) allocate(OBC%u(IsdB:IedB,jsd:jed,nz)) ; OBC%u(:,:,:) = 0.0 allocate(OBC%uh(IsdB:IedB,jsd:jed,nz)) ; OBC%uh(:,:,:) = 0.0 allocate(OBC%OBC_kind_u(IsdB:IedB,jsd:jed)) ; OBC%OBC_kind_u(:,:) = OBC_NONE + allocate(OBC%OBC_direction_u(IsdB:IedB,jsd:jed)) ; OBC%OBC_direction_u(:,:) = OBC_NONE do j=jsd,jed ; do I=IsdB,IedB if (OBC%OBC_mask_u(I,j)) OBC%OBC_kind_u(I,j) = OBC_SIMPLE enddo ; enddo @@ -335,6 +336,7 @@ subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg) allocate(OBC%v(isd:ied,JsdB:JedB,nz)) ; OBC%v(:,:,:) = 0.0 allocate(OBC%vh(isd:ied,JsdB:JedB,nz)) ; OBC%vh(:,:,:) = 0.0 allocate(OBC%OBC_kind_v(isd:ied,JsdB:JedB)) ; OBC%OBC_kind_v(:,:) = OBC_NONE + allocate(OBC%OBC_direction_v(isd:ied,JsdB:JedB)) ; OBC%OBC_direction_v(:,:) = OBC_NONE do J=JsdB,JedB ; do i=isd,ied if (OBC%OBC_mask_v(i,J)) OBC%OBC_kind_v(i,J) = OBC_SIMPLE enddo ; enddo diff --git a/src/user/user_initialization.F90 b/src/user/user_initialization.F90 index 33e1090617..1481c8845c 100644 --- a/src/user/user_initialization.F90 +++ b/src/user/user_initialization.F90 @@ -26,8 +26,8 @@ module user_initialization use MOM_io, only : close_file, fieldtype, file_exists use MOM_io, only : open_file, read_data, read_axis_data, SINGLE_FILE use MOM_io, only : write_field, slasher -use MOM_open_boundary, only : ocean_OBC_type, OBC_NONE, OBC_SIMPLE -use MOM_open_boundary, only : OBC_FLATHER_E, OBC_FLATHER_W, OBC_FLATHER_N, OBC_FLATHER_S +use MOM_open_boundary, only : ocean_OBC_type, OBC_NONE, OBC_SIMPLE, OBC_FLATHER +use MOM_open_boundary, only : OBC_DIRECTION_E, OBC_DIRECTION_W, OBC_DIRECTION_N, OBC_DIRECTION_S use MOM_sponge, only : set_up_sponge_field, initialize_sponge, sponge_CS use MOM_tracer_registry, only : tracer_registry_type, add_tracer_OBC_values use MOM_variables, only : thermo_var_ptrs