|
Intrepid
|
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 // This library is free software; you can redistribute it and/or modify 00011 // it under the terms of the GNU Lesser General Public License as 00012 // published by the Free Software Foundation; either version 2.1 of the 00013 // License, or (at your option) any later version. 00014 // 00015 // This library is distributed in the hope that it will be useful, but 00016 // WITHOUT ANY WARRANTY; without even the implied warranty of 00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00018 // Lesser General Public License for more details. 00019 // 00020 // You should have received a copy of the GNU Lesser General Public 00021 // License along with this library; if not, write to the Free Software 00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 00023 // USA 00024 // Questions? Contact Pavel Bochev (pbboche@sandia.gov) or 00025 // Denis Ridzal (dridzal@sandia.gov). 00026 // 00027 // ************************************************************************ 00028 // @HEADER 00029 00035 namespace Intrepid { 00036 00037 template<class Scalar, class ArrayScalar> 00038 Basis_HGRAD_WEDGE_C2_FEM<Scalar, ArrayScalar>::Basis_HGRAD_WEDGE_C2_FEM() 00039 { 00040 this -> basisCardinality_ = 18; 00041 this -> basisDegree_ = 2; 00042 this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Wedge<6> >() ); 00043 this -> basisType_ = BASIS_FEM_DEFAULT; 00044 this -> basisCoordinates_ = COORDINATES_CARTESIAN; 00045 this -> basisTagsAreSet_ = false; 00046 } 00047 00048 00049 template<class Scalar, class ArrayScalar> 00050 void Basis_HGRAD_WEDGE_C2_FEM<Scalar, ArrayScalar>::initializeTags() { 00051 00052 // Basis-dependent intializations 00053 int tagSize = 4; // size of DoF tag 00054 int posScDim = 0; // position in the tag, counting from 0, of the subcell dim 00055 int posScOrd = 1; // position in the tag, counting from 0, of the subcell ordinal 00056 int posDfOrd = 2; // position in the tag, counting from 0, of DoF ordinal relative to the subcell 00057 00058 // An array with local DoF tags assigned to basis functions, in the order of their local enumeration 00059 int tags[] = { 0, 0, 0, 1, 00060 0, 1, 0, 1, 00061 0, 2, 0, 1, 00062 0, 3, 0, 1, 00063 0, 4, 0, 1, 00064 0, 5, 0, 1, 00065 1, 0, 0, 1, 00066 1, 1, 0, 1, 00067 1, 2, 0, 1, 00068 1, 6, 0, 1, 00069 1, 7, 0, 1, 00070 1, 8, 0, 1, 00071 1, 3, 0, 1, 00072 1, 4, 0, 1, 00073 1, 5, 0, 1, 00074 2, 0, 0, 1, 00075 2, 1, 0, 1, 00076 2, 2, 0, 1 00077 }; 00078 00079 // Basis-independent function sets tag and enum data in tagToOrdinal_ and ordinalToTag_ arrays: 00080 Intrepid::setOrdinalTagData(this -> tagToOrdinal_, 00081 this -> ordinalToTag_, 00082 tags, 00083 this -> basisCardinality_, 00084 tagSize, 00085 posScDim, 00086 posScOrd, 00087 posDfOrd); 00088 } 00089 00090 00091 00092 template<class Scalar, class ArrayScalar> 00093 void Basis_HGRAD_WEDGE_C2_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues, 00094 const ArrayScalar & inputPoints, 00095 const EOperator operatorType) const { 00096 00097 // Verify arguments 00098 #ifdef HAVE_INTREPID_DEBUG 00099 Intrepid::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues, 00100 inputPoints, 00101 operatorType, 00102 this -> getBaseCellTopology(), 00103 this -> getCardinality() ); 00104 #endif 00105 00106 // Number of evaluation points = dim 0 of inputPoints 00107 int dim0 = inputPoints.dimension(0); 00108 00109 // Temporaries: (x,y,z) coordinates of the evaluation point 00110 Scalar x = 0.0; 00111 Scalar y = 0.0; 00112 Scalar z = 0.0; 00113 00114 switch (operatorType) { 00115 00116 case OPERATOR_VALUE: 00117 for (int i0 = 0; i0 < dim0; i0++) { 00118 x = inputPoints(i0, 0); 00119 y = inputPoints(i0, 1); 00120 z = inputPoints(i0, 2); 00121 00122 // outputValues is a rank-2 array with dimensions (basisCardinality_, dim0) 00123 outputValues(0, i0) = ((-1. + x + y)*(-1. + 2.*x + 2.*y)*(-1. + z)*z)/2.; 00124 outputValues(1, i0) = (x*(-1. + 2.*x)*(-1. + z)*z)/2.; 00125 outputValues(2, i0) = (y*(-1. + 2.*y)*(-1. + z)*z)/2.; 00126 outputValues(3, i0) = ((-1. + x + y)*(-1. + 2.*x + 2.*y)*z*(1. + z))/2.; 00127 outputValues(4, i0) = (x*(-1. + 2.*x)*z*(1. + z))/2.; 00128 outputValues(5, i0) = (y*(-1. + 2.*y)*z*(1. + z))/2.; 00129 00130 outputValues(6, i0) = -2.*x*(-1. + x + y)*(-1. + z)*z; 00131 outputValues(7, i0) = 2.*x*y*(-1. + z)*z; 00132 outputValues(8, i0) = -2.*y*(-1. + x + y)*(-1. + z)*z; 00133 outputValues(9, i0) = -((-1. + x + y)*(-1. + 2.*x + 2.*y)*(-1. + z)*(1. + z)); 00134 outputValues(10,i0) = -(x*(-1. + 2.*x)*(-1. + z)*(1. + z)); 00135 outputValues(11,i0) = -(y*(-1. + 2.*y)*(-1. + z)*(1. + z)); 00136 outputValues(12,i0) = -2.*x*(-1. + x + y)*z*(1. + z); 00137 outputValues(13,i0) = 2.*x*y*z*(1. + z); 00138 outputValues(14,i0) = -2.*y*(-1. + x + y)*z*(1. + z); 00139 outputValues(15,i0) = 4.*x*(-1. + x + y)*(-1. + z)*(1. + z); 00140 outputValues(16,i0) = -4.*x*y*(-1. + z)*(1. + z); 00141 outputValues(17,i0) = 4.*y*(-1. + x + y)*(-1. + z)*(1. + z); 00142 } 00143 break; 00144 00145 case OPERATOR_GRAD: 00146 case OPERATOR_D1: 00147 for (int i0 = 0; i0 < dim0; i0++) { 00148 x = inputPoints(i0,0); 00149 y = inputPoints(i0,1); 00150 z = inputPoints(i0,2); 00151 00152 // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim) 00153 outputValues(0, i0, 0) = ((-3 + 4*x + 4*y)*(-1 + z)*z)/2.; 00154 outputValues(0, i0, 1) = ((-3 + 4*x + 4*y)*(-1 + z)*z)/2.; 00155 outputValues(0, i0, 2) = ((-1 + x + y)*(-1 + 2*x + 2*y)*(-1 + 2*z))/2.; 00156 00157 outputValues(1, i0, 0) = ((-1 + 4*x)*(-1 + z)*z)/2.; 00158 outputValues(1, i0, 1) = 0.; 00159 outputValues(1, i0, 2) = (x*(-1 + 2*x)*(-1 + 2*z))/2.; 00160 00161 outputValues(2, i0, 0) = 0.; 00162 outputValues(2, i0, 1) = ((-1 + 4*y)*(-1 + z)*z)/2.; 00163 outputValues(2, i0, 2) = (y*(-1 + 2*y)*(-1 + 2*z))/2.; 00164 00165 outputValues(3, i0, 0) = ((-3 + 4*x + 4*y)*z*(1 + z))/2.; 00166 outputValues(3, i0, 1) = ((-3 + 4*x + 4*y)*z*(1 + z))/2.; 00167 outputValues(3, i0, 2) = ((-1 + x + y)*(-1 + 2*x + 2*y)*(1 + 2*z))/2.; 00168 00169 outputValues(4, i0, 0) = ((-1 + 4*x)*z*(1 + z))/2.; 00170 outputValues(4, i0, 1) = 0.; 00171 outputValues(4, i0, 2) = (x*(-1 + 2*x)*(1 + 2*z))/2.; 00172 00173 outputValues(5, i0, 0) = 0.; 00174 outputValues(5, i0, 1) = ((-1 + 4*y)*z*(1 + z))/2.; 00175 outputValues(5, i0, 2) = (y*(-1 + 2*y)*(1 + 2*z))/2.; 00176 00177 outputValues(6, i0, 0) = -2*(-1 + 2*x + y)*(-1 + z)*z; 00178 outputValues(6, i0, 1) = -2*x*(-1 + z)*z; 00179 outputValues(6, i0, 2) = 2*x*(-1 + x + y)*(1 - 2*z); 00180 00181 outputValues(7, i0, 0) = 2*y*(-1 + z)*z; 00182 outputValues(7, i0, 1) = 2*x*(-1 + z)*z; 00183 outputValues(7, i0, 2) = 2*x*y*(-1 + 2*z); 00184 00185 outputValues(8, i0, 0) = -2*y*(-1 + z)*z; 00186 outputValues(8, i0, 1) = -2*(-1 + x + 2*y)*(-1 + z)*z; 00187 outputValues(8, i0, 2) = 2*y*(-1 + x + y)*(1 - 2*z); 00188 00189 outputValues(9, i0, 0) = -(-3 + 4*x + 4*y)*(-1 + z*z); 00190 outputValues(9, i0, 1) = -(-3 + 4*x + 4*y)*(-1 + z*z); 00191 outputValues(9, i0, 2) = -2*(1 + 2*x*x - 3*y + 2*y*y + x*(-3 + 4*y))*z; 00192 00193 outputValues(10,i0, 0) = -(-1 + 4*x)*(-1 + z*z); 00194 outputValues(10,i0, 1) = 0; 00195 outputValues(10,i0, 2) = 2*(1 - 2*x)*x*z; 00196 00197 outputValues(11,i0, 0) = 0; 00198 outputValues(11,i0, 1) = -(-1 + 4*y)*(-1 + z*z); 00199 outputValues(11,i0, 2) = 2*(1 - 2*y)*y*z; 00200 00201 outputValues(12,i0, 0) = -2*(-1 + 2*x + y)*z*(1 + z); 00202 outputValues(12,i0, 1) = -2*x*z*(1 + z); 00203 outputValues(12,i0, 2) = -2*x*(-1 + x + y)*(1 + 2*z); 00204 00205 outputValues(13,i0, 0) = 2*y*z*(1 + z); 00206 outputValues(13,i0, 1) = 2*x*z*(1 + z); 00207 outputValues(13,i0, 2) = 2*x*y*(1 + 2*z); 00208 00209 outputValues(14,i0, 0) = -2*y*z*(1 + z); 00210 outputValues(14,i0, 1) = -2*(-1 + x + 2*y)*z*(1 + z); 00211 outputValues(14,i0, 2) = -2*y*(-1 + x + y)*(1 + 2*z); 00212 00213 outputValues(15,i0, 0) = 4*(-1 + 2*x + y)*(-1 + z*z); 00214 outputValues(15,i0, 1) = 4*x*(-1 + z)*(1 + z); 00215 outputValues(15,i0, 2) = 8*x*(-1 + x + y)*z; 00216 00217 outputValues(16,i0, 0) = -4*y*(-1 + z)*(1 + z); 00218 outputValues(16,i0, 1) = -4*x*(-1 + z)*(1 + z); 00219 outputValues(16,i0, 2) = -8*x*y*z; 00220 00221 outputValues(17,i0, 0) = 4*y*(-1 + z)*(1 + z); 00222 outputValues(17,i0, 1) = 4*(-1 + x + 2*y)*(-1 + z*z); 00223 outputValues(17,i0, 2) = 8*y*(-1 + x + y)*z; 00224 00225 } 00226 break; 00227 00228 case OPERATOR_CURL: 00229 TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument, 00230 ">>> ERROR (Basis_HGRAD_WEDGE_C2_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D"); 00231 break; 00232 00233 case OPERATOR_DIV: 00234 TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument, 00235 ">>> ERROR (Basis_HGRAD_WEDGE_C2_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D"); 00236 break; 00237 00238 case OPERATOR_D2: 00239 for (int i0 = 0; i0 < dim0; i0++) { 00240 x = inputPoints(i0,0); 00241 y = inputPoints(i0,1); 00242 z = inputPoints(i0,2); 00243 00244 outputValues(0, i0, 0) = 2.*(-1. + z)*z; 00245 outputValues(0, i0, 1) = 2.*(-1. + z)*z; 00246 outputValues(0, i0, 2) = ((-3. + 4.*x + 4.*y)*(-1. + 2.*z))/2.; 00247 outputValues(0, i0, 3) = 2.*(-1. + z)*z; 00248 outputValues(0, i0, 4) = ((-3. + 4.*x + 4.*y)*(-1. + 2.*z))/2.; 00249 outputValues(0, i0, 5) = (-1. + x + y)*(-1. + 2.*x + 2.*y); 00250 00251 outputValues(1, i0, 0) = 2.*(-1. + z)*z; 00252 outputValues(1, i0, 1) = 0.; 00253 outputValues(1, i0, 2) = ((-1. + 4.*x)*(-1. + 2.*z))/2.; 00254 outputValues(1, i0, 3) = 0.; 00255 outputValues(1, i0, 4) = 0.; 00256 outputValues(1, i0, 5) = x*(-1. + 2.*x); 00257 00258 outputValues(2, i0, 0) = 0.; 00259 outputValues(2, i0, 1) = 0.; 00260 outputValues(2, i0, 2) = 0.; 00261 outputValues(2, i0, 3) = 2.*(-1. + z)*z; 00262 outputValues(2, i0, 4) = ((-1. + 4.*y)*(-1. + 2.*z))/2.; 00263 outputValues(2, i0, 5) = y*(-1. + 2.*y); 00264 00265 outputValues(3, i0, 0) = 2.*z*(1. + z); 00266 outputValues(3, i0, 1) = 2.*z*(1. + z); 00267 outputValues(3, i0, 2) = ((-3. + 4.*x + 4.*y)*(1. + 2.*z))/2.; 00268 outputValues(3, i0, 3) = 2.*z*(1. + z); 00269 outputValues(3, i0, 4) = ((-3. + 4.*x + 4.*y)*(1. + 2.*z))/2.; 00270 outputValues(3, i0, 5) = (-1. + x + y)*(-1. + 2.*x + 2.*y); 00271 00272 outputValues(4, i0, 0) = 2.*z*(1. + z); 00273 outputValues(4, i0, 1) = 0.; 00274 outputValues(4, i0, 2) = ((-1. + 4.*x)*(1. + 2.*z))/2.;; 00275 outputValues(4, i0, 3) = 0.; 00276 outputValues(4, i0, 4) = 0.; 00277 outputValues(4, i0, 5) = x*(-1. + 2.*x); 00278 00279 outputValues(5, i0, 0) = 0.; 00280 outputValues(5, i0, 1) = 0.; 00281 outputValues(5, i0, 2) = 0.; 00282 outputValues(5, i0, 3) = 2.*z*(1. + z); 00283 outputValues(5, i0, 4) = ((-1. + 4.*y)*(1. + 2.*z))/2.; 00284 outputValues(5, i0, 5) = y*(-1. + 2.*y); 00285 00286 outputValues(6, i0, 0) = -4.*(-1. + z)*z; 00287 outputValues(6, i0, 1) = -2.*(-1. + z)*z; 00288 outputValues(6, i0, 2) = -2.*(-1. + 2.*x + y)*(-1. + 2.*z); 00289 outputValues(6, i0, 3) = 0.; 00290 outputValues(6, i0, 4) = x*(2. - 4.*z); 00291 outputValues(6, i0, 5) = -4.*x*(-1. + x + y); 00292 00293 outputValues(7, i0, 0) = 0.; 00294 outputValues(7, i0, 1) = 2.*(-1. + z)*z; 00295 outputValues(7, i0, 2) = 2.*y*(-1. + 2.*z); 00296 outputValues(7, i0, 3) = 0.; 00297 outputValues(7, i0, 4) = 2.*x*(-1. + 2.*z); 00298 outputValues(7, i0, 5) = 4.*x*y; 00299 00300 outputValues(8, i0, 0) = 0.; 00301 outputValues(8, i0, 1) = -2.*(-1. + z)*z; 00302 outputValues(8, i0, 2) = y*(2. - 4.*z); 00303 outputValues(8, i0, 3) = -4.*(-1. + z)*z; 00304 outputValues(8, i0, 4) = -2.*(-1. + x + 2.*y)*(-1. + 2.*z); 00305 outputValues(8, i0, 5) = -4.*y*(-1. + x + y); 00306 00307 outputValues(9, i0, 0) = 4. - 4.*z*z; 00308 outputValues(9, i0, 1) = 4. - 4.*z*z; 00309 outputValues(9, i0, 2) = -2.*(-3. + 4.*x + 4.*y)*z; 00310 outputValues(9, i0, 3) = 4. - 4.*z*z; 00311 outputValues(9, i0, 4) = -2.*(-3. + 4.*x + 4.*y)*z; 00312 outputValues(9, i0, 5) = -2.*(-1. + x + y)*(-1. + 2.*x + 2.*y); 00313 00314 outputValues(10,i0, 0) = 4. - 4.*z*z; 00315 outputValues(10,i0, 1) = 0.; 00316 outputValues(10,i0, 2) = (2. - 8.*x)*z; 00317 outputValues(10,i0, 3) = 0.; 00318 outputValues(10,i0, 4) = 0.; 00319 outputValues(10,i0, 5) = -2.*x*(-1. + 2.*x); 00320 00321 outputValues(11,i0, 0) = 0.; 00322 outputValues(11,i0, 1) = 0.; 00323 outputValues(11,i0, 2) = 0.; 00324 outputValues(11,i0, 3) = 4. - 4.*z*z; 00325 outputValues(11,i0, 4) = (2. - 8.*y)*z; 00326 outputValues(11,i0, 5) = -2.*y*(-1. + 2.*y); 00327 00328 outputValues(12,i0, 0) = -4.*z*(1. + z); 00329 outputValues(12,i0, 1) = -2.*z*(1. + z); 00330 outputValues(12,i0, 2) = -2.*(-1. + 2.*x + y)*(1. + 2.*z); 00331 outputValues(12,i0, 3) = 0.; 00332 outputValues(12,i0, 4) = -2.*(x + 2.*x*z); 00333 outputValues(12,i0, 5) = -4.*x*(-1. + x + y); 00334 00335 outputValues(13,i0, 0) = 0.; 00336 outputValues(13,i0, 1) = 2.*z*(1. + z); 00337 outputValues(13,i0, 2) = 2.*(y + 2.*y*z); 00338 outputValues(13,i0, 3) = 0.; 00339 outputValues(13,i0, 4) = 2.*(x + 2.*x*z); 00340 outputValues(13,i0, 5) = 4.*x*y; 00341 00342 outputValues(14,i0, 0) = 0.; 00343 outputValues(14,i0, 1) = -2.*z*(1. + z); 00344 outputValues(14,i0, 2) = -2.*(y + 2.*y*z); 00345 outputValues(14,i0, 3) = -4.*z*(1. + z); 00346 outputValues(14,i0, 4) = -2.*(-1. + x + 2.*y)*(1. + 2.*z); 00347 outputValues(14,i0, 5) = -4.*y*(-1. + x + y); 00348 00349 outputValues(15,i0, 0) = 8.*(-1. + z*z); 00350 outputValues(15,i0, 1) = 4.*(-1. + z*z); 00351 outputValues(15,i0, 2) = 8.*(-1. + 2.*x + y)*z; 00352 outputValues(15,i0, 3) = 0.; 00353 outputValues(15,i0, 4) = 8.*x*z; 00354 outputValues(15,i0, 5) = 8.*x*(-1. + x + y); 00355 00356 outputValues(16,i0, 0) = 0.; 00357 outputValues(16,i0, 1) = 4. - 4.*z*z; 00358 outputValues(16,i0, 2) = -8.*y*z; 00359 outputValues(16,i0, 3) = 0.; 00360 outputValues(16,i0, 4) = -8.*x*z; 00361 outputValues(16,i0, 5) = -8.*x*y; 00362 00363 00364 outputValues(17,i0, 0) = 0.; 00365 outputValues(17,i0, 1) = 4.*(-1. + z*z); 00366 outputValues(17,i0, 2) = 8.*y*z; 00367 outputValues(17,i0, 3) = 8.*(-1. + z*z); 00368 outputValues(17,i0, 4) = 8.*(-1. + x + 2.*y)*z; 00369 outputValues(17,i0, 5) = 8.*y*(-1. + x + y); 00370 } 00371 break; 00372 00373 case OPERATOR_D3: 00374 for (int i0 = 0; i0 < dim0; i0++) { 00375 x = inputPoints(i0,0); 00376 y = inputPoints(i0,1); 00377 z = inputPoints(i0,2); 00378 00379 outputValues(0, i0, 0) = 0.; 00380 outputValues(0, i0, 1) = 0.; 00381 outputValues(0, i0, 2) = -2. + 4.*z; 00382 outputValues(0, i0, 3) = 0.; 00383 outputValues(0, i0, 4) = -2. + 4.*z; 00384 outputValues(0, i0, 5) = -3. + 4.*x + 4.*y; 00385 outputValues(0, i0, 6) = 0.; 00386 outputValues(0, i0, 7) = -2. + 4.*z; 00387 outputValues(0, i0, 8) = -3. + 4.*x + 4.*y; 00388 outputValues(0, i0, 9) = 0.; 00389 00390 outputValues(1, i0, 0) = 0.; 00391 outputValues(1, i0, 1) = 0.; 00392 outputValues(1, i0, 2) = -2. + 4.*z; 00393 outputValues(1, i0, 3) = 0.; 00394 outputValues(1, i0, 4) = 0.; 00395 outputValues(1, i0, 5) = -1 + 4.*x; 00396 outputValues(1, i0, 6) = 0.; 00397 outputValues(1, i0, 7) = 0.; 00398 outputValues(1, i0, 8) = 0.; 00399 outputValues(1, i0, 9) = 0.; 00400 00401 outputValues(2, i0, 0) = 0.; 00402 outputValues(2, i0, 1) = 0.; 00403 outputValues(2, i0, 2) = 0.; 00404 outputValues(2, i0, 3) = 0.; 00405 outputValues(2, i0, 4) = 0.; 00406 outputValues(2, i0, 5) = 0.; 00407 outputValues(2, i0, 6) = 0.; 00408 outputValues(2, i0, 7) = -2. + 4.*z; 00409 outputValues(2, i0, 8) = -1 + 4.*y; 00410 outputValues(2, i0, 9) = 0.; 00411 00412 outputValues(3, i0, 0) = 0.; 00413 outputValues(3, i0, 1) = 0.; 00414 outputValues(3, i0, 2) = 2. + 4.*z; 00415 outputValues(3, i0, 3) = 0.; 00416 outputValues(3, i0, 4) = 2. + 4.*z; 00417 outputValues(3, i0, 5) = -3. + 4.*x + 4.*y; 00418 outputValues(3, i0, 6) = 0.; 00419 outputValues(3, i0, 7) = 2. + 4.*z; 00420 outputValues(3, i0, 8) = -3. + 4.*x + 4.*y; 00421 outputValues(3, i0, 9) = 0.; 00422 00423 outputValues(4, i0, 0) = 0.; 00424 outputValues(4, i0, 1) = 0.; 00425 outputValues(4, i0, 2) = 2. + 4.*z; 00426 outputValues(4, i0, 3) = 0.; 00427 outputValues(4, i0, 4) = 0.; 00428 outputValues(4, i0, 5) = -1 + 4.*x; 00429 outputValues(4, i0, 6) = 0.; 00430 outputValues(4, i0, 7) = 0.; 00431 outputValues(4, i0, 8) = 0.; 00432 outputValues(4, i0, 9) = 0.; 00433 00434 outputValues(5, i0, 0) = 0.; 00435 outputValues(5, i0, 1) = 0.; 00436 outputValues(5, i0, 2) = 0.; 00437 outputValues(5, i0, 3) = 0.; 00438 outputValues(5, i0, 4) = 0.; 00439 outputValues(5, i0, 5) = 0.; 00440 outputValues(5, i0, 6) = 0.; 00441 outputValues(5, i0, 7) = 2. + 4.*z; 00442 outputValues(5, i0, 8) = -1 + 4.*y; 00443 outputValues(5, i0, 9) = 0.; 00444 00445 outputValues(6, i0, 0) = 0.; 00446 outputValues(6, i0, 1) = 0.; 00447 outputValues(6, i0, 2) = 4. - 8.*z; 00448 outputValues(6, i0, 3) = 0.; 00449 outputValues(6, i0, 4) = 2. - 4.*z; 00450 outputValues(6, i0, 5) = -4.*(-1 + 2*x + y); 00451 outputValues(6, i0, 6) = 0.; 00452 outputValues(6, i0, 7) = 0.; 00453 outputValues(6, i0, 8) = -4.*x; 00454 outputValues(6, i0, 9) = 0.; 00455 00456 outputValues(7, i0, 0) = 0.; 00457 outputValues(7, i0, 1) = 0.; 00458 outputValues(7, i0, 2) = 0.; 00459 outputValues(7, i0, 3) = 0.; 00460 outputValues(7, i0, 4) = -2. + 4.*z; 00461 outputValues(7, i0, 5) = 4.*y; 00462 outputValues(7, i0, 6) = 0.; 00463 outputValues(7, i0, 7) = 0.; 00464 outputValues(7, i0, 8) = 4.*x; 00465 outputValues(7, i0, 9) = 0.; 00466 00467 outputValues(8, i0, 0) = 0.; 00468 outputValues(8, i0, 1) = 0.; 00469 outputValues(8, i0, 2) = 0.; 00470 outputValues(8, i0, 3) = 0.; 00471 outputValues(8, i0, 4) = 2. - 4.*z; 00472 outputValues(8, i0, 5) = -4.*y; 00473 outputValues(8, i0, 6) = 0.; 00474 outputValues(8, i0, 7) = 4. - 8.*z; 00475 outputValues(8, i0, 8) = -4.*(-1 + x + 2*y); 00476 outputValues(8, i0, 9) = 0.; 00477 00478 outputValues(9, i0, 0) = 0.; 00479 outputValues(9, i0, 1) = 0.; 00480 outputValues(9, i0, 2) = -8.*z; 00481 outputValues(9, i0, 3) = 0.; 00482 outputValues(9, i0, 4) = -8.*z; 00483 outputValues(9, i0, 5) = 6. - 8.*x - 8.*y; 00484 outputValues(9, i0, 6) = 0.; 00485 outputValues(9, i0, 7) = -8.*z; 00486 outputValues(9, i0, 8) = 6. - 8.*x - 8.*y; 00487 outputValues(9, i0, 9) = 0.; 00488 00489 outputValues(10,i0, 0) = 0.; 00490 outputValues(10,i0, 1) = 0.; 00491 outputValues(10,i0, 2) = -8.*z; 00492 outputValues(10,i0, 3) = 0.; 00493 outputValues(10,i0, 4) = 0.; 00494 outputValues(10,i0, 5) = 2. - 8.*x; 00495 outputValues(10,i0, 6) = 0.; 00496 outputValues(10,i0, 7) = 0.; 00497 outputValues(10,i0, 8) = 0.; 00498 outputValues(10,i0, 9) = 0.; 00499 00500 outputValues(11,i0, 0) = 0.; 00501 outputValues(11,i0, 1) = 0.; 00502 outputValues(11,i0, 2) = 0.; 00503 outputValues(11,i0, 3) = 0.; 00504 outputValues(11,i0, 4) = 0.; 00505 outputValues(11,i0, 5) = 0.; 00506 outputValues(11,i0, 6) = 0.; 00507 outputValues(11,i0, 7) = -8.*z; 00508 outputValues(11,i0, 8) = 2. - 8.*y; 00509 outputValues(11,i0, 9) = 0.; 00510 00511 outputValues(12,i0, 0) = 0.; 00512 outputValues(12,i0, 1) = 0.; 00513 outputValues(12,i0, 2) = -4. - 8.*z; 00514 outputValues(12,i0, 3) = 0.; 00515 outputValues(12,i0, 4) = -2. - 4.*z; 00516 outputValues(12,i0, 5) = -4.*(-1 + 2*x + y); 00517 outputValues(12,i0, 6) = 0.; 00518 outputValues(12,i0, 7) = 0.; 00519 outputValues(12,i0, 8) = -4.*x; 00520 outputValues(12,i0, 9) = 0.; 00521 00522 outputValues(13,i0, 0) = 0.; 00523 outputValues(13,i0, 1) = 0.; 00524 outputValues(13,i0, 2) = 0.; 00525 outputValues(13,i0, 3) = 0.; 00526 outputValues(13,i0, 4) = 2. + 4.*z; 00527 outputValues(13,i0, 5) = 4.*y; 00528 outputValues(13,i0, 6) = 0.; 00529 outputValues(13,i0, 7) = 0.; 00530 outputValues(13,i0, 8) = 4.*x; 00531 outputValues(13,i0, 9) = 0.; 00532 00533 outputValues(14,i0, 0) = 0.; 00534 outputValues(14,i0, 1) = 0.; 00535 outputValues(14,i0, 2) = 0.; 00536 outputValues(14,i0, 3) = 0.; 00537 outputValues(14,i0, 4) = -2. - 4.*z; 00538 outputValues(14,i0, 5) = -4.*y; 00539 outputValues(14,i0, 6) = 0.; 00540 outputValues(14,i0, 7) = -4. - 8.*z; 00541 outputValues(14,i0, 8) = -4.*(-1 + x + 2*y); 00542 outputValues(14,i0, 9) = 0.; 00543 00544 outputValues(15,i0, 0) = 0.; 00545 outputValues(15,i0, 1) = 0.; 00546 outputValues(15,i0, 2) = 16.*z; 00547 outputValues(15,i0, 3) = 0.; 00548 outputValues(15,i0, 4) = 8.*z; 00549 outputValues(15,i0, 5) = 8.*(-1 + 2*x + y); 00550 outputValues(15,i0, 6) = 0.; 00551 outputValues(15,i0, 7) = 0.; 00552 outputValues(15,i0, 8) = 8.*x; 00553 outputValues(15,i0, 9) = 0.; 00554 00555 outputValues(16,i0, 0) = 0.; 00556 outputValues(16,i0, 1) = 0.; 00557 outputValues(16,i0, 2) = 0.; 00558 outputValues(16,i0, 3) = 0.; 00559 outputValues(16,i0, 4) = -8.*z; 00560 outputValues(16,i0, 5) = -8.*y; 00561 outputValues(16,i0, 6) = 0.; 00562 outputValues(16,i0, 7) = 0.; 00563 outputValues(16,i0, 8) = -8.*x; 00564 outputValues(16,i0, 9) = 0.; 00565 00566 outputValues(17,i0, 0) = 0.; 00567 outputValues(17,i0, 1) = 0.; 00568 outputValues(17,i0, 2) = 0.; 00569 outputValues(17,i0, 3) = 0.; 00570 outputValues(17,i0, 4) = 8.*z; 00571 outputValues(17,i0, 5) = 8.*y; 00572 outputValues(17,i0, 6) = 0.; 00573 outputValues(17,i0, 7) = 16.*z; 00574 outputValues(17,i0, 8) = 8.*(-1 + x + 2*y); 00575 outputValues(17,i0, 9) = 0.; 00576 00577 } 00578 break; 00579 00580 case OPERATOR_D4: 00581 { 00582 // There are only few constant non-zero entries. Initialize by zero and then assign non-zero entries. 00583 int DkCardinality = Intrepid::getDkCardinality(operatorType, this -> basisCellTopology_.getDimension() ); 00584 for(int dofOrd = 0; dofOrd < this -> basisCardinality_; dofOrd++) { 00585 for (int i0 = 0; i0 < dim0; i0++) { 00586 for(int dkOrd = 0; dkOrd < DkCardinality; dkOrd++){ 00587 outputValues(dofOrd, i0, dkOrd) = 0.0; 00588 } 00589 } 00590 } 00591 00592 for (int i0 = 0; i0 < dim0; i0++) { 00593 00594 outputValues(0, i0, 5) = 4.; 00595 outputValues(0, i0, 8) = 4.; 00596 outputValues(0, i0,12) = 4.; 00597 00598 outputValues(1, i0, 5) = 4.; 00599 00600 outputValues(2, i0,12) = 4.; 00601 00602 outputValues(3, i0, 5) = 4.; 00603 outputValues(3, i0, 8) = 4.; 00604 outputValues(3, i0,12) = 4.; 00605 00606 outputValues(4, i0, 5) = 4.0; 00607 00608 outputValues(5, i0,12) = 4.0; 00609 00610 outputValues(6, i0, 5) =-8.; 00611 outputValues(6, i0, 8) =-4.; 00612 00613 outputValues(7, i0, 8) = 4.; 00614 00615 outputValues(8, i0, 8) =-4.; 00616 outputValues(8, i0,12) =-8.; 00617 00618 outputValues(9, i0, 5) =-8.; 00619 outputValues(9, i0, 8) =-8.; 00620 outputValues(9, i0,12) =-8.; 00621 00622 outputValues(10,i0, 5) =-8.; 00623 00624 outputValues(11,i0,12) =-8.; 00625 00626 outputValues(12,i0, 5) =-8.; 00627 outputValues(12,i0, 8) =-4.; 00628 00629 outputValues(13,i0, 8) = 4.; 00630 00631 outputValues(14,i0, 8) =-4; 00632 outputValues(14,i0,12) =-8.; 00633 00634 outputValues(15,i0, 5) =16.; 00635 outputValues(15,i0, 8) = 8.; 00636 00637 outputValues(16,i0, 8) =-8.; 00638 00639 00640 outputValues(17,i0, 8) = 8.; 00641 outputValues(17,i0,12) =16.; 00642 } 00643 } 00644 break; 00645 00646 case OPERATOR_D5: 00647 case OPERATOR_D6: 00648 case OPERATOR_D7: 00649 case OPERATOR_D8: 00650 case OPERATOR_D9: 00651 case OPERATOR_D10: 00652 { 00653 // outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, DkCardinality) 00654 int DkCardinality = Intrepid::getDkCardinality(operatorType, 00655 this -> basisCellTopology_.getDimension() ); 00656 for(int dofOrd = 0; dofOrd < this -> basisCardinality_; dofOrd++) { 00657 for (int i0 = 0; i0 < dim0; i0++) { 00658 for(int dkOrd = 0; dkOrd < DkCardinality; dkOrd++){ 00659 outputValues(dofOrd, i0, dkOrd) = 0.0; 00660 } 00661 } 00662 } 00663 } 00664 break; 00665 00666 default: 00667 TEST_FOR_EXCEPTION( !( Intrepid::isValidOperator(operatorType) ), std::invalid_argument, 00668 ">>> ERROR (Basis_HGRAD_WEDGE_C2_FEM): Invalid operator type"); 00669 } 00670 } 00671 00672 00673 00674 template<class Scalar, class ArrayScalar> 00675 void Basis_HGRAD_WEDGE_C2_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar& outputValues, 00676 const ArrayScalar & inputPoints, 00677 const ArrayScalar & cellVertices, 00678 const EOperator operatorType) const { 00679 TEST_FOR_EXCEPTION( (true), std::logic_error, 00680 ">>> ERROR (Basis_HGRAD_WEDGE_C2_FEM): FEM Basis calling an FVD member function"); 00681 } 00682 }// namespace Intrepid
1.7.4