diff --git a/AMDiS/src/AMDiS_fwd.h b/AMDiS/src/AMDiS_fwd.h
index dbbadad66b003488623e264fae976500740bd909..b2d5a33693a7c7fe2b4f3da91eacabab2cad474c 100644
--- a/AMDiS/src/AMDiS_fwd.h
+++ b/AMDiS/src/AMDiS_fwd.h
@@ -54,6 +54,7 @@ namespace AMDiS {
class Operator;
class ProblemInstat;
class ProblemIterationInterface;
+ class ProblemTimeInterface;
class ProblemVec;
class Projection;
class PreconditionerScal;
@@ -65,17 +66,20 @@ namespace AMDiS {
class VertexVector;
template<typename ReturnType, typename ArgumentType> class AbstractFunction;
+ template<typename T> class DOFIndexed;
template<typename T> class DOFVectorBase;
template<typename T> class DOFVector;
+ template<typename T> class DimVec;
+ template<typename T> class DimMat;
+ template<typename ITLSolver> class ITL_OEMSolver;
template<typename T> class ITL_Preconditioner;
template<typename T> class Matrix;
template<typename T> class MatVecMultiplier;
template<typename T> class SolverMatrix;
template<typename T> class Vector;
template<typename T> class WorldVector;
-
- template <typename ITLSolver> class ITL_OEMSolver;
-
+ template<typename T> class WorldMatrix;
+ template<typename T> class VectorOfFixVecs;
} // namespace AMDiS
#endif // AMDIS_AMDIS_FWD_INCLUDE
diff --git a/AMDiS/src/AbstractFunction.h b/AMDiS/src/AbstractFunction.h
index 5816de5b06b598f086772aac1fd1a4847716a26d..e1dcf88fd95031c7a4c78d6de4980777c2c626a5 100644
--- a/AMDiS/src/AbstractFunction.h
+++ b/AMDiS/src/AbstractFunction.h
@@ -27,10 +27,6 @@
namespace AMDiS {
- // ============================================================================
- // ===== class AbstractFunction ===============================================
- // ============================================================================
-
/**
* \ingroup Common
*
@@ -46,35 +42,25 @@ namespace AMDiS {
class AbstractFunction
{
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
AbstractFunction(int degree = 0) :
degree_(degree)
{}
virtual ~AbstractFunction() {}
- /** \brief
- * Returns \ref degree_.
- */
+ /// Returns \ref degree_.
inline int getDegree() const {
return degree_;
}
- /** \brief
- * Deligates the evaluation to overriden method f.
- */
+ /// Deligates the evaluation to overriden method f.
virtual ReturnType operator()(const ArgumentType& x) const = 0;
protected:
int degree_;
};
- // ============================================================================
- // ===== class BinaryAbstractFunction =========================================
- // ============================================================================
-
/**
* \ingroup Common
*
@@ -87,25 +73,19 @@ namespace AMDiS {
class BinaryAbstractFunction
{
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
BinaryAbstractFunction(int degree = 0) :
degree_(degree)
- {};
+ {}
- virtual ~BinaryAbstractFunction() {};
+ virtual ~BinaryAbstractFunction() {}
- /** \brief
- * Returns \ref degree_.
- */
+ /// Returns \ref degree_.
inline int getDegree() const {
return degree_;
- };
+ }
- /** \brief
- * Deligates the evaluation to overriden method f.
- */
+ /// Deligates the evaluation to overriden method f.
virtual ReturnType operator()(const ArgumentType1& x,
const ArgumentType2& y) const = 0;
@@ -113,10 +93,6 @@ namespace AMDiS {
int degree_;
};
- // ============================================================================
- // ===== class TertiaryAbstractFunction =======================================
- // ============================================================================
-
/**
* \ingroup Common
*
@@ -130,25 +106,19 @@ namespace AMDiS {
class TertiaryAbstractFunction
{
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
TertiaryAbstractFunction(int degree = 0) :
degree_(degree)
- {};
+ {}
- virtual ~TertiaryAbstractFunction() {};
+ virtual ~TertiaryAbstractFunction() {}
- /** \brief
- * Returns \ref degree_.
- */
+ /// Returns \ref degree_.
inline int getDegree() const {
return degree_;
- };
+ }
- /** \brief
- * function evaluation.
- */
+ /// function evaluation.
virtual ReturnType operator()(const ArgumentType1& x,
const ArgumentType2& y,
const ArgumentType3& z) const = 0;
@@ -157,11 +127,6 @@ namespace AMDiS {
int degree_;
};
-
- // ============================================================================
- // ===== class QuartAbstractFunction =======================================
- // ============================================================================
-
/**
* \ingroup Common
*
@@ -176,25 +141,19 @@ namespace AMDiS {
class QuartAbstractFunction
{
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
QuartAbstractFunction(int degree = 0) :
degree_(degree)
- {};
+ {}
- virtual ~QuartAbstractFunction() {};
+ virtual ~QuartAbstractFunction() {}
- /** \brief
- * Returns \ref degree_.
- */
+ /// Returns \ref degree_.
inline int getDegree() const {
return degree_;
- };
+ }
- /** \brief
- * function evaluation.
- */
+ /// function evaluation.
virtual ReturnType operator()(const ArgumentType1& x,
const ArgumentType2& y,
const ArgumentType3& z,
diff --git a/AMDiS/src/AdaptBase.h b/AMDiS/src/AdaptBase.h
index c7c862ca43d4905cba022ef746fbb6723e0c62b6..113d2e7530322fab74a82e291d3d5dbff32c4d9d 100644
--- a/AMDiS/src/AdaptBase.h
+++ b/AMDiS/src/AdaptBase.h
@@ -23,17 +23,10 @@
#define AMDIS_ADAPTBASE_H
#include <string>
+#include "AMDiS_fwd.h"
namespace AMDiS {
- class ProblemIterationInterface;
- class ProblemTimeInterface;
- class AdaptInfo;
-
- // ============================================================================
- // ===== class AdaptBase ======================================================
- // ============================================================================
-
/// Interface for adaption loops.
class AdaptBase
{
diff --git a/AMDiS/src/AdaptInfo.h b/AMDiS/src/AdaptInfo.h
index a6238309cbee24a4861d776856700129ed13f270..bd60c9279b27d0963e251723f5758999aa2bdac3 100644
--- a/AMDiS/src/AdaptInfo.h
+++ b/AMDiS/src/AdaptInfo.h
@@ -29,10 +29,6 @@
namespace AMDiS {
- // ===========================================================================
- // ===== class AdaptInfo =====================================================
- // ===========================================================================
-
/**
* \ingroup Adaption
*
@@ -49,9 +45,7 @@ namespace AMDiS {
*/
class ScalContent {
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
ScalContent(std::string prefix)
: est_sum(0.0),
est_t_sum(0.0),
@@ -86,49 +80,31 @@ namespace AMDiS {
timeErrLow = totalTol * relTimeErr * timeTheta2;
}
- /** \brief
- * Sum of all error estimates
- */
+ /// Sum of all error estimates
double est_sum;
- /** \brief
- * Sum of all time error estimates
- */
+ /// Sum of all time error estimates
double est_t_sum;
- /** \brief
- * maximal local error estimate.
- */
+ /// maximal local error estimate.
double est_max;
- /** \brief
- * Maximum of all time error estimates
- */
+ /// Maximum of all time error estimates
double est_t_max;
- /** \brief
- * Tolerance for the (absolute or relative) error
- */
+ /// Tolerance for the (absolute or relative) error
double spaceTolerance;
- /** \brief
- * Time tolerance.
- */
+ /// Time tolerance.
double timeTolerance;
- /** \brief
- * Lower bound for the time error.
- */
+ /// Lower bound for the time error.
double timeErrLow;
- /** \brief
- * true if coarsening is allowed, false otherwise.
- */
+ /// true if coarsening is allowed, false otherwise.
int coarsenAllowed;
- /** \brief
- * true if refinement is allowed, false otherwise.
- */
+ /// true if refinement is allowed, false otherwise.
int refinementAllowed;
/** \brief
@@ -136,20 +112,18 @@ namespace AMDiS {
* performed when an element is marked for refinement; usually the value is
* 1 or DIM
*/
- int refineBisections;
+ int refineBisections;
/** \brief
* parameter to tell the marking strategy how many bisections should
* be undone when an element is marked for coarsening; usually the value is
* 1 or DIM
*/
- int coarseBisections;
+ int coarseBisections;
};
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
AdaptInfo(std::string name_, int size = 1)
: name(name_),
spaceIteration(-1),
@@ -197,13 +171,10 @@ namespace AMDiS {
}
}
- /** \brief
- * Destructor.
- */
+ /// Destructor.
virtual ~AdaptInfo() {
- for (int i = 0; i < scalContents.getSize(); i++) {
+ for (int i = 0; i < scalContents.getSize(); i++)
delete scalContents[i];
- }
}
inline void reset()
@@ -220,185 +191,135 @@ namespace AMDiS {
GET_PARAMETER(0, name + "->timestep", "%f", ×tep);
}
- /** \brief
- * Returns whether space tolerance is reached.
- */
+ /// Returns whether space tolerance is reached.
virtual bool spaceToleranceReached() {
int size = scalContents.getSize();
- for (int i = 0; i < size; i++) {
- if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance)) {
+ for (int i = 0; i < size; i++)
+ if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance))
return false;
- }
- }
+
return true;
}
- /** \brief
- * Returns whether space tolerance of component i is reached.
- */
+ /// Returns whether space tolerance of component i is reached.
virtual bool spaceToleranceReached(int i) {
- if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance)) {
+ if (!(scalContents[i]->est_sum < scalContents[i]->spaceTolerance))
return false;
- } else {
+ else
return true;
- }
}
- /** \brief
- * Returns whether time tolerance is reached.
- */
+ /// Returns whether time tolerance is reached.
virtual bool timeToleranceReached() {
int size = scalContents.getSize();
- for (int i = 0; i < size; i++) {
- if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance)) {
+ for (int i = 0; i < size; i++)
+ if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance))
return false;
- }
- }
+
return true;
}
- /** \brief
- * Returns whether time tolerance of component i is reached.
- */
+ /// Returns whether time tolerance of component i is reached.
virtual bool timeToleranceReached(int i) {
- if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance)) {
+ if (!(scalContents[i]->est_t_sum < scalContents[i]->timeTolerance))
return false;
- } else {
+ else
return true;
- }
}
- /** \brief
- * Returns whether time error is under its lower bound.
- */
+ /// Returns whether time error is under its lower bound.
virtual bool timeErrorLow() {
int size = scalContents.getSize();
- for (int i = 0; i < size; i++) {
- if (!(scalContents[i]->est_t_sum < scalContents[i]->timeErrLow)) {
+ for (int i = 0; i < size; i++)
+ if (!(scalContents[i]->est_t_sum < scalContents[i]->timeErrLow))
return false;
- }
- }
+
return true;
}
- /** \brief
- * Print debug information about time error and its bound.
- */
+ /// Print debug information about time error and its bound.
void printTimeErrorLowInfo() {
- for (int i = 0; i < scalContents.getSize(); i++) {
+ for (int i = 0; i < scalContents.getSize(); i++)
std::cout << " Time error estimate = " << scalContents[i]->est_t_sum
<< " Time error bound = " << scalContents[i]->timeErrLow << "\n";
- }
}
- /** \brief
- * Returns \ref spaceIteration.
- */
+ /// Returns \ref spaceIteration.
inline int getSpaceIteration() {
return spaceIteration;
}
- /** \brief
- * Sets \ref spaceIteration.
- */
+ /// Sets \ref spaceIteration.
inline void setSpaceIteration(int it) {
spaceIteration = it;
}
- /** \brief
- * Returns \ref maxSpaceIteration.
- */
+ /// Returns \ref maxSpaceIteration.
inline int getMaxSpaceIteration() {
return maxSpaceIteration;
}
- /** \brief
- * Sets \ref maxSpaceIteration.
- */
+ /// Sets \ref maxSpaceIteration.
inline void setMaxSpaceIteration(int it) {
maxSpaceIteration = it;
}
- /** \brief
- * Increments \ref spaceIteration by 1;
- */
+ /// Increments \ref spaceIteration by 1;
inline void incSpaceIteration() {
spaceIteration++;
}
- /** \brief
- * Sets \ref timestepIteration.
- */
+ /// Sets \ref timestepIteration.
inline void setTimestepIteration(int it) {
timestepIteration = it;
}
- /** \brief
- * Returns \ref timestepIteration.
- */
+ /// Returns \ref timestepIteration.
inline int getTimestepIteration() {
return timestepIteration;
}
- /** \brief
- * Increments \ref timestepIteration by 1;
- */
+ /// Increments \ref timestepIteration by 1;
inline void incTimestepIteration() {
timestepIteration++;
}
- /** \brief
- * Returns \ref maxTimestepIteration.
- */
+ /// Returns \ref maxTimestepIteration.
inline int getMaxTimestepIteration() {
return maxTimestepIteration;
}
- /** \brief
- * Sets \ref maxTimestepIteration.
- */
+ /// Sets \ref maxTimestepIteration.
inline void setMaxTimestepIteration(int it) {
maxTimestepIteration = it;
}
- /** \brief
- * Sets \ref timeIteration.
- */
+ /// Sets \ref timeIteration.
inline void setTimeIteration(int it) {
timeIteration = it;
}
- /** \brief
- * Returns \ref timeIteration.
- */
+ /// Returns \ref timeIteration.
inline int getTimeIteration() {
return timeIteration;
}
- /** \brief
- * Increments \ref timesIteration by 1;
- */
+ /// Increments \ref timesIteration by 1;
inline void incTimeIteration() {
timeIteration++;
}
- /** \brief
- * Returns \ref maxTimeIteration.
- */
+ /// Returns \ref maxTimeIteration.
inline int getMaxTimeIteration() {
return maxTimeIteration;
}
- /** \brief
- * Sets \ref maxTimeIteration.
- */
+ /// Sets \ref maxTimeIteration.
inline void setMaxTimeIteration(int it) {
maxTimeIteration = it;
}
- /** \brief
- * Returns \ref timestepNumber.
- */
+ /// Returns \ref timestepNumber.
inline int getTimestepNumber() {
return timestepNumber;
}
@@ -408,100 +329,72 @@ namespace AMDiS {
return nTimesteps;
}
- /** \brief
- * Increments \ref timestepNumber by 1;
- */
+ /// Increments \ref timestepNumber by 1;
inline void incTimestepNumber() {
timestepNumber++;
}
- /** \brief
- * Sets \ref est_sum.
- */
+ /// Sets \ref est_sum.
inline void setEstSum(double e, int index) {
scalContents[index]->est_sum = e;
}
- /** \brief
- * Sets \ref est_max.
- */
+ /// Sets \ref est_max.
inline void setEstMax(double e, int index) {
scalContents[index]->est_max = e;
}
- /** \brief
- * Sets \ref est_max.
- */
+ /// Sets \ref est_max.
inline void setTimeEstMax(double e, int index) {
scalContents[index]->est_t_max = e;
}
- /** \brief
- * Sets \ref est_t_sum.
- */
+ /// Sets \ref est_t_sum.
inline void setTimeEstSum(double e, int index) {
scalContents[index]->est_t_sum = e;
}
- /** \brief
- * Returns \ref est_sum.
- */
+ /// Returns \ref est_sum.
inline double getEstSum(int index) {
return scalContents[index]->est_sum;
}
- /** \brief
- * Returns \ref est_t_sum.
- */
+ /// Returns \ref est_t_sum.
inline double getEstTSum(int index) {
return scalContents[index]->est_t_sum;
}
- /** \brief
- * Returns \ref est_max.
- */
+ /// Returns \ref est_max.
inline double getEstMax(int index) {
return scalContents[index]->est_max;
}
- /** \brief
- * Returns \ref est_max.
- */
+ /// Returns \ref est_max.
inline double getTimeEstMax(int index) {
return scalContents[index]->est_t_max;
}
- /** \brief
- * Returns \ref est_t_sum.
- */
+ /// Returns \ref est_t_sum.
inline double getTimeEstSum(int index) {
return scalContents[index]->est_t_sum;
}
- /** \brief
- * Returns \ref spaceTolerance.
- */
+ /// Returns \ref spaceTolerance.
inline double getSpaceTolerance(int index) {
return scalContents[index]->spaceTolerance;
}
- /** \brief
- * Sets \ref spaceTolerance.
- */
+ /// Sets \ref spaceTolerance.
inline void setSpaceTolerance(int index, double tol) {
scalContents[index]->spaceTolerance = tol;
}
- /** \brief
- * Returns \ref timeTolerance.
- */
+ /// Returns \ref timeTolerance.
inline double getTimeTolerance(int index) {
return scalContents[index]->timeTolerance;
}
- /** \brief
- * Sets \ref time
- */
+ /// Sets \ref time
inline double setTime(double t) {
time = t;
if (time > endTime)
@@ -512,23 +405,17 @@ namespace AMDiS {
return time;
}
- /** \brief
- * Gets \ref time
- */
+ /// Gets \ref time
inline double getTime() {
return time;
}
- /** \brief
- * Gets \ref &time
- */
+ /// Gets \ref &time
inline double *getTimePtr() {
return &time;
}
- /** \brief
- * Sets \ref timestep
- */
+ /// Sets \ref timestep
inline double setTimestep(double t) {
timestep = t;
if (timestep > maxTimestep) {
@@ -555,121 +442,87 @@ namespace AMDiS {
return !(time < endTime - DBL_TOL);
}
- /** \brief
- * Gets \ref timestep
- */
+ /// Gets \ref timestep
inline double getTimestep() {
return timestep;
}
- /** \brief
- * Sets \ref minTimestep
- */
+ /// Sets \ref minTimestep
inline void setMinTimestep(double t) {
minTimestep = t;
}
- /** \brief
- * Gets \ref minTimestep
- */
+ /// Gets \ref minTimestep
inline double getMinTimestep() {
return minTimestep;
}
- /** \brief
- * Sets \ref maxTimestep
- */
+ /// Sets \ref maxTimestep
inline void setMaxTimestep(double t) {
maxTimestep = t;
}
- /** \brief
- * Gets \ref maxTimestep
- */
+ /// Gets \ref maxTimestep
inline double getMaxTimestep() {
return maxTimestep;
}
- /** \brief
- * Gets \ref ×tep
- */
+ /// Gets \ref ×tep
inline double *getTimestepPtr() {
return ×tep;
}
- /** \brief
- * Sets \ref startTime = time
- */
+ /// Sets \ref startTime = time
inline void setStartTime(double time) {
startTime = time;
}
- /** \brief
- * Sets \ref endTime = time
- */
+ /// Sets \ref endTime = time
inline void setEndTime(double time) {
endTime = time;
}
- /** \brief
- * Returns \ref startTime
- */
+ /// Returns \ref startTime
inline double getStartTime() {
return startTime;
}
- /** \brief
- * Returns \ref endTime
- */
+ /// Returns \ref endTime
inline double getEndTime() {
return endTime;
}
- /** \brief
- * Returns \ref timeErrLow.
- */
+ /// Returns \ref timeErrLow.
inline double getTimeErrLow(int index) {
return scalContents[index]->timeErrLow;
}
- /** \brief
- * Returns whether coarsening is allowed or not.
- */
+ /// Returns whether coarsening is allowed or not.
inline bool isCoarseningAllowed(int index) {
return (scalContents[index]->coarsenAllowed == 1);
}
- /** \brief
- * Returns whether coarsening is allowed or not.
- */
+ /// Returns whether coarsening is allowed or not.
inline bool isRefinementAllowed(int index) {
return (scalContents[index]->refinementAllowed == 1);
}
- /** \brief
- *
- */
+ ///
inline void allowRefinement(bool allow, int index) {
scalContents[index]->refinementAllowed = allow;
}
- /** \brief
- *
- */
+ ///
inline void allowCoarsening(bool allow, int index) {
scalContents[index]->coarsenAllowed = allow;
}
- /** \brief
- * Returns \ref refineBisections
- */
+ /// Returns \ref refineBisections
inline const int getRefineBisections(int index) const {
return scalContents[index]->refineBisections;
}
- /** \brief
- * Returns \ref coarseBisections
- */
+ /// Returns \ref coarseBisections
inline const int getCoarseBisections(int index) const {
return scalContents[index]->coarseBisections;
}
@@ -710,9 +563,7 @@ namespace AMDiS {
return solverResidual;
}
- /** \brief
- * Returns true, if the adaptive procedure was deserialized from a file.
- */
+ /// Returns true, if the adaptive procedure was deserialized from a file.
const bool isDeserialized() const {
return isDeserialized_;
}
@@ -721,9 +572,7 @@ namespace AMDiS {
isDeserialized_ = b;
}
- /** \brief
- * Creates new scalContents with the given size.
- */
+ /// Creates new scalContents with the given size.
void setScalContents(int newSize);
/** \brief
@@ -741,20 +590,15 @@ namespace AMDiS {
timestepNumber = 0;
}
- // ===== Serialiazable implementation =====
void serialize(std::ostream& out);
void deserialize(std::istream& in);
protected:
- /** \brief
- * Name.
- */
+ /// Name.
std::string name;
- /** \brief
- * Current space iteration
- */
+ /// Current space iteration
int spaceIteration;
/** \brief
@@ -763,59 +607,37 @@ namespace AMDiS {
*/
int maxSpaceIteration;
- /** \brief
- * Current timestep iteration
- */
+ /// Current timestep iteration
int timestepIteration;
- /** \brief
- * Maximal number of iterations for choosing a timestep
- */
+ /// Maximal number of iterations for choosing a timestep
int maxTimestepIteration;
- /** \brief
- * Current time iteration
- */
+ /// Current time iteration
int timeIteration;
- /** \brief
- * Maximal number of time iterations
- */
+ /// Maximal number of time iterations
int maxTimeIteration;
- /** \brief
- * Actual time, end of time interval for current time step
- */
+ /// Actual time, end of time interval for current time step
double time;
- /** \brief
- * Initial time
- */
+ /// Initial time
double startTime;
- /** \brief
- * Final time
- */
+ /// Final time
double endTime;
- /** \brief
- * Current time step size
- */
+ /// Current time step size
double timestep;
- /** \brief
- * Minimal step size
- */
+ /// Minimal step size
double minTimestep;
- /** \brief
- * Maximal step size
- */
+ /// Maximal step size
double maxTimestep;
- /** \brief
- * Number of current time step
- */
+ /// Number of current time step
int timestepNumber;
/** \brief
@@ -825,35 +647,22 @@ namespace AMDiS {
*/
int nTimesteps;
- /** \brief
- * number of iterations needed of linear or nonlinear solver
- */
+ /// number of iterations needed of linear or nonlinear solver
int solverIterations;
-
- /** \brief
- * maximal number of iterations needed of linear or nonlinear solver
- */
+ /// maximal number of iterations needed of linear or nonlinear solver
int maxSolverIterations;
- /** \brief
- *
- */
+ ///
double solverTolerance;
- /** \brief
- *
- */
+ ///
double solverResidual;
- /** \brief
- * Scalar adapt infos.
- */
+ /// Scalar adapt infos.
Vector<ScalContent*> scalContents;
- /** \brief
- * Is true, if the adaptive procedure was deserialized from a file.
- */
+ /// Is true, if the adaptive procedure was deserialized from a file.
bool isDeserialized_;
};
diff --git a/AMDiS/src/AdaptInstationary.h b/AMDiS/src/AdaptInstationary.h
index 7f5803d9d5d79d2bdafb797de3ef2812cc26c96d..2d146d24bce2c0311dbfb6a4864f068ae63daa2d 100644
--- a/AMDiS/src/AdaptInstationary.h
+++ b/AMDiS/src/AdaptInstationary.h
@@ -29,16 +29,10 @@
#include "MemoryManager.h"
#include "AdaptInfo.h"
#include "AdaptBase.h"
+#include "AMDiS_fwd.h"
namespace AMDiS {
- class ProblemIterationInterface;
- class ProblemTimeInterface;
-
- // ============================================================================
- // ===== class AdaptInstationary ==============================================
- // ============================================================================
-
/** \ingroup Adaption
* \brief
* AdaptInstationary implements the adaptive procdure for time dependent
@@ -48,8 +42,6 @@ namespace AMDiS {
class AdaptInstationary : public AdaptBase
{
public:
- MEMORY_MANAGED(AdaptInstationary);
-
/** \brief
* Creates a AdaptInstationary object with the given name for the time
* dependent problem problemInstat.
diff --git a/AMDiS/src/AdaptStationary.h b/AMDiS/src/AdaptStationary.h
index c96204ee1b5f28a27425720a9513d5d533a7f33a..f5faa8c9d293e80c85f9e8019aa69425fcbab328 100644
--- a/AMDiS/src/AdaptStationary.h
+++ b/AMDiS/src/AdaptStationary.h
@@ -38,10 +38,6 @@
namespace AMDiS {
- // ============================================================================
- // ===== class AdaptStationary ================================================
- // ============================================================================
-
/** \ingroup Adaption
* \brief
* AdaptStationary contains information about the adaptive procedure and the
@@ -50,8 +46,6 @@ namespace AMDiS {
class AdaptStationary : public AdaptBase
{
public:
- MEMORY_MANAGED(AdaptStationary);
-
/// Creates a AdaptStationary object with given name.
AdaptStationary(const std::string &name,
ProblemIterationInterface *prob,
diff --git a/AMDiS/src/Assembler.cc b/AMDiS/src/Assembler.cc
index c6578a9bf7e5c7e41c6627c4806c7d142291da57..0df8b234bd2c4bb4f02d923a00a229740ab25a53 100644
--- a/AMDiS/src/Assembler.cc
+++ b/AMDiS/src/Assembler.cc
@@ -1,3 +1,7 @@
+#include <vector>
+#include <algorithm>
+#include <boost/numeric/mtl/mtl.hpp>
+
#include "Assembler.h"
#include "Operator.h"
#include "Element.h"
@@ -6,8 +10,6 @@
#include "ElementVector.h"
#include "DOFVector.h"
#include "OpenMP.h"
-#include <vector>
-#include <algorithm>
namespace AMDiS {
@@ -293,12 +295,12 @@ namespace AMDiS {
operat->uhOld->getLocalVector(auxEl, uhOldLoc);
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
+ const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
for (int i = 0; i < nBasFcts; i++) {
uhOldLoc2[i] = 0.0;
for (int j = 0; j < nBasFcts; j++) {
- uhOldLoc2[i] += (*m)[j][i] * uhOldLoc[i];
+ uhOldLoc2[i] += m[j][i] * uhOldLoc[i];
}
}
diff --git a/AMDiS/src/Assembler.h b/AMDiS/src/Assembler.h
index 649f38fa86c7d2b6d94286acfec52d689312be7a..1602ec5d1e6bbb475bfa81d505fd490e5867adc3 100644
--- a/AMDiS/src/Assembler.h
+++ b/AMDiS/src/Assembler.h
@@ -38,15 +38,10 @@
#include "SecondOrderAssembler.h"
#include "ElInfo.h"
#include "OpenMP.h"
+#include "AMDiS_fwd.h"
namespace AMDiS {
- class Element;
- class Quadrature;
- class ElementMatrix;
- class ElementVector;
- class Operator;
-
/**
* \ingroup Assembler
*
@@ -58,8 +53,6 @@ namespace AMDiS {
class Assembler
{
public:
- MEMORY_MANAGED(Assembler);
-
/// Constructor
Assembler(Operator *op,
const FiniteElemSpace *rowFESpace,
@@ -276,10 +269,6 @@ namespace AMDiS {
friend class SecondOrderAssembler;
};
- // ============================================================================
- // ===== class StandardAssembler =============================================
- // ============================================================================
-
/**
* \ingroup Assembler
*
@@ -289,8 +278,6 @@ namespace AMDiS {
class StandardAssembler : public Assembler
{
public:
- MEMORY_MANAGED(StandardAssembler);
-
/// Constructor
StandardAssembler(Operator *op,
Quadrature *quad2,
@@ -301,10 +288,6 @@ namespace AMDiS {
const FiniteElemSpace *colFESpace = NULL);
};
- // ============================================================================
- // ===== class OptimizedAssembler =============================================
- // ============================================================================
-
/**
* \ingroup Assembler
*
@@ -314,8 +297,6 @@ namespace AMDiS {
class OptimizedAssembler : public Assembler
{
public:
- MEMORY_MANAGED(OptimizedAssembler);
-
/// Constructor
OptimizedAssembler(Operator *op,
Quadrature *quad2,
diff --git a/AMDiS/src/BasisFunction.h b/AMDiS/src/BasisFunction.h
index 24206efe9d4ee658ef574473bc81c66430b8e563..b654d14bfb053b3621c6456de877beea1cc86902 100644
--- a/AMDiS/src/BasisFunction.h
+++ b/AMDiS/src/BasisFunction.h
@@ -30,26 +30,7 @@
namespace AMDiS {
- class DOFAdmin;
- class Element;
- class ElInfo;
- class RCNeighbourList;
- template<typename T> class WorldVector;
- template<typename T> class WorldMatrix;
- class Quadrature;
-
- template <typename ReturnType, typename ArgumentType> class AbstractFunction;
- // template <typename T> class DOFVector;
- template <typename T> class DOFIndexed;
- template <typename T> class DimVec;
- template <typename T> class DimMat;
- template <typename T, GeoIndex d> class FixVec;
- template <typename T> class VectorOfFixVecs;
-
-
- /** \brief
- * Function interface for evaluating basis functions.
- */
+ /// Function interface for evaluating basis functions.
class BasFctType
{
public:
@@ -60,10 +41,7 @@ namespace AMDiS {
virtual double operator()(const DimVec<double>&) const = 0;
};
-
- /** \brief
- * Function interface for evaluating gradients of basis functions.
- */
+ /// Function interface for evaluating gradients of basis functions.
class GrdBasFctType
{
public:
@@ -74,11 +52,8 @@ namespace AMDiS {
virtual void operator()(const DimVec<double>&,
DimVec<double>&) const = 0;
};
-
- /** \brief
- * Function interface for evaluating second derivative of basis functions.
- */
+ /// Function interface for evaluating second derivative of basis functions.
class D2BasFctType
{
public:
@@ -94,7 +69,6 @@ namespace AMDiS {
typedef GrdBasFctType *GBFptr;
typedef D2BasFctType *DBFptr;
-
/** \ingroup FEMSpace
* \brief
* Base class for finite element basis functions. In order to build up a
@@ -108,34 +82,24 @@ namespace AMDiS {
class BasisFunction
{
protected:
- /** \brief
- * Creates a BasisFunction object of given dim and degree
- */
+ /// Creates a BasisFunction object of given dim and degree
BasisFunction(const std::string& name, int dim, int degree);
- /** \brief
- * destructor
- */
+ /// destructor
virtual ~BasisFunction();
public:
- /** \brief
- * compares two BasisFunction objects.
- */
+ /// compares two BasisFunction objects.
virtual bool operator==(const BasisFunction& a) const {
return a.getName() == name;
}
- /** \brief
- * Returns !(*this == b)
- */
+ /// Returns !(*this == b)
inline bool operator!=(const BasisFunction& b) const {
return !(operator == (b));
}
- /** \brief
- * Used by \ref getDOFIndices and \ref getVec
- */
+ /// Used by \ref getDOFIndices and \ref getVec
virtual int* orderOfPositionIndices(const Element* el, GeoIndex position,
int positionIndex) const = 0;
@@ -225,33 +189,24 @@ namespace AMDiS {
AbstractFunction<double, WorldVector<double> > *f,
double *coeff) = 0;
-
- /** \brief
- * WorldVector<double> valued interpol function.
- */
+ /// WorldVector<double> valued interpol function.
virtual const WorldVector<double>*
interpol(const ElInfo *el_info, int no,
const int *b_no,
AbstractFunction<WorldVector<double>,WorldVector<double> > *f,
WorldVector<double> *vec) = 0;
- /** \brief
- * Returns the i-th local basis function
- */
+ /// Returns the i-th local basis function
inline BasFctType *getPhi(int i) const {
return (*phi)[i];
}
- /** \brief
- * Returns the gradient of the i-th local basis function
- */
+ /// Returns the gradient of the i-th local basis function
inline GrdBasFctType *getGrdPhi(int i) const {
- return (*grdPhi)[i];
+ return (*grdPhi)[i];
}
- /** \brief
- * Returns the second derivative of the i-th local basis function
- */
+ /// Returns the second derivative of the i-th local basis function
inline D2BasFctType *getD2Phi(int i) const {
return (*d2Phi)[i];
}
@@ -281,49 +236,35 @@ namespace AMDiS {
DOFVector<double>* /*fh*/)
{}
- /** \brief
- * WorldVector<double> valued l2ScpFctBas function
- */
+ /// WorldVector<double> valued l2ScpFctBas function
virtual void l2ScpFctBas(Quadrature* ,
AbstractFunction<WorldVector<double>, WorldVector<double> >* /*f*/,
DOFVector<WorldVector<double> >* /*fh*/)
{}
- /** \brief
- * Interpolates a DOFIndexed<double> after refinement
- */
- virtual void refineInter(DOFIndexed<double> *, RCNeighbourList*, int)
+ /// Interpolates a DOFIndexed<double> after refinement
+ virtual void refineInter(DOFIndexed<double> *, RCNeighbourList*, int)
{}
- /** \brief
- * Interpolates a DOFIndexed<double> after coarsening
- */
- virtual void coarseInter(DOFIndexed<double> *, RCNeighbourList*, int)
+ /// Interpolates a DOFIndexed<double> after coarsening
+ virtual void coarseInter(DOFIndexed<double> *, RCNeighbourList*, int)
{}
- /** \brief
- * Restricts a DOFIndexed<double> after coarsening
- */
- virtual void coarseRestr(DOFIndexed<double> *, RCNeighbourList*, int)
+ /// Restricts a DOFIndexed<double> after coarsening
+ virtual void coarseRestr(DOFIndexed<double> *, RCNeighbourList*, int)
{}
- /** \brief
- * Interpolates a DOFVector<WorldVector<double> > after refinement
- */
- virtual void refineInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
+ /// Interpolates a DOFVector<WorldVector<double> > after refinement
+ virtual void refineInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
{}
- /** \brief
- * Interpolates a DOFVector<WorldVector<double> > after coarsening
- */
- virtual void coarseInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
+ /// Interpolates a DOFVector<WorldVector<double> > after coarsening
+ virtual void coarseInter(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
{}
- /** \brief
- * Restricts a DOFVector<WorldVector<double> > after coarsening
- */
- virtual void coarseRestr(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
+ /// Restricts a DOFVector<WorldVector<double> > after coarsening
+ virtual void coarseRestr(DOFVector<WorldVector<double> >*, RCNeighbourList*, int)
{}
/// Returns local dof indices of the element for the given fe space.
@@ -379,52 +320,33 @@ namespace AMDiS {
const double* uh, WorldMatrix<double>* val) const;
protected:
- /** \brief
- * Textual description
- */
+ /// Textual description
std::string name;
- /** \brief
- * Number of basisfunctions on one Element
- */
+ /// Number of basisfunctions on one Element
int nBasFcts;
- /** \brief
- * Maximal degree of the basis functions
- */
+ /// Maximal degree of the basis functions
int degree;
- /** \brief
- * Dimension of the basis functions
- */
+ /// Dimension of the basis functions
int dim;
- /** \brief
- * Dimension of the world.
- */
+ /// Dimension of the world.
int dow;
- /** \brief
- * Number of DOFs at the different positions
- */
+ /// Number of DOFs at the different positions
DimVec<int> *nDOF;
- /** \brief
- * Vector of the local functions
- */
+ /// Vector of the local functions
std::vector<BasFctType*> *phi;
- /** \brief
- * Vector of gradients
- */
+ /// Vector of gradients
std::vector<GrdBasFctType*> *grdPhi;
- /** \brief
- * Vector of second derivatives
- */
+ /// Vector of second derivatives
std::vector<D2BasFctType*> *d2Phi;
-
/** \brief
* Is used by function evalGrdUh. To make it thread safe, we need a
* temporary DimVec for every thread.
diff --git a/AMDiS/src/Boundary.h b/AMDiS/src/Boundary.h
index c5ef7f33a17acc86dd0f1c4b68d781f6ae17b9d3..c18ee78dba90b4b591c5fa87fa424d236b5b7409 100644
--- a/AMDiS/src/Boundary.h
+++ b/AMDiS/src/Boundary.h
@@ -22,14 +22,12 @@
#ifndef AMDIS_BOUNDARY_H
#define AMDIS_BOUNDARY_H
+#include <map>
#include "Global.h"
#include "MemoryManager.h"
-#include <map>
namespace AMDiS {
- template<typename T> class WorldVector;
-
/// Contains boundary constants
typedef enum {
INTERIOR = 0, /**< interior => no boundary (b = 0) */
@@ -41,116 +39,6 @@ namespace AMDiS {
/// Type specifier for the different boundary types
typedef signed int BoundaryType;
- // /** \ingroup Triangulation
- // * \brief
- // * Holds information about the type of boundary associated to an edge/face,
- // * and how new vertices are projected to the boundary in the case of curved
- // * boundaries.
- //
- // class Boundary
- // {
- // public:
- // MEMORY_MANAGED(Boundary);
-
- // /** \brief
- // * constructor
- //
- // Boundary(BoundaryType type=0) {
- // bound = type;
- // };
-
- // /** \brief
- // * copy constructor
- //
- // Boundary(const Boundary& old) { bound = old.getBound(); };
-
- // /** \brief
- // * destructor
- //
- // virtual ~Boundary() {};
-
- // /** \brief
- // * assignment operator
- //
- // Boundary& operator=(const Boundary& old) {
- // if (this!=&old) bound = old.getBound();
- // return *this;
- // };
-
- // /** \brief
- // * Returns
- // * -true: if a new vertex should be projected to a curved boundary
- // * -false: otherwise
- //
- // virtual bool interpolateBoundary() { return false; };
-
- // /** \brief
- // * Projection to the curved boundary
- //
- // virtual void interpolateBoundary(WorldVector<double>& ) {};
-
- // /** \brief
- // * Returns \ref bound
- //
- // inline const BoundaryType getBound() const { return bound; };
-
- // /** \brief
- // * Returns
- // * -true: is \ref bound is INTERIOR
- // * -false: otherwise
- //
- // inline const bool isInterior() const {return (bound == INTERIOR);};
-
- // /** \brief
- // * Returns
- // * -true: is \ref bound is DIRICHLET
- // * -false: otherwise
- //
- // inline const bool isDirichlet() const {return (bound == DIRICHLET);};
-
- // /** \brief
- // * Returns
- // * -true: is \ref bound is NEUMANN
- // * -false: otherwise
- //
- // inline const bool isNeumann() const {return (bound == NEUMANN);};
-
- // /** \brief
- // * Returns
- // * -true: is \ref bound is ROBIN
- // * -false: otherwise
- //
- // inline const bool isRobin() const {return (bound == ROBIN);};
-
- // /** \brief
- // * Returns the new value of \ref bound with respect to its old value and
- // * the value of bound_.
- //
- // BoundaryType newVal(const BoundaryType bound_);
-
- // /** \brief
- // * Returns the Boundary with given type from \ref boundaryMap.
- //
- // static Boundary* getBoundary(BoundaryType type);
-
- // /** \brief
- // * Adds Boundary b to \ref boundaryMap.
- //
- // //static void addBoundary(Boundary *b);
-
- // protected:
- // /** \brief
- // * type of this boundary
- //
- // BoundaryType bound;
-
- // protected:
- // /** \brief
- // * stl map of all existing boundaries.
- //
- // static ::std::map<BoundaryType, Boundary*> boundaryMap;
- // };
-
BoundaryType newBound(BoundaryType oldBound, BoundaryType newBound);
}
diff --git a/AMDiS/src/BoundaryCondition.h b/AMDiS/src/BoundaryCondition.h
index 7f403143b4d97a83ffecf1e8da82e71f6d9a1028..33d62b1aabf652e5dd229922c759405cd615532a 100644
--- a/AMDiS/src/BoundaryCondition.h
+++ b/AMDiS/src/BoundaryCondition.h
@@ -28,14 +28,6 @@
namespace AMDiS {
- template<typename T> class DOFVectorBase;
- class Estimator;
- class ElInfo;
-
- // ============================================================================
- // ===== class BoundaryCondition ==============================================
- // ============================================================================
-
/**
* \ingroup Assembler
*
@@ -43,12 +35,10 @@ namespace AMDiS {
* Sub class of BoundaryCondition. Local boundary conditions are filled
* while mesh traversal.
*/
- class BoundaryCondition //: public BoundaryCondition
+ class BoundaryCondition
{
public:
- /** \brief
- * Constructor.
- */
+ /// Constructor.
BoundaryCondition(BoundaryType type,
const FiniteElemSpace *rowFESpace_,
const FiniteElemSpace *colFESpace_ = NULL)
@@ -57,35 +47,33 @@ namespace AMDiS {
colFESpace(colFESpace_)
{
if(!colFESpace) colFESpace = rowFESpace;
- };
+ }
- /** \brief
- * Returns \ref boundaryType.
- */
- inline BoundaryType getBoundaryType() { return boundaryType; };
+ /// Returns \ref boundaryType.
+ inline BoundaryType getBoundaryType() {
+ return boundaryType;
+ }
- /** \brief
- * Returns \ref rowFESpace.
- */
- inline const FiniteElemSpace *getRowFESpace() { return rowFESpace; };
+ /// Returns \ref rowFESpace.
+ inline const FiniteElemSpace *getRowFESpace() {
+ return rowFESpace;
+ }
- /** \brief
- * Returns \ref rowFESpace.
- */
- inline const FiniteElemSpace *getColFESpace() { return colFESpace; };
+ /// Returns \ref rowFESpace.
+ inline const FiniteElemSpace *getColFESpace() {
+ return colFESpace;
+ }
- virtual void initMatrix(DOFMatrix*) {};
+ virtual void initMatrix(DOFMatrix*) {}
- virtual void exitMatrix(DOFMatrix*) {};
+ virtual void exitMatrix(DOFMatrix*) {}
- virtual void initVector(DOFVectorBase<double>*) {};
+ virtual void initVector(DOFVectorBase<double>*) {}
- virtual void exitVector(DOFVectorBase<double>*) {};
+ virtual void exitVector(DOFVectorBase<double>*) {}
- /** \brief
- * Destructor.
- */
- virtual ~BoundaryCondition() {};
+ /// Destructor.
+ virtual ~BoundaryCondition() {}
/** \brief
* Adds the local boundary condition for elInfo to object.
@@ -96,7 +84,7 @@ namespace AMDiS {
ElInfo *elInfo,
const DegreeOfFreedom *dofIndices,
const BoundaryType *localBound,
- int nBasFcts) {};
+ int nBasFcts) {}
/** \brief
* Adds the local boundary condition for elInfo to vector.
@@ -107,20 +95,23 @@ namespace AMDiS {
ElInfo *elInfo,
const DegreeOfFreedom *dofIndices,
const BoundaryType *localBound,
- int nBasFcts) {};
+ int nBasFcts) {}
- /** \brief
- * Returns the boundary residual for the given element. Called by estimator.
- */
+ /// Returns the boundary residual for the given element. Called by estimator.
virtual double boundResidual(ElInfo *elInfo,
DOFMatrix *matrix,
- const DOFVectorBase<double> *dv) { return 0.0; };
+ const DOFVectorBase<double> *dv)
+ {
+ return 0.0;
+ }
/** \brief
* Returns whether the condition must be treated as dirichlet condition
* while assemblage.
*/
- virtual bool isDirichlet() { return false; };
+ virtual bool isDirichlet() {
+ return false;
+ }
protected:
/** \brief
@@ -130,14 +121,10 @@ namespace AMDiS {
*/
BoundaryType boundaryType;
- /** \brief
- * FiniteElemSpace for this BoundaryCondition.
- */
+ /// FiniteElemSpace for this BoundaryCondition.
const FiniteElemSpace *rowFESpace;
- /** \brief
- * FiniteElemSpace for this BoundaryCondition.
- */
+ /// FiniteElemSpace for this BoundaryCondition.
const FiniteElemSpace *colFESpace;
};
diff --git a/AMDiS/src/BoundaryManager.h b/AMDiS/src/BoundaryManager.h
index 07e61713b2ea17db6deb5081579aaca9445bb710..6aa72c70af3b7dd23862b7fbded2ad4f6cc5b35d 100644
--- a/AMDiS/src/BoundaryManager.h
+++ b/AMDiS/src/BoundaryManager.h
@@ -30,14 +30,6 @@
namespace AMDiS {
- class FiniteElemSpace;
- template<typename T> class Vector;
- template<typename T> class DOFVectorBase;
-
- // ============================================================================
- // ===== class BoundaryManager ================================================
- // ============================================================================
-
/**
* \ingroup Assembler
*
@@ -49,17 +41,13 @@ namespace AMDiS {
class BoundaryManager
{
public:
- MEMORY_MANAGED(BoundaryManager);
-
BoundaryManager(const FiniteElemSpace *feSpace);
BoundaryManager(BoundaryManager &bm);
~BoundaryManager();
- /** \brief
- * Adds a local boundary condition to the list of managed conditions.
- */
+ /// Adds a local boundary condition to the list of managed conditions.
void addBoundaryCondition(BoundaryCondition *localBC) {
BoundaryType type = localBC->getBoundaryType();
TEST_EXIT(localBCs[type] == NULL)
diff --git a/AMDiS/src/Cholesky.h b/AMDiS/src/Cholesky.h
index a2d5eae7f70d146a0d702af221509627c6555296..fdd3b0df9c082b71215326709fb6cdea7b4ed0a5 100644
--- a/AMDiS/src/Cholesky.h
+++ b/AMDiS/src/Cholesky.h
@@ -27,10 +27,6 @@
using namespace std;
using namespace AMDiS;
-// ============================================================================
-// ===== class Cholesky =======================================================
-// ============================================================================
-
/**
* \ingroup Solvers
*
@@ -42,8 +38,6 @@ using namespace AMDiS;
class Cholesky
{
public:
- MEMORY_MANAGED(Cholesky);
-
/** \brief
* Computes the Cholesky factor of a positive definite matrix A.
* On input, only the upper triangle of A need be given; it is not modified.
diff --git a/AMDiS/src/DOFVector.cc b/AMDiS/src/DOFVector.cc
index 51be662520ce080ed3bdd9ac20fdf86a08853a39..35fac1a3d9c00b1c339e04f91ecc3a95c09aba43 100644
--- a/AMDiS/src/DOFVector.cc
+++ b/AMDiS/src/DOFVector.cc
@@ -1,3 +1,5 @@
+#include <boost/numeric/mtl/mtl.hpp>
+
#include "DOFVector.h"
#include "Traverse.h"
#include "DualTraverse.h"
@@ -16,12 +18,12 @@ namespace AMDiS {
(const_cast<BasisFunction*>(feSpace->getBasisFcts()))->coarseRestr(this, &list, n);
break;
case COARSE_INTERPOL:
- if (feSpace == NULL || !feSpace) {
+ if (feSpace == NULL || !feSpace)
ERROR_EXIT("ERR1\n");
- }
- if (feSpace->getBasisFcts() == NULL || !(feSpace->getBasisFcts())) {
+
+ if (feSpace->getBasisFcts() == NULL || !(feSpace->getBasisFcts()))
ERROR_EXIT("ERR2\n");
- }
+
(const_cast<BasisFunction*>(feSpace->getBasisFcts()))->coarseInter(this, &list, n);
break;
default:
@@ -341,7 +343,7 @@ namespace AMDiS {
double *localVec = localVectors[myRank];
getLocalVector(largeElInfo->getElement(), localVec);
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
+ const mtl::dense2D<double> m = smallElInfo->getSubElemCoordsMat();
DimVec<double> &grd1 = *grdTmp[myRank];
int parts = Global::getGeo(PARTS, dim);
@@ -362,7 +364,7 @@ namespace AMDiS {
grdPhi->set(0.0);
for (int k = 0; k < nBasFcts; k++) {
(*(basFcts->getGrdPhi(k)))(quad->getLambda(i), tmp);
- tmp *= (*m)[j][k];
+ tmp *= m[j][k];
*grdPhi += tmp;
}
@@ -802,14 +804,14 @@ namespace AMDiS {
double *localVec = localVectors[omp_get_thread_num()];
getLocalVector(largeElInfo->getElement(), localVec);
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
+ const mtl::dense2D<double> m = smallElInfo->getSubElemCoordsMat();
for (int i = 0; i < nPoints; i++) {
result[i] = 0.0;
for (int j = 0; j < nBasFcts; j++) {
double val = 0.0;
for (int k = 0; k < nBasFcts; k++) {
- val += (*m)[j][k] * (*(basFcts->getPhi(k)))(quad->getLambda(i));
+ val += m[j][k] * (*(basFcts->getPhi(k)))(quad->getLambda(i));
}
result[i] += localVec[j] * val;
}
diff --git a/AMDiS/src/DualTraverse.cc b/AMDiS/src/DualTraverse.cc
index 6b8fc52e5b795688ef6fca7276b24275196cf7ee..aed16a83a7e13d5abeb8572a28d7a9fb4dd4755b 100644
--- a/AMDiS/src/DualTraverse.cc
+++ b/AMDiS/src/DualTraverse.cc
@@ -170,17 +170,12 @@ namespace AMDiS {
if (!fillSubElemMat)
return;
- DimMat<double> *subCoordsMat = elInfoSmall->getSubElemCoordsMat();
- if (!subCoordsMat) {
- subCoordsMat = NEW DimMat<double>(elInfo1->getMesh()->getDim());
- }
-
+ mtl::dense2D<double>& subCoordsMat =
+ const_cast<mtl::dense2D<double>&>(elInfoSmall->getSubElemCoordsMat());
if (elInfo1 == elInfoSmall) {
stack1.getCoordsInElem(elInfo2, basisFcts, subCoordsMat);
} else {
stack2.getCoordsInElem(elInfo1, basisFcts, subCoordsMat);
}
-
- elInfoSmall->setSubElemCoordsMat(subCoordsMat);
}
}
diff --git a/AMDiS/src/ElInfo.cc b/AMDiS/src/ElInfo.cc
index 9bbbcac9b165c5c4060c0ef0726af69e31faed42..7431e21ca3c430688cf1467df3b80759802162ae 100644
--- a/AMDiS/src/ElInfo.cc
+++ b/AMDiS/src/ElInfo.cc
@@ -29,8 +29,7 @@ namespace AMDiS {
neighbourCoord_(mesh_->getDim(), NO_INIT),
oppVertex_(mesh_->getDim(), NO_INIT),
grdLambda(mesh_->getDim(), NO_INIT),
- subElemCoordsMat(NULL),
- subElemCoordsMat_mtl(2, 2)
+ subElemCoordsMat(2, 2)
{
projection_.set(NULL);
@@ -41,11 +40,7 @@ namespace AMDiS {
}
ElInfo::~ElInfo()
- {
- if (subElemCoordsMat)
- DELETE subElemCoordsMat;
- }
-
+ {}
void ElInfo::coordToWorld(const DimVec<double>& l,
WorldVector<double>& w) const
@@ -157,9 +152,8 @@ namespace AMDiS {
if (fillFlag_.isSet(Mesh::FILL_GRD_LAMBDA)) {
det_ = calcGrdLambda(grdLambda);
} else {
- if (fillFlag_.isSet(Mesh::FILL_DET)) {
+ if (fillFlag_.isSet(Mesh::FILL_DET))
det_ = calcDet();
- }
}
}
diff --git a/AMDiS/src/ElInfo.h b/AMDiS/src/ElInfo.h
index ca786c39143f40b18a1a1015a50cbef27295f034..4a9fe8628c74b0a63aba969709385c2644469e9a 100644
--- a/AMDiS/src/ElInfo.h
+++ b/AMDiS/src/ElInfo.h
@@ -209,14 +209,10 @@ namespace AMDiS {
* Returns the element transformation matrix \ref subElemCoordsMat .
* This value is set only during dual traverse.
*/
- inline DimMat<double> *getSubElemCoordsMat() const {
+ inline const mtl::dense2D<double>& getSubElemCoordsMat() const {
return subElemCoordsMat;
}
- inline mtl::dense2D<double>& getSubElemCoordsMat_mtl4() {
- return subElemCoordsMat_mtl;
- }
-
/** \} */
/** \name setting methods
@@ -287,9 +283,9 @@ namespace AMDiS {
* Sets the element transformation matrix \ref subElemCoordsMat .
* This value is used only during dual traverse.
*/
- inline void setSubElemCoordsMat(DimMat<double> *mat) {
- subElemCoordsMat = mat;
- }
+// inline void setSubElemCoordsMat(DimMat<double> *mat) {
+// subElemCoordsMat = mat;
+// }
/** \} */
@@ -376,16 +372,9 @@ namespace AMDiS {
return(0.0);
}
- virtual void getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const = 0;
-
virtual void getRefSimplexCoords(const BasisFunction *basisFcts,
mtl::dense2D<double>& coords) const = 0;
- virtual void getSubElementCoords(const BasisFunction *basisFcts,
- int iChild,
- DimMat<double> *coords) const = 0;
-
virtual void getSubElementCoords(const BasisFunction *basisFcts,
int iChild,
mtl::dense2D<double>& coords) const = 0;
@@ -489,9 +478,7 @@ namespace AMDiS {
* coordinates on the larger element, to the barycentric coordinates of the smaller
* element.
*/
- DimMat<double> *subElemCoordsMat;
-
- mtl::dense2D<double> subElemCoordsMat_mtl;
+ mtl::dense2D<double> subElemCoordsMat;
public:
/** \brief
diff --git a/AMDiS/src/ElInfo1d.cc b/AMDiS/src/ElInfo1d.cc
index 910ac24bcd28994bce53524242f44915943550e5..c8bca3d0d3ce8734f889232bc8a642575471fca2 100644
--- a/AMDiS/src/ElInfo1d.cc
+++ b/AMDiS/src/ElInfo1d.cc
@@ -13,6 +13,18 @@
#include "DOFVector.h"
namespace AMDiS {
+
+ double ElInfo1d::mat_d1_val[2][2] = {{1.0, 0.0},
+ {0.0, 1.0}};
+ mtl::dense2D<double> ElInfo1d::mat_d1(mat_d1_val);
+
+ double ElInfo1d::mat_d1_left_val[2][2] = {{1.0, 0.5},
+ {0.0, 0.5}};
+ mtl::dense2D<double> ElInfo1d::mat_d1_left(mat_d1_left_val);
+
+ double ElInfo1d::mat_d1_right_val[2][2] = {{0.5, 1.0},
+ {0.5, 0.0}};
+ mtl::dense2D<double> ElInfo1d::mat_d1_right(mat_d1_right_val);
void ElInfo1d::fillMacroInfo(const MacroElement * mel)
{
@@ -268,58 +280,12 @@ namespace AMDiS {
}
}
- void ElInfo1d::getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const
- {
- FUNCNAME("ElInfo1d::getRefSimplexCoords()");
-
- switch (basisFcts->getDegree()) {
- case 1:
- (*coords)[0][0] = 1.0;
- (*coords)[1][0] = 0.0;
-
- (*coords)[0][1] = 0.0;
- (*coords)[1][1] = 1.0;
- break;
- case 2:
- ERROR_EXIT("DDA\n");
- break;
- default:
- ERROR_EXIT("Not yet implemented!\n");
- }
- }
-
void ElInfo1d::getRefSimplexCoords(const BasisFunction *basisFcts,
mtl::dense2D<double>& coords) const
{
TEST_EXIT(basisFcts->getDegree() == 1)("Wrong basis function degree!\n");
- double deg1[][2] = {{ 1.0, 0.0 },
- { 0.0, 1.0 }};
- coords = deg1;
- }
-
- void ElInfo1d::getSubElementCoords(const BasisFunction *basisFcts,
- int iChild,
- DimMat<double> *coords) const
- {
- FUNCNAME("ElInfo1d::getSubElementCoords()");
-
- switch (basisFcts->getDegree()) {
- case 1:
- if (iChild == 0) {
- (*coords)[0][1] = ((*coords)[0][0] + (*coords)[0][1]) * 0.5;
- (*coords)[1][1] = ((*coords)[1][0] + (*coords)[1][1]) * 0.5;
- } else {
- (*coords)[0][0] = ((*coords)[0][0] + (*coords)[0][1]) * 0.5;
- (*coords)[1][0] = ((*coords)[1][0] + (*coords)[1][1]) * 0.5;
- }
- break;
- case 2:
- break;
- default:
- ERROR_EXIT("Not yet implemented!\n");
- }
+ coords = mat_d1;
}
void ElInfo1d::getSubElementCoords(const BasisFunction *basisFcts,
@@ -328,20 +294,12 @@ namespace AMDiS {
{
FUNCNAME("ElInfo1d::getSubElementCoords()");
- double deg1_left[][2] = {{ 1.0, 0.5 },
- { 0.0, 0.5 }};
- mtl::dense2D<double> deg1_left_mat(deg1_left);
-
- double deg1_right[][2] = {{ 0.5, 1.0 },
- { 0.5, 0.0 }};
- mtl::dense2D<double> deg1_right_mat(deg1_right);
-
switch (basisFcts->getDegree()) {
case 1:
if (iChild == 0)
- coords *= deg1_left_mat;
+ coords *= mat_d1_left;
else
- coords *= deg1_right_mat;
+ coords *= mat_d1_right;
break;
case 2:
diff --git a/AMDiS/src/ElInfo1d.h b/AMDiS/src/ElInfo1d.h
index dee1d7ecca1c5e39f210722b27f011b7804778ce..d6ca2741a24e4e93d1c5f5c13f0f2767b0f70529 100644
--- a/AMDiS/src/ElInfo1d.h
+++ b/AMDiS/src/ElInfo1d.h
@@ -63,19 +63,25 @@ namespace AMDiS {
return (i + 1) % 2;
}
- void getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const;
-
void getRefSimplexCoords(const BasisFunction *basisFcts,
mtl::dense2D<double>& coords) const;
- void getSubElementCoords(const BasisFunction *basisFcts,
- int iChild,
- DimMat<double> *coords) const;
-
void getSubElementCoords(const BasisFunction *basisFcts,
int iChild,
mtl::dense2D<double>& coords) const;
+
+ protected:
+ static double mat_d1_val[2][2];
+
+ static mtl::dense2D<double> mat_d1;
+
+ static double mat_d1_left_val[2][2];
+
+ static mtl::dense2D<double> mat_d1_left;
+
+ static double mat_d1_right_val[2][2];
+
+ static mtl::dense2D<double> mat_d1_right;
};
}
diff --git a/AMDiS/src/ElInfo2d.cc b/AMDiS/src/ElInfo2d.cc
index 5b5814f246e0bdd9e208735c7df194df40091f9a..5524fa4b1231684a61d71c680b04fd37a9d1a0f5 100644
--- a/AMDiS/src/ElInfo2d.cc
+++ b/AMDiS/src/ElInfo2d.cc
@@ -14,6 +14,21 @@
namespace AMDiS {
+ double ElInfo2d::mat_d1_val[3][3] = {{1.0, 0.0, 0.0},
+ {0.0, 1.0, 0.0},
+ {0.0, 0.0, 1.0}};
+ mtl::dense2D<double> ElInfo2d::mat_d1(mat_d1_val);
+
+ double ElInfo2d::mat_d1_left_val[3][3] = {{0.0, 1.0, 0.5},
+ {0.0, 0.0, 0.5},
+ {1.0, 0.0, 0.0}};
+ mtl::dense2D<double> ElInfo2d::mat_d1_left(mat_d1_left_val);
+
+ double ElInfo2d::mat_d1_right_val[3][3] = {{0.0, 0.0, 0.5},
+ {1.0, 0.0, 0.5},
+ {0.0, 1.0, 0.0}};
+ mtl::dense2D<double> ElInfo2d::mat_d1_right(mat_d1_right_val);
+
ElInfo2d::ElInfo2d(Mesh *aMesh)
: ElInfo(aMesh)
{
@@ -578,7 +593,6 @@ namespace AMDiS {
return det;
}
-
/****************************************************************************/
/* calculate the normal of the element for dim of world = dim + 1 */
/* return the absulute value of the determinant from the */
@@ -606,50 +620,32 @@ namespace AMDiS {
}
void ElInfo2d::getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const
+ mtl::dense2D<double>& coords) const
{
- (*coords)[0][0] = 1.0;
- (*coords)[1][0] = 0.0;
- (*coords)[2][0] = 0.0;
-
- (*coords)[0][1] = 0.0;
- (*coords)[1][1] = 1.0;
- (*coords)[2][1] = 0.0;
+ TEST_EXIT(basisFcts->getDegree() == 1)("Wrong basis function degree!\n");
- (*coords)[0][2] = 0.0;
- (*coords)[1][2] = 0.0;
- (*coords)[2][2] = 1.0;
+ coords = mat_d1;
}
void ElInfo2d::getSubElementCoords(const BasisFunction *basisFcts,
int iChild,
- DimMat<double> *coords) const
+ mtl::dense2D<double>& coords) const
{
- double c0 = ((*coords)[0][0] + (*coords)[0][1]) * 0.5;
- double c1 = ((*coords)[1][0] + (*coords)[1][1]) * 0.5;
- double c2 = ((*coords)[2][0] + (*coords)[2][1]) * 0.5;
-
- if (iChild == 0) {
- (*coords)[0][1] = (*coords)[0][0];
- (*coords)[1][1] = (*coords)[1][0];
- (*coords)[2][1] = (*coords)[2][0];
-
- (*coords)[0][0] = (*coords)[0][2];
- (*coords)[1][0] = (*coords)[1][2];
- (*coords)[2][0] = (*coords)[2][2];
- } else {
- (*coords)[0][0] = (*coords)[0][1];
- (*coords)[1][0] = (*coords)[1][1];
- (*coords)[2][0] = (*coords)[2][1];
-
- (*coords)[0][1] = (*coords)[0][2];
- (*coords)[1][1] = (*coords)[1][2];
- (*coords)[2][1] = (*coords)[2][2];
+ FUNCNAME("ElInfo1d::getSubElementCoords()");
+
+ switch (basisFcts->getDegree()) {
+ case 1:
+ if (iChild == 0)
+ coords *= mat_d1_left;
+ else
+ coords *= mat_d1_right;
+
+ break;
+ case 2:
+ break;
+ default:
+ ERROR_EXIT("Not yet implemented!\n");
}
-
- (*coords)[0][2] = c0;
- (*coords)[1][2] = c1;
- (*coords)[2][2] = c2;
}
}
diff --git a/AMDiS/src/ElInfo2d.h b/AMDiS/src/ElInfo2d.h
index 5e5cf6161adc52a4b37d8ecbb1f758dab6e8c717..1e69dc9ccad61be8fdfb0e511f1a28708ca023e5 100644
--- a/AMDiS/src/ElInfo2d.h
+++ b/AMDiS/src/ElInfo2d.h
@@ -63,23 +63,27 @@ namespace AMDiS {
double getElementNormal(WorldVector<double> &normal) const;
void getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const;
-
- void getRefSimplexCoords(const BasisFunction *basisFcts,
- mtl::dense2D<double>& coords) const {}
-
- void getSubElementCoords(const BasisFunction *basisFcts,
- int iChild,
- DimMat<double> *coords) const;
+ mtl::dense2D<double>& coords) const;
void getSubElementCoords(const BasisFunction *basisFcts,
int iChild,
- mtl::dense2D<double>& coords) const {}
-
+ mtl::dense2D<double>& coords) const;
protected:
/// Temp vectors for function \ref calcGrdLambda.
WorldVector<double> *e1, *e2, *normal;
+
+ static double mat_d1_val[3][3];
+
+ static mtl::dense2D<double> mat_d1;
+
+ static double mat_d1_left_val[3][3];
+
+ static mtl::dense2D<double> mat_d1_left;
+
+ static double mat_d1_right_val[3][3];
+
+ static mtl::dense2D<double> mat_d1_right;
};
}
diff --git a/AMDiS/src/ElInfo3d.cc b/AMDiS/src/ElInfo3d.cc
index 196b45eb3267202de03f251c6aff317a41597d31..d8ae8dab2eb0a5a14e80da91618d00ef93afbaf9 100644
--- a/AMDiS/src/ElInfo3d.cc
+++ b/AMDiS/src/ElInfo3d.cc
@@ -610,64 +610,64 @@ namespace AMDiS {
}
void ElInfo3d::getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const
+ mtl::dense2D<double>& coords) const
{
- (*coords)[0][0] = 1.0;
- (*coords)[1][0] = 0.0;
- (*coords)[2][0] = 0.0;
- (*coords)[3][0] = 0.0;
-
- (*coords)[0][1] = 0.0;
- (*coords)[1][1] = 1.0;
- (*coords)[2][1] = 0.0;
- (*coords)[3][1] = 0.0;
-
- (*coords)[0][2] = 0.0;
- (*coords)[1][2] = 0.0;
- (*coords)[2][2] = 1.0;
- (*coords)[3][2] = 0.0;
-
- (*coords)[0][3] = 0.0;
- (*coords)[1][3] = 0.0;
- (*coords)[2][3] = 0.0;
- (*coords)[3][3] = 1.0;
+// (*coords)[0][0] = 1.0;
+// (*coords)[1][0] = 0.0;
+// (*coords)[2][0] = 0.0;
+// (*coords)[3][0] = 0.0;
+
+// (*coords)[0][1] = 0.0;
+// (*coords)[1][1] = 1.0;
+// (*coords)[2][1] = 0.0;
+// (*coords)[3][1] = 0.0;
+
+// (*coords)[0][2] = 0.0;
+// (*coords)[1][2] = 0.0;
+// (*coords)[2][2] = 1.0;
+// (*coords)[3][2] = 0.0;
+
+// (*coords)[0][3] = 0.0;
+// (*coords)[1][3] = 0.0;
+// (*coords)[2][3] = 0.0;
+// (*coords)[3][3] = 1.0;
}
void ElInfo3d::getSubElementCoords(const BasisFunction *basisFcts,
int iChild,
- DimMat<double> *coords) const
+ mtl::dense2D<double>& coords) const
{
FUNCNAME("ElInfo3d::getSubElementCoords()");
ERROR_EXIT("Not yet\n");
- double c0 = ((*coords)[0][0] + (*coords)[0][1]) * 0.5;
- double c1 = ((*coords)[1][0] + (*coords)[1][1]) * 0.5;
- double c2 = ((*coords)[2][0] + (*coords)[2][1]) * 0.5;
- double c3 = ((*coords)[3][0] + (*coords)[3][1]) * 0.5;
-
- if (iChild == 0) {
- (*coords)[0][1] = (*coords)[0][0];
- (*coords)[1][1] = (*coords)[1][0];
- (*coords)[2][1] = (*coords)[2][0];
-
- (*coords)[0][0] = (*coords)[0][2];
- (*coords)[1][0] = (*coords)[1][2];
- (*coords)[2][0] = (*coords)[2][2];
- } else {
- (*coords)[0][1] = (*coords)[0][2];
- (*coords)[1][1] = (*coords)[1][2];
- (*coords)[2][1] = (*coords)[2][2];
-
- (*coords)[0][0] = (*coords)[0][1];
- (*coords)[1][0] = (*coords)[1][1];
- (*coords)[2][0] = (*coords)[2][1];
- }
-
- (*coords)[0][3] = c0;
- (*coords)[1][3] = c1;
- (*coords)[2][3] = c2;
- (*coords)[3][3] = c3;
+// double c0 = ((*coords)[0][0] + (*coords)[0][1]) * 0.5;
+// double c1 = ((*coords)[1][0] + (*coords)[1][1]) * 0.5;
+// double c2 = ((*coords)[2][0] + (*coords)[2][1]) * 0.5;
+// double c3 = ((*coords)[3][0] + (*coords)[3][1]) * 0.5;
+
+// if (iChild == 0) {
+// (*coords)[0][1] = (*coords)[0][0];
+// (*coords)[1][1] = (*coords)[1][0];
+// (*coords)[2][1] = (*coords)[2][0];
+
+// (*coords)[0][0] = (*coords)[0][2];
+// (*coords)[1][0] = (*coords)[1][2];
+// (*coords)[2][0] = (*coords)[2][2];
+// } else {
+// (*coords)[0][1] = (*coords)[0][2];
+// (*coords)[1][1] = (*coords)[1][2];
+// (*coords)[2][1] = (*coords)[2][2];
+
+// (*coords)[0][0] = (*coords)[0][1];
+// (*coords)[1][0] = (*coords)[1][1];
+// (*coords)[2][0] = (*coords)[2][1];
+// }
+
+// (*coords)[0][3] = c0;
+// (*coords)[1][3] = c1;
+// (*coords)[2][3] = c2;
+// (*coords)[3][3] = c3;
}
}
diff --git a/AMDiS/src/ElInfo3d.h b/AMDiS/src/ElInfo3d.h
index 47c935bff46e53998b80c1b8acd787ad870bc20a..0ff51c7b4a1ff9b5b6efd0f857eb51415393e4d1 100644
--- a/AMDiS/src/ElInfo3d.h
+++ b/AMDiS/src/ElInfo3d.h
@@ -92,18 +92,11 @@ namespace AMDiS {
}
void getRefSimplexCoords(const BasisFunction *basisFcts,
- DimMat<double> *coords) const;
-
- void getRefSimplexCoords(const BasisFunction *basisFcts,
- mtl::dense2D<double>& coords) const {}
-
- void getSubElementCoords(const BasisFunction *basisFcts,
- int iChild,
- DimMat<double> *coords) const;
+ mtl::dense2D<double>& coords) const;
void getSubElementCoords(const BasisFunction *basisFcts,
int iChild,
- mtl::dense2D<double>& coords) const {}
+ mtl::dense2D<double>& coords) const;
protected:
/// \ref el 's type. Is Filled automatically by the traversal routines.
diff --git a/AMDiS/src/FirstOrderAssembler.cc b/AMDiS/src/FirstOrderAssembler.cc
index d14a984699e66f5e8cc700eb070f2f57a0f19614..f6a951b6c67b1fcaa26b116e46803668cce2ec8d 100644
--- a/AMDiS/src/FirstOrderAssembler.cc
+++ b/AMDiS/src/FirstOrderAssembler.cc
@@ -1,4 +1,6 @@
#include <vector>
+#include <boost/numeric/mtl/mtl.hpp>
+
#include "Assembler.h"
#include "FirstOrderAssembler.h"
#include "Operator.h"
@@ -53,9 +55,8 @@ namespace AMDiS {
op->firstOrderGrdPhi[myRank];
// check if a assembler is needed at all
- if (opTerms.size() == 0) {
+ if (opTerms.size() == 0)
return NULL;
- }
sort(opTerms.begin(), opTerms.end());
@@ -162,27 +163,23 @@ namespace AMDiS {
const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
-
- TEST_EXIT(m)("No subElemCoordsMat!\n");
+ const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
int nPoints = quadrature->getNumPoints();
int myRank = omp_get_thread_num();
VectorOfFixVecs<DimVec<double> > &Lb = tmpLb[myRank];
Lb.resize(nPoints);
- for (int iq = 0; iq < nPoints; iq++) {
+ for (int iq = 0; iq < nPoints; iq++)
Lb[iq].set(0.0);
- }
for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++) {
(static_cast<FirstOrderTerm*>((terms[myRank][i])))->
getLb(smallElInfo, nPoints, Lb);
}
- for (int iq = 0; iq < nPoints; iq++) {
+ for (int iq = 0; iq < nPoints; iq++)
Lb[iq] *= smallElInfo->getDet();
- }
DimVec<double> val1(dim, DEFAULT_VALUE, 0.0);
@@ -194,7 +191,7 @@ namespace AMDiS {
for (int k = 0; k < nCol; k++) {
(*(psi->getGrdPhi(k)))(quadrature->getLambda(iq), grdPsi);
- grdPsi *= (*m)[i][k];
+ grdPsi *= m[i][k];
val1 += grdPsi;
}
@@ -423,28 +420,22 @@ namespace AMDiS {
const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
-
- TEST_EXIT(m)("No subElemCoordsMat!\n");
+ const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
int nPoints = quadrature->getNumPoints();
int myRank = omp_get_thread_num();
VectorOfFixVecs<DimVec<double> > &Lb = tmpLb[myRank];
Lb.resize(nPoints);
- for (int iq = 0; iq < nPoints; iq++) {
+ for (int iq = 0; iq < nPoints; iq++)
Lb[iq].set(0.0);
- }
- for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++) {
+ for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++)
(static_cast<FirstOrderTerm*>((terms[myRank][i])))->
getLb(smallElInfo, nPoints, Lb);
- }
- for (int iq = 0; iq < nPoints; iq++) {
+ for (int iq = 0; iq < nPoints; iq++)
Lb[iq] *= smallElInfo->getDet();
- }
-
for (int i = 0; i < nRow; i++) {
for (int j = 0; j < nCol; j++) {
@@ -452,9 +443,8 @@ namespace AMDiS {
double val0 = 0.0;
- for (int k = 0; k < nCol; k++) {
- val0 = (*m)[i][k] * (*(psi->getPhi(k)))(quadrature->getLambda(iq));
- }
+ for (int k = 0; k < nCol; k++)
+ val0 = m[i][k] * (*(psi->getPhi(k)))(quadrature->getLambda(iq));
(*(phi->getGrdPhi(j)))(quadrature->getLambda(iq), grdPhi);
diff --git a/AMDiS/src/SecondOrderAssembler.cc b/AMDiS/src/SecondOrderAssembler.cc
index 727edcce5d136f1d7633404ba5233c31e660cc2e..fd932bb7d6665ef2974c448354f9824adf1c9a7b 100644
--- a/AMDiS/src/SecondOrderAssembler.cc
+++ b/AMDiS/src/SecondOrderAssembler.cc
@@ -1,4 +1,6 @@
#include <vector>
+#include <boost/numeric/mtl/mtl.hpp>
+
#include "Assembler.h"
#include "SecondOrderAssembler.h"
#include "Operator.h"
@@ -30,9 +32,8 @@ namespace AMDiS {
int myRank = omp_get_thread_num();
// check if a assembler is needed at all
- if (op->secondOrder[myRank].size() == 0) {
+ if (op->secondOrder[myRank].size() == 0)
return NULL;
- }
SecondOrderAssembler *newAssembler;
@@ -350,9 +351,7 @@ namespace AMDiS {
const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
-
- TEST_EXIT(m)("No subElemCoordsMat!\n");
+ const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
int nPoints = quadrature->getNumPoints();
@@ -398,7 +397,7 @@ namespace AMDiS {
for (int k = 0; k < nCol; k++) {
(*(phi->getGrdPhi(k)))(quadrature->getLambda(iq), grdPsi);
- grdPsi *= (*m)[j][k];
+ grdPsi *= m[j][k];
val1 += grdPsi;
}
diff --git a/AMDiS/src/Tetrahedron.h b/AMDiS/src/Tetrahedron.h
index 8a2c4de5b0affe6f5dd9283929f1a73f2c39a731..09c747e7de97e844dd2fc9ca2f6ad902bb7acef1 100644
--- a/AMDiS/src/Tetrahedron.h
+++ b/AMDiS/src/Tetrahedron.h
@@ -38,30 +38,23 @@ namespace AMDiS {
class Tetrahedron : public Element
{
public:
- MEMORY_MANAGED(Tetrahedron);
+
+ /// calls base class contructor.
+ Tetrahedron(Mesh* aMesh)
+ : Element(aMesh)
+ {}
- /** \brief
- * calls base class contructor.
- */
- Tetrahedron(Mesh* aMesh) : Element(aMesh) {};
-
- /** \brief
- * implements Element::clone
- */
+ /// implements Element::clone
inline Element *clone() {
return NEW Tetrahedron(mesh);
- };
+ }
- /** \brief
- * implements Element::getVertexOfEdge
- */
+ /// implements Element::getVertexOfEdge
inline int getVertexOfEdge(int i, int j) const {
return vertexOfEdge[i][j];
- };
+ }
- /** \brief
- * implements Element::getVertexOfPosition
- */
+ /// implements Element::getVertexOfPosition
int getVertexOfPosition(GeoIndex position,
int positionIndex,
int vertexIndex) const;
@@ -73,11 +66,9 @@ namespace AMDiS {
{0,1,-1,2},
{0,1,2,-1}};
return positionOfVertex[side][vertex];
- };
+ }
- /** \brief
- * implements Element::getGeo
- */
+ /// implements Element::getGeo
inline int getGeo(GeoIndex i) const {
switch(i) {
case VERTEX: case PARTS: case NEIGH:
@@ -103,26 +94,19 @@ namespace AMDiS {
ERROR_EXIT("invalid geo-index\n");
return 0;
}
- };
+ }
- /** \brief
- * implements Element::hasSide
- */
+ /// implements Element::hasSide
bool hasSide(Element* sideElem) const;
- /** \brief
- * implements Element::sortFaceIndices
- */
+ /// implements Element::sortFaceIndices
const FixVec<int,WORLD>& sortFaceIndices(int face,
FixVec<int,WORLD> *vec) const;
- /** \brief
- * implements Element::isLine. Returns false because this element is a
- * Tetrahedron
- */
+ /// implements Element::isLine. Returns false because this element is a Tetrahedron
inline bool isLine() const {
return false;
- };
+ }
/** \brief
* implements Element::isTriangle. Returns false because this element is a
@@ -130,7 +114,7 @@ namespace AMDiS {
*/
inline bool isTriangle() const {
return false;
- };
+ }
/** \brief
* implements Element::isTetrahedron. Returns true because this element is a
@@ -138,13 +122,11 @@ namespace AMDiS {
*/
inline bool isTetrahedron() const {
return true;
- };
-
- // ===== Serializable implementation =====
-
+ }
+
std::string getTypeName() const {
return "Tetrahedron";
- };
+ }
public:
/** \brief
@@ -192,61 +174,48 @@ namespace AMDiS {
*/
static const signed char childOrientation[3][2];
- /** \brief
- * edgeOfDOFs[i][j]: gives the local index of edge with vertices i and j
- */
+ /// edgeOfDOFs[i][j]: gives the local index of edge with vertices i and j
static const unsigned char edgeOfDOFs[4][4];
+ ///
static const int edgeOfFace[4][3];
- /** \brief
- * implements Element::getSideOfChild()
- */
+ /// implements Element::getSideOfChild()
virtual int getSideOfChild(int child, int side, int elType = 0) const {
FUNCNAME("Tetrahedron::getSideOfChild()");
TEST_EXIT_DBG(child==0 || child==1)("child must be in (0,1)\n");
TEST_EXIT_DBG(side >= 0 && side <= 3)("side must be between 0 and 3\n");
TEST_EXIT_DBG(elType >= 0 && elType <= 2)("elType must be between 0 and 2\n");
return sideOfChild[elType][child][side];
- };
+ }
- /** \brief
- * implements Element::getVertexOfParent()
- */
+ /// implements Element::getVertexOfParent()
virtual int getVertexOfParent(int child, int side, int elType = 0) const {
FUNCNAME("Tetrahedron::getVertexOfParent()");
TEST_EXIT_DBG(child==0 || child==1)("child must be in (0,1)\n");
TEST_EXIT_DBG(side >= 0 && side <= 3)("side must be between 0 and 3\n");
TEST_EXIT_DBG(elType >= 0 && elType <= 2)("elType must be between 0 and 2\n");
return vertexOfParent[elType][child][side];
- };
+ }
inline int getEdgeOfFace(int face, int edge) const {
TEST_EXIT_DBG(face >= 0 && face < 4)("invalid face\n");
TEST_EXIT_DBG(edge >= 0 && edge < 3)("invalid edge\n");
return edgeOfFace[face][edge];
- };
+ }
protected:
- /** \brief
- * vertexOfEdge[i][j] is the local number of the j-vertex of edge i
- */
+ /// vertexOfEdge[i][j] is the local number of the j-vertex of edge i
static const int vertexOfEdge[6][2];
- /** \brief
- * vertexOfFace[i][j] is the local number of the j-vertex of face i
- */
+ /// vertexOfFace[i][j] is the local number of the j-vertex of face i
static const int vertexOfFace[4][3];
- /** \brief
- *
- */
+ ///
static const int sideOfChild[3][2][4];
- /** \brief
- *
- */
+ ///
static const int vertexOfParent[3][2][4];
};
diff --git a/AMDiS/src/Traverse.cc b/AMDiS/src/Traverse.cc
index c878898f875be9587c5b77d700bd092d9efc2cb4..b972ec70b400e4edfc360ddca66c96f0d590df19 100644
--- a/AMDiS/src/Traverse.cc
+++ b/AMDiS/src/Traverse.cc
@@ -1068,20 +1068,6 @@ namespace AMDiS {
dynamic_cast<ElInfo3d*>(elinfo_stack[i])->update();
}
- void TraverseStack::getCoordsInElem(const ElInfo *upperElInfo,
- const BasisFunction *basisFcts,
- DimMat<double> *coords)
- {
- int levelDif = elinfo_stack[stack_used]->getLevel() - upperElInfo->getLevel();
-
- upperElInfo->getRefSimplexCoords(basisFcts, coords);
-
- for (int i = 1; i <= levelDif; i++)
- upperElInfo->getSubElementCoords(basisFcts,
- elinfo_stack[stack_used - levelDif + i]->getIChild(),
- coords);
- }
-
void TraverseStack::getCoordsInElem(const ElInfo *upperElInfo,
const BasisFunction *basisFcts,
mtl::dense2D<double>& coords)
diff --git a/AMDiS/src/Traverse.h b/AMDiS/src/Traverse.h
index 65d66ee6d36cdc1d46520c42c1f653fab0505ab7..2a16026885dc2b22ec795fa8ff6fcaccd362ca7d 100644
--- a/AMDiS/src/Traverse.h
+++ b/AMDiS/src/Traverse.h
@@ -117,9 +117,6 @@ namespace AMDiS {
/// Only for 3d: Calls update of all ElInfo3d onjects in \ref elinfo_stack
void update();
- void getCoordsInElem(const ElInfo *upperElInfo, const BasisFunction *basisFcts,
- DimMat<double> *coords);
-
void getCoordsInElem(const ElInfo *upperElInfo, const BasisFunction *basisFcts,
mtl::dense2D<double>& coords);
diff --git a/AMDiS/src/ZeroOrderAssembler.cc b/AMDiS/src/ZeroOrderAssembler.cc
index cdd5f726f5d6e44d26f0ab624ba4c96ab3d5d842..840c4fe3e7246d5eb96cbcd9c7be6882c8a69949 100644
--- a/AMDiS/src/ZeroOrderAssembler.cc
+++ b/AMDiS/src/ZeroOrderAssembler.cc
@@ -1,4 +1,6 @@
#include <vector>
+#include <boost/numeric/mtl/mtl.hpp>
+
#include "Assembler.h"
#include "ZeroOrderAssembler.h"
#include "Operator.h"
@@ -27,9 +29,8 @@ namespace AMDiS {
bool optimized)
{
// check if a assembler is needed at all
- if (op->zeroOrder.size() == 0) {
- return NULL;
- }
+ if (op->zeroOrder.size() == 0)
+ return NULL;
ZeroOrderAssembler *newAssembler;
@@ -151,18 +152,15 @@ namespace AMDiS {
const BasisFunction *psi = owner->getRowFESpace()->getBasisFcts();
const BasisFunction *phi = owner->getColFESpace()->getBasisFcts();
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
-
- TEST_EXIT(m)("No subElemCoordsMat!\n");
+ const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
int myRank = omp_get_thread_num();
int nPoints = quadrature->getNumPoints();
std::vector<double> c(nPoints, 0.0);
std::vector<OperatorTerm*>::iterator termIt;
- for (termIt = terms[myRank].begin(); termIt != terms[myRank].end(); ++termIt) {
+ for (termIt = terms[myRank].begin(); termIt != terms[myRank].end(); ++termIt)
(static_cast<ZeroOrderTerm*>((*termIt)))->getC(smallElInfo, nPoints, c);
- }
for (int i = 0; i < nRow; i++) {
for (int j = 0; j < nCol; j++) {
@@ -172,15 +170,13 @@ namespace AMDiS {
(*(phi->getPhi(i)))(quadrature->getLambda(iq));
double val1 = 0.0;
- for (int k = 0; k < nCol; k++) {
- val1 += (*m)[j][k] * (*(psi->getPhi(k)))(quadrature->getLambda(iq));
- }
+ for (int k = 0; k < nCol; k++)
+ val1 += m[j][k] * (*(psi->getPhi(k)))(quadrature->getLambda(iq));
val0 *= val1;
(*mat)[i][j] += val0;
}
-
}
}
}
@@ -279,7 +275,7 @@ namespace AMDiS {
int nPoints = quadrature->getNumPoints();
int myRank = omp_get_thread_num();
- DimMat<double> *m = smallElInfo->getSubElemCoordsMat();
+ const mtl::dense2D<double>& m = smallElInfo->getSubElemCoordsMat();
if (firstCall) {
#ifdef _OPENMP
@@ -312,7 +308,7 @@ namespace AMDiS {
double tmpval = 0.0;
for (int k = 0; k < nCol; k++) {
- tmpval += (*m)[j][k] * phi[k];
+ tmpval += m[j][k] * phi[k];
}
val *= tmpval;