From ec239f4ea7369a25118869fe7d3ccda60b77d78e Mon Sep 17 00:00:00 2001
From: Thomas Witkowski <thomas.witkowski@gmx.de>
Date: Thu, 4 Sep 2008 07:29:47 +0000
Subject: [PATCH] * This and that
---
AMDiS/compositeFEM/src/SubElementAssembler.cc | 4 +-
AMDiS/compositeFEM/src/SubElementAssembler.h | 4 +-
AMDiS/src/DOFMatrix.cc | 18 +-
AMDiS/src/DOFMatrix.h | 23 ++-
AMDiS/src/ProblemVec.cc | 157 +++++++++++++-----
AMDiS/src/ProblemVec.h | 11 ++
6 files changed, 162 insertions(+), 55 deletions(-)
diff --git a/AMDiS/compositeFEM/src/SubElementAssembler.cc b/AMDiS/compositeFEM/src/SubElementAssembler.cc
index 883f8c47..95f495f5 100644
--- a/AMDiS/compositeFEM/src/SubElementAssembler.cc
+++ b/AMDiS/compositeFEM/src/SubElementAssembler.cc
@@ -74,8 +74,6 @@ namespace AMDiS {
const ElInfo *elInfo,
ElementVector *userVec)
{
- double corrFactor;
-
/**
* Manipulate the quadratures of the SubAssemblers for subelement.
*/
@@ -88,7 +86,7 @@ namespace AMDiS {
* determinant of element has been used instead of the determinant of subelement. Thus
* the result must be corrected with respect to subelement.
*/
- corrFactor = subElInfo->getDet() / fabs(elInfo->getDet());
+ double corrFactor = subElInfo->getDet() / fabs(elInfo->getDet());
for (int i = 0; i < nRow; i++) {
(*userVec)[i] *= corrFactor;
}
diff --git a/AMDiS/compositeFEM/src/SubElementAssembler.h b/AMDiS/compositeFEM/src/SubElementAssembler.h
index e5697592..991b4b18 100644
--- a/AMDiS/compositeFEM/src/SubElementAssembler.h
+++ b/AMDiS/compositeFEM/src/SubElementAssembler.h
@@ -57,8 +57,8 @@ namespace AMDiS {
MEMORY_MANAGED(SubElementAssembler);
SubElementAssembler(Operator *op,
- const FiniteElemSpace *rowFESpace_,
- const FiniteElemSpace *colFESpace_=NULL);
+ const FiniteElemSpace *rowFESpace,
+ const FiniteElemSpace *colFESpace = NULL);
virtual ~SubElementAssembler()
{
diff --git a/AMDiS/src/DOFMatrix.cc b/AMDiS/src/DOFMatrix.cc
index 7a5c8f8a..aac17249 100644
--- a/AMDiS/src/DOFMatrix.cc
+++ b/AMDiS/src/DOFMatrix.cc
@@ -40,7 +40,7 @@ namespace AMDiS {
}
if (rowFESpace && rowFESpace->getAdmin())
- (const_cast<DOFAdmin*>( rowFESpace->getAdmin()))->addDOFIndexed(this);
+ (const_cast<DOFAdmin*>(rowFESpace->getAdmin()))->addDOFIndexed(this);
boundaryManager = NEW BoundaryManager;
elementMatrix = NEW ElementMatrix(rowFESpace->getBasisFcts()->getNumber(),
@@ -284,9 +284,12 @@ namespace AMDiS {
{
int j;
- for (j=0; j<static_cast<int>(matrix[a].size());j++)
- if (b==matrix[a][j].col) break;
- if (j!=static_cast<int>(matrix[a].size())) matrix[a][j].col=c;
+ for (j = 0; j<static_cast<int>(matrix[a].size());j++)
+ if (b == matrix[a][j].col)
+ break;
+
+ if (j != static_cast<int>(matrix[a].size()))
+ matrix[a][j].col = c;
}
void DOFMatrix::addSparseDOFEntry(double sign, int irow,
@@ -386,12 +389,13 @@ namespace AMDiS {
}
}
int usedSize = rowFESpace->getAdmin()->getUsedSize();
- for(i=0; i < usedSize; i++) {
+ for (i = 0; i < usedSize; i++) {
row = reinterpret_cast< std::vector<MatEntry>*>(&(matrix[i]));
int rowSize = static_cast<int>(row->size());
- for (j=0; j < rowSize; j++) {
+ for (j = 0; j < rowSize; j++) {
col = (*row)[j].col;
- if (entryUsed(i,j)) (*row)[j].col = newDOF[col];
+ if (entryUsed(i,j))
+ (*row)[j].col = newDOF[col];
}
}
}
diff --git a/AMDiS/src/DOFMatrix.h b/AMDiS/src/DOFMatrix.h
index c331b2c7..deb7a3ac 100644
--- a/AMDiS/src/DOFMatrix.h
+++ b/AMDiS/src/DOFMatrix.h
@@ -211,7 +211,7 @@ namespace AMDiS {
*/
DOFMatrix(const FiniteElemSpace* rowFESpace,
const FiniteElemSpace* colFESpace,
- std::string n="");
+ std::string n = "");
/** \brief
* Copy-Constructor
@@ -267,18 +267,22 @@ namespace AMDiS {
* Returns true if matrix[a][b] has entry \ref NO_MORE_ENTRIES
*/
inline bool noMoreEntries(DegreeOfFreedom a,int b) const {
- return (NO_MORE_ENTRIES==((matrix[a])[b]).col);
+ return (NO_MORE_ENTRIES == ((matrix[a])[b]).col);
};
/** \brief
* Returns \ref coupleMatrix.
*/
- inline bool isCoupleMatrix() { return coupleMatrix; };
+ inline bool isCoupleMatrix() {
+ return coupleMatrix;
+ };
/** \brief
* Returns \ref coupleMatrix.
*/
- inline void setCoupleMatrix(bool c) { coupleMatrix = c; };
+ inline void setCoupleMatrix(bool c) {
+ coupleMatrix = c;
+ };
/** \brief
* Matrix-matrix multiplication.
@@ -437,6 +441,14 @@ namespace AMDiS {
return matrix.size();
};
+ /** \brief
+ * Returns the number of used rows (equal to number of used DOFs in
+ * the row FE space).
+ */
+ inline int getUsedSize() const {
+ return rowFESpace->getAdmin()->getUsedSize();
+ };
+
/** \brief
* Returns number of cols. For that, the function iteratos over all
* rows and searchs for the entry with the highest col number.
@@ -453,7 +465,8 @@ namespace AMDiS {
}
}
- return max;
+ // Add one to the maximum value, because the indeces start with 0.
+ return max + 1;
};
/** \brief
diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc
index ea9ad47a..3ce3acd3 100644
--- a/AMDiS/src/ProblemVec.cc
+++ b/AMDiS/src/ProblemVec.cc
@@ -12,6 +12,7 @@
#include "FileWriter.h"
#include "CoarseningManager.h"
#include "RefinementManager.h"
+#include "DualTraverse.h"
#include "Mesh.h"
#include "OEMSolver.h"
#include "Preconditioner.h"
@@ -271,8 +272,6 @@ namespace AMDiS {
{
FUNCNAME("ProblemVec::createMatricesAndVectors()");
- int i;
-
// === create vectors and system matrix ===
systemMatrix_ = NEW Matrix<DOFMatrix*>(nComponents, nComponents);
@@ -282,7 +281,7 @@ namespace AMDiS {
char number[10];
std::string numberedName;
- for (i = 0; i < nComponents; i++) {
+ for (int i = 0; i < nComponents; i++) {
(*systemMatrix_)[i][i] = NEW DOFMatrix(componentSpaces[i],
componentSpaces[i], "A_ii");
(*systemMatrix_)[i][i]->setCoupleMatrix(false);
@@ -700,8 +699,6 @@ namespace AMDiS {
#pragma omp parallel for
#endif
for (i = 0; i < nComponents; i++) {
- const BasisFunction *basisFcts = componentSpaces[i]->getBasisFcts();
-
for (int j = 0; j < nComponents; j++) {
// Only if this variable is true, the current matrix will be assembled.
bool assembleMatrix = true;
@@ -729,43 +726,21 @@ namespace AMDiS {
if (assembleMatrix && matrix->getBoundaryManager())
matrix->getBoundaryManager()->initMatrix(matrix);
- BoundaryType *bound = NULL;
- if (useGetBound_) {
- bound = GET_MEMORY(BoundaryType, basisFcts->getNumber());
- }
-
- TraverseStack stack;
- ElInfo *elInfo = stack.traverseFirst(componentMeshes[i], -1, assembleFlag);
-
- while (elInfo) {
- if (useGetBound_) {
- basisFcts->getBound(elInfo, bound);
- }
-
- if (assembleMatrix) {
- matrix->assemble(1.0, elInfo, bound);
-
- if (matrix->getBoundaryManager()) {
- matrix->
- getBoundaryManager()->
- fillBoundaryConditions(elInfo, matrix);
- }
- }
-
- if (i == j) {
- rhs_->getDOFVector(i)->assemble(1.0, elInfo, bound);
- }
-
- elInfo = stack.traverseNext(elInfo);
+ if (componentSpaces[i] == componentSpaces[j]) {
+ assembleOnOneMesh(componentSpaces[i],
+ assembleFlag,
+ assembleMatrix ? matrix : NULL,
+ (i == j) ? rhs_->getDOFVector(i) : NULL);
+ } else {
+ assembleOnDifMeshes(componentSpaces[i], componentSpaces[j],
+ assembleFlag,
+ assembleMatrix ? matrix : NULL,
+ (i == j) ? rhs_->getDOFVector(i) : NULL);
}
-
+
if (assembleMatrix && matrix->getBoundaryManager())
matrix->getBoundaryManager()->exitMatrix(matrix);
- if (useGetBound_) {
- FREE_MEMORY(bound, BoundaryType, basisFcts->getNumber());
- }
-
assembledMatrix_[i][j] = true;
}
@@ -956,6 +931,112 @@ namespace AMDiS {
addBoundaryCondition(periodic);
}
+ void ProblemVec::assembleOnOneMesh(FiniteElemSpace *feSpace, Flag assembleFlag,
+ DOFMatrix *matrix, DOFVector<double> *vector)
+ {
+ Mesh *mesh = feSpace->getMesh();
+ const BasisFunction *basisFcts = feSpace->getBasisFcts();
+
+ BoundaryType *bound = NULL;
+ if (useGetBound_) {
+ bound = GET_MEMORY(BoundaryType, basisFcts->getNumber());
+ }
+
+ TraverseStack stack;
+ ElInfo *elInfo = stack.traverseFirst(mesh, -1, assembleFlag);
+
+ while (elInfo) {
+ if (useGetBound_) {
+ basisFcts->getBound(elInfo, bound);
+ }
+
+ if (matrix) {
+ matrix->assemble(1.0, elInfo, bound);
+
+ if (matrix->getBoundaryManager()) {
+ matrix->getBoundaryManager()->
+ fillBoundaryConditions(elInfo, matrix);
+ }
+ }
+
+ if (vector) {
+ vector->assemble(1.0, elInfo, bound);
+ }
+
+ elInfo = stack.traverseNext(elInfo);
+ }
+
+ if (useGetBound_) {
+ FREE_MEMORY(bound, BoundaryType, basisFcts->getNumber());
+ }
+ }
+
+ void ProblemVec::assembleOnDifMeshes(FiniteElemSpace *rowFeSpace, FiniteElemSpace *colFeSpace,
+ Flag assembleFlag,
+ DOFMatrix *matrix, DOFVector<double> *vector)
+ {
+ Mesh *rowMesh = rowFeSpace->getMesh();
+ Mesh *colMesh = colFeSpace->getMesh();
+
+ const BasisFunction *basisFcts = rowFeSpace->getBasisFcts();
+ BoundaryType *bound = NULL;
+ if (useGetBound_) {
+ bound = GET_MEMORY(BoundaryType, basisFcts->getNumber());
+ }
+
+ DualTraverse dualTraverse;
+ ElInfo *rowElInfo, *colElInfo;
+ ElInfo *largeElInfo, *smallElInfo;
+
+ bool cont = dualTraverse.traverseFirst(rowMesh, colMesh, -1, -1,
+ assembleFlag, assembleFlag,
+ &rowElInfo, &colElInfo,
+ &smallElInfo, &largeElInfo);
+ while (cont) {
+ Element *rowElem = rowElInfo->getElement();
+ Element *colElem = colElInfo->getElement();
+
+ if (rowElInfo->getLevel() != colElInfo->getLevel()) {
+ if (smallElInfo == rowElInfo)
+ std::cout << "Row = small\n";
+ else
+ std::cout << "Row = large\n";
+
+ if (largeElInfo == colElInfo)
+ std::cout << "Col = large\n";
+ else
+ std::cout << "Col = small\n";
+
+
+ ERROR_EXIT("KOMMT GLEICH!\n");
+ }
+
+ if (useGetBound_) {
+ basisFcts->getBound(rowElInfo, bound);
+ }
+
+ if (matrix) {
+ matrix->assemble(1.0, rowElInfo, bound);
+
+ if (matrix->getBoundaryManager()) {
+ matrix->getBoundaryManager()->
+ fillBoundaryConditions(rowElInfo, matrix);
+ }
+ }
+
+ if (vector) {
+ vector->assemble(1.0, rowElInfo, bound);
+ }
+
+ cont = dualTraverse.traverseNext(&rowElInfo, &colElInfo,
+ &smallElInfo, &largeElInfo);
+ }
+
+ if (useGetBound_) {
+ FREE_MEMORY(bound, BoundaryType, basisFcts->getNumber());
+ }
+ }
+
void ProblemVec::writeResidualMesh(AdaptInfo *adaptInfo, const std::string name)
{
FUNCNAME("ProblemVec::writeResidualMesh()");
diff --git a/AMDiS/src/ProblemVec.h b/AMDiS/src/ProblemVec.h
index 6dca6598..6d510825 100644
--- a/AMDiS/src/ProblemVec.h
+++ b/AMDiS/src/ProblemVec.h
@@ -289,6 +289,17 @@ namespace AMDiS {
inline void allowFirstRefinement() {
allowFirstRef_ = true;
};
+
+ /** \brief
+ * This function assembles a DOFMatrix and a DOFVector for the case,
+ * the meshes from row and col FE-space are equal.
+ */
+ void assembleOnOneMesh(FiniteElemSpace *feSpace, Flag assembleFlag,
+ DOFMatrix *matrix, DOFVector<double> *vector);
+
+ void assembleOnDifMeshes(FiniteElemSpace *rowFeSpace, FiniteElemSpace *colFeSpace,
+ Flag assembleFlag,
+ DOFMatrix *matrix, DOFVector<double> *vector);
// ===== getting-methods ======================================================
--
GitLab