Intrepid2
Intrepid2_OrientationToolsDefModifyPoints.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid2 Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Kyungjoo Kim (kyukim@sandia.gov), or
38 // Mauro Perego (mperego@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 
43 
48 #ifndef __INTREPID2_ORIENTATIONTOOLS_DEF_MODIFY_POINTS_HPP__
49 #define __INTREPID2_ORIENTATIONTOOLS_DEF_MODIFY_POINTS_HPP__
50 
51 // disable clang warnings
52 #if defined (__clang__) && !defined (__INTEL_COMPILER)
53 #pragma clang system_header
54 #endif
55 
56 namespace Intrepid2 {
57 
58  namespace Impl {
59 
60  // ------------------------------------------------------------------------------------
61  // Modified points according to orientations
62  //
63  //
64  template<typename VT>
65  KOKKOS_INLINE_FUNCTION
66  void
69  const VT pt,
70  const ordinal_type ort) {
71 #ifdef HAVE_INTREPID2_DEBUG
72  INTREPID2_TEST_FOR_ABORT( !( -1.0 <= pt && pt <= 1.0 ),
73  ">>> ERROR (Intrepid::OrientationTools::getModifiedLinePoint): "
74  "Input point is out of range [-1, 1].");
75 #endif
76 
77  switch (ort) {
78  case 0: ot = pt; break;
79  case 1: ot = -pt; break;
80  default:
81  INTREPID2_TEST_FOR_ABORT( true,
82  ">>> ERROR (Intrepid2::OrientationTools::getModifiedLinePoint): "
83  "Orientation is invalid (0--1)." );
84  }
85  }
86 
87  template<typename JacobianViewType>
88  KOKKOS_INLINE_FUNCTION
89  void
91  getLineJacobian(JacobianViewType jacobian, const ordinal_type ort) {
92 #ifdef HAVE_INTREPID2_DEBUG
93  INTREPID2_TEST_FOR_ABORT( ort <0 || ort >1,
94  ">>> ERROR (Intrepid2::OrientationTools::getLineJacobian): " \
95  "Orientation is invalid (0--1)." );
96 #endif
97 
98  ordinal_type jac[2] = { 1, -1 };
99 
100  jacobian(0,0) = jac[ort];
101  }
102 
103  template<typename VT>
104  KOKKOS_INLINE_FUNCTION
105  void
108  VT &ot1,
109  const VT pt0,
110  const VT pt1,
111  const ordinal_type ort) {
112  const VT lambda[3] = { 1.0 - pt0 - pt1,
113  pt0,
114  pt1 };
115 
116 #ifdef HAVE_INTREPID2_DEBUG
117  INTREPID2_TEST_FOR_ABORT( !( 0.0 <= lambda[0] && lambda[0] <= 1.0 ),
118  ">>> ERROR (Intrepid::OrientationTools::getModifiedTrianglePoint): " \
119  "Computed bicentric coordinate (lamba[0]) is out of range [0, 1].");
120 
121  INTREPID2_TEST_FOR_ABORT( !( 0.0 <= lambda[1] && lambda[1] <= 1.0 ),
122  ">>> ERROR (Intrepid::OrientationTools::getModifiedTrianglePoint): " \
123  "Computed bicentric coordinate (lamba[1]) is out of range [0, 1].");
124 
125  INTREPID2_TEST_FOR_ABORT( !( 0.0 <= lambda[2] && lambda[2] <= 1.0 ),
126  ">>> ERROR (Intrepid::OrientationTools::getModifiedTrianglePoint): "
127  "Computed bicentric coordinate (lamba[2]) is out of range [0, 1].");
128 #endif
129 
130  switch (ort) {
131  case 0: ot0 = lambda[1]; ot1 = lambda[2]; break;
132  case 1: ot0 = lambda[0]; ot1 = lambda[1]; break;
133  case 2: ot0 = lambda[2]; ot1 = lambda[0]; break;
134 
135  case 3: ot0 = lambda[2]; ot1 = lambda[1]; break;
136  case 4: ot0 = lambda[0]; ot1 = lambda[2]; break;
137  case 5: ot0 = lambda[1]; ot1 = lambda[0]; break;
138  default:
139  INTREPID2_TEST_FOR_ABORT( true,
140  ">>> ERROR (Intrepid2::OrientationTools::getModifiedTrianglePoint): " \
141  "Orientation is invalid (0--5)." );
142  }
143  }
144 
145  template<typename JacobianViewType>
146  KOKKOS_INLINE_FUNCTION
147  void
149  getTriangleJacobian(JacobianViewType jacobian, const ordinal_type ort) {
150 #ifdef HAVE_INTREPID2_DEBUG
151  INTREPID2_TEST_FOR_ABORT( ort <0 || ort >5,
152  ">>> ERROR (Intrepid2::OrientationTools::getTriangleJacobian): " \
153  "Orientation is invalid (0--5)." );
154 #endif
155 
156  ordinal_type jac[6][2][2] = { { { 1, 0 },
157  { 0, 1 } }, // 0
158  { { -1, -1 },
159  { 1, 0 } }, // 1
160  { { 0, 1 },
161  { -1, -1 } }, // 2
162  { { 0, 1 },
163  { 1, 0 } }, // 3
164  { { -1, -1 },
165  { 0, 1 } }, // 4
166  { { 1, 0 },
167  { -1, -1 } } }; // 5
168 
169  jacobian(0,0) = jac[ort][0][0];
170  jacobian(0,1) = jac[ort][0][1];
171  jacobian(1,0) = jac[ort][1][0];
172  jacobian(1,1) = jac[ort][1][1];
173  }
174 
175  template<typename VT>
176  KOKKOS_INLINE_FUNCTION
177  void
180  VT &ot1,
181  const VT pt0,
182  const VT pt1,
183  const ordinal_type ort) {
184 #ifdef HAVE_INTREPID2_DEBUG
185  INTREPID2_TEST_FOR_ABORT( !( -1.0 <= pt0 && pt0 <= 1.0 ),
186  ">>> ERROR (Intrepid::OrientationTools::getModifiedQuadrilateralPoint): " \
187  "Input point(0) is out of range [-1, 1].");
188 
189  INTREPID2_TEST_FOR_ABORT( !( -1.0 <= pt1 && pt1 <= 1.0 ),
190  ">>> ERROR (Intrepid::OrientationTools::getModifiedQuadrilateralPoint): " \
191  "Input point(1) is out of range [-1, 1].");
192 #endif
193 
194  const VT lambda[2][2] = { { pt0, -pt0 },
195  { pt1, -pt1 } };
196 
197  switch (ort) {
198  case 0: ot0 = lambda[0][0]; ot1 = lambda[1][0]; break;
199  case 1: ot0 = lambda[1][1]; ot1 = lambda[0][0]; break;
200  case 2: ot0 = lambda[0][1]; ot1 = lambda[1][1]; break;
201  case 3: ot0 = lambda[1][0]; ot1 = lambda[0][1]; break;
202  // flip
203  case 4: ot0 = lambda[1][0]; ot1 = lambda[0][0]; break;
204  case 5: ot0 = lambda[0][1]; ot1 = lambda[1][0]; break;
205  case 6: ot0 = lambda[1][1]; ot1 = lambda[0][1]; break;
206  case 7: ot0 = lambda[0][0]; ot1 = lambda[1][1]; break;
207  default:
208  INTREPID2_TEST_FOR_ABORT( true,
209  ">>> ERROR (Intrepid2::OrientationTools::getModifiedQuadrilateralPoint): " \
210  "Orientation is invalid (0--7)." );
211  }
212  }
213 
214  template<typename JacobianViewType>
215  KOKKOS_INLINE_FUNCTION
216  void
218  getQuadrilateralJacobian(JacobianViewType jacobian, const ordinal_type ort) {
219 #ifdef HAVE_INTREPID2_DEBUG
220  INTREPID2_TEST_FOR_ABORT( ort <0 || ort >7,
221  ">>> ERROR (Intrepid2::OrientationTools::getQuadrilateralJacobian): " \
222  "Orientation is invalid (0--7)." );
223 #endif
224 
225  ordinal_type jac[8][2][2] = { { { 1, 0 },
226  { 0, 1 } }, // 0
227  { { 0, -1 },
228  { 1, 0 } }, // 1
229  { { -1, 0 },
230  { 0, -1 } }, // 2
231  { { 0, 1 },
232  { -1, 0 } }, // 3
233  { { 0, 1 },
234  { 1, 0 } }, // 4
235  { { -1, 0 },
236  { 0, 1 } }, // 5
237  { { 0, -1 },
238  { -1, 0 } }, // 6
239  { { 1, 0 },
240  { 0, -1 } } }; // 7
241 
242  jacobian(0,0) = jac[ort][0][0];
243  jacobian(0,1) = jac[ort][0][1];
244  jacobian(1,0) = jac[ort][1][0];
245  jacobian(1,1) = jac[ort][1][1];
246  }
247 
248  template<typename outPointViewType,
249  typename refPointViewType>
250  inline
251  void
253  mapToModifiedReference(outPointViewType outPoints,
254  const refPointViewType refPoints,
255  const shards::CellTopology cellTopo,
256  const ordinal_type cellOrt) {
257 #ifdef HAVE_INTREPID2_DEBUG
258  {
259  const auto cellDim = cellTopo.getDimension();
260  INTREPID2_TEST_FOR_EXCEPTION( !( (1 <= cellDim) && (cellDim <= 2 ) ), std::invalid_argument,
261  ">>> ERROR (Intrepid::OrientationTools::mapToModifiedReference): " \
262  "Method defined only for 1 and 2-dimensional subcells.");
263 
264  INTREPID2_TEST_FOR_EXCEPTION( !( outPoints.extent(0) == refPoints.extent(0) ), std::invalid_argument,
265  ">>> ERROR (Intrepid::OrientationTools::mapToModifiedReference): " \
266  "Size of input and output point arrays does not match each other.");
267  }
268 #endif
269 
270  // Apply the parametrization map to every point in parameter domain
271  const ordinal_type numPts = outPoints.extent(0);
272  const auto key = cellTopo.getBaseCellTopologyData()->key;
273  switch (key) {
274  case shards::Line<>::key : {
275  for (ordinal_type pt=0;pt<numPts;++pt)
276  getModifiedLinePoint(outPoints(pt, 0),
277  refPoints(pt, 0),
278  cellOrt);
279  break;
280  }
281  case shards::Triangle<>::key : {
282  for (ordinal_type pt=0;pt<numPts;++pt)
283  getModifiedTrianglePoint(outPoints(pt, 0), outPoints(pt, 1),
284  refPoints(pt, 0), refPoints(pt, 1),
285  cellOrt);
286  break;
287  }
288  case shards::Quadrilateral<>::key : {
289  for (ordinal_type pt=0;pt<numPts;++pt)
290  getModifiedQuadrilateralPoint(outPoints(pt, 0), outPoints(pt, 1),
291  refPoints(pt, 0), refPoints(pt, 1),
292  cellOrt);
293  break;
294  }
295  default: {
296  INTREPID2_TEST_FOR_WARNING( true,
297  ">>> ERROR (Intrepid2::OrientationTools::mapToModifiedReference): " \
298  "Invalid cell topology." );
299  break;
300  }
301  }
302  }
303 
304 
305  template<typename outPointViewType>
306  inline
307  void
309  getJacobianOfOrientationMap(outPointViewType jacobian,
310  const shards::CellTopology cellTopo,
311  const ordinal_type cellOrt) {
312 #ifdef HAVE_INTREPID2_DEBUG
313  {
314  const auto cellDim = cellTopo.getDimension();
315  INTREPID2_TEST_FOR_EXCEPTION( !( (1 <= cellDim) && (cellDim <= 2 ) ), std::invalid_argument,
316  ">>> ERROR (Intrepid::OrientationTools::getJacobianOfOrientationMap): " \
317  "Method defined only for 1 and 2-dimensional subcells.");
318 
319  INTREPID2_TEST_FOR_ABORT( jacobian.rank() != 2,
320  ">>> ERROR (Intrepid2::OrientationTools::getJacobianOfOrientationMap): " \
321  "Jacobian should have rank 2" );
322 
323  INTREPID2_TEST_FOR_EXCEPTION( ((jacobian.extent(0) != cellDim) || (jacobian.extent(1) != cellDim)), std::invalid_argument,
324  ">>> ERROR (Intrepid::OrientationTools::getJacobianOfOrientationMap): " \
325  "Size of jacobian is not compatible with cell topology.");
326  }
327 #endif
328 
329  const auto key = cellTopo.getBaseCellTopologyData()->key;
330  switch (key) {
331  case shards::Line<>::key :
332  getLineJacobian(jacobian, cellOrt);
333  break;
334  case shards::Triangle<>::key :
335  getTriangleJacobian(jacobian, cellOrt);
336  break;
337  case shards::Quadrilateral<>::key :
338  getQuadrilateralJacobian(jacobian, cellOrt);
339  break;
340  default: {
341  INTREPID2_TEST_FOR_WARNING( true,
342  ">>> ERROR (Intrepid2::OrientationTools::mapToModifiedReference): " \
343  "Invalid cell topology." );
344  break;
345  }
346  }
347  }
348 
349 
350 
351  }
352 }
353 
354 #endif
static KOKKOS_INLINE_FUNCTION void getModifiedLinePoint(ValueType &ot, const ValueType pt, const ordinal_type ort)
Computes modified point for line segment.
static KOKKOS_INLINE_FUNCTION void getModifiedQuadrilateralPoint(ValueType &ot0, ValueType &ot1, const ValueType pt0, const ValueType pt1, const ordinal_type ort)
Computes modified point for quadrilateral.
static void getJacobianOfOrientationMap(JacobianViewType jacobian, const shards::CellTopology cellTopo, const ordinal_type cellOrt)
Computes jacobian of the parameterization maps of 1- and 2-subcells with orientation.
static KOKKOS_INLINE_FUNCTION void getModifiedTrianglePoint(ValueType &ot0, ValueType &ot1, const ValueType pt0, const ValueType pt1, const ordinal_type ort)
Computes modified point for triangle.
static KOKKOS_INLINE_FUNCTION void getTriangleJacobian(JacobianViewType jacobian, const ordinal_type ort)
Computes Jacobian of orientation map for triangle.
static KOKKOS_INLINE_FUNCTION void getLineJacobian(JacobianViewType jacobian, const ordinal_type ort)
Computes Jacobian of orientation map for line segment.
static KOKKOS_INLINE_FUNCTION void getQuadrilateralJacobian(JacobianViewType jacobian, const ordinal_type ort)
Computes Jacobian of orientation map for quadrilateral.
static void mapToModifiedReference(outPointViewType outPoints, const refPointViewType refPoints, const shards::CellTopology cellTopo, const ordinal_type cellOrt=0)
Computes modified parameterization maps of 1- and 2-subcells with orientation.