ALL 0.9.3
A Loadbalacing Library
Loading...
Searching...
No Matches
ALL_Staggered_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!Copyright 2020-2020 Stephan Schulz, Forschungszentrum Juelich GmbH, 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
31 use all_module
32 use iso_c_binding
33 use, intrinsic :: iso_fortran_env, only: stdin=>input_unit&
34 , stdout=>output_unit&
35 , stderr=>error_unit&
36 , file_storage_size
37#ifndef ALL_USE_F08
38 use mpi
39#else
40 use mpi_f08
41#endif
42 implicit none
43 integer, parameter :: real64 = selected_real_kind(15)
44
45 call main()
46contains
47 subroutine print_loc(my_rank, maximum_rank, my_location, number_of_processes)
48 integer, intent(in) :: my_rank, maximum_rank
49 integer, dimension(3), intent(in) :: my_location, number_of_processes
50 integer :: rank, error
51 do rank=0, maximum_rank-1
52 if(rank == my_rank) then
53 write(stdout,'(a,i3.3,a,i3,",",i3,",",i3,a,i3,",",i3,",",i3,a)')&
54 "[", my_rank, "] Location: (", my_location, ")/(", number_of_processes, ")"
55 flush(stdout)
56 call mpi_barrier(mpi_comm_world, error)
57 else
58 call mpi_barrier(mpi_comm_world, error)
59 endif
60 enddo
61 end subroutine
62
63 subroutine print_domain(my_rank, maximum_rank, domain_vertices)
64 integer, intent(in) :: my_rank, maximum_rank
65 real(real64), dimension(3,2), intent(in) :: domain_vertices
66 integer :: rank, error
67 do rank=0, maximum_rank-1
68 if(rank == my_rank) then
69 write(stdout,'("[",i3.3,"]",a,es12.4,a,es12.4,a,es12.4)')&
70 my_rank, " Lower: ",&
71 domain_vertices(1,1), achar(9),&
72 domain_vertices(2,1), achar(9),&
73 domain_vertices(3,1)
74 write(stdout,'("[",i3.3,"]",a,es12.4,a,es12.4,a,es12.4)')&
75 my_rank, " Upper: ",&
76 domain_vertices(1,2), achar(9),&
77 domain_vertices(2,2), achar(9),&
78 domain_vertices(3,2)
79 flush(stdout)
80 call mpi_barrier(mpi_comm_world,error)
81 else
82 call mpi_barrier(mpi_comm_world,error)
83 endif
84 enddo
85 end subroutine
86
87 subroutine print_work(my_rank, maximum_rank, my_work)
88 integer, intent(in) :: my_rank, maximum_rank
89 real(real64), intent(in) :: my_work
90 integer :: rank, error
91 do rank=0, maximum_rank-1
92 if(rank == my_rank) then
93 write(stdout,'(a,i3.3,a,es12.4)')&
94 "[", my_rank, "] Work: ", my_work
95 flush(stdout)
96 call mpi_barrier(mpi_comm_world, error)
97 else
98 call mpi_barrier(mpi_comm_world, error)
99 endif
100 enddo
101 end subroutine
102
103 subroutine print_testing_output(my_rank, maximum_rank, new_vertices, timestep)
104 integer, intent(in) :: my_rank, maximum_rank, timestep
105 real(real64), dimension(3,2), intent(in) :: new_vertices
106 integer :: rank, error
107 do rank=0, maximum_rank-1
108 if(rank == my_rank) then
109 !write(stdout,'(a,i4,",",i3.3,a,3f11.6)')&
110 ! "[", timestep, my_rank, "] Result Width: ",&
111 ! new_vertices(:,2)-new_vertices(:,1)
112 !flush(stdout)
113 write(stdout,'(a,i4,",",i3.3,",",i2.2,a,3f11.6)')&
114 "[", timestep, my_rank, 0, "] Result Vertex: ",&
115 new_vertices(:,1)
116 flush(stdout)
117 write(stdout,'(a,i4,",",i3.3,",",i2.2,a,3f11.6)')&
118 "[", timestep, my_rank, 1, "] Result Vertex: ",&
119 new_vertices(:,2)
120 flush(stdout)
121 call mpi_barrier(mpi_comm_world, error)
122 else
123 call mpi_barrier(mpi_comm_world, error)
124 endif
125 enddo
126 end subroutine
127
128 subroutine print_binary(my_rank, my_work, new_vertices, my_location, number_of_processes, fh, timestep)
129 use iso_fortran_env, only: error_unit
130
131#ifndef ALL_USE_F08
132 integer, intent(in) :: fh
133 integer, dimension(MPI_STATUS_SIZE) :: state
134#else
135 type(mpi_file), intent(in) :: fh
136 type(mpi_status) :: state
137#endif
138 integer, intent(in) :: my_rank, timestep
139 integer, dimension(3), intent(in) :: my_location, number_of_processes
140 real(real64), intent(in) :: my_work
141 real(real64), dimension(3,2), intent(in) :: new_vertices
142 integer(MPI_OFFSET_KIND) :: offset
143 real(real64), dimension(11) :: buffer
144 integer :: error
145
146 offset = ((timestep-1) * number_of_processes(3) * 11 + my_rank * 11) * 8
147
148 buffer(1) = real(my_rank,real64)
149 buffer(2) = my_work
150 buffer(3) = new_vertices(1,1)
151 buffer(4) = new_vertices(2,1)
152 buffer(5) = new_vertices(3,1)
153 buffer(6) = new_vertices(1,2)
154 buffer(7) = new_vertices(2,2)
155 buffer(8) = new_vertices(3,2)
156 buffer(9) = my_location(1)
157 buffer(10) = my_location(2)
158 buffer(11) = my_location(3)
159
160 call mpi_file_write_at(fh, offset, buffer, 11, mpi_double_precision, state, error);
161 !write(ERROR_UNIT,"(11f10.3,a,i0,a,i0)") buffer," ",offset, " ", timestep
162
163 end subroutine
164
165 subroutine main()
166 integer :: error
167 integer :: current_step
168 integer, parameter :: number_of_steps = 50
169 integer, parameter :: dimensions = 3
170 real(real64), parameter :: loadbalancer_gamma = 0. ! ignored for staggered method
171 integer, dimension(dimensions) :: my_location, number_of_processes
172 real(real64), dimension(dimensions) :: minimum_domain_size
173 real(real64), dimension(dimensions,2) :: domain_vertices, new_vertices
174 real(real64), parameter :: domain_size = 1.0
175 real(real64) :: my_work
176 integer :: my_rank, maximum_rank
177 integer :: i,j
178 type(all_t) :: jall
179
180 ! number of command line arguments
181 integer :: n_clargs
182 ! command line arguments
183 character(256), dimension(:), allocatable :: clargs
184 character(256) :: filename, filename2
185 integer(MPI_OFFSET_KIND) :: offset
186 integer :: test_file
187 real(real64) :: d
188#ifndef ALL_USE_F08
189 integer :: fh
190 integer, dimension(MPI_STATUS_SIZE) :: state
191#else
192 type(mpi_file) :: fh
193 type(mpi_status) :: state
194#endif
195
196#ifdef ALL_VTK_OUTPUT_EXAMPLE
197 character (len=ALL_ERROR_LENGTH) :: error_msg
198#endif
199
200 call mpi_init(error)
201
202 ! The ALL_TENSOR method can be used just like the staggered method.
203 ! Just exchange ALL_STAGGERED with ALL_TENSOR in the init call.
204 call jall%init(all_staggered, dimensions, loadbalancer_gamma)
205 my_location(:) = 0
206 ! All domains are placed along a line in z direction, evne though they are three dimensional
207 call mpi_comm_rank(mpi_comm_world, my_location(3), error)
208 my_rank = my_location(3)
209
210 number_of_processes(:) = 1
211 call mpi_comm_size(mpi_comm_world, number_of_processes(3), error)
212 maximum_rank = number_of_processes(3)
213
214 n_clargs = command_argument_count()
215 if (n_clargs > 1) then
216 call get_command_argument(2,clargs(2))
217 write(filename,"(a,i0,a)") trim(clargs(2)),number_of_processes,trim(".bin")
218 else
219 write(filename,"(a,i0,a)") trim("./ALL_Staggered_f_"),number_of_processes(3),trim(".bin")
220 endif
221 call mpi_file_open(mpi_comm_world,filename,ior(mpi_mode_create, mpi_mode_wronly), mpi_info_null, fh, error)
222
223 if(my_rank == 0) then
224 write(stdout,'(a,i0)') "Ranks: ", maximum_rank
225 write(stdout,'(a,i0)') "Number of Steps: ", number_of_steps
226 flush(stdout)
227 endif
228
229 call mpi_barrier(mpi_comm_world, error)
230 call print_loc(my_rank, maximum_rank, my_location, number_of_processes)
231
232 ! For a cartesian communicator this is not required, but we are using
233 ! MPI_COMM_WORLD here.
234 call jall%set_proc_grid_params(my_location, number_of_processes)
235
236 ! If we want ot set a minimum domain size:
237 minimum_domain_size(:) = 0.1
238 call jall%set_min_domain_size(minimum_domain_size)
239
240 call jall%set_communicator(mpi_comm_world)
241
242 ! We also set the optional process tag for the output.
243 ! This can be useful if we want to know which of 'our'
244 ! ranks is which in the output produces by the library.
245 ! The ranks used inside the library do not necessarily
246 ! match our own numbering.
247 call jall%set_proc_tag(my_rank)
248
249 call jall%setup()
250
251 ! A first domain distribution must be given to the balancer.
252 ! We use the provided ALL::Point class to define the vertices,
253 ! but a simple double array can also be used. We need 2 vertices
254 ! which correspond to the two opposing corners.
255 ! We create a cubic domain initially
256 do i=1, 2
257 do j=1,dimensions
258 domain_vertices(j,i) = (my_location(j)+i-1) * domain_size
259 enddo
260 enddo
261 call print_domain(my_rank, maximum_rank, domain_vertices)
262 call jall%set_vertices(domain_vertices)
263
264 ! Calculate the work fo our domain. Here we just use
265 my_work = my_rank + 1.
266 call jall%set_work(my_work)
267 call print_work(my_rank, maximum_rank, my_work)
268 do current_step=1, number_of_steps
269 ! In a real code we need to set the updated work in each
270 ! iteration before calling balance()
271 if(my_rank == 0) then
272 write(stdout,'(a,i0,"/",i0)') "Starting step: ", current_step, number_of_steps
273 flush(stdout)
274 endif
275#ifdef ALL_VTK_OUTPUT_EXAMPLE
276 call all_reset_error()
277 call jall%print_vtk_outlines(current_step)
278 if(all_error() /= 0) then
279 error_msg = all_error_description()
280 write(stdout, '(a)') error_msg
281 ! Maybe also abort if the output is necesssary or handle this
282 ! some other way.
283 call mpi_abort(mpi_comm_world, 2, error)
284 stop
285 endif
286#endif
287 call jall%balance()
288
289 call jall%get_vertices(new_vertices)
290
291 !call print_testing_output(my_rank, maximum_rank, new_vertices, current_step)
292 call print_binary(my_rank, my_work, new_vertices, my_location, number_of_processes, fh, current_step)
293
294 call mpi_barrier(mpi_comm_world, error)
295 enddo
296#ifdef ALL_VTK_OUTPUT_EXAMPLE
297 call all_reset_error()
298 call jall%print_vtk_outlines(current_step)
299 if(all_error() /= 0) then
300 error_msg = all_error_description()
301 write(stdout, '(a)') error_msg
302 ! Maybe also abort if the output is necesssary or handle this
303 ! some other way.
304 call mpi_abort(mpi_comm_world, 2, error)
305 stop
306 endif
307#endif
308
309 call jall%finalize()
310 call mpi_file_close(fh, error)
311
312 call mpi_finalize(error)
313 end subroutine
314end program
void print_loc(int rank, int *loc, int *size)
int main(int argc, char **argv)
void print_binary(int step, int rank, double work, std::vector< ALL::Point< double > > &vertices, int *loc, int *size, MPI_File fh)
void print_testing_output(int rank, std::vector< ALL::Point< double > > &vertices, int timestep)
void print_domain(int rank, double *verts)
void print_work(int rank, double work)
program all_staggered_f
The function API for ALL.
integer(c_int), parameter, public all_staggered
integer function, public all_error()
subroutine, public all_reset_error()
character(len=all_error_length) function, public all_error_description()
The object oriented API is this object. It contains all relevant functions.