libMesh::H1FETransformation< OutputShape > Class Template Reference

#include <fe_transformation_base.h>

Inheritance diagram for libMesh::H1FETransformation< OutputShape >:

Public Member Functions

 H1FETransformation ()
 
virtual ~H1FETransformation ()
 
virtual void init_map_phi (const FEGenericBase< OutputShape > &fe) const libmesh_override
 
virtual void init_map_dphi (const FEGenericBase< OutputShape > &fe) const libmesh_override
 
virtual void init_map_d2phi (const FEGenericBase< OutputShape > &fe) const libmesh_override
 
virtual void map_phi (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< OutputShape > &, std::vector< std::vector< OutputShape > > &) const libmesh_override
 
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 libmesh_override
 
virtual void map_d2phi (const unsigned int dim, 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 libmesh_override
 
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 libmesh_override
 
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 libmesh_override
 
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 UniquePtr< 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
Date
2012

Definition at line 28 of file fe_transformation_base.h.

Constructor & Destructor Documentation

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

Definition at line 44 of file h1_fe_transformation.h.

45  : FETransformationBase<OutputShape>(){}

Member Function Documentation

template<typename OutputShape >
UniquePtr< FETransformationBase< OutputShape > > libMesh::FETransformationBase< OutputShape >::build ( const FEType type)
staticinherited

Builds an FETransformation object based on the finite element type

Definition at line 27 of file fe_transformation_base.C.

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

Referenced by libMesh::FETransformationBase< OutputShape >::~FETransformationBase().

28 {
29  switch (fe_type.family)
30  {
31  // H1 Conforming Elements
32  case LAGRANGE:
33  case HIERARCHIC:
34  case BERNSTEIN:
35  case SZABAB:
36  case CLOUGH: // PB: Really H2
37  case HERMITE: // PB: Really H2
38  case SUBDIVISION:
39  case LAGRANGE_VEC:
40  case MONOMIAL: // PB: Shouldn't this be L2 conforming?
41  case XYZ: // PB: Shouldn't this be L2 conforming?
42  case L2_HIERARCHIC: // PB: Shouldn't this be L2 conforming?
43  case L2_LAGRANGE: // PB: Shouldn't this be L2 conforming?
44  case JACOBI_20_00: // PB: For infinite elements...
45  case JACOBI_30_00: // PB: For infinite elements...
46  return UniquePtr<FETransformationBase<OutputShape> >(new H1FETransformation<OutputShape>);
47 
48  // HCurl Conforming Elements
49  case NEDELEC_ONE:
50  return UniquePtr<FETransformationBase<OutputShape> >(new HCurlFETransformation<OutputShape>);
51 
52  // HDiv Conforming Elements
53  // L2 Conforming Elements
54 
55  // Other...
56  case SCALAR:
57  // Should never need this for SCALARs
58  return UniquePtr<FETransformationBase<OutputShape> >(new H1FETransformation<OutputShape>);
59 
60  default:
61  libmesh_error_msg("Unknown family = " << fe_type.family);
62  }
63 
64  libmesh_error_msg("We'll never get here!");
65  return UniquePtr<FETransformationBase<OutputShape> >();
66 }
template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::init_map_d2phi ( const FEGenericBase< OutputShape > &  fe) const
virtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 42 of file h1_fe_transformation.C.

References libMesh::FEMap::get_d2xidxyz2(), libMesh::FEMap::get_dxidx(), and libMesh::FEAbstract::get_fe_map().

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

43 {
44  fe.get_fe_map().get_dxidx();
45 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
46  fe.get_fe_map().get_d2xidxyz2();
47 #endif
48 }
template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::init_map_dphi ( const FEGenericBase< OutputShape > &  fe) const
virtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 34 of file h1_fe_transformation.C.

References libMesh::FEMap::get_dxidx(), and libMesh::FEAbstract::get_fe_map().

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

35 {
36  fe.get_fe_map().get_dxidx();
37 }
template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::init_map_phi ( const FEGenericBase< OutputShape > &  fe) const
virtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 26 of file h1_fe_transformation.C.

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

27 {
28  // Nothing needed
29 }
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 >.

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

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

Definition at line 526 of file h1_fe_transformation.C.

531 {
532  libmesh_error_msg("Computing the curl of a shape function only \nmakes sense for vector-valued elements.");
533 }
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

Definition at line 536 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.

541 {
542  switch (dim)
543  {
544  case 0:
545  case 1:
546  libmesh_error_msg("The curl only make sense in 2D and 3D");
547 
548  case 2:
549  {
550  const std::vector<std::vector<RealGradient> > & dphidxi = fe.get_dphidxi();
551  const std::vector<std::vector<RealGradient> > & dphideta = fe.get_dphideta();
552 
553  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
554  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
555 #if LIBMESH_DIM > 2
556  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
557 #endif
558 
559  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
560  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
561 #if LIBMESH_DIM > 2
562  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
563 #endif
564 
565  /*
566  For 2D elements in 3D space:
567 
568  curl = ( -dphi_y/dz, dphi_x/dz, dphi_y/dx - dphi_x/dy )
569  */
570  for (std::size_t i=0; i<curl_phi.size(); i++)
571  for (std::size_t p=0; p<curl_phi[i].size(); p++)
572  {
573 
574  Real dphiy_dx = (dphidxi[i][p].slice(1))*dxidx_map[p]
575  + (dphideta[i][p].slice(1))*detadx_map[p];
576 
577  Real dphix_dy = (dphidxi[i][p].slice(0))*dxidy_map[p]
578  + (dphideta[i][p].slice(0))*detady_map[p];
579 
580  curl_phi[i][p].slice(2) = dphiy_dx - dphix_dy;
581 
582 #if LIBMESH_DIM > 2
583  Real dphiy_dz = (dphidxi[i][p].slice(1))*dxidz_map[p]
584  + (dphideta[i][p].slice(1))*detadz_map[p];
585 
586  Real dphix_dz = (dphidxi[i][p].slice(0))*dxidz_map[p]
587  + (dphideta[i][p].slice(0))*detadz_map[p];
588 
589  curl_phi[i][p].slice(0) = -dphiy_dz;
590  curl_phi[i][p].slice(1) = dphix_dz;
591 #endif
592  }
593 
594  break;
595  }
596  case 3:
597  {
598  const std::vector<std::vector<RealGradient> > & dphidxi = fe.get_dphidxi();
599  const std::vector<std::vector<RealGradient> > & dphideta = fe.get_dphideta();
600  const std::vector<std::vector<RealGradient> > & dphidzeta = fe.get_dphidzeta();
601 
602  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
603  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
604  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
605 
606  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
607  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
608  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
609 
610  const std::vector<Real> & dzetadx_map = fe.get_fe_map().get_dzetadx();
611  const std::vector<Real> & dzetady_map = fe.get_fe_map().get_dzetady();
612  const std::vector<Real> & dzetadz_map = fe.get_fe_map().get_dzetadz();
613 
614  /*
615  In 3D: curl = ( dphi_z/dy - dphi_y/dz, dphi_x/dz - dphi_z/dx, dphi_y/dx - dphi_x/dy )
616  */
617  for (std::size_t i=0; i<curl_phi.size(); i++)
618  for (std::size_t p=0; p<curl_phi[i].size(); p++)
619  {
620  Real dphiz_dy = (dphidxi[i][p].slice(2))*dxidy_map[p]
621  + (dphideta[i][p].slice(2))*detady_map[p]
622  + (dphidzeta[i][p].slice(2))*dzetady_map[p];
623 
624  Real dphiy_dz = (dphidxi[i][p].slice(1))*dxidz_map[p]
625  + (dphideta[i][p].slice(1))*detadz_map[p]
626  + (dphidzeta[i][p].slice(1))*dzetadz_map[p];
627 
628  Real dphix_dz = (dphidxi[i][p].slice(0))*dxidz_map[p]
629  + (dphideta[i][p].slice(0))*detadz_map[p]
630  + (dphidzeta[i][p].slice(0))*dzetadz_map[p];
631 
632  Real dphiz_dx = (dphidxi[i][p].slice(2))*dxidx_map[p]
633  + (dphideta[i][p].slice(2))*detadx_map[p]
634  + (dphidzeta[i][p].slice(2))*dzetadx_map[p];
635 
636  Real dphiy_dx = (dphidxi[i][p].slice(1))*dxidx_map[p]
637  + (dphideta[i][p].slice(1))*detadx_map[p]
638  + (dphidzeta[i][p].slice(1))*dzetadx_map[p];
639 
640  Real dphix_dy = (dphidxi[i][p].slice(0))*dxidy_map[p]
641  + (dphideta[i][p].slice(0))*detady_map[p]
642  + (dphidzeta[i][p].slice(0))*dzetady_map[p];
643 
644  curl_phi[i][p].slice(0) = dphiz_dy - dphiy_dz;
645 
646  curl_phi[i][p].slice(1) = dphix_dz - dphiz_dx;
647 
648  curl_phi[i][p].slice(2) = dphiy_dx - dphix_dy;
649  }
650 
651  break;
652  }
653  default:
654  libmesh_error_msg("Invalid dim = " << dim);
655  }
656 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::map_d2phi ( const unsigned int  dim,
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

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

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 239 of file h1_fe_transformation.C.

References libMesh::FEMap::get_d2etadxyz2(), 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_d2xidxyz2(), libMesh::FEMap::get_d2zetadxyz2(), 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().

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

249 {
250  switch (dim)
251  {
252  case 0: // No derivatives in 0D
253  {
254  for (std::size_t i=0; i<d2phi.size(); i++)
255  for (std::size_t p=0; p<d2phi[i].size(); p++)
256  d2phi[i][p] = 0.;
257 
258  break;
259  }
260 
261  case 1:
262  {
263  const std::vector<std::vector<OutputShape> > & d2phidxi2 = fe.get_d2phidxi2();
264 
265  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
266 #if LIBMESH_DIM>1
267  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
268 #endif
269 #if LIBMESH_DIM>2
270  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
271 #endif
272 
273  // Shape function derivatives in reference space
274  const std::vector<std::vector<OutputShape> > & dphidxi = fe.get_dphidxi();
275 
276  // Inverse map second derivatives
277  const std::vector<std::vector<Real> > & d2xidxyz2 = fe.get_fe_map().get_d2xidxyz2();
278 
279  for (std::size_t i=0; i<d2phi.size(); i++)
280  for (std::size_t p=0; p<d2phi[i].size(); p++)
281  {
282  // phi_{x x}
283  d2phi[i][p].slice(0).slice(0) = d2phidx2[i][p] =
284  d2phidxi2[i][p]*dxidx_map[p]*dxidx_map[p] + // (xi_x)^2 * phi_{xi xi}
285  d2xidxyz2[p][0]*dphidxi[i][p]; // xi_{x x} * phi_{xi}
286 
287 #if LIBMESH_DIM>1
288  // phi_{x y}
289  d2phi[i][p].slice(0).slice(1) = d2phi[i][p].slice(1).slice(0) = d2phidxdy[i][p] =
290  d2phidxi2[i][p]*dxidx_map[p]*dxidy_map[p] + // xi_x * xi_y * phi_{xi xi}
291  d2xidxyz2[p][1]*dphidxi[i][p]; // xi_{x y} * phi_{xi}
292 #endif
293 
294 #if LIBMESH_DIM>2
295  // phi_{x z}
296  d2phi[i][p].slice(0).slice(2) = d2phi[i][p].slice(2).slice(0) = d2phidxdz[i][p] =
297  d2phidxi2[i][p]*dxidx_map[p]*dxidz_map[p] + // xi_x * xi_z * phi_{xi xi}
298  d2xidxyz2[p][2]*dphidxi[i][p]; // xi_{x z} * phi_{xi}
299 #endif
300 
301 
302 #if LIBMESH_DIM>1
303  // phi_{y y}
304  d2phi[i][p].slice(1).slice(1) = d2phidy2[i][p] =
305  d2phidxi2[i][p]*dxidy_map[p]*dxidy_map[p] + // (xi_y)^2 * phi_{xi xi}
306  d2xidxyz2[p][3]*dphidxi[i][p]; // xi_{y y} * phi_{xi}
307 #endif
308 
309 #if LIBMESH_DIM>2
310  // phi_{y z}
311  d2phi[i][p].slice(1).slice(2) = d2phi[i][p].slice(2).slice(1) = d2phidydz[i][p] =
312  d2phidxi2[i][p]*dxidy_map[p]*dxidz_map[p] + // xi_y * xi_z * phi_{xi xi}
313  d2xidxyz2[p][4]*dphidxi[i][p]; // xi_{y z} * phi_{xi}
314 
315  // phi_{z z}
316  d2phi[i][p].slice(2).slice(2) = d2phidz2[i][p] =
317  d2phidxi2[i][p]*dxidz_map[p]*dxidz_map[p] + // (xi_z)^2 * phi_{xi xi}
318  d2xidxyz2[p][5]*dphidxi[i][p]; // xi_{z z} * phi_{xi}
319 #endif
320  }
321  break;
322  }
323 
324  case 2:
325  {
326  const std::vector<std::vector<OutputShape> > & d2phidxi2 = fe.get_d2phidxi2();
327  const std::vector<std::vector<OutputShape> > & d2phidxideta = fe.get_d2phidxideta();
328  const std::vector<std::vector<OutputShape> > & d2phideta2 = fe.get_d2phideta2();
329 
330  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
331  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
332 #if LIBMESH_DIM > 2
333  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
334 #endif
335 
336  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
337  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
338 #if LIBMESH_DIM > 2
339  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
340 #endif
341 
342  // Shape function derivatives in reference space
343  const std::vector<std::vector<OutputShape> > & dphidxi = fe.get_dphidxi();
344  const std::vector<std::vector<OutputShape> > & dphideta = fe.get_dphideta();
345 
346  // Inverse map second derivatives
347  const std::vector<std::vector<Real> > & d2xidxyz2 = fe.get_fe_map().get_d2xidxyz2();
348  const std::vector<std::vector<Real> > & d2etadxyz2 = fe.get_fe_map().get_d2etadxyz2();
349 
350  for (std::size_t i=0; i<d2phi.size(); i++)
351  for (std::size_t p=0; p<d2phi[i].size(); p++)
352  {
353  // phi_{x x}
354  d2phi[i][p].slice(0).slice(0) = d2phidx2[i][p] =
355  d2phidxi2[i][p]*dxidx_map[p]*dxidx_map[p] + // (xi_x)^2 * phi_{xi xi}
356  d2phideta2[i][p]*detadx_map[p]*detadx_map[p] + // (eta_x)^2 * phi_{eta eta}
357  2*d2phidxideta[i][p]*dxidx_map[p]*detadx_map[p] + // 2 * xi_x * eta_x * phi_{xi eta}
358  d2xidxyz2[p][0]*dphidxi[i][p] + // xi_{x x} * phi_{xi}
359  d2etadxyz2[p][0]*dphideta[i][p]; // eta_{x x} * phi_{eta}
360 
361  // phi_{x y}
362  d2phi[i][p].slice(0).slice(1) = d2phi[i][p].slice(1).slice(0) = d2phidxdy[i][p] =
363  d2phidxi2[i][p]*dxidx_map[p]*dxidy_map[p] + // xi_x * xi_y * phi_{xi xi}
364  d2phideta2[i][p]*detadx_map[p]*detady_map[p] + // eta_x * eta_y * phi_{eta eta}
365  d2phidxideta[i][p]*(dxidx_map[p]*detady_map[p] + detadx_map[p]*dxidy_map[p]) + // (xi_x*eta_y + eta_x*xi_y) * phi_{xi eta}
366  d2xidxyz2[p][1]*dphidxi[i][p] + // xi_{x y} * phi_{xi}
367  d2etadxyz2[p][1]*dphideta[i][p]; // eta_{x y} * phi_{eta}
368 
369 #if LIBMESH_DIM > 2
370  // phi_{x z}
371  d2phi[i][p].slice(0).slice(2) = d2phi[i][p].slice(2).slice(0) = d2phidxdz[i][p] =
372  d2phidxi2[i][p]*dxidx_map[p]*dxidz_map[p] + // xi_x * xi_z * phi_{xi xi}
373  d2phideta2[i][p]*detadx_map[p]*detadz_map[p] + // eta_x * eta_z * phi_{eta eta}
374  d2phidxideta[i][p]*(dxidx_map[p]*detadz_map[p] + detadx_map[p]*dxidz_map[p]) + // (xi_x*eta_z + eta_x*xi_z) * phi_{xi eta}
375  d2xidxyz2[p][2]*dphidxi[i][p] + // xi_{x z} * phi_{xi}
376  d2etadxyz2[p][2]*dphideta[i][p]; // eta_{x z} * phi_{eta}
377 #endif
378 
379  // phi_{y y}
380  d2phi[i][p].slice(1).slice(1) = d2phidy2[i][p] =
381  d2phidxi2[i][p]*dxidy_map[p]*dxidy_map[p] + // (xi_y)^2 * phi_{xi xi}
382  d2phideta2[i][p]*detady_map[p]*detady_map[p] + // (eta_y)^2 * phi_{eta eta}
383  2*d2phidxideta[i][p]*dxidy_map[p]*detady_map[p] + // 2 * xi_y * eta_y * phi_{xi eta}
384  d2xidxyz2[p][3]*dphidxi[i][p] + // xi_{y y} * phi_{xi}
385  d2etadxyz2[p][3]*dphideta[i][p]; // eta_{y y} * phi_{eta}
386 
387 #if LIBMESH_DIM > 2
388  // phi_{y z}
389  d2phi[i][p].slice(1).slice(2) = d2phi[i][p].slice(2).slice(1) = d2phidydz[i][p] =
390  d2phidxi2[i][p]*dxidy_map[p]*dxidz_map[p] + // xi_y * xi_z * phi_{xi xi}
391  d2phideta2[i][p]*detady_map[p]*detadz_map[p] + // eta_y * eta_z * phi_{eta eta}
392  d2phidxideta[i][p]*(dxidy_map[p]*detadz_map[p] + detady_map[p]*dxidz_map[p]) + // (xi_y*eta_z + eta_y*xi_z) * phi_{xi eta}
393  d2xidxyz2[p][4]*dphidxi[i][p] + // xi_{y z} * phi_{xi}
394  d2etadxyz2[p][4]*dphideta[i][p]; // eta_{y z} * phi_{eta}
395 
396  // phi_{z z}
397  d2phi[i][p].slice(2).slice(2) = d2phidz2[i][p] =
398  d2phidxi2[i][p]*dxidz_map[p]*dxidz_map[p] + // (xi_z)^2 * phi_{xi xi}
399  d2phideta2[i][p]*detadz_map[p]*detadz_map[p] + // (eta_z)^2 * phi_{eta eta}
400  2*d2phidxideta[i][p]*dxidz_map[p]*detadz_map[p] + // 2 * xi_z * eta_z * phi_{xi eta}
401  d2xidxyz2[p][5]*dphidxi[i][p] + // xi_{z z} * phi_{xi}
402  d2etadxyz2[p][5]*dphideta[i][p]; // eta_{z z} * phi_{eta}
403 #endif
404  }
405 
406  break;
407  }
408 
409  case 3:
410  {
411  const std::vector<std::vector<OutputShape> > & d2phidxi2 = fe.get_d2phidxi2();
412  const std::vector<std::vector<OutputShape> > & d2phidxideta = fe.get_d2phidxideta();
413  const std::vector<std::vector<OutputShape> > & d2phideta2 = fe.get_d2phideta2();
414  const std::vector<std::vector<OutputShape> > & d2phidxidzeta = fe.get_d2phidxidzeta();
415  const std::vector<std::vector<OutputShape> > & d2phidetadzeta = fe.get_d2phidetadzeta();
416  const std::vector<std::vector<OutputShape> > & d2phidzeta2 = fe.get_d2phidzeta2();
417 
418  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
419  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
420  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
421 
422  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
423  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
424  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
425 
426  const std::vector<Real> & dzetadx_map = fe.get_fe_map().get_dzetadx();
427  const std::vector<Real> & dzetady_map = fe.get_fe_map().get_dzetady();
428  const std::vector<Real> & dzetadz_map = fe.get_fe_map().get_dzetadz();
429 
430  // Shape function derivatives in reference space
431  const std::vector<std::vector<OutputShape> > & dphidxi = fe.get_dphidxi();
432  const std::vector<std::vector<OutputShape> > & dphideta = fe.get_dphideta();
433  const std::vector<std::vector<OutputShape> > & dphidzeta = fe.get_dphidzeta();
434 
435  // Inverse map second derivatives
436  const std::vector<std::vector<Real> > & d2xidxyz2 = fe.get_fe_map().get_d2xidxyz2();
437  const std::vector<std::vector<Real> > & d2etadxyz2 = fe.get_fe_map().get_d2etadxyz2();
438  const std::vector<std::vector<Real> > & d2zetadxyz2 = fe.get_fe_map().get_d2zetadxyz2();
439 
440  for (std::size_t i=0; i<d2phi.size(); i++)
441  for (std::size_t p=0; p<d2phi[i].size(); p++)
442  {
443  // phi_{x x}
444  d2phi[i][p].slice(0).slice(0) = d2phidx2[i][p] =
445  d2phidxi2[i][p]*dxidx_map[p]*dxidx_map[p] + // (xi_x)^2 * phi_{xi xi}
446  d2phideta2[i][p]*detadx_map[p]*detadx_map[p] + // (eta_x)^2 * phi_{eta eta}
447  d2phidzeta2[i][p]*dzetadx_map[p]*dzetadx_map[p] + // (zeta_x)^2 * phi_{zeta zeta}
448  2*d2phidxideta[i][p]*dxidx_map[p]*detadx_map[p] + // 2 * xi_x * eta_x * phi_{xi eta}
449  2*d2phidxidzeta[i][p]*dxidx_map[p]*dzetadx_map[p] + // 2 * xi_x * zeta_x * phi_{xi zeta}
450  2*d2phidetadzeta[i][p]*detadx_map[p]*dzetadx_map[p] + // 2 * eta_x * zeta_x * phi_{eta zeta}
451  d2xidxyz2[p][0]*dphidxi[i][p] + // xi_{x x} * phi_{xi}
452  d2etadxyz2[p][0]*dphideta[i][p] + // eta_{x x} * phi_{eta}
453  d2zetadxyz2[p][0]*dphidzeta[i][p]; // zeta_{x x} * phi_{zeta}
454 
455  // phi_{x y}
456  d2phi[i][p].slice(0).slice(1) = d2phi[i][p].slice(1).slice(0) = d2phidxdy[i][p] =
457  d2phidxi2[i][p]*dxidx_map[p]*dxidy_map[p] + // xi_x * xi_y * phi_{xi xi}
458  d2phideta2[i][p]*detadx_map[p]*detady_map[p] + // eta_x * eta_y * phi_{eta eta}
459  d2phidzeta2[i][p]*dzetadx_map[p]*dzetady_map[p] + // zeta_x * zeta_y * phi_{zeta zeta}
460  d2phidxideta[i][p]*(dxidx_map[p]*detady_map[p] + detadx_map[p]*dxidy_map[p]) + // (xi_x*eta_y + eta_x*xi_y) * phi_{xi eta}
461  d2phidxidzeta[i][p]*(dxidx_map[p]*dzetady_map[p] + dzetadx_map[p]*dxidy_map[p]) + // (zeta_x*xi_y + xi_x*zeta_y) * phi_{xi zeta}
462  d2phidetadzeta[i][p]*(detadx_map[p]*dzetady_map[p] + dzetadx_map[p]*detady_map[p]) + // (zeta_x*eta_y + eta_x*zeta_y) * phi_{eta zeta}
463  d2xidxyz2[p][1]*dphidxi[i][p] + // xi_{x y} * phi_{xi}
464  d2etadxyz2[p][1]*dphideta[i][p] + // eta_{x y} * phi_{eta}
465  d2zetadxyz2[p][1]*dphidzeta[i][p]; // zeta_{x y} * phi_{zeta}
466 
467  // phi_{x z}
468  d2phi[i][p].slice(0).slice(2) = d2phi[i][p].slice(2).slice(0) = d2phidy2[i][p] =
469  d2phidxi2[i][p]*dxidx_map[p]*dxidz_map[p] + // xi_x * xi_z * phi_{xi xi}
470  d2phideta2[i][p]*detadx_map[p]*detadz_map[p] + // eta_x * eta_z * phi_{eta eta}
471  d2phidzeta2[i][p]*dzetadx_map[p]*dzetadz_map[p] + // zeta_x * zeta_z * phi_{zeta zeta}
472  d2phidxideta[i][p]*(dxidx_map[p]*detadz_map[p] + detadx_map[p]*dxidz_map[p]) + // (xi_x*eta_z + eta_x*xi_z) * phi_{xi eta}
473  d2phidxidzeta[i][p]*(dxidx_map[p]*dzetadz_map[p] + dzetadx_map[p]*dxidz_map[p]) + // (zeta_x*xi_z + xi_x*zeta_z) * phi_{xi zeta}
474  d2phidetadzeta[i][p]*(detadx_map[p]*dzetadz_map[p] + dzetadx_map[p]*detadz_map[p]) + // (zeta_x*eta_z + eta_x*zeta_z) * phi_{eta zeta}
475  d2xidxyz2[p][2]*dphidxi[i][p] + // xi_{x z} * phi_{xi}
476  d2etadxyz2[p][2]*dphideta[i][p] + // eta_{x z} * phi_{eta}
477  d2zetadxyz2[p][2]*dphidzeta[i][p]; // zeta_{x z} * phi_{zeta}
478 
479  // phi_{y y}
480  d2phi[i][p].slice(1).slice(1) = d2phidxdz[i][p] =
481  d2phidxi2[i][p]*dxidy_map[p]*dxidy_map[p] + // (xi_y)^2 * phi_{xi xi}
482  d2phideta2[i][p]*detady_map[p]*detady_map[p] + // (eta_y)^2 * phi_{eta eta}
483  d2phidzeta2[i][p]*dzetady_map[p]*dzetady_map[p] + // (zeta_y)^2 * phi_{zeta zeta}
484  2*d2phidxideta[i][p]*dxidy_map[p]*detady_map[p] + // 2 * xi_y * eta_y * phi_{xi eta}
485  2*d2phidxidzeta[i][p]*dxidy_map[p]*dzetady_map[p] + // 2 * xi_y * zeta_y * phi_{xi zeta}
486  2*d2phidetadzeta[i][p]*detady_map[p]*dzetady_map[p] + // 2 * eta_y * zeta_y * phi_{eta zeta}
487  d2xidxyz2[p][3]*dphidxi[i][p] + // xi_{y y} * phi_{xi}
488  d2etadxyz2[p][3]*dphideta[i][p] + // eta_{y y} * phi_{eta}
489  d2zetadxyz2[p][3]*dphidzeta[i][p]; // zeta_{y y} * phi_{zeta}
490 
491  // phi_{y z}
492  d2phi[i][p].slice(1).slice(2) = d2phi[i][p].slice(2).slice(1) = d2phidydz[i][p] =
493  d2phidxi2[i][p]*dxidy_map[p]*dxidz_map[p] + // xi_y * xi_z * phi_{xi xi}
494  d2phideta2[i][p]*detady_map[p]*detadz_map[p] + // eta_y * eta_z * phi_{eta eta}
495  d2phidzeta2[i][p]*dzetady_map[p]*dzetadz_map[p] + // zeta_y * zeta_z * phi_{zeta zeta}
496  d2phidxideta[i][p]*(dxidy_map[p]*detadz_map[p] + detady_map[p]*dxidz_map[p]) + // (xi_y*eta_z + eta_y*xi_z) * phi_{xi eta}
497  d2phidxidzeta[i][p]*(dxidy_map[p]*dzetadz_map[p] + dzetady_map[p]*dxidz_map[p]) + // (zeta_y*xi_z + xi_y*zeta_z) * phi_{xi zeta}
498  d2phidetadzeta[i][p]*(detady_map[p]*dzetadz_map[p] + dzetady_map[p]*detadz_map[p]) + // (zeta_y*eta_z + eta_y*zeta_z) * phi_{eta zeta}
499  d2xidxyz2[p][4]*dphidxi[i][p] + // xi_{y z} * phi_{xi}
500  d2etadxyz2[p][4]*dphideta[i][p] + // eta_{y z} * phi_{eta}
501  d2zetadxyz2[p][4]*dphidzeta[i][p]; // zeta_{y z} * phi_{zeta}
502 
503  // phi_{z z}
504  d2phi[i][p].slice(2).slice(2) = d2phidz2[i][p] =
505  d2phidxi2[i][p]*dxidz_map[p]*dxidz_map[p] + // (xi_z)^2 * phi_{xi xi}
506  d2phideta2[i][p]*detadz_map[p]*detadz_map[p] + // (eta_z)^2 * phi_{eta eta}
507  d2phidzeta2[i][p]*dzetadz_map[p]*dzetadz_map[p] + // (zeta_z)^2 * phi_{zeta zeta}
508  2*d2phidxideta[i][p]*dxidz_map[p]*detadz_map[p] + // 2 * xi_z * eta_z * phi_{xi eta}
509  2*d2phidxidzeta[i][p]*dxidz_map[p]*dzetadz_map[p] + // 2 * xi_z * zeta_z * phi_{xi zeta}
510  2*d2phidetadzeta[i][p]*detadz_map[p]*dzetadz_map[p] + // 2 * eta_z * zeta_z * phi_{eta zeta}
511  d2xidxyz2[p][5]*dphidxi[i][p] + // xi_{z z} * phi_{xi}
512  d2etadxyz2[p][5]*dphideta[i][p] + // eta_{z z} * phi_{eta}
513  d2zetadxyz2[p][5]*dphidzeta[i][p]; // zeta_{z z} * phi_{zeta}
514  }
515 
516  break;
517  }
518 
519  default:
520  libmesh_error_msg("Invalid dim = " << dim);
521  } // switch(dim)
522 }
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 >.

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

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

Definition at line 660 of file h1_fe_transformation.C.

665 {
666  libmesh_error_msg("Computing the divergence of a shape function only \nmakes sense for vector-valued elements.");
667 }
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

Definition at line 671 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.

676 {
677  switch (dim)
678  {
679  case 0:
680  case 1:
681  libmesh_error_msg("The divergence only make sense in 2D and 3D");
682 
683  case 2:
684  {
685  const std::vector<std::vector<RealGradient> > & dphidxi = fe.get_dphidxi();
686  const std::vector<std::vector<RealGradient> > & dphideta = fe.get_dphideta();
687 
688  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
689  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
690 
691  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
692  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
693 
694  // In 2D: div = dphi_x/dx + dphi_y/dy
695  for (std::size_t i=0; i<div_phi.size(); i++)
696  for (std::size_t p=0; p<div_phi[i].size(); p++)
697  {
698 
699  Real dphix_dx = (dphidxi[i][p].slice(0))*dxidx_map[p]
700  + (dphideta[i][p].slice(0))*detadx_map[p];
701 
702  Real dphiy_dy = (dphidxi[i][p].slice(1))*dxidy_map[p]
703  + (dphideta[i][p].slice(1))*detady_map[p];
704 
705  div_phi[i][p] = dphix_dx + dphiy_dy;
706  }
707  break;
708  }
709  case 3:
710  {
711  const std::vector<std::vector<RealGradient> > & dphidxi = fe.get_dphidxi();
712  const std::vector<std::vector<RealGradient> > & dphideta = fe.get_dphideta();
713  const std::vector<std::vector<RealGradient> > & dphidzeta = fe.get_dphidzeta();
714 
715  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
716  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
717  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
718 
719  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
720  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
721  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
722 
723  const std::vector<Real> & dzetadx_map = fe.get_fe_map().get_dzetadx();
724  const std::vector<Real> & dzetady_map = fe.get_fe_map().get_dzetady();
725  const std::vector<Real> & dzetadz_map = fe.get_fe_map().get_dzetadz();
726 
727  /*
728  In 3D: div = dphi_x/dx + dphi_y/dy + dphi_z/dz
729  */
730  for (std::size_t i=0; i<div_phi.size(); i++)
731  for (std::size_t p=0; p<div_phi[i].size(); p++)
732  {
733  Real dphix_dx = (dphidxi[i][p].slice(0))*dxidx_map[p]
734  + (dphideta[i][p].slice(0))*detadx_map[p]
735  + (dphidzeta[i][p].slice(0))*dzetadx_map[p];
736 
737  Real dphiy_dy = (dphidxi[i][p].slice(1))*dxidy_map[p]
738  + (dphideta[i][p].slice(1))*detady_map[p]
739  + (dphidzeta[i][p].slice(1))*dzetady_map[p];
740 
741  Real dphiz_dz = (dphidxi[i][p].slice(2))*dxidz_map[p]
742  + (dphideta[i][p].slice(2))*detadz_map[p]
743  + (dphidzeta[i][p].slice(2))*dzetadz_map[p];
744 
745  div_phi[i][p] = dphix_dx + dphiy_dy + dphiz_dz;
746  }
747 
748  break;
749  }
750 
751  default: \
752  libmesh_error_msg("Invalid dim = " << dim); \
753  } // switch(dim)
754 
755  return;
756 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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
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 107 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().

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

115 {
116  switch(dim)
117  {
118  case 0: // No derivatives in 0D
119  {
120  for (std::size_t i=0; i<dphi.size(); i++)
121  for (std::size_t p=0; p<dphi[i].size(); p++)
122  dphi[i][p] = 0.;
123  break;
124  }
125 
126  case 1:
127  {
128  const std::vector<std::vector<OutputShape> > & dphidxi = fe.get_dphidxi();
129 
130  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
131 #if LIBMESH_DIM>1
132  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
133 #endif
134 #if LIBMESH_DIM>2
135  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
136 #endif
137 
138  for (std::size_t i=0; i<dphi.size(); i++)
139  for (std::size_t p=0; p<dphi[i].size(); p++)
140  {
141  // dphi/dx = (dphi/dxi)*(dxi/dx)
142  dphi[i][p].slice(0) = dphidx[i][p] = dphidxi[i][p]*dxidx_map[p];
143 
144 #if LIBMESH_DIM>1
145  dphi[i][p].slice(1) = dphidy[i][p] = dphidxi[i][p]*dxidy_map[p];
146 #endif
147 #if LIBMESH_DIM>2
148  dphi[i][p].slice(2) = dphidz[i][p] = dphidxi[i][p]*dxidz_map[p];
149 #endif
150  }
151 
152  break;
153  }
154 
155  case 2:
156  {
157  const std::vector<std::vector<OutputShape> > & dphidxi = fe.get_dphidxi();
158  const std::vector<std::vector<OutputShape> > & dphideta = fe.get_dphideta();
159 
160  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
161  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
162 #if LIBMESH_DIM > 2
163  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
164 #endif
165 
166  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
167  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
168 #if LIBMESH_DIM > 2
169  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
170 #endif
171 
172  for (std::size_t i=0; i<dphi.size(); i++)
173  for (std::size_t p=0; p<dphi[i].size(); p++)
174  {
175  // dphi/dx = (dphi/dxi)*(dxi/dx) + (dphi/deta)*(deta/dx)
176  dphi[i][p].slice(0) = dphidx[i][p] = (dphidxi[i][p]*dxidx_map[p] +
177  dphideta[i][p]*detadx_map[p]);
178 
179  // dphi/dy = (dphi/dxi)*(dxi/dy) + (dphi/deta)*(deta/dy)
180  dphi[i][p].slice(1) = dphidy[i][p] = (dphidxi[i][p]*dxidy_map[p] +
181  dphideta[i][p]*detady_map[p]);
182 
183 #if LIBMESH_DIM > 2
184  // dphi/dz = (dphi/dxi)*(dxi/dz) + (dphi/deta)*(deta/dz)
185  dphi[i][p].slice(2) = dphidz[i][p] = (dphidxi[i][p]*dxidz_map[p] +
186  dphideta[i][p]*detadz_map[p]);
187 #endif
188  }
189 
190  break;
191  }
192 
193  case 3:
194  {
195  const std::vector<std::vector<OutputShape> > & dphidxi = fe.get_dphidxi();
196  const std::vector<std::vector<OutputShape> > & dphideta = fe.get_dphideta();
197  const std::vector<std::vector<OutputShape> > & dphidzeta = fe.get_dphidzeta();
198 
199  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
200  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
201  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
202 
203  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
204  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
205  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
206 
207  const std::vector<Real> & dzetadx_map = fe.get_fe_map().get_dzetadx();
208  const std::vector<Real> & dzetady_map = fe.get_fe_map().get_dzetady();
209  const std::vector<Real> & dzetadz_map = fe.get_fe_map().get_dzetadz();
210 
211  for (std::size_t i=0; i<dphi.size(); i++)
212  for (std::size_t p=0; p<dphi[i].size(); p++)
213  {
214  // dphi/dx = (dphi/dxi)*(dxi/dx) + (dphi/deta)*(deta/dx) + (dphi/dzeta)*(dzeta/dx);
215  dphi[i][p].slice(0) = dphidx[i][p] = (dphidxi[i][p]*dxidx_map[p] +
216  dphideta[i][p]*detadx_map[p] +
217  dphidzeta[i][p]*dzetadx_map[p]);
218 
219  // dphi/dy = (dphi/dxi)*(dxi/dy) + (dphi/deta)*(deta/dy) + (dphi/dzeta)*(dzeta/dy);
220  dphi[i][p].slice(1) = dphidy[i][p] = (dphidxi[i][p]*dxidy_map[p] +
221  dphideta[i][p]*detady_map[p] +
222  dphidzeta[i][p]*dzetady_map[p]);
223 
224  // dphi/dz = (dphi/dxi)*(dxi/dz) + (dphi/deta)*(deta/dz) + (dphi/dzeta)*(dzeta/dz);
225  dphi[i][p].slice(2) = dphidz[i][p] = (dphidxi[i][p]*dxidz_map[p] +
226  dphideta[i][p]*detadz_map[p] +
227  dphidzeta[i][p]*dzetadz_map[p]);
228  }
229  break;
230  }
231 
232  default:
233  libmesh_error_msg("Invalid dim = " << dim);
234  } // switch(dim)
235 }
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
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 51 of file h1_fe_transformation.C.

References libMesh::FEAbstract::get_fe_type().

Referenced by libMesh::H1FETransformation< OutputShape >::~H1FETransformation().

56 {
57  switch(dim)
58  {
59  case 0:
60  {
61  for (std::size_t i=0; i<phi.size(); i++)
62  {
63  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
64  for (std::size_t p=0; p<phi[i].size(); p++)
65  FEInterface::shape<OutputShape>(0, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
66  }
67  break;
68  }
69  case 1:
70  {
71  for (std::size_t i=0; i<phi.size(); i++)
72  {
73  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
74  for (std::size_t p=0; p<phi[i].size(); p++)
75  FEInterface::shape<OutputShape>(1, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
76  }
77  break;
78  }
79  case 2:
80  {
81  for (std::size_t i=0; i<phi.size(); i++)
82  {
83  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
84  for (std::size_t p=0; p<phi[i].size(); p++)
85  FEInterface::shape<OutputShape>(2, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
86  }
87  break;
88  }
89  case 3:
90  {
91  for (std::size_t i=0; i<phi.size(); i++)
92  {
93  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
94  for (std::size_t p=0; p<phi[i].size(); p++)
95  FEInterface::shape<OutputShape>(3, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
96  }
97  break;
98  }
99 
100  default:
101  libmesh_error_msg("Invalid dim = " << dim);
102  }
103 }

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