libMesh::H1FETransformation< OutputShape > Class Template Reference

#include <h1_fe_transformation.h>

Inheritance diagram for libMesh::H1FETransformation< OutputShape >:

List of all members.

Public Member Functions

 H1FETransformation ()
virtual ~H1FETransformation ()
virtual void map_phi (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< OutputShape > &, std::vector< std::vector< OutputShape > > &) const
virtual void map_dphi (const unsigned int dim, const Elem *const elem, const std::vector< Point > &qp, const FEGenericBase< OutputShape > &fe, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputGradient > > &dphi, std::vector< std::vector< OutputShape > > &dphidx, std::vector< std::vector< OutputShape > > &dphidy, std::vector< std::vector< OutputShape > > &dphidz) const
virtual void map_d2phi (const unsigned int dim, const Elem *const elem, const std::vector< Point > &qp, const FEGenericBase< OutputShape > &fe, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputTensor > > &d2phi, std::vector< std::vector< OutputShape > > &d2phidx2, std::vector< std::vector< OutputShape > > &d2phidxdy, std::vector< std::vector< OutputShape > > &d2phidxdz, std::vector< std::vector< OutputShape > > &d2phidy2, std::vector< std::vector< OutputShape > > &d2phidydz, std::vector< std::vector< OutputShape > > &d2phidz2) const
virtual void map_curl (const unsigned int dim, const Elem *const elem, const std::vector< Point > &qp, const FEGenericBase< OutputShape > &fe, std::vector< std::vector< OutputShape > > &curl_phi) const
virtual void map_div (const unsigned int dim, const Elem *const elem, const std::vector< Point > &qp, const FEGenericBase< OutputShape > &fe, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputDivergence > > &div_phi) const
template<>
void map_curl (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< Real > &, std::vector< std::vector< Real > > &) const
template<>
void map_curl (const unsigned int dim, const Elem *const, const std::vector< Point > &, const FEGenericBase< RealGradient > &fe, std::vector< std::vector< RealGradient > > &curl_phi) const
template<>
void map_div (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< Real > &, std::vector< std::vector< FEGenericBase< Real >::OutputDivergence > > &) const
template<>
void map_div (const unsigned int dim, const Elem *const, const std::vector< Point > &, const FEGenericBase< RealGradient > &fe, std::vector< std::vector< FEGenericBase< RealGradient >::OutputDivergence > > &div_phi) const

Static Public Member Functions

static AutoPtr
< FETransformationBase
< OutputShape > > 
build (const FEType &type)

Detailed Description

template<typename OutputShape>
class libMesh::H1FETransformation< OutputShape >

This class handles the computation of the shape functions in the physical domain for H1 conforming elements. This class assumes the FEGenericBase object has been initialized in the reference domain (i.e. init_shape_functions has been called).

Author:
Paul T. Bauman, 2012

Definition at line 35 of file h1_fe_transformation.h.


Constructor & Destructor Documentation

template<typename OutputShape >
libMesh::H1FETransformation< OutputShape >::H1FETransformation (  )  [inline]

Definition at line 39 of file h1_fe_transformation.h.

00040       : FETransformationBase<OutputShape>(){};

template<typename OutputShape >
virtual libMesh::H1FETransformation< OutputShape >::~H1FETransformation (  )  [inline, virtual]

Definition at line 42 of file h1_fe_transformation.h.

00042 {};


Member Function Documentation

template<typename OutputShape >
AutoPtr< FETransformationBase< OutputShape > > libMesh::FETransformationBase< OutputShape >::build ( const FEType type  )  [inline, static, inherited]

Builds an FETransformation object based on the finite element type

Definition at line 26 of file fe_transformation_base.C.

References libMeshEnums::BERNSTEIN, libMeshEnums::CLOUGH, libMesh::FEType::family, libMeshEnums::HERMITE, libMeshEnums::HIERARCHIC, libMeshEnums::JACOBI_20_00, libMeshEnums::JACOBI_30_00, libMeshEnums::L2_HIERARCHIC, libMeshEnums::L2_LAGRANGE, libMeshEnums::LAGRANGE, libMeshEnums::LAGRANGE_VEC, libMeshEnums::MONOMIAL, libMeshEnums::NEDELEC_ONE, libMeshEnums::SCALAR, libMeshEnums::SZABAB, and libMeshEnums::XYZ.

00027   {
00028     switch( fe_type.family )
00029       {
00030         /* H1 Conforming Elements */
00031       case LAGRANGE:
00032       case HIERARCHIC:
00033       case BERNSTEIN:
00034       case SZABAB:
00035       case CLOUGH: /* PB: Really H2 */
00036       case HERMITE: /* PB: Really H2 */
00037       case LAGRANGE_VEC:
00038       case MONOMIAL: /* PB: Shouldn't this be L2 conforming? */
00039       case XYZ: /* PB: Shouldn't this be L2 conforming? */
00040       case L2_HIERARCHIC: /* PB: Shouldn't this be L2 conforming? */
00041       case L2_LAGRANGE: /* PB: Shouldn't this be L2 conforming? */
00042       case JACOBI_20_00: /* PB: For infinite elements... */
00043       case JACOBI_30_00: /* PB: For infinite elements... */
00044         {
00045           AutoPtr<FETransformationBase<OutputShape> > ap( new H1FETransformation<OutputShape> );
00046           return ap;
00047         }
00048         /* HCurl Conforming Elements */
00049       case NEDELEC_ONE:
00050         {
00051           AutoPtr<FETransformationBase<OutputShape> > ap( new HCurlFETransformation<OutputShape> );
00052           return ap;
00053         }
00054 
00055         /* HDiv Conforming Elements */
00056         /* L2 Conforming Elements */
00057 
00058         /* Other... */
00059       case SCALAR:
00060         {
00061           // Should never need this for SCALARs
00062           AutoPtr<FETransformationBase<OutputShape> > ap( new H1FETransformation<OutputShape> );
00063           return ap;
00064         }
00065         
00066       default:
00067         libmesh_error();
00068       }
00069 
00070     // Should never get here...
00071     libmesh_error();
00072     AutoPtr<FETransformationBase<OutputShape> > ap( NULL );
00073     return ap;
00074   }

template<>
void libMesh::H1FETransformation< RealGradient >::map_curl ( const unsigned int  dim,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< RealGradient > &  fe,
std::vector< std::vector< RealGradient > > &  curl_phi 
) const [inline]

Definition at line 468 of file h1_fe_transformation.C.

References libMesh::FEMap::get_detadx(), libMesh::FEMap::get_detady(), libMesh::FEMap::get_detadz(), libMesh::FEGenericBase< OutputType >::get_dphideta(), libMesh::FEGenericBase< OutputType >::get_dphidxi(), libMesh::FEGenericBase< OutputType >::get_dphidzeta(), libMesh::FEMap::get_dxidx(), libMesh::FEMap::get_dxidy(), libMesh::FEMap::get_dxidz(), libMesh::FEMap::get_dzetadx(), libMesh::FEMap::get_dzetady(), libMesh::FEMap::get_dzetadz(), libMesh::FEAbstract::get_fe_map(), and libMesh::Real.

00473   {
00474     switch(dim)
00475       {
00476         // The curl only make sense in 2D and 3D
00477       case 0: 
00478       case 1:
00479         {
00480           libmesh_error();
00481         }
00482       case 2:
00483         {
00484           const std::vector<std::vector<RealGradient> >& dphidxi = fe.get_dphidxi();
00485           const std::vector<std::vector<RealGradient> >& dphideta = fe.get_dphideta();
00486 
00487           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00488           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00489 #if LIBMESH_DIM > 2
00490           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00491 #endif
00492 
00493           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00494           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00495 #if LIBMESH_DIM > 2
00496           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00497 #endif
00498 
00499           /*
00500             For 2D elements in 3D space:
00501 
00502             curl = ( -dphi_y/dz, dphi_x/dz, dphi_y/dx - dphi_x/dy )
00503            */
00504           for (unsigned int i=0; i<curl_phi.size(); i++)
00505             for (unsigned int p=0; p<curl_phi[i].size(); p++)
00506               {
00507 
00508                 Real dphiy_dx = (dphidxi[i][p].slice(1))*dxidx_map[p]
00509                                 + (dphideta[i][p].slice(1))*detadx_map[p];
00510 
00511                 Real dphix_dy = (dphidxi[i][p].slice(0))*dxidy_map[p]
00512                                 + (dphideta[i][p].slice(0))*detady_map[p];
00513 
00514                 curl_phi[i][p].slice(2) = dphiy_dx - dphix_dy;
00515 
00516 #if LIBMESH_DIM > 2
00517                 Real dphiy_dz = (dphidxi[i][p].slice(1))*dxidz_map[p]
00518                                 + (dphideta[i][p].slice(1))*detadz_map[p];
00519 
00520                 Real dphix_dz = (dphidxi[i][p].slice(0))*dxidz_map[p]
00521                                 + (dphideta[i][p].slice(0))*detadz_map[p];
00522                 
00523                 curl_phi[i][p].slice(0) = -dphiy_dz;
00524                 curl_phi[i][p].slice(1) = dphix_dz;
00525 #endif
00526               }
00527 
00528           break;
00529         }
00530       case 3:
00531         {
00532           const std::vector<std::vector<RealGradient> >& dphidxi = fe.get_dphidxi();
00533           const std::vector<std::vector<RealGradient> >& dphideta = fe.get_dphideta();
00534           const std::vector<std::vector<RealGradient> >& dphidzeta = fe.get_dphidzeta();
00535 
00536           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00537           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00538           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00539           
00540           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00541           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00542           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00543           
00544           const std::vector<Real>& dzetadx_map = fe.get_fe_map().get_dzetadx();
00545           const std::vector<Real>& dzetady_map = fe.get_fe_map().get_dzetady();
00546           const std::vector<Real>& dzetadz_map = fe.get_fe_map().get_dzetadz();
00547 
00548           /*
00549             In 3D: curl = ( dphi_z/dy - dphi_y/dz, dphi_x/dz - dphi_z/dx, dphi_y/dx - dphi_x/dy )
00550            */
00551           for (unsigned int i=0; i<curl_phi.size(); i++)
00552             for (unsigned int p=0; p<curl_phi[i].size(); p++)
00553               {
00554                 Real dphiz_dy = (dphidxi[i][p].slice(2))*dxidy_map[p]
00555                                 + (dphideta[i][p].slice(2))*detady_map[p]
00556                                 + (dphidzeta[i][p].slice(2))*dzetady_map[p];
00557 
00558                 Real dphiy_dz = (dphidxi[i][p].slice(1))*dxidz_map[p]
00559                                 + (dphideta[i][p].slice(1))*detadz_map[p]
00560                                 + (dphidzeta[i][p].slice(1))*dzetadz_map[p];
00561 
00562                 Real dphix_dz = (dphidxi[i][p].slice(0))*dxidz_map[p]
00563                                 + (dphideta[i][p].slice(0))*detadz_map[p]
00564                                 + (dphidzeta[i][p].slice(0))*dzetadz_map[p];
00565 
00566                 Real dphiz_dx = (dphidxi[i][p].slice(2))*dxidx_map[p]
00567                                 + (dphideta[i][p].slice(2))*detadx_map[p]
00568                                 + (dphidzeta[i][p].slice(2))*dzetadx_map[p];
00569                 
00570                 Real dphiy_dx = (dphidxi[i][p].slice(1))*dxidx_map[p]
00571                                 + (dphideta[i][p].slice(1))*detadx_map[p]
00572                                 + (dphidzeta[i][p].slice(1))*dzetadx_map[p];
00573 
00574                 Real dphix_dy = (dphidxi[i][p].slice(0))*dxidy_map[p]
00575                                 + (dphideta[i][p].slice(0))*detady_map[p]
00576                                 + (dphidzeta[i][p].slice(0))*dzetady_map[p];
00577 
00578                 curl_phi[i][p].slice(0) = dphiz_dy - dphiy_dz;
00579 
00580                 curl_phi[i][p].slice(1) = dphix_dz - dphiz_dx;
00581 
00582                 curl_phi[i][p].slice(2) = dphiy_dx - dphix_dy;
00583               }
00584 
00585           break;
00586         }
00587       default:
00588         libmesh_error();
00589       }
00590 
00591     return;
00592   }

