MOAB: Mesh Oriented datABase  (version 5.4.1)
HexLagrangeShape.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2006 Lawrence Livermore National Laboratory.  Under
00005     the terms of Contract B545069 with the University of Wisconsin --
00006     Madison, Lawrence Livermore National Laboratory retains certain
00007     rights in this software.
00008 
00009     This library is free software; you can redistribute it and/or
00010     modify it under the terms of the GNU Lesser General Public
00011     License as published by the Free Software Foundation; either
00012     version 2.1 of the License, or (at your option) any later version.
00013 
00014     This library is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017     Lesser General Public License for more details.
00018 
00019     You should have received a copy of the GNU Lesser General Public License
00020     (lgpl.txt) along with this library; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 
00023     [email protected]
00024 
00025   ***************************************************************** */
00026 /** \file HexLagrangeShape.cpp
00027  *  \author Nicholas Voshell
00028  */
00029 
00030 #include "HexLagrangeShape.hpp"
00031 #include "MsqError.hpp"
00032 
00033 namespace MBMesquite
00034 {
00035 
00036 EntityTopology HexLagrangeShape::element_topology() const
00037 {
00038     return HEXAHEDRON;
00039 }
00040 
00041 int HexLagrangeShape::num_nodes() const
00042 {
00043     return 27;
00044 }
00045 
00046 void HexLagrangeShape::coefficients( Sample loc,
00047                                      NodeSet nodeset,
00048                                      double* coeff_out,
00049                                      size_t* indices_out,
00050                                      size_t& num_coeff,
00051                                      MsqError& err ) const
00052 {
00053     if( nodeset.num_nodes() != 27 )
00054     {
00055         MSQ_SETERR( err )
00056         ( "Mapping function supports only 27-node hexahedra with no slaved nodes.", MsqError::UNSUPPORTED_ELEMENT );
00057         return;
00058     }
00059 
00060     switch( loc.dimension )
00061     {
00062         case 0:  // Corner sample point - assume that it is there always
00063             num_coeff      = 1;
00064             indices_out[0] = loc.number;
00065             coeff_out[0]   = 1.0;
00066             break;
00067         case 1:  // Line sample point - check if it is there,
00068             num_coeff      = 1;
00069             indices_out[0] = loc.number + 8;
00070             coeff_out[0]   = 1.0;
00071             break;
00072         case 2:  // Face sample point - check if it is there,
00073             num_coeff      = 1;
00074             indices_out[0] = loc.number + 20;
00075             coeff_out[0]   = 1.0;
00076             break;
00077         case 3:  // Center sample point - check if it is there,
00078             num_coeff      = 1;
00079             indices_out[0] = 26;
00080             coeff_out[0]   = 1.0;
00081             break;
00082         default:
00083             MSQ_SETERR( err )
00084             ( MsqError::UNSUPPORTED_ELEMENT,
00085               "Request for dimension %d mapping function value"
00086               "for a quadratic hexahedral element",
00087               loc.dimension );
00088     }
00089 }
00090 
00091 void HexLagrangeShape::derivatives( Sample loc,
00092                                     NodeSet nodeset,
00093                                     size_t* vertices,
00094                                     MsqVector< 3 >* derivs,
00095                                     size_t& num_vtx,
00096                                     MsqError& err ) const
00097 {
00098     if( nodeset.num_nodes() != 27 )
00099     {
00100         MSQ_SETERR( err )
00101         ( "Mapping function supports only 27-node hexahedra with no slaved nodes.", MsqError::UNSUPPORTED_ELEMENT );
00102         return;
00103     }
00104 
00105     // r coordinate
00106     const int HLS1[] = { 1, 3, 3, 1, 1, 3, 3, 1, 2, 3, 2, 1, 1, 3, 3, 1, 2, 3, 2, 1, 2, 3, 2, 1, 2, 2, 2 };
00107     // s coordinate
00108     const int HLS2[] = { 1, 1, 3, 3, 1, 1, 3, 3, 1, 2, 3, 2, 1, 1, 3, 3, 1, 2, 3, 2, 1, 2, 3, 2, 2, 2, 2 };
00109     // t coordinate
00110     const int HLS3[] = { 1, 1, 1, 1, 3, 3, 3, 3, 1, 1, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 2, 2, 2, 2, 1, 3, 2 };
00111 
00112     const int HLSup1[] = { 12, 13, 14, 15, -1, -1, -1, -1, 20, 21, 22, 23, 4, 5,
00113                            6,  7,  -1, -1, -1, -1, 16, 17, 18, 19, 26, -1, 25 };
00114     const int HLSup2[] = { 4,  5,  6,  7,  -1, -1, -1, -1, 16, 17, 18, 19, -1, -1,
00115                            -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, 25, -1, -1 };
00116     const int HLSdn1[] = { -1, -1, -1, -1, 12, 13, 14, 15, -1, -1, -1, -1, 0, 1,
00117                            2,  3,  20, 21, 22, 23, 8,  9,  10, 11, -1, 26, 24 };
00118     const int HLSdn2[] = { -1, -1, -1, -1, 0,  1,  2,  3,  -1, -1, -1, -1, -1, -1,
00119                            -1, -1, 8,  9,  10, 11, -1, -1, -1, -1, -1, 24, -1 };
00120 
00121     const int HLSlt1[] = { -1, 8,  10, -1, -1, 16, 18, -1, 0,  24, 3,  -1, -1, 20,
00122                            22, -1, 4,  25, 7,  -1, 12, 26, 15, -1, 11, 19, 23 };
00123     const int HLSlt2[] = { -1, 0,  3,  -1, -1, 4,  7,  -1, -1, 11, -1, -1, -1, 12,
00124                            15, -1, -1, 19, -1, -1, -1, 23, -1, -1, -1, -1, -1 };
00125     const int HLSrt1[] = { 8,  -1, -1, 10, 16, -1, -1, 18, 1,  -1, 2, 24, 20, -1,
00126                            -1, 22, 5,  -1, 6,  25, 13, -1, 14, 26, 9, 17, 21 };
00127     const int HLSrt2[] = { 1,  -1, -1, 2,  5,  -1, -1, 6,  -1, -1, -1, 9,  13, -1,
00128                            -1, 14, -1, -1, -1, 17, -1, -1, -1, 21, -1, -1, -1 };
00129 
00130     const int HLSft1[] = { -1, -1, 9,  11, -1, -1, 17, 19, -1, 1,  24, 0,  -1, -1,
00131                            21, 23, -1, 5,  25, 4,  -1, 13, 26, 12, 8,  16, 20 };
00132     const int HLSft2[] = { -1, -1, 1,  0,  -1, -1, 5,  4,  -1, -1, 8,  -1, -1, -1,
00133                            13, 12, -1, -1, 16, -1, -1, -1, 20, -1, -1, -1, -1 };
00134     const int HLSbk1[] = { 11, 9,  -1, -1, 19, 17, -1, -1, 24, 2,  -1, 3,  23, 21,
00135                            -1, -1, 25, 6,  -1, 7,  26, 14, -1, 15, 10, 18, 22 };
00136     const int HLSbk2[] = { 3,  2,  -1, -1, 7,  6,  -1, -1, 10, -1, -1, -1, 15, 14,
00137                            -1, -1, 18, -1, -1, -1, 22, -1, -1, -1, -1, -1, -1 };
00138 
00139     double entries[] = { 0, -3, 0, 3 };
00140 
00141     int location = loc.number;
00142 
00143     switch( loc.dimension )
00144     {
00145         case 1:
00146             location += 8;
00147             break;
00148         case 2:
00149             location += 20;
00150             break;
00151         case 3:
00152             location += 26;
00153             break;
00154     }
00155 
00156     num_vtx = 0;
00157 
00158     if( location != 26 )
00159     {
00160         vertices[num_vtx]  = location;
00161         derivs[num_vtx][0] = entries[HLS1[location]];
00162         derivs[num_vtx][1] = entries[HLS2[location]];
00163         derivs[num_vtx][2] = entries[HLS3[location]];
00164         num_vtx++;
00165     }
00166 
00167     if( HLSup1[location] != -1 )
00168     {
00169         vertices[num_vtx] = HLSup1[location];
00170         if( HLS3[location] == 2 )
00171         {
00172             derivs[num_vtx][0] = 0;
00173             derivs[num_vtx][1] = 0;
00174             derivs[num_vtx][2] = 1;
00175         }
00176         else
00177         {
00178             derivs[num_vtx][0] = 0;
00179             derivs[num_vtx][1] = 0;
00180             derivs[num_vtx][2] = 4;
00181         }
00182         num_vtx++;
00183     }
00184 
00185     if( HLSdn1[location] != -1 )
00186     {
00187         vertices[num_vtx] = HLSdn1[location];
00188         if( HLS3[location] == 2 )
00189         {
00190             derivs[num_vtx][0] = 0;
00191             derivs[num_vtx][1] = 0;
00192             derivs[num_vtx][2] = -1;
00193         }
00194         else
00195         {
00196             derivs[num_vtx][0] = 0;
00197             derivs[num_vtx][1] = 0;
00198             derivs[num_vtx][2] = -4;
00199         }
00200         num_vtx++;
00201     }
00202 
00203     if( HLSup2[location] != -1 )
00204     {
00205         vertices[num_vtx]  = HLSup2[location];
00206         derivs[num_vtx][0] = 0;
00207         derivs[num_vtx][1] = 0;
00208         derivs[num_vtx][2] = -1;
00209         num_vtx++;
00210     }
00211 
00212     if( HLSdn2[location] != -1 )
00213     {
00214         vertices[num_vtx]  = HLSdn2[location];
00215         derivs[num_vtx][0] = 0;
00216         derivs[num_vtx][1] = 0;
00217         derivs[num_vtx][2] = 1;
00218         num_vtx++;
00219     }
00220 
00221     if( HLSlt1[location] != -1 )
00222     {
00223         vertices[num_vtx] = HLSlt1[location];
00224         if( HLS1[location] == 2 )
00225         {
00226             derivs[num_vtx][0] = -1;
00227             derivs[num_vtx][1] = 0;
00228             derivs[num_vtx][2] = 0;
00229         }
00230         else
00231         {
00232             derivs[num_vtx][0] = -4;
00233             derivs[num_vtx][1] = 0;
00234             derivs[num_vtx][2] = 0;
00235         }
00236         num_vtx++;
00237     }
00238 
00239     if( HLSrt1[location] != -1 )
00240     {
00241         vertices[num_vtx] = HLSrt1[location];
00242         if( HLS1[location] == 2 )
00243         {
00244             derivs[num_vtx][0] = 1;
00245             derivs[num_vtx][1] = 0;
00246             derivs[num_vtx][2] = 0;
00247         }
00248         else
00249         {
00250             derivs[num_vtx][0] = 4;
00251             derivs[num_vtx][1] = 0;
00252             derivs[num_vtx][2] = 0;
00253         }
00254         num_vtx++;
00255     }
00256 
00257     if( HLSlt2[location] != -1 )
00258     {
00259         vertices[num_vtx]  = HLSlt2[location];
00260         derivs[num_vtx][0] = 1;
00261         derivs[num_vtx][1] = 0;
00262         derivs[num_vtx][2] = 0;
00263         num_vtx++;
00264     }
00265 
00266     if( HLSrt2[location] != -1 )
00267     {
00268         vertices[num_vtx]  = HLSrt2[location];
00269         derivs[num_vtx][0] = -1;
00270         derivs[num_vtx][1] = 0;
00271         derivs[num_vtx][2] = 0;
00272         num_vtx++;
00273     }
00274 
00275     if( HLSft1[location] != -1 )
00276     {
00277         vertices[num_vtx] = HLSft1[location];
00278         if( HLS2[location] == 2 )
00279         {
00280             derivs[num_vtx][0] = 0;
00281             derivs[num_vtx][1] = -1;
00282             derivs[num_vtx][2] = 0;
00283         }
00284         else
00285         {
00286             derivs[num_vtx][0] = 0;
00287             derivs[num_vtx][1] = -4;
00288             derivs[num_vtx][2] = 0;
00289         }
00290         num_vtx++;
00291     }
00292 
00293     if( HLSbk1[location] != -1 )
00294     {
00295         vertices[num_vtx] = HLSbk1[location];
00296         if( HLS2[location] == 2 )
00297         {
00298             derivs[num_vtx][0] = 0;
00299             derivs[num_vtx][1] = 1;
00300             derivs[num_vtx][2] = 0;
00301         }
00302         else
00303         {
00304             derivs[num_vtx][0] = 0;
00305             derivs[num_vtx][1] = 4;
00306             derivs[num_vtx][2] = 0;
00307         }
00308         num_vtx++;
00309     }
00310 
00311     if( HLSft2[location] != -1 )
00312     {
00313         vertices[num_vtx]  = HLSft2[location];
00314         derivs[num_vtx][0] = 0;
00315         derivs[num_vtx][1] = 1;
00316         derivs[num_vtx][2] = 0;
00317         num_vtx++;
00318     }
00319 
00320     if( HLSbk2[location] != -1 )
00321     {
00322         vertices[num_vtx]  = HLSbk2[location];
00323         derivs[num_vtx][0] = 0;
00324         derivs[num_vtx][1] = -1;
00325         derivs[num_vtx][2] = 0;
00326         num_vtx++;
00327     }
00328 }
00329 
00330 void HexLagrangeShape::ideal( Sample, MsqMatrix< 3, 3 >& J, MsqError& ) const
00331 {
00332     J = MsqMatrix< 3, 3 >( 1.0 );
00333 }
00334 
00335 }  // namespace MBMesquite
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines