← Back to branch summary

~chaffra/+junk/trilinos

« back to all changes in this revision

Viewing changes to packages/intrepid/test/Discretization/Basis/HCURL_WEDGE_I1_FEM/test_01.cpp

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme, Christophe Prud'homme, Johannes Ring
  • Date: 2009-12-13 12:53:22 UTC
  • mfrom: (5.1.2 sid)
  • Revision ID: james.westby@ubuntu.com-20091213125322-in0nrdjc55deqsw9
Tags: 10.0.3.dfsg-1
[Christophe Prud'homme]
* New upstream release

[Johannes Ring]
* debian/patches/libname.patch: Add prefix 'libtrilinos_' to all
  libraries. 
* debian/patches/soname.patch: Add soversion to libraries.
* debian/watch: Update download URL.
* debian/control:
  - Remove python-numeric from Build-Depends (virtual package).
  - Remove automake and autotools from Build-Depends and add cmake to
    reflect switch to CMake.
  - Add python-support to Build-Depends.
* debian/rules: 
  - Cleanup and updates for switch to CMake.

Show diffs side-by-side

added added

removed removed

Lines of Context:
packages/intrepid/test/Discretization/Basis/HCURL_WEDGE_I1_FEM/test_01.cpp
 
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
// This library is free software; you can redistribute it and/or modify
 
11
// it under the terms of the GNU Lesser General Public License as
 
12
// published by the Free Software Foundation; either version 2.1 of the
 
13
// License, or (at your option) any later version.
 
14
//
 
15
// This library is distributed in the hope that it will be useful, but
 
16
// WITHOUT ANY WARRANTY; without even the implied warranty of
 
17
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 
18
// Lesser General Public License for more details.
 
19
//
 
20
// You should have received a copy of the GNU Lesser General Public
 
21
// License along with this library; if not, write to the Free Software
 
22
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
 
23
// USA
 
24
// Questions? Contact Pavel Bochev  (pbboche@sandia.gov),
 
25
//                    Denis Ridzal  (dridzal@sandia.gov),
 
26
//                    Kara Peterson (kjpeter@sandia.gov).
 
27
//
 
28
// ************************************************************************
 
29
// @HEADER
 
30
 
 
31
/** \file   test_01.cpp
 
32
    \brief  Unit tests for the Intrepid::C_WEDGE_I1_FEM class.
 
33
    \author Created by P. Bochev, D. Ridzal, and K. Peterson.
 
34
*/
 
35
#include "Intrepid_FieldContainer.hpp"
 
36
#include "Intrepid_HCURL_WEDGE_I1_FEM.hpp"
 
37
#include "Teuchos_oblackholestream.hpp"
 
38
#include "Teuchos_RCP.hpp"
 
39
#include "Teuchos_GlobalMPISession.hpp"
 
40
 
 
41
using namespace std;
 
42
using namespace Intrepid;
 
43
 
 
44
#define INTREPID_TEST_COMMAND( S , throwCounter, nException )                                                              \
 
45
{                                                                                                                          \
 
46
  ++nException;                                                                                                            \
 
47
  try {                                                                                                                    \
 
48
    S ;                                                                                                                    \
 
49
  }                                                                                                                        \
 
50
  catch (std::logic_error err) {                                                                                           \
 
51
      ++throwCounter;                                                                                                      \
 
52
      *outStream << "Expected Error " << nException << " -------------------------------------------------------------\n"; \
 
53
      *outStream << err.what() << '\n';                                                                                    \
 
54
      *outStream << "-------------------------------------------------------------------------------" << "\n\n";           \
 
55
  };                                                                                                                       \
 
56
}
 
57
 
 
58
int main(int argc, char *argv[]) {
 
59
 
 
60
  Teuchos::GlobalMPISession mpiSession(&argc, &argv);
 
61
  
 
62
  // This little trick lets us print to std::cout only if
 
63
  // a (dummy) command-line argument is provided.
 
64
  int iprint     = argc - 1;
 
65
  Teuchos::RCP<std::ostream> outStream;
 
66
  Teuchos::oblackholestream bhs; // outputs nothing
 
67
  if (iprint > 0)
 
68
    outStream = Teuchos::rcp(&std::cout, false);
 
69
  else
 
70
    outStream = Teuchos::rcp(&bhs, false);
 
71
  
 
72
  // Save the format state of the original std::cout.
 
73
  Teuchos::oblackholestream oldFormatState;
 
74
  oldFormatState.copyfmt(std::cout);
 
75
  
 
76
  *outStream \
 
77
    << "===============================================================================\n" \
 
78
    << "|                                                                             |\n" \
 
79
    << "|                 Unit Test (Basis_HCURL_WEDGE_I1_FEM)                        |\n" \
 
80
    << "|                                                                             |\n" \
 
81
    << "|     1) Conversion of Dof tags into Dof ordinals and back                    |\n" \
 
82
    << "|     2) Basis values for VALUE and CURL operators                            |\n" \
 
83
    << "|                                                                             |\n" \
 
84
    << "|  Questions? Contact  Pavel Bochev  (pbboche@sandia.gov),                    |\n" \
 
85
    << "|                      Denis Ridzal  (dridzal@sandia.gov),                    |\n" \
 
86
    << "|                      Kara Peterson (kjpeter@sandia.gov).                    |\n" \
 
87
    << "|                                                                             |\n" \
 
88
    << "|  Intrepid's website: http://trilinos.sandia.gov/packages/intrepid           |\n" \
 
89
    << "|  Trilinos website:   http://trilinos.sandia.gov                             |\n" \
 
90
    << "|                                                                             |\n" \
 
91
    << "===============================================================================\n"\
 
92
    << "| TEST 1: Basis creation, exception testing                                   |\n"\
 
93
    << "===============================================================================\n";
 
94
  
 
95
  // Define basis and error flag
 
96
  Basis_HCURL_WEDGE_I1_FEM<double, FieldContainer<double> > wedgeBasis;
 
97
  int errorFlag = 0;
 
98
 
 
99
  // Initialize throw counter for exception testing
 
100
  int nException     = 0;
 
101
  int throwCounter   = 0;
 
102
 
 
103
  // Define array containing the 6 vertices of the reference WEDGE and 6 other points.
 
104
  FieldContainer<double> wedgeNodes(12, 3);
 
105
  wedgeNodes(0,0) =  0.0;  wedgeNodes(0,1) =  0.0;  wedgeNodes(0,2) = -1.0;
 
106
  wedgeNodes(1,0) =  1.0;  wedgeNodes(1,1) =  0.0;  wedgeNodes(1,2) = -1.0;
 
107
  wedgeNodes(2,0) =  0.0;  wedgeNodes(2,1) =  1.0;  wedgeNodes(2,2) = -1.0;
 
108
  wedgeNodes(3,0) =  0.0;  wedgeNodes(3,1) =  0.0;  wedgeNodes(3,2) =  1.0;
 
109
  wedgeNodes(4,0) =  1.0;  wedgeNodes(4,1) =  0.0;  wedgeNodes(4,2) =  1.0;
 
110
  wedgeNodes(5,0) =  0.0;  wedgeNodes(5,1) =  1.0;  wedgeNodes(5,2) =  1.0;
 
111
 
 
112
  wedgeNodes(6,0) =  0.25; wedgeNodes(6,1) =  0.5;  wedgeNodes(6,2) = -1.0;
 
113
  wedgeNodes(7,0) =  0.5;  wedgeNodes(7,1) =  0.25; wedgeNodes(7,2) =  0.0;
 
114
  wedgeNodes(8,0) =  0.25; wedgeNodes(8,1) =  0.25; wedgeNodes(8,2) =  1.0;
 
115
  wedgeNodes(9,0) =  0.25; wedgeNodes(9,1) =  0.0;  wedgeNodes(9,2) =  0.75;
 
116
  wedgeNodes(10,0)=  0.0;  wedgeNodes(10,1)=  0.5;  wedgeNodes(10,2)= -0.25;
 
117
  wedgeNodes(11,0)=  0.5;  wedgeNodes(11,1)=  0.5;  wedgeNodes(11,2)=  0.0;
 
118
 
 
119
 
 
120
 
 
121
  // Generic array for the output values; needs to be properly resized depending on the operator type
 
122
  FieldContainer<double> vals;
 
123
 
 
124
  try{
 
125
    // exception #1: GRAD cannot be applied to HCURL functions 
 
126
    // resize vals to rank-3 container with dimensions (num. basis functions, num. points, arbitrary)
 
127
    vals.resize(wedgeBasis.getCardinality(), wedgeNodes.dimension(0), 3 );
 
128
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(vals, wedgeNodes, OPERATOR_GRAD), throwCounter, nException );
 
129
 
 
130
    // exception #2: DIV cannot be applied to HCURL functions
 
131
    // resize vals to rank-2 container with dimensions (num. basis functions, num. points)
 
132
    vals.resize(wedgeBasis.getCardinality(), wedgeNodes.dimension(0));
 
133
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(vals, wedgeNodes, OPERATOR_DIV), throwCounter, nException );
 
134
        
 
135
    // Exceptions 3-7: all bf tags/bf Ids below are wrong and should cause getDofOrdinal() and 
 
136
    // getDofTag() to access invalid array elements thereby causing bounds check exception
 
137
    // exception #3
 
138
    INTREPID_TEST_COMMAND( wedgeBasis.getDofOrdinal(3,0,0), throwCounter, nException );
 
139
    // exception #4
 
140
    INTREPID_TEST_COMMAND( wedgeBasis.getDofOrdinal(1,1,1), throwCounter, nException );
 
141
    // exception #5
 
142
    INTREPID_TEST_COMMAND( wedgeBasis.getDofOrdinal(0,4,1), throwCounter, nException );
 
143
    // exception #6
 
144
    INTREPID_TEST_COMMAND( wedgeBasis.getDofTag(10), throwCounter, nException );
 
145
    // exception #7
 
146
    INTREPID_TEST_COMMAND( wedgeBasis.getDofTag(-1), throwCounter, nException );
 
147
    
 
148
#ifdef HAVE_INTREPID_DEBUG
 
149
    // Exceptions 8- test exception handling with incorrectly dimensioned input/output arrays
 
150
    // exception #8: input points array must be of rank-2
 
151
    FieldContainer<double> badPoints1(4, 5, 3);
 
152
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
 
153
    
 
154
    // exception #9 dimension 1 in the input point array must equal space dimension of the cell
 
155
    FieldContainer<double> badPoints2(4, 2);
 
156
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(vals, badPoints1, OPERATOR_VALUE), throwCounter, nException );
 