template<>
void libMesh::H1FETransformation< Real >::map_curl ( const unsigned int  ,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< Real > &  ,
std::vector< std::vector< Real > > &   
) const [inline]

Definition at line 456 of file h1_fe_transformation.C.

References libMesh::err.

00461   {
00462     libMesh::err << "Computing the curl of a shape function only\n"
00463                  << "makes sense for vector-valued elements." << std::endl;
00464     libmesh_error();
00465   }

template<typename OutputShape >
virtual void libMesh::H1FETransformation< OutputShape >::map_curl ( const unsigned int  dim,
const Elem *const   elem,
const std::vector< Point > &  qp,
const FEGenericBase< OutputShape > &  fe,
std::vector< std::vector< OutputShape > > &  curl_phi 
) const [virtual]

Evaluates the shape function curl in physical coordinates based on H1 conforming finite element transformation.

Implements libMesh::FETransformationBase< OutputShape >.

template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::map_d2phi ( const unsigned int  dim,
const Elem *const   elem,
const std::vector< Point > &  qp,
const FEGenericBase< OutputShape > &  fe,
std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputTensor > > &  d2phi,
std::vector< std::vector< OutputShape > > &  d2phidx2,
std::vector< std::vector< OutputShape > > &  d2phidxdy,
std::vector< std::vector< OutputShape > > &  d2phidxdz,
std::vector< std::vector< OutputShape > > &  d2phidy2,
std::vector< std::vector< OutputShape > > &  d2phidydz,
std::vector< std::vector< OutputShape > > &  d2phidz2 
) const [inline, virtual]

Evaluates shape function Hessians in physical coordinates based on H1 conforming finite element transformation. FIXME: These are currently not calculated correctly for non-affine elements. The second derivative terms of the FEMap are not implemented.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 219 of file h1_fe_transformation.C.

References libMesh::err, libMesh::FEGenericBase< OutputType >::get_d2phideta2(), libMesh::FEGenericBase< OutputType >::get_d2phidetadzeta(), libMesh::FEGenericBase< OutputType >::get_d2phidxi2(), libMesh::FEGenericBase< OutputType >::get_d2phidxideta(), libMesh::FEGenericBase< OutputType >::get_d2phidxidzeta(), libMesh::FEGenericBase< OutputType >::get_d2phidzeta2(), libMesh::FEMap::get_detadx(), libMesh::FEMap::get_detady(), libMesh::FEMap::get_detadz(), libMesh::FEMap::get_dxidx(), libMesh::FEMap::get_dxidy(), libMesh::FEMap::get_dxidz(), libMesh::FEMap::get_dzetadx(), libMesh::FEMap::get_dzetady(), libMesh::FEMap::get_dzetadz(), libMesh::FEAbstract::get_fe_map(), and libMesh::Elem::has_affine_map().

00230   {
00231     libmesh_do_once(
00232                     if (!elem->has_affine_map())
00233                       {
00234                         libMesh::err << "WARNING: Second derivatives are not currently "
00235                                      << "correctly calculated on non-affine elements!"
00236                                      << std::endl;
00237                       }
00238                     );
00239       
00240 
00241     switch(dim)
00242       {
00243       case 0: // No derivatives in 0D
00244         {
00245           for (unsigned int i=0; i<d2phi.size(); i++)
00246             for (unsigned int p=0; p<d2phi[i].size(); p++)
00247               {
00248                 d2phi[i][p] = 0.;
00249               }
00250 
00251           break;
00252         }
00253 
00254       case 1:
00255         {
00256           const std::vector<std::vector<OutputShape> >& d2phidxi2 = fe.get_d2phidxi2();
00257 
00258           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00259 #if LIBMESH_DIM>1
00260           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00261 #endif
00262 #if LIBMESH_DIM>2
00263           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00264 #endif
00265           
00266           for (unsigned int i=0; i<d2phi.size(); i++)
00267             for (unsigned int p=0; p<d2phi[i].size(); p++)
00268               {
00269                 d2phi[i][p].slice(0).slice(0) = d2phidx2[i][p] =
00270                   d2phidxi2[i][p]*dxidx_map[p]*dxidx_map[p];
00271 #if LIBMESH_DIM>1
00272                 d2phi[i][p].slice(0).slice(1) = 
00273                   d2phi[i][p].slice(1).slice(0) = d2phidxdy[i][p] =
00274                   d2phidxi2[i][p]*dxidx_map[p]*dxidy_map[p];
00275 
00276                 d2phi[i][p].slice(1).slice(1) = d2phidy2[i][p] =
00277                   d2phidxi2[i][p]*dxidy_map[p]*dxidy_map[p];
00278 #endif
00279 #if LIBMESH_DIM>2
00280                 d2phi[i][p].slice(0).slice(2) =
00281                   d2phi[i][p].slice(2).slice(0) = d2phidxdz[i][p] =
00282                   d2phidxi2[i][p]*dxidx_map[p]*dxidz_map[p];
00283 
00284                 d2phi[i][p].slice(1).slice(2) =
00285                   d2phi[i][p].slice(2).slice(1) = d2phidydz[i][p] =
00286                   d2phidxi2[i][p]*dxidy_map[p]*dxidz_map[p];
00287 
00288                 d2phi[i][p].slice(2).slice(2) = d2phidz2[i][p] =
00289                   d2phidxi2[i][p]*dxidz_map[p]*dxidz_map[p];
00290 #endif
00291                 }
00292           break;
00293         }
00294 
00295       case 2:
00296         {
00297           const std::vector<std::vector<OutputShape> >& d2phidxi2 = fe.get_d2phidxi2();
00298           const std::vector<std::vector<OutputShape> >& d2phidxideta = fe.get_d2phidxideta();
00299           const std::vector<std::vector<OutputShape> >& d2phideta2 = fe.get_d2phideta2();
00300 
00301           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00302           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00303 #if LIBMESH_DIM > 2
00304           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00305 #endif
00306           
00307           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00308           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00309 #if LIBMESH_DIM > 2
00310           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00311 #endif
00312           
00313           for (unsigned int i=0; i<d2phi.size(); i++)
00314             for (unsigned int p=0; p<d2phi[i].size(); p++)
00315               {
00316                 d2phi[i][p].slice(0).slice(0) = d2phidx2[i][p] =
00317                   d2phidxi2[i][p]*dxidx_map[p]*dxidx_map[p] +
00318                   2*d2phidxideta[i][p]*dxidx_map[p]*detadx_map[p] +
00319                   d2phideta2[i][p]*detadx_map[p]*detadx_map[p];
00320 
00321                 d2phi[i][p].slice(0).slice(1) =
00322                   d2phi[i][p].slice(1).slice(0) = d2phidxdy[i][p] =
00323                   d2phidxi2[i][p]*dxidx_map[p]*dxidy_map[p] +
00324                   d2phidxideta[i][p]*dxidx_map[p]*detady_map[p] +
00325                   d2phideta2[i][p]*detadx_map[p]*detady_map[p] +
00326                   d2phidxideta[i][p]*detadx_map[p]*dxidy_map[p];
00327 
00328                 d2phi[i][p].slice(1).slice(1) = d2phidy2[i][p] =
00329                   d2phidxi2[i][p]*dxidy_map[p]*dxidy_map[p] +
00330                   2*d2phidxideta[i][p]*dxidy_map[p]*detady_map[p] +
00331                   d2phideta2[i][p]*detady_map[p]*detady_map[p];
00332 
00333 #if LIBMESH_DIM > 2
00334                 d2phi[i][p].slice(0).slice(2) =
00335                   d2phi[i][p].slice(2).slice(0) = d2phidxdz[i][p] =
00336                   d2phidxi2[i][p]*dxidx_map[p]*dxidz_map[p] +
00337                   d2phidxideta[i][p]*dxidx_map[p]*detadz_map[p] +
00338                   d2phideta2[i][p]*detadx_map[p]*detadz_map[p] +
00339                   d2phidxideta[i][p]*detadx_map[p]*dxidz_map[p];
00340 
00341                 d2phi[i][p].slice(1).slice(2) =
00342                   d2phi[i][p].slice(2).slice(1) = d2phidydz[i][p] =
00343                   d2phidxi2[i][p]*dxidy_map[p]*dxidz_map[p] +
00344                   d2phidxideta[i][p]*dxidy_map[p]*detadz_map[p] +
00345                   d2phideta2[i][p]*detady_map[p]*detadz_map[p] +
00346                   d2phidxideta[i][p]*detady_map[p]*dxidz_map[p];
00347 
00348                 d2phi[i][p].slice(2).slice(2) = d2phidz2[i][p] =
00349                   d2phidxi2[i][p]*dxidz_map[p]*dxidz_map[p] +
00350                   2*d2phidxideta[i][p]*dxidz_map[p]*detadz_map[p] +
00351                   d2phideta2[i][p]*detadz_map[p]*detadz_map[p];
00352 #endif
00353               }
00354           
00355           break;
00356         }
00357 
00358       case 3:
00359         {
00360           const std::vector<std::vector<OutputShape> >& d2phidxi2 = fe.get_d2phidxi2();
00361           const std::vector<std::vector<OutputShape> >& d2phidxideta = fe.get_d2phidxideta();
00362           const std::vector<std::vector<OutputShape> >& d2phideta2 = fe.get_d2phideta2();
00363           const std::vector<std::vector<OutputShape> >& d2phidxidzeta = fe.get_d2phidxidzeta();
00364           const std::vector<std::vector<OutputShape> >& d2phidetadzeta = fe.get_d2phidetadzeta();
00365           const std::vector<std::vector<OutputShape> >& d2phidzeta2 = fe.get_d2phidzeta2();
00366 
00367           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00368           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00369           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00370           
00371           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00372           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00373           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00374           
00375           const std::vector<Real>& dzetadx_map = fe.get_fe_map().get_dzetadx();
00376           const std::vector<Real>& dzetady_map = fe.get_fe_map().get_dzetady();
00377           const std::vector<Real>& dzetadz_map = fe.get_fe_map().get_dzetadz();
00378           
00379           for (unsigned int i=0; i<d2phi.size(); i++)
00380             for (unsigned int p=0; p<d2phi[i].size(); p++)
00381               {
00382                 d2phi[i][p].slice(0).slice(0) = d2phidx2[i][p] =
00383                   d2phidxi2[i][p]*dxidx_map[p]*dxidx_map[p] +
00384                   2*d2phidxideta[i][p]*dxidx_map[p]*detadx_map[p] +
00385                   2*d2phidxidzeta[i][p]*dxidx_map[p]*dzetadx_map[p] +
00386                   2*d2phidetadzeta[i][p]*detadx_map[p]*dzetadx_map[p] +
00387                   d2phideta2[i][p]*detadx_map[p]*detadx_map[p] +
00388                   d2phidzeta2[i][p]*dzetadx_map[p]*dzetadx_map[p];
00389 
00390                 d2phi[i][p].slice(0).slice(1) =
00391                   d2phi[i][p].slice(1).slice(0) = d2phidxdy[i][p] =
00392                   d2phidxi2[i][p]*dxidx_map[p]*dxidy_map[p] +
00393                   d2phidxideta[i][p]*dxidx_map[p]*detady_map[p] +
00394                   d2phidxidzeta[i][p]*dxidx_map[p]*dzetady_map[p] +
00395                   d2phideta2[i][p]*detadx_map[p]*detady_map[p] +
00396                   d2phidxideta[i][p]*detadx_map[p]*dxidy_map[p] +
00397                   d2phidetadzeta[i][p]*detadx_map[p]*dzetady_map[p] +
00398                   d2phidzeta2[i][p]*dzetadx_map[p]*dzetady_map[p] +
00399                   d2phidxidzeta[i][p]*dzetadx_map[p]*dxidy_map[p] +
00400                   d2phidetadzeta[i][p]*dzetadx_map[p]*detady_map[p];
00401 
00402                 d2phi[i][p].slice(0).slice(2) =
00403                   d2phi[i][p].slice(2).slice(0) = d2phidy2[i][p] =
00404                   d2phidxi2[i][p]*dxidx_map[p]*dxidz_map[p] +
00405                   d2phidxideta[i][p]*dxidx_map[p]*detadz_map[p] +
00406                   d2phidxidzeta[i][p]*dxidx_map[p]*dzetadz_map[p] +
00407                   d2phideta2[i][p]*detadx_map[p]*detadz_map[p] +
00408                   d2phidxideta[i][p]*detadx_map[p]*dxidz_map[p] +
00409                   d2phidetadzeta[i][p]*detadx_map[p]*dzetadz_map[p] +
00410                   d2phidzeta2[i][p]*dzetadx_map[p]*dzetadz_map[p] +
00411                   d2phidxidzeta[i][p]*dzetadx_map[p]*dxidz_map[p] +
00412                   d2phidetadzeta[i][p]*dzetadx_map[p]*detadz_map[p];
00413 
00414                 d2phi[i][p].slice(1).slice(1) = d2phidxdz[i][p] =
00415                   d2phidxi2[i][p]*dxidy_map[p]*dxidy_map[p] +
00416                   2*d2phidxideta[i][p]*dxidy_map[p]*detady_map[p] +
00417                   2*d2phidxidzeta[i][p]*dxidy_map[p]*dzetady_map[p] +
00418                   2*d2phidetadzeta[i][p]*detady_map[p]*dzetady_map[p] +
00419                   d2phideta2[i][p]*detady_map[p]*detady_map[p] +
00420                   d2phidzeta2[i][p]*dzetady_map[p]*dzetady_map[p];
00421 
00422                 d2phi[i][p].slice(1).slice(2) =
00423                   d2phi[i][p].slice(2).slice(1) = d2phidydz[i][p] =
00424                   d2phidxi2[i][p]*dxidy_map[p]*dxidz_map[p] +
00425                   d2phidxideta[i][p]*dxidy_map[p]*detadz_map[p] +
00426                   d2phidxidzeta[i][p]*dxidy_map[p]*dzetadz_map[p] +
00427                   d2phideta2[i][p]*detady_map[p]*detadz_map[p] +
00428                   d2phidxideta[i][p]*detady_map[p]*dxidz_map[p] +
00429                   d2phidetadzeta[i][p]*detady_map[p]*dzetadz_map[p] +
00430                   d2phidzeta2[i][p]*dzetady_map[p]*dzetadz_map[p] +
00431                   d2phidxidzeta[i][p]*dzetady_map[p]*dxidz_map[p] +
00432                   d2phidetadzeta[i][p]*dzetady_map[p]*detadz_map[p];
00433 
00434                 d2phi[i][p].slice(2).slice(2) = d2phidz2[i][p] =
00435                   d2phidxi2[i][p]*dxidz_map[p]*dxidz_map[p] +
00436                   2*d2phidxideta[i][p]*dxidz_map[p]*detadz_map[p] +
00437                   2*d2phidxidzeta[i][p]*dxidz_map[p]*dzetadz_map[p] +
00438                   2*d2phidetadzeta[i][p]*detadz_map[p]*dzetadz_map[p] +
00439                   d2phideta2[i][p]*detadz_map[p]*detadz_map[p] +
00440                   d2phidzeta2[i][p]*dzetadz_map[p]*dzetadz_map[p];
00441               }
00442 
00443           break;
00444         }
00445 
00446       default:
00447         libmesh_error();
00448 
00449       } // switch(dim)
00450 
00451     return;
00452   }

template<>
void libMesh::H1FETransformation< RealGradient >::map_div ( const unsigned int  dim,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< RealGradient > &  fe,
std::vector< std::vector< FEGenericBase< RealGradient >::OutputDivergence > > &  div_phi 
) const [inline]

Definition at line 611 of file h1_fe_transformation.C.

References libMesh::FEMap::get_detadx(), libMesh::FEMap::get_detady(), libMesh::FEMap::get_detadz(), libMesh::FEGenericBase< OutputType >::get_dphideta(), libMesh::FEGenericBase< OutputType >::get_dphidxi(), libMesh::FEGenericBase< OutputType >::get_dphidzeta(), libMesh::FEMap::get_dxidx(), libMesh::FEMap::get_dxidy(), libMesh::FEMap::get_dxidz(), libMesh::FEMap::get_dzetadx(), libMesh::FEMap::get_dzetady(), libMesh::FEMap::get_dzetadz(), libMesh::FEAbstract::get_fe_map(), and libMesh::Real.

00616   {
00617     switch(dim)
00618       {
00619         // The divergence only make sense in 2D and 3D
00620       case 0: 
00621       case 1:
00622         {
00623           libmesh_error();
00624         }
00625       case 2:
00626         {
00627           const std::vector<std::vector<RealGradient> >& dphidxi = fe.get_dphidxi();
00628           const std::vector<std::vector<RealGradient> >& dphideta = fe.get_dphideta();
00629           
00630           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00631           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00632 
00633           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00634           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00635 
00636           /*
00637             In 2D: div = dphi_x/dx + dphi_y/dy
00638            */
00639           for (unsigned int i=0; i<div_phi.size(); i++)
00640             for (unsigned int p=0; p<div_phi[i].size(); p++)
00641               {
00642 
00643                 Real dphix_dx = (dphidxi[i][p].slice(0))*dxidx_map[p]
00644                                 + (dphideta[i][p].slice(0))*detadx_map[p];
00645 
00646                 Real dphiy_dy = (dphidxi[i][p].slice(1))*dxidy_map[p]
00647                                 + (dphideta[i][p].slice(1))*detady_map[p];
00648 
00649                 div_phi[i][p] = dphix_dx + dphiy_dy;
00650               }
00651           break;
00652         }
00653       case 3:
00654         {
00655           const std::vector<std::vector<RealGradient> >& dphidxi = fe.get_dphidxi();
00656           const std::vector<std::vector<RealGradient> >& dphideta = fe.get_dphideta();
00657           const std::vector<std::vector<RealGradient> >& dphidzeta = fe.get_dphidzeta();
00658 
00659           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00660           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00661           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00662           
00663           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00664           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00665           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00666           
00667           const std::vector<Real>& dzetadx_map = fe.get_fe_map().get_dzetadx();
00668           const std::vector<Real>& dzetady_map = fe.get_fe_map().get_dzetady();
00669           const std::vector<Real>& dzetadz_map = fe.get_fe_map().get_dzetadz();
00670 
00671           /*
00672             In 3D: div = dphi_x/dx + dphi_y/dy + dphi_z/dz
00673            */
00674           for (unsigned int i=0; i<div_phi.size(); i++)
00675             for (unsigned int p=0; p<div_phi[i].size(); p++)
00676               {
00677                 Real dphix_dx = (dphidxi[i][p].slice(0))*dxidx_map[p]
00678                                 + (dphideta[i][p].slice(0))*detadx_map[p]
00679                                 + (dphidzeta[i][p].slice(0))*dzetadx_map[p];
00680 
00681                 Real dphiy_dy = (dphidxi[i][p].slice(1))*dxidy_map[p]
00682                                 + (dphideta[i][p].slice(1))*detady_map[p]
00683                                 + (dphidzeta[i][p].slice(1))*dzetady_map[p];
00684 
00685                 Real dphiz_dz = (dphidxi[i][p].slice(2))*dxidz_map[p]
00686                                 + (dphideta[i][p].slice(2))*detadz_map[p]
00687                                 + (dphidzeta[i][p].slice(2))*dzetadz_map[p];
00688 
00689                 div_phi[i][p] = dphix_dx + dphiy_dy + dphiz_dz;
00690               }
00691 
00692           break;
00693         }
00694       } // switch(dim)
00695 
00696     return;
00697   }

template<>
void libMesh::H1FETransformation< Real >::map_div ( const unsigned int  ,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< Real > &  ,
std::vector< std::vector< FEGenericBase< Real >::OutputDivergence > > &   
) const [inline]

Definition at line 597 of file h1_fe_transformation.C.

References libMesh::err.

00602   {
00603     libMesh::err << "Computing the divergence of a shape function only\n"
00604                  << "makes sense for vector-valued elements." << std::endl;
00605     libmesh_error();
00606   }

template<typename OutputShape >
virtual void libMesh::H1FETransformation< OutputShape >::map_div ( const unsigned int  dim,
const Elem *const   elem,
const std::vector< Point > &  qp,
const FEGenericBase< OutputShape > &  fe,
std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputDivergence > > &  div_phi 
) const [virtual]

Evaluates the shape function divergence in physical coordinates based on H1 conforming finite element transformation.

Implements libMesh::FETransformationBase< OutputShape >.

template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::map_dphi ( const unsigned int  dim,
const Elem *const  elem,
const std::vector< Point > &  qp,
const FEGenericBase< OutputShape > &  fe,
std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputGradient > > &  dphi,
std::vector< std::vector< OutputShape > > &  dphidx,
std::vector< std::vector< OutputShape > > &  dphidy,
std::vector< std::vector< OutputShape > > &  dphidz 
) const [inline, virtual]

Evaluates shape function gradients in physical coordinates for H1 conforming elements. dphi/dx = dphi/dxi * dxi/dx, etc.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 82 of file h1_fe_transformation.C.

References libMesh::FEMap::get_detadx(), libMesh::FEMap::get_detady(), libMesh::FEMap::get_detadz(), libMesh::FEGenericBase< OutputType >::get_dphideta(), libMesh::FEGenericBase< OutputType >::get_dphidxi(), libMesh::FEGenericBase< OutputType >::get_dphidzeta(), libMesh::FEMap::get_dxidx(), libMesh::FEMap::get_dxidy(), libMesh::FEMap::get_dxidz(), libMesh::FEMap::get_dzetadx(), libMesh::FEMap::get_dzetady(), libMesh::FEMap::get_dzetadz(), and libMesh::FEAbstract::get_fe_map().

00090   {
00091     switch(dim)
00092       {
00093       case 0: // No derivatives in 0D
00094         {
00095           for (unsigned int i=0; i<dphi.size(); i++)
00096             for (unsigned int p=0; p<dphi[i].size(); p++)
00097               {
00098                 dphi[i][p] = 0.;
00099               }
00100           break;
00101         }
00102 
00103       case 1:
00104         {
00105           const std::vector<std::vector<OutputShape> >& dphidxi = fe.get_dphidxi();
00106 
00107           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00108 #if LIBMESH_DIM>1
00109           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00110 #endif
00111 #if LIBMESH_DIM>2
00112           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00113 #endif
00114           
00115           for (unsigned int i=0; i<dphi.size(); i++)
00116             for (unsigned int p=0; p<dphi[i].size(); p++)
00117               {
00118                 // dphi/dx    = (dphi/dxi)*(dxi/dx)
00119                 dphi[i][p].slice(0) = dphidx[i][p] = dphidxi[i][p]*dxidx_map[p];
00120                 
00121 #if LIBMESH_DIM>1
00122                 dphi[i][p].slice(1)  = dphidy[i][p] = dphidxi[i][p]*dxidy_map[p];
00123 #endif
00124 #if LIBMESH_DIM>2
00125                 dphi[i][p].slice(2) = dphidz[i][p] = dphidxi[i][p]*dxidz_map[p];
00126 #endif
00127               }
00128           
00129           break;
00130         }
00131         
00132       case 2:
00133         {
00134           const std::vector<std::vector<OutputShape> >& dphidxi = fe.get_dphidxi();
00135           const std::vector<std::vector<OutputShape> >& dphideta = fe.get_dphideta();
00136 
00137           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00138           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00139 #if LIBMESH_DIM > 2
00140           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00141 #endif
00142 
00143           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00144           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00145 #if LIBMESH_DIM > 2
00146           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00147 #endif
00148 
00149           for (unsigned int i=0; i<dphi.size(); i++)
00150             for (unsigned int p=0; p<dphi[i].size(); p++)
00151               {
00152                 // dphi/dx    = (dphi/dxi)*(dxi/dx) + (dphi/deta)*(deta/dx)
00153                 dphi[i][p].slice(0) = dphidx[i][p] = (dphidxi[i][p]*dxidx_map[p] +
00154                                                       dphideta[i][p]*detadx_map[p]);
00155                 
00156                 // dphi/dy    = (dphi/dxi)*(dxi/dy) + (dphi/deta)*(deta/dy)
00157                 dphi[i][p].slice(1) = dphidy[i][p] = (dphidxi[i][p]*dxidy_map[p] +
00158                                                       dphideta[i][p]*detady_map[p]);
00159 
00160 #if LIBMESH_DIM > 2
00161                 // dphi/dz    = (dphi/dxi)*(dxi/dz) + (dphi/deta)*(deta/dz)
00162                 dphi[i][p].slice(2) = dphidz[i][p] = (dphidxi[i][p]*dxidz_map[p] +
00163                                                       dphideta[i][p]*detadz_map[p]);
00164 #endif
00165               }
00166 
00167           break;
00168         }
00169 
00170       case 3: 
00171         {
00172           const std::vector<std::vector<OutputShape> >& dphidxi = fe.get_dphidxi();
00173           const std::vector<std::vector<OutputShape> >& dphideta = fe.get_dphideta();
00174           const std::vector<std::vector<OutputShape> >& dphidzeta = fe.get_dphidzeta();
00175 
00176           const std::vector<Real>& dxidx_map = fe.get_fe_map().get_dxidx();
00177           const std::vector<Real>& dxidy_map = fe.get_fe_map().get_dxidy();
00178           const std::vector<Real>& dxidz_map = fe.get_fe_map().get_dxidz();
00179           
00180           const std::vector<Real>& detadx_map = fe.get_fe_map().get_detadx();
00181           const std::vector<Real>& detady_map = fe.get_fe_map().get_detady();
00182           const std::vector<Real>& detadz_map = fe.get_fe_map().get_detadz();
00183           
00184           const std::vector<Real>& dzetadx_map = fe.get_fe_map().get_dzetadx();
00185           const std::vector<Real>& dzetady_map = fe.get_fe_map().get_dzetady();
00186           const std::vector<Real>& dzetadz_map = fe.get_fe_map().get_dzetadz();
00187           
00188           for (unsigned int i=0; i<dphi.size(); i++)
00189             for (unsigned int p=0; p<dphi[i].size(); p++)
00190               {
00191                 // dphi/dx    = (dphi/dxi)*(dxi/dx) + (dphi/deta)*(deta/dx) + (dphi/dzeta)*(dzeta/dx);
00192                 dphi[i][p].slice(0) = dphidx[i][p] = (dphidxi[i][p]*dxidx_map[p] +
00193                                                       dphideta[i][p]*detadx_map[p] +
00194                                                       dphidzeta[i][p]*dzetadx_map[p]);
00195                   
00196                 // dphi/dy    = (dphi/dxi)*(dxi/dy) + (dphi/deta)*(deta/dy) + (dphi/dzeta)*(dzeta/dy);
00197                 dphi[i][p].slice(1) = dphidy[i][p] = (dphidxi[i][p]*dxidy_map[p] +
00198                                                       dphideta[i][p]*detady_map[p] +
00199                                                       dphidzeta[i][p]*dzetady_map[p]);
00200                 
00201                 // dphi/dz    = (dphi/dxi)*(dxi/dz) + (dphi/deta)*(deta/dz) + (dphi/dzeta)*(dzeta/dz);
00202                 dphi[i][p].slice(2) = dphidz[i][p] = (dphidxi[i][p]*dxidz_map[p] +
00203                                                       dphideta[i][p]*detadz_map[p] +
00204                                                       dphidzeta[i][p]*dzetadz_map[p]);
00205               }
00206           break;
00207         }
00208 
00209       default:
00210         libmesh_error();
00211 
00212       } // switch(dim)
00213 
00214     return;
00215   }

