Commit 75731942 authored by Dennis Gläser's avatar Dennis Gläser Committed by Timo Koch
Browse files

[test][multidomain][facet] Add fracture grid creator test

parent 370dc816
add_subdirectory(embedded)
add_subdirectory(poromechanics)
add_subdirectory(boundary)
add_subdirectory(facet)
dune_symlink_to_source_files(FILES "grid.msh" "test_gridcreator.input")
dune_add_test(NAME test_facetgridcreator_alu
CMAKE_GUARD "( dune-foamgrid_FOUND AND dune-alugrid_FOUND )"
COMPILE_DEFINITIONS BULKGRIDTYPE=Dune::ALUGrid<3,3,Dune::simplex,Dune::conforming>
SOURCES test_gridcreator.cc
COMMAND ./test_facetgridcreator_alu
CMD_ARGS test_gridcreator.input)
dune_add_test(NAME test_facetgridcreator_ug
CMAKE_GUARD "( dune-foamgrid_FOUND AND dune-uggrid_FOUND )"
COMPILE_DEFINITIONS BULKGRIDTYPE=Dune::UGGrid<3>
SOURCES test_gridcreator.cc
COMMAND ./test_facetgridcreator_ug
CMD_ARGS test_gridcreator.input)
#install sources
install(FILES
test_gridcreator.cc
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dumux/test/multidomain/facet)
This diff is collapsed.
// -*- 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 Tests the grid creator class for models using facet coupling.
*/
#include <config.h>
#include <iostream>
#include <dune/common/exceptions.hh>
#include <dune/common/float_cmp.hh>
#include <dune/grid/uggrid.hh>
#include <dune/alugrid/grid.hh>
#include <dune/foamgrid/foamgrid.hh>
#include <dune/grid/io/file/vtk/vtkwriter.hh>
#include <dumux/common/parameters.hh>
#include <dumux/multidomain/facet/gridcreator.hh>
#ifndef BULKGRIDTYPE // default to ug grid if not provided by CMake
#define BULKGRIDTYPE Dune::UGGrid<3>
#endif
int main (int argc, char *argv[]) try
{
// maybe initialize mpi
Dune::MPIHelper::instance(argc, argv);
// parse command line argument parameters
Dumux::Parameters::init(argc, argv);
using BulkGrid = BULKGRIDTYPE;
using FacetGrid = Dune::FoamGrid<2, 3>;
using EdgeGrid = Dune::FoamGrid<1, 3>;
using GridCreator = Dumux::FacetCouplingGridCreator<BulkGrid, FacetGrid, EdgeGrid>;
GridCreator gridCreator;
gridCreator.makeGrids("grid.msh");
// check grid sizes
const auto& bulkGridView = gridCreator.grid<0>().leafGridView();
const auto& facetGridView = gridCreator.grid<1>().leafGridView();
const auto& edgeGridView = gridCreator.grid<2>().leafGridView();
if (bulkGridView.size(0) != 507)
DUNE_THROW(Dune::InvalidStateException, "Bulk grid has " << bulkGridView.size(0) << " instead of 507 elements");
if (bulkGridView.size(3) != 153)
DUNE_THROW(Dune::InvalidStateException, "Bulk grid has " << bulkGridView.size(3) << " instead of 153 vertices");
if (facetGridView.size(0) != 32)
DUNE_THROW(Dune::InvalidStateException, "Facet grid has " << facetGridView.size(0) << " instead of 32 elements");
if (facetGridView.size(2) != 23)
DUNE_THROW(Dune::InvalidStateException, "Facet grid has " << facetGridView.size(2) << " instead of 23 vertices");
if (edgeGridView.size(0) != 2)
DUNE_THROW(Dune::InvalidStateException, "Edge grid has " << edgeGridView.size(0) << " instead of 2 elements");
if (edgeGridView.size(1) != 3)
DUNE_THROW(Dune::InvalidStateException, "Edge grid has " << edgeGridView.size(1) << " instead of 3 vertices");
// all entities should have the domain index 10
for (const auto& e : elements(bulkGridView))
{
const auto domainMarker = gridCreator.getElementDomainMarker<0>(e);
if (domainMarker != 10)
DUNE_THROW(Dune::InvalidStateException, "Bulk element has domain marker " << domainMarker << " instead of 10");
}
for (const auto& e : elements(facetGridView))
{
const auto domainMarker = gridCreator.getElementDomainMarker<1>(e);
if (domainMarker != 10)
DUNE_THROW(Dune::InvalidStateException, "Facet element has domain marker " << domainMarker << " instead of 10");
}
for (const auto& e : elements(edgeGridView))
{
const auto domainMarker = gridCreator.getElementDomainMarker<2>(e);
if (domainMarker != 10)
DUNE_THROW(Dune::InvalidStateException, "Edge element has domain marker " << domainMarker << " instead of 10");
}
// Check boundary boundary segments. There should be mesh
// file-defined segments on all sides except the positive x-direction
std::size_t bSegmentCount = 0;
std::size_t negXCount, posYCount, negYCount, posZCount, negZCount;
negXCount = posYCount = negYCount = posZCount = negZCount = 0;
for (const auto& e : elements(bulkGridView))
{
for (const auto& is : intersections(bulkGridView, e))
{
if (is.boundary())
{
const auto n = is.centerUnitOuterNormal();
if (gridCreator.gridFactory<0>().wasInserted(is))
{
bSegmentCount++;
const auto marker = gridCreator.getBoundaryDomainMarker<0>(is);
if (Dune::FloatCmp::eq(n[0], -1.0, 1e-6)) // neg x-dir
{
if (marker != 4)
DUNE_THROW(Dune::InvalidStateException, "negative x-face should have marker 4, but has " << marker);
negXCount++;
}
else if (Dune::FloatCmp::eq(n[1], 1.0, 1e-6)) // pos y-dir
{
if (marker != 5)
DUNE_THROW(Dune::InvalidStateException, "positive y-face should have marker 5, but has " << marker);
posYCount++;
}
else if (Dune::FloatCmp::eq(n[1], -1.0, 1e-6)) // neg y-dir
{
if (marker != 6)
DUNE_THROW(Dune::InvalidStateException, "negative y-face should have marker 6, but has " << marker);
negYCount++;
}
else if (Dune::FloatCmp::eq(n[2], 1.0, 1e-6)) // pos z-dir
{
if (marker != 1)
DUNE_THROW(Dune::InvalidStateException, "positive z-face should have marker 1, but has " << marker);
posZCount++;
}
else if (Dune::FloatCmp::eq(n[2], -1.0, 1e-6)) // neg z-dir
{
if (marker != 2)
DUNE_THROW(Dune::InvalidStateException, "negative z-face should have marker 2, but has " << marker);
negZCount++;
}
}
else
{
if (!Dune::FloatCmp::eq(n[0], 1.0, 1e-6)) // pos x-dir
DUNE_THROW(Dune::InvalidStateException, "Found a boundary segment not specified in .msh file that is not on pos-x face");
}
}
}
}
//! make sure we found the right number of boundary segments
if (bSegmentCount != 200) DUNE_THROW(Dune::InvalidStateException, "Found " << bSegmentCount << " instead of 200 boundary segments");
if (negXCount != 40) DUNE_THROW(Dune::InvalidStateException, "Found " << negXCount << " instead of 40 boundary segments in neg x-direction");
if (posYCount != 40) DUNE_THROW(Dune::InvalidStateException, "Found " << posYCount << " instead of 40 boundary segments in pos y-direction");
if (negYCount != 40) DUNE_THROW(Dune::InvalidStateException, "Found " << negYCount << " instead of 40 boundary segments in neg y-direction");
if (posZCount != 40) DUNE_THROW(Dune::InvalidStateException, "Found " << posZCount << " instead of 40 boundary segments in pos z-direction");
if (negZCount != 40) DUNE_THROW(Dune::InvalidStateException, "Found " << negZCount << " instead of 40 boundary segments in neg z-direction");
//! check if we found the right number of embeddings
const auto& edgeEmbeddings = gridCreator.embeddedEntityMap(2);
const auto& edgeEmbedments = gridCreator.embedmentMap(2);
if (edgeEmbeddings.size() != 0) DUNE_THROW(Dune::InvalidStateException, "The grid with lowest dimension can't have embedded entities");
if (edgeEmbedments.size() != 2) DUNE_THROW(Dune::InvalidStateException, "Found " << edgeEmbedments.size() << " instead of 2 edge element embedments");
for (const auto& embedments : edgeEmbedments)
if (embedments.second.size() != 4)
DUNE_THROW(Dune::InvalidStateException, "edge element is embedded in " << embedments.second.size() << " facet elements instead of 4");
const auto& facetEmbeddings = gridCreator.embeddedEntityMap(1);
const auto& facetEmbedments = gridCreator.embedmentMap(1);
if (facetEmbeddings.size() != 8) DUNE_THROW(Dune::InvalidStateException, "Found " << facetEmbeddings.size() << " instead of 8 embeddings in facet grid");
if (facetEmbedments.size() != 32) DUNE_THROW(Dune::InvalidStateException, "Found " << facetEmbedments.size() << " instead of 32 facet element embedments");
for (const auto& embedments : facetEmbedments)
if (embedments.second.size() != 2)
DUNE_THROW(Dune::InvalidStateException, "facet element is embedded in " << embedments.second.size() << " bulk elements instead of 2");
for (const auto& embeddings : facetEmbeddings)
if (embeddings.second.size() != 1)
DUNE_THROW(Dune::InvalidStateException, "facet element has " << embeddings.second.size() << " embedded entities instead of 1");
const auto& bulkEmbeddings = gridCreator.embeddedEntityMap(0);
const auto& bulkEmbedments = gridCreator.embedmentMap(0);
if (bulkEmbeddings.size() != 56) DUNE_THROW(Dune::InvalidStateException, "Found " << bulkEmbeddings.size() << " instead of 56 embeddings in bulk grid");
if (bulkEmbedments.size() != 0) DUNE_THROW(Dune::InvalidStateException, "The grid with highest dimension can't have embedments");
std::size_t singleEmbeddings = 0;
std::size_t doubleEmbeddings = 0;
for (const auto& embeddings : bulkEmbeddings)
{
if (embeddings.second.size() != 1 && embeddings.second.size() != 2)
DUNE_THROW(Dune::InvalidStateException, "bulk element has " << embeddings.second.size() << " embedded entities instead of 1 or 2");
if (embeddings.second.size() == 1) singleEmbeddings++;
if (embeddings.second.size() == 2) doubleEmbeddings++;
}
if (singleEmbeddings != 48) DUNE_THROW(Dune::InvalidStateException, "Found " << singleEmbeddings << " instead of 48 bulk elements with 1 embedding");
if (doubleEmbeddings != 8) DUNE_THROW(Dune::InvalidStateException, "Found " << doubleEmbeddings << " instead of 8 bulk elements with 2 embeddings");
// test access to embeddings for single elements
singleEmbeddings = 0;
doubleEmbeddings = 0;
for (const auto& e : elements(bulkGridView))
{
if (gridCreator.embeddedEntityIndices<0>(e).size() == 1)
singleEmbeddings++;
else if (gridCreator.embeddedEntityIndices<0>(e).size() == 2)
doubleEmbeddings++;
else if (gridCreator.embeddedEntityIndices<0>(e).size() != 0)
DUNE_THROW(Dune::InvalidStateException, "wrong number of embeddings!");
if (gridCreator.embedmentEntityIndices<0>(e).size() != 0)
DUNE_THROW(Dune::InvalidStateException, "bulk grid can't be embedded anywhere!");
}
if (singleEmbeddings != 48) DUNE_THROW(Dune::InvalidStateException, "Found " << singleEmbeddings << " instead of 48 bulk elements with 1 embedding");
if (doubleEmbeddings != 8) DUNE_THROW(Dune::InvalidStateException, "Found " << doubleEmbeddings << " instead of 8 bulk elements with 2 embeddings");
std::size_t embeddings = 0;
std::size_t embedments = 0;
for (const auto& e : elements(facetGridView))
{
if (gridCreator.embeddedEntityIndices<1>(e).size() == 1)
embeddings++;
else if (gridCreator.embeddedEntityIndices<1>(e).size() != 0)
DUNE_THROW(Dune::InvalidStateException, "wrong number of embeddings!");
if (gridCreator.embedmentEntityIndices<1>(e).size() == 2)
embedments++;
else if (gridCreator.embedmentEntityIndices<1>(e).size() != 0)
DUNE_THROW(Dune::InvalidStateException, "wrong number of embedments!");
}
if (embeddings != 8) DUNE_THROW(Dune::InvalidStateException, "Found " << embeddings << " instead of 8 embeddings in facet grid");
if (embedments != 32) DUNE_THROW(Dune::InvalidStateException, "Found " << embedments << " instead of 32 facet element embedments");
embeddings = 0;
embedments = 0;
for (const auto& e : elements(edgeGridView))
{
if (gridCreator.embeddedEntityIndices<2>(e).size() != 0)
DUNE_THROW(Dune::InvalidStateException, "wrong number of embeddings!");
if (gridCreator.embedmentEntityIndices<2>(e).size() == 4)
embedments++;
else if (gridCreator.embedmentEntityIndices<2>(e).size() != 0)
DUNE_THROW(Dune::InvalidStateException, "wrong number of embedments!");
}
if (embeddings != 0) DUNE_THROW(Dune::InvalidStateException, "The grid with lowest dimension can't have embedded entities");
if (embedments!= 2) DUNE_THROW(Dune::InvalidStateException, "Found " << edgeEmbedments.size() << " instead of 2 edge element embedments");
// write .vtk file for each grid
using BulkWriter = Dune::VTKWriter<typename BulkGrid::LeafGridView>;
BulkWriter bulkWriter(bulkGridView); bulkWriter.write("bulkgrid");
using FacetWriter = Dune::VTKWriter<typename FacetGrid::LeafGridView>;
FacetWriter facetWriter(facetGridView); facetWriter.write("facetgrid");
using EdgeWriter = Dune::VTKWriter<typename EdgeGrid::LeafGridView>;
EdgeWriter edgeWriter(edgeGridView); edgeWriter.write("edgegrid");
}
// //////////////////////////////////
// Error handler
// /////////////////////////////////
catch (Dune::Exception e) {
std::cout << e << std::endl;
return 1;
}
[Grid]
File = ./grid.msh
Verbosity = true
DomainMarkers = true
BoundarySegments = true
GmshPhysicalEntityThreshold = 10
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment