MOM_boundary_update.F90

1! This file is part of MOM6, the Modular Ocean Model version 6.
2! See the LICENSE file for licensing information.
3! SPDX-License-Identifier: Apache-2.0
4
5!> Controls where open boundary conditions are applied
7
8use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end, clock_routine
9use mom_diag_mediator, only : time_type
10use mom_error_handler, only : mom_mesg, mom_error, fatal, warning
11use mom_file_parser, only : get_param, log_version, param_file_type, log_param
12use mom_grid, only : ocean_grid_type
19use mom_tracer_registry, only : tracer_registry_type
31
32implicit none ; private
33
34#include <MOM_memory.h>
35
37public update_obc_data
38
39!> The control structure for the MOM_boundary_update module
40type, public :: update_obc_cs ; private
41 logical :: use_files = .false. !< If true, use external files for the open boundary.
42 logical :: use_kelvin = .false. !< If true, use the Kelvin wave open boundary.
43 logical :: use_tidal_bay = .false. !< If true, use the tidal_bay open boundary.
44 logical :: use_shelfwave = .false. !< If true, use the shelfwave open boundary.
45 logical :: use_dyed_channel = .false. !< If true, use the dyed channel open boundary.
46 logical :: debug_obcs = .false. !< If true, write verbose OBC values for debugging purposes.
47 integer :: nk_obc_debug = 0 !< The number of layers of OBC segment data to write out
48 !! in full when DEBUG_OBCS is true.
49 !>@{ Pointers to the control structures for named OBC specifications
50 type(file_obc_cs), pointer :: file_obc_csp => null()
51 type(kelvin_obc_cs), pointer :: kelvin_obc_csp => null()
52 type(tidal_bay_obc_cs) :: tidal_bay_obc
53 type(shelfwave_obc_cs), pointer :: shelfwave_obc_csp => null()
54 type(dyed_channel_obc_cs), pointer :: dyed_channel_obc_csp => null()
55 !>@}
56end type update_obc_cs
57
58integer :: id_clock_pass !< A CPU time clock ID
59
60! character(len=40) :: mdl = "MOM_boundary_update" ! This module's name.
61
62contains
63
64!> The following subroutines and associated definitions provide the
65!! machinery to register and call the subroutines that initialize
66!! open boundary conditions.
67subroutine call_obc_register(G, GV, US, param_file, CS, OBC, tr_Reg)
68 type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
69 type(verticalgrid_type), intent(in) :: gv !< Ocean vertical grid structure
70 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
71 type(param_file_type), intent(in) :: param_file !< Parameter file to parse
72 type(update_obc_cs), pointer :: cs !< Control structure for OBCs
73 type(ocean_obc_type), pointer :: obc !< Open boundary structure
74 type(tracer_registry_type), pointer :: tr_reg !< Tracer registry.
75
76 ! Local variables
77 logical :: debug
78 character(len=200) :: config
79 character(len=40) :: mdl = "MOM_boundary_update" ! This module's name.
80 ! This include declares and sets the variable "version".
81# include "version_variable.h"
82 if (associated(cs)) then
83 call mom_error(warning, "call_OBC_register called with an associated "// &
84 "control structure.")
85 return
86 else ; allocate(cs) ; endif
87
88 call log_version(param_file, mdl, version, "")
89
90 call get_param(param_file, mdl, "USE_FILE_OBC", cs%use_files, &
91 "If true, use external files for the open boundary.", &
92 default=.false.)
93 call get_param(param_file, mdl, "USE_TIDAL_BAY_OBC", cs%use_tidal_bay, &
94 "If true, use the tidal_bay open boundary.", &
95 default=.false.)
96 call get_param(param_file, mdl, "USE_KELVIN_WAVE_OBC", cs%use_Kelvin, &
97 "If true, use the Kelvin wave open boundary.", &
98 default=.false.)
99 call get_param(param_file, mdl, "USE_SHELFWAVE_OBC", cs%use_shelfwave, &
100 "If true, use the shelfwave open boundary.", &
101 default=.false.)
102 call get_param(param_file, mdl, "USE_DYED_CHANNEL_OBC", cs%use_dyed_channel, &
103 "If true, use the dyed channel open boundary.", &
104 default=.false.)
105 call get_param(param_file, mdl, "OBC_USER_CONFIG", config, &
106 "A string that sets how the user code is invoked to set open boundary data: \n"//&
107 " DOME - specified inflow on northern boundary\n"//&
108 " dyed_channel - supercritical with dye on the inflow boundary\n"//&
109 " dyed_obcs - circle_obcs with dyes on the open boundaries\n"//&
110 " Kelvin - barotropic Kelvin wave forcing on the western boundary\n"//&
111 " shelfwave - Flather with shelf wave forcing on western boundary\n"//&
112 " supercritical - now only needed here for the allocations\n"//&
113 " tidal_bay - Flather with tidal forcing on eastern boundary\n"//&
114 " USER - user specified", default="none", do_not_log=.true.)
115 call get_param(param_file, mdl, "DEBUG", debug, &
116 "If true, write out verbose debugging data.", &
117 default=.false., debuggingparam=.true.)
118 call get_param(param_file, mdl, "DEBUG_OBCS", cs%debug_OBCs, &
119 "If true, write out verbose debugging data about OBCs.", &
120 default=.false., debuggingparam=.true.)
121 call get_param(param_file, mdl, "NK_OBC_DEBUG", cs%nk_OBC_debug, &
122 "The number of layers of OBC segment data to write out in full "//&
123 "when DEBUG_OBCS is true.", &
124 default=0, debuggingparam=.true., do_not_log=.not.cs%debug_OBCs)
125
126 if (cs%use_files) cs%use_files = &
127 register_file_obc(param_file, cs%file_OBC_CSp, us, &
128 obc%OBC_Reg)
129
130 if (trim(config) == "DOME") then
131 call register_dome_obc(param_file, us, obc, tr_reg)
132! elseif (trim(config) == "tidal_bay") then
133! elseif (trim(config) == "Kelvin") then
134! elseif (trim(config) == "shelfwave") then
135! elseif (trim(config) == "dyed_channel") then
136 endif
137
138 if (cs%use_tidal_bay) cs%use_tidal_bay = &
139 register_tidal_bay_obc(param_file, cs%tidal_bay_OBC, us, &
140 obc%OBC_Reg)
141 if (cs%use_Kelvin) cs%use_Kelvin = &
142 register_kelvin_obc(param_file, cs%Kelvin_OBC_CSp, us, &
143 obc%OBC_Reg)
144 if (cs%use_shelfwave) cs%use_shelfwave = &
145 register_shelfwave_obc(param_file, cs%shelfwave_OBC_CSp, g, us, &
146 obc%OBC_Reg)
147 if (cs%use_dyed_channel) cs%use_dyed_channel = &
148 register_dyed_channel_obc(param_file, cs%dyed_channel_OBC_CSp, us, &
149 obc%OBC_Reg)
150
151end subroutine call_obc_register
152
153!> Calls appropriate routine to update the open boundary conditions.
154subroutine update_obc_data(OBC, G, GV, US, tv, h, CS, Time)
155 type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
156 type(verticalgrid_type), intent(in) :: gv !< Ocean vertical grid structure
157 type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
158 type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamics structure
159 real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: h !< layer thicknesses [H ~> m or kg m-2]
160 type(ocean_obc_type), pointer :: obc !< Open boundary structure
161 type(update_obc_cs), pointer :: cs !< Control structure for OBCs
162 type(time_type), intent(in) :: time !< Model time
163
164 if (cs%use_tidal_bay) &
165 call tidal_bay_set_obc_data(obc, cs%tidal_bay_OBC, g, gv, us, h, time)
166 if (cs%use_Kelvin) &
167 call kelvin_set_obc_data(obc, cs%Kelvin_OBC_CSp, g, gv, us, h, time)
168 if (cs%use_shelfwave) &
169 call shelfwave_set_obc_data(obc, cs%shelfwave_OBC_CSp, g, gv, us, h, time)
170 if (cs%use_dyed_channel) &
171 call dyed_channel_update_flow(obc, cs%dyed_channel_OBC_CSp, g, gv, us, h, time)
172 if (obc%any_needs_IO_for_data .or. obc%add_tide_constituents) then
173 call read_obc_segment_data(g, gv, us, obc, tv, h, time)
174 call update_obc_segment_data(g, gv, us, obc, h, time)
175 endif
176 if (cs%debug_OBCs) call chksum_obc_segments(obc, g, gv, us, cs%nk_OBC_debug)
177
178end subroutine update_obc_data
179
180!> Clean up the OBC registry.
181subroutine obc_register_end(CS)
182 type(update_obc_cs), pointer :: cs !< Control structure for OBCs
183
184 if (cs%use_files) call file_obc_end(cs%file_OBC_CSp)
185 if (cs%use_Kelvin) call kelvin_obc_end(cs%Kelvin_OBC_CSp)
186
187 if (associated(cs)) deallocate(cs)
188end subroutine obc_register_end
189
190!> \namespace mom_boundary_update
191!! This module updates the open boundary arrays when time-varying.
192!! It caused a circular dependency with the tidal_bay and other setups when in
193!! MOM_open_boundary.
194!!
195!! A small fragment of the grid is shown below:
196!!
197!! j+1 x ^ x ^ x At x: q, CoriolisBu
198!! j+1 > o > o > At ^: v, tauy
199!! j x ^ x ^ x At >: u, taux
200!! j > o > o > At o: h, bathyT, buoy, tr, T, S, Rml, ustar
201!! j-1 x ^ x ^ x
202!! i-1 i i+1 At x & ^:
203!! i i+1 At > & o:
204!!
205!! The boundaries always run through q grid points (x).
206
207end module mom_boundary_update