template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::map_phi ( const unsigned int  dim,
const Elem * const  elem,
const std::vector< Point > &  qp,
const FEGenericBase< OutputShape > &  fe,
std::vector< std::vector< OutputShape > > &  phi 
) const [inline, virtual]

Evaluates shape functions in physical coordinates for H1 conforming elements. In this case $ \phi(x) = \phi(\xi) $

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 25 of file h1_fe_transformation.C.

References libMesh::FEAbstract::get_fe_type().

00030   {
00031     switch(dim)
00032       {
00033       case 0:
00034         {
00035           for (unsigned int i=0; i<phi.size(); i++)
00036             {
00037               libmesh_assert_equal_to ( qp.size(), phi[i].size() );
00038               for (unsigned int p=0; p<phi[i].size(); p++)
00039                 FEInterface::shape<OutputShape>(0, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
00040             }
00041           break;
00042         }
00043       case 1:
00044         {
00045           for (unsigned int i=0; i<phi.size(); i++)
00046             {
00047               libmesh_assert_equal_to ( qp.size(), phi[i].size() );
00048               for (unsigned int p=0; p<phi[i].size(); p++)
00049                 FEInterface::shape<OutputShape>(1, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
00050             }
00051           break;
00052         }
00053       case 2:
00054         {
00055           for (unsigned int i=0; i<phi.size(); i++)
00056             {
00057               libmesh_assert_equal_to ( qp.size(), phi[i].size() );
00058               for (unsigned int p=0; p<phi[i].size(); p++)
00059                 FEInterface::shape<OutputShape>(2, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
00060             }
00061           break;
00062         }
00063       case 3:
00064         {
00065           for (unsigned int i=0; i<phi.size(); i++)
00066             {
00067               libmesh_assert_equal_to ( qp.size(), phi[i].size() );
00068               for (unsigned int p=0; p<phi[i].size(); p++)
00069                 FEInterface::shape<OutputShape>(3, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
00070             }
00071           break;
00072         }
00073       default:
00074         libmesh_error();
00075       }
00076 
00077     return;
00078   }


The documentation for this class was generated from the following files:

Site Created By: libMesh Developers
Last modified: February 05 2013 19:55:26 UTC

Hosted By:
SourceForge.net Logo