14 #define _XOPEN_SOURCE 700
25 #include <sys/types.h>
32 #include "sion_error_handler.h"
44 int _sion_errorprint_mpi(
int rc,
int level,
const char *format, ...)
46 int rank=-1, thread=-1;
50 MPI_Comm_rank(MPI_COMM_WORLD, &rank);
54 rc=__sion_errorprint_vargs(rc, level, rank, thread, format, ap);
76 int rc = SION_SUCCESS;
80 MPI_Comm_size(gComm, >asks);
81 MPI_Comm_rank(gComm, &gRank);
83 if (gtasks < numFiles) {
84 return(_sion_errorprint_mpi(SION_NOT_SUCCESS,_SION_ERROR_RETURN,
85 "_sion_gen_info_from_gcomm_mpi: Number of files(%d) cannot be bigger the the number of tasks(%d), aborting...\n",
89 task_per_file = gtasks / numFiles;
92 if (gRank >= (numFiles * task_per_file)) {
93 *filenumber = numFiles - 1;
96 *filenumber = gRank / task_per_file;
99 if(*filenumber == numFiles - 1) {
100 *lrank=gRank-(numFiles - 1)*task_per_file;
101 *lsize=gtasks-(numFiles - 1)*task_per_file;
103 *lrank=gRank%task_per_file;
104 *lsize=task_per_file;
107 DPRINTFP((1,
"_sion_get_info_from_splitted_comm_mpi", gRank,
"Global communicator divided in %d local communicators (%f: %d of %d)\n",
108 numFiles,*filenumber,*lrank,*lsize));
113 int _sion_get_info_from_splitted_comm_mpi(MPI_Comm gComm, MPI_Comm lComm,
int *numComm,
int *CommNumber,
int *lrank,
int *lsize) {
114 int gSize, gRank, lSize, lRank;
115 int *allRanks=NULL, *allSizes=NULL, i, ntasks;
116 int ncomms, color = 0;
117 int rc = SION_SUCCESS;
120 MPI_Comm_size(gComm, &gSize);
121 MPI_Comm_rank(gComm, &gRank);
123 if (lComm != MPI_COMM_NULL) {
124 MPI_Comm_size(lComm, &lSize);
125 MPI_Comm_rank(lComm, &lRank);
132 DPRINTFP((32,
"_sion_get_info_from_splitted_comm_mpi", gRank,
"lRank: %d\n", lRank));
136 allRanks = (
int *) malloc(gSize *
sizeof(
int));
137 allSizes = (
int *) malloc(gSize *
sizeof(
int));
138 if ((allRanks == NULL) || (allSizes == NULL)) {
141 return(_sion_errorprint_mpi(SION_NOT_SUCCESS,_SION_ERROR_RETURN,
"_sion_get_info_from_splitted_comm_mpi: Cannot allocate temp arrays of size %lu\n", (
unsigned long) gSize *
sizeof(
int)));
146 MPI_Gather(&lRank, 1, MPI_INT, allRanks, 1, MPI_INT, 0, gComm);
147 MPI_Gather(&lSize, 1, MPI_INT, allSizes, 1, MPI_INT, 0, gComm);
152 for (i = 0; i < gSize; i++) {
153 if (allRanks[i] == 0) {
156 MPI_Send(&ncomms, 1, MPI_INT, i, i, gComm);
157 DPRINTFP((32,
"sion_paropen_comms_mpi", gRank,
"Sent ncomms=%05d to %d with TAG %d\n", ncomms, i, i));
164 ntasks += allSizes[i];
169 if ((lRank == 0) && (gRank != 0)) {
170 MPI_Recv(&color, 1, MPI_INT, 0, gRank, gComm, &status);
171 DPRINTFP((32,
"sion_paropen_comms_mpi", gRank,
"Received ncomms=%05d from %d with TAG %d\n", color, status.MPI_SOURCE, status.MPI_TAG));
175 MPI_Bcast(&ncomms, 1, MPI_INT, 0, gComm);
176 MPI_Bcast(&ntasks, 1, MPI_INT, 0, gComm);
179 DPRINTFP((32,
"_sion_get_info_from_splitted_comm_mpi", gRank,
"#Comms=%05d #Tasks=%05d #TotalTasks=%05d\n", ncomms, ntasks, gSize));
181 if (lComm != MPI_COMM_NULL) {
182 MPI_Bcast(&color, 1, MPI_INT, 0, lComm);
186 if(allRanks) free(allRanks);
187 if(allSizes) free(allSizes);
198 int _sion_get_numfiles_from_file_mpi(
char *fname) {
203 sion_int64 chunksize;
204 sion_int32 fsblksize;
207 sid = _sion_open_rank(fname,
"br", &chunksize, &fsblksize, &rank, &fileptr);
209 numfiles = sion_get_number_of_files(sid);
211 _sion_close_sid(sid);
213 DPRINTFP((1,
"_sion_get_numfiles_from_file_mpi", _SION_DEFAULT_RANK,
"sion file %s has %d files\n", fname, numfiles));
218 int _sion_get_filenumber_from_files_mpi(
char *fname, MPI_Comm gComm,
int *numfiles,
int *filenumber,
int *lRank) {
220 int sid, gSize, gRank, ntasks, nfiles;
221 int rc = SION_SUCCESS;
223 sion_int64 *chunksizes;
224 sion_int32 fsblksize;
226 int mapping_size = 0;
227 sion_int32 *mapping = NULL;
230 MPI_Comm_size(gComm, &gSize);
231 MPI_Comm_rank(gComm, &gRank);
233 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"before open\n"));
236 chunksizes=NULL;globalranks=NULL;
237 sid=
_sion_open_read(fname,_SION_FMODE_READ|_SION_FMODE_ANSI,_SION_READ_MASTER_ONLY_OF_MULTI_FILES,
238 &ntasks,&nfiles,&chunksizes,&fsblksize,&globalranks,&fileptr);
241 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"after open\n"));
243 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"sion file %d files rc=%d\n", *numfiles, rc));
250 MPI_Bcast(numfiles, 1, MPI_INT, 0, gComm);
252 if((gRank == 0) && (*numfiles>1)) {
253 if(mapping_size!=gSize) {
254 return(_sion_errorprint_mpi(SION_NOT_SUCCESS,_SION_ERROR_RETURN,
"_sion_get_filenumber_from_files_mpi: Incorrect sum of ntasks of files %d <> %d\n", mapping_size, gSize));
259 return(_sion_errorprint_mpi(SION_NOT_SUCCESS,_SION_ERROR_RETURN,
"_sion_get_filenumber_from_files_mpi: could not get numfiles from sion file\n"));
263 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"before scatter\n"));
264 MPI_Scatter(mapping, 2, SION_MPI_INT32, lpos, 2, SION_MPI_INT32, 0, gComm);
267 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"after scatter\n"));
271 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"only one file -> filenumber=%d lRank=%d\n",*filenumber,*lRank));
276 if (sid>=0) _sion_close_sid(sid);
279 DPRINTFP((1,
"_sion_get_filenumber_from_files_mpi", gRank,
"close rc=%d\n",rc));
int sion_get_mapping(int sid, int *mapping_size, sion_int32 **mapping, int *numfiles)
Returns pointers to the internal field mapping.
int _sion_open_read(const char *fname, sion_int64 file_mode_flags, int read_all, int *ntasks, int *nfiles, sion_int64 **chunksizes, sion_int32 *fsblksize, int **globalranks, FILE **fileptr)
internal sion serial open function for reading on one or more files
int _sion_gen_info_from_gcomm_mpi(int numFiles, MPI_Comm gComm, int *filenumber, int *lrank, int *lsize)
Splits a Communicator in numfiles different communicators.