157
    
 
158
    // exception #10 output values must be of rank-3 for OPERATOR_VALUE
 
159
    FieldContainer<double> badVals1(4, 3);
 
160
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(badVals1, wedgeNodes, OPERATOR_VALUE), throwCounter, nException );
 
161
 
 
162
    // exception #11 output values must be of rank-3 for OPERATOR_CURL
 
163
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(badVals1, wedgeNodes, OPERATOR_CURL), throwCounter, nException );
 
164
    
 
165
    // exception #12 incorrect 0th dimension of output array (must equal number of basis functions)
 
166
    FieldContainer<double> badVals2(wedgeBasis.getCardinality() + 1, wedgeNodes.dimension(0), 3);
 
167
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(badVals2, wedgeNodes, OPERATOR_VALUE), throwCounter, nException );
 
168
    
 
169
    // exception #13 incorrect 1st dimension of output array (must equal number of points)
 
170
    FieldContainer<double> badVals3(wedgeBasis.getCardinality(), wedgeNodes.dimension(0) + 1, 3);
 
171
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(badVals3, wedgeNodes, OPERATOR_VALUE), throwCounter, nException );
 
172
 
 
173
    // exception #14: incorrect 2nd dimension of output array (must equal the space dimension)
 
174
    FieldContainer<double> badVals4(wedgeBasis.getCardinality(), wedgeNodes.dimension(0), 4);
 
