Intrepid
Intrepid_HGRAD_PYR_I2_FEMDef.hpp
Go to the documentation of this file.
1 #ifndef INTREPID_HGRAD_PYR_I2_FEMDEF_HPP
2 #define INTREPID_HGRAD_PYR_I2_FEMDEF_HPP
3 // @HEADER
4 // ************************************************************************
5 //
6 // Intrepid Package
7 // Copyright (2007) Sandia Corporation
8 //
9 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
10 // license for use of this work by or on behalf of the U.S. Government.
11 //
12 // Redistribution and use in source and binary forms, with or without
13 // modification, are permitted provided that the following conditions are
14 // met:
15 //
16 // 1. Redistributions of source code must retain the above copyright
17 // notice, this list of conditions and the following disclaimer.
18 //
19 // 2. Redistributions in binary form must reproduce the above copyright
20 // notice, this list of conditions and the following disclaimer in the
21 // documentation and/or other materials provided with the distribution.
22 //
23 // 3. Neither the name of the Corporation nor the names of the
24 // contributors may be used to endorse or promote products derived from
25 // this software without specific prior written permission.
26 //
27 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38 //
39 // Questions? Contact Pavel Bochev (pbboche@sandia.gov)
40 // Denis Ridzal (dridzal@sandia.gov), or
41 // Kara Peterson (kjpeter@sandia.gov)
42 //
43 // ************************************************************************
44 // @HEADER
45 
51 namespace Intrepid {
52 
53  template<class Scalar, class ArrayScalar>
54  Basis_HGRAD_PYR_I2_FEM<Scalar, ArrayScalar>::Basis_HGRAD_PYR_I2_FEM()
55  {
56  this -> basisCardinality_ = 13;
57  this -> basisDegree_ = 2;
58  this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Pyramid<5> >() );
59  this -> basisType_ = BASIS_FEM_DEFAULT;
60  this -> basisCoordinates_ = COORDINATES_CARTESIAN;
61  this -> basisTagsAreSet_ = false;
62  }
63 
64 
65 template<class Scalar, class ArrayScalar>
66 void Basis_HGRAD_PYR_I2_FEM<Scalar, ArrayScalar>::initializeTags() {
67 
68  // Basis-dependent intializations
69  int tagSize = 4; // size of DoF tag
70  int posScDim = 0; // position in the tag, counting from 0, of the subcell dim
71  int posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
72  int posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
73 
74  // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
75  int tags[] = { 0, 0, 0, 1,
76  0, 1, 0, 1,
77  0, 2, 0, 1,
78  0, 3, 0, 1,
79  0, 4, 0, 1,
80  1, 0, 0, 1,
81  1, 1, 0, 1,
82  1, 2, 0, 1,
83  1, 3, 0, 1,
84  1, 4, 0, 1,
85  1, 5, 0, 1,
86  1, 6, 0, 1,
87  1, 7, 0, 1,
88  };
89 
90  // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
91  Intrepid::setOrdinalTagData(this -> tagToOrdinal_,
92  this -> ordinalToTag_,
93  tags,
94  this -> basisCardinality_,
95  tagSize,
96  posScDim,
97  posScOrd,
98  posDfOrd);
99 }
100 
101 
102 
103 template<class Scalar, class ArrayScalar>
104 void Basis_HGRAD_PYR_I2_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
105  const ArrayScalar & inputPoints,
106  const EOperator operatorType) const {
107 
108  // Verify arguments
109 #ifdef HAVE_INTREPID_DEBUG
110  Intrepid::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
111  inputPoints,
112  operatorType,
113  this -> getBaseCellTopology(),
114  this -> getCardinality() );
115 #endif
116 
117  // Number of evaluation points = dim 0 of inputPoints
118  int dim0 = inputPoints.dimension(0);
119 
120  // Temporaries: (x,y,z) coordinates of the evaluation point
121  Scalar x = 0.0;
122  Scalar y = 0.0;
123  Scalar z = 0.0;
124  const Scalar eps = std::numeric_limits<Scalar>::epsilon( );
125 
126  switch (operatorType) {
127 
128  case OPERATOR_VALUE:
129  for (int i0 = 0; i0 < dim0; i0++) {
130  x = inputPoints(i0, 0);
131  y = inputPoints(i0, 1);
132  z = inputPoints(i0, 2);
133 
134  //be sure that the basis functions are defined when z is very close to 1.
135 
136  if(fabs(z-1.0) < eps) {
137  if(z <= 1.0) z = 1.0-eps;
138  else z = 1.0+eps;
139  //z = 1.0-eps;
140  }
141 
142  Scalar w = 1.0/(1.0 - z);
143 
144  // outputValues is a rank-2 array with dimensions (basisCardinality_, dim0)
145  outputValues(0, i0) = 0.25 * (-x - y - 1.0)*((1.0-x)*(1.0-y) - z + x*y*z*w);
146  outputValues(1, i0) = 0.25 * ( x - y - 1.0)*((1.0+x)*(1.0-y) - z - x*y*z*w);
147  outputValues(2, i0) = 0.25 * ( x + y - 1.0)*((1.0+x)*(1.0+y) - z + x*y*z*w);
148  outputValues(3, i0) = 0.25 * (-x + y - 1.0)*((1.0-x)*(1.0+y) - z - x*y*z*w);
149 
150  outputValues(4, i0) = z * (2.0*z - 1.0);
151 
152  outputValues(5, i0) = 0.5 * (1.0 + x - z)*(1.0 - x - z)*(1.0 - y - z)*w;
153  outputValues(6, i0) = 0.5 * (1.0 + y - z)*(1.0 - y - z)*(1.0 + x - z)*w;
154  outputValues(7, i0) = 0.5 * (1.0 + x - z)*(1.0 - x - z)*(1.0 + y - z)*w;
155  outputValues(8, i0) = 0.5 * (1.0 + y - z)*(1.0 - y - z)*(1.0 - x - z)*w;
156 
157  outputValues(9, i0) = z*(1.0 - x - z)*(1.0 - y - z)*w;
158  outputValues(10,i0) = z*(1.0 + x - z)*(1.0 - y - z)*w;
159  outputValues(11,i0) = z*(1.0 + x - z)*(1.0 + y - z)*w;
160  outputValues(12,i0) = z*(1.0 - x - z)*(1.0 + y - z)*w;
161  }
162  break;
163 
164  case OPERATOR_GRAD:
165  case OPERATOR_D1:
166  for (int i0 = 0; i0 < dim0; i0++) {
167  x = inputPoints(i0,0);
168  y = inputPoints(i0,1);
169  z = inputPoints(i0,2);
170 
171  //be sure that the basis functions are defined when z is very close to 1.
172 
173  if(fabs(z-1.0) < eps) {
174  if(z <= 1.0) z = 1.0-eps;
175  else z = 1.0+eps;
176  }
177 
178  Scalar w = 1.0/(1.0 - z);
179 
180  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
181  outputValues(0, i0, 0) = 0.25*(-1.0-x-y)*(-1.0+y + y*z*w) - 0.25*((1.0-x)*(1.0-y)-z + x*y*z*w);
182  outputValues(0, i0, 1) = 0.25*(-1.0-x-y)*(-1.0+x + x*z*w) - 0.25*((1.0-x)*(1.0-y)-z + x*y*z*w);
183  outputValues(0, i0, 2) = 0.25*(-1.0-x-y)*(-1.0 + x*y*w + x*y*z*w*w);
184 
185  outputValues(1, i0, 0) = 0.25*(-1.0+x-y)*( 1.0-y - y*z*w) + 0.25*((1.0+x)*(1.0-y)-z - x*y*z*w);
186  outputValues(1, i0, 1) = 0.25*(-1.0+x-y)*(-1.0-x - x*z*w) - 0.25*((1.0+x)*(1.0-y)-z - x*y*z*w);
187  outputValues(1, i0, 2) = 0.25*(-1.0+x-y)*(-1.0 - x*y*w - x*y*z*w*w);
188 
189  outputValues(2, i0, 0) = 0.25*(-1.0+x+y)*(1.0+y + y*z*w) + 0.25*((1.0+x)*(1.0+y)-z + x*y*z*w);
190  outputValues(2, i0, 1) = 0.25*(-1.0+x+y)*(1.0+x + x*z*w) + 0.25*((1.0+x)*(1.0+y)-z + x*y*z*w);
191  outputValues(2, i0, 2) = 0.25*(-1.0+x+y)*(-1.0 + x*y*w + x*y*z*w*w);
192 
193  outputValues(3, i0, 0) = 0.25*(-1.0-x+y)*(-1.0-y - y*z*w) - 0.25*((1.0-x)*(1.0+y)-z - x*y*z*w);
194  outputValues(3, i0, 1) = 0.25*(-1.0-x+y)*( 1.0-x - x*z*w) + 0.25*((1.0-x)*(1.0+y)-z - x*y*z*w);
195  outputValues(3, i0, 2) = 0.25*(-1.0-x+y)*(-1.0 - x*y*w - x*y*z*w*w);
196 
197  outputValues(4, i0, 0) = 0.0;
198  outputValues(4, i0, 1) = 0.0;
199  outputValues(4, i0, 2) = -1.0 + 4.0*z;
200 
201  outputValues(5, i0, 0) = -x*w*(1.0-y-z);
202  outputValues(5, i0, 1) = -0.5*(1.0-x-z)*(1.0+x-z)*w;
203  outputValues(5, i0, 2) = 0.5*y*x*x*w*w + 0.5*y - 1.0+z;
204 
205  outputValues(6, i0, 0) = 0.5*(1.0-y-z)*(1.0+y-z)*w;
206  outputValues(6, i0, 1) = -y*w*(1.0+x-z);
207  outputValues(6, i0, 2) = -0.5*x*y*y*w*w - 0.5*x - 1.0+z;
208 
209  outputValues(7, i0, 0) = -x*w*(1.0+y-z);
210  outputValues(7, i0, 1) = 0.5*(1.0-x-z)*(1.0+x-z)*w;
211  outputValues(7, i0, 2) = -0.5*y*x*x*w*w - 0.5*y - 1.0+z;
212 
213  outputValues(8, i0, 0) = -0.5*(1.0-y-z)*(1.0+y-z)*w;
214  outputValues(8, i0, 1) = -y*w*(1.0-x-z);
215  outputValues(8, i0, 2) = 0.5*x*y*y*w*w + 0.5*x - 1.0+z;
216 
217  outputValues(9, i0, 0) = -(1.0-y-z)*z*w;
218  outputValues(9, i0, 1) = -(1.0-x-z)*z*w;
219  outputValues(9, i0, 2) = x*y*w*w + 1.0 - x - y - 2.0*z;
220 
221  outputValues(10,i0, 0) = (1.0-y-z)*z*w;
222  outputValues(10,i0, 1) = -(1.0+x-z)*z*w;
223  outputValues(10,i0, 2) = -x*y*w*w + 1.0 + x - y - 2.0*z;
224 
225  outputValues(11,i0, 0) = (1.0+y-z)*z*w;
226  outputValues(11,i0, 1) = (1.0+x-z)*z*w;
227  outputValues(11,i0, 2) = x*y*w*w + 1.0 + x + y - 2.0*z;
228 
229  outputValues(12,i0, 0) = -(1.0+y-z)*z*w;
230  outputValues(12,i0, 1) = (1.0-x-z)*z*w;
231  outputValues(12,i0, 2) = -x*y*w*w + 1.0 - x + y - 2.0*z;
232 
233  }
234  break;
235 
236  case OPERATOR_CURL:
237  TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
238  ">>> ERROR (Basis_HGRAD_PYR_I2_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
239  break;
240 
241  case OPERATOR_DIV:
242  TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
243  ">>> ERROR (Basis_HGRAD_PYR_I2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
244  break;
245 
246  case OPERATOR_D2:
247  for (int i0 = 0; i0 < dim0; i0++) {
248  x = inputPoints(i0,0);
249  y = inputPoints(i0,1);
250  z = inputPoints(i0,2);
251 
252  if(fabs(z-1.0) < eps) {
253  if(z <= 1.0) z = 1.0-eps;
254  else z = 1.0+eps;
255  }
256  Scalar w = 1.0/(1.0 - z);
257 
258  outputValues(0, i0, 0) = -0.5*(-1.0+y+y*z*w);
259  outputValues(0, i0, 1) = -(-0.25 + 0.5*x + 0.5*y + 0.5*z)*w;
260  outputValues(0, i0, 2) = 0.25 + (-0.25*y-0.5*x*y-0.25*y*y)*w*w;
261  outputValues(0, i0, 3) = -0.5*(-1.0+x+x*z*w);
262  outputValues(0, i0, 4) = 0.25 + (-0.25*x*x-0.25*x-0.5*x*y)*w*w;
263  outputValues(0, i0, 5) = 0.5*x*y*(-1.0-x-y)*w*w*w;
264 
265  outputValues(1, i0, 0) = 0.5*(1.0-y-y*z*w);
266  outputValues(1, i0, 1) =-(0.25 + 0.5*x - 0.5*y - 0.5*z)*w;
267  outputValues(1, i0, 2) = -0.25 + (0.25*y-0.5*x*y+0.25*y*y)*w*w;
268  outputValues(1, i0, 3) = -0.5*(-1.0-x-x*z*w);
269  outputValues(1, i0, 4) = 0.25 + (-0.25*x*x + 0.25*x + 0.5*x*y)*w*w;
270  outputValues(1, i0, 5) = -0.5*x*y*(-1.0+x-y)*w*w*w;
271 
272  outputValues(2, i0, 0) = 0.5*(1.0+y+y*z*w);
273  outputValues(2, i0, 1) =-(-0.25 - 0.5*x - 0.5*y + 0.5*z)*w;
274  outputValues(2, i0, 2) = -0.25 + (-0.25*y+0.5*x*y+0.25*y*y)*w*w;
275  outputValues(2, i0, 3) = 0.5*(1.0+x+x*z*w);
276  outputValues(2, i0, 4) = -0.25 + (0.25*x*x -0.25*x + 0.5*x*y)*w*w;
277  outputValues(2, i0, 5) = 0.5*x*y*(-1.0+x+y)*w*w*w;
278 
279  outputValues(3, i0, 0) = -0.5*(-1.0-y-y*z*w);
280  outputValues(3, i0, 1) =-(0.25 - 0.5*x + 0.5*y - 0.5*z)*w;
281  outputValues(3, i0, 2) = 0.25 + (0.25*y+0.5*x*y-0.25*y*y)*w*w;
282  outputValues(3, i0, 3) = 0.5*(1.0-x-x*z*w);
283  outputValues(3, i0, 4) = -0.25 + (0.25*x + 0.25*x*x - 0.5*x*y)*w*w;
284  outputValues(3, i0, 5) = -0.5*x*y*(-1.0-x+y)*w*w*w;
285 
286  outputValues(4, i0, 0) = 0.0;
287  outputValues(4, i0, 1) = 0.0;
288  outputValues(4, i0, 2) = 0.0;
289  outputValues(4, i0, 3) = 0.0;
290  outputValues(4, i0, 4) = 0.0;
291  outputValues(4, i0, 5) = 4.0;
292 
293  outputValues(5, i0, 0) = -(1.0-y-z)*w;
294  outputValues(5, i0, 1) = x*w;
295  outputValues(5, i0, 2) = x*y*w*w;
296  outputValues(5, i0, 3) = 0.0;
297  outputValues(5, i0, 4) = 0.5*x*x*w*w + 0.5;
298  outputValues(5, i0, 5) = x*x*y*w*w*w + 1.0;
299 
300  outputValues(6, i0, 0) = 0.0;
301  outputValues(6, i0, 1) = -y*w;
302  outputValues(6, i0, 2) = -0.5*y*y*w*w - 0.5;
303  outputValues(6, i0, 3) =-(1.0+x-z)*w;
304  outputValues(6, i0, 4) = -x*y*w*w;
305  outputValues(6, i0, 5) = -x*y*y*w*w*w + 1.0;
306 
307  outputValues(7, i0, 0) = -(1.0+y-z)*w;
308  outputValues(7, i0, 1) = -x*w;
309  outputValues(7, i0, 2) = -x*y*w*w;
310  outputValues(7, i0, 3) = 0.0;
311  outputValues(7, i0, 4) = -0.5*x*x*w*w - 0.5;
312  outputValues(7, i0, 5) = -x*x*y*w*w*w + 1.0;
313 
314  outputValues(8, i0, 0) = 0.0;
315  outputValues(8, i0, 1) = y*w;
316  outputValues(8, i0, 2) = 0.5*y*y*w*w + 0.5;
317  outputValues(8, i0, 3) = -(1.0-x-z)*w;
318  outputValues(8, i0, 4) = x*y*w*w;
319  outputValues(8, i0, 5) = x*y*y*w*w*w + 1.0;
320 
321  outputValues(9, i0, 0) = 0.0;
322  outputValues(9, i0, 1) = z*w;
323  outputValues(9, i0, 2) = y*w*w - 1.0;
324  outputValues(9, i0, 3) = 0.0;
325  outputValues(9, i0, 4) = x*w*w - 1.0;
326  outputValues(9, i0, 5) = 2.0*x*y*w*w*w - 2.0;
327 
328  outputValues(10,i0, 0) = 0.0;
329  outputValues(10,i0, 1) = -z*w;
330  outputValues(10,i0, 2) = -y*w*w + 1.0;
331  outputValues(10,i0, 3) = 0.0;
332  outputValues(10,i0, 4) = -x*w*w - 1.0;
333  outputValues(10,i0, 5) = -2.0*x*y*w*w*w - 2.0;
334 
335  outputValues(11,i0, 0) = 0.0;
336  outputValues(11,i0, 1) = z*w;
337  outputValues(11,i0, 2) = y*w*w + 1.0;
338  outputValues(11,i0, 3) = 0.0;
339  outputValues(11,i0, 4) = x*w*w + 1.0;
340  outputValues(11,i0, 5) = 2.0*x*y*w*w*w - 2.0;
341 
342  outputValues(12,i0, 0) = 0.0;
343  outputValues(12,i0, 1) = -z*w;
344  outputValues(12,i0, 2) = -y*w*w - 1.0;
345  outputValues(12,i0, 3) = 0.0;
346  outputValues(12,i0, 4) = -x*w*w + 1.0;
347  outputValues(12,i0, 5) = -2.0*x*y*w*w*w - 2.0;
348 
349  }
350  break;
351 
352  case OPERATOR_D3:
353  for (int i0 = 0; i0 < dim0; i0++) {
354  x = inputPoints(i0,0);
355  y = inputPoints(i0,1);
356  z = inputPoints(i0,2);
357 
358  if(fabs(z-1.0) < eps) {
359  if(z <= 1.0) z = 1.0-eps;
360  else z = 1.0+eps;
361  }
362 
363  Scalar w = 1.0/(1.0 - z);
364 
365  outputValues(0, i0, 0) = 0.0;
366  outputValues(0, i0, 1) = -0.5*w;
367  outputValues(0, i0, 2) = -0.5*y*w*w;
368  outputValues(0, i0, 3) = -0.5*w;
369  outputValues(0, i0, 4) =(-0.25 - 0.5*x - 0.5*y)*w*w;
370  outputValues(0, i0, 5) =-(0.5*y + x*y + 0.5*y*y)*w*w*w;
371  outputValues(0, i0, 6) = 0.0;
372  outputValues(0, i0, 7) = -0.5*x*w*w;
373  outputValues(0, i0, 8) =-(0.5*x + x*y + 0.5*x*x)*w*w*w;
374  outputValues(0, i0, 9) = 1.5*x*y*(-1.0 - x - y)*w*w*w*w;
375 
376  outputValues(1, i0, 0) = 0.0;
377  outputValues(1, i0, 1) = -0.5*w;
378  outputValues(1, i0, 2) = -0.5*y*w*w;
379  outputValues(1, i0, 3) = 0.5*w;
380  outputValues(1, i0, 4) = (0.25 - 0.5*x - 0.5*y)*w*w;
381  outputValues(1, i0, 5) =-(-0.5*y + x*y - 0.5*y*y)*w*w*w;
382  outputValues(1, i0, 6) = 0.0;
383  outputValues(1, i0, 7) = 0.5*x*w*w;
384  outputValues(1, i0, 8) =-(-0.5*x - x*y + 0.5*x*x)*w*w*w;
385  outputValues(1, i0, 9) = -1.5*x*y*(-1.0 + x - y)*w*w*w*w;
386 
387  outputValues(2, i0, 0) = 0.0;
388  outputValues(2, i0, 1) = 0.5*w;
389  outputValues(2, i0, 2) = 0.5*y*w*w;
390  outputValues(2, i0, 3) = 0.5*w;
391  outputValues(2, i0, 4) = (-0.25 + 0.5*x + 0.5*y)*w*w;
392  outputValues(2, i0, 5) =-(0.5*y - x*y - 0.5*y*y)*w*w*w;
393  outputValues(2, i0, 6) = 0.0;
394  outputValues(2, i0, 7) = 0.5*x*w*w;
395  outputValues(2, i0, 8) =-(0.5*x - x*y - 0.5*x*x)*w*w*w;
396  outputValues(2, i0, 9) = 1.5*x*y*(-1.0 + x + y)*w*w*w*w;
397 
398  outputValues(3, i0, 0) = 0.0;
399  outputValues(3, i0, 1) = 0.5*w;
400  outputValues(3, i0, 2) = 0.5*y*w*w;
401  outputValues(3, i0, 3) = -0.5*w;
402  outputValues(3, i0, 4) = (0.25 + 0.5*x - 0.5*y)*w*w;
403  outputValues(3, i0, 5) =-(-0.5*y - x*y + 0.5*y*y)*w*w*w;
404  outputValues(3, i0, 6) = 0.0;
405  outputValues(3, i0, 7) = -0.5*x*w*w;
406  outputValues(3, i0, 8) =-(-0.5*x + x*y - 0.5*x*x)*w*w*w;
407  outputValues(3, i0, 9) = -1.5*x*y*(-1.0 - x + y)*w*w*w*w;
408 
409  outputValues(4, i0, 0) = 0.0;
410  outputValues(4, i0, 1) = 0.0;
411  outputValues(4, i0, 2) = 0.0;
412  outputValues(4, i0, 3) = 0.0;
413  outputValues(4, i0, 4) = 0.0;
414  outputValues(4, i0, 5) = 0.0;
415  outputValues(4, i0, 6) = 0.0;
416  outputValues(4, i0, 7) = 0.0;
417  outputValues(4, i0, 8) = 0.0;
418  outputValues(4, i0, 9) = 0.0;
419 
420  outputValues(5, i0, 0) = 0.0;
421  outputValues(5, i0, 1) = w;
422  outputValues(5, i0, 2) = y*w*w;
423  outputValues(5, i0, 3) = 0.0;
424  outputValues(5, i0, 4) = x*w*w;
425  outputValues(5, i0, 5) = 2.0*y*x*w*w*w;
426  outputValues(5, i0, 6) = 0.0;
427  outputValues(5, i0, 7) = 0.0;
428  outputValues(5, i0, 8) = x*x*w*w*w;
429  outputValues(5, i0, 9) = 3.0*x*x*y*w*w*w*w;
430 
431  outputValues(6, i0, 0) = 0.0;
432  outputValues(6, i0, 1) = 0.0;
433  outputValues(6, i0, 2) = 0.0;
434  outputValues(6, i0, 3) = -w;
435  outputValues(6, i0, 4) = -y*w*w;
436  outputValues(6, i0, 5) = -y*y*w*w*w;
437  outputValues(6, i0, 6) = 0.0;
438  outputValues(6, i0, 7) = -x*w*w ;
439  outputValues(6, i0, 8) = -2.0*x*y*w*w*w;
440  outputValues(6, i0, 9) = -3.0*x*y*y*w*w*w*w;
441 
442  outputValues(7, i0, 0) = 0.0;
443  outputValues(7, i0, 1) = -w;
444  outputValues(7, i0, 2) = -y*w*w;
445  outputValues(7, i0, 3) = 0.0;
446  outputValues(7, i0, 4) = -x*w*w;
447  outputValues(7, i0, 5) = -2.0*x*y*w*w*w;
448  outputValues(7, i0, 6) = 0.0;
449  outputValues(7, i0, 7) = 0.0;
450  outputValues(7, i0, 8) = -x*x*w*w*w;
451  outputValues(7, i0, 9) = -3.0*x*x*y*w*w*w*w;
452 
453  outputValues(8, i0, 0) = 0.0;
454  outputValues(8, i0, 1) = 0.0;
455  outputValues(8, i0, 2) = 0.0;
456  outputValues(8, i0, 3) = w;
457  outputValues(8, i0, 4) = y*w*w;
458  outputValues(8, i0, 5) = y*y*w*w*w;
459  outputValues(8, i0, 6) = 0.0;
460  outputValues(8, i0, 7) = x*w*w;
461  outputValues(8, i0, 8) = 2.0*x*y*w*w*w;
462  outputValues(8, i0, 9) = 3.0*x*y*y*w*w*w*w ;
463 
464  outputValues(9, i0, 0) = 0.0;
465  outputValues(9, i0, 1) = 0.0;
466  outputValues(9, i0, 2) = 0.0;
467  outputValues(9, i0, 3) = 0.0;
468  outputValues(9, i0, 4) = w*w;
469  outputValues(9, i0, 5) = 2.0*y*w*w*w;
470  outputValues(9, i0, 6) = 0.0;
471  outputValues(9, i0, 7) = 0.0;
472  outputValues(9, i0, 8) = 2.0*x*w*w*w;
473  outputValues(9, i0, 9) = 6.0*x*y*w*w*w*w;
474 
475  outputValues(10,i0, 0) = 0.0;
476  outputValues(10,i0, 1) = 0.0;
477  outputValues(10,i0, 2) = 0.0;
478  outputValues(10,i0, 3) = 0.0;
479  outputValues(10,i0, 4) = -w*w;
480  outputValues(10,i0, 5) = -2.0*y*w*w*w;
481  outputValues(10,i0, 6) = 0.0;
482  outputValues(10,i0, 7) = 0.0;
483  outputValues(10,i0, 8) = -2.0*x*w*w*w;
484  outputValues(10,i0, 9) = -6.0*x*y*w*w*w*w;
485 
486  outputValues(11,i0, 0) = 0.0;
487  outputValues(11,i0, 1) = 0.0;
488  outputValues(11,i0, 2) = 0.0;
489  outputValues(11,i0, 3) = 0.0;
490  outputValues(11,i0, 4) = w*w;
491  outputValues(11,i0, 5) = 2.0*y*w*w*w;
492  outputValues(11,i0, 6) = 0.0;
493  outputValues(11,i0, 7) = 0.0;
494  outputValues(11,i0, 8) = 2.0*x*w*w*w;
495  outputValues(11,i0, 9) = 6.0*x*y*w*w*w*w;
496 
497  outputValues(12,i0, 0) = 0.0;
498  outputValues(12,i0, 1) = 0.0;
499  outputValues(12,i0, 2) = 0.0;
500  outputValues(12,i0, 3) = 0.0;
501  outputValues(12,i0, 4) = -w*w;
502  outputValues(12,i0, 5) = -2.0*y*w*w*w;
503  outputValues(12,i0, 6) = 0.0;
504  outputValues(12,i0, 7) = 0.0;
505  outputValues(12,i0, 8) = -2.0*x*w*w*w;
506  outputValues(12,i0, 9) = -6.0*x*y*w*w*w*w;
507 
508  }
509  break;
510 
511  case OPERATOR_D4:
512  case OPERATOR_D5:
513  case OPERATOR_D6:
514  case OPERATOR_D7:
515  case OPERATOR_D8:
516  case OPERATOR_D9:
517  case OPERATOR_D10:
518  {
519  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, DkCardinality)
520  int DkCardinality = Intrepid::getDkCardinality(operatorType,
521  this -> basisCellTopology_.getDimension() );
522  for(int dofOrd = 0; dofOrd < this -> basisCardinality_; dofOrd++) {
523  for (int i0 = 0; i0 < dim0; i0++) {
524  for(int dkOrd = 0; dkOrd < DkCardinality; dkOrd++){
525  outputValues(dofOrd, i0, dkOrd) = 0.0;
526  }
527  }
528  }
529  }
530  break;
531 
532  default:
533  TEUCHOS_TEST_FOR_EXCEPTION( !( Intrepid::isValidOperator(operatorType) ), std::invalid_argument,
534  ">>> ERROR (Basis_HGRAD_PYR_I2_FEM): Invalid operator type");
535  }
536 }
537 
538 
539 
540 template<class Scalar, class ArrayScalar>
541 void Basis_HGRAD_PYR_I2_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar& outputValues,
542  const ArrayScalar & inputPoints,
543  const ArrayScalar & cellVertices,
544  const EOperator operatorType) const {
545  TEUCHOS_TEST_FOR_EXCEPTION( (true), std::logic_error,
546  ">>> ERROR (Basis_HGRAD_PYR_I2_FEM): FEM Basis calling an FVD member function");
547 }
548 }// namespace Intrepid
549 #endif
Intrepid::Basis_HGRAD_PYR_I2_FEM::getValues
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
FEM basis evaluation on a reference Pyramid cell.
Definition: Intrepid_HGRAD_PYR_I2_FEMDef.hpp:104
Intrepid::Basis_HGRAD_PYR_I2_FEM::initializeTags
void initializeTags()
Initializes tagToOrdinal_ and ordinalToTag_ lookup arrays.
Definition: Intrepid_HGRAD_PYR_I2_FEMDef.hpp:66
Intrepid::Basis_HGRAD_PYR_I2_FEM::Basis_HGRAD_PYR_I2_FEM
Basis_HGRAD_PYR_I2_FEM()
Constructor.
Definition: Intrepid_HGRAD_PYR_I2_FEMDef.hpp:54
Generated by   doxygen 1.8.5