ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_module.F90
Go to the documentation of this file.
1!Copyright 2018-2020 Rene Halver, Forschungszentrum Juelich GmbH, Germany
2!Copyright 2018-2020 Godehard Sutmann, Forschungszentrum Juelich GmbH, Germany
3!Copyright 2019-2020 Stephan Schulz, Ruhr-Universität Bochum, Germany
4!
5!Redistribution and use in source and binary forms, with or without modification, are
6!permitted provided that the following conditions are met:
7!
8!1. Redistributions of source code must retain the above copyright notice, this list
9! of conditions and the following disclaimer.
10!
11!2. Redistributions in binary form must reproduce the above copyright notice, this
12! list of conditions and the following disclaimer in the documentation and/or
13! other materials provided with the distribution.
14!
15!3. Neither the name of the copyright holder nor the names of its contributors
16! may be used to endorse or promote products derived from this software without
17! specific prior written permission.
18!
19!THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
20!EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
21!OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
22!SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
23!INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
24!TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
25!BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
26!CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
27!ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
28!SUCH DAMAGE.
29! -DCMAKE_Fortran_FLAGS='-std=f2003 -Wall -Wextra -fbacktrace' -DCM_ALL_FORTRAN=ON -DCMAKE_BUILD_TYPE=Debug -DCM_ALL_DEBUG=ON
30
31!> The function API for ALL
32!!
33!! The ``*_c`` functions should never be called directly, but their non suffixed
34!! counterparts should. Safety checks and conveniences provided by Fortran are
35!! available there. For a simple example on how to use this see ``examples/``
36!! directory.
38 use iso_c_binding
39 use iso_fortran_env
40 implicit none
41 private
42 ! The different loadbalancing methods. Must be kept in sync with C side.
43 integer(c_int), public, parameter :: all_staggered = 0
44 integer(c_int), public, parameter :: all_tensor = 1
45 integer(c_int), public, parameter :: all_forcebased = 2
46#ifdef ALL_VORONOI_ACTIVE
47 integer(c_int), public, parameter :: all_voronoi = 3
48#endif
49 integer(c_int), public, parameter :: all_histogram = 4
50 integer(c_int), public, parameter :: all_tensor_max = 5
51 integer(c_int), public, parameter :: all_unimplemented = 128
52 ! The different error IDs. Must be kept in sync with CustomException class
53 integer(c_int), public, parameter :: all_error_generic = 1
54 integer(c_int), public, parameter :: all_error_pointdimensionmissmatch = 2
55 integer(c_int), public, parameter :: all_error_invalidcommtype = 3
56 integer(c_int), public, parameter :: all_error_invalidargument = 4
57 integer(c_int), public, parameter :: all_error_outofbounds = 5
58 integer(c_int), public, parameter :: all_error_internal = 6
59 integer(c_int), public, parameter :: all_error_filesystem = 6
60 ! Maximum character length of error descriptions
61 integer, public, parameter :: all_error_length = 1024
62 ! Direct interface with C wrapper
63 ! TODO add intent() to dummy arguments
64 interface
65 function all_init_c(method, dim, gamma) result(this) bind(C)
66 use iso_c_binding
67 integer(c_int), value :: method
68 integer(c_int), value :: dim
69 real(c_double), value :: gamma
70 type(c_ptr) :: this
71 end function
72 subroutine all_finalize_c(obj) bind(C)
73 use iso_c_binding
74 type(c_ptr), value :: obj
75 end subroutine
76 subroutine all_set_gamma_c(obj, gamma) bind(C)
77 use iso_c_binding
78 type(c_ptr), value :: obj
79 real(c_double), value :: gamma
80 end subroutine
81 subroutine all_set_proc_grid_params_c(obj, nloc, loc, nsize, size) bind(C)
82 use iso_c_binding
83 type(c_ptr), value :: obj
84 integer(c_int), value :: nloc
85 integer(c_int), dimension(nloc) :: loc
86 integer(c_int), value :: nsize
87 integer(c_int), dimension(nsize) :: size
88 end subroutine
89 subroutine all_set_proc_tag_c(obj, tag) bind(C)
90 use iso_c_binding
91 type(c_ptr), value :: obj
92 integer(c_int), value :: tag
93 end subroutine
94 subroutine all_set_min_domain_size_c(obj, dim, domain_size) bind(C)
95 use iso_c_binding
96 type(c_ptr), value :: obj
97 integer(c_int), value :: dim
98 real(c_double), dimension(dim) :: domain_size
99 end subroutine
100 subroutine all_set_work_c(obj, work) bind(C)
101 use iso_c_binding
102 type(c_ptr), value :: obj
103 real(c_double), value :: work
104 end subroutine
105 subroutine all_set_work_multi_c(obj, work, dim) bind(C)
106 use iso_c_binding
107 type(c_ptr), value :: obj
108 integer(c_int), value :: dim
109 real(c_double), dimension(dim) :: work
110 end subroutine
111 subroutine all_set_vertices_c(obj, n, dim, vertices) bind(C)
112 use iso_c_binding
113 type(c_ptr), value :: obj
114 integer(c_int), value :: n
115 integer(c_int), value :: dim
116 real(c_double), dimension(n*dim) :: vertices
117 end subroutine
118 subroutine all_set_communicator_c(obj, comm) bind(C)
119 use iso_c_binding
120 type(c_ptr), value :: obj
121 integer, value :: comm
122 end subroutine
123 subroutine all_set_sys_size_c(obj, size, dim) bind(C)
124 use iso_c_binding
125 type(c_ptr), value :: obj
126 integer(c_int), value :: dim
127 real(c_double), dimension(dim) :: size
128 end subroutine
129 subroutine all_set_method_data_histogram_c(obj, nbins) bind(C)
130 use iso_c_binding
131 type(c_ptr), value :: obj
132 integer(c_int), dimension(3) :: nbins
133 end subroutine
134 subroutine all_setup_c(obj) bind(C)
135 use iso_c_binding
136 type(c_ptr), value :: obj
137 end subroutine
138 subroutine all_balance_c(obj) bind(C)
139 use iso_c_binding
140 type(c_ptr), value :: obj
141 end subroutine
142 subroutine all_get_gamma_c(obj, gamma) bind(C)
143 use iso_c_binding
144 type(c_ptr), value :: obj
145 real(c_double) :: gamma
146 end subroutine
147 subroutine all_get_number_of_vertices_c(obj, n) bind(C)
148 use iso_c_binding
149 type(c_ptr), value :: obj
150 integer(c_int) :: n
151 end subroutine
152 subroutine all_get_vertices_c(obj, n, vertices) bind(C)
153 use iso_c_binding
154 type(c_ptr), value :: obj
155 integer(c_int), value :: n
156 real(c_double), dimension(*) :: vertices
157 end subroutine
158 subroutine all_get_prev_vertices_c(obj, n, vertices) bind(C)
159 use iso_c_binding
160 type(c_ptr), value :: obj
161 integer(c_int), value :: n
162 real(c_double), dimension(*) :: vertices
163 end subroutine
164 subroutine all_get_dimension_c(obj, dim) bind(C)
165 use iso_c_binding
166 type(c_ptr), value :: obj
167 integer(c_int) :: dim
168 end subroutine
169 subroutine all_get_length_of_work_c(obj, length) bind(C)
170 use iso_c_binding
171 type(c_ptr), value :: obj
172 integer(c_int) :: length
173 end subroutine
174 subroutine all_get_work_c(obj, work) bind(C)
175 use iso_c_binding
176 type(c_ptr), value :: obj
177 real(c_double) :: work
178 end subroutine
179 subroutine all_get_work_array_c(obj, work, length) bind(C)
180 use iso_c_binding
181 type(c_ptr), value :: obj
182 integer(c_int), value :: length
183 real(c_double), dimension(length) :: work
184 end subroutine
185 subroutine all_get_number_of_neighbors_c(obj, count) bind(C)
186 use iso_c_binding
187 type(c_ptr), value :: obj
188 integer(c_int) :: count
189 end subroutine
190 subroutine all_get_neighbors_c(obj, neighbors, count) bind(C)
191 use iso_c_binding
192 type(c_ptr), value :: obj
193 integer(c_int), value :: count
194 integer(c_int), dimension(count) :: neighbors
195 end subroutine
196#ifdef ALL_VTK_OUTPUT
197 subroutine all_print_vtk_outlines_c(obj, step) bind(C)
198 use iso_c_binding
199 type(c_ptr), value :: obj
200 integer(c_int), value :: step
201 end subroutine
202 subroutine all_print_vtk_vertices_c(obj, step) bind(C)
203 use iso_c_binding
204 type(c_ptr), value :: obj
205 integer(c_int), value :: step
206 end subroutine
207#endif
208 function all_errno_c() result(errno) bind(C)
209 use iso_c_binding
210 integer(c_int) :: errno
211 end function
212 subroutine all_reset_errno_c() bind(C)
213 end subroutine
214 subroutine all_errdesc_c(text, length) bind(C)
215 use iso_c_binding
216 character(len=1, kind=c_char) :: text(*)
217 integer(c_size_t), value :: length
218 end subroutine
219 end interface
220
221 !> The object oriented API is this object. It contains all relevant functions
222 type, public :: all_t
223 type(c_ptr), private :: object = c_null_ptr !< pointer to C++ object used on C side
224 integer, private :: dim = 0 !< dimensionality of system, set during init
225 contains
226 procedure :: init => all_init
227 procedure :: finalize => all_finalize
228 procedure :: set_gamma => all_set_gamma
229 procedure :: set_proc_grid_params => all_set_proc_grid_params
230 procedure :: set_proc_tag => all_set_proc_tag
231 procedure :: set_min_domain_size => all_set_min_domain_size
232 procedure :: set_work => all_set_work
233 procedure :: set_work_multi => all_set_work_multi
234 procedure :: set_vertices => all_set_vertices
235#ifdef ALL_USE_F08
236 procedure :: set_communicator => all_set_communicator_f08
237#else
238 procedure :: set_communicator => all_set_communicator_int
239#endif
240 procedure :: set_sys_size => all_set_sys_size
241 procedure :: set_method_data_histgram => all_set_method_data_histgram
242 procedure :: setup => all_setup
243 procedure :: balance => all_balance
244 procedure :: get_gamma => all_get_gamma
245 procedure :: get_number_of_vertices => all_get_number_of_vertices
246 procedure :: get_vertices => all_get_vertices
247 procedure :: get_vertices_alloc => all_get_vertices_alloc
248 procedure :: get_prev_vertices => all_get_prev_vertices
249 procedure :: get_dimension => all_get_dimension
250 procedure :: get_length_of_work => all_get_length_of_work
251 procedure :: get_work => all_get_work
252 procedure :: get_work_array => all_get_work_array
253 procedure :: get_number_of_neighbors => all_get_number_of_neighbors
254 procedure :: get_neighbors => all_get_neighbors
255#ifdef ALL_VTK_OUTPUT
256 procedure :: print_vtk_outlines => all_print_vtk_outlines
257 procedure :: print_vtk_vertices => all_print_vtk_vertices
258#endif
259 end type
260
261 ! This is the old, but still supported API of separate functions
262 public :: all_init
263 public :: all_finalize
264 public :: all_set_gamma
266 public :: all_set_proc_tag
268 public :: all_set_work
269 public :: all_set_work_multi
270 public :: all_set_vertices
271 public :: all_set_communicator
272 public :: all_set_sys_size
274 public :: all_setup
275 public :: all_balance
276 public :: all_get_gamma
278 public :: all_get_vertices
279 public :: all_get_vertices_alloc
280 public :: all_get_prev_vertices
281 public :: all_get_dimension
282 public :: all_get_length_of_work
283 public :: all_get_work
284 public :: all_get_work_array
286 public :: all_get_neighbors
287#ifdef ALL_VTK_OUTPUT
288 public :: all_print_vtk_outlines
289 public :: all_print_vtk_vertices
290#endif
291 public :: all_error
292 public :: all_reset_error
293 public :: all_error_description
294
296 module procedure all_set_communicator_int
297#ifdef ALL_USE_F08
298 module procedure all_set_communicator_f08
299#endif
300 end interface
301contains
302 !> Initialises the loadbalancer
303 subroutine all_init(this, method, dim, gamma)
304 class(all_t), intent(out) :: this !< teh ALL object is returned
305 integer, intent(in) :: method !< Must be one of the ALL_* method values
306 integer, intent(in) :: dim !< dimensionality of system
307 real(c_double), intent(in) :: gamma !< gamma value for load balancer (ignored for TENSOR and STAGGERED)
308 this%object = all_init_c(method, dim, gamma)
309 this%dim = dim
310 end subroutine
311 !> Delete the loadbalance object
312 subroutine all_finalize(this)
313 class(all_t), intent(in) :: this
314 call all_finalize_c(this%object)
315 end subroutine
316 !> Set gamma value for load balancer
317 subroutine all_set_gamma(this, gamma)
318 class(all_t), intent(in) :: this
319 real(c_double), intent(in) :: gamma
320 call all_set_gamma_c(this%object, gamma)
321 end subroutine
322 !> Set the grid parameters for this process
323 subroutine all_set_proc_grid_params(this, loc, ranks)
324 class(all_t), intent(in) :: this
325 integer, dimension(this%dim), intent(in) :: loc !< index of domain in `dim` directions (0-indexed!)
326 integer, dimension(this%dim), intent(in) :: ranks !< total number of domains in `dim` directions
327 call all_set_proc_grid_params_c(this%object, this%dim, loc, this%dim, ranks)
328 end subroutine
329 !> Set process identifying tag for output
330 subroutine all_set_proc_tag(this, tag)
331 class(all_t), intent(in) :: this
332 integer, intent(in) :: tag !< tag of local process, only output in VTK outlines
333 call all_set_proc_tag_c(this%object, tag)
334 end subroutine
335 !> Set a minimum domain size
336 subroutine all_set_min_domain_size(this, domain_size)
337 class(all_t), intent(in) :: this
338 real(c_double), dimension(this%dim), intent(in) :: domain_size !< minimum domain size
339 call all_set_min_domain_size_c(this%object, this%dim, domain_size)
340 end subroutine
341 !> Set the work of this process
342 subroutine all_set_work(this, work)
343 class(all_t), intent(in) :: this
344 real(c_double), intent(in) :: work !< work of this domain
345 call all_set_work_c(this%object, work)
346 end subroutine
347 !> Set multi dimensional work of this process
348 subroutine all_set_work_multi(this, work)
349 class(all_t), intent(in) :: this
350 real(c_double), dimension(:), intent(in) :: work !< multi dimensional work of this domain
351#ifndef NDEBUG
352 if (size(work) /= this%dim) then
353 write (error_unit, '(a)') "ALL_set_work_multi: Wrong dimension for work"
354 stop
355 end if
356#endif
357 call all_set_work_multi_c(this%object, work, size(work))
358 end subroutine
359 !> Set new vertices
360 !!
361 !! @todo write interface for single rank array of vertices
362 subroutine all_set_vertices(this, vertices)
363 class(all_t), intent(in) :: this
364 real(c_double), dimension(:, :), intent(in) :: vertices !< vertices of domain, for `n` domains must have shape: (dim,n)
365#ifndef NDEBUG
366 if (size(vertices, 1) /= this%dim) then
367 write (error_unit, '(a)') "ALL_set_vertices: Wrong dimension for vertices"
368 stop
369 end if
370#endif
371 call all_set_vertices_c(this%object, size(vertices, 2), this%dim, vertices)
372 end subroutine
373 !> Set the MPI communicator with an mpi_f08 object
374#ifdef ALL_USE_F08
375 subroutine all_set_communicator_f08(this, comm)
376 use mpi_f08
377 class(all_t), intent(in) :: this
378 type(mpi_comm), intent(in) :: comm !< MPI Communicator
379 call all_set_communicator_c(this%object, comm%mpi_val)
380 end subroutine
381#endif
382 !> Set the MPI communicator with an ``mpi`` oder ``mpif.h`` communicator
383 subroutine all_set_communicator_int(this, comm)
384 class(all_t), intent(in) :: this
385 integer, intent(in) :: comm !< MPI Communicator, not type checked!
386 call all_set_communicator_c(this%object, comm)
387 end subroutine
388 !> Set size of system, which is required for the histogram method
389 subroutine all_set_sys_size(this, syssize)
390 class(all_t), intent(in) :: this
391 real(c_double), dimension(:), intent(in) :: syssize !< system size
392#ifndef NDEBUG
393 if (size(syssize) /= this%dim) then
394 write (error_unit, '(a)') "ALL_set_size: Wrong dimension for size"
395 stop
396 end if
397#endif
398 call all_set_sys_size_c(this%object, syssize, size(syssize))
399 end subroutine
400 !> Set number of bins for histogram method
401 subroutine all_set_method_data_histgram(this, nbins)
402 class(all_t), intent(in) :: this
403 integer(c_int), dimension(3), intent(in) :: nbins !< Number of bins per dimension
404 call all_set_method_data_histogram_c(this%object, nbins)
405 end subroutine
406 !> Set up the loadbalancer
407 subroutine all_setup(this)
408 class(all_t), intent(in) :: this
409 call all_setup_c(this%object)
410 end subroutine
411 !> Run loadbalancer and calculate new vertices
412 subroutine all_balance(this)
413 class(all_t), intent(in) :: this
414 call all_balance_c(this%object)
415 end subroutine
416 !> Retrieve currently set gamma value of balancer
417 subroutine all_get_gamma(this, gamma)
418 class(all_t), intent(in) :: this
419 real(c_double), intent(out) :: gamma
420 call all_get_gamma_c(this%object, gamma)
421 end subroutine
422 !> Retrieve number of vertices for current vertices
423 subroutine all_get_number_of_vertices(this, n)
424 class(all_t), intent(in) :: this
425 integer(c_int), intent(out) :: n !< set to number of new vertices
426 call all_get_number_of_vertices_c(this%object, n)
427 end subroutine
428 !> Retrieve current vertices
429 subroutine all_get_vertices(this, vertices)
430 class(all_t), intent(in) :: this
431 real(c_double), dimension(:, :), intent(out) :: vertices !< set to new vertices, must be exact size (dim,n)
432 integer(c_int) :: n_verts
433 call this%get_number_of_vertices(n_verts)
434#ifndef NDEBUG
435 if (size(vertices, 1) /= this%dim .or. size(vertices, 2) < n_verts) then
436 write (error_unit, '(a)') "ALL_get_vertices: vertices array not large enough!"
437 stop
438 end if
439#endif
440 call all_get_vertices_c(this%object, n_verts, vertices)
441 end subroutine
442 !> Same function as get_vertices, but takes an allocatable array, and resizes automatically
443 subroutine all_get_vertices_alloc(this, vertices)
444 class(all_t), intent(in) :: this
445 real(c_double), dimension(:, :), allocatable, intent(inout) :: vertices !< set to new vertices, may be reallocated to fit
446 integer(c_int) :: n_verts
447 call this%get_number_of_vertices(n_verts)
448 if (allocated(vertices)) then
449 if (size(vertices, 1) /= this%dim .or. size(vertices, 2) < n_verts) then
450 deallocate (vertices)
451 end if
452 end if
453 if (.not. allocated(vertices)) allocate (vertices(this%dim, n_verts))
454 call all_get_vertices_c(this%object, n_verts, vertices)
455 end subroutine
456 !> Retrieve vertices from before loadbalancing
457 subroutine all_get_prev_vertices(this, vertices)
458 class(all_t), intent(in) :: this
459 real(c_double), dimension(:, :), intent(out) :: vertices !< set to prev vertices, must be exact size (dim,n), unchecked!
460#ifndef NDEBUG
461 ! TODO Check against number of vertices not only dimensionality
462 if (size(vertices, 1) /= this%dim) then
463 write (error_unit, '(a)') "ALL_get_prev_vertices: vertices array not large enough!"
464 stop
465 end if
466#endif
467 call all_get_prev_vertices_c(this%object, size(vertices, 2), vertices)
468 end subroutine
469 !> Get current dimension from loadbalancer
470 subroutine all_get_dimension(this, dim)
471 class(all_t), intent(in) :: this
472 integer(c_int), intent(out) :: dim
473 call all_get_dimension_c(this%object, dim)
474 end subroutine
475 !> Retrieve length of work array
476 subroutine all_get_length_of_work(this, length)
477 class(all_t), intent(in) :: this
478 integer(c_int), intent(out) :: length
479 call all_get_length_of_work_c(this%object, length)
480 end subroutine
481 !> Retrieve first element of work array
482 subroutine all_get_work(this, work)
483 class(all_t), intent(in) :: this
484 real(c_double), intent(out) :: work
485 call all_get_work_c(this%object, work)
486 end subroutine
487 !> Retrieve work array, which must already be the correct size
488 subroutine all_get_work_array(this, work)
489 class(all_t), intent(in) :: this
490 real(c_double), dimension(:), intent(out) :: work
491 integer :: length
492 call all_get_length_of_work(this, length)
493 if (size(work) /= length) then
494 write (error_unit, '(a)') "ALL_get_work_array: work has wrong length!"
495 stop
496 end if
497 call all_get_work_array_c(this%object, work, size(work))
498 end subroutine
499 !> Retrieve number of neighbors (i.e. length of neighbors list)
500 subroutine all_get_number_of_neighbors(this, count)
501 class(all_t), intent(in) :: this
502 integer(c_int), intent(out) :: count
503 call all_get_number_of_vertices_c(this%object, count)
504 end subroutine
505 !> Retrieve list of neighboring ranks (must be correct size already)
506 subroutine all_get_neighbors(this, neighbors)
507 class(all_t), intent(in) :: this
508 integer(c_int), dimension(:), intent(out) :: neighbors
509 integer :: count
510 call all_get_number_of_neighbors(this, count)
511 if (size(neighbors) /= count) then
512 write (error_unit, '(a)') "ALL_get_neighbors: neighbors has wrong length!"
513 stop
514 end if
515 call all_get_neighbors_c(this%object, neighbors, size(neighbors))
516 end subroutine
517
518#ifdef ALL_VTK_OUTPUT
519 !> Print VTK outlines (must be enabled in build step)
520 subroutine all_print_vtk_outlines(this, step)
521 class(all_t), intent(in) :: this
522 integer(c_int), intent(in) :: step !< current step, used for filename
523 call all_print_vtk_outlines_c(this%object, step)
524 end subroutine
525 !> Print VTK domain vertices (must be enabled in build step)
526 subroutine all_print_vtk_vertices(this, step)
527 class(all_t), intent(in) :: this
528 integer(c_int), intent(in) :: step !< current step, used for filename
529 call all_print_vtk_vertices_c(this%object, step)
530 end subroutine
531#endif
532 function all_error() result(error)
533 integer :: error
534 error = all_errno_c()
535 end function
536 subroutine all_reset_error()
537 call all_reset_errno_c()
538 end subroutine
539 function all_error_description() result(desc)
540 character(len=ALL_ERROR_LENGTH) :: desc
541 call all_errdesc_c(desc, int(all_error_length, c_size_t))
542 end function
543end module
544
545! VIM modeline for indentation
546! vim: sw=4 et
int all_errno_c(void)
void all_get_prev_vertices_c(ALL_t *all_obj, int n_vertices, double *prevVertices)
void all_set_gamma_c(ALL_t *all_obj, double gamma)
void all_setup_c(ALL_t *all_obj)
void all_balance_c(ALL_t *all_obj)
void all_set_work_c(ALL_t *all_obj, double work)
void all_get_length_of_work_c(ALL_t *all_obj, int *length)
ALL_t * all_init_c(ALL::LB_t method, const int dim, double gamma)
void all_errdesc_c(char *description, size_t len)
void all_get_vertices_c(ALL_t *all_obj, int n_vertices, double *vertices)
void all_get_work_array_c(ALL_t *all_obj, double *work, int length)
void all_finalize_c(ALL_t *all_obj)
void all_set_sys_size_c(ALL_t *all_obj, double *size, int dim)
void all_set_communicator_c(ALL_t *all_obj, MPI_Fint fcomm)
void all_set_min_domain_size_c(ALL_t *all_obj, int dim, double *domain_size)
void all_get_number_of_vertices_c(ALL_t *all_obj, int *n_vertices)
void all_set_work_multi_c(ALL_t *all_obj, double *work, int dim)
void all_set_vertices_c(ALL_t *all_obj, const int n, const int dim, const double *vertices)
void all_get_number_of_neighbors_c(ALL_t *all_obj, int *count)
void all_reset_errno_c()
void all_print_vtk_outlines_c(ALL_t *all_obj known_unused, int known_unused step)
void all_print_vtk_vertices_c(ALL_t *all_obj known_unused, int known_unused step)
void all_get_gamma_c(ALL_t *all_obj, double *gamma)
void all_get_neighbors_c(ALL_t *all_obj, int *neighbors, int count)
void all_set_proc_tag_c(ALL_t *all_obj, int tag)
void all_set_proc_grid_params_c(ALL_t *all_obj, int nloc, int *loc, int nsize, int *size)
void all_get_work_c(ALL_t *all_obj, double *work)
void all_set_method_data_histogram_c(ALL_t *all_obj, int *nbins)
void all_get_dimension_c(ALL_t *all_obj, int *dim)
The function API for ALL.
subroutine, public all_get_neighbors(this, neighbors)
Retrieve list of neighboring ranks (must be correct size already)
subroutine, public all_get_prev_vertices(this, vertices)
Retrieve vertices from before loadbalancing.
subroutine, public all_set_proc_grid_params(this, loc, ranks)
Set the grid parameters for this process.
integer(c_int), parameter, public all_error_outofbounds
integer(c_int), parameter, public all_staggered
subroutine, public all_get_work_array(this, work)
Retrieve work array, which must already be the correct size.
integer function, public all_error()
subroutine, public all_get_dimension(this, dim)
Get current dimension from loadbalancer.
subroutine, public all_set_min_domain_size(this, domain_size)
Set a minimum domain size.
subroutine, public all_reset_error()
subroutine, public all_get_work(this, work)
Retrieve first element of work array.
subroutine, public all_init(this, method, dim, gamma)
Initialises the loadbalancer.
subroutine, public all_set_work_multi(this, work)
Set multi dimensional work of this process.
subroutine, public all_setup(this)
Set up the loadbalancer.
integer(c_int), parameter, public all_forcebased
subroutine, public all_set_proc_tag(this, tag)
Set process identifying tag for output.
subroutine, public all_set_vertices(this, vertices)
Set new vertices.
subroutine, public all_print_vtk_vertices(this, step)
Print VTK domain vertices (must be enabled in build step)
subroutine, public all_get_length_of_work(this, length)
Retrieve length of work array.
subroutine, public all_print_vtk_outlines(this, step)
Print VTK outlines (must be enabled in build step)
subroutine, public all_get_gamma(this, gamma)
Retrieve currently set gamma value of balancer.
integer(c_int), parameter, public all_voronoi
subroutine, public all_get_number_of_vertices(this, n)
Retrieve number of vertices for current vertices.
integer, parameter, public all_error_length
integer(c_int), parameter, public all_unimplemented
integer(c_int), parameter, public all_error_invalidcommtype
subroutine, public all_finalize(this)
Delete the loadbalance object.
subroutine, public all_get_number_of_neighbors(this, count)
Retrieve number of neighbors (i.e. length of neighbors list)
subroutine, public all_set_work(this, work)
Set the work of this process.
subroutine, public all_balance(this)
Run loadbalancer and calculate new vertices.
subroutine, public all_set_sys_size(this, syssize)
Set size of system, which is required for the histogram method.
integer(c_int), parameter, public all_error_pointdimensionmissmatch
subroutine, public all_set_gamma(this, gamma)
Set gamma value for load balancer.
integer(c_int), parameter, public all_error_filesystem
integer(c_int), parameter, public all_tensor
integer(c_int), parameter, public all_error_generic
integer(c_int), parameter, public all_error_invalidargument
integer(c_int), parameter, public all_histogram
subroutine, public all_get_vertices(this, vertices)
Retrieve current vertices.
subroutine, public all_get_vertices_alloc(this, vertices)
Same function as get_vertices, but takes an allocatable array, and resizes automatically.
integer(c_int), parameter, public all_error_internal
integer(c_int), parameter, public all_tensor_max
character(len=all_error_length) function, public all_error_description()
subroutine, public all_set_method_data_histgram(this, nbins)
Set number of bins for histogram method.
The object oriented API is this object. It contains all relevant functions.