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

This class handles the computation of the shape functions in the physical domain for HCurl conforming elements. More...

#include <fe_transformation_base.h>

Inheritance diagram for libMesh::HCurlFETransformation< OutputShape >:
[legend]

Public Member Functions

 HCurlFETransformation ()
 
virtual ~HCurlFETransformation ()=default
 
virtual void init_map_phi (const FEGenericBase< OutputShape > &fe) const override
 Pre-requests any necessary data from FEMap. More...
 
virtual void init_map_dphi (const FEGenericBase< OutputShape > &fe) const override
 Pre-requests any necessary data from FEMap. More...
 
virtual void init_map_d2phi (const FEGenericBase< OutputShape > &fe) const override
 Pre-requests any necessary data from FEMap. More...
 
virtual void 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, bool add_p_level=true) const override
 Evaluates shape functions in physical coordinates for \( H(curl) \) conforming elements. More...
 
virtual void map_dphi (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< OutputShape > &, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputGradient >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &) const override
 Evaluates shape function gradients in physical coordinates for \( H(curl) \) conforming elements. More...
 
virtual void map_d2phi (const unsigned int, const std::vector< Point > &, const FEGenericBase< OutputShape > &, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputTensor >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &, std::vector< std::vector< OutputShape >> &) const override
 Evaluates shape function Hessians in physical coordinates based on \( H(curl) \) conforming finite element transformation. More...
 
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
 Evaluates the shape function curl in physical coordinates based on \( H(curl) \) conforming finite element transformation. More...
 
virtual void map_div (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< OutputShape > &, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputDivergence >> &) const override
 Evaluates the shape function divergence in physical coordinates based on \( H(curl) \) conforming finite element transformation. More...
 
template<>
void init_map_phi (const FEGenericBase< Real > &) const
 
template<>
void init_map_dphi (const FEGenericBase< Real > &) const
 
template<>
void init_map_d2phi (const FEGenericBase< Real > &) const
 
template<>
void map_phi (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< Real > &, std::vector< std::vector< Real >> &, bool) const
 
template<>
void map_curl (const unsigned int, const Elem *const, const std::vector< Point > &, const FEGenericBase< Real > &, std::vector< std::vector< Real >> &) const
 

Static Public Member Functions

static std::unique_ptr< FETransformationBase< OutputShape > > build (const FEType &type)
 Builds an FETransformation object based on the finite element type. More...
 

Detailed Description

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

This class handles the computation of the shape functions in the physical domain for HCurl 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 29 of file fe_transformation_base.h.

Constructor & Destructor Documentation

◆ HCurlFETransformation()

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

Definition at line 40 of file hcurl_fe_transformation.h.

41  : FETransformationBase<OutputShape>() {}

◆ ~HCurlFETransformation()

template<typename OutputShape >
virtual libMesh::HCurlFETransformation< OutputShape >::~HCurlFETransformation ( )
virtualdefault

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 32 of file fe_transformation_base.C.

References libMesh::BERNSTEIN, libMesh::CLOUGH, libMesh::Utility::enum_to_string(), libMesh::FEType::family, libMesh::HERMITE, libMesh::HIERARCHIC, libMesh::HIERARCHIC_VEC, libMesh::JACOBI_20_00, libMesh::JACOBI_30_00, libMesh::L2_HIERARCHIC, libMesh::L2_HIERARCHIC_VEC, libMesh::L2_LAGRANGE, libMesh::L2_LAGRANGE_VEC, libMesh::L2_RAVIART_THOMAS, libMesh::LAGRANGE, libMesh::LAGRANGE_VEC, libMesh::MONOMIAL, libMesh::MONOMIAL_VEC, libMesh::NEDELEC_ONE, libMesh::RATIONAL_BERNSTEIN, libMesh::RAVIART_THOMAS, libMesh::SCALAR, libMesh::SIDE_HIERARCHIC, libMesh::SUBDIVISION, libMesh::SZABAB, and libMesh::XYZ.