175
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(badVals4, wedgeNodes, OPERATOR_VALUE), throwCounter, nException );
 
176
    
 
177
    // exception #15: incorrect 2nd dimension of output array (must equal the space dimension)
 
178
    INTREPID_TEST_COMMAND( wedgeBasis.getValues(badVals4, wedgeNodes, OPERATOR_CURL), throwCounter, nException );
 
179
#endif
 
180
    
 
181
  }
 
182
  catch (std::logic_error err) {
 
183
    *outStream << "UNEXPECTED ERROR !!! ----------------------------------------------------------\n";
 
184
    *outStream << err.what() << '\n';
 
185
    *outStream << "-------------------------------------------------------------------------------" << "\n\n";
 
186
    errorFlag = -1000;
 
187
  };
 
188
  
 
189
  // Check if number of thrown exceptions matches the one we expect 
 
190
  // Note Teuchos throw number will not pick up exceptions 3-7 and therefore will not match.
 
191
  if (throwCounter != nException) {
 
192
    errorFlag++;
 
193
    *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
 
194
  }
 
195
  
 
196
  *outStream \
 
197
    << "\n"
 
198
    << "===============================================================================\n"\
 
199
    << "| TEST 2: correctness of tag to enum and enum to tag lookups                  |\n"\
 
200
    << "===============================================================================\n";
 
