diff --git a/AMDiS/CMakeLists.txt b/AMDiS/CMakeLists.txt
index d7822b4b887157d57048d1a9168195c8c629d64d..c98889a9c23c622b9d72f3c0a83dc160c856ebc3 100644
--- a/AMDiS/CMakeLists.txt
+++ b/AMDiS/CMakeLists.txt
@@ -492,8 +492,7 @@ if(ENABLE_EXTENSIONS)
${EXTENSIONS_DIR}/base_problems/NavierStokes_TaylorHood_RB.cc
${EXTENSIONS_DIR}/base_problems/NavierStokes_TH_MultiPhase.cc
${EXTENSIONS_DIR}/base_problems/NavierStokes_TH_MultiPhase_RB.cc
- ${EXTENSIONS_DIR}/base_problems/PhaseFieldCrystal_Base.cc
- # ${EXTENSIONS_DIR}/base_problems/PhaseFieldCrystal.cc
+ ${EXTENSIONS_DIR}/base_problems/PhaseFieldCrystal.cc
${EXTENSIONS_DIR}/base_problems/PhaseFieldCrystal_Phase.cc
${EXTENSIONS_DIR}/base_problems/PhaseFieldCrystal_RB.cc)
list(APPEND COMPILEFLAGS "-DHAVE_BASE_PROBLEMS=1")
diff --git a/AMDiS/src/AbstractFunction.h b/AMDiS/src/AbstractFunction.h
index dfbc2f297b3fd0476d62ff164e28e2ac06f3080b..233fc2a3df822a8fd4bbd5063e7c6da55f5baa7a 100644
--- a/AMDiS/src/AbstractFunction.h
+++ b/AMDiS/src/AbstractFunction.h
@@ -27,6 +27,10 @@
#include "Global.h"
+#include <boost/preprocessor/arithmetic/inc.hpp>
+#include <boost/preprocessor/repetition/enum_params.hpp>
+#include <boost/preprocessor/repetition/repeat.hpp>
+
namespace AMDiS {
/**
@@ -170,6 +174,53 @@ namespace AMDiS {
};
+ ///////////////////////////////////////////////////////////////
+ // test of AbstractFunction with arbitrary number of arguments
+
+ #define ABSTRACT_FUNCTION_MACRO(z, n, _) \
+ template< typename ReturnType, \
+ BOOST_PP_ENUM_PARAMS_Z(z, n, typename T) \
+ > class AbstractFunction ## n { \
+ AbstractFunction ## n (int degree = 0) : \
+ degree_(degree) \
+ {} \
+ virtual ~AbstractFunction ## n () {} \
+ inline int getDegree() const \
+ { \
+ return degree_; \
+ } \
+ virtual ReturnType operator()(BOOST_PP_ENUM_BINARY_PARAMS_Z(z, n, const T, &t)) const = 0; \
+ \
+ protected: \
+ int degree_; \
+ }; \
+
+// template< typename ReturnType,
+// BOOST_PP_ENUM_PARAMS(N, typename T) // expands to typename T0, typename T1, typename T2...
+// >
+// class CONCAT_STR(AbstractFunction,N)
+// {
+// public:
+// CONCAT_STR(AbstractFunction,N)(int degree = 0) :
+// degree_(degree)
+// {}
+//
+// virtual ~CONCAT_STR(AbstractFunction,N)() {}
+//
+// /// Returns \ref degree_.
+// inline int getDegree() const
+// {
+// return degree_;
+// }
+//
+// /// function evaluation.
+// virtual ReturnType operator()(BOOST_PP_ENUM_BINARY_PARAMS(N, const T, &t)) const = 0;
+//
+// protected:
+// int degree_;
+// };
+
+ BOOST_PP_REPEAT_FROM_TO(1, 11, ABSTRACT_FUNCTION_MACRO, nil)
}
diff --git a/AMDiS/src/AdaptInstationary.cc b/AMDiS/src/AdaptInstationary.cc
index 6b6372dcdd3240aeedf8ff3088dd120934620b36..ef1d7ea00376f3d317c57e5de6c16e4f61f1590e 100644
--- a/AMDiS/src/AdaptInstationary.cc
+++ b/AMDiS/src/AdaptInstationary.cc
@@ -192,7 +192,7 @@ namespace AMDiS {
adaptInfo->setSpaceIteration(0);
- // === Do only space iterations only if the maximum is higher than 0. ===
+ // === Do space iterations only if the maximum is higher than 0. ===
if (adaptInfo->getMaxSpaceIteration() > 0) {
diff --git a/AMDiS/src/Global.h b/AMDiS/src/Global.h
index c0af4d7b31fb1256bec86cc3148a8a902b3a80cb..920b9ccf52f815b689fac661f9e45e8dbaba96db 100644
--- a/AMDiS/src/Global.h
+++ b/AMDiS/src/Global.h
@@ -65,6 +65,7 @@
#include <boost/algorithm/string.hpp>
#include <boost/algorithm/string/trim.hpp>
+#include <boost/math/special_functions/fpclassify.hpp>
#include "boost/tuple/tuple.hpp"
#include "AMDiS_fwd.h"
#include "OpenMP.h"
@@ -147,6 +148,12 @@ namespace AMDiS {
return access(filename.c_str(), F_OK) == 0;
#endif
}
+
+ /// check for inf and nan values
+ inline bool isNumber(double val)
+ {
+ return !boost::math::isnan(val) && !boost::math::isinf(val);
+ }
/// trim std::string
diff --git a/AMDiS/src/ProblemStat.cc b/AMDiS/src/ProblemStat.cc
index 6915927b684547d7e4242f1a8f139bd5679c86ee..618dd4e4dda2a688256c9bf358a5c7f51aca49eb 100644
--- a/AMDiS/src/ProblemStat.cc
+++ b/AMDiS/src/ProblemStat.cc
@@ -514,6 +514,7 @@ namespace AMDiS {
{
FUNCNAME("ProblemStat::createSolver()");
+ // definition of standard-backends
#if defined HAVE_PARALLEL_PETSC
string backend("p_petsc");
#elif defined HAVE_PARALLEL_MTL
@@ -576,7 +577,7 @@ namespace AMDiS {
for (int j = 0; j < nComponents; j++)
estimator[i]->addSystem((*systemMatrix)[i][j],
solution->getDOFVector(j),
- rhs->getDOFVector(j));
+ rhs->getDOFVector(j)); // TODO: hier eventuell (i) statt (j) ???
}
}
@@ -625,7 +626,7 @@ namespace AMDiS {
solutionList));
}
- // Create own filewriters for each components of the problem
+ // Create own filewriters for each component of the problem
for (int i = 0; i < nComponents; i++) {
numberedName = name + "->output[" + boost::lexical_cast<string>(i) + "]";
filename = "";
@@ -659,13 +660,12 @@ namespace AMDiS {
return;
}
- clock_t first = clock();
-
+ Timer t;
solver->solveSystem(solverMatrix, *solution, *rhs,
createMatrixData, storeMatrixData);
INFO(info, 8)("solution of discrete system needed %.5f seconds\n",
- TIME_USED(first, clock()));
+ t.elapsed());
adaptInfo->setSolverIterations(solver->getIterations());
adaptInfo->setMaxSolverIterations(solver->getMaxIterations());
@@ -678,11 +678,7 @@ namespace AMDiS {
{
FUNCNAME("ProblemStat::estimate()");
-#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
- double first = MPI::Wtime();
-#else
- clock_t first = clock();
-#endif
+ Timer t;
if (computeExactError) {
computeError(adaptInfo);
@@ -708,12 +704,9 @@ namespace AMDiS {
#ifdef HAVE_PARALLEL_DOMAIN_AMDIS
MPI::COMM_WORLD.Barrier();
- INFO(info, 8)("estimation of the error needed %.5f seconds\n",
- MPI::Wtime() - first);
-#else
- INFO(info, 8)("estimation of the error needed %.5f seconds\n",
- TIME_USED(first, clock()));
#endif
+ INFO(info, 8)("estimation of the error needed %.5f seconds\n",
+ t.elapsed());
}
@@ -1506,8 +1499,7 @@ namespace AMDiS {
if (matrix) {
matrix->assemble(1.0, elInfo, bound);
- // Take the matrix boundary manager from the public matrix,
- // but assemble the boundary conditions on the thread private matrix.
+ // assemble the boundary conditions on the matrix.
if (matrix->getBoundaryManager())
matrix->getBoundaryManager()->fillBoundaryConditions(elInfo, matrix);
}
@@ -1518,9 +1510,6 @@ namespace AMDiS {
elInfo = stack.traverseNext(elInfo);
}
- // == Finally, if we have assembled in parallel, we have to add the thread ==
- // == private matrix and vector to the global one. ==
-
if (matrix) {
matrix->clearDirichletRows();
matrix->finishAssembling();
diff --git a/AMDiS/src/solver/LinearSolver.h b/AMDiS/src/solver/LinearSolver.h
index 2c3b5060ce13695985e3d5820c0468a9f51eee99..555126ce79851b3b723377e62c3d786949bd93c1 100644
--- a/AMDiS/src/solver/LinearSolver.h
+++ b/AMDiS/src/solver/LinearSolver.h
@@ -95,7 +95,7 @@ namespace AMDiS {
print_cycle(100),
iterations(-1),
error(-1),
- breakTolNotReached(false),
+ breakTolNotReached(true),
calculateResidual(false)
{
Parameters::get(name + "->tolerance", tolerance);
@@ -136,7 +136,7 @@ namespace AMDiS {
MSG("Residual norm: ||b-Ax|| = %e\n", residual);
- TEST_EXIT(residual <= tolerance || !breakTolNotReached)
+ TEST_EXIT((isNumber(residual) && residual <= tolerance) || !breakTolNotReached)
("Tolerance tol = %e could not be reached!\n Set tolerance by '->solver->tolerance:' \n", tolerance);
}
return error_code;
diff --git a/AMDiS/src/solver/UmfPackSolver.h b/AMDiS/src/solver/UmfPackSolver.h
index f674180375661be49c9d4edbf5b67123d4e52bd3..02f4accafd8ae444de367bd6e0618154c192e0b4 100644
--- a/AMDiS/src/solver/UmfPackSolver.h
+++ b/AMDiS/src/solver/UmfPackSolver.h
@@ -67,13 +67,12 @@ namespace AMDiS {
int code = (*solver)(x, b);
- if (oem.getInfo() > 0) {
- VectorType r(b);
- r -= A * x;
- double residual = two_norm(r);
- oem.setResidual(residual);
- oem.setErrorCode(code);
- }
+ VectorType r(b);
+ r -= A * x;
+ double residual = two_norm(r);
+ oem.setResidual(residual);
+ oem.setErrorCode(code);
+
return code;
}