diff --git a/AMDiS/src/DOFMatrix.cc b/AMDiS/src/DOFMatrix.cc
index 843e42e5826c3fb41391820ba931b0b3866a2eab..d0f9edfc06ec2dda45bc9673330437f430b384cc 100644
--- a/AMDiS/src/DOFMatrix.cc
+++ b/AMDiS/src/DOFMatrix.cc
@@ -147,8 +147,6 @@ namespace AMDiS {
}
}
-
- // ok
DOFMatrix& DOFMatrix::operator=(const DOFMatrix& rhs)
{
TEST_EXIT(rhs.inserter == 0 && inserter == 0)("Cannot copy during insertion");
@@ -158,6 +156,7 @@ namespace AMDiS {
operators = rhs.operators;
operatorFactor = rhs.operatorFactor;
matrix = rhs.matrix;
+ coupleMatrix = rhs.coupleMatrix;
if (rhs.boundaryManager) {
boundaryManager = NEW BoundaryManager(*rhs.boundaryManager);
} else {
@@ -173,7 +172,6 @@ namespace AMDiS {
return *this;
}
-
void DOFMatrix::addElementMatrix(double sign,
const ElementMatrix &elMat,
const BoundaryType *bound,
@@ -198,7 +196,7 @@ namespace AMDiS {
if (condition && condition->isDirichlet()) {
if (!coupleMatrix)
- ins[row][row]= 1.0;
+ applyDBCs.insert(static_cast<int>(row));
} else
for (int j = 0; j < n_col; j++) { // for all columns
col = elMat.colIndices[j];
@@ -211,7 +209,7 @@ namespace AMDiS {
}
}
- double DOFMatrix::logAcc(DegreeOfFreedom a,DegreeOfFreedom b) const
+ double DOFMatrix::logAcc(DegreeOfFreedom a, DegreeOfFreedom b) const
{
return matrix[a][b];
}
@@ -366,20 +364,9 @@ namespace AMDiS {
return fillFlag;
}
- void DOFMatrix::axpy(double a,
- const DOFMatrix& x,
- const DOFMatrix& y)
+ void DOFMatrix::axpy(double a, const DOFMatrix& x, const DOFMatrix& y)
{
- FUNCNAME("DOFMatrix::axpy");
- TEST_EXIT(x.getRowFESpace() == y.getRowFESpace() &&
- rowFESpace == x.getRowFESpace())
- ("row fe-spaces not equal\n");
- TEST_EXIT(x.getColFESpace() == y.getColFESpace() &&
- colFESpace == x.getColFESpace())
- ("col fe-spaces not equal\n");
-
matrix+= a * x.matrix + y.matrix;
-
}
@@ -402,21 +389,12 @@ namespace AMDiS {
void DOFMatrix::removeRowsWithDBC(std::set<int> *rows)
{
- ERROR_EXIT("TODO: removeRowsWithDBC");
-
-// for (std::set<int>::iterator it = rows->begin();
-// it != rows->end();
-// ++it) {
-// if (coupleMatrix) {
-// matrix[*it].resize(0);
-// } else {
-// matrix[*it].resize(1);
-// matrix[*it][0].col = *it;
-// matrix[*it][0].entry = 1.0;
-// }
-// }
-
-// rows->clear();
+ inserter_type &ins= *inserter;
+
+ for (std::set<int>::iterator it = rows->begin(); it != rows->end(); ++it)
+ ins[*it][*it] = 1.0;
+
+ rows->clear();
}
void DOFMatrix::createPictureFile(const char* filename, int dim)
diff --git a/AMDiS/src/DOFMatrix.h b/AMDiS/src/DOFMatrix.h
index 5321efc2e824c36a909c44ce995a12e63cd7164e..c4c5da5c1098439103cdc3f083c26a2b028b0d03 100644
--- a/AMDiS/src/DOFMatrix.h
+++ b/AMDiS/src/DOFMatrix.h
@@ -256,38 +256,26 @@ namespace AMDiS {
const bool& operator[](int i) const {ERROR_EXIT("Shouldn't be used, only fake."); return dummy;}
#endif
- /** \brief
- * DOFMatrix does not need to be compressed before assembling, when using MTL4.
- */
+ /// DOFMatrix does not need to be compressed before assembling, when using MTL4.
void compressDOFIndexed(int first, int last, std::vector<DegreeOfFreedom> &newDOF) {}
- /** \brief
- * Implements DOFIndexedBase::freeDOFContent()
- */
+ /// Implements DOFIndexedBase::freeDOFContent()
virtual void freeDOFContent(int index);
- /** \brief
- * Returns \ref coupleMatrix.
- */
+ /// Returns \ref coupleMatrix.
inline bool isCoupleMatrix() {
return coupleMatrix;
}
- /** \brief
- * Returns \ref coupleMatrix.
- */
+ /// Returns \ref coupleMatrix.
inline void setCoupleMatrix(bool c) {
coupleMatrix = c;
}
- /** \brief
- * a*x + y
- */
+ /// a * x + y
void axpy(double a, const DOFMatrix& x, const DOFMatrix& y);
- /** \brief
- * Multiplication with a scalar.
- */
+ /// Multiplication with a scalar.
void scal(double s);
/** \brief
@@ -427,7 +415,7 @@ namespace AMDiS {
/// Returns number of rows (\ref matrix.size())
inline int getSize() const {
- return num_rows(matrix);
+ return num_rows(matrix);
}
/** \brief
@@ -438,7 +426,6 @@ namespace AMDiS {
return rowFESpace->getAdmin()->getUsedSize();
}
-
// Only fake, shouldn't be called
/** \brief
* Returns number of cols. For that, the function iteratos over all
@@ -510,6 +497,7 @@ namespace AMDiS {
return boundaryManager;
}
+ /// Returns a pointer to \ref applyDBCs.
std::set<int>* getApplyDBCs() {
return &applyDBCs;
}
@@ -621,7 +609,7 @@ namespace AMDiS {
std::string name;
/// Sparse matrix, type is a template parameter by default compressed2D<double>
- base_matrix_type matrix;
+ base_matrix_type matrix;
/// Used while mesh traversal
static DOFMatrix *traversePtr;
@@ -644,13 +632,22 @@ namespace AMDiS {
///
BoundaryManager *boundaryManager;
- ///
+ /** \brief
+ * If false, the matrix is a diagonal matrix within a matrix of DOF matrices.
+ * Otherwise the value is true, and the matrix is an off-diagonal matrix.
+ */
bool coupleMatrix;
/// Temporary variable used in assemble()
ElementMatrix *elementMatrix;
- ///
+ /* \brief
+ * A set of row indices. When assembling the DOF matrix, all rows, that
+ * correspond to a dof at a dirichlet boundary, are ignored and the row is
+ * left blank. After assembling, the diagonal element of the matrix must be
+ * set to 1. The indices of all rows, where this should be done, are stored
+ * in this set.
+ */
std::set<int> applyDBCs;
#ifdef HAVE_PARALLEL_AMDIS
@@ -659,7 +656,7 @@ namespace AMDiS {
#endif
/// Inserter object: implemented as pointer, allocated and deallocated as needed
- inserter_type *inserter;
+ inserter_type *inserter;
friend class DOFAdmin;
friend class DOFVector<double>;