dune-grid-glue  2.7.0
contactmerge.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
8 #ifndef DUNE_GRIDGLUE_MERGING_CONTACTMERGE_HH
9 #define DUNE_GRIDGLUE_MERGING_CONTACTMERGE_HH
10 
11 
12 #include <iostream>
13 #include <fstream>
14 #include <iomanip>
15 #include <vector>
16 #include <algorithm>
17 #include <limits>
18 #include <memory>
19 #include <functional>
20 
21 #include <dune/common/fvector.hh>
22 #include <dune/common/exceptions.hh>
23 #include <dune/common/bitsetvector.hh>
24 #include <dune/common/deprecated.hh>
25 
26 #include <dune/grid/common/grid.hh>
27 
30 
31 namespace Dune {
32 namespace GridGlue {
33 
39 template<int dimworld, typename T = double>
41 : public StandardMerge<T,dimworld-1,dimworld-1,dimworld>
42 {
43  static constexpr int dim = dimworld-1;
44 
45  static_assert( dim==1 || dim==2,
46  "ContactMerge yet only handles the cases dim==1 and dim==2!");
47 
49 public:
50 
51  /* E X P O R T E D T Y P E S A N D C O N S T A N T S */
52 
54  typedef T ctype;
55 
57  typedef Dune::FieldVector<T, dimworld> WorldCoords;
58 
60  typedef Dune::FieldVector<T, dim> LocalCoords;
61 
71  ContactMerge(const T allowedOverlap=T(0),
72  std::function<WorldCoords(WorldCoords)> domainDirections=nullptr,
73  std::function<WorldCoords(WorldCoords)> targetDirections=nullptr,
75  : domainDirections_(domainDirections), targetDirections_(targetDirections),
76  overlap_(allowedOverlap), type_(type)
77  {}
78 
84  ContactMerge(const T allowedOverlap, ProjectionType type)
85  : overlap_(allowedOverlap),
86  type_(type)
87  {}
88 
97  inline
98  void setSurfaceDirections(std::function<WorldCoords(WorldCoords)> domainDirections,
99  std::function<WorldCoords(WorldCoords)> targetDirections)
100  {
101  domainDirections_ = domainDirections;
102  targetDirections_ = targetDirections;
103  this->valid = false;
104  }
105 
107  void setOverlap(T overlap)
108  {
109  overlap_ = overlap;
110  }
111 
113  T getOverlap() const
114  {
115  return overlap_;
116  }
117 
121  void minNormalAngle(T angle)
122  {
123  using std::cos;
124  maxNormalProduct_ = cos(angle);
125  }
126 
130  T minNormalAngle() const
131  {
132  using std::acos;
133  return acos(maxNormalProduct_);
134  }
135 
136 protected:
137  typedef typename StandardMerge<T,dimworld-1,dimworld-1,dimworld>::SimplicialIntersection SimplicialIntersection;
138 
139 private:
143  std::function<WorldCoords(WorldCoords)> domainDirections_;
144  std::vector<WorldCoords> nodalDomainDirections_;
145 
154  std::function<WorldCoords(WorldCoords)> targetDirections_;
155  std::vector<WorldCoords> nodalTargetDirections_;
156 
158  T overlap_;
159 
161  ProjectionType type_;
162 
166  T maxNormalProduct_ = T(-0.1);
167 
172  void computeIntersections(const Dune::GeometryType& grid1ElementType,
173  const std::vector<Dune::FieldVector<T,dimworld> >& grid1ElementCorners,
174  std::bitset<(1<<dim)>& neighborIntersects1,
175  unsigned int grid1Index,
176  const Dune::GeometryType& grid2ElementType,
177  const std::vector<Dune::FieldVector<T,dimworld> >& grid2ElementCorners,
178  std::bitset<(1<<dim)>& neighborIntersects2,
179  unsigned int grid2Index,
180  std::vector<SimplicialIntersection>& intersections) override;
181 
185 protected:
186  void build(const std::vector<Dune::FieldVector<T,dimworld> >& grid1Coords,
187  const std::vector<unsigned int>& grid1Elements,
188  const std::vector<Dune::GeometryType>& grid1ElementTypes,
189  const std::vector<Dune::FieldVector<T,dimworld> >& grid2Coords,
190  const std::vector<unsigned int>& grid2Elements,
191  const std::vector<Dune::GeometryType>& grid2ElementTypes) override
192  {
193  std::cout<<"ContactMerge building grid!\n";
194  // setup the nodal direction vectors
195  setupNodalDirections(grid1Coords, grid1Elements, grid1ElementTypes,
196  grid2Coords, grid2Elements, grid2ElementTypes);
197 
198  Base::build(grid1Coords, grid1Elements, grid1ElementTypes,
199  grid2Coords, grid2Elements, grid2ElementTypes);
200 
201  }
202 
203 private:
204 
206  static LocalCoords localCornerCoords(int i, const Dune::GeometryType& gt)
207  {
208  const auto& ref = Dune::ReferenceElements<T,dim>::general(gt);
209  return ref.position(i,dim);
210  }
211 
212 protected:
213 
215  void computeCyclicOrder(const std::vector<std::array<LocalCoords,2> >& polytopeCorners,
216  const LocalCoords& center, std::vector<int>& ordering) const;
217 
219  void setupNodalDirections(const std::vector<WorldCoords>& coords1,
220  const std::vector<unsigned int>& elements1,
221  const std::vector<Dune::GeometryType>& elementTypes1,
222  const std::vector<WorldCoords>& coords2,
223  const std::vector<unsigned int>& elements2,
224  const std::vector<Dune::GeometryType>& elementTypes2);
225 
227  void computeOuterNormalField(const std::vector<WorldCoords>& coords,
228  const std::vector<unsigned int>& elements,
229  const std::vector<Dune::GeometryType>& elementTypes,
230  std::vector<WorldCoords>& normals);
231 
233  void removeDoubles(std::vector<std::array<LocalCoords,2> >& polytopeCorners);
234 };
235 
236 } /* namespace GridGlue */
237 } /* namespace Dune */
238 
239 #include "contactmerge.cc"
240 
241 #endif // DUNE_GRIDGLUE_MERGING_CONTACTMERGE_HH
Dune
Definition: gridglue.hh:35
Dune::GridGlue::ContactMerge::minNormalAngle
void minNormalAngle(T angle)
set minimum angle in radians between normals at x and Φ(x)
Definition: contactmerge.hh:121
contactmerge.cc
Dune::GridGlue::ContactMerge::removeDoubles
void removeDoubles(std::vector< std::array< LocalCoords, 2 > > &polytopeCorners)
Remove all multiples.
Definition: contactmerge.cc:333
Dune::GridGlue::StandardMerge
Common base class for many merger implementations: produce pairs of entities that may intersect.
Definition: standardmerge.hh:54
Dune::GridGlue::ContactMerge
Merge two codimension-1 surfaces that may be a positive distance apart.
Definition: contactmerge.hh:40
Dune::GridGlue::ContactMerge::setupNodalDirections
void setupNodalDirections(const std::vector< WorldCoords > &coords1, const std::vector< unsigned int > &elements1, const std::vector< Dune::GeometryType > &elementTypes1, const std::vector< WorldCoords > &coords2, const std::vector< unsigned int > &elements2, const std::vector< Dune::GeometryType > &elementTypes2)
Setup the direction vectors containing the directions for each vertex.
Definition: contactmerge.cc:267
gridglue.hh
Central component of the module implementing the coupling of two grids.
Dune::GridGlue::ContactMerge::minNormalAngle
T minNormalAngle() const
get minimum angle in radians between normals at x and Φ(x)
Definition: contactmerge.hh:130
Dune::GridGlue::ContactMerge::setSurfaceDirections
void setSurfaceDirections(std::function< WorldCoords(WorldCoords)> domainDirections, std::function< WorldCoords(WorldCoords)> targetDirections)
Set surface direction functions.
Definition: contactmerge.hh:98
Dune::GridGlue::ContactMerge::WorldCoords
Dune::FieldVector< T, dimworld > WorldCoords
the coordinate type used in this interface
Definition: contactmerge.hh:57
Dune::GridGlue::ContactMerge::getOverlap
T getOverlap() const
Get the allowed overlap of the surfaces.
Definition: contactmerge.hh:113
Dune::GridGlue::ContactMerge::build
void build(const std::vector< Dune::FieldVector< T, dimworld > > &grid1Coords, const std::vector< unsigned int > &grid1Elements, const std::vector< Dune::GeometryType > &grid1ElementTypes, const std::vector< Dune::FieldVector< T, dimworld > > &grid2Coords, const std::vector< unsigned int > &grid2Elements, const std::vector< Dune::GeometryType > &grid2ElementTypes) override
builds the merged grid
Definition: contactmerge.hh:186
Dune::GridGlue::ContactMerge::CLOSEST_POINT
@ CLOSEST_POINT
Definition: contactmerge.hh:63
Dune::GridGlue::ContactMerge::LocalCoords
Dune::FieldVector< T, dim > LocalCoords
the coordinate type used in this interface
Definition: contactmerge.hh:60
Dune::GridGlue::ContactMerge::OUTER_NORMAL
@ OUTER_NORMAL
Definition: contactmerge.hh:63
Dune::GridGlue::ContactMerge::ProjectionType
ProjectionType
Type of the projection, closest point or outer normal projection.
Definition: contactmerge.hh:63
Dune::GridGlue::StandardMerge< double, dimworld-1, dimworld-1, dimworld >::valid
bool valid
Definition: standardmerge.hh:84
Dune::GridGlue::ContactMerge::computeOuterNormalField
void computeOuterNormalField(const std::vector< WorldCoords > &coords, const std::vector< unsigned int > &elements, const std::vector< Dune::GeometryType > &elementTypes, std::vector< WorldCoords > &normals)
If no direction field was specified compute the outer normal field.
Definition: contactmerge.cc:294
Dune::GridGlue::ContactMerge::setOverlap
void setOverlap(T overlap)
Set the allowed overlap of the surfaces.
Definition: contactmerge.hh:107
Dune::GridGlue::ContactMerge::ContactMerge
ContactMerge(const T allowedOverlap=T(0), std::function< WorldCoords(WorldCoords)> domainDirections=nullptr, std::function< WorldCoords(WorldCoords)> targetDirections=nullptr, ProjectionType type=OUTER_NORMAL)
Construct merger given overlap and possible projection directions.
Definition: contactmerge.hh:71
Dune::GridGlue::ContactMerge::SimplicialIntersection
StandardMerge< T, dimworld-1, dimworld-1, dimworld >::SimplicialIntersection SimplicialIntersection
Definition: contactmerge.hh:137
standardmerge.hh
Common base class for many merger implementations: produce pairs of entities that may intersect.
Dune::GridGlue::ContactMerge::computeCyclicOrder
void computeCyclicOrder(const std::vector< std::array< LocalCoords, 2 > > &polytopeCorners, const LocalCoords &center, std::vector< int > &ordering) const
Order the corners of the intersection polytope in cyclic order.
Definition: contactmerge.cc:212
Dune::GridGlue::ContactMerge::ctype
T ctype
the numeric type used in this interface
Definition: contactmerge.hh:54
Dune::GridGlue::StandardMerge::build
void build(const std::vector< Dune::FieldVector< T, dimworld > > &grid1_Coords, const std::vector< unsigned int > &grid1_elements, const std::vector< Dune::GeometryType > &grid1_element_types, const std::vector< Dune::FieldVector< T, dimworld > > &grid2_coords, const std::vector< unsigned int > &grid2_elements, const std::vector< Dune::GeometryType > &grid2_element_types) override
builds the merged grid
Definition: standardmerge.hh:390
Dune::GridGlue::ContactMerge::ContactMerge
ContactMerge(const T allowedOverlap, ProjectionType type)
Construct merger given overlap and type of the projection.
Definition: contactmerge.hh:84