Intrepid2
Intrepid2_PyramidCoords.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 Mauro Perego (mperego@sandia.gov) or
38 // Nate Roberts (nvrober@sandia.gov)
39 //
40 // ************************************************************************
41 // @HEADER
42 
48 #ifndef Intrepid2_Intrepid2_PyramidCoords_h
49 #define Intrepid2_Intrepid2_PyramidCoords_h
50 
51 #include <Kokkos_DynRankView.hpp>
52 
53 #include <Intrepid2_config.h>
54 
55 namespace Intrepid2
56 {
58  template<class PointScalar>
59  KOKKOS_INLINE_FUNCTION
60  void affinePyramid(Kokkos::Array<PointScalar,5> &lambda,
61  Kokkos::Array<Kokkos::Array<PointScalar,3>,5> &lambdaGrad,
62  Kokkos::Array<Kokkos::Array<PointScalar,3>,2> &mu,
63  Kokkos::Array<Kokkos::Array<Kokkos::Array<PointScalar,3>,3>,2> &muGrad,
64  Kokkos::Array<Kokkos::Array<PointScalar,2>,3> &nu,
65  Kokkos::Array<Kokkos::Array<Kokkos::Array<PointScalar,3>,2>,3> &nuGrad,
66  Kokkos::Array<PointScalar,3> &coords)
67  {
68  const auto & x = coords[0];
69  const auto & y = coords[1];
70  const auto & z = coords[2];
71  nu[0][0] = 1. - x - z; // nu_0^{\zeta,\xi_1}
72  nu[0][1] = 1. - y - z; // nu_0^{\zeta,\xi_2}
73  nu[1][0] = x ; // nu_1^{\zeta,\xi_1}
74  nu[1][1] = y ; // nu_1^{\zeta,\xi_2}
75  nu[2][0] = z ; // nu_2^{\zeta,\xi_1}
76  nu[2][1] = z ; // nu_2^{\zeta,\xi_2}
77 
78  nuGrad[0][0][0] = -1. ; // nu_0^{\zeta,\xi_1}_dxi_1
79  nuGrad[0][0][1] = 0. ; // nu_0^{\zeta,\xi_1}_dxi_2
80  nuGrad[0][0][2] = -1. ; // nu_0^{\zeta,\xi_1}_dzeta
81 
82  nuGrad[0][1][0] = 0. ; // nu_0^{\zeta,\xi_2}_dxi_1
83  nuGrad[0][1][1] = -1. ; // nu_0^{\zeta,\xi_2}_dxi_2
84  nuGrad[0][1][2] = -1. ; // nu_0^{\zeta,\xi_2}_dzeta
85 
86  nuGrad[1][0][0] = 1. ; // nu_1^{\zeta,\xi_1}_dxi_1
87  nuGrad[1][0][1] = 0. ; // nu_1^{\zeta,\xi_1}_dxi_2
88  nuGrad[1][0][2] = 0. ; // nu_1^{\zeta,\xi_1}_dzeta
89 
90  nuGrad[1][1][0] = 0. ; // nu_1^{\zeta,\xi_2}_dxi_1
91  nuGrad[1][1][1] = 1. ; // nu_1^{\zeta,\xi_2}_dxi_2
92  nuGrad[1][1][2] = 0. ; // nu_1^{\zeta,\xi_2}_dzeta
93 
94  nuGrad[2][0][0] = 0. ; // nu_2^{\zeta,\xi_1}_dxi_1
95  nuGrad[2][0][1] = 0. ; // nu_2^{\zeta,\xi_1}_dxi_2
96  nuGrad[2][0][2] = 1. ; // nu_2^{\zeta,\xi_1}_dzeta
97 
98  nuGrad[2][1][0] = 0. ; // nu_2^{\zeta,\xi_2}_dxi_1
99  nuGrad[2][1][1] = 0. ; // nu_2^{\zeta,\xi_2}_dxi_2
100  nuGrad[2][1][2] = 1. ; // nu_2^{\zeta,\xi_2}_dzeta
101 
102  // (1 - z) goes in denominator -- so we check for 1-z=0
103  auto & muZ_0 = mu[0][2];
104  auto & muZ_1 = mu[1][2];
105  const double epsilon = 1e-12;
106  muZ_0 = (fabs(1.-z) > epsilon) ? 1. - z : epsilon;
107  muZ_1 = (fabs(1.-z) > epsilon) ? z : 1. - epsilon;
108  PointScalar scaling = 1. / muZ_0;
109  mu[0][0] = 1. - x * scaling;
110  mu[0][1] = 1. - y * scaling;
111  mu[1][0] = x * scaling;
112  mu[1][1] = y * scaling;
113 
114  PointScalar scaling2 = scaling * scaling;
115  muGrad[0][0][0] = -scaling ;
116  muGrad[0][0][1] = 0. ;
117  muGrad[0][0][2] = - x * scaling2;
118 
119  muGrad[0][1][0] = 0. ;
120  muGrad[0][1][1] = -scaling ;
121  muGrad[0][1][2] = -y * scaling2;
122 
123  muGrad[0][2][0] = 0. ;
124  muGrad[0][2][1] = 0. ;
125  muGrad[0][2][2] = -1. ;
126 
127  muGrad[1][0][0] = scaling ;
128  muGrad[1][0][1] = 0. ;
129  muGrad[1][0][2] = x * scaling2;
130 
131  muGrad[1][1][0] = 0. ;
132  muGrad[1][1][1] = scaling ;
133  muGrad[1][1][2] = y * scaling2;
134 
135  muGrad[1][2][0] = 0. ;
136  muGrad[1][2][1] = 0. ;
137  muGrad[1][2][2] = 1. ;
138 
139  lambda[0] = nu[0][0] * mu[0][1];
140  lambda[1] = nu[0][1] * mu[1][0];
141  lambda[2] = nu[1][0] * mu[1][1];
142  lambda[3] = nu[1][1] * mu[0][0];
143  lambda[4] = z;
144 
145  for (int d=0; d<3; d++)
146  {
147  lambdaGrad[0][d] = nu[0][0] * muGrad[0][1][d] + nuGrad[0][0][d] * mu[0][1];
148  lambdaGrad[1][d] = nu[0][1] * muGrad[1][0][d] + nuGrad[0][1][d] * mu[1][0];
149  lambdaGrad[2][d] = nu[1][0] * muGrad[1][1][d] + nuGrad[1][0][d] * mu[1][1];
150  lambdaGrad[3][d] = nu[1][1] * muGrad[0][0][d] + nuGrad[1][1][d] * mu[0][0];
151  }
152  lambdaGrad[4][0] = 0;
153  lambdaGrad[4][1] = 0;
154  lambdaGrad[4][2] = 1;
155  }
156 
158  template<class PointScalar>
159  KOKKOS_INLINE_FUNCTION
160  void transformToESEASPyramid( PointScalar &x_eseas, PointScalar &y_eseas, PointScalar &z_eseas,
161  const PointScalar &x_int2, const PointScalar &y_int2, const PointScalar &z_int2)
162  {
163  x_eseas = (x_int2 + 1. - z_int2) / 2.;
164  y_eseas = (y_int2 + 1. - z_int2) / 2.;
165  z_eseas = z_int2;
166  }
167 
169  template<class OutputScalar>
170  KOKKOS_INLINE_FUNCTION
171  void transformFromESEASPyramidGradient( OutputScalar &dx_int2, OutputScalar &dy_int2, OutputScalar &dz_int2,
172  const OutputScalar &dx_eseas, const OutputScalar &dy_eseas, const OutputScalar &dz_eseas)
173  {
174  dx_int2 = dx_eseas / 2.;
175  dy_int2 = dy_eseas / 2.;
176  dz_int2 = dz_eseas - dx_int2 - dy_int2;
177  }
178 } // end namespace Intrepid2
179 
180 #endif /* Intrepid2_Intrepid2_PyramidCoords_h */
Generated by   doxygen 1.8.5