51 template<
class Scalar,
class ArrayScalar>
53 const ArrayScalar & ptsClosed ,
54 const ArrayScalar & ptsOpen):
55 closedBasis_( order , ptsClosed ),
56 openBasis_( order-1 , ptsOpen ),
57 closedPts_( ptsClosed ),
61 this ->
basisCardinality_ = 3 * closedBasis_.getCardinality() * openBasis_.getCardinality() * openBasis_.getCardinality();
62 this ->
basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
67 Array<Array<RCP<Basis<Scalar,ArrayScalar > > > > bases(3);
68 bases[0].resize(3); bases[1].resize(3); bases[2].resize(3);
69 bases[0][0] = rcp( &closedBasis_ ,
false );
70 bases[0][1] = rcp( &openBasis_ ,
false );
71 bases[0][2] = rcp( &openBasis_ ,
false );
72 bases[1][0] = rcp( &openBasis_ ,
false );
73 bases[1][1] = rcp( &closedBasis_ ,
false );
74 bases[1][2] = rcp( &openBasis_ ,
false );
75 bases[2][0] = rcp( &openBasis_ ,
false );
76 bases[2][1] = rcp( &openBasis_ ,
false );
77 bases[2][2] = rcp( &closedBasis_ ,
false );
78 this->setBases( bases );
82 template<
class Scalar,
class ArrayScalar>
84 closedBasis_( order , pointType==POINTTYPE_SPECTRAL?POINTTYPE_SPECTRAL:POINTTYPE_EQUISPACED ),
85 openBasis_( order-1 , pointType==POINTTYPE_SPECTRAL?POINTTYPE_SPECTRAL_OPEN:POINTTYPE_EQUISPACED ),
86 closedPts_( order+1 , 1 ),
90 this ->
basisCardinality_ = 3 * closedBasis_.getCardinality() * openBasis_.getCardinality() * openBasis_.getCardinality();
91 this ->
basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
96 PointTools::getLattice<Scalar,FieldContainer<Scalar> >( closedPts_ ,
97 shards::CellTopology(shards::getCellTopologyData<shards::Line<2> >()) ,
100 pointType==POINTTYPE_SPECTRAL?POINTTYPE_WARPBLEND:POINTTYPE_EQUISPACED );
102 if (pointType == POINTTYPE_SPECTRAL)
104 PointTools::getGaussPoints<Scalar,FieldContainer<Scalar> >( openPts_ ,
109 PointTools::getLattice<Scalar,FieldContainer<Scalar> >( openPts_ ,
110 shards::CellTopology(shards::getCellTopologyData<shards::Line<2> >()) ,
113 POINTTYPE_EQUISPACED );
118 Array<Array<RCP<Basis<Scalar,ArrayScalar > > > > bases(3);
119 bases[0].resize(3); bases[1].resize(3); bases[2].resize(3);
120 bases[0][0] = rcp( &closedBasis_ ,
false );
121 bases[0][1] = rcp( &openBasis_ ,
false );
122 bases[0][2] = rcp( &openBasis_ ,
false );
123 bases[1][0] = rcp( &openBasis_ ,
false );
124 bases[1][1] = rcp( &closedBasis_ ,
false );
125 bases[1][2] = rcp( &openBasis_ ,
false );
126 bases[2][0] = rcp( &openBasis_ ,
false );
127 bases[2][1] = rcp( &openBasis_ ,
false );
128 bases[2][2] = rcp( &closedBasis_ ,
false );
129 this->setBases( bases );
135 template<
class Scalar,
class ArrayScalar>
144 std::vector<int> tags( tagSize * this->getCardinality() );
146 const std::vector<std::vector<int> >& closedDofTags = closedBasis_.getAllDofTags();
147 const std::vector<std::vector<int> >& openDofTags = openBasis_.getAllDofTags();
149 std::map<int,std::map<int,int> > total_dof_per_entity;
150 std::map<int,std::map<int,int> > current_dof_per_entity;
153 for (
int i=0;i<4;i++) {
154 total_dof_per_entity[0][i] = 0;
155 current_dof_per_entity[0][i] = 0;
158 for (
int i=0;i<12;i++) {
159 total_dof_per_entity[1][i] = 0;
160 current_dof_per_entity[1][i] = 0;
163 for (
int i=0;i<6;i++) {
164 total_dof_per_entity[2][i] = 0;
165 current_dof_per_entity[2][i] = 0;
167 total_dof_per_entity[3][0] = 0;
168 current_dof_per_entity[3][0] = 0;
171 for (
int i=0;i<6;i++) {
172 total_dof_per_entity[2][i] = openBasis_.getCardinality() * openBasis_.getCardinality();
175 total_dof_per_entity[2][0] = this->getCardinality() - 6 * openBasis_.getCardinality()* openBasis_.getCardinality();
180 for (
int k=0;k<openBasis_.getCardinality();k++) {
181 const int odimk = openDofTags[k][0];
182 const int oentk = openDofTags[k][1];
183 for (
int j=0;j<openBasis_.getCardinality();j++) {
184 const int odimj = openDofTags[j][0];
185 const int oentj = openDofTags[j][1];
186 for (
int i=0;i<closedBasis_.getCardinality();i++) {
187 const int cdim = closedDofTags[i][0];
188 const int cent = closedDofTags[i][1];
192 tags[4*tagcur] = dofdim;
193 tags[4*tagcur+1] = dofent;
194 tags[4*tagcur+2] = current_dof_per_entity[dofdim][dofent];
195 current_dof_per_entity[dofdim][dofent]++;
196 tags[4*tagcur+3] = total_dof_per_entity[dofdim][dofent];
204 for (
int k=0;k<openBasis_.getCardinality();k++) {
205 const int odimk = openDofTags[k][0];
206 const int oentk = openDofTags[k][1];
207 for (
int j=0;j<closedBasis_.getCardinality();j++) {
208 const int cdim = closedDofTags[j][0];
209 const int cent = closedDofTags[j][1];
210 for (
int i=0;i<openBasis_.getCardinality();i++) {
211 const int odimi = openDofTags[i][0];
212 const int oenti = openDofTags[i][1];
216 tags[4*tagcur] = dofdim;
217 tags[4*tagcur+1] = dofent;
218 tags[4*tagcur+2] = current_dof_per_entity[dofdim][dofent];
219 current_dof_per_entity[dofdim][dofent]++;
220 tags[4*tagcur+3] = total_dof_per_entity[dofdim][dofent];
228 for (
int k=0;k<closedBasis_.getCardinality();k++) {
229 const int cdim = closedDofTags[k][0];
230 const int cent = closedDofTags[k][1];
231 for (
int j=0;j<openBasis_.getCardinality();j++) {
232 const int odimj = openDofTags[j][0];
233 const int oentj = openDofTags[j][1];
234 for (
int i=0;i<openBasis_.getCardinality();i++) {
235 const int odimi = openDofTags[i][0];
236 const int oenti = openDofTags[i][1];
240 tags[4*tagcur] = dofdim;
241 tags[4*tagcur+1] = dofent;
242 tags[4*tagcur+2] = current_dof_per_entity[dofdim][dofent];
243 current_dof_per_entity[dofdim][dofent]++;
244 tags[4*tagcur+3] = total_dof_per_entity[dofdim][dofent];
259 this -> ordinalToTag_,
261 this -> basisCardinality_,
269 template<
class Scalar,
class ArrayScalar>
271 const ArrayScalar & inputPoints,
275 #ifdef HAVE_INTREPID_DEBUG 276 Intrepid::getValues_HDIV_Args<Scalar, ArrayScalar>(outputValues,
279 this -> getBaseCellTopology(),
280 this -> getCardinality() );
284 int dim0 = inputPoints.dimension(0);
292 for (
int i=0;i<dim0;i++) {
293 xPoints(i,0) = inputPoints(i,0);
294 yPoints(i,0) = inputPoints(i,1);
295 zPoints(i,0) = inputPoints(i,2);
298 switch (operatorType) {
308 closedBasis_.getValues( closedBasisValsXPts , xPoints , OPERATOR_VALUE );
309 closedBasis_.getValues( closedBasisValsYPts , yPoints , OPERATOR_VALUE );
310 closedBasis_.getValues( closedBasisValsZPts , zPoints , OPERATOR_VALUE );
311 openBasis_.getValues( openBasisValsXPts , xPoints , OPERATOR_VALUE );
312 openBasis_.getValues( openBasisValsYPts , yPoints , OPERATOR_VALUE );
313 openBasis_.getValues( openBasisValsZPts , zPoints , OPERATOR_VALUE );
317 for (
int k=0;k<openBasis_.getCardinality();k++) {
318 for (
int j=0;j<openBasis_.getCardinality();j++) {
319 for (
int i=0;i<closedBasis_.getCardinality();i++) {
320 for (
int l=0;l<dim0;l++) {
321 outputValues(bfcur,l,0) = closedBasisValsXPts(i,l) * openBasisValsYPts(j,l) * openBasisValsZPts(k,l);
322 outputValues(bfcur,l,1) = 0.0;
323 outputValues(bfcur,l,2) = 0.0;
331 for (
int k=0;k<openBasis_.getCardinality();k++) {
332 for (
int j=0;j<closedBasis_.getCardinality();j++) {
333 for (
int i=0;i<openBasis_.getCardinality();i++) {
334 for (
int l=0;l<dim0;l++) {
335 outputValues(bfcur,l,0) = 0.0;
336 outputValues(bfcur,l,1) = openBasisValsXPts(i,l) * closedBasisValsYPts(j,l) * openBasisValsZPts(k,l);
337 outputValues(bfcur,l,2) = 0.0;
345 for (
int k=0;k<closedBasis_.getCardinality();k++) {
346 for (
int j=0;j<openBasis_.getCardinality();j++) {
347 for (
int i=0;i<openBasis_.getCardinality();i++) {
348 for (
int l=0;l<dim0;l++) {
349 outputValues(bfcur,l,0) = 0.0;
350 outputValues(bfcur,l,1) = 0.0;
351 outputValues(bfcur,l,2) = openBasisValsXPts(i,l) * openBasisValsYPts(j,l) * closedBasisValsZPts(k,l);
370 closedBasis_.getValues( closedBasisDerivsXPts , xPoints , OPERATOR_D1 );
371 closedBasis_.getValues( closedBasisDerivsYPts , yPoints , OPERATOR_D1 );
372 closedBasis_.getValues( closedBasisDerivsZPts , zPoints , OPERATOR_D1 );
373 openBasis_.getValues( openBasisValsXPts , xPoints , OPERATOR_VALUE );
374 openBasis_.getValues( openBasisValsYPts , yPoints , OPERATOR_VALUE );
375 openBasis_.getValues( openBasisValsZPts , zPoints , OPERATOR_VALUE );
380 for (
int k=0;k<openBasis_.getCardinality();k++) {
381 for (
int j=0;j<openBasis_.getCardinality();j++) {
382 for (
int i=0;i<closedBasis_.getCardinality();i++) {
383 for (
int l=0;l<dim0;l++) {
384 outputValues(bfcur,l) = closedBasisDerivsXPts(i,l,0) * openBasisValsYPts(j,l) * openBasisValsZPts(k,l);
392 for (
int k=0;k<openBasis_.getCardinality();k++) {
393 for (
int j=0;j<closedBasis_.getCardinality();j++) {
394 for (
int i=0;i<openBasis_.getCardinality();i++) {
395 for (
int l=0;l<dim0;l++) {
396 outputValues(bfcur,l) = openBasisValsXPts(i,l) * closedBasisDerivsYPts(j,l,0) * openBasisValsZPts(k,l);
404 for (
int k=0;k<closedBasis_.getCardinality();k++) {
405 for (
int j=0;j<openBasis_.getCardinality();j++) {
406 for (
int i=0;i<openBasis_.getCardinality();i++) {
407 for (
int l=0;l<dim0;l++) {
408 outputValues(bfcur,l) = openBasisValsXPts(i,l) * openBasisValsYPts(j,l) * closedBasisDerivsZPts(k,l,0);
417 TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
418 ">>> ERROR (Basis_HDIV_HEX_In_FEM): CURL is invalid operator for HDIV Basis Functions");
422 TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_GRAD), std::invalid_argument,
423 ">>> ERROR (Basis_HDIV_HEX_In_FEM): GRAD is invalid operator for HDIV Basis Functions");
436 TEUCHOS_TEST_FOR_EXCEPTION( ( (operatorType == OPERATOR_D1) ||
437 (operatorType == OPERATOR_D2) ||
438 (operatorType == OPERATOR_D3) ||
439 (operatorType == OPERATOR_D4) ||
440 (operatorType == OPERATOR_D5) ||
441 (operatorType == OPERATOR_D6) ||
442 (operatorType == OPERATOR_D7) ||
443 (operatorType == OPERATOR_D8) ||
444 (operatorType == OPERATOR_D9) ||
445 (operatorType == OPERATOR_D10) ),
446 std::invalid_argument,
447 ">>> ERROR (Basis_HDIV_HEX_In_FEM): Invalid operator type");
451 TEUCHOS_TEST_FOR_EXCEPTION( ( (operatorType != OPERATOR_VALUE) &&
452 (operatorType != OPERATOR_GRAD) &&
453 (operatorType != OPERATOR_CURL) &&
454 (operatorType != OPERATOR_DIV) &&
455 (operatorType != OPERATOR_D1) &&
456 (operatorType != OPERATOR_D2) &&
457 (operatorType != OPERATOR_D3) &&
458 (operatorType != OPERATOR_D4) &&
459 (operatorType != OPERATOR_D5) &&
460 (operatorType != OPERATOR_D6) &&
461 (operatorType != OPERATOR_D7) &&
462 (operatorType != OPERATOR_D8) &&
463 (operatorType != OPERATOR_D9) &&
464 (operatorType != OPERATOR_D10) ),
465 std::invalid_argument,
466 ">>> ERROR (Basis_HDIV_HEX_In_FEM): Invalid operator type");
472 template<
class Scalar,
class ArrayScalar>
474 const ArrayScalar & inputPoints,
475 const ArrayScalar & cellVertices,
477 TEUCHOS_TEST_FOR_EXCEPTION( (
true), std::logic_error,
478 ">>> ERROR (Basis_HDIV_HEX_In_FEM): FEM Basis calling an FVD member function");
481 template<
class Scalar,
class ArrayScalar>
487 for (
int k=0;k<openPts_.dimension(0);k++)
489 for (
int j=0;j<openPts_.dimension(0);j++)
491 for (
int i=0;i<closedPts_.dimension(0);i++)
493 DofCoords(cur,0) = closedPts_(i,0);
494 DofCoords(cur,1) = openPts_(j,0);
495 DofCoords(cur,2) = openPts_(k,0);
502 for (
int k=0;k<openPts_.dimension(0);k++)
504 for (
int j=0;j<closedPts_.dimension(0);j++)
506 for (
int i=0;i<openPts_.dimension(0);i++)
508 DofCoords(cur,0) = openPts_(i,0);
509 DofCoords(cur,1) = closedPts_(j,0);
510 DofCoords(cur,2) = openPts_(k,0);
517 for (
int k=0;k<closedPts_.dimension(0);k++)
519 for (
int j=0;j<openPts_.dimension(0);j++)
521 for (
int i=0;i<openPts_.dimension(0);i++)
523 DofCoords(cur,0) = openPts_(i,0);
524 DofCoords(cur,1) = openPts_(j,0);
525 DofCoords(cur,2) = closedPts_(k,0);
virtual void getDofCoords(ArrayScalar &DofCoords) const
Returns spatial locations (coordinates) of degrees of freedom on a reference cell; defined for interp...
EBasis basisType_
Type of the basis.
void setOrdinalTagData(std::vector< std::vector< std::vector< int > > > &tagToOrdinal, std::vector< std::vector< int > > &ordinalToTag, const int *tags, const int basisCard, const int tagSize, const int posScDim, const int posScOrd, const int posDfOrd)
Fills ordinalToTag_ and tagToOrdinal_ by basis-specific tag data.
Basis_HDIV_HEX_In_FEM(int order, const ArrayScalar &ptsClosed, const ArrayScalar &ptsOpen)
Constructor.
bool basisTagsAreSet_
"true" if tagToOrdinal_ and ordinalToTag_ have been initialized
EOperator
Enumeration of primitive operators available in Intrepid. Primitive operators act on reconstructed fu...
ECoordinates basisCoordinates_
The coordinate system for which the basis is defined.
shards::CellTopology basisCellTopology_
Base topology of the cells for which the basis is defined. See the Shards package http://trilinos...
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
Evaluation of a FEM basis on a reference Hexahedral cell.
EPointType
Enumeration of types of point distributions in Intrepid.
static void lineProduct3d(const int dim0, const int entity0, const int dim1, const int entity1, const int dim2, const int entity2, int &resultdim, int &resultentity)
void initializeTags()
Initializes tagToOrdinal_ and ordinalToTag_ lookup arrays.
int basisDegree_
Degree of the largest complete polynomial space that can be represented by the basis.
int basisCardinality_
Cardinality of the basis, i.e., the number of basis functions/degrees-of-freedom. ...