33 {
34  switch (fe_type.family)
35  {
36  // H1 Conforming Elements
37  case LAGRANGE:
38  case HIERARCHIC:
39  case BERNSTEIN:
40  case SZABAB:
41  case CLOUGH: // PB: Really H2
42  case HERMITE: // PB: Really H2
43  case SUBDIVISION:
44  case LAGRANGE_VEC:
45  case HIERARCHIC_VEC:
46  case MONOMIAL: // PB: Shouldn't this be L2 conforming?
47  case MONOMIAL_VEC: // PB: Shouldn't this be L2 conforming?
48  case XYZ: // PB: Shouldn't this be L2 conforming?
49  case RATIONAL_BERNSTEIN:
50  case L2_HIERARCHIC:
51  case L2_HIERARCHIC_VEC:
52  case SIDE_HIERARCHIC:
53  case L2_LAGRANGE: // PB: Shouldn't this be L2 conforming?
54  case L2_LAGRANGE_VEC: // PB: Shouldn't this be L2 conforming?
55  case JACOBI_20_00: // PB: For infinite elements...
56  case JACOBI_30_00: // PB: For infinite elements...
57  return std::make_unique<H1FETransformation<OutputShape>>();
58 
59  // HCurl Conforming Elements
60  case NEDELEC_ONE:
61  return std::make_unique<HCurlFETransformation<OutputShape>>();
62 
63  // HDiv Conforming Elements
64  case RAVIART_THOMAS:
65  case L2_RAVIART_THOMAS:
66  return std::make_unique<HDivFETransformation<OutputShape>>();
67 
68  // L2 Conforming Elements
69 
70  // Other...
71  case SCALAR:
72  // Should never need this for SCALARs
73  return std::make_unique<H1FETransformation<OutputShape>>();
74 
75  default:
76  libmesh_error_msg("Unknown family = " << Utility::enum_to_string(fe_type.family));
77  }
78 }
std::string enum_to_string(const T e)

◆ init_map_d2phi() [1/2]

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

Pre-requests any necessary data from FEMap.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 42 of file hcurl_fe_transformation.C.

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

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

◆ init_map_d2phi() [2/2]

template<>
void libMesh::HCurlFETransformation< Real >::init_map_d2phi ( const FEGenericBase< Real > &  ) const

Definition at line 285 of file hcurl_fe_transformation.C.

286 {
287  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
288 }

◆ init_map_dphi() [1/2]

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

Pre-requests any necessary data from FEMap.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 34 of file hcurl_fe_transformation.C.

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

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

◆ init_map_dphi() [2/2]

template<>
void libMesh::HCurlFETransformation< Real >::init_map_dphi ( const FEGenericBase< Real > &  ) const

Definition at line 279 of file hcurl_fe_transformation.C.

280 {
281  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
282 }

◆ init_map_phi() [1/2]

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

Pre-requests any necessary data from FEMap.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 26 of file hcurl_fe_transformation.C.

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

27 {
28  fe.get_fe_map().get_dxidx();
29 }

◆ init_map_phi() [2/2]

template<>
void libMesh::HCurlFETransformation< Real >::init_map_phi ( const FEGenericBase< Real > &  ) const

Definition at line 273 of file hcurl_fe_transformation.C.

274 {
275  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
276 }

◆ map_curl() [1/2]

template<typename OutputShape >
void libMesh::HCurlFETransformation< 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 \( H(curl) \) conforming finite element transformation.

In 2-D, the transformation is \( \nabla \times \phi = J^{-1} * \nabla \times \hat{\phi} \) where \( J = \det( dx/d\xi ) \) In 3-D, the transformation is \( \nabla \times \phi = J^{-1} dx/d\xi \nabla \times \hat{\phi} \)

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 161 of file hcurl_fe_transformation.C.

References dim, 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_dxyzdeta(), libMesh::FEMap::get_dxyzdxi(), libMesh::FEMap::get_dxyzdzeta(), libMesh::FEAbstract::get_fe_map(), libMesh::FEMap::get_jacobian(), libMesh::index_range(), and libMesh::Real.

166 {
167  switch (dim)
168  {
169  case 0:
170  case 1:
171  libmesh_error_msg("These element transformations only make sense in 2D and 3D.");
172 
173  case 2:
174  {
175  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
176  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
177  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
178 
179  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
180  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
181  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
182 
183  const std::vector<std::vector<OutputShape>> & dphi_dxi = fe.get_dphidxi();
184  const std::vector<std::vector<OutputShape>> & dphi_deta = fe.get_dphideta();
185 
186  // In 2D (embedded in 3D): curl(phi)_j = C(j) * curl(\hat{phi}) where C(j)
187  // is the jth column cofactor of (dx/dxi)^{-1} = [ dxi/dx dxi/dy dxi/dz
188  // deta/dx deta/dy deta/dz ]
189  for (auto i : index_range(curl_phi))
190  for (auto p : index_range(curl_phi[i]))
191  {
192  const Real curl_ref = dphi_dxi[i][p].slice(1) - dphi_deta[i][p].slice(0);
193  curl_phi[i][p].slice(0) = curl_ref * (dxidy_map[p]*detadz_map[p] - dxidz_map[p]*detady_map[p]);
194  curl_phi[i][p].slice(1) = curl_ref * -(dxidx_map[p]*detadz_map[p] - dxidz_map[p]*detadx_map[p]);
195  curl_phi[i][p].slice(2) = curl_ref * (dxidx_map[p]*detady_map[p] - dxidy_map[p]*detadx_map[p]);
196  }
197 
198  break;
199  }
200 
201  case 3:
202  {
203  const std::vector<std::vector<OutputShape>> & dphi_dxi = fe.get_dphidxi();
204  const std::vector<std::vector<OutputShape>> & dphi_deta = fe.get_dphideta();
205  const std::vector<std::vector<OutputShape>> & dphi_dzeta = fe.get_dphidzeta();
206 
207  const std::vector<RealGradient> & dxyz_dxi = fe.get_fe_map().get_dxyzdxi();
208  const std::vector<RealGradient> & dxyz_deta = fe.get_fe_map().get_dxyzdeta();
209  const std::vector<RealGradient> & dxyz_dzeta = fe.get_fe_map().get_dxyzdzeta();
210 
211  const std::vector<Real> & J = fe.get_fe_map().get_jacobian();
212 
213  for (auto i : index_range(curl_phi))
214  for (auto p : index_range(curl_phi[i]))
215  {
216  Real dx_dxi = dxyz_dxi[p](0);
217  Real dx_deta = dxyz_deta[p](0);
218  Real dx_dzeta = dxyz_dzeta[p](0);
219 
220  Real dy_dxi = dxyz_dxi[p](1);
221  Real dy_deta = dxyz_deta[p](1);
222  Real dy_dzeta = dxyz_dzeta[p](1);
223 
224  Real dz_dxi = dxyz_dxi[p](2);
225  Real dz_deta = dxyz_deta[p](2);
226  Real dz_dzeta = dxyz_dzeta[p](2);
227 
228  const Real inv_jac = 1.0/J[p];
229 
230  /* In 3D: curl(phi) = J^{-1} dx/dxi * curl(\hat{phi})
231 
232  dx/dxi = [ dx/dxi dx/deta dx/dzeta
233  dy/dxi dy/deta dy/dzeta
234  dz/dxi dz/deta dz/dzeta ]
235 
236  curl(u) = [ du_z/deta - du_y/dzeta
237  du_x/dzeta - du_z/dxi
238  du_y/dxi - du_x/deta ]
239  */
240  curl_phi[i][p].slice(0) = inv_jac*( dx_dxi*( dphi_deta[i][p].slice(2) -
241  dphi_dzeta[i][p].slice(1) ) +
242  dx_deta*( dphi_dzeta[i][p].slice(0) -
243  dphi_dxi[i][p].slice(2) ) +
244  dx_dzeta*( dphi_dxi[i][p].slice(1) -
245  dphi_deta[i][p].slice(0) ) );
246 
247  curl_phi[i][p].slice(1) = inv_jac*( dy_dxi*( dphi_deta[i][p].slice(2) -
248  dphi_dzeta[i][p].slice(1) ) +
249  dy_deta*( dphi_dzeta[i][p].slice(0)-
250  dphi_dxi[i][p].slice(2) ) +
251  dy_dzeta*( dphi_dxi[i][p].slice(1) -
252  dphi_deta[i][p].slice(0) ) );
253 
254  curl_phi[i][p].slice(2) = inv_jac*( dz_dxi*( dphi_deta[i][p].slice(2) -
255  dphi_dzeta[i][p].slice(1) ) +
256  dz_deta*( dphi_dzeta[i][p].slice(0) -
257  dphi_dxi[i][p].slice(2) ) +
258  dz_dzeta*( dphi_dxi[i][p].slice(1) -
259  dphi_deta[i][p].slice(0) ) );
260  }
261 
262  break;
263  }
264 
265  default:
266  libmesh_error_msg("Invalid dim = " << dim);
267  } // switch(dim)
268 }
unsigned int dim
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117

