3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
gridmanager_sp.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_IO_GRID_MANAGER_SP_HH
25#define DUMUX_IO_GRID_MANAGER_SP_HH
26
27// SPGrid specific includes
28#if HAVE_DUNE_SPGRID
29#include <dune/grid/spgrid.hh>
30#include <dune/grid/spgrid/dgfparser.hh>
31#endif
32
33#ifndef DUMUX_IO_GRID_MANAGER_BASE_HH
35#endif
36
37namespace Dumux {
38
39#if HAVE_DUNE_SPGRID
40
49template<class ct, int dim, template< int > class Ref, class Comm>
50class GridManager<Dune::SPGrid<ct, dim, Ref, Comm>>
51: public GridManagerBase<Dune::SPGrid<ct, dim, Ref, Comm>>
52{
53public:
55 using ParentType = GridManagerBase<Grid>;
56
60 void init(const std::string& paramGroup = "")
61 {
62 const auto overlap = getParamFromGroup<int>(paramGroup, "Grid.Overlap", 1);
63 if (overlap == 0)
64 DUNE_THROW(Dune::NotImplemented, "dune-spgrid does currently not support zero overlap!");
65
66 // try to create it from file
67 if (hasParamInGroup(paramGroup, "Grid.File"))
68 {
69 ParentType::makeGridFromDgfFile(getParamFromGroup<std::string>(paramGroup, "Grid.File"));
70 ParentType::maybeRefineGrid(paramGroup);
71 ParentType::loadBalance();
72 return;
73 }
74 // Didn't find a way to construct the grid
75 else if (hasParamInGroup(paramGroup, "Grid.UpperRight"))
76 {
77 using GlobalPosition = Dune::FieldVector<ct, dim>;
78 const auto lowerLeft = getParamFromGroup<GlobalPosition>(paramGroup, "Grid.LowerLeft", GlobalPosition(0.0));
79 const auto upperRight = getParamFromGroup<GlobalPosition>(paramGroup, "Grid.UpperRight");
80
81 using IntArray = std::array<int, dim>;
82 IntArray cells; cells.fill(1);
83 cells = getParamFromGroup<IntArray>(paramGroup, "Grid.Cells", cells);
84
85 const auto periodic = getParamFromGroup<std::bitset<dim>>(paramGroup, "Grid.Periodic", std::bitset<dim>{});
86 IntArray spOverlap; spOverlap.fill(overlap);
87
88 using Domain = typename Grid::Domain;
89 std::vector< typename Domain::Cube > cubes;
90 cubes.push_back( typename Domain::Cube( lowerLeft, upperRight ) );
91 Domain domain( cubes, typename Domain::Topology( static_cast<unsigned int>(periodic.to_ulong()) ) );
92 ParentType::gridPtr() = std::make_shared<Grid>( domain, cells, spOverlap );
93 ParentType::maybeRefineGrid(paramGroup);
94 ParentType::loadBalance();
95 }
96 else
97 {
98 const auto prefix = paramGroup.empty() ? paramGroup : paramGroup + ".";
99 DUNE_THROW(ParameterException, "Please supply a grid file in " << prefix << "Grid.File or " << prefix << "Grid.UpperRight/Cells.");
100 }
101 }
102};
103
104#endif // HAVE_DUNE_SPGRID
105
106} // end namespace Dumux
107
108#endif
Provides a grid manager for all supported grid managers with input file interfaces....
bool hasParamInGroup(const std::string &paramGroup, const std::string &param)
Check whether a key exists in the parameter tree with a model group prefix.
Definition: parameters.hh:177
Definition: adapt.hh:29
Definition: common/pdesolver.hh:36
Definition: momentum/velocitygradients.hh:38
Grid Grid
Definition: gridmanager_base.hh:69
void init(const std::string &modelParamGroup="")
Make the grid. Implement this method in the specialization of this class for a grid type.
Definition: gridmanager_base.hh:75