Skip to content

Commit

Permalink
Merge branch 'user/ksh/open_bc' into project/obc_segment_fixes
Browse files Browse the repository at this point in the history
  • Loading branch information
adcroft committed Sep 14, 2016
2 parents 36275ae + cb553eb commit 4a23711
Show file tree
Hide file tree
Showing 20 changed files with 1,372 additions and 458 deletions.
106 changes: 51 additions & 55 deletions src/core/MOM_barotropic.F90
Original file line number Diff line number Diff line change
Expand Up @@ -363,9 +363,7 @@ module MOM_barotropic
OBC_mask_v => NULL()
integer, dimension(:,:), pointer :: &
OBC_direction_u => NULL(), &
OBC_direction_v => NULL(), &
OBC_kind_u => NULL(), &
OBC_kind_v => NULL()
OBC_direction_v => NULL()
real, dimension(:,:), pointer :: &
Cg_u => NULL(), & ! The external wave speed at u-points, in m s-1.
Cg_v => NULL(), & ! The external wave speed at u-points, in m s-1.
Expand Down Expand Up @@ -2376,11 +2374,11 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans,

if (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_SIMPLE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%specified) then
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) then
elseif (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) 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
Expand Down Expand Up @@ -2422,7 +2420,7 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans,
endif
endif

if (BT_OBC%OBC_kind_u(I,j) /= OBC_SIMPLE) then
if (.not. OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%specified) then
if (use_BT_cont) then
uhbt(I,j) = find_uhbt(vel_trans,BTCL_u(I,j)) + uhbt0(I,j)
else
Expand All @@ -2436,11 +2434,11 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans,

if (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_SIMPLE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%specified) then
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) then
elseif (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) 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
Expand Down Expand Up @@ -2490,7 +2488,7 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans,
endif
endif

if (BT_OBC%OBC_kind_v(i,J) /= OBC_SIMPLE) then
if (.not. OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%specified) then
if (use_BT_cont) then
vhbt(i,J) = find_vhbt(vel_trans,BTCL_v(i,J)) + vhbt0(i,J)
else
Expand Down Expand Up @@ -2538,50 +2536,54 @@ 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_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)
do j=js,je ; do I=is-1,ie ; if (OBC%OBC_segment_u(I,j) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) 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

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_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)
do J=js-1,je ; do i=is,ie ; if (OBC%OBC_segment_v(i,J) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) 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
Expand Down Expand Up @@ -2640,7 +2642,6 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D
allocate(BT_OBC%ubt_outer(isdw-1:iedw,jsdw:jedw)) ; BT_OBC%ubt_outer(:,:) = 0.0
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
Expand All @@ -2649,13 +2650,11 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D
allocate(BT_OBC%vbt_outer(isdw:iedw,jsdw-1:jedw)) ; BT_OBC%vbt_outer(:,:) = 0.0
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
Expand All @@ -2664,7 +2663,7 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D
enddo ; enddo ; enddo
endif
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_SIMPLE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%specified) then
if (use_BT_cont) then
BT_OBC%ubt_outer(I,j) = uhbt_to_ubt(BT_OBC%uhbt(I,j),BTCL_u(I,j))
else
Expand All @@ -2691,7 +2690,6 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D
if (associated(OBC%OBC_mask_v)) then
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
Expand All @@ -2701,7 +2699,7 @@ subroutine set_up_BT_OBC(OBC, eta, BT_OBC, G, GV, MS, halo, use_BT_cont, Datu, D
endif

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_SIMPLE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%specified) then
if (use_BT_cont) then
BT_OBC%vbt_outer(i,J) = vhbt_to_vbt(BT_OBC%vhbt(i,J),BTCL_v(i,J))
else
Expand Down Expand Up @@ -2741,7 +2739,6 @@ subroutine destroy_BT_OBC(BT_OBC)
type(BT_OBC_type), intent(inout) :: 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)
Expand All @@ -2750,7 +2747,6 @@ subroutine destroy_BT_OBC(BT_OBC)
deallocate(BT_OBC%eta_outer_u)

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)
Expand Down
88 changes: 88 additions & 0 deletions src/core/MOM_boundary_update.F90
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@
! This file is part of MOM6. See LICENSE.md for the license.
!> Controls where open boundary conditions are applied
module MOM_boundary_update

