Intrepid
test_01.cpp
Go to the documentation of this file.
1// @HEADER
2// ************************************************************************
3//
4// Intrepid 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 Pavel Bochev (pbboche@sandia.gov)
38// Denis Ridzal (dridzal@sandia.gov), or
39// Kara Peterson (kjpeter@sandia.gov)
40//
41// ************************************************************************
42// @HEADER
43
50#include "Teuchos_oblackholestream.hpp"
51#include "Teuchos_RCP.hpp"
52#include "Teuchos_GlobalMPISession.hpp"
53
54using namespace std;
55using namespace Intrepid;
56
57#define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \
58{ \
59 ++nException; \
60 try { \
61 S ; \
62 } \
63 catch (const std::logic_error & err) { \
64 ++throwCounter; \
65 *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
66 *outStream << err.what() << '\n'; \
67 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \
68 }; \
69}
70
71int main(int argc, char *argv[]) {
72
73 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
74
75 // This little trick lets us print to std::cout only if
76 // a (dummy) command-line argument is provided.
77 int iprint = argc - 1;
78 Teuchos::RCP<std::ostream> outStream;
79 Teuchos::oblackholestream bhs; // outputs nothing
80 if (iprint > 0)
81 outStream = Teuchos::rcp(&std::cout, false);
82 else
83 outStream = Teuchos::rcp(&bhs, false);
84
85 // Save the format state of the original std::cout.
86 Teuchos::oblackholestream oldFormatState;
87 oldFormatState.copyfmt(std::cout);
88
89 *outStream \
90 << "===============================================================================\n" \
91 << "| |\n" \
92 << "| Unit Test (Basis_HCURL_TET_I1_FEM) |\n" \
93 << "| |\n" \
94 << "| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
95 << "| 2) Basis values for VALUE and CURL operators |\n" \
96 << "| |\n" \
97 << "| Questions? Contact Pavel Bochev (pbboche@sandia.gov) or |\n" \
98 << "| Denis Ridzal (dridzal@sandia.gov). |\n" \
99 << "| |\n" \
100 << "| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
101 << "| Trilinos website: http://trilinos.sandia.gov |\n" \
102 << "| |\n" \
103 << "===============================================================================\n"\
104 << "| TEST 1: Basis creation, exception testing |\n"\
105 << "===============================================================================\n";
106
107 // Define basis and error flag
109 int errorFlag = 0;
110
111// Define throw number for exception testing
112 int nException = 0;
113 int throwCounter = 0;
114
115 // Define array containing the 4 vertices of the reference TET and its 6 edge centers.
116 FieldContainer<double> tetNodes(10, 3);
117 tetNodes(0,0) = 0.0; tetNodes(0,1) = 0.0; tetNodes(0,2) = 0.0;
118 tetNodes(1,0) = 1.0; tetNodes(1,1) = 0.0; tetNodes(1,2) = 0.0;
119 tetNodes(2,0) = 0.0; tetNodes(2,1) = 1.0; tetNodes(2,2) = 0.0;
120 tetNodes(3,0) = 0.0; tetNodes(3,1) = 0.0; tetNodes(3,2) = 1.0;
121 tetNodes(4,0) = 0.5; tetNodes(4,1) = 0.0; tetNodes(4,2) = 0.0;
122 tetNodes(5,0) = 0.5; tetNodes(5,1) = 0.5; tetNodes(5,2) = 0.0;
123 tetNodes(6,0) = 0.0; tetNodes(6,1) = 0.5; tetNodes(6,2) = 0.0;
124 tetNodes(7,0) = 0.0; tetNodes(7,1) = 0.0; tetNodes(7,2) = 0.5;
125 tetNodes(8,0) = 0.5; tetNodes(8,1) = 0.0; tetNodes(8,2) = 0.5;
126 tetNodes(9,0) = 0.0; tetNodes(9,1) = 0.5; tetNodes(9,2) = 0.5;
127
128
129 // Generic array for the output values; needs to be properly resized depending on the operator type
131
132
133 try{
134 // exception #1: GRAD cannot be applied to HCURL functions
135 // resize vals to rank-3 container with dimensions (num. basis functions, num. points, arbitrary)
136 vals.resize(tetBasis.getCardinality(), tetNodes.dimension(0), 3 );
137 INTREPID_TEST_COMMAND( tetBasis.getValues(vals, tetNodes, OPERATOR_GRAD), throwCounter, nException );
138
139 // exception #2: DIV cannot be applied to HCURL functions
140 // resize vals to rank-2 container with dimensions (num. basis functions, num. points)
141 vals.resize(tetBasis.getCardinality(), tetNodes.dimension(0) );
142 INTREPID_TEST_COMMAND( tetBasis.getValues(vals, tetNodes, OPERATOR_DIV), throwCounter, nException );
143
144 // Exceptions 3-7: all bf tags/bf Ids below are wrong and should cause getDofOrdinal() and
145 // getDofTag() to access invalid array elements thereby causing bounds check exception
146 // exception #3
147 INTREPID_TEST_COMMAND( tetBasis.getDofOrdinal(3,0,0), throwCounter, nException );
148 // exception #4
149 INTREPID_TEST_COMMAND( tetBasis.getDofOrdinal(1,1,1), throwCounter, nException );
150 // exception #5
151 INTREPID_TEST_COMMAND( tetBasis.getDofOrdinal(0,4,1), throwCounter, nException );
152 // exception #6
153 INTREPID_TEST_COMMAND( tetBasis.getDofTag(7), throwCounter, nException );
154 // exception #7
155 INTREPID_TEST_COMMAND( tetBasis.getDofTag(-1), throwCounter, nException );
156
157#ifdef HAVE_INTREPID_DEBUG
158 // Exceptions 8-15 test exception handling with incorrectly dimensioned input/output arrays
159 // exception #8: input points array must be of rank-2
160 FieldContainer<double> badPoints1(4, 5, 3);
161 INTREPID_TEST_COMMAND( tetBasis.getValues(vals,badPoints1,OPERATOR_VALUE), throwCounter, nException );
162
163 // exception #9 dimension 1 in the input point array must equal space dimension of the cell
164 FieldContainer<double> badPoints2(4, 2);
165 INTREPID_TEST_COMMAND( tetBasis.getValues(vals,badPoints2,OPERATOR_VALUE), throwCounter, nException );
166
167 // exception #10 output values must be of rank-3 for OPERATOR_VALUE
168 FieldContainer<double> badVals1(4, 3);
169 INTREPID_TEST_COMMAND( tetBasis.getValues(badVals1,tetNodes,OPERATOR_VALUE), throwCounter, nException );
170
171 // exception #11 output values must be of rank-3 for OPERATOR_CURL
172 INTREPID_TEST_COMMAND( tetBasis.getValues(badVals1,tetNodes,OPERATOR_CURL), throwCounter, nException );
173
174 // exception #12 incorrect 0th dimension of output array (must equal number of basis functions)
175 FieldContainer<double> badVals2(tetBasis.getCardinality() + 1, tetNodes.dimension(0), 3);
176 INTREPID_TEST_COMMAND( tetBasis.getValues(badVals2,tetNodes,OPERATOR_VALUE), throwCounter, nException );
177
178 // exception #13 incorrect 1st dimension of output array (must equal number of points)
179 FieldContainer<double> badVals3(tetBasis.getCardinality(), tetNodes.dimension(0) + 1, 3);
180 INTREPID_TEST_COMMAND( tetBasis.getValues(badVals3,tetNodes,OPERATOR_VALUE), throwCounter, nException );
181
182 // exception #14: incorrect 2nd dimension of output array (must equal the space dimension)
183 FieldContainer<double> badVals4(tetBasis.getCardinality(), tetNodes.dimension(0), 4);
184 INTREPID_TEST_COMMAND( tetBasis.getValues(badVals4,tetNodes,OPERATOR_VALUE), throwCounter, nException );
185
186 // exception #15: incorrect 2nd dimension of output array (must equal the space dimension)
187 INTREPID_TEST_COMMAND( tetBasis.getValues(badVals4,tetNodes,OPERATOR_CURL), throwCounter, nException );
188#endif
189
190 }
191 catch (const std::logic_error & err) {
192 *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
193 *outStream << err.what() << '\n';
194 *outStream << "-------------------------------------------------------------------------------" << "\n\n";
195 errorFlag = -1000;
196 };
197
198 // Check if number of thrown exceptions matches the one we expect
199 // Note Teuchos throw number will not pick up exceptions 3-7 and therefore will not match.
200 if (throwCounter != nException) {
201 errorFlag++;
202 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
203 }
204
205 *outStream \
206 << "\n"
207 << "===============================================================================\n"\
208 << "| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
209 << "===============================================================================\n";
210
211 try{
212 std::vector<std::vector<int> > allTags = tetBasis.getAllDofTags();
213
214 // Loop over all tags, lookup the associated dof enumeration and then lookup the tag again
215 for (unsigned i = 0; i < allTags.size(); i++) {
216 int bfOrd = tetBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
217
218 std::vector<int> myTag = tetBasis.getDofTag(bfOrd);
219 if( !( (myTag[0] == allTags[i][0]) &&
220 (myTag[1] == allTags[i][1]) &&
221 (myTag[2] == allTags[i][2]) &&
222 (myTag[3] == allTags[i][3]) ) ) {
223 errorFlag++;
224 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
225 *outStream << " getDofOrdinal( {"
226 << allTags[i][0] << ", "
227 << allTags[i][1] << ", "
228 << allTags[i][2] << ", "
229 << allTags[i][3] << "}) = " << bfOrd <<" but \n";
230 *outStream << " getDofTag(" << bfOrd << ") = { "
231 << myTag[0] << ", "
232 << myTag[1] << ", "
233 << myTag[2] << ", "
234 << myTag[3] << "}\n";
235 }
236 }
237
238 // Now do the same but loop over basis functions
239 for( int bfOrd = 0; bfOrd < tetBasis.getCardinality(); bfOrd++) {
240 std::vector<int> myTag = tetBasis.getDofTag(bfOrd);
241 int myBfOrd = tetBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
242 if( bfOrd != myBfOrd) {
243 errorFlag++;
244 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
245 *outStream << " getDofTag(" << bfOrd << ") = { "
246 << myTag[0] << ", "
247 << myTag[1] << ", "
248 << myTag[2] << ", "
249 << myTag[3] << "} but getDofOrdinal({"
250 << myTag[0] << ", "
251 << myTag[1] << ", "
252 << myTag[2] << ", "
253 << myTag[3] << "} ) = " << myBfOrd << "\n";
254 }
255 }
256 }
257 catch (const std::logic_error & err){
258 *outStream << err.what() << "\n\n";
259 errorFlag = -1000;
260 };
261
262 *outStream \
263 << "\n"
264 << "===============================================================================\n"\
265 << "| TEST 3: correctness of basis function values |\n"\
266 << "===============================================================================\n";
267
268 outStream -> precision(20);
269
270 // VALUE: Each row pair gives the 6x3 correct basis set values at an evaluation point: (P,F,D) layout
271 double basisValues[] = {
272 // 4 vertices
273 1.0,0.,0., 0.,0.,0., 0.,-1.0,0., 0.,0.,1.0,
274 0.,0.,0., 0.,0.,0.,
275
276 1.0,1.0,1.0, 0.,1.,0., 0.,0.,0., 0.,0.,0.,
277 0.,0.,1., 0.,0.,0.,
278
279 0.,0.,0., -1.,0.,0., -1.0,-1.0,-1.0,
280 0.,0.,0., 0.,0.,0., 0.,0.,1.,
281
282 0.,0.,0., 0.,0.,0., 0.,0.,0., 1.0,1.0,1.0,
283 -1.,0.,0., 0.,-1.,0.,
284
285 // 6 edge centers
286 1.0,0.5,0.5, 0.,0.5,0., 0.,-0.5,0.,
287 0.,0.,0.5, 0.,0.,0.5, 0.,0.,0.,
288
289 0.5,0.5,0.5, -0.5,0.5,0.,
290 -0.5,-0.5,-0.5, 0.,0.,0., 0.,0.,0.5, 0.,0.,0.5,
291
292 0.5,0.,0., -0.5,0.,0., -0.5,-1.0,-0.5,
293 0.,0.,0.5, 0.,0.,0., 0.,0.,0.5,
294
295 0.5,0.,0., 0.,0.,0., 0.,-0.5,0., 0.5,0.5,1.0,
296 -0.5,0.,0., 0.,-0.5,0.,
297
298 0.5,0.5,0.5, 0.,0.5,0., 0.,0.,0., 0.5,0.5,0.5,
299 -0.5,0.,0.5, 0.,-0.5,0.,
300
301 0.,0.,0., -0.5,0.,0., -0.5,-0.5,-0.5, 0.5,0.5,0.5,
302 -0.5,0.,0., 0.,-0.5,0.5
303 };
304
305 // CURL: each row pair gives the 3x12 correct values of the curls of the 12 basis functions: (P,F,D) layout
306 double basisCurls[] = {
307 // 4 vertices
308 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
309 0.,-2.0,0., 2.0,0.,0.,
310
311 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
312 0.,-2.0,0., 2.0,0.,0.,
313
314 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
315 0.,-2.0,0., 2.0,0.,0.,
316
317 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
318 0.,-2.0,0., 2.0,0.,0.,
319
320 // 6 edge centers
321 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
322 0.,-2.0,0., 2.0,0.,0.,
323
324 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
325 0.,-2.0,0., 2.0,0.,0.,
326
327 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
328 0.,-2.0,0., 2.0,0.,0.,
329
330 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
331 0.,-2.0,0., 2.0,0.,0.,
332
333 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
334 0.,-2.0,0., 2.0,0.,0.,
335
336 0.,-2.0,2.0, 0.,0.,2.0, -2.0,0.,2.0, -2.0,2.0,0.,
337 0.,-2.0,0., 2.0,0.,0.,
338 };
339
340 try{
341
342 // Dimensions for the output arrays:
343 int numFields = tetBasis.getCardinality();
344 int numPoints = tetNodes.dimension(0);
345 int spaceDim = tetBasis.getBaseCellTopology().getDimension();
346
347 // Generic array for values and curls that will be properly sized before each call
349
350 // Check VALUE of basis functions: resize vals to rank-3 container:
351 vals.resize(numFields, numPoints, spaceDim);
352 tetBasis.getValues(vals, tetNodes, OPERATOR_VALUE);
353 for (int i = 0; i < numFields; i++) {
354 for (int j = 0; j < numPoints; j++) {
355 for (int k = 0; k < spaceDim; k++) {
356
357 // compute offset for (P,F,D) data layout: indices are P->j, F->i, D->k
358 int l = k + i * spaceDim + j * spaceDim * numFields;
359 if (std::abs(vals(i,j,k) - basisValues[l]) > INTREPID_TOL) {
360 errorFlag++;
361 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
362
363 // Output the multi-index of the value where the error is:
364 *outStream << " At multi-index { ";
365 *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
366 *outStream << "} computed value: " << vals(i,j,k)
367 << " but reference value: " << basisValues[l] << "\n";
368 }
369 }
370 }
371 }
372
373 // Check CURL of basis function: resize vals to rank-3 container
374 vals.resize(numFields, numPoints, spaceDim);
375 tetBasis.getValues(vals, tetNodes, OPERATOR_CURL);
376 for (int i = 0; i < numFields; i++) {
377 for (int j = 0; j < numPoints; j++) {
378 for (int k = 0; k < spaceDim; k++) {
379
380 // compute offset for (P,F,D) data layout: indices are P->j, F->i, D->k
381 int l = k + i * spaceDim + j * spaceDim * numFields;
382 if (std::abs(vals(i,j,k) - basisCurls[l]) > INTREPID_TOL) {
383 errorFlag++;
384 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
385
386 // Output the multi-index of the value where the error is:
387 *outStream << " At multi-index { ";
388 *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
389 *outStream << "} computed curl component: " << vals(i,j,k)
390 << " but reference curl component: " << basisCurls[l] << "\n";
391 }
392 }
393 }
394 }
395
396 }
397
398 // Catch unexpected errors
399 catch (const std::logic_error & err) {
400 *outStream << err.what() << "\n\n";
401 errorFlag = -1000;
402 };
403
404 *outStream \
405 << "\n"
406 << "===============================================================================\n"\
407 << "| TEST 4: correctness of DoF locations |\n"\
408 << "===============================================================================\n";
409
410 try{
411 Teuchos::RCP<Basis<double, FieldContainer<double> > > basis =
412 Teuchos::rcp(new Basis_HCURL_TET_I1_FEM<double, FieldContainer<double> >);
413 Teuchos::RCP<DofCoordsInterface<FieldContainer<double> > > coord_iface =
414 Teuchos::rcp_dynamic_cast<DofCoordsInterface<FieldContainer<double> > >(basis);
415
416 int spaceDim = 3;
418 FieldContainer<double> bvals(basis->getCardinality(), basis->getCardinality(),spaceDim); // last dimension is spatial dim
419
420 // Check exceptions.
421#ifdef HAVE_INTREPID_DEBUG
422 cvals.resize(1,2,3);
423 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
424 cvals.resize(3,2);
425 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
426 cvals.resize(4,2);
427 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException );
428#endif
429 cvals.resize(6,spaceDim);
430 INTREPID_TEST_COMMAND( coord_iface->getDofCoords(cvals), throwCounter, nException ); nException--;
431 // Check if number of thrown exceptions matches the one we expect
432 if (throwCounter != nException) {
433 errorFlag++;
434 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
435 }
436
437 // Check mathematical correctness
438 FieldContainer<double> tangents(basis->getCardinality(),spaceDim); // tangents at each point basis point
439 tangents(0,0) = 1.0; tangents(0,1) = 0.0; tangents(0,2) = 0.0;
440 tangents(1,0) = -1.0; tangents(1,1) = 1.0; tangents(1,2) = 0.0;
441 tangents(2,0) = 0.0; tangents(2,1) = -1.0; tangents(2,2) = 0.0;
442 tangents(3,0) = 0.0; tangents(3,1) = 0.0; tangents(3,2) = 1.0;
443 tangents(4,0) = -1.0; tangents(4,1) = 0.0; tangents(4,2) = 1.0;
444 tangents(5,0) = 0.0; tangents(5,1) = -1.0; tangents(5,2) = 1.0;
445
446 basis->getValues(bvals, cvals, OPERATOR_VALUE);
447 char buffer[120];
448 for (int i=0; i<bvals.dimension(0); i++) {
449 for (int j=0; j<bvals.dimension(1); j++) {
450
451 double tangent = 0.0;
452 for(int d=0;d<spaceDim;d++)
453 tangent += bvals(i,j,d)*tangents(j,d);
454
455 if ((i != j) && (std::abs(tangent - 0.0) > INTREPID_TOL)) {
456 errorFlag++;
457 sprintf(buffer, "\nValue of basis function %d at (%6.4e, %6.4e, %6.4e) is %6.4e but should be %6.4e!\n", i, cvals(i,0), cvals(i,1), cvals(i,2), tangent, 0.0);
458 *outStream << buffer;
459 }
460 else if ((i == j) && (std::abs(tangent - 1.0) > INTREPID_TOL)) {
461 errorFlag++;
462 sprintf(buffer, "\nValue of basis function %d at (%6.4e, %6.4e, %6.4e) is %6.4e but should be %6.4e!\n", i, cvals(i,0), cvals(i,1), cvals(i,2), tangent, 1.0);
463 *outStream << buffer;
464 }
465 }
466 }
467
468 }
469 catch (const std::logic_error & err){
470 *outStream << err.what() << "\n\n";
471 errorFlag = -1000;
472 };
473
474 if (errorFlag != 0)
475 std::cout << "End Result: TEST FAILED\n";
476 else
477 std::cout << "End Result: TEST PASSED\n";
478
479 // reset format state of std::cout
480 std::cout.copyfmt(oldFormatState);
481
482 return errorFlag;
483}
Header file for utility class to provide multidimensional containers.
Header file for the Intrepid::HCURL_TET_I1_FEM class.
Implementation of the default H(curl)-compatible FEM basis of degree 1 on Tetrahedron cell.
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
Evaluation of a FEM basis on a reference Tetrahedron cell.
virtual int getDofOrdinal(const int subcDim, const int subcOrd, const int subcDofOrd)
DoF tag to ordinal lookup.
virtual const std::vector< std::vector< int > > & getAllDofTags()
Retrieves all DoF tags.
virtual const shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation http://trilin...
virtual const std::vector< int > & getDofTag(const int dofOrd)
DoF ordinal to DoF tag lookup.
virtual int getCardinality() const
Returns cardinality of the basis.
Implementation of a templated lexicographical container for a multi-indexed scalar quantity....
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0.