ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_test_f_obj.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
30! module for ALL access from Fortran
32 use all_module
33 use iso_c_binding
34#ifndef ALL_USE_F08
35 use mpi
36#else
37 use mpi_f08
38#endif
39 implicit none
40
41#ifndef ALL_USE_F08
42 integer :: cart_comm
43#else
44 type(mpi_comm) :: cart_comm
45#endif
46 integer :: dims(3)
47 integer :: coords(3)
48 logical :: period(3)
49 real(8) :: length(3)
50 real(8) :: vertices(3,2)
51 integer :: error
52 integer :: rank, n_ranks
53 real(8),dimension(:,:),allocatable :: new_vertices
54 integer :: i
55
56 type(all_t) :: lb
57
58 call mpi_init(error)
59
60 call mpi_comm_rank(mpi_comm_world,rank,error)
61 call mpi_comm_size(mpi_comm_world,n_ranks,error)
62
63 ! create cartesian communicator
64 dims = 0
65 period = .true.
66 call mpi_dims_create(n_ranks,3,dims, error)
67 call mpi_cart_create(mpi_comm_world, 3, dims, period, .true., cart_comm, error)
68
69 ! compute cell length
70 length = 1.0d0 / real(dims,8)
71
72 ! compute vertices of local domain
73 call mpi_cart_coords(cart_comm, rank, 3, coords, error)
74 vertices(:,1) = real(coords,8) * length
75 vertices(:,2) = real(coords+1,8) * length
76
77 do i = 0, n_ranks-1
78 if (i == rank) then
79 write(*,"(a,i7,2(a,3es14.7))") "rank: ", rank, " old vertices: ", vertices(:,1), ", ", vertices(:,2)
80 call mpi_barrier(cart_comm, error)
81 else
82 call mpi_barrier(cart_comm, error)
83 end if
84 end do
85
86 call lb%init(all_staggered, 3, 4.0d0)
87 call lb%set_proc_tag(rank)
88 call lb%set_work(real( product(coords,1)*64,8))
89 call lb%set_vertices(vertices)
90 ! if using a non cartesian communicator, this call would be required
91 ! call lb%set_proc_grid_params(coords,dims)
92 call lb%set_communicator(cart_comm)
93
94 call lb%setup()
95 call lb%balance()
96
97 ! If not allocatable as target, must use
98 ! lb%get_number_of_vertices and lb%get_vertices instead!
99 call lb%get_vertices_alloc(new_vertices)
100
101 do i = 0, n_ranks-1
102 if (i == rank) then
103 write(*,"(2(a,i7),2(a,3es14.7))") "rank: ", rank, " ", size(new_vertices,2), " new vertices: ", new_vertices(:,1), &
104 ", ", new_vertices(:,2)
105 call mpi_barrier(cart_comm, error)
106 else
107 call mpi_barrier(cart_comm, error)
108 end if
109 end do
110 deallocate(new_vertices)
111 call lb%finalize()
112 call mpi_finalize(error);
113
114END PROGRAM
115
program all_test_f
The function API for ALL.
integer(c_int), parameter, public all_staggered
The object oriented API is this object. It contains all relevant functions.