Commit 5055d8c6 authored by Dion Häfner's avatar Dion Häfner

Merge branch 'debug/fem-order-adaptivity' into 'master'

use grid adaptivity in all orders and introduce order=3

This fixes #13. After all, it was just an issue of calling `adapt_grid()` with the correct arguments.

I introduced `FEorder=3` as new supported order and removed `FEorder=0`. As reported to me at the Duner User Meeting, a higher than linear FEM order is desireable for DG problems. So I see no point in retaining a 0th order, i.e. piecewise constant functions.

I tested adaptive grid refinement for all orders, all dimensions, and rectangular and simplex grids. Everything seems to work fine. Notably, the performance difference between the structured YaspGrid and the unstructured UG grid is quite strong for the third order - at least in Dune v2.4.

See merge request !5
parents f147acad 942292e1
......@@ -81,8 +81,7 @@ adding an empty line, make text **bold** or ``monospaced``.
<definition> Order of the finite element method used. Values &gt; 1 are not
thoroughly tested. </definition>
<suggestion> 1 </suggestion>
<values> 0, 1, 2 </values>
<comment> Values &gt; 1 are not thoroughly tested </comment>
<values> 1, 2, 3 </values>
</parameter>
<parameter name="gridFile">
......@@ -182,7 +181,7 @@ adding an empty line, make text **bold** or ``monospaced``.
<parameter name="startTimestep">
<definition> Value of the first time step in seconds. </definition>
<values> float </values>
<suggestion> 1000 </suggestion>
<suggestion> 100 </suggestion>
</parameter>
<parameter name="maxTimestep">
......
......@@ -122,15 +122,15 @@ int main(int argc, char** argv)
auto grid = Dune::Dorie::build_grid_gmsh<Dune::UGGrid<2>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
case 0:
RichardsSolverSimplex<RF,0>(grid,gv,inifile,helper);
break;
case 1:
RichardsSolverSimplex<RF,1>(grid,gv,inifile,helper);
break;
case 2:
RichardsSolverSimplex<RF,2>(grid,gv,inifile,helper);
break;
case 3:
RichardsSolverSimplex<RF,3>(grid,gv,inifile,helper);
break;
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......@@ -142,15 +142,15 @@ int main(int argc, char** argv)
// Avoid the creation of triangles upon refinement
grid->setClosureType(Dune::UGGrid<2>::ClosureType::NONE);
switch(FEorder){
case 0:
RichardsSolverRectangular<RF,0>(grid,gv,inifile,helper);
break;
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
break;
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......@@ -159,15 +159,15 @@ int main(int argc, char** argv)
auto grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<2>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
case 0:
RichardsSolverRectangular<RF,0>(grid,gv,inifile,helper);
break;
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
break;
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......@@ -183,15 +183,15 @@ int main(int argc, char** argv)
auto grid = Dune::Dorie::build_grid_gmsh<Dune::UGGrid<3>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
case 0:
RichardsSolverSimplex<RF,0>(grid,gv,inifile,helper);
break;
case 1:
RichardsSolverSimplex<RF,1>(grid,gv,inifile,helper);
break;
case 2:
RichardsSolverSimplex<RF,2>(grid,gv,inifile,helper);
break;
case 3:
RichardsSolverSimplex<RF,3>(grid,gv,inifile,helper);
break;
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......@@ -203,15 +203,15 @@ int main(int argc, char** argv)
// Avoid the creation of triangles upon refinement
grid->setClosureType(Dune::UGGrid<3>::ClosureType::NONE);
switch(FEorder){
case 0:
RichardsSolverRectangular<RF,0>(grid,gv,inifile,helper);
break;
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
break;
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......@@ -220,15 +220,15 @@ int main(int argc, char** argv)
auto grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<3>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
case 0:
RichardsSolverRectangular<RF,0>(grid,gv,inifile,helper);
break;
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
break;
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......
......@@ -20,6 +20,7 @@ class AdaptivityController
const double alpha; //!< Refinement factor
const double beta; //!< Coarsening factor
const float adaptivityThreshold; //!< Global error threshold below which no refinement is applied
const int intorder; //!< Polynomial order of the local basis functions
int verbose; //!< Verbosity level
const int dim; //!< Spatial dimensions
......@@ -39,6 +40,7 @@ class AdaptivityController
alpha(inifile_.get<float>("adaptivity.refinementFraction")),
beta(inifile_.get<float>("adaptivity.coarseningFraction")),
adaptivityThreshold(inifile_.get<float>("adaptivity.threshold")),
intorder(inifile_.get<int>("grid.FEorder")),
verbose(inifile_.get<int>("output.verbose")),
dim(dim_), param(param_), inifile(inifile_), gv(gv_), fem(fem_)
{
......@@ -134,7 +136,7 @@ class AdaptivityController
t_mark = timer3.elapsed();
timer3.reset();
Dune::PDELab::adapt_grid( *grid, gfs, u, 2 );
Dune::PDELab::adapt_grid( *grid, gfs, u, intorder*2 );
t_adapt = timer3.elapsed();
timer3.reset();
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment