ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_test_f.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!
4!Redistribution and use in source and binary forms, with or without modification, are
5!permitted provided that the following conditions are met:
6!
7!1. Redistributions of source code must retain the above copyright notice, this list
8! of conditions and the following disclaimer.
9!
10!2. Redistributions in binary form must reproduce the above copyright notice, this
11! list of conditions and the following disclaimer in the documentation and/or
12! other materials provided with the distribution.
13!
14!3. Neither the name of the copyright holder nor the names of its contributors
15! may be used to endorse or promote products derived from this software without
16! specific prior written permission.
17!
18!THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
19!EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
20!OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
21!SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
22!INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
23!TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
24!BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25!CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
26!ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
27!SUCH DAMAGE.
28
29! module for ALL access from Fortran
31 USE all_module
32 USE iso_c_binding
33 use mpi
34 IMPLICIT NONE
35
36 INTEGER :: cart_comm
37 INTEGER :: dims(3)
38 INTEGER :: coords(3)
39 LOGICAL :: period(3)
40 REAL(8) :: length(3)
41 REAL(8) :: vertices(3,2)
42 INTEGER :: error
43 INTEGER :: rank, n_ranks
44 INTEGER :: n_vertices
45 REAL(8),DIMENSION(:,:),ALLOCATABLE :: new_vertices
46 INTEGER :: i
47 CHARACTER(LEN=ALL_ERROR_LENGTH) :: all_error_desc
48
49 TYPE(all_t) :: obj
50
51 CALL mpi_init(error)
52
53 CALL mpi_comm_rank(mpi_comm_world,rank,error)
54 CALL mpi_comm_size(mpi_comm_world,n_ranks,error)
55
56 ! create cartesian communicator
57 dims = 0
58 period = .true.
59 CALL mpi_dims_create(n_ranks,3,dims, error)
60 CALL mpi_cart_create(mpi_comm_world, 3, dims, period, .true., cart_comm, error)
61
62 ! compute cell length
63 length = 1.0d0 / real(dims,8)
64
65 ! compute vertices of local domain
66 CALL mpi_cart_coords(cart_comm, rank, 3, coords, error)
67 vertices(:,1) = real(coords,8) * length
68 vertices(:,2) = real(coords+1,8) * length
69
70 DO i = 0, n_ranks-1
71 IF (i == rank) THEN
72 WRITE(*,"(a,i7,2(a,3es14.7))") "rank: ", rank, " old vertices: ", vertices(:,1), ", ", vertices(:,2)
73 call mpi_barrier(cart_comm, error)
74 ELSE
75 call mpi_barrier(cart_comm, error)
76 END IF
77 END DO
78
79 CALL all_reset_error()
80 CALL all_init(obj,all_staggered,3,4.0d0)
81 IF (all_error() /= 0) THEN
82 print*, "Error: ", all_error()
83 all_error_desc= all_error_description()
84 print*, "Message: ", all_error_desc
85 call mpi_abort(mpi_comm_world, 1, error)
86 stop
87 END IF
88 CALL all_set_proc_tag(obj, rank)
89 CALL all_set_work(obj, real( product(coords,1)*64,8) )
90 CALL all_set_vertices(obj,vertices)
91 ! if using a non cartesian communicator, this call would be required
92 !CALL ALL_set_proc_grid_params(obj,coords,dims)
93 CALL all_set_communicator(obj,cart_comm)
94
95 CALL all_setup(obj)
96 CALL all_balance(obj)
97
98 CALL all_get_number_of_vertices(obj,n_vertices)
99 ALLOCATE(new_vertices(3,n_vertices))
100 CALL all_get_vertices(obj,new_vertices)
101
102 DO i = 0, n_ranks-1
103 IF (i == rank) THEN
104 WRITE(*,"(2(a,i7),2(a,3es14.7))") "rank: ", rank, " ", n_vertices, " new vertices: ", new_vertices(:,1), &
105 ", ", new_vertices(:,2)
106 call mpi_barrier(cart_comm, error)
107 ELSE
108 call mpi_barrier(cart_comm, error)
109 END IF
110 END DO
111 DEALLOCATE(new_vertices)
112 CALL all_finalize(obj);
113 CALL mpi_finalize(error);
114
115END PROGRAM
116
program all_test_f
The function API for ALL.
integer(c_int), parameter, public all_staggered
integer function, public all_error()
subroutine, public all_reset_error()
subroutine, public all_init(this, method, dim, gamma)
Initialises the loadbalancer.
subroutine, public all_setup(this)
Set up the loadbalancer.
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_get_number_of_vertices(this, n)
Retrieve number of vertices for current vertices.
subroutine, public all_finalize(this)
Delete the loadbalance object.
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_get_vertices(this, vertices)
Retrieve current vertices.
character(len=all_error_length) function, public all_error_description()
The object oriented API is this object. It contains all relevant functions.