◆ map_curl() [2/2]

template<>
void libMesh::HCurlFETransformation< 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 302 of file hcurl_fe_transformation.C.

307 {
308  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
309 }

◆ map_d2phi()

template<typename OutputShape >
virtual void libMesh::HCurlFETransformation< OutputShape >::map_d2phi ( const unsigned int  ,
const std::vector< Point > &  ,
const FEGenericBase< OutputShape > &  ,
std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputTensor >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &   
) const
inlineoverridevirtual

Evaluates shape function Hessians in physical coordinates based on \( H(curl) \) conforming finite element transformation.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 96 of file hcurl_fe_transformation.h.

106  {
107  libmesh_warning("WARNING: Shape function Hessians for HCurl elements are not currently being computed!");
108  }

◆ map_div()

template<typename OutputShape >
virtual void libMesh::HCurlFETransformation< OutputShape >::map_div ( const unsigned int  ,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< OutputShape > &  ,
std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputDivergence >> &   
) const
inlineoverridevirtual

Evaluates the shape function divergence in physical coordinates based on \( H(curl) \) conforming finite element transformation.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 128 of file hcurl_fe_transformation.h.

133  {
134  libmesh_warning("WARNING: Shape function divergences for HCurl elements are not currently being computed!");
135  }

◆ map_dphi()

template<typename OutputShape >
virtual void libMesh::HCurlFETransformation< OutputShape >::map_dphi ( const unsigned int  ,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< OutputShape > &  ,
std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputGradient >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &  ,
std::vector< std::vector< OutputShape >> &   
) const
inlineoverridevirtual

Evaluates shape function gradients in physical coordinates for \( H(curl) \) conforming elements.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 79 of file hcurl_fe_transformation.h.

87  {
88  libmesh_warning("WARNING: Shape function gradients for HCurl elements are not currently being computed!");
89  }

◆ map_phi() [1/2]

template<typename OutputShape >
void libMesh::HCurlFETransformation< 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,
bool  add_p_level = true 
) const
overridevirtual

Evaluates shape functions in physical coordinates for \( H(curl) \) conforming elements.

In this case \( \phi = (dx/d\xi)^{-T} \hat{\phi} \), where \( (dx/d\xi)^{-T} \) is the inverse-transpose of the Jacobian matrix of the element map.

Note
Here \( x, \xi \) are vectors.

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 53 of file hcurl_fe_transformation.C.

References dim, 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(), libMesh::FEAbstract::get_fe_type(), libMesh::index_range(), and libMesh::FEInterface::shape().

