From 26262d3622d601b2808d65829471a501da4f8300 Mon Sep 17 00:00:00 2001 From: Arne Morten Kvarving Date: Mon, 27 Aug 2018 15:03:50 +0200 Subject: [PATCH] added support for K refinement cycles --- CMakeLists.txt | 4 ++++ SIMPoisson.h | 26 ++++++++++++++++++++---- Test/Cube-kref.reg | 44 +++++++++++++++++++++++++++++++++++++++++ Test/Cube-kref.xinp | 40 +++++++++++++++++++++++++++++++++++++ Test/Square-kref.reg | 44 +++++++++++++++++++++++++++++++++++++++++ Test/Square-kref.xinp | 46 +++++++++++++++++++++++++++++++++++++++++++ Test/Square.vreg | 4 ++++ main_Poisson.C | 46 +++++++++++++++++++++++++++++++++++++++++-- 8 files changed, 248 insertions(+), 6 deletions(-) create mode 100644 Test/Cube-kref.reg create mode 100644 Test/Cube-kref.xinp create mode 100644 Test/Square-kref.reg create mode 100644 Test/Square-kref.xinp diff --git a/CMakeLists.txt b/CMakeLists.txt index dd5b464..50b12e2 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -75,6 +75,10 @@ IF(NOT MPI_FOUND OR IFEM_SERIAL_TESTS_IN_PARALLEL) ifem_add_vtf_test(round_stair_3patch_3D_LR_adap.vreg Poisson) endif() ifem_add_vtf_test(Square.vreg Poisson) + if(PETSC_FOUND) + ifem_add_test(Cube-kref.reg Poisson) + ifem_add_test(Square-kref.reg Poisson) + endif() endif() list(APPEND TEST_APPS Poisson) diff --git a/SIMPoisson.h b/SIMPoisson.h index 8abb6a6..4f6a435 100644 --- a/SIMPoisson.h +++ b/SIMPoisson.h @@ -16,15 +16,22 @@ #include "SIMMultiPatchModelGen.h" #include "Poisson.h" +#include "AlgEqSystem.h" #include "AnaSol.h" +#include "ASMbase.h" #include "SIM1D.h" #include "SIM2D.h" #include "SIM3D.h" #include "Utilities.h" #include "DataExporter.h" #include "IFEM.h" +#include "LinSolParams.h" #include "tinyxml.h" +#include "TimeStep.h" #include +#ifdef HAS_PETSC +#include "PETScMatrix.h" +#endif /*! @@ -34,6 +41,8 @@ template class SIMPoisson : public SIMMultiPatchModelGen { public: + using SetupProps = bool; //!< A bool is the only setup property + //! \brief Default constructor. explicit SIMPoisson(bool checkRHS = false) : SIMMultiPatchModelGen(1,checkRHS), @@ -137,8 +146,8 @@ template class SIMPoisson : public SIMMultiPatchModelGen return this->writeGlvBC(nBlock); } - //! \brief Assembles and solves the linear system. - bool solveStep(TimeStep&) + //! \brief Initialize simulator. + bool init() { if (!this->setMode(Dim::opt.eig == 0 ? SIM::STATIC : SIM::VIBRATION)) return false; @@ -149,6 +158,12 @@ template class SIMPoisson : public SIMMultiPatchModelGen this->initSystem(Dim::opt.solver,1,0); this->setQuadratureRule(Dim::opt.nGauss[0],true); + return true; + } + + //! \brief Assembles and solves the linear system. + bool solveStep(TimeStep&) + { if (!this->assembleSystem()) return false; else if (vizRHS && Dim::opt.eig == 0) @@ -187,7 +202,7 @@ template class SIMPoisson : public SIMMultiPatchModelGen } //! \brief Saves solution-dependent quantities to file for postprocessing. - bool saveStep(TimeStep&, int& nBlock) + bool saveStep(TimeStep& tp, int& nBlock) { if (!asciiFile.empty()) { @@ -253,7 +268,10 @@ template class SIMPoisson : public SIMMultiPatchModelGen if (!this->writeGlvN(myNorm,1,nBlock,prefix.data())) return false; - return this->writeGlvStep(1,0.0,1); + if (tp.iter == 0) + return this->writeGlvStep(1,0.0,1); + else + return this->writeGlvStep(tp.iter,tp.iter,1); } //! \brief Prints a norm group to the log stream. diff --git a/Test/Cube-kref.reg b/Test/Cube-kref.reg new file mode 100644 index 0000000..d6d3c91 --- /dev/null +++ b/Test/Cube-kref.reg @@ -0,0 +1,44 @@ +Cube-kref.xinp -kref + + Doing 3 k-refinement cycles. + K-solver tolerance 1e-07 +Number of elements 64 +Number of nodes 125 +Number of dofs 125 +Number of unknowns 27 +K-refinement cycle 1 +Number of elements 64 +Number of nodes 216 +Number of dofs 216 +Number of unknowns 64 +L2-norm : 0.242676 +Max displacement : 0.994383 +Energy norm a(u^h,u^h)^0.5 : 1.92322 +External energy (h,u^h)^0.5 : 1.92322 +Exact norm a(u,u)^0.5 : 1.92382 +Exact error a(e,e)^0.5, e=u-u^h : 0.0483303 +Exact relative error (%) : 2.5122 +K-refinement cycle 2 +Number of elements 64 +Number of nodes 343 +Number of dofs 343 +Number of unknowns 125 +L2-norm : 0.225482 +Max displacement : 1.38805 +Energy norm a(u^h,u^h)^0.5 : 1.92374 +External energy (h,u^h)^0.5 : 1.92374 +Exact norm a(u,u)^0.5 : 1.92382 +Exact error a(e,e)^0.5, e=u-u^h : 0.0180016 +Exact relative error (%) : 0.935719 +K-refinement cycle 3 +Number of elements 64 +Number of nodes 512 +Number of dofs 512 +Number of unknowns 216 +L2-norm : 0.227619 +Max displacement : 1.16167 +Energy norm a(u^h,u^h)^0.5 : 1.92382 +External energy (h,u^h)^0.5 : 1.92382 +Exact norm a(u,u)^0.5 : 1.92382 +Exact error a(e,e)^0.5, e=u-u^h : 0.000660832 +Exact relative error (%) : 0.0343499 diff --git a/Test/Cube-kref.xinp b/Test/Cube-kref.xinp new file mode 100644 index 0000000..fd33ecc --- /dev/null +++ b/Test/Cube-kref.xinp @@ -0,0 +1,40 @@ + + + + + + + + + + + 1 2 3 4 5 6 + + + + + + + + + + + + 3*PI*PI*sin(x*PI)*sin(y*PI)*sin(z*PI) + + sin(PI*x)*sin(PI*y)*sin(z*PI) + -PI*cos(PI*x)*sin(PI*y)*sin(z*PI) | + -PI*sin(PI*x)*cos(PI*y)*sin(PI*z) | + -PI*sin(PI*x)*sin(PI*y)*cos(PI*z) + + + + + + + + cg + 1e-9 + + + diff --git a/Test/Square-kref.reg b/Test/Square-kref.reg new file mode 100644 index 0000000..8d62ed9 --- /dev/null +++ b/Test/Square-kref.reg @@ -0,0 +1,44 @@ +Square-kref.xinp -kref + + Doing 3 k-refinement cycles. + K-solver tolerance 1e-07 +Number of elements 64 +Number of nodes 81 +Number of dofs 81 +Number of unknowns 72 +K-refinement cycle 1 +Number of elements 64 +Number of nodes 100 +Number of dofs 100 +Number of unknowns 90 +L2-norm : 0.986916 +Max displacement : 2.00012 +Energy norm a(u^h,u^h)^0.5 : 5.32009 +External energy (h,u^h)^0.5 : 5.32009 +Exact norm a(u,u)^0.5 : 5.32155 +Exact error a(e,e)^0.5, e=u-u^h : 0.125387 +Exact relative error (%) : 2.3562 +K-refinement cycle 2 +Number of elements 64 +Number of nodes 121 +Number of dofs 121 +Number of unknowns 110 +L2-norm : 1.02204 +Max displacement : 2.2147 +Energy norm a(u^h,u^h)^0.5 : 5.32152 +External energy (h,u^h)^0.5 : 5.32152 +Exact norm a(u,u)^0.5 : 5.32155 +Exact error a(e,e)^0.5, e=u-u^h : 0.0180157 +Exact relative error (%) : 0.338542 +K-refinement cycle 3 +Number of elements 64 +Number of nodes 144 +Number of dofs 144 +Number of unknowns 132 +L2-norm : 1.04048 +Max displacement : 2.10254 +Energy norm a(u^h,u^h)^0.5 : 5.32155 +External energy (h,u^h)^0.5 : 5.32155 +Exact norm a(u,u)^0.5 : 5.32155 +Exact error a(e,e)^0.5, e=u-u^h : 0.00224356 +Exact relative error (%) : 0.0421598 diff --git a/Test/Square-kref.xinp b/Test/Square-kref.xinp new file mode 100644 index 0000000..19b119e --- /dev/null +++ b/Test/Square-kref.xinp @@ -0,0 +1,46 @@ + + + + + + + + + + + 4 + + + 3 + + + + + + + + + + + + + PI*PI*cos(PI*x)*(2-y) + + cos(PI*x)*(2-y) + PI*sin(PI*x)*(2-y)|cos(PI*x) + + 1.0 + + + + + + cg + 1e-9 + + + + + + + diff --git a/Test/Square.vreg b/Test/Square.vreg index ae04c87..fd07adf 100644 --- a/Test/Square.vreg +++ b/Test/Square.vreg @@ -16,3 +16,7 @@ Step 1 Scalar: 'FE q_y' Scalar: 'Exact q_x' Scalar: 'Exact q_y' + Element block 2 + 32 nodes + 32 elements + Vector: 'Heat flux' diff --git a/main_Poisson.C b/main_Poisson.C index 78ea814..78ee89d 100644 --- a/main_Poisson.C +++ b/main_Poisson.C @@ -15,6 +15,7 @@ #include "SIMPoisson.h" #include "SIMargsBase.h" #include "SIMSolverAdap.h" +#include "SIMSolverKRef.h" #include "Utilities.h" #include "Profiler.h" @@ -67,10 +68,31 @@ int runSimulator(char* infile, bool checkRHS, if (model.opt.dumpHDF5(infile)) solver.handleDataOutput(model.opt.hdf5); + model.init(); return solver.solveProblem(infile,"Solving the Poisson problem"); } +#ifdef HAS_PETSC +template +int runSimulatorK(char* infile, bool checkRHS, + const std::vector& ignoredPatches, + bool fixDup, bool vizRHS, bool dumpASCII) +{ + SIMKRefWrap> model(checkRHS), model2(checkRHS); + SIMSolverKRef> solver(model, model2); + + if (!solver.read(infile)) + return 1; + + if (model.opt.dumpHDF5(infile)) + solver.handleDataOutput(model.opt.hdf5); + + return solver.solveProblem(infile,"Solving the Poisson problem"); +} +#endif + + /*! \brief Main program for the NURBS-based isogeometric Poisson equation solver. @@ -119,6 +141,7 @@ int main (int argc, char** argv) bool vizRHS = false; bool fixDup = false; bool dumpASCII = false; + bool kref = false; char* infile = nullptr; SIMargsBase args("poisson"); @@ -137,6 +160,8 @@ int main (int argc, char** argv) SIMbase::ignoreDirichlet = true; else if (!strcmp(argv[i],"-checkRHS")) checkRHS = true; + else if (!strcmp(argv[i],"-kref")) + kref = true; else if (!strcmp(argv[i],"-vizRHS")) vizRHS = true; else if (!strcmp(argv[i],"-fixDup")) @@ -197,8 +222,25 @@ int main (int argc, char** argv) case 3: return runSimulator(infile, checkRHS, ignoredPatches, fixDup, vizRHS, dumpASCII); - } - else + } + else if (kref) { +#ifndef HAS_PETSC + std::cerr << "K-refinement requested, but built without PETSc support. Bailing.." << std::endl; + return 1; +#else + switch (args.dim) { + case 1: + return runSimulatorK(infile, checkRHS, ignoredPatches, + fixDup, vizRHS, dumpASCII); + case 2: + return runSimulatorK(infile, checkRHS, ignoredPatches, + fixDup, vizRHS, dumpASCII); + case 3: + return runSimulatorK(infile, checkRHS, ignoredPatches, + fixDup, vizRHS, dumpASCII); + } +#endif + } else switch (args.dim) { case 1: return runSimulator(infile, checkRHS, ignoredPatches,