Skip to content

Commit

Permalink
Fixes needed to "not" extend OBCs beyond domain.
Browse files Browse the repository at this point in the history 
- Added %Flather logical to segment type to differentiate Flather (for
  barotropic) from Orlanski (baroclinic) boundary conditions.
- Bugfix: Barotropic-direction south had wrong sign for setting the
  tangential ubt. Does not appear to change answers?
- Changed syntax of OBC_SEGMENT_%%% string to allow multiple boundary
  conditions type per segment, e.g. FLATHER,ORLANSKI
- zonal_mass_flux() and meridional_mass_flux() in continuity_PPM now
  use %Flather instead of %radiation just before reconciling fluxes with
  the barotropic fluxes.
- Fixed need to extend boundaries beyond domain in run-time parameters.
  Extended domains used to be needed to avoid masking the corner vorticity
  points. Added code to detect corner joins of segments and handle explicitly.
- Removed code to interpret I=-1:N+1 (extend boundaries no longer needed).
- Use "hack" in setup_*_point_obc() to extend segments by one point (in code)
  to recover "corner" behavior of old OBCs.
- Added back setting of tangential velocity components in continuity_PPM.
- set_Flather_data() has commented out code to set "h" across corners when
  segments are not extended.
- Cleaned up setting of "direction" in setup_u_point_obc() and
  setup_v_point_obc().
- Simplied using of directions/flags in continuity_PPM.
- No answer changes with consistent switch to SEGMENT parameters.

TODO:
- Allow coexistence of specified and radiative BCs within zonal_mass_flux()
  and meridional_mass_flux().
- Understand why the above bugfix did not change answers?
- Figure out why extending segments to recover corner behavior is necessary:
  LAPLACIAN=BIHARMONIC=False still failed, suggesting PV term needs examination.
  • Loading branch information
@adcroft
adcroft committed Sep 14, 2016
1 parent 0aef287 commit fcdb84e
Show file tree
Hide file tree
Showing 4 changed files with 283 additions and 217 deletions.
10 changes: 5 additions & 5 deletions src/core/MOM_barotropic.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2378,7 +2378,7 @@ 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 (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) then
elseif (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%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
Expand Down Expand Up @@ -2411,7 +2411,7 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans,
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
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)
Expand All @@ -2438,7 +2438,7 @@ 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 (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) then
elseif (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%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
Expand Down Expand Up @@ -2537,7 +2537,7 @@ 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 (OBC%OBC_segment_u(I,j) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%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
Expand All @@ -2564,7 +2564,7 @@ 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 (OBC%OBC_segment_v(i,J) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%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
Expand Down
72 changes: 46 additions & 26 deletions src/core/MOM_continuity_PPM.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -191,22 +191,26 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC,
call cpu_clock_end(id_clock_update)

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
do k=1,nz
do j=LB%jsh,LB%jeh ; do I=LB%ish,LB%ieh+1
if (OBC%OBC_segment_u(I-1,j) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I-1,j))%radiation &
.and. (OBC%OBC_direction_u(I-1,j) == OBC_DIRECTION_E)) &
if (OBC%OBC_direction_u(I-1,j) == OBC_DIRECTION_E) &
h(i,j,k) = h_input(i-1,j,k)
endif
enddo
do i=LB%ish-1,LB%ieh
if (OBC%OBC_segment_u(I,j) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation &
.and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W)) &
if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) &
h(i,j,k) = h_input(i+1,j,k)
endif
enddo
enddo ; enddo
enddo ; enddo
do J=LB%jsh-1,LB%jeh ; do i=LB%ish-1,LB%ieh+1
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%Flather .and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_E)) &
v(i,J,k) = v(i-1,J,k)
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%Flather .and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_W)) &
v(i,J,k) = v(i+1,J,k)
enddo ; enddo
enddo
endif
LB%ish = G%isc ; LB%ieh = G%iec ; LB%jsh = G%jsc ; LB%jeh = G%jec

