diff --git a/test/discretization/staggered/CMakeLists.txt b/test/discretization/staggered/CMakeLists.txt index 6f026f178477d077ab4794ad95fd444ba29a28d5..ca77cf7288376bf9f6b0554b8d5beca0b3403115 100644 --- a/test/discretization/staggered/CMakeLists.txt +++ b/test/discretization/staggered/CMakeLists.txt @@ -1,12 +1,11 @@ dune_add_test(NAME test_staggeredfvgeometry - SOURCES test_staggeredfvgeometry.cc - COMPILE_DEFINITIONS ENABLE_CACHING=false) + SOURCES test_staggeredfvgeometry.cc) -dune_add_test(NAME test_staggeredfvgeometry_caching - SOURCES test_staggeredfvgeometry.cc - COMPILE_DEFINITIONS ENABLE_CACHING=true) +dune_add_test(NAME test_staggered_free_flow_geometry + SOURCES test_staggered_free_flow_geometry.cc) #install sources install(FILES test_staggeredfvgeometry.cc +test_staggered_free_flow_geometry.cc DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/test/discretization/staggered) diff --git a/test/discretization/staggered/test_staggered_free_flow_geometry.cc b/test/discretization/staggered/test_staggered_free_flow_geometry.cc new file mode 100644 index 0000000000000000000000000000000000000000..101b088e7f39a895d84c720e43ac4462425d1499 --- /dev/null +++ b/test/discretization/staggered/test_staggered_free_flow_geometry.cc @@ -0,0 +1,164 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +/***************************************************************************** + * See the file COPYING for full copying permissions. * + * * + * This program is free software: you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation, either version 2 of the License, or * + * (at your option) any later version. * + * * + * This program is distributed in the hope that it will be useful, * + * but WITHOUT ANY WARRANTY; without even the implied warranty of * + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * + * GNU General Public License for more details. * + * * + * You should have received a copy of the GNU General Public License * + * along with this program. If not, see <http://www.gnu.org/licenses/>. * + *****************************************************************************/ +/*! + * \file + * \brief Test for finite volume element geometry, sub control volume, and sub + control volume faces + */ +#include <config.h> + +#include <iostream> +#include <utility> +#include <iomanip> + +#include <dune/common/test/iteratortest.hh> +#include <dune/grid/utility/structuredgridfactory.hh> +#include <dune/grid/yaspgrid.hh> +#include <dune/grid/common/mcmgmapper.hh> + +#include <dumux/common/properties.hh> +#include <dumux/common/properties/model.hh> +#include <dumux/discretization/staggered/properties.hh> +#include <dumux/discretization/staggered/fvgridgeometry.hh> +#include <dumux/discretization/staggered/fvelementgeometry.hh> + +#include <dumux/discretization/staggered/properties.hh> +#include <dumux/freeflow/staggered/properties.hh> + +namespace Dumux +{ + +template<class TypeTag> +class MockProblem +{ }; + + +namespace Properties +{ +NEW_TYPE_TAG(TestStaggeredFreeFlowGeometry, INHERITS_FROM(StaggeredModel, NavierStokes)); + +SET_TYPE_PROP(TestStaggeredFreeFlowGeometry, Grid, Dune::YaspGrid<2>); + +// Set the problem property +SET_TYPE_PROP(TestStaggeredFreeFlowGeometry, Problem, MockProblem<TypeTag> ); + +SET_BOOL_PROP(TestStaggeredFreeFlowGeometry, EnableFVGridGeometryCache, true); +} + +} + +template<class T> +class NoopFunctor { +public: + NoopFunctor() {} + void operator()(const T& t){} +}; + +int main (int argc, char *argv[]) try +{ + // maybe initialize mpi + Dune::MPIHelper::instance(argc, argv); + + std::cout << "Checking the FVGeometries, SCVs and SCV faces" << std::endl; + + // aliases + using TypeTag = TTAG(TestStaggeredFreeFlowGeometry); + using Grid = typename GET_PROP_TYPE(TypeTag, Grid); + using GridView = typename Grid::LeafGridView; + + constexpr int dim = GridView::dimension; + constexpr int dimworld = GridView::dimensionworld; + + using GlobalPosition = typename Dune::FieldVector<Grid::ctype, dimworld>; + using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume); + using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace); + using FVGridGeometry = typename GET_PROP_TYPE(TypeTag, FVGridGeometry); + + // make a grid + GlobalPosition lower(0.0); + GlobalPosition upper(1.0); + std::array<unsigned int, dim> els{{2, 4}}; + std::shared_ptr<Grid> grid = Dune::StructuredGridFactory<Grid>::createCubeGrid(lower, upper, els); + auto leafGridView = grid->leafGridView(); + + FVGridGeometry fvGridGeometry(leafGridView); + fvGridGeometry.update(); + + std::cout << "Abbreviatons:\n" + << "pos - postition of face center\n" + << "fIdx - face index\n" + << "dofIdx (self/oppo.) - dofIdx on face (self/opposite)\n" + << "norm in/out - dofIdx on face normal to own face (within own element / in adjacent element)\n" + << "par - dofIdx on face parallel to own face\n" + << "normDist - distance bewteen the dofs on the faces normal the own face\n" + << "parDist - distance bewteen the dof on the parallel face and the one on the own face\n" + << "norm in/out - dofIdx on face normal to own face (within own element / in adjacent element)" << std::endl; + + // iterate over elements. For every element get fv geometry and loop over scvs and scvfaces + for (const auto& element : elements(leafGridView)) + { + auto eIdx = fvGridGeometry.elementMapper().index(element); + std::cout << std::endl << "Checking fvGeometry of element " << eIdx << std::endl; + auto fvGeometry = localView(fvGridGeometry); + fvGeometry.bind(element); + + auto range = scvs(fvGeometry); + NoopFunctor<SubControlVolume> op; + if(0 != testForwardIterator(range.begin(), range.end(), op)) + DUNE_THROW(Dune::Exception, "Iterator does not fulfill the forward iterator concept"); + + for (auto&& scv : scvs(fvGeometry)) + { + std::cout << "-- scv " << scv.dofIndex() << " center at: " << scv.center() << std::endl; + } + + auto range2 = scvfs(fvGeometry); + NoopFunctor<SubControlVolumeFace> op2; + if(0 != testForwardIterator(range2.begin(), range2.end(), op2)) + DUNE_THROW(Dune::Exception, "Iterator does not fulfill the forward iterator concept"); + + + for (auto&& scvf : scvfs(fvGeometry)) + { + std::cout << std::fixed << std::left << std::setprecision(2) + << "pos "<< scvf.ipGlobal() + << "; fIdx " << std::setw(3) << scvf.index() + << "; dofIdx (self/oppo.) " << std::setw(3) << scvf.dofIndex() << "/" << std::setw(3) <<scvf.dofIndexOpposingFace() + << "; dist (self/oppo.) " << std::setw(3) << scvf.selfToOppositeDistance() + << ", norm1 in/out " << std::setw(3) << scvf.pairData(0).normalPair.first << "/" << std::setw(3) << scvf.pairData(0).normalPair.second + << ", norm2 in/out " << std::setw(3) << scvf.pairData(1).normalPair.first << "/" << std::setw(3) << scvf.pairData(1).normalPair.second + << ", par1 " << std::setw(3) << std::setw(3) << scvf.pairData(0).outerParallelFaceDofIdx + << ", par2 " << std::setw(3) << std::setw(3) << scvf.pairData(1).outerParallelFaceDofIdx + << ", normDist1 " << std::setw(3) << scvf.pairData(0).normalDistance + << ", normDist2 " << std::setw(3) << scvf.pairData(1).normalDistance + << ", parDist1 " << std::setw(3) << scvf.pairData(0).parallelDistance + << ", parDist2 " << std::setw(3) << scvf.pairData(1).parallelDistance; + if (scvf.boundary()) std::cout << " (on boundary)"; + std::cout << std::endl; + } + } +} +// ////////////////////////////////// +// Error handler +// ///////////////////////////////// +catch (Dune::Exception e) { + + std::cout << e << std::endl; + return 1; +} diff --git a/test/discretization/staggered/test_staggeredfvgeometry.cc b/test/discretization/staggered/test_staggeredfvgeometry.cc index 48504ff483eca461496598f1117f02c59d5db671..255cced8efda4385c49d7f2a44028f348502c113 100644 --- a/test/discretization/staggered/test_staggeredfvgeometry.cc +++ b/test/discretization/staggered/test_staggeredfvgeometry.cc @@ -32,7 +32,8 @@ #include <dune/grid/yaspgrid.hh> #include <dune/grid/common/mcmgmapper.hh> -#include <dumux/common/basicproperties.hh> +#include <dumux/common/properties.hh> +#include <dumux/common/properties/model.hh> #include <dumux/discretization/staggered/properties.hh> #include <dumux/discretization/staggered/fvgridgeometry.hh> #include <dumux/discretization/staggered/fvelementgeometry.hh> @@ -40,13 +41,48 @@ namespace Dumux { +//! Dummy flux variables class so that we can update the connectivity map +class MockFluxVariables +{ +public: + + template<class Map, class Element, class FvGeometry, class Scvf> + void computeCellCenterToCellCenterStencil(Map& map, + const Element& element, + const FvGeometry& fvGeometry, + const Scvf& scvf) + {} + + template<class Map, class Element, class FvGeometry, class Scvf> + void computeCellCenterToFaceStencil(Map& map, + const Element& element, + const FvGeometry& fvGeometry, + const Scvf& scvf) + {} + + template<class Map, class FvGeometry, class Scvf> + void computeFaceToCellCenterStencil(Map& map, + const FvGeometry& fvGeometry, + const Scvf& scvf) + {} + + template<class Map, class FvGeometry, class Scvf> + void computeFaceToFaceStencil(Map& map, + const FvGeometry& fvGeometry, + const Scvf& scvf) + {} + +}; + namespace Properties { -NEW_TYPE_TAG(TestFVGeometry, INHERITS_FROM(StaggeredModel, NumericModel)); +NEW_TYPE_TAG(TestFVGeometry, INHERITS_FROM(StaggeredModel, ModelProperties)); SET_TYPE_PROP(TestFVGeometry, Grid, Dune::YaspGrid<2>); -SET_BOOL_PROP(TestFVGeometry, EnableFVGridGeometryCache, ENABLE_CACHING); +SET_TYPE_PROP(TestFVGeometry, FluxVariables, MockFluxVariables); + +SET_BOOL_PROP(TestFVGeometry, EnableFVGridGeometryCache, true); } } @@ -88,12 +124,6 @@ int main (int argc, char *argv[]) try FVGridGeometry fvGridGeometry(leafGridView); fvGridGeometry.update(); - std::cout << "Abbreviatons:\n" - << "ip - postition of face center\n" - << "face - face index\n" - << "self/oppo - dofIdx on intersection (self/opposite)\n" - << "norm in/out - dofIdx on side normal to intersection (within own element / in adjacent element)" << std::endl; - // iterate over elements. For every element get fv geometry and loop over scvs and scvfaces for (const auto& element : elements(leafGridView)) { @@ -104,36 +134,23 @@ int main (int argc, char *argv[]) try auto range = scvs(fvGeometry); NoopFunctor<SubControlVolume> op; - if(0 != testForwardIterator(range.begin(), range.end(), [](const SubControlVolumeFace& scvf){})) + if(0 != testForwardIterator(range.begin(), range.end(), op)) DUNE_THROW(Dune::Exception, "Iterator does not fulfill the forward iterator concept"); for (auto&& scv : scvs(fvGeometry)) { - std::cout << "-- scv " << scv.dofIndex() << " center at: " << scv.center() << std::endl; + std::cout << "-- scv " << scv.indexInElement() << " center at: " << scv.center() << " , volume: " << scv.volume() << std::endl; } auto range2 = scvfs(fvGeometry); NoopFunctor<SubControlVolumeFace> op2; - if(0 != testForwardIterator(range2.begin(), range2.end(), [](const SubControlVolumeFace& scvf){})) + if(0 != testForwardIterator(range2.begin(), range2.end(), op2)) DUNE_THROW(Dune::Exception, "Iterator does not fulfill the forward iterator concept"); - for (auto&& scvf : scvfs(fvGeometry)) { - std::cout << std::fixed << std::left << std::setprecision(2) - << "ip "<< scvf.ipGlobal() - << "; face " << std::setw(3) << scvf.index() - << "; self/oppo " << std::setw(3) << scvf.dofIndex() << "/" << std::setw(3) <<scvf.dofIndexOpposingFace() - << "; dist self/oppo " << std::setw(3) << scvf.selfToOppositeDistance() - << ", norm1 in/out " << std::setw(3) << scvf.pairData(0).normalPair.first << "/" << std::setw(3) << scvf.pairData(0).normalPair.second - << ", norm2 in/out " << std::setw(3) << scvf.pairData(1).normalPair.first << "/" << std::setw(3) << scvf.pairData(1).normalPair.second - << ", par1 in/out " << std::setw(3) << scvf.dofIndex() << "/" << std::setw(3) << scvf.pairData(0).outerParallelFaceDofIdx - << ", par2 in/out " << std::setw(3) << scvf.dofIndex() << "/" << std::setw(3) << scvf.pairData(1).outerParallelFaceDofIdx - << ", normDist1 " << std::setw(3) << scvf.pairData(0).normalDistance - << ", normDist2 " << std::setw(3) << scvf.pairData(1).normalDistance - << ", parDist1 " << std::setw(3) << scvf.pairData(0).parallelDistance - << ", parDist2 " << std::setw(3) << scvf.pairData(1).parallelDistance; - if (scvf.boundary()) std::cout << " (on boundary)"; + std::cout << "-- scvf " << scvf.index() << " ip at: " << scvf.ipGlobal() << " normal: " << scvf.unitOuterNormal(); + if (scvf.boundary()) std::cout << " (on boundary)."; std::cout << std::endl; } }