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