Expand All @@ -228,18 +232,22 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC,
do k=1,nz
do J=LB%jsh,LB%jeh+1 ; do i=LB%ish-1,LB%ieh+1
if (OBC%OBC_segment_v(i,J-1) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J-1))%radiation &
.and. (OBC%OBC_direction_v(i,J-1) == OBC_DIRECTION_N)) &
if (OBC%OBC_direction_v(i,J-1) == OBC_DIRECTION_N) &
h(i,j,k) = h_input(i,j-1,k)
endif
enddo ; enddo
do J=LB%jsh-1,LB%jeh ; do i=LB%ish-1,LB%ieh+1
if (OBC%OBC_segment_v(i,J) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation &
.and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S)) &
if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) &
h(i,j,k) = h_input(i,j+1,k)
endif
enddo ; enddo
do j=LB%jsh-1,LB%jeh+1 ; do I=LB%ish-1,LB%ieh
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%Flather .and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_N)) &
u(I,j,k) = u(I,j-1,k)
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%Flather .and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_S)) &
u(I,j,k) = u(I,j+1,k)
enddo ; enddo
enddo
endif
else ! .not. x_first
Expand All @@ -261,18 +269,22 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC,
do k=1,nz
do J=LB%jsh,LB%jeh+1 ; do i=LB%ish-1,LB%ieh+1
if (OBC%OBC_segment_v(i,J-1) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J-1))%radiation &
.and. (OBC%OBC_direction_v(i,J-1) == OBC_DIRECTION_N)) &
if (OBC%OBC_direction_v(i,J-1) == OBC_DIRECTION_N) &
h(i,j,k) = h_input(i,j-1,k)
endif
enddo ; enddo
do J=LB%jsh-1,LB%jeh ; do i=LB%ish-1,LB%ieh+1
if (OBC%OBC_segment_v(i,J) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation &
.and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S)) &
if (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_S) &
h(i,j,k) = h_input(i,j+1,k)
endif
enddo ; enddo
do j=LB%jsh-1,LB%jeh+1 ; do I=LB%ish-1,LB%ieh
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%Flather .and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_N)) &
u(I,j,k) = u(I,j-1,k)
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%Flather .and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_S)) &
u(I,j,k) = u(I,j+1,k)
enddo ; enddo
enddo
endif

Expand All @@ -292,22 +304,26 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, CS, uhbt, vhbt, OBC,
call cpu_clock_end(id_clock_update)

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
do k=1,nz
do j=LB%jsh,LB%jeh ; do I=LB%ish,LB%ieh+1
if (OBC%OBC_segment_u(I-1,j) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I-1,j))%radiation &
.and. (OBC%OBC_direction_u(I-1,j) == OBC_DIRECTION_E)) &
if (OBC%OBC_direction_u(I-1,j) == OBC_DIRECTION_E) &
h(i,j,k) = h_input(i-1,j,k)
endif
enddo
do i=LB%ish-1,LB%ieh
if (OBC%OBC_segment_u(I,j) /= OBC_NONE) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation &
.and. (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W)) &
if (OBC%OBC_direction_u(I,j) == OBC_DIRECTION_W) &
h(i,j,k) = h_input(i+1,j,k)
endif
enddo
enddo ; enddo
enddo ; enddo
do J=LB%jsh-1,LB%jeh ; do i=LB%ish-1,LB%ieh+1
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%Flather .and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_E)) &
v(i,J,k) = v(i-1,J,k)
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%Flather .and. (OBC%OBC_direction_v(i,J) == OBC_DIRECTION_W)) &
v(i,J,k) = v(i+1,J,k)
enddo ; enddo
enddo
endif
endif

Expand Down Expand Up @@ -524,8 +540,10 @@ subroutine zonal_mass_flux(u, h_in, uh, dt, G, GV, CS, LB, uhbt, OBC, &
(OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%specified))
if (.not.do_i(I)) any_simple_OBC = .true.
enddo ; else if (apply_OBC_flather) then ; do I=ish-1,ieh
! This is a tangential condition and is needed for unknown reasons and
! probably implies that we made a calculation elsewhere that we should not have.
do_i(I) = .not.(OBC%OBC_mask_u(I,j) .and. &
(OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) .and. &
(OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%Flather) .and. &
((OBC%OBC_direction_u(I,j) == OBC_DIRECTION_N) .or. &
(OBC%OBC_direction_u(I,j) == OBC_DIRECTION_S)))
enddo ; else ; do I=ish-1,ieh
Expand Down Expand Up @@ -1290,8 +1308,10 @@ subroutine meridional_mass_flux(v, h_in, vh, dt, G, GV, CS, LB, vhbt, OBC, &
(OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%specified))
if (.not.do_i(i)) any_simple_OBC = .true.
enddo ; else if (apply_OBC_flather) then ; do i=ish,ieh
! This is a tangential condition and is needed for unknown reasons and
! probably implies that we made a calculation elsewhere that we should not have.
do_i(i) = .not.(OBC%OBC_mask_v(i,J) .and. &
(OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) .and. &
(OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%Flather) .and. &
((OBC%OBC_direction_v(i,J) == OBC_DIRECTION_E) .or. &
(OBC%OBC_direction_v(i,J) == OBC_DIRECTION_W)))
enddo ; else ; do i=ish,ieh
Expand Down
10 changes: 5 additions & 5 deletions src/core/MOM_legacy_barotropic.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2236,7 +2236,7 @@ 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 (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) then
elseif (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%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
Expand Down Expand Up @@ -2269,7 +2269,7 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans,
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
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)
Expand All @@ -2296,7 +2296,7 @@ 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 (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) then
elseif (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%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
Expand Down Expand Up @@ -2395,7 +2395,7 @@ 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 (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%radiation) then
if (OBC%OBC_segment_list(OBC%OBC_segment_u(I,j))%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
Expand All @@ -2422,7 +2422,7 @@ 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 (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%radiation) then
if (OBC%OBC_segment_list(OBC%OBC_segment_v(i,J))%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
Expand Down
Loading

0 comments on commit fcdb84e

Please sign in to comment.