201
  
 
202
  try{
 
203
    std::vector<std::vector<int> > allTags = wedgeBasis.getAllDofTags();
 
204
    
 
205
    // Loop over all tags, lookup the associated dof enumeration and then lookup the tag again
 
206
    for (unsigned i = 0; i < allTags.size(); i++) {
 
207
      int bfOrd  = wedgeBasis.getDofOrdinal(allTags[i][0], allTags[i][1], allTags[i][2]);
 
208
      
 
209
      std::vector<int> myTag = wedgeBasis.getDofTag(bfOrd);
 
210
       if( !( (myTag[0] == allTags[i][0]) &&
 
211
              (myTag[1] == allTags[i][1]) &&
 
212
              (myTag[2] == allTags[i][2]) &&
 
213
              (myTag[3] == allTags[i][3]) ) ) {
 
214
        errorFlag++;
 
215
        *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
 
216
        *outStream << " getDofOrdinal( {" 
 
217
          << allTags[i][0] << ", " 
 
218
          << allTags[i][1] << ", " 
 
219
          << allTags[i][2] << ", " 
 
220
          << allTags[i][3] << "}) = " << bfOrd <<" but \n";   
 
221
        *outStream << " getDofTag(" << bfOrd << ") = { "
 
222
          << myTag[0] << ", " 
 
223
          << myTag[1] << ", " 
 
224
          << myTag[2] << ", " 
 
225
          << myTag[3] << "}\n";        
 
226
      }
 
227
    }
 
228
    
 
229
    // Now do the same but loop over basis functions
 
230
    for( int bfOrd = 0; bfOrd < wedgeBasis.getCardinality(); bfOrd++) {
 
231
      std::vector<int> myTag  = wedgeBasis.getDofTag(bfOrd);
 
232
      int myBfOrd = wedgeBasis.getDofOrdinal(myTag[0], myTag[1], myTag[2]);
 
233
      if( bfOrd != myBfOrd) {
 
234
        errorFlag++;
 
235
        *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
 
236
        *outStream << " getDofTag(" << bfOrd << ") = { "
 
237
          << myTag[0] << ", " 
 
238
          << myTag[1] << ", " 
 
239
          << myTag[2] << ", " 
 
240
          << myTag[3] << "} but getDofOrdinal({" 
 
241
          << myTag[0] << ", " 
 
242
          << myTag[1] << ", " 
 
243
          << myTag[2] << ", " 
 
244
          << myTag[3] << "} ) = " << myBfOrd << "\n";
 
245
      }
 
246
    }
 
247
  }
 
