Skip to content
Snippets Groups Projects
Commit cd9137e6 authored by graeser's avatar graeser
Browse files

Extract reusable functionality to header

parent 319d3275
Branches
No related tags found
No related merge requests found
......@@ -9,7 +9,6 @@
// included standard library headers
#include <iostream>
#include <array>
// included dune-common headers
#include <dune/common/parallel/mpihelper.hh>
......@@ -17,64 +16,12 @@
#include <dune/common/fvector.hh>
#include <dune/common/stringutility.hh>
// included dune-geometry headers
#include <dune/geometry/quadraturerules.hh>
// included dune-grid headers
#include <dune/grid/io/file/vtk/vtkwriter.hh>
#include <dune/grid/utility/structuredgridfactory.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/uggrid.hh>
// included dune-fu-tutorial headers
template<class Grid>
auto createCubeGrid()
{
static const int dim = Grid::dimension;
auto lowerLeft = Dune::FieldVector<double,dim>(0);
auto upperRight = Dune::FieldVector<double,dim>(1);
auto elementsPerDirection = std::array<unsigned int,dim>();
for(auto& e : elementsPerDirection)
e = 2;
return Dune::StructuredGridFactory<Grid>::createCubeGrid(lowerLeft, upperRight, elementsPerDirection);
}
template<class Grid>
auto createSimplexGrid()
{
static const int dim = Grid::dimension;
auto lowerLeft = Dune::FieldVector<double,dim>(0);
auto upperRight = Dune::FieldVector<double,dim>(1);
auto elementsPerDirection = std::array<unsigned int,dim>();
for(auto& e : elementsPerDirection)
e = 2;
return Dune::StructuredGridFactory<Grid>::createSimplexGrid(lowerLeft, upperRight, elementsPerDirection);
}
template<class GridView>
void writeGridView(const GridView& gridView, std::string postFix)
{
Dune::VTKWriter<GridView> vtkWriter(gridView);
vtkWriter.write(std::string("04-gridviews-")+postFix);
}
template<class Grid>
void writeAllGridViews(const Grid& grid, std::string gridName)
{
for(int level = 0; level <= grid.maxLevel(); ++level)
writeGridView(grid.levelGridView(level), Dune::formatString(gridName+"-level-%02d", level));
writeGridView(grid.leafGridView(), gridName+"-leaf");
}
#include "04-gridviews.hh"
......
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_FUTUTORIAL_SRC_04GRIDVIEWS_HH
#define DUNE_FUTUTORIAL_SRC_04GRIDVIEWS_HH
// included standard library headers
#include <string>
#include <array>
// included dune-common headers
#include <dune/common/fvector.hh>
#include <dune/common/stringutility.hh>
// included dune-grid headers
#include <dune/grid/io/file/vtk/vtkwriter.hh>
#include <dune/grid/utility/structuredgridfactory.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/uggrid.hh>
template<class Grid>
auto createCubeGrid()
{
static const int dim = Grid::dimension;
auto lowerLeft = Dune::FieldVector<double,dim>(0);
auto upperRight = Dune::FieldVector<double,dim>(1);
auto elementsPerDirection = std::array<unsigned int,dim>();
for(auto& e : elementsPerDirection)
e = 2;
return Dune::StructuredGridFactory<Grid>::createCubeGrid(lowerLeft, upperRight, elementsPerDirection);
}
template<class Grid>
auto createSimplexGrid()
{
static const int dim = Grid::dimension;
auto lowerLeft = Dune::FieldVector<double,dim>(0);
auto upperRight = Dune::FieldVector<double,dim>(1);
auto elementsPerDirection = std::array<unsigned int,dim>();
for(auto& e : elementsPerDirection)
e = 2;
return Dune::StructuredGridFactory<Grid>::createSimplexGrid(lowerLeft, upperRight, elementsPerDirection);
}
template<class GridView>
void writeGridView(const GridView& gridView, std::string postFix)
{
Dune::VTKWriter<GridView> vtkWriter(gridView);
vtkWriter.write(std::string("04-gridviews-")+postFix);
}
template<class Grid>
void writeAllGridViews(const Grid& grid, std::string gridName)
{
for(int level = 0; level <= grid.maxLevel(); ++level)
writeGridView(grid.levelGridView(level), Dune::formatString(gridName+"-level-%02d", level));
writeGridView(grid.leafGridView(), gridName+"-leaf");
}
#endif // DUNE_FUTUTORIAL_SRC_04GRIDVIEWS_HH
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment