8 #ifndef DUNE_GRIDGLUE_MERGING_CONTACTMERGE_HH 9 #define DUNE_GRIDGLUE_MERGING_CONTACTMERGE_HH 21 #include <dune/common/fvector.hh> 22 #include <dune/common/function.hh> 23 #include <dune/common/exceptions.hh> 24 #include <dune/common/bitsetvector.hh> 25 #include <dune/common/deprecated.hh> 27 #include <dune/grid/common/grid.hh> 40 template<
int dimworld,
typename T =
double>
44 enum {dim = dimworld-1};
46 static_assert( dim==1 || dim==2,
47 "ContactMerge yet only handles the cases dim==1 and dim==2!");
70 DUNE_DEPRECATED_MSG(
"Please use a std::function<FieldVector(FieldVector)> to prescribe non-default projections")
72 const
Dune::VirtualFunction<WorldCoords,WorldCoords>* domainDirections,
73 const
Dune::VirtualFunction<WorldCoords,WorldCoords>* targetDirections)
74 : overlap_(allowedOverlap)
76 if (domainDirections) {
77 domainDirections_ = [domainDirections](
const WorldCoords& in) {
79 domainDirections->evaluate(in,out);
83 if (targetDirections) {
84 targetDirections_ = [targetDirections](
const WorldCoords& in) {
86 targetDirections->evaluate(in,out);
100 std::function<
WorldCoords(WorldCoords)> domainDirections=
nullptr,
101 std::function<
WorldCoords(WorldCoords)> targetDirections=
nullptr)
102 : domainDirections_(domainDirections), targetDirections_(targetDirections),
103 overlap_(allowedOverlap)
116 std::function<
WorldCoords(WorldCoords)> targetDirections)
118 domainDirections_ = domainDirections;
119 targetDirections_ = targetDirections;
131 DUNE_DEPRECATED_MSG(
"Please use a std::function<FieldVector(FieldVector)> to prescribe non-default projections")
134 const Dune::VirtualFunction<WorldCoords,WorldCoords>* targetDirections)
136 domainDirections_ = [domainDirections](
const WorldCoords& in) {
138 domainDirections->evaluate(in,out);
142 targetDirections_ = [targetDirections](
const WorldCoords& in) {
144 targetDirections->evaluate(in,out);
168 maxNormalProduct_ = cos(angle);
177 return acos(maxNormalProduct_);
187 std::function<WorldCoords(WorldCoords)> domainDirections_;
188 std::vector<WorldCoords> nodalDomainDirections_;
198 std::function<WorldCoords(WorldCoords)> targetDirections_;
199 std::vector<WorldCoords> nodalTargetDirections_;
207 T maxNormalProduct_ = T(-0.1);
213 void computeIntersections(
const Dune::GeometryType& grid1ElementType,
214 const std::vector<Dune::FieldVector<T,dimworld> >& grid1ElementCorners,
215 std::bitset<(1<<dim)>& neighborIntersects1,
216 unsigned int grid1Index,
217 const Dune::GeometryType& grid2ElementType,
218 const std::vector<Dune::FieldVector<T,dimworld> >& grid2ElementCorners,
219 std::bitset<(1<<dim)>& neighborIntersects2,
220 unsigned int grid2Index,
227 void build(
const std::vector<Dune::FieldVector<T,dimworld> >& grid1Coords,
228 const std::vector<unsigned int>& grid1Elements,
229 const std::vector<Dune::GeometryType>& grid1ElementTypes,
230 const std::vector<Dune::FieldVector<T,dimworld> >& grid2Coords,
231 const std::vector<unsigned int>& grid2Elements,
232 const std::vector<Dune::GeometryType>& grid2ElementTypes)
234 std::cout<<
"ContactMerge building grid!\n";
237 grid2Coords, grid2Elements, grid2ElementTypes);
239 Base::build(grid1Coords, grid1Elements, grid1ElementTypes,
240 grid2Coords, grid2Elements, grid2ElementTypes);
247 static LocalCoords localCornerCoords(
int i,
const Dune::GeometryType& gt)
249 const Dune::ReferenceElement<T,dim>& ref = Dune::ReferenceElements<T,dim>::general(gt);
250 return ref.position(i,dim);
256 void computeCyclicOrder(
const std::vector<std::array<LocalCoords,2> >& polytopeCorners,
257 const LocalCoords& center, std::vector<int>& ordering)
const;
261 const std::vector<unsigned int>& elements1,
262 const std::vector<Dune::GeometryType>& elementTypes1,
263 const std::vector<WorldCoords>& coords2,
264 const std::vector<unsigned int>& elements2,
265 const std::vector<Dune::GeometryType>& elementTypes2);
269 const std::vector<unsigned int>& elements,
270 const std::vector<Dune::GeometryType>& elementTypes,
271 std::vector<WorldCoords>& normals);
274 void removeDoubles(std::vector<std::array<LocalCoords,2> >& polytopeCorners);
282 #endif // DUNE_GRIDGLUE_MERGING_CONTACTMERGE_HH Common base class for many merger implementations: produce pairs of entities that may intersect...
Definition: standardmerge.hh:53
void setSurfaceDirections(std::function< WorldCoords(WorldCoords)> domainDirections, std::function< WorldCoords(WorldCoords)> targetDirections)
Set surface direction functions.
Definition: contactmerge.hh:115
StandardMerge< T, dimworld-1, dimworld-1, dimworld >::RemoteSimplicialIntersection RemoteSimplicialIntersection
Definition: contactmerge.hh:181
void removeDoubles(std::vector< std::array< LocalCoords, 2 > > &polytopeCorners)
Remove all multiples.
Definition: contactmerge.cc:326
Dune::FieldVector< T, dimworld > WorldCoords
the coordinate type used in this interface
Definition: contactmerge.hh:58
Merge two codimension-1 surfaces that may be a positive distance apart.
Definition: contactmerge.hh:41
IteratorRange<... > intersections(const GridGlue<... > &glue, const Reverse<... > &reverse=!reversed)
Iterate over all intersections of a GridGlue.
void setOverlap(T overlap)
Set the allowed overlap of the surfaces.
Definition: contactmerge.hh:151
T getOverlap() const
Get the allowed overlap of the surfaces.
Definition: contactmerge.hh:157
bool valid
Definition: standardmerge.hh:75
Dune::FieldVector< T, dim > LocalCoords
the coordinate type used in this interface
Definition: contactmerge.hh:61
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)
builds the merged grid
Definition: contactmerge.hh:227
void setSurfaceDirections(const Dune::VirtualFunction< WorldCoords, WorldCoords > *domainDirections, const Dune::VirtualFunction< WorldCoords, WorldCoords > *targetDirections)
Set surface direction functions.
Definition: contactmerge.hh:133
void computeCyclicOrder(const std::vector< std::array< LocalCoords, 2 > > &polytopeCorners, const LocalCoords ¢er, std::vector< int > &ordering) const
Order the corners of the intersection polytope in cyclic order.
Definition: contactmerge.cc:205
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:287
Common base class for many merger implementations: produce pairs of entities that may intersect...
virtual 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)
builds the merged grid
Definition: standardmerge.hh:465
void minNormalAngle(T angle)
set minimum angle in radians between normals at x and Φ(x)
Definition: contactmerge.hh:165
T ctype
the numeric type used in this interface
Definition: contactmerge.hh:55
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:260
T minNormalAngle() const
get minimum angle in radians between normals at x and Φ(x)
Definition: contactmerge.hh:174
Definition: gridglue.hh:30
Central component of the module implementing the coupling of two grids.
ContactMerge(const T allowedOverlap=T(0), std::function< WorldCoords(WorldCoords)> domainDirections=nullptr, std::function< WorldCoords(WorldCoords)> targetDirections=nullptr)
Construct merger given overlap and possible projection directions.
Definition: contactmerge.hh:99