libMesh::HCurlFETransformation< OutputShape > Class Template Reference

#include <fe_transformation_base.h>

Inheritance diagram for libMesh::HCurlFETransformation< OutputShape >:

Public Member Functions

 HCurlFETransformation ()
 
virtual ~HCurlFETransformation ()
 
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 dim, const Elem *const elem, const std::vector< Point > &qp, const FEGenericBase< OutputShape > &fe, std::vector< std::vector< OutputShape > > &phi) const libmesh_override
 
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 libmesh_override
 
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 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, const Elem *const , const std::vector< Point > &, const FEGenericBase< OutputShape > &, std::vector< std::vector< typename FEGenericBase< OutputShape >::OutputDivergence > > &) const libmesh_override
 
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 > > &) 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 UniquePtr< FETransformationBase< OutputShape > > build (const FEType &type)
 

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

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

Definition at line 40 of file hcurl_fe_transformation.h.

41  : 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::HCurlFETransformation< OutputShape >::init_map_d2phi ( const FEGenericBase< OutputShape > &  fe) const
virtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 41 of file hcurl_fe_transformation.C.

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

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

42 {
43  fe.get_fe_map().get_dxidx();
44 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
45  fe.get_fe_map().get_d2xidxyz2();
46 #endif
47 }
template<>
void libMesh::HCurlFETransformation< Real >::init_map_d2phi ( const FEGenericBase< Real > &  ) const

Definition at line 269 of file hcurl_fe_transformation.C.

270 {
271  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
272 }
template<typename OutputShape >
void libMesh::HCurlFETransformation< OutputShape >::init_map_dphi ( const FEGenericBase< OutputShape > &  fe) const
virtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 33 of file hcurl_fe_transformation.C.

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

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

34 {
35  fe.get_fe_map().get_dxidx();
36 }
template<>
void libMesh::HCurlFETransformation< Real >::init_map_dphi ( const FEGenericBase< Real > &  ) const

Definition at line 263 of file hcurl_fe_transformation.C.

264 {
265  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
266 }
template<typename OutputShape >
void libMesh::HCurlFETransformation< OutputShape >::init_map_phi ( const FEGenericBase< OutputShape > &  fe) const
virtual

Pre-requests any necessary data from FEMap

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 25 of file hcurl_fe_transformation.C.

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

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

26 {
27  fe.get_fe_map().get_dxidx();
28 }
template<>
void libMesh::HCurlFETransformation< Real >::init_map_phi ( const FEGenericBase< Real > &  ) const

Definition at line 257 of file hcurl_fe_transformation.C.

258 {
259  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
260 }
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
virtual

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 153 of file hcurl_fe_transformation.C.

References libMesh::FEGenericBase< OutputType >::get_dphideta(), libMesh::FEGenericBase< OutputType >::get_dphidxi(), libMesh::FEGenericBase< OutputType >::get_dphidzeta(), libMesh::FEMap::get_dxyzdeta(), libMesh::FEMap::get_dxyzdxi(), libMesh::FEMap::get_dxyzdzeta(), libMesh::FEAbstract::get_fe_map(), libMesh::FEMap::get_jacobian(), and libMesh::Real.

Referenced by libMesh::HCurlFETransformation< OutputShape >::map_d2phi().