248
  catch (std::logic_error err){
 
249
    *outStream << err.what() << "\n\n";
 
250
    errorFlag = -1000;
 
251
  };
 
252
  
 
253
  *outStream \
 
254
    << "\n"
 
255
    << "===============================================================================\n"\
 
256
    << "| TEST 3: correctness of basis function values                                |\n"\
 
257
    << "===============================================================================\n";
 
258
  
 
259
  outStream -> precision(20);
 
260
  
 
261
  // VALUE: Each row pair gives the 9x3 correct basis set values at an evaluation point: (P,F,D) layout
 
262
  double basisValues[] = {
 
263
    1.00000, 0, 0, 0, 0, 0, 0, -1.00000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
 
264
    0, 0.500000, 0, 0, 0, 0, 0, 0, 1.00000, 1.00000, 0, 0, 1.00000, 0, 0, \
 
265
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.500000, 0, 0, 0, 0, \
 
266
    0, 0, -1.00000, 0, 0, -1.00000, -1.00000, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
 
267
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0.500000, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
 
268
    1.00000, 0, 0, 0, 0, 0, 0, -1.00000, 0, 0, 0, 0.500000, 0, 0, 0, 0, \
 
269
    0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1.00000, 1.00000, 0, 0, 1.00000, 0, \
 
270
    0, 0, 0, 0, 0, 0, 0, 0, 0.500000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, \
 
271
    0, 0, 0, -1.00000, 0, 0, -1.00000, -1.00000, 0, 0, 0, 0, 0, 0, 0, 0, \
 
272
    0, 0.500000, 0.500000, 0.250000, 0, -0.500000, 0.250000, 0, \
 
273
    -0.500000, -0.750000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.125000, \
 
274
    0, 0, 0.125000, 0, 0, 0.250000, 0.375000, 0.250000, 0, -0.125000, \
 
275
    0.250000, 0, -0.125000, -0.250000, 0, 0.375000, 0.250000, 0, \
 
276
    -0.125000, 0.250000, 0, -0.125000, -0.250000, 0, 0, 0, 0.125000, 0, \
 
277
    0, 0.250000, 0, 0, 0.125000, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.750000, \
 
278
    0.250000, 0, -0.250000, 0.250000, 0, -0.250000, -0.750000, 0, 0, 0, \
 
279
    0.250000, 0, 0, 0.125000, 0, 0, 0.125000, 0.125000, 0.0312500, 0, 0, \
 
280
    0.0312500, 0, 0, -0.0937500, 0, 0.875000, 0.218750, 0, 0, 0.218750, \
 
281
    0, 0, -0.656250, 0, 0, 0, 0.375000, 0, 0, 0.125000, 0, 0, 0, \
 
282
    0.312500, 0, 0, -0.312500, 0, 0, -0.312500, -0.625000, 0, 0.187500, \
 
283
    0, 0, -0.187500, 0, 0, -0.187500, -0.375000, 0, 0, 0, 0.250000, 0, 0, \
 
284
    0, 0, 0, 0.250000, 0.250000, 0.250000, 0, -0.250000, 0.250000, 0, \
 
285
    -0.250000, -0.250000, 0, 0.250000, 0.250000, 0, -0.250000, 0.250000, \
 
286
    0, -0.250000, -0.250000, 0, 0, 0, 0, 0, 0, 0.250000, 0, 0, 0.250000    
 
287
  };
 
288
  
 
289
  // CURL: each row pair gives the 9x3 correct values of the curls of the 9 basis functions: (P,F,D) layout
 