! This file is part of MOM6. See LICENSE.md for the license.

use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end, CLOCK_ROUTINE
use MOM_diag_mediator, only : diag_ctrl, time_type
use MOM_domains, only : pass_var, pass_vector
use MOM_domains, only : To_All, SCALAR_PAIR, CGRID_NE
use MOM_error_handler, only : MOM_mesg, MOM_error, FATAL, WARNING
use MOM_file_parser, only : get_param, log_version, param_file_type, log_param
use MOM_grid, only : ocean_grid_type
use MOM_dyn_horgrid, only : dyn_horgrid_type
use MOM_io, only : EAST_FACE, NORTH_FACE
use MOM_io, only : slasher, read_data
use MOM_open_boundary, only : ocean_obc_type
use MOM_tracer_registry, only : add_tracer_OBC_values, tracer_registry_type
use MOM_variables, only : thermo_var_ptrs
use tidal_bay_initialization, only : tidal_bay_set_OBC_data

implicit none ; private

#include <MOM_memory.h>

public update_OBC_data

integer :: id_clock_pass

character(len=40) :: mod = "MOM_boundary_update" ! This module's name.
! This include declares and sets the variable "version".
#include "version_variable.h"

contains

!> Calls appropriate routine to update the open boundary conditions.
subroutine update_OBC_data(OBC, G, h, Time)
type(ocean_grid_type), intent(in) :: G !< Ocean grid structure
real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(inout) :: h !< layer thickness
type(ocean_OBC_type), pointer :: OBC !< Open boundary structure
type(time_type), intent(in) :: Time !< Model time
! Local variables
logical :: read_OBC_eta = .false.
logical :: read_OBC_uv = .false.
logical :: read_OBC_TS = .false.
integer :: i, j, k, itt, is, ie, js, je, isd, ied, jsd, jed, nz
integer :: isd_off, jsd_off
integer :: IsdB, IedB, JsdB, JedB
character(len=40) :: mod = "set_Flather_Bdry_Conds" ! This subroutine's name.
character(len=200) :: filename, OBC_file, inputdir ! Strings for file/path

real :: temp_u(G%domain%niglobal+1,G%domain%njglobal)
real :: temp_v(G%domain%niglobal,G%domain%njglobal+1)

real, pointer, dimension(:,:,:) :: &
OBC_T_u => NULL(), & ! These arrays should be allocated and set to
OBC_T_v => NULL(), & ! specify the values of T and S that should come
OBC_S_u => NULL(), &
OBC_S_v => NULL()

is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = G%ke
isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed
IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB

if (OBC%OBC_values_config == "tidal_bay") then
call tidal_bay_set_OBC_data(OBC, G, h, Time)
endif

end subroutine update_OBC_data

!> \namespace mom_boundary_update
!! This module updates the open boundary arrays when time-varying.
!! It caused a circular dependency with the tidal_bay setup when
!! MOM_open_boundary.
!!
!! A small fragment of the grid is shown below:
!!
!! j+1 x ^ x ^ x At x: q, CoriolisBu
!! j+1 > o > o > At ^: v, tauy
!! j x ^ x ^ x At >: u, taux
!! j > o > o > At o: h, bathyT, buoy, tr, T, S, Rml, ustar
!! j-1 x ^ x ^ x
!! i-1 i i+1 At x & ^:
!! i i+1 At > & o:
!!
!! The boundaries always run through q grid points (x).

end module MOM_boundary_update
Loading

0 comments on commit 4a23711

Please sign in to comment.