158 {
159  switch (dim)
160  {
161  case 0:
162  case 1:
163  libmesh_error_msg("These element transformations only make sense in 2D and 3D.");
164 
165  case 2:
166  {
167  const std::vector<std::vector<OutputShape> > & dphi_dxi = fe.get_dphidxi();
168  const std::vector<std::vector<OutputShape> > & dphi_deta = fe.get_dphideta();
169 
170  const std::vector<Real> & J = fe.get_fe_map().get_jacobian();
171 
172  // FIXME: I don't think this is valid for 2D elements in 3D space
173  /* In 2D: curl(phi) = J^{-1} * curl(\hat{phi}) */
174  for (std::size_t i=0; i<curl_phi.size(); i++)
175  for (std::size_t p=0; p<curl_phi[i].size(); p++)
176  {
177  curl_phi[i][p].slice(0) = curl_phi[i][p].slice(1) = 0.0;
178 
179  curl_phi[i][p].slice(2) = ( dphi_dxi[i][p].slice(1) - dphi_deta[i][p].slice(0) )/J[p];
180  }
181 
182  break;
183  }
184 
185  case 3:
186  {
187  const std::vector<std::vector<OutputShape> > & dphi_dxi = fe.get_dphidxi();
188  const std::vector<std::vector<OutputShape> > & dphi_deta = fe.get_dphideta();
189  const std::vector<std::vector<OutputShape> > & dphi_dzeta = fe.get_dphidzeta();
190 
191  const std::vector<RealGradient> & dxyz_dxi = fe.get_fe_map().get_dxyzdxi();
192  const std::vector<RealGradient> & dxyz_deta = fe.get_fe_map().get_dxyzdeta();
193  const std::vector<RealGradient> & dxyz_dzeta = fe.get_fe_map().get_dxyzdzeta();
194 
195  const std::vector<Real> & J = fe.get_fe_map().get_jacobian();
196 
197  for (std::size_t i=0; i<curl_phi.size(); i++)
198  for (std::size_t p=0; p<curl_phi[i].size(); p++)
199  {
200  Real dx_dxi = dxyz_dxi[p](0);
201  Real dx_deta = dxyz_deta[p](0);
202  Real dx_dzeta = dxyz_dzeta[p](0);
203 
204  Real dy_dxi = dxyz_dxi[p](1);
205  Real dy_deta = dxyz_deta[p](1);
206  Real dy_dzeta = dxyz_dzeta[p](1);
207 
208  Real dz_dxi = dxyz_dxi[p](2);
209  Real dz_deta = dxyz_deta[p](2);
210  Real dz_dzeta = dxyz_dzeta[p](2);
211 
212  const Real inv_jac = 1.0/J[p];
213 
214  /* In 3D: curl(phi) = J^{-1} dx/dxi * curl(\hat{phi})
215 
216  dx/dxi = [ dx/dxi dx/deta dx/dzeta
217  dy/dxi dy/deta dy/dzeta
218  dz/dxi dz/deta dz/dzeta ]
219 
220  curl(u) = [ du_z/deta - du_y/dzeta
221  du_x/dzeta - du_z/dxi
222  du_y/dxi - du_x/deta ]
223  */
224  curl_phi[i][p].slice(0) = inv_jac*( dx_dxi*( dphi_deta[i][p].slice(2) -
225  dphi_dzeta[i][p].slice(1) ) +
226  dx_deta*( dphi_dzeta[i][p].slice(0) -
227  dphi_dxi[i][p].slice(2) ) +
228  dx_dzeta*( dphi_dxi[i][p].slice(1) -
229  dphi_deta[i][p].slice(0) ) );
230 
231  curl_phi[i][p].slice(1) = inv_jac*( dy_dxi*( dphi_deta[i][p].slice(2) -
232  dphi_dzeta[i][p].slice(1) ) +
233  dy_deta*( dphi_dzeta[i][p].slice(0)-
234  dphi_dxi[i][p].slice(2) ) +
235  dy_dzeta*( dphi_dxi[i][p].slice(1) -
236  dphi_deta[i][p].slice(0) ) );
237 
238  curl_phi[i][p].slice(2) = inv_jac*( dz_dxi*( dphi_deta[i][p].slice(2) -
239  dphi_dzeta[i][p].slice(1) ) +
240  dz_deta*( dphi_dzeta[i][p].slice(0) -
241  dphi_dxi[i][p].slice(2) ) +
242  dz_dzeta*( dphi_dxi[i][p].slice(1) -
243  dphi_deta[i][p].slice(0) ) );
244  }
245 
246  break;
247  }
248 
249  default:
250  libmesh_error_msg("Invalid dim = " << dim);
251  } // switch(dim)
252 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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 285 of file hcurl_fe_transformation.C.

290 {
291  libmesh_error_msg("HCurl transformations only make sense for vector-valued elements.");
292 }
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
inlinevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 95 of file hcurl_fe_transformation.h.

References libMesh::HCurlFETransformation< OutputShape >::map_curl().

105  {
106  libmesh_warning("WARNING: Shape function Hessians for HCurl elements are not currently being computed!");
107  }
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
inlinevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 127 of file hcurl_fe_transformation.h.

132  {
133  libmesh_warning("WARNING: Shape function divergences for HCurl elements are not currently being computed!");
134  }
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
inlinevirtual

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

Implements libMesh::FETransformationBase< OutputShape >.

Definition at line 78 of file hcurl_fe_transformation.h.

86  {
87  libmesh_warning("WARNING: Shape function gradients for HCurl elements are not currently being computed!");
88  }
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 
) const
virtual

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 52 of file hcurl_fe_transformation.C.

References libMesh::FEMap::get_detadx(), libMesh::FEMap::get_detady(), libMesh::FEMap::get_detadz(), libMesh::FEMap::get_dxidx(), libMesh::FEMap::get_dxidy(), libMesh::FEMap::get_dxidz(), libMesh::FEMap::get_dzetadx(), libMesh::FEMap::get_dzetady(), libMesh::FEMap::get_dzetadz(), libMesh::FEAbstract::get_fe_map(), and libMesh::FEAbstract::get_fe_type().

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