290
  double basisCurls[] = {  
 
291
    0, -0.500000, 2.00000, 0, 0, 2.00000, -0.500000, 0, 2.00000, 0, \
 
292
    0.500000, 0, 0, 0, 0, 0.500000, 0, 0, -0.500000, 0.500000, 0, 0, \
 
293
    -0.500000, 0, 0.500000, 0, 0, 0.500000, -0.500000, 2.00000, 0.500000, \
 
294
    0, 2.00000, 0, 0, 2.00000, -0.500000, 0.500000, 0, -0.500000, 0, 0, \
 
295
    0, 0, 0, -0.500000, 0.500000, 0, 0, -0.500000, 0, 0.500000, 0, 0, 0, \
 
296
    0, 2.00000, 0, 0.500000, 2.00000, -0.500000, 0.500000, 2.00000, 0, 0, \
 
297
    0, 0, -0.500000, 0, 0.500000, -0.500000, 0, -0.500000, 0.500000, 0, \
 
298
    0, -0.500000, 0, 0.500000, 0, 0, 0, -0.500000, 0, 0, 0, 0, -0.500000, \
 
299
    0, 0, 0, 0.500000, 2.00000, 0, 0, 2.00000, 0.500000, 0, 2.00000, \
 
300
    -0.500000, 0.500000, 0, 0, -0.500000, 0, 0.500000, 0, 0, 0.500000, \
 
301
    -0.500000, 0, 0.500000, 0, 0, 0, 0, 0, -0.500000, 0.500000, 2.00000, \
 
302
    -0.500000, 0, 2.00000, 0, 0, 2.00000, -0.500000, 0.500000, 0, 0, \
 
303
    -0.500000, 0, 0.500000, 0, 0, 0, 0, 0, 0, 0.500000, 0, -0.500000, \
 
304
    0.500000, 0, 0, 0, 2.00000, 0, -0.500000, 2.00000, 0.500000, \
 
305
    -0.500000, 2.00000, -0.500000, 0.500000, 0, 0, -0.500000, 0, \
 
306
    0.500000, 0, 0, 0.125000, -0.250000, 2.00000, 0.125000, 0.250000, \
 
307
    2.00000, -0.375000, 0.250000, 2.00000, -0.125000, 0.250000, 0, \
 
308
    -0.125000, -0.250000, 0, 0.375000, -0.250000, 0, -0.500000, 0.500000, \
 
309
    0, 0, -0.500000, 0, 0.500000, 0, 0, 0.250000, -0.375000, 1.00000, \
 
310
    0.250000, 0.125000, 1.00000, -0.250000, 0.125000, 1.00000, -0.250000, \
 
311
    0.375000, 1.00000, -0.250000, -0.125000, 1.00000, 0.250000, \
 
312
    -0.125000, 1.00000, -0.500000, 0.500000, 0, 0, -0.500000, 0, \
 
313
    0.500000, 0, 0, 0.125000, -0.375000, 0, 0.125000, 0.125000, 0, \
 
314
    -0.375000, 0.125000, 0, -0.125000, 0.375000, 2.00000, -0.125000, \
 
315
    -0.125000, 2.00000, 0.375000, -0.125000, 2.00000, -0.500000, \
 
316
    0.500000, 0, 0, -0.500000, 0, 0.500000, 0, 0, 0.125000, -0.500000, \
 
317
    0.250000, 0.125000, 0, 0.250000, -0.375000, 0, 0.250000, -0.125000, \
 
318
    0.500000, 1.75000, -0.125000, 0, 1.75000, 0.375000, 0, 1.75000, \
 
319
    -0.500000, 0.500000, 0, 0, -0.500000, 0, 0.500000, 0, 0, 0, \
 
320
    -0.250000, 1.25000, 0, 0.250000, 1.25000, -0.500000, 0.250000, \
 
321
    1.25000, 0, 0.250000, 0.750000, 0, -0.250000, 0.750000, 0.500000, \
 
322
    -0.250000, 0.750000, -0.500000, 0.500000, 0, 0, -0.500000, 0, \
 
323
    0.500000, 0, 0, 0.250000, -0.250000, 1.00000, 0.250000, 0.250000, \
 
324
    1.00000, -0.250000, 0.250000, 1.00000, -0.250000, 0.250000, 1.00000, \
 
325
    -0.250000, -0.250000, 1.00000, 0.250000, -0.250000, 1.00000, \
 
326
    -0.500000, 0.500000, 0, 0, -0.500000, 0, 0.500000, 0, 0
 
327
  };
 
