diff --git a/AMDiS/src/parallel/MatrixNnzStructure.cc b/AMDiS/src/parallel/MatrixNnzStructure.cc
index 82cc355a3c3036e173c54f403ec4583d97441ab5..382267c84edb6cf92f395e0ededda45cf642aa03 100644
--- a/AMDiS/src/parallel/MatrixNnzStructure.cc
+++ b/AMDiS/src/parallel/MatrixNnzStructure.cc
@@ -10,6 +10,7 @@
// See also license.opensource.txt in the distribution.
+#include <algorithm>
#include "Global.h"
#include "DOFMatrix.h"
#include "parallel/MatrixNnzStructure.h"
@@ -286,8 +287,8 @@ namespace AMDiS {
nSum = 0;
for (int i = 0; i < nRankRows; i++) {
- nMax = std::max(nMax, nnz.onnz[i]);
- nSum += nnz.onnz[i];
+ nMax = std::max(nMax, onnz[i]);
+ nSum += onnz[i];
}
MSG("NNZ in offdiag block: max = %d, avrg = %.0f\n",
diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc
index 418fb714a173cc9c3afcb55d1cc46f5fd46ea30e..aa4a7847cbd6f602ba443cd7dccdd2068588db3f 100644
--- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc
+++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc
@@ -1016,189 +1016,4 @@ namespace AMDiS {
}
}
-#if 0
- void PetscSolverGlobalMatrix::createPetscNnzStructure(Matrix<DOFMatrix*> *mat)
- {
- FUNCNAME("PetscSolverGlobalMatrix::createPetscNnzStructure()");
-
- vector<const FiniteElemSpace*> feSpaces = getFeSpaces(mat);
- int nRankRows = interiorMap->getRankDofs();
- int rankStartIndex = interiorMap->getStartDofs();
-
- nnzInterior.create(nRankRows, nRankRows);
-
- 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 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.
- map<int, MatrixNnzEntry> sendMatrixEntry;
-
- // Maps to each DOF that must be send to another rank the rank number of the
- // receiving rank.
- map<pair<DegreeOfFreedom, int>, int> sendDofToRank;
-
-
- // First, create for all ranks, to which we send data to, MatrixNnzEntry
- // object with 0 entries.
- for (unsigned int i = 0; i < feSpaces.size(); i++) {
- for (DofComm::Iterator it(meshDistributor->getDofComm().getRecvDofs(), feSpaces[i]);
- !it.end(); it.nextRank()) {
- sendMatrixEntry[it.getRank()].resize(0);
-
- for (; !it.endDofIter(); it.nextDof())
- sendDofToRank[make_pair(it.getDofIndex(), i)] = it.getRank();
- }
- }
-
- // Create list of ranks from which we receive data from.
- std::set<int> recvFromRank;
- for (unsigned int i = 0; i < feSpaces.size(); i++)
- for (DofComm::Iterator it(meshDistributor->getDofComm().getSendDofs(), feSpaces[i]);
- !it.end(); it.nextRank())
- recvFromRank.insert(it.getRank());
-
-
- // === Traverse matrices to create nnz data. ===
-
- int nComponents = mat->getNumRows();
- for (int i = 0; i < nComponents; i++) {
- for (int j = 0; j < nComponents; j++) {
- if (!(*mat)[i][j])
- continue;
-
- TEST_EXIT_DBG((*mat)[i][j]->getRowFeSpace() == feSpaces[i])
- ("Should not happen!\n");
- TEST_EXIT_DBG((*mat)[i][j]->getColFeSpace() == feSpaces[j])
- ("Should not happen!\n");
-
- Matrix bmat = (*mat)[i][j]->getBaseMatrix();
-
- traits::col<Matrix>::type col(bmat);
- traits::const_value<Matrix>::type value(bmat);
-
- typedef traits::range_generator<row, Matrix>::type cursor_type;
- typedef traits::range_generator<nz, cursor_type>::type icursor_type;
-
- for (cursor_type cursor = begin<row>(bmat),
- cend = end<row>(bmat); cursor != cend; ++cursor) {
- int globalRowDof = (*interiorMap)[feSpaces[i]][*cursor].global;
-
- // The corresponding global matrix row index of the current row DOF.
- int petscRowIdx = interiorMap->getMatIndex(i, globalRowDof);
- if ((*interiorMap)[feSpaces[i]].isRankDof(*cursor)) {
-
- // === The current row DOF is a rank DOF, so create the ===
- // === corresponding nnz values directly on rank's nnz data. ===
-
- // This is the local row index of the local PETSc matrix.
- int localPetscRowIdx = petscRowIdx - rankStartIndex;
-
- TEST_EXIT_DBG(localPetscRowIdx >= 0 && localPetscRowIdx < nRankRows)
- ("Should not happen! \n Debug info: DOF = %d globalRowIndx = %d petscRowIdx = %d localPetscRowIdx = %d rStart = %d compontens = %d from %d nRankRows = %d\n",
- *cursor,
- (*interiorMap)[feSpaces[i]][*cursor].global,
- petscRowIdx,
- localPetscRowIdx,
- rankStartIndex,
- i,
- nComponents,
- nRankRows);
-
-
- // Traverse all non zero entries in this row.
- for (icursor_type icursor = begin<nz>(cursor),
- icend = end<nz>(cursor); icursor != icend; ++icursor) {
- int globalColDof = (*interiorMap)[feSpaces[j]][col(*icursor)].global;
- int petscColIdx = interiorMap->getMatIndex(j, globalColDof);
-
- if (value(*icursor) != 0.0 || petscRowIdx == petscColIdx) {
- // The row DOF is a rank DOF, if also the column is a rank DOF,
- // increment the matInteriorDiagNnz values for this row,
- // otherwise the matInteriorDiagNnz value.
- if (petscColIdx >= rankStartIndex &&
- petscColIdx < rankStartIndex + nRankRows)
- nnzInterior.dnnz[localPetscRowIdx]++;
- else
- nnzInterior.onnz[localPetscRowIdx]++;
- }
- }
- } else {
- // === The current row DOF is not a rank DOF, i.e., its values ===
- // === are also created on this rank, but afterthere they will ===
- // === be send to another rank. So we need to send also the ===
- // === corresponding nnz structure of this row to the corres- ===
- // === ponding rank. ===
-
- // Send all non zero entries to the member of the row DOF.
- int sendToRank = sendDofToRank[make_pair(*cursor, i)];
-
- for (icursor_type icursor = begin<nz>(cursor),
- icend = end<nz>(cursor); icursor != icend; ++icursor) {
- if (value(*icursor) != 0.0) {
- int globalColDof =
- (*interiorMap)[feSpaces[j]][col(*icursor)].global;
- int petscColIdx = interiorMap->getMatIndex(j, globalColDof);
-
- sendMatrixEntry[sendToRank].
- push_back(make_pair(petscRowIdx, petscColIdx));
- }
- }
-
- } // if (isRankDof[*cursor]) ... else ...
- } // for each row in mat[i][j]
- }
- }
-
- // === Send and recv the nnz row structure to/from other ranks. ===
-
- StdMpi<MatrixNnzEntry> stdMpi(mpiCommGlobal, true);
- stdMpi.send(sendMatrixEntry);
- for (std::set<int>::iterator it = recvFromRank.begin();
- it != recvFromRank.end(); ++it)
- stdMpi.recv(*it);
- stdMpi.startCommunication();
-
-
- // === Evaluate the nnz structure this rank got from other ranks and add ===
- // === it to the PETSc nnz data structure. ===
-
- 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++) {
- int r = it->second[i].first;
- int c = it->second[i].second;
-
- int localRowIdx = r - rankStartIndex;
-
- TEST_EXIT_DBG(localRowIdx >= 0 && localRowIdx < nRankRows)
- ("Got row index %d/%d (nRankRows = %d) from rank %d. Should not happen!\n",
- r, localRowIdx, nRankRows, it->first);
-
- if (c < rankStartIndex || c >= rankStartIndex + nRankRows)
- nnzInterior.onnz[localRowIdx]++;
- else
- nnzInterior.dnnz[localRowIdx]++;
- }
- }
- }
-
- // The above algorithm for calculating the number of nnz per row over-
- // approximates the value, i.e., the number is always equal or larger to
- // the real number of nnz values in the global parallel matrix. For small
- // matrices, the problem may arise, that the result is larger than the
- // number of elements in a row. This is fixed in the following.
-
- if (nRankRows < 100)
- for (int i = 0; i < nRankRows; i++)
- nnzInterior.dnnz[i] = std::min(nnzInterior.dnnz[i], nRankRows);
- }
-#endif
-
}
diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.h b/AMDiS/src/parallel/PetscSolverGlobalMatrix.h
index c8624b53c71fc4c678c920f955ad48517a1a4aa6..a70d5590e2dff66624357faffaf59714508cd54e 100644
--- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.h
+++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.h
@@ -75,9 +75,6 @@ namespace AMDiS {
virtual void exitPreconditioner(PC pc);
- /// Creates a new non zero pattern structure for the PETSc matrix.
- // void createPetscNnzStructure(Matrix<DOFMatrix*> *mat);
-
/// Takes a DOF matrix and sends the values to the global PETSc matrix.
void setDofMatrix(DOFMatrix* mat, int nRowMat = 0, int nColMat = 0);