diff --git a/AMDiS/src/BiCGSolver.h b/AMDiS/src/BiCGSolver.h index 2094d71780597ec2926b60cbd0cee33eca000714..e428ad4a1abe6546d778affb6b5b39bd6d277b04 100644 --- a/AMDiS/src/BiCGSolver.h +++ b/AMDiS/src/BiCGSolver.h @@ -77,7 +77,7 @@ namespace AMDiS { /** \brief * realisation of OEMSolver::solveSystem */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/BiCGSolver.hh b/AMDiS/src/BiCGSolver.hh index b540eb3f612200ce857973389734a5bf14172703..14965297745c106903a597f87a95e3423366d5d9 100644 --- a/AMDiS/src/BiCGSolver.hh +++ b/AMDiS/src/BiCGSolver.hh @@ -35,7 +35,7 @@ namespace AMDiS { template<typename VectorType> int BiCGSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("BiCGSolver::solveSystem"); VectorType *t = h; diff --git a/AMDiS/src/BiCGStab.h b/AMDiS/src/BiCGStab.h index fc8e6ce6b0523b4bb6eb3b0251d172a41826a3ee..e965c80d73cbbf0ab991d4fbca0098bc6eef0ec7 100644 --- a/AMDiS/src/BiCGStab.h +++ b/AMDiS/src/BiCGStab.h @@ -80,7 +80,7 @@ namespace AMDiS * realisation of OEMSolver::solveSystem */ int solveSystem(MatVecMultiplier<VectorType> *mv, - VectorType *x, VectorType *b); + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/BiCGStab.hh b/AMDiS/src/BiCGStab.hh index 1058576d4cef6cdb571596b7fd5ee05ee7702fa5..7479cfdfdedede1daac01ad3a980151a51012aed 100644 --- a/AMDiS/src/BiCGStab.hh +++ b/AMDiS/src/BiCGStab.hh @@ -42,7 +42,7 @@ namespace AMDiS template<typename VectorType> int BiCGStab<VectorType>::solveSystem(MatVecMultiplier<VectorType> *mv, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("BiCGStab::solveSystem"); diff --git a/AMDiS/src/BiCGStab2.h b/AMDiS/src/BiCGStab2.h index 8767dd964b8463698cfb5adf38a4d033a824e35b..054a5bfcc4ce149ba5513223704c7550a7d2f4d6 100644 --- a/AMDiS/src/BiCGStab2.h +++ b/AMDiS/src/BiCGStab2.h @@ -79,7 +79,7 @@ namespace AMDiS * realisation of OEMSolver::solveSystem */ int solveSystem(MatVecMultiplier<VectorType> *mv, - VectorType *x, VectorType *b); + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/BiCGStab2.hh b/AMDiS/src/BiCGStab2.hh index 139fff73e4d294dc150007c56b5b611408f542ce..e0bdedbc5b461b0267ab33c3aed2208e419241e0 100644 --- a/AMDiS/src/BiCGStab2.hh +++ b/AMDiS/src/BiCGStab2.hh @@ -42,7 +42,7 @@ namespace AMDiS template<typename VectorType> int BiCGStab2<VectorType>::solveSystem(MatVecMultiplier<VectorType> *mv, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("BiCGStab2::solveSystem()"); diff --git a/AMDiS/src/BiCGStab_M.h b/AMDiS/src/BiCGStab_M.h index 1cada652dcc65282299e0eba2b97bc918900bda9..c9f36d0f82a31acae12879e801efa76fadd19817 100644 --- a/AMDiS/src/BiCGStab_M.h +++ b/AMDiS/src/BiCGStab_M.h @@ -84,7 +84,7 @@ namespace AMDiS { /** \brief * Implements OEMSolver<VectorType>::solve(). */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b, bool reuseMatrix); private: // Pointers to memory needed for solveSystem diff --git a/AMDiS/src/BiCGStab_M.hh b/AMDiS/src/BiCGStab_M.hh index aee65d69d1cfe5f4eb5651b3e2aab42ec171a461..57789557398fef54be5a818079b1145465cfe4c2 100644 --- a/AMDiS/src/BiCGStab_M.hh +++ b/AMDiS/src/BiCGStab_M.hh @@ -34,7 +34,7 @@ namespace AMDiS { template<typename VectorType> int BiCGStab_M<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("BiCGStab_M::solveSystem"); double old_res = -1.0; diff --git a/AMDiS/src/CGSolver.h b/AMDiS/src/CGSolver.h index 81ece30489b6081c5734c5c9475c42502b9e3041..f39a96907cb58a737854e8f5a1cab9289c1c9be8 100644 --- a/AMDiS/src/CGSolver.h +++ b/AMDiS/src/CGSolver.h @@ -102,7 +102,7 @@ namespace AMDiS { /** \brief * Implements OEMSolver<VectorType>::solve(). */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b, bool reuseMatrix); private: VectorType *p, *r, *v, *z; diff --git a/AMDiS/src/CGSolver.hh b/AMDiS/src/CGSolver.hh index ed43b1df1305add1fa49edc1ef98041162f139eb..7b1fa1e71ea145dafc2cc274d92a8da4ce67a252 100644 --- a/AMDiS/src/CGSolver.hh +++ b/AMDiS/src/CGSolver.hh @@ -16,7 +16,7 @@ namespace AMDiS { template<typename VectorType> int CGSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("CGSolver::solve()"); diff --git a/AMDiS/src/Flag.h b/AMDiS/src/Flag.h index 404643818eed061897a062b2bf166940b1512c6d..d393a94c702eebec78aad7c480861a10d5fa6353 100644 --- a/AMDiS/src/Flag.h +++ b/AMDiS/src/Flag.h @@ -41,36 +41,36 @@ namespace AMDiS { /** \brief * Constructs a unset Flag */ - inline Flag() : flags(0) {}; + inline Flag() : flags(0) {} /** \brief * Constructs a Flag initialized by f */ - inline Flag(const unsigned long f) : flags(f) {}; + inline Flag(const unsigned long f) : flags(f) {} /** \brief * Copy constructor */ - inline Flag(const Flag& f) : flags(f.flags) {}; + inline Flag(const Flag& f) : flags(f.flags) {} /** \brief * Destructor */ - inline ~Flag() {}; + inline ~Flag() {} /** \brief * Compares two Flags */ inline bool operator==(const Flag& f) const { return (flags == f.flags); - }; + } /** \brief * Compares two Flags */ inline bool operator!=(const Flag& f) const { return !(f == *this); - }; + } /** \brief * Assignment operator @@ -79,60 +79,60 @@ namespace AMDiS { if (this != &f) flags = f.flags; return *this; - }; + } /** \brief * Typecast */ inline operator bool() const { return isAnySet(); - }; + } /** \brief * Set \ref flags */ inline void setFlags(const unsigned long f) { flags = f; - }; + } /** \brief * Set \ref flags */ inline void setFlags(const Flag& f) { flags = f.flags; - }; + } /** \brief * Sets \ref flags to \ref flags | f */ inline void setFlag(const unsigned long f) { flags |= f; - }; + } /** \brief * Sets \ref flags to \ref flags | f.flags */ inline void setFlag(const Flag& f) { flags |= f.flags; - }; + } /** \brief * Sets \ref flags to \ref flags & ~f */ inline void unsetFlag(const unsigned long f) { flags &= ~f; - }; + } /** \brief * Sets \ref flags to \ref flags & ~f.flags */ inline void unsetFlag(const Flag& f) { flags &= ~f.flags; - }; + } inline const unsigned long getFlags() const { return flags; - }; + } /** \brief * Returns \ref flags | f.flags @@ -141,7 +141,7 @@ namespace AMDiS { Flag r(flags); r.setFlag(f); return r; - }; + } /** \brief * Returns \ref flags & ~f.flags @@ -150,7 +150,7 @@ namespace AMDiS { Flag r(flags); r.unsetFlag(f); return r; - }; + } /** \brief * Returns \ref flags | f.flags @@ -159,7 +159,7 @@ namespace AMDiS { Flag r(flags); r.setFlag(f); return r; - }; + } /** \brief * Returns \ref flags & f.flags @@ -168,7 +168,7 @@ namespace AMDiS { Flag r(flags); r.flags &= f.flags; return r; - }; + } /** \brief * Sets \ref flags to \ref flags &= f.flags @@ -176,7 +176,7 @@ namespace AMDiS { inline Flag operator&=(const Flag& f) { flags &= f.flags; return *this; - }; + } /** \brief * Returns \ref flags ^ f.flags @@ -185,7 +185,7 @@ namespace AMDiS { Flag r(flags); r.flags ^=f.flags; return r; - }; + } /** \brief * Sets \ref flags to \ref flags & f.flags @@ -195,7 +195,7 @@ namespace AMDiS { flags |= f.flags; }; return *this; - }; + } /** \brief * Returns ~\ref flags @@ -204,28 +204,28 @@ namespace AMDiS { Flag r; r.flags = ~flags; return r; - }; + } /** \brief * Checks whether all set bits of f.flags are set in \ref flags too. */ inline bool isSet(const Flag& f) const { return ((flags&f.flags) == f.flags); - }; + } /** \brief * Returns !\ref isSet(f) */ inline bool isUnset(const Flag& f) const { return !isSet(f); - }; + } /** \brief * Returns true if \ref flags != 0 */ inline bool isAnySet() const { return (flags != 0); - }; + } protected: /** \brief diff --git a/AMDiS/src/GMResSolver.h b/AMDiS/src/GMResSolver.h index 6cbf0e185f238c5de345d8813e1acec84c53279f..9a03ed3487e5df831bc4d06a3623c0f0c6c51295 100644 --- a/AMDiS/src/GMResSolver.h +++ b/AMDiS/src/GMResSolver.h @@ -86,7 +86,8 @@ namespace AMDiS { /** \brief * realisation of OEMSolver::solveSystem */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/GMResSolver.hh b/AMDiS/src/GMResSolver.hh index e35fc959c3d0615ef9c7ab3ec89ed117fda30ee9..cf442a3eabaf86782220cc961eb12ae3803d8471 100644 --- a/AMDiS/src/GMResSolver.hh +++ b/AMDiS/src/GMResSolver.hh @@ -157,8 +157,7 @@ namespace AMDiS { template<typename VectorType> int GMResSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, - VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("GMResSolver::solveSystem()"); diff --git a/AMDiS/src/GMResSolver2.h b/AMDiS/src/GMResSolver2.h index ace8097247e998c0c18a2e35bfaa82ff8f41089c..6e820643311c510918b23db2edeacced1d0ea7fc 100644 --- a/AMDiS/src/GMResSolver2.h +++ b/AMDiS/src/GMResSolver2.h @@ -85,7 +85,8 @@ namespace AMDiS { /** \brief * realisation of OEMSolver::solveSystem */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/GMResSolver2.hh b/AMDiS/src/GMResSolver2.hh index 4ffee64efcf12e16e8ac16736d8e49521b30d6da..63fb90c24c6e0eb971682398e311d77f2c73c4f1 100644 --- a/AMDiS/src/GMResSolver2.hh +++ b/AMDiS/src/GMResSolver2.hh @@ -157,8 +157,8 @@ namespace AMDiS { template<typename VectorType> int GMResSolver2<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, - VectorType *b) + VectorType *x, VectorType *b, + bool reuseMatrix) { FUNCNAME("GMResSolver2::solveSystem()"); diff --git a/AMDiS/src/MultiGridWrapper.h b/AMDiS/src/MultiGridWrapper.h index 841d52ac1b26790b769452816f85cd69fdf9a231..a0f6dcf158d22b38c8d515238f00757a81d622ec 100644 --- a/AMDiS/src/MultiGridWrapper.h +++ b/AMDiS/src/MultiGridWrapper.h @@ -73,11 +73,12 @@ namespace AMDiS { int solveSystem(MatVecMultiplier<DOFVector<double> > *mv, DOFVector<double> *x, - DOFVector<double> *b) + DOFVector<double> *b, + bool reuseMatrix) { DOFMatrix *matrix = dynamic_cast<StandardMatVec<DOFMatrix, DOFVector<double> >*>(mv)->getMatrix(); - return mgSolver_->solve(*matrix, *x, *b); + return mgSolver_->solve(*matrix, *x, *b, reuseMatrix); }; protected: @@ -128,11 +129,12 @@ namespace AMDiS { int solveSystem(MatVecMultiplier<SystemVector> *mv, SystemVector *x, - SystemVector *b) + SystemVector *b, + bool reuseMatrix) { Matrix<DOFMatrix*> *matrix = dynamic_cast<StandardMatVec<Matrix<DOFMatrix*>, SystemVector>*>(mv)->getMatrix(); - return mgSolver_->solve(*matrix, *x, *b); + return mgSolver_->solve(*matrix, *x, *b, reuseMatrix); }; protected: diff --git a/AMDiS/src/ODirSolver.h b/AMDiS/src/ODirSolver.h index 9fab35ad277e96a89661f399f9c3fb22a920d24f..15ba9c97cb160cee64c3af2eac5c5c07bde9ecf0 100644 --- a/AMDiS/src/ODirSolver.h +++ b/AMDiS/src/ODirSolver.h @@ -76,7 +76,8 @@ namespace AMDiS { /** \brief * realisation of OEMSolver::solveSystem */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/ODirSolver.hh b/AMDiS/src/ODirSolver.hh index a7142a31e41a57202d1035e9908a8af473e601f5..11fcdaa2e7e6a89e42d59c04d3dd27bdb20c453b 100644 --- a/AMDiS/src/ODirSolver.hh +++ b/AMDiS/src/ODirSolver.hh @@ -35,7 +35,7 @@ namespace AMDiS { template<typename VectorType> int ODirSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("ODirSolver::solveSystem"); diff --git a/AMDiS/src/OEMSolver.h b/AMDiS/src/OEMSolver.h index f16053b7a42af7ddce05cbb816210c545b58ccce..19a158fd86e8b006233c6f8c9c08526e937c5ac8 100644 --- a/AMDiS/src/OEMSolver.h +++ b/AMDiS/src/OEMSolver.h @@ -96,7 +96,8 @@ namespace AMDiS { int solve(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b, Preconditioner<VectorType> *lPrecon = NULL, - Preconditioner<VectorType> *rPrecon = NULL) + Preconditioner<VectorType> *rPrecon = NULL, + bool reuseMatrix = false) { FUNCNAME("OEMSolver::solve()"); TEST_EXIT(mv)("no matVec\n"); @@ -113,7 +114,7 @@ namespace AMDiS { if (rightPrecon) rightPrecon->init(); - int iter = solveSystem(mv, x, b); + int iter = solveSystem(mv, x, b, reuseMatrix); if (leftPrecon) leftPrecon->exit(); @@ -124,7 +125,7 @@ namespace AMDiS { exit(); return iter; - }; + } // ===== getting-methods ====================================================== @@ -137,28 +138,28 @@ namespace AMDiS { */ inline const std::string& getName() { return name; - }; + } /** \brief * Returns \ref tolerance */ inline double getTolerance() { return tolerance; - }; + } /** \brief * Returns \ref max_iter */ inline int getMaxIterations() { return max_iter; - }; + } /** \brief * Returns \ref residual */ inline double getResidual() { return residual; - }; + } /** \} */ @@ -173,37 +174,36 @@ namespace AMDiS { */ inline void setTolerance(double tol) { tolerance = tol; - }; + } /** \brief * Sets \ref relative */ inline void setRelative(bool rel) { relative = rel; - }; + } /** \brief * Sets \ref max_iter */ inline void setMaxIterations(int i) { max_iter = i; - }; + } /** \brief * Sets \ref info */ inline void setInfo(int i) { info = i; - }; + } inline void setVectorCreator(CreatorInterface<VectorType> *creator) { vectorCreator = creator; - }; + } - inline CreatorInterface<VectorType>* getVectorCreator() - { + inline CreatorInterface<VectorType>* getVectorCreator() { return vectorCreator; - }; + } /** \} */ @@ -222,8 +222,9 @@ namespace AMDiS { * To be overloaded by concrete solvers. */ virtual int solveSystem(MatVecMultiplier<VectorType> *mv, - VectorType *x, - VectorType *b) = 0; + VectorType *x, + VectorType *b, + bool reuseMatrix) = 0; /** \brief @@ -294,14 +295,14 @@ namespace AMDiS { class OEMSolverCreator : public CreatorInterface<OEMSolver<T> > { public: - virtual ~OEMSolverCreator() {}; + virtual ~OEMSolverCreator() {} /** \brief * Sets \ref problem */ void setName(std::string name_) { name = name_; - }; + } protected: /** \brief diff --git a/AMDiS/src/OResSolver.h b/AMDiS/src/OResSolver.h index ddeadce05f2bfd675bcb4d1e07b01eadc7faeefb..8d731b37e688426db31c1c2a5a69be6d72afd998 100644 --- a/AMDiS/src/OResSolver.h +++ b/AMDiS/src/OResSolver.h @@ -76,7 +76,8 @@ namespace AMDiS { /** \brief * realisation of OEMSolver::solveSystem */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/OResSolver.hh b/AMDiS/src/OResSolver.hh index e3fca6d6cba6b41fa35a7f1f22282642c598c65e..74e95865f31552ba193cfa8910d8c30dd0e1df80 100644 --- a/AMDiS/src/OResSolver.hh +++ b/AMDiS/src/OResSolver.hh @@ -39,7 +39,8 @@ namespace AMDiS { template<typename VectorType> int OResSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, + bool reuseMatrix) { FUNCNAME("OResSolver::solveSystem"); diff --git a/AMDiS/src/PardisoSolver.cc b/AMDiS/src/PardisoSolver.cc index f4ee3a8ffd916cb7b5192299cdb40eb993c15129..bf395892465375f5992196acba6b0d0644c576cb 100644 --- a/AMDiS/src/PardisoSolver.cc +++ b/AMDiS/src/PardisoSolver.cc @@ -11,7 +11,8 @@ namespace AMDiS { template<> int PardisoSolver<SystemVector>::solveSystem(MatVecMultiplier<SystemVector> *matVec, - SystemVector *x, SystemVector *b) + SystemVector *x, SystemVector *b, + bool reuseMatrix) { FUNCNAME("PardisoSolver::solveSystem()"); diff --git a/AMDiS/src/PardisoSolver.h b/AMDiS/src/PardisoSolver.h index 11af655dbc67bd97fc4cb6cbd9fc51b58c5d55d6..e41415461d4c3a89e13207082767dcc113a69307 100644 --- a/AMDiS/src/PardisoSolver.h +++ b/AMDiS/src/PardisoSolver.h @@ -102,7 +102,8 @@ namespace AMDiS { /** \brief * Implements OEMSolver<VectorType>::solve(). */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); private: /** \brief diff --git a/AMDiS/src/ProblemInstat.h b/AMDiS/src/ProblemInstat.h index e3bd1d1a3fd6dc1f702b240a8fa1ee984189ab9f..e1621aaa2f27b3de8335b6585e5e37bf8d67575e 100644 --- a/AMDiS/src/ProblemInstat.h +++ b/AMDiS/src/ProblemInstat.h @@ -51,7 +51,7 @@ namespace AMDiS { ProblemStatBase *initialProb) : name(name_), initialProblem(initialProb ? initialProb : this) - {}; + {} /** \brief * Destructor. @@ -68,53 +68,53 @@ namespace AMDiS { /** \brief */ - virtual void solve(AdaptInfo *adaptInfo) {}; + virtual void solve(AdaptInfo *adaptInfo, bool fixedMatrix) {} /** \brief */ - virtual void estimate(AdaptInfo *adaptInfo) {}; + virtual void estimate(AdaptInfo *adaptInfo) {} /** \brief */ - virtual void buildBeforeRefine(AdaptInfo *adaptInfo, Flag) {}; + virtual void buildBeforeRefine(AdaptInfo *adaptInfo, Flag) {} /** \brief */ - virtual void buildBeforeCoarsen(AdaptInfo *adaptInfo, Flag) {}; + virtual void buildBeforeCoarsen(AdaptInfo *adaptInfo, Flag) {} /** \brief */ - virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag) {}; + virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag, bool, bool) {} /** \brief */ virtual Flag markElements(AdaptInfo *adaptInfo) { return 0; - }; + } /** \brief */ virtual Flag refineMesh(AdaptInfo *adaptInfo) { return 0; - }; + } /** \brief */ virtual Flag coarsenMesh(AdaptInfo *adaptInfo) { return 0; - }; + } /** \brief * Implementation of ProblemTimeInterface::closeTimestep(). */ - virtual void closeTimestep() {}; + virtual void closeTimestep() {} /** \brief * Returns \ref name. */ inline const std::string& getName() { return name; - }; + } /** \brief * Used by \ref problemInitial diff --git a/AMDiS/src/ProblemIterationInterface.h b/AMDiS/src/ProblemIterationInterface.h index d2f6e7e316cdaccb0404373a692ca89dda10b54a..742f9be822480e6a385396c30f6fbe58f7e543f6 100644 --- a/AMDiS/src/ProblemIterationInterface.h +++ b/AMDiS/src/ProblemIterationInterface.h @@ -29,11 +29,13 @@ namespace AMDiS { class AdaptInfo; class ProblemStatBase; - const Flag BUILD = 0x01L; - const Flag ADAPT = 0x02L; - const Flag SOLVE = 0x04L; - const Flag ESTIMATE = 0x08L; - const Flag MARK = 0x10L; + const Flag BUILD = 1; // Assemble vectors and matrices + const Flag BUILD_RHS = 2; // Assemble rhs vectors only + const Flag ADAPT = 4; // Run adaption procedure + const Flag SOLVE = 8; // Solve system + const Flag SOLVE_RHS = 16; // Solve system, where only rhs vectors have changed + const Flag ESTIMATE = 32; // Estimate error + const Flag MARK = 64; // Mark elements const Flag FULL_ITERATION = BUILD | ADAPT | SOLVE | ESTIMATE | MARK; const Flag NO_ADAPTION = BUILD | SOLVE | ESTIMATE; diff --git a/AMDiS/src/ProblemScal.cc b/AMDiS/src/ProblemScal.cc index cf959792931fb0dc4199cdc8086781244254081d..ede6b31c47c364a52b792de4342bd4e1fb9c387c 100644 --- a/AMDiS/src/ProblemScal.cc +++ b/AMDiS/src/ProblemScal.cc @@ -192,7 +192,7 @@ namespace AMDiS { } } - void ProblemScal::solve(AdaptInfo *adaptInfo) + void ProblemScal::solve(AdaptInfo *adaptInfo, bool fixedMatrix) { FUNCNAME("Problem::solve()"); if (!solver_) { @@ -542,7 +542,8 @@ namespace AMDiS { } } - void ProblemScal::buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag) + void ProblemScal::buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag, + bool assembleMatrix, bool assembleVector) { FUNCNAME("ProblemScal::buildAfterCoarsen()"); diff --git a/AMDiS/src/ProblemScal.h b/AMDiS/src/ProblemScal.h index d60bf76e9c0ea106b91940e4943ee13f6edc4e62..fb6e52cc13b1aec7fb8f6976d66b8d925fa28235 100644 --- a/AMDiS/src/ProblemScal.h +++ b/AMDiS/src/ProblemScal.h @@ -130,7 +130,7 @@ namespace AMDiS { * Implementation of ProblemStatBase::solve(). Deligates the solving * of problems system to \ref solver. */ - virtual void solve(AdaptInfo *adaptInfo); + virtual void solve(AdaptInfo *adaptInfo, bool fixedMatrix); /** \brief * Implementation of ProblemStatBase::estimate(). Deligates the estimation @@ -172,21 +172,23 @@ namespace AMDiS { * Implementation of ProblemStatBase::buildAfterCoarsen(). * Assembles \ref A and \ref rhs. */ - virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag); + virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag, + bool assembleMatrix = true, + bool assembleVector = true); /** \brief * Returns number of managed problems */ virtual int getNumProblems() { return 1; - }; + } /** \brief * Implementation of ProblemStatBase::getNumComponents() */ virtual int getNumComponents() { return 1; - }; + } /** \brief * Returns the problem with the given number. If only one problem @@ -194,7 +196,7 @@ namespace AMDiS { */ virtual ProblemStatBase *getProblem(int number = 0) { return this; - }; + } /** \brief * Writes output files. diff --git a/AMDiS/src/ProblemStatBase.h b/AMDiS/src/ProblemStatBase.h index 264de49fe4f8e392f48f7f77da64c04f948061a9..12d9633c47b79155f80b08a69daba84e956089e5 100644 --- a/AMDiS/src/ProblemStatBase.h +++ b/AMDiS/src/ProblemStatBase.h @@ -102,8 +102,11 @@ namespace AMDiS { /** \brief * Assembling of system matrices and vectors after coarsening. + * By the last two parameters, assembling can be restricted to either + * matrices or vectors only. */ - virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag) = 0; + virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag, + bool assembleMatrix, bool assembleVector) = 0; /** \brief * Refinement of the mesh. @@ -118,7 +121,7 @@ namespace AMDiS { /** \brief * Solves the assembled system. The result is an approximative solution. */ - virtual void solve(AdaptInfo *adaptInfo) = 0; + virtual void solve(AdaptInfo *adaptInfo, bool fixedMatrix) = 0; /** \brief * A posteriori error estimation of the calculated solution. Should store diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc index ea6c7e702557885041d1063dfe6fad7b1b189f5a..89f7b4da61f213014449ecd5a2cab817bdc7cc32 100644 --- a/AMDiS/src/ProblemVec.cc +++ b/AMDiS/src/ProblemVec.cc @@ -506,7 +506,7 @@ namespace AMDiS { { } - void ProblemVec::solve(AdaptInfo *adaptInfo) + void ProblemVec::solve(AdaptInfo *adaptInfo, bool fixedMatrix) { FUNCNAME("Problem::solve()"); @@ -520,7 +520,8 @@ namespace AMDiS { #endif clock_t first = clock(); - int iter = solver_->solve(matVec_, solution_, rhs_, leftPrecon_, rightPrecon_); + int iter = solver_->solve(matVec_, solution_, rhs_, + leftPrecon_, rightPrecon_, fixedMatrix); #ifdef _OPENMP INFO(info_, 8)("solution of discrete system needed %.5f seconds system time / %.5f seconds wallclock time\n", @@ -648,7 +649,8 @@ namespace AMDiS { return StandardProblemIteration::oneIteration(adaptInfo, toDo); } - void ProblemVec::buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag) + void ProblemVec::buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag, + bool asmMatrix, bool asmVector) { FUNCNAME("ProblemVec::buildAfterCoarsen()"); @@ -687,7 +689,7 @@ namespace AMDiS { if ((*systemMatrix_)[i][j]) { // The matrix should not be deleted, if it was assembled before // and it is marked to be assembled only once. - if (!(assembleMatrixOnlyOnce_[i][j] && assembledMatrix_[i][j])) { + if (!(assembleMatrixOnlyOnce_[i][j] && assembledMatrix_[i][j]) && asmMatrix) { (*systemMatrix_)[i][j]->clear(); } } @@ -712,6 +714,10 @@ namespace AMDiS { assembleMatrix = false; } + if (!asmMatrix) { + assembleMatrix = false; + } + // If the matrix should not be assembled, the rhs vector has to be considered. // This will be only done, if i == j. So, if both is not true, we can jump // to the next matrix. @@ -727,7 +733,7 @@ namespace AMDiS { assembleOnOneMesh(componentSpaces[i], assembleFlag, assembleMatrix ? matrix : NULL, - (i == j) ? rhs_->getDOFVector(i) : NULL); + ((i == j) && asmVector) ? rhs_->getDOFVector(i) : NULL); // std::cout << "--- Finished ---\n"; // std::cout << "--- Dim: " << matrix->getUsedSize() << "x" << matrix->getNumCols() << "---\n"; } else { @@ -735,7 +741,7 @@ namespace AMDiS { assembleOnDifMeshes(componentSpaces[i], componentSpaces[j], assembleFlag, assembleMatrix ? matrix : NULL, - (i == j) ? rhs_->getDOFVector(i) : NULL); + ((i == j) && asmVector) ? rhs_->getDOFVector(i) : NULL); // std::cout << "--- Finished ---\n"; // std::cout << "--- Dim: " << matrix->getUsedSize() << "x" << matrix->getNumCols() << "---\n"; } diff --git a/AMDiS/src/ProblemVec.h b/AMDiS/src/ProblemVec.h index e3c6b6500d17b8dd9a80fd2285991f38b20abb6d..e840b33dae3458869526b59970d480c01c6dbf35 100644 --- a/AMDiS/src/ProblemVec.h +++ b/AMDiS/src/ProblemVec.h @@ -151,7 +151,7 @@ namespace AMDiS { * Implementation of ProblemStatBase::solve(). Deligates the solving * of problems system to \ref solver. */ - virtual void solve(AdaptInfo *adaptInfo); + virtual void solve(AdaptInfo *adaptInfo, bool fixedMatrix = false); /** \brief * Implementation of ProblemStatBase::estimate(). Deligates the estimation @@ -191,9 +191,12 @@ namespace AMDiS { /** \brief * Implementation of ProblemStatBase::buildAfterCoarsen(). - * Assembles \ref A and \ref rhs. + * Assembles \ref A and \ref rhs. With the last two parameters, assembling + * can be restricted to matrices or vectors only. */ - virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag); + virtual void buildAfterCoarsen(AdaptInfo *adaptInfo, Flag flag, + bool assembleMatrix = true, + bool assembleVector = true); /** \brief * Determines the execution order of the single adaption steps. If adapt is diff --git a/AMDiS/src/SolutionDataStorage.h b/AMDiS/src/SolutionDataStorage.h index 15319ac1a13c54f635284ea8d809cc43d3a808fb..4cd0b0ebaea72e37aec5516c04a2bc83a8166433 100644 --- a/AMDiS/src/SolutionDataStorage.h +++ b/AMDiS/src/SolutionDataStorage.h @@ -82,8 +82,9 @@ namespace AMDiS { } void deleteFeSpace(std::vector<FiniteElemSpace*> feSpaces) { - for (int i = 0; i < static_cast<int>(feSpaces.size()); i++) + for (int i = 0; i < static_cast<int>(feSpaces.size()); i++) { DELETE feSpaces[i]; + } feSpaces.empty(); } diff --git a/AMDiS/src/SolutionDataStorage.hh b/AMDiS/src/SolutionDataStorage.hh index 62e810f1481433b37f4682dd24c08b662a70e9a7..1bf73e5cb37fc41a7c7082f0c71566ea9df1b4df 100644 --- a/AMDiS/src/SolutionDataStorage.hh +++ b/AMDiS/src/SolutionDataStorage.hh @@ -61,9 +61,13 @@ namespace AMDiS { double timestamp, typename SolutionHelper<T>::type feSpace) { + FUNCNAME("SolutionDataStorage<T>::push()"); + push(solution, timestamp); // Store fe space only, if we do not have a fixed fe space. + TEST_EXIT(!fixedFESpace)("push wit fe space not possible!\n"); + if (!fixedFESpace) { feSpaces.push_back(feSpace); } @@ -111,8 +115,10 @@ namespace AMDiS { DELETE solutions[i]; } - for (int i = 0; i < static_cast<int>(feSpaces.size()); i++) { - deleteFeSpace(feSpaces[i]); + if (!fixedFESpace) { + for (int i = 0; i < static_cast<int>(feSpaces.size()); i++) { + deleteFeSpace(feSpaces[i]); + } } solutions.empty(); diff --git a/AMDiS/src/StandardProblemIteration.cc b/AMDiS/src/StandardProblemIteration.cc index 706b4beabdc4801f65a8e95531ecaed53c7b4d6d..ac9742c9a6ff1e9ae4106a99141910655dd27786 100644 --- a/AMDiS/src/StandardProblemIteration.cc +++ b/AMDiS/src/StandardProblemIteration.cc @@ -24,7 +24,10 @@ namespace AMDiS { Flag flag = buildAndAdapt(adaptInfo, toDo); if (toDo.isSet(SOLVE)) - problem_->solve(adaptInfo); + problem_->solve(adaptInfo, false); + + if (toDo.isSet(SOLVE_RHS)) + problem_->solve(adaptInfo, true); if (toDo.isSet(ESTIMATE)) problem_->estimate(adaptInfo); @@ -72,7 +75,11 @@ namespace AMDiS { } if (toDo.isSet(BUILD)) { - problem_->buildAfterCoarsen(adaptInfo, markFlag); + problem_->buildAfterCoarsen(adaptInfo, markFlag, true, true); + } + + if (toDo.isSet(BUILD_RHS)) { + problem_->buildAfterCoarsen(adaptInfo, markFlag, false, true); } return flag; diff --git a/AMDiS/src/StandardProblemIteration.h b/AMDiS/src/StandardProblemIteration.h index bb9e8d58c2be4e5641a8eab3697f34b3bff45f36..b1b6faf9b481e6e8ee34f4882c028dcd6cb27be1 100644 --- a/AMDiS/src/StandardProblemIteration.h +++ b/AMDiS/src/StandardProblemIteration.h @@ -71,14 +71,14 @@ namespace AMDiS { */ int getNumProblems() { return 1; - }; + } /** \brief * implementation of \ref ProblemIterationInterface::getProblem(int) */ ProblemStatBase *getProblem(int number = 0) { return problem_; - }; + } /** \brief * Returns the name of the problem. diff --git a/AMDiS/src/TFQMR.h b/AMDiS/src/TFQMR.h index 67087e1d874cb8f1dc09a6e71bba6c93ce09f437..e127e26719132cea19152f3ed2ae70956dc72498 100644 --- a/AMDiS/src/TFQMR.h +++ b/AMDiS/src/TFQMR.h @@ -72,7 +72,8 @@ namespace AMDiS { /** \brief * realisation of OEMSolver::solveSystem */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); /** \brief * realisation of OEMSolver::init diff --git a/AMDiS/src/TFQMR.hh b/AMDiS/src/TFQMR.hh index a3773e7514d3909198b803f37a8e4b057237b8fa..d2ae5e926868added5c6083471eee5b2b600e230 100644 --- a/AMDiS/src/TFQMR.hh +++ b/AMDiS/src/TFQMR.hh @@ -43,8 +43,9 @@ namespace AMDiS { template<typename VectorType> int TFQMR<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, - VectorType *b) + VectorType *x, + VectorType *b, + bool reuseMatrix) { DOFMatrix *m = (dynamic_cast<StandardMatVec<DOFMatrix, VectorType> *>(matVec))->getMatrix(); diff --git a/AMDiS/src/UmfPackSolver.h b/AMDiS/src/UmfPackSolver.h index fcbc07c274bd057c78657d40e332c49c659be8f3..6bc57c890e569c8588579f27c3a3043dc1902795 100644 --- a/AMDiS/src/UmfPackSolver.h +++ b/AMDiS/src/UmfPackSolver.h @@ -82,7 +82,7 @@ namespace AMDiS { void init() { p = this->vectorCreator->create(); r = this->vectorCreator->create(); - }; + } /** \brief * Implements OEMSolver<VectorType>::exit(). @@ -90,12 +90,12 @@ namespace AMDiS { void exit() { this->vectorCreator->free(p); this->vectorCreator->free(r); - }; + } /** \brief * Implements OEMSolver<VectorType>::solve(). */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b, bool reuseMatrix); private: /** \brief @@ -109,6 +109,11 @@ namespace AMDiS { */ void *symbolic_; + /** \brief + * Stores the result of umfpack_di_numeric, if multiple right hand sides will occure. + */ + void *numeric; + /** \brief * If the symbolic analysis should be done only once (for example, if the matrices * to solve have all the same pattern because of no adaptivity), this variable @@ -116,6 +121,12 @@ namespace AMDiS { * of solveSystem(). */ int storeSymbolic_; + + /** \brief + * If not zero, Umfpack is prepared to solve multiple following systems with equal + * system matrix but different right hand side vectors. + */ + int multipleRhs; }; } diff --git a/AMDiS/src/UmfPackSolver.hh b/AMDiS/src/UmfPackSolver.hh index 42e47c8b651a5560968748b9661496a7287a7f02..2d7655db296b76f8713f3414623d446bbf0c196c 100644 --- a/AMDiS/src/UmfPackSolver.hh +++ b/AMDiS/src/UmfPackSolver.hh @@ -13,9 +13,11 @@ namespace AMDiS { UmfPackSolver<VectorType>::UmfPackSolver(std::string name) : OEMSolver<VectorType>(name), symbolic_(NULL), - storeSymbolic_(0) + storeSymbolic_(0), + multipleRhs(0) { GET_PARAMETER(0, name + "->store symbolic", "%d", &storeSymbolic_); + GET_PARAMETER(0, name + "->multiple rhs", "%d", &multipleRhs); } template<typename VectorType> @@ -28,7 +30,8 @@ namespace AMDiS { template<typename VectorType> int UmfPackSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, + bool reuseMatrix) { FUNCNAME("UmfPackSolver::solveSystem()"); @@ -114,7 +117,6 @@ namespace AMDiS { } - void *numeric; double Control[UMFPACK_CONTROL]; double Info[UMFPACK_INFO]; int status = 0; @@ -128,15 +130,24 @@ namespace AMDiS { if (!(storeSymbolic_) || (storeSymbolic_ && (symbolic_ == NULL))) { - status = umfpack_di_symbolic(newMatrixSize, newMatrixSize, Ap, Ai, Ax, &symbolic_, Control, Info); + status = umfpack_di_symbolic(newMatrixSize, newMatrixSize, Ap, Ai, Ax, + &symbolic_, Control, Info); if (!status == UMFPACK_OK) { ERROR_EXIT("UMFPACK Error in function umfpack_di_symbolic"); } } - - status = umfpack_di_numeric(Ap, Ai, Ax, symbolic_, &numeric, Control, Info); - if (!status == UMFPACK_OK) { - ERROR_EXIT("UMFPACK Error in function umfpack_di_numeric"); + + if (!((multipleRhs != 0) && (reuseMatrix))) { + // With multiple rhs are allowed, but not this time, free the matrix + // factorization of last solution step, + if ((multipleRhs != 0) && !reuseMatrix) { + umfpack_di_free_numeric(&numeric); + } + + status = umfpack_di_numeric(Ap, Ai, Ax, symbolic_, &numeric, Control, Info); + if (!status == UMFPACK_OK) { + ERROR_EXIT("UMFPACK Error in function umfpack_di_numeric"); + } } status = umfpack_di_solve(UMFPACK_A, Ap, Ai, Ax, xvec, bvec, numeric, Control, Info); @@ -147,7 +158,10 @@ namespace AMDiS { if (!storeSymbolic_) { umfpack_di_free_symbolic(&symbolic_); } - umfpack_di_free_numeric(&numeric); + if (multipleRhs == 0) { + std::cout << "FREE THE STUFF" << std::endl; + umfpack_di_free_numeric(&numeric); + } // Copy and resort solution. diff --git a/AMDiS/src/VecSymSolver.h b/AMDiS/src/VecSymSolver.h index 632719ef0da29754032b72e02a95e0b7d1a94a48..130d50b3a5dfffc67c2fb8aba948eb0986119dfe 100644 --- a/AMDiS/src/VecSymSolver.h +++ b/AMDiS/src/VecSymSolver.h @@ -89,7 +89,8 @@ namespace AMDiS { /** \brief * Implements OEMSolver<VectorType>::solve(). */ - int solveSystem(MatVecMultiplier<VectorType> *mv, VectorType *x, VectorType *b); + int solveSystem(MatVecMultiplier<VectorType> *mv, + VectorType *x, VectorType *b, bool reuseMatrix); private: VectorType *r, *p; diff --git a/AMDiS/src/VecSymSolver.hh b/AMDiS/src/VecSymSolver.hh index 3e5f522d287469f8ee3ac0e778a19540dd7d7b56..76a1cd6fbaaf08559a341a83c2af19f25cad2438 100644 --- a/AMDiS/src/VecSymSolver.hh +++ b/AMDiS/src/VecSymSolver.hh @@ -227,7 +227,7 @@ namespace AMDiS { template<typename VectorType> int VecSymSolver<VectorType>::solveSystem(MatVecMultiplier<VectorType> *matVec, - VectorType *x, VectorType *b) + VectorType *x, VectorType *b, bool reuseMatrix) { FUNCNAME("VecSymSolver::solveSystem()");