diff --git a/AMDiS/libtool b/AMDiS/libtool index 1d1978368d0efde93663faf1df97bfe1eea1ae55..b70a25ce066c61dee0221d0f2befd34eb4875f21 100755 --- a/AMDiS/libtool +++ b/AMDiS/libtool @@ -44,7 +44,7 @@ available_tags=" CXX F77" # ### BEGIN LIBTOOL CONFIG -# Libtool was configured on host deimos101: +# Libtool was configured on host p2d108: # Shell to use when invoking shell scripts. SHELL="/bin/sh" @@ -82,13 +82,13 @@ AR="ar" AR_FLAGS="cru" # A C compiler. -LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc" +LTCC="gcc" # LTCC compiler flags. LTCFLAGS="-g -O2" # A language-specific compiler. -CC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc" +CC="gcc" # Is the compiler the GNU C compiler? with_gcc=yes @@ -174,7 +174,7 @@ dlopen_self=unknown dlopen_self_static=unknown # Compiler flag to prevent dynamic linking. -link_static_flag="" +link_static_flag="-static" # Compiler flag to turn off builtin functions. no_builtin_flag=" -fno-builtin" @@ -6763,7 +6763,7 @@ build_old_libs=`case $build_libtool_libs in yes) $echo no;; *) $echo yes;; esac` # End: # ### BEGIN LIBTOOL TAG CONFIG: CXX -# Libtool was configured on host deimos101: +# Libtool was configured on host p2d108: # Shell to use when invoking shell scripts. SHELL="/bin/sh" @@ -6801,13 +6801,13 @@ AR="ar" AR_FLAGS="cru" # A C compiler. -LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc" +LTCC="gcc" # LTCC compiler flags. LTCFLAGS="-g -O2" # A language-specific compiler. -CC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpiCC" +CC="g++" # Is the compiler the GNU C compiler? with_gcc=yes @@ -6893,7 +6893,7 @@ dlopen_self=unknown dlopen_self_static=unknown # Compiler flag to prevent dynamic linking. -link_static_flag="" +link_static_flag="-static" # Compiler flag to turn off builtin functions. no_builtin_flag=" -fno-builtin" @@ -6960,11 +6960,11 @@ predeps="" # Dependencies to place after the objects being linked to create a # shared library. -postdeps="-lmpi_cxx -lmpi -lopen-rte -lopen-pal -libverbs -lrt -lnuma -ldl -lnsl -lutil -ldl -lstdc++ -lm -lgcc_s -lpthread -lc -lgcc_s" +postdeps="-lstdc++ -lm -lgcc_s -lc -lgcc_s" # The library search path used internally by the compiler when linking # a shared library. -compiler_lib_search_path=`echo "-L/usr/lib64 -L/licsoft/libraries/openmpi/1.2.6/64bit/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2 -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64 -L/lib/../lib64 -L/usr/lib/../lib64 -L/fastfs/wir/local/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../.." | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"` +compiler_lib_search_path=`echo "-L/usr/lib64/gcc/x86_64-suse-linux/4.1.2 -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64 -L/lib/../lib64 -L/usr/lib/../lib64 -L/fastfs/wir/local/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../.." | $SED -e "s@${gcc_dir}@\${gcc_dir}@g;s@${gcc_ver}@\${gcc_ver}@g"` # Method to check whether dependent libraries are shared objects. deplibs_check_method="pass_all" @@ -7071,7 +7071,7 @@ include_expsyms="" # ### BEGIN LIBTOOL TAG CONFIG: F77 -# Libtool was configured on host deimos101: +# Libtool was configured on host p2d108: # Shell to use when invoking shell scripts. SHELL="/bin/sh" @@ -7109,7 +7109,7 @@ AR="ar" AR_FLAGS="cru" # A C compiler. -LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc" +LTCC="gcc" # LTCC compiler flags. LTCFLAGS="-g -O2" diff --git a/AMDiS/src/AdaptInfo.h b/AMDiS/src/AdaptInfo.h index d493bc1a7c67a30045f3404e5579ba4618e3703d..ab86c4581391a5368efcc881c58e4e578e5b02da 100644 --- a/AMDiS/src/AdaptInfo.h +++ b/AMDiS/src/AdaptInfo.h @@ -47,19 +47,18 @@ namespace AMDiS { /// Constructor. ScalContent(std::string prefix) : est_sum(0.0), - est_t_sum(0.0), - est_max(0.0), - est_t_max(0.0), - fac_max(0.0), - fac_sum(1.0), - spaceTolerance(0.0), - timeTolerance(0.0), - timeErrLow(0.0), - coarsenAllowed(0), - refinementAllowed(1), - refineBisections(1), - coarseBisections(1) - + est_t_sum(0.0), + est_max(0.0), + est_t_max(0.0), + fac_max(0.0), + fac_sum(1.0), + spaceTolerance(0.0), + timeTolerance(0.0), + timeErrLow(0.0), + coarsenAllowed(0), + refinementAllowed(1), + refineBisections(1), + coarseBisections(1) { double timeTheta2 = 0.3; @@ -147,7 +146,8 @@ namespace AMDiS { solverTolerance(1e-8), solverResidual(0.0), scalContents(size), - isDeserialized_(false) + deserialized(false), + rosenbrockMode(false) { GET_PARAMETER(0, name_ + "->start time", "%f", &startTime); time = startTime; @@ -615,13 +615,18 @@ namespace AMDiS { inline const int getCoarseBisections(int index) const { return scalContents[index]->coarseBisections; - } + } inline int getSize() { return scalContents.size(); } + inline bool getRosenbrockMode() + { + return rosenbrockMode; + } + inline void setSolverIterations(int it) { solverIterations = it; @@ -665,12 +670,17 @@ namespace AMDiS { /// Returns true, if the adaptive procedure was deserialized from a file. const bool isDeserialized() const { - return isDeserialized_; + return deserialized; } inline void setIsDeserialized(bool b) { - isDeserialized_ = b; + deserialized = b; + } + + inline void setRosenbrockMode(bool b) + { + rosenbrockMode = b; } /// Creates new scalContents with the given size. @@ -767,7 +777,10 @@ namespace AMDiS { std::vector<ScalContent*> scalContents; /// Is true, if the adaptive procedure was deserialized from a file. - bool isDeserialized_; + bool deserialized; + + /// Is true, if the time adaption is controlled by a Rosenbrock Method. + bool rosenbrockMode; }; } diff --git a/AMDiS/src/DOFVector.h b/AMDiS/src/DOFVector.h index f2771e079e2295eb8730221283c8926231dcedec..bbb2b0ad5094dc89ac4cdc5f6c3aa3bb72509b4d 100644 --- a/AMDiS/src/DOFVector.h +++ b/AMDiS/src/DOFVector.h @@ -523,19 +523,22 @@ namespace AMDiS { /// vec[i] = v.vec[i] void copy(const DOFVector<T>& v); - /// Returns minimum of DOFVector + /// Returns minimum of DOFVector. T min() const; - /// Returns maximum of DOFVector + /// Returns maximum of DOFVector. T max() const; - /// Returns absolute maximum of DOFVector + /// Returns absolute maximum of DOFVector. T absMax() const; - /// Used by interpol while mesh traversal + /// Returns the average value of the DOFVector. + T average() const; + + /// Used by interpol while mesh traversal. void interpolFct(ElInfo* elinfo); - /// Prints \ref vec to stdout + /// Prints \ref vec to stdout. void print() const; /// diff --git a/AMDiS/src/DOFVector.hh b/AMDiS/src/DOFVector.hh index 1a92055050a09aaccc2d62ca0224054bdf9fc696..b1412a09b2d278f32fcf8f88903316c8239eab4e 100644 --- a/AMDiS/src/DOFVector.hh +++ b/AMDiS/src/DOFVector.hh @@ -289,6 +289,7 @@ namespace AMDiS { return m; } + template<typename T> T DOFVector<T>::max() const { @@ -302,6 +303,7 @@ namespace AMDiS { m = std::max(m, *vecIterator); return m; + } template<typename T> @@ -310,6 +312,26 @@ namespace AMDiS { return std::max(abs(max()), abs(min())); } + + template<typename T> + T DOFVector<T>::average() const + { + FUNCNAME("DOFVector<T>::average()"); + + checkFeSpace(this->feSpace, vec); + + int count = 0; + T m = 0; + Iterator vecIterator(const_cast<DOFIndexed<T>*>(dynamic_cast<const DOFIndexed<T>*>(this)), USED_DOFS); + for (vecIterator.reset(); !vecIterator.end(); ++vecIterator) { + m += *vecIterator; + count++; + } + + return m / count; + } + + template<typename T> void DOFVector<T>::print() const { diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc index 52d210534852ba737396e7ad2b753bc3b226be32..6e56072cb3abecb4805d35acc51b35a7fb89ccd1 100644 --- a/AMDiS/src/ProblemVec.cc +++ b/AMDiS/src/ProblemVec.cc @@ -538,10 +538,11 @@ namespace AMDiS { adaptInfo->setEstSum(scalEstimator->getErrorSum(), i); adaptInfo->setEstMax(scalEstimator->getErrorMax(), i); - adaptInfo->setTimeEstSum(scalEstimator->getTimeEst(), i); - adaptInfo->setTimeEstMax(scalEstimator->getTimeEstMax(), i); - } else { - WARNING("No estimator for component %d\n" , i); + + if (adaptInfo->getRosenbrockMode() == false) { + adaptInfo->setTimeEstSum(scalEstimator->getTimeEst(), i); + adaptInfo->setTimeEstMax(scalEstimator->getTimeEstMax(), i); + } } } } @@ -623,7 +624,7 @@ namespace AMDiS { FUNCNAME("ProblemVec::buildAfterCoarsen()"); if (dualMeshTraverseRequired()) { - dualAssemble(adaptInfo, flag); + dualAssemble(adaptInfo, flag, asmMatrix, asmVector); return; } @@ -779,10 +780,13 @@ namespace AMDiS { } - void ProblemVec::dualAssemble(AdaptInfo *adaptInfo, Flag flag) + void ProblemVec::dualAssemble(AdaptInfo *adaptInfo, Flag flag, + bool asmMatrix, bool asmVector) { FUNCNAME("ProblemVec::dualAssemble()"); + TEST_EXIT(asmVector)("Not yet implemented!\n"); + for (unsigned int i = 0; i < meshes.size(); i++) meshes[i]->dofCompress(); @@ -828,7 +832,7 @@ namespace AMDiS { bool assembleMatrix = true; // The DOFMatrix which should be assembled (or not, if assembleMatrix // will be set to false). - DOFMatrix *matrix = (*systemMatrix)[i][j]; + DOFMatrix *matrix = (asmMatrix ? (*systemMatrix)[i][j] : NULL); if (writeAsmInfo && matrix) { for (std::vector<Operator*>::iterator it = matrix->getOperatorsBegin(); @@ -874,6 +878,10 @@ namespace AMDiS { if (!assembleMatrix && i != j) if (matrix) nnz += matrix->getBaseMatrix().nnz(); + + if (matrix && !assembleMatrix) { + ERROR_EXIT("Not yet implemented!\n"); + } } } @@ -899,7 +907,7 @@ namespace AMDiS { for (int i = 0; i < nComponents; i++) { for (int j = 0; j < nComponents; j++) { - DOFMatrix *matrix = (*systemMatrix)[i][j]; + DOFMatrix *matrix = (asmMatrix ? (*systemMatrix)[i][j] : NULL); if (!matrix) continue; @@ -983,19 +991,18 @@ namespace AMDiS { for (int i = 0; i < nComponents; i++) { for (int j = 0; j < nComponents; j++) { - DOFMatrix *matrix = (*systemMatrix)[i][j]; + DOFMatrix *matrix = (asmMatrix ? (*systemMatrix)[i][j] : NULL); - if (matrix) - matrix->removeRowsWithDBC(matrix->getApplyDBCs()); + if (!matrix) + continue; - if (matrix) - matrix->finishInsertion(); + matrix->removeRowsWithDBC(matrix->getApplyDBCs()); + matrix->finishInsertion(); - if (matrix && matrix->getBoundaryManager()) + if (matrix->getBoundaryManager()) matrix->getBoundaryManager()->exitMatrix(matrix); - if (matrix) - nnz += matrix->getBaseMatrix().nnz(); + nnz += matrix->getBaseMatrix().nnz(); } // And now assemble boundary conditions on the vectors diff --git a/AMDiS/src/ProblemVec.h b/AMDiS/src/ProblemVec.h index 800ba3ab2750e1677705dae97b8bbecfdf1e8cf2..862574a2f34dd769eb2d8fbfcff00a4a4793dee0 100644 --- a/AMDiS/src/ProblemVec.h +++ b/AMDiS/src/ProblemVec.h @@ -182,7 +182,8 @@ namespace AMDiS { bool dualMeshTraverseRequired(); - void dualAssemble(AdaptInfo *adaptInfo, Flag flag); + void dualAssemble(AdaptInfo *adaptInfo, Flag flag, + bool asmMatrix = true, bool asmVector = true); void createPrecon(); diff --git a/AMDiS/src/ResidualEstimator.cc b/AMDiS/src/ResidualEstimator.cc index 49aa4efed8b25423411857c74e51452d69da25c9..f8a741fd862eb719c458d8bd54b806cffd828bad 100644 --- a/AMDiS/src/ResidualEstimator.cc +++ b/AMDiS/src/ResidualEstimator.cc @@ -179,9 +179,9 @@ namespace AMDiS { } if (output) { - MSG("estimate = %.8e\n", est_sum); + MSG("estimate for component %d = %.8e\n", row, est_sum); if (C3) - MSG("time estimate = %.8e\n", est_t_sum); + MSG("time estimate for component %d = %.8e\n", row, est_t_sum); } delete [] riq; diff --git a/AMDiS/src/VtkWriter.hh b/AMDiS/src/VtkWriter.hh index ff0045fd1a351d73374e82e9e711b0e05c31c037..8a711b0524682dc0274915df65a1c6cb766ba41f 100644 --- a/AMDiS/src/VtkWriter.hh +++ b/AMDiS/src/VtkWriter.hh @@ -126,7 +126,8 @@ namespace AMDiS { DOFVector<double>::Iterator valueIt(values, USED_DOFS); DOFVector< std::list<WorldVector<double> > >::Iterator coordIt(dofCoords, USED_DOFS); - file << std::fixed; + file << std::scientific; + file.precision(15); // Write the values for all vertex DOFs. for (intPointIt.reset(), valueIt.reset(), coordIt.reset(); diff --git a/AMDiS/src/time/RosenbrockAdaptInstationary.cc b/AMDiS/src/time/RosenbrockAdaptInstationary.cc index ccf07696d6dcff3f1f43bafefbb83e5ef2d05cb3..523235ffbe3c9f154287a8ec3aa4730a7d6108a4 100644 --- a/AMDiS/src/time/RosenbrockAdaptInstationary.cc +++ b/AMDiS/src/time/RosenbrockAdaptInstationary.cc @@ -15,15 +15,23 @@ namespace AMDiS { firstTimestep(true), lastTimestepRejected(false), succRejection(false), + fixFirstTimesteps(0), tauGamma(1.0) { + FUNCNAME("RosenbrockAdaptInstationary::RosenbrockAdaptInstationary()"); + std::string str; GET_PARAMETER(0, name + "->rosenbrock method", &str); RosenbrockMethodCreator *creator = dynamic_cast<RosenbrockMethodCreator*>(CreatorMap<RosenbrockMethod>::getCreator(str)); rosenbrockMethod = creator->create(); + TEST_EXIT_DBG(rosenbrockMethod)("This should not happen!\n"); + + GET_PARAMETER(0, name + "->fix first timesteps", "%d", &fixFirstTimesteps); problemStat->setRosenbrockMethod(rosenbrockMethod); + + adaptInfo->setRosenbrockMode(true); } @@ -58,61 +66,63 @@ namespace AMDiS { double newTimestep = 0.0; double order = rosenbrockMethod->getOrder(); - if (errorEst < timeTol) { + if (errorEst < timeTol || fixFirstTimesteps > 0 || firstTimestep) { double fac = 1.0; - - if (firstTimestep || succRejection) { - fac = pow((timeTol / errorEst), 1.0 / order); + + if (fixFirstTimesteps > 0) { + newTimestep = adaptInfo->getTimestep(); } else { - fac = adaptInfo->getTimestep() / tauAcc * - pow((timeTol * estAcc / (errorEst * errorEst)), 1.0 / order); + if (firstTimestep || succRejection) { + fac = pow((timeTol / errorEst), 1.0 / order); + } else { + fac = adaptInfo->getTimestep() / tauAcc * + pow((timeTol * estAcc / (errorEst * errorEst)), 1.0 / order); + } + fac = std::min(fac, 3.0); + newTimestep = fac * adaptInfo->getTimestep(); + newTimestep *= 0.95; } - fac = std::min(fac, 3.0); - newTimestep = fac * adaptInfo->getTimestep(); tauAcc = adaptInfo->getTimestep(); estAcc = errorEst; - + lastTimestepRejected = false; succRejection = false; } else { if (lastTimestepRejected) { succRejection = true; - + double reducedP = log(errorEst / estRej) / log(adaptInfo->getTimestep() / tauRej); if (reducedP < order && reducedP > 0.0) order = reducedP; } - - newTimestep = pow((timeTol / errorEst), 1.0 / order) * adaptInfo->getTimestep(); + newTimestep = pow((timeTol / errorEst), 1.0 / order) * adaptInfo->getTimestep(); + newTimestep *= 0.95; + tauRej = adaptInfo->getTimestep(); estRej = errorEst; - + lastTimestepRejected = true; } - - newTimestep *= 0.95; + - if (errorEst < timeTol || firstTimestep) { - INFO(info, 6)("Accepted timestep at time = %f with timestep = %f\n", + if (errorEst < timeTol || fixFirstTimesteps) { + INFO(info, 6)("Accepted timestep at time = %e with timestep = %e\n", adaptInfo->getTime() - adaptInfo->getTimestep(), adaptInfo->getTimestep()); adaptInfo->setTimestep(newTimestep); - if (firstTimestep) - firstTimestep = false; - rejected = false; for (int i = 0; i < adaptInfo->getSize(); i++) - INFO(info, 6)("time estimator for component %d = %f\n", + INFO(info, 6)("time estimate for component %d = %e\n", i, adaptInfo->getTimeEstSum(i)); } else { - INFO(info, 6)("Rejected timestep at time = %f with timestep = %f\n", + INFO(info, 6)("Rejected timestep at time = %e with timestep = %e\n", adaptInfo->getTime() - adaptInfo->getTimestep(), adaptInfo->getTimestep()); @@ -122,6 +132,12 @@ namespace AMDiS { rejected = true; } + if (firstTimestep) + firstTimestep = false; + + if (fixFirstTimesteps > 0) + fixFirstTimesteps--; + } while (rejected); rosenbrockStat->acceptTimestep(); diff --git a/AMDiS/src/time/RosenbrockAdaptInstationary.h b/AMDiS/src/time/RosenbrockAdaptInstationary.h index dc7548ac2df4b4587913d95dbb3384829c91a372..b93a236455efe922756d486e7db6812f2baac7f6 100644 --- a/AMDiS/src/time/RosenbrockAdaptInstationary.h +++ b/AMDiS/src/time/RosenbrockAdaptInstationary.h @@ -40,21 +40,43 @@ namespace AMDiS { AdaptInfo *initialInfo, time_t initialTimestamp = 0); + /// Runs the Rosenbrock loop until one timestep is accepted. void oneTimestep(); + /// Returns \ref tauGamma. double* getTauGamma() { return &tauGamma; } protected: + /// Pointer to the Rosenbrock method that should be used. RosenbrockMethod *rosenbrockMethod; + /// Pointer to the stationary problem; RosenbrockStationary *rosenbrockStat; - bool firstTimestep, lastTimestepRejected, succRejection; + /// Indicates, if this is the very first timestep. + bool firstTimestep; - double tauAcc, estAcc, tauRej, estRej; + /// If true, the last timestep was rejected. + bool lastTimestepRejected; + + + /// If true, more than one of the last timesteps were rejected. + bool succRejection; + + /// If greater than 0, than for the first given number of timesteps the timestep + /// will be not changed and is set to the very first one. + int fixFirstTimesteps; + + double tauAcc; + + double estAcc; + + double tauRej; + + double estRej; double tauGamma; };