Skip to content

Commit

Permalink
Split out DOME OBC positions
Browse files Browse the repository at this point in the history 
- Previously, s/r DOME_set_Open_Bdry_Conds() set both the OBC masks
  and the OBC data (inflow conditions).
- DOME_set_OBC_positions() now sets the OBC masks.
- DOME_set_OBC_data() replaces DOME_set_Open_Bdry_Conds() and only sets
  inflow data.
- This is one step in a larger re-factor of OBC code.
  • Loading branch information
@adcroft
adcroft committed Jun 14, 2016
1 parent 89b858f commit 83e50eb
Show file tree
Hide file tree
Showing 2 changed files with 67 additions and 76 deletions.
5 changes: 3 additions & 2 deletions src/initialization/MOM_state_initialization.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@ module MOM_state_initialization
use user_initialization, only : user_init_temperature_salinity
use user_initialization, only : user_set_Open_Bdry_Conds, user_initialize_sponges
use DOME_initialization, only : DOME_initialize_thickness
use DOME_initialization, only : DOME_set_Open_Bdry_Conds
use DOME_initialization, only : DOME_set_OBC_positions, DOME_set_OBC_data
use DOME_initialization, only : DOME_initialize_sponges
use ISOMIP_initialization, only : ISOMIP_initialize_thickness
use ISOMIP_initialization, only : ISOMIP_initialize_sponges
Expand Down Expand Up @@ -440,7 +440,8 @@ subroutine MOM_initialize_state(u, v, h, tv, Time, G, GV, PF, dirs, &
" \t USER - call a user modified routine.", default="file", &
fail_if_missing=.true.)
if (trim(config) == "DOME") then
call DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, PF, tracer_Reg)
call DOME_set_OBC_positions(G, PF, OBC)
call DOME_set_OBC_data(OBC, tv, G, GV, PF, tracer_Reg)
elseif (trim(config) == "USER") then
call user_set_Open_Bdry_Conds(OBC, tv, G, PF, tracer_Reg)
else
Expand Down
138 changes: 64 additions & 74 deletions src/user/DOME_initialization.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,8 @@ module DOME_initialization
public DOME_initialize_topography
public DOME_initialize_thickness
public DOME_initialize_sponges
public DOME_set_Open_Bdry_Conds
public DOME_set_OBC_positions
public DOME_set_OBC_data

contains

Expand Down Expand Up @@ -230,31 +231,69 @@ subroutine DOME_initialize_sponges(G, GV, tv, PF, CSp)
endif

end subroutine DOME_initialize_sponges
! -----------------------------------------------------------------------------

! -----------------------------------------------------------------------------
!> Set the positions of the open boundary needed for the DOME experiment.
subroutine DOME_set_OBC_positions(G, param_file, OBC)
type(ocean_grid_type), intent(in) :: G !< Grid structure.
type(param_file_type), intent(in) :: param_file !< Parameter file handle.
type(ocean_OBC_type), pointer :: OBC !< Open boundary control structure.
! Local variables
character(len=40) :: mod = "DOME_set_OBC_positions" ! This subroutine's name.
integer :: i, j
logical :: any_OBC ! Set to true if any points in this subdomain use OBCs

if (.not.associated(OBC)) allocate(OBC)

call get_param(param_file, mod, "APPLY_OBC_U", OBC%apply_OBC_u, &
"If true, open boundary conditions may be set at some \n"//&
"u-points, with the configuration controlled by OBC_CONFIG", &
default=.false.)
call get_param(param_file, mod, "APPLY_OBC_V", OBC%apply_OBC_v, &
"If true, open boundary conditions may be set at some \n"//&
"v-points, with the configuration controlled by OBC_CONFIG", &
default=.false.)

