Commit 5a7ff74f authored by Thomas Witkowski's avatar Thomas Witkowski

Removed parallel dof mapping access in MeshDistributor.

parent 90357cd5
...@@ -145,70 +145,16 @@ namespace AMDiS { ...@@ -145,70 +145,16 @@ namespace AMDiS {
return feSpaces; return feSpaces;
} }
/// Returns the number of DOFs in rank's domain for a given FE space. /// Returns the DOF mapping object, \ref dofMap.
inline int getNumberRankDofs(const FiniteElemSpace *feSpace) inline ParallelDofMapping& getDofMap()
{
return dofMap[feSpace].nRankDofs;
}
/// Returns the number of DOFs in rank's domain for a set of FE spaces.
inline int getNumberRankDofs(vector<const FiniteElemSpace*>& feSpaces)
{
FUNCNAME("MeshDistributor::getNumberRankDofs()");
int result = 0;
for (unsigned int i = 0; i < feSpaces.size(); i++)
result += dofMap[feSpaces[i]].nRankDofs;
return result;
}
/// Returns the first global DOF index of an FE space, owned by rank.
inline int getStartDofs(const FiniteElemSpace *feSpace)
{
return dofMap[feSpace].rStartDofs;
}
/// Returns the first global DOF index for a set of FE spaces, owned by rank.
inline int getStartDofs(vector<const FiniteElemSpace*>& feSpaces)
{
FUNCNAME("MeshDistributor::getStartDofs()");
int result = 0;
for (unsigned int i = 0; i < feSpaces.size(); i++)
result += dofMap[feSpaces[i]].rStartDofs;
return result;
}
/// Returns the global number of DOFs for a given FE space.
inline int getNumberOverallDofs(const FiniteElemSpace *feSpace)
{
return dofMap[feSpace].nOverallDofs;
}
/// Returns the global number of DOFs for a set of FE spaces.
inline int getNumberOverallDofs(vector<const FiniteElemSpace*>& feSpaces)
{
FUNCNAME("MeshDistributor::getNumberOverallDofs()");
int result = 0;
for (unsigned int i = 0; i < feSpaces.size(); i++)
result += dofMap[feSpaces[i]].nOverallDofs;
return result;
}
inline map<DegreeOfFreedom, MultiIndex>& getMapDofToGlobal(const FiniteElemSpace *feSpace)
{ {
return dofMap[feSpace].getMap(); return dofMap;
} }
/// Maps a local DOF to its global index. /// Returns the periodic mapping handler, \ref periodicMap.
inline DegreeOfFreedom mapDofToGlobal(const FiniteElemSpace *feSpace, inline PeriodicMap& getPeriodicMap()
DegreeOfFreedom dof)
{ {
return dofMap[feSpace][dof].global; return periodicMap;
} }
/// Returns for a global index the DOF index in rank's subdomain. As there /// Returns for a global index the DOF index in rank's subdomain. As there
...@@ -218,25 +164,6 @@ namespace AMDiS { ...@@ -218,25 +164,6 @@ namespace AMDiS {
DegreeOfFreedom mapGlobalToLocal(const FiniteElemSpace *feSpace, DegreeOfFreedom mapGlobalToLocal(const FiniteElemSpace *feSpace,
DegreeOfFreedom dof); DegreeOfFreedom dof);
/// Maps a local DOF to its local index.
inline DegreeOfFreedom mapLocalToDof(const FiniteElemSpace *feSpace,
DegreeOfFreedom dof)
{
return dofMap[feSpace][dof].local;
}
/// Returns the DOF mapping object, \ref dofMap.
inline ParallelDofMapping& getDofMap()
{
return dofMap;
}
/// Returns the periodic mapping handler, \ref periodicMap.
inline PeriodicMap& getPeriodicMap()
{
return periodicMap;
}
DofComm& getSendDofs() DofComm& getSendDofs()
{ {
return sendDofs; return sendDofs;
...@@ -252,13 +179,6 @@ namespace AMDiS { ...@@ -252,13 +179,6 @@ namespace AMDiS {
return periodicDofs; return periodicDofs;
} }
/// Return true, if the given DOF is owned by the rank. If false, the DOF
/// is in rank's partition, but is owned by some other rank.
inline bool getIsRankDof(const FiniteElemSpace *feSpace, DegreeOfFreedom dof)
{
return dofMap[feSpace].isRankDof(dof);
}
inline long getLastMeshChangeIndex() inline long getLastMeshChangeIndex()
{ {
return lastMeshChangeIndex; return lastMeshChangeIndex;
......
...@@ -775,7 +775,7 @@ namespace AMDiS { ...@@ -775,7 +775,7 @@ namespace AMDiS {
for (it.reset(); !it.end(); ++it) { for (it.reset(); !it.end(); ++it) {
file << it.getDOFIndex() << " " file << it.getDOFIndex() << " "
<< pdb.dofMap[feSpace][it.getDOFIndex()].global << " " << pdb.dofMap[feSpace][it.getDOFIndex()].global << " "
<< pdb.getIsRankDof(feSpace, it.getDOFIndex()); << pdb.dofMap[feSpace].isRankDof(it.getDOFIndex());
for (int i = 0; i < pdb.mesh->getDim(); i++) for (int i = 0; i < pdb.mesh->getDim(); i++)
file << " " << (*it)[i]; file << " " << (*it)[i];
file << "\n"; file << "\n";
......
...@@ -336,7 +336,7 @@ namespace AMDiS { ...@@ -336,7 +336,7 @@ namespace AMDiS {
// === create local indices of the primals starting at zero. === // === create local indices of the primals starting at zero. ===
for (DofIndexSet::iterator it = primals.begin(); it != primals.end(); ++it) for (DofIndexSet::iterator it = primals.begin(); it != primals.end(); ++it)
if (meshDistributor->getIsRankDof(feSpace, *it)) if (meshDistributor->getDofMap()[feSpace].isRankDof(*it))
primalDofMap[feSpace].insertRankDof(*it); primalDofMap[feSpace].insertRankDof(*it);
else else
primalDofMap[feSpace].insert(*it); primalDofMap[feSpace].insert(*it);
...@@ -422,7 +422,8 @@ namespace AMDiS { ...@@ -422,7 +422,8 @@ namespace AMDiS {
int nRankLagrange = 0; int nRankLagrange = 0;
DofMap& dualMap = dualDofMap[feSpace].getMap(); DofMap& dualMap = dualDofMap[feSpace].getMap();
for (DofMap::iterator it = dualMap.begin(); it != dualMap.end(); ++it) { for (DofMap::iterator it = dualMap.begin(); it != dualMap.end(); ++it) {
if (meshDistributor->getIsRankDof(feSpace, it->first)) {
if (meshDistributor->getDofMap()[feSpace].isRankDof(it->first)) {
lagrangeMap[feSpace].insertRankDof(it->first, nRankLagrange); lagrangeMap[feSpace].insertRankDof(it->first, nRankLagrange);
int degree = boundaryDofRanks[feSpace][it->first].size(); int degree = boundaryDofRanks[feSpace][it->first].size();
nRankLagrange += (degree * (degree - 1)) / 2; nRankLagrange += (degree * (degree - 1)) / 2;
......
...@@ -21,13 +21,14 @@ namespace AMDiS { ...@@ -21,13 +21,14 @@ namespace AMDiS {
FUNCNAME("PetscSolverGlobalBlockMatrix::fillPetscMatrix()"); FUNCNAME("PetscSolverGlobalBlockMatrix::fillPetscMatrix()");
TEST_EXIT_DBG(meshDistributor)("No mesh distributor object defined!\n"); TEST_EXIT_DBG(meshDistributor)("No mesh distributor object defined!\n");
TEST_EXIT_DBG(dofMap)("No parallel mapping object defined!\n");
TEST_EXIT_DBG(mat)("No DOF matrix defined!\n"); TEST_EXIT_DBG(mat)("No DOF matrix defined!\n");
double wtime = MPI::Wtime(); double wtime = MPI::Wtime();
const FiniteElemSpace *feSpace = meshDistributor->getFeSpace(0); const FiniteElemSpace *feSpace = meshDistributor->getFeSpace(0);
nComponents = mat->getNumRows(); nComponents = mat->getNumRows();
int nRankRows = meshDistributor->getNumberRankDofs(feSpace); int nRankRows = (*dofMap)[feSpace].nRankDofs;
int nOverallRows = meshDistributor->getNumberOverallDofs(feSpace); int nOverallRows = (*dofMap)[feSpace].nOverallDofs;
#if (DEBUG != 0) #if (DEBUG != 0)
MSG("Fill petsc matrix 1 needed %.5f seconds\n", MPI::Wtime() - wtime); MSG("Fill petsc matrix 1 needed %.5f seconds\n", MPI::Wtime() - wtime);
...@@ -107,8 +108,8 @@ namespace AMDiS { ...@@ -107,8 +108,8 @@ namespace AMDiS {
nComponents = vec->getSize(); nComponents = vec->getSize();
const FiniteElemSpace *feSpace = meshDistributor->getFeSpace(0); const FiniteElemSpace *feSpace = meshDistributor->getFeSpace(0);
int nRankRows = meshDistributor->getNumberRankDofs(feSpace); int nRankRows = (*dofMap)[feSpace].nRankDofs;
int nOverallRows = meshDistributor->getNumberOverallDofs(feSpace); int nOverallRows = (*dofMap)[feSpace].nOverallDofs;
nestVec.resize(nComponents); nestVec.resize(nComponents);
...@@ -150,12 +151,14 @@ namespace AMDiS { ...@@ -150,12 +151,14 @@ namespace AMDiS {
Vec tmp; Vec tmp;
VecNestGetSubVec(petscSolVec, i, &tmp); VecNestGetSubVec(petscSolVec, i, &tmp);
int nRankDofs = meshDistributor->getNumberRankDofs(feSpace); int nRankDofs = (*dofMap)[feSpace].nRankDofs;
PetscScalar *vecPointer; PetscScalar *vecPointer;
VecGetArray(tmp, &vecPointer); VecGetArray(tmp, &vecPointer);
for (int j = 0; j < nRankDofs; j++) DofMap& d = (*dofMap)[feSpace].getMap();
dofvec[meshDistributor->mapLocalToDof(feSpace, j)] = vecPointer[j]; for (DofMap::iterator it = d.begin(); it != d.end(); ++it)
if (it->second.local != -1)
dofvec[it->first] = vecPointer[it->second.local];
VecRestoreArray(tmp, &vecPointer); VecRestoreArray(tmp, &vecPointer);
} }
...@@ -214,8 +217,8 @@ namespace AMDiS { ...@@ -214,8 +217,8 @@ namespace AMDiS {
typedef traits::range_generator<row, Matrix>::type cursor_type; typedef traits::range_generator<row, Matrix>::type cursor_type;
typedef traits::range_generator<nz, cursor_type>::type icursor_type; typedef traits::range_generator<nz, cursor_type>::type icursor_type;
int dispRowIndex = meshDistributor->getNumberRankDofs(feSpace) * dispRowBlock; int dispRowIndex = (*dofMap)[feSpace].nRankDofs * dispRowBlock;
int dispColIndex = meshDistributor->getNumberRankDofs(feSpace) * dispColBlock; int dispColIndex = (*dofMap)[feSpace].nRankDofs * dispColBlock;
vector<int> cols; vector<int> cols;
vector<double> values; vector<double> values;
...@@ -229,8 +232,7 @@ namespace AMDiS { ...@@ -229,8 +232,7 @@ namespace AMDiS {
cend = end<row>(mat->getBaseMatrix()); cursor != cend; ++cursor) { cend = end<row>(mat->getBaseMatrix()); cursor != cend; ++cursor) {
// Global index of the current row DOF. // Global index of the current row DOF.
int rowIndex = int rowIndex = (*dofMap)[feSpace][*cursor].global + dispRowIndex;
meshDistributor->mapDofToGlobal(feSpace, *cursor) + dispRowIndex;
cols.clear(); cols.clear();
values.clear(); values.clear();
...@@ -238,8 +240,7 @@ namespace AMDiS { ...@@ -238,8 +240,7 @@ namespace AMDiS {
for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor);
icursor != icend; ++icursor) { icursor != icend; ++icursor) {
// Global index of the current column index. // Global index of the current column index.
int colIndex = int colIndex = (*dofMap)[feSpace][col(*icursor)].global + dispColIndex;
meshDistributor->mapDofToGlobal(feSpace, col(*icursor)) + dispColIndex;
// Ignore all zero entries, expect it is a diagonal entry. // Ignore all zero entries, expect it is a diagonal entry.
if (value(*icursor) == 0.0 && rowIndex != colIndex) if (value(*icursor) == 0.0 && rowIndex != colIndex)
...@@ -266,7 +267,7 @@ namespace AMDiS { ...@@ -266,7 +267,7 @@ namespace AMDiS {
// Traverse all used DOFs in the dof vector. // Traverse all used DOFs in the dof vector.
DOFVector<double>::Iterator dofIt(vec, USED_DOFS); DOFVector<double>::Iterator dofIt(vec, USED_DOFS);
for (dofIt.reset(); !dofIt.end(); ++dofIt) { for (dofIt.reset(); !dofIt.end(); ++dofIt) {
int index = meshDistributor->mapDofToGlobal(feSpace, dofIt.getDOFIndex()); int index = (*dofMap)[feSpace][dofIt.getDOFIndex()].global;
double value = *dofIt; double value = *dofIt;
VecSetValues(petscVec, 1, &index, &value, ADD_VALUES); VecSetValues(petscVec, 1, &index, &value, ADD_VALUES);
......
...@@ -24,11 +24,13 @@ namespace AMDiS { ...@@ -24,11 +24,13 @@ namespace AMDiS {
TEST_EXIT_DBG(meshDistributor)("No mesh distributor object defined!\n"); TEST_EXIT_DBG(meshDistributor)("No mesh distributor object defined!\n");
TEST_EXIT_DBG(dofMap)("No parallel mapping object defined!\n"); TEST_EXIT_DBG(dofMap)("No parallel mapping object defined!\n");
TEST_EXIT_DBG(mat)("No DOF matrix defined!\n"); TEST_EXIT_DBG(mat)("No DOF matrix defined!\n");
double wtime = MPI::Wtime(); double wtime = MPI::Wtime();
vector<const FiniteElemSpace*> feSpaces = getFeSpaces(mat); vector<const FiniteElemSpace*> feSpaces = getFeSpaces(mat);
int nRankRows = meshDistributor->getNumberRankDofs(feSpaces); dofMap->update(feSpaces);
int nOverallRows = meshDistributor->getNumberOverallDofs(feSpaces);
int nRankRows = dofMap->getRankDofs();
int nOverallRows = dofMap->getOverallDofs();
// === Create PETSc vector (solution and a temporary vector). === // === Create PETSc vector (solution and a temporary vector). ===
...@@ -86,9 +88,8 @@ namespace AMDiS { ...@@ -86,9 +88,8 @@ namespace AMDiS {
#if (DEBUG != 0) #if (DEBUG != 0)
int a, b; int a, b;
MatGetOwnershipRange(petscMatrix, &a, &b); MatGetOwnershipRange(petscMatrix, &a, &b);
TEST_EXIT(a == meshDistributor->getStartDofs(feSpaces)) TEST_EXIT(a == dofMap->getStartDofs())("Wrong matrix ownership range!\n");
("Wrong matrix ownership range!\n"); TEST_EXIT(b == dofMap->getStartDofs() + nRankRows)
TEST_EXIT(b == meshDistributor->getStartDofs(feSpaces) + nRankRows)
("Wrong matrix ownership range!\n"); ("Wrong matrix ownership range!\n");
#endif #endif
...@@ -131,11 +132,11 @@ namespace AMDiS { ...@@ -131,11 +132,11 @@ namespace AMDiS {
FUNCNAME("PetscSolverGlobalMatrix::fillPetscRhs()"); FUNCNAME("PetscSolverGlobalMatrix::fillPetscRhs()");
TEST_EXIT_DBG(vec)("No DOF vector defined!\n"); TEST_EXIT_DBG(vec)("No DOF vector defined!\n");
TEST_EXIT_DBG(meshDistributor)("No mesh distributor defined!\n"); TEST_EXIT_DBG(dofMap)("No parallel DOF map defined!\n");
vector<const FiniteElemSpace*> feSpaces = getFeSpaces(vec); vector<const FiniteElemSpace*> feSpaces = getFeSpaces(vec);
int nRankRows = meshDistributor->getNumberRankDofs(feSpaces); int nRankRows = dofMap->getRankDofs();
int nOverallRows = meshDistributor->getNumberOverallDofs(feSpaces); int nOverallRows = dofMap->getOverallDofs();
VecCreateMPI(mpiComm, nRankRows, nOverallRows, &petscRhsVec); VecCreateMPI(mpiComm, nRankRows, nOverallRows, &petscRhsVec);
...@@ -194,11 +195,10 @@ namespace AMDiS { ...@@ -194,11 +195,10 @@ namespace AMDiS {
int c = 0; int c = 0;
for (int i = 0; i < nComponents; i++) { for (int i = 0; i < nComponents; i++) {
DOFVector<double> &dv = *(vec.getDOFVector(i)); DOFVector<double> &dv = *(vec.getDOFVector(i));
const FiniteElemSpace *feSpace = dv.getFeSpace(); DofMap& d = (*dofMap)[dv.getFeSpace()].getMap();
int nRankDofs = meshDistributor->getNumberRankDofs(feSpace); for (DofMap::iterator it = d.begin(); it != d.end(); ++it)
if (it->second.local != -1)
for (int j = 0; j < nRankDofs; j++) dv[it->first] = vecPointer[c++];
dv[meshDistributor->mapLocalToDof(feSpace, j)] = vecPointer[c++];
} }
VecRestoreArray(petscSolVec, &vecPointer); VecRestoreArray(petscSolVec, &vecPointer);
...@@ -270,8 +270,8 @@ namespace AMDiS { ...@@ -270,8 +270,8 @@ namespace AMDiS {
const FiniteElemSpace *colFe = mat->getColFeSpace(); const FiniteElemSpace *colFe = mat->getColFeSpace();
// Global index of the current row DOF. // Global index of the current row DOF.
int globalRowDof = int globalRowDof = (*dofMap)[rowFe][*cursor].global;
meshDistributor->mapDofToGlobal(rowFe, *cursor);
// Test if the current row DOF is a periodic DOF. // Test if the current row DOF is a periodic DOF.
bool periodicRow = perMap.isPeriodic(rowFe, globalRowDof); bool periodicRow = perMap.isPeriodic(rowFe, globalRowDof);
...@@ -288,8 +288,7 @@ namespace AMDiS { ...@@ -288,8 +288,7 @@ namespace AMDiS {
icursor != icend; ++icursor) { icursor != icend; ++icursor) {
// Global index of the current column index. // Global index of the current column index.
int globalColDof = int globalColDof = (*dofMap)[colFe][col(*icursor)].global;
meshDistributor->mapDofToGlobal(colFe, col(*icursor));
// Test if the current col dof is a periodic dof. // Test if the current col dof is a periodic dof.
bool periodicCol = perMap.isPeriodic(colFe, globalColDof); bool periodicCol = perMap.isPeriodic(colFe, globalColDof);
// Get PETSc's mat col index. // Get PETSc's mat col index.
...@@ -364,14 +363,12 @@ namespace AMDiS { ...@@ -364,14 +363,12 @@ namespace AMDiS {
for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor);
icursor != icend; ++icursor) { icursor != icend; ++icursor) {
// Global index of the current column index. // Global index of the current column index.
int globalColDof = int globalColDof = (*dofMap)[colFe][col(*icursor)].global;
meshDistributor->mapDofToGlobal(colFe, col(*icursor));
// Ignore all zero entries, expect it is a diagonal entry. // Ignore all zero entries, expect it is a diagonal entry.
if (value(*icursor) == 0.0 && globalRowDof != globalColDof) if (value(*icursor) == 0.0 && globalRowDof != globalColDof)
continue; continue;
// === Add all periodic associations of both, the row and the column === // === Add all periodic associations of both, the row and the column ===
// === indices to the set perAsc. === // === indices to the set perAsc. ===
...@@ -470,12 +467,12 @@ namespace AMDiS { ...@@ -470,12 +467,12 @@ namespace AMDiS {
// Traverse all used DOFs in the dof vector. // Traverse all used DOFs in the dof vector.
DOFVector<double>::Iterator dofIt(vec, USED_DOFS); DOFVector<double>::Iterator dofIt(vec, USED_DOFS);
for (dofIt.reset(); !dofIt.end(); ++dofIt) { for (dofIt.reset(); !dofIt.end(); ++dofIt) {
if (rankOnly && !meshDistributor->getIsRankDof(feSpace, dofIt.getDOFIndex())) if (rankOnly && !(*dofMap)[feSpace].isRankDof(dofIt.getDOFIndex()))
continue; continue;
// Calculate global row index of the DOF. // Calculate global row index of the DOF.
DegreeOfFreedom globalRowDof = DegreeOfFreedom globalRowDof =
meshDistributor->mapDofToGlobal(feSpace, dofIt.getDOFIndex()); (*dofMap)[feSpace][dofIt.getDOFIndex()].global;
// Get PETSc's mat index of the row DOF. // Get PETSc's mat index of the row DOF.
int index = dofToMatIndex.get(nRowVec, globalRowDof); int index = dofToMatIndex.get(nRowVec, globalRowDof);
...@@ -507,8 +504,9 @@ namespace AMDiS { ...@@ -507,8 +504,9 @@ namespace AMDiS {
TEST_EXIT_DBG(!o_nnz)("There is something wrong!\n"); TEST_EXIT_DBG(!o_nnz)("There is something wrong!\n");
vector<const FiniteElemSpace*> feSpaces = getFeSpaces(mat); vector<const FiniteElemSpace*> feSpaces = getFeSpaces(mat);
int nRankRows = meshDistributor->getNumberRankDofs(feSpaces); int nRankRows = dofMap->getRankDofs();
int rankStartIndex = meshDistributor->getStartDofs(feSpaces); int rankStartIndex = dofMap->getStartDofs();
d_nnz = new int[nRankRows]; d_nnz = new int[nRankRows];
o_nnz = new int[nRankRows]; o_nnz = new int[nRankRows];
for (int i = 0; i < nRankRows; i++) { for (int i = 0; i < nRankRows; i++) {
...@@ -576,15 +574,12 @@ namespace AMDiS { ...@@ -576,15 +574,12 @@ namespace AMDiS {
for (cursor_type cursor = begin<row>(bmat), for (cursor_type cursor = begin<row>(bmat),
cend = end<row>(bmat); cursor != cend; ++cursor) { cend = end<row>(bmat); cursor != cend; ++cursor) {
int globalRowDof = (*dofMap)[feSpaces[i]][*cursor].global;
int globalRowDof =
meshDistributor->mapDofToGlobal(feSpaces[i], *cursor);
// The corresponding global matrix row index of the current row DOF. // The corresponding global matrix row index of the current row DOF.
int petscRowIdx = dofToMatIndex.get(i, globalRowDof); int petscRowIdx = dofToMatIndex.get(i, globalRowDof);
if ((*dofMap)[feSpaces[i]].isRankDof(*cursor)) {
if (meshDistributor->getIsRankDof(feSpaces[i], *cursor)) {
// === The current row DOF is a rank DOF, so create the === // === The current row DOF is a rank DOF, so create the ===
// === corresponding nnz values directly on rank's nnz data. === // === corresponding nnz values directly on rank's nnz data. ===
...@@ -593,8 +588,8 @@ namespace AMDiS { ...@@ -593,8 +588,8 @@ namespace AMDiS {
TEST_EXIT_DBG(localPetscRowIdx >= 0 && localPetscRowIdx < nRankRows) TEST_EXIT_DBG(localPetscRowIdx >= 0 && localPetscRowIdx < nRankRows)
("Should not happen! \n Debug info: localRowIdx = %d globalRowIndx = %d petscRowIdx = %d localPetscRowIdx = %d rStart = %d nCompontens = %d nRankRows = %d\n", ("Should not happen! \n Debug info: localRowIdx = %d globalRowIndx = %d petscRowIdx = %d localPetscRowIdx = %d rStart = %d nCompontens = %d nRankRows = %d\n",
*cursor, *cursor,
meshDistributor->mapDofToGlobal(feSpaces[i], *cursor), (*dofMap)[feSpaces[i]][*cursor].global,
petscRowIdx, petscRowIdx,
localPetscRowIdx, localPetscRowIdx,
rankStartIndex, rankStartIndex,
...@@ -605,8 +600,7 @@ namespace AMDiS { ...@@ -605,8 +600,7 @@ namespace AMDiS {
// Traverse all non zero entries in this row. // Traverse all non zero entries in this row.
for (icursor_type icursor = begin<nz>(cursor), for (icursor_type icursor = begin<nz>(cursor),
icend = end<nz>(cursor); icursor != icend; ++icursor) { icend = end<nz>(cursor); icursor != icend; ++icursor) {
int globalColDof = int globalColDof = (*dofMap)[feSpaces[j]][col(*icursor)].global;
meshDistributor->mapDofToGlobal(feSpaces[j], col(*icursor));
int petscColIdx = dofToMatIndex.get(j, globalColDof); int petscColIdx = dofToMatIndex.get(j, globalColDof);
if (value(*icursor) != 0.0 || petscRowIdx == petscColIdx) { if (value(*icursor) != 0.0 || petscRowIdx == petscColIdx) {
...@@ -633,8 +627,8 @@ namespace AMDiS { ...@@ -633,8 +627,8 @@ namespace AMDiS {
for (icursor_type icursor = begin<nz>(cursor), for (icursor_type icursor = begin<nz>(cursor),
icend = end<nz>(cursor); icursor != icend; ++icursor) { icend = end<nz>(cursor); icursor != icend; ++icursor) {
if (value(*icursor) != 0.0) { if (value(*icursor) != 0.0) {
int globalColDof = int globalColDof =
meshDistributor->mapDofToGlobal(feSpaces[j], col(*icursor)); (*dofMap)[feSpaces[j]][col(*icursor)].global;
int petscColIdx = dofToMatIndex.get(j, globalColDof); int petscColIdx = dofToMatIndex.get(j, globalColDof);
sendMatrixEntry[sendToRank]. sendMatrixEntry[sendToRank].
...@@ -649,7 +643,7 @@ namespace AMDiS { ...@@ -649,7 +643,7 @@ namespace AMDiS {
// === Send and recv the nnz row structure to/from other ranks. === // === Send and recv the nnz row structure to/from other ranks. ===
StdMpi<MatrixNnzEntry> stdMpi(meshDistributor->getMpiComm(), true); StdMpi<MatrixNnzEntry> stdMpi(mpiComm, true);
stdMpi.send(sendMatrixEntry); stdMpi.send(sendMatrixEntry);
for (std::set<int>::iterator it = recvFromRank.begin(); for (std::set<int>::iterator it = recvFromRank.begin();
it != recvFromRank.end(); ++it) it != recvFromRank.end(); ++it)
...@@ -697,8 +691,8 @@ namespace AMDiS { ...@@ -697,8 +691,8 @@ namespace AMDiS {
{ {
FUNCNAME("PetscSolverGlobalMatrix::createGlobalDofMapping()"); FUNCNAME("PetscSolverGlobalMatrix::createGlobalDofMapping()");
int offset = meshDistributor->getStartDofs(feSpaces); int offset = dofMap->getStartDofs();
Mesh *mesh = meshDistributor->getMesh(); Mesh *mesh = feSpaces[0]->getMesh();
dofToMatIndex.clear(); dofToMatIndex.clear();
...@@ -709,28 +703,25 @@ namespace AMDiS { ...@@ -709,28 +703,25 @@ namespace AMDiS {
mesh->getAllDofs(feSpaces[i], rankDofSet); mesh->getAllDofs(feSpaces[i], rankDofSet);
for (std::set<const DegreeOfFreedom*>::iterator it = rankDofSet.begin(); for (std::set<const DegreeOfFreedom*>::iterator it = rankDofSet.begin();
it != rankDofSet.end(); ++it) it != rankDofSet.end(); ++it)
if (meshDistributor->getIsRankDof(feSpaces[i], **it)) { if ((*dofMap)[feSpaces[i]].isRankDof(**it)) {
int globalIndex = int globalIndex = (*dofMap)[feSpaces[i]][**it].global;
meshDistributor->mapDofToGlobal(feSpaces[i], **it);
int globalMatIndex = int globalMatIndex =
globalIndex - meshDistributor->getStartDofs(feSpaces[i]) + offset; globalIndex - (*dofMap)[feSpaces[i]].rStartDofs + offset;
dofToMatIndex.add(i, globalIndex, globalMatIndex); dofToMatIndex.add(i, globalIndex, globalMatIndex);
} }
// === Communicate interior boundary DOFs between domains. === // === Communicate interior boundary DOFs between domains. ===
StdMpi<vector<int> > stdMpi(meshDistributor->getMpiComm()); StdMpi<vector<int> > stdMpi(mpiComm);
for (DofComm::Iterator it(meshDistributor->getSendDofs(), feSpaces[i]); for (DofComm::Iterator it(meshDistributor->getSendDofs(), feSpaces[i]);
!it.end(); it.nextRank()) { !it.end(); it.nextRank()) {
vector<DegreeOfFreedom> sendGlobalDofs; vector<DegreeOfFreedom> sendGlobalDofs;
for (; !it.endDofIter(); it.nextDof()) {