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
51#include "Teuchos_oblackholestream.hpp"
52#include "Teuchos_RCP.hpp"
53#include "Teuchos_GlobalMPISession.hpp"
54
55using namespace std;
56using namespace Intrepid;
57
58#define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \
59{ \
60 ++nException; \
61 try { \
62 S ; \
63 } \
64 catch (const std::logic_error & err) { \
65 ++throwCounter; \
66 *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
67 *outStream << err.what() << '\n'; \
68 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \
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_HGRAD_TRI_C2_FEM) |\n" \
93 << "| |\n" \
94 << "| 1) Conversion of Dof tags into Dof ordinals and back |\n" \
95 << "| 2) Basis values for VALUE, GRAD, CURL, and Dk operators |\n" \
96 << "| |\n" \
97 << "| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \
98 << "| Denis Ridzal (dridzal@sandia.gov), |\n" \
99 << "| Kara Peterson (dridzal@sandia.gov). |\n" \
100 << "| |\n" \
101 << "| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
102 << "| Trilinos website: http://trilinos.sandia.gov |\n" \
103 << "| |\n" \
104 << "===============================================================================\n"\
105 << "| TEST 1: Basis creation, exception testing |\n"\
106 << "===============================================================================\n";
107
108
109 // Define basis and error flag
111 int errorFlag = 0;
112
113 // Initialize throw counter for exception testing
114 int nException = 0;
115 int throwCounter = 0;
116
117 // The unisolvent set for the default Lagrangian basis of degree 2 uses the standard equispaced
118 // lattice on Triangle consisting of the 3 vertices and the 3 edge midpoints. To check basis
119 // correctness it suffices to evaluate each basis function at the lattice points, defined below,
120 // but we throw in an additional random point.
121 FieldContainer<double> triNodes(7, 2);
122 triNodes(0,0) = 0.0; triNodes(0,1) = 0.0;
123 triNodes(1,0) = 1.0; triNodes(1,1) = 0.0;
124 triNodes(2,0) = 0.0; triNodes(2,1) = 1.0;
125 triNodes(3,0) = 0.5; triNodes(3,1) = 0.0;
126 triNodes(4,0) = 0.5; triNodes(4,1) = 0.5;
127 triNodes(5,0) = 0.0; triNodes(5,1) = 0.5;
128 triNodes(6,0) = 0.1732; triNodes(6,1) = 0.6714;
129
130 // Generic array for the output values; needs to be properly resized depending on the operator type
132
133 try{
134 // exception #1: DIV cannot be applied to scalar functions
135 // resize vals to rank-2 container with dimensions (num. points, num. basis functions)
136 vals.resize(triBasis.getCardinality(), triNodes.dimension(0) );
137 INTREPID_TEST_COMMAND( triBasis.getValues(vals, triNodes, OPERATOR_DIV), throwCounter, nException );
138
139 // Exceptions 2-6: all bf tags/bf Ids below are wrong and should cause getDofOrdinal() and
140 // getDofTag() to access invalid array elements thereby causing bounds check exception
141 // exception #2
142 INTREPID_TEST_COMMAND( triBasis.getDofOrdinal(2,0,0), throwCounter, nException );
143 // exception #3
144 INTREPID_TEST_COMMAND( triBasis.getDofOrdinal(1,1,1), throwCounter, nException );
145 // exception #4
146 INTREPID_TEST_COMMAND( triBasis.getDofOrdinal(0,4,0), throwCounter, nException );
147 // exception #5
148 INTREPID_TEST_COMMAND( triBasis.getDofTag(6), throwCounter, nException );
149 // exception #6
150 INTREPID_TEST_COMMAND( triBasis.getDofTag(-1), throwCounter, nException );
151
152#ifdef HAVE_INTREPID_DEBUG
153 // Exceptions 7-17 test exception handling with incorrectly dimensioned input/output arrays
154 // exception #7: input points array must be of rank-2
155 FieldContainer<double> badPoints1(4, 5, 3);
156 INTREPID_TEST_COMMAND( triBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
157
158 // exception #8 dimension 1 in the input point array must equal space dimension of the cell
159 FieldContainer<double> badPoints2(4, triBasis.getBaseCellTopology().getDimension() + 1);
160 INTREPID_TEST_COMMAND( triBasis.getValues(vals, badPoints2, OPERATOR_VALUE), throwCounter, nException );
161
162 // exception #9 output values must be of rank-2 for OPERATOR_VALUE
163 FieldContainer<double> badVals1(4, 3, 1);
164 INTREPID_TEST_COMMAND( triBasis.getValues(badVals1, triNodes, OPERATOR_VALUE), throwCounter, nException );
165
166 // exception #10 output values must be of rank-3 for OPERATOR_GRAD
167 FieldContainer<double> badVals2(4, 3);
168 INTREPID_TEST_COMMAND( triBasis.getValues(badVals2, triNodes, OPERATOR_GRAD), throwCounter, nException );
169
170 // exception #11 output values must be of rank-3 for OPERATOR_CURL
171 INTREPID_TEST_COMMAND( triBasis.getValues(badVals2, triNodes, OPERATOR_CURL), throwCounter, nException );
172
173 // exception #12 output values must be of rank-3 for OPERATOR_D2
174 INTREPID_TEST_COMMAND( triBasis.getValues(badVals2, triNodes, OPERATOR_D2), throwCounter, nException );
175
176 // exception #13 incorrect 1st dimension of output array (must equal number of basis functions)
177 FieldContainer<double> badVals3(triBasis.getCardinality() + 1, triNodes.dimension(0));
178 INTREPID_TEST_COMMAND( triBasis.getValues(badVals3, triNodes, OPERATOR_VALUE), throwCounter, nException );
179
180 // exception #14 incorrect 0th dimension of output array (must equal number of points)
181 FieldContainer<double> badVals4(triBasis.getCardinality(), triNodes.dimension(0) + 1);
182 INTREPID_TEST_COMMAND( triBasis.getValues(badVals4, triNodes, OPERATOR_VALUE), throwCounter, nException );
183
184 // exception #15: incorrect 2nd dimension of output array (must equal the space dimension)
185 FieldContainer<double> badVals5(triBasis.getCardinality(), triNodes.dimension(0), triBasis.getBaseCellTopology().getDimension() + 1);
186 INTREPID_TEST_COMMAND( triBasis.getValues(badVals5, triNodes, OPERATOR_GRAD), throwCounter, nException );
187
188 // exception #16: incorrect 2nd dimension of output array (must equal D2 cardinality in 2D)
189 FieldContainer<double> badVals6(triBasis.getCardinality(), triNodes.dimension(0), 40);
190 INTREPID_TEST_COMMAND( triBasis.getValues(badVals6, triNodes, OPERATOR_D2), throwCounter, nException );
191
192 // exception #17: incorrect 2nd dimension of output array (must equal D3 cardinality in 2D)
193 INTREPID_TEST_COMMAND( triBasis.getValues(badVals6, triNodes, OPERATOR_D3), throwCounter, nException );
194#endif
195
196 }
197 catch (const std::logic_error & err) {
198 *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
199 *outStream << err.what() << '\n';
200 *outStream << "-------------------------------------------------------------------------------" << "\n\n";
201 errorFlag = -1000;
202 };
203
204 // Check if number of thrown exceptions matches the one we expect
205 if (throwCounter != nException) {
206 errorFlag++;
207 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
208 }
209
210 *outStream \
211 << "\n"
212 << "===============================================================================\n"\
213 << "| TEST 2: correctness of tag to enum and enum to tag lookups |\n"\
214 << "===============================================================================\n";
215
216 try{
217 std::vector<std::vector<int> > allTags = triBasis.getAllDofTags();
218
219 // Loop over all tags, lookup the associated dof enumeration and then lookup the tag again
220 for (unsigned i = 0; i < allTags.size(); i++) {
221 int bfOrd = triBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
222
223 std::vector<int> myTag = triBasis.getDofTag(bfOrd);
224 if( !( (myTag[0] == allTags[i][0]) &&
225 (myTag[1] == allTags[i][1]) &&
226 (myTag[2] == allTags[i][2]) &&
227 (myTag[3] == allTags[i][3]) ) ) {
228 errorFlag++;
229 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
230 *outStream << " getDofOrdinal( {"
231 << allTags[i][0] << ", "
232 << allTags[i][1] << ", "
233 << allTags[i][2] << ", "
234 << allTags[i][3] << "}) = " << bfOrd <<" but \n";
235 *outStream << " getDofTag(" << bfOrd << ") = { "
236 << myTag[0] << ", "
237 << myTag[1] << ", "
238 << myTag[2] << ", "
239 << myTag[3] << "}\n";
240 }
241 }
242
243 // Now do the same but loop over basis functions
244 for( int bfOrd = 0; bfOrd < triBasis.getCardinality(); bfOrd++) {
245 std::vector<int> myTag = triBasis.getDofTag(bfOrd);
246 int myBfOrd = triBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
247 if( bfOrd != myBfOrd) {
248 errorFlag++;
249 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
250 *outStream << " getDofTag(" << bfOrd << ") = { "
251 << myTag[0] << ", "
252 << myTag[1] << ", "
253 << myTag[2] << ", "
254 << myTag[3] << "} but getDofOrdinal({"
255 << myTag[0] << ", "
256 << myTag[1] << ", "
257 << myTag[2] << ", "
258 << myTag[3] << "} ) = " << myBfOrd << "\n";
259 }
260 }
261 }
262 catch (const std::logic_error & err){
263 *outStream << err.what() << "\n\n";
264 errorFlag = -1000;
265 };
266
267 *outStream \
268 << "\n"
269 << "===============================================================================\n"\
270 << "| TEST 3: correctness of basis function values |\n"\
271 << "===============================================================================\n";
272
273 outStream -> precision(20);
274
275 // VALUE: Correct basis values in (F,P) format: each row gives the 6 correct basis values at
276 // the lattice points (3 vertices and 3 edge midpoints)
277 double basisValues[] = {
278 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, -0.10710168,
279 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, -0.11320352,
280 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.23015592,
281 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.10766112,
282 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.46514592,
283 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.41734224
284 };
285
286 // GRAD and D1: Correct gradients and D1 in (F,P,D) format: each row contains 6x2 values of
287 // gradients of basis functions. Same array is used to check correctness of CURL.
288 double basisGrads[] = {
289 // V0 V1 V2 M0 M1 M2 RP
290 -3.0, -3.0, 1.0, 1.0, 1.0, 1.0, -1.0, -1.0, 1.0, 1.0, -1.0,-1.0, 0.3784, 0.3784,
291 -1.0, 0.0, 3.0, 0.0, -1.0, 0.0, 1.0, 0.0, 1.0, 0.0, -1.0, 0.0, -0.3072, 0.0,
292 0.0, -1.0, 0.0,-1.0, 0.0, 3.0, 0.0, -1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.6856,
293 4.0, 0.0, -4.0,-4.0, 0.0, 0.0, 0.0, -2.0, -2.0,-2.0, 2.0, 0.0, -0.0712, -0.6928,
294 0.0, 0.0, 0.0, 4.0, 4.0, 0.0, 0.0, 2.0, 2.0, 2.0, 2.0, 0.0, 2.6856, 0.6928,
295 0.0, 4.0, 0.0, 0.0, -4.0,-4.0, 0.0, 2.0, -2.0,-2.0, -2.0, 0.0, -2.6856, -2.064
296 };
297
298 // D2: Correct multiset of second order partials in (F,P,Dk) format. D2 cardinality = 3 for 2D
299 double basisD2[] = {
300 // V0 V1 V2 M0 M1 M2 RP
301 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0,
302 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0,
303 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0,
304 -8.0,-4.0, 0.0, -8.0,-4.0, 0.0, -8.0,-4.0, 0.0, -8.0,-4.0, 0.0, -8.0,-4.0, 0.0, -8.0,-4.0, 0.0, -8.0,-4.0, 0.0,
305 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0, 0.0, 4.0, 0.0,
306 0.0,-4.0,-8.0, 0.0,-4.0,-8.0, 0.0,-4.0,-8.0, 0.0,-4.0,-8.0, 0.0,-4.0,-8.0, 0.0,-4.0,-8.0, 0.0,-4.0,-8.0
307 };
308 try{
309
310 // Dimensions for the output arrays:
311 int numFields = triBasis.getCardinality();
312 int numPoints = triNodes.dimension(0);
313 int spaceDim = triBasis.getBaseCellTopology().getDimension();
314
315 // Generic array for values, grads, curls, etc. that will be properly sized before each call
317
318 // Check VALUE of basis functions: resize vals to rank-2 container:
319 vals.resize(numFields, numPoints);
320 triBasis.getValues(vals, triNodes, OPERATOR_VALUE);
321 for (int i = 0; i < numFields; i++) {
322 for (int j = 0; j < numPoints; j++) {
323
324 // Compute offset for (F,P) container
325 int l = j + i * numPoints;
326 if (std::abs(vals(i,j) - basisValues[l]) > INTREPID_TOL) {
327 errorFlag++;
328 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
329
330 // Output the multi-index of the value where the error is:
331 *outStream << " At multi-index { ";
332 *outStream << i << " ";*outStream << j << " ";
333 *outStream << "} computed value: " << vals(i,j)
334 << " but reference value: " << basisValues[l] << "\n";
335 }
336 }
337 }
338
339 // Check GRAD of basis function: resize vals to rank-3 container
340 vals.resize(numFields, numPoints, spaceDim);
341 triBasis.getValues(vals, triNodes, OPERATOR_GRAD);
342 for (int i = 0; i < numFields; i++) {
343 for (int j = 0; j < numPoints; j++) {
344 for (int k = 0; k < spaceDim; k++) {
345 // basisGrads is (F,P,D), compute offset:
346 int l = k + j * spaceDim + i * spaceDim * numPoints;
347 if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
348 errorFlag++;
349 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
350
351 // Output the multi-index of the value where the error is:
352 *outStream << " At multi-index { ";
353 *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
354 *outStream << "} computed grad component: " << vals(i,j,k)
355 << " but reference grad component: " << basisGrads[l] << "\n";
356 }
357 }
358 }
359 }
360
361 // Check D1 of basis function (do not resize vals because it has the correct size: D1 = GRAD)
362 triBasis.getValues(vals, triNodes, OPERATOR_D1);
363 for (int i = 0; i < numFields; i++) {
364 for (int j = 0; j < numPoints; j++) {
365 for (int k = 0; k < spaceDim; k++) {
366 // basisGrads is (F,P,D), compute offset:
367 int l = k + j * spaceDim + i * spaceDim * numPoints;
368 if (std::abs(vals(i,j,k) - basisGrads[l]) > INTREPID_TOL) {
369 errorFlag++;
370 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
371
372 // Output the multi-index of the value where the error is:
373 *outStream << " At multi-index { ";
374 *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
375 *outStream << "} computed D1 component: " << vals(i,j,k)
376 << " but reference D1 component: " << basisGrads[l] << "\n";
377 }
378 }
379 }
380 }
381
382 // Check CURL of basis function: resize vals just for illustration!
383 vals.resize(numFields, numPoints, spaceDim);
384 triBasis.getValues(vals, triNodes, OPERATOR_CURL);
385 for (int i = 0; i < numFields; i++) {
386 for (int j = 0; j < numPoints; j++) {
387 // We will use "rotated" basisGrads to check CURL: get offsets to extract (u_y, -u_x)
388 int curl_0 = 1 + j * spaceDim + i * spaceDim * numPoints; // position of y-derivative
389 int curl_1 = 0 + j * spaceDim + i * spaceDim * numPoints; // position of x-derivative
390
391 double curl_value_0 = basisGrads[curl_0];
392 double curl_value_1 =-basisGrads[curl_1];
393 if (std::abs(vals(i,j,0) - curl_value_0) > INTREPID_TOL) {
394 errorFlag++;
395 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
396 // Output the multi-index of the value where the error is:
397 *outStream << " At multi-index { ";
398 *outStream << i << " ";*outStream << j << " ";*outStream << 0 << " ";
399 *outStream << "} computed curl component: " << vals(i,j,0)
400 << " but reference curl component: " << curl_value_0 << "\n";
401 }
402 if (std::abs(vals(i,j,1) - curl_value_1) > INTREPID_TOL) {
403 errorFlag++;
404 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
405 // Output the multi-index of the value where the error is:
406 *outStream << " At multi-index { ";
407 *outStream << i << " ";*outStream << j << " ";*outStream << 1 << " ";
408 *outStream << "} computed curl component: " << vals(i,j,1)
409 << " but reference curl component: " << curl_value_1 << "\n";
410 }
411 }
412 }
413
414 // Check D2 of basis function
415 int D2cardinality = Intrepid::getDkCardinality(OPERATOR_D2, spaceDim);
416 vals.resize(numFields, numPoints, D2cardinality);
417 triBasis.getValues(vals, triNodes, OPERATOR_D2);
418 for (int i = 0; i < numFields; i++) {
419 for (int j = 0; j < numPoints; j++) {
420 for (int k = 0; k < D2cardinality; k++) {
421
422 // basisD2 is (F,P,Dk), compute offset:
423 int l = k + j * D2cardinality + i * D2cardinality * numPoints;
424 if (std::abs(vals(i,j,k) - basisD2[l]) > INTREPID_TOL) {
425 errorFlag++;
426 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
427
428 // Output the multi-index of the value where the error is:
429 *outStream << " At multi-index { ";
430 *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
431 *outStream << "} computed D2 component: " << vals(i,j,k)
432 << " but reference D2 component: " << basisD2[l] << "\n";
433 }
434 }
435 }
436 }
437
438 // Check all higher derivatives - must be zero.
439 for(EOperator op = OPERATOR_D3; op < OPERATOR_MAX; op++) {
440
441 // The last dimension is the number of kth derivatives and needs to be resized for every Dk
442 int DkCardin = Intrepid::getDkCardinality(op, spaceDim);
443 vals.resize(numFields, numPoints, DkCardin);
444
445 triBasis.getValues(vals, triNodes, op);
446 for (int i = 0; i < vals.size(); i++) {
447 if (std::abs(vals[i]) > INTREPID_TOL) {
448 errorFlag++;
449 *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
450
451 // Get the multi-index of the value where the error is and the operator order
452 std::vector<int> myIndex;
453 vals.getMultiIndex(myIndex,i);
454 int ord = Intrepid::getOperatorOrder(op);
455 *outStream << " At multi-index { ";
456 for(int j = 0; j < vals.rank(); j++) {
457 *outStream << myIndex[j] << " ";
458 }
459 *outStream << "} computed D"<< ord <<" component: " << vals[i]
460 << " but reference D" << ord << " component: 0 \n";
461 }
462 }
463 }
464 }
465
466 // Catch unexpected errors
467 catch (const std::logic_error & err) {
468 *outStream << err.what() << "\n\n";
469 errorFlag = -1000;
470 };
471
472 if (errorFlag != 0)
473 std::cout << "End Result: TEST FAILED\n";
474 else
475 std::cout << "End Result: TEST PASSED\n";
476
477 // reset format state of std::cout
478 std::cout.copyfmt(oldFormatState);
479
480 return errorFlag;
481}
Header file for utility class to provide multidimensional containers.
Header file for the Intrepid::HGRAD_TRI_C2_FEM class.
int getOperatorOrder(const EOperator operatorType)
Returns order of an operator.
int getDkCardinality(const EOperator operatorType, const int spaceDim)
Returns cardinality of Dk, i.e., the number of all derivatives of order k.
Implementation of the default H(grad)-compatible FEM basis of degree 2 on Triangle cell.
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
Evaluation of a FEM basis on a reference Triangle 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....
int size() const
Returns size of the FieldContainer defined as the product of its dimensions.
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0.
void getMultiIndex(int &i0, const int valueEnum) const
Returns the multi-index of a value, based on its enumeration, as a list, for rank-1 containers.
int rank() const
Return rank of the FieldContainer = number of indices used to tag the multi-indexed value.