Public Member Functions | Static Public Member Functions | List of all members
libMesh::H1FETransformation< OutputShape > Class Template Reference

#include <fe_transformation_base.h>

Inheritance diagram for libMesh::H1FETransformation< OutputShape >:
Inheritance graph
[legend]

Public Member Functions

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

◆ H1FETransformation()

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

Definition at line 44 of file h1_fe_transformation.h.

45  : FETransformationBase<OutputShape>(){}

◆ ~H1FETransformation()

template<typename OutputShape >
virtual libMesh::H1FETransformation< OutputShape >::~H1FETransformation ( )
inlinevirtual

Definition at line 47 of file h1_fe_transformation.h.

47 {}

Member Function Documentation

◆ build()

template<typename OutputShape >
std::unique_ptr< FETransformationBase< OutputShape > > libMesh::FETransformationBase< OutputShape >::build ( const FEType &  type)
staticinherited

Builds an FETransformation object based on the finite element type

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

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

◆ init_map_d2phi()

template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::init_map_d2phi ( const FEGenericBase< OutputShape > &  fe) const
overridevirtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 43 of file h1_fe_transformation.C.

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

44 {
45  fe.get_fe_map().get_dxidx();
46 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
47  fe.get_fe_map().get_d2xidxyz2();
48 #endif
49 }

◆ init_map_dphi()

template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::init_map_dphi ( const FEGenericBase< OutputShape > &  fe) const
overridevirtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 35 of file h1_fe_transformation.C.

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

36 {
37  fe.get_fe_map().get_dxidx();
38 }

◆ init_map_phi()

template<typename OutputShape >
void libMesh::H1FETransformation< OutputShape >::init_map_phi ( const FEGenericBase< OutputShape > &  fe) const
overridevirtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 27 of file h1_fe_transformation.C.

28 {
29  // Nothing needed
30 }

◆ map_curl() [1/3]

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
overridevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

◆ map_curl() [2/3]

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 527 of file h1_fe_transformation.C.

532 {
533  libmesh_error_msg("Computing the curl of a shape function only \nmakes sense for vector-valued elements.");
534 }

◆ map_curl() [3/3]

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 537 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(), libMesh::index_range(), and libMesh::Real.

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

◆ map_d2phi()

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
overridevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 240 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(), libMesh::FEAbstract::get_fe_map(), and libMesh::index_range().

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

◆ map_div() [1/3]

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
overridevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

◆ map_div() [2/3]

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 656 of file h1_fe_transformation.C.

661 {
662  libmesh_error_msg("Computing the divergence of a shape function only \nmakes sense for vector-valued elements.");
663 }

◆ map_div() [3/3]

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 667 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(), libMesh::index_range(), and libMesh::Real.

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

◆ map_dphi()

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
overridevirtual

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 108 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::index_range().

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

◆ map_phi()

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
overridevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 52 of file h1_fe_transformation.C.

References libMesh::FEAbstract::get_fe_type(), and libMesh::index_range().

57 {
58  switch(dim)
59  {
60  case 0:
61  {
62  for (auto i : index_range(phi))
63  {
64  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
65  for (auto p : index_range(phi[i]))
66  FEInterface::shape<OutputShape>(0, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
67  }
68  break;
69  }
70  case 1:
71  {
72  for (auto i : index_range(phi))
73  {
74  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
75  for (auto p : index_range(phi[i]))
76  FEInterface::shape<OutputShape>(1, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
77  }
78  break;
79  }
80  case 2:
81  {
82  for (auto i : index_range(phi))
83  {
84  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
85  for (auto p : index_range(phi[i]))
86  FEInterface::shape<OutputShape>(2, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
87  }
88  break;
89  }
90  case 3:
91  {
92  for (auto i : index_range(phi))
93  {
94  libmesh_assert_equal_to ( qp.size(), phi[i].size() );
95  for (auto p : index_range(phi[i]))
96  FEInterface::shape<OutputShape>(3, fe.get_fe_type(), elem, i, qp[p], phi[i][p]);
97  }
98  break;
99  }
100 
101  default:
102  libmesh_error_msg("Invalid dim = " << dim);
103  }
104 }
libMesh::index_range
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Definition: int_range.h:104
The documentation for this class was generated from the following files:
generated by doxygen