any_OBC = .false.
if (OBC%apply_OBC_u) then
! Set where u points are determined by OBCs.
!allocate(OBC_mask_u(IsdB:IedB,jsd:jed)) ; OBC_mask_u(:,:) = .false.
call MOM_error(FATAL,"DOME_initialization, DOME_set_OBC_positions: "//&
"APPLY_OBC_U=True is not coded for the DOME experiment")
endif
if (OBC%apply_OBC_v) then
! Set where v points are determined by OBCs.
allocate(OBC%OBC_mask_v(G%isd:G%ied,G%JsdB:G%JedB)) ; OBC%OBC_mask_v(:,:) = .false.
do J=G%JsdB,G%JedB ; do i=G%isd,G%ied
if ((G%geoLonCv(i,J) > 1000.0) .and. (G%geoLonCv(i,J) < 1100.0) .and. &
(abs(G%geoLatCv(i,J) - G%gridLatB(G%JegB)) < 0.1)) then
OBC%OBC_mask_v(i,J) = .true.
any_OBC = .true.
endif
enddo ; enddo
endif
if (.not.any_OBC) then
! If this PE does not have any OBC points then we do not need the mask
OBC%apply_OBC_v = .false.
deallocate(OBC%OBC_mask_v)
endif
end subroutine DOME_set_OBC_positions

!> This subroutine sets the properties of flow at open boundary conditions.
!! This particular example is for the DOME inflow describe in Legg et al. 2006.
subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg)
type(ocean_OBC_type), pointer :: OBC !< This open boundary condition type specifies
!! whether, where, and what open boundary
!! conditions are used.
type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers to any
!! available thermodynamic fields, including potential
!! temperature and salinity or mixed layer density. Absent
!! fields have NULL ptrs.
!! temperature and salinity or mixed layer density. Absent
!! fields have NULL ptrs.
type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure.
type(param_file_type), intent(in) :: param_file !< A structure indicating the open file
!! to parse for model parameter values.
type(tracer_registry_type), pointer :: tr_Reg !< Tracer registry.

logical :: any_OBC ! Set to true if any points in this subdomain use
! open boundary conditions.
logical, pointer, dimension(:,:) :: &
OBC_mask_u => NULL(), & ! These arrays are true at zonal or meridional
OBC_mask_v => NULL() ! velocity points that have prescribed open boundary
! conditions.
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
Expand All @@ -276,7 +315,7 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
! thickness D_edge, in the same units as lat.
real :: Ri_trans ! The shear Richardson number in the transition
! region of the specified shear profile.
character(len=40) :: mod = "DOME_set_Open_Bdry_Conds" ! This subroutine's name.
character(len=40) :: mod = "DOME_set_OBC_data" ! This subroutine's name.
integer :: i, j, k, itt, is, ie, js, je, isd, ied, jsd, jed, nz
integer :: IsdB, IedB, JsdB, JedB

Expand All @@ -289,66 +328,18 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
Ri_trans = 1.0/3.0 ! The shear Richardson number in the transition region
! region of the specified shear profile.

call get_param(param_file, mod, "APPLY_OBC_U", apply_OBC_u, &
"If true, open boundary conditions may be set at some \n"//&
"u-points, with the configuration controlled by OBC_CONFIG", &
default=.false.)
call get_param(param_file, mod, "APPLY_OBC_V", apply_OBC_v, &
"If true, open boundary conditions may be set at some \n"//&
"v-points, with the configuration controlled by OBC_CONFIG", &
default=.false.)

if (apply_OBC_u) then
! Determine where u points are applied.
allocate(OBC_mask_u(IsdB:IedB,jsd:jed)) ; OBC_mask_u(:,:) = .false.
any_OBC = .false.
do j=jsd,jed ; do I=IsdB,IedB
! if (SOME_TEST_FOR_U_OPEN_BCS) then
! OBC_mask_u(I,j) = .true. ; any_OBC = .true.
! endif
enddo ; enddo
if (.not.any_OBC) then
! This processor has no u points at which open boundary conditions are
! to be applied.
apply_OBC_u = .false.
deallocate(OBC_mask_u)
endif
endif
if (apply_OBC_v) then
! Determine where v points are applied.
allocate(OBC_mask_v(isd:ied,JsdB:JedB)) ; OBC_mask_v(:,:) = .false.
any_OBC = .false.
do J=JsdB,JedB ; do i=isd,ied
if ((G%geoLonCv(i,J) > 1000.0) .and. (G%geoLonCv(i,J) < 1100.0) .and. &
(abs(G%geoLatCv(i,J) - G%gridLatB(G%JegB)) < 0.1)) then
OBC_mask_v(i,J) = .true. ; any_OBC = .true.
endif
enddo ; enddo
if (.not.any_OBC) then
! This processor has no v points at which open boundary conditions are
! to be applied.
apply_OBC_v = .false.
deallocate(OBC_mask_v)
endif
endif

