MOM6
MOM_dyn_horgrid.F90
Go to the documentation of this file.
1 !> Contains a shareable dynamic type for describing horizontal grids and metric data
2 !! and utilty routines that work on this type.
4 
5 ! This file is part of MOM6. See LICENSE.md for the license.
6 
9 use mom_error_handler, only : mom_error, mom_mesg, fatal, warning
11 
12 implicit none ; private
13 
16 
17 ! A note on unit descriptions in comments: MOM6 uses units that can be rescaled for dimensional
18 ! consistency testing. These are noted in comments with units like Z, H, L, and T, along with
19 ! their mks counterparts with notation like "a velocity [Z T-1 ~> m s-1]". If the units
20 ! vary with the Boussinesq approximation, the Boussinesq variant is given first.
21 
22 !> Describes the horizontal ocean grid with only dynamic memory arrays
23 type, public :: dyn_horgrid_type
24  type(mom_domain_type), pointer :: domain => null() !< Ocean model domain
25  type(mom_domain_type), pointer :: domain_aux => null() !< A non-symmetric auxiliary domain type.
26  type(hor_index_type) :: hi !< Horizontal index ranges
27 
28  integer :: isc !< The start i-index of cell centers within the computational domain
29  integer :: iec !< The end i-index of cell centers within the computational domain
30  integer :: jsc !< The start j-index of cell centers within the computational domain
31  integer :: jec !< The end j-index of cell centers within the computational domain
32 
33  integer :: isd !< The start i-index of cell centers within the data domain
34  integer :: ied !< The end i-index of cell centers within the data domain
35  integer :: jsd !< The start j-index of cell centers within the data domain
36  integer :: jed !< The end j-index of cell centers within the data domain
37 
38  integer :: isg !< The start i-index of cell centers within the global domain
39  integer :: ieg !< The end i-index of cell centers within the global domain
40  integer :: jsg !< The start j-index of cell centers within the global domain
41  integer :: jeg !< The end j-index of cell centers within the global domain
42 
43  integer :: iscb !< The start i-index of cell vertices within the computational domain
44  integer :: iecb !< The end i-index of cell vertices within the computational domain
45  integer :: jscb !< The start j-index of cell vertices within the computational domain
46  integer :: jecb !< The end j-index of cell vertices within the computational domain
47 
48  integer :: isdb !< The start i-index of cell vertices within the data domain
49  integer :: iedb !< The end i-index of cell vertices within the data domain
50  integer :: jsdb !< The start j-index of cell vertices within the data domain
51  integer :: jedb !< The end j-index of cell vertices within the data domain
52 
53  integer :: isgb !< The start i-index of cell vertices within the global domain
54  integer :: iegb !< The end i-index of cell vertices within the global domain
55  integer :: jsgb !< The start j-index of cell vertices within the global domain
56  integer :: jegb !< The end j-index of cell vertices within the global domain
57 
58  integer :: isd_global !< The value of isd in the global index space (decompoistion invariant).
59  integer :: jsd_global !< The value of isd in the global index space (decompoistion invariant).
60  integer :: idg_offset !< The offset between the corresponding global and local i-indices.
61  integer :: jdg_offset !< The offset between the corresponding global and local j-indices.
62  logical :: symmetric !< True if symmetric memory is used.
63 
64  logical :: nonblocking_updates !< If true, non-blocking halo updates are
65  !! allowed. The default is .false. (for now).
66  integer :: first_direction !< An integer that indicates which direction is to be updated first in
67  !! directionally split parts of the calculation. This can be altered
68  !! during the course of the run via calls to set_first_direction.
69 
70  real, allocatable, dimension(:,:) :: &
71  mask2dt, & !< 0 for land points and 1 for ocean points on the h-grid [nondim].
72  geolatt, & !< The geographic latitude at q points [degrees of latitude] or [m].
73  geolont, & !< The geographic longitude at q points [degrees of longitude] or [m].
74  dxt, & !< dxT is delta x at h points [L ~> m].
75  idxt, & !< 1/dxT [L-1 ~> m-1].
76  dyt, & !< dyT is delta y at h points [L ~> m].
77  idyt, & !< IdyT is 1/dyT [L-1 ~> m-1].
78  areat, & !< The area of an h-cell [L2 ~> m2].
79  iareat !< 1/areaT [L-2 ~> m-2].
80  real, allocatable, dimension(:,:) :: sin_rot
81  !< The sine of the angular rotation between the local model grid's northward
82  !! and the true northward directions [nondim].
83  real, allocatable, dimension(:,:) :: cos_rot
84  !< The cosine of the angular rotation between the local model grid's northward
85  !! and the true northward directions [nondim].
86 
87  real, allocatable, dimension(:,:) :: &
88  mask2dcu, & !< 0 for boundary points and 1 for ocean points on the u grid [nondim].
89  geolatcu, & !< The geographic latitude at u points [degrees of latitude] or [m].
90  geoloncu, & !< The geographic longitude at u points [degrees of longitude] or [m].
91  dxcu, & !< dxCu is delta x at u points [L ~> m].
92  idxcu, & !< 1/dxCu [L-1 ~> m-1].
93  dycu, & !< dyCu is delta y at u points [L ~> m].
94  idycu, & !< 1/dyCu [L-1 ~> m-1].
95  dy_cu, & !< The unblocked lengths of the u-faces of the h-cell [L ~> m].
96  iareacu, & !< The masked inverse areas of u-grid cells [L-2 ~> m-2].
97  areacu !< The areas of the u-grid cells [L2 ~> m2].
98 
99  real, allocatable, dimension(:,:) :: &
100  mask2dcv, & !< 0 for boundary points and 1 for ocean points on the v grid [nondim].
101  geolatcv, & !< The geographic latitude at v points [degrees of latitude] or [m].
102  geoloncv, & !< The geographic longitude at v points [degrees of longitude] or [m].
103  dxcv, & !< dxCv is delta x at v points [L ~> m].
104  idxcv, & !< 1/dxCv [L-1 ~> m-1].
105  dycv, & !< dyCv is delta y at v points [L ~> m].
106  idycv, & !< 1/dyCv [L-1 ~> m-1].
107  dx_cv, & !< The unblocked lengths of the v-faces of the h-cell [L ~> m].
108  iareacv, & !< The masked inverse areas of v-grid cells [L-2 ~> m-2].
109  areacv !< The areas of the v-grid cells [L2 ~> m2].
110 
111  real, allocatable, dimension(:,:) :: &
112  mask2dbu, & !< 0 for boundary points and 1 for ocean points on the q grid [nondim].
113  geolatbu, & !< The geographic latitude at q points [degrees of latitude] or [m].
114  geolonbu, & !< The geographic longitude at q points [degrees of longitude] or [m].
115  dxbu, & !< dxBu is delta x at q points [L ~> m].
116  idxbu, & !< 1/dxBu [L-1 ~> m-1].
117  dybu, & !< dyBu is delta y at q points [L ~> m].
118  idybu, & !< 1/dyBu [L-1 ~> m-1].
119  areabu, & !< areaBu is the area of a q-cell [L ~> m]
120  iareabu !< IareaBu = 1/areaBu [L-2 ~> m-2].
121 
122  real, pointer, dimension(:) :: gridlatt => null()
123  !< The latitude of T points for the purpose of labeling the output axes.
124  !! On many grids this is the same as geoLatT.
125  real, pointer, dimension(:) :: gridlatb => null()
126  !< The latitude of B points for the purpose of labeling the output axes.
127  !! On many grids this is the same as geoLatBu.
128  real, pointer, dimension(:) :: gridlont => null()
129  !< The longitude of T points for the purpose of labeling the output axes.
130  !! On many grids this is the same as geoLonT.
131  real, pointer, dimension(:) :: gridlonb => null()
132  !< The longitude of B points for the purpose of labeling the output axes.
133  !! On many grids this is the same as geoLonBu.
134  character(len=40) :: &
135  x_axis_units, & !< The units that are used in labeling the x coordinate axes.
136  y_axis_units !< The units that are used in labeling the y coordinate axes.
137  ! Except on a Cartesian grid, these are usually some variant of "degrees".
138 
139  real, allocatable, dimension(:,:) :: &
140  bathyt !< Ocean bottom depth at tracer points, in depth units [Z ~> m].
141 
142  logical :: bathymetry_at_vel !< If true, there are separate values for the
143  !! basin depths at velocity points. Otherwise the effects of
144  !! of topography are entirely determined from thickness points.
145  real, allocatable, dimension(:,:) :: &
146  dblock_u, & !< Topographic depths at u-points at which the flow is blocked [Z ~> m].
147  dopen_u !< Topographic depths at u-points at which the flow is open at width dy_Cu [Z ~> m].
148  real, allocatable, dimension(:,:) :: &
149  dblock_v, & !< Topographic depths at v-points at which the flow is blocked [Z ~> m].
150  dopen_v !< Topographic depths at v-points at which the flow is open at width dx_Cv [Z ~> m].
151  real, allocatable, dimension(:,:) :: &
152  coriolisbu !< The Coriolis parameter at corner points [T-1 ~> s-1].
153  real, allocatable, dimension(:,:) :: &
154  df_dx, & !< Derivative d/dx f (Coriolis parameter) at h-points [T-1 L-1 ~> s-1 m-1].
155  df_dy !< Derivative d/dy f (Coriolis parameter) at h-points [T-1 L-1 ~> s-1 m-1].
156 
157  ! These variables are global sums that are useful for 1-d diagnostics and should not be rescaled.
158  real :: areat_global !< Global sum of h-cell area [m2]
159  real :: iareat_global !< Global sum of inverse h-cell area (1/areaT_global) [m-2]
160 
161  ! These parameters are run-time parameters that are used during some
162  ! initialization routines (but not all)
163  real :: south_lat !< The latitude (or y-coordinate) of the first v-line
164  real :: west_lon !< The longitude (or x-coordinate) of the first u-line
165  real :: len_lat = 0. !< The latitudinal (or y-coord) extent of physical domain
166  real :: len_lon = 0. !< The longitudinal (or x-coord) extent of physical domain
167  real :: rad_earth = 6.378e6 !< The radius of the planet [m].
168  real :: max_depth !< The maximum depth of the ocean [Z ~> m].
169 end type dyn_horgrid_type
170 
171 contains
172 
173 !---------------------------------------------------------------------
174 !> Allocate memory used by the dyn_horgrid_type and related structures.
175 subroutine create_dyn_horgrid(G, HI, bathymetry_at_vel)
176  type(dyn_horgrid_type), pointer :: g !< A pointer to the dynamic horizontal grid type
177  type(hor_index_type), intent(in) :: hi !< A hor_index_type for array extents
178  logical, optional, intent(in) :: bathymetry_at_vel !< If true, there are
179  !! separate values for the basin depths at velocity
180  !! points. Otherwise the effects of topography are
181  !! entirely determined from thickness points.
182  integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb, isg, ieg, jsg, jeg
183 
184  ! This subroutine allocates the lateral elements of the dyn_horgrid_type that
185  ! are always used and zeros them out.
186 
187  if (associated(g)) then
188  call mom_error(warning, "create_dyn_horgrid called with an associated horgrid_type.")
189  else
190  allocate(g)
191  endif
192 
193  g%HI = hi
194 
195  g%isc = hi%isc ; g%iec = hi%iec ; g%jsc = hi%jsc ; g%jec = hi%jec
196  g%isd = hi%isd ; g%ied = hi%ied ; g%jsd = hi%jsd ; g%jed = hi%jed
197  g%isg = hi%isg ; g%ieg = hi%ieg ; g%jsg = hi%jsg ; g%jeg = hi%jeg
198 
199  g%IscB = hi%IscB ; g%IecB = hi%IecB ; g%JscB = hi%JscB ; g%JecB = hi%JecB
200  g%IsdB = hi%IsdB ; g%IedB = hi%IedB ; g%JsdB = hi%JsdB ; g%JedB = hi%JedB
201  g%IsgB = hi%IsgB ; g%IegB = hi%IegB ; g%JsgB = hi%JsgB ; g%JegB = hi%JegB
202 
203  g%idg_offset = hi%idg_offset ; g%jdg_offset = hi%jdg_offset
204  g%isd_global = g%isd + hi%idg_offset ; g%jsd_global = g%jsd + hi%jdg_offset
205  g%symmetric = hi%symmetric
206 
207  g%bathymetry_at_vel = .false.
208  if (present(bathymetry_at_vel)) g%bathymetry_at_vel = bathymetry_at_vel
209 
210  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
211  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
212  isg = g%isg ; ieg = g%ieg ; jsg = g%jsg ; jeg = g%jeg
213 
214  allocate(g%dxT(isd:ied,jsd:jed)) ; g%dxT(:,:) = 0.0
215  allocate(g%dxCu(isdb:iedb,jsd:jed)) ; g%dxCu(:,:) = 0.0
216  allocate(g%dxCv(isd:ied,jsdb:jedb)) ; g%dxCv(:,:) = 0.0
217  allocate(g%dxBu(isdb:iedb,jsdb:jedb)) ; g%dxBu(:,:) = 0.0
218  allocate(g%IdxT(isd:ied,jsd:jed)) ; g%IdxT(:,:) = 0.0
219  allocate(g%IdxCu(isdb:iedb,jsd:jed)) ; g%IdxCu(:,:) = 0.0
220  allocate(g%IdxCv(isd:ied,jsdb:jedb)) ; g%IdxCv(:,:) = 0.0
221  allocate(g%IdxBu(isdb:iedb,jsdb:jedb)) ; g%IdxBu(:,:) = 0.0
222 
223  allocate(g%dyT(isd:ied,jsd:jed)) ; g%dyT(:,:) = 0.0
224  allocate(g%dyCu(isdb:iedb,jsd:jed)) ; g%dyCu(:,:) = 0.0
225  allocate(g%dyCv(isd:ied,jsdb:jedb)) ; g%dyCv(:,:) = 0.0
226  allocate(g%dyBu(isdb:iedb,jsdb:jedb)) ; g%dyBu(:,:) = 0.0
227  allocate(g%IdyT(isd:ied,jsd:jed)) ; g%IdyT(:,:) = 0.0
228  allocate(g%IdyCu(isdb:iedb,jsd:jed)) ; g%IdyCu(:,:) = 0.0
229  allocate(g%IdyCv(isd:ied,jsdb:jedb)) ; g%IdyCv(:,:) = 0.0
230  allocate(g%IdyBu(isdb:iedb,jsdb:jedb)) ; g%IdyBu(:,:) = 0.0
231 
232  allocate(g%areaT(isd:ied,jsd:jed)) ; g%areaT(:,:) = 0.0
233  allocate(g%IareaT(isd:ied,jsd:jed)) ; g%IareaT(:,:) = 0.0
234  allocate(g%areaBu(isdb:iedb,jsdb:jedb)) ; g%areaBu(:,:) = 0.0
235  allocate(g%IareaBu(isdb:iedb,jsdb:jedb)) ; g%IareaBu(:,:) = 0.0
236 
237  allocate(g%mask2dT(isd:ied,jsd:jed)) ; g%mask2dT(:,:) = 0.0
238  allocate(g%mask2dCu(isdb:iedb,jsd:jed)) ; g%mask2dCu(:,:) = 0.0
239  allocate(g%mask2dCv(isd:ied,jsdb:jedb)) ; g%mask2dCv(:,:) = 0.0
240  allocate(g%mask2dBu(isdb:iedb,jsdb:jedb)) ; g%mask2dBu(:,:) = 0.0
241  allocate(g%geoLatT(isd:ied,jsd:jed)) ; g%geoLatT(:,:) = 0.0
242  allocate(g%geoLatCu(isdb:iedb,jsd:jed)) ; g%geoLatCu(:,:) = 0.0
243  allocate(g%geoLatCv(isd:ied,jsdb:jedb)) ; g%geoLatCv(:,:) = 0.0
244  allocate(g%geoLatBu(isdb:iedb,jsdb:jedb)) ; g%geoLatBu(:,:) = 0.0
245  allocate(g%geoLonT(isd:ied,jsd:jed)) ; g%geoLonT(:,:) = 0.0
246  allocate(g%geoLonCu(isdb:iedb,jsd:jed)) ; g%geoLonCu(:,:) = 0.0
247  allocate(g%geoLonCv(isd:ied,jsdb:jedb)) ; g%geoLonCv(:,:) = 0.0
248  allocate(g%geoLonBu(isdb:iedb,jsdb:jedb)) ; g%geoLonBu(:,:) = 0.0
249 
250  allocate(g%dx_Cv(isd:ied,jsdb:jedb)) ; g%dx_Cv(:,:) = 0.0
251  allocate(g%dy_Cu(isdb:iedb,jsd:jed)) ; g%dy_Cu(:,:) = 0.0
252 
253  allocate(g%areaCu(isdb:iedb,jsd:jed)) ; g%areaCu(:,:) = 0.0
254  allocate(g%areaCv(isd:ied,jsdb:jedb)) ; g%areaCv(:,:) = 0.0
255  allocate(g%IareaCu(isdb:iedb,jsd:jed)) ; g%IareaCu(:,:) = 0.0
256  allocate(g%IareaCv(isd:ied,jsdb:jedb)) ; g%IareaCv(:,:) = 0.0
257 
258  allocate(g%bathyT(isd:ied, jsd:jed)) ; g%bathyT(:,:) = 0.0
259  allocate(g%CoriolisBu(isdb:iedb, jsdb:jedb)) ; g%CoriolisBu(:,:) = 0.0
260  allocate(g%dF_dx(isd:ied, jsd:jed)) ; g%dF_dx(:,:) = 0.0
261  allocate(g%dF_dy(isd:ied, jsd:jed)) ; g%dF_dy(:,:) = 0.0
262 
263  allocate(g%sin_rot(isd:ied,jsd:jed)) ; g%sin_rot(:,:) = 0.0
264  allocate(g%cos_rot(isd:ied,jsd:jed)) ; g%cos_rot(:,:) = 1.0
265 
266  if (g%bathymetry_at_vel) then
267  allocate(g%Dblock_u(isdb:iedb, jsd:jed)) ; g%Dblock_u(:,:) = 0.0
268  allocate(g%Dopen_u(isdb:iedb, jsd:jed)) ; g%Dopen_u(:,:) = 0.0
269  allocate(g%Dblock_v(isd:ied, jsdb:jedb)) ; g%Dblock_v(:,:) = 0.0
270  allocate(g%Dopen_v(isd:ied, jsdb:jedb)) ; g%Dopen_v(:,:) = 0.0
271  endif
272 
273  ! gridLonB and gridLatB are used as edge values in some cases, so they
274  ! always need to use symmetric memory allcoations.
275  allocate(g%gridLonT(isg:ieg)) ; g%gridLonT(:) = 0.0
276  allocate(g%gridLonB(isg-1:ieg)) ; g%gridLonB(:) = 0.0
277  allocate(g%gridLatT(jsg:jeg)) ; g%gridLatT(:) = 0.0
278  allocate(g%gridLatB(jsg-1:jeg)) ; g%gridLatB(:) = 0.0
279 
280 end subroutine create_dyn_horgrid
281 
282 !> rescale_dyn_horgrid_bathymetry permits a change in the internal units for the bathymetry on the
283 !! grid, both rescaling the depths and recording the new internal depth units.
284 subroutine rescale_dyn_horgrid_bathymetry(G, m_in_new_units)
285  type(dyn_horgrid_type), intent(inout) :: g !< The dynamic horizontal grid type
286  real, intent(in) :: m_in_new_units !< The new internal representation of 1 m depth.
287 
288  ! Local variables
289  real :: rescale
290  integer :: i, j, isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
291 
292  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
293  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
294 
295  if (m_in_new_units == 1.0) return
296  if (m_in_new_units < 0.0) &
297  call mom_error(fatal, "rescale_grid_bathymetry: Negative depth units are not permitted.")
298  if (m_in_new_units == 0.0) &
299  call mom_error(fatal, "rescale_grid_bathymetry: Zero depth units are not permitted.")
300 
301  rescale = 1.0 / m_in_new_units
302  do j=jsd,jed ; do i=isd,ied
303  g%bathyT(i,j) = rescale*g%bathyT(i,j)
304  enddo ; enddo
305  if (g%bathymetry_at_vel) then ; do j=jsd,jed ; do i=isdb,iedb
306  g%Dblock_u(i,j) = rescale*g%Dblock_u(i,j) ; g%Dopen_u(i,j) = rescale*g%Dopen_u(i,j)
307  enddo ; enddo ; endif
308  if (g%bathymetry_at_vel) then ; do j=jsdb,jedb ; do i=isd,ied
309  g%Dblock_v(i,j) = rescale*g%Dblock_v(i,j) ; g%Dopen_v(i,j) = rescale*g%Dopen_v(i,j)
310  enddo ; enddo ; endif
311  g%max_depth = rescale*g%max_depth
312 
313 end subroutine rescale_dyn_horgrid_bathymetry
314 
315 !> set_derived_dyn_horgrid calculates metric terms that are derived from other metrics.
316 subroutine set_derived_dyn_horgrid(G, US)
317  type(dyn_horgrid_type), intent(inout) :: g !< The dynamic horizontal grid type
318  type(unit_scale_type), optional, intent(in) :: us !< A dimensional unit scaling type
319 ! Various inverse grid spacings and derived areas are calculated within this
320 ! subroutine.
321  real :: m_to_l ! A unit conversion factor [L m-1 ~> nondim]
322  real :: l_to_m ! A unit conversion factor [L m-1 ~> nondim]
323  integer :: i, j, isd, ied, jsd, jed
324  integer :: isdb, iedb, jsdb, jedb
325  m_to_l = 1.0 ; if (present(us)) m_to_l = us%m_to_L
326  l_to_m = 1.0 ; if (present(us)) l_to_m = us%L_to_m
327 
328  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
329  isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
330 
331  do j=jsd,jed ; do i=isd,ied
332  if (g%dxT(i,j) < 0.0) g%dxT(i,j) = 0.0
333  if (g%dyT(i,j) < 0.0) g%dyT(i,j) = 0.0
334  g%IdxT(i,j) = adcroft_reciprocal(g%dxT(i,j))
335  g%IdyT(i,j) = adcroft_reciprocal(g%dyT(i,j))
336  g%IareaT(i,j) = adcroft_reciprocal(g%areaT(i,j))
337  enddo ; enddo
338 
339  do j=jsd,jed ; do i=isdb,iedb
340  if (g%dxCu(i,j) < 0.0) g%dxCu(i,j) = 0.0
341  if (g%dyCu(i,j) < 0.0) g%dyCu(i,j) = 0.0
342  g%IdxCu(i,j) = adcroft_reciprocal(g%dxCu(i,j))
343  g%IdyCu(i,j) = adcroft_reciprocal(g%dyCu(i,j))
344  enddo ; enddo
345 
346  do j=jsdb,jedb ; do i=isd,ied
347  if (g%dxCv(i,j) < 0.0) g%dxCv(i,j) = 0.0
348  if (g%dyCv(i,j) < 0.0) g%dyCv(i,j) = 0.0
349  g%IdxCv(i,j) = adcroft_reciprocal(g%dxCv(i,j))
350  g%IdyCv(i,j) = adcroft_reciprocal(g%dyCv(i,j))
351  enddo ; enddo
352 
353  do j=jsdb,jedb ; do i=isdb,iedb
354  if (g%dxBu(i,j) < 0.0) g%dxBu(i,j) = 0.0
355  if (g%dyBu(i,j) < 0.0) g%dyBu(i,j) = 0.0
356 
357  g%IdxBu(i,j) = adcroft_reciprocal(g%dxBu(i,j))
358  g%IdyBu(i,j) = adcroft_reciprocal(g%dyBu(i,j))
359  ! areaBu has usually been set to a positive area elsewhere.
360  if (g%areaBu(i,j) <= 0.0) g%areaBu(i,j) = g%dxBu(i,j) * g%dyBu(i,j)
361  g%IareaBu(i,j) = adcroft_reciprocal(g%areaBu(i,j))
362  enddo ; enddo
363 
364 end subroutine set_derived_dyn_horgrid
365 
366 !> Adcroft_reciprocal(x) = 1/x for |x|>0 or 0 for x=0.
367 function adcroft_reciprocal(val) result(I_val)
368  real, intent(in) :: val !< The value being inverted.
369  real :: i_val !< The Adcroft reciprocal of val.
370 
371  i_val = 0.0 ; if (val /= 0.0) i_val = 1.0/val
372 end function adcroft_reciprocal
373 
374 !---------------------------------------------------------------------
375 !> Release memory used by the dyn_horgrid_type and related structures.
376 subroutine destroy_dyn_horgrid(G)
377  type(dyn_horgrid_type), pointer :: g !< The dynamic horizontal grid type
378 
379  if (.not.associated(g)) then
380  call mom_error(fatal, "destroy_dyn_horgrid called with an unassociated horgrid_type.")
381  endif
382 
383  deallocate(g%dxT) ; deallocate(g%dxCu) ; deallocate(g%dxCv) ; deallocate(g%dxBu)
384  deallocate(g%IdxT) ; deallocate(g%IdxCu) ; deallocate(g%IdxCv) ; deallocate(g%IdxBu)
385 
386  deallocate(g%dyT) ; deallocate(g%dyCu) ; deallocate(g%dyCv) ; deallocate(g%dyBu)
387  deallocate(g%IdyT) ; deallocate(g%IdyCu) ; deallocate(g%IdyCv) ; deallocate(g%IdyBu)
388 
389  deallocate(g%areaT) ; deallocate(g%IareaT)
390  deallocate(g%areaBu) ; deallocate(g%IareaBu)
391  deallocate(g%areaCu) ; deallocate(g%IareaCu)
392  deallocate(g%areaCv) ; deallocate(g%IareaCv)
393 
394  deallocate(g%mask2dT) ; deallocate(g%mask2dCu)
395  deallocate(g%mask2dCv) ; deallocate(g%mask2dBu)
396 
397  deallocate(g%geoLatT) ; deallocate(g%geoLatCu)
398  deallocate(g%geoLatCv) ; deallocate(g%geoLatBu)
399  deallocate(g%geoLonT) ; deallocate(g%geoLonCu)
400  deallocate(g%geoLonCv) ; deallocate(g%geoLonBu)
401 
402  deallocate(g%dx_Cv) ; deallocate(g%dy_Cu)
403 
404  deallocate(g%bathyT) ; deallocate(g%CoriolisBu)
405  deallocate(g%dF_dx) ; deallocate(g%dF_dy)
406  deallocate(g%sin_rot) ; deallocate(g%cos_rot)
407 
408  if (allocated(g%Dblock_u)) deallocate(g%Dblock_u)
409  if (allocated(g%Dopen_u)) deallocate(g%Dopen_u)
410  if (allocated(g%Dblock_v)) deallocate(g%Dblock_v)
411  if (allocated(g%Dopen_v)) deallocate(g%Dopen_v)
412 
413  deallocate(g%gridLonT) ; deallocate(g%gridLatT)
414  deallocate(g%gridLonB) ; deallocate(g%gridLatB)
415 
416  deallocate(g%Domain%mpp_domain)
417  deallocate(g%Domain)
418 
419  deallocate(g)
420 
421 end subroutine destroy_dyn_horgrid
422 
423 end module mom_dyn_horgrid
mom_dyn_horgrid::adcroft_reciprocal
real function adcroft_reciprocal(val)
Adcroft_reciprocal(x) = 1/x for |x|>0 or 0 for x=0.
Definition: MOM_dyn_horgrid.F90:368
mom_dyn_horgrid::destroy_dyn_horgrid
subroutine, public destroy_dyn_horgrid(G)
Release memory used by the dyn_horgrid_type and related structures.
Definition: MOM_dyn_horgrid.F90:377
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_hor_index
Defines the horizontal index type (hor_index_type) used for providing index ranges.
Definition: MOM_hor_index.F90:2
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_dyn_horgrid::set_derived_dyn_horgrid
subroutine, public set_derived_dyn_horgrid(G, US)
set_derived_dyn_horgrid calculates metric terms that are derived from other metrics.
Definition: MOM_dyn_horgrid.F90:317
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_dyn_horgrid::create_dyn_horgrid
subroutine, public create_dyn_horgrid(G, HI, bathymetry_at_vel)
Allocate memory used by the dyn_horgrid_type and related structures.
Definition: MOM_dyn_horgrid.F90:176
mom_hor_index::hor_index_type
Container for horizontal index ranges for data, computational and global domains.
Definition: MOM_hor_index.F90:15
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
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_dyn_horgrid::rescale_dyn_horgrid_bathymetry
subroutine, public rescale_dyn_horgrid_bathymetry(G, m_in_new_units)
rescale_dyn_horgrid_bathymetry permits a change in the internal units for the bathymetry on the grid,...
Definition: MOM_dyn_horgrid.F90:285
mom_domains::mom_domain_type
The MOM_domain_type contains information about the domain decompositoin.
Definition: MOM_domains.F90:99
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