da_cv_to_global.inc
References to this file elsewhere.
1 subroutine da_cv_to_global(cv_size, cv_size_global, x, xp, mzs, xg)
2
3
4 !-----------------------------------------------------------------------
5 ! Purpose: Gathers local cv-array x into domain cv-array xg(:).
6 ! Global cv-array xg will only be valid on the "monitor" task.
7 !
8 ! Must be called by all MPI tasks.
9 !-----------------------------------------------------------------------
10
11 implicit none
12
13 integer, intent(in) :: cv_size ! Size of local cv-array
14 integer, intent(in) :: cv_size_global ! Size of domain cv-array
15 real, intent(in) :: x(1:cv_size) ! local cv-array
16 type (xpose_type), intent(in) :: xp ! decomposed dimensions
17 integer, intent(in) :: mzs(:) ! mz for each variable
18 ! (to identify empty or 2D arrays)
19 real, intent(inout) :: xg(1:cv_size_global) ! global cv-array
20
21 #ifdef DM_PARALLEL
22
23 ! Local declarations
24 type (vp_type) :: vv_x ! Grdipt/EOF cv_array (local)
25 type (vp_type) :: vv_xg ! Grdipt/EOF cv_array (global)
26 type (xpose_type) :: xpg ! global dimensions
27 integer :: ids, ide, jds, jde, kds, &
28 ims, ime, jms, jme, kms, &
29 ips, ipe, jps, jpe, kps
30
31 !
32 ! Gather to mimic serial summation order.
33 !
34
35 ! k?e varies with variable v1 - v5
36 ids = xp%ids; ide = xp%ide; jds = xp%jds; jde = xp%jde; kds = xp%kds
37 ims = xp%ims; ime = xp%ime; jms = xp%jms; jme = xp%jme; kms = xp%kms
38 ips = xp%ips; ipe = xp%ipe; jps = xp%jps; jpe = xp%jpe; kps = xp%kps
39
40 ! TOdo: encapsulate this crap!
41 ! allocate vv_x
42 allocate(vv_x%v1(xp%ims:xp%ime,xp%jms:xp%jme,mzs(1)))
43 allocate(vv_x%v2(xp%ims:xp%ime,xp%jms:xp%jme,mzs(2)))
44 allocate(vv_x%v3(xp%ims:xp%ime,xp%jms:xp%jme,mzs(3)))
45 allocate(vv_x%v4(xp%ims:xp%ime,xp%jms:xp%jme,mzs(4)))
46 allocate(vv_x%v5(xp%ims:xp%ime,xp%jms:xp%jme,mzs(5)))
47 allocate(vv_x%alpha(xp%ims:xp%ime,xp%jms:xp%jme,mzs(6)))
48
49 call da_cv_to_vv (cv_size, x, xp, mzs, vv_x)
50
51 if (rootproc) then
52 ! allocate vv_xg
53 allocate(vv_xg%v1(xp%ids:xp%ide,xp%jds:xp%jde,mzs(1)))
54 allocate(vv_xg%v2(xp%ids:xp%ide,xp%jds:xp%jde,mzs(2)))
55 allocate(vv_xg%v3(xp%ids:xp%ide,xp%jds:xp%jde,mzs(3)))
56 allocate(vv_xg%v4(xp%ids:xp%ide,xp%jds:xp%jde,mzs(4)))
57 allocate(vv_xg%v5(xp%ids:xp%ide,xp%jds:xp%jde,mzs(5)))
58 allocate(vv_xg%alpha(xp%ids:xp%ide,xp%jds:xp%jde,mzs(6)))
59 else
60 ! Allocate dummy array for non-monitor process to keep Fortran
61 ! CICO happy...
62 allocate(vv_xg%v1(1,1,1))
63 allocate(vv_xg%v2(1,1,1))
64 allocate(vv_xg%v3(1,1,1))
65 allocate(vv_xg%v4(1,1,1))
66 allocate(vv_xg%v5(1,1,1))
67 allocate(vv_xg%alpha(1,1,1))
68 end if
69
70 ! TOdo: encapsulate this crap!
71 ! gather to global data structures
72 ! it is possible to gather straight into a globally-sized cv-array
73 ! TOdo: add this optimization later
74 call da_wv_patch_to_global(vv_x%v1, vv_xg%v1, &
75 xp%domdesc, mzs(1), &
76 ids, ide, jds, jde, kds, &
77 ims, ime, jms, jme, kms, &
78 ips, ipe, jps, jpe, kps)
79 call da_wv_patch_to_global(vv_x%v2, vv_xg%v2, &
80 xp%domdesc, mzs(2), &
81 ids, ide, jds, jde, kds, &
82 ims, ime, jms, jme, kms, &
83 ips, ipe, jps, jpe, kps)
84 call da_wv_patch_to_global(vv_x%v3, vv_xg%v3, &
85 xp%domdesc, mzs(3), &
86 ids, ide, jds, jde, kds, &
87 ims, ime, jms, jme, kms, &
88 ips, ipe, jps, jpe, kps)
89 call da_wv_patch_to_global(vv_x%v4, vv_xg%v4, &
90 xp%domdesc, mzs(4), &
91 ids, ide, jds, jde, kds, &
92 ims, ime, jms, jme, kms, &
93 ips, ipe, jps, jpe, kps)
94 call da_wv_patch_to_global(vv_x%v5, vv_xg%v5, &
95 xp%domdesc, mzs(5), &
96 ids, ide, jds, jde, kds, &
97 ims, ime, jms, jme, kms, &
98 ips, ipe, jps, jpe, kps)
99 call da_wv_patch_to_global(vv_x%alpha, vv_xg%alpha, &
100 xp%domdesc, mzs(6), &
101 ids, ide, jds, jde, kds, &
102 ims, ime, jms, jme, kms, &
103 ips, ipe, jps, jpe, kps)
104
105 ! deallocate vv_x
106 deallocate (vv_x%v1, vv_x%v2, vv_x%v3, vv_x%v4, vv_x%v5)
107
108 if (rootproc) then
109 ! finally, collapse data back into a globally-sized cv-array
110 xpg%ids = xp%ids; xpg%ide = xp%ide
111 xpg%ims = xp%ids; xpg%ime = xp%ide
112 xpg%its = xp%ids; xpg%ite = xp%ide
113 xpg%jds = xp%jds; xpg%jde = xp%jde
114 xpg%jms = xp%jds; xpg%jme = xp%jde
115 xpg%jts = xp%jds; xpg%jte = xp%jde
116 xpg%kds = xp%kds; xpg%kde = xp%kde
117 xpg%kms = xp%kds; xpg%kme = xp%kde
118 xpg%kts = xp%kds; xpg%kte = xp%kde
119 call da_vv_to_cv(vv_xg, xpg, mzs, cv_size_global, xg)
120 end if
121
122 ! deallocate vv_xg
123 deallocate(vv_xg%v1, vv_xg%v2, vv_xg%v3, vv_xg%v4, vv_xg%v5, vv_xg%alpha)
124
125 #endif
126
127 end subroutine da_cv_to_global
128
129