Intrepid2
Intrepid2_CellToolsDefInclusion.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Intrepid2 Package
4//
5// Copyright 2007 NTESS and the Intrepid2 contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
10
16#ifndef __INTREPID2_CELLTOOLS_DEF_INCLUSION_HPP__
17#define __INTREPID2_CELLTOOLS_DEF_INCLUSION_HPP__
18
19// disable clang warnings
20#if defined (__clang__) && !defined (__INTEL_COMPILER)
21#pragma clang system_header
22#endif
23
24namespace Intrepid2 {
25
26 //============================================================================================//
27 // //
28 // Inclusion tests //
29 // //
30 //============================================================================================//
31
32
33 template<typename DeviceType>
34 template<typename PointViewType>
35 bool
36 CellTools<DeviceType>::
37 checkPointInclusion( const PointViewType point,
38 const shards::CellTopology cellTopo,
39 const typename ScalarTraits<typename PointViewType::value_type>::scalar_type threshold) {
40#ifdef HAVE_INTREPID2_DEBUG
41 INTREPID2_TEST_FOR_EXCEPTION( point.rank() != 1, std::invalid_argument,
42 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): Point must have rank 1. ");
43 INTREPID2_TEST_FOR_EXCEPTION( point.extent(0) != cellTopo.getDimension(), std::invalid_argument,
44 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): Point and cell dimensions do not match. ");
45#endif
46 bool testResult = true;
47
48 // A cell with extended topology has the same reference cell as a cell with base topology.
49 // => testing for inclusion in a reference Triangle<> and a reference Triangle<6> relies on
50 // on the same set of inequalities. To eliminate unnecessary cases we switch on the base topology
51 const auto key = cellTopo.getBaseKey();
52 switch (key) {
53
54 case shards::Line<>::key :
55 testResult = PointInclusion<shards::Line<>::key>::check(point, threshold);
56 break;
57 case shards::Triangle<>::key :
58 testResult = PointInclusion<shards::Triangle<>::key>::check(point, threshold);
59 break;
60 case shards::Quadrilateral<>::key :
61 testResult = PointInclusion<shards::Quadrilateral<>::key>::check(point, threshold);
62 break;
63 case shards::Tetrahedron<>::key :
64 testResult = PointInclusion<shards::Tetrahedron<>::key>::check(point, threshold);
65 break;
66 case shards::Hexahedron<>::key :
67 testResult = PointInclusion<shards::Hexahedron<>::key>::check(point, threshold);
68 break;
69 case shards::Wedge<>::key :
70 testResult = PointInclusion<shards::Wedge<>::key>::check(point, threshold);
71 break;
72 case shards::Pyramid<>::key :
73 testResult = PointInclusion<shards::Pyramid<>::key>::check(point, threshold);
74 break;
75 default:
76 INTREPID2_TEST_FOR_EXCEPTION( !( (key == shards::Line<>::key ) ||
77 (key == shards::Triangle<>::key) ||
78 (key == shards::Quadrilateral<>::key) ||
79 (key == shards::Tetrahedron<>::key) ||
80 (key == shards::Hexahedron<>::key) ||
81 (key == shards::Wedge<>::key) ||
82 (key == shards::Pyramid<>::key) ),
83 std::invalid_argument,
84 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): Invalid cell topology. ");
85 }
86 return testResult;
87 }
88
89
90
91 template<unsigned cellTopologyKey,
92 typename OutputViewType,
93 typename InputViewType>
94 struct checkPointInclusionFunctor {
95 OutputViewType output_;
96 InputViewType input_;
97 using ScalarType = typename ScalarTraits<typename InputViewType::value_type>::scalar_type;
98 ScalarType threshold_;
99
100 KOKKOS_INLINE_FUNCTION
101 checkPointInclusionFunctor( OutputViewType output,
102 const InputViewType input,
103 const ScalarType threshold)
104 : output_(output),
105 input_(input),
106 threshold_(threshold) {}
107
108 KOKKOS_INLINE_FUNCTION
109 void
110 operator()(const ordinal_type i) const {
111 const auto in = Kokkos::subview(input_,i,Kokkos::ALL());
112 const auto check = PointInclusion<cellTopologyKey>::check(in, threshold_);
113 output_(i) = check;
114 }
115
116 KOKKOS_INLINE_FUNCTION
117 void
118 operator()(const ordinal_type i, const ordinal_type j) const {
119 const auto in = Kokkos::subview(input_,i,j,Kokkos::ALL());
120 const auto check = PointInclusion<cellTopologyKey>::check(in, threshold_);
121 output_(i,j) = check;
122 }
123 };
124
125
126 template<typename DeviceType>
127 template<unsigned cellTopologyKey,
128 typename OutputViewType,
129 typename InputViewType>
130 void CellTools<DeviceType>::
131 checkPointwiseInclusion( OutputViewType inCell,
132 const InputViewType points,
133 const typename ScalarTraits<typename InputViewType::value_type>::scalar_type threshold) {
134
135 using FunctorType = checkPointInclusionFunctor<cellTopologyKey,decltype(inCell),decltype(points)>;
136 if (points.rank() == 2) { // inCell.rank() == 1
137 Kokkos::RangePolicy<typename DeviceType::execution_space> policy(0, points.extent(0));
138 Kokkos::parallel_for(policy, FunctorType(inCell, points, threshold));
139 } else { //points.rank() == 3, inCell.rank() == 2
140 Kokkos::MDRangePolicy<typename DeviceType::execution_space,Kokkos::Rank<2>> policy({0,0},{points.extent(0),points.extent(1)});
141 Kokkos::parallel_for(policy, FunctorType(inCell, points, threshold));
142 }
143 }
144
145
146 template<typename DeviceType>
147 template<typename InCellViewType,
148 typename InputViewType>
149 void
150 CellTools<DeviceType>::
151 checkPointwiseInclusion( InCellViewType inCell,
152 const InputViewType points,
153 const shards::CellTopology cellTopo,
154 const typename ScalarTraits<typename InputViewType::value_type>::scalar_type threshold ) {
155#ifdef HAVE_INTREPID2_DEBUG
156 {
157 INTREPID2_TEST_FOR_EXCEPTION( (inCell.rank() != 1) && (inCell.rank() != 2), std::invalid_argument,
158 ">>> ERROR (Intrepid2::CellTools::checkPointwiseInclusion): InCell must have rank 1 or 2. ");
159 INTREPID2_TEST_FOR_EXCEPTION( points.extent(points.rank()-1) != cellTopo.getDimension(), std::invalid_argument,
160 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): Points and cell dimensions do not match. ");
161 INTREPID2_TEST_FOR_EXCEPTION( inCell.rank() != (points.rank()-1), std::invalid_argument,
162 ">>> ERROR (Intrepid2::CellTools::checkPointwiseInclusion): rank difference between inCell and points is 1.");
163 const ordinal_type iend = inCell.rank();
164
165 for (ordinal_type i=0;i<iend;++i) {
166 INTREPID2_TEST_FOR_EXCEPTION( inCell.extent(i) != points.extent(i), std::invalid_argument,
167 ">>> ERROR (Intrepid2::CellTools::checkPointwiseInclusion): dimension mismatch between inCell and points. " );
168 }
169 }
170#endif
171
172 const auto key = cellTopo.getBaseKey();
173 switch (key) {
174
175 case shards::Line<>::key :
176 checkPointwiseInclusion<shards::Line<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
177 break;
178
179 case shards::Triangle<>::key :
180 checkPointwiseInclusion<shards::Triangle<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
181 break;
182
183 case shards::Quadrilateral<>::key :
184 checkPointwiseInclusion<shards::Quadrilateral<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
185 break;
186
187 case shards::Tetrahedron<>::key :
188 checkPointwiseInclusion<shards::Tetrahedron<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
189 break;
190
191 case shards::Hexahedron<>::key :
192 checkPointwiseInclusion<shards::Hexahedron<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
193 break;
194
195 case shards::Wedge<>::key :
196 checkPointwiseInclusion<shards::Wedge<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
197 break;
198
199 case shards::Pyramid<>::key :
200 checkPointwiseInclusion<shards::Pyramid<>::key,decltype(inCell),decltype(points)>(inCell, points, threshold);
201 break;
202
203 default:
204 INTREPID2_TEST_FOR_EXCEPTION( false,
205 std::invalid_argument,
206 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): Invalid cell topology. ");
207 }
208 }
209
210 template<typename DeviceType>
211 template<typename inCellValueType, class ...inCellProperties,
212 typename pointValueType, class ...pointProperties,
213 typename cellWorksetValueType, class ...cellWorksetProperties>
214 void
215 CellTools<DeviceType>::
216 checkPointwiseInclusion( Kokkos::DynRankView<inCellValueType,inCellProperties...> inCell,
217 const Kokkos::DynRankView<pointValueType,pointProperties...> points,
218 const Kokkos::DynRankView<cellWorksetValueType,cellWorksetProperties...> cellWorkset,
219 const shards::CellTopology cellTopo,
220 const typename ScalarTraits<pointValueType>::scalar_type threshold ) {
221#ifdef HAVE_INTREPID2_DEBUG
222 {
223 const auto key = cellTopo.getBaseKey();
224 INTREPID2_TEST_FOR_EXCEPTION( key != shards::Line<>::key &&
225 key != shards::Triangle<>::key &&
226 key != shards::Quadrilateral<>::key &&
227 key != shards::Tetrahedron<>::key &&
228 key != shards::Hexahedron<>::key &&
229 key != shards::Wedge<>::key &&
230 key != shards::Pyramid<>::key,
231 std::invalid_argument,
232 ">>> ERROR (Intrepid2::CellTools::checkPointwiseInclusion): cell topology not supported");
233 INTREPID2_TEST_FOR_EXCEPTION( inCell.rank() != 2, std::invalid_argument,
234 ">>> ERROR (Intrepid2::CellTools::checkPointwiseInclusion): InCell must have rank 2. ");
235 INTREPID2_TEST_FOR_EXCEPTION( cellWorkset.rank() != 3, std::invalid_argument,
236 ">>> ERROR (Intrepid2::CellTools::checkPointwiseInclusion): cellWorkset must have rank 3. ");
237 INTREPID2_TEST_FOR_EXCEPTION( cellWorkset.extent(2) != cellTopo.getDimension(), std::invalid_argument,
238 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): cellWorkset and points have incompatible dimensions. ");
239 INTREPID2_TEST_FOR_EXCEPTION( (points.rank() == 3) && (points.extent(0) != cellWorkset.extent(0)) , std::invalid_argument,
240 ">>> ERROR (Intrepid2::CellTools::checkPointInclusion): cellWorkset and points have incompatible dimensions. ");
241 }
242#endif
243 const ordinal_type
244 numCells = cellWorkset.extent(0),
245 numPoints = points.extent(points.rank()-2),
246 spaceDim = cellTopo.getDimension();
247
248 using result_layout = typename DeduceLayout< decltype(points) >::result_layout;
249 auto vcprop = Kokkos::common_view_alloc_prop(points);
250 using common_value_type = typename decltype(vcprop)::value_type;
251 Kokkos::DynRankView< common_value_type, result_layout, DeviceType > refPoints ( Kokkos::view_alloc("CellTools::checkPointwiseInclusion::refPoints", vcprop), numCells, numPoints, spaceDim);
252
253 // expect refPoints(CPD), points (CPD or PD), cellWorkset(CND)
254 if(points.rank() == 3)
255 mapToReferenceFrame(refPoints, points, cellWorkset, cellTopo);
256 else { //points.rank() == 2
257 Kokkos::DynRankView< common_value_type, result_layout, DeviceType > cellPoints ( Kokkos::view_alloc("CellTools::checkPointwiseInclusion::physCellPoints", vcprop), numCells, numPoints, spaceDim);
258 RealSpaceTools<DeviceType>::clone(cellPoints,points);
259 mapToReferenceFrame(refPoints, cellPoints, cellWorkset, cellTopo);
260 }
261 checkPointwiseInclusion(inCell, refPoints, cellTopo, threshold);
262 }
263
264}
265
266#endif
Intrepid2::CellTools::checkPointInclusion
static bool checkPointInclusion(const PointViewType point, const shards::CellTopology cellTopo, const typename ScalarTraits< typename PointViewType::value_type >::scalar_type thres=threshold< typename ScalarTraits< typename PointViewType::value_type >::scalar_type >())
Checks if a point belongs to a reference cell.
Definition Intrepid2_CellToolsDefInclusion.hpp:37
Intrepid2::CellTools::checkPointwiseInclusion
static void checkPointwiseInclusion(OutputViewType inCell, const InputViewType points, const typename ScalarTraits< typename InputViewType::value_type >::scalar_type thresh=threshold< typename ScalarTraits< typename InputViewType::value_type >::scalar_type >())
Checks every point for inclusion in the reference cell of a given topology. The points can belong to ...
Definition Intrepid2_CellToolsDefInclusion.hpp:131
Intrepid2::RealSpaceTools::clone
static void clone(Kokkos::DynRankView< outputValueType, outputProperties... > output, const Kokkos::DynRankView< inputValueType, inputProperties... > input)
Clone input array.
Definition Intrepid2_RealSpaceToolsDef.hpp:301
Intrepid2::DeduceLayout
layout deduction (temporary meta-function)
Definition Intrepid2_Utils.hpp:185
Intrepid2::PointInclusion
This class implements a check function that determines whether a given point is inside or outside the...
Definition Intrepid2_CellData.hpp:333
Generated by doxygen 1.9.8