Intrepid
Intrepid_HGRAD_QUAD_C1_FEMDef.hpp
Go to the documentation of this file.
1 #ifndef INTREPID_HGRAD_QUAD_C1_FEMDEF_HPP
2 #define INTREPID_HGRAD_QUAD_C1_FEMDEF_HPP
3 // @HEADER
4 // ************************************************************************
5 //
6 // Intrepid Package
7 // Copyright (2007) Sandia Corporation
8 //
9 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
10 // license for use of this work by or on behalf of the U.S. Government.
11 //
12 // Redistribution and use in source and binary forms, with or without
13 // modification, are permitted provided that the following conditions are
14 // met:
15 //
16 // 1. Redistributions of source code must retain the above copyright
17 // notice, this list of conditions and the following disclaimer.
18 //
19 // 2. Redistributions in binary form must reproduce the above copyright
20 // notice, this list of conditions and the following disclaimer in the
21 // documentation and/or other materials provided with the distribution.
22 //
23 // 3. Neither the name of the Corporation nor the names of the
24 // contributors may be used to endorse or promote products derived from
25 // this software without specific prior written permission.
26 //
27 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38 //
39 // Questions? Contact Pavel Bochev (pbboche@sandia.gov)
40 // Denis Ridzal (dridzal@sandia.gov), or
41 // Kara Peterson (kjpeter@sandia.gov)
42 //
43 // ************************************************************************
44 // @HEADER
45 
51 namespace Intrepid {
52 
53  template<class Scalar, class ArrayScalar>
55  {
56  this -> basisCardinality_ = 4;
57  this -> basisDegree_ = 1;
58  this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Quadrilateral<4> >() );
59  this -> basisType_ = BASIS_FEM_DEFAULT;
60  this -> basisCoordinates_ = COORDINATES_CARTESIAN;
61  this -> basisTagsAreSet_ = false;
62  }
63 
64 
65 template<class Scalar, class ArrayScalar>
67 
68  // Basis-dependent intializations
69  int tagSize = 4; // size of DoF tag, i.e., number of fields in the tag
70  int posScDim = 0; // position in the tag, counting from 0, of the subcell dim
71  int posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal
72  int posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell
73 
74  // An array with local DoF tags assigned to basis functions, in the order of their local enumeration
75  int tags[] = { 0, 0, 0, 1,
76  0, 1, 0, 1,
77  0, 2, 0, 1,
78  0, 3, 0, 1};
79 
80  // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
81  Intrepid::setOrdinalTagData(this -> tagToOrdinal_,
82  this -> ordinalToTag_,
83  tags,
84  this -> basisCardinality_,
85  tagSize,
86  posScDim,
87  posScOrd,
88  posDfOrd);
89 }
90 
91 
92 
93 template<class Scalar, class ArrayScalar>
95  const ArrayScalar & inputPoints,
96  const EOperator operatorType) const {
97 
98  // Verify arguments
99 #ifdef HAVE_INTREPID_DEBUG
100  Intrepid::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
101  inputPoints,
102  operatorType,
103  this -> getBaseCellTopology(),
104  this -> getCardinality() );
105 #endif
106 
107  // Number of evaluation points = dim 0 of inputPoints
108  int dim0 = inputPoints.dimension(0);
109 
110  // Temporaries: (x,y) coordinates of the evaluation point
111  Scalar x = 0.0;
112  Scalar y = 0.0;
113 
114  switch (operatorType) {
115 
116  case OPERATOR_VALUE:
117  for (int i0 = 0; i0 < dim0; i0++) {
118  x = inputPoints(i0, 0);
119  y = inputPoints(i0, 1);
120 
121  // outputValues is a rank-2 array with dimensions (basisCardinality_, dim0)
122  outputValues(0, i0) = (1.0 - x)*(1.0 - y)/4.0;
123  outputValues(1, i0) = (1.0 + x)*(1.0 - y)/4.0;
124  outputValues(2, i0) = (1.0 + x)*(1.0 + y)/4.0;
125  outputValues(3, i0) = (1.0 - x)*(1.0 + y)/4.0;
126  }
127  break;
128 
129  case OPERATOR_GRAD:
130  case OPERATOR_D1:
131  for (int i0 = 0; i0 < dim0; i0++) {
132  x = inputPoints(i0,0);
133  y = inputPoints(i0,1);
134 
135  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
136  outputValues(0, i0, 0) = -(1.0 - y)/4.0;
137  outputValues(0, i0, 1) = -(1.0 - x)/4.0;
138 
139  outputValues(1, i0, 0) = (1.0 - y)/4.0;
140  outputValues(1, i0, 1) = -(1.0 + x)/4.0;
141 
142  outputValues(2, i0, 0) = (1.0 + y)/4.0;
143  outputValues(2, i0, 1) = (1.0 + x)/4.0;
144 
145  outputValues(3, i0, 0) = -(1.0 + y)/4.0;
146  outputValues(3, i0, 1) = (1.0 - x)/4.0;
147  }
148  break;
149 
150  case OPERATOR_CURL:
151  for (int i0 = 0; i0 < dim0; i0++) {
152  x = inputPoints(i0,0);
153  y = inputPoints(i0,1);
154 
155  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
156  outputValues(0, i0, 0) = -(1.0 - x)/4.0;
157  outputValues(0, i0, 1) = (1.0 - y)/4.0;
158 
159  outputValues(1, i0, 0) = -(1.0 + x)/4.0;
160  outputValues(1, i0, 1) = -(1.0 - y)/4.0;
161 
162  outputValues(2, i0, 0) = (1.0 + x)/4.0;
163  outputValues(2, i0, 1) = -(1.0 + y)/4.0;
164 
165  outputValues(3, i0, 0) = (1.0 - x)/4.0;
166  outputValues(3, i0, 1) = (1.0 + y)/4.0;
167  }
168  break;
169 
170  case OPERATOR_DIV:
171  TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
172  ">>> ERROR (Basis_HGRAD_QUAD_C1_FEM): DIV is invalid operator for rank-0 (scalar) functions in 2D");
173  break;
174 
175  case OPERATOR_D2:
176  for (int i0 = 0; i0 < dim0; i0++) {
177 
178  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, D2Cardinality=3)
179  outputValues(0, i0, 0) = 0.0;
180  outputValues(0, i0, 1) = 0.25;
181  outputValues(0, i0, 2) = 0.0;
182 
183  outputValues(1, i0, 0) = 0.0;
184  outputValues(1, i0, 1) = -0.25;
185  outputValues(1, i0, 2) = 0.0;
186 
187  outputValues(2, i0, 0) = 0.0;
188  outputValues(2, i0, 1) = 0.25;
189  outputValues(2, i0, 2) = 0.0;
190 
191  outputValues(3, i0, 0) = 0.0;
192  outputValues(3, i0, 1) = -0.25;
193  outputValues(3, i0, 2) = 0.0;
194  }
195  break;
196 
197  case OPERATOR_D3:
198  case OPERATOR_D4:
199  case OPERATOR_D5:
200  case OPERATOR_D6:
201  case OPERATOR_D7:
202  case OPERATOR_D8:
203  case OPERATOR_D9:
204  case OPERATOR_D10:
205  {
206  // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, DkCardinality)
207  int DkCardinality = Intrepid::getDkCardinality(operatorType,
208  this -> basisCellTopology_.getDimension() );
209  for(int dofOrd = 0; dofOrd < this -> basisCardinality_; dofOrd++) {
210  for (int i0 = 0; i0 < dim0; i0++) {
211  for(int dkOrd = 0; dkOrd < DkCardinality; dkOrd++){
212  outputValues(dofOrd, i0, dkOrd) = 0.0;
213  }
214  }
215  }
216  }
217  break;
218 
219  default:
220  TEUCHOS_TEST_FOR_EXCEPTION( !( Intrepid::isValidOperator(operatorType) ), std::invalid_argument,
221  ">>> ERROR (Basis_HGRAD_QUAD_C1_FEM): Invalid operator type");
222  }
223 }
224 
225 
226 
227 template<class Scalar, class ArrayScalar>
229  const ArrayScalar & inputPoints,
230  const ArrayScalar & cellVertices,
231  const EOperator operatorType) const {
232  TEUCHOS_TEST_FOR_EXCEPTION( (true), std::logic_error,
233  ">>> ERROR (Basis_HGRAD_QUAD_C1_FEM): FEM Basis calling an FVD member function");
234 }
235 
236 
237 
238 template<class Scalar, class ArrayScalar>
240 #ifdef HAVE_INTREPID_DEBUG
241  // Verify rank of output array.
242  TEUCHOS_TEST_FOR_EXCEPTION( !(DofCoords.rank() == 2), std::invalid_argument,
243  ">>> ERROR: (Intrepid::Basis_HGRAD_QUAD_C1_FEM::getDofCoords) rank = 2 required for DofCoords array");
244  // Verify 0th dimension of output array.
245  TEUCHOS_TEST_FOR_EXCEPTION( !( DofCoords.dimension(0) == this -> basisCardinality_ ), std::invalid_argument,
246  ">>> ERROR: (Intrepid::Basis_HGRAD_QUAD_C1_FEM::getDofCoords) mismatch in number of DoF and 0th dimension of DofCoords array");
247  // Verify 1st dimension of output array.
248  TEUCHOS_TEST_FOR_EXCEPTION( !( DofCoords.dimension(1) == (int)(this -> basisCellTopology_.getDimension()) ), std::invalid_argument,
249  ">>> ERROR: (Intrepid::Basis_HGRAD_QUAD_C1_FEM::getDofCoords) incorrect reference cell (1st) dimension in DofCoords array");
250 #endif
251 
252  DofCoords(0,0) = -1.0; DofCoords(0,1) = -1.0;
253  DofCoords(1,0) = 1.0; DofCoords(1,1) = -1.0;
254  DofCoords(2,0) = 1.0; DofCoords(2,1) = 1.0;
255  DofCoords(3,0) = -1.0; DofCoords(3,1) = 1.0;
256 }
257 
258 }// namespace Intrepid
259 #endif
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.
EOperator
Enumeration of primitive operators available in Intrepid. Primitive operators act on reconstructed fu...
int getDkCardinality(const EOperator operatorType, const int spaceDim)
Returns cardinality of Dk, i.e., the number of all derivatives of order k.
void getDofCoords(ArrayScalar &DofCoords) const
Returns spatial locations (coordinates) of degrees of freedom on a reference Quadrilateral.
int isValidOperator(const EOperator operatorType)
Verifies validity of an operator enum.
void getValues(ArrayScalar &outputValues, const ArrayScalar &inputPoints, const EOperator operatorType) const
FEM basis evaluation on a reference Quadrilateral cell.
void initializeTags()
Initializes tagToOrdinal_ and ordinalToTag_ lookup arrays.