26 #ifndef _SF_PARALLEL_UTILS_H 27 #define _SF_PARALLEL_UTILS_H 55 MPI_Comm_size(comm, &size); MPI_Comm_rank(comm, &rank);
62 T bsize = (gmax - gmin) / size + 1;
69 for(
size_t i=0; i<dest.size(); i++)
70 dest[i] = (idx[i] - gmin) / bsize;
77 for(
size_t i=0; i<perm.size(); i++)
78 snd_idx[i] = idx[perm[i]];
105 template<
class T,
class V>
110 MPI_Comm_size(comm, &size); MPI_Comm_rank(comm, &rank);
117 T bsize = (gmax - gmin) / size + 1;
124 for(
size_t i=0; i<dest.size(); i++)
125 dest[i] = (idx[i] - gmin) / bsize;
133 for(
size_t i=0; i<perm.size(); i++) {
134 snd_idx[i] = idx[perm[i]];
135 snd_val[i] = val[perm[i]];
153 template<
class T,
class V>
158 MPI_Comm_size(comm, &size); MPI_Comm_rank(comm, &rank);
165 T bsize = (gmax - gmin) / size + 1;
173 for(
size_t i=0; i<dest.size(); i++)
174 dest[i] = (idx[i] - gmin) / bsize;
184 for(
size_t i=0; i<perm.size(); i++) {
185 snd_idx[i] = idx[perm[i]];
186 snd_cnt[i] = cnt[perm[i]];
190 for(
size_t i=0; i<perm.size(); i++) {
191 const V* read = val.
data() + dsp[perm[i]];
192 V* write = snd_val.data() + snd_dsp[i];
194 for(T j=0; j<snd_cnt[i]; j++)
202 grph_entr.
configure(dest, snd_cnt, comm);
205 vector<T> rec_cnt(rsize), rec_dsp, out_dsp;
208 out_val.
resize(rsize_entr);
221 for(
size_t i=0; i<perm.size(); i++)
222 out_cnt[i] = rec_cnt[perm[i]];
226 for(
size_t i=0; i<perm.size(); i++) {
227 const V* read = rec_val.
data() + rec_dsp[perm[i]];
228 V* write = out_val.
data() + out_dsp[i];
230 for(T j=0; j<out_cnt[i]; j++)
253 MPI_Comm_rank(comm, &rank); MPI_Comm_size(comm, &size);
254 long int lsize = vec.
size();
261 nwr += fwrite(vec.
data(),
sizeof(V), vec.
size(), fd);
267 for(
int pid=1; pid < size; pid++)
270 MPI_Send(&lsize, 1, MPI_LONG, 0,
SF_MPITAG, comm);
271 MPI_Send(vec.
data(), lsize*
sizeof(V), MPI_BYTE, 0,
SF_MPITAG, comm);
273 else if (rank == 0) {
277 MPI_Recv(&rsize, 1, MPI_LONG, pid,
SF_MPITAG, comm, &stat);
280 MPI_Recv(wbuff.data(), rsize*
sizeof(V), MPI_BYTE, pid,
SF_MPITAG, comm, &stat);
282 nwr += fwrite(wbuff.data(),
sizeof(V), rsize, fd);
288 MPI_Bcast(&nwr, 1, MPI_LONG, 0, comm);
296 size_t root_write(FILE* fd, V* vec,
const size_t vec_size, MPI_Comm comm)
299 vecbuff.
assign(vec_size, vec,
false);
303 vecbuff.
assign(0, NULL,
false);
325 MPI_Comm_rank(comm, &rank); MPI_Comm_size(comm, &size);
326 long int lsize = vec.
size();
334 nrd += fread(vec.
data(),
sizeof(V), vec.
size(), fd);
338 for(
int pid=1; pid < size; pid++)
342 MPI_Recv(&rsize, 1, MPI_LONG, pid,
SF_MPITAG, comm, &stat);
345 nrd += fread(rbuff.data(),
sizeof(V), rsize, fd);
347 MPI_Send(rbuff.data(), rsize*
sizeof(V), MPI_BYTE, pid,
SF_MPITAG, comm);
351 MPI_Send(&lsize, 1, MPI_LONG, 0,
SF_MPITAG, comm);
354 MPI_Recv(vec.
data(), lsize*
sizeof(V), MPI_BYTE, 0,
SF_MPITAG, comm, &stat);
357 MPI_Bcast(&nrd, 1, MPI_LONG, 0, comm);
379 MPI_Comm_rank(comm, &rank); MPI_Comm_size(comm, &size);
380 long int lsize = vec.
size();
392 for(
size_t i=0; i<vec.
size(); i++) {
393 nrd += fscanf(fd,
"%ld", &ibuff);
398 for(
size_t i=0; i<vec.
size(); i++) {
399 nrd += fscanf(fd,
"%lf", &fbuff);
407 for(
int pid=1; pid < size; pid++)
411 MPI_Recv(&rsize, 1, MPI_LONG, pid,
SF_MPITAG, comm, &stat);
414 for(
long int i=0; i<rsize; i++) {
416 nrd += fscanf(fd,
"%ld", &ibuff);
420 nrd += fscanf(fd,
"%lf", &fbuff);
425 MPI_Send(rbuff.
data(), rsize*
sizeof(V), MPI_BYTE, pid,
SF_MPITAG, comm);
429 MPI_Send(&lsize, 1, MPI_LONG, 0,
SF_MPITAG, comm);
432 MPI_Recv(vec.
data(), lsize*
sizeof(V), MPI_BYTE, 0,
SF_MPITAG, comm, &stat);
435 MPI_Bcast(&nrd, 1, MPI_LONG, 0, comm);
443 MPI_Comm_size(comm, &size); MPI_Comm_rank(comm, &rank);
445 size_t line_count = 0;
447 FILE* fd = fopen(file.c_str(),
"r");
453 divide(fsize, size, chunk_sizes);
456 for(
int chunk_size : chunk_sizes) {
458 fread(chunk.
data(), chunk_size, 1, fd);
460 for(
unsigned char c : chunk)
461 if(c ==
'\n') line_count++;
466 fprintf(stderr,
"%s error: Cannot open file %s!\n", __func__, file.c_str());
470 MPI_Bcast(&line_count,
sizeof(
size_t), MPI_BYTE, 0, comm);
478 size_t root_read(FILE* fd, V* vec,
const size_t vec_size, MPI_Comm comm)
481 vecbuff.
assign(vec_size, vec,
false);
483 size_t nrd =
root_read(fd, vecbuff, comm);
485 vecbuff.
assign(0, NULL,
false);
494 size_t root_read_ascii(FILE* fd, V* vec,
const size_t vec_size, MPI_Comm comm,
bool int_type)
497 vecbuff.
assign(vec_size, vec,
false);
501 vecbuff.
assign(0, NULL,
false);
521 template<
class T,
class V>
547 template<
class T,
class V>
554 sort_parallel(comm, idx, cnt, vec, srt_idx, srt_cnt, srt_vec);
562 template<
class T,
class V>
568 idxbuff.
assign(vec_size, idx,
false);
569 vecbuff.
assign(vec_size, vec,
false);
573 idxbuff.
assign(0, NULL,
false);
574 vecbuff.
assign(0, NULL,
false);
581 template<
class T,
class V>
583 const size_t idx_size,
const size_t vec_size, MPI_Comm comm)
588 idxbuff.
assign(idx_size, idx,
false);
589 cntbuff.
assign(idx_size, cnt,
false);
590 vecbuff.
assign(vec_size, vec,
false);
594 idxbuff.
assign(0, NULL,
false);
595 cntbuff.
assign(0, NULL,
false);
596 vecbuff.
assign(0, NULL,
false);
606 MPI_Comm_size(comm, &size), MPI_Comm_rank(comm, &rank);
609 for (
size_t i=0; i<vec.
size() / dpn; i++ ) {
610 for(
short j=0; j<dpn-1; j++)
611 fprintf(fd,
"%g ",
double(vec[i*dpn + j]) );
612 fprintf(fd,
"%g\n",
double(vec[i*dpn + (dpn-1)]) );
617 for(
int pid=1; pid < size; pid++)
621 MPI_Recv(&rsize, 1, MPI_LONG, pid,
SF_MPITAG, comm, &stat);
623 MPI_Recv(wbuff.
data(), rsize*
sizeof(T), MPI_BYTE, pid,
SF_MPITAG, comm, &stat);
625 for (
size_t i=0; i<wbuff.
size() / dpn; i++ ) {
626 for(
short j=0; j<dpn-1; j++)
627 fprintf(fd,
"%g ",
double(wbuff[i*dpn + j]) );
629 fprintf(fd,
"%g\n",
double(wbuff[i*dpn + (dpn-1)]) );
633 long int lsize = vec.
size();
634 MPI_Send(&lsize, 1, MPI_LONG, 0,
SF_MPITAG, comm);
635 MPI_Send(vec.
data(), lsize*
sizeof(T), MPI_BYTE, 0,
SF_MPITAG, comm);
640 template<
class T,
class S>
642 std::string file,
short dpn = 1)
647 MPI_Comm_rank(comm, &rank);
655 fd = fopen(file.c_str(),
"w");
657 fprintf(stderr,
"%s error: Cannot open file %s for writing! Aborting!\n", __func__, file.c_str());
T global_min(const vector< T > &vec, MPI_Comm comm)
Compute the global minimum of a distributed vector.
The vector class and related algorithms.
size_t root_read_ascii(FILE *fd, vector< V > &vec, MPI_Comm comm, bool int_data)
Read binary data into a vector.
void MPI_Exchange(commgraph< T > &grph, vector< S > &send, vector< S > &recv, MPI_Comm comm)
Exchange data in parallel over MPI.
void sort_parallel(MPI_Comm comm, const vector< T > &idx, vector< T > &out_idx)
Sort index values parallel ascending across the ranks.
void interval(vector< T > &vec, size_t start, size_t end)
Create an integer interval between start and end.
size_t file_size(FILE *fd)
return file size from a file descriptor
void write_data_ascii(const MPI_Comm comm, const vector< T > &idx, const vector< S > &data, std::string file, short dpn=1)
void divide(const size_t gsize, const size_t num_parts, vector< T > &loc_sizes)
divide gsize into num_parts local parts with even distribution of the remainder
#define SF_MPITAG
the MPI tag when communicating
T * data()
Pointer to the vector's start.
size_t root_write_ordered(FILE *fd, const vector< T > &idx, const vector< V > &vec, MPI_Comm comm)
Write index value pairs to disk in ordered permutation.
Functions related to mesh IO.
size_t root_read(FILE *fd, vector< V > &vec, MPI_Comm comm)
Read binary data into a vector.
void configure(const vector< V > &dest, MPI_Comm comm)
Set up the communication graph.
T sum(const vector< T > &vec)
Compute sum of a vector's entries.
size_t root_count_ascii_lines(std::string file, MPI_Comm comm)
count the lines in a ascii file
void binary_sort_copy(vector< T > &_V, vector< S > &_W)
T global_max(const vector< T > &vec, MPI_Comm comm)
Compute the global maximum of a distributed vector.
size_t size() const
The current size of the vector.
A vector storing arbitrary data.
void assign(InputIterator s, InputIterator e)
Assign a memory range.
vector< T > rcnt
Number of elements received from each rank.
void binary_sort(vector< T > &_V)
size_t root_write(FILE *fd, const vector< V > &vec, MPI_Comm comm)
Write vector data binary to disk.
void resize(size_t n)
Resize a vector.
void print_vector(MPI_Comm comm, const vector< T > &vec, const short dpn, FILE *fd)
void dsp_from_cnt(const vector< T > &cnt, vector< T > &dsp)
Compute displacements from counts.