8 #include <dune/common/filledarray.hh> 9 #include <dune/common/fvector.hh> 10 #include <dune/common/typeutilities.hh> 13 #include <dune/grid/albertagrid/albertareader.hh> 16 #include <dune/grid/io/file/gmshreader.hh> 17 #include <dune/grid/utility/structuredgridfactory.hh> 20 #include <dune/vtk/vtkreader.hh> 21 #include <dune/vtk/gridcreators/lagrangegridcreator.hh> 25 #include <dune/gmsh4/gmsh4reader.hh> 26 #include <dune/gmsh4/gridcreators/lagrangegridcreator.hh> 27 #include <dune/gmsh4/utility/version.hh> 30 #include <amdis/AdaptiveGrid.hpp> 31 #include <amdis/Initfile.hpp> 32 #include <amdis/Output.hpp> 33 #include <amdis/common/Concepts.hpp> 34 #include <amdis/common/Filesystem.hpp> 35 #include <amdis/common/TypeTraits.hpp> 36 #include <amdis/utility/MacroGridFactory.hpp> 41 template <
int dim,
class Coordinates>
class YaspGrid;
42 template <
class ct,
int dim>
class EquidistantCoordinates;
43 template <
class ct,
int dim>
class EquidistantOffsetCoordinates;
44 template <
class ct,
int dim>
class TensorProductCoordinates;
54 enum { dimension = G::dimension };
55 enum { dimworld = G::dimensionworld };
58 using HostGrid =
typename Grid::HostGrid;
60 using ctype =
typename Grid::ctype;
71 static std::shared_ptr<Grid>
create(std::string name)
90 auto filename = Parameters::get<std::string>(name_ +
"->macro file name");
91 auto structured = Parameters::get<std::string>(name_ +
"->structured");
92 std::unique_ptr<HostGrid> gridPtr;
95 gridPtr = create_unstructured_grid(filename.value());
96 }
else if (structured) {
98 if (structured.value() ==
"cube") {
99 gridPtr = create_cube_grid();
100 }
else if (structured && structured.value() ==
"simplex") {
101 gridPtr = create_simplex_grid();
103 error_exit(
"Unknown structured grid type. Must be either `cube` or `simplex` in parameter [meshName]->structured.");
107 gridPtr = create_by_gridtype<HostGrid>(Dune::PriorityTag<42>{});
111 auto globalRefinements = Parameters::get<int>(name_ +
"->global refinements");
112 if (globalRefinements)
113 gridPtr->globalRefine(globalRefinements.value());
114 auto loadBalance = Parameters::get<bool>(name_ +
"->load balance");
115 if (loadBalance && loadBalance.value())
116 gridPtr->loadBalance();
118 return construct(std::move(gridPtr));
124 return create_structured_grid([](
auto&& lower,
auto&& upper,
auto&& numCells)
126 return Dune::StructuredGridFactory<HostGrid>::createCubeGrid(lower, upper, numCells);
133 return create_structured_grid([](
auto&& lower,
auto&& upper,
auto&& numCells)
135 return Dune::StructuredGridFactory<HostGrid>::createSimplexGrid(lower, upper, numCells);
152 static std::unique_ptr<Grid> construct(std::unique_ptr<Grid> hostGrid)
154 return std::move(hostGrid);
158 static std::unique_ptr<Grid> construct(std::unique_ptr<HG> hostGrid)
160 return std::make_unique<Grid>(std::move(hostGrid));
164 template <
class Size =
unsigned int,
class Factory>
165 std::unique_ptr<HostGrid> create_structured_grid(Factory factory)
const 168 Dune::FieldVector<ctype,int(dimworld)> lower(0);
170 Dune::FieldVector<ctype,int(dimworld)> upper(1);
172 auto numCells = Dune::filledArray<std::size_t(dimension),Size>(1);
174 Parameters::get(name_ +
"->min corner", lower);
175 Parameters::get(name_ +
"->max corner", upper);
176 Parameters::get(name_ +
"->num cells", numCells);
177 return factory(lower, upper, numCells);
181 std::unique_ptr<HostGrid> create_unstructured_grid(std::string
const& filename)
const 188 if (Dune::Gmsh4::fileVersion(filename)[0] >= 4)
189 return Dune::Gmsh4Reader<HostGrid, Dune::Gmsh4::LagrangeGridCreator<HostGrid>>::createGridFromFile(filename);
192 return read_gmsh_file<HostGrid>(filename, Dune::PriorityTag<42>{});
195 else if (ext ==
".vtu") {
196 return Dune::VtkReader<HostGrid, Dune::Vtk::LagrangeGridCreator<HostGrid>>::createGridFromFile(filename);
199 else if (ext ==
".1d" || ext ==
".2d" || ext ==
".3d" || ext ==
".amc") {
200 return read_alberta_file<HostGrid>(filename, Dune::PriorityTag<42>{});
203 error_exit(
"Cannot read grid file. Unknown file extension.");
208 template <
class Gr
idType>
209 using SupportsGmshReader = decltype(std::declval<Dune::GridFactory<GridType>>().insertBoundarySegment(
210 std::declval<std::vector<unsigned int>>(),
211 std::declval<std::shared_ptr<Dune::BoundarySegment<GridType::dimension, GridType::dimensionworld> >>()) );
213 template <
class Gr
idType,
class LC>
214 using SupportsInsertionIndex = decltype(std::declval<Dune::GridFactory<GridType>>().insertionIndex(std::declval<LC>()));
217 template <
class GridType = HostGrid,
218 REQUIRES(Dune::Std::is_detected<SupportsGmshReader, GridType>::value)>
219 std::unique_ptr<GridType> read_gmsh_file(std::string
const& filename, Dune::PriorityTag<1>)
const 221 Dune::GmshReader<GridType> reader;
222 Dune::GridFactory<GridType> factory;
223 std::vector<int> boundaryIds, elementIds;
224 reader.read(factory, filename, boundaryIds, elementIds);
226 using HasInsertionIndexEntity
227 = Dune::Std::is_detected<SupportsInsertionIndex, GridType, typename GridType::template Codim<0>::Entity>;
228 using HasInsertionIndexIntersection
229 = Dune::Std::is_detected<SupportsInsertionIndex, GridType, typename GridType::LeafIntersection>;
231 auto gridPtr = factory.createGrid();
232 if ((boundaryIds.empty() && elementIds.empty()) ||
233 (!HasInsertionIndexEntity::value && !HasInsertionIndexIntersection::value))
234 return std::unique_ptr<GridType>(std::move(gridPtr));
237 if (!boundaryIds.empty() && HasInsertionIndexIntersection::value)
238 boundaryIds_.resize(gridPtr->numBoundarySegments());
239 if (!elementIds.empty() && HasInsertionIndexEntity::value)
240 elementIds_.resize(gridPtr->size(0));
242 auto const& indexSet = gridPtr->leafIndexSet();
243 for (
auto const& e : elements(gridPtr->leafGridView())) {
244 if constexpr(HasInsertionIndexEntity::value) {
245 if (!elementIds.empty())
246 elementIds_[indexSet.index(e)] = elementIds[factory.insertionIndex(e)];
249 if (boundaryIds.empty() || !e.hasBoundaryIntersections())
252 for (
auto const& it : intersections(gridPtr->leafGridView(), e)) {
253 if constexpr(HasInsertionIndexIntersection::value) {
255 boundaryIds_[it.boundarySegmentIndex()]
256 = boundaryIds[factory.insertionIndex(it)];
261 return std::unique_ptr<GridType>(std::move(gridPtr));
265 template <
class Gr
idType = HostGr
id>
266 std::unique_ptr<GridType> read_gmsh_file(std::string
const&, Dune::PriorityTag<0>)
const 268 error_exit(
"Gmsh reader not supported for this grid type!");
275 template <
class Gr
idType>
276 using IsAlbertaGrid = decltype(std::declval<GridType>().alberta2dune(0,0));
279 template <
class GridType = HostGrid,
280 REQUIRES(GridType::dimensionworld == DIM_OF_WORLD),
281 REQUIRES(Dune::Std::is_detected<IsAlbertaGrid, GridType>::value)>
282 std::unique_ptr<GridType> read_alberta_file(std::string
const& filename, Dune::PriorityTag<3>)
const 284 return std::make_unique<GridType>(filename);
288 template <
class GridType = HostGrid,
289 REQUIRES(GridType::dimensionworld == DIM_OF_WORLD)>
290 std::unique_ptr<GridType> read_alberta_file(std::string
const& filename, Dune::PriorityTag<2>)
const 292 Dune::GridFactory<GridType> factory;
293 if (Environment::mpiRank() == 0) {
294 Dune::AlbertaReader<GridType> reader;
295 reader.readGrid(filename, factory);
297 return std::unique_ptr<GridType>{factory.createGrid()};
301 template <
class GridType = HostGrid,
302 REQUIRES(GridType::dimensionworld != DIM_OF_WORLD)>
303 std::unique_ptr<GridType> read_alberta_file(std::string
const& filename, Dune::PriorityTag<1>)
const 305 error_exit(
"AlbertaGrid (and AlbertaReader) require WORLDDIM == Grid::dimensionworld. Change the cmake parameters!");
310 template <
class Gr
idType>
311 std::unique_ptr<GridType> read_alberta_file(std::string
const&, Dune::PriorityTag<0>)
const 313 error_exit(
"AlbertaGrid (and AlbertaReader) not available. Set AlbertaFlags to your executable in cmake!");
321 template <
class GridType = HostGrid,
322 REQUIRES(Dune::Std::is_detected<IsAlbertaGrid, GridType>::value)>
323 std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<3>)
const 325 return create_structured_grid([](
auto&& lower,
auto&& upper,
auto&& numCells)
333 template <
class GridType = HostGrid,
334 class =
typename GridType::YGrid>
335 std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<2>)
const 338 Parameters::get(name_ +
"->overlap", overlap);
339 std::string periodic(dimension,
'0');
340 Parameters::get(name_ +
"->periodic", periodic);
342 return create_yaspgrid(
Types<GridType>{}, overlap, std::bitset<dimension>(periodic));
345 template <
int dim,
class ct>
346 std::unique_ptr<HostGrid> create_yaspgrid(
Types<Dune::YaspGrid<dim,Dune::EquidistantCoordinates<ct,dim>>>,
347 int overlap, std::bitset<dimension>
const& per)
const 349 return create_structured_grid<int>([&](
auto&& ,
auto&& upper, std::array<int,dimension>
const& numCells)
351 return std::make_unique<HostGrid>(upper, numCells, per, overlap);
355 template <
int dim,
class ct>
356 std::unique_ptr<HostGrid> create_yaspgrid(
Types<Dune::YaspGrid<dim,Dune::EquidistantOffsetCoordinates<ct,dim>>>,
357 int overlap, std::bitset<dimension>
const& per)
const 359 return create_structured_grid<int>([&](
auto&& lower,
auto&& upper, std::array<int,dimension>
const& numCells)
361 return std::make_unique<HostGrid>(lower, upper, numCells, per, overlap);
365 template <
int dim,
class ct>
366 std::unique_ptr<HostGrid> create_yaspgrid(
Types<Dune::YaspGrid<dim,Dune::TensorProductCoordinates<ct,dim>>>,
367 int, std::bitset<dimension>
const&)
const 369 error_exit(
"MeshCreator cannot create YaspGrid with TensorProductCoordinates.");
376 template <
class GridType = HostGrid,
377 class =
typename GridType::ReferenceCube,
379 std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<1>)
const 381 return create_structured_grid<int>([](
auto&& lower,
auto&& upper, std::array<int,dimension>
const& numCells)
384 return std::make_unique<GridType>(lower, upper, multiIndex);
390 template <
class Gr
idType = HostGr
id>
391 std::unique_ptr<GridType> create_by_gridtype(Dune::PriorityTag<0>)
const 393 error_exit(
"Don't know how to create the grid.");
399 mutable std::vector<int> boundaryIds_;
400 mutable std::vector<int> elementIds_;
constexpr bool MultiIndex
A multi-index type.
Definition: Concepts.hpp:150
path extension() const
Returns the file extension path component.
Definition: Filesystem.cpp:83
Definition: AdaptiveGrid.hpp:373
A variadic type list.
Definition: TypeTraits.hpp:88
void error_exit(std::string const &str, Args &&... args)
print a message and exit
Definition: Output.hpp:142
typename Impl::AdaptiveGridImpl< HostGrid >::type AdaptiveGrid_t
Definition: AdaptiveGrid.hpp:367
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
std::unique_ptr< HostGrid > create_simplex_grid() const
Create a structured simplex grid.
Definition: MeshCreator.hpp:131
std::vector< int > const & elementIds() const
Return the filled vector of element ids. NOTE: not all creators support reading this.
Definition: MeshCreator.hpp:146
A creator class for dune grids.
Definition: MeshCreator.hpp:52
static std::shared_ptr< Grid > create(std::string name)
Static create mthod. See create()
Definition: MeshCreator.hpp:71
std::unique_ptr< HostGrid > create_cube_grid() const
Create a structured cube grid.
Definition: MeshCreator.hpp:122
Definition: Filesystem.hpp:11
MeshCreator(std::string const &name)
Construct a new MeshCreator.
Definition: MeshCreator.hpp:66
std::unique_ptr< Grid > create() const
Create a new grid by inspecting the intifile parameter group [meshName]
Definition: MeshCreator.hpp:88
std::vector< int > const & boundaryIds() const
Return the filled vector of boundary ids. NOTE: not all creators support reading this.
Definition: MeshCreator.hpp:140
Definition: MacroGridFactory.hpp:105