PartMC  2.3.0
integer_rmap.F90
Go to the documentation of this file.
1 ! Copyright (C) 2011-2012 Matthew West
2 ! Licensed under the GNU General Public License version 2 or (at your
3 ! option) any later version. See the file COPYING for details.
4 
5 !> \file
6 !> The pmc_integer_rmap module.
7 
8 !> The integer_rmap_t structure and assocated subroutines.
10 
12  use pmc_util
13  use pmc_mpi
14 
15  !> A map from integers to integers, together with its multi-valued
16  !> inverse.
17  !!
18  !! The forward map takes integer \f$i\f$ in the domain
19  !! 1,...,n_domain to an integer \f$j\f$ in the range
20  !! 1,...,n_range. This is stored with <tt>j =
21  !! integer_rmap%%forward%%entry(i)</tt>. This map will generally not be
22  !! one-to-one or onto.
23  !!
24  !! The inverse map is multi-valued, with
25  !! <tt>integer_rmap%%inverse(j)</tt> containing all the inverses of
26  !! \f$j\f$. The entries in the inverse map are given by
27  !! <tt>inverse_rmap%%index</tt>. The relationships between
28  !! the forward and reverse maps are as follows.
29  !!
30  !! Given \f$i\f$, let:
31  !! <pre>
32  !! j = integer_rmap%%forward%%entry(i)
33  !! k = integer_rmap%%index%%entry(i)
34  !! </pre>
35  !! Then:
36  !! <pre>
37  !! integer_rmap%%inverse(j)%%entry(k) == i
38  !! </pre>
39  !!
40  !! Alternatively, given \f$j\f$ and \f$k\f$, let:
41  !! <pre>
42  !! i = integer_rmap%%inverse(j)%%entry(k)
43  !! </pre>
44  !! Then:
45  !! <pre>
46  !! integer_rmap%%forward%%entry(i) == j
47  !! integer_rmap%%index%%entry(i) == k
48  !! </pre>
50  !> Forward map (single valued).
51  type(integer_varray_t) :: forward
52  !> Inverse map (multi-valued).
53  type(integer_varray_t), allocatable :: inverse(:)
54  !> Forward map to inverse map entries (single valued).
55  type(integer_varray_t) :: index
56  end type integer_rmap_t
57 
58 contains
59 
60 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
61 
62  !> Allocates an empty structure.
63  elemental subroutine integer_rmap_allocate(integer_rmap)
64 
65  !> Structure to initialize.
66  type(integer_rmap_t), intent(out) :: integer_rmap
67 
68  call integer_varray_allocate(integer_rmap%forward)
69  allocate(integer_rmap%inverse(0))
70  call integer_varray_allocate(integer_rmap%index)
71 
72  end subroutine integer_rmap_allocate
73 
74 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
75 
76  !> Allocates a structure with the given size.
77  elemental subroutine integer_rmap_allocate_size(integer_rmap, n_range)
78 
79  !> Structure to initialize.
80  type(integer_rmap_t), intent(out) :: integer_rmap
81  !> Size of range space.
82  integer, intent(in) :: n_range
83 
84  call integer_varray_allocate(integer_rmap%forward)
85  allocate(integer_rmap%inverse(n_range))
86  call integer_varray_allocate(integer_rmap%inverse)
87  call integer_varray_allocate(integer_rmap%index)
88 
89  end subroutine integer_rmap_allocate_size
90 
91 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
92 
93  !> Deallocates a previously allocated structure.
94  elemental subroutine integer_rmap_deallocate(integer_rmap)
95 
96  !> Structure to deallocate.
97  type(integer_rmap_t), intent(inout) :: integer_rmap
98 
99  call integer_varray_deallocate(integer_rmap%forward)
100  call integer_varray_deallocate(integer_rmap%inverse)
101  deallocate(integer_rmap%inverse)
102  call integer_varray_deallocate(integer_rmap%index)
103 
104  end subroutine integer_rmap_deallocate
105 
106 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
107 
108  !> Resets an integer_rmap to have zero particles per bin.
109  elemental subroutine integer_rmap_zero(integer_rmap)
110 
111  !> Structure to zero.
112  type(integer_rmap_t), intent(inout) :: integer_rmap
113 
114  call integer_varray_zero(integer_rmap%forward)
115  call integer_varray_zero(integer_rmap%inverse)
116  call integer_varray_zero(integer_rmap%index)
117 
118  end subroutine integer_rmap_zero
119 
120 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
121 
122  !> Set the map value of the next free domain value to \c i_range.
123  subroutine integer_rmap_append(integer_rmap, i_range)
124 
125  !> Map to append to.
126  type(integer_rmap_t), intent(inout) :: integer_rmap
127  !> Range value.
128  integer, intent(in) :: i_range
129 
130  call assert(549740445, i_range >= 1)
131  call assert(145872613, i_range <= size(integer_rmap%inverse))
132 
133  ! grow map by one element
134  call integer_varray_append(integer_rmap%forward, i_range)
135  call integer_varray_append(integer_rmap%inverse(i_range), &
136  integer_rmap%forward%n_entry)
137  call integer_varray_append(integer_rmap%index, &
138  integer_rmap%inverse(i_range)%n_entry)
139 
140  end subroutine integer_rmap_append
141 
142 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
143 
144  !> Change the map value of \c i_domain to \c i_range.
145  subroutine integer_rmap_change(integer_rmap, i_domain, i_range)
146 
147  !> Map to change.
148  type(integer_rmap_t), intent(inout) :: integer_rmap
149  !> Domain value.
150  integer, intent(in) :: i_domain
151  !> Range value.
152  integer, intent(in) :: i_range
153 
154  integer :: i_range_old, i_index_old, i_domain_shifted
155 
156  call assert(709581778, i_domain >= 1)
157  call assert(494892311, i_domain <= integer_rmap%forward%n_entry)
158 
159  call assert(590911054, i_range >= 1)
160  call assert(859774512, i_range <= size(integer_rmap%inverse))
161 
162  i_range_old = integer_rmap%forward%entry(i_domain)
163  if (i_range_old == i_range) return
164  i_index_old = integer_rmap%index%entry(i_domain)
165 
166  ! remove the old inverse map
167  call integer_varray_remove_entry(integer_rmap%inverse(i_range_old), &
168  i_index_old)
169  if (i_index_old <= integer_rmap%inverse(i_range_old)%n_entry) then
170  ! the removed entry wasn't the last one, so the last entry
171  ! was moved and needs fixing
172  i_domain_shifted = integer_rmap%inverse(i_range_old)%entry(i_index_old)
173  integer_rmap%index%entry(i_domain_shifted) = i_index_old
174  end if
175 
176  ! set the new map and inverse
177  integer_rmap%forward%entry(i_domain) = i_range
178  call integer_varray_append(integer_rmap%inverse(i_range), i_domain)
179  integer_rmap%index%entry(i_domain) = integer_rmap%inverse(i_range)%n_entry
180 
181  end subroutine integer_rmap_change
182 
183 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
184 
185  !> Replace the map at the given \c i_domain with the map value of
186  !> the last entry, and delete the last entry.
187  subroutine integer_rmap_remove(integer_rmap, i_domain)
188 
189  !> Map to remove from.
190  type(integer_rmap_t), intent(inout) :: integer_rmap
191  !> Domain value to replace.
192  integer, intent(in) :: i_domain
193 
194  integer :: i_range_old, i_index_old, i_domain_shifted, i_range_fix
195  integer :: i_index_fix, i_domain_fix
196 
197  call assert(745161821, i_domain >= 1)
198  call assert(143043782, i_domain <= integer_rmap%forward%n_entry)
199 
200  ! Deleting particles shifts the end particles into the empty slots
201  ! in the aero_particle_array and the aero_sorted forward and
202  ! reverse indexes. All must be fixed in the right order to
203  ! maintain consistency.
204 
205  i_range_old = integer_rmap%forward%entry(i_domain)
206  i_index_old = integer_rmap%index%entry(i_domain)
207 
208  i_domain_shifted = integer_rmap%forward%n_entry ! old loc of shifted value
209  if (i_domain_shifted /= i_domain) then
210  i_range_fix = integer_rmap%forward%entry(i_domain_shifted)
211  i_index_fix = integer_rmap%index%entry(i_domain_shifted)
212  integer_rmap%inverse(i_range_fix)%entry(i_index_fix) = i_domain
213  end if
214 
215  ! remove the particle from the forward map (with the side effect
216  ! of fixing the forward map for the shifted value)
217  call integer_varray_remove_entry(integer_rmap%forward, i_domain)
218  call integer_varray_remove_entry(integer_rmap%index, i_domain)
219 
220  ! remove the inverse map
221  i_index_fix = integer_rmap%inverse(i_range_old)%n_entry
222  i_domain_fix = integer_rmap%inverse(i_range_old)%entry(i_index_fix)
223  call integer_varray_remove_entry(integer_rmap%inverse(i_range_old), &
224  i_index_old)
225 
226  if (i_index_fix /= i_index_old) then
227  ! fix index map
228  integer_rmap%index%entry(i_domain_fix) = i_index_old
229  end if
230 
231  end subroutine integer_rmap_remove
232 
233 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
234 
235  !> Check that the data is consistent.
236  subroutine integer_rmap_check(integer_rmap, name, n_domain, n_range, &
237  continue_on_error)
238 
239  !> Structure to check.
240  type(integer_rmap_t) :: integer_rmap
241  !> Check name.
242  character(len=*), intent(in) :: name
243  !> Number of domain items.
244  integer, intent(in) :: n_domain
245  !> Number of image items.
246  integer, intent(in) :: n_range
247  !> Whether to continue despite error.
248  logical, intent(in) :: continue_on_error
249 
250  integer :: i_domain, i_range, i_index
251 
252  if ((n_domain /= integer_rmap%forward%n_entry) &
253  .or. (n_domain /= integer_rmap%index%n_entry) &
254  .or. (n_range /= size(integer_rmap%inverse))) then
255  write(0,*) 'ERROR integer_rmap A:', name
256  write(0,*) 'n_domain', n_domain
257  write(0,*) 'n_range', n_range
258  write(0,*) 'integer_rmap%forward%n_entry', integer_rmap%forward%n_entry
259  write(0,*) 'integer_rmap%index%n_entry', integer_rmap%index%n_entry
260  write(0,*) 'size(integer_rmap%inverse)', size(integer_rmap%inverse)
261  call assert(973643016, continue_on_error)
262  end if
263 
264  do i_domain = 1,n_domain
265  i_range = integer_rmap%forward%entry(i_domain)
266  if ((i_range < 1) .or. (i_range > n_range)) then
267  write(0,*) 'ERROR integer_rmap B:', name
268  write(0,*) 'i_domain', i_domain
269  write(0,*) 'i_range', i_range
270  write(0,*) 'n_range', n_range
271  call assert(798857945, continue_on_error)
272  end if
273 
274  i_index = integer_rmap%index%entry(i_domain)
275  if ((i_index < 1) &
276  .or. (i_index > integer_rmap%inverse(i_range)%n_entry)) then
277  write(0,*) 'ERROR integer_rmap C:', name
278  write(0,*) 'i_domain', i_domain
279  write(0,*) 'i_range', i_range
280  write(0,*) 'i_index', i_index
281  write(0,*) 'integer_rmap%inverse(i_range)%n_entry', &
282  integer_rmap%inverse(i_range)%n_entry
283  call assert(823748734, continue_on_error)
284  end if
285  if (i_domain /= integer_rmap%inverse(i_range)%entry(i_index)) then
286  write(0,*) 'ERROR integer_rmap D:', name
287  write(0,*) 'i_domain', i_domain
288  write(0,*) 'i_range', i_range
289  write(0,*) 'i_index', i_index
290  write(0,*) 'integer_rmap%inverse(i_range)%entry(i_index)', &
291  integer_rmap%inverse(i_range)%entry(i_index)
292  call assert(735205557, continue_on_error)
293  end if
294  end do
295 
296  do i_range = 1,n_range
297  do i_index = 1,integer_rmap%inverse(i_range)%n_entry
298  i_domain = integer_rmap%inverse(i_range)%entry(i_index)
299  if ((i_domain < 1) .or. (i_domain > n_domain)) then
300  write(0,*) 'ERROR integer_rmap E:', name
301  write(0,*) 'i_range', i_range
302  write(0,*) 'i_index', i_index
303  write(0,*) 'i_domain', i_domain
304  write(0,*) 'n_domain', n_domain
305  call assert(502643520, continue_on_error)
306  end if
307  if ((i_range /= integer_rmap%forward%entry(i_domain)) &
308  .or. (i_index /= integer_rmap%index%entry(i_domain))) then
309  write(0,*) 'ERROR integer_rmap F:', name
310  write(0,*) 'i_domain', i_domain
311  write(0,*) 'i_range', i_range
312  write(0,*) 'integer_rmap%forward%entry(i_domain)', &
313  integer_rmap%forward%entry(i_domain)
314  write(0,*) 'i_index', i_index
315  write(0,*) 'integer_rmap%index%entry(i_domain)', &
316  integer_rmap%index%entry(i_domain)
317  call assert(544747928, continue_on_error)
318  end if
319  end do
320  end do
321 
322  end subroutine integer_rmap_check
323 
324 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
325 
326  !> Determines the number of bytes required to pack the given value.
327  integer function pmc_mpi_pack_size_integer_rmap(val)
328 
329  !> Value to pack.
330  type(integer_rmap_t), intent(in) :: val
331 
332  integer :: i, total_size
333 
334  total_size = 0
335  total_size = total_size + pmc_mpi_pack_size_integer(size(val%inverse))
336  total_size = total_size + pmc_mpi_pack_size_integer_varray(val%forward)
337  do i = 1,size(val%inverse)
338  total_size = total_size &
339  + pmc_mpi_pack_size_integer_varray(val%inverse(i))
340  end do
341  total_size = total_size + pmc_mpi_pack_size_integer_varray(val%index)
342  pmc_mpi_pack_size_integer_rmap = total_size
343 
344  end function pmc_mpi_pack_size_integer_rmap
345 
346 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
347 
348  !> Packs the given value into the buffer, advancing position.
349  subroutine pmc_mpi_pack_integer_rmap(buffer, position, val)
350 
351  !> Memory buffer.
352  character, intent(inout) :: buffer(:)
353  !> Current buffer position.
354  integer, intent(inout) :: position
355  !> Value to pack.
356  type(integer_rmap_t), intent(in) :: val
357 
358 #ifdef PMC_USE_MPI
359  integer :: prev_position, i
360 
361  prev_position = position
362  call pmc_mpi_pack_integer(buffer, position, size(val%inverse))
363  call pmc_mpi_pack_integer_varray(buffer, position, val%forward)
364  do i = 1,size(val%inverse)
365  call pmc_mpi_pack_integer_varray(buffer, position, val%inverse(i))
366  end do
367  call pmc_mpi_pack_integer_varray(buffer, position, val%index)
368  call assert(533568488, &
369  position - prev_position <= pmc_mpi_pack_size_integer_rmap(val))
370 #endif
371 
372  end subroutine pmc_mpi_pack_integer_rmap
373 
374 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
375 
376  !> Unpacks the given value from the buffer, advancing position.
377  subroutine pmc_mpi_unpack_integer_rmap(buffer, position, val)
378 
379  !> Memory buffer.
380  character, intent(inout) :: buffer(:)
381  !> Current buffer position.
382  integer, intent(inout) :: position
383  !> Value to pack.
384  type(integer_rmap_t), intent(inout) :: val
385 
386 #ifdef PMC_USE_MPI
387  integer :: prev_position, i, n
388 
389  prev_position = position
390  call pmc_mpi_unpack_integer(buffer, position, n)
391  call integer_rmap_deallocate(val)
392  call integer_rmap_allocate_size(val, n)
393  call pmc_mpi_unpack_integer_varray(buffer, position, val%forward)
394  do i = 1,size(val%inverse)
395  call pmc_mpi_unpack_integer_varray(buffer, position, val%inverse(i))
396  end do
397  call pmc_mpi_unpack_integer_varray(buffer, position, val%index)
398  call assert(663161025, &
399  position - prev_position <= pmc_mpi_pack_size_integer_rmap(val))
400 #endif
401 
402  end subroutine pmc_mpi_unpack_integer_rmap
403 
404 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
405 
406 end module pmc_integer_rmap
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:347
elemental subroutine integer_rmap_allocate_size(integer_rmap, n_range)
Allocates a structure with the given size.
The integer_varray_t structure and assocated subroutines.
integer function pmc_mpi_pack_size_integer_rmap(val)
Determines the number of bytes required to pack the given value.
subroutine integer_rmap_check(integer_rmap, name, n_domain, n_range, continue_on_error)
Check that the data is consistent.
subroutine integer_varray_append(integer_varray, val)
Adds the given number to the end of the array.
subroutine integer_varray_remove_entry(integer_varray, index)
Removes the entry at the given index, repacking values to maintain contiguous data.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:534
A map from integers to integers, together with its multi-valued inverse.
Common utility subroutines.
Definition: util.F90:9
elemental subroutine integer_varray_deallocate(integer_varray)
Deallocates a previously allocated structure.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:771
elemental subroutine integer_varray_zero(integer_varray)
Resets an integer_varray to have zero particles per bin.
Wrapper functions for MPI.
Definition: mpi.F90:13
elemental subroutine integer_rmap_deallocate(integer_rmap)
Deallocates a previously allocated structure.
elemental subroutine integer_varray_allocate(integer_varray)
Allocates an empty structure.
subroutine pmc_mpi_unpack_integer_rmap(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine integer_rmap_remove(integer_rmap, i_domain)
Replace the map at the given i_domain with the map value of the last entry, and delete the last entry...
subroutine pmc_mpi_pack_integer_varray(buffer, position, val)
Packs the given value into the buffer, advancing position.
A variable-length 1D array of integers.
The integer_rmap_t structure and assocated subroutines.
Definition: integer_rmap.F90:9
subroutine pmc_mpi_pack_integer_rmap(buffer, position, val)
Packs the given value into the buffer, advancing position.
elemental subroutine integer_rmap_zero(integer_rmap)
Resets an integer_rmap to have zero particles per bin.
subroutine pmc_mpi_unpack_integer_varray(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function pmc_mpi_pack_size_integer_varray(val)
Determines the number of bytes required to pack the given value.
elemental subroutine integer_rmap_allocate(integer_rmap)
Allocates an empty structure.
subroutine integer_rmap_change(integer_rmap, i_domain, i_range)
Change the map value of i_domain to i_range.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:102
subroutine integer_rmap_append(integer_rmap, i_range)
Set the map value of the next free domain value to i_range.