mpi_assembly_begin.h
Go to the documentation of this file.
#ifndef _RHEO_MPI_ASSEMBLY_BEGIN_H
#define _RHEO_MPI_ASSEMBLY_BEGIN_H
namespace rheolef {
/*F:
NAME: mpi_assembly_begin -- for array or matrix (@PACKAGE@ @VERSION@)
DESCRIPTION:
Start a dense array or a sparse matrix assembly.
COMPLEXITY:
Assume that stash has indexes in increasing order.
cpu complexity : O(stash.size + nproc)
memory complexity : O(stash.size + nproc)
msg size complexity : O(stash.size + nproc)
When assembling a finite element matrix, the stash.size is at boundaries
of the mesh partition, and stash.size = O(nrow^alpha), where alpha=(d-1)/d
and d=2,3. Using nproc <= O(nrow^alpha) is performant.
NOTE:
The stash may be sorted by increasing nows and column.
Inspirated from Petsc (petsc/src/vec/vec/impls/mpi/pdvec.c), here with
a pre-sorted stash, thanks to the stl::map data structure.
AUTHORS:
LMC-IMAG, 38041 Grenoble cedex 9, France
| Pierre.Saramito@imag.fr
DATE: 23 march 1999
END:
*/
//<mpi_assembly_begin:
template <
class Stash,
class Message,
class InputIterator>
typename Stash::size_type
mpi_assembly_begin (
// input:
const Stash& stash,
InputIterator first_stash_idx, // wrapper in shash
InputIterator last_stash_idx,
const distributor& ownership,
// ouput:
Message& receive, // buffer
Message& send) // buffer
{
typedef typename Stash::size_type size_type;
mpi::communicator comm = ownership.comm();
size_type my_proc = ownership.process();
size_type nproc = ownership.n_process();
distributor::tag_type tag = distributor::get_new_tag();
// ----------------------------------------------------------------
// 1) count the messages contained in stash by process id
// ----------------------------------------------------------------
// assume that stash elements are sorted by increasing stash_idx (e.g. stash = stl::map)
// cpu complexity = O(stash.size + nproc)
// mem complexity = O(stash.size + nproc)
std::vector<size_type> msg_size(nproc, 0);
std::vector<size_type> msg_mark(nproc, 0);
size_type send_nproc = 0;
{
size_type iproc = 0;
size_type i = 0;
for (InputIterator iter_idx = first_stash_idx; iter_idx != last_stash_idx; iter_idx++, i++) {
for (; iproc < nproc; iproc++) {
if (*iter_idx >= ownership[iproc] && *iter_idx < ownership[iproc+1]) {
msg_size [iproc]++;
if (!msg_mark[iproc]) {
msg_mark[iproc] = 1;
send_nproc++;
}
break;
}
}
assert_macro (iproc != nproc, "bad stash data: index "<<*iter_idx<<" out of range [0:"<<ownership[nproc]<<"[");
}
} // end block
// ----------------------------------------------------------------
// 2) avoid to send message to my-proc in counting
// ----------------------------------------------------------------
if (msg_size [my_proc] != 0) {
msg_size [my_proc] = 0;
msg_mark [my_proc] = 0;
send_nproc--;
}
// ----------------------------------------------------------------
// 3) compute number of messages to be send to my_proc
// ----------------------------------------------------------------
// msg complexity : O(nproc) or O(log(nproc)), depending on reduce implementation
std::vector<size_type> work (nproc);
mpi::all_reduce (
comm,
msg_mark.begin().operator->(),
nproc,
work.begin().operator->(),
std::plus<size_type>());
size_type receive_nproc = work [my_proc];
// ----------------------------------------------------------------
// 4) compute messages max size to be send to my_proc
// ----------------------------------------------------------------
// msg complexity : O(nproc) or O(log(nproc)), depending on reduce implementation
mpi::all_reduce (
comm,
msg_size.begin().operator->(),
nproc,
work.begin().operator->(),
mpi::maximum<size_type>());
size_type receive_max_size = work [my_proc];
// ----------------------------------------------------------------
// 5) post receive: exchange the buffer adresses between processes
// ----------------------------------------------------------------
// Each message will consist of ordered pairs (global index,value).
// since we don't know how long each indiidual message is,
// we allocate the largest : receive_nproc*receive_max_size
// potentially, this is a lot of wasted space
// TODO: how to optimize the receive.data buffer ?
// cpu complexity : O(nproc)
// mem complexity : O(nproc*(stash.size/nproc)) = O(stash.size), worst case ?
// msg complexity : O(nproc)
receive.waits.clear();
receive.data.resize (receive_nproc*receive_max_size);
for (size_t i_receive = 0; i_receive < receive_nproc; i_receive++) {
mpi::request i_req = comm.irecv (
mpi::any_source,
tag,
receive.data.begin().operator->() + i_receive*receive_max_size,
receive_max_size);
receive.waits.push_back (std::make_pair(i_receive, i_req));
}
// ----------------------------------------------------------------
// 6) copy stash in send buffer
// ----------------------------------------------------------------
// since the stash is sorted by increasing order => simple copy
// cpu complexity : O(stash.size)
// mem complexity : O(stash.size)
send.data.resize (stash.size());
copy (stash.begin(), stash.end(), send.data.begin());
// ---------------------------------------------------------------------------
// 7) do send
// ---------------------------------------------------------------------------
// cpu complexity : O(nproc)
// mem complexity : O(send_nproc) \approx O(nproc), worst case
// msg complexity : O(stash.size)
send.waits.clear();
{
size_type i_send = 0;
size_type i_start = 0;
for (size_type iproc = 0; iproc < nproc; iproc++) {
size_type i_msg_size = msg_size[iproc];
if (i_msg_size == 0) continue;
mpi::request i_req = comm.isend (
iproc,
tag,
send.data.begin().operator->() + i_start,
i_msg_size);
send.waits.push_back(std::make_pair(i_send,i_req));
i_send++;
i_start += i_msg_size;
}
} // end block
return receive_max_size;
}
//>mpi_assembly_begin:
} // namespace rheolef
#endif //_RHEO_MPI_ASSEMBLY_BEGIN_H
size_type
field::size_type size_type
Definition: branch.cc:425
rheolef::distributor
see the distributor page for the full documentation
Definition: distributor.h:62
rheolef::distributor::n_process
size_type n_process() const
number of processes
Definition: distributor.h:169
rheolef::distributor::get_new_tag
static tag_type get_new_tag()
returns a new tag
Definition: distributor.cc:133
rheolef::distributor::process
size_type process() const
current process id
Definition: distributor.h:179
rheolef::distributor::tag_type
int tag_type
Definition: distributor.h:71
rheolef::distributor::comm
const communicator_type & comm() const
Definition: distributor.h:145
rheolef::size_type
size_t size_type
Definition: basis_get.cc:76
assert_macro
#define assert_macro(ok_condition, message)
Definition: dis_macros.h:113
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
rheolef::mpi_assembly_begin
Stash::size_type mpi_assembly_begin(const Stash &stash, InputIterator first_stash_idx, InputIterator last_stash_idx, const distributor &ownership, Message &receive, Message &send)
Definition: mpi_assembly_begin.h:57