Skip to content

Commit

Permalink
Some updates, not compiling.
Browse files Browse the repository at this point in the history
  • Loading branch information
kshedstrom committed Jan 12, 2017
1 parent e3dab80 commit 7d28845
Show file tree
Hide file tree
Showing 4 changed files with 73 additions and 27 deletions.
66 changes: 39 additions & 27 deletions src/core/MOM_open_boundary.F90
Original file line number Diff line number Diff line change
Expand Up @@ -81,11 +81,11 @@ module MOM_open_boundary
real :: Tnudge_in !< Inverse nudging timescale on inflow (1/s).
real :: Tnudge_out !< Inverse nudging timescale on outflow (1/s).
logical :: on_pe !< true if segment is located in the computational domain
integer, dimension(4) :: global_indices ! is,ie,js,je list of global indices for segment
real, pointer, dimension(:,:) :: Cg=>NULL() !<The external gravity
!<wave speed (m -s) at OBC-points.
real, pointer, dimension(:,:) :: Htot=>NULL() !<The total column thickness (m) at OBC-points.
real, pointer, dimension(:,:,:) :: h=>NULL() !<The cell thickness (m) at OBC-points.
real, pointer, dimension(:,:,:) :: un=>NULL() !<The grid-oriented normal layer velocity (m s-1).
real, pointer, dimension(:,:,:) :: unh=>NULL() !<The grid-oriented normal layer transports (m3 s-1).
real, pointer, dimension(:,:) :: unhbt=>NULL() !<The vertically summed normal layer transports (m3 s-1).
real, pointer, dimension(:,:) :: unbt=>NULL() !<The vertically-averaged normal velocity (m s-1).
Expand Down Expand Up @@ -390,7 +390,6 @@ subroutine initialize_segment_data(G, OBC, PF)

return


end subroutine initialize_segment_data

!< Define indices for segment and store in hor_index_type
Expand Down Expand Up @@ -469,8 +468,7 @@ subroutine setup_u_point_obc(OBC, G, segment_str, l_seg)
call parse_segment_str(G%ieg, G%jeg, segment_str, I_obc, Js_obc, Je_obc, action_str )

call setup_segment_indices(G, OBC%OBC_segment_number(l_seg),I_obc,I_obc,Js_obc,Je_obc)

OBC%OBC_segment_number(l_seg)%global_indices = (/I_obc,I_obc,Js_obc,Je_obc/)
call allocate_OBC_segment_data(G, OBC%OBC_segment_number(l_seg))

I_obc = I_obc - G%idg_offset ! Convert to local tile indices on this tile
Js_obc = Js_obc - G%jdg_offset ! Convert to local tile indices on this tile
Expand All @@ -492,8 +490,6 @@ subroutine setup_u_point_obc(OBC, G, segment_str, l_seg)
j=js_obc;js_obc=je_obc;je_obc=j
endif

OBC%OBC_segment_number(l_seg)%global_indices = (/I_obc+G%idg_offset,I_obc+G%idg_offset,Js_obc+G%jdg_offset,Je_obc+G%jdg_offset/)

OBC%OBC_segment_number(l_seg)%on_pe = .false.

do a_loop = 1,5 ! up to 5 options available
Expand Down Expand Up @@ -532,7 +528,7 @@ subroutine setup_u_point_obc(OBC, G, segment_str, l_seg)
elseif (trim(action_str(a_loop)) == 'SIMPLE') then
OBC%OBC_segment_number(l_seg)%specified = .true.
OBC%OBC_segment_number(l_seg)%values_needed = .true.
! OBC%update_OBC = .true.
OBC%update_OBC = .true.
OBC%specified_u_BCs_exist_globally = .true. ! This avoids deallocation
! Hack to undo the hack above for SIMPLE BCs
Js_obc = Js_obc + 1 ; Je_obc = Je_obc - 1
Expand Down Expand Up @@ -598,6 +594,7 @@ subroutine setup_v_point_obc(OBC, G, segment_str, l_seg)
call parse_segment_str(G%ieg, G%jeg, segment_str, J_obc, Is_obc, Ie_obc, action_str )

call setup_segment_indices(G, OBC%OBC_segment_number(l_seg),Is_obc,Ie_obc,J_obc,J_obc)
call allocate_OBC_segment_data(G, OBC%OBC_segment_number(l_seg))

J_obc = J_obc - G%jdg_offset ! Convert to local tile indices on this tile
Is_obc = Is_obc - G%idg_offset ! Convert to local tile indices on this tile
Expand All @@ -618,8 +615,6 @@ subroutine setup_v_point_obc(OBC, G, segment_str, l_seg)
i=is_obc;is_obc=ie_obc;ie_obc=i
endif

OBC%OBC_segment_number(l_seg)%global_indices = (/Is_obc+G%idg_offset,Ie_obc+G%idg_offset,J_obc+G%jdg_offset,J_obc+G%jdg_offset/)

OBC%OBC_segment_number(l_seg)%on_pe = .false.

