19 #ifndef DUNE_GRIDGLUE_MERGING_CONFORMINGMERGE_HH
20 #define DUNE_GRIDGLUE_MERGING_CONFORMINGMERGE_HH
27 #include <dune/common/fmatrix.hh>
28 #include <dune/common/fvector.hh>
30 #include <dune/geometry/referenceelements.hh>
44 template<
int dim,
int dimworld,
typename T =
double>
75 void computeIntersections(
const Dune::GeometryType& grid1ElementType,
76 const std::vector<Dune::FieldVector<T,dimworld> >& grid1ElementCorners,
77 std::bitset<(1<<dim)>& neighborIntersects1,
78 unsigned int grid1Index,
79 const Dune::GeometryType& grid2ElementType,
80 const std::vector<Dune::FieldVector<T,dimworld> >& grid2ElementCorners,
81 std::bitset<(1<<dim)>& neighborIntersects2,
82 unsigned int grid2Index,
83 std::vector<SimplicialIntersection>& intersections);
94 template<
int dim,
int dimworld,
typename T>
97 template<
int dim,
int dimworld,
typename T>
98 void ConformingMerge<dim, dimworld, T>::computeIntersections(
const Dune::GeometryType& grid1ElementType,
99 const std::vector<Dune::FieldVector<T,dimworld> >& grid1ElementCorners,
100 std::bitset<(1<<dim)>& neighborIntersects1,
101 unsigned int grid1Index,
102 const Dune::GeometryType& grid2ElementType,
103 const std::vector<Dune::FieldVector<T,dimworld> >& grid2ElementCorners,
104 std::bitset<(1<<dim)>& neighborIntersects2,
105 unsigned int grid2Index,
111 assert((
unsigned int)(Dune::ReferenceElements<T,dim>::general(grid1ElementType).size(dim)) == grid1ElementCorners.size());
112 assert((
unsigned int)(Dune::ReferenceElements<T,dim>::general(grid2ElementType).size(dim)) == grid2ElementCorners.size());
114 neighborIntersects1.reset();
115 neighborIntersects2.reset();
118 if (grid1ElementType != grid2ElementType)
124 std::vector<int> other(grid1ElementCorners.size(), -1);
126 for (
unsigned int i=0; i<grid1ElementCorners.size(); i++) {
128 for (
unsigned int j=0; j<grid2ElementCorners.size(); j++) {
130 if ( (grid1ElementCorners[i]-grid2ElementCorners[j]).two_norm() < tolerance_ ) {
149 const auto& refElement = Dune::ReferenceElements<T,dim>::general(grid1ElementType);
152 if (grid1ElementType.isSimplex()) {
156 for (
int i=0; i<refElement.size(dim); i++) {
157 intersections.back().corners0[0][i] = refElement.position(i,dim);
158 intersections.back().corners1[0][i] = refElement.position(other[i],dim);
161 }
else if (dim == 2 && grid1ElementType.isQuadrilateral()) {
164 const unsigned int subVertices[2][3] = {{0,1,3}, {0,3,2}};
166 for (
int i=0; i<2; i++) {
168 SimplicialIntersection newSimplicialIntersection(grid1Index, grid2Index);
170 for (
int j=0; j<dim+1; j++) {
171 newSimplicialIntersection.corners0[0][j] = refElement.position(subVertices[i][j],dim);
172 newSimplicialIntersection.corners1[0][j] = refElement.position(subVertices[i][other[j]],dim);
179 }
else if (grid1ElementType.isHexahedron()) {
183 const unsigned int subVertices[5][4] = {{0,1,3,5}, {0,3,2,6}, {4,5,0,6}, {6,7,6,3}, {6,0,5,3}};
185 for (
int i=0; i<5; i++) {
187 SimplicialIntersection newSimplicialIntersection(grid1Index, grid2Index);
189 for (
int j=0; j<dim+1; j++) {
190 newSimplicialIntersection.corners0[0][j] = refElement.position(subVertices[i][j],dim);
191 newSimplicialIntersection.corners1[0][j] = refElement.position(subVertices[i][other[j]],dim);
199 DUNE_THROW(Dune::GridError,
"Unsupported element type");
207 #endif // DUNE_GRIDGLUE_MERGING_CONFORMINGMERGE_HH