ROMS
Loading...
Searching...
No Matches
mod_extract.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#ifdef GRID_EXTRACT
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! Grid extration module for I/O purposes. !
14! !
15! Hz Thicknesses (m) of vertical RHO-points. !
16! CosAngler Cosine of curvilinear angle, angler. !
17! SinAngler Sine of curvilinear angle, angler. !
18! angler Angle (radians) between XI-axis and true EAST at !
19! RHO-points. !
20! dmde ETA-derivative of inverse metric factor pm, !
21! d(1/pm)/d(ETA). !
22! dndx XI-derivative of inverse metric factor pn, !
23! d(1/pn)/d(XI). !
24! f Coriolis parameter (1/s). !
25! h Bottom depth (m) at RHO-points. !
26! latp Latitude (degrees_north) at PSI-points. !
27! latr Latitude (degrees_north) at RHO-points. !
28! latu Latitude (degrees_north) at U-points. !
29! latv Latitude (degrees_north) at V-points. !
30! lonp Longitude (degrees_east) at PSI-points. !
31! lonr Longitude (degrees_east) at RHO-points. !
32! lonu Longitude (degrees_east) at U-points. !
33! lonv Longitude (degrees_east) at V-points. !
34! pm Coordinate transformation metric "m" (1/meters) !
35! associated with the differential distances in XI. !
36! pn Coordinate transformation metric "n" (1/meters) !
37! associated with the differential distances in ETA. !
38# ifdef MASKING
39! pmask Slipperiness time-independent mask at PSI-points: !
40! (0=Land, 1=Sea, 2=no-slip). !
41! rmask Time-independent mask at RHO-points (0=Land, 1=Sea). !
42! umask Time-independent mask at U-points (0=Land, 1=Sea). !
43! vmask Time-independent mask at V-points (0=Land, 1=Sea). !
44# endif
45! xp XI-coordinates (m) at PSI-points. !
46! xr XI-coordinates (m) at RHO-points. !
47! xu XI-coordinates (m) at U-points. !
48! xv XI-coordinates (m) at V-points. !
49! yp ETA-coordinates (m) at PSI-points. !
50! yr ETA-coordinates (m) at RHO-points. !
51! yu ETA-coordinates (m) at U-points. !
52! yv ETA-coordinates (m) at V-points. !
53! z_r Actual depths (m) at horizontal RHO-points and !
54! vertical RHO-points. !
55! !
56! Packed global arrays needed for field extraction by interpolation: !
57!
58! Gx_psi Global X-coordinate at PSI-points. !
59! Gx_rho Global X-coordinate at RHO-points. !
60! Gx_u Global X-coordinate at U-points. !
61! Gx_v Global X-coordinate at V-points. !
62! Gy_psi Global Y-coordinate at PSI-points. !
63! Gy_rho Global Y-coordinate at RHO-points. !
64! Gy_u Global Y-coordinate at U-points. !
65! Gy_v Global Y-coordinate at V-points. !
66
67! !
68!=======================================================================
69!
70 USE mod_kinds
71!
72 implicit none
73!
74 PUBLIC :: allocate_extract
75 PUBLIC :: deallocate_extract
76 PUBLIC :: initialize_extract
77!
78!-----------------------------------------------------------------------
79! Define T_EXTRACT structure.
80!-----------------------------------------------------------------------
81!
82 TYPE t_extract
83!
84! Parameters.
85!
86 real(r8), allocatable :: Hmin(:)
87 real(r8), allocatable :: Hmax(:)
88 real(r8), allocatable :: LonMin(:)
89 real(r8), allocatable :: LonMax(:)
90 real(r8), allocatable :: LatMin(:)
91 real(r8), allocatable :: LatMax(:)
92!
93! Nonlinear model state.
94!
95 real(r8), pointer :: angler(:,:)
96 real(r8), pointer :: CosAngler(:,:)
97 real(r8), pointer :: SinAngler(:,:)
98
99# if defined CURVGRID && defined UV_ADV
100 real(r8), pointer :: dmde(:,:)
101 real(r8), pointer :: dndx(:,:)
102# endif
103 real(r8), pointer :: f(:,:)
104 real(r8), pointer :: h(:,:)
105 real(r8), pointer :: latp(:,:)
106 real(r8), pointer :: latr(:,:)
107 real(r8), pointer :: latu(:,:)
108 real(r8), pointer :: latv(:,:)
109 real(r8), pointer :: lonp(:,:)
110 real(r8), pointer :: lonr(:,:)
111 real(r8), pointer :: lonu(:,:)
112 real(r8), pointer :: lonv(:,:)
113 real(r8), pointer :: MyLon(:,:)
114 real(r8), pointer :: omn(:,:)
115 real(r8), pointer :: om_p(:,:)
116 real(r8), pointer :: om_r(:,:)
117 real(r8), pointer :: om_u(:,:)
118 real(r8), pointer :: om_v(:,:)
119 real(r8), pointer :: on_p(:,:)
120 real(r8), pointer :: on_r(:,:)
121 real(r8), pointer :: on_u(:,:)
122 real(r8), pointer :: on_v(:,:)
123 real(r8), pointer :: pm(:,:)
124 real(r8), pointer :: pn(:,:)
125 real(r8), pointer :: xp(:,:)
126 real(r8), pointer :: xr(:,:)
127 real(r8), pointer :: xu(:,:)
128 real(r8), pointer :: xv(:,:)
129 real(r8), pointer :: yp(:,:)
130 real(r8), pointer :: yr(:,:)
131 real(r8), pointer :: yu(:,:)
132 real(r8), pointer :: yv(:,:)
133# ifdef SOLVE3D
134 real(r8), pointer :: Hz(:,:,:)
135 real(r8), pointer :: z_r(:,:,:)
136 real(r8), pointer :: z_v(:,:,:)
137 real(r8), pointer :: z_w(:,:,:)
138# endif
139# ifdef MASKING
140 real(r8), pointer :: pmask(:,:)
141 real(r8), pointer :: rmask(:,:)
142 real(r8), pointer :: umask(:,:)
143 real(r8), pointer :: vmask(:,:)
144# endif
145!
146! Set global geometry 1D arrays, packed in column-major order,
147! fractional I- and J-coordinates of the extracted data with
148! respect the donor grid ti facilitate bilinear or bicubic
149! interpolations. The 3D fields are interpolated level by level,
150! so only horizontal fractional coordinates are needed for each
151! staggered location.
152!
153 real(r8), pointer :: Gx_psi(:)
154 real(r8), pointer :: Gx_rho(:)
155 real(r8), pointer :: Gx_u(:)
156 real(r8), pointer :: Gx_v(:)
157
158 real(r8), pointer :: Gy_psi(:)
159 real(r8), pointer :: Gy_rho(:)
160 real(r8), pointer :: Gy_u(:)
161 real(r8), pointer :: Gy_v(:)
162
163# ifdef MASKING
164 real(r8), pointer :: Gmask_psi(:)
165 real(r8), pointer :: Gmask_rho(:)
166 real(r8), pointer :: Gmask_u(:)
167 real(r8), pointer :: Gmask_v(:)
168# endif
169
170 real(r8), pointer :: Iout_psi(:)
171 real(r8), pointer :: Iout_rho(:)
172 real(r8), pointer :: Iout_u(:)
173 real(r8), pointer :: Iout_v(:)
174
175 real(r8), pointer :: Jout_psi(:)
176 real(r8), pointer :: Jout_rho(:)
177 real(r8), pointer :: Jout_u(:)
178 real(r8), pointer :: Jout_v(:)
179
180 END TYPE t_extract
181!
182 TYPE (T_EXTRACT), allocatable :: EXTRACT(:)
183!
184 CONTAINS
185!
186 SUBROUTINE allocate_extract (ng, Extract_Flag, LBi, UBi, LBj, UBj)
187!
188!=======================================================================
189! !
190! This routine allocates all variables in the module for all nested !
191! grids. !
192! !
193!=======================================================================
194!
195 USE mod_param
196!
197! Local variable declarations.
198!
199 integer, intent(in) :: ng, Extract_Flag
200 integer, intent(in) :: LBi, UBi, LBj, UBj
201!
202! Local variable declarations.
203!
204 integer :: my_size
205!
206 real(r8) :: size2d
207!
208!-----------------------------------------------------------------------
209! Allocate and initialize module variables.
210!-----------------------------------------------------------------------
211!
212 IF (ng.eq.1) allocate ( extract(ngrids) )
213!
214! Parameters.
215!
216 IF (.not.allocated(extract(ng) % Hmin)) &
217 & allocate ( extract(ng) % Hmin(ngrids) )
218
219 IF (.not.allocated(extract(ng) % Hmax)) &
220 & allocate ( extract(ng) % Hmax(ngrids) )
221
222 IF (.not.allocated(extract(ng) % LonMin)) &
223 & allocate ( extract(ng) % LonMin(ngrids) )
224
225 IF (.not.allocated(extract(ng) % LonMax)) &
226 & allocate ( extract(ng) % LonMax(ngrids) )
227
228 IF (.not.allocated(extract(ng) % LatMin)) &
229 & allocate ( extract(ng) % LatMin(ngrids) )
230
231 IF (.not.allocated(extract(ng) % LatMax)) &
232 & allocate ( extract(ng) % LatMax(ngrids) )
233!
234! Set horizontal array size.
235!
236 size2d=real((ubi-lbi+1)*(ubj-lbj+1),r8)
237!
238! Nonlinear model state.
239!
240 allocate ( extract(ng) % angler(lbi:ubi,lbj:ubj) )
241 dmem(ng)=dmem(ng)+size2d
242
243 allocate ( extract(ng) % CosAngler(lbi:ubi,lbj:ubj) )
244 dmem(ng)=dmem(ng)+size2d
245
246 allocate ( extract(ng) % SinAngler(lbi:ubi,lbj:ubj) )
247 dmem(ng)=dmem(ng)+size2d
248
249# if defined CURVGRID && defined UV_ADV
250 allocate ( extract(ng) % dmde(lbi:ubi,lbj:ubj) )
251 dmem(ng)=dmem(ng)+size2d
252
253 allocate ( extract(ng) % dndx(lbi:ubi,lbj:ubj) )
254 dmem(ng)=dmem(ng)+size2d
255# endif
256
257 allocate ( extract(ng) % f(lbi:ubi,lbj:ubj) )
258 dmem(ng)=dmem(ng)+size2d
259
260 allocate ( extract(ng) % h(lbi:ubi,lbj:ubj) )
261 dmem(ng)=dmem(ng)+size2d
262
263 allocate ( extract(ng) % latp(lbi:ubi,lbj:ubj) )
264 dmem(ng)=dmem(ng)+size2d
265
266 allocate ( extract(ng) % latr(lbi:ubi,lbj:ubj) )
267 dmem(ng)=dmem(ng)+size2d
268
269 allocate ( extract(ng) % latu(lbi:ubi,lbj:ubj) )
270 dmem(ng)=dmem(ng)+size2d
271
272 allocate ( extract(ng) % latv(lbi:ubi,lbj:ubj) )
273 dmem(ng)=dmem(ng)+size2d
274
275 allocate ( extract(ng) % lonp(lbi:ubi,lbj:ubj))
276 dmem(ng)=dmem(ng)+size2d
277
278 allocate ( extract(ng) % lonr(lbi:ubi,lbj:ubj))
279 dmem(ng)=dmem(ng)+size2d
280
281 allocate ( extract(ng) % lonu(lbi:ubi,lbj:ubj))
282 dmem(ng)=dmem(ng)+size2d
283
284 allocate ( extract(ng) % lonv(lbi:ubi,lbj:ubj))
285 dmem(ng)=dmem(ng)+size2d
286
287 allocate ( extract(ng) % Mylon(lbi:ubi,lbj:ubj))
288 dmem(ng)=dmem(ng)+size2d
289
290 allocate ( extract(ng) % pm(lbi:ubi,lbj:ubj) )
291 dmem(ng)=dmem(ng)+size2d
292
293 allocate ( extract(ng) % pn(lbi:ubi,lbj:ubj) )
294 dmem(ng)=dmem(ng)+size2d
295
296 allocate ( extract(ng) % xp(lbi:ubi,lbj:ubj) )
297 dmem(ng)=dmem(ng)+size2d
298
299 allocate ( extract(ng) % xr(lbi:ubi,lbj:ubj) )
300 dmem(ng)=dmem(ng)+size2d
301
302 allocate ( extract(ng) % xu(lbi:ubi,lbj:ubj) )
303 dmem(ng)=dmem(ng)+size2d
304
305 allocate ( extract(ng) % xv(lbi:ubi,lbj:ubj) )
306 dmem(ng)=dmem(ng)+size2d
307
308 allocate ( extract(ng) % yp(lbi:ubi,lbj:ubj) )
309 dmem(ng)=dmem(ng)+size2d
310
311 allocate ( extract(ng) % yr(lbi:ubi,lbj:ubj) )
312 dmem(ng)=dmem(ng)+size2d
313
314 allocate ( extract(ng) % yu(lbi:ubi,lbj:ubj) )
315 dmem(ng)=dmem(ng)+size2d
316
317 allocate ( extract(ng) % yv(lbi:ubi,lbj:ubj) )
318 dmem(ng)=dmem(ng)+size2d
319
320# ifdef SOLVE3D
321 allocate ( extract(ng) % Hz(lbi:ubi,lbj:ubj,n(ng)) )
322 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
323
324 allocate ( extract(ng) % z_r(lbi:ubi,lbj:ubj,n(ng)) )
325 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
326
327 allocate ( extract(ng) % z_v(lbi:ubi,lbj:ubj,n(ng)) )
328 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
329
330 allocate ( extract(ng) % z_w(lbi:ubi,lbj:ubj,0:n(ng)) )
331 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
332# endif
333
334# ifdef MASKING
335 allocate ( extract(ng) % pmask(lbi:ubi,lbj:ubj) )
336 dmem(ng)=dmem(ng)+size2d
337
338 allocate ( extract(ng) % rmask(lbi:ubi,lbj:ubj) )
339 dmem(ng)=dmem(ng)+size2d
340
341 allocate ( extract(ng) % umask(lbi:ubi,lbj:ubj) )
342 dmem(ng)=dmem(ng)+size2d
343
344 allocate ( extract(ng) % vmask(lbi:ubi,lbj:ubj) )
345 dmem(ng)=dmem(ng)+size2d
346# endif
347!
348! Allocate global arrays for interpolation
349!
350 IF (extract_flag.ge.1) THEN
351 my_size=(xtr_iobounds(ng)%IUB_psi-xtr_iobounds(ng)%ILB_psi+1)* &
352 & (xtr_iobounds(ng)%JUB_psi-xtr_iobounds(ng)%JLB_psi+1)
353
354 allocate ( extract(ng) % Gx_psi(my_size) )
355 dmem(ng)=dmem(ng)+real(my_size,r8)
356
357 allocate ( extract(ng) % Gy_psi(my_size) )
358 dmem(ng)=dmem(ng)+real(my_size,r8)
359
360# ifdef MASKING
361 allocate ( extract(ng) % Gmask_psi(my_size) )
362 dmem(ng)=dmem(ng)+real(my_size,r8)
363# endif
364
365 allocate ( extract(ng) % Iout_psi(my_size) )
366 dmem(ng)=dmem(ng)+real(my_size,r8)
367
368 allocate ( extract(ng) % Jout_psi(my_size) )
369 dmem(ng)=dmem(ng)+real(my_size,r8)
370!
371 my_size=(xtr_iobounds(ng)%IUB_rho-xtr_iobounds(ng)%ILB_rho+1)* &
372 & (xtr_iobounds(ng)%JUB_rho-xtr_iobounds(ng)%JLB_rho+1)
373
374 allocate ( extract(ng) % Gx_rho(my_size) )
375 dmem(ng)=dmem(ng)+real(my_size,r8)
376
377 allocate ( extract(ng) % Gy_rho(my_size) )
378 dmem(ng)=dmem(ng)+real(my_size,r8)
379
380# ifdef MASKING
381 allocate ( extract(ng) % Gmask_rho(my_size) )
382 dmem(ng)=dmem(ng)+real(my_size,r8)
383# endif
384
385 allocate ( extract(ng) % Iout_rho(my_size) )
386 dmem(ng)=dmem(ng)+real(my_size,r8)
387
388 allocate ( extract(ng) % Jout_rho(my_size) )
389 dmem(ng)=dmem(ng)+real(my_size,r8)
390!
391 my_size=(xtr_iobounds(ng)%IUB_u-xtr_iobounds(ng)%ILB_u+1)* &
392 & (xtr_iobounds(ng)%JUB_u-xtr_iobounds(ng)%JLB_u+1)
393
394 allocate ( extract(ng) % Gx_u(my_size) )
395 dmem(ng)=dmem(ng)+real(my_size,r8)
396
397 allocate ( extract(ng) % Gy_u(my_size) )
398 dmem(ng)=dmem(ng)+real(my_size,r8)
399
400# ifdef MASKING
401 allocate ( extract(ng) % Gmask_u(my_size) )
402 dmem(ng)=dmem(ng)+real(my_size,r8)
403# endif
404
405 allocate ( extract(ng) % Iout_u(my_size) )
406 dmem(ng)=dmem(ng)+real(my_size,r8)
407
408 allocate ( extract(ng) % Jout_u(my_size) )
409 dmem(ng)=dmem(ng)+real(my_size,r8)
410!
411 my_size=(xtr_iobounds(ng)%IUB_v-xtr_iobounds(ng)%ILB_v+1)* &
412 & (xtr_iobounds(ng)%JUB_v-xtr_iobounds(ng)%JLB_v+1)
413
414 allocate ( extract(ng) % Gx_v(my_size) )
415 dmem(ng)=dmem(ng)+real(my_size,r8)
416
417 allocate ( extract(ng) % Gy_v(my_size) )
418 dmem(ng)=dmem(ng)+real(my_size,r8)
419
420# ifdef MASKING
421 allocate ( extract(ng) % Gmask_v(my_size) )
422 dmem(ng)=dmem(ng)+real(my_size,r8)
423# endif
424
425 allocate ( extract(ng) % Iout_v(my_size) )
426 dmem(ng)=dmem(ng)+real(my_size,r8)
427
428 allocate ( extract(ng) % Jout_v(my_size) )
429 dmem(ng)=dmem(ng)+real(my_size,r8)
430 END IF
431!
432 RETURN
433 END SUBROUTINE allocate_extract
434!
435 SUBROUTINE deallocate_extract (ng)
436!
437!=======================================================================
438! !
439! This routine deallocates all variables in the module for all nested !
440! grids. !
441! !
442!=======================================================================
443!
444 USE mod_param, ONLY : ngrids
445# ifdef SUBOBJECT_DEALLOCATION
446 USE destroy_mod, ONLY : destroy
447# endif
448!
449! Imported variable declarations.
450!
451 integer, intent(in) :: ng
452!
453! Local variable declarations.
454!
455 character (len=*), parameter :: MyFile = &
456 & __FILE__//", deallocate_extract"
457
458# ifdef SUBOBJECT_DEALLOCATION
459!
460!-----------------------------------------------------------------------
461! Deallocate each variable in the derived-type T_GRID structure
462! separately.
463!-----------------------------------------------------------------------
464!
465! Nonlinear model state.
466!
467 IF (.not.destroy(ng, extract(ng)%angler, myfile, &
468 & __line__, 'EXTRACT(ng)%angler')) RETURN
469
470 IF (.not.destroy(ng, extract(ng)%CosAngler, myfile, &
471 & __line__, 'EXTRACT(ng)%CosAngler')) RETURN
472
473 IF (.not.destroy(ng, extract(ng)%SinAngler, myfile, &
474 & __line__, 'EXTRACT(ng)%SinAngler')) RETURN
475
476# if defined CURVGRID && defined UV_ADV
477 IF (.not.destroy(ng, extract(ng)%dmde, myfile, &
478 & __line__, 'EXTRACT(ng)%dmde')) RETURN
479
480 IF (.not.destroy(ng, extract(ng)%dndx, myfile, &
481 & __line__, 'EXTRACT(ng)%dndx')) RETURN
482# endif
483
484 IF (.not.destroy(ng, extract(ng)%f, myfile, &
485 & __line__, 'EXTRACT(ng)%f')) RETURN
486
487 IF (.not.destroy(ng, extract(ng)%h, myfile, &
488 & __line__, 'EXTRACT(ng)%h')) RETURN
489
490 IF (.not.destroy(ng, extract(ng)%latp, myfile, &
491 & __line__, 'EXTRACT(ng)%latp')) RETURN
492
493 IF (.not.destroy(ng, extract(ng)%latr, myfile, &
494 & __line__, 'EXTRACT(ng)%latr')) RETURN
495
496 IF (.not.destroy(ng, extract(ng)%latu, myfile, &
497 & __line__, 'EXTRACT(ng)%latu')) RETURN
498
499 IF (.not.destroy(ng, extract(ng)%latv, myfile, &
500 & __line__, 'EXTRACT(ng)%latv')) RETURN
501
502 IF (.not.destroy(ng, extract(ng)%lonp, myfile, &
503 & __line__, 'EXTRACT(ng)%lonp')) RETURN
504
505 IF (.not.destroy(ng, extract(ng)%lonr, myfile, &
506 & __line__, 'EXTRACT(ng)%lonr')) RETURN
507
508 IF (.not.destroy(ng, extract(ng)%lonu, myfile, &
509 & __line__, 'EXTRACT(ng)%lonu')) RETURN
510
511 IF (.not.destroy(ng, extract(ng)%lonv, myfile, &
512 & __line__, 'EXTRACT(ng)%lonv')) RETURN
513
514 IF (.not.destroy(ng, extract(ng)%Mylon, myfile, &
515 & __line__, 'EXTRACT(ng)%Mylon')) RETURN
516
517 IF (.not.destroy(ng, extract(ng)%pm, myfile, &
518 & __line__, 'EXTRACT(ng)%pm')) RETURN
519
520 IF (.not.destroy(ng, extract(ng)%pn, myfile, &
521 & __line__, 'EXTRACT(ng)%pn')) RETURN
522
523
524 IF (.not.destroy(ng, extract(ng)%xp, myfile, &
525 & __line__, 'EXTRACT(ng)%xp')) RETURN
526
527 IF (.not.destroy(ng, extract(ng)%xr, myfile, &
528 & __line__, 'EXTRACT(ng)%xr')) RETURN
529
530 IF (.not.destroy(ng, extract(ng)%xu, myfile, &
531 & __line__, 'EXTRACT(ng)%xu')) RETURN
532
533 IF (.not.destroy(ng, extract(ng)%xv, myfile, &
534 & __line__, 'EXTRACT(ng)%xv')) RETURN
535
536 IF (.not.destroy(ng, extract(ng)%yp, myfile, &
537 & __line__, 'EXTRACT(ng)%yp')) RETURN
538
539 IF (.not.destroy(ng, extract(ng)%yr, myfile, &
540 & __line__, 'EXTRACT(ng)%yr')) RETURN
541
542 IF (.not.destroy(ng, extract(ng)%yu, myfile, &
543 & __line__, 'EXTRACT(ng)%yu')) RETURN
544
545 IF (.not.destroy(ng, extract(ng)%yv, myfile, &
546 & __line__, 'EXTRACT(ng)%yv')) RETURN
547
548# ifdef SOLVE3D
549 IF (.not.destroy(ng, extract(ng)%Hz, myfile, &
550 & __line__, 'EXTRACT(ng)%Hz')) RETURN
551
552 IF (.not.destroy(ng, extract(ng)%z_r, myfile, &
553 & __line__, 'EXTRACT(ng)%z_r')) RETURN
554
555 IF (.not.destroy(ng, extract(ng)%z_v, myfile, &
556 & __line__, 'GRID(ng)%z_v')) RETURN
557
558 IF (.not.destroy(ng, extract(ng)%z_w, myfile, &
559 & __line__, 'GRID(ng)%z_w')) RETURN
560# endif
561
562# ifdef MASKING
563 IF (.not.destroy(ng, extract(ng)%pmask, myfile, &
564 & __line__, 'EXTRACT(ng)%pmask')) RETURN
565
566 IF (.not.destroy(ng, extract(ng)%rmask, myfile, &
567 & __line__, 'EXTRACT(ng)%rmask')) RETURN
568
569 IF (.not.destroy(ng, extract(ng)%umask, myfile, &
570 & __line__, 'EXTRACT(ng)%umask')) RETURN
571
572 IF (.not.destroy(ng, extract(ng)%vmask, myfile, &
573 & __line__, 'EXTRACT(ng)%vmask')) RETURN
574# endif
575# endif
576!
577!-----------------------------------------------------------------------
578! Deallocate derived-type GRID structure.
579!-----------------------------------------------------------------------
580!
581 IF (ng.eq.ngrids) THEN
582 IF (allocated(extract)) deallocate ( extract )
583 END IF
584!
585 RETURN
586 END SUBROUTINE deallocate_extract
587!
588 SUBROUTINE initialize_extract (ng, tile, model)
589!
590!=======================================================================
591! !
592! This routine initialize all variables in the module using first !
593! touch distribution policy. In shared-memory configuration, this !
594! operation actually performs propagation of the "shared arrays" !
595! across the cluster, unless another policy is specified to !
596! override the default. !
597! !
598!=======================================================================
599!
600 USE mod_param
601 USE mod_scalars
602!
603! Imported variable declarations.
604!
605 integer, intent(in) :: ng, tile, model
606!
607! Local variable declarations.
608!
609 integer :: Imin, Imax, Jmin, Jmax
610 integer :: i, j
611# ifdef SOLVE3D
612 integer :: k
613# endif
614
615 real(r8), parameter :: IniVal = 0.0_r8
616 real(r8) :: IniMetricVal
617
618# include "set_bounds.h"
619!
620! Set array initialization range.
621!
622# ifdef DISTRIBUTE
623 imin=xtr_bounds(ng)%LBi(tile)
624 imax=xtr_bounds(ng)%UBi(tile)
625 jmin=xtr_bounds(ng)%LBj(tile)
626 jmax=xtr_bounds(ng)%UBj(tile)
627# else
628 IF (xtr_domain(ng)%Western_Edge(tile)) THEN
629 imin=xtr_bounds(ng)%LBi(tile)
630 ELSE
631 imin=xtr_bounds(ng)%Istr(tile)
632 END IF
633 IF (xtr_domain(ng)%Eastern_Edge(tile)) THEN
634 imax=xtr_bounds(ng)%UBi(tile)
635 ELSE
636 imax=xtr_bounds(ng)%Iend(tile)
637 END IF
638 IF (domain(ng)%Southern_Edge(tile)) THEN
639 jmin=xtr_bounds(ng)%LBj(tile)
640 ELSE
641 jmin=xtr_bounds(ng)%Jstr(tile)
642 END IF
643 IF (domain(ng)%Northern_Edge(tile)) THEN
644 jmax=xtr_bounds(ng)%UBj(tile)
645 ELSE
646 jmax=xtr_bounds(ng)%Jend(tile)
647 END IF
648# endif
649!
650! Set initialization value that it is special in nexting to just
651! load contact points that have not been initialized from the
652! regular physical grid. This is done to make sure that all these
653! important metric values have been set-up correctly.
654!
655# ifdef NESTING
656 inimetricval=spval ! very large value
657# else
658 inimetricval=inival
659# endif
660!
661!-----------------------------------------------------------------------
662! Initialize module variables.
663!-----------------------------------------------------------------------
664!
665! Nonlinear model state.
666!
667 IF ((model.eq.0).or.(model.eq.inlm)) THEN
668 DO j=jmin,jmax
669 DO i=imin,imax
670 extract(ng) % angler(i,j) = inimetricval
671 extract(ng) % CosAngler(i,j) = inival
672 extract(ng) % SinAngler(i,j) = inival
673
674# if defined CURVGRID && defined UV_ADV
675 extract(ng) % dmde(i,j) = inimetricval
676 extract(ng) % dndx(i,j) = inimetricval
677# endif
678 extract(ng) % f(i,j) = inimetricval
679 extract(ng) % h(i,j) = inimetricval
680
681 extract(ng) % latp(i,j) = inival
682 extract(ng) % latr(i,j) = inimetricval
683 extract(ng) % latu(i,j) = inimetricval
684 extract(ng) % latv(i,j) = inimetricval
685 extract(ng) % lonp(i,j) = inival
686 extract(ng) % lonr(i,j) = inimetricval
687 extract(ng) % lonu(i,j) = inimetricval
688 extract(ng) % lonv(i,j) = inimetricval
689 extract(ng) % MyLon(i,j) = inimetricval
690
691 extract(ng) % pm(i,j) = inimetricval
692 extract(ng) % pn(i,j) = inimetricval
693
694 extract(ng) % xp(i,j) = inival
695 extract(ng) % xr(i,j) = inimetricval
696 extract(ng) % xu(i,j) = inimetricval
697 extract(ng) % xv(i,j) = inimetricval
698 extract(ng) % yp(i,j) = inival
699 extract(ng) % yr(i,j) = inimetricval
700 extract(ng) % yu(i,j) = inimetricval
701 extract(ng) % yv(i,j) = inimetricval
702
703# ifdef MASKING
704 extract(ng) % pmask(i,j) = inival
705 extract(ng) % rmask(i,j) = inimetricval
706 extract(ng) % umask(i,j) = inimetricval
707 extract(ng) % vmask(i,j) = inimetricval
708# endif
709 END DO
710
711# ifdef SOLVE3D
712 DO k=1,n(ng)
713 DO i=imin,imax
714 extract(ng) % Hz(i,j,k) = inival
715 extract(ng) % z_r(i,j,k) = inival
716 extract(ng) % z_v(i,j,k) = inival
717 END DO
718 END DO
719 DO k=0,n(ng)
720 DO i=imin,imax
721 extract(ng) % z_w(i,j,k) = inival
722 END DO
723 END DO
724# endif
725 END DO
726 END IF
727!
728 RETURN
729 END SUBROUTINE initialize_extract
730#endif
731 END MODULE mod_extract
integer, parameter inlm
Definition mod_param.F:662
integer, dimension(:), allocatable n
Definition mod_param.F:479
real(r8), dimension(:), allocatable dmem
Definition mod_param.F:137
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer ngrids
Definition mod_param.F:113
real(dp), parameter spval