Skip to content
Snippets Groups Projects
Commit 5c80160b authored by Stephan Schulz's avatar Stephan Schulz
Browse files

add Tensor example

parent ba148928
Branches
Tags
No related merge requests found
/*
Copyright 2020-2020 Stephan Schulz, Forschungszentrum Juelich GmbH, Germany
Copyright 2018-2020 Rene Halver, Forschungszentrum Juelich GmbH, Germany
Copyright 2018-2020 Godehard Sutmann, Forschungszentrum Juelich GmbH, Germany
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:
1. Redistributions of source code must retain the above copyright notice,
this list of conditions and the following disclaimer.
2. Redistributions in binary form must reproduce the above copyright notice,
this list of conditions and the following disclaimer in the documentation
and/or other materials provided with the distribution.
3. Neither the name of the copyright holder nor the names of its contributors
may be used to endorse or promote products derived from this software without
specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
POSSIBILITY OF SUCH DAMAGE.
*/
#include <stdlib.h>
#include <stdio.h>
//#include <time.h>
#include <mpi.h>
//#define ALL_VTK_OUTPUT
#ifdef USE_AMALGAMATED
#include "ALL_Amalgamated.hpp"
#else
#include <ALL.hpp>
#endif
// Run fun in order of ranks
// Todo(s.schulz): This seems to only work roughly with the result width an 32 ranks, with up to 16 it seems to work correctly.
// Adding sleep(1) also orders everything correctly. So this is probably a flushing problem.
// It also exists for the cout stream with endl.
#define MPI_RUN_ORDER(comm, rank, max_ranks, fun) {int MPI_RO_IT;\
for(MPI_RO_IT=0;MPI_RO_IT<max_ranks;MPI_RO_IT++)\
{\
if(MPI_RO_IT==rank)\
{\
fun;\
MPI_Barrier(comm);\
} else {\
MPI_Barrier(comm);\
}\
}\
}
// Quick and dirty helper function. Assumes comm, rank and max_ranks
// CAVEAT: the function call must be wrapped in () if it contains a comma
#define MPI_RUN_ORDER_DEF(fun) MPI_RUN_ORDER(MPI_COMM_WORLD, MyRank, MaximumRank, fun)
//void millisleep(unsigned long ms)
//{
// struct timespec SleepTime;
// SleepTime.tv_sec = ms/1000;
// SleepTime.tv_nsec = (ms%1000)*1000000L;
// struct timespec RemainingTime;
// while(nanosleep(&SleepTime, &RemainingTime))
// {
// SleepTime.tv_sec = RemainingTime.tv_sec;
// SleepTime.tv_nsec = RemainingTime.tv_nsec;
// }
//}
void print_width(int rank, double width, double bottom, double top)
{
printf("[%03d] Result Width: %10.6f (%10.6f -- %10.6f)\n", rank, width, bottom, top);
fflush(stdout);
}
void print_testing_output(int rank, std::vector<ALL::Point<double>>& vertices, int timestep)
{
// printf("[%4d,%03d] Result Width: %10.6f %10.6f %10.6f\n",
// timestep,
// rank,
// vertices.at(1)[0]-vertices.at(0)[0],
// vertices.at(1)[1]-vertices.at(0)[1],
// vertices.at(1)[2]-vertices.at(0)[2]);
// fflush(stdout);
printf("[%4d,%03d,%02d] Result Vertex: %10.6f %10.6f %10.6f\n",
timestep,
rank,
0,
vertices.at(0)[0],
vertices.at(0)[1],
vertices.at(0)[2]);
fflush(stdout);
printf("[%4d,%03d,%02d] Result Vertex: %10.6f %10.6f %10.6f\n",
timestep,
rank,
1,
vertices.at(1)[0],
vertices.at(1)[1],
vertices.at(1)[2]);
fflush(stdout);
}
void print_loc(int rank, int* loc, int* size)
{
printf("[%03d] Location: (%3d,%3d,%3d)/(%3d,%3d,%3d)\n", rank, loc[0], loc[1], loc[2], size[0], size[1], size[2]);
fflush(stdout);
}
void print_domain(int rank, double* verts)
{
printf("[%03d] Lower: %g\t%g\t%g\n", rank, verts[0], verts[1], verts[2]);
printf("[%03d] Upper: %g\t%g\t%g\n", rank, verts[3], verts[4], verts[5]);
fflush(stdout);
}
void print_work(int rank, double work)
{
printf("[%03d] Work: %g\n", rank, work);
fflush(stdout);
}
void convert_verts(std::vector<ALL::Point<double>>* vv, double* verts)
{
verts[0] = vv->at(0)[0];
verts[1] = vv->at(0)[1];
verts[2] = vv->at(0)[2];
verts[3] = vv->at(1)[0];
verts[4] = vv->at(1)[1];
verts[5] = vv->at(1)[2];
}
int main(int argc, char** argv)
{
MPI_Init(&argc,&argv);
int CurrentStep = 0;
const int NumberOfSteps = 50;
const int Dimensions = 3;
const int LoadbalancerGamma = 0; // ignored for tensor method
// Todo(s.schulz): maybe change to automatic destruction on scope end
ALL::ALL<double, double> *jall = new ALL::ALL<double, double>(ALL::TENSOR, Dimensions, LoadbalancerGamma);
int MyLocation[3] = {0};
// All domains are placed along a line in z direction, even though they are three dimensional
MPI_Comm_rank(MPI_COMM_WORLD,&MyLocation[2]);
int MyRank = MyLocation[2];
int NumberOfProcesses[3] = {1,1,1};
MPI_Comm_size(MPI_COMM_WORLD, &NumberOfProcesses[2]);
int MaximumRank = NumberOfProcesses[2];
if(MyRank==0)
{
printf("Ranks: %d\nNumber of Steps: %d\n", MaximumRank, NumberOfSteps);
fflush(stdout);
}
MPI_Barrier(MPI_COMM_WORLD);
MPI_RUN_ORDER_DEF((print_loc(MyRank, MyLocation, NumberOfProcesses)));
// For a cartesian communicator this is not required, but we are using
// MPI_COMM_WORLD here.
std::vector<int> MyLocationVector(MyLocation, MyLocation+3);
std::vector<int> NumberOfProcessesVector(NumberOfProcesses, NumberOfProcesses+3);
jall->setProcGridParams(MyLocationVector, NumberOfProcessesVector);
// If we want to set a minimum domain size:
std::vector<double> MinimumDomainSize{0.1, 0.1, 0.1};
jall->setMinDomainSize(MinimumDomainSize);
jall->setCommunicator(MPI_COMM_WORLD);
// We also set the optional process tag for the output.
// This can be useful if we want to know which of 'our'
// ranks is which in the output produces by the library.
// The ranks used inside the library do not necessarily
// match our own numbering.
jall->setProcTag(MyRank);
jall->setup();
// Todo(s.schulz): document: what exactly must be set before setup()?
// A first domain distribution must be given to the balancer.
// We use the provided ALL::Point class to define the vertices,
// but a simple double array can also be used. We need 2 vertices
// which correspond to the two opposing corners.
std::vector<ALL::Point<double>> DomainVertices(2, ALL::Point<double>(3));
const double DomainSize = 1.0; // Domain size
// We create a cubic domain initially
for(int VertexIndex=0; VertexIndex<2; VertexIndex++)
{
for(int DimensionIndex=0; DimensionIndex<Dimensions; DimensionIndex++)
{
DomainVertices.at(VertexIndex)[DimensionIndex] = (MyLocation[DimensionIndex]+VertexIndex) * DomainSize;
}
}
double VertexArray[6];
convert_verts(&DomainVertices, VertexArray);
MPI_RUN_ORDER_DEF((print_domain(MyRank, VertexArray)));
jall->setVertices(DomainVertices);
// Calculate the work of our domain. Here we just use
double MyWork = (double) MyRank + 1.;
jall->setWork(MyWork);
MPI_RUN_ORDER_DEF((print_work(MyRank,MyWork)));
for(CurrentStep=0; CurrentStep<NumberOfSteps; CurrentStep++)
{
// In a real code we need to set the updated work in each
// iteration before calling balance()
if(MyRank==0)
{
printf("Starting step: %d/%d\n", CurrentStep+1, NumberOfSteps);
fflush(stdout);
}
#ifdef ALL_VTK_OUTPUT_EXAMPLE
try {
jall->printVTKoutlines(CurrentStep);
} catch (ALL::FilesystemErrorException &e) {
std::cout << e.what() << std::endl;
std::exit(2);
// Maybe also treat this error in some way, e.g. whether the simulation
// should abort if no output could be written or not.
}
#endif
jall->balance();
std::vector<ALL::Point<double>> NewVertices = jall->getVertices();
//MPI_RUN_ORDER(MPI_COMM_WORLD, MyRank, MaximumRank, (print_width(MyRank, NewVertices.at(1)[2]-NewVertices.at(0)[2], NewVertices.at(0)[2], NewVertices.at(1)[2])));
MPI_RUN_ORDER_DEF((print_testing_output(MyRank, NewVertices, CurrentStep+1)));
// Maybe print our new domain? Or not..
//convert_verts(&NewVertices, VertexArray);
//MPI_RUN_ORDER_DEF((print_domain(MyRank, VertexArray)));
//jall->getWork(MyWork);
//MPI_RUN_ORDER_DEF((print_work(MyRank,MyWork)));
MPI_Barrier(MPI_COMM_WORLD);
}
#ifdef ALL_VTK_OUTPUT_EXAMPLE
try {
jall->printVTKoutlines(CurrentStep);
} catch (ALL::FilesystemErrorException &e) {
std::cout << e.what() << std::endl;
}
#endif
delete jall;
MPI_Finalize();
return EXIT_SUCCESS;
}
!Copyright 2018-2020 Rene Halver, Forschungszentrum Juelich GmbH, Germany
!Copyright 2018-2020 Godehard Sutmann, Forschungszentrum Juelich GmbH, Germany
!Copyright 2020-2020 Stephan Schulz, Forschungszentrum Juelich GmbH, Germany
!
!Redistribution and use in source and binary forms, with or without modification, are
!permitted provided that the following conditions are met:
!
!1. Redistributions of source code must retain the above copyright notice, this list
! of conditions and the following disclaimer.
!
!2. Redistributions in binary form must reproduce the above copyright notice, this
! list of conditions and the following disclaimer in the documentation and/or
! other materials provided with the distribution.
!
!3. Neither the name of the copyright holder nor the names of its contributors
! may be used to endorse or promote products derived from this software without
! specific prior written permission.
!
!THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
!EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
!OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
!SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
!INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
!TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
!BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
!CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
!ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
!SUCH DAMAGE.
program ALL_Staggered_f
use ALL_module
use iso_c_binding
use, intrinsic :: iso_fortran_env, only: stdin=>input_unit&
, stdout=>output_unit&
, stderr=>error_unit&
, file_storage_size
#ifndef ALL_USE_F08
use mpi
#else
use mpi_f08
#endif
implicit none
integer, parameter :: real64 = selected_real_kind(15)
call main()
contains
subroutine print_loc(my_rank, maximum_rank, my_location, number_of_processes)
integer, intent(in) :: my_rank, maximum_rank
integer, dimension(3), intent(in) :: my_location, number_of_processes
integer :: rank, error
do rank=0, maximum_rank-1
if(rank == my_rank) then
write(stdout,'(a,i3.3,a,i3,",",i3,",",i3,a,i3,",",i3,",",i3,a)')&
"[", my_rank, "] Location: (", my_location, ")/(", number_of_processes, ")"
flush(stdout)
call MPI_Barrier(MPI_COMM_WORLD, error)
else
call MPI_Barrier(MPI_COMM_WORLD, error)
endif
enddo
end subroutine
subroutine print_domain(my_rank, maximum_rank, domain_vertices)
integer, intent(in) :: my_rank, maximum_rank
real(real64), dimension(3,2), intent(in) :: domain_vertices
integer :: rank, error
do rank=0, maximum_rank-1
if(rank == my_rank) then
write(stdout,'("[",i3.3,"]",a,es12.4,a,es12.4,a,es12.4)')&
my_rank, " Lower: ",&
domain_vertices(1,1), achar(9),&
domain_vertices(2,1), achar(9),&
domain_vertices(3,1)
write(stdout,'("[",i3.3,"]",a,es12.4,a,es12.4,a,es12.4)')&
my_rank, " Upper: ",&
domain_vertices(1,2), achar(9),&
domain_vertices(2,2), achar(9),&
domain_vertices(3,2)
flush(stdout)
call MPI_Barrier(MPI_COMM_WORLD,error)
else
call MPI_Barrier(MPI_COMM_WORLD,error)
endif
enddo
end subroutine
subroutine print_work(my_rank, maximum_rank, my_work)
integer, intent(in) :: my_rank, maximum_rank
real(real64), intent(in) :: my_work
integer :: rank, error
do rank=0, maximum_rank-1
if(rank == my_rank) then
write(stdout,'(a,i3.3,a,es12.4)')&
"[", my_rank, "] Work: ", my_work
flush(stdout)
call MPI_Barrier(MPI_COMM_WORLD, error)
else
call MPI_Barrier(MPI_COMM_WORLD, error)
endif
enddo
end subroutine
subroutine print_testing_output(my_rank, maximum_rank, new_vertices, timestep)
integer, intent(in) :: my_rank, maximum_rank, timestep
real(real64), dimension(3,2), intent(in) :: new_vertices
integer :: rank, error
do rank=0, maximum_rank-1
if(rank == my_rank) then
!write(stdout,'(a,i4,",",i3.3,a,3f11.6)')&
! "[", timestep, my_rank, "] Result Width: ",&
! new_vertices(:,2)-new_vertices(:,1)
!flush(stdout)
write(stdout,'(a,i4,",",i3.3,",",i2.2,a,3f11.6)')&
"[", timestep, my_rank, 0, "] Result Vertex: ",&
new_vertices(:,1)
flush(stdout)
write(stdout,'(a,i4,",",i3.3,",",i2.2,a,3f11.6)')&
"[", timestep, my_rank, 1, "] Result Vertex: ",&
new_vertices(:,2)
flush(stdout)
call MPI_Barrier(MPI_COMM_WORLD, error)
else
call MPI_Barrier(MPI_COMM_WORLD, error)
endif
enddo
end subroutine
subroutine main()
integer :: error
integer :: current_step
integer, parameter :: number_of_steps = 50
integer, parameter :: dimensions = 3
real(real64), parameter :: loadbalancer_gamma = 0. ! ignored for tensor method
integer, dimension(dimensions) :: my_location, number_of_processes
real(real64), dimension(dimensions) :: minimum_domain_size
real(real64), dimension(dimensions,2) :: domain_vertices, new_vertices
real(real64), parameter :: domain_size = 1.0
real(real64) :: my_work
integer :: my_rank, maximum_rank
integer :: i,j
type(ALL_t) :: jall
#ifdef ALL_VTK_OUTPUT_EXAMPLE
character (len=ALL_ERROR_LENGTH) :: error_msg
#endif
call MPI_Init(error)
call jall%init(ALL_TENSOR, dimensions, loadbalancer_gamma)
my_location(:) = 0
! All domains are placed along a line in z direction, evne though they are three dimensional
call MPI_Comm_rank(MPI_COMM_WORLD, my_location(3), error)
my_rank = my_location(3)
number_of_processes(:) = 1
call MPI_Comm_size(MPI_COMM_WORLD, number_of_processes(3), error)
maximum_rank = number_of_processes(3)
if(my_rank == 0) then
write(stdout,'(a,i0)') "Ranks: ", maximum_rank
write(stdout,'(a,i0)') "Number of Steps: ", number_of_steps
flush(stdout)
endif
call MPI_Barrier(MPI_COMM_WORLD, error)
call print_loc(my_rank, maximum_rank, my_location, number_of_processes)
! For a cartesian communicator this is not required, but we are using
! MPI_COMM_WORLD here.
call jall%set_proc_grid_params(my_location, number_of_processes)
! If we want ot set a minimum domain size:
minimum_domain_size(:) = 0.1
call jall%set_min_domain_size(minimum_domain_size)
call jall%set_communicator(MPI_COMM_WORLD)
! We also set the optional process tag for the output.
! This can be useful if we want to know which of 'our'
! ranks is which in the output produces by the library.
! The ranks used inside the library do not necessarily
! match our own numbering.
call jall%set_proc_tag(my_rank)
call jall%setup()
! A first domain distribution must be given to the balancer.
! We use the provided ALL::Point class to define the vertices,
! but a simple double array can also be used. We need 2 vertices
! which correspond to the two opposing corners.
! We create a cubic domain initially
do i=1, 2
do j=1,dimensions
domain_vertices(j,i) = (my_location(j)+i-1) * domain_size
enddo
enddo
call print_domain(my_rank, maximum_rank, domain_vertices)
call jall%set_vertices(domain_vertices)
! Calculate the work fo our domain. Here we just use
my_work = my_rank + 1.
call jall%set_work(my_work)
call print_work(my_rank, maximum_rank, my_work)
do current_step=1, number_of_steps
! In a real code we need to set the updated work in each
! iteration before calling balance()
if(my_rank == 0) then
write(stdout,'(a,i0,"/",i0)') "Starting step: ", current_step, number_of_steps
flush(stdout)
endif
#ifdef ALL_VTK_OUTPUT_EXAMPLE
call ALL_reset_error()
call jall%print_vtk_outlines(current_step)
if(ALL_error() /= 0) then
error_msg = ALL_error_description()
write(stdout, '(a)') error_msg
! Maybe also abort if the output is necesssary or handle this
! some other way.
call MPI_Abort(MPI_COMM_WORLD, 2, error)
stop
endif
#endif
call jall%balance()
call jall%get_vertices(new_vertices)
call print_testing_output(my_rank, maximum_rank, new_vertices, current_step)
call MPI_Barrier(MPI_COMM_WORLD, error)
enddo
#ifdef ALL_VTK_OUTPUT_EXAMPLE
call ALL_reset_error()
call jall%print_vtk_outlines(current_step)
if(ALL_error() /= 0) then
error_msg = ALL_error_description()
write(stdout, '(a)') error_msg
! Maybe also abort if the output is necesssary or handle this
! some other way.
call MPI_Abort(MPI_COMM_WORLD, 2, error)
stop
endif
#endif
call jall%finalize()
call MPI_Finalize(error)
end subroutine
end program
...@@ -54,6 +54,18 @@ install(TARGETS ALL_Staggered ...@@ -54,6 +54,18 @@ install(TARGETS ALL_Staggered
LIBRARY DESTINATION lib LIBRARY DESTINATION lib
INCLUDES DESTINATION include) INCLUDES DESTINATION include)
add_executable(ALL_Tensor ALL_Tensor.cpp)
if(CM_ALL_VTK_OUTPUT)
target_compile_definitions(ALL_Tensor PRIVATE ALL_VTK_OUTPUT_EXAMPLE)
target_include_directories(ALL_Tensor PUBLIC ${VTK_INCLUDE_DIRS})
target_link_libraries(ALL_Tensor PUBLIC ${VTK_LIBRARY_DIRS})
endif(CM_ALL_VTK_OUTPUT)
target_link_libraries (ALL_Tensor LINK_PUBLIC ALL)
install(TARGETS ALL_Tensor
RUNTIME DESTINATION bin
LIBRARY DESTINATION lib
INCLUDES DESTINATION include)
if(CM_ALL_VORONOI) if(CM_ALL_VORONOI)
add_executable(ALL_Voronoi ALL_Voronoi.cpp) add_executable(ALL_Voronoi ALL_Voronoi.cpp)
if(CM_ALL_VTK_OUTPUT) if(CM_ALL_VTK_OUTPUT)
...@@ -125,4 +137,20 @@ if(CM_ALL_FORTRAN) ...@@ -125,4 +137,20 @@ if(CM_ALL_FORTRAN)
RUNTIME DESTINATION bin RUNTIME DESTINATION bin
LIBRARY DESTINATION lib LIBRARY DESTINATION lib
INCLUDES DESTINATION include) INCLUDES DESTINATION include)
add_executable(ALL_Tensor_f ALL_Tensor_f.F90)
if(CM_ALL_VTK_OUTPUT)
target_compile_definitions(ALL_Tensor_f PRIVATE ALL_VTK_OUTPUT_EXAMPLE)
target_include_directories(ALL_Tensor_f PUBLIC ${VTK_INCLUDE_DIRS})
target_link_libraries(ALL_Tensor_f PUBLIC ${VTK_LIBRARY_DIRS})
endif(CM_ALL_VTK_OUTPUT)
target_link_libraries(ALL_Tensor_f LINK_PUBLIC ALL)
target_link_libraries(ALL_Tensor_f LINK_PUBLIC ALL_fortran)
target_include_directories(ALL_Tensor_f PUBLIC ${CMAKE_BINARY_DIR}/modules)
list(APPEND CMAKE_MODULE_PATH ${CMAKE_BINARY_DIR}/modules)
set_property(TARGET ALL_Tensor_f PROPERTY LINKER_LANGUAGE Fortran)
install(TARGETS ALL_Tensor_f
RUNTIME DESTINATION bin
LIBRARY DESTINATION lib
INCLUDES DESTINATION include)
endif(CM_ALL_FORTRAN) endif(CM_ALL_FORTRAN)
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment