From edd91679b193bc5d0622325b16ea39a792c46b96 Mon Sep 17 00:00:00 2001 From: Elias Pipping <elias.pipping@fu-berlin.de> Date: Sat, 14 Dec 2013 13:36:58 +0100 Subject: [PATCH] [Cleanup] Construct boundary from faces --- src/one-body-sample.cc | 62 ++++++++++++++++++++++++------------------ 1 file changed, 36 insertions(+), 26 deletions(-) diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc index 504e5e7a..2d3a781a 100644 --- a/src/one-body-sample.cc +++ b/src/one-body-sample.cc @@ -137,42 +137,52 @@ int main(int argc, char *argv[]) { GridView const leafView = grid->leafView(); // }}} - // Set up the boundary - Dune::BitSetVector<dims> dirichletNodes(fineVertexCount); - Dune::BitSetVector<dims> noNodes(fineVertexCount); - BoundaryPatch<GridView> const neumannBoundary(leafView); + // Set up faces + BoundaryPatch<GridView> lowerFace(leafView); + BoundaryPatch<GridView> upperFace(leafView); + { + auto const isClose = [](double a, + double b) { return std::abs(a - b) < 1e-14; }; + lowerFace.insertFacesByProperty([&]( + typename GridView::Intersection const &in) { + return isClose(lowerLeft[1], in.geometry().center()[1]); + }); + + upperFace.insertFacesByProperty([&]( + typename GridView::Intersection const &in) { + return isClose(upperRight[1], in.geometry().center()[1]); + }); + } Vector vertexCoordinates(fineVertexCount); - Dune::BitSetVector<1> frictionalNodes(fineVertexCount, false); { Dune::MultipleCodimMultipleGeomTypeMapper< GridView, Dune::MCMGVertexLayout> const vertexMapper(leafView); for (auto it = leafView.begin<dims>(); it != leafView.end<dims>(); ++it) { - assert(it->geometry().corners() == 1); - size_t const id = vertexMapper.map(*it); - vertexCoordinates[id] = it->geometry().corner(0); - auto const &localCoordinates = vertexCoordinates[id]; - - // lower face - if (localCoordinates[1] == lowerLeft[1]) { - frictionalNodes[id] = true; - velocityDirichletNodes[id][1] = true; - } + auto const geometry = it->geometry(); + assert(geometry.corners() == 1); + vertexCoordinates[vertexMapper.map(*it)] = geometry.corner(0); + } + } - // upper face - else if (localCoordinates[1] == upperRight[1]) - velocityDirichletNodes[id] = true; + // Neumann boundary + BoundaryPatch<GridView> const neumannBoundary(leafView); - // right face (except for both corners) - else if (localCoordinates[0] == upperRight[0]) - ; + // Frictional Boundary + BoundaryPatch<GridView> const &frictionalBoundary = lowerFace; + Dune::BitSetVector<1> frictionalNodes(fineVertexCount); + frictionalBoundary.getVertices(frictionalNodes); - // left face (except for both corners) - else if (localCoordinates[0] == lowerLeft[0]) - ; - } + // Dirichlet Boundary + Dune::BitSetVector<dims> noNodes(fineVertexCount); + Dune::BitSetVector<dims> dirichletNodes(fineVertexCount); + for (size_t i = 0; i < fineVertexCount; ++i) { + if (lowerFace.containsVertex(i)) + dirichletNodes[i][1] = true; + + if (upperFace.containsVertex(i)) + dirichletNodes[i] = true; } - BoundaryPatch<GridView> const frictionalBoundary(leafView, frictionalNodes); // Set up functions for time-dependent boundary conditions using Function = Dune::VirtualFunction<double, double>; -- GitLab