46 #ifndef MUELU_PERFUTILS_DEF_HPP
47 #define MUELU_PERFUTILS_DEF_HPP
53 #include <Teuchos_CommHelpers.hpp>
59 #include <Xpetra_CrsMatrixWrap.hpp>
68 void calculateStats(Type& minVal, Type& maxVal,
double& avgVal,
double& devVal,
int& minProc,
int& maxProc,
const RCP<
const Teuchos::Comm<int> >& comm,
int numActiveProcs,
const Type& v) {
69 Type sumVal, sum2Val, v2 = v * v;
79 w = (minVal == v) ? comm->getRank() : -1;
81 w = (maxVal == v) ? comm->getRank() : -1;
85 MT avgVal_MT = Teuchos::as<MT>(avgVal);
92 double avgVal, devVal;
94 calculateStats<Type>(minVal, maxVal, avgVal, devVal, minProc, maxProc, comm, numActiveProcs, v);
98 std::ostringstream buf;
100 if ((avgVal != zero) && (paramList.is_null() || !paramList->isParameter(
"print abs") || paramList->get<
bool>(
"print abs") ==
false)) {
101 double relDev = (devVal / avgVal) * 100;
104 buf <<
"avg = " << std::scientific << std::setw(10) << std::setprecision(2) << avgVal <<
", "
105 <<
"dev = " << std::fixed << std::setw(6) << std::setprecision(1) << relDev <<
"%, "
106 <<
"min = " << std::fixed << std::setw(7) << std::setprecision(1) << std::setw(7) << relMin <<
"%"
107 <<
" (" << std::scientific << std::setw(10) << std::setprecision(2) << minVal <<
" on " << std::fixed << std::setw(4) << minProc <<
"), "
108 <<
"max = " << std::fixed << std::setw(7) << std::setprecision(1) << relMax <<
"%"
109 <<
" (" << std::scientific << std::setw(10) << std::setprecision(2) << maxVal <<
" on " << std::fixed << std::setw(4) << maxProc <<
")";
111 double relDev = (avgVal != zero ? (devVal / avgVal) * 100 : zero);
112 buf <<
"avg = " << std::scientific << std::setw(10) << std::setprecision(2) << avgVal <<
", "
113 <<
"dev = " << std::fixed << std::setw(6) << std::setprecision(1) << relDev <<
"%, "
114 <<
"min = " << std::scientific << std::setw(10) << std::setprecision(2) << minVal
115 <<
" (on " << std::fixed << std::setw(4) << minProc <<
"), "
116 <<
"max = " << std::scientific << std::setw(10) << std::setprecision(2) << maxVal
117 <<
" (on " << std::fixed << std::setw(4) << maxProc <<
")";
123 bool cmp_less(
typename Map::value_type& v1,
typename Map::value_type& v2) {
124 return v1.second < v2.second;
127 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
134 std::ostringstream ss;
136 ss << msgTag <<
" size = " << A.getGlobalNumRows() <<
" x " << A.getGlobalNumCols();
137 if (A.haveGlobalConstants())
138 ss <<
", nnz = " << A.getGlobalNumEntries();
144 bool printLoadBalanceInfo =
false, printCommInfo =
false, printEntryStats =
false;
145 if (params->
isParameter(
"printLoadBalancingInfo") && params->
get<
bool>(
"printLoadBalancingInfo"))
146 printLoadBalanceInfo =
true;
147 if (params->
isParameter(
"printCommInfo") && params->
get<
bool>(
"printCommInfo"))
148 printCommInfo =
true;
149 if (params->
isParameter(
"printEntryStats") && params->
get<
bool>(
"printEntryStats"))
150 printEntryStats =
true;
152 if (!printLoadBalanceInfo && !printCommInfo && !printEntryStats)
158 size_t numMyNnz = A.getLocalNumEntries(), numMyRows = A.getLocalNumRows();
162 bool activeProc =
true;
163 int numProc = origComm->getSize();
164 int numActiveProcs = 0;
167 MPI_Comm rawComm = (*mpiComm->getRawMpiComm())();
169 std::vector<size_t> numRowsPerProc(numProc);
170 Teuchos::gatherAll(*origComm, 1, &numMyRows, numProc, &numRowsPerProc[0]);
173 bool rootFlag =
true;
174 for (
int i = 0; i < numProc; i++) {
175 if (numRowsPerProc[i]) {
184 if (numMyRows == 0) {
189 if (numMyRows == 0) {
202 absList.
set(
"print abs",
true);
208 crsA = crsWrapA->getCrsMatrix();
209 if (printEntryStats && !crsA.is_null()) {
211 typedef typename STS::magnitudeType magnitudeType;
221 crsA->getAllValues(rowptr_RCP, colind_RCP, vals_RCP);
222 crsA->getLocalDiagOffsets(offsets_RCP);
223 rowptr = rowptr_RCP();
225 offsets = offsets_RCP();
227 Scalar val, minVal, maxVal;
228 magnitudeType absVal, minAbsVal, maxAbsVal;
230 minVal = STS::rmax();
231 maxVal = STS::rmin();
232 minAbsVal = MTS::rmax();
233 maxAbsVal = MTS::zero();
235 for (
int i = 0; i < offsets.
size(); i++) {
236 val = vals[rowptr[i] + offsets[i]];
237 if (STS::real(val) < STS::real(minVal))
239 if (STS::real(val) > STS::real(maxVal))
241 absVal = STS::magnitude(val);
242 minAbsVal = std::min(minAbsVal, absVal);
243 maxAbsVal = std::max(maxAbsVal, absVal);
246 ss << msgTag <<
" diag min : " << stringStats<Scalar>(origComm, numActiveProcs, minVal) << std::endl;
247 ss << msgTag <<
" diag max : " << stringStats<Scalar>(origComm, numActiveProcs, maxVal) << std::endl;
248 ss << msgTag <<
" abs(diag) min : " << stringStats<Scalar>(origComm, numActiveProcs, minAbsVal) << std::endl;
249 ss << msgTag <<
" abs(diag) max : " << stringStats<Scalar>(origComm, numActiveProcs, maxAbsVal) << std::endl;
253 minVal = STS::rmax();
254 maxVal = STS::rmin();
255 minAbsVal = MTS::rmax();
256 maxAbsVal = MTS::zero();
258 for (
int i = 0; i < vals.
size(); i++) {
260 if (STS::real(val) < STS::real(minVal))
262 if (STS::real(val) > STS::real(maxVal))
264 absVal = STS::magnitude(val);
265 minAbsVal = std::min(minAbsVal, absVal);
266 maxAbsVal = std::max(maxAbsVal, absVal);
269 ss << msgTag <<
" entry min : " << stringStats<Scalar>(origComm, numActiveProcs, minVal) << std::endl;
270 ss << msgTag <<
" entry max : " << stringStats<Scalar>(origComm, numActiveProcs, maxVal) << std::endl;
271 ss << msgTag <<
" abs(entry) min : " << stringStats<Scalar>(origComm, numActiveProcs, minAbsVal) << std::endl;
272 ss << msgTag <<
" abs(entry) max : " << stringStats<Scalar>(origComm, numActiveProcs, maxAbsVal) << std::endl;
276 if (printLoadBalanceInfo) {
277 ss << msgTag <<
" Load balancing info" << std::endl;
278 ss << msgTag <<
" # active processes: " << numActiveProcs <<
"/" << numProc << std::endl;
279 ss << msgTag <<
" # rows per proc : " << stringStats<global_size_t>(origComm, numActiveProcs, numMyRows) << std::endl;
280 ss << msgTag <<
" # nnz per proc : " << stringStats<global_size_t>(origComm, numActiveProcs, numMyNnz) << std::endl;
283 if (printCommInfo && numActiveProcs != 1) {
284 typedef std::map<int, size_t> map_type;
288 if (exportPIDs.
size())
289 for (
int i = 0; i < exportPIDs.
size(); i++)
290 neighMap[exportPIDs[i]]++;
294 size_t numExportSend = 0;
295 size_t numImportSend = 0;
301 numExportSend = (!exporter.
is_null() ? exporter->getNumExportIDs() : 0);
302 numImportSend = (!importer.
is_null() ? importer->getNumExportIDs() : 0);
303 numMsgs = neighMap.size();
304 map_type::const_iterator it = std::min_element(neighMap.begin(), neighMap.end(), cmp_less<map_type>);
305 minMsg = (it != neighMap.end() ? it->second : 0);
306 it = std::max_element(neighMap.begin(), neighMap.end(), cmp_less<map_type>);
307 maxMsg = (it != neighMap.end() ? it->second : 0);
310 ss << msgTag <<
" Communication info" << std::endl;
311 ss << msgTag <<
" # num export send : " << stringStats<global_size_t>(origComm, numActiveProcs, numExportSend) << std::endl;
312 ss << msgTag <<
" # num import send : " << stringStats<global_size_t>(origComm, numActiveProcs, numImportSend) << std::endl;
313 ss << msgTag <<
" # num msgs : " << stringStats<global_size_t>(origComm, numActiveProcs, numMsgs, rcpFromRef(absList)) << std::endl;
314 ss << msgTag <<
" # min msg size : " << stringStats<global_size_t>(origComm, numActiveProcs, minMsg) << std::endl;
315 ss << msgTag <<
" # max msg size : " << stringStats<global_size_t>(origComm, numActiveProcs, maxMsg) << std::endl;
321 int strLength = outstr.size();
322 MPI_Bcast(&strLength, 1, MPI_INT, root, rawComm);
323 if (origComm->getRank() != root)
324 outstr.resize(strLength);
325 MPI_Bcast(&outstr[0], strLength, MPI_CHAR, root, rawComm);
331 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
335 std::ostringstream ss;
339 bool activeProc =
true;
340 int numActiveProcs = origComm->getSize();
343 MPI_Comm rawComm = (*mpiComm->getRawMpiComm())();
349 absList.
set(
"print abs",
true);
351 typedef std::map<int, size_t> map_type;
354 if (exportPIDs.
size())
355 for (
int i = 0; i < exportPIDs.
size(); i++)
356 neighMap[exportPIDs[i]]++;
359 size_t numImportSend = 0;
365 numImportSend = importer->getNumExportIDs();
366 numMsgs = neighMap.size();
367 map_type::const_iterator it = std::min_element(neighMap.begin(), neighMap.end(), cmp_less<map_type>);
368 minMsg = (it != neighMap.end() ? it->second : 0);
369 it = std::max_element(neighMap.begin(), neighMap.end(), cmp_less<map_type>);
370 maxMsg = (it != neighMap.end() ? it->second : 0);
373 ss << msgTag <<
" Communication info" << std::endl;
374 ss << msgTag <<
" # num import send : " << stringStats<global_size_t>(origComm, numActiveProcs, numImportSend) << std::endl;
375 ss << msgTag <<
" # num msgs : " << stringStats<global_size_t>(origComm, numActiveProcs, numMsgs, rcpFromRef(absList)) << std::endl;
376 ss << msgTag <<
" # min msg size : " << stringStats<global_size_t>(origComm, numActiveProcs, minMsg) << std::endl;
377 ss << msgTag <<
" # max msg size : " << stringStats<global_size_t>(origComm, numActiveProcs, maxMsg) << std::endl;
382 int strLength = outstr.size();
383 MPI_Bcast(&strLength, 1, MPI_INT, root, rawComm);
384 if (origComm->getRank() != root)
385 outstr.resize(strLength);
386 MPI_Bcast(&outstr[0], strLength, MPI_CHAR, root, rawComm);
392 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
397 std::ostringstream out;
400 int myRank = comm->getRank();
402 out << msgTag <<
" " << myRank <<
":";
404 RCP<const Import> importer = (A.getCrsGraph() != Teuchos::null ? A.getCrsGraph()->getImporter() : Teuchos::null);
412 if (exportPIDs.
size()) {
414 int neigh = exportPIDs[0];
416 for (
int i = 1; i < exportPIDs.
size(); i++) {
417 if (exportPIDs[i] != exportPIDs[i - 1]) {
418 out <<
" " << neigh <<
"(" << weight <<
")";
420 neigh = exportPIDs[i];
427 out <<
" " << neigh <<
"(" << weight <<
")" << std::endl;
433 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
440 bool hasCrsGraph =
true;
453 #endif // MUELU_PERFUTILS_DEF_HPP
static bool CheckMatrix(const Matrix &A)
std::string stringStats(const RCP< const Teuchos::Comm< int > > &comm, int numActiveProcs, const Type &v, RCP< ParameterList > paramList=Teuchos::null)
#define MueLu_sumAll(rcpComm, in, out)
#define MueLu_maxAll(rcpComm, in, out)
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
static magnitudeType real(T a)
#define MueLu_minAll(rcpComm, in, out)
static std::string CommPattern(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
static std::string PrintImporterInfo(RCP< const Import > importer, const std::string &msgTag)
bool isParameter(const std::string &name) const
MueLu::DefaultScalar Scalar
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
void calculateStats(Type &minVal, Type &maxVal, double &avgVal, double &devVal, int &minProc, int &maxProc, const RCP< const Teuchos::Comm< int > > &comm, int numActiveProcs, const Type &v)
bool cmp_less(typename Map::value_type &v1, typename Map::value_type &v2)