diff --git a/AMDiS/libtool b/AMDiS/libtool
index 85b81df55dafa7bf065db39a64a411063aebf89f..e8aeacc60e813eddedb5385b394097a7cd2f25b2 100755
--- a/AMDiS/libtool
+++ b/AMDiS/libtool
@@ -44,7 +44,7 @@ available_tags=" CXX F77"
# ### BEGIN LIBTOOL CONFIG
-# Libtool was configured on host deimos104:
+# Libtool was configured on host p1d098:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
@@ -82,13 +82,13 @@ AR="ar"
AR_FLAGS="cru"
# A C compiler.
-LTCC="gcc"
+LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
# LTCC compiler flags.
LTCFLAGS="-g -O3"
# A language-specific compiler.
-CC="gcc"
+CC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
# Is the compiler the GNU C compiler?
with_gcc=yes
@@ -174,7 +174,7 @@ dlopen_self=unknown
dlopen_self_static=unknown
# Compiler flag to prevent dynamic linking.
-link_static_flag="-static"
+link_static_flag=""
# Compiler flag to turn off builtin functions.
no_builtin_flag=" -fno-builtin"
@@ -6763,7 +6763,7 @@ build_old_libs=`case $build_libtool_libs in yes) $echo no;; *) $echo yes;; esac`
# End:
# ### BEGIN LIBTOOL TAG CONFIG: CXX
-# Libtool was configured on host deimos104:
+# Libtool was configured on host p1d098:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
@@ -6801,13 +6801,13 @@ AR="ar"
AR_FLAGS="cru"
# A C compiler.
-LTCC="gcc"
+LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
# LTCC compiler flags.
LTCFLAGS="-g -O3"
# A language-specific compiler.
-CC="g++"
+CC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicxx"
# Is the compiler the GNU C compiler?
with_gcc=yes
@@ -6893,7 +6893,7 @@ dlopen_self=unknown
dlopen_self_static=unknown
# Compiler flag to prevent dynamic linking.
-link_static_flag="-static"
+link_static_flag=""
# Compiler flag to turn off builtin functions.
no_builtin_flag=" -fno-builtin"
@@ -6960,11 +6960,11 @@ predeps=""
# Dependencies to place after the objects being linked to create a
# shared library.
-postdeps="-lstdc++ -lm -lgcc_s -lc -lgcc_s"
+postdeps="-lmpi_cxx -lmpi -lopen-rte -lopen-pal -libverbs -lrt -lnuma -ldl -lnsl -lutil -ldl -lstdc++ -lm -lgcc_s -lpthread -lc -lgcc_s"
# The library search path used internally by the compiler when linking
# a shared library.
-compiler_lib_search_path=`echo "-L/usr/lib64/gcc/x86_64-suse-linux/4.1.2 -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64 -L/lib/../lib64 -L/usr/lib/../lib64 -L/fastfs/wir/local/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../.." | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"`
+compiler_lib_search_path=`echo "-L/usr/lib64 -L/licsoft/libraries/openmpi/1.2.6/64bit/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2 -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64 -L/lib/../lib64 -L/usr/lib/../lib64 -L/fastfs/wir/local/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../.." | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"`
# Method to check whether dependent libraries are shared objects.
deplibs_check_method="pass_all"
@@ -7071,7 +7071,7 @@ include_expsyms=""
# ### BEGIN LIBTOOL TAG CONFIG: F77
-# Libtool was configured on host deimos104:
+# Libtool was configured on host p1d098:
# Shell to use when invoking shell scripts.
SHELL="/bin/sh"
@@ -7109,7 +7109,7 @@ AR="ar"
AR_FLAGS="cru"
# A C compiler.
-LTCC="gcc"
+LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
# LTCC compiler flags.
LTCFLAGS="-g -O3"
diff --git a/AMDiS/src/parallel/PetscSolver.cc b/AMDiS/src/parallel/PetscSolver.cc
index 8c6d4245cb78d7eff8857425382fd5f9b3829806..5353ac6f17d0fc1ea312262ca33a954a75d45b91 100644
--- a/AMDiS/src/parallel/PetscSolver.cc
+++ b/AMDiS/src/parallel/PetscSolver.cc
@@ -24,10 +24,12 @@
namespace AMDiS {
+ using namespace std;
+
PetscErrorCode myKSPMonitor(KSP ksp, PetscInt iter, PetscReal rnorm, void *)
{
if (iter % 100 == 0 && MPI::COMM_WORLD.Get_rank() == 0)
- std::cout << "[0] Petsc-Iteration " << iter << ": " << rnorm << std::endl;
+ cout << "[0] Petsc-Iteration " << iter << ": " << rnorm << endl;
return 0;
}
@@ -66,12 +68,12 @@ namespace AMDiS {
typedef traits::range_generator<row, Matrix>::type cursor_type;
typedef traits::range_generator<nz, cursor_type>::type icursor_type;
- std::vector<int> cols;
- std::vector<double> values;
+ vector<int> cols;
+ vector<double> values;
cols.reserve(300);
values.reserve(300);
- std::vector<int> globalCols;
+ vector<int> globalCols;
// === Traverse all rows of the dof matrix and insert row wise the values ===
@@ -86,7 +88,9 @@ namespace AMDiS {
bool periodicRow = meshDistributor->isPeriodicDof(globalRowDof);
if (!periodicRow) {
- // Calculate petsc row index.
+ // === Row DOF index is not periodic. ===
+
+ // Calculate PETSc row index.
int rowIndex = globalRowDof * dispMult + dispAddRow;
cols.clear();
@@ -99,35 +103,54 @@ namespace AMDiS {
int globalColDof = meshDistributor->mapLocalToGlobal(col(*icursor));
// Test if the current col dof is a periodic dof.
bool periodicCol = meshDistributor->isPeriodicDof(globalColDof);
+ // Calculate PETSc col index.
+ int colIndex = globalColDof * dispMult + dispAddCol;
- if (value(*icursor) == 0.0 && globalRowDof != globalColDof)
+ // Ignore all zero entries, expect it is a diagonal entry.
+ // if (value(*icursor) == 0.0 && globalRowDof != globalColDof)
+ if (value(*icursor) == 0.0 && rowIndex != colIndex)
continue;
if (!periodicCol) {
- // Calculate the exact position of the column index in the petsc matrix.
+ // Calculate the exact position of the column index in the PETSc matrix.
cols.push_back(globalColDof * dispMult + dispAddCol);
values.push_back(value(*icursor));
} else {
- std::set<int>& perColAsc =
- meshDistributor->getPerDofAssociations(globalColDof);
+ // === Row index is not periodic, but column index is. ===
+ // Create set of all periodic associations of the column index.
std::set<int> perAsc;
+ std::set<int>& perColAsc =
+ meshDistributor->getPerDofAssociations(globalColDof);
for (std::set<int>::iterator it = perColAsc.begin();
it != perColAsc.end(); ++it)
if (*it >= -3)
perAsc.insert(*it);
+ // Scale value to the number of periodic associations of the column index.
double scaledValue =
- value(*icursor) * std::pow(0.5, static_cast<double>(perAsc.size()));
- std::vector<int> newCols;
+ value(*icursor) * pow(0.5, static_cast<double>(perAsc.size()));
+
+
+ // === Create set of all matrix column indices due to the periodic ===
+ // === associations of the column DOF index. ===
+
+ vector<int> newCols;
+
+ // First, add the original matrix index.
newCols.push_back(globalColDof);
+ // And add all periodic matrix indices.
for (std::set<int>::iterator it = perAsc.begin();
it != perAsc.end(); ++it) {
- int nCols = static_cast<int>(newCols.size());
- for (int i = 0; i < nCols; i++) {
+ for (unsigned int i = 0; i < newCols.size(); i++) {
TEST_EXIT(meshDistributor->isPeriodicDof(newCols[i], *it))
- ("Should not happen: %d %d\n", *it, newCols[i]);
+ ("Wrong periodic DOF associations at boundary %d with DOF %d!\n",
+ *it, newCols[i]);
+
+ MSG("MAP DOF %d -> %d\n",
+ newCols[i],
+ meshDistributor->getPeriodicMapping(newCols[i], *it));
newCols.push_back(meshDistributor->getPeriodicMapping(newCols[i], *it));
}
@@ -143,18 +166,30 @@ namespace AMDiS {
MatSetValues(petscMatrix, 1, &rowIndex, cols.size(),
&(cols[0]), &(values[0]), ADD_VALUES);
} else {
- std::map<int, std::vector<int> > colsMap;
- std::map<int, std::vector<double> > valsMap;
+ // === Row DOF index is periodic. ===
+ // Because this row is periodic, we will have to add the entries of this
+ // matrix row to multiple rows. The following maps store to each row an
+ // array of column indices and values of the entries that must be added to
+ // the PETSc matrix.
+ map<int, vector<int> > colsMap;
+ map<int, vector<double> > valsMap;
+
+ // Traverse all column entries.
for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor);
icursor != icend; ++icursor) {
// Global index of the current column index.
int globalColDof = meshDistributor->mapLocalToGlobal(col(*icursor));
+ // Ignore all zero entries, expect it is a diagonal entry.
if (value(*icursor) == 0.0 && globalRowDof != globalColDof)
continue;
+
+ // === Add all periodic associations of both, the row and the column ===
+ // === indices to the set perAsc. ===
+
std::set<int> perAsc;
if (meshDistributor->isPeriodicDof(globalColDof)) {
@@ -173,11 +208,21 @@ namespace AMDiS {
if (*it >= -3)
perAsc.insert(*it);
+ // Scale the value with respect to the number of periodic associations.
double scaledValue =
- value(*icursor) * std::pow(0.5, static_cast<double>(perAsc.size()));
- std::vector<std::pair<int, int> > entry;
- entry.push_back(std::make_pair(globalRowDof, globalColDof));
+ value(*icursor) * pow(0.5, static_cast<double>(perAsc.size()));
+
+
+ // === Create all matrix entries with respect to the periodic ===
+ // === associations of the row and column indices. ===
+
+ vector<pair<int, int> > entry;
+
+ // First, add the original entry.
+ entry.push_back(make_pair(globalRowDof, globalColDof));
+ // Then, traverse the periodic associations of the row and column indices
+ // and create the corresponding entries.
for (std::set<int>::iterator it = perAsc.begin(); it != perAsc.end(); ++it) {
int nEntry = static_cast<int>(entry.size());
for (int i = 0; i < nEntry; i++) {
@@ -194,23 +239,28 @@ namespace AMDiS {
perColDof = entry[i].second;
- entry.push_back(std::make_pair(perRowDof, perColDof));
+ entry.push_back(make_pair(perRowDof, perColDof));
}
}
- for (std::vector<std::pair<int, int> >::iterator eIt = entry.begin();
+
+ // === Translate the matrix entries to PETSc's matrix.
+
+ for (vector<pair<int, int> >::iterator eIt = entry.begin();
eIt != entry.end(); ++eIt) {
- // Calculate petsc row index.
+ // Calculate row and column indices of the PETSc matrix.
int rowIndex = eIt->first * dispMult + dispAddRow;
-
int colIndex = eIt->second * dispMult + dispAddCol;
+
colsMap[rowIndex].push_back(colIndex);
valsMap[rowIndex].push_back(scaledValue);
}
}
- for (std::map<int, std::vector<int> >::iterator rowIt = colsMap.begin();
+ // === Finally, add all periodic rows to the PETSc matrix. ===
+
+ for (map<int, vector<int> >::iterator rowIt = colsMap.begin();
rowIt != colsMap.end(); ++rowIt) {
TEST_EXIT_DBG(rowIt->second.size() == valsMap[rowIt->first].size())
("Should not happen!\n");
@@ -275,18 +325,18 @@ namespace AMDiS {
using mtl::tag::row; using mtl::tag::nz; using mtl::begin; using mtl::end;
namespace traits = mtl::traits;
typedef DOFMatrix::base_matrix_type Matrix;
- typedef std::vector<std::pair<int, int> > MatrixNnzEntry;
- typedef std::map<int, DofContainer> RankToDofContainer;
+ typedef vector<pair<int, int> > MatrixNnzEntry;
+ typedef map<int, DofContainer> RankToDofContainer;
// Stores to each rank a list of nnz entries (i.e. pairs of row and column index)
// that this rank will send to. These nnz entries will be assembled on this rank,
// but because the row DOFs are not DOFs of this rank they will be send to the
// owner of the row DOFs.
- std::map<int, MatrixNnzEntry> sendMatrixEntry;
+ map<int, MatrixNnzEntry> sendMatrixEntry;
// Maps to each DOF that must be send to another rank the rank number of the
// receiving rank.
- std::map<DegreeOfFreedom, int> sendDofToRank;
+ map<DegreeOfFreedom, int> sendDofToRank;
// First, create for all ranks we send data to MatrixNnzEntry object with 0 entries.
@@ -324,9 +374,9 @@ namespace AMDiS {
int globalRowDof = meshDistributor->mapLocalToGlobal(*cursor);
- std::vector<int> rows;
+ vector<int> rows;
rows.push_back(globalRowDof);
- std::vector<int> rowRank;
+ vector<int> rowRank;
if (meshDistributor->getIsRankDof(*cursor)) {
rowRank.push_back(meshDistributor->getMpiRank());
} else {
@@ -387,7 +437,7 @@ namespace AMDiS {
meshDistributor->mapLocalToGlobal(col(*icursor)) * nComponents + j;
sendMatrixEntry[sendToRank].
- push_back(std::make_pair(petscRowIdx, petscColIdx));
+ push_back(make_pair(petscRowIdx, petscColIdx));
}
}
@@ -409,7 +459,7 @@ namespace AMDiS {
// === Evaluate the nnz structure this rank got from other ranks and add it to ===
// === the PETSc nnz data structure. ===
- for (std::map<int, MatrixNnzEntry>::iterator it = stdMpi.getRecvData().begin();
+ for (map<int, MatrixNnzEntry>::iterator it = stdMpi.getRecvData().begin();
it != stdMpi.getRecvData().end(); ++it) {
if (it->second.size() > 0) {
for (unsigned int i = 0; i < it->second.size(); i++) {