da_cv_to_global.inc
References to this file elsewhere.
1 subroutine da_cv_to_global(cv_size, cv_size_global, x, grid, 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(domain), intent(in) :: grid
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 if (trace_use) call da_trace_entry("da_cv_to_global")
32
33 !
34 ! Gather to mimic serial summation order.
35 !
36
37 ! k?e varies with variable v1 - v5
38 ids = ids; ide = ide; jds = jds; jde = jde; kds = kds
39 ims = ims; ime = ime; jms = jms; jme = jme; kms = grid%xp%kms
40 ips = grid%xp%ips; ipe = grid%xp%ipe; jps = grid%xp%jps; jpe = grid%xp%jpe; kps = grid%xp%kps
41
42 ! TOdo: encapsulate this crap!
43 ! allocate vv_x
44 allocate(vv_x%v1(ims:ime,jms:jme,mzs(1)))
45 allocate(vv_x%v2(ims:ime,jms:jme,mzs(2)))
46 allocate(vv_x%v3(ims:ime,jms:jme,mzs(3)))
47 allocate(vv_x%v4(ims:ime,jms:jme,mzs(4)))
48 allocate(vv_x%v5(ims:ime,jms:jme,mzs(5)))
49 allocate(vv_x%alpha(ims:ime,jms:jme,mzs(6)))
50
51 call da_cv_to_vv (cv_size, x, mzs, vv_x)
52
53 if (rootproc) then
54 ! allocate vv_xg
55 allocate(vv_xg%v1(ids:ide,jds:jde,mzs(1)))
56 allocate(vv_xg%v2(ids:ide,jds:jde,mzs(2)))
57 allocate(vv_xg%v3(ids:ide,jds:jde,mzs(3)))
58 allocate(vv_xg%v4(ids:ide,jds:jde,mzs(4)))
59 allocate(vv_xg%v5(ids:ide,jds:jde,mzs(5)))
60 allocate(vv_xg%alpha(ids:ide,jds:jde,mzs(6)))
61 else
62 ! Allocate dummy array for non-monitor process to keep Fortran
63 ! CICO happy...
64 allocate(vv_xg%v1(1,1,1))
65 allocate(vv_xg%v2(1,1,1))
66 allocate(vv_xg%v3(1,1,1))
67 allocate(vv_xg%v4(1,1,1))
68 allocate(vv_xg%v5(1,1,1))
69 allocate(vv_xg%alpha(1,1,1))
70 end if
71
72 ! TOdo: encapsulate this crap!
73 ! gather to global data structures
74 ! it is possible to gather straight into a globally-sized cv-array
75 ! TOdo: add this optimization later
76 call da_patch_to_global(grid, vv_x%v1, vv_xg%v1, mzs(1))
77 call da_patch_to_global(grid, vv_x%v2, vv_xg%v2, mzs(2))
78 call da_patch_to_global(grid, vv_x%v3, vv_xg%v3, mzs(3))
79 call da_patch_to_global(grid, vv_x%v4, vv_xg%v4, mzs(4))
80 call da_patch_to_global(grid, vv_x%v5, vv_xg%v5, mzs(5))
81 call da_patch_to_global(grid, vv_x%alpha, vv_xg%alpha, mzs(6))
82
83 ! deallocate vv_x
84 deallocate (vv_x%v1, vv_x%v2, vv_x%v3, vv_x%v4, vv_x%v5)
85
86 if (rootproc) then
87 ! finally, collapse data back into a globally-sized cv-array
88 xpg%ids = ids; xpg%ide = ide
89 xpg%ims = ids; xpg%ime = ide
90 xpg%its = ids; xpg%ite = ide
91 xpg%jds = jds; xpg%jde = jde
92 xpg%jms = jds; xpg%jme = jde
93 xpg%jts = jds; xpg%jte = jde
94 xpg%kds = kds; xpg%kde = kde
95 xpg%kms = kds; xpg%kme = kde
96 xpg%kts = kds; xpg%kte = kde
97 call da_vv_to_cv(vv_xg, xpg, mzs, cv_size_global, xg)
98 end if
99
100 ! deallocate vv_xg
101 deallocate(vv_xg%v1, vv_xg%v2, vv_xg%v3, vv_xg%v4, vv_xg%v5, vv_xg%alpha)
102
103 if (trace_use) call da_trace_exit("da_cv_to_global")
104
105 #endif
106
107 end subroutine da_cv_to_global
108
109