46#ifdef ALL_VORONOI_ACTIVE
65 function all_init_c(method, dim, gamma)
result(this)
bind(C)
67 integer(c_int),
value :: method
68 integer(c_int),
value :: dim
69 real(c_double),
value :: gamma
74 type(c_ptr),
value :: obj
78 type(c_ptr),
value :: obj
79 real(c_double),
value :: gamma
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
91 type(c_ptr),
value :: obj
92 integer(c_int),
value :: tag
96 type(c_ptr),
value :: obj
97 integer(c_int),
value :: dim
98 real(c_double),
dimension(dim) :: domain_size
102 type(c_ptr),
value :: obj
103 real(c_double),
value :: work
107 type(c_ptr),
value :: obj
108 integer(c_int),
value :: dim
109 real(c_double),
dimension(dim) :: work
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
120 type(c_ptr),
value :: obj
121 integer,
value :: comm
125 type(c_ptr),
value :: obj
126 integer(c_int),
value :: dim
127 real(c_double),
dimension(dim) :: size
131 type(c_ptr),
value :: obj
132 integer(c_int),
dimension(3) :: nbins
136 type(c_ptr),
value :: obj
140 type(c_ptr),
value :: obj
144 type(c_ptr),
value :: obj
145 real(c_double) :: gamma
149 type(c_ptr),
value :: obj
154 type(c_ptr),
value :: obj
155 integer(c_int),
value :: n
156 real(c_double),
dimension(*) :: vertices
160 type(c_ptr),
value :: obj
161 integer(c_int),
value :: n
162 real(c_double),
dimension(*) :: vertices
166 type(c_ptr),
value :: obj
167 integer(c_int) :: dim
171 type(c_ptr),
value :: obj
172 integer(c_int) :: length
176 type(c_ptr),
value :: obj
177 real(c_double) :: work
181 type(c_ptr),
value :: obj
182 integer(c_int),
value :: length
183 real(c_double),
dimension(length) :: work
187 type(c_ptr),
value :: obj
188 integer(c_int) :: count
192 type(c_ptr),
value :: obj
193 integer(c_int),
value :: count
194 integer(c_int),
dimension(count) :: neighbors
199 type(c_ptr),
value :: obj
200 integer(c_int),
value :: step
204 type(c_ptr),
value :: obj
205 integer(c_int),
value :: step
210 integer(c_int) :: errno
216 character(len=1, kind=c_char) :: text(*)
217 integer(c_size_t),
value :: length
223 type(c_ptr),
private :: object = c_null_ptr
224 integer,
private :: dim = 0
236 procedure :: set_communicator => all_set_communicator_f08
238 procedure :: set_communicator => all_set_communicator_int
296 module procedure all_set_communicator_int
298 module procedure all_set_communicator_f08
304 class(
all_t),
intent(out) :: this
305 integer,
intent(in) :: method
306 integer,
intent(in) :: dim
307 real(c_double),
intent(in) :: gamma
313 class(
all_t),
intent(in) :: this
318 class(
all_t),
intent(in) :: this
319 real(c_double),
intent(in) :: gamma
324 class(
all_t),
intent(in) :: this
325 integer,
dimension(this%dim),
intent(in) :: loc
326 integer,
dimension(this%dim),
intent(in) :: ranks
331 class(
all_t),
intent(in) :: this
332 integer,
intent(in) :: tag
337 class(
all_t),
intent(in) :: this
338 real(c_double),
dimension(this%dim),
intent(in) :: domain_size
343 class(
all_t),
intent(in) :: this
344 real(c_double),
intent(in) :: work
349 class(
all_t),
intent(in) :: this
350 real(c_double),
dimension(:),
intent(in) :: work
352 if (
size(work) /= this%dim)
then
353 write (error_unit,
'(a)')
"ALL_set_work_multi: Wrong dimension for work"
363 class(
all_t),
intent(in) :: this
364 real(c_double),
dimension(:, :),
intent(in) :: vertices
366 if (
size(vertices, 1) /= this%dim)
then
367 write (error_unit,
'(a)')
"ALL_set_vertices: Wrong dimension for vertices"
375 subroutine all_set_communicator_f08(this, comm)
377 class(
all_t),
intent(in) :: this
378 type(mpi_comm),
intent(in) :: comm
383 subroutine all_set_communicator_int(this, comm)
384 class(
all_t),
intent(in) :: this
385 integer,
intent(in) :: comm
390 class(
all_t),
intent(in) :: this
391 real(c_double),
dimension(:),
intent(in) :: syssize
393 if (
size(syssize) /= this%dim)
then
394 write (error_unit,
'(a)')
"ALL_set_size: Wrong dimension for size"
402 class(
all_t),
intent(in) :: this
403 integer(c_int),
dimension(3),
intent(in) :: nbins
408 class(
all_t),
intent(in) :: this
413 class(
all_t),
intent(in) :: this
418 class(
all_t),
intent(in) :: this
419 real(c_double),
intent(out) :: gamma
424 class(
all_t),
intent(in) :: this
425 integer(c_int),
intent(out) :: n
430 class(
all_t),
intent(in) :: this
431 real(c_double),
dimension(:, :),
intent(out) :: vertices
432 integer(c_int) :: n_verts
433 call this%get_number_of_vertices(n_verts)
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!"
444 class(
all_t),
intent(in) :: this
445 real(c_double),
dimension(:, :),
allocatable,
intent(inout) :: vertices
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)
453 if (.not.
allocated(vertices))
allocate (vertices(this%dim, n_verts))
458 class(
all_t),
intent(in) :: this
459 real(c_double),
dimension(:, :),
intent(out) :: vertices
462 if (
size(vertices, 1) /= this%dim)
then
463 write (error_unit,
'(a)')
"ALL_get_prev_vertices: vertices array not large enough!"
471 class(
all_t),
intent(in) :: this
472 integer(c_int),
intent(out) :: dim
477 class(
all_t),
intent(in) :: this
478 integer(c_int),
intent(out) :: length
483 class(
all_t),
intent(in) :: this
484 real(c_double),
intent(out) :: work
489 class(
all_t),
intent(in) :: this
490 real(c_double),
dimension(:),
intent(out) :: work
493 if (
size(work) /= length)
then
494 write (error_unit,
'(a)')
"ALL_get_work_array: work has wrong length!"
501 class(
all_t),
intent(in) :: this
502 integer(c_int),
intent(out) :: count
507 class(
all_t),
intent(in) :: this
508 integer(c_int),
dimension(:),
intent(out) :: neighbors
511 if (
size(neighbors) /= count)
then
512 write (error_unit,
'(a)')
"ALL_get_neighbors: neighbors has wrong length!"
521 class(
all_t),
intent(in) :: this
522 integer(c_int),
intent(in) :: step
527 class(
all_t),
intent(in) :: this
528 integer(c_int),
intent(in) :: step
540 character(len=ALL_ERROR_LENGTH) :: desc
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_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.