PartMC  2.3.0
aero_particle_array.F90
Go to the documentation of this file.
1 ! Copyright (C) 2005-2011 Nicole Riemer and 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_aero_particle_array module.
7 
8 !> The aero_particle_array_t structure and assoicated subroutines.
10 
12  use pmc_util
13  use pmc_spec_file
14  use pmc_mpi
15 #ifdef PMC_USE_MPI
16  use mpi
17 #endif
18 
19  !> 1-D arrays of particles, used by aero_state to build a ragged
20  !> array.
21  !!
22  !! One aero_particle_array is generally a list of particles in a
23  !! single size bin, but the basic type can be used for any list of
24  !! particles.
25  !!
26  !! To give a reasonable tradeoff between frequent re-allocs and
27  !! memory usage, the length of an aero_particle_array is generally a
28  !! bit longer than the number of particles stored in it. When the
29  !! array is full then a larger array is allocated with new extra
30  !! space. As a balance between memory usage and frequency of
31  !! re-allocs the length of the array is currently doubled when
32  !! necessary and halved when possible.
33  !!
34  !! The true allocated length of the aero_particle_array can be
35  !! obtained by size(aero_particle_array%%particle), while the number
36  !! of used particle slots in it is given by
37  !! aero_particle_array%%n_part. It must be that
38  !! aero_particle_array%%n_part is less than or equal to
39  !! size(aero_particle_array%%particle).
41  !> Number of particles.
42  integer :: n_part
43  !> Particle array.
44  type(aero_particle_t), pointer :: particle(:)
45  end type aero_particle_array_t
46 
47 contains
48 
49 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
50 
51  !> Allocates and initializes.
52  subroutine aero_particle_array_allocate(aero_particle_array)
53 
54  !> Result.
55  type(aero_particle_array_t), intent(out) :: aero_particle_array
56 
57  aero_particle_array%n_part = 0
58  allocate(aero_particle_array%particle(0))
59 
60  end subroutine aero_particle_array_allocate
61 
62 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
63 
64  !> Allocates and initializes to the given size.
65  subroutine aero_particle_array_allocate_size(aero_particle_array, n_part)
66 
67  !> Result.
68  type(aero_particle_array_t), intent(out) :: aero_particle_array
69  !> Number of particles.
70  integer, intent(in) :: n_part
71 
72  integer :: i
73 
74  aero_particle_array%n_part = n_part
75  allocate(aero_particle_array%particle(n_part))
76  do i = 1,n_part
77  call aero_particle_allocate(aero_particle_array%particle(i))
78  end do
79 
81 
82 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
83 
84  !> Deallocates.
85  subroutine aero_particle_array_deallocate(aero_particle_array)
86 
87  !> Structure to deallocate.
88  type(aero_particle_array_t), intent(inout) :: aero_particle_array
89 
90  integer :: i
91 
92  do i = 1,aero_particle_array%n_part
93  call aero_particle_deallocate(aero_particle_array%particle(i))
94  end do
95  deallocate(aero_particle_array%particle)
96 
97  end subroutine aero_particle_array_deallocate
98 
99 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
100 
101  !> Copies aero_particle_array_from to aero_particle_array_to, both
102  !> of which must already be allocated.
103  subroutine aero_particle_array_copy(aero_particle_array_from, &
104  aero_particle_array_to)
105 
106  !> Origin structure.
107  type(aero_particle_array_t), intent(in) :: aero_particle_array_from
108  !> Destination structure.
109  type(aero_particle_array_t), intent(inout) :: aero_particle_array_to
110 
111  integer :: i
112 
113  call aero_particle_array_deallocate(aero_particle_array_to)
114  call aero_particle_array_allocate_size(aero_particle_array_to, &
115  aero_particle_array_from%n_part)
116  do i = 1,aero_particle_array_from%n_part
117  call aero_particle_copy(aero_particle_array_from%particle(i), &
118  aero_particle_array_to%particle(i))
119  end do
120 
121  end subroutine aero_particle_array_copy
122 
123 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
124 
125  !> Resets an aero_particle_array to contain zero particles.
126  subroutine aero_particle_array_zero(aero_particle_array)
127 
128  !> Structure to reset.
129  type(aero_particle_array_t), intent(inout) :: aero_particle_array
130 
131  call aero_particle_array_deallocate(aero_particle_array)
132  allocate(aero_particle_array%particle(0))
133  aero_particle_array%n_part = 0
134 
135  end subroutine aero_particle_array_zero
136 
137 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
138 
139  !> Changes the given aero_particle_array to exactly the given
140  !> new_length.
141  !!
142  !! This function should not be called directly, but rather use
143  !! aero_particle_array_enlarge(), aero_particle_array_enlarge_to()
144  !! or aero_particle_array_shrink().
145  subroutine aero_particle_array_realloc(aero_particle_array, new_length)
146 
147  !> Array to reallocate (must already be allocated on entry).
148  type(aero_particle_array_t), intent(inout) :: aero_particle_array
149  !> New length of the array.
150  integer, intent(in) :: new_length
151 
152  integer :: n_part, i
153  type(aero_particle_t), pointer :: new_particles(:)
154 
155  n_part = aero_particle_array%n_part
156  call assert(867444847, new_length >= n_part)
157  allocate(new_particles(new_length))
158  do i = 1,aero_particle_array%n_part
159  call aero_particle_shift(aero_particle_array%particle(i), &
160  new_particles(i))
161  end do
162  deallocate(aero_particle_array%particle)
163  aero_particle_array%particle => new_particles
164 
165  end subroutine aero_particle_array_realloc
166 
167 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
168 
169  !> Enlarges the given aero_particle_array by at least one element
170  !!
171  !! Currently this doubles the length.
172  subroutine aero_particle_array_enlarge(aero_particle_array)
173 
174  !> Array to enlarge.
175  type(aero_particle_array_t), intent(inout) :: aero_particle_array
176 
177  integer :: length, new_length
178 
179  length = size(aero_particle_array%particle)
180  new_length = max(length * 2, length + 1)
181  call aero_particle_array_realloc(aero_particle_array, new_length)
182 
183  end subroutine aero_particle_array_enlarge
184 
185 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
186 
187  !> Enlarges the given array so that it is at least of size n.
188  subroutine aero_particle_array_enlarge_to(aero_particle_array, n)
189 
190  !> Array to enlarge.
191  type(aero_particle_array_t), intent(inout) :: aero_particle_array
192  !> Minimum new size of array.
193  integer, intent(in) :: n
194 
195  do while (size(aero_particle_array%particle) < n)
196  call aero_particle_array_enlarge(aero_particle_array)
197  end do
198 
199  end subroutine aero_particle_array_enlarge_to
200 
201 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
202 
203  !> Possibly shrinks the storage of the given array, ensuring that
204  !> it can still store the allocated particles.
205  subroutine aero_particle_array_shrink(aero_particle_array)
206 
207  !> Array to shrink.
208  type(aero_particle_array_t), intent(inout) :: aero_particle_array
209 
210  integer :: n_part, length, new_length
211 
212  n_part = aero_particle_array%n_part
213  length = size(aero_particle_array%particle)
214  new_length = length / 2
215  do while ((n_part <= new_length) .and. (length > 0))
216  call aero_particle_array_realloc(aero_particle_array, new_length)
217  length = size(aero_particle_array%particle)
218  new_length = length / 2
219  end do
220 
221  end subroutine aero_particle_array_shrink
222 
223 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
224 
225  !> Adds the given particle to the end of the array.
226  subroutine aero_particle_array_add_particle(aero_particle_array, &
227  aero_particle)
228 
229  !> Array to add to.
230  type(aero_particle_array_t), intent(inout) :: aero_particle_array
231  !> Particle to add.
232  type(aero_particle_t), intent(in) :: aero_particle
233 
234  integer :: n
235 
236  n = aero_particle_array%n_part + 1
237  call aero_particle_array_enlarge_to(aero_particle_array, n)
238  call aero_particle_allocate(aero_particle_array%particle(n))
239  call aero_particle_copy(aero_particle, &
240  aero_particle_array%particle(n))
241  aero_particle_array%n_part = aero_particle_array%n_part + 1
242 
243  end subroutine aero_particle_array_add_particle
244 
245 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
246 
247  !> Removes the particle at the given index.
248  subroutine aero_particle_array_remove_particle(aero_particle_array, &
249  index)
250 
251  !> Array to remove from.
252  type(aero_particle_array_t), intent(inout) :: aero_particle_array
253  !> Index of particle to remove.
254  integer, intent(in) :: index
255 
256  call assert(992946227, index >= 1)
257  call assert(711246139, index <= aero_particle_array%n_part)
258  call aero_particle_deallocate(aero_particle_array%particle(index))
259  if (index < aero_particle_array%n_part) then
260  ! shift last particle into empty slot to preserve dense packing
261  call aero_particle_shift( &
262  aero_particle_array%particle(aero_particle_array%n_part), &
263  aero_particle_array%particle(index))
264  end if
265  aero_particle_array%n_part = aero_particle_array%n_part - 1
266  call aero_particle_array_shrink(aero_particle_array)
267 
269 
270 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
271 
272  !> Determines the number of bytes required to pack the given value.
273  integer function pmc_mpi_pack_size_apa(val)
274 
275  !> Value to pack.
276  type(aero_particle_array_t), intent(in) :: val
277 
278  integer :: i, total_size
279 
280  total_size = 0
281  total_size = total_size + pmc_mpi_pack_size_integer(val%n_part)
282  do i = 1,val%n_part
283  total_size = total_size &
284  + pmc_mpi_pack_size_aero_particle(val%particle(i))
285  end do
286  pmc_mpi_pack_size_apa = total_size
287 
288  end function pmc_mpi_pack_size_apa
289 
290 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
291 
292  !> Packs the given value into the buffer, advancing position.
293  subroutine pmc_mpi_pack_aero_particle_array(buffer, position, val)
294 
295  !> Memory buffer.
296  character, intent(inout) :: buffer(:)
297  !> Current buffer position.
298  integer, intent(inout) :: position
299  !> Value to pack.
300  type(aero_particle_array_t), intent(in) :: val
301 
302 #ifdef PMC_USE_MPI
303  integer :: prev_position, i
304 
305  prev_position = position
306  call pmc_mpi_pack_integer(buffer, position, val%n_part)
307  do i = 1,val%n_part
308  call pmc_mpi_pack_aero_particle(buffer, position, val%particle(i))
309  end do
310  call assert(803856329, &
311  position - prev_position <= pmc_mpi_pack_size_apa(val))
312 #endif
313 
314  end subroutine pmc_mpi_pack_aero_particle_array
315 
316 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
317 
318  !> Unpacks the given value from the buffer, advancing position.
319  subroutine pmc_mpi_unpack_aero_particle_array(buffer, position, val)
320 
321  !> Memory buffer.
322  character, intent(inout) :: buffer(:)
323  !> Current buffer position.
324  integer, intent(inout) :: position
325  !> Value to pack.
326  type(aero_particle_array_t), intent(inout) :: val
327 
328 #ifdef PMC_USE_MPI
329  integer :: prev_position, i
330 
332  prev_position = position
333  call pmc_mpi_unpack_integer(buffer, position, val%n_part)
334  allocate(val%particle(val%n_part))
335  do i = 1,val%n_part
336  call aero_particle_allocate(val%particle(i))
337  call pmc_mpi_unpack_aero_particle(buffer, position, val%particle(i))
338  end do
339  call assert(138783294, &
340  position - prev_position <= pmc_mpi_pack_size_apa(val))
341 #endif
342 
344 
345 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
346 
347 end module pmc_aero_particle_array
subroutine pmc_mpi_pack_aero_particle(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_particle_array_copy(aero_particle_array_from, aero_particle_array_to)
Copies aero_particle_array_from to aero_particle_array_to, both of which must already be allocated...
subroutine aero_particle_array_deallocate(aero_particle_array)
Deallocates.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
Definition: mpi.F90:347
1-D arrays of particles, used by aero_state to build a ragged array.
subroutine pmc_mpi_unpack_aero_particle(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine aero_particle_array_zero(aero_particle_array)
Resets an aero_particle_array to contain zero particles.
The aero_particle_t structure and associated subroutines.
subroutine aero_particle_deallocate(aero_particle)
Deallocates memory associated with an aero_particle_t.
subroutine aero_particle_array_allocate_size(aero_particle_array, n_part)
Allocates and initializes to the given size.
subroutine aero_particle_shift(aero_particle_from, aero_particle_to)
Shift data from one aero_particle_t to another and free the first one.
subroutine aero_particle_array_realloc(aero_particle_array, new_length)
Changes the given aero_particle_array to exactly the given new_length.
The aero_particle_array_t structure and assoicated subroutines.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
Definition: mpi.F90:534
subroutine aero_particle_array_add_particle(aero_particle_array, aero_particle)
Adds the given particle to the end of the array.
Common utility subroutines.
Definition: util.F90:9
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Definition: mpi.F90:771
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine aero_particle_array_shrink(aero_particle_array)
Possibly shrinks the storage of the given array, ensuring that it can still store the allocated parti...
subroutine aero_particle_allocate(aero_particle)
Allocates memory in an aero_particle_t.
Single aerosol particle data structure.
Reading formatted text input.
Definition: spec_file.F90:43
subroutine aero_particle_array_remove_particle(aero_particle_array, index)
Removes the particle at the given index.
subroutine aero_particle_array_enlarge_to(aero_particle_array, n)
Enlarges the given array so that it is at least of size n.
subroutine pmc_mpi_unpack_aero_particle_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine aero_particle_array_allocate(aero_particle_array)
Allocates and initializes.
subroutine aero_particle_array_enlarge(aero_particle_array)
Enlarges the given aero_particle_array by at least one element.
integer function pmc_mpi_pack_size_aero_particle(val)
Determines the number of bytes required to pack the given value.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
Definition: util.F90:102
subroutine pmc_mpi_pack_aero_particle_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine aero_particle_copy(aero_particle_from, aero_particle_to)
Copies a particle.
integer function pmc_mpi_pack_size_apa(val)
Determines the number of bytes required to pack the given value.