if (.not.(apply_OBC_u .or. apply_OBC_v)) return
if (.not.associated(OBC)) return
if (.not.(OBC%apply_OBC_u .or. OBC%apply_OBC_v)) return

if (.not.associated(OBC)) allocate(OBC)

if (apply_OBC_u) then
OBC%apply_OBC_u = .true.
OBC%OBC_mask_u => OBC_mask_u
if (OBC%apply_OBC_u) then
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
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
endif
if (apply_OBC_v) then
OBC%apply_OBC_v = .true.
OBC%OBC_mask_v => OBC_mask_v
if (OBC%apply_OBC_v) then
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
Expand All @@ -357,7 +348,7 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
enddo ; enddo
endif

if (apply_OBC_v) then
if (OBC%apply_OBC_v) then
g_prime_tot = (G%g_Earth/GV%Rho0)*2.0
Def_Rad = sqrt(D_edge*g_prime_tot) / (1.0e-4*1000.0)
tr_0 = (-D_edge*sqrt(D_edge*g_prime_tot)*0.5e3*Def_Rad) * GV%m_to_H
Expand All @@ -380,7 +371,7 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
if (k == nz) tr_k = tr_k + tr_0 * (2.0/(Ri_trans*(2.0+Ri_trans))) * &
log((2.0+Ri_trans)/(2.0-Ri_trans))
do J=JsdB,JedB ; do i=isd,ied
if (OBC_mask_v(i,J)) then
if (OBC%OBC_mask_v(i,J)) then
! This needs to be unneccesarily complicated without symmetric memory.
lon_im1 = 2.0*G%geoLonCv(i,J) - G%geoLonBu(I,J)
! if (isd > IsdB) lon_im1 = G%geoLonBu(I-1,J)
Expand All @@ -394,9 +385,9 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
enddo
endif

if (apply_OBC_u) then
if (OBC%apply_OBC_u) then
do k=1,nz ; do j=jsd,jed ; do I=IsdB,IedB
if (OBC_mask_u(I,j)) then
if (OBC%OBC_mask_u(I,j)) then
! An appropriate expression for the zonal inflow velocities and
! transports should go here.
OBC%uh(I,j,k) = 0.0 * GV%m_to_H ; OBC%u(I,j,k) = 0.0
Expand All @@ -408,7 +399,7 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)

! The inflow values of temperature and salinity also need to be set here if
! these variables are used. The following code is just a naive example.
if (apply_OBC_u .or. apply_OBC_v) then
if (OBC%apply_OBC_u .or. OBC%apply_OBC_v) then
if (associated(tv%S)) then
! In this example, all S inflows have values of 35 psu.
call add_tracer_OBC_values("S", tr_Reg, OBC_inflow=35.0)
Expand All @@ -428,13 +419,13 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
do k=1,nz ; T0(k) = T0(k) + (GV%Rlay(k)-rho_guess(k)) / drho_dT(k) ; enddo
enddo

if (apply_OBC_u) then
if (OBC%apply_OBC_u) then
allocate(OBC_T_u(IsdB:IedB,jsd:jed,nz))
do k=1,nz ; do j=jsd,jed ; do I=IsdB,IedB
OBC_T_u(I,j,k) = T0(k)
enddo ; enddo ; enddo
endif
if (apply_OBC_v) then
if (OBC%apply_OBC_v) then
allocate(OBC_T_v(isd:ied,JsdB:JedB,nz))
do k=1,nz ; do J=JsdB,JedB ; do i=isd,ied
OBC_T_v(i,J,k) = T0(k)
Expand All @@ -445,8 +436,7 @@ subroutine DOME_set_Open_Bdry_Conds(OBC, tv, G, GV, param_file, tr_Reg)
endif
endif

end subroutine DOME_set_Open_Bdry_Conds
! -----------------------------------------------------------------------------
end subroutine DOME_set_OBC_data

!> \class DOME_initialization
!!
Expand Down

0 comments on commit 83e50eb

Please sign in to comment.