Intrepid
http://trilinos.sandia.gov/packages/docs/dev/packages/intrepid/src/Discretization/Basis/Intrepid_HCURL_HEX_I1_FEMDef.hpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ************************************************************************
00003 //
00004 //                           Intrepid Package
00005 //                 Copyright (2007) Sandia Corporation
00006 //
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Pavel Bochev  (pbboche@sandia.gov)
00038 //                    Denis Ridzal  (dridzal@sandia.gov), or
00039 //                    Kara Peterson (kjpeter@sandia.gov)
00040 //
00041 // ************************************************************************
00042 // @HEADER
00043 
00049 namespace Intrepid {
00050 
00051 template<class Scalar, class ArrayScalar>
00052 Basis_HCURL_HEX_I1_FEM<Scalar,ArrayScalar>::Basis_HCURL_HEX_I1_FEM()
00053   {
00054     this -> basisCardinality_  = 12;
00055     this -> basisDegree_       = 1;
00056     this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
00057     this -> basisType_         = BASIS_FEM_DEFAULT;
00058     this -> basisCoordinates_  = COORDINATES_CARTESIAN;
00059     this -> basisTagsAreSet_   = false;
00060   }
00061   
00062 template<class Scalar, class ArrayScalar>
00063 void Basis_HCURL_HEX_I1_FEM<Scalar, ArrayScalar>::initializeTags() {
00064   
00065   // Basis-dependent intializations
00066   int tagSize  = 4;        // size of DoF tag
00067   int posScDim = 0;        // position in the tag, counting from 0, of the subcell dim 
00068   int posScOrd = 1;        // position in the tag, counting from 0, of the subcell ordinal
00069   int posDfOrd = 2;        // position in the tag, counting from 0, of DoF ordinal relative to the subcell
00070 
00071   // An array with local DoF tags assigned to basis functions, in the order of their local enumeration 
00072   int tags[]  = {
00073                   1, 0, 0, 1,
00074                   1, 1, 0, 1,
00075                   1, 2, 0, 1,
00076                   1, 3, 0, 1,
00077                   1, 4, 0, 1,
00078                   1, 5, 0, 1,
00079                   1, 6, 0, 1,
00080                   1, 7, 0, 1,
00081                   1, 8, 0, 1,
00082                   1, 9, 0, 1,
00083                   1, 10, 0, 1,
00084                   1, 11, 0, 1 };
00085   
00086   // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays:
00087   Intrepid::setOrdinalTagData(this -> tagToOrdinal_,
00088                               this -> ordinalToTag_,
00089                               tags,
00090                               this -> basisCardinality_,
00091                               tagSize,
00092                               posScDim,
00093                               posScOrd,
00094                               posDfOrd);
00095 }
00096 
00097 
00098 
00099 template<class Scalar, class ArrayScalar>
00100 void Basis_HCURL_HEX_I1_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar &        outputValues,
00101                                                             const ArrayScalar &  inputPoints,
00102                                                             const EOperator      operatorType) const {
00103   
00104 // Verify arguments
00105 #ifdef HAVE_INTREPID_DEBUG
00106   Intrepid::getValues_HCURL_Args<Scalar, ArrayScalar>(outputValues,
00107                                                       inputPoints,
00108                                                       operatorType,
00109                                                       this -> getBaseCellTopology(),
00110                                                       this -> getCardinality() );
00111 #endif
00112   
00113  // Number of evaluation points = dim 0 of inputPoints
00114   int dim0 = inputPoints.dimension(0);
00115 
00116   // Temporaries: (x,y,z) coordinates of the evaluation point
00117   Scalar x = 0.0;                                    
00118   Scalar y = 0.0;                                    
00119   Scalar z = 0.0;                                    
00120   
00121   switch (operatorType) {
00122     case OPERATOR_VALUE:
00123       for (int i0 = 0; i0 < dim0; i0++) {
00124         x = inputPoints(i0, 0);
00125         y = inputPoints(i0, 1);
00126         z = inputPoints(i0, 2);
00127         
00128         // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
00129         outputValues(0, i0, 0) = (1.0 - y)*(1.0 - z)/8.0;
00130         outputValues(0, i0, 1) = 0.0;
00131         outputValues(0, i0, 2) = 0.0;
00132 
00133         outputValues(1, i0, 0) = 0.0;
00134         outputValues(1, i0, 1) = (1.0 + x)*(1.0 - z)/8.0;
00135         outputValues(1, i0, 2) = 0.0;
00136 
00137         outputValues(2, i0, 0) = -(1.0 + y)*(1.0 - z)/8.0;
00138         outputValues(2, i0, 1) = 0.0;
00139         outputValues(2, i0, 2) = 0.0;
00140 
00141         outputValues(3, i0, 0) = 0.0;
00142         outputValues(3, i0, 1) = -(1.0 - x)*(1.0 - z)/8.0;
00143         outputValues(3, i0, 2) = 0.0;
00144 
00145         outputValues(4, i0, 0) = (1.0 - y)*(1.0 + z)/8.0;
00146         outputValues(4, i0, 1) = 0.0;
00147         outputValues(4, i0, 2) = 0.0;
00148 
00149         outputValues(5, i0, 0) = 0.0;
00150         outputValues(5, i0, 1) = (1.0 + x)*(1.0 + z)/8.0;
00151         outputValues(5, i0, 2) = 0.0;
00152 
00153         outputValues(6, i0, 0) = -(1.0 + y)*(1.0 + z)/8.0;
00154         outputValues(6, i0, 1) = 0.0;
00155         outputValues(6, i0, 2) = 0.0;
00156 
00157         outputValues(7, i0, 0) = 0.0;
00158         outputValues(7, i0, 1) = -(1.0 - x)*(1.0 + z)/8.0;
00159         outputValues(7, i0, 2) = 0.0;
00160 
00161         outputValues(8, i0, 0) = 0.0;
00162         outputValues(8, i0, 1) = 0.0;
00163         outputValues(8, i0, 2) = (1.0 - x)*(1.0 - y)/8.0;
00164 
00165         outputValues(9, i0, 0) = 0.0;
00166         outputValues(9, i0, 1) = 0.0;
00167         outputValues(9, i0, 2) = (1.0 + x)*(1.0 - y)/8.0;
00168 
00169         outputValues(10, i0, 0) = 0.0;
00170         outputValues(10, i0, 1) = 0.0;
00171         outputValues(10, i0, 2) = (1.0 + x)*(1.0 + y)/8.0;
00172 
00173         outputValues(11, i0, 0) = 0.0;
00174         outputValues(11, i0, 1) = 0.0;
00175         outputValues(11, i0, 2) = (1.0 - x)*(1.0 + y)/8.0;
00176       }
00177       break;
00178       
00179     case OPERATOR_CURL:
00180       for (int i0 = 0; i0 < dim0; i0++) {
00181         x = inputPoints(i0, 0);
00182         y = inputPoints(i0, 1);
00183         z = inputPoints(i0, 2);
00184         
00185         // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
00186         outputValues(0, i0, 0) = 0.0;
00187         outputValues(0, i0, 1) = -(1.0 - y)/8.0;
00188         outputValues(0, i0, 2) = (1.0 - z)/8.0;
00189 
00190         outputValues(1, i0, 0) = (1.0 + x)/8.0;
00191         outputValues(1, i0, 1) = 0.0;
00192         outputValues(1, i0, 2) = (1.0 - z)/8.0;
00193 
00194         outputValues(2, i0, 0) = 0.0;
00195         outputValues(2, i0, 1) = (1.0 + y)/8.0;
00196         outputValues(2, i0, 2) = (1.0 - z)/8.0;
00197 
00198         outputValues(3, i0, 0) = -(1.0 - x)/8.0;
00199         outputValues(3, i0, 1) = 0.0;
00200         outputValues(3, i0, 2) = (1.0 - z)/8.0;
00201 
00202         outputValues(4, i0, 0) = 0.0;
00203         outputValues(4, i0, 1) = (1.0 - y)/8.0;
00204         outputValues(4, i0, 2) = (1.0 + z)/8.0;
00205 
00206         outputValues(5, i0, 0) = -(1.0 + x)/8.0;
00207         outputValues(5, i0, 1) = 0.0;
00208         outputValues(5, i0, 2) = (1.0 + z)/8.0;
00209 
00210         outputValues(6, i0, 0) = 0.0;
00211         outputValues(6, i0, 1) = -(1.0 + y)/8.0;
00212         outputValues(6, i0, 2) = (1.0 + z)/8.0;
00213 
00214         outputValues(7, i0, 0) = (1.0 - x)/8.0;
00215         outputValues(7, i0, 1) = 0.0;
00216         outputValues(7, i0, 2) = (1.0 + z)/8.0;
00217 
00218         outputValues(8, i0, 0) = -(1.0 - x)/8.0;
00219         outputValues(8, i0, 1) = (1.0 - y)/8.0;
00220         outputValues(8, i0, 2) = 0.0;
00221 
00222         outputValues(9, i0, 0) = -(1.0 + x)/8.0;
00223         outputValues(9, i0, 1) = -(1.0 - y)/8.0;
00224         outputValues(9, i0, 2) = 0.0;
00225 
00226         outputValues(10, i0, 0) = (1.0 + x)/8.0;
00227         outputValues(10, i0, 1) = -(1.0 + y)/8.0;
00228         outputValues(10, i0, 2) = 0.0;
00229 
00230         outputValues(11, i0, 0) = (1.0 - x)/8.0;
00231         outputValues(11, i0, 1) = (1.0 + y)/8.0;
00232         outputValues(11, i0, 2) = 0.0;
00233       }
00234       break;
00235       
00236     case OPERATOR_DIV:
00237        TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
00238                           ">>> ERROR (Basis_HCURL_HEX_I1_FEM): DIV is invalid operator for HCURL Basis Functions");
00239       break;
00240       
00241     case OPERATOR_GRAD:
00242        TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_GRAD), std::invalid_argument,
00243                           ">>> ERROR (Basis_HCURL_HEX_I1_FEM): GRAD is invalid operator for HCURL Basis Functions");
00244       break;
00245 
00246     case OPERATOR_D1:
00247     case OPERATOR_D2:
00248     case OPERATOR_D3:
00249     case OPERATOR_D4:
00250     case OPERATOR_D5:
00251     case OPERATOR_D6:
00252     case OPERATOR_D7:
00253     case OPERATOR_D8:
00254     case OPERATOR_D9:
00255     case OPERATOR_D10:
00256       TEUCHOS_TEST_FOR_EXCEPTION( ( (operatorType == OPERATOR_D1)    ||
00257                             (operatorType == OPERATOR_D2)    ||
00258                             (operatorType == OPERATOR_D3)    ||
00259                             (operatorType == OPERATOR_D4)    ||
00260                             (operatorType == OPERATOR_D5)    ||
00261                             (operatorType == OPERATOR_D6)    ||
00262                             (operatorType == OPERATOR_D7)    ||
00263                             (operatorType == OPERATOR_D8)    ||
00264                             (operatorType == OPERATOR_D9)    ||
00265                             (operatorType == OPERATOR_D10) ),
00266                           std::invalid_argument,
00267                           ">>> ERROR (Basis_HCURL_HEX_I1_FEM): Invalid operator type");
00268       break;
00269 
00270     default:
00271       TEUCHOS_TEST_FOR_EXCEPTION( ( (operatorType != OPERATOR_VALUE) &&
00272                             (operatorType != OPERATOR_GRAD)  &&
00273                             (operatorType != OPERATOR_CURL)  &&
00274                             (operatorType != OPERATOR_DIV)   &&
00275                             (operatorType != OPERATOR_D1)    &&
00276                             (operatorType != OPERATOR_D2)    &&
00277                             (operatorType != OPERATOR_D3)    &&
00278                             (operatorType != OPERATOR_D4)    &&
00279                             (operatorType != OPERATOR_D5)    &&
00280                             (operatorType != OPERATOR_D6)    &&
00281                             (operatorType != OPERATOR_D7)    &&
00282                             (operatorType != OPERATOR_D8)    &&
00283                             (operatorType != OPERATOR_D9)    &&
00284                             (operatorType != OPERATOR_D10) ),
00285                           std::invalid_argument,
00286                           ">>> ERROR (Basis_HCURL_HEX_I1_FEM): Invalid operator type");
00287   }
00288 }
00289 
00290 
00291   
00292 template<class Scalar, class ArrayScalar>
00293 void Basis_HCURL_HEX_I1_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar&           outputValues,
00294                                                             const ArrayScalar &    inputPoints,
00295                                                             const ArrayScalar &    cellVertices,
00296                                                             const EOperator        operatorType) const {
00297   TEUCHOS_TEST_FOR_EXCEPTION( (true), std::logic_error,
00298                       ">>> ERROR (Basis_HCURL_HEX_I1_FEM): FEM Basis calling an FVD member function");
00299                                                              }
00300 
00301 }// namespace Intrepid