From 8a79f6861c2a2bbc0ce5c9cfb5c2cdc6bb1ed1e3 Mon Sep 17 00:00:00 2001 From: Andreas Naumann <andreas.naumann@tu-dresden.de> Date: Tue, 30 Nov 2010 08:59:25 +0000 Subject: [PATCH] added test data creation programs --- AMDiS/CMakeLists.txt | 12 +- AMDiS/test/CMakeLists.txt | 2 + AMDiS/test/datacreation/CMakeLists.txt | 13 ++ AMDiS/test/datacreation/src/BallProject.cpp | 111 +++++++++++++++ AMDiS/test/datacreation/src/BallProject.h | 24 ++++ AMDiS/test/datacreation/src/BunnyProject.cpp | 91 ++++++++++++ AMDiS/test/datacreation/src/BunnyProject.h | 23 +++ AMDiS/test/datacreation/src/ElliptProject.cpp | 94 +++++++++++++ AMDiS/test/datacreation/src/ElliptProject.h | 21 +++ AMDiS/test/datacreation/src/HeatProject.cpp | 132 ++++++++++++++++++ AMDiS/test/datacreation/src/HeatProject.h | 28 ++++ AMDiS/test/datacreation/src/HeatProject.hh | 118 ++++++++++++++++ .../test/datacreation/src/NeumannProject.cpp | 99 +++++++++++++ AMDiS/test/datacreation/src/NeumannProject.h | 20 +++ .../datacreation/src/ParametricProject.cpp | 72 ++++++++++ .../test/datacreation/src/ParametricProject.h | 23 +++ .../datacreation/src/ParametricProject.hh | 94 +++++++++++++ .../test/datacreation/src/PeriodicProject.cpp | 90 ++++++++++++ AMDiS/test/datacreation/src/PeriodicProject.h | 20 +++ AMDiS/test/datacreation/src/Project.h | 63 +++++++++ AMDiS/test/datacreation/src/ProjectList.cpp | 38 +++++ AMDiS/test/datacreation/src/SphereProject.cpp | 81 +++++++++++ AMDiS/test/datacreation/src/SphereProject.h | 23 +++ AMDiS/test/datacreation/src/TorusProject.h | 21 +++ .../datacreation/src/VecelliptProject.cpp | 105 ++++++++++++++ .../test/datacreation/src/VecelliptProject.h | 24 ++++ AMDiS/test/datacreation/src/creator_base.cpp | 13 ++ 27 files changed, 1448 insertions(+), 7 deletions(-) create mode 100644 AMDiS/test/CMakeLists.txt create mode 100644 AMDiS/test/datacreation/CMakeLists.txt create mode 100644 AMDiS/test/datacreation/src/BallProject.cpp create mode 100644 AMDiS/test/datacreation/src/BallProject.h create mode 100644 AMDiS/test/datacreation/src/BunnyProject.cpp create mode 100644 AMDiS/test/datacreation/src/BunnyProject.h create mode 100644 AMDiS/test/datacreation/src/ElliptProject.cpp create mode 100644 AMDiS/test/datacreation/src/ElliptProject.h create mode 100644 AMDiS/test/datacreation/src/HeatProject.cpp create mode 100644 AMDiS/test/datacreation/src/HeatProject.h create mode 100644 AMDiS/test/datacreation/src/HeatProject.hh create mode 100644 AMDiS/test/datacreation/src/NeumannProject.cpp create mode 100644 AMDiS/test/datacreation/src/NeumannProject.h create mode 100644 AMDiS/test/datacreation/src/ParametricProject.cpp create mode 100644 AMDiS/test/datacreation/src/ParametricProject.h create mode 100644 AMDiS/test/datacreation/src/ParametricProject.hh create mode 100644 AMDiS/test/datacreation/src/PeriodicProject.cpp create mode 100644 AMDiS/test/datacreation/src/PeriodicProject.h create mode 100644 AMDiS/test/datacreation/src/Project.h create mode 100644 AMDiS/test/datacreation/src/ProjectList.cpp create mode 100644 AMDiS/test/datacreation/src/SphereProject.cpp create mode 100644 AMDiS/test/datacreation/src/SphereProject.h create mode 100644 AMDiS/test/datacreation/src/TorusProject.h create mode 100644 AMDiS/test/datacreation/src/VecelliptProject.cpp create mode 100644 AMDiS/test/datacreation/src/VecelliptProject.h create mode 100644 AMDiS/test/datacreation/src/creator_base.cpp diff --git a/AMDiS/CMakeLists.txt b/AMDiS/CMakeLists.txt index 81447bb1..889e6d90 100644 --- a/AMDiS/CMakeLists.txt +++ b/AMDiS/CMakeLists.txt @@ -28,8 +28,7 @@ set(ENABLE_PARMETIS off) option(ENABLE_UMFPACK "use umfpack solver" false) option(ENABLE_MKL "use the mkl" false) SET(MKL_DIR "" CACHE PATH "MKL directory") -#option(ENABLE_TESTS "compile the tests" false) -SET(ENABLE_TESTS false) +option(ENABLE_TESTS "compile the tests" false) if(ENABLE_INTEL) Message("please set the icc manually") @@ -233,6 +232,7 @@ if(ENABLE_UMFPACK) INSTALL(FILES ${LIB_DIR}/AMD/Lib/libamd.a DESTINATION lib/amdis/amd/) SET(RPM_DEPEND_STR "blas") + LIST(APPEND AMDiS_LIBS amdis blas amd umfpack) endif(ENABLE_UMFPACK) if(ENABLE_MKL) @@ -254,7 +254,7 @@ include_directories(${SOURCE_DIR}) add_library(amdis SHARED ${AMDIS_SRC} ${PARALLEL_DOMAIN_AMDIS_SRC}) add_library(compositeFEM SHARED ${COMPOSITE_FEM_SRC}) target_link_libraries(compositeFEM amdis) - +LIST(APPEND AMDiS_LIBS amdis boost_system boost_iostreams) if(WIN32) SET(COMPILEFLAGS "${COMPILEFLAGS} -D_SCL_SECURE_NO_WARNINGS -D_CRT_SECURE_NO_WARNINGS") endif(WIN32) @@ -329,8 +329,6 @@ set(CPACK_RPM_PACKAGE_REQUIRES "boost-devel >= 1.42, ${RPM_DEPEND_STR}") set(CPACK_DEBIAN_PACKAGE_DEPENDS "libboost-dev (>= 1.42)") include(CPack) if(ENABLE_TESTS) -#ENABLE_TESTING() -#add_test(demo_test run_test.sh) - -#add_test(demo_test ${CMAKE_CTEST_COMMAND} --build-and-test "${CMAKE_SOURCE_DIR}/test/demo_test" "${CMAKE_BINARY_DIR}/test/demo_test" --build-generator "${CMAKE_GENERATOR}" --build-makeprogram "${CMAKE_MAKE_PROGRAM}" "${CMAKE_BINARY_DIR}/test/demo_test") +ENABLE_TESTING() +add_subdirectory(test) endif(ENABLE_TESTS) diff --git a/AMDiS/test/CMakeLists.txt b/AMDiS/test/CMakeLists.txt new file mode 100644 index 00000000..e5de2886 --- /dev/null +++ b/AMDiS/test/CMakeLists.txt @@ -0,0 +1,2 @@ +project(tests) +add_subdirectory(datacreation) diff --git a/AMDiS/test/datacreation/CMakeLists.txt b/AMDiS/test/datacreation/CMakeLists.txt new file mode 100644 index 00000000..45277e22 --- /dev/null +++ b/AMDiS/test/datacreation/CMakeLists.txt @@ -0,0 +1,13 @@ +project(testdatacreation) + message("source-dir: ${AMDiS_SOURCE_DIR}") + include_directories(${AMDiS_SOURCE_DIR}) + file(GLOB PROJECTFILES src/*Project.cpp) + foreach(projectfile ${PROJECTFILES}) + #create creatorname + get_filename_component(CppName ${projectfile}, NAME_WE) + string(REPLACE "Project" "" creatorsuffix ${CppName}) + set(PROJECTINCLUDE "${CppName}.h") + configure_file(src/creator_base.cpp src/creator${creatorsuffix}.cpp @ONLY) + add_executable(creator${creatorsuffix} src/creator${creatorsuffix}.cpp src/ProjectList.cpp ${projectfile}) + target_link_libraries(creator${creatorsuffix} ${AMDiS_LIBS}) + endforeach(projectfile) diff --git a/AMDiS/test/datacreation/src/BallProject.cpp b/AMDiS/test/datacreation/src/BallProject.cpp new file mode 100644 index 00000000..33a45c47 --- /dev/null +++ b/AMDiS/test/datacreation/src/BallProject.cpp @@ -0,0 +1,111 @@ +#include "BallProject.h" +#include "AdaptStationary.h" +/// Dirichlet boundary function +class G : public AbstractFunction<double, WorldVector<double> > +{ +public: + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return exp(-10.0 * (x * x)); + } +}; + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + int dim = x.getSize(); + double r2 = x * x; + double ux = exp(-10.0*r2); + return -(400.0 * r2 - 20.0 * dim) * ux; + } +}; + +BallDemo::BallDemo(): + ball("ball"), + ballCenter(NULL), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{ +} + +BallDemo::~BallDemo() { + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; + if(ballCenter != NULL) + delete ballCenter; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void BallDemo::create(std::string& filename) { + // ===== init parameters ===== + Parameters::init(false, filename); + ballCenter = new WorldVector< double >(); + ballCenter->set(0.0); + + // ===== create projection ===== + new AMDiS::BallProject(1, + BOUNDARY_PROJECTION, + *ballCenter, + 1.0); + + // ===== create and init the scalar problem ===== + ball.initialize(INIT_ALL); + + // === create adapt info === + adaptInfo = new AdaptInfo("ball->adapt"); + + // === create adapt === + adapt = new AdaptStationary("ball->adapt", &ball, adaptInfo); + + // ===== create matrix operator ===== + matrixOperator=new Operator(ball.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + ball.addMatrixOperator(matrixOperator); + + // ===== create rhs operator ===== + int degree = ball.getFeSpace()->getBasisFcts()->getDegree(); + rhsOperator = new Operator(ball.getFeSpace()); + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + ball.addVectorOperator(rhsOperator); + + // ===== add boundary conditions ===== + ball.addDirichletBC(1, new G); + + +} + +int BallDemo::solve(SolutionInformation& solinfo) { + assert(adaptInfo!=NULL); + assert(adapt!=NULL); + assert(matrixOperator!=NULL); + assert(rhsOperator!=NULL); + + int retCode = adapt->adapt(); + solinfo.dofVec = ball.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { + list.clear(); + + Project* demo = new BallDemo(); + ProjectInfo ballInfo(demo, "init/ball.dat.2d", "../testdata/balldata_2d"); + list.push_back(ballInfo); + ballInfo = ProjectInfo(demo, "init/ball.dat.3d", "../testdata/balldata_3d"); + list.push_back(ballInfo); +} diff --git a/AMDiS/test/datacreation/src/BallProject.h b/AMDiS/test/datacreation/src/BallProject.h new file mode 100644 index 00000000..424b626e --- /dev/null +++ b/AMDiS/test/datacreation/src/BallProject.h @@ -0,0 +1,24 @@ +#ifndef BALLPROJECT_H +#define BALLPROJECT_H + +#include "Project.h" +#include "FixVec.h" +#include "ProblemScal.h" + +using namespace AMDiS; + +class BallDemo : public Project { + ProblemScal ball; + WorldVector< double >* ballCenter; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + public: + BallDemo() ; + ~BallDemo() ; + void create(std::string& ) ; + int solve(SolutionInformation&); + +}; +#endif diff --git a/AMDiS/test/datacreation/src/BunnyProject.cpp b/AMDiS/test/datacreation/src/BunnyProject.cpp new file mode 100644 index 00000000..501f7902 --- /dev/null +++ b/AMDiS/test/datacreation/src/BunnyProject.cpp @@ -0,0 +1,91 @@ +#include "BunnyProject.h" +#include "AdaptStationary.h" + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return -2 * x[0]; + } +}; + +/// boundary +class G : public AbstractFunction<double, WorldVector<double> > +{ +public: + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return 10000.0; + } +}; + +Bunnydemo::Bunnydemo(): + bunny("bunny"), + ballCenter(NULL), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{} + +Bunnydemo::~Bunnydemo() { + if(ballCenter != NULL) + delete ballCenter; + if(adaptInfo != NULL) + delete adaptInfo; + if(adapt != NULL) + delete adapt; + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; +} + +void Bunnydemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(false, filename); + + // ===== create projection ===== + ballCenter = new WorldVector< double > (); + ballCenter->set(0.0); + new BallProject(1, VOLUME_PROJECTION, *ballCenter, 1.0); + + // ===== create and init the scalar problem ===== + bunny.initialize(INIT_ALL); + + // === create adapt info === + adaptInfo = new AdaptInfo("bunny->adapt"); + + // === create adapt === + adapt = new AdaptStationary("bunny->adapt", &bunny, adaptInfo); + + // ===== create matrix operator ===== + matrixOperator = new Operator(bunny.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + bunny.addMatrixOperator(matrixOperator); + + // ===== create rhs operator ===== + rhsOperator = new Operator(bunny.getFeSpace()); + + int degree = bunny.getFeSpace()->getBasisFcts()->getDegree(); + + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + bunny.addVectorOperator(rhsOperator); + +} + +int Bunnydemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.dofVec = bunny.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/BunnyProject.h b/AMDiS/test/datacreation/src/BunnyProject.h new file mode 100644 index 00000000..59b897a7 --- /dev/null +++ b/AMDiS/test/datacreation/src/BunnyProject.h @@ -0,0 +1,23 @@ +#ifndef BUNNYPROJECT_H +#define BUNNYPROJECT_H + +#include "Project.h" +#include "FixVec.h" +#include "ProblemScal.h" + +class Bunnydemo : public Project { + ProblemScal bunny; + WorldVector< double >* ballCenter; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + + public: + Bunnydemo(); + ~Bunnydemo(); + + void create(string filename); + int solve(SolutionInformation&); +}; +#endif diff --git a/AMDiS/test/datacreation/src/ElliptProject.cpp b/AMDiS/test/datacreation/src/ElliptProject.cpp new file mode 100644 index 00000000..685364c5 --- /dev/null +++ b/AMDiS/test/datacreation/src/ElliptProject.cpp @@ -0,0 +1,94 @@ +#include "ElliptProject.h" +#include "AdaptStationary.h" +/// Dirichlet boundary function +class G : public AbstractFunction<double, WorldVector<double> > +{ +public: + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return exp(-10.0 * (x * x)); + } +}; + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + int dow = x.getSize(); + double r2 = (x * x); + double ux = exp(-10.0 * r2); + return -(400.0 * r2 - 20.0 * dow) * ux; + } +}; + +Elliptdemo::Elliptdemo(): + ellipt("ellipt"), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{} + +Elliptdemo::~Elliptdemo() { + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void Elliptdemo::create(std::string filename) { + // ===== init parameters ===== + Parameters::init(true, filename); + + + // ===== create and init the scalar problem ===== + ellipt.initialize(INIT_ALL); + + + // === create adapt info === + adaptInfo = new AdaptInfo("ellipt->adapt"); + + + // === create adapt === + adapt = new AdaptStationary("ellipt->adapt", &ellipt, adaptInfo); + + + // ===== create matrix operator ===== + matrixOperator = new Operator(ellipt.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + ellipt.addMatrixOperator(matrixOperator); + + + // ===== create rhs operator ===== + int degree = ellipt.getFeSpace()->getBasisFcts()->getDegree(); + rhsOperator = new Operator(ellipt.getFeSpace()); + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + ellipt.addVectorOperator(rhsOperator); + + + // ===== add boundary conditions ===== + ellipt.addDirichletBC(1, new G); + + +} + +int Elliptdemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.dofVec = ellipt.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/ElliptProject.h b/AMDiS/test/datacreation/src/ElliptProject.h new file mode 100644 index 00000000..50f49690 --- /dev/null +++ b/AMDiS/test/datacreation/src/ElliptProject.h @@ -0,0 +1,21 @@ +#ifndef ELLIPTPROJECT_H +#define ELLIPTPROJECT_H + +#include "Project.h" +#include "ProblemScal.h" + +class Elliptdemo : public Project { + ProblemScal ellipt; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + + public: + Elliptdemo(); + ~Elliptdemo(); + + void create(std::string filename); + int solve(SolutionInformation& info); +}; +#endif diff --git a/AMDiS/test/datacreation/src/HeatProject.cpp b/AMDiS/test/datacreation/src/HeatProject.cpp new file mode 100644 index 00000000..db758b30 --- /dev/null +++ b/AMDiS/test/datacreation/src/HeatProject.cpp @@ -0,0 +1,132 @@ +#include "HeatProject.h" +#include "TimedObject.h" + +/// Dirichlet boundary function +class G : public AbstractFunction<double, WorldVector<double> >, + public TimedObject +{ +public: + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return sin(M_PI * (*timePtr)) * exp(-10.0 * (x * x)); + } +}; + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> >, + public TimedObject +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + int dim = x.getSize(); + double r2 = x * x; + double ux = sin(M_PI * (*timePtr)) * exp(-10.0 * r2); + double ut = M_PI * cos(M_PI * (*timePtr)) * exp(-10.0 * r2); + return ut -(400.0 * r2 - 20.0 * dim) * ux; + } +}; + +Heatdemo::Heatdemo(): + heatSpace("heat->space"), + heat(heatSpace), + adaptInfo(NULL), + adaptInfoInitial(NULL), + A(NULL), + C(NULL), + Fop(NULL) +{} + +Heatdemo::~Heatdemo() { + if(adaptInfo != NULL) + delete adaptInfo; + if(adaptInfoInitial != NULL) + delete adaptInfoInitial; + + if(A != NULL) + delete A; + if(C != NULL) + delete C; + if(Fop != NULL) + delete Fop; +} + +void Heatdemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(false, filename); +// Parameters::readArgv(argc, argv); + + + // ===== create and init stationary problem ===== + heatSpace.initialize(INIT_ALL); + + + // ===== create instationary problem ===== + heat.initialize(INIT_ALL); + + // create adapt info + adaptInfo = new AdaptInfo("heat->adapt"); + + // create initial adapt info + adaptInfoInitial = new AdaptInfo("heat->initial->adapt"); + + // create instationary adapt + adaptInstat = new AdaptInstationary("heat->adapt", + &heatSpace, + adaptInfo, + &heat, + adaptInfoInitial); + + + // ===== create rhs functions ===== + int degree = heatSpace.getFeSpace()->getBasisFcts()->getDegree(); + F *rhsFct = new F(degree); + rhsFct->setTimePtr(heat.getRHSTimePtr()); + + + // ===== create operators ===== + double one = 1.0; + double zero = 0.0; + + // create laplace + A = new Operator(heatSpace.getFeSpace()); + A->addSecondOrderTerm(new Laplace_SOT); + A->setUhOld(heat.getOldSolution()); + if (*(heat.getThetaPtr()) != 0.0) + heatSpace.addMatrixOperator(A, heat.getThetaPtr(), &one); + if (*(heat.getTheta1Ptr()) != 0.0) + heatSpace.addVectorOperator(A, heat.getTheta1Ptr(), &zero); + + // create zero order operator + C = new Operator(heatSpace.getFeSpace()); + C->addZeroOrderTerm(new Simple_ZOT); + C->setUhOld(heat.getOldSolution()); + heatSpace.addMatrixOperator(C, heat.getTau1Ptr(), heat.getTau1Ptr()); + heatSpace.addVectorOperator(C, heat.getTau1Ptr(), heat.getTau1Ptr()); + + // create RHS operator + Fop = new Operator(heatSpace.getFeSpace()); + Fop->addZeroOrderTerm(new CoordsAtQP_ZOT(rhsFct)); + heatSpace.addVectorOperator(Fop); + + + // ===== create boundary functions ===== + G *boundaryFct = new G; + boundaryFct->setTimePtr(heat.getBoundaryTimePtr()); + heat.setExactSolution(boundaryFct); + +} + +int Heatdemo::solve(SolutionInformation& info) { + int retCode = adaptInstat->adapt(); + info.dofVec = heatSpace.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/HeatProject.h b/AMDiS/test/datacreation/src/HeatProject.h new file mode 100644 index 00000000..407b77a7 --- /dev/null +++ b/AMDiS/test/datacreation/src/HeatProject.h @@ -0,0 +1,28 @@ +#ifndef HEATPROJECT_H +#define HEATPROJECT_H + +#include "Project.h" + +#include "ProblemScal.h" +#include "ProblemInstat.h" + +#include "HeatProject.hh" +class Heatdemo : public Project { + ProblemScal heatSpace; + Heat heat; + AdaptInfo* adaptInfo; + AdaptInfo* adaptInfoInitial; + AdaptInstationary* adaptInstat; + Operator* A; + Operator* C; + Operator* Fop; + + public: + Heatdemo(); + ~Heatdemo(); + + void create(string filename); + int solve(SolutionInformation& info); + +}; +#endif diff --git a/AMDiS/test/datacreation/src/HeatProject.hh b/AMDiS/test/datacreation/src/HeatProject.hh new file mode 100644 index 00000000..a68199b8 --- /dev/null +++ b/AMDiS/test/datacreation/src/HeatProject.hh @@ -0,0 +1,118 @@ +// =========================================================================== +// ===== instationary problem ================================================ +// =========================================================================== + +/// Instationary problem +class Heat : public ProblemInstatScal +{ +public: + Heat(ProblemScal &heatSpace) + : ProblemInstatScal("heat", heatSpace) + { + // init theta scheme + theta = -1.0; + GET_PARAMETER(0, name + "->theta", "%f", &theta); + TEST_EXIT(theta >= 0)("theta not set!\n"); + if (theta == 0.0) { + WARNING("You are using the explicit Euler scheme\n"); + WARNING("Use a sufficiently small time step size!!!\n"); + } + MSG("theta = %f\n", theta); + theta1 = theta - 1; + } + + // ===== ProblemInstatBase methods =================================== + + /// set the time in all needed functions! + void setTime(AdaptInfo *adaptInfo) + { + rhsTime = adaptInfo->getTime() - (1 - theta) * adaptInfo->getTimestep(); + boundaryTime = adaptInfo->getTime(); + tau1 = 1.0 / adaptInfo->getTimestep(); + } + + void closeTimestep(AdaptInfo *adaptInfo) + { + ProblemInstatScal::closeTimestep(adaptInfo); + WAIT; + } + + // ===== initial problem methods ===================================== + + /// Used by \ref problemInitial to solve the system of the initial problem + void solve(AdaptInfo *adaptInfo) + { + problemStat->getSolution()->interpol(exactSolution); + } + + /// Used by \ref problemInitial to do error estimation for the initial problem. + void estimate(AdaptInfo *adaptInfo) + { + double errMax, errSum; + + errSum = Error<double>::L2Err(*exactSolution, + *(problemStat->getSolution()), 0, &errMax, false); + adaptInfo->setEstSum(errSum, 0); + adaptInfo->setEstMax(errMax, 0); + } + + // ===== setting methods =============================================== + + /// Sets \ref exactSolution; + void setExactSolution(AbstractFunction<double, WorldVector<double> > *fct) + { + exactSolution = fct; + } + + // ===== getting methods =============================================== + + /// Returns pointer to \ref theta. + double *getThetaPtr() + { + return θ + } + + /// Returns pointer to \ref theta1. + double *getTheta1Ptr() + { + return &theta1; + } + + /// Returns pointer to \ref tau1 + double *getTau1Ptr() + { + return &tau1; + } + + /// Returns pointer to \ref rhsTime. + double *getRHSTimePtr() + { + return &rhsTime; + } + + /// Returns pointer to \ref theta1. + double *getBoundaryTimePtr() + { + return &boundaryTime; + } + +private: + /// Used for theta scheme. + double theta; + + /// theta - 1 + double theta1; + + /// 1.0 / timestep + double tau1; + + /// time for right hand side functions. + double rhsTime; + + /// time for boundary functions. + double boundaryTime; + + /// Pointer to boundary function. Needed for initial problem. + AbstractFunction<double, WorldVector<double> > *exactSolution; +}; + diff --git a/AMDiS/test/datacreation/src/NeumannProject.cpp b/AMDiS/test/datacreation/src/NeumannProject.cpp new file mode 100644 index 00000000..a76479aa --- /dev/null +++ b/AMDiS/test/datacreation/src/NeumannProject.cpp @@ -0,0 +1,99 @@ +#include "NeumannProject.h" +#include "AdaptInfo.h" +#include "AdaptStationary.h" + +class N : public AbstractFunction<double, WorldVector<double> > +{ +public: + double operator()(const WorldVector<double>& x) const + { + return 1.0; + } +}; + +/// Dirichlet boundary function +class G : public AbstractFunction<double, WorldVector<double> > +{ +public: + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return exp(-10.0 * (x * x)); + } +}; + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + int dow = x.getSize(); + double r2 = (x*x); + double ux = exp(-10.0*r2); + return -(400.0*r2 - 20.0*dow)*ux; + } +}; + +Neumanndemo::Neumanndemo(): + neumann("neumann"), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{} + +Neumanndemo::~Neumanndemo() { + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void Neumanndemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(true, filename); + + // ===== create and init the scalar problem ===== + neumann.initialize(INIT_ALL); + + // === create adapt info === + adaptInfo = new AdaptInfo("neumann->adapt", 1); + + // === create adapt === + adapt = new AdaptStationary("neumann->adapt", + &neumann, + adaptInfo); + + // ===== create matrix operator ===== + matrixOperator = new Operator(neumann.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + neumann.addMatrixOperator(matrixOperator); + + // ===== create rhs operator ===== + int degree = neumann.getFeSpace()->getBasisFcts()->getDegree(); + rhsOperator = new Operator(neumann.getFeSpace()); + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + neumann.addVectorOperator(rhsOperator); + + // ===== add boundary conditions ===== + neumann.addNeumannBC(1, new N); + neumann.addDirichletBC(2, new G); + +} + +int Neumanndemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.dofVec = neumann.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/NeumannProject.h b/AMDiS/test/datacreation/src/NeumannProject.h new file mode 100644 index 00000000..495c44f5 --- /dev/null +++ b/AMDiS/test/datacreation/src/NeumannProject.h @@ -0,0 +1,20 @@ +#ifndef NEUMANNPROJECT_H +#define NEUMANNPROJECT_H + +#include "Project.h" +#include "ProblemScal.h" +class Neumanndemo : public Project { + ProblemScal neumann; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + + public: + Neumanndemo(); + ~Neumanndemo(); + + void create(string filename); + int solve(SolutionInformation&); +}; +#endif diff --git a/AMDiS/test/datacreation/src/ParametricProject.cpp b/AMDiS/test/datacreation/src/ParametricProject.cpp new file mode 100644 index 00000000..8eff746c --- /dev/null +++ b/AMDiS/test/datacreation/src/ParametricProject.cpp @@ -0,0 +1,72 @@ +#include "ParametricProject.h" +#include "AdaptStationary.h" +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return -2.0 * x[0]; + } +}; + +Parametricdemo::Parametricdemo(): + parametric("parametric"), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{} + +Parametricdemo::~Parametricdemo() { + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void Parametricdemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(false, filename); + + // ===== create and init the scalar problem ===== + parametric.initialize(INIT_ALL); + + // === create adapt info === + adaptInfo = new AdaptInfo("parametric->adapt", 1); + + // === create adapt === + adapt = new AdaptStationary("parametric->adapt", + ¶metric, + adaptInfo); + + // ===== create matrix operator ===== + matrixOperator = new Operator(parametric.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + parametric.addMatrixOperator(matrixOperator); + + // ===== create rhs operator ===== + rhsOperator = new Operator( parametric.getFeSpace()); + + int degree = parametric.getFeSpace()->getBasisFcts()->getDegree(); + + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + parametric.addVectorOperator(rhsOperator); + +} + +int Parametricdemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.dofVec = parametric.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/ParametricProject.h b/AMDiS/test/datacreation/src/ParametricProject.h new file mode 100644 index 00000000..2e0861c6 --- /dev/null +++ b/AMDiS/test/datacreation/src/ParametricProject.h @@ -0,0 +1,23 @@ +#ifndef PARAMETRICPROJECT_H +#define PARAMETRICPROJECT_H + +#include "Project.h" +#include "ProblemScal.h" + +#include "ParametricProject.hh" +class Parametricdemo : public Project { + ParametricSphere parametric; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + + public: + Parametricdemo(); + ~Parametricdemo(); + + void create(string filename); + int solve(SolutionInformation&); +}; + +#endif diff --git a/AMDiS/test/datacreation/src/ParametricProject.hh b/AMDiS/test/datacreation/src/ParametricProject.hh new file mode 100644 index 00000000..75b1e5ad --- /dev/null +++ b/AMDiS/test/datacreation/src/ParametricProject.hh @@ -0,0 +1,94 @@ +/** \brief + * Implementation of a scalar problem. In \ref buildBeforeCoarsen() parametric + * coordinates for the vertices are filled in a DOFVector. This DOFVector is + * used in \ref parametric to parametrize elements while mesh traversal. + */ +class ParametricSphere : public ProblemScal +{ +public: + /// constructor + ParametricSphere(const char *name) + : ProblemScal(name), + parametric(NULL) + {} + + /// destructor + ~ParametricSphere() + { + for (int i = 0; i < Global::getGeo(WORLD); i++) + delete parametricCoords[i]; + + delete parametric; + } + + /// initialization of the base class and creation of \ref parametric. + void initialize(Flag initFlag, + ProblemScal *adoptProblem = NULL, + Flag adoptFlag = INIT_NOTHING) + { + ProblemScal::initialize(initFlag, adoptProblem, adoptFlag); + + // ===== create projection ===== + WorldVector<double> ballCenter; + ballCenter.set(0.0); + new BallProject(1, + VOLUME_PROJECTION, + ballCenter, + 1.0); + + // ===== create coords vector ===== + for (int i = 0; i < Global::getGeo(WORLD); i++) + parametricCoords[i] = new DOFVector<double>(this->getFeSpace(), + "parametric coords"); + + // ===== create parametric object ===== + parametric = new ParametricFirstOrder(¶metricCoords); + + // ===== enable parametric traverse ===== + this->getMesh()->setParametric(parametric); + }; + + /** \brief + * Implementation of ProblemStatBase::buildBeforeCoarsen(). + */ + void buildBeforeCoarsen(AdaptInfo *adaptInfo, Flag flag) { + FUNCNAME("ParametricSphere::buildAfterCoarsen()"); + MSG("calculation of parametric coordinates\n"); + int preDOFs = this->getFeSpace()->getAdmin()->getNumberOfPreDOFs(VERTEX); + int dim = this->getMesh()->getDim(); + int dow = Global::getGeo(WORLD); + WorldVector<double> vertexCoords; + const DegreeOfFreedom **dof; + DegreeOfFreedom vertexIndex; + + // ===== disable parametric traverse ===== + this->getMesh()->setParametric(NULL); + + TraverseStack stack; + ElInfo *elInfo = NULL; + elInfo = stack.traverseFirst(this->getMesh(), -1, + Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS); + while (elInfo) { + dof = elInfo->getElement()->getDof(); + for (int i = 0; i < dim + 1; i++) { + vertexCoords = elInfo->getCoord(i); + Projection::getProjection(1)->project(vertexCoords); + for (int j = 0; j < dow; j++) + (*(parametricCoords[j]))[dof[i][preDOFs]] = vertexCoords[j]; + } + elInfo = stack.traverseNext(elInfo); + } + + + // ===== enable parametric traverse ===== + this->getMesh()->setParametric(parametric); + } + +protected: + /// DOFVector storing parametric coordinates. + WorldVector<DOFVector<double>*> parametricCoords; + + /// Parametrizes vertex coordinates while mesh traversal. + ParametricFirstOrder *parametric; +}; + diff --git a/AMDiS/test/datacreation/src/PeriodicProject.cpp b/AMDiS/test/datacreation/src/PeriodicProject.cpp new file mode 100644 index 00000000..900cb396 --- /dev/null +++ b/AMDiS/test/datacreation/src/PeriodicProject.cpp @@ -0,0 +1,90 @@ +#include "PeriodicProject.h" +#include "AdaptStationary.h" +/// Dirichlet boundary function +class G : public AbstractFunction<double, WorldVector<double> > +{ +public: + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return 0.0; + } +}; + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + int dim = x.getSize(); + double r2 = (x*x); + double ux = exp(-10.0*r2); + return -(400.0*r2 - 20.0*dim)*ux; + } +}; + +Periodicdemo::Periodicdemo(): + periodic("periodic"), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{} + +Periodicdemo::~Periodicdemo() { + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void Periodicdemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(false, filename); + + // ===== create and init the scalar problem ===== + periodic.initialize(INIT_ALL); + + // === create adapt info === + adaptInfo = new AdaptInfo("periodic->adapt", 1); + + // === create adapt === + adapt = new AdaptStationary("periodic->adapt", + &periodic, + adaptInfo); + + // ===== create matrix operator ===== + matrixOperator = new Operator(periodic.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + periodic.addMatrixOperator(matrixOperator); + + // ===== create rhs operator ===== + int degree = periodic.getFeSpace()->getBasisFcts()->getDegree(); + rhsOperator = new Operator(periodic.getFeSpace()); + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + periodic.addVectorOperator(rhsOperator); + + // ===== add boundary conditions ===== + periodic.addDirichletBC(1, new G); + + // ===== add periodic conditions ===== + periodic.addPeriodicBC(-1); + +} + +int Periodicdemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.dofVec = periodic.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/PeriodicProject.h b/AMDiS/test/datacreation/src/PeriodicProject.h new file mode 100644 index 00000000..420a3bc4 --- /dev/null +++ b/AMDiS/test/datacreation/src/PeriodicProject.h @@ -0,0 +1,20 @@ +#ifndef PERIODICPROJECT_H +#define PERIODICPROJECT_H + +#include "Project.h" +#include "ProblemScal.h" + +class Periodicdemo : public Project { + ProblemScal periodic; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + public: + Periodicdemo(); + ~Periodicdemo(); + + void create(string filename); + int solve(SolutionInformation&); +}; +#endif diff --git a/AMDiS/test/datacreation/src/Project.h b/AMDiS/test/datacreation/src/Project.h new file mode 100644 index 00000000..91fb92ea --- /dev/null +++ b/AMDiS/test/datacreation/src/Project.h @@ -0,0 +1,63 @@ +#ifndef PROJECT_H +#define PROJECT_H + +#include <string> +#include <list> + +#include "DOFVector.h" +#include "SystemVector.h" + +using namespace AMDiS; +using namespace std; + +struct SolutionInformation { + SystemVector* sysVec; + DOFVector< double >* dofVec; + SolutionInformation(): + sysVec(NULL), + dofVec(NULL) + {} + +}; + +void write(SolutionInformation&, std::string filename); +bool compare(SolutionInformation&, std::string filename); +class Project { + public: + virtual void create(std::string& ) = 0; + virtual int solve(SolutionInformation&) = 0; + +}; + +class ProjectInfo { + string initfilename; + Project* project; + string outputfilename; + + public: + ProjectInfo(Project* project, + string initfilename, + string outputfilename): + initfilename(initfilename), + project(project), + outputfilename(outputfilename) {} + + void create() { + assert(initfilename.length() > 0); + assert(project!=NULL); + project->create(initfilename); + } + + int solve(SolutionInformation& info) { + assert(project!=NULL); + return project->solve(info); + } + + inline string getFilename() const { + return outputfilename; + } +}; + +typedef list< ProjectInfo > ProjectList; +void createProjectList(ProjectList& ); +#endif diff --git a/AMDiS/test/datacreation/src/ProjectList.cpp b/AMDiS/test/datacreation/src/ProjectList.cpp new file mode 100644 index 00000000..f42076d3 --- /dev/null +++ b/AMDiS/test/datacreation/src/ProjectList.cpp @@ -0,0 +1,38 @@ +#include "Project.h" + +#include "BallProject.h" + +#include "io/ArhWriter.h" +#include "io/ArhReader.h" +#include <vector> + +using namespace AMDiS; + +void write(SolutionInformation& info, std::string filename) { + if(info.sysVec != NULL) { + std::vector< DOFVector< double >* > data(info.sysVec->getNumVectors()); + Mesh* mesh = info.sysVec->getDOFVector(0)->getFeSpace()->getMesh(); + for(int i=0; i < info.sysVec->getNumVectors(); ++i) { + data[i] = info.sysVec->getDOFVector(i); + } + ArhWriter::write(filename, mesh, data); + }else if(info.dofVec != NULL) { + Mesh* mesh = info.dofVec->getFeSpace()->getMesh(); + ArhWriter::write(filename, mesh, info.dofVec); + }else + assert(false); +} + +bool compare(SolutionInformation& info, std::string filename) { + if(info.sysVec != NULL) { + assert(false); + }else if(info.dofVec != NULL) { + DOFVector< double > fileVec(info.dofVec->getFeSpace(), "fileVec"); + Mesh fileMesh("",2); + fileMesh = *(info.dofVec->getFeSpace()->getMesh()); + ArhReader::read(filename, &fileMesh, &fileVec); + }else + assert(false); + return false; +} + diff --git a/AMDiS/test/datacreation/src/SphereProject.cpp b/AMDiS/test/datacreation/src/SphereProject.cpp new file mode 100644 index 00000000..3ef2e3f5 --- /dev/null +++ b/AMDiS/test/datacreation/src/SphereProject.cpp @@ -0,0 +1,81 @@ +#include "SphereProject.h" +#include "AdaptStationary.h" + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return -2.0 * x[0]; + } +}; + +Spheredemo::Spheredemo(): + sphere("sphere"), + ballCenter(NULL), + adaptInfo(NULL), + adapt(NULL), + matrixOperator(NULL), + rhsOperator(NULL) +{} + +Spheredemo::~Spheredemo() { + if(matrixOperator != NULL) + delete matrixOperator; + if(rhsOperator != NULL) + delete rhsOperator; + if(ballCenter != NULL) + delete ballCenter; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void Spheredemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(false, filename); + + // ===== create projection ===== + ballCenter = new WorldVector< double >; + ballCenter->set(0.0); + new BallProject(1, BOUNDARY_PROJECTION, *ballCenter, 1.0); + + // ===== create and init the scalar problem ===== + sphere.initialize(INIT_ALL); + + // === create adapt info === + adaptInfo = new AdaptInfo("sphere->adapt"); + + // === create adapt === + adapt = new AdaptStationary("sphere->adapt", + &sphere, + adaptInfo); + + // ===== create matrix operator ===== + matrixOperator = new Operator(sphere.getFeSpace()); + matrixOperator->addSecondOrderTerm(new Laplace_SOT); + sphere.addMatrixOperator(matrixOperator); + + // ===== create rhs operator ===== + rhsOperator = new Operator(sphere.getFeSpace()); + + int degree = sphere.getFeSpace()->getBasisFcts()->getDegree(); + + rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + sphere.addVectorOperator(rhsOperator); + +} + +int Spheredemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.dofVec = sphere.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/SphereProject.h b/AMDiS/test/datacreation/src/SphereProject.h new file mode 100644 index 00000000..02bb071a --- /dev/null +++ b/AMDiS/test/datacreation/src/SphereProject.h @@ -0,0 +1,23 @@ +#ifndef SPHEREPROJECT_H +#define SPHEREPROJECT_H + +#include "Project.h" +#include "ProblemScal.h" +#include "FixVec.h" + +class Spheredemo : public Project { + ProblemScal sphere; + WorldVector< double >* ballCenter; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + + public: + Spheredemo(); + ~Spheredemo(); + + void create(string filename); + int solve(SolutionInformation& ); +}; +#endif diff --git a/AMDiS/test/datacreation/src/TorusProject.h b/AMDiS/test/datacreation/src/TorusProject.h new file mode 100644 index 00000000..58126e67 --- /dev/null +++ b/AMDiS/test/datacreation/src/TorusProject.h @@ -0,0 +1,21 @@ +#ifndef TORUSPROJECT_H +#define TORUSPROJECT_H + +#include "Project.h" +#include "ProblemScal.h" + +class Torusdemo : public Project { + ProblemScal torus; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator; + Operator* rhsOperator; + public: + Torusdemo(); + ~Torusdemo(); + + void create(string filename); + int solve(SolutionInformation&); +}; + +#endif diff --git a/AMDiS/test/datacreation/src/VecelliptProject.cpp b/AMDiS/test/datacreation/src/VecelliptProject.cpp new file mode 100644 index 00000000..c7ad0a32 --- /dev/null +++ b/AMDiS/test/datacreation/src/VecelliptProject.cpp @@ -0,0 +1,105 @@ +#include "VecelliptProject.h" +#include "AdaptStationary.h" + +/// Dirichlet boundary function +class G : public AbstractFunction<double, WorldVector<double> > +{ +public: + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + return exp(-10.0 * (x * x)); + } +}; + +/// RHS function +class F : public AbstractFunction<double, WorldVector<double> > +{ +public: + F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {} + + /// Implementation of AbstractFunction::operator(). + double operator()(const WorldVector<double>& x) const + { + int dim = x.getSize(); + double r2 = (x * x); + double ux = exp(-10.0 * r2); + return -(400.0 * r2 - 20.0 * dim) * ux; + }; +}; + +Vecelliptdemo::Vecelliptdemo(): + vecellipt("vecellipt"), + adaptInfo(NULL), + adapt(NULL), + matrixOperator00(NULL), + matrixOperator10(NULL), + matrixOperator11(NULL), + rhsOperator0(NULL) +{} + +Vecelliptdemo::~Vecelliptdemo() { + if(matrixOperator00 != NULL) + delete matrixOperator00; + if(matrixOperator10 != NULL) + delete matrixOperator10; + if(matrixOperator11 != NULL) + delete matrixOperator11; + if(rhsOperator0 != NULL) + delete rhsOperator0; + if(adapt != NULL) + delete adapt; + if(adaptInfo != NULL) + delete adaptInfo; +} + +void Vecelliptdemo::create(string filename) { + // ===== init parameters ===== + Parameters::init(false, filename); + + + // ===== create and init the scalar problem ===== + vecellipt.initialize(INIT_ALL); + + + // === create adapt info === + adaptInfo = new AdaptInfo("vecellipt->adapt", vecellipt.getNumComponents()); + + + // === create adapt === + adapt = new AdaptStationary("vecellipt->adapt", &vecellipt, adaptInfo); + + + // ===== create matrix operators ===== + matrixOperator00 = new Operator(vecellipt.getFeSpace(0), vecellipt.getFeSpace(0)); + matrixOperator00->addSecondOrderTerm(new Laplace_SOT); + vecellipt.addMatrixOperator(matrixOperator00, 0, 0); + + matrixOperator10 = new Operator(vecellipt.getFeSpace(1), vecellipt.getFeSpace(0)); + matrixOperator10->addZeroOrderTerm(new Simple_ZOT); + vecellipt.addMatrixOperator(matrixOperator10, 1, 0); + + matrixOperator11 = new Operator(vecellipt.getFeSpace(1), vecellipt.getFeSpace(1)); + matrixOperator11->addZeroOrderTerm(new Simple_ZOT(-1.0)); + vecellipt.addMatrixOperator(matrixOperator11, 1, 1); + + + // ===== create rhs operator ===== + int degree = vecellipt.getFeSpace(0)->getBasisFcts()->getDegree(); + rhsOperator0 = new Operator(vecellipt.getFeSpace(0)); + rhsOperator0->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree))); + vecellipt.addVectorOperator(rhsOperator0, 0); + + // ===== add boundary conditions ===== + vecellipt.addDirichletBC(1, 0, 0, new G); + +} + +int Vecelliptdemo::solve(SolutionInformation& info) { + int retCode = adapt->adapt(); + info.sysVec = vecellipt.getSolution(); + return retCode; +} + +void createProjectList(ProjectList& list) { +} diff --git a/AMDiS/test/datacreation/src/VecelliptProject.h b/AMDiS/test/datacreation/src/VecelliptProject.h new file mode 100644 index 00000000..c76ca4a4 --- /dev/null +++ b/AMDiS/test/datacreation/src/VecelliptProject.h @@ -0,0 +1,24 @@ +#ifndef VECELLIPTPROJECT_H +#define VECELLIPTPROJECT_H + +#include "Project.h" +#include "ProblemVec.h" + +class Vecelliptdemo : public Project { + ProblemVec vecellipt; + AdaptInfo* adaptInfo; + AdaptStationary* adapt; + Operator* matrixOperator00; + Operator* matrixOperator10; + Operator* matrixOperator11; + Operator* rhsOperator0; + + public: + Vecelliptdemo(); + ~Vecelliptdemo(); + + void create(string filename); + int solve(SolutionInformation& ); +}; + +#endif diff --git a/AMDiS/test/datacreation/src/creator_base.cpp b/AMDiS/test/datacreation/src/creator_base.cpp new file mode 100644 index 00000000..bb776f98 --- /dev/null +++ b/AMDiS/test/datacreation/src/creator_base.cpp @@ -0,0 +1,13 @@ +#include "Project.h" +#include "@PROJECTINCLUDE@" +int main(int argc, char** argv) { + ProjectList list; + createProjectList(list); + for(ProjectList::iterator it=list.begin() ; it != list.end(); ++it) { + it->create(); + SolutionInformation info; + it->solve(info); + write(info, it->getFilename()); + } + return 0; +} -- GitLab