Skip to content

GitLab

Commit e53e8cf3 authored by Thomas Witkowski's avatar Thomas Witkowski
Browse files

* Some more small changes

parent 69dadcaa
Hide whitespace changes
Inline Side-by-side
Showing with 167 additions and 254 deletions
AMDiS/compositeFEM/src/CFE_Integration.cc
View file @ e53e8cf3
......@@ -288,7 +288,7 @@ namespace AMDiS {
}
// calculate determinant for surface
surfDet = ElInfo::calcDet(worldCoords);
surfDet = loc_elInfo->calcDet(worldCoords);
return surfDet;
}
......
AMDiS/src/Assembler.cc
View file @ e53e8cf3
......@@ -11,12 +11,12 @@
namespace AMDiS {
Assembler::Assembler(Operator *op,
const FiniteElemSpace *rowFESpace_,
const FiniteElemSpace *colFESpace_)
Assembler::Assembler(Operator *op,
const FiniteElemSpace *row,
const FiniteElemSpace *col)
: operat(op),
rowFESpace(rowFESpace_),
colFESpace(colFESpace_ ? colFESpace_ : rowFESpace_),
rowFESpace(row),
colFESpace(col ? col : row),
nRow(rowFESpace->getBasisFcts()->getNumber()),
nCol(colFESpace->getBasisFcts()->getNumber()),
remember(true),
......@@ -28,7 +28,19 @@ namespace AMDiS {
lastVecEl(NULL),
lastTraverseId(-1)
{}
{
elementMatrix = NEW ElementMatrix(nRow, nCol);
elementVector = NEW ElementVector(nRow);
}
Assembler::~Assembler()
{
if (elementMatrix)
DELETE elementMatrix;
if (elementVector)
DELETE elementVector;
}
void Assembler::calculateElementMatrix(const ElInfo *elInfo,
ElementMatrix *userMat,
......@@ -40,11 +52,6 @@ namespace AMDiS {
rememberElMat = true;
}
if (rememberElMat && !elementMatrix)
elementMatrix = NEW ElementMatrix(nRow, nCol);
checkQuadratures();
Element *el = elInfo->getElement();
if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) {
......@@ -93,15 +100,9 @@ namespace AMDiS {
rememberElMat = true;
}
if (rememberElMat && !elementMatrix)
elementMatrix = NEW ElementMatrix(nRow, nCol);
Element *el = smallElInfo->getElement();
Element *el = smallElInfo->getElement();
lastVecEl = lastMatEl = NULL;
checkQuadratures();
if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) {
initElement(smallElInfo, largeElInfo);
}
......@@ -148,12 +149,7 @@ namespace AMDiS {
rememberElVec = true;
}
if (rememberElVec && !elementVector) {
elementVector = NEW ElementVector(nRow);
}
Element *el = elInfo->getElement();
checkQuadratures();
if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) {
initElement(elInfo);
......@@ -202,12 +198,7 @@ namespace AMDiS {
rememberElVec = true;
}
if (rememberElVec && !elementVector) {
elementVector = NEW ElementVector(nRow);
}
Element *el = mainElInfo->getElement();
checkQuadratures();
if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) {
initElement(auxElInfo);
......@@ -333,8 +324,6 @@ namespace AMDiS {
const ElInfo *largeElInfo,
Quadrature *quad)
{
checkQuadratures();
if (secondOrderAssembler)
secondOrderAssembler->initElement(smallElInfo, largeElInfo, quad);
if (firstOrderAssemblerGrdPsi)
......@@ -345,50 +334,6 @@ namespace AMDiS {
zeroOrderAssembler->initElement(smallElInfo, largeElInfo, quad);
}
OptimizedAssembler::OptimizedAssembler(Operator *op,
Quadrature *quad2,
Quadrature *quad1GrdPsi,
Quadrature *quad1GrdPhi,
Quadrature *quad0,
const FiniteElemSpace *rowFESpace_,
const FiniteElemSpace *colFESpace_)
: Assembler(op, rowFESpace_, colFESpace_)
{
bool opt = (rowFESpace_ == colFESpace_);
// create sub assemblers
secondOrderAssembler =
SecondOrderAssembler::getSubAssembler(op, this, quad2, opt);
firstOrderAssemblerGrdPsi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPsi, GRD_PSI, opt);
firstOrderAssemblerGrdPhi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPhi, GRD_PHI, opt);
zeroOrderAssembler =
ZeroOrderAssembler::getSubAssembler(op, this, quad0, opt);
}
StandardAssembler::StandardAssembler(Operator *op,
Quadrature *quad2,
Quadrature *quad1GrdPsi,
Quadrature *quad1GrdPhi,
Quadrature *quad0,
const FiniteElemSpace *rowFESpace_,
const FiniteElemSpace *colFESpace_)
: Assembler(op, rowFESpace_, colFESpace_)
{
remember = false;
// create sub assemblers
secondOrderAssembler =
SecondOrderAssembler::getSubAssembler(op, this, quad2, false);
firstOrderAssemblerGrdPsi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPsi, GRD_PSI, false);
firstOrderAssemblerGrdPhi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPhi, GRD_PHI, false);
zeroOrderAssembler =
ZeroOrderAssembler::getSubAssembler(op, this, quad0, false);
}
void Assembler::initElementMatrix(ElementMatrix *elMat,
const ElInfo *rowElInfo,
const ElInfo *colElInfo)
......@@ -473,4 +418,53 @@ namespace AMDiS {
lastVecEl = NULL;
lastMatEl = NULL;
}
OptimizedAssembler::OptimizedAssembler(Operator *op,
Quadrature *quad2,
Quadrature *quad1GrdPsi,
Quadrature *quad1GrdPhi,
Quadrature *quad0,
const FiniteElemSpace *row,
const FiniteElemSpace *col)
: Assembler(op, row, col)
{
bool opt = (row == col);
// create sub assemblers
secondOrderAssembler =
SecondOrderAssembler::getSubAssembler(op, this, quad2, opt);
firstOrderAssemblerGrdPsi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPsi, GRD_PSI, opt);
firstOrderAssemblerGrdPhi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPhi, GRD_PHI, opt);
zeroOrderAssembler =
ZeroOrderAssembler::getSubAssembler(op, this, quad0, opt);
checkQuadratures();
}
StandardAssembler::StandardAssembler(Operator *op,
Quadrature *quad2,
Quadrature *quad1GrdPsi,
Quadrature *quad1GrdPhi,
Quadrature *quad0,
const FiniteElemSpace *row,
const FiniteElemSpace *col)
: Assembler(op, row, col)
{
remember = false;
// create sub assemblers
secondOrderAssembler =
SecondOrderAssembler::getSubAssembler(op, this, quad2, false);
firstOrderAssemblerGrdPsi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPsi, GRD_PSI, false);
firstOrderAssemblerGrdPhi =
FirstOrderAssembler::getSubAssembler(op, this, quad1GrdPhi, GRD_PHI, false);
zeroOrderAssembler =
ZeroOrderAssembler::getSubAssembler(op, this, quad0, false);
checkQuadratures();
}
}
AMDiS/src/Assembler.h
View file @ e53e8cf3
......@@ -60,12 +60,13 @@ namespace AMDiS {
public:
MEMORY_MANAGED(Assembler);
/// Constructor.
/// Constructor
Assembler(Operator *op,
const FiniteElemSpace *rowFESpace,
const FiniteElemSpace *colFESpace = NULL);
virtual ~Assembler() {}
/// Destructor
~Assembler();
void initElementMatrix(ElementMatrix *elMat,
const ElInfo *rowElInfo,
......@@ -200,12 +201,13 @@ namespace AMDiS {
void matVecAssemble(const ElInfo *elInfo, ElementVector *vec);
///
void matVecAssemble(const ElInfo *mainElInfo, const ElInfo *auxElInfo,
const ElInfo *smallElInfo, const ElInfo *largeElInfo,
ElementVector *vec);
/** \brief
* Checks whether quadratures for sub assemblers are already set.
* Checks whether quadratures for subassemblers are already set.
* If not they will be created.
*/
void checkQuadratures();
......
AMDiS/src/DOFMatrix.cc
View file @ e53e8cf3
......@@ -443,7 +443,7 @@ namespace AMDiS {
std::vector<double*>::iterator factorIt = operatorFactor.begin();
for (; it != operators.end(); ++it, ++factorIt) {
if ((*it)->getNeedDualTraverse() == false) {
(*it)->getElementMatrix(elInfo,
(*it)->getElementMatrix(elInfo,
elementMatrix,
*factorIt ? **factorIt : 1.0);
}
......@@ -460,7 +460,7 @@ namespace AMDiS {
TEST_EXIT_DBG(op)("No operator!\n");
op->getAssembler(omp_get_thread_num())->initElementMatrix(elementMatrix, elInfo);
op->getElementMatrix(elInfo, elementMatrix);
op->getElementMatrix(elInfo, elementMatrix, factor);
addElementMatrix(factor, *elementMatrix, bound);
}
......
AMDiS/src/ElInfo.cc
View file @ e53e8cf3
......@@ -78,7 +78,7 @@ namespace AMDiS {
return calcDet(coord_);
}
double ElInfo::calcDet(const FixVec<WorldVector<double>, VERTEX> &coords)
double ElInfo::calcDet(const FixVec<WorldVector<double>, VERTEX> &coords) const
{
FUNCNAME("ElInfo::calcDet()");
......@@ -92,34 +92,38 @@ namespace AMDiS {
if (dim == 0)
return 1.0;
// dim = dim of world
WorldVector<double> e1, e2, e3, v0;
v0 = coords[0];
for (int i = 0; i < dow; i++) {
e1[i] = coords[1][i] - v0[i];
if (dim > 1)
e2[i] = coords[2][i] - v0[i];
if (dim > 2)
e3[i] = coords[3][i] - v0[i];
}
switch (dow) {
case 1:
det = coords[1][0] - coords[0][0];
break;
case 2:
if (dim == 1) {
WorldVector<double> e1;
e1[0] = coords[1][0] - coords[0][0];
e1[1] = coords[1][1] - coords[0][1];
det = norm(&e1);
} else {
det = e1[0] * e2[1] - e1[1] * e2[0];
det = (coords[1][0] - coords[0][0]) * (coords[2][1] - coords[0][1]) -
(coords[1][1] - coords[0][1]) * (coords[2][0] - coords[0][0]);
}
break;
case 3:
{
WorldVector<double> n;
WorldVector<double> e1, e2, e3, n;
for (int i = 0; i < dow; i++) {
e1[i] = coords[1][i] - coords[0][i];
if (dim > 1)
e2[i] = coords[2][i] - coords[0][i];
if (dim > 2)
e3[i] = coords[3][i] - coords[0][i];
}
if (dim > 1) {
n[0] = e1[1] * e2[2] - e1[2] * e2[1];
n[1] = e1[2] * e2[0] - e1[0] * e2[2];
......
AMDiS/src/ElInfo.h
View file @ e53e8cf3
......@@ -383,7 +383,7 @@ namespace AMDiS {
* Used by non static method \ref calcDet(). Calculates the determinant
* for a given vector of vertex coordinates.
*/
static double calcDet(const FixVec<WorldVector<double>, VERTEX> &coords);
double calcDet(const FixVec<WorldVector<double>, VERTEX> &coords) const;
/** \brief
* Checks whether flag is set in ElInfo's \ref fillFlag_. If not, the program
......
AMDiS/src/Operator.cc
View file @ e53e8cf3
......@@ -124,7 +124,7 @@ namespace AMDiS {
{
FUNCNAME("OperatorTerm::getGradientsAtQPs()");
TEST_EXIT(elInfo->getMesh() == vec->getFESpace()->getMesh())
TEST_EXIT_DBG(elInfo->getMesh() == vec->getFESpace()->getMesh())
("There is something wrong!\n");
return subAssembler->getGradientsAtQPs(vec, elInfo, quad);
......
AMDiS/src/ProblemVec.cc
View file @ e53e8cf3
......@@ -1193,9 +1193,6 @@ namespace AMDiS {
&mainElInfo, &auxElInfo,
&smallElInfo, &largeElInfo);
while (cont) {
Element *mainElem = mainElInfo->getElement();
Element *auxElem = auxElInfo->getElement();
if (useGetBound) {
basisFcts->getBound(mainElInfo, bound);
}
......
AMDiS/src/SubAssembler.h
View file @ e53e8cf3
......@@ -54,9 +54,7 @@ namespace AMDiS {
bool optimized,
FirstOrderType type = GRD_PHI);
/** \brief
* Destructor
*/
/// Destructor
virtual ~SubAssembler() {}
/** \brief
......@@ -79,23 +77,17 @@ namespace AMDiS {
virtual void calculateElementVector(const ElInfo *elInfo,
ElementVector *vec) = 0;
/** \brief
* Returns \ref terms
*/
/// Returns \ref terms
inline std::vector<OperatorTerm*> *getTerms() {
return &terms[omp_get_thread_num()];
}
/** \brief
* Returns \ref quadrature.
*/
/// Returns \ref quadrature.
inline Quadrature *getQuadrature() {
return quadrature;
}
/** \brief
* Sets \ref quadrature to q.
*/
/// Sets \ref quadrature to q.
inline void setQuadrature(Quadrature* q) {
quadrature = q;
}
......@@ -116,9 +108,7 @@ namespace AMDiS {
const ElInfo* elInfo,
Quadrature *quad = NULL);
/** \brief
*
*/
///
double* getVectorAtQPs(DOFVectorBase<double>* dv,
const ElInfo* smallElInfo,
const ElInfo* largeElInfo,
......@@ -132,9 +122,7 @@ namespace AMDiS {
const ElInfo* elInfo,
Quadrature *quad = NULL);
/** \brief
*
*/
///
WorldVector<double>* getGradientsAtQPs(DOFVectorBase<double>* dv,
const ElInfo* smallElInfo,
const ElInfo* largeElInfo,
......@@ -149,32 +137,24 @@ namespace AMDiS {
const ElInfo *largeElInfo = NULL,
Quadrature *quad = NULL);
/** \brief
* Returns \ref psiFast.
*/
/// Returns \ref psiFast.
const FastQuadrature *getPsiFast() const {
return psiFast;
}
/** \brief
* Returns \ref phiFast.
*/
// Returns \ref phiFast.
const FastQuadrature *getPhiFast() const {
return phiFast;
}
protected:
/** \brief
* Updates \ref psiFast and \ref phiFast.
*/
/// Updates \ref psiFast and \ref phiFast.
FastQuadrature *updateFastQuadrature(FastQuadrature *quadFast,
const BasisFunction *psi,
Flag updateFlag);
protected:
/** \brief
* Problem dimension
*/
/// Problem dimension
int dim;
/** \brief
......@@ -189,32 +169,24 @@ namespace AMDiS {
*/
int nCol;
/** \brief
* Used for \ref getVectorAtQPs().
*/
/// Used for \ref getVectorAtQPs().
class ValuesAtQPs {
public:
Vector<double> values;
bool valid;
};
/** \brief
* Used for \ref getGradientsAtQPs().
*/
/// Used for \ref getGradientsAtQPs().
class GradientsAtQPs {
public:
Vector<WorldVector<double> > values;
bool valid;
};
/** \brief
* Used for \ref getVectorAtQPs().
*/
/// Used for \ref getVectorAtQPs().
std::map<const DOFVectorBase<double>*, ValuesAtQPs*> valuesAtQPs;
/** \brief
* Used for \ref getGradientsAtQPs().
*/
/// Used for \ref getGradientsAtQPs().
std::map<const DOFVectorBase<double>*, GradientsAtQPs*> gradientsAtQPs;
/** \brief
......@@ -223,54 +195,34 @@ namespace AMDiS {
*/
WorldVector<double> *coordsAtQPs;
/** \brief
* Used for \ref getCoordsAtQPs().
*/
/// Used for \ref getCoordsAtQPs().
bool coordsValid;
/** \brief
* Used for \ref getCoordsAtQP(). Stores the number of allocated WorldVectors.
*/
/// Used for \ref getCoordsAtQP(). Stores the number of allocated WorldVectors.
int coordsNumAllocated;
/** \brief
* Needed Quadrature. Constructed in the constructor of SubAssembler
*/
/// Needed Quadrature. Constructed in the constructor of SubAssembler
Quadrature *quadrature;
/** \brief
* FastQuadrature for row basis functions
*/
/// FastQuadrature for row basis functions
FastQuadrature *psiFast;
/** \brief
* FastQuadrature for column basis functions
*/
/// FastQuadrature for column basis functions
FastQuadrature *phiFast;
/** \brief
* Corresponding Assembler
*/
/// Corresponding Assembler
Assembler* owner;
/** \brief
* Flag that specifies whether the element matrix is symmetric.
*/
/// Flag that specifies whether the element matrix is symmetric.
bool symmetric;
/** \brief
* List of all terms with a contribution to this SubAssembler
*/
/// List of all terms with a contribution to this SubAssembler
std::vector< std::vector<OperatorTerm*> > terms;
/** \brief
*
*/
///
bool opt;
/** \brief
*
*/
///
bool firstCall;
friend class Assembler;
......
AMDiS/src/SubElInfo.cc
View file @ e53e8cf3
......@@ -32,7 +32,7 @@ namespace AMDiS {
worldCoords[i]);
}
det = ElInfo::calcDet(worldCoords);
det = elInfo->calcDet(worldCoords);
}
}
AMDiS/src/SubElInfo.h
View file @ e53e8cf3
......@@ -30,23 +30,15 @@ namespace AMDiS {
public:
MEMORY_MANAGED(SubElInfo);
/**
* Constructor
*/
SubElInfo(VectorOfFixVecs<DimVec<double> > *lambda,
const ElInfo *elInfo);
/// Constructor
SubElInfo(VectorOfFixVecs<DimVec<double> > *lambda, const ElInfo *elInfo);
/**
* Destructor.
*/
/// Destructor
~SubElInfo() {
DELETE lambda;
}
/**
* Get b-th coordinate of the a-th vertex of subelement (barycentric
* coordinates).
*/
/// Get b-th coordinate of the a-th vertex of subelement (barycentriccoordinates).
inline const double getLambda(int a, int b) const {
if (lambda)
return (*lambda)[a][b];
......@@ -54,42 +46,30 @@ namespace AMDiS {
return 0.0;
}
/**
* Get coordinates of a-th vertex of subelement (barycentric coordinates).
*/
/// Get coordinates of a-th vertex of subelement (barycentric coordinates).
inline const DimVec<double>& getLambda(int a) const {
return (*lambda)[a];
}
/**
* Get coordinates of all vertices of subelement (barycentric coordinates).
*/
/// Get coordinates of all vertices of subelement (barycentric coordinates).
inline VectorOfFixVecs<DimVec<double> > *getLambda() const {
return lambda;
}
/**
* Get determinant corresponding to subelement.
*/
/// Get determinant corresponding to subelement.
inline const double getDet() const {
return det;
}
protected:
/**