Commit 97c8809d authored by Lukas Riedel's avatar Lukas Riedel

Make Simulation build GridCreator by itself. Tests broken.

* Remove grid creation from executables.
* RichardsSimulation now builds its own Grid(Creator).
* Add new input key `parameters.file` to `parameters.xml` and remove all
    other keys of this section.
* Adapt other test executables.
* Add test for grid creation.
    * Add mesh files for grid creation test.
    * Add Python script for writing mapping files to be read by the
        test.
* Remove `grid-test.cc`.

The testing suite is generally not adapted to the new input scheme.
All tests of the main executable tests are broken.
parent 9c6599fa
......@@ -159,12 +159,12 @@ test:mass_conserve:
- $DUNECONTROL --only=dorie configure
- $DUNECONTROL --only=dorie make test_mass_conservation
test:new-parameters:
test:grid-creator:
<<: *test
script:
- $DUNECONTROL --only=dorie configure
- $DUNECONTROL --only=dorie make test_param
- python3 test/maps/create_param_maps.py test/maps/cell_ids.h5
- $DUNECONTROL --only=dorie make test_grid_creator
# --- Deploy jobs ---
deploy:dockerhub-devel: &deploy
......
......@@ -135,6 +135,21 @@ adding an empty line, make text **bold** or ``monospaced``.
<values> int </values>
<suggestion> 0 </suggestion>
</parameter>
<parameter name="mappingFile">
<definition> The H5 file containing the mapping from cell to medium index.
Specify the dataset inside the file with the next key.
</definition>
<values> path </values>
</parameter>
<parameter name="mappingFileDataset">
<definition> The H5 dataset inside ``grid.mappingFile`` containing the
mapping from cell to medium index.
</definition>
<values> path </values>
</parameter>
</category>
<category name="boundary">
......@@ -248,33 +263,10 @@ adding an empty line, make text **bold** or ``monospaced``.
</category>
<category name="parameters">
<parameter name="arrayFile">
<definition> Path to an HDF5 array file containing all parameters used in the chosen
parameterization. Can be created using the ``dorie pfg`` command. </definition>
<values> path </values>
<comment> Create with 'dorie pfg' </comment>
</parameter>
<parameter name="scale">
<definition> Scaling factor of the parameter field. A value &gt; 1 zooms into
the field; a value &lt; 1 zooms out. This may of course change the statistical
properties of the field (e.g. correlation lengths).
<parameter name="file">
<definition> YAML file containing the parameter definitions.
</definition>
<values> float &times; float (&times; float) </values>
<suggestion> 1 1 </suggestion>
</parameter>
<parameter name="offset">
<definition> The parameter field is shifted by this value (in physical units) </definition>
<values> float &times; float (&times; float) </values>
<suggestion> 0 0 </suggestion>
</parameter>
<parameter name="interpolation">
<definition> DEPRECATED: Sets the interpolation behavior when querying parameter field values
at a certain grid cell. </definition>
<values> nearest </values>
<suggestion> nearest </suggestion>
<values> path </values>
</parameter>
</category>
......
......@@ -132,21 +132,20 @@ int main(int argc, char** argv)
if (dim==2)
{
if (gtype == "gmsh"){
auto grid = Dune::Dorie::build_grid_gmsh<Dune::UGGrid<2>>(inifile,helper);
if(adaptivity){
switch(FEorder){
case 1:{
Sim<SimplexAdaptive<2,1>> sim(helper,grid,inifile);
Sim<SimplexAdaptive<2,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<SimplexAdaptive<2,2>> sim(helper,grid,inifile);
Sim<SimplexAdaptive<2,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<SimplexAdaptive<2,3>> sim(helper,grid,inifile);
Sim<SimplexAdaptive<2,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -157,17 +156,17 @@ int main(int argc, char** argv)
else{ // no adaptivity
switch(FEorder){
case 1:{
Sim<Simplex<2,1>> sim(helper,grid,inifile);
Sim<Simplex<2,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<Simplex<2,2>> sim(helper,grid,inifile);
Sim<Simplex<2,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<Simplex<2,3>> sim(helper,grid,inifile);
Sim<Simplex<2,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -178,20 +177,19 @@ int main(int argc, char** argv)
}
else if (gtype == "rectangular"){
if(adaptivity){
auto grid = Dune::Dorie::build_grid_cube<Dune::UGGrid<2>>(inifile,helper);
switch(FEorder){
case 1:{
Sim<CubeAdaptive<2,1>> sim(helper,grid,inifile);
Sim<CubeAdaptive<2,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<CubeAdaptive<2,2>> sim(helper,grid,inifile);
Sim<CubeAdaptive<2,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<CubeAdaptive<2,3>> sim(helper,grid,inifile);
Sim<CubeAdaptive<2,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -200,20 +198,19 @@ int main(int argc, char** argv)
}
}
else{ // no adaptivity
auto grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<2>>(inifile,helper);
switch(FEorder){
case 1:{
Sim<Cube<2,1>> sim(helper,grid,inifile);
Sim<Cube<2,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<Cube<2,2>> sim(helper,grid,inifile);
Sim<Cube<2,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<Cube<2,3>> sim(helper,grid,inifile);
Sim<Cube<2,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -229,21 +226,20 @@ int main(int argc, char** argv)
else if (dim==3)
{
if (gtype == "gmsh"){
auto grid = Dune::Dorie::build_grid_gmsh<Dune::UGGrid<3>>(inifile,helper);
if(adaptivity){
switch(FEorder){
case 1:{
Sim<SimplexAdaptive<3,1>> sim(helper,grid,inifile);
Sim<SimplexAdaptive<3,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<SimplexAdaptive<3,2>> sim(helper,grid,inifile);
Sim<SimplexAdaptive<3,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<SimplexAdaptive<3,3>> sim(helper,grid,inifile);
Sim<SimplexAdaptive<3,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -254,17 +250,17 @@ int main(int argc, char** argv)
else{ // no adaptivity
switch(FEorder){
case 1:{
Sim<Simplex<3,1>> sim(helper,grid,inifile);
Sim<Simplex<3,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<Simplex<3,2>> sim(helper,grid,inifile);
Sim<Simplex<3,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<Simplex<3,3>> sim(helper,grid,inifile);
Sim<Simplex<3,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -275,20 +271,19 @@ int main(int argc, char** argv)
}
else if (gtype == "rectangular"){
if(adaptivity){
auto grid = Dune::Dorie::build_grid_cube<Dune::UGGrid<3>>(inifile,helper);
switch(FEorder){
case 1:{
Sim<CubeAdaptive<3,1>> sim(helper,grid,inifile);
Sim<CubeAdaptive<3,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<CubeAdaptive<3,2>> sim(helper,grid,inifile);
Sim<CubeAdaptive<3,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<CubeAdaptive<3,3>> sim(helper,grid,inifile);
Sim<CubeAdaptive<3,3>> sim(inifile, helper);
sim.run();
break;
}
......@@ -297,20 +292,19 @@ int main(int argc, char** argv)
}
}
else{ // no adaptivity
auto grid = Dune::Dorie::build_grid_cube<Dune::YaspGrid<3>>(inifile,helper);
switch(FEorder){
case 1:{
Sim<Cube<3,1>> sim(helper,grid,inifile);
Sim<Cube<3,1>> sim(inifile, helper);
sim.run();
break;
}
case 2:{
Sim<Cube<3,2>> sim(helper,grid,inifile);
Sim<Cube<3,2>> sim(inifile, helper);
sim.run();
break;
}
case 3:{
Sim<Cube<3,3>> sim(helper,grid,inifile);
Sim<Cube<3,3>> sim(inifile, helper);
sim.run();
break;
}
......
add_executable(grid-test grid-test.cc)
dune_target_link_libraries(grid-test ${DUNE_LIBS})
\ No newline at end of file
# add_executable(grid-test grid-test.cc)
# dune_target_link_libraries(grid-test ${DUNE_LIBS})
\ No newline at end of file
#ifdef HAVE_CONFIG_H
#include "config.h"
#endif
#include <string>
#include <iostream>
#include <cstdio>
#include <thread>
#include <dune/common/exceptions.hh>
#include <dune/common/parametertree.hh>
#include <dune/common/parametertreeparser.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/uggrid.hh>
#include <dune/grid/io/file/vtk/vtkwriter.hh>
#include <dune/dorie/solver/util_grid_creator.hh>
#include <dune/dorie/solver/param_interface.hh>
#include <dune/dorie/interface/util.hh>
using SimplexTraits = Dune::Dorie::BaseTraits<Dune::UGGrid<2>,
Dune::GeometryType::BasicType::simplex,true,false>;
// using CubeTraits = Dune::Dorie::BaseTraits<Dune::YaspGrid<dim>,
// Dune::GeometryType::BasicType::cube,true,false>;
int main(int argc, char** argv)
{
try{
//Initialize Mpi
Dune::MPIHelper& helper = Dune::MPIHelper::instance(argc, argv);
if (argc!=2)
DUNE_THROW(Dune::IOError,"No parameter file specified!");
const std::string inifilename = argv[1];
// Read ini file
Dune::ParameterTree inifile;
Dune::ParameterTreeParser ptreeparser;
ptreeparser.readINITree(inifilename, inifile);
// maybe attach debugger
if (helper.size() > 1 && inifile.get<bool>("misc.debugMode", false))
{
int i = 0;
std::cout << "Attach now" << std::endl;
while (0 == i)
std::this_thread::sleep_for(std::chrono::milliseconds(100));
}
// build the grid creator
using GridType = Dune::UGGrid<2>;
Dune::Dorie::GridCreator<GridType> gc(inifile, helper);
auto [grid, mapper] = gc.get_grid_and_mapper();
auto element_map = mapper.get_element_index_map();
auto boundary_map = mapper.get_boundary_index_map();
std::cout << "Process " << helper.rank() << std::endl;
std::cout << "Element index map (" << element_map.size() << "): ";
for (auto&& value : element_map) {
std::cout << value << " ";
}
std::cout << std::endl;
// auto gv = grid->levelGridView(0);
// Dune::VTKWriter writer(gv);
// writer.addCellData(element_map, "mapping");
// writer.pwrite(
// inifile.get<std::string>("output.fileName"),
// inifile.get<std::string>("output.outputPath"),
// "./");
std::cout << "Adding FlowParameters" << std::endl;
Dune::Dorie::FlowParameters<SimplexTraits> fparam(inifile, grid,
element_map);
grid->globalRefine(1);
auto gv = grid->leafGridView();
for (auto&& elem : elements(gv)) {
fparam.bind(elem);
auto& cache = fparam.cache();
auto param = std::get<std::shared_ptr<Dune::Dorie::RichardsParameterization<SimplexTraits>>>(cache);
auto values = param->parameters();
assert(values.at("k0") == 1E-5 || values.at("k0") == 2.2E-5);
}
std::cout << "Finished!" << std::endl;
return 0;
}
catch (Dune::Exception &e){
std::cerr << "Dune reported error: " << e << std::endl;
return 1;
}
catch (...){
std::cerr << "Unknown exception thrown!" << std::endl;
throw;
return 1;
}
}
\ No newline at end of file
......@@ -4,9 +4,14 @@ namespace Dune{
namespace Dorie{
template<typename Traits>
RichardsSimulation<Traits>::RichardsSimulation (Dune::MPIHelper& _helper, std::shared_ptr<Grid> _grid, Dune::ParameterTree& _inifile) :
helper(_helper), grid(_grid), inifile(_inifile), gv(grid->leafGridView()),
mbe_slop(estimate_mbe_entries<typename MBE::size_type>(Traits::dim,Traits::GridGeometryType)),
RichardsSimulation<Traits>::RichardsSimulation (
Dune::ParameterTree& _inifile,
Dune::MPIHelper& _helper
):
helper(_helper),
inifile(_inifile),
mbe_slop(estimate_mbe_entries<typename MBE::size_type>(
Traits::dim,Traits::GridGeometryType)),
mbe_tlop(1),
verbose(inifile.get<int>("output.verbose")),
output_type(_inifile.get<bool>("output.asciiVtk") ? Dune::VTK::OutputType::ascii : Dune::VTK::OutputType::base64)
......@@ -14,26 +19,26 @@ RichardsSimulation<Traits>::RichardsSimulation (Dune::MPIHelper& _helper, std::s
Dune::Timer timer;
const int verbose_low = helper.rank() == 0 ? verbose : 0;
// create the grid
GridCreator<Grid> grid_creator(inifile, helper);
auto mapper = grid_creator.get_mapper();
grid = mapper.get_grid();
gv = std::make_shared<GV>(grid->leafGridView());
// --- Grid Function Space ---
gfs = std::make_shared<GFS>(GFSHelper::create(gv));
gfs = std::make_unique<GFS>(GFSHelper::create(*gv));
gfs->name("matric_head");
cc = std::make_unique<CC>();
Dune::PDELab::constraints(*gfs,*cc,false);
// --- Parameter class ---
InterpolatorFactory intfac(inifile,helper,verbose_low);
interp = intfac.create();
ParameterizationFactory pfac(inifile,helper,*interp,verbose_low);
param = pfac.create();
// --- Create the new parameter class
fparam = std::make_unique<FlowParameters>(inifile, grid);
fparam->copy_from_old_parameters(*param);
auto element_map = mapper.get_element_index_map();
fparam = std::make_shared<FlowParameters>(inifile, grid, element_map);
// --- Operator Helper Classes ---
fboundary = std::make_unique<FlowBoundary>(inifile);
fsource = std::make_unique<FlowSource>(inifile);
finitial = std::make_unique<FlowInitial>(inifile, gv);
finitial = std::make_unique<FlowInitial>(inifile, *gv);
// --- Local Operators ---
#ifdef EXPERIMENTAL_DG_FEATURES
......@@ -44,13 +49,13 @@ RichardsSimulation<Traits>::RichardsSimulation (Dune::MPIHelper& _helper, std::s
const auto upwinding = std::get<OP::RichardsDGUpwinding::Type>(settings);
const auto weights = std::get<OP::RichardsDGWeights::Type>(settings);
slop = std::make_unique<SLOP>(inifile, *fparam, gv, *fboundary, *fsource,
slop = std::make_unique<SLOP>(inifile, *fparam, *gv, *fboundary, *fsource,
method, upwinding, weights);
#else
slop = std::make_unique<SLOP>(inifile, *fparam, gv, *fboundary, *fsource);
slop = std::make_unique<SLOP>(inifile, *fparam, *gv, *fboundary, *fsource);
#endif // EXPERIMENTAL_DG_FEATURES
tlop = std::make_unique<TLOP>(inifile, *fparam, gv);
tlop = std::make_unique<TLOP>(inifile, *fparam, *gv);
controller = std::make_unique<CalculationController>(inifile,*fboundary,helper);
// --- Solution Vectors and Initial Condition ---
......@@ -63,10 +68,12 @@ RichardsSimulation<Traits>::RichardsSimulation (Dune::MPIHelper& _helper, std::s
// --- Utility Class Setup --- //
if constexpr (Traits::write_output)
{
const int subsamling_lvl = _inifile.get<int>("output.subsamplingLevel", 0);
const int subsamling_lvl =
_inifile.get<int>("output.subsamplingLevel", 0);
const auto subsamling_intervals = Dune::refinementLevels(subsamling_lvl);
auto subsamling_vtk = std::make_shared<Dune::SubsamplingVTKWriter<GV>>(gv,
subsamling_intervals);
auto subsamling_vtk =
std::make_shared<Dune::SubsamplingVTKWriter<GV>>(*gv,
subsamling_intervals);
vtkwriter = std::make_unique<Writer>(subsamling_vtk,
inifile.get<std::string>("output.fileName"),
......@@ -91,7 +98,7 @@ void RichardsSimulation<Traits>::operator_setup ()
igo = std::make_unique<IGO>(*go0,*go1);
// --- Solvers ---
lsgfs = std::make_unique<LSGFS>(LSGFSHelper::create(gv));
lsgfs = std::make_unique<LSGFS>(LSGFSHelper::create(*gv));
lscc = std::make_unique<LSCC>();
ls = std::make_unique<LS>(*igo,*cc,*lsgfs,*lscc,1000,0,true,true);
......@@ -121,7 +128,7 @@ void RichardsSimulation<Traits>::operator_setup ()
auto DOFnumber = gfs->globalSize();
if(verbose>1)
std::cout << " Process " << helper.rank() << ": number of DOF: " << DOFnumber << std::endl;
DOFnumber = gv.comm().sum(DOFnumber);
DOFnumber = gv->comm().sum(DOFnumber);
if (helper.rank()==0)
std::cout << " Total number of DOF: " << DOFnumber << std::endl;
}
......@@ -168,10 +175,6 @@ bool RichardsSimulation<Traits>::compute_time_step ()
template<typename Traits>
void RichardsSimulation<Traits>::run ()
{
// delete old paramterization storage
param.reset();
interp.reset();
if constexpr (Traits::write_output)
write_data(controller->getTime());
......@@ -185,7 +188,7 @@ void RichardsSimulation<Traits>::run ()
write_data(controller->getTime());
if(controller->doStep()
&& adaptivity->adapt_grid(*grid, gv, *gfs, *fparam, *fboundary,
&& adaptivity->adapt_grid(*grid, *gv, *gfs, *fparam, *fboundary,
time+dt, // need "old" boundary condition
*u))
{
......@@ -199,9 +202,9 @@ void RichardsSimulation<Traits>::update_adapters ()
{
udgf = std::make_shared<GFMatricHead>(gfs, u);
fluxdgf = std::make_shared<GFWaterFlux>(gfs, u, fparam);
condgf = std::make_shared<GFConductivity>(gv, fparam);
waterdgf = std::make_shared<GFWaterContent>(udgf, gv, fparam);
satdgf = std::make_shared<GFSaturation>(udgf, gv, fparam);
condgf = std::make_shared<GFConductivity>(*gv, fparam);
waterdgf = std::make_shared<GFWaterContent>(udgf, *gv, fparam);
satdgf = std::make_shared<GFSaturation>(udgf, *gv, fparam);
}
template<typename Traits>
......
......@@ -16,8 +16,6 @@
#include "util.hh"
#include "adaptivity.hh"
#include "../solver/util_interpolator.hh"
#include "../solver/param_base.hh"
#include "../solver/param_factory.hh"
#include "../solver/param_interface.hh"
#include "../solver/richards_initial.hh"
#include "../solver/richards_boundary.hh"
......@@ -25,6 +23,7 @@
#include "../solver/operator_DG.hh"
#include "../solver/util_controller.hh"
#include "../solver/util_writer.hh"
#include "../solver/util_grid_creator.hh"
#include "../solver/adapters/water_content.hh"
#include "../solver/adapters/saturation.hh"
#include "../solver/adapters/water_flux.hh"
......@@ -69,15 +68,7 @@ struct RichardsSimulationTraits : public BaseTraits
using LSCC = typename LSGFSHelper::CC;
// -- DORiE Class Definitions -- //
/// Parameter Interpolator
using Interpolator = Dune::Dorie::InterpolatorBase<BaseTraits,BaseTraits::dim>;
/// Factory class for Interpolator
using InterpolatorFactory = Dune::Dorie::InterpolatorFactory<BaseTraits,BaseTraits::dim>;
/// Class handling soil parameters
using Parameters = Dune::Dorie::ParameterizationBase<BaseTraits,Interpolator>;
/// Factory class for Parmeter class
using ParameterizationFactory = Dune::Dorie::ParameterizationFactory<BaseTraits,Interpolator>;
/// Class handling soil parameters
/// Class defining the soil parameterization
using FlowParameters = Dune::Dorie::FlowParameters<BaseTraits>;
/// Class defining boundary conditions
using FlowBoundary = Dune::Dorie::FlowBoundary<BaseTraits>;
......@@ -173,15 +164,7 @@ private:
using LSCC = typename Traits::LSCC;
// -- DORiE Class Definitions -- //
/// Parameter Interpolator
using Interpolator = typename Traits::Interpolator;
/// Factory class for Interpolator
using InterpolatorFactory = typename Traits::InterpolatorFactory;
/// Class handling soil parameters
using Parameters = typename Traits::Parameters;
/// Factory class for Parmeter class
using ParameterizationFactory = typename Traits::ParameterizationFactory;
/// Class handling soil parameters
/// Class defining the soil parameterization
using FlowParameters = typename Traits::FlowParameters;
/// Class defining boundary conditions
using FlowBoundary = typename Traits::FlowBoundary;
......@@ -241,17 +224,16 @@ private:
protected:
Dune::MPIHelper& helper;
std::shared_ptr<Grid> grid;
Dune::ParameterTree& inifile;
GV gv;
std::shared_ptr<GFS> gfs;
std::shared_ptr<Grid> grid;
std::shared_ptr<GV> gv;
std::unique_ptr<CC> cc;
std::unique_ptr<LSGFS> lsgfs;
std::unique_ptr<LSCC> lscc;
std::unique_ptr<Interpolator> interp;
std::unique_ptr<Parameters> param;
std::shared_ptr<FlowParameters> fparam;
std::unique_ptr<FlowBoundary> fboundary;
std::unique_ptr<FlowSource> fsource;
......@@ -297,10 +279,10 @@ public:
* simulation.
*
* @param _helper Dune MPI instance controller
* @param[in] _grid Shared pointer to the grid
* @param _inifile Dune parameter file parser
*/
RichardsSimulation (Dune::MPIHelper& _helper, std::shared_ptr<Grid> _grid, Dune::ParameterTree& _inifile);
RichardsSimulation (Dune::ParameterTree& _inifile,
Dune::MPIHelper& _helper);
/*------------------------------------------------------------------------*//**
* @brief Create the Grid Operators, Solvers, and Time Step Operator.
......
......@@ -65,7 +65,7 @@ public:
_config(config),
_verbose(_config.get<int>("output.verbose"))
{
const auto grid_type = _config.get<std::string>("grid.type");
const auto grid_type = _config.get<std::string>("grid.gridType");
if (grid_type == "rectangular")
_grid_mode = GridMode::regular;
else if (grid_type == "gmsh")
......
#ifndef TEST_GRID_CREATION_HH
#define TEST_GRID_CREATION_HH
#include <string>
#include <dune/common/exceptions.hh>
#include <dune/common/parametertree.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/grid/common/gridview.hh>
#include <dune/grid/io/file/vtk/vtkwriter.hh>
#include <dune/grid/common/mcmgmapper.hh>
#include "../solver/util_grid_creator.hh"
/// Create the grid and check if the mapping is as expected
template<typename Grid>
void create_grid_and_test_mapping (const Dune::ParameterTree& config,
const Dune::MPIHelper& helper,
const int count=0)
{
// create the grid mapper