! if (Ie_obc>Is_obc) OBC%OBC_segment_number(l_seg)%direction = OBC_DIRECTION_S
Expand Down Expand Up @@ -660,7 +655,7 @@ subroutine setup_v_point_obc(OBC, G, segment_str, l_seg)
elseif (trim(action_str(a_loop)) == 'SIMPLE') then
OBC%OBC_segment_number(l_seg)%specified = .true.
OBC%OBC_segment_number(l_seg)%values_needed = .true.
! OBC%update_OBC = .true.
OBC%update_OBC = .true.
OBC%specified_v_BCs_exist_globally = .true. ! This avoids deallocation
! Hack to undo the hack above for SIMPLE BCs
Is_obc = Is_obc + 1 ; Ie_obc = Ie_obc - 1
Expand Down Expand Up @@ -1721,7 +1716,7 @@ subroutine fill_OBC_halos(G, GV, OBC, tv, h, tracer_Reg)
isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed
nz=G%ke

if (.not.associated(tv%T)) return
if (.not. associated(OBC) .or. .not. associated(tv%T)) return

call pass_var(tv%T,G%Domain)
call pass_var(tv%S,G%Domain)
Expand Down Expand Up @@ -1828,6 +1823,38 @@ subroutine fill_OBC_halos(G, GV, OBC, tv, h, tracer_Reg)
endif
end subroutine fill_OBC_halos

!> Allocate segment data fields
subroutine allocate_OBC_segment_data(G, segment)
type(dyn_horgrid_type), intent(in) :: G !< Ocean grid structure
type(OBC_segment_type), intent(inout) :: segment !< Open boundary segment
! Local variables
integer :: isd, ied, jsd, jed
integer :: IsdB, IedB, JsdB, JedB
character(len=40) :: mod = "allocate_OBC_segment_data" ! This subroutine's name.

isd = segment%HI%isd ; ied = segment%HI%ied
jsd = segment%HI%jsd ; jed = segment%HI%jed
IsdB = segment%HI%IsdB ; IedB = segment%HI%IedB
JsdB = segment%HI%JsdB ; JedB = segment%HI%JedB

if (.not. associated(seg)) return

if (.not. segment%on_pe) return ! continue to next segment if not in computational domain

if (.not. ASSOCIATED(segment%Cg)) then ! finishing allocating storage for segments
allocate(segment%Cg(isd:ied,jsd:jed)); segment%Cg(:,:)=0.
allocate(segment%Htot(isd:ied,jsd:jed)); segment%Htot(:,:)=0.0
allocate(segment%h(isd:ied,jsd:jed,G%ke)); segment%h(:,:,:)=0.0
allocate(segment%un(isd:ied,jsd:jed,G%ke)); segment%un(:,:,:)=0.0
allocate(segment%unh(isd:ied,jsd:jed,G%ke)); segment%unh(:,:,:)=0.0
allocate(segment%unhbt(isd:ied,jsd:jed)); segment%unhbt(:,:)=0.0
allocate(segment%unbt(isd:ied,jsd:jed)); segment%unbt(:,:)=0.0
allocate(segment%tangent_vel(IsdB:IedB,JsdB:JedB,G%ke)); segment%tangent_vel(:,:,:)=0.0
allocate(segment%tangent_vel_bt(IsdB:IedB,JsdB:JedB)); segment%tangent_vel_bt(:,:)=0.0
allocate(segment%eta(isd:ied,jsd:jed)); segment%eta(:,:)=0.0
endif
enddo
end subroutine allocate_OBC_segment_data

subroutine update_OBC_segment_data(G, GV, OBC, tv, h, Time)

Expand Down Expand Up @@ -1860,10 +1887,8 @@ subroutine update_OBC_segment_data(G, GV, OBC, tv, h, Time)
IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB
nz=G%ke


if (.not. associated(OBC)) return


do n = 1, OBC%number_of_segments
segment => OBC%OBC_segment_number(n)

Expand All @@ -1885,20 +1910,6 @@ subroutine update_OBC_segment_data(G, GV, OBC, tv, h, Time)
is_obc=is_obc+1
endif


! ni_seg = max(segment%Ie_obc-segment%Is_obc,1)
! nj_seg = max(segment%Je_obc-segment%Js_obc,1)

if (.not. ASSOCIATED(segment%Cg)) then ! finishing allocating storage for segments
allocate(segment%Cg(is_obc:ie_obc,js_obc:je_obc));segment%Cg(:,:)=0.
allocate(segment%Htot(is_obc:ie_obc,js_obc:je_obc));segment%Htot(:,:)=0.0
allocate(segment%h(is_obc:ie_obc,js_obc:je_obc,G%ke));segment%h(:,:,:)=0.0
allocate(segment%unh(is_obc:ie_obc,js_obc:je_obc,G%ke));segment%unh(:,:,:)=0.0
allocate(segment%unhbt(is_obc:ie_obc,js_obc:je_obc));segment%unhbt(:,:)=0.0
allocate(segment%unbt(is_obc:ie_obc,js_obc:je_obc));segment%unbt(:,:)=0.0
allocate(segment%eta(is_obc:ie_obc,js_obc:je_obc));segment%eta(:,:)=0.0
endif

