8 use mom_domains,
only : get_global_shape, get_domain_extent_dsamp2
12 implicit none ;
private
14 #include <MOM_memory.h>
16 public mom_grid_init, mom_grid_end, set_derived_metrics, set_first_direction
17 public ispointincell,
hor_index_type, get_global_grid_size, rescale_grid_bathymetry
67 logical :: nonblocking_updates
69 integer :: first_direction
75 real allocable_,
dimension(NIMEM_,NJMEM_) :: &
76 mask2dt, & !< 0 for land points and 1 for ocean points on the h-grid. Nd.
77 geolatt, & !< The geographic latitude at q points in degrees of latitude or m.
78 geolont, & !< The geographic longitude at q points in degrees of longitude or m.
79 dxt, & !< dxT is delta x at h points [m].
80 idxt, & !< 1/dxT [m-1].
81 dyt, & !< dyT is delta y at h points [m].
82 idyt, & !< IdyT is 1/dyT [m-1].
83 areat, & !< The area of an h-cell [m2].
84 iareat, & !< 1/areaT [m-2].
85 sin_rot, & !< The sine of the angular rotation between the local model grid
's northward
86 !! and the true northward directions.
87 cos_rot !< The cosine of the angular rotation between the local model grid's northward
90 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_) :: &
91 mask2dcu, & !< 0 for boundary points and 1 for ocean points on the u grid. Nondim.
92 geolatcu, & !< The geographic latitude at u points in degrees of latitude or m.
93 geoloncu, & !< The geographic longitude at u points in degrees of longitude or m.
94 dxcu, & !< dxCu is delta x at u points [m].
95 idxcu, & !< 1/dxCu [m-1].
96 dycu, & !< dyCu is delta y at u points [m].
97 idycu, & !< 1/dyCu [m-1].
98 dy_cu, & !< The unblocked lengths of the u-faces of the h-cell [m].
99 iareacu, & !< The masked inverse areas of u-grid cells [m2].
102 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_) :: &
103 mask2dcv, & !< 0 for boundary points and 1 for ocean points on the v grid. Nondim.
104 geolatcv, & !< The geographic latitude at v points in degrees of latitude or m.
105 geoloncv, & !< The geographic longitude at v points in degrees of longitude or m.
106 dxcv, & !< dxCv is delta x at v points [m].
107 idxcv, & !< 1/dxCv [m-1].
108 dycv, & !< dyCv is delta y at v points [m].
109 idycv, & !< 1/dyCv [m-1].
110 dx_cv, & !< The unblocked lengths of the v-faces of the h-cell [m].
111 iareacv, & !< The masked inverse areas of v-grid cells [m2].
114 real allocable_,
dimension(NIMEMB_PTR_,NJMEMB_PTR_) :: &
115 mask2dbu, & !< 0 for boundary points and 1 for ocean points on the q grid. Nondim.
116 geolatbu, & !< The geographic latitude at q points in degrees of latitude or m.
117 geolonbu, & !< The geographic longitude at q points in degrees of longitude or m.
118 dxbu, & !< dxBu is delta x at q points [m].
119 idxbu, & !< 1/dxBu [m-1].
120 dybu, & !< dyBu is delta y at q points [m].
121 idybu, & !< 1/dyBu [m-1].
122 areabu, & !< areaBu is the area of a q-cell [m2]
125 real,
pointer,
dimension(:) :: &
126 gridlatt => null(), &
130 real,
pointer,
dimension(:) :: &
131 gridlont => null(), &
135 character(len=40) :: &
136 x_axis_units, & !< The units that are used in labeling the x coordinate axes.
139 real allocable_,
dimension(NIMEM_,NJMEM_) :: &
142 logical :: bathymetry_at_vel
145 real allocable_,
dimension(NIMEMB_PTR_,NJMEM_) :: &
146 dblock_u, & !< Topographic depths at u-points at which the flow is blocked [Z ~> m].
148 real allocable_,
dimension(NIMEM_,NJMEMB_PTR_) :: &
149 dblock_v, & !< Topographic depths at v-points at which the flow is blocked [Z ~> m].
151 real allocable_,
dimension(NIMEMB_PTR_,NJMEMB_PTR_) :: &
153 real allocable_,
dimension(NIMEM_,NJMEM_) :: &
154 df_dx, & !< Derivative d/dx f (Coriolis parameter) at h-points [T-1 m-1 ~> s-1 m-1].
160 real :: iareat_global
172 real :: rad_earth = 6.378e6
179 subroutine mom_grid_init(G, param_file, HI, global_indexing, bathymetry_at_vel)
183 optional,
intent(in) :: hi
184 logical,
optional,
intent(in) :: global_indexing
187 logical,
optional,
intent(in) :: bathymetry_at_vel
193 #include "version_variable.h"
194 integer :: isd, ied, jsd, jed, nk
195 integer :: isdb, iedb, jsdb, jedb
196 integer :: ied_max, jed_max
197 integer :: niblock, njblock, nihalo, njhalo, nblocks, n, i, j
198 logical :: local_indexing
201 integer,
allocatable,
dimension(:) :: ibegin, iend, jbegin, jend
202 character(len=40) :: mod_nm =
"MOM_grid"
207 "Parameters providing information about the lateral grid.")
210 call get_param(param_file, mod_nm,
"NIBLOCK", niblock,
"The number of blocks "// &
211 "in the x-direction on each processor (for openmp).", default=1, &
213 call get_param(param_file, mod_nm,
"NJBLOCK", njblock,
"The number of blocks "// &
214 "in the y-direction on each processor (for openmp).", default=1, &
217 if (
present(hi))
then
220 g%isc = hi%isc ; g%iec = hi%iec ; g%jsc = hi%jsc ; g%jec = hi%jec
221 g%isd = hi%isd ; g%ied = hi%ied ; g%jsd = hi%jsd ; g%jed = hi%jed
222 g%isg = hi%isg ; g%ieg = hi%ieg ; g%jsg = hi%jsg ; g%jeg = hi%jeg
224 g%IscB = hi%IscB ; g%IecB = hi%IecB ; g%JscB = hi%JscB ; g%JecB = hi%JecB
225 g%IsdB = hi%IsdB ; g%IedB = hi%IedB ; g%JsdB = hi%JsdB ; g%JedB = hi%JedB
226 g%IsgB = hi%IsgB ; g%IegB = hi%IegB ; g%JsgB = hi%JsgB ; g%JegB = hi%JegB
228 g%idg_offset = hi%idg_offset ; g%jdg_offset = hi%jdg_offset
229 g%isd_global = g%isd + hi%idg_offset ; g%jsd_global = g%jsd + hi%jdg_offset
230 g%symmetric = hi%symmetric
232 local_indexing = .true.
233 if (
present(global_indexing)) local_indexing = .not.global_indexing
234 call hor_index_init(g%Domain, g%HI, param_file, &
235 local_indexing=local_indexing)
239 call get_domain_extent(g%Domain, g%isc, g%iec, g%jsc, g%jec, &
240 g%isd, g%ied, g%jsd, g%jed, &
241 g%isg, g%ieg, g%jsg, g%jeg, &
242 g%idg_offset, g%jdg_offset, g%symmetric, &
243 local_indexing=local_indexing)
244 g%isd_global = g%isd+g%idg_offset ; g%jsd_global = g%jsd+g%jdg_offset
247 g%nonblocking_updates = g%Domain%nonblocking_updates
250 g%IscB = g%isc ; g%JscB = g%jsc
251 g%IsdB = g%isd ; g%JsdB = g%jsd
252 g%IsgB = g%isg ; g%JsgB = g%jsg
253 if (g%symmetric)
then
254 g%IscB = g%isc-1 ; g%JscB = g%jsc-1
255 g%IsdB = g%isd-1 ; g%JsdB = g%jsd-1
256 g%IsgB = g%isg-1 ; g%JsgB = g%jsg-1
258 g%IecB = g%iec ; g%JecB = g%jec
259 g%IedB = g%ied ; g%JedB = g%jed
260 g%IegB = g%ieg ; g%JegB = g%jeg
262 call mom_mesg(
" MOM_grid.F90, MOM_grid_init: allocating metrics", 5)
264 call allocate_metrics(g)
266 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
267 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
269 g%bathymetry_at_vel = .false.
270 if (
present(bathymetry_at_vel)) g%bathymetry_at_vel = bathymetry_at_vel
271 if (g%bathymetry_at_vel)
then
272 alloc_(g%Dblock_u(isdb:iedb, jsd:jed)) ; g%Dblock_u(:,:) = 0.0
273 alloc_(g%Dopen_u(isdb:iedb, jsd:jed)) ; g%Dopen_u(:,:) = 0.0
274 alloc_(g%Dblock_v(isd:ied, jsdb:jedb)) ; g%Dblock_v(:,:) = 0.0
275 alloc_(g%Dopen_v(isd:ied, jsdb:jedb)) ; g%Dopen_v(:,:) = 0.0
279 nihalo = g%Domain%nihalo
280 njhalo = g%Domain%njhalo
281 nblocks = niblock * njblock
282 if (nblocks < 1)
call mom_error(fatal,
"MOM_grid_init: " // &
283 "nblocks(=NI_BLOCK*NJ_BLOCK) must be no less than 1")
285 allocate(ibegin(niblock), iend(niblock), jbegin(njblock), jend(njblock))
286 call compute_block_extent(g%HI%isc,g%HI%iec,niblock,ibegin,iend)
287 call compute_block_extent(g%HI%jsc,g%HI%jec,njblock,jbegin,jend)
290 if (iend(i)-ibegin(i) > iend(niblock)-ibegin(niblock) )
call mom_error(fatal, &
291 "MOM_grid_init: the last block size in x-direction is not the largest")
294 if (jend(j)-jbegin(j) > jend(njblock)-jbegin(njblock) )
call mom_error(fatal, &
295 "MOM_grid_init: the last block size in y-direction is not the largest")
299 allocate(g%Block(nblocks))
300 ied_max = 1 ; jed_max = 1
305 i = mod((n-1), niblock) + 1
306 j = (n-1)/niblock + 1
308 g%Block(n)%isd = 1 ; g%Block(n)%jsd = 1
309 g%Block(n)%isc = g%Block(n)%isd+nihalo
310 g%Block(n)%jsc = g%Block(n)%jsd+njhalo
311 g%Block(n)%iec = g%Block(n)%isc + iend(i) - ibegin(i)
312 g%Block(n)%jec = g%Block(n)%jsc + jend(j) - jbegin(j)
313 g%Block(n)%ied = g%Block(n)%iec + nihalo
314 g%Block(n)%jed = g%Block(n)%jec + njhalo
315 g%Block(n)%IscB = g%Block(n)%isc; g%Block(n)%IecB = g%Block(n)%iec
316 g%Block(n)%JscB = g%Block(n)%jsc; g%Block(n)%JecB = g%Block(n)%jec
319 if (g%symmetric)
then
320 if (i==1) g%Block(n)%IscB = g%Block(n)%IscB-1
321 if (j==1) g%Block(n)%JscB = g%Block(n)%JscB-1
323 g%Block(n)%IsdB = g%Block(n)%isd; g%Block(n)%IedB = g%Block(n)%ied
324 g%Block(n)%JsdB = g%Block(n)%jsd; g%Block(n)%JedB = g%Block(n)%jed
327 if (g%symmetric)
then
328 g%Block(n)%IsdB = g%Block(n)%IsdB-1
329 g%Block(n)%JsdB = g%Block(n)%JsdB-1
331 g%Block(n)%idg_offset = (ibegin(i) - g%Block(n)%isc) + g%HI%idg_offset
332 g%Block(n)%jdg_offset = (jbegin(j) - g%Block(n)%jsc) + g%HI%jdg_offset
335 ied_max = max(ied_max, g%Block(n)%ied)
336 jed_max = max(jed_max, g%Block(n)%jed)
342 g%Block(n)%ied = ied_max ; g%Block(n)%IedB = ied_max
343 g%Block(n)%jed = jed_max ; g%Block(n)%JedB = jed_max
347 if ( g%block(nblocks)%ied+g%block(nblocks)%idg_offset > g%HI%ied + g%HI%idg_offset ) &
348 call mom_error(fatal,
"MOM_grid_init: G%ied_bk > G%ied")
349 if ( g%block(nblocks)%jed+g%block(nblocks)%jdg_offset > g%HI%jed + g%HI%jdg_offset ) &
350 call mom_error(fatal,
"MOM_grid_init: G%jed_bk > G%jed")
352 call get_domain_extent_dsamp2(g%Domain, g%HId2%isc, g%HId2%iec, g%HId2%jsc, g%HId2%jec,&
353 g%HId2%isd, g%HId2%ied, g%HId2%jsd, g%HId2%jed,&
354 g%HId2%isg, g%HId2%ieg, g%HId2%jsg, g%HId2%jeg)
357 g%HId2%IscB = g%HId2%isc ; g%HId2%JscB = g%HId2%jsc
358 g%HId2%IsdB = g%HId2%isd ; g%HId2%JsdB = g%HId2%jsd
359 g%HId2%IsgB = g%HId2%isg ; g%HId2%JsgB = g%HId2%jsg
360 if (g%symmetric)
then
361 g%HId2%IscB = g%HId2%isc-1 ; g%HId2%JscB = g%HId2%jsc-1
362 g%HId2%IsdB = g%HId2%isd-1 ; g%HId2%JsdB = g%HId2%jsd-1
363 g%HId2%IsgB = g%HId2%isg-1 ; g%HId2%JsgB = g%HId2%jsg-1
365 g%HId2%IecB = g%HId2%iec ; g%HId2%JecB = g%HId2%jec
366 g%HId2%IedB = g%HId2%ied ; g%HId2%JedB = g%HId2%jed
367 g%HId2%IegB = g%HId2%ieg ; g%HId2%JegB = g%HId2%jeg
369 end subroutine mom_grid_init
373 subroutine rescale_grid_bathymetry(G, m_in_new_units)
375 real,
intent(in) :: m_in_new_units
379 integer :: i, j, isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
381 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
382 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
384 if (m_in_new_units == 1.0)
return
385 if (m_in_new_units < 0.0) &
386 call mom_error(fatal,
"rescale_grid_bathymetry: Negative depth units are not permitted.")
387 if (m_in_new_units == 0.0) &
388 call mom_error(fatal,
"rescale_grid_bathymetry: Zero depth units are not permitted.")
390 rescale = 1.0 / m_in_new_units
391 do j=jsd,jed ;
do i=isd,ied
392 g%bathyT(i,j) = rescale*g%bathyT(i,j)
394 if (g%bathymetry_at_vel)
then ;
do j=jsd,jed ;
do i=isdb,iedb
395 g%Dblock_u(i,j) = rescale*g%Dblock_u(i,j) ; g%Dopen_u(i,j) = rescale*g%Dopen_u(i,j)
396 enddo ;
enddo ;
endif
397 if (g%bathymetry_at_vel)
then ;
do j=jsdb,jedb ;
do i=isd,ied
398 g%Dblock_v(i,j) = rescale*g%Dblock_v(i,j) ; g%Dopen_v(i,j) = rescale*g%Dopen_v(i,j)
399 enddo ;
enddo ;
endif
400 g%max_depth = rescale*g%max_depth
402 end subroutine rescale_grid_bathymetry
405 subroutine set_derived_metrics(G)
409 integer :: i, j, isd, ied, jsd, jed
410 integer :: isdb, iedb, jsdb, jedb
412 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
413 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
415 do j=jsd,jed ;
do i=isd,ied
416 if (g%dxT(i,j) < 0.0) g%dxT(i,j) = 0.0
417 if (g%dyT(i,j) < 0.0) g%dyT(i,j) = 0.0
418 g%IdxT(i,j) = adcroft_reciprocal(g%dxT(i,j))
419 g%IdyT(i,j) = adcroft_reciprocal(g%dyT(i,j))
420 g%IareaT(i,j) = adcroft_reciprocal(g%areaT(i,j))
423 do j=jsd,jed ;
do i=isdb,iedb
424 if (g%dxCu(i,j) < 0.0) g%dxCu(i,j) = 0.0
425 if (g%dyCu(i,j) < 0.0) g%dyCu(i,j) = 0.0
426 g%IdxCu(i,j) = adcroft_reciprocal(g%dxCu(i,j))
427 g%IdyCu(i,j) = adcroft_reciprocal(g%dyCu(i,j))
430 do j=jsdb,jedb ;
do i=isd,ied
431 if (g%dxCv(i,j) < 0.0) g%dxCv(i,j) = 0.0
432 if (g%dyCv(i,j) < 0.0) g%dyCv(i,j) = 0.0
433 g%IdxCv(i,j) = adcroft_reciprocal(g%dxCv(i,j))
434 g%IdyCv(i,j) = adcroft_reciprocal(g%dyCv(i,j))
437 do j=jsdb,jedb ;
do i=isdb,iedb
438 if (g%dxBu(i,j) < 0.0) g%dxBu(i,j) = 0.0
439 if (g%dyBu(i,j) < 0.0) g%dyBu(i,j) = 0.0
441 g%IdxBu(i,j) = adcroft_reciprocal(g%dxBu(i,j))
442 g%IdyBu(i,j) = adcroft_reciprocal(g%dyBu(i,j))
444 if (g%areaBu(i,j) <= 0.0) g%areaBu(i,j) = g%dxBu(i,j) * g%dyBu(i,j)
445 g%IareaBu(i,j) = adcroft_reciprocal(g%areaBu(i,j))
447 end subroutine set_derived_metrics
450 function adcroft_reciprocal(val)
result(I_val)
451 real,
intent(in) :: val
454 i_val = 0.0 ;
if (val /= 0.0) i_val = 1.0/val
455 end function adcroft_reciprocal
458 logical function ispointincell(G, i, j, x, y)
460 integer,
intent(in) :: i
461 integer,
intent(in) :: j
462 real,
intent(in) :: x
463 real,
intent(in) :: y
465 real :: xne, xnw, xse, xsw, yne, ynw, yse, ysw
466 real :: p0, p1, p2, p3, l0, l1, l2, l3
467 ispointincell = .false.
468 xne = g%geoLonBu(i ,j ) ; yne = g%geoLatBu(i ,j )
469 xnw = g%geoLonBu(i-1,j ) ; ynw = g%geoLatBu(i-1,j )
470 xse = g%geoLonBu(i ,j-1) ; yse = g%geoLatBu(i ,j-1)
471 xsw = g%geoLonBu(i-1,j-1) ; ysw = g%geoLatBu(i-1,j-1)
473 if (x<min(xne,xnw,xse,xsw) .or. x>max(xne,xnw,xse,xsw) .or. &
474 y<min(yne,ynw,yse,ysw) .or. y>max(yne,ynw,yse,ysw) )
then
477 l0 = (x-xsw)*(yse-ysw) - (y-ysw)*(xse-xsw)
478 l1 = (x-xse)*(yne-yse) - (y-yse)*(xne-xse)
479 l2 = (x-xne)*(ynw-yne) - (y-yne)*(xnw-xne)
480 l3 = (x-xnw)*(ysw-ynw) - (y-ynw)*(xsw-xnw)
482 p0 = sign(1., l0) ;
if (l0 == 0.) p0=0.
483 p1 = sign(1., l1) ;
if (l1 == 0.) p1=0.
484 p2 = sign(1., l2) ;
if (l2 == 0.) p2=0.
485 p3 = sign(1., l3) ;
if (l3 == 0.) p3=0.
487 if ( (abs(p0)+abs(p2)) + (abs(p1)+abs(p3)) == abs((p0+p2) + (p1+p3)) )
then
490 end function ispointincell
493 subroutine set_first_direction(G, y_first)
495 integer,
intent(in) :: y_first
497 g%first_direction = y_first
498 end subroutine set_first_direction
501 subroutine get_global_grid_size(G, niglobal, njglobal)
503 integer,
intent(out) :: niglobal
504 integer,
intent(out) :: njglobal
506 call get_global_shape(g%domain, niglobal, njglobal)
508 end subroutine get_global_grid_size
511 subroutine allocate_metrics(G)
513 integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, isg, ieg, jsg, jeg
518 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
519 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
520 isg = g%isg ; ieg = g%ieg ; jsg = g%jsg ; jeg = g%jeg
522 alloc_(g%dxT(isd:ied,jsd:jed)) ; g%dxT(:,:) = 0.0
523 alloc_(g%dxCu(isdb:iedb,jsd:jed)) ; g%dxCu(:,:) = 0.0
524 alloc_(g%dxCv(isd:ied,jsdb:jedb)) ; g%dxCv(:,:) = 0.0
525 alloc_(g%dxBu(isdb:iedb,jsdb:jedb)) ; g%dxBu(:,:) = 0.0
526 alloc_(g%IdxT(isd:ied,jsd:jed)) ; g%IdxT(:,:) = 0.0
527 alloc_(g%IdxCu(isdb:iedb,jsd:jed)) ; g%IdxCu(:,:) = 0.0
528 alloc_(g%IdxCv(isd:ied,jsdb:jedb)) ; g%IdxCv(:,:) = 0.0
529 alloc_(g%IdxBu(isdb:iedb,jsdb:jedb)) ; g%IdxBu(:,:) = 0.0
531 alloc_(g%dyT(isd:ied,jsd:jed)) ; g%dyT(:,:) = 0.0
532 alloc_(g%dyCu(isdb:iedb,jsd:jed)) ; g%dyCu(:,:) = 0.0
533 alloc_(g%dyCv(isd:ied,jsdb:jedb)) ; g%dyCv(:,:) = 0.0
534 alloc_(g%dyBu(isdb:iedb,jsdb:jedb)) ; g%dyBu(:,:) = 0.0
535 alloc_(g%IdyT(isd:ied,jsd:jed)) ; g%IdyT(:,:) = 0.0
536 alloc_(g%IdyCu(isdb:iedb,jsd:jed)) ; g%IdyCu(:,:) = 0.0
537 alloc_(g%IdyCv(isd:ied,jsdb:jedb)) ; g%IdyCv(:,:) = 0.0
538 alloc_(g%IdyBu(isdb:iedb,jsdb:jedb)) ; g%IdyBu(:,:) = 0.0
540 alloc_(g%areaT(isd:ied,jsd:jed)) ; g%areaT(:,:) = 0.0
541 alloc_(g%IareaT(isd:ied,jsd:jed)) ; g%IareaT(:,:) = 0.0
542 alloc_(g%areaBu(isdb:iedb,jsdb:jedb)) ; g%areaBu(:,:) = 0.0
543 alloc_(g%IareaBu(isdb:iedb,jsdb:jedb)) ; g%IareaBu(:,:) = 0.0
545 alloc_(g%mask2dT(isd:ied,jsd:jed)) ; g%mask2dT(:,:) = 0.0
546 alloc_(g%mask2dCu(isdb:iedb,jsd:jed)) ; g%mask2dCu(:,:) = 0.0
547 alloc_(g%mask2dCv(isd:ied,jsdb:jedb)) ; g%mask2dCv(:,:) = 0.0
548 alloc_(g%mask2dBu(isdb:iedb,jsdb:jedb)) ; g%mask2dBu(:,:) = 0.0
549 alloc_(g%geoLatT(isd:ied,jsd:jed)) ; g%geoLatT(:,:) = 0.0
550 alloc_(g%geoLatCu(isdb:iedb,jsd:jed)) ; g%geoLatCu(:,:) = 0.0
551 alloc_(g%geoLatCv(isd:ied,jsdb:jedb)) ; g%geoLatCv(:,:) = 0.0
552 alloc_(g%geoLatBu(isdb:iedb,jsdb:jedb)) ; g%geoLatBu(:,:) = 0.0
553 alloc_(g%geoLonT(isd:ied,jsd:jed)) ; g%geoLonT(:,:) = 0.0
554 alloc_(g%geoLonCu(isdb:iedb,jsd:jed)) ; g%geoLonCu(:,:) = 0.0
555 alloc_(g%geoLonCv(isd:ied,jsdb:jedb)) ; g%geoLonCv(:,:) = 0.0
556 alloc_(g%geoLonBu(isdb:iedb,jsdb:jedb)) ; g%geoLonBu(:,:) = 0.0
558 alloc_(g%dx_Cv(isd:ied,jsdb:jedb)) ; g%dx_Cv(:,:) = 0.0
559 alloc_(g%dy_Cu(isdb:iedb,jsd:jed)) ; g%dy_Cu(:,:) = 0.0
561 alloc_(g%areaCu(isdb:iedb,jsd:jed)) ; g%areaCu(:,:) = 0.0
562 alloc_(g%areaCv(isd:ied,jsdb:jedb)) ; g%areaCv(:,:) = 0.0
563 alloc_(g%IareaCu(isdb:iedb,jsd:jed)) ; g%IareaCu(:,:) = 0.0
564 alloc_(g%IareaCv(isd:ied,jsdb:jedb)) ; g%IareaCv(:,:) = 0.0
566 alloc_(g%bathyT(isd:ied, jsd:jed)) ; g%bathyT(:,:) = 0.0
567 alloc_(g%CoriolisBu(isdb:iedb, jsdb:jedb)) ; g%CoriolisBu(:,:) = 0.0
568 alloc_(g%dF_dx(isd:ied, jsd:jed)) ; g%dF_dx(:,:) = 0.0
569 alloc_(g%dF_dy(isd:ied, jsd:jed)) ; g%dF_dy(:,:) = 0.0
571 alloc_(g%sin_rot(isd:ied,jsd:jed)) ; g%sin_rot(:,:) = 0.0
572 alloc_(g%cos_rot(isd:ied,jsd:jed)) ; g%cos_rot(:,:) = 1.0
574 allocate(g%gridLonT(isg:ieg)) ; g%gridLonT(:) = 0.0
575 allocate(g%gridLonB(g%IsgB:g%IegB)) ; g%gridLonB(:) = 0.0
576 allocate(g%gridLatT(jsg:jeg)) ; g%gridLatT(:) = 0.0
577 allocate(g%gridLatB(g%JsgB:g%JegB)) ; g%gridLatB(:) = 0.0
579 end subroutine allocate_metrics
582 subroutine mom_grid_end(G)
585 dealloc_(g%dxT) ; dealloc_(g%dxCu) ; dealloc_(g%dxCv) ; dealloc_(g%dxBu)
586 dealloc_(g%IdxT) ; dealloc_(g%IdxCu) ; dealloc_(g%IdxCv) ; dealloc_(g%IdxBu)
588 dealloc_(g%dyT) ; dealloc_(g%dyCu) ; dealloc_(g%dyCv) ; dealloc_(g%dyBu)
589 dealloc_(g%IdyT) ; dealloc_(g%IdyCu) ; dealloc_(g%IdyCv) ; dealloc_(g%IdyBu)
591 dealloc_(g%areaT) ; dealloc_(g%IareaT)
592 dealloc_(g%areaBu) ; dealloc_(g%IareaBu)
593 dealloc_(g%areaCu) ; dealloc_(g%IareaCu)
594 dealloc_(g%areaCv) ; dealloc_(g%IareaCv)
596 dealloc_(g%mask2dT) ; dealloc_(g%mask2dCu)
597 dealloc_(g%mask2dCv) ; dealloc_(g%mask2dBu)
599 dealloc_(g%geoLatT) ; dealloc_(g%geoLatCu)
600 dealloc_(g%geoLatCv) ; dealloc_(g%geoLatBu)
601 dealloc_(g%geoLonT) ; dealloc_(g%geoLonCu)
602 dealloc_(g%geoLonCv) ; dealloc_(g%geoLonBu)
604 dealloc_(g%dx_Cv) ; dealloc_(g%dy_Cu)
606 dealloc_(g%bathyT) ; dealloc_(g%CoriolisBu)
607 dealloc_(g%dF_dx) ; dealloc_(g%dF_dy)
608 dealloc_(g%sin_rot) ; dealloc_(g%cos_rot)
610 if (g%bathymetry_at_vel)
then
611 dealloc_(g%Dblock_u) ; dealloc_(g%Dopen_u)
612 dealloc_(g%Dblock_v) ; dealloc_(g%Dopen_v)
615 deallocate(g%gridLonT) ;
deallocate(g%gridLatT)
616 deallocate(g%gridLonB) ;
deallocate(g%gridLatB)
618 deallocate(g%Domain%mpp_domain)
621 end subroutine mom_grid_end