Commit 68275e39 authored by Thomas Witkowski's avatar Thomas Witkowski
Browse files

* independent meshes now also for 2d possible

parent 8bf583ea
......@@ -84,6 +84,8 @@ namespace AMDiS {
void Assembler::calculateElementMatrix(const ElInfo *rowElInfo,
const ElInfo *colElInfo,
const ElInfo *smallElInfo,
const ElInfo *largeElInfo,
ElementMatrix *userMat,
double factor)
{
......@@ -122,14 +124,24 @@ namespace AMDiS {
ElementMatrix *mat = rememberElMat ? elementMatrix : userMat;
if (secondOrderAssembler)
secondOrderAssembler->calculateElementMatrix(rowElInfo, colElInfo, mat);
if (firstOrderAssemblerGrdPsi)
firstOrderAssemblerGrdPsi->calculateElementMatrix(rowElInfo, colElInfo, mat);
if (firstOrderAssemblerGrdPhi)
firstOrderAssemblerGrdPhi->calculateElementMatrix(rowElInfo, colElInfo, mat);
if (secondOrderAssembler) {
secondOrderAssembler->calculateElementMatrix(rowElInfo, colElInfo,
smallElInfo, largeElInfo, mat);
}
if (firstOrderAssemblerGrdPsi) {
firstOrderAssemblerGrdPsi->calculateElementMatrix(rowElInfo, colElInfo,
smallElInfo, largeElInfo, mat);
}
if (firstOrderAssemblerGrdPhi) {
firstOrderAssemblerGrdPhi->calculateElementMatrix(rowElInfo, colElInfo,
smallElInfo, largeElInfo, mat);
}
if (zeroOrderAssembler) {
zeroOrderAssembler->calculateElementMatrix(rowElInfo, colElInfo, mat);
zeroOrderAssembler->calculateElementMatrix(rowElInfo, colElInfo,
smallElInfo, largeElInfo, mat);
}
if (rememberElMat && userMat) {
......
......@@ -68,7 +68,8 @@ namespace AMDiS {
const FiniteElemSpace *rowFESpace,
const FiniteElemSpace *colFESpace = NULL);
virtual ~Assembler() {};
virtual ~Assembler()
{}
ElementMatrix *initElementMatrix(ElementMatrix *elMat,
const ElInfo *rowElInfo,
......@@ -87,6 +88,8 @@ namespace AMDiS {
void calculateElementMatrix(const ElInfo *rowElInfo,
const ElInfo *colElInfo,
const ElInfo *smallElInfo,
const ElInfo *largeElInfo,
ElementMatrix *userMat,
double factor = 1.0);
......@@ -102,49 +105,49 @@ namespace AMDiS {
*/
inline const FiniteElemSpace* getRowFESpace() {
return rowFESpace;
};
}
/** \brief
* Returns \ref colFESpace.
*/
inline const FiniteElemSpace* getColFESpace() {
return colFESpace;
};
}
/** \brief
* Returns \ref nRow.
*/
inline int getNRow() {
return nRow;
};
}
/** \brief
* Returns \ref nCol.
*/
inline int getNCol() {
return nCol;
};
}
/** \brief
* Sets \ref rememberElMat.
*/
inline void rememberElementMatrix(bool rem) {
rememberElMat = rem;
};
}
/** \brief
* Sets \ref rememberElVec.
*/
inline void rememberElementVector(bool rem) {
rememberElVec = rem;
};
}
/** \brief
* Returns \ref zeroOrderAssembler.
*/
inline ZeroOrderAssembler* getZeroOrderAssembler() {
return zeroOrderAssembler;
};
}
/** \brief
* Returns \ref firstOrderAssemblerGrdPsi or \ref firstOrderAssemblerGrdPhi
......@@ -155,21 +158,21 @@ namespace AMDiS {
return (type == GRD_PSI) ?
firstOrderAssemblerGrdPsi :
firstOrderAssemblerGrdPhi;
};
}
/** \brief
* Returns \ref secondOrderAssembler.
*/
inline SecondOrderAssembler* getSecondOrderAssembler() {
return secondOrderAssembler;
};
}
/** \brief
* Returns \ref operat;
*/
inline Operator* getOperator() {
return operat;
};
}
/** \brief
* Initialisation for the given ElInfo. The call is deligated to
......@@ -223,7 +226,7 @@ namespace AMDiS {
lastVecEl = lastMatEl = NULL;
lastTraverseId = ElInfo::traverseId[omp_get_thread_num()];
}
};
}
/** \brief
* Checks whether quadratures for sub assemblers are already set.
......
......@@ -276,21 +276,21 @@ namespace AMDiS {
*/
inline BasFctType *getPhi(int i) const {
return (*phi)[i];
};
}
/** \brief
* Returns the gradient of the i-th local basis function
*/
inline GrdBasFctType *getGrdPhi(int i) const {
return (*grdPhi)[i];
};
}
/** \brief
* Returns the second derivative of the i-th local basis function
*/
inline D2BasFctType *getD2Phi(int i) const {
return (*d2Phi)[i];
};
}
/** \brief
* Approximates the L2 scalar products of a given function with all basis
......@@ -314,45 +314,53 @@ namespace AMDiS {
*/
virtual void l2ScpFctBas(Quadrature*,
AbstractFunction<double, WorldVector<double> >* /*f*/,
DOFVector<double>* /*fh*/) {};
DOFVector<double>* /*fh*/)
{}
/** \brief
* WorldVector<double> valued l2ScpFctBas function
*/
virtual void l2ScpFctBas(Quadrature* ,
AbstractFunction<WorldVector<double>, WorldVector<double> >* /*f*/,
DOFVector<WorldVector<double> >* /*fh*/) {};
DOFVector<WorldVector<double> >* /*fh*/)
{}
/** \brief
* Interpolates a DOFIndexed<double> after refinement
*/
virtual void refineInter(DOFIndexed<double> *, RCNeighbourList*, int){};
virtual void refineInter(DOFIndexed<double> *, RCNeighbourList*, int)
{}
/** \brief
* Interpolates a DOFIndexed<double> after coarsening
*/
virtual void coarseInter(DOFIndexed<double> *, RCNeighbourList*, int){};
virtual void coarseInter(DOFIndexed<double> *, RCNeighbourList*, int)
{}
/** \brief
* Restricts a DOFIndexed<double> after coarsening
*/
virtual void coarseRestr(DOFIndexed<double> *, RCNeighbourList*, int){};
virtual void coarseRestr(DOFIndexed<double> *, RCNeighbourList*, int)
{}
/** \brief
* Interpolates a DOFVector<WorldVector<double> > after refinement
*/
virtual void refineInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int){};
virtual void refineInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
{}
/** \brief
* Interpolates a DOFVector<WorldVector<double> > after coarsening
*/
virtual void coarseInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int){};
virtual void coarseInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
{}
/** \brief
* Restricts a DOFVector<WorldVector<double> > after coarsening
*/
virtual void coarseRestr(DOFVector<WorldVector<double> >*, RCNeighbourList*, int){};
virtual void coarseRestr(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
{}
/** \brief
* Returns local dof indices of the element for the given fe space.
......@@ -362,7 +370,7 @@ namespace AMDiS {
DegreeOfFreedom*) const
{
return NULL;
};
}
/** \brief
* Returns local dof indices of the element for the given fe space.
......@@ -370,7 +378,7 @@ namespace AMDiS {
virtual void getLocalIndicesVec(const Element*,
const DOFAdmin*,
Vector<DegreeOfFreedom>*) const
{};
{}
/** \brief
......
......@@ -250,14 +250,6 @@ namespace AMDiS {
bound ? boundaryManager->getBoundaryCondition(bound[i]) : NULL;
if (condition && condition->isDirichlet()) {
/* MatrixRow *matrixRow = &(matrix[row]);
if (coupleMatrix) {
matrixRow->resize(0);
} else {
matrixRow->resize(1);
((*matrixRow)[0]).col = row;
((*matrixRow)[0]).entry = 1.0;
} */
applyDBCs.insert(static_cast<int>(row));
} else {
for (int j = 0; j < nCol; j++) { // for all columns
......@@ -312,7 +304,7 @@ namespace AMDiS {
("Column index %d out of range 0-%d\n", jcol, colFESpace->getAdmin()->getUsedSize() - 1);
// first entry is diagonal entry
if ((rowSize == 0) && (rowFESpace == colFESpace)) {
if ((rowSize == 0) && (colFESpace == rowFESpace)) {
MatEntry newEntry = {irow, 0.0};
row->push_back(newEntry);
rowSize = 1;
......@@ -458,11 +450,13 @@ namespace AMDiS {
*factorIt ? **factorIt : 1.0);
}
}
addElementMatrix(factor, *elementMatrix, bound);
}
void DOFMatrix::assemble(double factor, ElInfo *rowElInfo, ElInfo *colElInfo,
void DOFMatrix::assemble(double factor,
ElInfo *rowElInfo, ElInfo *colElInfo,
ElInfo *smallElInfo, ElInfo *largeElInfo,
const BoundaryType *bound, Operator *op)
{
FUNCNAME("DOFMatrix::assemble()");
......@@ -476,19 +470,30 @@ namespace AMDiS {
initElementMatrix(elementMatrix, rowElInfo, colElInfo);
if (op) {
op->getElementMatrix(rowElInfo, colElInfo, elementMatrix);
op->getElementMatrix(rowElInfo, colElInfo,
smallElInfo, largeElInfo,
elementMatrix);
} else {
std::vector<Operator*>::iterator it;
std::vector<double*>::iterator factorIt;
for (it = operators.begin(), factorIt = operatorFactor.begin();
it != operators.end();
++it, ++factorIt) {
(*it)->getElementMatrix(rowElInfo, colElInfo,
elementMatrix,
*factorIt ? **factorIt : 1.0);
// If both levels are equal, we may use the standard assembler for
// one element.
if (rowElInfo->getLevel() == colElInfo->getLevel()) {
(*it)->getElementMatrix(rowElInfo,
elementMatrix,
*factorIt ? **factorIt : 1.0);
} else {
(*it)->getElementMatrix(rowElInfo, colElInfo,
smallElInfo, largeElInfo,
elementMatrix,
*factorIt ? **factorIt : 1.0);
}
}
}
addElementMatrix(factor, *elementMatrix, bound);
}
......@@ -671,6 +676,23 @@ namespace AMDiS {
}
}
int DOFMatrix::getNumCols() const {
int max = 0;
int rowSize = static_cast<int>(matrix.size());
for (int i = 0; i < rowSize; i++) {
int colSize = static_cast<int>(matrix[i].size());
for (int j = 0; j < colSize; j++) {
if (matrix[i][j].col > max) {
max = matrix[i][j].col;
}
}
}
// Add one to the maximum value, because the indeces start with 0.
return max + 1;
}
void DOFMatrix::createPictureFile(const char* filename, int dim)
{
png_structp png_ptr = png_create_write_struct(PNG_LIBPNG_VER_STRING, NULL, NULL, NULL);
......
......@@ -88,11 +88,11 @@ namespace AMDiS {
public:
MatEntryValueLess(const double& value)
: value_(value)
{};
{}
bool operator()(const MatEntry& m) {
return (fabs(m.entry) < value_);
};
}
};
......@@ -106,11 +106,11 @@ namespace AMDiS {
public:
MatEntryColLess(const int& col)
: col_(col)
{};
{}
bool operator()(const MatEntry& m) {
return (fabs(m.col) < col_);
};
}
};
......@@ -121,7 +121,7 @@ namespace AMDiS {
struct CmpMatEntryCol : public std::binary_function<MatEntry, MatEntry, bool> {
bool operator()(const MatEntry& m1, const MatEntry m2) const {
return m1.col < m2.col;
};
}
};
......@@ -132,7 +132,7 @@ namespace AMDiS {
struct CmpMatEntryValue : public std::binary_function<MatEntry, MatEntry, bool> {
bool operator()(const MatEntry& m1, const MatEntry m2) const {
return m1.entry < m2.entry;
};
}
};
......@@ -143,7 +143,7 @@ namespace AMDiS {
struct CmpMatEntryAbsValueLess : public std::binary_function<MatEntry, MatEntry, bool> {
bool operator()(const MatEntry& m1, const MatEntry m2) const {
return fabs(m1.entry) < fabs(m2.entry);
};
}
};
/** \brief
......@@ -153,7 +153,7 @@ namespace AMDiS {
struct CmpMatEntryAbsValueGreater : public std::binary_function<MatEntry, MatEntry, bool> {
bool operator()(const MatEntry& m1, const MatEntry m2) const {
return fabs(m1.entry) > fabs(m2.entry);
};
}
};
// ============================================================================
......@@ -183,7 +183,7 @@ namespace AMDiS {
Iterator(DOFIndexed< std::vector<MatEntry> > *c,
DOFIteratorType type)
: DOFIterator< std::vector<MatEntry> >(c, type)
{};
{}
};
/** \brief
......@@ -237,14 +237,14 @@ namespace AMDiS {
*/
std::vector< std::vector<MatEntry> >::iterator begin() {
return matrix.begin();
};
}
/** \brief
* Returns an iterator to the end of matrix rows (\ref matrix.end())
*/
std::vector< std::vector<MatEntry> >::iterator end() {
return matrix.end();
};
}
/** \brief
* Used by DOFAdmin to compress the DOFMatrix
......@@ -262,28 +262,28 @@ namespace AMDiS {
*/
inline bool entryUsed(DegreeOfFreedom a,int b) const {
return (((matrix[a])[b]).col >= 0);
};
}
/** \brief
* 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);
};
}
/** \brief
* Returns \ref coupleMatrix.
*/
inline bool isCoupleMatrix() {
return coupleMatrix;
};
}
/** \brief
* Returns \ref coupleMatrix.
*/
inline void setCoupleMatrix(bool c) {
coupleMatrix = c;
};
}
/** \brief
* Matrix-matrix multiplication.
......@@ -310,31 +310,31 @@ namespace AMDiS {
operators.push_back(op);
operatorFactor.push_back(factor);
operatorEstFactor.push_back(estFactor);
};
}
inline std::vector<double*>::iterator getOperatorFactorBegin() {
return operatorFactor.begin();
};
}
inline std::vector<double*>::iterator getOperatorFactorEnd() {
return operatorFactor.end();
};
}
inline std::vector<double*>::iterator getOperatorEstFactorBegin() {
return operatorEstFactor.begin();
};
}
inline std::vector<double*>::iterator getOperatorEstFactorEnd() {
return operatorEstFactor.end();
};
}
inline std::vector<Operator*>::iterator getOperatorsBegin() {
return operators.begin();
};
}
inline std::vector<Operator*>::iterator getOperatorsEnd() {
return operators.end();
};
}
Flag getAssembleFlag();
......@@ -365,7 +365,9 @@ namespace AMDiS {
const BoundaryType *bound, Operator *op = NULL);
void assemble(double factor, ElInfo *rowElInfo, ElInfo *colElInfo,
void assemble(double factor,
ElInfo *rowElInfo, ElInfo *colElInfo,
ElInfo *smallElInfo, ElInfo *largeElInfo,
const BoundaryType *bound, Operator *op = NULL);
......@@ -382,25 +384,25 @@ namespace AMDiS {
*/
std::vector< std::vector<MatEntry> >& getMatrix() {
return matrix;
};
}
void setMatrix(std::vector< std::vector<MatEntry> > m) {
matrix = m;
};
}
/** \brief
* Returns \ref matrix[n]
*/
const std::vector<MatEntry>& getRow(int n) const {
return matrix[n];
};
}
/** \brief
* Returns \ref matrix[n]
*/
std::vector<MatEntry>& getRow(int n) {
return matrix[n];
};
}
/** \brief
* Returns whether restriction should be performed after coarsening
......@@ -408,35 +410,35 @@ namespace AMDiS {
*/
virtual bool coarseRestrict() {
return false;
};
}
/** \brief
* Returns const \ref rowFESpace
*/
const FiniteElemSpace* getRowFESpace() const {
return rowFESpace;
};
}
/** \brief
* Returns const \ref colFESpace
*/
const FiniteElemSpace* getColFESpace() const {
return colFESpace;
};
}
/** \brief
* Returns const \ref rowFESpace
*/
const FiniteElemSpace* getFESpace() const {
return rowFESpace;
};
}
/** \brief
* Returns number of rows (\ref matrix.size())
*/
inline int getSize() const {
return matrix.size();
};
}
/** \brief
* Returns the number of used rows (equal to number of used DOFs in
......@@ -444,34 +446,20 @@ namespace AMDiS {
*/
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.
*/
int getNumCols() const {
int max = 0;
int rowSize = static_cast<int>(matrix.size());
for (int i = 0; i < rowSize; i++) {
int colSize = static_cast<int>(matrix[i].size());
for (int j = 0; j < colSize; j++) {
if (matrix[i][j].col > max) {
max = matrix[i][j].col;
}
}
}
// Add one to the maximum value, because the indeces start with 0.
return max + 1;
};
int getNumCols() const;
/** \brief
* Returns \ref name
*/
inline const std::string& getName() const {
return name;
};
}
/** \brief
* Resizes \ref matrix to n rows
......@@ -479,7 +467,7 @@ namespace AMDiS {
inline void resize(int n) {
TEST_EXIT_DBG(n >= 0)("Can't resize DOFMatrix to negative size\n");