Go to the documentation of this file.
18 implicit none ;
private
20 #include <MOM_memory.h>
27 real :: tide_flow = 3.0e6
38 character(len=32) :: casename =
"tidal bay"
40 if (
associated(cs))
then
41 call mom_error(warning,
"register_tidal_bay_OBC called with an "// &
42 "associated control structure.")
48 call register_obc(casename, param_file, obc_reg)
57 if (
associated(cs))
then
69 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
70 type(time_type),
intent(in) :: time
74 real :: my_flux, total_area
76 real,
allocatable :: my_area(:,:)
77 character(len=40) :: mdl =
"tidal_bay_set_OBC_data"
78 integer :: i, j, k, itt, is, ie, js, je, isd, ied, jsd, jed, nz, n
79 integer :: isdb, iedb, jsdb, jedb
82 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
83 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
84 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
88 if (.not.
associated(obc))
return
90 allocate(my_area(1:1,js:je))
92 time_sec = time_type_to_real(time)
93 cff = 0.1*sin(2.0*pi*time_sec/(12.0*3600.0))
96 segment => obc%segment(1)
98 do j=segment%HI%jsc,segment%HI%jec ;
do i=segment%HI%IscB,segment%HI%IecB
99 if (obc%segnum_u(i,j) /= obc_none)
then
101 my_area(1,j) = my_area(1,j) + h(i,j,k)*g%US%L_to_m*g%dyCu(i,j)
106 my_flux = - cs%tide_flow*sin(2.0*pi*time_sec/(12.0*3600.0))
108 do n = 1, obc%number_of_segments
109 segment => obc%segment(n)
111 if (.not. segment%on_pe) cycle
113 segment%normal_vel_bt(:,:) = g%US%m_s_to_L_T*my_flux/total_area
114 segment%eta(:,:) = cff
subroutine, public tidal_bay_set_obc_data(OBC, CS, G, h, Time)
This subroutine sets the properties of flow at open boundary conditions.
Wraps the FMS time manager functions.
Provides a transparent vertical ocean grid type and supporting routines.
An overloaded interface to log version information about modules.
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
A structure that can be parsed to read and document run-time parameters.
An overloaded interface to read and log the values of various types of parameters.
logical function, public register_tidal_bay_obc(param_file, CS, OBC_Reg)
Add tidal bay to OBC registry.
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Type to carry basic OBC information needed for updating values.
Describes the vertical ocean grid, including unit conversion factors.
subroutine, public register_obc(name, param_file, Reg)
register open boundary objects for boundary updates.
Controls where open boundary conditions are applied.
subroutine, public tidal_bay_obc_end(CS)
Clean up the tidal bay OBC from registry.
The MOM6 facility to parse input files for runtime parameters.
Provides the ocean grid type.
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Control structure for tidal bay open boundaries.
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Find an accurate and order-invariant sum of distributed 2d or 3d fields.
integer, parameter, public obc_none
Indicates the use of no open boundary.
Configures the model for the "tidal_bay" experiment. tidal_bay = Tidally resonant bay from Zygmunt Ko...
Open boundary segment data structure.
Routines for error handling and I/O management.
Describes the horizontal ocean grid with only dynamic memory arrays.
Ocean grid type. See mom_grid for details.