57 {
58  switch (dim)
59  {
60  case 0:
61  case 1:
62  libmesh_error_msg("These element transformations only make sense in 2D and 3D.");
63 
64  case 2:
65  {
66  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
67  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
68 
69  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
70  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
71 
72  // FIXME: Need to update for 2D elements in 3D space
73  /* phi = (dx/dxi)^-T * \hat{phi}
74  In 2D:
75  (dx/dxi)^{-1} = [ dxi/dx dxi/dy
76  deta/dx deta/dy ]
77 
78  so: dxi/dx^{-T} * \hat{phi} = [ dxi/dx deta/dx [ \hat{phi}_xi
79  dxi/dy deta/dy ] \hat{phi}_eta ]
80 
81  or in indicial notation: phi_j = xi_{i,j}*\hat{phi}_i */
82 
83  for (std::size_t i=0; i<phi.size(); i++)
84  for (std::size_t p=0; p<phi[i].size(); p++)
85  {
86  // Need to temporarily cache reference shape functions
87  // TODO: PB: Might be worth trying to build phi_ref separately to see
88  // if we can get vectorization
89  OutputShape phi_ref;
90  FEInterface::shape<OutputShape>(2, fe.get_fe_type(), elem, i, qp[p], phi_ref);
91 
92  phi[i][p](0) = dxidx_map[p]*phi_ref.slice(0) + detadx_map[p]*phi_ref.slice(1);
93 
94  phi[i][p](1) = dxidy_map[p]*phi_ref.slice(0) + detady_map[p]*phi_ref.slice(1);
95  }
96 
97  break;
98  }
99 
100  case 3:
101  {
102  const std::vector<Real> & dxidx_map = fe.get_fe_map().get_dxidx();
103  const std::vector<Real> & dxidy_map = fe.get_fe_map().get_dxidy();
104  const std::vector<Real> & dxidz_map = fe.get_fe_map().get_dxidz();
105 
106  const std::vector<Real> & detadx_map = fe.get_fe_map().get_detadx();
107  const std::vector<Real> & detady_map = fe.get_fe_map().get_detady();
108  const std::vector<Real> & detadz_map = fe.get_fe_map().get_detadz();
109 
110  const std::vector<Real> & dzetadx_map = fe.get_fe_map().get_dzetadx();
111  const std::vector<Real> & dzetady_map = fe.get_fe_map().get_dzetady();
112  const std::vector<Real> & dzetadz_map = fe.get_fe_map().get_dzetadz();
113 
114  /* phi = (dx/dxi)^-T * \hat{phi}
115  In 3D:
116  dx/dxi^-1 = [ dxi/dx dxi/dy dxi/dz
117  deta/dx deta/dy deta/dz
118  dzeta/dx dzeta/dy dzeta/dz]
119 
120  so: dxi/dx^-T * \hat{phi} = [ dxi/dx deta/dx dzeta/dx [ \hat{phi}_xi
121  dxi/dy deta/dy dzeta/dy \hat{phi}_eta
122  dxi/dz deta/dz dzeta/dz ] \hat{phi}_zeta ]
123 
124  or in indicial notation: phi_j = xi_{i,j}*\hat{phi}_i */
125 
126  for (std::size_t i=0; i<phi.size(); i++)
127  for (std::size_t p=0; p<phi[i].size(); p++)
128  {
129  // Need to temporarily cache reference shape functions
130  // TODO: PB: Might be worth trying to build phi_ref separately to see
131  // if we can get vectorization
132  OutputShape phi_ref;
133  FEInterface::shape<OutputShape>(3, fe.get_fe_type(), elem, i, qp[p], phi_ref);
134 
135  phi[i][p].slice(0) = dxidx_map[p]*phi_ref.slice(0) + detadx_map[p]*phi_ref.slice(1)
136  + dzetadx_map[p]*phi_ref.slice(2);
137 
138  phi[i][p].slice(1) = dxidy_map[p]*phi_ref.slice(0) + detady_map[p]*phi_ref.slice(1)
139  + dzetady_map[p]*phi_ref.slice(2);
140 
141  phi[i][p].slice(2) = dxidz_map[p]*phi_ref.slice(0) + detadz_map[p]*phi_ref.slice(1)
142  + dzetadz_map[p]*phi_ref.slice(2);
143  }
144  break;
145  }
146 
147  default:
148  libmesh_error_msg("Invalid dim = " << dim);
149  } // switch(dim)
150 }
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 > > &   
) const

Definition at line 275 of file hcurl_fe_transformation.C.

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

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