Commit 441579ca authored by Lukas Riedel's avatar Lukas Riedel

adapted dorie.cc to interface changes. added headers for subfolder includes....

adapted dorie.cc to interface changes. added headers for subfolder includes. changed factory enable_if
parent b655a279
......@@ -2,4 +2,7 @@ add_subdirectory(interface)
add_subdirectory(solver)
add_executable("dorie" dorie.cc)
target_link_dune_default_libraries("dorie")
\ No newline at end of file
target_link_dune_default_libraries("dorie")
add_executable("sim" EXCLUDE_FROM_ALL sim.cc)
target_link_dune_default_libraries("sim")
\ No newline at end of file
// -*- tab-width: 2; indent-tabs-mode: nil -*-
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <dune/dorie/dorie.hh> // System and DUNE Headers
#include <dune/dorie/solver/util_typedefs.hh>
#include <dune/dorie/solver/util_vtk_adapters.hh>
#include <dune/dorie/solver/util_grid_creator.hh>
#include <dune/dorie/solver/util_boundary_condition.hh>
#include <dune/dorie/solver/util_controller.hh>
#include <dune/dorie/solver/util_adaptivity.hh>
#include <dune/dorie/solver/util_h5tools.hh>
#include <dune/dorie/solver/richards_initial.hh>
#include <dune/dorie/solver/richards_boundary.hh>
#include <dune/dorie/solver/richards_source.hh>
#include <dune/dorie/solver/util_interpolator.hh>
#include <dune/dorie/solver/param_base.hh>
#include <dune/dorie/solver/param_van_genuchten.hh>
#include <dune/dorie/solver/param_factory.hh>
#include <dune/dorie/solver/operator_DG.hh>
#include <dune/dorie/solver/operator_error_indicator.hh>
#include <dune/dorie/solver/dune_implicitonestep.hh>
#include <dune/dorie/solver/solver.hh>
#include <dune/dorie/solver/solver_prep.hh>
//===============================================================
......@@ -38,6 +22,29 @@
* \see Dune::Dorie::FlowSource Class handling Source Term query functions
*/
/// Resolve the second (default) template parameter of YaspGrid
template<int dim>
using YaspGrid = Dune::YaspGrid<dim>;
template<typename Traits>
using Sim = Dune::Dorie::Simulation<Traits>;
template<int dim, int order>
using Simplex = Dune::Dorie::BaseTraits<Dune::UGGrid,
Dune::GeometryType::BasicType::simplex,dim,order,true,false>;
template<int dim, int order>
using SimplexAdaptive = Dune::Dorie::BaseTraits<Dune::UGGrid,
Dune::GeometryType::BasicType::simplex,dim,order,true,true>;
template<int dim, int order>
using Cube = Dune::Dorie::BaseTraits<YaspGrid,
Dune::GeometryType::BasicType::cube,dim,order,true,false>;
template<int dim, int order>
using CubeAdaptive = Dune::Dorie::BaseTraits<Dune::UGGrid,
Dune::GeometryType::BasicType::cube,dim,order,true,true>;
/// Main Program Function: Initialize parameter file reader, build grid and call Richards Solver.
/** As simplex and rectangular grids demand different FiniteElementMaps, the program calls different functions for these tasks.
* The objects and types are then passed to the generic solver function.\n
......@@ -100,9 +107,6 @@ int main(int argc, char** argv)
if (result != 0)
DUNE_THROW(Dune::IOError,"Output folder " << outputPath << " not writable");
// Define precision (RangeField)
typedef double RF;
if (helper.rank()==0){
std::cout << "INPUT FILE: " << inifilename << std::endl;
std::cout << "BC FILE: " << inifile.get<std::string>("boundary.file") << std::endl;
......@@ -117,54 +121,90 @@ int main(int argc, char** argv)
{
if (gtype == "gmsh"){
auto grid = Dune::Dorie::build_grid_gmsh<Dune::UGGrid<2>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
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!");
if(adaptivity){
switch(FEorder){
case 1:{
Sim<SimplexAdaptive<2,1>> sim(helper,grid,inifile);
sim.run();
break;
}
case 2:{
Sim<SimplexAdaptive<2,2>> sim(helper,grid,inifile);
sim.run();
break;
}
case 3:{
Sim<SimplexAdaptive<2,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
else{ // no adaptivity
switch(FEorder){
case 1:{
Sim<Simplex<2,1>> sim(helper,grid,inifile);
sim.run();
break;
}
case 2:{
Sim<Simplex<2,2>> sim(helper,grid,inifile);
sim.run();
break;
}
case 3:{
Sim<Simplex<2,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
}
else if (gtype == "rectangular"){
if(adaptivity){
auto grid = Dune::Dorie::build_grid_cube<Dune::UGGrid<2>>(inifile,helper);
auto gv = grid->leafGridView();
// Avoid the creation of triangles upon refinement
grid->setClosureType(Dune::UGGrid<2>::ClosureType::NONE);
switch(FEorder){
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
case 1:{
Sim<CubeAdaptive<2,1>> sim(helper,grid,inifile);
sim.run();
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
}
case 2:{
Sim<CubeAdaptive<2,2>> sim(helper,grid,inifile);
sim.run();
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
}
case 3:{
Sim<CubeAdaptive<2,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
else{ // no adaptivity
auto grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<2>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
case 1:{
Sim<Cube<2,1>> sim(helper,grid,inifile);
sim.run();
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
}
case 2:{
Sim<Cube<2,2>> sim(helper,grid,inifile);
sim.run();
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
}
case 3:{
Sim<Cube<2,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
......@@ -178,61 +218,97 @@ int main(int argc, char** argv)
{
if (gtype == "gmsh"){
auto grid = Dune::Dorie::build_grid_gmsh<Dune::UGGrid<3>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
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!");
if(adaptivity){
switch(FEorder){
case 1:{
Sim<SimplexAdaptive<3,1>> sim(helper,grid,inifile);
sim.run();
break;
}
case 2:{
Sim<SimplexAdaptive<3,2>> sim(helper,grid,inifile);
sim.run();
break;
}
case 3:{
Sim<SimplexAdaptive<3,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
else{ // no adaptivity
switch(FEorder){
case 1:{
Sim<Simplex<3,1>> sim(helper,grid,inifile);
sim.run();
break;
}
case 2:{
Sim<Simplex<3,2>> sim(helper,grid,inifile);
sim.run();
break;
}
case 3:{
Sim<Simplex<3,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
}
else if (gtype == "rectangular"){
if(adaptivity){
auto grid = Dune::Dorie::build_grid_cube<Dune::UGGrid<3>>(inifile,helper);
auto gv = grid->leafGridView();
// Avoid the creation of triangles upon refinement
grid->setClosureType(Dune::UGGrid<3>::ClosureType::NONE);
switch(FEorder){
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
case 1:{
Sim<CubeAdaptive<3,1>> sim(helper,grid,inifile);
sim.run();
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
}
case 2:{
Sim<CubeAdaptive<3,2>> sim(helper,grid,inifile);
sim.run();
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
}
case 3:{
Sim<CubeAdaptive<3,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
else{ // no adaptivity
auto grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<3>>(inifile,helper);
auto gv = grid->leafGridView();
switch(FEorder){
case 1:
RichardsSolverRectangular<RF,1>(grid,gv,inifile,helper);
case 1:{
Sim<Cube<3,1>> sim(helper,grid,inifile);
sim.run();
break;
case 2:
RichardsSolverRectangular<RF,2>(grid,gv,inifile,helper);
}
case 2:{
Sim<Cube<3,2>> sim(helper,grid,inifile);
sim.run();
break;
case 3:
RichardsSolverRectangular<RF,3>(grid,gv,inifile,helper);
}
case 3:{
Sim<Cube<3,3>> sim(helper,grid,inifile);
sim.run();
break;
}
default:
DUNE_THROW(Dune::NotImplemented,"Finite Element Order (grid.FEorder) not supported!");
}
}
}
else
DUNE_THROW(Dune::NotImplemented,"Grid Type (grid.gridType) not supported!");
DUNE_THROW(Dune::NotImplemented,"Grid Type not supported!");
}
// grid_dim != 2,3
......
#ifndef DORIE_HH
#define DORIE_HH
//#pragma GCC system_header
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <math.h>
#include <iostream>
#include <vector>
#include <map>
......@@ -22,47 +15,39 @@
#include <unistd.h>
#include <hdf5.h>
#include <dune/grid/uggrid.hh>
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <dune/common/ios_state.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/common/exceptions.hh>
#include <dune/common/fvector.hh>
#include <dune/common/timer.hh>
#include <dune/common/parametertreeparser.hh>
#include <dune/common/std/memory.hh>
#include <dune/grid/uggrid.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
#include <dune/grid/io/file/gmshreader.hh>
#include <dune/grid/utility/structuredgridfactory.hh>
#include <dune/istl/bvector.hh>
#include <dune/istl/operators.hh>
#include <dune/istl/solvers.hh>
#include <dune/istl/preconditioners.hh>
#include <dune/istl/io.hh>
#include <dune/istl/superlu.hh>
#include <dune/istl/owneroverlapcopy.hh>
#include <dune/pdelab/newton/newton.hh>
#include <dune/geometry/type.hh>
#include <dune/pdelab/adaptivity/adaptivity.hh>
#include <dune/pdelab/backend/istl.hh>
#include <dune/pdelab/backend/interface.hh>
#include <dune/pdelab/common/function.hh>
#include <dune/pdelab/common/vtkexport.hh>
#include <dune/pdelab/finiteelementmap/qkfem.hh>
#include <dune/pdelab/finiteelementmap/p0fem.hh>
#include <dune/pdelab/finiteelementmap/pkfem.hh>
#include <dune/pdelab/finiteelementmap/qkdg.hh>
#include <dune/pdelab/constraints/p0.hh>
#include <dune/pdelab/constraints/p0ghost.hh>
#include <dune/pdelab/constraints/common/constraints.hh>
#include <dune/pdelab/constraints/conforming.hh>
#include <dune/pdelab/finiteelementmap/qkdg.hh>
#include <dune/pdelab/finiteelementmap/pkfem.hh>
#include <dune/pdelab/gridfunctionspace/gridfunctionspace.hh>
#include <dune/pdelab/gridfunctionspace/gridfunctionspaceutilities.hh>
#include <dune/pdelab/gridfunctionspace/genericdatahandle.hh>
#include <dune/pdelab/gridfunctionspace/interpolate.hh>
#include <dune/pdelab/gridfunctionspace/vtk.hh>
#include <dune/pdelab/gridfunctionspace/subspace.hh>
#include <dune/pdelab/gridoperator/gridoperator.hh>
#include <dune/pdelab/gridoperator/onestep.hh>
#include <dune/pdelab/backend/istl.hh>
#include <dune/pdelab/backend/interface.hh>
#include <dune/pdelab/stationary/linearproblem.hh>
#include <dune/pdelab/instationary/onestep.hh>
#include <dune/pdelab/adaptivity/adaptivity.hh>
#include <dune/pdelab/newton/newton.hh>
#include <dune/pdelab/stationary/linearproblem.hh>
#include <dune/dorie/solver.hh>
#include <dune/dorie/interface.hh>
#endif // DORIE_HH
#ifndef DUNE_DORIE_INTERFACE_HH
#define DUNE_DORIE_INTERFACE_HH
#include <dune/dorie/interface/util.hh>
#include <dune/dorie/interface/adaptivity.hh>
#include <dune/dorie/interface/output.hh>
#include <dune/dorie/interface/simulation.hh>
#include <dune/dorie/interface/knofu.hh>
#endif // DUNE_DORIE_INTERFACE_HH
add_executable("sim" EXCLUDE_FROM_ALL sim.cc)
target_link_dune_default_libraries("sim")
\ No newline at end of file
......@@ -229,30 +229,31 @@ private:
public:
/// Create the factory, taking the parameters for building an AdaptivityHandler
/** Assert that adaptive grid (UG) is used.
/**
* \param _inifile Parameter file parser
* \param _grid Grid to adapt (reference is not saved)
*/
AdaptivityHandlerFactory (Dune::ParameterTree& _inifile, Grid& _grid) :
inifile(_inifile),
grid(_grid)
{
static_assert(std::is_same<Grid,Dune::UGGrid<2>>::value ||
std::is_same<Grid,Dune::UGGrid<3>>::value,
"Adaptivity is only implemented for UGGrid!");
}
{ }
/// Create a dummy AdaptivityHandler
/** Assert that adaptive grid (UG) is used.
*/
template<bool active = enabled>
std::unique_ptr<AHB> create (typename std::enable_if_t<!active,int> = 0)
std::enable_if_t<!active,std::unique_ptr<AHB>> create ()
{
return std::make_unique<AHB>();
}
/// Create a true AdaptivityHandler
template<bool active = enabled>
std::unique_ptr<AHB> create (typename std::enable_if_t<active,int> = 0)
std::enable_if_t<active,std::unique_ptr<AHB>> create ()
{
static_assert(std::is_same<Grid,Dune::UGGrid<2>>::value ||
std::is_same<Grid,Dune::UGGrid<3>>::value,
"Adaptivity is only implemented for UGGrid!");
std::unique_ptr<AHB> p;
p = std::make_unique<AdaptivityHandler<Traits,GFS,Param,U>>(inifile,grid);
return p;
......
......@@ -2,75 +2,7 @@
#include "config.h"
#endif
#include <iostream>
#include <vector>
#include <map>
#include <string>
#include <stdlib.h>
#include <fstream>
#include <cmath>
#include <memory>
#include <utility>
#include <sys/stat.h>
#include <sys/types.h>
#include <unistd.h>
#include <hdf5.h>
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <dune/common/ios_state.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/common/exceptions.hh>
#include <dune/common/timer.hh>
#include <dune/common/parametertreeparser.hh>
#include <dune/grid/uggrid.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
#include <dune/grid/io/file/gmshreader.hh>
#include <dune/grid/utility/structuredgridfactory.hh>
#include <dune/geometry/type.hh>
#include <dune/pdelab/adaptivity/adaptivity.hh>
#include <dune/pdelab/backend/istl.hh>
#include <dune/pdelab/backend/interface.hh>
#include <dune/pdelab/common/function.hh>
#include <dune/pdelab/constraints/p0.hh>
#include <dune/pdelab/constraints/p0ghost.hh>
#include <dune/pdelab/constraints/common/constraints.hh>
#include <dune/pdelab/finiteelementmap/qkdg.hh>
#include <dune/pdelab/finiteelementmap/pkfem.hh>
#include <dune/pdelab/gridfunctionspace/gridfunctionspace.hh>
#include <dune/pdelab/gridoperator/gridoperator.hh>
#include <dune/pdelab/gridoperator/onestep.hh>
#include <dune/pdelab/instationary/onestep.hh>
#include <dune/pdelab/newton/newton.hh>
#include <dune/pdelab/stationary/linearproblem.hh>
#include <dune/dorie/solver/util_typedefs.hh>
#include <dune/dorie/solver/util_vtk_adapters.hh>
#include <dune/dorie/solver/util_grid_creator.hh>
#include <dune/dorie/solver/util_boundary_condition.hh>
#include <dune/dorie/solver/util_controller.hh>
#include <dune/dorie/solver/util_adaptivity.hh>
#include <dune/dorie/solver/util_h5tools.hh>
#include <dune/dorie/solver/richards_initial.hh>
#include <dune/dorie/solver/richards_boundary.hh>
#include <dune/dorie/solver/richards_source.hh>
#include <dune/dorie/solver/util_interpolator.hh>
#include <dune/dorie/solver/param_base.hh>
#include <dune/dorie/solver/param_van_genuchten.hh>
#include <dune/dorie/solver/param_factory.hh>
#include <dune/dorie/solver/operator_DG.hh>
#include <dune/dorie/solver/operator_error_indicator.hh>
#include <dune/dorie/solver/dune_implicitonestep.hh>
#include <dune/dorie/interface/util.hh>
#include <dune/dorie/interface/output.hh>
#include <dune/dorie/interface/adaptivity.hh>
#include <dune/dorie/interface/simulation.hh>
#include <dune/dorie/interface/knofu.hh>
#include <dune/dorie/dorie.hh>
template<int dim>
using YaspGrid = Dune::YaspGrid<dim>;
......@@ -100,7 +32,7 @@ int main(int argc, char **argv)
auto yasp_grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<2>>(inifile,helper);
auto ug_grid = Dune::Dorie::build_grid_cube<Dune::UGGrid<2>>(inifile,helper);
Dune::Dorie::Simulation<Dune::Dorie::BaseTraits<Dune::UGGrid,Dune::GeometryType::BasicType::cube,2,1,true,true>> sim(helper,ug_grid,inifile);
Dune::Dorie::Simulation<Dune::Dorie::BaseTraits<YaspGrid,Dune::GeometryType::BasicType::cube,2,1,true,false>> sim(helper,yasp_grid,inifile);
//sim.run();
......
#ifndef DUNE_DORIE_SOLVER_HH
#define DUNE_DORIE_SOLVER_HH
#include <dune/dorie/solver/util_typedefs.hh>
#include <dune/dorie/solver/util_vtk_adapters.hh>
#include <dune/dorie/solver/util_grid_creator.hh>
#include <dune/dorie/solver/util_boundary_condition.hh>
#include <dune/dorie/solver/util_controller.hh>
#include <dune/dorie/solver/util_adaptivity.hh>
#include <dune/dorie/solver/util_h5tools.hh>
#include <dune/dorie/solver/richards_initial.hh>
#include <dune/dorie/solver/richards_boundary.hh>
#include <dune/dorie/solver/richards_source.hh>
#include <dune/dorie/solver/util_interpolator.hh>
#include <dune/dorie/solver/param_base.hh>
#include <dune/dorie/solver/param_van_genuchten.hh>
#include <dune/dorie/solver/param_factory.hh>
#include <dune/dorie/solver/operator_DG.hh>
#include <dune/dorie/solver/operator_error_indicator.hh>
#include <dune/dorie/solver/dune_implicitonestep.hh>
#endif // DUNE_DORIE_SOLVER_HH
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