50 #include "Teuchos_oblackholestream.hpp"
51 #include "Teuchos_RCP.hpp"
52 #include "Teuchos_GlobalMPISession.hpp"
57 using namespace Intrepid;
59 #define INTREPID_TEST_COMMAND( S , throwCounter, nException ) \
65 catch (const std::logic_error & err) { \
67 *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
68 *outStream << err.what() << '\n'; \
69 *outStream << "-------------------------------------------------------------------------------" << "\n\n"; \
73 int main(
int argc,
char *argv[]) {
75 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
79 int iprint = argc - 1;
80 Teuchos::RCP<std::ostream> outStream;
81 Teuchos::oblackholestream bhs;
83 outStream = Teuchos::rcp(&std::cout,
false);
85 outStream = Teuchos::rcp(&bhs,
false);
88 Teuchos::oblackholestream oldFormatState;
89 oldFormatState.copyfmt(std::cout);
92 <<
"===============================================================================\n" \
94 <<
"| Unit Test (Basis_HGRAD_TET_COMP12_FEM) |\n" \
96 <<
"| 1) Evaluation of Basis Function Values |\n" \
98 <<
"| Questions? Contact Pavel Bochev (pbboche@sandia.gov), |\n" \
99 <<
"| Denis Ridzal (dridzal@sandia.gov), |\n" \
100 <<
"| Kara Peterson (kjpeter@sandia.gov), |\n" \
101 <<
"| Jake Ostien (jtostie@sandia.gov). |\n" \
103 <<
"| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
104 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" \
106 <<
"===============================================================================\n"\
107 <<
"| TEST 1: Basis creation, values |\n"\
108 <<
"===============================================================================\n";
120 tetNodes(0,0) = 0.0; tetNodes(0,1) = 0.0; tetNodes(0,2) = 0.0;
121 tetNodes(1,0) = 1.0; tetNodes(1,1) = 0.0; tetNodes(1,2) = 0.0;
122 tetNodes(2,0) = 0.0; tetNodes(2,1) = 1.0; tetNodes(2,2) = 0.0;
123 tetNodes(3,0) = 0.0; tetNodes(3,1) = 0.0; tetNodes(3,2) = 1.0;
124 tetNodes(4,0) = 0.5; tetNodes(4,1) = 0.0; tetNodes(4,2) = 0.0;
125 tetNodes(5,0) = 0.5; tetNodes(5,1) = 0.5; tetNodes(5,2) = 0.0;
126 tetNodes(6,0) = 0.0; tetNodes(6,1) = 0.5; tetNodes(6,2) = 0.0;
127 tetNodes(7,0) = 0.0; tetNodes(7,1) = 0.0; tetNodes(7,2) = 0.5;
128 tetNodes(8,0) = 0.5; tetNodes(8,1) = 0.0; tetNodes(8,2) = 0.5;
129 tetNodes(9,0) = 0.0; tetNodes(9,1) = 0.5; tetNodes(9,2) = 0.5;
133 tetPoints(0,0) = 0.25; tetPoints(0,1) = 0.25; tetPoints(0,2) = 0.25;
134 tetPoints(1,0) = 0.5; tetPoints(1,1) = (1./6.); tetPoints(1,2) = (1./6.);
135 tetPoints(2,0) = (1./6.); tetPoints(2,1) = 0.5; tetPoints(2,2) = (1./6.);
136 tetPoints(3,0) = (1./6.); tetPoints(3,1) = (1./6.); tetPoints(3,2) = 0.5;
137 tetPoints(4,0) = (1./6.); tetPoints(4,1) = (1./6.); tetPoints(4,2) = (1./6.);
139 tetPoints(5,0) = 0.1381966011250105151795413165634361882280;
140 tetPoints(5,1) = 0.1381966011250105151795413165634361882280;
141 tetPoints(5,2) = 0.1381966011250105151795413165634361882280;
143 tetPoints(6,0) = 0.5854101966249684544613760503096914353161;
144 tetPoints(6,1) = 0.1381966011250105151795413165634361882280;
145 tetPoints(6,2) = 0.1381966011250105151795413165634361882280;
147 tetPoints(7,0) = 0.1381966011250105151795413165634361882280;
148 tetPoints(7,1) = 0.5854101966249684544613760503096914353161;
149 tetPoints(7,2) = 0.1381966011250105151795413165634361882280;
151 tetPoints(8,0) = 0.1381966011250105151795413165634361882280;
152 tetPoints(8,1) = 0.1381966011250105151795413165634361882280;
153 tetPoints(8,2) = 0.5854101966249684544613760503096914353161;
156 outStream -> precision(20);
159 double nodalBasisValues[] = {
161 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
162 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
163 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
164 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
166 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0,
167 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0,
168 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0,
169 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0,
170 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0,
171 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0
173 double pointBasisValues[] = {
175 0.0, 0.0, 0.0, 0.0, 1./6., 1./6., 1./6., 1./6., 1./6., 1./6.,
177 0.0, 0.0, 0.0, 0.0, 1./3., 1./3., 0.0, 0.0, 1./3., 0.0,
179 0.0, 0.0, 0.0, 0.0, 0.0, 1./3., 1./3., 0.0, 0.0, 1./3.,
181 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1./3., 1./3., 1./3.,
183 0.0, 0.0, 0.0, 0.0, 1./3., 0.0, 1./3., 1./3., 0.0, 0.0,
185 0.170820393249936908922752100619382870632, 0.0, 0.0, 0.0, 0.276393202250021030359082633126872376456, 0.0, 0.276393202250021030359082633126872376456, 0.276393202250021030359082633126872376456, 0.0, 0.0,
187 0.0, 0.170820393249936908922752100619382870632, 0.0, 0.0, 0.276393202250021030359082633126872376456, 0.276393202250021030359082633126872376456, 0.0, 0.0, 0.276393202250021030359082633126872376456, 0.0,
189 0.0, 0.0, 0.170820393249936908922752100619382870632, 0.0, 0.0, 0.276393202250021030359082633126872376456, 0.276393202250021030359082633126872376456, 0.0, 0.0, 0.276393202250021030359082633126872376456,
191 0.0, 0.0, 0.0, 0.170820393249936908922752100619382870632, 0.0, 0.0, 0.0, 0.276393202250021030359082633126872376456, 0.276393202250021030359082633126872376456, 0.276393202250021030359082633126872376456,
195 double pointBasisGrads[] = {
197 -1./4., -1./4., -1./4., \
201 0.0, -3./4., -3./4., \
203 -3./4., 0.0, -3./4., \
204 -3./4., -3./4., 0.0, \
209 -1./24., -1./24., -1./24., \
213 -35./36., -19./12., -19./12., \
214 11./18., 19./12., 0.0, \
215 -17./36., 0.0, -1./3., \
216 -17./36., -1./3., 0.0, \
217 11./18., 0.0, 19./12., \
218 -5./36., 1./3., 1./3., \
221 -1./24., -1./24., -1./24., \
225 0.0, -17./36., -1./3., \
226 19./12., 11./18., 0.0, \
227 -19./12., -35./36., -19./12., \
228 -1./3., -17./36., 0.0, \
229 1./3., -5./36., 1./3., \
230 0.0, 11./18., 19./12., \
233 -1./24., -1./24., -1./24., \
237 0.0, -1./3., -17./36., \
238 1./3., 1./3., -5./36., \
239 -1./3., 0.0, -17./36., \
240 -19./12., -19./12., -35./36., \
241 19./12., 0.0, 11./18., \
242 0.0, 19./12., 11./18., \
245 -7./8., -7./8., -7./8., \
249 35./36., -11./18., -11./18., \
250 17./36., 17./36., 5./36., \
251 -11./18., 35./36., -11./18., \
252 -11./18., -11./18., 35./36., \
253 17./36., 5./36., 17./36., \
254 5./36., 17./36., 17./36., \
257 -1.088525491562421136153440125774228588290, -1.088525491562421136153440125774228588290, -1.088525491562421136153440125774228588290, \
258 -0.029508497187473712051146708591409529430, 0.0, 0.0, \
259 0.0, -0.029508497187473712051146708591409529430, 0.0, \
260 0.0, 0.0, -0.029508497187473712051146708591409529430, \
261 1.30437298687487732290535130675991113734, -0.563661001875017525299235527605726980380, -0.563661001875017525299235527605726980380, \
262 0.377322003750035050598471055211453960760, 0.377322003750035050598471055211453960760, 0.186338998124982474700764472394273019620, \
263 -0.563661001875017525299235527605726980380, 1.30437298687487732290535130675991113734, -0.563661001875017525299235527605726980380, \
264 -0.563661001875017525299235527605726980380, -0.563661001875017525299235527605726980380, 1.30437298687487732290535130675991113734, \
265 0.377322003750035050598471055211453960760, 0.186338998124982474700764472394273019620, 0.377322003750035050598471055211453960760, \
266 0.186338998124982474700764472394273019620, 0.377322003750035050598471055211453960760, 0.377322003750035050598471055211453960760, \
269 0.029508497187473712051146708591409529430, 0.029508497187473712051146708591409529430, 0.029508497187473712051146708591409529430, \
270 1.088525491562421136153440125774228588290, 0.0, 0.0, \
271 0.0, -0.029508497187473712051146708591409529430, 0.0, \
272 0.0, 0.0, -0.029508497187473712051146708591409529430, \
273 -1.30437298687487732290535130675991113734, -1.868033988749894848204586834365638117720, -1.868033988749894848204586834365638117720, \
274 0.563661001875017525299235527605726980380, 1.868033988749894848204586834365638117720, 0.0, \
275 -0.377322003750035050598471055211453960760, 0.0, -0.190983005625052575897706582817180941140, \
276 -0.377322003750035050598471055211453960760, -0.190983005625052575897706582817180941140, 0.0, \
277 0.563661001875017525299235527605726980380, 0.0, 1.868033988749894848204586834365638117720, \
278 -0.186338998124982474700764472394273019620, 0.190983005625052575897706582817180941140, 0.19098300562505257589770658281718094114, \
281 0.029508497187473712051146708591409529430, 0.029508497187473712051146708591409529430, 0.029508497187473712051146708591409529430, \
282 -0.029508497187473712051146708591409529430, 0.0, 0.0, \
283 0.0, 1.088525491562421136153440125774228588290, 0.0, \
284 0.0, 0.0, -0.029508497187473712051146708591409529430, \
285 0.0, -0.377322003750035050598471055211453960760, -0.190983005625052575897706582817180941140, \
286 1.868033988749894848204586834365638117720, 0.563661001875017525299235527605726980380, 0.0, \
287 -1.868033988749894848204586834365638117720, -1.30437298687487732290535130675991113734, -1.868033988749894848204586834365638117720, \
288 -0.190983005625052575897706582817180941140, -0.377322003750035050598471055211453960760, 0.0, \
289 0.190983005625052575897706582817180941140, -0.186338998124982474700764472394273019620, 0.190983005625052575897706582817180941140, \
290 0.0, 0.563661001875017525299235527605726980380, 1.868033988749894848204586834365638117720, \
293 0.029508497187473712051146708591409529430, 0.029508497187473712051146708591409529430, 0.029508497187473712051146708591409529430, \
294 -0.029508497187473712051146708591409529430, 0.0, 0.0, \
295 0.0, -0.029508497187473712051146708591409529430, 0.0, \
296 0.0, 0.0, 1.088525491562421136153440125774228588290, \
297 0.0, -0.190983005625052575897706582817180941140, -0.377322003750035050598471055211453960760, \
298 0.190983005625052575897706582817180941140, 0.190983005625052575897706582817180941140, -0.186338998124982474700764472394273019620, \
299 -0.190983005625052575897706582817180941140, 0.0, -0.377322003750035050598471055211453960760, \
300 -1.868033988749894848204586834365638117720, -1.868033988749894848204586834365638117720, -1.30437298687487732290535130675991113734,
301 1.868033988749894848204586834365638117720, 0.0, 0.563661001875017525299235527605726980380, \
302 0.0, 1.868033988749894848204586834365638117720, 0.563661001875017525299235527605726980380, \
309 int numNodes = tetNodes.dimension(0);
316 *outStream <<
" check VALUE of basis functions at nodes\n";
317 vals.
resize(numFields, numNodes);
318 tetBasis.
getValues(vals, tetNodes, OPERATOR_VALUE);
320 for (
int i = 0; i < numFields; i++) {
321 for (
int j = 0; j < numNodes; j++) {
322 int l = i + j * numFields;
323 if (std::abs(vals(i,j) - nodalBasisValues[l]) > INTREPID_TOL) {
325 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
328 *outStream <<
" At multi-index { ";
329 *outStream << i <<
" ";*outStream << j <<
" ";
330 *outStream <<
"} computed value: " << vals(i,j)
331 <<
" but reference value: " << nodalBasisValues[l] <<
"\n";
337 *outStream <<
" check VALUE of basis functions at points\n";
339 vals.
resize(numFields, numPoints);
341 tetBasis.
getValues(vals, tetPoints, OPERATOR_VALUE);
343 for (
int i = 0; i < numFields; i++) {
344 for (
int j = 0; j < numPoints; j++) {
345 int l = i + j * numFields;
346 if (std::abs(vals(i,j) - pointBasisValues[l]) > INTREPID_TOL) {
348 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
351 *outStream <<
" At multi-index { ";
352 *outStream << i <<
" ";*outStream << j <<
" ";
353 *outStream <<
"} computed value: " << vals(i,j)
354 <<
" but reference value: " << pointBasisValues[l] <<
"\n";
360 *outStream <<
" check VALUE of basis functions at random points\n";
361 int numRandomPoints = 16384;
363 vals.
resize(numFields, numRandomPoints);
368 std::random_device rd;
369 std::mt19937 gen(rd());
370 std::uniform_real_distribution<> dis(0, 1);
371 while (point < numRandomPoints) {
377 if (r + s + t > 1.0)
continue;
379 tetRandomPoints(point, 0) = r;
380 tetRandomPoints(point, 1) = s;
381 tetRandomPoints(point, 2) = t;
386 tetBasis.
getValues(vals, tetRandomPoints, OPERATOR_VALUE);
388 for (
int j = 0; j < numRandomPoints; j++) {
390 for (
int i = 0; i < numFields; i++) {
393 if (std::abs(sum - 1.0) > INTREPID_TOL) {
395 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
398 *outStream <<
" Composite tet basis functions";
399 *outStream <<
" are not summing to 1.0\n";
400 *outStream <<
" sum : " << sum <<
"\n";
405 numPoints = tetPoints.dimension(0);
406 vals.
resize(numFields, numPoints, spaceDim);
408 tetBasis.
getValues(vals, tetPoints, OPERATOR_GRAD);
410 for (
int i = 0; i < numFields; i++) {
411 for (
int j = 0; j < numPoints; j++) {
412 for (
int k = 0; k < spaceDim; k++) {
413 int l = k + i * spaceDim + j * spaceDim * numFields;
414 if (std::abs(vals(i,j,k) - pointBasisGrads[l]) > INTREPID_TOL) {
416 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
419 *outStream <<
" At multi-index { ";
420 *outStream << i <<
" ";*outStream << j <<
" ";*outStream << k <<
" ";
421 *outStream <<
"} computed grad component: " << vals(i,j,k)
422 <<
" but reference grad component: " << pointBasisGrads[l] <<
"\n";
429 catch (
const std::logic_error & err) {
430 *outStream << err.what() <<
"\n\n";
436 std::cout <<
"End Result: TEST FAILED\n";
438 std::cout <<
"End Result: TEST PASSED\n";
441 std::cout.copyfmt(oldFormatState);
virtual int getCardinality() const
Returns cardinality of the basis.
int dimension(const int whichDim) const
Returns the specified dimension.
Header file for utility class to provide multidimensional containers.
void resize(const int dim0)
Resizes FieldContainer to a rank-1 container with the specified dimension, initialized by 0...
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
FEM basis evaluation on a reference Tetrahedron cell.
void initialize(const Scalar value=0)
Initializes a field container by assigning value to all its elements.
Header file for the Intrepid::HGRAD_TET_COMP12_FEM class.
virtual const shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation http://trilin...