Line data Source code
1 : module dynamics_vars
2 :
3 : !-----------------------------------------------------------------------
4 : ! CAM fvcore internal variables
5 : !
6 : ! !REVISION HISTORY:
7 : ! 01.06.06 Sawyer Consolidated from various code snippets
8 : ! 03.06.25 Sawyer Cleaned up, used ParPatternCopy (Create)
9 : ! 03.08.05 Sawyer Removed rayf_init and hswf_init, related vars
10 : ! 03.10.22 Sawyer pmgrid removed (now spmd_dyn)
11 : ! 03.11.18 Sawyer Removed set_eta (ak, bk, now read from restart)
12 : ! 03.12.04 Sawyer Moved T_FVDYCORE_GRID here (removed some vars)
13 : ! 04.08.25 Sawyer Removed all module data members, now GRID only
14 : ! 04.10.06 Sawyer Added spmd_dyn vars here; ESMF transpose vars
15 : ! 05.04.12 Sawyer Added support for r4/r8 tracers
16 : ! 05.05.24 Sawyer CAM/GEOS5 merge (removed GEOS_mod dependencies)
17 : ! 05.06.10 Sawyer Scaled down version for CAM (no ESMF)
18 : ! 05.11.10 Sawyer Removed dyn_interface (now in dyn_comp)
19 : ! 06.03.01 Sawyer Removed m_ttrans, q_to_qxy, qxy_to_q, etc.
20 : ! 06.05.09 Sawyer Added CONSV to dyn_state (conserve energy)
21 : ! 06.08.27 Sawyer Removed unused ESMF code for RouteHandle
22 : !-----------------------------------------------------------------------
23 :
24 : use shr_kind_mod, only: r8=>shr_kind_r8
25 : use pmgrid, only: plon, plat, plev
26 :
27 : use decompmodule, only: decomptype
28 : use ghostmodule, only: ghosttype
29 :
30 : use cam_logfile, only: iulog
31 : use cam_abortutils, only: endrun
32 :
33 : #if defined(SPMD)
34 : use parutilitiesmodule, only: parpatterntype, REAL4, INT4
35 : #endif
36 :
37 : implicit none
38 : private
39 : save
40 :
41 : public :: &
42 : t_fvdycore_vars, &
43 : t_fvdycore_grid, &
44 : t_fvdycore_constants, &
45 : t_fvdycore_state, &
46 : grid_vars_init, &
47 : dynamics_clean
48 :
49 : #ifdef SPMD
50 : public :: spmd_vars_init
51 : #endif
52 :
53 : ! T_FVDYCORE_VARS contains the prognostic variables for FVdycore
54 :
55 : type T_FVDYCORE_VARS
56 : real(r8), dimension(:,:,: ), pointer :: U ! U winds (D-grid)
57 : real(r8), dimension(:,:,: ), pointer :: V ! V winds (D-grid)
58 : real(r8), dimension(:,:,: ), pointer :: PT ! scaled virtual pot. temp.
59 : real(r8), dimension(:,:,: ), pointer :: PE ! Pressure at layer edges
60 : real(r8), dimension(:,:,: ), pointer :: PKZ ! P^kappa mean
61 : real(r8), dimension(:,:,:,:), pointer :: tracer ! Tracers
62 : end type T_FVDYCORE_VARS
63 :
64 : ! T_FVDYCORE_GRID contains information about the horizontal and vertical
65 : ! discretization and decompositions.
66 :
67 : type T_FVDYCORE_GRID
68 :
69 : ! PILGRIM communication information
70 :
71 : integer :: twod_decomp = 0 ! 1 for multi-2D decompositions, 0 otherwise
72 : ! To assure that the latitudinal decomposition operates
73 : ! as efficiently as before, a separate parameter "twod_decomp" has
74 : ! been defined; a value of 1 refers to the multi-2D decomposition with
75 : ! transposes; a value of 0 means that the decomposition is effectively
76 : ! one-dimensional, thereby enabling the transpose logic to be skipped;
77 : ! there is an option to force computation of transposes even for case
78 : ! where decomposition is effectively 1-D.
79 :
80 : integer :: npes_xy= 1 ! number of PEs for XY decomposition
81 : integer :: npes_yz= 1 ! number of PEs for YZ decomposition
82 : integer :: myid_y = 0 ! subdomain index (0-based) in latitude (y)
83 : integer :: myid_z = 0 ! subdomain index (0 based) in level (z)
84 : integer :: npr_y = 1 ! number of subdomains in y
85 : integer :: npr_z = 1 ! number of subdomains in z
86 :
87 : integer :: myidxy_x = 0 ! subdomain index (0-based) in longitude (x) (second. decomp.)
88 : integer :: myidxy_y = 0 ! subdomain index (0 based) in latitude (y) (second. decomp.)
89 : integer :: nprxy_x = 1 ! number of subdomains in x (second. decomp.)
90 : integer :: nprxy_y = 1 ! number of subdomains in y (second. decomp.)
91 : integer :: iam = 0 !
92 :
93 : integer :: mod_method = 0 ! 1 for mpi derived types with transposes, 0 for contiguous buffers
94 : integer :: mod_geopk = 0 ! 1 for mpi derived types with transposes, 0 for contiguous buffers
95 : integer :: mod_gatscat = 0 ! 1 for mpi derived types with transposes, 0 for contiguous buffers
96 :
97 : type(decomptype) :: strip2d, strip2dx, strip3dxyz, strip3dxzy, &
98 : strip3dxyzp, strip3zaty, strip3dxzyp, &
99 : strip3yatz, strip3yatzp, strip3zatypt, &
100 : strip3kxyz, strip3kxzy, strip3kxyzp, strip3kxzyp, &
101 : strip3dyz, checker3kxy
102 :
103 : integer :: commdyn ! communicator for all dynamics
104 : integer :: commxy ! communicator for XY decomposition
105 : integer :: commyz ! communicator for YZ decomposition
106 : integer :: commnyz ! communicator for multiple YZ decomposition
107 :
108 : integer :: comm_y ! communicator in latitude
109 : integer :: comm_z ! communicator in vertical
110 : integer :: commxy_x ! communicator in longitude (xy second. decomp.)
111 : integer :: commxy_y ! communicator in latitude (xy second. decomp.)
112 : logical :: geopkdist ! use distributed method for geopotential calculation
113 : ! with 2D decomp.
114 : logical :: geopk16byte ! use Z-parallel distributed method for geopotential
115 : ! calculation with 2D decomp.; otherwise use Z-serial
116 : ! pipeline algorithm when using distributed algoritm
117 : integer :: geopkblocks ! number of stages to use in Z-serial pipeline
118 : ! (non-transpose) geopotential algorithm
119 : integer :: modc_dynrun(4) ! 1: mod_comm irregular underlying communication method for dyn_run/misc
120 : ! 2: mod_comm irregular communication handshaking for dyn_run/misc
121 : ! 3: mod_comm irregular communication send protocol for dyn_run/misc
122 : ! 4: mod_comm irregular communication nonblocking request throttle for dyn_run/misc
123 : integer :: modc_cdcore(4) ! 1: mod_comm irregular underlying communication method for cd_core/geopk
124 : ! 2: mod_comm irregular communication handshaking for cd_core/geopk
125 : ! 3: geopk_d and mod_comm irregular communication send protocol for cd_core/geopk
126 : ! 4: mod_comm irregular communication nonblocking request throttle for cd_core/geopk
127 : integer :: modc_gather(4) ! 1: mod_comm irregular underlying communication method for gather
128 : ! 2: mod_comm irregular communication handshaking for gather
129 : ! 3: mod_comm irregular communication send protocol for gather
130 : ! 4: mod_comm irregular communication nonblocking request throttle for gather
131 : integer :: modc_scatter(4) ! 1: mod_comm irregular underlying communication method for scatter
132 : ! 2: mod_comm irregular communication handshaking for scatter
133 : ! 3: mod_comm irregular communication send protocol for scatter
134 : ! 4: mod_comm irregular communication nonblocking request throttle for scatter
135 : integer :: modc_tracer(4) ! 1: mod_comm irregular underlying communication method for multiple tracers
136 : ! 2: mod_comm irregular communication handshaking for multiple tracers
137 : ! 3: mod_comm irregular communication send protocol for multiple tracers
138 : ! 4: mod_comm irregular communication nonblocking request throttle for multiple tracers
139 : integer :: modc_onetwo ! one or two simultaneous mod_comm irregular communications (excl. tracers)
140 : integer :: modc_tracers ! max number of tracers for simultaneous mod_comm irregular communications
141 :
142 : #if defined(SPMD)
143 : type (ghosttype) :: ghostu_yz, ghostv_yz, ghostpt_yz, &
144 : ghostpe_yz, ghostpkc_yz
145 : type (parpatterntype) :: u_to_uxy, uxy_to_u, v_to_vxy, vxy_to_v, &
146 : ikj_yz_to_xy, ikj_xy_to_yz, &
147 : ijk_yz_to_xy, ijk_xy_to_yz, &
148 : pe_to_pexy, pexy_to_pe, &
149 : pt_to_ptxy, ptxy_to_pt, pkxy_to_pkc, &
150 : r4_xy_to_yz, r4_yz_to_xy, q3_to_qxy3, qxy3_to_q3, &
151 : xy2d_to_yz2d, yz2d_to_xy2d, scatter_3d, gather_3d, &
152 : g_2dxy_r8, g_2dxy_r4, g_2dxy_i4, &
153 : s_2dxy_r8, s_2dxy_r4, s_2dxy_i4, &
154 : g_3dxyz_r8, g_3dxyz_r4, g_3dxyzp_r8, g_3dxyzp_r4, &
155 : s_3dxyz_r8, s_3dxyz_r4, s_3dxyzp_r8, s_3dxyzp_r4
156 : #endif
157 :
158 : ! END PILGRIM communication information
159 :
160 :
161 : integer :: JFIRST = 1 ! Start latitude (exclusive)
162 : integer :: JLAST = plat ! End latitude (exclusive)
163 :
164 : integer :: ng_c = 0 ! Ccore ghosting
165 : integer :: ng_d = 0 ! Dcore ghosting
166 : integer :: ng_s = 0 ! Staggered grid ghosting for
167 : ! certain arrays, max(ng_c+1,ng_d)
168 : ! For 2D decomposition
169 :
170 : integer :: IFIRSTXY = 1 ! Start longitude (exclusive)
171 : integer :: ILASTXY = plon ! End longitude (exclusive)
172 : integer :: JFIRSTXY = 1 ! Start latitude (exclusive)
173 : integer :: JLASTXY = plat ! End latitude (exclusive)
174 :
175 : integer :: IM ! Full longitude dim
176 : integer :: JM ! Full latitude dim (including poles)
177 :
178 : real(r8) :: DL
179 : real(r8) :: DP
180 : real(r8) :: ACAP
181 : real(r8) :: RCAP
182 :
183 : real(r8), dimension(:), pointer :: COSP ! Cosine of lat angle -- volume mean
184 : real(r8), dimension(:), pointer :: SINP ! Sine of lat angle -- volume mean
185 : real(r8), dimension(:), pointer :: COSE ! Cosine at finite volume edge
186 : real(r8), dimension(:), pointer :: SINE ! Sine at finite volume edge
187 : real(r8), dimension(:), pointer :: ACOSP ! Reciprocal of cosine of lat angle
188 :
189 : real(r8), dimension(:), pointer :: ACOSU ! Reciprocal of cosine of lat angle (staggered)
190 :
191 : real(r8), dimension(:), pointer :: COSLON ! Cosine of longitudes - volume center
192 : real(r8), dimension(:), pointer :: SINLON ! Sine of longitudes - volume center
193 : real(r8), dimension(:), pointer :: COSL5 ! Cosine of longitudes - volume center
194 : real(r8), dimension(:), pointer :: SINL5 ! Sine of longitudes - volume center
195 :
196 : ! Variables which are used repeatedly in CD_CORE
197 :
198 : integer :: js2g0
199 : integer :: jn2g0
200 : integer :: jn1g1
201 :
202 : real(r8), pointer :: trigs(:)
203 : real(r8), pointer :: fc(:), f0(:)
204 : real(r8), pointer :: dc(:,:), de(:,:), sc(:), se(:)
205 : real(r8), pointer :: cdx(:,:), cdy(:,:)
206 : real(r8), pointer :: cdx4(:,:), cdy4(:,:) ! for div4 damping
207 : real(r8), pointer :: cdxde(:,:), cdxdp(:,:),cdyde(:,:),cdydp(:,:) ! for del2 damping
208 : real(r8), pointer :: cdxdiv(:,:), cdydiv(:,:), cdtau4(:,:) ! for del2 damping
209 :
210 : real(r8), pointer :: dcdiv4(:,:), dediv4(:,:), scdiv4(:), sediv4(:) ! for div4 damping
211 :
212 : real(r8), pointer :: dtdx(:), dtdxe(:), txe5(:), dtxe5(:)
213 : real(r8), pointer :: dyce(:), dx(:) , rdx(:), cy(:)
214 : real(r8), pointer :: dtdx2(:), dtdx4(:), dxdt(:), dxe(:)
215 : real(r8), pointer :: cye(:), dycp(:), rdxe(:)
216 :
217 : real(r8) :: rdy, dtdy, dydt, dtdy5, tdy5
218 : real(r8) :: dt0 = 0
219 :
220 : integer :: ifax(13)
221 :
222 : real(r8) :: zt_c
223 : real(r8) :: zt_d
224 :
225 : ! This part refers to the vertical grid
226 :
227 : integer :: KM ! Numer of levels
228 : integer :: KMAX ! KM+1 (?)
229 :
230 : ! For 2D decomposition
231 :
232 : integer :: KFIRST = 1 ! Start level (exclusive)
233 : integer :: KLAST = plev ! End level (exclusive)
234 : integer :: KLASTP = plev+1 ! klast+1, except km+1 when klastp=km+1
235 :
236 : integer :: KORD ! monotonicity order for mapping (te_map)
237 : integer :: KS ! Number of true pressure levels (out of KM+1)
238 : real(r8) :: PTOP ! pressure at top (ak(1))
239 : real(r8) :: PINT ! initial pressure (ak(km+1))
240 : real(r8), dimension(:), pointer :: AK ! Sigma mapping
241 : real(r8), dimension(:), pointer :: BK ! Sigma mapping
242 :
243 : ! Tracers
244 :
245 : integer :: NQ ! Number of advected tracers
246 : integer :: NTOTQ ! Total number of tracers (NQ <= NC)
247 :
248 : integer :: ct_overlap ! nonzero for overlap of cd_core and trac2d, 0 otherwise
249 : integer :: trac_decomp ! size of tracer domain decomposition for trac2d
250 :
251 : ! Extra subdomain bounds for cd_core/trac2d overlap and trac2d decomposition
252 : ! Relevant for secondary yz decomposition only; refers back to primary yz decomposition
253 : integer :: JFIRSTCT ! jfirst
254 : integer :: JLASTCT ! jlast
255 : integer :: KFIRSTCT ! kfirst
256 : integer :: KLASTCT ! klast
257 :
258 : ! Bounds for tracer decomposition
259 : integer, dimension(:), pointer :: ktloa ! lower tracer index (global map)
260 : integer, dimension(:), pointer :: kthia ! upper tracer index (global map)
261 : integer :: ktlo ! lower tracer index (local)
262 : integer :: kthi ! upper tracer index (local)
263 :
264 : logical :: high_alt ! high-altitude physics parameters switch
265 :
266 : end type T_FVDYCORE_GRID
267 :
268 : ! Constants used by fvcore
269 : type T_FVDYCORE_CONSTANTS
270 : real(r8) :: pi
271 : real(r8) :: omega ! angular velocity of earth's rotation
272 : real(r8) :: cp ! heat capacity of air at constant pressure
273 : real(r8) :: ae ! radius of the earth (m)
274 : real(r8) :: rair ! Gas constant of the air
275 : real(r8) :: cappa ! Cappa?
276 : real(r8) :: zvir ! RWV/RAIR-1
277 : end type T_FVDYCORE_CONSTANTS
278 :
279 : type t_fvdycore_state
280 : type (t_fvdycore_vars) :: vars
281 : type (t_fvdycore_grid ) :: grid
282 : type (t_fvdycore_constants) :: constants
283 : real(r8) :: DT ! Large time step
284 : real(r8) :: CHECK_DT ! Time step to check maxmin
285 : integer :: ICD, JCD ! Algorithm orders (C Grid)
286 : integer :: IORD, JORD ! Algorithm orders (D Grid)
287 : integer :: KORD ! Vertical order
288 : integer :: TE_METHOD ! method for total energy mapping (te_map)
289 : logical :: CONSV ! dycore conserves tot. en.
290 : integer :: NSPLIT
291 : integer :: NSPLTRAC
292 : integer :: NSPLTVRM
293 : integer :: FILTCW ! filter c-grid winds if positive
294 : integer :: fft_flt ! 0 => FFT/algebraic filter; 1 => FFT filter
295 : integer :: div24del2flag ! 2 for 2nd order div damping, 4 for 4th order div damping,
296 : ! 42 for 4th order div damping plus 2nd order velocity damping
297 : real(r8) :: del2coef ! strength of 2nd order velocity damping
298 : logical :: high_order_top! use normal 4-order PPM calculation near the model top
299 : logical :: am_geom_crrct ! apply correction for angular momentum (AM) conservation in geometry
300 : logical :: am_correction ! apply correction for angular momentum (AM) conservation in SW eqns
301 : logical :: am_fixer ! apply global fixer to conserve AM
302 : logical :: am_fix_lbl ! apply global AM fixer level by level
303 : logical :: am_diag ! turns on an AM diagnostic calculations
304 : end type t_fvdycore_state
305 :
306 : !========================================================================================
307 : contains
308 : !========================================================================================
309 :
310 : #if defined(SPMD)
311 :
312 1536 : subroutine spmd_vars_init(imxy, jmxy, jmyz, kmyz, grid)
313 :
314 : ! Initialize SPMD related variables.
315 : ! !REVISION HISTORY:
316 : ! 02.11.08 Sawyer Creation
317 : ! 03.05.07 Sawyer Use ParPatternCopy for q_to_qxy, etc.
318 : ! 03.07.23 Sawyer Removed dependency on constituents module
319 : ! 03.09.10 Sawyer Reactivated u_to_uxy, etc, redefined pe2pexy
320 : ! 03.11.19 Sawyer Merged in CAM code with mod_method
321 : ! 04.08.25 Sawyer Added GRID as argument
322 :
323 : use decompmodule, only: decompcreate, decompfree
324 : use ghostmodule, only : ghostcreate, ghostfree
325 : use parutilitiesmodule, only : gid, parpatterncreate, parsplit
326 : use mpishorthand, only: mpiint
327 :
328 : ! Arguments
329 :
330 : integer, dimension(:), intent(in) :: imxy
331 : integer, dimension(:), intent(in) :: jmxy
332 : integer, dimension(:), intent(in) :: jmyz
333 : integer, dimension(:), intent(in) :: kmyz
334 :
335 : type( t_fvdycore_grid ), intent(inout) :: grid
336 :
337 : ! local variables:
338 :
339 : type(decomptype) :: global2d, local2d
340 :
341 : integer :: im, jm, km ! Global dims
342 : integer :: nq
343 :
344 : integer :: nprxy_x ! XY decomp - Nr in X
345 : integer :: nprxy_y ! XY decomp - Nr in Y
346 : integer :: npryz_y ! YZ decomp - Nr in Y
347 : integer :: npryz_z ! YZ decomp - Nr in Z
348 : integer :: npes_xy ! XY decomp - Total nr
349 : integer :: npes_yz ! YZ decomp - Total nr
350 :
351 : integer :: commxy ! Communicator for XY decomp
352 : integer :: commyz ! Communicator for YZ decomp
353 : integer :: commnyz ! Communicator for multiple YZ decomp
354 :
355 : integer :: jfirstxy, jlastxy
356 : integer :: jfirst, jlast
357 : integer :: kfirst, klast
358 :
359 : integer :: ng_s, ng_d ! Ghost widths
360 :
361 : integer :: rank_y, rank_z, rankxy_x, rankxy_y ! Currently not used
362 : integer :: size_y, size_z, sizexy_x, sizexy_y ! Currently not used
363 :
364 : integer :: xdist(1), ydistk(1), zdist1(1), zdistxy(1) ! non-distributed dims
365 1536 : integer, allocatable :: xdist_global(:), ydist_global(:)
366 1536 : integer, allocatable :: zdist(:) ! number of levels per subdomain
367 :
368 : integer :: ig1, ig2, jg1, jg2
369 : integer :: jg1d, jg2d, jg1s, jg2s
370 : integer :: kg1, kg2, kg2p
371 :
372 : integer :: ct_overlap
373 : integer :: trac_decomp
374 : integer :: ktmod, ml
375 : integer :: myidmod
376 : integer :: ictstuff(4)
377 : integer :: kquot, krem, krun, kt, mlt
378 : !---------------------------------------------------------------------------
379 :
380 1536 : im = grid%im
381 1536 : jm = grid%jm
382 1536 : km = grid%km
383 1536 : nq = grid%nq
384 :
385 1536 : nprxy_x = grid%nprxy_x
386 1536 : nprxy_y = grid%nprxy_y
387 1536 : npryz_y = grid%npr_y
388 1536 : npryz_z = grid%npr_z
389 1536 : npes_xy = grid%npes_xy
390 1536 : npes_yz = grid%npes_yz
391 :
392 1536 : commxy = grid%commxy
393 1536 : commyz = grid%commyz
394 1536 : commnyz = grid%commnyz
395 :
396 1536 : jfirstxy = grid%jfirstxy
397 1536 : jlastxy = grid%jlastxy
398 :
399 1536 : jfirst = grid%jfirst
400 1536 : jlast = grid%jlast
401 1536 : kfirst = grid%kfirst
402 1536 : klast = grid%klast
403 :
404 1536 : ng_s = grid%ng_s
405 1536 : ng_d = grid%ng_d
406 :
407 : ! Split communicators
408 1536 : call parsplit(commyz, grid%myid_z, gid, grid%comm_y, rank_y, size_y)
409 1536 : call parsplit(commyz, grid%myid_y, gid, grid%comm_z, rank_z, size_z)
410 1536 : call parsplit(commxy, grid%myidxy_y, gid, grid%commxy_x, rankxy_x, sizexy_x)
411 1536 : call parsplit(commxy, grid%myidxy_x, gid, grid%commxy_y, rankxy_y, sizexy_y)
412 :
413 :
414 : ! create decompositions for CAM data structures
415 :
416 4608 : allocate(xdist_global(nprxy_x))
417 4608 : allocate(ydist_global(nprxy_y))
418 4608 : allocate(zdist (npryz_z))
419 1536 : xdist(1) = im
420 :
421 : ! Create PILGRIM decompositions (see decompmodule)
422 :
423 1536 : if (gid < npes_xy) then
424 19968 : xdist_global = 0
425 99840 : ydist_global = 0
426 1536 : xdist_global(1) = im
427 1536 : ydist_global(1) = jm
428 : call decompcreate(nprxy_x, nprxy_y, xdist_global, &
429 1536 : ydist_global, global2d)
430 1536 : call decompcreate(nprxy_x, nprxy_y, imxy, jmxy, local2d )
431 :
432 : ! Decompositions needed on xy decomposition for parpatterncreate
433 :
434 1536 : call decompcreate( 1, npryz_y, xdist, jmyz, grid%strip2d )
435 : call decompcreate( 1, npryz_y, npryz_z, xdist, &
436 1536 : jmyz, kmyz, grid%strip3dxyz )
437 : call decompcreate( "xzy", 1, npryz_z, grid%npr_y, xdist, &
438 1536 : kmyz, jmyz, grid%strip3dxzy )
439 :
440 : ! For y communication within z subdomain (klast version)
441 : ! Use myidmod to have valid index for inactive processes
442 : ! for smaller yz decomposition
443 1536 : myidmod = mod(grid%myid_z, grid%npr_z) ! = myid_z for active yz process
444 1536 : zdist1(1) = kmyz(myidmod+1)
445 : call decompcreate( 1, npryz_y, 1, xdist, jmyz, zdist1, &
446 1536 : grid%strip3yatz )
447 :
448 : ! For z communication within y subdomain
449 :
450 1536 : ydistk(1) = jmyz(grid%myid_y+1)
451 : call decompcreate( 1, 1, npryz_z, xdist, ydistk, kmyz, &
452 1536 : grid%strip3zaty )
453 :
454 : ! Arrays dimensioned plev+1
455 :
456 19968 : zdist(:) = kmyz(:)
457 1536 : zdist(npryz_z) = kmyz(npryz_z) + 1
458 : call decompcreate( 1, npryz_y, npryz_z, xdist, jmyz, zdist,&
459 1536 : grid%strip3dxyzp )
460 : call decompcreate( "xzy", 1, npryz_z, npryz_y, &
461 1536 : xdist, zdist, jmyz, grid%strip3dxzyp )
462 :
463 : ! Arrays dimensioned plev+1, within y subdomain
464 :
465 1536 : ydistk(1) = jmyz(grid%myid_y+1)
466 : call decompcreate( "xzy", 1, npryz_z, 1, xdist, zdist, ydistk, &
467 1536 : grid%strip3zatypt )
468 :
469 : ! For y communication within z subdomain (klast+1 version)
470 : ! Use myidmod to have valid index for inactive processes
471 : ! for smaller yz decomposition
472 1536 : myidmod = mod(grid%myid_z, grid%npr_z) ! = myid_z for active yz process
473 1536 : zdist1(1) = kmyz(myidmod+1)
474 : call decompcreate( 1, npryz_y, 1, xdist, jmyz, zdist1, &
475 1536 : grid%strip3yatzp )
476 :
477 : ! For the 2D XY-YZ data transfer, we need a short 3D array
478 19968 : zdist(:) = 1 ! One copy on each z PE set
479 : call decompcreate( 1, npryz_y, npryz_z, &
480 1536 : xdist, jmyz, zdist, grid%strip3dyz )
481 : end if
482 :
483 : ! Secondary xy decomposition
484 :
485 1536 : if (grid%twod_decomp == 1) then
486 1536 : if (gid < npes_xy) then
487 1536 : zdistxy(1) = npryz_z ! All npr_z copies on 1 PE
488 : call decompcreate( nprxy_x, nprxy_y, 1, &
489 1536 : imxy, jmxy, zdistxy, grid%checker3kxy )
490 1536 : zdistxy(1) = km
491 : call decompcreate( nprxy_x, nprxy_y, 1, &
492 1536 : imxy, jmxy, zdistxy, grid%strip3kxyz )
493 : call decompcreate( "xzy", nprxy_x, 1, nprxy_y, &
494 1536 : imxy, zdistxy, jmxy, grid%strip3kxzy )
495 :
496 1536 : zdistxy(1) = zdistxy(1) + 1
497 : call decompcreate( nprxy_x, nprxy_y, 1, &
498 1536 : imxy, jmxy, zdistxy, grid%strip3kxyzp )
499 : call decompcreate( "xzy", nprxy_x, 1, nprxy_y, &
500 1536 : imxy, zdistxy, jmxy, grid%strip3kxzyp )
501 1536 : zdistxy(1) = jlastxy - jfirstxy + 1
502 1536 : call decompcreate( nprxy_x, 1, imxy, zdistxy, grid%strip2dx )
503 : end if
504 : end if
505 :
506 1536 : deallocate(zdist)
507 1536 : deallocate(ydist_global)
508 1536 : deallocate(xdist_global)
509 :
510 1536 : if ( grid%twod_decomp == 1 ) then
511 :
512 : ! Initialize ghost regions
513 :
514 : ! Set limits for ghostcreate
515 1536 : ig1 = 1
516 1536 : ig2 = im
517 1536 : jg1 = jfirst
518 1536 : jg2 = jlast
519 1536 : jg1d = jfirst-ng_d
520 1536 : jg1s = jfirst-ng_s
521 1536 : jg2d = jlast+ng_d
522 1536 : jg2s = jlast+ng_s
523 1536 : kg1 = kfirst
524 1536 : kg2 = klast
525 1536 : kg2p = klast+1
526 :
527 : ! Call ghostcreate with null ranges for non-yz processes
528 1536 : if (gid >= npes_yz) then
529 0 : ig1 = im/2
530 0 : ig2 = ig1 - 1
531 0 : jg1 = (jfirst+jlast)/2
532 0 : jg2 = jg1 - 1
533 0 : jg1d = jg1
534 0 : jg1s = jg1
535 0 : jg2d = jg2
536 0 : jg2s = jg2
537 0 : kg1 = (kfirst+klast)/2
538 0 : kg2 = kg1 - 1
539 0 : kg2p = kg2
540 : end if
541 :
542 : ! Ghosted decompositions needed on xy decomposition for parpatterncreate
543 1536 : if (gid < npes_xy) then
544 : call ghostcreate( grid%strip3dxyz, gid, im, ig1, ig2, .true., &
545 : jm, jg1d, jg2s, .false., &
546 1536 : km, kg1, kg2, .false., grid%ghostu_yz )
547 : call ghostcreate( grid%strip3dxyz, gid, im, ig1, ig2, .true., &
548 : jm, jg1s, jg2d, .false., &
549 1536 : km, kg1, kg2, .false., grid%ghostv_yz )
550 : call ghostcreate( grid%strip3dxyz, gid, im, ig1, ig2, .true., &
551 : jm, jg1d, jg2d, .false., &
552 1536 : km, kg1, kg2, .false., grid%ghostpt_yz )
553 : call ghostcreate( grid%strip3dxzyp, gid, im, ig1, ig2, .true., &
554 : km+1, kg1, kg2p, .false., &
555 1536 : jm, jg1, jg2, .false., grid%ghostpe_yz)
556 : call ghostcreate( grid%strip3dxyzp, gid, im, ig1, ig2, .true., &
557 : jm, jg1, jg2, .false., &
558 1536 : km+1, kg1, kg2p, .false., grid%ghostpkc_yz)
559 : end if
560 :
561 : ! Initialize transposes
562 :
563 1536 : if (gid < npes_xy) then
564 : call parpatterncreate(commxy, grid%ghostu_yz, grid%strip3kxyz, &
565 1536 : grid%u_to_uxy, mod_method=grid%mod_method)
566 : call parpatterncreate(commxy, grid%strip3kxyz,grid%ghostu_yz, &
567 1536 : grid%uxy_to_u, mod_method=grid%mod_method)
568 : call parpatterncreate(commxy, grid%ghostv_yz, grid%strip3kxyz, &
569 1536 : grid%v_to_vxy, mod_method=grid%mod_method)
570 : call parpatterncreate(commxy, grid%strip3kxyz, grid%ghostv_yz, &
571 1536 : grid%vxy_to_v, mod_method=grid%mod_method)
572 : call parpatterncreate(commxy, grid%strip3dxyz, grid%strip3kxyz,&
573 1536 : grid%ijk_yz_to_xy, mod_method=grid%mod_method)
574 : call parpatterncreate(commxy, grid%strip3kxyz, grid%strip3dxyz,&
575 1536 : grid%ijk_xy_to_yz, mod_method=grid%mod_method)
576 : call parpatterncreate(commxy, grid%strip3dxzy, grid%strip3kxzy,&
577 1536 : grid%ikj_yz_to_xy, mod_method=grid%mod_method)
578 : call parpatterncreate(commxy, grid%strip3kxzy, grid%strip3dxzy,&
579 1536 : grid%ikj_xy_to_yz, mod_method=grid%mod_method)
580 :
581 : ! Note PE <-> PEXY has been redefined for PEXY ijk, but PE ikj
582 :
583 : call parpatterncreate(commxy, grid%ghostpe_yz, grid%strip3kxzyp, &
584 1536 : grid%pe_to_pexy, mod_method=grid%mod_method)
585 : call parpatterncreate(commxy, grid%strip3kxzyp, grid%ghostpe_yz, &
586 1536 : grid%pexy_to_pe, mod_method=grid%mod_method)
587 : call parpatterncreate(commxy, grid%ghostpt_yz, grid%strip3kxyz, &
588 1536 : grid%pt_to_ptxy, mod_method=grid%mod_method)
589 : call parpatterncreate(commxy, grid%strip3kxyz, grid%ghostpt_yz, &
590 1536 : grid%ptxy_to_pt, mod_method=grid%mod_method)
591 : call parpatterncreate(commxy, grid%strip3dxyz, grid%strip3kxyz, &
592 : grid%r4_yz_to_xy, mod_method=grid%mod_method, &
593 1536 : T = REAL4 )
594 : call parpatterncreate(commxy, grid%strip3kxyz, grid%strip3dxyz, &
595 : grid%r4_xy_to_yz, mod_method=grid%mod_method, &
596 1536 : T = REAL4 )
597 : call parpatterncreate(commxy, grid%strip3kxyzp, grid%ghostpkc_yz, &
598 1536 : grid%pkxy_to_pkc, mod_method=grid%mod_method)
599 :
600 : ! These are for 'transposing' 2D arrays from XY YZ
601 : call parpatterncreate(commxy, grid%checker3kxy, grid%strip3dyz, &
602 1536 : grid%xy2d_to_yz2d, mod_method=grid%mod_method)
603 : call parpatterncreate(commxy, grid%strip3dyz, grid%checker3kxy, &
604 1536 : grid%yz2d_to_xy2d, mod_method=grid%mod_method)
605 : end if
606 :
607 : ! Free unneeded decompositions
608 :
609 1536 : call decompfree(grid%strip3dxzyp)
610 1536 : call decompfree(grid%strip3dyz)
611 1536 : call decompfree(grid%strip3yatz)
612 1536 : call decompfree(grid%strip3yatzp)
613 1536 : call decompfree(grid%strip3zaty)
614 1536 : call decompfree(grid%strip3zatypt)
615 1536 : call decompfree(grid%strip3kxyz)
616 1536 : call decompfree(grid%strip3kxzy)
617 1536 : call decompfree(grid%strip3kxyzp)
618 1536 : call decompfree(grid%strip3kxzyp)
619 1536 : call decompfree(grid%checker3kxy)
620 :
621 1536 : call ghostfree(grid%ghostu_yz)
622 1536 : call ghostfree(grid%ghostv_yz)
623 1536 : call ghostfree(grid%ghostpt_yz)
624 1536 : call ghostfree(grid%ghostpe_yz)
625 1536 : call ghostfree(grid%ghostpkc_yz)
626 :
627 : end if
628 :
629 : ! Define scatter and gather patterns for 2D and 3D unghosted arrays
630 :
631 1536 : if (gid < npes_xy) then
632 : call parpatterncreate( commxy, global2d, local2d, grid%s_2dxy_r8, &
633 1536 : mod_method=grid%mod_gatscat )
634 : call parpatterncreate( commxy, local2d, global2d, grid%g_2dxy_r8, &
635 1536 : mod_method=grid%mod_gatscat )
636 :
637 : call parpatterncreate( commxy, global2d, local2d, grid%s_2dxy_r4, &
638 1536 : mod_method=grid%mod_gatscat, T = REAL4 )
639 : call parpatterncreate( commxy, local2d, global2d, grid%g_2dxy_r4, &
640 1536 : mod_method=grid%mod_gatscat, T = REAL4 )
641 :
642 : call parpatterncreate( commxy, global2d, local2d, grid%s_2dxy_i4, &
643 1536 : mod_method=grid%mod_gatscat, T = INT4 )
644 : call parpatterncreate( commxy, local2d, global2d, grid%g_2dxy_i4, &
645 1536 : mod_method=grid%mod_gatscat, T = INT4 )
646 :
647 : ! 3D XYZ patterns, will replace XZY patterns eventually
648 :
649 1536 : call parpatterncreate( commxy, grid%s_2dxy_r8, grid%s_3dxyz_r8, km )
650 1536 : call parpatterncreate( commxy, grid%g_2dxy_r8, grid%g_3dxyz_r8, km )
651 1536 : call parpatterncreate( commxy, grid%s_2dxy_r8, grid%s_3dxyzp_r8, km+1 )
652 1536 : call parpatterncreate( commxy, grid%g_2dxy_r8, grid%g_3dxyzp_r8, km+1 )
653 :
654 1536 : call parpatterncreate( commxy, grid%s_2dxy_r4, grid%s_3dxyz_r4, km )
655 1536 : call parpatterncreate( commxy, grid%g_2dxy_r4, grid%g_3dxyz_r4, km )
656 1536 : call parpatterncreate( commxy, grid%s_2dxy_r4, grid%s_3dxyzp_r4, km+1 )
657 1536 : call parpatterncreate( commxy, grid%g_2dxy_r4, grid%g_3dxyzp_r4, km+1 )
658 :
659 1536 : call decompfree( global2d )
660 1536 : call decompfree( local2d )
661 : end if
662 :
663 : ! Secondary subdomain limits for cd_core/trac2d overlap and trac2d decomposition
664 :
665 1536 : ct_overlap = grid%ct_overlap
666 1536 : trac_decomp = grid%trac_decomp
667 1536 : grid%jfirstct = grid%jfirst
668 1536 : grid%jlastct = grid%jlast
669 1536 : grid%kfirstct = grid%kfirst
670 1536 : grid%klastct = grid%klast
671 1536 : if (ct_overlap > 0) then
672 : mlt = 2
673 1536 : elseif (trac_decomp .gt. 1) then
674 : mlt = trac_decomp
675 : else
676 : mlt = 1
677 : end if
678 :
679 : if (mlt > 1) then
680 0 : if (gid < npes_yz) then
681 0 : ictstuff(1) = grid%jfirstct
682 0 : ictstuff(2) = grid%jlastct
683 0 : ictstuff(3) = grid%kfirstct
684 0 : ictstuff(4) = grid%klastct
685 0 : do ml = 2, mlt
686 0 : call mpisend(ictstuff, 4, mpiint, gid+(ml-1)*npes_yz, gid+(ml-1)*npes_yz, commnyz)
687 : enddo
688 0 : elseif (gid < mlt*npes_yz) then
689 0 : ktmod = gid/npes_yz
690 0 : call mpirecv(ictstuff, 4, mpiint, gid-ktmod*npes_yz, gid, commnyz)
691 0 : grid%jfirstct = ictstuff(1)
692 0 : grid%jlastct = ictstuff(2)
693 0 : grid%kfirstct = ictstuff(3)
694 0 : grid%klastct = ictstuff(4)
695 : end if
696 : end if
697 :
698 1536 : if (trac_decomp .gt. 1) then
699 0 : kquot = nq / trac_decomp
700 0 : krem = nq - kquot * trac_decomp
701 0 : krun = 0
702 0 : do kt = 1, krem
703 0 : grid%ktloa(kt) = krun + 1
704 0 : krun = krun + kquot + 1
705 0 : grid%kthia(kt) = krun
706 : enddo
707 0 : do kt = krem+1, trac_decomp
708 0 : grid%ktloa(kt) = krun + 1
709 0 : krun = krun + kquot
710 0 : grid%kthia(kt) = krun
711 : enddo
712 0 : ktmod = gid/npes_yz + 1
713 0 : ktmod = min(ktmod, trac_decomp)
714 0 : grid%ktlo = grid%ktloa(ktmod)
715 0 : grid%kthi = grid%kthia(ktmod)
716 : endif
717 :
718 3072 : end subroutine spmd_vars_init
719 :
720 : #endif
721 :
722 : !========================================================================================
723 :
724 1536 : subroutine grid_vars_init(pi, ae, om, dt, fft_flt, &
725 : am_geom_crrct, grid)
726 :
727 : ! Initialize FV specific GRID vars
728 : !
729 : ! !REVISION HISTORY:
730 : ! 00.01.10 Grant Creation using code from SJ Lin
731 : ! 01.06.06 Sawyer Modified for dynamics_vars
732 : ! 04.08.25 Sawyer Now updates GRID
733 : ! 05.06.30 Sawyer Added initializations from cd_core
734 : ! 06.09.15 Sawyer PI now passed as argument
735 :
736 : use pft_module, only: pftinit, pft_cf
737 :
738 : ! Arguments
739 : real(r8), intent(in) :: pi
740 : real(r8), intent(in) :: ae ! radius of the earth (m)
741 : real(r8), intent(in) :: om ! angular velocity of earth's rotation
742 : real(r8), intent(in) :: dt
743 : integer, intent(in) :: fft_flt
744 : logical, intent(in) :: am_geom_crrct
745 :
746 : type( T_FVDYCORE_GRID ), intent(inout) :: grid
747 :
748 : ! local variables:
749 : integer :: im
750 : integer :: jm
751 :
752 : real(r8) :: ph5 ! This is to ensure 64-bit for any choice of r8
753 :
754 : integer :: i, j, imh
755 : real(r8) :: zam5, zamda
756 : integer :: js2g0, jn2g0, jn1g1, js2gc, jn1gc
757 : integer :: js2gs, jn2gd, jn1gs
758 :
759 1536 : real(r8), pointer :: cosp(:), sinp(:), cose(:), sine(:), acosp(:), acosu(:)
760 1536 : real(r8), pointer :: coslon(:), sinlon(:), cosl5(:), sinl5(:)
761 :
762 : real(r8) :: rat, ycrit
763 : real(r8) :: dt5
764 :
765 : character(len=*), parameter :: sub='grid_vars_init'
766 : !---------------------------------------------------------------------------
767 :
768 1536 : im = grid%im
769 1536 : jm = grid%jm
770 :
771 1536 : grid%dl = (pi+pi)/im
772 1536 : grid%dp = pi/(jm-1)
773 :
774 4608 : allocate(grid%cosp(jm))
775 3072 : allocate(grid%sinp(jm))
776 3072 : allocate(grid%cose(jm))
777 3072 : allocate(grid%sine(jm))
778 3072 : allocate(grid%acosp(jm))
779 3072 : allocate(grid%acosu(jm))
780 :
781 4608 : allocate(grid%coslon(im))
782 3072 : allocate(grid%sinlon(im))
783 3072 : allocate(grid%cosl5(im))
784 3072 : allocate(grid%sinl5(im))
785 :
786 1536 : cosp => grid%cosp
787 1536 : sinp => grid%sinp
788 1536 : cose => grid%cose
789 1536 : sine => grid%sine
790 1536 : acosp => grid%acosp
791 1536 : acosu => grid%acosu
792 :
793 1536 : coslon => grid%coslon
794 1536 : sinlon => grid%sinlon
795 1536 : cosl5 => grid%cosl5
796 1536 : sinl5 => grid%sinl5
797 :
798 : ! philosophy below: edge values = local true values; centred values = area averages
799 :
800 294912 : do j = 2, jm
801 293376 : ph5 = -0.5_r8*pi + ((j-1)-0.5_r8)*(pi/(jm-1))
802 294912 : sine(j) = sin(ph5)
803 : end do
804 :
805 : ! cos(theta) at cell center distretized as
806 : !
807 : ! cos(theta) = d(sin(theta))/d(theta)
808 :
809 1536 : cosp( 1) = 0._r8
810 1536 : cosp(jm) = 0._r8
811 293376 : do j = 2, jm-1
812 293376 : cosp(j) = (sine(j+1)-sine(j)) / grid%dp
813 : end do
814 :
815 : ! Define cosine at edges..
816 :
817 1536 : if (am_geom_crrct) then
818 0 : do j = 2, jm
819 0 : ph5 = -0.5_r8*pi + ((j-1)-0.5_r8)*(pi/(jm-1._r8))
820 0 : cose(j) = cos(ph5)
821 : end do
822 : else
823 294912 : do j = 2, jm
824 294912 : cose(j) = 0.5_r8 * (cosp(j-1) + cosp(j)) ! dsine/dpe between j+1 and j-1
825 : end do
826 : end if
827 1536 : cose(1) = cose(2)
828 :
829 293376 : do j = 2, jm-1
830 293376 : acosu(j) = 2._r8 / (cose(j) + cose(j+1))
831 : end do
832 :
833 1536 : sinp( 1) = -1._r8
834 1536 : sinp(jm) = 1._r8
835 1536 : if (am_geom_crrct) then
836 0 : do j = 2, jm-1
837 0 : sinp(j) = (cose(j) - cose(j+1))/grid%dp ! sqrt(cosp^2+sinp^2)=1
838 : end do
839 : else
840 293376 : do j = 2, jm-1
841 293376 : sinp(j) = 0.5_r8 * (sine(j) + sine(j+1)) ! 2*sinp*cosp*dp=d(cose^2)
842 : end do
843 : end if
844 :
845 : ! Pole cap area and inverse
846 1536 : grid%acap = im*(1._r8+sine(2)) / grid%dp
847 1536 : grid%rcap = 1._r8 / grid%acap
848 :
849 1536 : imh = im/2
850 1536 : if (im /= 2*imh) then
851 0 : write(iulog,*) sub//': ERROR: im must be an even integer'
852 0 : call endrun(sub//': ERROR: im must be an even integer')
853 : end if
854 :
855 : ! Define logitude at the center of the volume
856 : ! i=1, Zamda = -pi
857 :
858 222720 : do i = 1, imh
859 221184 : zam5 = ((i-1)-0.5_r8) * grid%dl
860 221184 : cosl5(i) = cos(zam5)
861 221184 : cosl5(i+imh) = -cosl5(i)
862 221184 : sinl5(i) = sin(zam5)
863 221184 : sinl5(i+imh) = -sinl5(i)
864 221184 : zamda = (i-1)*grid%dl
865 221184 : coslon(i) = cos(zamda)
866 221184 : coslon(i+imh) = -coslon(i)
867 221184 : sinlon(i) = sin(zamda)
868 222720 : sinlon(i+imh) = -sinlon(i)
869 : end do
870 :
871 1536 : acosp( 1) = grid%rcap * im
872 1536 : acosp(jm) = grid%rcap * im
873 293376 : do j = 2, jm-1
874 293376 : acosp(j) = 1._r8 / cosp(j)
875 : enddo
876 :
877 : ! cd_core initializations
878 :
879 4608 : allocate(grid%dtdx(jm))
880 3072 : allocate(grid%dtdx2(jm))
881 3072 : allocate(grid%dtdx4(jm))
882 3072 : allocate(grid%dtdxe(jm))
883 3072 : allocate(grid%dxdt(jm))
884 3072 : allocate(grid%dxe(jm))
885 3072 : allocate(grid%cye(jm))
886 3072 : allocate(grid%dycp(jm))
887 3072 : allocate(grid%rdxe(jm))
888 3072 : allocate(grid%txe5(jm))
889 3072 : allocate(grid%dtxe5(jm))
890 3072 : allocate(grid%dyce(jm))
891 3072 : allocate(grid%dx(jm))
892 3072 : allocate(grid%rdx(jm))
893 3072 : allocate(grid%cy(jm))
894 :
895 1536 : js2g0 = max(2,grid%jfirst)
896 1536 : jn2g0 = min(jm-1,grid%jlast)
897 1536 : jn1g1 = min(jm,grid%jlast+1)
898 1536 : js2gc = max(2,grid%jfirst-grid%ng_c) ! NG lats on S (starting at 2)
899 1536 : jn1gc = min(jm,grid%jlast+grid%ng_c) ! ng_c lats on N (ending at jm)
900 :
901 1536 : grid%js2g0 = js2g0
902 1536 : grid%jn2g0 = jn2g0
903 1536 : grid%jn1g1 = jn1g1
904 :
905 1536 : js2gs = max(2,grid%jfirst-grid%ng_s)
906 1536 : jn2gd = min(jm-1,grid%jlast+grid%ng_d)
907 1536 : jn1gs = min(jm,grid%jlast+grid%ng_s)
908 :
909 4608 : allocate(grid%sc(js2g0:jn2g0))
910 4608 : allocate(grid%se(js2g0:jn1g1))
911 6144 : allocate(grid%dc(im,js2g0:jn2g0))
912 6144 : allocate(grid%de(im,js2g0:jn1g1))
913 :
914 4608 : allocate(grid%scdiv4(js2gs:jn2gd)) !for filtering of u and v in div4 damping
915 4608 : allocate(grid%sediv4(js2gs:jn1gs)) !for filtering of u and v in div4 damping
916 6144 : allocate(grid%dcdiv4(im,js2gs:jn2gd))!for filtering of u and v in div4 damping
917 6144 : allocate(grid%dediv4(im,js2gs:jn1gs))!for filtering of u and v in div4 damping
918 :
919 1536 : call pftinit(im, fft_flt)
920 :
921 : ! Determine ycrit such that effective DX >= DY
922 1536 : rat = real(im,r8)/real(2*(jm-1),r8)
923 1536 : ycrit = acos( min(0.81_r8, rat) ) * (180._r8/pi)
924 :
925 : call pft_cf(im, jm, js2g0, jn2g0, jn1g1, &
926 : grid%sc, grid%se, grid%dc, grid%de, &
927 1536 : grid%cosp, grid%cose, ycrit)
928 :
929 : !for filtering of u and v in div4 damping
930 : call pft_cf(im, jm, js2gs, jn2gd, jn1gs, &
931 : grid%scdiv4, grid%sediv4, grid%dcdiv4, grid%dediv4, &
932 1536 : grid%cosp, grid%cose, ycrit)
933 :
934 6144 : allocate( grid%cdx (js2g0:jn1g1,grid%kfirst:grid%klast) )
935 4608 : allocate( grid%cdy (js2g0:jn1g1,grid%kfirst:grid%klast) )
936 :
937 4608 : allocate( grid%cdx4 (js2g0:jn1g1,grid%kfirst:grid%klast) )!for div4 damping
938 4608 : allocate( grid%cdy4 (js2g0:jn1g1,grid%kfirst:grid%klast) )!for div4 damping
939 :
940 4608 : allocate( grid%cdxde (js2g0:jn1g1,grid%kfirst:grid%klast) )!for del2 damping
941 4608 : allocate( grid%cdxdp (js2g0:jn1g1,grid%kfirst:grid%klast) )!for del2 damping
942 4608 : allocate( grid%cdyde (js2g0:jn1g1,grid%kfirst:grid%klast) )!for del2 damping
943 4608 : allocate( grid%cdydp (js2g0:jn1g1,grid%kfirst:grid%klast) )!for del2 damping
944 :
945 6144 : allocate( grid%cdxdiv(jm,grid%kfirst:grid%klast) )!for div4 damping
946 4608 : allocate( grid%cdydiv(jm,grid%kfirst:grid%klast) )!for div4 damping
947 4608 : allocate( grid%cdtau4(js2g0:jn1g1,grid%kfirst:grid%klast) )!for div4 damping
948 :
949 4608 : allocate( grid%f0(grid%jfirst-grid%ng_s-1:grid%jlast+grid%ng_d) )
950 4608 : allocate( grid%fc(js2gc:jn1gc) )
951 :
952 16704 : do j = max(1,grid%jfirst-grid%ng_s-1), min(jm,grid%jlast+grid%ng_d)
953 16704 : grid%f0(j) = (om+om)*grid%sinp(j)
954 : end do
955 :
956 : ! Compute coriolis parameter at cell corners.
957 :
958 1536 : if (am_geom_crrct) then
959 0 : do j = js2gc, jn1gc
960 0 : grid%fc(j) = (om+om)*grid%sine(j)
961 : end do
962 : else
963 12168 : do j = js2gc, jn1gc ! Not the issue with ng_c = ng_d
964 12168 : grid%fc(j) = 0.5_r8*(grid%f0(j) + grid%f0(j-1))
965 : end do
966 : end if
967 :
968 1536 : grid%dt0 = 0._r8
969 1536 : dt5 = 0.5_r8*dt
970 :
971 1536 : grid%rdy = 1._r8/(ae*grid%dp)
972 1536 : grid%dtdy = dt *grid%rdy
973 1536 : grid%dtdy5 = dt5*grid%rdy
974 1536 : grid%dydt = (ae*grid%dp) / dt
975 1536 : grid%tdy5 = 0.5_r8/grid%dtdy
976 :
977 1536 : end subroutine grid_vars_init
978 :
979 : !========================================================================================
980 :
981 0 : subroutine dynamics_clean(grid)
982 :
983 : ! Arguments
984 : type(t_fvdycore_grid), intent(inout) :: grid
985 :
986 : ! Temporary data structures
987 :
988 0 : if(associated(GRID%SINLON )) deallocate(GRID%SINLON)
989 0 : if(associated(GRID%COSLON )) deallocate(GRID%COSLON)
990 0 : if(associated(GRID%SINL5 )) deallocate(GRID%SINL5)
991 0 : if(associated(GRID%COSL5 )) deallocate(GRID%COSL5)
992 :
993 0 : if(associated(GRID%ACOSP )) deallocate(GRID%ACOSP)
994 0 : if(associated(GRID%ACOSU )) deallocate(GRID%ACOSU)
995 0 : if(associated(GRID%SINP )) deallocate(GRID%SINP)
996 0 : if(associated(GRID%COSP )) deallocate(GRID%COSP)
997 0 : if(associated(GRID%SINE )) deallocate(GRID%SINE)
998 0 : if(associated(GRID%COSE )) deallocate(GRID%COSE)
999 0 : if(associated(GRID%AK )) deallocate(GRID%AK)
1000 0 : if(associated(GRID%BK )) deallocate(GRID%BK)
1001 :
1002 : ! cd_core variables
1003 :
1004 0 : if(associated( grid%dtdx )) deallocate(grid%dtdx)
1005 0 : if(associated( grid%dtdx2 )) deallocate(grid%dtdx2)
1006 0 : if(associated( grid%dtdx4 )) deallocate(grid%dtdx4)
1007 0 : if(associated( grid%dtdxe )) deallocate(grid%dtdxe)
1008 0 : if(associated( grid%dxdt )) deallocate(grid%dxdt)
1009 0 : if(associated( grid%dxe )) deallocate(grid%dxe)
1010 0 : if(associated( grid%cye )) deallocate(grid%cye)
1011 0 : if(associated( grid%dycp )) deallocate(grid%dycp)
1012 0 : if(associated( grid%rdxe )) deallocate(grid%rdxe)
1013 0 : if(associated( grid%txe5 )) deallocate(grid%txe5)
1014 0 : if(associated( grid%dtxe5 )) deallocate(grid%dtxe5)
1015 0 : if(associated( grid%dyce )) deallocate(grid%dyce)
1016 0 : if(associated( grid%dx )) deallocate(grid%dx)
1017 0 : if(associated( grid%rdx )) deallocate(grid%rdx)
1018 0 : if(associated( grid%cy )) deallocate(grid%cy)
1019 :
1020 0 : if(associated( grid%sc )) deallocate(grid%sc)
1021 0 : if(associated( grid%se )) deallocate(grid%se)
1022 0 : if(associated( grid%dc )) deallocate(grid%dc)
1023 0 : if(associated( grid%de )) deallocate(grid%de)
1024 :
1025 0 : if(associated( grid%cdx )) deallocate(grid%cdx)
1026 0 : if(associated( grid%cdy )) deallocate(grid%cdy)
1027 0 : if(associated( grid%cdx4 )) deallocate(grid%cdx4)
1028 0 : if(associated( grid%cdy4 )) deallocate(grid%cdy4)
1029 0 : if(associated( grid%cdxde )) deallocate(grid%cdxde)
1030 0 : if(associated( grid%cdxdp )) deallocate(grid%cdxdp)
1031 0 : if(associated( grid%cdydp )) deallocate(grid%cdydp)
1032 0 : if(associated( grid%cdyde )) deallocate(grid%cdyde)
1033 0 : if(associated( grid%cdxdiv)) deallocate(grid%cdxdiv)
1034 0 : if(associated( grid%cdydiv)) deallocate(grid%cdydiv)
1035 0 : if(associated( grid%cdtau4)) deallocate(grid%cdtau4)
1036 :
1037 0 : if(associated( grid%scdiv4)) deallocate(grid%scdiv4)
1038 0 : if(associated( grid%sediv4)) deallocate(grid%sediv4)
1039 0 : if(associated( grid%dcdiv4)) deallocate(grid%dcdiv4)
1040 0 : if(associated( grid%dediv4)) deallocate(grid%dediv4)
1041 :
1042 0 : if(associated( grid%f0 )) deallocate(grid%f0)
1043 0 : if(associated( grid%fc )) deallocate(grid%fc)
1044 :
1045 : #if defined(SPMD)
1046 0 : call spmd_vars_clean(grid)
1047 : #endif
1048 :
1049 0 : end subroutine dynamics_clean
1050 :
1051 : !========================================================================================
1052 :
1053 : #if defined(SPMD)
1054 0 : subroutine spmd_vars_clean(grid)
1055 :
1056 : use parutilitiesmodule, only : parpatternfree
1057 :
1058 : ! args
1059 : type (T_FVDYCORE_GRID), intent(inout) :: grid
1060 : !-----------------------------------------------------------------------
1061 :
1062 0 : if ( grid%twod_decomp == 1 ) then
1063 :
1064 : ! Clean transposes
1065 :
1066 0 : call parpatternfree(grid%commxy, grid%u_to_uxy)
1067 0 : call parpatternfree(grid%commxy, grid%uxy_to_u)
1068 0 : call parpatternfree(grid%commxy, grid%v_to_vxy)
1069 0 : call parpatternfree(grid%commxy, grid%vxy_to_v)
1070 0 : call parpatternfree(grid%commxy, grid%ijk_yz_to_xy)
1071 0 : call parpatternfree(grid%commxy, grid%ijk_xy_to_yz)
1072 0 : call parpatternfree(grid%commxy, grid%ikj_xy_to_yz)
1073 0 : call parpatternfree(grid%commxy, grid%ikj_yz_to_xy)
1074 0 : call parpatternfree(grid%commxy, grid%pe_to_pexy)
1075 0 : call parpatternfree(grid%commxy, grid%pexy_to_pe)
1076 0 : call parpatternfree(grid%commxy, grid%pt_to_ptxy)
1077 0 : call parpatternfree(grid%commxy, grid%ptxy_to_pt)
1078 0 : call parpatternfree(grid%commxy, grid%r4_xy_to_yz)
1079 0 : call parpatternfree(grid%commxy, grid%r4_yz_to_xy)
1080 0 : call parpatternfree(grid%commxy, grid%pkxy_to_pkc)
1081 0 : call parpatternfree(grid%commxy, grid%xy2d_to_yz2d)
1082 0 : call parpatternfree(grid%commxy, grid%yz2d_to_xy2d)
1083 : endif
1084 :
1085 0 : end subroutine spmd_vars_clean
1086 : #endif
1087 :
1088 : !========================================================================================
1089 :
1090 0 : end module dynamics_vars
1091 :
|