328
  
 
329
  try{
 
330
        
 
331
    // Dimensions for the output arrays:
 
332
    int numFields = wedgeBasis.getCardinality();
 
333
    int numPoints = wedgeNodes.dimension(0);
 
334
    int spaceDim  = wedgeBasis.getBaseCellTopology().getDimension();
 
335
    
 
336
    // Generic array for values and curls that will be properly sized before each call
 
337
    FieldContainer<double> vals;
 
338
    
 
339
    // Check VALUE of basis functions: resize vals to rank-3 container:
 
340
    vals.resize(numFields, numPoints, spaceDim);
 
341
    wedgeBasis.getValues(vals, wedgeNodes, OPERATOR_VALUE);
 
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
          
 
346
          // compute offset for (P,F,D) data layout: indices are P->j, F->i, D->k
 
347
           int l = k + i * spaceDim + j * spaceDim * numFields;
 
348
           if (std::abs(vals(i,j,k) - basisValues[l]) > INTREPID_TOL) {
 
349
             errorFlag++;
 
350
             *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
 
351
 
 
352
             // Output the multi-index of the value where the error is:
 
353
             *outStream << " At multi-index { ";
 
354
             *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
 
355
             *outStream << "}  computed value: " << vals(i,j,k)
 
356
               << " but reference value: " << basisValues[l] << "\n";
 
357
            }
 
358
         }
 
359
      }
 
360
    }
 
361
 
 
362
    // Check CURL of basis function: resize vals to rank-3 container
 
363
    vals.resize(numFields, numPoints, spaceDim);
 
364
    wedgeBasis.getValues(vals, wedgeNodes, OPERATOR_CURL);
 
365
    for (int i = 0; i < numFields; i++) {
 
366
      for (int j = 0; j < numPoints; j++) {
 
367
        for (int k = 0; k < spaceDim; k++) {
 
368
          
 
369
          // compute offset for (P,F,D) data layout: indices are P->j, F->i, D->k
 
370
           int l = k + i * spaceDim + j * spaceDim * numFields;
 
371
           if (std::abs(vals(i,j,k) - basisCurls[l]) > INTREPID_TOL) {
 
372
             errorFlag++;
 
373
             *outStream << std::setw(70) << "^^^^----FAILURE!" << "\n";
 
374
 
 
375
             // Output the multi-index of the value where the error is:
 
376
             *outStream << " At multi-index { ";
 
377
             *outStream << i << " ";*outStream << j << " ";*outStream << k << " ";
 
378
             *outStream << "}  computed curl component: " << vals(i,j,k)
 
379
               << " but reference curl component: " << basisCurls[l] << "\n";
 
380
            }
 
381
         }
 
382
      }
 
383
    }
 
384
    
 
385
   }    
 
386
  
 
387
  // Catch unexpected errors
 
388
  catch (std::logic_error err) {
 
389
    *outStream << err.what() << "\n\n";
 
390
    errorFlag = -1000;
 
391
  };
 
392
  
 
393
  if (errorFlag != 0)
 
394
    std::cout << "End Result: TEST FAILED\n";
 
395
  else
 
396
    std::cout << "End Result: TEST PASSED\n";
 
397
  
 
398
  // reset format state of std::cout
 
399
  std::cout.copyfmt(oldFormatState);
 
400
  
 
401
  return errorFlag;
 
402
}