! do j=jsd,jed ; do I=isd,ied-1
! if (segment%direction == OBC_DIRECTION_E .and. OBC%OBC_segment_u(I,j) /= OBC_NONE ) then
! do k=1,nz
Expand Down Expand Up @@ -2043,6 +2054,7 @@ subroutine update_OBC_segment_data(G, GV, OBC, tv, h, Time)
do i=is_obc,ie_obc
segment%unhbt(i,j) = 0.0
do k=1,G%ke
segment%un(i,j,k) = segment%field(m)%buffer_dst(i,j,k)
segment%unh(i,j,k) = segment%field(m)%buffer_dst(i,j,k)*segment%h(i,j,k)
segment%unhbt(i,j)= segment%unhbt(i,j)+segment%unh(i,j,k)
enddo
Expand Down
1 change: 1 addition & 0 deletions src/initialization/MOM_state_initialization.F90
Original file line number Diff line number Diff line change
Expand Up @@ -449,6 +449,7 @@ subroutine MOM_initialize_state(u, v, h, tv, Time, G, GV, PF, dirs, &
if (open_boundary_query(OBC, apply_specified_OBC=.true.)) then
if (trim(config) == "DOME") then
call DOME_set_OBC_data(OBC, tv, G, GV, PF, tracer_Reg)
OBC%update_OBC = .false.
elseif (trim(config) == "USER") then
call user_set_OBC_data(OBC, tv, G, PF, tracer_Reg)
elseif (.not. trim(config) == "none") then
Expand Down
21 changes: 21 additions & 0 deletions src/user/DOME_initialization.F90
Original file line number Diff line number Diff line change
Expand Up @@ -269,6 +269,10 @@ subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg)
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
type(OBC_segment_type), pointer :: segment
integer :: ni_seg, nj_seg ! number of src gridpoints along the segments
integer :: i2, j2 ! indices for referencing local domain array
integer :: ishift, jshift ! offsets for staggered locations

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
Expand All @@ -288,6 +292,13 @@ subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg)
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

if (OBC%number_of_segments .ne. 1) then
print *, 'Error in DOME OBC segment setup'
return !!! Need a better error message here
endif
segment => OBC%OBC_segment_number(1)
if (.not. segment%on_pe) cycle ! continue to next segment if not in computational domain

do k=1,nz
rst = -1.0
if (k>1) rst = -1.0 + (real(k-1)-0.5)/real(nz-1)
Expand All @@ -305,6 +316,7 @@ subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg)
(2.0 - Ri_trans))
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))
! Old way
do J=JsdB,JedB ; do i=isd,ied
if (OBC%OBC_segment_v(i,J) /= OBC_NONE) then
! This needs to be unneccesarily complicated without symmetric memory.
Expand All @@ -317,6 +329,15 @@ subroutine DOME_set_OBC_data(OBC, tv, G, GV, param_file, tr_Reg)
OBC%vh(i,J,k) = 0.0 ; OBC%v(i,J,k) = 0.0
endif
enddo ; enddo
! New way
isd = segment%HI%isd ; ied = segment%HI%ied
JsdB = segment%HI%JsdB ; JedB = segment%HI%JedB
do J=JsdB,JedB ; do i=isd,ied
lon_im1 = 2.0*G%geoLonCv(i,J) - G%geoLonBu(I,J)
segment%unh(i,J,k) = tr_k * (exp(-2.0*(lon_im1 - 1000.0)/Def_Rad) -&
exp(-2.0*(G%geoLonBu(I,J) - 1000.0)/Def_Rad))
segment%un(i,J,k) = v_k * exp(-2.0*(G%geoLonCv(i,J) - 1000.0)/Def_Rad)
enddo ; enddo
enddo

! The inflow values of temperature and salinity also need to be set here if
Expand Down
12 changes: 12 additions & 0 deletions src/user/tidal_bay_initialization.F90
Original file line number Diff line number Diff line change
Expand Up @@ -75,6 +75,7 @@ subroutine tidal_bay_set_OBC_data(OBC, G, h, Time)
enddo ; enddo
my_flux = -tide_flow*SIN(2.0*PI*time_sec/(12.0*3600.0))

! Old way
do j=jsd,jed ; do I=IsdB,IedB
if (OBC%OBC_segment_u(I,j) /= OBC_NONE) then
OBC%eta_outer_u(I,j) = cff
Expand All @@ -88,6 +89,17 @@ subroutine tidal_bay_set_OBC_data(OBC, G, h, Time)
endif
enddo ; enddo

! New way
do n = 1, OBC%number_of_segments
segment => OBC%OBC_segment_number(n)

if (.not. segment%on_pe) cycle ! continue to next segment if not in computational domain

segment%unbt(:,:) = my_flux/my_area
segment%eta(:,:) = cff

enddo ! end segment loop

end subroutine tidal_bay_set_OBC_data

!> \class tidal_bay_Initialization
Expand Down

0 comments on commit 7d28845

Please sign in to comment.