Intrepid
Intrepid_HGRAD_HEX_C2_FEMDef.hpp
Go to the documentation of this file.
1#ifndef INTREPID_HGRAD_HEX_C2_FEMDEF_HPP
2#define INTREPID_HGRAD_HEX_C2_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
51namespace Intrepid {
52
53
54template<class Scalar, class ArrayScalar>
56 {
57 this -> basisCardinality_ = 27;
58 this -> basisDegree_ = 2;
59 this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
60 this -> basisType_ = BASIS_FEM_DEFAULT;
61 this -> basisCoordinates_ = COORDINATES_CARTESIAN;
62 this -> basisTagsAreSet_ = false;
63 }
64
65
66
67template<class Scalar, class ArrayScalar>
69
70 // Basis-dependent intializations
71 int tagSize = 4; // size of DoF tag, i.e., number of fields in the tag
72 int posScDim = 0; // position in the tag, counting from 0, of the subcell dim
73 int posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
74 int posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
75
76 // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
77 int tags[] = { 0, 0, 0, 1, // Nodes 0 to 7 follow vertex order of the topology
78 0, 1, 0, 1,
79 0, 2, 0, 1,
80 0, 3, 0, 1,
81 0, 4, 0, 1,
82 0, 5, 0, 1,
83 0, 6, 0, 1,
84 0, 7, 0, 1,
85 1, 0, 0, 1, // Node 8 -> edge 0
86 1, 1, 0, 1, // Node 9 -> edge 1
87 1, 2, 0, 1, // Node 10 -> edge 2
88 1, 3, 0, 1, // Node 11 -> edge 3
89 1, 8, 0, 1, // Node 12 -> edge 8
90 1, 9, 0, 1, // Node 13 -> edge 9
91 1,10, 0, 1, // Node 14 -> edge 10
92 1,11, 0, 1, // Node 15 -> edge 11
93 1, 4, 0, 1, // Node 16 -> edge 4
94 1, 5, 0, 1, // Node 17 -> edge 5
95 1, 6, 0, 1, // Node 18 -> edge 6
96 1, 7, 0, 1, // Node 19 -> edge 7
97 3, 0, 0, 1, // Node 20 -> Hexahedron
98 2, 4, 0, 1, // Node 21 -> face 4
99 2, 5, 0, 1, // Node 22 -> face 5
100 2, 3, 0, 1, // Node 23 -> face 3
101 2, 1, 0, 1, // Node 24 -> face 1
102 2, 0, 0, 1, // Node 25 -> face 0
103 2, 2, 0, 1, // Node 26 -> face 2
104 };
105
106 // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
107 Intrepid::setOrdinalTagData(this -> tagToOrdinal_,
108 this -> ordinalToTag_,
109 tags,
110 this -> basisCardinality_,
111 tagSize,
112 posScDim,
113 posScOrd,
114 posDfOrd);
115}
116
117
118
119template<class Scalar, class ArrayScalar>
121 const ArrayScalar & inputPoints,
122 const EOperator operatorType) const {
123
124 // Verify arguments
125#ifdef HAVE_INTREPID_DEBUG
126 Intrepid::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
127 inputPoints,
128 operatorType,
129 this -> getBaseCellTopology(),
130 this -> getCardinality() );
131#endif
132
133 // Number of evaluation points = dim 0 of inputPoints
134 int dim0 = inputPoints.dimension(0);
135
136 // Temporaries: (x,y,z) coordinates of the evaluation point
137 Scalar x = 0.0;
138 Scalar y = 0.0;
139 Scalar z = 0.0;
140
141 switch (operatorType) {
142
143 case OPERATOR_VALUE:
144 for (int i0 = 0; i0 < dim0; i0++) {
145 x = inputPoints(i0, 0);
146 y = inputPoints(i0, 1);
147 z = inputPoints(i0, 2);
148
149 // outputValues is a rank-2 array with dimensions (basisCardinality_, dim0)
150 outputValues( 0, i0) = 0.125*(-1. + x)*x*(-1. + y)*y*(-1. + z)*z;
151 outputValues( 1, i0) = 0.125*x*(1.+ x)*(-1. + y)*y*(-1. + z)*z;
152 outputValues( 2, i0) = 0.125*x*(1.+ x)*y*(1.+ y)*(-1. + z)*z;
153 outputValues( 3, i0) = 0.125*(-1. + x)*x*y*(1.+ y)*(-1. + z)*z;
154 outputValues( 4, i0) = 0.125*(-1. + x)*x*(-1. + y)*y*z*(1.+ z);
155 outputValues( 5, i0) = 0.125*x*(1.+ x)*(-1. + y)*y*z*(1.+ z);
156 outputValues( 6, i0) = 0.125*x*(1.+ x)*y*(1.+ y)*z*(1.+ z);
157 outputValues( 7, i0) = 0.125*(-1. + x)*x*y*(1.+ y)*z*(1.+ z);
158 outputValues( 8, i0) = 0.25*(1. - x)*(1. + x)*(-1. + y)*y*(-1. + z)*z;
159 outputValues( 9, i0) = 0.25*x*(1.+ x)*(1. - y)*(1. + y)*(-1. + z)*z;
160 outputValues(10, i0) = 0.25*(1. - x)*(1. + x)*y*(1.+ y)*(-1. + z)*z;
161 outputValues(11, i0) = 0.25*(-1. + x)*x*(1. - y)*(1. + y)*(-1. + z)*z;
162 outputValues(12, i0) = 0.25*(-1. + x)*x*(-1. + y)*y*(1. - z)*(1. + z);
163 outputValues(13, i0) = 0.25*x*(1.+ x)*(-1. + y)*y*(1. - z)*(1. + z);
164 outputValues(14, i0) = 0.25*x*(1.+ x)*y*(1.+ y)*(1. - z)*(1. + z);
165 outputValues(15, i0) = 0.25*(-1. + x)*x*y*(1.+ y)*(1. - z)*(1. + z);
166 outputValues(16, i0) = 0.25*(1. - x)*(1. + x)*(-1. + y)*y*z*(1.+ z);
167 outputValues(17, i0) = 0.25*x*(1.+ x)*(1. - y)*(1. + y)*z*(1.+ z);
168 outputValues(18, i0) = 0.25*(1. - x)*(1. + x)*y*(1.+ y)*z*(1.+ z);
169 outputValues(19, i0) = 0.25*(-1. + x)*x*(1. - y)*(1. + y)*z*(1.+ z);
170 outputValues(20, i0) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
171 outputValues(21, i0) = 0.5*(1. - x)*(1. + x)*(1. - y)*(1. + y)*(-1. + z)*z;
172 outputValues(22, i0) = 0.5*(1. - x)*(1. + x)*(1. - y)*(1. + y)*z*(1.+ z);
173 outputValues(23, i0) = 0.5*(-1. + x)*x*(1. - y)*(1. + y)*(1. - z)*(1. + z);
174 outputValues(24, i0) = 0.5*x*(1.+ x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
175 outputValues(25, i0) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y*(1. - z)*(1. + z);
176 outputValues(26, i0) = 0.5*(1. - x)*(1. + x)*y*(1.+ y)*(1. - z)*(1. + z);
177 }
178 break;
179
180 case OPERATOR_GRAD:
181 case OPERATOR_D1:
182 for (int i0 = 0; i0 < dim0; i0++) {
183 x = inputPoints(i0,0);
184 y = inputPoints(i0,1);
185 z = inputPoints(i0,2);
186
187 // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
188 outputValues(0, i0, 0) = (-0.125 + 0.25*x)*(-1. + y)*y*(-1. + z)*z;
189 outputValues(0, i0, 1) = (-1. + x)*x*(-0.125 + 0.25*y)*(-1. + z)*z;
190 outputValues(0, i0, 2) = (-1. + x)*x*(-1. + y)*y*(-0.125 + 0.25*z);
191
192 outputValues(1, i0, 0) = (0.125 + 0.25*x)*(-1. + y)*y*(-1. + z)*z;
193 outputValues(1, i0, 1) = x*(1. + x)*(-0.125 + 0.25*y)*(-1. + z)*z;
194 outputValues(1, i0, 2) = x*(1. + x)*(-1. + y)*y*(-0.125 + 0.25*z);
195
196 outputValues(2, i0, 0) = (0.125 + 0.25*x)*y*(1. + y)*(-1. + z)*z;
197 outputValues(2, i0, 1) = x*(1. + x)*(0.125 + 0.25*y)*(-1. + z)*z;
198 outputValues(2, i0, 2) = x*(1. + x)*y*(1. + y)*(-0.125 + 0.25*z);
199
200 outputValues(3, i0, 0) = (-0.125 + 0.25*x)*y*(1. + y)*(-1. + z)*z;
201 outputValues(3, i0, 1) = (-1. + x)*x*(0.125 + 0.25*y)*(-1. + z)*z;
202 outputValues(3, i0, 2) = (-1. + x)*x*y*(1. + y)*(-0.125 + 0.25*z);
203
204 outputValues(4, i0, 0) = (-0.125 + 0.25*x)*(-1. + y)*y*z*(1. + z);
205 outputValues(4, i0, 1) = (-1. + x)*x*(-0.125 + 0.25*y)*z*(1. + z);
206 outputValues(4, i0, 2) = (-1. + x)*x*(-1. + y)*y*(0.125 + 0.25*z);
207
208 outputValues(5, i0, 0) = (0.125 + 0.25*x)*(-1. + y)*y*z*(1. + z);
209 outputValues(5, i0, 1) = x*(1. + x)*(-0.125 + 0.25*y)*z*(1. + z);
210 outputValues(5, i0, 2) = x*(1. + x)*(-1. + y)*y*(0.125 + 0.25*z);
211
212 outputValues(6, i0, 0) = (0.125 + 0.25*x)*y*(1. + y)*z*(1. + z);
213 outputValues(6, i0, 1) = x*(1. + x)*(0.125 + 0.25*y)*z*(1. + z);
214 outputValues(6, i0, 2) = x*(1. + x)*y*(1. + y)*(0.125 + 0.25*z);
215
216 outputValues(7, i0, 0) = (-0.125 + 0.25*x)*y*(1. + y)*z*(1. + z);
217 outputValues(7, i0, 1) = (-1. + x)*x*(0.125 + 0.25*y)*z*(1. + z);
218 outputValues(7, i0, 2) = (-1. + x)*x*y*(1. + y)*(0.125 + 0.25*z);
219
220 outputValues(8, i0, 0) = -0.5*x*(-1. + y)*y*(-1. + z)*z;
221 outputValues(8, i0, 1) = (1. - x)*(1. + x)*(-0.25 + 0.5*y)*(-1. + z)*z;
222 outputValues(8, i0, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(-0.25 + 0.5*z);
223
224 outputValues(9, i0, 0) = (0.25 + 0.5*x)*(1. - y)*(1. + y)*(-1. + z)*z;
225 outputValues(9, i0, 1) = x*(1. + x)*(-0.5*y)*(-1. + z)*z;
226 outputValues(9, i0, 2) = x*(1. + x)*(1. - y)*(1. + y)*(-0.25 + 0.5*z);
227
228 outputValues(10,i0, 0) = -0.5*x*y*(1. + y)*(-1. + z)*z;
229 outputValues(10,i0, 1) = (1. - x)*(1. + x)*(0.25 + 0.5*y)*(-1. + z)*z;
230 outputValues(10,i0, 2) = (1. - x)*(1. + x)*y*(1. + y)*(-0.25 + 0.5*z);
231
232 outputValues(11,i0, 0) = (-0.25 + 0.5*x)*(1. - y)*(1. + y)*(-1. + z)*z;
233 outputValues(11,i0, 1) = (-1. + x)*x*(-0.5*y)*(-1. + z)*z;
234 outputValues(11,i0, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(-0.25 + 0.5*z);
235
236 outputValues(12,i0, 0) = (-0.25 + 0.5*x)*(-1. + y)*y*(1. - z)*(1. + z);
237 outputValues(12,i0, 1) = (-1. + x)*x*(-0.25 + 0.5*y)*(1. - z)*(1. + z);
238 outputValues(12,i0, 2) = (-1. + x)*x*(-1. + y)*y*(-0.5*z);
239
240 outputValues(13,i0, 0) = (0.25 + 0.5*x)*(-1. + y)*y*(1. - z)*(1. + z);
241 outputValues(13,i0, 1) = x*(1. + x)*(-0.25 + 0.5*y)*(1. - z)*(1. + z);
242 outputValues(13,i0, 2) = x*(1. + x)*(-1. + y)*y*(-0.5*z);
243
244 outputValues(14,i0, 0) = (0.25 + 0.5*x)*y*(1. + y)*(1. - z)*(1. + z);
245 outputValues(14,i0, 1) = x*(1. + x)*(0.25 + 0.5*y)*(1. - z)*(1. + z);
246 outputValues(14,i0, 2) = x*(1. + x)*y*(1. + y)*(-0.5*z);
247
248 outputValues(15,i0, 0) = (-0.25 + 0.5*x)*y*(1. + y)*(1. - z)*(1. + z);
249 outputValues(15,i0, 1) = (-1. + x)*x*(0.25 + 0.5*y)*(1. - z)*(1. + z);
250 outputValues(15,i0, 2) = (-1. + x)*x*y*(1. + y)*(-0.5*z);
251
252 outputValues(16,i0, 0) = -0.5*x*(-1. + y)*y*z*(1. + z);
253 outputValues(16,i0, 1) = (1. - x)*(1. + x)*(-0.25 + 0.5*y)*z*(1. + z);
254 outputValues(16,i0, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(0.25 + 0.5*z);
255
256 outputValues(17,i0, 0) = (0.25 + 0.5*x)*(1. - y)*(1. + y)*z*(1. + z);
257 outputValues(17,i0, 1) = x*(1. + x)*(-0.5*y)*z*(1. + z);
258 outputValues(17,i0, 2) = x*(1. + x)*(1. - y)*(1. + y)*(0.25 + 0.5*z);
259
260 outputValues(18,i0, 0) = -0.5*x*y*(1. + y)*z*(1. + z);
261 outputValues(18,i0, 1) = (1. - x)*(1. + x)*(0.25 + 0.5*y)*z*(1. + z);
262 outputValues(18,i0, 2) = (1. - x)*(1. + x)*y*(1. + y)*(0.25 + 0.5*z);
263
264 outputValues(19,i0, 0) = (-0.25 + 0.5*x)*(1. - y)*(1. + y)*z*(1. + z);
265 outputValues(19,i0, 1) = (-1. + x)*x*(-0.5*y)*z*(1. + z);
266 outputValues(19,i0, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(0.25 + 0.5*z);
267
268 outputValues(20,i0, 0) = -2.*x*(1. - y)*(1. + y)*(1. - z)*(1. + z);
269 outputValues(20,i0, 1) = (1. - x)*(1. + x)*(-2.*y)*(1. - z)*(1. + z);
270 outputValues(20,i0, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(-2.*z);
271
272 outputValues(21,i0, 0) = -x*(1. - y)*(1. + y)*(-1. + z)*z;
273 outputValues(21,i0, 1) = (1. - x)*(1. + x)*(-y)*(-1. + z)*z;
274 outputValues(21,i0, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(-0.5 + z);
275
276 outputValues(22,i0, 0) = -x*(1. - y)*(1. + y)*z*(1. + z);
277 outputValues(22,i0, 1) = (1. - x)*(1. + x)*(-y)*z*(1. + z);
278 outputValues(22,i0, 2) = (1. - x)*(1. + x)*(1. - y)*(1. + y)*(0.5 + z);
279
280 outputValues(23,i0, 0) = (-0.5 + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
281 outputValues(23,i0, 1) = (-1. + x)*x*(-y)*(1. - z)*(1. + z);
282 outputValues(23,i0, 2) = (-1. + x)*x*(1. - y)*(1. + y)*(-z);
283
284 outputValues(24,i0, 0) = (0.5 + x)*(1. - y)*(1. + y)*(1. - z)*(1. + z);
285 outputValues(24,i0, 1) = x*(1. + x)*(-y)*(1. - z)*(1. + z);
286 outputValues(24,i0, 2) = x*(1. + x)*(1. - y)*(1. + y)*(-z);
287
288 outputValues(25,i0, 0) = -x*(-1. + y)*y*(1. - z)*(1. + z);
289 outputValues(25,i0, 1) = (1. - x)*(1. + x)*(-0.5 + y)*(1. - z)*(1. + z);
290 outputValues(25,i0, 2) = (1. - x)*(1. + x)*(-1. + y)*y*(-z);
291
292 outputValues(26,i0, 0) = -x*y*(1. + y)*(1. - z)*(1. + z);
293 outputValues(26,i0, 1) = (1. - x)*(1. + x)*(0.5 + y)*(1. - z)*(1. + z);
294 outputValues(26,i0, 2) = (1. - x)*(1. + x)*y*(1. + y)*(-z);
295 }
296 break;
297
298 case OPERATOR_CURL:
299 TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
300 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
301 break;
302
303 case OPERATOR_DIV:
304 TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
305 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
306 break;
307
308 case OPERATOR_D2:
309 for (int i0 = 0; i0 < dim0; i0++) {
310 x = inputPoints(i0,0);
311 y = inputPoints(i0,1);
312 z = inputPoints(i0,2);
313
314 // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, D2Cardinality = 6)
315 outputValues(0, i0, 0) = 0.25*(-1. + y)*y*(-1. + z)*z;
316 outputValues(0, i0, 1) = (-0.125 + y*(0.25 - 0.25*z) + x*(0.25 + y*(-0.5 + 0.5*z) - 0.25*z) + 0.125*z)*z;
317 outputValues(0, i0, 2) = y*(-0.125 + x*(0.25 + y*(-0.25 + 0.5*z) - 0.5*z) + y*(0.125 - 0.25*z) + 0.25*z);
318 outputValues(0, i0, 3) = 0.25*(-1. + x)*x*(-1. + z)*z;
319 outputValues(0, i0, 4) = x*(-0.125 + y*(0.25 - 0.5*z) + x*(0.125 + y*(-0.25 + 0.5*z) - 0.25*z) + 0.25*z);
320 outputValues(0, i0, 5) = 0.25*(-1. + x)*x*(-1. + y)*y;
321
322 outputValues(1, i0, 0) = 0.25*(-1. + y)*y*(-1. + z)*z;
323 outputValues(1, i0, 1) = (0.125 + x*(0.25 + y*(-0.5 + 0.5*z) - 0.25*z) + y*(-0.25 + 0.25*z) - 0.125*z)*z;
324 outputValues(1, i0, 2) = y*(0.125 + x*(0.25 + y*(-0.25 + 0.5*z) - 0.5*z) + y*(-0.125 + 0.25*z) - 0.25*z);
325 outputValues(1, i0, 3) = 0.25*x*(1 + x)*(-1. + z)*z;
326 outputValues(1, i0, 4) = x*(1. + x)*(0.125 + y*(-0.25 + 0.5*z) - 0.25*z);
327 outputValues(1, i0, 5) = 0.25*x*(1 + x)*(-1. + y)*y;
328
329 outputValues(2, i0, 0) = 0.25*y*(1 + y)*(-1. + z)*z;
330 outputValues(2, i0, 1) = (0.125 + x*(0.25 + 0.5*y) + 0.25*y)*(-1. + z)*z;
331 outputValues(2, i0, 2) = y*(1. + y)*(-0.125 + x*(-0.25 + 0.5*z) + 0.25*z);
332 outputValues(2, i0, 3) = 0.25*x*(1 + x)*(-1. + z)*z;
333 outputValues(2, i0, 4) = x*(1. + x)*(-0.125 + y*(-0.25 + 0.5*z) + 0.25*z);
334 outputValues(2, i0, 5) = 0.25*x*(1 + x)*y*(1 + y);
335
336 outputValues(3, i0, 0) = 0.25*y*(1 + y)*(-1. + z)*z;
337 outputValues(3, i0, 1) = (0.125 + y*(0.25 - 0.25*z) + x*(-0.25 + y*(-0.5 + 0.5*z) + 0.25*z) - 0.125*z)*z;
338 outputValues(3, i0, 2) = y*(1. + y)*(0.125 + x*(-0.25 + 0.5*z) - 0.25*z);
339 outputValues(3, i0, 3) = 0.25*(-1. + x)*x*(-1. + z)*z;
340 outputValues(3, i0, 4) = x*(0.125 + y*(0.25 - 0.5*z) + x*(-0.125 + y*(-0.25 + 0.5*z) + 0.25*z) - 0.25*z);
341 outputValues(3, i0, 5) = 0.25*(-1. + x)*x*y*(1 + y);
342
343 outputValues(4, i0, 0) = 0.25*(-1. + y)*y*z*(1 + z);
344 outputValues(4, i0, 1) = (0.125 + x*(-0.25 + 0.5*y) - 0.25*y)*z*(1. + z);
345 outputValues(4, i0, 2) = y*(0.125 + x*(-0.25 + y*(0.25 + 0.5*z) - 0.5*z) + y*(-0.125 - 0.25*z) + 0.25*z);
346 outputValues(4, i0, 3) = 0.25*(-1. + x)*x*z*(1 + z);
347 outputValues(4, i0, 4) = x*(0.125 + y*(-0.25 - 0.5*z) + x*(-0.125 + y*(0.25 + 0.5*z) - 0.25*z) + 0.25*z);
348 outputValues(4, i0, 5) = 0.25*(-1. + x)*x*(-1. + y)*y;
349
350 outputValues(5, i0, 0) = 0.25*(-1. + y)*y*z*(1 + z);
351 outputValues(5, i0, 1) = (-0.125 + x*(-0.25 + 0.5*y) + 0.25*y)*z*(1. + z);
352 outputValues(5, i0, 2) = (-1. + y)*y*(0.125 + x*(0.25 + 0.5*z) + 0.25*z);
353 outputValues(5, i0, 3) = 0.25*x*(1 + x)*z*(1 + z);
354 outputValues(5, i0, 4) = x*(1. + x)*(-0.125 + y*(0.25 + 0.5*z) - 0.25*z);
355 outputValues(5, i0, 5) = 0.25*x*(1 + x)*(-1. + y)*y;
356
357 outputValues(6, i0, 0) = 0.25*y*(1 + y)*z*(1 + z);
358 outputValues(6, i0, 1) = (0.125 + x*(0.25 + 0.5*y) + 0.25*y)*z*(1. + z);
359 outputValues(6, i0, 2) = y*(1. + y)*(0.125 + x*(0.25 + 0.5*z) + 0.25*z);
360 outputValues(6, i0, 3) = 0.25*x*(1 + x)*z*(1 + z);
361 outputValues(6, i0, 4) = x*(1. + x)*(0.125 + y*(0.25 + 0.5*z) + 0.25*z);
362 outputValues(6, i0, 5) = 0.25*x*(1 + x)*y*(1 + y);
363
364 outputValues(7, i0, 0) = 0.25*y*(1 + y)*z*(1 + z);
365 outputValues(7, i0, 1) = (-0.125 + x*(0.25 + 0.5*y) - 0.25*y)*z*(1. + z);
366 outputValues(7, i0, 2) = y*(1. + y)*(-0.125 + x*(0.25 + 0.5*z) - 0.25*z);
367 outputValues(7, i0, 3) = 0.25*(-1. + x)*x*z*(1 + z);
368 outputValues(7, i0, 4) = (-1. + x)*x*(0.125 + y*(0.25 + 0.5*z) + 0.25*z);
369 outputValues(7, i0, 5) = 0.25*(-1. + x)*x*y*(1 + y);
370
371 outputValues(8, i0, 0) = -0.5*(-1. + y)*y*(-1. + z)*z;
372 outputValues(8, i0, 1) = (0. + x*(-0.5 + y))*z + (x*(0.5 - y) )*(z*z);
373 outputValues(8, i0, 2) = (y*y)*(x*(0.5 - z) ) + y*(x*(-0.5 + z));
374 outputValues(8, i0, 3) = 0.5*(1. - x)*(1. + x)*(-1. + z)*z;
375 outputValues(8, i0, 4) = 0.25 + (x*x)*(-0.25 + y*(0.5 - z) + 0.5*z) - 0.5*z + y*(-0.5 + z);
376 outputValues(8, i0, 5) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y;
377
378 outputValues(9, i0, 0) = 0.5*(1. - y)*(1. + y)*(-1. + z)*z;
379 outputValues(9, i0, 1) = (0.5*y + x*(y))*z + (x*(-y) - 0.5*y)*(z*z);
380 outputValues(9, i0, 2) = -0.25 + (y*y)*(0.25 - 0.5*z) + 0.5*z + x*(-0.5 + (y*y)*(0.5 - z) + z);
381 outputValues(9, i0, 3) = -0.5*x*(1 + x)*(-1. + z)*z;
382 outputValues(9, i0, 4) = x*(y*(0.5 - z) ) + (x*x)*(y*(0.5 - z) );
383 outputValues(9, i0, 5) = 0.5*x*(1 + x)*(1. - y)*(1. + y);
384
385 outputValues(10,i0, 0) = -0.5*y*(1 + y)*(-1. + z)*z;
386 outputValues(10,i0, 1) = (0. + x*(0.5 + y))*z + (x*(-0.5 - y) )*(z*z);
387 outputValues(10,i0, 2) = y*(x*(0.5 - z) ) + (y*y)*(x*(0.5 - z) );
388 outputValues(10,i0, 3) = 0.5*(1. - x)*(1. + x)*(-1. + z)*z;
389 outputValues(10,i0, 4) = -0.25 + (x*x)*(0.25 + y*(0.5 - z) - 0.5*z) + 0.5*z + y*(-0.5 + z);
390 outputValues(10,i0, 5) = 0.5*(1. - x)*(1. + x)*y*(1 + y);
391
392 outputValues(11,i0, 0) = 0.5*(1. - y)*(1. + y)*(-1. + z)*z;
393 outputValues(11,i0, 1) = (-0.5*y + x*(y))*z + (x*(-y) + 0.5*y)*(z*z);
394 outputValues(11,i0, 2) = 0.25 + (y*y)*(-0.25 + 0.5*z) - 0.5*z + x*(-0.5 + (y*y)*(0.5 - z) + z);
395 outputValues(11,i0, 3) = -0.5*(-1. + x)*x*(-1. + z)*z;
396 outputValues(11,i0, 4) = (x*x)*(y*(0.5 - z) ) + x*(y*(-0.5 + z));
397 outputValues(11,i0, 5) = 0.5*(-1. + x)*x*(1. - y)*(1. + y);
398
399 outputValues(12,i0, 0) = 0.5*(-1. + y)*y*(1. - z)*(1. + z);
400 outputValues(12,i0, 1) = 0.25 - 0.25*(z*z) + y*(-0.5 + 0.5*(z*z)) + x*(-0.5 + 0.5*(z*z) + y*(1. - (z*z)));
401 outputValues(12,i0, 2) = (y*y)*(x*(-z) + 0.5*z) + y*(-0.5*z + x*(z));
402 outputValues(12,i0, 3) = 0.5*(-1. + x)*x*(1. - z)*(1. + z);
403 outputValues(12,i0, 4) = (x*x)*(y*(-z) + 0.5*z) + x*(-0.5*z + y*(z));
404 outputValues(12,i0, 5) = -0.5*(-1. + x)*x*(-1. + y)*y;
405
406 outputValues(13,i0, 0) = 0.5*(-1. + y)*y*(1. - z)*(1. + z);
407 outputValues(13,i0, 1) = -0.25 + 0.25*(z*z) + y*(0.5 - 0.5*(z*z)) + x*(-0.5 + 0.5*(z*z) + y*(1. - (z*z)));
408 outputValues(13,i0, 2) = (y*y)*(x*(-z) - 0.5*z) + y*(0.5*z + x*(z));
409 outputValues(13,i0, 3) = 0.5*x*(1 + x)*(1. - z)*(1. + z);
410 outputValues(13,i0, 4) = x*(y*(-z) + 0.5*z) + (x*x)*(y*(-z) + 0.5*z);
411 outputValues(13,i0, 5) = -0.5*x*(1 + x)*(-1. + y)*y;
412
413 outputValues(14,i0, 0) = 0.5*y*(1 + y)*(1. - z)*(1. + z);
414 outputValues(14,i0, 1) = 0.25 - 0.25*(z*z) + y*(0.5 - 0.5*(z*z)) + x*(0.5 - 0.5*(z*z) + y*(1. - (z*z)));
415 outputValues(14,i0, 2) = y*(x*(-z) - 0.5*z) + (y*y)*(x*(-z) - 0.5*z);
416 outputValues(14,i0, 3) = 0.5*x*(1 + x)*(1. - z)*(1. + z);
417 outputValues(14,i0, 4) = x*(y*(-z) - 0.5*z) + (x*x)*(y*(-z) - 0.5*z);
418 outputValues(14,i0, 5) = -0.5*x*(1 + x)*y*(1 + y);
419
420 outputValues(15,i0, 0) = 0.5*y*(1 + y)*(1. - z)*(1. + z);
421 outputValues(15,i0, 1) = -0.25 + 0.25*(z*z) + y*(-0.5 + 0.5*(z*z)) + x*(0.5 - 0.5*(z*z) + y*(1. - (z*z)));
422 outputValues(15,i0, 2) = y*(x*(-z) + 0.5*z) + (y*y)*(x*(-z) + 0.5*z);
423 outputValues(15,i0, 3) = 0.5*(-1. + x)*x*(1. - z)*(1. + z);
424 outputValues(15,i0, 4) = (x*x)*(y*(-z) - 0.5*z) + x*(0.5*z + y*(z));
425 outputValues(15,i0, 5) = -0.5*(-1. + x)*x*y*(1 + y);
426
427 outputValues(16,i0, 0) = -0.5*(-1. + y)*y*z*(1 + z);
428 outputValues(16,i0, 1) = (x*(0.5 - y) )*z + (x*(0.5 - y) )*(z*z);
429 outputValues(16,i0, 2) = (y*y)*(x*(-0.5 - z) ) + y*(x*(0.5 + z));
430 outputValues(16,i0, 3) = 0.5*(1. - x)*(1. + x)*z*(1 + z);
431 outputValues(16,i0, 4) = -0.25 + (x*x)*(0.25 + y*(-0.5 - z) + 0.5*z) - 0.5*z + y*(0.5 + z);
432 outputValues(16,i0, 5) = 0.5*(1. - x)*(1. + x)*(-1. + y)*y;
433
434 outputValues(17,i0, 0) = 0.5*(1. - y)*(1. + y)*z*(1 + z);
435 outputValues(17,i0, 1) = (x*(-y) - 0.5*y)*z + (x*(-y) - 0.5*y)*(z*z);
436 outputValues(17,i0, 2) = 0.25 + (y*y)*(-0.25 - 0.5*z) + 0.5*z + x*(0.5 + (y*y)*(-0.5 - z) + z);
437 outputValues(17,i0, 3) = -0.5*x*(1 + x)*z*(1 + z);
438 outputValues(17,i0, 4) = x*(y*(-0.5 - z) ) + (x*x)*(y*(-0.5 - z) );
439 outputValues(17,i0, 5) = 0.5*x*(1 + x)*(1. - y)*(1. + y);
440
441 outputValues(18,i0, 0) = -0.5*y*(1 + y)*z*(1 + z);
442 outputValues(18,i0, 1) = (x*(-0.5 - y) )*z + (x*(-0.5 - y) )*(z*z);
443 outputValues(18,i0, 2) = y*(x*(-0.5 - z) ) + (y*y)*(x*(-0.5 - z) );
444 outputValues(18,i0, 3) = 0.5*(1. - x)*(1. + x)*z*(1 + z);
445 outputValues(18,i0, 4) = 0.25 + (x*x)*(-0.25 + y*(-0.5 - z) - 0.5*z) + 0.5*z + y*(0.5 + z);
446 outputValues(18,i0, 5) = 0.5*(1. - x)*(1. + x)*y*(1 + y);
447
448 outputValues(19,i0, 0) = 0.5*(1. - y)*(1. + y)*z*(1 + z);
449 outputValues(19,i0, 1) = (x*(-y) + 0.5*y)*z + (x*(-y) + 0.5*y)*(z*z);
450 outputValues(19,i0, 2) = -0.25 + (y*y)*(0.25 + 0.5*z) - 0.5*z + x*(0.5 + (y*y)*(-0.5 - z) + z);
451 outputValues(19,i0, 3) = -0.5*(-1. + x)*x*z*(1 + z);
452 outputValues(19,i0, 4) = (x*x)*(y*(-0.5 - z) ) + x*(y*(0.5 + z));
453 outputValues(19,i0, 5) = 0.5*(-1. + x)*x*(1. - y)*(1. + y);
454
455 outputValues(20,i0, 0) = -2.*(1. - y)*(1. + y)*(1. - z)*(1. + z);
456 outputValues(20,i0, 1) = -4.*x*y*(-1. + z*z);
457 outputValues(20,i0, 2) = x*((y*y)*(-4.*z) + 4.*z);
458 outputValues(20,i0, 3) = -2.*(1. - x)*(1. + x)*(1. - z)*(1. + z);
459 outputValues(20,i0, 4) = (x*x)*(y*(-4.*z) ) + y*(4.*z);
460 outputValues(20,i0, 5) = -2.*(1. - x)*(1. + x)*(1. - y)*(1. + y);
461
462 outputValues(21,i0, 0) = -(1. - y)*(1. + y)*(-1. + z)*z;
463 outputValues(21,i0, 1) = (x*(-2.*y) )*z + (0. + x*(2.*y))*(z*z);
464 outputValues(21,i0, 2) = x*(1. - 2.*z + (y*y)*(-1. + 2.*z));
465 outputValues(21,i0, 3) = -(1. - x)*(1. + x)*(-1. + z)*z;
466 outputValues(21,i0, 4) = y*(1. - 2.*z) + (x*x)*(y*(-1. + 2.*z));
467 outputValues(21,i0, 5) = (1. - x)*(1. + x)*(1. - y)*(1. + y);
468
469 outputValues(22,i0, 0) = -(1. - y)*(1. + y)*z*(1 + z);
470 outputValues(22,i0, 1) = (0. + x*(2.*y))*z + (0. + x*(2.*y))*(z*z);
471 outputValues(22,i0, 2) = x*(-1. - 2.*z + (y*y)*(1. + 2.*z));
472 outputValues(22,i0, 3) = -(1. - x)*(1. + x)*z*(1 + z);
473 outputValues(22,i0, 4) = y*(-1. - 2.*z) + (x*x)*(y*(1. + 2.*z));
474 outputValues(22,i0, 5) = (1. - x)*(1. + x)*(1. - y)*(1. + y);
475
476 outputValues(23,i0, 0) = (1. - y)*(1. + y)*(1. - z)*(1. + z);
477 outputValues(23,i0, 1) = (-1. + 2.*x)*y*(-1. + z*z);
478 outputValues(23,i0, 2) = (-1. + 2.*x)*(-1. + y*y)*z;
479 outputValues(23,i0, 3) =-(-1. + x)*x*(1. - z)*(1. + z);
480 outputValues(23,i0, 4) = 2.*(-1. + x)*x*y*z;
481 outputValues(23,i0, 5) =-(-1. + x)*x*(1. - y)*(1. + y);
482
483 outputValues(24,i0, 0) = (1. - y)*(1. + y)*(1. - z)*(1. + z);
484 outputValues(24,i0, 1) = (1. + 2.*x)*y*(-1. + z*z);
485 outputValues(24,i0, 2) = (1. + 2.*x)*(-1. + y*y)*z;
486 outputValues(24,i0, 3) = x*(1. + x)*(-1. + z)*(1. + z);
487 outputValues(24,i0, 4) = 2.*x*(1. + x)*y*z;
488 outputValues(24,i0, 5) = x*(1. + x)*(-1. + y)*(1. + y);
489
490 outputValues(25,i0, 0) = -(-1. + y)*y*(1. - z)*(1. + z);
491 outputValues(25,i0, 1) = x*(-1. + 2.*y)*(-1. + z*z);
492 outputValues(25,i0, 2) = 2.*x*(-1. + y)*y*z;
493 outputValues(25,i0, 3) = (1. - x)*(1. + x)*(1. - z)*(1. + z);
494 outputValues(25,i0, 4) = (-1. + x*x)*(-1. + 2.*y)*z;
495 outputValues(25,i0, 5) =-(1. - x)*(1. + x)*(-1. + y)*y;
496
497 outputValues(26,i0, 0) = y*(1. + y)*(-1. + z)*(1. + z);
498 outputValues(26,i0, 1) = x*(1. + 2.*y)*(-1. + z*z);
499 outputValues(26,i0, 2) = 2.*x*y*(1. + y)*z;
500 outputValues(26,i0, 3) = (-1. + x)*(1. + x)*(-1. + z)*(1. + z);
501 outputValues(26,i0, 4) = (-1. + x*x)*(1. + 2.*y)*z;
502 outputValues(26,i0, 5) = (-1. + x)*(1. + x)*y*(1. + y);
503 }
504 break;
505
506 case OPERATOR_D3:
507 for (int i0 = 0; i0 < dim0; i0++) {
508 x = inputPoints(i0,0);
509 y = inputPoints(i0,1);
510 z = inputPoints(i0,2);
511
512 outputValues(0,i0, 0) = 0.;
513 outputValues(0,i0, 1) = ((-1.+ 2.*y)*(-1.+ z)*z)/4.;
514 outputValues(0,i0, 2) = ((-1.+ y)*y*(-1.+ 2.*z))/4.;
515 outputValues(0,i0, 3) = ((-1.+ 2.*x)*(-1.+ z)*z)/4.;
516 outputValues(0,i0, 4) = ((-1.+ 2.*x)*(-1.+ 2.*y)*(-1.+ 2.*z))/8.;
517 outputValues(0,i0, 5) = ((-1.+ 2.*x)*(-1.+ y)*y)/4.;
518 outputValues(0,i0, 6) = 0.;
519 outputValues(0,i0, 7) = ((-1.+ x)*x*(-1.+ 2.*z))/4.;
520 outputValues(0,i0, 8) = ((-1.+ x)*x*(-1.+ 2.*y))/4.;
521 outputValues(0,i0, 9) = 0.;
522
523 outputValues(1, i0, 0) = 0.;
524 outputValues(1, i0, 1) = ((-1.+ 2.*y)*(-1.+ z)*z)/4.;
525 outputValues(1, i0, 2) = ((-1.+ y)*y*(-1.+ 2.*z))/4.;
526 outputValues(1, i0, 3) = ((1.+ 2.*x)*(-1.+ z)*z)/4.;
527 outputValues(1, i0, 4) = ((1.+ 2.*x)*(-1.+ 2.*y)*(-1.+ 2.*z))/8.;
528 outputValues(1, i0, 5) = ((1.+ 2.*x)*(-1.+ y)*y)/4.;
529 outputValues(1, i0, 6) = 0.;
530 outputValues(1, i0, 7) = (x*(1.+ x)*(-1.+ 2.*z))/4.;
531 outputValues(1, i0, 8) = (x*(1.+ x)*(-1.+ 2.*y))/4.;
532 outputValues(1, i0, 9) = 0.;
533
534 outputValues(2, i0, 0) = 0.;
535 outputValues(2, i0, 1) = ((1.+ 2.*y)*(-1.+ z)*z)/4.;
536 outputValues(2, i0, 2) = (y*(1.+ y)*(-1.+ 2.*z))/4.;
537 outputValues(2, i0, 3) = ((1.+ 2.*x)*(-1.+ z)*z)/4.;
538 outputValues(2, i0, 4) = ((1.+ 2.*x)*(1.+ 2.*y)*(-1.+ 2.*z))/8.;
539 outputValues(2, i0, 5) = ((1.+ 2.*x)*y*(1.+ y))/4.;
540 outputValues(2, i0, 6) = 0.;
541 outputValues(2, i0, 7) = (x*(1.+ x)*(-1.+ 2.*z))/4.;
542 outputValues(2, i0, 8) = (x*(1.+ x)*(1.+ 2.*y))/4.;
543 outputValues(2, i0, 9) = 0.;
544
545 outputValues(3, i0, 0) = 0.;
546 outputValues(3, i0, 1) = ((1.+ 2.*y)*(-1.+ z)*z)/4.;
547 outputValues(3, i0, 2) = (y*(1.+ y)*(-1.+ 2.*z))/4.;
548 outputValues(3, i0, 3) = ((-1.+ 2.*x)*(-1.+ z)*z)/4.;
549 outputValues(3, i0, 4) = ((-1.+ 2.*x)*(1.+ 2.*y)*(-1.+ 2.*z))/8.;
550 outputValues(3, i0, 5) = ((-1.+ 2.*x)*y*(1.+ y))/4.;
551 outputValues(3, i0, 6) = 0.;
552 outputValues(3, i0, 7) = ((-1.+ x)*x*(-1.+ 2.*z))/4.;
553 outputValues(3, i0, 8) = ((-1.+ x)*x*(1.+ 2.*y))/4.;
554 outputValues(3, i0, 9) = 0.;
555
556 outputValues(4, i0, 0) = 0.;
557 outputValues(4, i0, 1) = ((-1.+ 2.*y)*z*(1.+ z))/4.;
558 outputValues(4, i0, 2) = ((-1.+ y)*y*(1.+ 2.*z))/4.;
559 outputValues(4, i0, 3) = ((-1.+ 2.*x)*z*(1.+ z))/4.;
560 outputValues(4, i0, 4) = ((-1.+ 2.*x)*(-1.+ 2.*y)*(1.+ 2.*z))/8.;
561 outputValues(4, i0, 5) = ((-1.+ 2.*x)*(-1.+ y)*y)/4.;
562 outputValues(4, i0, 6) = 0.;
563 outputValues(4, i0, 7) = ((-1.+ x)*x*(1.+ 2.*z))/4.;
564 outputValues(4, i0, 8) = ((-1.+ x)*x*(-1.+ 2.*y))/4.;
565 outputValues(4, i0, 9) = 0.;
566
567 outputValues(5, i0, 0) = 0.;
568 outputValues(5, i0, 1) = ((-1.+ 2.*y)*z*(1.+ z))/4.;
569 outputValues(5, i0, 2) = ((-1.+ y)*y*(1.+ 2.*z))/4.;
570 outputValues(5, i0, 3) = ((1.+ 2.*x)*z*(1.+ z))/4.;
571 outputValues(5, i0, 4) = ((1.+ 2.*x)*(-1.+ 2.*y)*(1.+ 2.*z))/8.;
572 outputValues(5, i0, 5) = ((1.+ 2.*x)*(-1.+ y)*y)/4.;
573 outputValues(5, i0, 6) = 0.;
574 outputValues(5, i0, 7) = (x*(1.+ x)*(1.+ 2.*z))/4.;
575 outputValues(5, i0, 8) = (x*(1.+ x)*(-1.+ 2.*y))/4.;
576 outputValues(5, i0, 9) = 0.;
577
578 outputValues(6, i0, 0) = 0.;
579 outputValues(6, i0, 1) = ((1.+ 2.*y)*z*(1.+ z))/4.;
580 outputValues(6, i0, 2) = (y*(1.+ y)*(1.+ 2.*z))/4.;
581 outputValues(6, i0, 3) = ((1.+ 2.*x)*z*(1.+ z))/4.;
582 outputValues(6, i0, 4) = ((1.+ 2.*x)*(1.+ 2.*y)*(1.+ 2.*z))/8.;
583 outputValues(6, i0, 5) = ((1.+ 2.*x)*y*(1.+ y))/4.;
584 outputValues(6, i0, 6) = 0.;
585 outputValues(6, i0, 7) = (x*(1.+ x)*(1.+ 2.*z))/4.;
586 outputValues(6, i0, 8) = (x*(1.+ x)*(1.+ 2.*y))/4.;
587 outputValues(6, i0, 9) = 0.;
588
589 outputValues(7, i0, 0) = 0.;
590 outputValues(7, i0, 1) = ((1.+ 2.*y)*z*(1.+ z))/4.;
591 outputValues(7, i0, 2) = (y*(1.+ y)*(1.+ 2.*z))/4.;
592 outputValues(7, i0, 3) = ((-1.+ 2.*x)*z*(1.+ z))/4.;
593 outputValues(7, i0, 4) = ((-1.+ 2.*x)*(1.+ 2.*y)*(1.+ 2.*z))/8.;
594 outputValues(7, i0, 5) = ((-1.+ 2.*x)*y*(1.+ y))/4.;
595 outputValues(7, i0, 6) = 0.;
596 outputValues(7, i0, 7) = ((-1.+ x)*x*(1.+ 2.*z))/4.;
597 outputValues(7, i0, 8) = ((-1.+ x)*x*(1.+ 2.*y))/4.;
598 outputValues(7, i0, 9) = 0.;
599
600 outputValues(8, i0, 0) = 0.;
601 outputValues(8, i0, 1) = -((-1.+ 2.*y)*(-1.+ z)*z)/2.;
602 outputValues(8, i0, 2) = -((-1.+ y)*y*(-1.+ 2.*z))/2.;
603 outputValues(8, i0, 3) = -(x*(-1.+ z)*z);
604 outputValues(8, i0, 4) = -(x*(-1.+ 2.*y)*(-1.+ 2.*z))/2.;
605 outputValues(8, i0, 5) = -(x*(-1.+ y)*y);
606 outputValues(8, i0, 6) = 0.;
607 outputValues(8, i0, 7) = -((-1.+ (x*x))*(-1.+ 2.*z))/2.;
608 outputValues(8, i0, 8) = -((-1.+ (x*x))*(-1.+ 2.*y))/2.;
609 outputValues(8, i0, 9) = 0.;
610
611 outputValues(9, i0, 0) = 0.;
612 outputValues(9, i0, 1) = -(y*(-1.+ z)*z);
613 outputValues(9, i0, 2) = -((-1.+ (y*y))*(-1.+ 2.*z))/2.;
614 outputValues(9, i0, 3) = -((1.+ 2.*x)*(-1.+ z)*z)/2.;
615 outputValues(9, i0, 4) = -((1.+ 2.*x)*y*(-1.+ 2.*z))/2.;
616 outputValues(9, i0, 5) = -((1.+ 2.*x)*(-1.+ (y*y)))/2.;
617 outputValues(9, i0, 6) = 0.;
618 outputValues(9, i0, 7) = -(x*(1.+ x)*(-1.+ 2.*z))/2.;
619 outputValues(9, i0, 8) = -(x*(1.+ x)*y);
620 outputValues(9, i0, 9) = 0.;
621
622 outputValues(10,i0, 0) = 0.;
623 outputValues(10,i0, 1) = -((1.+ 2.*y)*(-1.+ z)*z)/2.;
624 outputValues(10,i0, 2) = -(y*(1.+ y)*(-1.+ 2.*z))/2.;
625 outputValues(10,i0, 3) = -(x*(-1.+ z)*z);
626 outputValues(10,i0, 4) = -(x*(1.+ 2.*y)*(-1.+ 2.*z))/2.;
627 outputValues(10,i0, 5) = -(x*y*(1.+ y));
628 outputValues(10,i0, 6) = 0.;
629 outputValues(10,i0, 7) = -((-1.+ (x*x))*(-1.+ 2.*z))/2.;
630 outputValues(10,i0, 8) = -((-1.+ (x*x))*(1.+ 2.*y))/2.;
631 outputValues(10,i0, 9) = 0.;
632
633 outputValues(11,i0, 0) = 0.;
634 outputValues(11,i0, 1) = -(y*(-1.+ z)*z);
635 outputValues(11,i0, 2) = -((-1.+ (y*y))*(-1.+ 2.*z))/2.;
636 outputValues(11,i0, 3) = -((-1.+ 2.*x)*(-1.+ z)*z)/2.;
637 outputValues(11,i0, 4) = -((-1.+ 2.*x)*y*(-1.+ 2.*z))/2.;
638 outputValues(11,i0, 5) = -((-1.+ 2.*x)*(-1.+ (y*y)))/2.;
639 outputValues(11,i0, 6) = 0.;
640 outputValues(11,i0, 7) = -((-1.+ x)*x*(-1.+ 2.*z))/2.;
641 outputValues(11,i0, 8) = -((-1.+ x)*x*y);
642 outputValues(11,i0, 9) = 0.;
643
644 outputValues(12,i0, 0) = 0.;
645 outputValues(12,i0, 1) = -((-1.+ 2.*y)*(-1.+ (z*z)))/2.;
646 outputValues(12,i0, 2) = -((-1.+ y)*y*z);
647 outputValues(12,i0, 3) = -((-1.+ 2.*x)*(-1.+ (z*z)))/2.;
648 outputValues(12,i0, 4) = -((-1.+ 2.*x)*(-1.+ 2.*y)*z)/2.;
649 outputValues(12,i0, 5) = -((-1.+ 2.*x)*(-1.+ y)*y)/2.;
650 outputValues(12,i0, 6) = 0.;
651 outputValues(12,i0, 7) = -((-1.+ x)*x*z);
652 outputValues(12,i0, 8) = -((-1.+ x)*x*(-1.+ 2.*y))/2.;
653 outputValues(12,i0, 9) = 0.;
654
655 outputValues(13,i0, 0) = 0.;
656 outputValues(13,i0, 1) = -((-1.+ 2.*y)*(-1.+ (z*z)))/2.;
657 outputValues(13,i0, 2) = -((-1.+ y)*y*z);
658 outputValues(13,i0, 3) = -((1.+ 2.*x)*(-1.+ (z*z)))/2.;
659 outputValues(13,i0, 4) = -((1.+ 2.*x)*(-1.+ 2.*y)*z)/2.;
660 outputValues(13,i0, 5) = -((1.+ 2.*x)*(-1.+ y)*y)/2.;
661 outputValues(13,i0, 6) = 0.;
662 outputValues(13,i0, 7) = -(x*(1.+ x)*z);
663 outputValues(13,i0, 8) = -(x*(1.+ x)*(-1.+ 2.*y))/2.;
664 outputValues(13,i0, 9) = 0.;
665
666 outputValues(14,i0, 0) = 0.;
667 outputValues(14,i0, 1) = -((1.+ 2.*y)*(-1.+ (z*z)))/2.;
668 outputValues(14,i0, 2) = -(y*(1.+ y)*z);
669 outputValues(14,i0, 3) = -((1.+ 2.*x)*(-1.+ (z*z)))/2.;
670 outputValues(14,i0, 4) = -((1.+ 2.*x)*(1.+ 2.*y)*z)/2.;
671 outputValues(14,i0, 5) = -((1.+ 2.*x)*y*(1.+ y))/2.;
672 outputValues(14,i0, 6) = 0.;
673 outputValues(14,i0, 7) = -(x*(1.+ x)*z);
674 outputValues(14,i0, 8) = -(x*(1.+ x)*(1.+ 2.*y))/2.;
675 outputValues(14,i0, 9) = 0.;
676
677 outputValues(15,i0, 0) = 0.;
678 outputValues(15,i0, 1) = -((1.+ 2.*y)*(-1.+ (z*z)))/2.;
679 outputValues(15,i0, 2) = -(y*(1.+ y)*z);
680 outputValues(15,i0, 3) = -((-1.+ 2.*x)*(-1.+ (z*z)))/2.;
681 outputValues(15,i0, 4) = -((-1.+ 2.*x)*(1.+ 2.*y)*z)/2.;
682 outputValues(15,i0, 5) = -((-1.+ 2.*x)*y*(1.+ y))/2.;
683 outputValues(15,i0, 6) = 0.;
684 outputValues(15,i0, 7) = -((-1.+ x)*x*z);
685 outputValues(15,i0, 8) = -((-1.+ x)*x*(1.+ 2.*y))/2.;
686 outputValues(15,i0, 9) = 0.;
687
688 outputValues(16,i0, 0) = 0.;
689 outputValues(16,i0, 1) = -((-1.+ 2.*y)*z*(1.+ z))/2.;
690 outputValues(16,i0, 2) = -((-1.+ y)*y*(1.+ 2.*z))/2.;
691 outputValues(16,i0, 3) = -(x*z*(1.+ z));
692 outputValues(16,i0, 4) = -(x*(-1.+ 2.*y)*(1.+ 2.*z))/2.;
693 outputValues(16,i0, 5) = -(x*(-1.+ y)*y);
694 outputValues(16,i0, 6) = 0.;
695 outputValues(16,i0, 7) = -((-1.+ (x*x))*(1.+ 2.*z))/2.;
696 outputValues(16,i0, 8) = -((-1.+ (x*x))*(-1.+ 2.*y))/2.;
697 outputValues(16,i0, 9) = 0.;
698
699 outputValues(17,i0, 0) = 0.;
700 outputValues(17,i0, 1) = -(y*z*(1.+ z));
701 outputValues(17,i0, 2) = -((-1.+ (y*y))*(1.+ 2.*z))/2.;
702 outputValues(17,i0, 3) = -((1.+ 2.*x)*z*(1.+ z))/2.;
703 outputValues(17,i0, 4) = -((1.+ 2.*x)*y*(1.+ 2.*z))/2.;
704 outputValues(17,i0, 5) = -((1.+ 2.*x)*(-1.+ (y*y)))/2.;
705 outputValues(17,i0, 6) = 0.;
706 outputValues(17,i0, 7) = -(x*(1.+ x)*(1.+ 2.*z))/2.;
707 outputValues(17,i0, 8) = -(x*(1.+ x)*y);
708 outputValues(17,i0, 9) = 0.;
709
710 outputValues(18,i0, 0) = 0.;
711 outputValues(18,i0, 1) = -((1.+ 2.*y)*z*(1.+ z))/2.;
712 outputValues(18,i0, 2) = -(y*(1.+ y)*(1.+ 2.*z))/2.;
713 outputValues(18,i0, 3) = -(x*z*(1.+ z));
714 outputValues(18,i0, 4) = -(x*(1.+ 2.*y)*(1.+ 2.*z))/2.;
715 outputValues(18,i0, 5) = -(x*y*(1.+ y));
716 outputValues(18,i0, 6) = 0.;
717 outputValues(18,i0, 7) = -((-1.+ (x*x))*(1.+ 2.*z))/2.;
718 outputValues(18,i0, 8) = -((-1.+ (x*x))*(1.+ 2.*y))/2.;
719 outputValues(18,i0, 9) = 0.;
720
721 outputValues(19,i0, 0) = 0.;
722 outputValues(19,i0, 1) = -(y*z*(1.+ z));
723 outputValues(19,i0, 2) = -((-1.+ (y*y))*(1.+ 2.*z))/2.;
724 outputValues(19,i0, 3) = -((-1.+ 2.*x)*z*(1.+ z))/2.;
725 outputValues(19,i0, 4) = -((-1.+ 2.*x)*y*(1.+ 2.*z))/2.;
726 outputValues(19,i0, 5) = -((-1.+ 2.*x)*(-1.+ (y*y)))/2.;
727 outputValues(19,i0, 6) = 0.;
728 outputValues(19,i0, 7) = -((-1.+ x)*x*(1.+ 2.*z))/2.;
729 outputValues(19,i0, 8) = -((-1.+ x)*x*y);
730 outputValues(19,i0, 9) = 0.;
731
732 outputValues(20,i0, 0) = 0.;
733 outputValues(20,i0, 1) = -4*y*(-1.+ (z*z));
734 outputValues(20,i0, 2) = -4*(-1.+ (y*y))*z;
735 outputValues(20,i0, 3) = -4*x*(-1.+ (z*z));
736 outputValues(20,i0, 4) = -8*x*y*z;
737 outputValues(20,i0, 5) = -4*x*(-1.+ (y*y));
738 outputValues(20,i0, 6) = 0.;
739 outputValues(20,i0, 7) = -4*(-1.+ (x*x))*z;
740 outputValues(20,i0, 8) = -4*(-1.+ (x*x))*y;
741 outputValues(20,i0, 9) = 0.;
742
743 outputValues(21,i0, 0) = 0.;
744 outputValues(21,i0, 1) = 2.*y*(-1.+ z)*z;
745 outputValues(21,i0, 2) = (-1.+ (y*y))*(-1.+ 2.*z);
746 outputValues(21,i0, 3) = 2.*x*(-1.+ z)*z;
747 outputValues(21,i0, 4) = 2.*x*y*(-1.+ 2.*z);
748 outputValues(21,i0, 5) = 2.*x*(-1.+ (y*y));
749 outputValues(21,i0, 6) = 0.;
750 outputValues(21,i0, 7) = (-1.+ (x*x))*(-1.+ 2.*z);
751 outputValues(21,i0, 8) = 2.*(-1.+ (x*x))*y;
752 outputValues(21,i0, 9) = 0.;
753
754 outputValues(22,i0, 0) = 0.;
755 outputValues(22,i0, 1) = 2.*y*z*(1.+ z);
756 outputValues(22,i0, 2) = (-1.+ (y*y))*(1.+ 2.*z);
757 outputValues(22,i0, 3) = 2.*x*z*(1.+ z);
758 outputValues(22,i0, 4) = 2.*x*y*(1.+ 2.*z);
759 outputValues(22,i0, 5) = 2.*x*(-1.+ (y*y));
760 outputValues(22,i0, 6) = 0.;
761 outputValues(22,i0, 7) = (-1.+ (x*x))*(1.+ 2.*z);
762 outputValues(22,i0, 8) = 2.*(-1.+ (x*x))*y;
763 outputValues(22,i0, 9) = 0.;
764
765 outputValues(23,i0, 0) = 0.;
766 outputValues(23,i0, 1) = 2.*y*(-1.+ (z*z));
767 outputValues(23,i0, 2) = 2.*(-1.+ (y*y))*z;
768 outputValues(23,i0, 3) = (-1.+ 2.*x)*(-1.+ (z*z));
769 outputValues(23,i0, 4) = 2.*(-1.+ 2.*x)*y*z;
770 outputValues(23,i0, 5) = (-1.+ 2.*x)*(-1.+ (y*y));
771 outputValues(23,i0, 6) = 0.;
772 outputValues(23,i0, 7) = 2.*(-1.+ x)*x*z;
773 outputValues(23,i0, 8) = 2.*(-1.+ x)*x*y;
774 outputValues(23,i0, 9) = 0.;
775
776 outputValues(24,i0, 0) = 0.;
777 outputValues(24,i0, 1) = 2.*y*(-1.+ (z*z));
778 outputValues(24,i0, 2) = 2.*(-1.+ (y*y))*z;
779 outputValues(24,i0, 3) = (1.+ 2.*x)*(-1.+ (z*z));
780 outputValues(24,i0, 4) = 2.*(1.+ 2.*x)*y*z;
781 outputValues(24,i0, 5) = (1.+ 2.*x)*(-1.+ (y*y));
782 outputValues(24,i0, 6) = 0.;
783 outputValues(24,i0, 7) = 2.*x*(1.+ x)*z;
784 outputValues(24,i0, 8) = 2.*x*(1.+ x)*y;
785 outputValues(24,i0, 9) = 0.;
786
787 outputValues(25,i0, 0) = 0.;
788 outputValues(25,i0, 1) = (-1.+ 2.*y)*(-1.+ (z*z));
789 outputValues(25,i0, 2) = 2.*(-1.+ y)*y*z;
790 outputValues(25,i0, 3) = 2.*x*(-1.+ (z*z));
791 outputValues(25,i0, 4) = 2.*x*(-1.+ 2.*y)*z;
792 outputValues(25,i0, 5) = 2.*x*(-1.+ y)*y;
793 outputValues(25,i0, 6) = 0.;
794 outputValues(25,i0, 7) = 2.*(-1.+ (x*x))*z;
795 outputValues(25,i0, 8) = (-1.+ (x*x))*(-1.+ 2.*y);
796 outputValues(25,i0, 9) = 0.;
797
798 outputValues(26,i0, 0) = 0.;
799 outputValues(26,i0, 1) = (1.+ 2.*y)*(-1.+ (z*z));
800 outputValues(26,i0, 2) = 2.*y*(1.+ y)*z;
801 outputValues(26,i0, 3) = 2.*x*(-1.+ (z*z));
802 outputValues(26,i0, 4) = 2.*x*(1.+ 2.*y)*z;
803 outputValues(26,i0, 5) = 2.*x*y*(1.+ y);
804 outputValues(26,i0, 6) = 0.;
805 outputValues(26,i0, 7) = 2.*(-1.+ (x*x))*z;
806 outputValues(26,i0, 8) = (-1.+ (x*x))*(1.+ 2.*y);
807 outputValues(26,i0, 9) = 0.;
808 }
809 break;
810
811 case OPERATOR_D4:
812 {
813 // Non-zero entries have Dk (derivative cardinality) indices {3,4,5,7,8,12}, all other entries are 0.
814 // Intitialize array by zero and then fill only non-zero entries.
815 int DkCardinality = Intrepid::getDkCardinality(operatorType, this -> basisCellTopology_.getDimension() );
816 for(int dofOrd = 0; dofOrd < this -> basisCardinality_; dofOrd++) {
817 for (int i0 = 0; i0 < dim0; i0++) {
818 for(int dkOrd = 0; dkOrd < DkCardinality; dkOrd++){
819 outputValues(dofOrd, i0, dkOrd) = 0.0;
820 }
821 }
822 }
823
824 for (int i0 = 0; i0 < dim0; i0++) {
825 x = inputPoints(i0,0);
826 y = inputPoints(i0,1);
827 z = inputPoints(i0,2);
828
829 outputValues(0, i0, 3) = ((-1.+ z)*z)/2.;
830 outputValues(0, i0, 4) = ((-1.+ 2.*y)*(-1.+ 2.*z))/4.;
831 outputValues(0, i0, 5) = ((-1.+ y)*y)/2.;
832 outputValues(0, i0, 7) = ((-1.+ 2.*x)*(-1.+ 2.*z))/4.;
833 outputValues(0, i0, 8) = ((-1.+ 2.*x)*(-1.+ 2.*y))/4.;
834 outputValues(0, i0, 12)= ((-1.+ x)*x)/2.;
835
836 outputValues(1, i0, 3) = ((-1.+ z)*z)/2.;
837 outputValues(1, i0, 4) = ((-1.+ 2.*y)*(-1.+ 2.*z))/4.;
838 outputValues(1, i0, 5) = ((-1.+ y)*y)/2.;
839 outputValues(1, i0, 7) = ((1. + 2.*x)*(-1.+ 2.*z))/4.;
840 outputValues(1, i0, 8) = ((1. + 2.*x)*(-1.+ 2.*y))/4.;
841 outputValues(1, i0, 12)= (x*(1. + x))/2.;
842
843 outputValues(2, i0, 3) = ((-1.+ z)*z)/2.;
844 outputValues(2, i0, 4) = ((1. + 2.*y)*(-1.+ 2.*z))/4.;
845 outputValues(2, i0, 5) = (y*(1. + y))/2.;
846 outputValues(2, i0, 7) = ((1. + 2.*x)*(-1.+ 2.*z))/4.;
847 outputValues(2, i0, 8) = ((1. + 2.*x)*(1. + 2.*y))/4.;
848 outputValues(2, i0, 12)= (x*(1. + x))/2.;
849
850 outputValues(3, i0, 3) = ((-1.+ z)*z)/2.;
851 outputValues(3, i0, 4) = ((1. + 2.*y)*(-1.+ 2.*z))/4.;
852 outputValues(3, i0, 5) = (y*(1. + y))/2.;
853 outputValues(3, i0, 7) = ((-1.+ 2.*x)*(-1.+ 2.*z))/4.;
854 outputValues(3, i0, 8) = ((-1.+ 2.*x)*(1. + 2.*y))/4.;
855 outputValues(3, i0, 12)= ((-1.+ x)*x)/2.;
856
857 outputValues(4, i0, 3) = (z*(1. + z))/2.;
858 outputValues(4, i0, 4) = ((-1.+ 2.*y)*(1. + 2.*z))/4.;
859 outputValues(4, i0, 5) = ((-1.+ y)*y)/2.;
860 outputValues(4, i0, 7) = ((-1.+ 2.*x)*(1. + 2.*z))/4.;
861 outputValues(4, i0, 8) = ((-1.+ 2.*x)*(-1.+ 2.*y))/4.;
862 outputValues(4, i0, 12)= ((-1.+ x)*x)/2.;
863
864 outputValues(5, i0, 3) = (z*(1. + z))/2.;
865 outputValues(5, i0, 4) = ((-1.+ 2.*y)*(1. + 2.*z))/4.;
866 outputValues(5, i0, 5) = ((-1.+ y)*y)/2.;
867 outputValues(5, i0, 7) = ((1. + 2.*x)*(1. + 2.*z))/4.;
868 outputValues(5, i0, 8) = ((1. + 2.*x)*(-1.+ 2.*y))/4.;
869 outputValues(5, i0, 12)= (x*(1. + x))/2.;
870
871 outputValues(6, i0, 3) = (z*(1. + z))/2.;
872 outputValues(6, i0, 4) = ((1. + 2.*y)*(1. + 2.*z))/4.;
873 outputValues(6, i0, 5) = (y*(1. + y))/2.;
874 outputValues(6, i0, 7) = ((1. + 2.*x)*(1. + 2.*z))/4.;
875 outputValues(6, i0, 8) = ((1. + 2.*x)*(1. + 2.*y))/4.;
876 outputValues(6, i0, 12)= (x*(1. + x))/2.;
877
878 outputValues(7, i0, 3) = (z*(1. + z))/2.;
879 outputValues(7, i0, 4) = ((1. + 2.*y)*(1. + 2.*z))/4.;
880 outputValues(7, i0, 5) = (y*(1. + y))/2.;
881 outputValues(7, i0, 7) = ((-1.+ 2.*x)*(1. + 2.*z))/4.;
882 outputValues(7, i0, 8) = ((-1.+ 2.*x)*(1. + 2.*y))/4.;
883 outputValues(7, i0, 12)= ((-1.+ x)*x)/2.;
884
885 outputValues(8, i0, 3) = -((-1.+ z)*z);
886 outputValues(8, i0, 4) = -0.5 + y + z - 2.*y*z;
887 outputValues(8, i0, 5) = -((-1.+ y)*y);
888 outputValues(8, i0, 7) = x - 2.*x*z;
889 outputValues(8, i0, 8) = x - 2.*x*y;
890 outputValues(8, i0, 12)= 1. - x*x;
891
892 outputValues(9, i0, 3) = -((-1.+ z)*z);
893 outputValues(9, i0, 4) = y - 2.*y*z;
894 outputValues(9, i0, 5) = 1 - y*y;
895 outputValues(9, i0, 7) = 0.5 + x - z - 2.*x*z;
896 outputValues(9, i0, 8) = -((1. + 2.*x)*y);
897 outputValues(9, i0, 12)= -(x*(1. + x));
898
899 outputValues(10,i0, 3) = -((-1.+ z)*z);
900 outputValues(10,i0, 4) = 0.5 + y - z - 2.*y*z;
901 outputValues(10,i0, 5) = -(y*(1. + y));
902 outputValues(10,i0, 7) = x - 2.*x*z;
903 outputValues(10,i0, 8) = -(x*(1. + 2.*y));
904 outputValues(10,i0, 12)= 1. - x*x;
905
906 outputValues(11,i0, 3) = -((-1.+ z)*z);
907 outputValues(11,i0, 4) = y - 2.*y*z;
908 outputValues(11,i0, 5) = 1. - y*y;
909 outputValues(11,i0, 7) = -0.5 + x + z - 2.*x*z;
910 outputValues(11,i0, 8) = y - 2.*x*y;
911 outputValues(11,i0, 12)= -((-1.+ x)*x);
912
913 outputValues(12,i0, 3) = 1. - z*z;
914 outputValues(12,i0, 4) = z - 2.*y*z;
915 outputValues(12,i0, 5) = -((-1.+ y)*y);
916 outputValues(12,i0, 7) = z - 2.*x*z;
917 outputValues(12,i0, 8) = -0.5 + x + y - 2.*x*y;
918 outputValues(12,i0, 12)= -((-1.+ x)*x);
919
920 outputValues(13,i0, 3) = 1. - z*z;
921 outputValues(13,i0, 4) = z - 2.*y*z;
922 outputValues(13,i0, 5) = -((-1.+ y)*y);
923 outputValues(13,i0, 7) = -((1. + 2.*x)*z);
924 outputValues(13,i0, 8) = 0.5 + x - y - 2.*x*y;
925 outputValues(13,i0, 12)= -(x*(1. + x));
926
927 outputValues(14,i0, 3) = 1. - z*z;
928 outputValues(14,i0, 4) = -((1. + 2.*y)*z);
929 outputValues(14,i0, 5) = -(y*(1. + y));
930 outputValues(14,i0, 7) = -((1. + 2.*x)*z);
931 outputValues(14,i0, 8) = -((1. + 2.*x)*(1. + 2.*y))/2.;
932 outputValues(14,i0, 12)= -(x*(1. + x));
933
934 outputValues(15,i0, 3) = 1. - z*z;
935 outputValues(15,i0, 4) = -((1. + 2.*y)*z);
936 outputValues(15,i0, 5) = -(y*(1. + y));
937 outputValues(15,i0, 7) = z - 2.*x*z;
938 outputValues(15,i0, 8) = 0.5 + y - x*(1. + 2.*y);
939 outputValues(15,i0, 12)= -((-1.+ x)*x);
940
941 outputValues(16,i0, 3) = -(z*(1. + z));
942 outputValues(16,i0, 4) = 0.5 + z - y*(1. + 2.*z);
943 outputValues(16,i0, 5) = -((-1.+ y)*y);
944 outputValues(16,i0, 7) = -(x*(1. + 2.*z));
945 outputValues(16,i0, 8) = x - 2.*x*y;
946 outputValues(16,i0, 12)= 1. - x*x;
947
948 outputValues(17,i0, 3) = -(z*(1. + z));
949 outputValues(17,i0, 4) = -(y*(1. + 2.*z));
950 outputValues(17,i0, 5) = 1. - y*y;
951 outputValues(17,i0, 7) = -((1. + 2.*x)*(1. + 2.*z))/2.;
952 outputValues(17,i0, 8) = -((1. + 2.*x)*y);
953 outputValues(17,i0, 12)= -(x*(1. + x));
954
955 outputValues(18,i0, 3) = -(z*(1. + z));
956 outputValues(18,i0, 4) = -((1. + 2.*y)*(1. + 2.*z))/2.;
957 outputValues(18,i0, 5) = -(y*(1. + y));
958 outputValues(18,i0, 7) = -(x*(1. + 2.*z));
959 outputValues(18,i0, 8) = -(x*(1. + 2.*y));
960 outputValues(18,i0, 12)= 1. - x*x;
961
962 outputValues(19,i0, 3) = -(z*(1. + z));
963 outputValues(19,i0, 4) = -(y*(1. + 2.*z));
964 outputValues(19,i0, 5) = 1. - y*y;
965 outputValues(19,i0, 7) = 0.5 + z - x*(1. + 2.*z);
966 outputValues(19,i0, 8) = y - 2.*x*y;
967 outputValues(19,i0, 12)= -((-1.+ x)*x);
968
969 outputValues(20,i0, 3) = 4. - 4.*z*z;
970 outputValues(20,i0, 4) = -8.*y*z;
971 outputValues(20,i0, 5) = 4. - 4.*y*y;
972 outputValues(20,i0, 7) = -8.*x*z;
973 outputValues(20,i0, 8) = -8.*x*y;
974 outputValues(20,i0, 12)= 4. - 4.*x*x;
975
976 outputValues(21,i0, 3) = 2.*(-1.+ z)*z;
977 outputValues(21,i0, 4) = 2.*y*(-1.+ 2.*z);
978 outputValues(21,i0, 5) = 2.*(-1.+ y*y);
979 outputValues(21,i0, 7) = 2.*x*(-1.+ 2.*z);
980 outputValues(21,i0, 8) = 4.*x*y;
981 outputValues(21,i0, 12)= 2.*(-1.+ x*x);
982
983 outputValues(22,i0, 3) = 2.*z*(1. + z);
984 outputValues(22,i0, 4) = 2.*(y + 2.*y*z);
985 outputValues(22,i0, 5) = 2.*(-1.+ y*y);
986 outputValues(22,i0, 7) = 2.*(x + 2.*x*z);
987 outputValues(22,i0, 8) = 4.*x*y;
988 outputValues(22,i0, 12)= 2.*(-1.+ x*x);
989
990 outputValues(23,i0, 3) = 2.*(-1.+ z*z);
991 outputValues(23,i0, 4) = 4.*y*z;
992 outputValues(23,i0, 5) = 2.*(-1.+ y*y);
993 outputValues(23,i0, 7) = 2.*(-1.+ 2.*x)*z;
994 outputValues(23,i0, 8) = 2.*(-1.+ 2.*x)*y;
995 outputValues(23,i0, 12)= 2.*(-1.+ x)*x;
996
997 outputValues(24,i0, 3) = 2.*(-1.+ z*z);
998 outputValues(24,i0, 4) = 4.*y*z;
999 outputValues(24,i0, 5) = 2.*(-1.+ y*y);
1000 outputValues(24,i0, 7) = 2.*(z + 2.*x*z);
1001 outputValues(24,i0, 8) = 2.*(y + 2.*x*y);
1002 outputValues(24,i0, 12)= 2.*x*(1. + x);
1003
1004 outputValues(25,i0, 3) = 2.*(-1.+ z*z);
1005 outputValues(25,i0, 4) = 2.*(-1.+ 2.*y)*z;
1006 outputValues(25,i0, 5) = 2.*(-1.+ y)*y;
1007 outputValues(25,i0, 7) = 4.*x*z;
1008 outputValues(25,i0, 8) = 2.*x*(-1.+ 2.*y);
1009 outputValues(25,i0, 12)= 2.*(-1.+ x*x);
1010
1011 outputValues(26,i0, 3) = 2.*(-1.+ z*z);
1012 outputValues(26,i0, 4) = 2.*(z + 2.*y*z);
1013 outputValues(26,i0, 5) = 2.*y*(1. + y);
1014 outputValues(26,i0, 7) = 4.*x*z;
1015 outputValues(26,i0, 8) = 2.*(x + 2.*x*y);
1016 outputValues(26,i0, 12)= 2.*(-1.+ x*x);
1017 }
1018 }
1019 break;
1020
1021 case OPERATOR_D5:
1022 case OPERATOR_D6:
1023 TEUCHOS_TEST_FOR_EXCEPTION( true, std::invalid_argument,
1024 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): operator not supported");
1025 break;
1026
1027 case OPERATOR_D7:
1028 case OPERATOR_D8:
1029 case OPERATOR_D9:
1030 case OPERATOR_D10:
1031 {
1032 // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, DkCardinality)
1033 int DkCardinality = Intrepid::getDkCardinality(operatorType,
1034 this -> basisCellTopology_.getDimension() );
1035 for(int dofOrd = 0; dofOrd < this -> basisCardinality_; dofOrd++) {
1036 for (int i0 = 0; i0 < dim0; i0++) {
1037 for(int dkOrd = 0; dkOrd < DkCardinality; dkOrd++){
1038 outputValues(dofOrd, i0, dkOrd) = 0.0;
1039 }
1040 }
1041 }
1042 }
1043 break;
1044
1045 default:
1046 TEUCHOS_TEST_FOR_EXCEPTION( !( Intrepid::isValidOperator(operatorType) ), std::invalid_argument,
1047 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): Invalid operator type");
1048 }
1049}
1050
1051
1052
1053template<class Scalar, class ArrayScalar>
1055 const ArrayScalar & inputPoints,
1056 const ArrayScalar & cellVertices,
1057 const EOperator operatorType) const {
1058 TEUCHOS_TEST_FOR_EXCEPTION( (true), std::logic_error,
1059 ">>> ERROR (Basis_HGRAD_HEX_C2_FEM): FEM Basis calling an FVD member function");
1060 }
1061
1062template<class Scalar, class ArrayScalar>
1064#ifdef HAVE_INTREPID_DEBUG
1065 // Verify rank of output array.
1066 TEUCHOS_TEST_FOR_EXCEPTION( !(DofCoords.rank() == 2), std::invalid_argument,
1067 ">>> ERROR: (Intrepid::Basis_HGRAD_HEX_C2_FEM::getDofCoords) rank = 2 required for DofCoords array");
1068 // Verify 0th dimension of output array.
1069 TEUCHOS_TEST_FOR_EXCEPTION( !( DofCoords.dimension(0) == this -> basisCardinality_ ), std::invalid_argument,
1070 ">>> ERROR: (Intrepid::Basis_HGRAD_HEX_C2_FEM::getDofCoords) mismatch in number of DoF and 0th dimension of DofCoords array");
1071 // Verify 1st dimension of output array.
1072 TEUCHOS_TEST_FOR_EXCEPTION( !( DofCoords.dimension(1) == (int)(this -> basisCellTopology_.getDimension()) ), std::invalid_argument,
1073 ">>> ERROR: (Intrepid::Basis_HGRAD_HEX_C2_FEM::getDofCoords) incorrect reference cell (1st) dimension in DofCoords array");
1074#endif
1075
1076 DofCoords(0,0) = -1.0; DofCoords(0,1) = -1.0; DofCoords(0,2) = -1.0;
1077 DofCoords(1,0) = 1.0; DofCoords(1,1) = -1.0; DofCoords(1,2) = -1.0;
1078 DofCoords(2,0) = 1.0; DofCoords(2,1) = 1.0; DofCoords(2,2) = -1.0;
1079 DofCoords(3,0) = -1.0; DofCoords(3,1) = 1.0; DofCoords(3,2) = -1.0;
1080 DofCoords(4,0) = -1.0; DofCoords(4,1) = -1.0; DofCoords(4,2) = 1.0;
1081 DofCoords(5,0) = 1.0; DofCoords(5,1) = -1.0; DofCoords(5,2) = 1.0;
1082 DofCoords(6,0) = 1.0; DofCoords(6,1) = 1.0; DofCoords(6,2) = 1.0;
1083 DofCoords(7,0) = -1.0; DofCoords(7,1) = 1.0; DofCoords(7,2) = 1.0;
1084
1085 DofCoords(8,0) = 0.0; DofCoords(8,1) = -1.0; DofCoords(8,2) = -1.0;
1086 DofCoords(9,0) = 1.0; DofCoords(9,1) = 0.0; DofCoords(9,2) = -1.0;
1087 DofCoords(10,0) = 0.0; DofCoords(10,1) = 1.0; DofCoords(10,2) = -1.0;
1088 DofCoords(11,0) = -1.0; DofCoords(11,1) = 0.0; DofCoords(11,2) = -1.0;
1089 DofCoords(12,0) = -1.0; DofCoords(12,1) = -1.0; DofCoords(12,2) = 0.0;
1090 DofCoords(13,0) = 1.0; DofCoords(13,1) = -1.0; DofCoords(13,2) = 0.0;
1091 DofCoords(14,0) = 1.0; DofCoords(14,1) = 1.0; DofCoords(14,2) = 0.0;
1092 DofCoords(15,0) = -1.0; DofCoords(15,1) = 1.0; DofCoords(15,2) = 0.0;
1093 DofCoords(16,0) = 0.0; DofCoords(16,1) = -1.0; DofCoords(16,2) = 1.0;
1094 DofCoords(17,0) = 1.0; DofCoords(17,1) = 0.0; DofCoords(17,2) = 1.0;
1095 DofCoords(18,0) = 0.0; DofCoords(18,1) = 1.0; DofCoords(18,2) = 1.0;
1096 DofCoords(19,0) = -1.0; DofCoords(19,1) = 0.0; DofCoords(19,2) = 1.0;
1097
1098 DofCoords(20,0) = 0.0; DofCoords(20,1) = 0.0; DofCoords(20,2) = 0.0;
1099
1100 DofCoords(21,0) = 0.0; DofCoords(21,1) = 0.0; DofCoords(21,2) = -1.0;
1101 DofCoords(22,0) = 0.0; DofCoords(22,1) = 0.0; DofCoords(22,2) = 1.0;
1102 DofCoords(23,0) = -1.0; DofCoords(23,1) = 0.0; DofCoords(23,2) = 0.0;
1103 DofCoords(24,0) = 1.0; DofCoords(24,1) = 0.0; DofCoords(24,2) = 0.0;
1104 DofCoords(25,0) = 0.0; DofCoords(25,1) = -1.0; DofCoords(25,2) = 0.0;
1105 DofCoords(26,0) = 0.0; DofCoords(26,1) = 1.0; DofCoords(26,2) = 0.0;
1106}
1107
1108}// namespace Intrepid
1109#endif
int isValidOperator(const EOperator operatorType)
Verifies validity of an operator enum.
void setOrdinalTagData(std::vector< std::vector< std::vector< int > > > &tagToOrdinal, std::vector< std::vector< int > > &ordinalToTag, const int *tags, const int basisCard, const int tagSize, const int posScDim, const int posScOrd, const int posDfOrd)
Fills ordinalToTag_ and tagToOrdinal_ by basis-specific tag data.
int getDkCardinality(const EOperator operatorType, const int spaceDim)
Returns cardinality of Dk, i.e., the number of all derivatives of order k.
void getDofCoords(ArrayScalar &DofCoords) const
Returns spatial locations (coordinates) of degrees of freedom on a reference Quadrilateral.
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
Evaluation of a FEM basis on a reference Hexahedron cell.
void initializeTags()
Initializes tagToOrdinal_ and ordinalToTag_ lookup arrays.