51 #include "Teuchos_oblackholestream.hpp"
52 #include "Teuchos_RCP.hpp"
53 #include "Teuchos_GlobalMPISession.hpp"
56 #include "Shards_CellTopology.hpp"
59 using namespace Intrepid;
65 int main(
int argc,
char *argv[]) {
67 Teuchos::GlobalMPISession mpiSession(&argc, &argv);
70 int iprint = argc - 1;
72 Teuchos::RCP<std::ostream> outStream;
73 Teuchos::oblackholestream bhs;
76 outStream = Teuchos::rcp(&std::cout,
false);
78 outStream = Teuchos::rcp(&bhs,
false);
81 Teuchos::oblackholestream oldFormatState;
82 oldFormatState.copyfmt(std::cout);
85 <<
"===============================================================================\n" \
87 <<
"| Unit Test HDIV_TRI_In_FEM |\n" \
89 <<
"| 1) Tests triangular Raviart-Thomas basis |\n" \
91 <<
"| Questions? Contact Pavel Bochev (pbboche@sandia.gov) or |\n" \
92 <<
"| Denis Ridzal (dridzal@sandia.gov) or |\n" \
93 <<
"| Robert Kirby (robert.c.kirby@ttu.edu) |\n" \
95 <<
"| Intrepid's website: http://trilinos.sandia.gov/packages/intrepid |\n" \
96 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" \
98 <<
"===============================================================================\n";
111 PointTools::getLattice<double,FieldContainer<double> >( lattice ,
115 POINTTYPE_EQUISPACED );
117 myBasis.
getValues( myBasisValues , lattice , OPERATOR_VALUE );
119 const double fiat_vals[] = {
120 0.000000000000000e+00, -2.000000000000000e+00,
121 2.500000000000000e-01, -5.000000000000000e-01,
122 -1.000000000000000e+00, 1.000000000000000e+00,
123 0.000000000000000e+00, -2.500000000000000e-01,
124 -5.000000000000000e-01, 5.000000000000000e-01,
125 0.000000000000000e+00, 0.000000000000000e+00,
126 0.000000000000000e+00, 1.000000000000000e+00,
127 2.500000000000000e-01, -5.000000000000000e-01,
128 2.000000000000000e+00, -2.000000000000000e+00,
129 0.000000000000000e+00, 5.000000000000000e-01,
130 2.500000000000000e-01, -2.500000000000000e-01,
131 0.000000000000000e+00, 0.000000000000000e+00,
132 0.000000000000000e+00, 0.000000000000000e+00,
133 2.500000000000000e-01, 0.000000000000000e+00,
134 2.000000000000000e+00, 0.000000000000000e+00,
135 0.000000000000000e+00, -5.000000000000000e-01,
136 2.500000000000000e-01, 2.500000000000000e-01,
137 0.000000000000000e+00, -1.000000000000000e+00,
138 0.000000000000000e+00, 0.000000000000000e+00,
139 -5.000000000000000e-01, 0.000000000000000e+00,
140 -1.000000000000000e+00, 0.000000000000000e+00,
141 0.000000000000000e+00, 2.500000000000000e-01,
142 2.500000000000000e-01, 2.500000000000000e-01,
143 0.000000000000000e+00, 2.000000000000000e+00,
144 1.000000000000000e+00, 0.000000000000000e+00,
145 5.000000000000000e-01, 0.000000000000000e+00,
146 0.000000000000000e+00, 0.000000000000000e+00,
147 -5.000000000000000e-01, 2.500000000000000e-01,
148 -2.500000000000000e-01, 2.500000000000000e-01,
149 -2.000000000000000e+00, 2.000000000000000e+00,
150 -2.000000000000000e+00, 0.000000000000000e+00,
151 -2.500000000000000e-01, 0.000000000000000e+00,
152 0.000000000000000e+00, 0.000000000000000e+00,
153 -5.000000000000000e-01, 2.500000000000000e-01,
154 5.000000000000000e-01, -5.000000000000000e-01,
155 1.000000000000000e+00, -1.000000000000000e+00,
156 0.000000000000000e+00, 0.000000000000000e+00,
157 1.500000000000000e+00, 0.000000000000000e+00,
158 0.000000000000000e+00, 0.000000000000000e+00,
159 0.000000000000000e+00, 7.500000000000000e-01,
160 7.500000000000000e-01, -7.500000000000000e-01,
161 0.000000000000000e+00, 0.000000000000000e+00,
162 0.000000000000000e+00, 0.000000000000000e+00,
163 7.500000000000000e-01, 0.000000000000000e+00,
164 0.000000000000000e+00, 0.000000000000000e+00,
165 0.000000000000000e+00, 1.500000000000000e+00,
166 -7.500000000000000e-01, 7.500000000000000e-01,
167 0.000000000000000e+00, 0.000000000000000e+00
171 for (
int i=0;i<myBasisValues.dimension(0);i++) {
172 for (
int j=0;j<myBasisValues.dimension(1);j++) {
173 for (
int k=0;k<myBasisValues.dimension(2);k++) {
174 if (std::abs( myBasisValues(i,j,k) - fiat_vals[cur] ) > INTREPID_TOL ) {
176 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
179 *outStream <<
" At multi-index { ";
180 *outStream << i <<
" " << j <<
" " << k;
181 *outStream <<
"} computed value: " << myBasisValues(i,j,k)
182 <<
" but correct value: " << fiat_vals[cur] <<
"\n";
183 *outStream <<
"Difference: " << std::abs( myBasisValues(i,j,k) - fiat_vals[cur] ) <<
"\n";
190 catch (
const std::exception & err) {
191 *outStream << err.what() <<
"\n\n";
202 PointTools::getLattice<double,FieldContainer<double> >( lattice ,
206 POINTTYPE_EQUISPACED );
208 myBasis.
getValues( myBasisDivs , lattice , OPERATOR_DIV );
211 const double fiat_divs[] = {
212 7.000000000000000e+00,
213 2.500000000000000e+00,
214 -2.000000000000000e+00,
215 2.500000000000000e+00,
216 -2.000000000000000e+00,
217 -2.000000000000000e+00,
218 -2.000000000000000e+00,
219 2.500000000000000e+00,
220 7.000000000000000e+00,
221 -2.000000000000000e+00,
222 2.500000000000000e+00,
223 -2.000000000000000e+00,
224 -2.000000000000000e+00,
225 2.500000000000000e+00,
226 7.000000000000000e+00,
227 -2.000000000000000e+00,
228 2.500000000000000e+00,
229 -2.000000000000000e+00,
230 -2.000000000000000e+00,
231 -2.000000000000000e+00,
232 -2.000000000000000e+00,
233 2.500000000000000e+00,
234 2.500000000000000e+00,
235 7.000000000000000e+00,
236 -2.000000000000000e+00,
237 -2.000000000000000e+00,
238 -2.000000000000000e+00,
239 2.500000000000000e+00,
240 2.500000000000000e+00,
241 7.000000000000000e+00,
242 7.000000000000000e+00,
243 2.500000000000000e+00,
244 -2.000000000000000e+00,
245 2.500000000000000e+00,
246 -2.000000000000000e+00,
247 -2.000000000000000e+00,
248 9.000000000000000e+00,
249 0.000000000000000e+00,
250 -9.000000000000000e+00,
251 4.500000000000000e+00,
252 -4.500000000000000e+00,
253 0.000000000000000e+00,
254 9.000000000000000e+00,
255 4.500000000000000e+00,
256 0.000000000000000e+00,
257 0.000000000000000e+00,
258 -4.500000000000000e+00,
259 -9.000000000000000e+00
263 for (
int i=0;i<myBasisDivs.dimension(0);i++) {
264 for (
int j=0;j<myBasisDivs.dimension(1);j++) {
265 if (std::abs( myBasisDivs(i,j) - fiat_divs[cur] ) > INTREPID_TOL ) {
267 *outStream << std::setw(70) <<
"^^^^----FAILURE!" <<
"\n";
270 *outStream <<
" At multi-index { ";
271 *outStream << i <<
" " << j;
272 *outStream <<
"} computed value: " << myBasisDivs(i,j)
273 <<
" but correct value: " << fiat_divs[cur] <<
"\n";
274 *outStream <<
"Difference: " << std::abs( myBasisDivs(i,j) - fiat_divs[cur] ) <<
"\n";
280 catch (
const std::exception & err) {
281 *outStream << err.what() <<
"\n\n";
287 std::cout <<
"End Result: TEST FAILED\n";
289 std::cout <<
"End Result: TEST PASSED\n";
292 std::cout.copyfmt(oldFormatState);
virtual int getCardinality() const
Returns cardinality of the basis.
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
Evaluation of a FEM basis on a reference Triangle cell.
Header file for the Intrepid::HDIV_TRI_In_FEM class.
Header file for utility class to provide multidimensional containers.
Implementation of the default H(div)-compatible Raviart-Thomas basis of arbitrary degree on Triangle ...
virtual const shards::CellTopology getBaseCellTopology() const
Returns the base cell topology for which the basis is defined. See Shards documentation http://trilin...