MOM6
tidal_bay_initialization.F90
Go to the documentation of this file.
1 !> Configures the model for the "tidal_bay" experiment.
2 !! tidal_bay = Tidally resonant bay from Zygmunt Kowalik's class on tides.
4 
5 ! This file is part of MOM6. See LICENSE.md for the license.
6 
7 use mom_coms, only : reproducing_sum
9 use mom_error_handler, only : mom_mesg, mom_error, fatal, warning, is_root_pe
11 use mom_grid, only : ocean_grid_type
16 use mom_time_manager, only : time_type, time_type_to_real
17 
18 implicit none ; private
19 
20 #include <MOM_memory.h>
21 
24 
25 !> Control structure for tidal bay open boundaries.
26 type, public :: tidal_bay_obc_cs ; private
27  real :: tide_flow = 3.0e6 !< Maximum tidal flux.
28 end type tidal_bay_obc_cs
29 
30 contains
31 
32 !> Add tidal bay to OBC registry.
33 function register_tidal_bay_obc(param_file, CS, OBC_Reg)
34  type(param_file_type), intent(in) :: param_file !< parameter file.
35  type(tidal_bay_obc_cs), pointer :: cs !< tidal bay control structure.
36  type(obc_registry_type), pointer :: obc_reg !< OBC registry.
37  logical :: register_tidal_bay_obc
38  character(len=32) :: casename = "tidal bay" !< This case's name.
39 
40  if (associated(cs)) then
41  call mom_error(warning, "register_tidal_bay_OBC called with an "// &
42  "associated control structure.")
43  return
44  endif
45  allocate(cs)
46 
47  ! Register the open boundaries.
48  call register_obc(casename, param_file, obc_reg)
49  register_tidal_bay_obc = .true.
50 
51 end function register_tidal_bay_obc
52 
53 !> Clean up the tidal bay OBC from registry.
54 subroutine tidal_bay_obc_end(CS)
55  type(tidal_bay_obc_cs), pointer :: cs !< tidal bay control structure.
56 
57  if (associated(cs)) then
58  deallocate(cs)
59  endif
60 end subroutine tidal_bay_obc_end
61 
62 !> This subroutine sets the properties of flow at open boundary conditions.
63 subroutine tidal_bay_set_obc_data(OBC, CS, G, h, Time)
64  type(ocean_obc_type), pointer :: obc !< This open boundary condition type specifies
65  !! whether, where, and what open boundary
66  !! conditions are used.
67  type(tidal_bay_obc_cs), pointer :: cs !< tidal bay control structure.
68  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure.
69  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< layer thickness.
70  type(time_type), intent(in) :: time !< model time.
71 
72  ! The following variables are used to set up the transport in the tidal_bay example.
73  real :: time_sec, cff
74  real :: my_flux, total_area
75  real :: pi
76  real, allocatable :: my_area(:,:)
77  character(len=40) :: mdl = "tidal_bay_set_OBC_data" ! This subroutine's name.
78  integer :: i, j, k, itt, is, ie, js, je, isd, ied, jsd, jed, nz, n
79  integer :: isdb, iedb, jsdb, jedb
80  type(obc_segment_type), pointer :: segment => null()
81 
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
85 
86  pi = 4.0*atan(1.0)
87 
88  if (.not.associated(obc)) return
89 
90  allocate(my_area(1:1,js:je))
91 
92  time_sec = time_type_to_real(time)
93  cff = 0.1*sin(2.0*pi*time_sec/(12.0*3600.0))
94  my_area=0.0
95  my_flux=0.0
96  segment => obc%segment(1)
97 
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
100  do k=1,nz
101  my_area(1,j) = my_area(1,j) + h(i,j,k)*g%US%L_to_m*g%dyCu(i,j)
102  enddo
103  endif
104  enddo ; enddo
105  total_area = reproducing_sum(my_area)
106  my_flux = - cs%tide_flow*sin(2.0*pi*time_sec/(12.0*3600.0))
107 
108  do n = 1, obc%number_of_segments
109  segment => obc%segment(n)
110 
111  if (.not. segment%on_pe) cycle
112 
113  segment%normal_vel_bt(:,:) = g%US%m_s_to_L_T*my_flux/total_area
114  segment%eta(:,:) = cff
115 
116  enddo ! end segment loop
117 
118 end subroutine tidal_bay_set_obc_data
119 
120 end module tidal_bay_initialization
tidal_bay_initialization::tidal_bay_set_obc_data
subroutine, public tidal_bay_set_obc_data(OBC, CS, G, h, Time)
This subroutine sets the properties of flow at open boundary conditions.
Definition: tidal_bay_initialization.F90:64
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_dyn_horgrid
Contains a shareable dynamic type for describing horizontal grids and metric data and utilty routines...
Definition: MOM_dyn_horgrid.F90:3
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
tidal_bay_initialization::register_tidal_bay_obc
logical function, public register_tidal_bay_obc(param_file, CS, OBC_Reg)
Add tidal bay to OBC registry.
Definition: tidal_bay_initialization.F90:34
mom_coms
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Definition: MOM_coms.F90:3
mom_open_boundary::obc_registry_type
Type to carry basic OBC information needed for updating values.
Definition: MOM_open_boundary.F90:294
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
mom_open_boundary::register_obc
subroutine, public register_obc(name, param_file, Reg)
register open boundary objects for boundary updates.
Definition: MOM_open_boundary.F90:3899
mom_open_boundary
Controls where open boundary conditions are applied.
Definition: MOM_open_boundary.F90:2
tidal_bay_initialization::tidal_bay_obc_end
subroutine, public tidal_bay_obc_end(CS)
Clean up the tidal bay OBC from registry.
Definition: tidal_bay_initialization.F90:55
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_open_boundary::ocean_obc_type
Open-boundary data.
Definition: MOM_open_boundary.F90:195
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
tidal_bay_initialization::tidal_bay_obc_cs
Control structure for tidal bay open boundaries.
Definition: tidal_bay_initialization.F90:26
mom_error_handler::mom_error
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...
Definition: MOM_error_handler.F90:72
mom_coms::reproducing_sum
Find an accurate and order-invariant sum of distributed 2d or 3d fields.
Definition: MOM_coms.F90:54
mom_open_boundary::obc_none
integer, parameter, public obc_none
Indicates the use of no open boundary.
Definition: MOM_open_boundary.F90:58
tidal_bay_initialization
Configures the model for the "tidal_bay" experiment. tidal_bay = Tidally resonant bay from Zygmunt Ko...
Definition: tidal_bay_initialization.F90:3
mom_open_boundary::obc_segment_type
Open boundary segment data structure.
Definition: MOM_open_boundary.F90:103
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_dyn_horgrid::dyn_horgrid_type
Describes the horizontal ocean grid with only dynamic memory arrays.
Definition: MOM_dyn_horgrid.F90:23
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26