59 {
60  switch (dim)
61  {
62  case 0:
63  case 1:
64  libmesh_error_msg("These element transformations only make sense in 2D and 3D.");
65 
66  case 2:
67  {
68  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
69  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
70  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
71 
72  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
73  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
74  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
75 
76  // phi = (dx/dxi)^-T * \hat{phi}
77  // In 2D (embedded in 3D):
78  // (dx/dxi)^{-1} = [ dxi/dx dxi/dy dxi/dz
79  // deta/dx deta/dy deta/dz ]
80  //
81  // so: dxi/dx^{-T} * \hat{phi} = [ dxi/dx deta/dx [ \hat{phi}_xi
82  // dxi/dy deta/dy \hat{phi}_eta ]
83  // dxi/dz deta/dz ]
84  //
85  // or in indicial notation: phi_j = xi_{i,j}*\hat{phi}_i
86  for (auto i : index_range(phi))
87  for (auto p : index_range(phi[i]))
88  {
89  // Need to temporarily cache reference shape functions
90  // TODO: PB: Might be worth trying to build phi_ref separately to see
91  // if we can get vectorization
92  // We are computing mapping basis functions, so we explicitly ignore
93  // any non-zero p_level() the Elem might have.
94  OutputShape phi_ref;
95  FEInterface::shape(fe.get_fe_type(), /*extra_order=*/0, elem, i, qp[p], phi_ref);
96 
97  phi[i][p](0) = dxidx_map[p]*phi_ref.slice(0) + detadx_map[p]*phi_ref.slice(1);
98 
99  phi[i][p](1) = dxidy_map[p]*phi_ref.slice(0) + detady_map[p]*phi_ref.slice(1);
100 
101  phi[i][p](2) = dxidz_map[p]*phi_ref.slice(0) + detadz_map[p]*phi_ref.slice(1);
102  }
103 
104  break;
105  }
106 
107  case 3:
108  {
109  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
110  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
111  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
112 
113  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
114  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
115  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
116 
117  const std::vector<Real> & dzetadx_map = fe.get_fe_map().get_dzetadx();
118  const std::vector<Real> & dzetady_map = fe.get_fe_map().get_dzetady();
119  const std::vector<Real> & dzetadz_map = fe.get_fe_map().get_dzetadz();
120 
121  // phi = (dx/dxi)^-T * \hat{phi}
122  // In 3D:
123  // dx/dxi^-1 = [ dxi/dx dxi/dy dxi/dz
124  // deta/dx deta/dy deta/dz
125  // dzeta/dx dzeta/dy dzeta/dz]
126  //
127  // so: dxi/dx^-T * \hat{phi} = [ dxi/dx deta/dx dzeta/dx [ \hat{phi}_xi
128  // dxi/dy deta/dy dzeta/dy \hat{phi}_eta
129  // dxi/dz deta/dz dzeta/dz ] \hat{phi}_zeta ]
130  //
131  // or in indicial notation: phi_j = xi_{i,j}*\hat{phi}_i
132  for (auto i : index_range(phi))
133  for (auto p : index_range(phi[i]))
134  {
135  // Need to temporarily cache reference shape functions
136  // TODO: PB: Might be worth trying to build phi_ref separately to see
137  // if we can get vectorization
138  // We are computing mapping basis functions, so we explicitly ignore
139  // any non-zero p_level() the Elem might have.
140  OutputShape phi_ref;
141  FEInterface::shape(fe.get_fe_type(), /*extra_order=*/0, elem, i, qp[p], phi_ref);
142 
143  phi[i][p].slice(0) = dxidx_map[p]*phi_ref.slice(0) + detadx_map[p]*phi_ref.slice(1)
144  + dzetadx_map[p]*phi_ref.slice(2);
145 
146  phi[i][p].slice(1) = dxidy_map[p]*phi_ref.slice(0) + detady_map[p]*phi_ref.slice(1)
147  + dzetady_map[p]*phi_ref.slice(2);
148 
149  phi[i][p].slice(2) = dxidz_map[p]*phi_ref.slice(0) + detadz_map[p]*phi_ref.slice(1)
150  + dzetadz_map[p]*phi_ref.slice(2);
151  }
152  break;
153  }
154 
155  default:
156  libmesh_error_msg("Invalid dim = " << dim);
157  } // switch(dim)
158 }
unsigned int dim
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:760
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117

◆ map_phi() [2/2]

template<>
void libMesh::HCurlFETransformation< Real >::map_phi ( const unsigned int  ,
const Elem const,
const std::vector< Point > &  ,
const FEGenericBase< Real > &  ,
std::vector< std::vector< Real >> &  ,
bool   
) const

Definition at line 291 of file hcurl_fe_transformation.C.

297 {
298  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
299 }

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