libMesh::FEMap Class Reference

Computes finite element mapping function values, gradients, etc. More...

#include <fe_map.h>

Inheritance diagram for libMesh::FEMap:

Public Member Functions

 FEMap ()
 
virtual ~FEMap ()
 
template<unsigned int Dim>
void init_reference_to_physical_map (const std::vector< Point > &qp, const Elem *elem)
 
void compute_single_point_map (const unsigned int dim, const std::vector< Real > &qw, const Elem *elem, unsigned int p, const std::vector< const Node * > &elem_nodes, bool compute_second_derivatives)
 
virtual void compute_affine_map (const unsigned int dim, const std::vector< Real > &qw, const Elem *elem)
 
virtual void compute_null_map (const unsigned int dim, const std::vector< Real > &qw)
 
virtual void compute_map (const unsigned int dim, const std::vector< Real > &qw, const Elem *elem, bool calculate_d2phi)
 
virtual void compute_face_map (int dim, const std::vector< Real > &qw, const Elem *side)
 
void compute_edge_map (int dim, const std::vector< Real > &qw, const Elem *side)
 
template<unsigned int Dim>
void init_face_shape_functions (const std::vector< Point > &qp, const Elem *side)
 
template<unsigned int Dim>
void init_edge_shape_functions (const std::vector< Point > &qp, const Elem *edge)
 
const std::vector< Point > & get_xyz () const
 
const std::vector< Real > & get_jacobian () const
 
const std::vector< Real > & get_JxW () const
 
const std::vector< RealGradient > & get_dxyzdxi () const
 
const std::vector< RealGradient > & get_dxyzdeta () const
 
const std::vector< RealGradient > & get_dxyzdzeta () const
 
const std::vector< RealGradient > & get_d2xyzdxi2 () const
 
const std::vector< RealGradient > & get_d2xyzdeta2 () const
 
const std::vector< RealGradient > & get_d2xyzdzeta2 () const
 
const std::vector< RealGradient > & get_d2xyzdxideta () const
 
const std::vector< RealGradient > & get_d2xyzdxidzeta () const
 
const std::vector< RealGradient > & get_d2xyzdetadzeta () const
 
const std::vector< Real > & get_dxidx () const
 
const std::vector< Real > & get_dxidy () const
 
const std::vector< Real > & get_dxidz () const
 
const std::vector< Real > & get_detadx () const
 
const std::vector< Real > & get_detady () const
 
const std::vector< Real > & get_detadz () const
 
const std::vector< Real > & get_dzetadx () const
 
const std::vector< Real > & get_dzetady () const
 
const std::vector< Real > & get_dzetadz () const
 
const std::vector< std::vector< Real > > & get_d2xidxyz2 () const
 
const std::vector< std::vector< Real > > & get_d2etadxyz2 () const
 
const std::vector< std::vector< Real > > & get_d2zetadxyz2 () const
 
const std::vector< std::vector< Real > > & get_psi () const
 
const std::vector< std::vector< Real > > & get_phi_map () const
 
const std::vector< std::vector< Real > > & get_dphidxi_map () const
 
const std::vector< std::vector< Real > > & get_dphideta_map () const
 
const std::vector< std::vector< Real > > & get_dphidzeta_map () const
 
const std::vector< std::vector< Point > > & get_tangents () const
 
const std::vector< Point > & get_normals () const
 
const std::vector< Real > & get_curvatures () const
 
void print_JxW (std::ostream &os) const
 
void print_xyz (std::ostream &os) const
 
std::vector< std::vector< Real > > & get_psi ()
 
std::vector< std::vector< Real > > & get_dpsidxi ()
 
std::vector< std::vector< Real > > & get_dpsideta ()
 
std::vector< std::vector< Real > > & get_d2psidxi2 ()
 
std::vector< std::vector< Real > > & get_d2psidxideta ()
 
std::vector< std::vector< Real > > & get_d2psideta2 ()
 
std::vector< std::vector< Real > > & get_phi_map ()
 
std::vector< std::vector< Real > > & get_dphidxi_map ()
 
std::vector< std::vector< Real > > & get_dphideta_map ()
 
std::vector< std::vector< Real > > & get_dphidzeta_map ()
 
std::vector< std::vector< Real > > & get_d2phidxi2_map ()
 
std::vector< std::vector< Real > > & get_d2phidxideta_map ()
 
std::vector< std::vector< Real > > & get_d2phidxidzeta_map ()
 
std::vector< std::vector< Real > > & get_d2phideta2_map ()
 
std::vector< std::vector< Real > > & get_d2phidetadzeta_map ()
 
std::vector< std::vector< Real > > & get_d2phidzeta2_map ()
 
std::vector< Real > & get_JxW ()
 

Static Public Member Functions

static UniquePtr< FEMapbuild (FEType fe_type)
 

Protected Member Functions

void determine_calculations ()
 
void resize_quadrature_map_vectors (const unsigned int dim, unsigned int n_qp)
 
Real dxdxi_map (const unsigned int p) const
 
Real dydxi_map (const unsigned int p) const
 
Real dzdxi_map (const unsigned int p) const
 
Real dxdeta_map (const unsigned int p) const
 
Real dydeta_map (const unsigned int p) const
 
Real dzdeta_map (const unsigned int p) const
 
Real dxdzeta_map (const unsigned int p) const
 
Real dydzeta_map (const unsigned int p) const
 
Real dzdzeta_map (const unsigned int p) const
 

Protected Attributes

std::vector< Pointxyz
 
std::vector< RealGradientdxyzdxi_map
 
std::vector< RealGradientdxyzdeta_map
 
std::vector< RealGradientdxyzdzeta_map
 
std::vector< RealGradientd2xyzdxi2_map
 
std::vector< RealGradientd2xyzdxideta_map
 
std::vector< RealGradientd2xyzdeta2_map
 
std::vector< RealGradientd2xyzdxidzeta_map
 
std::vector< RealGradientd2xyzdetadzeta_map
 
std::vector< RealGradientd2xyzdzeta2_map
 
std::vector< Realdxidx_map
 
std::vector< Realdxidy_map
 
std::vector< Realdxidz_map
 
std::vector< Realdetadx_map
 
std::vector< Realdetady_map
 
std::vector< Realdetadz_map
 
std::vector< Realdzetadx_map
 
std::vector< Realdzetady_map
 
std::vector< Realdzetadz_map
 
std::vector< std::vector< Real > > d2xidxyz2_map
 
std::vector< std::vector< Real > > d2etadxyz2_map
 
std::vector< std::vector< Real > > d2zetadxyz2_map
 
std::vector< std::vector< Real > > phi_map
 
std::vector< std::vector< Real > > dphidxi_map
 
std::vector< std::vector< Real > > dphideta_map
 
std::vector< std::vector< Real > > dphidzeta_map
 
std::vector< std::vector< Real > > d2phidxi2_map
 
std::vector< std::vector< Real > > d2phidxideta_map
 
std::vector< std::vector< Real > > d2phidxidzeta_map
 
std::vector< std::vector< Real > > d2phideta2_map
 
std::vector< std::vector< Real > > d2phidetadzeta_map
 
std::vector< std::vector< Real > > d2phidzeta2_map
 
std::vector< std::vector< Real > > psi_map
 
std::vector< std::vector< Real > > dpsidxi_map
 
std::vector< std::vector< Real > > dpsideta_map
 
std::vector< std::vector< Real > > d2psidxi2_map
 
std::vector< std::vector< Real > > d2psidxideta_map
 
std::vector< std::vector< Real > > d2psideta2_map
 
std::vector< std::vector< Point > > tangents
 
std::vector< Pointnormals
 
std::vector< Realcurvatures
 
std::vector< Realjac
 
std::vector< RealJxW
 
bool calculations_started
 
bool calculate_xyz
 
bool calculate_dxyz
 
bool calculate_d2xyz
 

Private Member Functions

void compute_inverse_map_second_derivs (unsigned p)
 

Private Attributes

std::vector< const Node * > elem_nodes
 

Friends

template<unsigned int Dim, FEFamily T>
class FE
 

Detailed Description

Computes finite element mapping function values, gradients, etc.

Class contained in FE that encapsulates mapping (i.e. from physical space to reference space and vice-versa) quantities and operations.

Author
Paul Bauman
Date
2012

Definition at line 45 of file fe_map.h.

Constructor & Destructor Documentation

libMesh::FEMap::FEMap ( )

Definition at line 41 of file fe_map.C.

Referenced by build().

41  :
42  calculations_started(false),
43  calculate_xyz(false),
44  calculate_dxyz(false),
45  calculate_d2xyz(false)
46 {}
bool calculate_dxyz
Definition: fe_map.h:884
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
bool calculate_xyz
Definition: fe_map.h:879

Member Function Documentation

UniquePtr< FEMap > libMesh::FEMap::build ( FEType  fe_type)
static

Definition at line 50 of file fe_map.C.

References libMesh::FEType::family, FEMap(), and libMesh::XYZ.

Referenced by ~FEMap().

51 {
52  switch( fe_type.family )
53  {
54  case XYZ:
55  return UniquePtr<FEMap>(new FEXYZMap);
56 
57  default:
58  return UniquePtr<FEMap>(new FEMap);
59  }
60 
61  libmesh_error_msg("We'll never get here!");
62  return UniquePtr<FEMap>();
63 }
void libMesh::FEMap::compute_affine_map ( const unsigned int  dim,
const std::vector< Real > &  qw,
const Elem elem 
)
virtual

Compute the jacobian and some other additional data fields. Takes the integration weights as input, along with a pointer to the element. The element is assumed to have a constant Jacobian

Definition at line 1223 of file fe_map.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, compute_single_point_map(), d2xyzdeta2_map, d2xyzdetadzeta_map, d2xyzdxi2_map, d2xyzdxideta_map, d2xyzdxidzeta_map, d2xyzdzeta2_map, detadx_map, detady_map, detadz_map, dxidx_map, dxidy_map, dxidz_map, dxyzdeta_map, dxyzdxi_map, dxyzdzeta_map, dzetadx_map, dzetady_map, dzetadz_map, elem_nodes, jac, JxW, libMesh::libmesh_assert(), n_nodes, libMesh::Elem::n_nodes(), libMesh::Elem::node_ptr(), phi_map, resize_quadrature_map_vectors(), and xyz.

Referenced by compute_map(), and ~FEMap().

1226 {
1227  // Start logging the map computation.
1228  LOG_SCOPE("compute_affine_map()", "FEMap");
1229 
1230  libmesh_assert(elem);
1231 
1232  const unsigned int n_qp = cast_int<unsigned int>(qw.size());
1233 
1234  // Resize the vectors to hold data at the quadrature points
1235  this->resize_quadrature_map_vectors(dim, n_qp);
1236 
1237  // Determine the nodes contributing to element elem
1238  unsigned int n_nodes = elem->n_nodes();
1239  elem_nodes.resize(elem->n_nodes());
1240  for (unsigned int i=0; i<n_nodes; i++)
1241  elem_nodes[i] = elem->node_ptr(i);
1242 
1243  // Compute map at quadrature point 0
1244  this->compute_single_point_map(dim, qw, elem, 0, elem_nodes, /*compute_second_derivatives=*/false);
1245 
1246  // Compute xyz at all other quadrature points
1247  if (calculate_xyz)
1248  for (unsigned int p=1; p<n_qp; p++)
1249  {
1250  xyz[p].zero();
1251  for (std::size_t i=0; i<phi_map.size(); i++) // sum over the nodes
1252  xyz[p].add_scaled (*elem_nodes[i], phi_map[i][p] );
1253  }
1254 
1255  // Copy other map data from quadrature point 0
1256  if (calculate_dxyz)
1257  for (unsigned int p=1; p<n_qp; p++) // for each extra quadrature point
1258  {
1259  dxyzdxi_map[p] = dxyzdxi_map[0];
1260  dxidx_map[p] = dxidx_map[0];
1261  dxidy_map[p] = dxidy_map[0];
1262  dxidz_map[p] = dxidz_map[0];
1263 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1264  // The map should be affine, so second derivatives are zero
1265  if (calculate_d2xyz)
1266  d2xyzdxi2_map[p] = 0.;
1267 #endif
1268  if (dim > 1)
1269  {
1270  dxyzdeta_map[p] = dxyzdeta_map[0];
1271  detadx_map[p] = detadx_map[0];
1272  detady_map[p] = detady_map[0];
1273  detadz_map[p] = detadz_map[0];
1274 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1275  if (calculate_d2xyz)
1276  {
1277  d2xyzdxideta_map[p] = 0.;
1278  d2xyzdeta2_map[p] = 0.;
1279  }
1280 #endif
1281  if (dim > 2)
1282  {
1283  dxyzdzeta_map[p] = dxyzdzeta_map[0];
1284  dzetadx_map[p] = dzetadx_map[0];
1285  dzetady_map[p] = dzetady_map[0];
1286  dzetadz_map[p] = dzetadz_map[0];
1287 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1288  if (calculate_d2xyz)
1289  {
1290  d2xyzdxidzeta_map[p] = 0.;
1291  d2xyzdetadzeta_map[p] = 0.;
1292  d2xyzdzeta2_map[p] = 0.;
1293  }
1294 #endif
1295  }
1296  }
1297  jac[p] = jac[0];
1298  JxW[p] = JxW[0] / qw[0] * qw[p];
1299  }
1300 }
bool calculate_dxyz
Definition: fe_map.h:884
std::vector< RealGradient > d2xyzdzeta2_map
Definition: fe_map.h:671
std::vector< Real > dzetady_map
Definition: fe_map.h:723
std::vector< Real > dxidz_map
Definition: fe_map.h:691
std::vector< std::vector< Real > > phi_map
Definition: fe_map.h:754
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
libmesh_assert(j)
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
std::vector< Real > dzetadx_map
Definition: fe_map.h:717
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
const dof_id_type n_nodes
Definition: tecplot_io.C:67
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
std::vector< Real > dzetadz_map
Definition: fe_map.h:729
bool calculate_d2xyz
Definition: fe_map.h:889
std::vector< RealGradient > d2xyzdetadzeta_map
Definition: fe_map.h:665
std::vector< Real > dxidx_map
Definition: fe_map.h:679
bool calculate_xyz
Definition: fe_map.h:879
std::vector< Real > dxidy_map
Definition: fe_map.h:685
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector< Point > xyz
Definition: fe_map.h:615
std::vector< const Node * > elem_nodes
Definition: fe_map.h:908
std::vector< Real > detady_map
Definition: fe_map.h:704
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
void resize_quadrature_map_vectors(const unsigned int dim, unsigned int n_qp)
Definition: fe_map.C:1142
std::vector< Real > jac
Definition: fe_map.h:863
std::vector< Real > detadz_map
Definition: fe_map.h:710
void compute_single_point_map(const unsigned int dim, const std::vector< Real > &qw, const Elem *elem, unsigned int p, const std::vector< const Node * > &elem_nodes, bool compute_second_derivatives)
Definition: fe_map.C:403
std::vector< Real > detadx_map
Definition: fe_map.h:698
std::vector< RealGradient > d2xyzdxidzeta_map
Definition: fe_map.h:659
void libMesh::FEMap::compute_edge_map ( int  dim,
const std::vector< Real > &  qw,
const Elem side 
)

Same as before, but for an edge. Useful for some projections.

Definition at line 921 of file fe_boundary.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, compute_face_map(), curvatures, d2psidxi2_map, d2xyzdeta2_map, d2xyzdxi2_map, d2xyzdxideta_map, determine_calculations(), dpsidxi_map, dxdxi_map(), dxyzdeta_map, dxyzdxi_map, dydxi_map(), dzdxi_map(), libMesh::FE< Dim, T >::edge_reinit(), init_edge_shape_functions(), init_face_shape_functions(), JxW, libMesh::libmesh_assert(), normals, libMesh::Elem::point(), psi_map, libMesh::Real, libMesh::FE< Dim, T >::reinit(), libMesh::FE< Dim, T >::side_map(), tangents, and xyz.

Referenced by ~FEMap().

924 {
925  libmesh_assert(edge);
926 
927  if (dim == 2)
928  {
929  // A 2D finite element living in either 2D or 3D space.
930  // The edges here are the sides of the element, so the
931  // (misnamed) compute_face_map function does what we want
932  this->compute_face_map(dim, qw, edge);
933  return;
934  }
935 
936  libmesh_assert_equal_to (dim, 3); // 1D is unnecessary and currently unsupported
937 
938  LOG_SCOPE("compute_edge_map()", "FEMap");
939 
940  // We're calculating now!
941  this->determine_calculations();
942 
943  // The number of quadrature points.
944  const unsigned int n_qp = cast_int<unsigned int>(qw.size());
945 
946  // Resize the vectors to hold data at the quadrature points
947  if (calculate_xyz)
948  this->xyz.resize(n_qp);
949  if (calculate_dxyz)
950  {
951  this->dxyzdxi_map.resize(n_qp);
952  this->dxyzdeta_map.resize(n_qp);
953  this->tangents.resize(n_qp);
954  this->normals.resize(n_qp);
955  this->JxW.resize(n_qp);
956  }
957  if (calculate_d2xyz)
958  {
959  this->d2xyzdxi2_map.resize(n_qp);
960  this->d2xyzdxideta_map.resize(n_qp);
961  this->d2xyzdeta2_map.resize(n_qp);
962  this->curvatures.resize(n_qp);
963  }
964 
965  // Clear the entities that will be summed
966  for (unsigned int p=0; p<n_qp; p++)
967  {
968  if (calculate_xyz)
969  this->xyz[p].zero();
970  if (calculate_dxyz)
971  {
972  this->tangents[p].resize(1);
973  this->dxyzdxi_map[p].zero();
974  this->dxyzdeta_map[p].zero();
975  }
976  if (calculate_d2xyz)
977  {
978  this->d2xyzdxi2_map[p].zero();
979  this->d2xyzdxideta_map[p].zero();
980  this->d2xyzdeta2_map[p].zero();
981  }
982  }
983 
984  // compute x, dxdxi at the quadrature points
985  for (std::size_t i=0; i<this->psi_map.size(); i++) // sum over the nodes
986  {
987  const Point & edge_point = edge->point(i);
988 
989  for (unsigned int p=0; p<n_qp; p++) // for each quadrature point...
990  {
991  if (calculate_xyz)
992  this->xyz[p].add_scaled (edge_point, this->psi_map[i][p]);
993  if (calculate_dxyz)
994  this->dxyzdxi_map[p].add_scaled (edge_point, this->dpsidxi_map[i][p]);
995  if (calculate_d2xyz)
996  this->d2xyzdxi2_map[p].add_scaled (edge_point, this->d2psidxi2_map[i][p]);
997  }
998  }
999 
1000  // Compute the tangents at the quadrature point
1001  // FIXME: normals (plural!) and curvatures are uncalculated
1002  if (calculate_dxyz)
1003  for (unsigned int p=0; p<n_qp; p++)
1004  {
1005  const Point n = this->dxyzdxi_map[p].cross(this->dxyzdeta_map[p]);
1006  this->tangents[p][0] = this->dxyzdxi_map[p].unit();
1007 
1008  // compute the jacobian at the quadrature points
1009  const Real the_jac = std::sqrt(this->dxdxi_map(p)*this->dxdxi_map(p) +
1010  this->dydxi_map(p)*this->dydxi_map(p) +
1011  this->dzdxi_map(p)*this->dzdxi_map(p));
1012 
1013  libmesh_assert_greater (the_jac, 0.);
1014 
1015  this->JxW[p] = the_jac*qw[p];
1016  }
1017 }
std::vector< std::vector< Real > > psi_map
Definition: fe_map.h:808
bool calculate_dxyz
Definition: fe_map.h:884
Real dxdxi_map(const unsigned int p) const
Definition: fe_map.h:546
Real dydxi_map(const unsigned int p) const
Definition: fe_map.h:554
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
libmesh_assert(j)
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
bool calculate_d2xyz
Definition: fe_map.h:889
Real dzdxi_map(const unsigned int p) const
Definition: fe_map.h:562
bool calculate_xyz
Definition: fe_map.h:879
std::vector< Real > curvatures
Definition: fe_map.h:858
std::vector< std::vector< Real > > d2psidxi2_map
Definition: fe_map.h:827
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector< Point > normals
Definition: fe_map.h:851
std::vector< Point > xyz
Definition: fe_map.h:615
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< std::vector< Real > > dpsidxi_map
Definition: fe_map.h:814
std::vector< std::vector< Point > > tangents
Definition: fe_map.h:846
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
void determine_calculations()
Definition: fe_map.h:524
virtual void compute_face_map(int dim, const std::vector< Real > &qw, const Elem *side)
Definition: fe_boundary.C:573
void libMesh::FEMap::compute_face_map ( int  dim,
const std::vector< Real > &  qw,
const Elem side 
)
virtual

Same as compute_map, but for a side. Useful for boundary integration.

Reimplemented in libMesh::FEXYZMap.

Definition at line 573 of file fe_boundary.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, libMesh::TypeVector< T >::cross(), curvatures, d2psideta2_map, d2psidxi2_map, d2psidxideta_map, d2xyzdeta2_map, d2xyzdxi2_map, d2xyzdxideta_map, libMesh::Elem::default_order(), determine_calculations(), dpsideta_map, dpsidxi_map, dxdeta_map(), dxdxi_map(), dxyzdeta_map, dxyzdxi_map, dydeta_map(), dydxi_map(), dzdeta_map(), dzdxi_map(), libMesh::FE< Dim, T >::inverse_map(), JxW, libMesh::libmesh_assert(), libMesh::FE< Dim, T >::map_eta(), libMesh::FE< Dim, T >::map_xi(), libMesh::FE< Dim, T >::n_shape_functions(), libMesh::Elem::node_id(), normals, libMesh::Elem::parent(), libMesh::Elem::point(), psi_map, libMesh::Real, tangents, libMesh::Elem::type(), libMesh::TypeVector< T >::unit(), and xyz.

Referenced by compute_edge_map(), and ~FEMap().

575 {
577 
578  // We're calculating now!
579  this->determine_calculations();
580 
581  LOG_SCOPE("compute_face_map()", "FEMap");
582 
583  // The number of quadrature points.
584  const unsigned int n_qp = cast_int<unsigned int>(qw.size());
585 
586  switch (dim)
587  {
588  case 1:
589  {
590  // A 1D finite element, currently assumed to be in 1D space
591  // This means the boundary is a "0D finite element", a
592  // NODEELEM.
593 
594  // Resize the vectors to hold data at the quadrature points
595  {
596  if (calculate_xyz)
597  this->xyz.resize(n_qp);
598  if (calculate_dxyz)
599  normals.resize(n_qp);
600 
601  if (calculate_dxyz)
602  this->JxW.resize(n_qp);
603  }
604 
605  // If we have no quadrature points, there's nothing else to do
606  if (!n_qp)
607  break;
608 
609  // We need to look back at the full edge to figure out the normal
610  // vector
611  const Elem * elem = side->parent();
612  libmesh_assert (elem);
613  if (calculate_dxyz)
614  {
615  if (side->node_id(0) == elem->node_id(0))
616  normals[0] = Point(-1.);
617  else
618  {
619  libmesh_assert_equal_to (side->node_id(0),
620  elem->node_id(1));
621  normals[0] = Point(1.);
622  }
623  }
624 
625  // Calculate x at the point
626  if (calculate_xyz)
627  libmesh_assert_equal_to (this->psi_map.size(), 1);
628  // In the unlikely event we have multiple quadrature
629  // points, they'll be in the same place
630  for (unsigned int p=0; p<n_qp; p++)
631  {
632  if (calculate_xyz)
633  {
634  this->xyz[p].zero();
635  this->xyz[p].add_scaled (side->point(0), this->psi_map[0][p]);
636  }
637  if (calculate_dxyz)
638  {
639  normals[p] = normals[0];
640  this->JxW[p] = 1.0*qw[p];
641  }
642  }
643 
644  // done computing the map
645  break;
646  }
647 
648  case 2:
649  {
650  // A 2D finite element living in either 2D or 3D space.
651  // This means the boundary is a 1D finite element, i.e.
652  // and EDGE2 or EDGE3.
653  // Resize the vectors to hold data at the quadrature points
654  {
655  if (calculate_xyz)
656  this->xyz.resize(n_qp);
657  if (calculate_dxyz)
658  {
659  this->dxyzdxi_map.resize(n_qp);
660  this->tangents.resize(n_qp);
661  this->normals.resize(n_qp);
662 
663  this->JxW.resize(n_qp);
664  }
665 
666  if (calculate_d2xyz)
667  {
668  this->d2xyzdxi2_map.resize(n_qp);
669  this->curvatures.resize(n_qp);
670  }
671  }
672 
673  // Clear the entities that will be summed
674  // Compute the tangent & normal at the quadrature point
675  for (unsigned int p=0; p<n_qp; p++)
676  {
677  if (calculate_xyz)
678  this->xyz[p].zero();
679  if (calculate_dxyz)
680  {
681  this->tangents[p].resize(LIBMESH_DIM-1); // 1 Tangent in 2D, 2 in 3D
682  this->dxyzdxi_map[p].zero();
683  }
684  if (calculate_d2xyz)
685  this->d2xyzdxi2_map[p].zero();
686  }
687 
688  const unsigned int n_mapping_shape_functions =
690  side->default_order());
691 
692  // compute x, dxdxi at the quadrature points
693  for (unsigned int i=0; i<n_mapping_shape_functions; i++) // sum over the nodes
694  {
695  const Point & side_point = side->point(i);
696 
697  for (unsigned int p=0; p<n_qp; p++) // for each quadrature point...
698  {
699  if (calculate_xyz)
700  this->xyz[p].add_scaled (side_point, this->psi_map[i][p]);
701  if (calculate_dxyz)
702  this->dxyzdxi_map[p].add_scaled (side_point, this->dpsidxi_map[i][p]);
703  if (calculate_d2xyz)
704  this->d2xyzdxi2_map[p].add_scaled(side_point, this->d2psidxi2_map[i][p]);
705  }
706  }
707 
708  // Compute the tangent & normal at the quadrature point
709  if (calculate_dxyz)
710  {
711  for (unsigned int p=0; p<n_qp; p++)
712  {
713  // The first tangent comes from just the edge's Jacobian
714  this->tangents[p][0] = this->dxyzdxi_map[p].unit();
715 
716 #if LIBMESH_DIM == 2
717  // For a 2D element living in 2D, the normal is given directly
718  // from the entries in the edge Jacobian.
719  this->normals[p] = (Point(this->dxyzdxi_map[p](1), -this->dxyzdxi_map[p](0), 0.)).unit();
720 
721 #elif LIBMESH_DIM == 3
722  // For a 2D element living in 3D, there is a second tangent.
723  // For the second tangent, we need to refer to the full
724  // element's (not just the edge's) Jacobian.
725  const Elem * elem = side->parent();
726  libmesh_assert(elem);
727 
728  // Inverse map xyz[p] to a reference point on the parent...
729  Point reference_point = FE<2,LAGRANGE>::inverse_map(elem, this->xyz[p]);
730 
731  // Get dxyz/dxi and dxyz/deta from the parent map.
732  Point dx_dxi = FE<2,LAGRANGE>::map_xi (elem, reference_point);
733  Point dx_deta = FE<2,LAGRANGE>::map_eta(elem, reference_point);
734 
735  // The second tangent vector is formed by crossing these vectors.
736  tangents[p][1] = dx_dxi.cross(dx_deta).unit();
737 
738  // Finally, the normal in this case is given by crossing these
739  // two tangents.
740  normals[p] = tangents[p][0].cross(tangents[p][1]).unit();
741 #endif
742 
743 
744  // The curvature is computed via the familiar Frenet formula:
745  // curvature = [d^2(x) / d (xi)^2] dot [normal]
746  // For a reference, see:
747  // F.S. Merritt, Mathematics Manual, 1962, McGraw-Hill, p. 310
748  //
749  // Note: The sign convention here is different from the
750  // 3D case. Concave-upward curves (smiles) have a positive
751  // curvature. Concave-downward curves (frowns) have a
752  // negative curvature. Be sure to take that into account!
753  if (calculate_d2xyz)
754  {
755  const Real numerator = this->d2xyzdxi2_map[p] * this->normals[p];
756  const Real denominator = this->dxyzdxi_map[p].norm_sq();
757  libmesh_assert_not_equal_to (denominator, 0);
758  curvatures[p] = numerator / denominator;
759  }
760  }
761 
762  // compute the jacobian at the quadrature points
763  for (unsigned int p=0; p<n_qp; p++)
764  {
765  const Real the_jac = this->dxyzdxi_map[p].norm();
766 
767  libmesh_assert_greater (the_jac, 0.);
768 
769  this->JxW[p] = the_jac*qw[p];
770  }
771  }
772 
773  // done computing the map
774  break;
775  }
776 
777 
778 
779  case 3:
780  {
781  // A 3D finite element living in 3D space.
782  // Resize the vectors to hold data at the quadrature points
783  {
784  if (calculate_xyz)
785  this->xyz.resize(n_qp);
786  if (calculate_dxyz)
787  {
788  this->dxyzdxi_map.resize(n_qp);
789  this->dxyzdeta_map.resize(n_qp);
790  this->tangents.resize(n_qp);
791  this->normals.resize(n_qp);
792  this->JxW.resize(n_qp);
793  }
794  if (calculate_d2xyz)
795  {
796  this->d2xyzdxi2_map.resize(n_qp);
797  this->d2xyzdxideta_map.resize(n_qp);
798  this->d2xyzdeta2_map.resize(n_qp);
799  this->curvatures.resize(n_qp);
800  }
801  }
802 
803  // Clear the entities that will be summed
804  for (unsigned int p=0; p<n_qp; p++)
805  {
806  if (calculate_xyz)
807  this->xyz[p].zero();
808  if (calculate_dxyz)
809  {
810  this->tangents[p].resize(LIBMESH_DIM-1); // 1 Tangent in 2D, 2 in 3D
811  this->dxyzdxi_map[p].zero();
812  this->dxyzdeta_map[p].zero();
813  }
814  if (calculate_d2xyz)
815  {
816  this->d2xyzdxi2_map[p].zero();
817  this->d2xyzdxideta_map[p].zero();
818  this->d2xyzdeta2_map[p].zero();
819  }
820  }
821 
822  const unsigned int n_mapping_shape_functions =
824  side->default_order());
825 
826  // compute x, dxdxi at the quadrature points
827  for (unsigned int i=0; i<n_mapping_shape_functions; i++) // sum over the nodes
828  {
829  const Point & side_point = side->point(i);
830 
831  for (unsigned int p=0; p<n_qp; p++) // for each quadrature point...
832  {
833  if (calculate_xyz)
834  this->xyz[p].add_scaled (side_point, this->psi_map[i][p]);
835  if (calculate_dxyz)
836  {
837  this->dxyzdxi_map[p].add_scaled (side_point, this->dpsidxi_map[i][p]);
838  this->dxyzdeta_map[p].add_scaled(side_point, this->dpsideta_map[i][p]);
839  }
840  if (calculate_d2xyz)
841  {
842  this->d2xyzdxi2_map[p].add_scaled (side_point, this->d2psidxi2_map[i][p]);
843  this->d2xyzdxideta_map[p].add_scaled(side_point, this->d2psidxideta_map[i][p]);
844  this->d2xyzdeta2_map[p].add_scaled (side_point, this->d2psideta2_map[i][p]);
845  }
846  }
847  }
848 
849  // Compute the tangents, normal, and curvature at the quadrature point
850  if (calculate_dxyz)
851  {
852  for (unsigned int p=0; p<n_qp; p++)
853  {
854  const Point n = this->dxyzdxi_map[p].cross(this->dxyzdeta_map[p]);
855  this->normals[p] = n.unit();
856  this->tangents[p][0] = this->dxyzdxi_map[p].unit();
857  this->tangents[p][1] = n.cross(this->dxyzdxi_map[p]).unit();
858 
859  if (calculate_d2xyz)
860  {
861  // Compute curvature using the typical nomenclature
862  // of the first and second fundamental forms.
863  // For reference, see:
864  // 1) http://mathworld.wolfram.com/MeanCurvature.html
865  // (note -- they are using inward normal)
866  // 2) F.S. Merritt, Mathematics Manual, 1962, McGraw-Hill
867  const Real L = -this->d2xyzdxi2_map[p] * this->normals[p];
868  const Real M = -this->d2xyzdxideta_map[p] * this->normals[p];
869  const Real N = -this->d2xyzdeta2_map[p] * this->normals[p];
870  const Real E = this->dxyzdxi_map[p].norm_sq();
871  const Real F = this->dxyzdxi_map[p] * this->dxyzdeta_map[p];
872  const Real G = this->dxyzdeta_map[p].norm_sq();
873 
874  const Real numerator = E*N -2.*F*M + G*L;
875  const Real denominator = E*G - F*F;
876  libmesh_assert_not_equal_to (denominator, 0.);
877  curvatures[p] = 0.5*numerator/denominator;
878  }
879  }
880 
881  // compute the jacobian at the quadrature points, see
882  // http://sp81.msi.umn.edu:999/fluent/fidap/help/theory/thtoc.htm
883  for (unsigned int p=0; p<n_qp; p++)
884  {
885  const Real g11 = (dxdxi_map(p)*dxdxi_map(p) +
886  dydxi_map(p)*dydxi_map(p) +
887  dzdxi_map(p)*dzdxi_map(p));
888 
889  const Real g12 = (dxdxi_map(p)*dxdeta_map(p) +
890  dydxi_map(p)*dydeta_map(p) +
891  dzdxi_map(p)*dzdeta_map(p));
892 
893  const Real g21 = g12;
894 
895  const Real g22 = (dxdeta_map(p)*dxdeta_map(p) +
896  dydeta_map(p)*dydeta_map(p) +
897  dzdeta_map(p)*dzdeta_map(p));
898 
899 
900  const Real the_jac = std::sqrt(g11*g22 - g12*g21);
901 
902  libmesh_assert_greater (the_jac, 0.);
903 
904  this->JxW[p] = the_jac*qw[p];
905  }
906  }
907 
908  // done computing the map
909  break;
910  }
911 
912 
913  default:
914  libmesh_error_msg("Invalid dimension dim = " << dim);
915  }
916 }
std::vector< std::vector< Real > > dpsideta_map
Definition: fe_map.h:820
std::vector< std::vector< Real > > psi_map
Definition: fe_map.h:808
std::vector< std::vector< Real > > d2psideta2_map
Definition: fe_map.h:841
bool calculate_dxyz
Definition: fe_map.h:884
static Point map_xi(const Elem *elem, const Point &reference_point)
Definition: fe_map.C:2002
Real dxdxi_map(const unsigned int p) const
Definition: fe_map.h:546
unsigned short int side
Definition: xdr_io.C:49
Real dydxi_map(const unsigned int p) const
Definition: fe_map.h:554
std::vector< std::vector< Real > > d2psidxideta_map
Definition: fe_map.h:834
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
libmesh_assert(j)
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
bool calculate_d2xyz
Definition: fe_map.h:889
Real dzdxi_map(const unsigned int p) const
Definition: fe_map.h:562
bool calculate_xyz
Definition: fe_map.h:879
Real dxdeta_map(const unsigned int p) const
Definition: fe_map.h:570
Real dydeta_map(const unsigned int p) const
Definition: fe_map.h:578
std::vector< Real > curvatures
Definition: fe_map.h:858
std::vector< std::vector< Real > > d2psidxi2_map
Definition: fe_map.h:827
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector< Point > normals
Definition: fe_map.h:851
std::vector< Point > xyz
Definition: fe_map.h:615
static Point map_eta(const Elem *elem, const Point &reference_point)
Definition: fe_map.C:2030
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< std::vector< Real > > dpsidxi_map
Definition: fe_map.h:814
virtual unsigned int n_shape_functions() const libmesh_override
Definition: fe.C:36
std::vector< std::vector< Point > > tangents
Definition: fe_map.h:846
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
void determine_calculations()
Definition: fe_map.h:524
Real dzdeta_map(const unsigned int p) const
Definition: fe_map.h:586
static Point inverse_map(const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_map.C:1566
void libMesh::FEMap::compute_inverse_map_second_derivs ( unsigned  p)
private

A helper function used by FEMap::compute_single_point_map() to compute second derivatives of the inverse map.

Definition at line 1446 of file fe_map.C.

References A, d2etadxyz2_map, d2xidxyz2_map, d2xyzdeta2_map, d2xyzdetadzeta_map, d2xyzdxi2_map, d2xyzdxideta_map, d2xyzdxidzeta_map, d2xyzdzeta2_map, d2zetadxyz2_map, detadx_map, detady_map, detadz_map, dxidx_map, dxidy_map, dxidz_map, dzetadx_map, dzetady_map, and dzetadz_map.

Referenced by compute_single_point_map().

1447 {
1448 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1449  // Only certain second derivatives are valid depending on the
1450  // dimension...
1451  std::set<unsigned> valid_indices;
1452 
1453  // Construct J^{-1}, A, and B matrices (see JWP's notes for details)
1454  // for cases in which the element dimension matches LIBMESH_DIM.
1455 #if LIBMESH_DIM==1
1456  RealTensor
1457  Jinv(dxidx_map[p], 0., 0.,
1458  0., 0., 0.,
1459  0., 0., 0.),
1460 
1461  A(d2xyzdxi2_map[p](0), 0., 0.,
1462  0., 0., 0.,
1463  0., 0., 0.),
1464 
1465  B(0., 0., 0.,
1466  0., 0., 0.,
1467  0., 0., 0.);
1468 
1470  dxi (dxidx_map[p], 0., 0.),
1471  deta (0., 0., 0.),
1472  dzeta(0., 0., 0.);
1473 
1474  // In 1D, we have only the xx second derivative
1475  valid_indices.insert(0);
1476 
1477 #elif LIBMESH_DIM==2
1478  RealTensor
1479  Jinv(dxidx_map[p], dxidy_map[p], 0.,
1480  detadx_map[p], detady_map[p], 0.,
1481  0., 0., 0.),
1482 
1483  A(d2xyzdxi2_map[p](0), d2xyzdeta2_map[p](0), 0.,
1484  d2xyzdxi2_map[p](1), d2xyzdeta2_map[p](1), 0.,
1485  0., 0., 0.),
1486 
1487  B(d2xyzdxideta_map[p](0), 0., 0.,
1488  d2xyzdxideta_map[p](1), 0., 0.,
1489  0., 0., 0.);
1490 
1492  dxi (dxidx_map[p], dxidy_map[p], 0.),
1493  deta (detadx_map[p], detady_map[p], 0.),
1494  dzeta(0., 0., 0.);
1495 
1496  // In 2D, we have xx, xy, and yy second derivatives
1497  const unsigned tmp[3] = {0,1,3};
1498  valid_indices.insert(tmp, tmp+3);
1499 
1500 #elif LIBMESH_DIM==3
1501  RealTensor
1502  Jinv(dxidx_map[p], dxidy_map[p], dxidz_map[p],
1503  detadx_map[p], detady_map[p], detadz_map[p],
1504  dzetadx_map[p], dzetady_map[p], dzetadz_map[p]),
1505 
1506  A(d2xyzdxi2_map[p](0), d2xyzdeta2_map[p](0), d2xyzdzeta2_map[p](0),
1507  d2xyzdxi2_map[p](1), d2xyzdeta2_map[p](1), d2xyzdzeta2_map[p](1),
1508  d2xyzdxi2_map[p](2), d2xyzdeta2_map[p](2), d2xyzdzeta2_map[p](2)),
1509 
1513 
1515  dxi (dxidx_map[p], dxidy_map[p], dxidz_map[p]),
1516  deta (detadx_map[p], detady_map[p], detadz_map[p]),
1517  dzeta(dzetadx_map[p], dzetady_map[p], dzetadz_map[p]);
1518 
1519  // In 3D, we have xx, xy, xz, yy, yz, and zz second derivatives
1520  const unsigned tmp[6] = {0,1,2,3,4,5};
1521  valid_indices.insert(tmp, tmp+6);
1522 
1523 #endif
1524 
1525  // For (s,t) in {(x,x), (x,y), (x,z), (y,y), (y,z), (z,z)}, compute the
1526  // vector of inverse map second derivatives [xi_{s t}, eta_{s t}, zeta_{s t}]
1527  unsigned ctr=0;
1528  for (unsigned s=0; s<3; ++s)
1529  for (unsigned t=s; t<3; ++t)
1530  {
1531  if (valid_indices.count(ctr))
1532  {
1534  v1(dxi(s)*dxi(t),
1535  deta(s)*deta(t),
1536  dzeta(s)*dzeta(t)),
1537 
1538  v2(dxi(s)*deta(t) + deta(s)*dxi(t),
1539  dxi(s)*dzeta(t) + dzeta(s)*dxi(t),
1540  deta(s)*dzeta(t) + dzeta(s)*deta(t));
1541 
1542  // Compute the inverse map second derivatives
1543  RealVectorValue v3 = -Jinv*(A*v1 + B*v2);
1544 
1545  // Store them in the appropriate locations in the class data structures
1546  d2xidxyz2_map[p][ctr] = v3(0);
1547 
1548  if (LIBMESH_DIM > 1)
1549  d2etadxyz2_map[p][ctr] = v3(1);
1550 
1551  if (LIBMESH_DIM > 2)
1552  d2zetadxyz2_map[p][ctr] = v3(2);
1553  }
1554 
1555  // Increment the counter
1556  ctr++;
1557  }
1558 
1559 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
1560 }
std::vector< std::vector< Real > > d2etadxyz2_map
Definition: fe_map.h:742
VectorValue< Real > RealVectorValue
std::vector< RealGradient > d2xyzdzeta2_map
Definition: fe_map.h:671
std::vector< Real > dzetady_map
Definition: fe_map.h:723
std::vector< std::vector< Real > > d2xidxyz2_map
Definition: fe_map.h:736
RealTensorValue RealTensor
std::vector< Real > dxidz_map
Definition: fe_map.h:691
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
std::vector< Real > dzetadx_map
Definition: fe_map.h:717
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
std::vector< Real > dzetadz_map
Definition: fe_map.h:729
std::vector< std::vector< Real > > d2zetadxyz2_map
Definition: fe_map.h:748
std::vector< RealGradient > d2xyzdetadzeta_map
Definition: fe_map.h:665
std::vector< Real > dxidx_map
Definition: fe_map.h:679
std::vector< Real > dxidy_map
Definition: fe_map.h:685
static PetscErrorCode Mat * A
std::vector< Real > detady_map
Definition: fe_map.h:704
std::vector< Real > detadz_map
Definition: fe_map.h:710
std::vector< Real > detadx_map
Definition: fe_map.h:698
std::vector< RealGradient > d2xyzdxidzeta_map
Definition: fe_map.h:659
void libMesh::FEMap::compute_map ( const unsigned int  dim,
const std::vector< Real > &  qw,
const Elem elem,
bool  calculate_d2phi 
)
virtual

Compute the jacobian and some other additional data fields. Takes the integration weights as input, along with a pointer to the element. Also takes a boolean parameter indicating whether second derivatives need to be calculated, allowing us to potentially skip unnecessary, expensive computations.

Definition at line 1379 of file fe_map.C.

References compute_affine_map(), compute_null_map(), compute_single_point_map(), elem_nodes, libMesh::MeshTools::Subdivision::find_one_ring(), libMesh::Elem::has_affine_map(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::Elem::n_nodes(), libMesh::Elem::node_ptr(), resize_quadrature_map_vectors(), libMesh::TRI3SUBDIVISION, and libMesh::Elem::type().

Referenced by ~FEMap().

1383 {
1384  if (!elem)
1385  {
1386  compute_null_map(dim, qw);
1387  return;
1388  }
1389 
1390  if (elem->has_affine_map())
1391  {
1392  compute_affine_map(dim, qw, elem);
1393  return;
1394  }
1395 
1396  // Start logging the map computation.
1397  LOG_SCOPE("compute_map()", "FEMap");
1398 
1399  libmesh_assert(elem);
1400 
1401  const unsigned int n_qp = cast_int<unsigned int>(qw.size());
1402 
1403  // Resize the vectors to hold data at the quadrature points
1404  this->resize_quadrature_map_vectors(dim, n_qp);
1405 
1406  // Determine the nodes contributing to element elem
1407  std::vector<const Node *> elem_nodes;
1408  if (elem->type() == TRI3SUBDIVISION)
1409  {
1410  // Subdivision surface FE require the 1-ring around elem
1411  libmesh_assert_equal_to (dim, 2);
1412  const Tri3Subdivision * sd_elem = static_cast<const Tri3Subdivision *>(elem);
1413  MeshTools::Subdivision::find_one_ring(sd_elem, elem_nodes);
1414  }
1415  else
1416  {
1417  // All other FE use only the nodes of elem itself
1418  elem_nodes.resize(elem->n_nodes(), libmesh_nullptr);
1419  for (unsigned int i=0; i<elem->n_nodes(); i++)
1420  elem_nodes[i] = elem->node_ptr(i);
1421  }
1422 
1423  // Compute map at all quadrature points
1424  for (unsigned int p=0; p!=n_qp; p++)
1425  this->compute_single_point_map(dim, qw, elem, p, elem_nodes, calculate_d2phi);
1426 }
virtual void compute_affine_map(const unsigned int dim, const std::vector< Real > &qw, const Elem *elem)
Definition: fe_map.C:1223
const class libmesh_nullptr_t libmesh_nullptr
libmesh_assert(j)
void find_one_ring(const Tri3Subdivision *elem, std::vector< const Node * > &nodes)
std::vector< const Node * > elem_nodes
Definition: fe_map.h:908
virtual void compute_null_map(const unsigned int dim, const std::vector< Real > &qw)
Definition: fe_map.C:1304
void resize_quadrature_map_vectors(const unsigned int dim, unsigned int n_qp)
Definition: fe_map.C:1142
void compute_single_point_map(const unsigned int dim, const std::vector< Real > &qw, const Elem *elem, unsigned int p, const std::vector< const Node * > &elem_nodes, bool compute_second_derivatives)
Definition: fe_map.C:403
void libMesh::FEMap::compute_null_map ( const unsigned int  dim,
const std::vector< Real > &  qw 
)
virtual

Assign a fake jacobian and some other additional data fields. Takes the integration weights as input. For use on non-element evaluations.

Definition at line 1304 of file fe_map.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, d2xyzdeta2_map, d2xyzdetadzeta_map, d2xyzdxi2_map, d2xyzdxideta_map, d2xyzdxidzeta_map, d2xyzdzeta2_map, detadx_map, detady_map, detadz_map, dxidx_map, dxidy_map, dxidz_map, dxyzdeta_map, dxyzdxi_map, dxyzdzeta_map, dzetadx_map, dzetady_map, dzetadz_map, jac, JxW, resize_quadrature_map_vectors(), and xyz.

Referenced by compute_map(), and ~FEMap().

1306 {
1307  // Start logging the map computation.
1308  LOG_SCOPE("compute_null_map()", "FEMap");
1309 
1310  const unsigned int n_qp = cast_int<unsigned int>(qw.size());
1311 
1312  // Resize the vectors to hold data at the quadrature points
1313  this->resize_quadrature_map_vectors(dim, n_qp);
1314 
1315  // Compute "fake" xyz
1316  for (unsigned int p=1; p<n_qp; p++)
1317  {
1318  if (calculate_xyz)
1319  xyz[p].zero();
1320 
1321  if (calculate_dxyz)
1322  {
1323  dxyzdxi_map[p] = 0;
1324  dxidx_map[p] = 0;
1325  dxidy_map[p] = 0;
1326  dxidz_map[p] = 0;
1327  }
1328 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1329  if (calculate_d2xyz)
1330  {
1331  d2xyzdxi2_map[p] = 0;
1332  }
1333 #endif
1334  if (dim > 1)
1335  {
1336  if (calculate_dxyz)
1337  {
1338  dxyzdeta_map[p] = 0;
1339  detadx_map[p] = 0;
1340  detady_map[p] = 0;
1341  detadz_map[p] = 0;
1342  }
1343 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1344  if (calculate_d2xyz)
1345  {
1346  d2xyzdxideta_map[p] = 0.;
1347  d2xyzdeta2_map[p] = 0.;
1348  }
1349 #endif
1350  if (dim > 2)
1351  {
1352  if (calculate_dxyz)
1353  {
1354  dxyzdzeta_map[p] = 0;
1355  dzetadx_map[p] = 0;
1356  dzetady_map[p] = 0;
1357  dzetadz_map[p] = 0;
1358  }
1359 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1360  if (calculate_d2xyz)
1361  {
1362  d2xyzdxidzeta_map[p] = 0;
1363  d2xyzdetadzeta_map[p] = 0;
1364  d2xyzdzeta2_map[p] = 0;
1365  }
1366 #endif
1367  }
1368  }
1369  if (calculate_dxyz)
1370  {
1371  jac[p] = 1;
1372  JxW[p] = qw[p];
1373  }
1374  }
1375 }
bool calculate_dxyz
Definition: fe_map.h:884
std::vector< RealGradient > d2xyzdzeta2_map
Definition: fe_map.h:671
std::vector< Real > dzetady_map
Definition: fe_map.h:723
std::vector< Real > dxidz_map
Definition: fe_map.h:691
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
std::vector< Real > dzetadx_map
Definition: fe_map.h:717
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
std::vector< Real > dzetadz_map
Definition: fe_map.h:729
bool calculate_d2xyz
Definition: fe_map.h:889
std::vector< RealGradient > d2xyzdetadzeta_map
Definition: fe_map.h:665
std::vector< Real > dxidx_map
Definition: fe_map.h:679
bool calculate_xyz
Definition: fe_map.h:879
std::vector< Real > dxidy_map
Definition: fe_map.h:685
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector< Point > xyz
Definition: fe_map.h:615
std::vector< Real > detady_map
Definition: fe_map.h:704
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
void resize_quadrature_map_vectors(const unsigned int dim, unsigned int n_qp)
Definition: fe_map.C:1142
std::vector< Real > jac
Definition: fe_map.h:863
std::vector< Real > detadz_map
Definition: fe_map.h:710
std::vector< Real > detadx_map
Definition: fe_map.h:698
std::vector< RealGradient > d2xyzdxidzeta_map
Definition: fe_map.h:659
void libMesh::FEMap::compute_single_point_map ( const unsigned int  dim,
const std::vector< Real > &  qw,
const Elem elem,
unsigned int  p,
const std::vector< const Node * > &  elem_nodes,
bool  compute_second_derivatives 
)

Compute the jacobian and some other additional data fields at the single point with index p. Takes the integration weights as input, along with a pointer to the element and a list of points that contribute to the element. Also takes a boolean flag telling whether second derivatives should actually be computed.

Definition at line 403 of file fe_map.C.

References A, calculate_d2xyz, calculate_dxyz, calculate_xyz, calculations_started, compute_inverse_map_second_derivs(), d2etadxyz2_map, d2phideta2_map, d2phidetadzeta_map, d2phidxi2_map, d2phidxideta_map, d2phidxidzeta_map, d2phidzeta2_map, d2xidxyz2_map, d2xyzdeta2_map, d2xyzdetadzeta_map, d2xyzdxi2_map, d2xyzdxideta_map, d2xyzdxidzeta_map, d2xyzdzeta2_map, d2zetadxyz2_map, detadx_map, detady_map, detadz_map, dphideta_map, dphidxi_map, dphidzeta_map, dxdeta_map(), dxdxi_map(), dxdzeta_map(), dxidx_map, dxidy_map, dxidz_map, dxyzdeta_map, dxyzdxi_map, dxyzdzeta_map, dydeta_map(), dydxi_map(), dydzeta_map(), dzdeta_map(), dzdxi_map(), dzdzeta_map(), dzetadx_map, dzetady_map, dzetadz_map, libMesh::err, libMesh::DofObject::id(), jac, JxW, libMesh::libmesh_assert(), phi_map, libMesh::Elem::print_info(), libMesh::Real, libMesh::DenseMatrix< T >::vector_mult(), and xyz.

Referenced by compute_affine_map(), compute_map(), and ~FEMap().

409 {
410  libmesh_assert(elem);
412 
413  if (calculate_xyz)
414  libmesh_assert_equal_to(phi_map.size(), elem_nodes.size());
415 
416  switch (dim)
417  {
418  //--------------------------------------------------------------------
419  // 0D
420  case 0:
421  {
423  if (calculate_xyz)
424  xyz[p] = *elem_nodes[0];
425  if (calculate_dxyz)
426  {
427  jac[p] = 1.0;
428  JxW[p] = qw[p];
429  }
430  break;
431  }
432 
433  //--------------------------------------------------------------------
434  // 1D
435  case 1:
436  {
437  // Clear the entities that will be summed
438  if (calculate_xyz)
439  xyz[p].zero();
440  if (calculate_dxyz)
441  dxyzdxi_map[p].zero();
442 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
443  if (calculate_d2xyz)
444  {
445  d2xyzdxi2_map[p].zero();
446  // Inverse map second derivatives
447  d2xidxyz2_map[p].assign(6, 0.);
448  }
449 #endif
450 
451  // compute x, dx, d2x at the quadrature point
452  for (std::size_t i=0; i<elem_nodes.size(); i++) // sum over the nodes
453  {
454  // Reference to the point, helps eliminate
455  // exessive temporaries in the inner loop
457  const Point & elem_point = *elem_nodes[i];
458 
459  if (calculate_xyz)
460  xyz[p].add_scaled (elem_point, phi_map[i][p] );
461  if (calculate_dxyz)
462  dxyzdxi_map[p].add_scaled (elem_point, dphidxi_map[i][p]);
463 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
464  if (calculate_d2xyz)
465  d2xyzdxi2_map[p].add_scaled(elem_point, d2phidxi2_map[i][p]);
466 #endif
467  }
468 
469  // Compute the jacobian
470  //
471  // 1D elements can live in 2D or 3D space.
472  // The transformation matrix from local->global
473  // coordinates is
474  //
475  // T = | dx/dxi |
476  // | dy/dxi |
477  // | dz/dxi |
478  //
479  // The generalized determinant of T (from the
480  // so-called "normal" eqns.) is
481  // jac = "det(T)" = sqrt(det(T'T))
482  //
483  // where T'= transpose of T, so
484  //
485  // jac = sqrt( (dx/dxi)^2 + (dy/dxi)^2 + (dz/dxi)^2 )
486 
487  if (calculate_dxyz)
488  {
489  jac[p] = dxyzdxi_map[p].norm();
490 
491  if (jac[p] <= 0.)
492  {
493  // Don't call print_info() recursively if we're already
494  // failing. print_info() calls Elem::volume() which may
495  // call FE::reinit() and trigger the same failure again.
496  static bool failing = false;
497  if (!failing)
498  {
499  failing = true;
500  elem->print_info(libMesh::err);
501  if (calculate_xyz)
502  {
503  libmesh_error_msg("ERROR: negative Jacobian " \
504  << jac[p] \
505  << " at point " \
506  << xyz[p] \
507  << " in element " \
508  << elem->id());
509  }
510  else
511  {
512  // In this case xyz[p] is not defined, so don't
513  // try to print it out.
514  libmesh_error_msg("ERROR: negative Jacobian " \
515  << jac[p] \
516  << " at point index " \
517  << p \
518  << " in element " \
519  << elem->id());
520  }
521  }
522  else
523  {
524  // We were already failing when we called this, so just
525  // stop the current computation and return with
526  // incomplete results.
527  return;
528  }
529  }
530 
531  // The inverse Jacobian entries also come from the
532  // generalized inverse of T (see also the 2D element
533  // living in 3D code).
534  const Real jacm2 = 1./jac[p]/jac[p];
535  dxidx_map[p] = jacm2*dxdxi_map(p);
536 #if LIBMESH_DIM > 1
537  dxidy_map[p] = jacm2*dydxi_map(p);
538 #endif
539 #if LIBMESH_DIM > 2
540  dxidz_map[p] = jacm2*dzdxi_map(p);
541 #endif
542 
543  JxW[p] = jac[p]*qw[p];
544  }
545 
546 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
547 
548  if (calculate_d2xyz)
549  {
550 #if LIBMESH_DIM == 1
551  // Compute inverse map second derivatives for 1D-element-living-in-1D case
553 #elif LIBMESH_DIM == 2
554  // Compute inverse map second derivatives for 1D-element-living-in-2D case
555  // See JWP notes for details
556 
557  // numer = x_xi*x_{xi xi} + y_xi*y_{xi xi}
558  Real numer =
559  dxyzdxi_map[p](0)*d2xyzdxi2_map[p](0) +
560  dxyzdxi_map[p](1)*d2xyzdxi2_map[p](1);
561 
562  // denom = (x_xi)^2 + (y_xi)^2 must be >= 0.0
563  Real denom =
564  dxyzdxi_map[p](0)*dxyzdxi_map[p](0) +
565  dxyzdxi_map[p](1)*dxyzdxi_map[p](1);
566 
567  if (denom <= 0.0)
568  {
569  // Don't call print_info() recursively if we're already
570  // failing. print_info() calls Elem::volume() which may
571  // call FE::reinit() and trigger the same failure again.
572  static bool failing = false;
573  if (!failing)
574  {
575  failing = true;
576  elem->print_info(libMesh::err);
577  libmesh_error_msg("Encountered invalid 1D element!");
578  }
579  else
580  {
581  // We were already failing when we called this, so just
582  // stop the current computation and return with
583  // incomplete results.
584  return;
585  }
586  }
587 
588  // xi_{x x}
589  d2xidxyz2_map[p][0] = -numer * dxidx_map[p]*dxidx_map[p] / denom;
590 
591  // xi_{x y}
592  d2xidxyz2_map[p][1] = -numer * dxidx_map[p]*dxidy_map[p] / denom;
593 
594  // xi_{y y}
595  d2xidxyz2_map[p][3] = -numer * dxidy_map[p]*dxidy_map[p] / denom;
596 
597 #elif LIBMESH_DIM == 3
598  // Compute inverse map second derivatives for 1D-element-living-in-3D case
599  // See JWP notes for details
600 
601  // numer = x_xi*x_{xi xi} + y_xi*y_{xi xi} + z_xi*z_{xi xi}
602  Real numer =
603  dxyzdxi_map[p](0)*d2xyzdxi2_map[p](0) +
604  dxyzdxi_map[p](1)*d2xyzdxi2_map[p](1) +
605  dxyzdxi_map[p](2)*d2xyzdxi2_map[p](2);
606 
607  // denom = (x_xi)^2 + (y_xi)^2 + (z_xi)^2 must be >= 0.0
608  Real denom =
609  dxyzdxi_map[p](0)*dxyzdxi_map[p](0) +
610  dxyzdxi_map[p](1)*dxyzdxi_map[p](1) +
611  dxyzdxi_map[p](2)*dxyzdxi_map[p](2);
612 
613  if (denom <= 0.0)
614  {
615  // Don't call print_info() recursively if we're already
616  // failing. print_info() calls Elem::volume() which may
617  // call FE::reinit() and trigger the same failure again.
618  static bool failing = false;
619  if (!failing)
620  {
621  failing = true;
622  elem->print_info(libMesh::err);
623  libmesh_error_msg("Encountered invalid 1D element!");
624  }
625  else
626  {
627  // We were already failing when we called this, so just
628  // stop the current computation and return with
629  // incomplete results.
630  return;
631  }
632  }
633 
634  // xi_{x x}
635  d2xidxyz2_map[p][0] = -numer * dxidx_map[p]*dxidx_map[p] / denom;
636 
637  // xi_{x y}
638  d2xidxyz2_map[p][1] = -numer * dxidx_map[p]*dxidy_map[p] / denom;
639 
640  // xi_{x z}
641  d2xidxyz2_map[p][2] = -numer * dxidx_map[p]*dxidz_map[p] / denom;
642 
643  // xi_{y y}
644  d2xidxyz2_map[p][3] = -numer * dxidy_map[p]*dxidy_map[p] / denom;
645 
646  // xi_{y z}
647  d2xidxyz2_map[p][4] = -numer * dxidy_map[p]*dxidz_map[p] / denom;
648 
649  // xi_{z z}
650  d2xidxyz2_map[p][5] = -numer * dxidz_map[p]*dxidz_map[p] / denom;
651 #endif
652  } // calculate_d2xyz
653 
654 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
655 
656  // done computing the map
657  break;
658  }
659 
660 
661  //--------------------------------------------------------------------
662  // 2D
663  case 2:
664  {
665  //------------------------------------------------------------------
666  // Compute the (x,y) values at the quadrature points,
667  // the Jacobian at the quadrature points
668 
669  if (calculate_xyz)
670  xyz[p].zero();
671 
672  if (calculate_dxyz)
673  {
674  dxyzdxi_map[p].zero();
675  dxyzdeta_map[p].zero();
676  }
677 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
678  if (calculate_d2xyz)
679  {
680  d2xyzdxi2_map[p].zero();
681  d2xyzdxideta_map[p].zero();
682  d2xyzdeta2_map[p].zero();
683  // Inverse map second derivatives
684  d2xidxyz2_map[p].assign(6, 0.);
685  d2etadxyz2_map[p].assign(6, 0.);
686  }
687 #endif
688 
689 
690  // compute (x,y) at the quadrature points, derivatives once
691  for (std::size_t i=0; i<elem_nodes.size(); i++) // sum over the nodes
692  {
693  // Reference to the point, helps eliminate
694  // exessive temporaries in the inner loop
696  const Point & elem_point = *elem_nodes[i];
697 
698  if (calculate_xyz)
699  xyz[p].add_scaled (elem_point, phi_map[i][p] );
700 
701  if (calculate_dxyz)
702  {
703  dxyzdxi_map[p].add_scaled (elem_point, dphidxi_map[i][p] );
704  dxyzdeta_map[p].add_scaled (elem_point, dphideta_map[i][p]);
705  }
706 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
707  if (calculate_d2xyz)
708  {
709  d2xyzdxi2_map[p].add_scaled (elem_point, d2phidxi2_map[i][p]);
710  d2xyzdxideta_map[p].add_scaled (elem_point, d2phidxideta_map[i][p]);
711  d2xyzdeta2_map[p].add_scaled (elem_point, d2phideta2_map[i][p]);
712  }
713 #endif
714  }
715 
716  if (calculate_dxyz)
717  {
718  // compute the jacobian once
719  const Real dx_dxi = dxdxi_map(p),
720  dx_deta = dxdeta_map(p),
721  dy_dxi = dydxi_map(p),
722  dy_deta = dydeta_map(p);
723 
724 #if LIBMESH_DIM == 2
725  // Compute the Jacobian. This assumes the 2D face
726  // lives in 2D space
727  //
728  // Symbolically, the matrix determinant is
729  //
730  // | dx/dxi dx/deta |
731  // jac = | dy/dxi dy/deta |
732  //
733  // jac = dx/dxi*dy/deta - dx/deta*dy/dxi
734  jac[p] = (dx_dxi*dy_deta - dx_deta*dy_dxi);
735 
736  if (jac[p] <= 0.)
737  {
738  // Don't call print_info() recursively if we're already
739  // failing. print_info() calls Elem::volume() which may
740  // call FE::reinit() and trigger the same failure again.
741  static bool failing = false;
742  if (!failing)
743  {
744  failing = true;
745  elem->print_info(libMesh::err);
746  if (calculate_xyz)
747  {
748  libmesh_error_msg("ERROR: negative Jacobian " \
749  << jac[p] \
750  << " at point " \
751  << xyz[p] \
752  << " in element " \
753  << elem->id());
754  }
755  else
756  {
757  // In this case xyz[p] is not defined, so don't
758  // try to print it out.
759  libmesh_error_msg("ERROR: negative Jacobian " \
760  << jac[p] \
761  << " at point index " \
762  << p \
763  << " in element " \
764  << elem->id());
765  }
766  }
767  else
768  {
769  // We were already failing when we called this, so just
770  // stop the current computation and return with
771  // incomplete results.
772  return;
773  }
774  }
775 
776  JxW[p] = jac[p]*qw[p];
777 
778  // Compute the shape function derivatives wrt x,y at the
779  // quadrature points
780  const Real inv_jac = 1./jac[p];
781 
782  dxidx_map[p] = dy_deta*inv_jac; //dxi/dx = (1/J)*dy/deta
783  dxidy_map[p] = -dx_deta*inv_jac; //dxi/dy = -(1/J)*dx/deta
784  detadx_map[p] = -dy_dxi* inv_jac; //deta/dx = -(1/J)*dy/dxi
785  detady_map[p] = dx_dxi* inv_jac; //deta/dy = (1/J)*dx/dxi
786 
787  dxidz_map[p] = detadz_map[p] = 0.;
788 
789  if (compute_second_derivatives)
791 
792 #else // LIBMESH_DIM == 3
793 
794  const Real dz_dxi = dzdxi_map(p),
795  dz_deta = dzdeta_map(p);
796 
797  // Compute the Jacobian. This assumes a 2D face in
798  // 3D space.
799  //
800  // The transformation matrix T from local to global
801  // coordinates is
802  //
803  // | dx/dxi dx/deta |
804  // T = | dy/dxi dy/deta |
805  // | dz/dxi dz/deta |
806  // note det(T' T) = det(T')det(T) = det(T)det(T)
807  // so det(T) = std::sqrt(det(T' T))
808  //
809  //----------------------------------------------
810  // Notes:
811  //
812  // dX = R dXi -> R'dX = R'R dXi
813  // (R^-1)dX = dXi [(R'R)^-1 R']dX = dXi
814  //
815  // so R^-1 = (R'R)^-1 R'
816  //
817  // and R^-1 R = (R'R)^-1 R'R = I.
818  //
819  const Real g11 = (dx_dxi*dx_dxi +
820  dy_dxi*dy_dxi +
821  dz_dxi*dz_dxi);
822 
823  const Real g12 = (dx_dxi*dx_deta +
824  dy_dxi*dy_deta +
825  dz_dxi*dz_deta);
826 
827  const Real g21 = g12;
828 
829  const Real g22 = (dx_deta*dx_deta +
830  dy_deta*dy_deta +
831  dz_deta*dz_deta);
832 
833  const Real det = (g11*g22 - g12*g21);
834 
835  if (det <= 0.)
836  {
837  // Don't call print_info() recursively if we're already
838  // failing. print_info() calls Elem::volume() which may
839  // call FE::reinit() and trigger the same failure again.
840  static bool failing = false;
841  if (!failing)
842  {
843  failing = true;
844  elem->print_info(libMesh::err);
845  if (calculate_xyz)
846  {
847  libmesh_error_msg("ERROR: negative Jacobian " \
848  << det \
849  << " at point " \
850  << xyz[p] \
851  << " in element " \
852  << elem->id());
853  }
854  else
855  {
856  // In this case xyz[p] is not defined, so don't
857  // try to print it out.
858  libmesh_error_msg("ERROR: negative Jacobian " \
859  << det \
860  << " at point index " \
861  << p \
862  << " in element " \
863  << elem->id());
864  }
865  }
866  else
867  {
868  // We were already failing when we called this, so just
869  // stop the current computation and return with
870  // incomplete results.
871  return;
872  }
873  }
874 
875  const Real inv_det = 1./det;
876  jac[p] = std::sqrt(det);
877 
878  JxW[p] = jac[p]*qw[p];
879 
880  const Real g11inv = g22*inv_det;
881  const Real g12inv = -g12*inv_det;
882  const Real g21inv = -g21*inv_det;
883  const Real g22inv = g11*inv_det;
884 
885  dxidx_map[p] = g11inv*dx_dxi + g12inv*dx_deta;
886  dxidy_map[p] = g11inv*dy_dxi + g12inv*dy_deta;
887  dxidz_map[p] = g11inv*dz_dxi + g12inv*dz_deta;
888 
889  detadx_map[p] = g21inv*dx_dxi + g22inv*dx_deta;
890  detady_map[p] = g21inv*dy_dxi + g22inv*dy_deta;
891  detadz_map[p] = g21inv*dz_dxi + g22inv*dz_deta;
892 
893 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
894 
895  if (calculate_d2xyz)
896  {
897  // Compute inverse map second derivative values for
898  // 2D-element-living-in-3D case. We pursue a least-squares
899  // solution approach for this "non-square" case, see JWP notes
900  // for details.
901 
902  // A = [ x_{xi xi} x_{eta eta} ]
903  // [ y_{xi xi} y_{eta eta} ]
904  // [ z_{xi xi} z_{eta eta} ]
905  DenseMatrix<Real> A(3,2);
906  A(0,0) = d2xyzdxi2_map[p](0); A(0,1) = d2xyzdeta2_map[p](0);
907  A(1,0) = d2xyzdxi2_map[p](1); A(1,1) = d2xyzdeta2_map[p](1);
908  A(2,0) = d2xyzdxi2_map[p](2); A(2,1) = d2xyzdeta2_map[p](2);
909 
910  // J^T, the transpose of the Jacobian matrix
911  DenseMatrix<Real> JT(2,3);
912  JT(0,0) = dx_dxi; JT(0,1) = dy_dxi; JT(0,2) = dz_dxi;
913  JT(1,0) = dx_deta; JT(1,1) = dy_deta; JT(1,2) = dz_deta;
914 
915  // (J^T J)^(-1), this has already been computed for us above...
916  DenseMatrix<Real> JTJinv(2,2);
917  JTJinv(0,0) = g11inv; JTJinv(0,1) = g12inv;
918  JTJinv(1,0) = g21inv; JTJinv(1,1) = g22inv;
919 
920  // Some helper variables
922  dxi (dxidx_map[p], dxidy_map[p], dxidz_map[p]),
923  deta (detadx_map[p], detady_map[p], detadz_map[p]);
924 
925  // To be filled in below
926  DenseVector<Real> tmp1(2);
927  DenseVector<Real> tmp2(3);
928  DenseVector<Real> tmp3(2);
929 
930  // For (s,t) in {(x,x), (x,y), (x,z), (y,y), (y,z), (z,z)}, compute the
931  // vector of inverse map second derivatives [xi_{s t}, eta_{s t}]
932  unsigned ctr=0;
933  for (unsigned s=0; s<3; ++s)
934  for (unsigned t=s; t<3; ++t)
935  {
936  // Construct tmp1 = [xi_s*xi_t, eta_s*eta_t]
937  tmp1(0) = dxi(s)*dxi(t);
938  tmp1(1) = deta(s)*deta(t);
939 
940  // Compute tmp2 = A * tmp1
941  A.vector_mult(tmp2, tmp1);
942 
943  // Compute scalar value "alpha"
944  Real alpha = dxi(s)*deta(t) + deta(s)*dxi(t);
945 
946  // Compute tmp2 <- tmp2 + alpha * x_{xi eta}
947  for (unsigned i=0; i<3; ++i)
948  tmp2(i) += alpha*d2xyzdxideta_map[p](i);
949 
950  // Compute tmp3 = J^T * tmp2
951  JT.vector_mult(tmp3, tmp2);
952 
953  // Compute tmp1 = (J^T J)^(-1) * tmp3. tmp1 is available for us to reuse.
954  JTJinv.vector_mult(tmp1, tmp3);
955 
956  // Fill in appropriate entries, don't forget to multiply by -1!
957  d2xidxyz2_map[p][ctr] = -tmp1(0);
958  d2etadxyz2_map[p][ctr] = -tmp1(1);
959 
960  // Increment the counter
961  ctr++;
962  }
963  }
964 
965 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
966 
967 #endif // LIBMESH_DIM == 3
968  }
969  // done computing the map
970  break;
971  }
972 
973 
974 
975  //--------------------------------------------------------------------
976  // 3D
977  case 3:
978  {
979  //------------------------------------------------------------------
980  // Compute the (x,y,z) values at the quadrature points,
981  // the Jacobian at the quadrature point
982 
983  // Clear the entities that will be summed
984  if (calculate_xyz)
985  xyz[p].zero ();
986  if (calculate_dxyz)
987  {
988  dxyzdxi_map[p].zero ();
989  dxyzdeta_map[p].zero ();
990  dxyzdzeta_map[p].zero ();
991  }
992 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
993  if (calculate_d2xyz)
994  {
995  d2xyzdxi2_map[p].zero();
996  d2xyzdxideta_map[p].zero();
997  d2xyzdxidzeta_map[p].zero();
998  d2xyzdeta2_map[p].zero();
999  d2xyzdetadzeta_map[p].zero();
1000  d2xyzdzeta2_map[p].zero();
1001  // Inverse map second derivatives
1002  d2xidxyz2_map[p].assign(6, 0.);
1003  d2etadxyz2_map[p].assign(6, 0.);
1004  d2zetadxyz2_map[p].assign(6, 0.);
1005  }
1006 #endif
1007 
1008 
1009  // compute (x,y,z) at the quadrature points,
1010  // dxdxi, dydxi, dzdxi,
1011  // dxdeta, dydeta, dzdeta,
1012  // dxdzeta, dydzeta, dzdzeta all once
1013  for (std::size_t i=0; i<elem_nodes.size(); i++) // sum over the nodes
1014  {
1015  // Reference to the point, helps eliminate
1016  // exessive temporaries in the inner loop
1018  const Point & elem_point = *elem_nodes[i];
1019 
1020  if (calculate_xyz)
1021  xyz[p].add_scaled (elem_point, phi_map[i][p] );
1022  if (calculate_dxyz)
1023  {
1024  dxyzdxi_map[p].add_scaled (elem_point, dphidxi_map[i][p] );
1025  dxyzdeta_map[p].add_scaled (elem_point, dphideta_map[i][p] );
1026  dxyzdzeta_map[p].add_scaled (elem_point, dphidzeta_map[i][p]);
1027  }
1028 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1029  if (calculate_d2xyz)
1030  {
1031  d2xyzdxi2_map[p].add_scaled (elem_point,
1032  d2phidxi2_map[i][p]);
1033  d2xyzdxideta_map[p].add_scaled (elem_point,
1034  d2phidxideta_map[i][p]);
1035  d2xyzdxidzeta_map[p].add_scaled (elem_point,
1036  d2phidxidzeta_map[i][p]);
1037  d2xyzdeta2_map[p].add_scaled (elem_point,
1038  d2phideta2_map[i][p]);
1039  d2xyzdetadzeta_map[p].add_scaled (elem_point,
1040  d2phidetadzeta_map[i][p]);
1041  d2xyzdzeta2_map[p].add_scaled (elem_point,
1042  d2phidzeta2_map[i][p]);
1043  }
1044 #endif
1045  }
1046 
1047  if (calculate_dxyz)
1048  {
1049  // compute the jacobian
1050  const Real
1051  dx_dxi = dxdxi_map(p), dy_dxi = dydxi_map(p), dz_dxi = dzdxi_map(p),
1052  dx_deta = dxdeta_map(p), dy_deta = dydeta_map(p), dz_deta = dzdeta_map(p),
1053  dx_dzeta = dxdzeta_map(p), dy_dzeta = dydzeta_map(p), dz_dzeta = dzdzeta_map(p);
1054 
1055  // Symbolically, the matrix determinant is
1056  //
1057  // | dx/dxi dy/dxi dz/dxi |
1058  // jac = | dx/deta dy/deta dz/deta |
1059  // | dx/dzeta dy/dzeta dz/dzeta |
1060  //
1061  // jac = dx/dxi*(dy/deta*dz/dzeta - dz/deta*dy/dzeta) +
1062  // dy/dxi*(dz/deta*dx/dzeta - dx/deta*dz/dzeta) +
1063  // dz/dxi*(dx/deta*dy/dzeta - dy/deta*dx/dzeta)
1064 
1065  jac[p] = (dx_dxi*(dy_deta*dz_dzeta - dz_deta*dy_dzeta) +
1066  dy_dxi*(dz_deta*dx_dzeta - dx_deta*dz_dzeta) +
1067  dz_dxi*(dx_deta*dy_dzeta - dy_deta*dx_dzeta));
1068 
1069  if (jac[p] <= 0.)
1070  {
1071  // Don't call print_info() recursively if we're already
1072  // failing. print_info() calls Elem::volume() which may
1073  // call FE::reinit() and trigger the same failure again.
1074  static bool failing = false;
1075  if (!failing)
1076  {
1077  failing = true;
1078  elem->print_info(libMesh::err);
1079  if (calculate_xyz)
1080  {
1081  libmesh_error_msg("ERROR: negative Jacobian " \
1082  << jac[p] \
1083  << " at point " \
1084  << xyz[p] \
1085  << " in element " \
1086  << elem->id());
1087  }
1088  else
1089  {
1090  // In this case xyz[p] is not defined, so don't
1091  // try to print it out.
1092  libmesh_error_msg("ERROR: negative Jacobian " \
1093  << jac[p] \
1094  << " at point index " \
1095  << p \
1096  << " in element " \
1097  << elem->id());
1098  }
1099  }
1100  else
1101  {
1102  // We were already failing when we called this, so just
1103  // stop the current computation and return with
1104  // incomplete results.
1105  return;
1106  }
1107  }
1108 
1109  JxW[p] = jac[p]*qw[p];
1110 
1111  // Compute the shape function derivatives wrt x,y at the
1112  // quadrature points
1113  const Real inv_jac = 1./jac[p];
1114 
1115  dxidx_map[p] = (dy_deta*dz_dzeta - dz_deta*dy_dzeta)*inv_jac;
1116  dxidy_map[p] = (dz_deta*dx_dzeta - dx_deta*dz_dzeta)*inv_jac;
1117  dxidz_map[p] = (dx_deta*dy_dzeta - dy_deta*dx_dzeta)*inv_jac;
1118 
1119  detadx_map[p] = (dz_dxi*dy_dzeta - dy_dxi*dz_dzeta )*inv_jac;
1120  detady_map[p] = (dx_dxi*dz_dzeta - dz_dxi*dx_dzeta )*inv_jac;
1121  detadz_map[p] = (dy_dxi*dx_dzeta - dx_dxi*dy_dzeta )*inv_jac;
1122 
1123  dzetadx_map[p] = (dy_dxi*dz_deta - dz_dxi*dy_deta )*inv_jac;
1124  dzetady_map[p] = (dz_dxi*dx_deta - dx_dxi*dz_deta )*inv_jac;
1125  dzetadz_map[p] = (dx_dxi*dy_deta - dy_dxi*dx_deta )*inv_jac;
1126  }
1127 
1128  if (compute_second_derivatives)
1130 
1131  // done computing the map
1132  break;
1133  }
1134 
1135  default:
1136  libmesh_error_msg("Invalid dim = " << dim);
1137  }
1138 }
std::vector< std::vector< Real > > dphidzeta_map
Definition: fe_map.h:769
void compute_inverse_map_second_derivs(unsigned p)
Definition: fe_map.C:1446
Real dxdzeta_map(const unsigned int p) const
Definition: fe_map.h:594
std::vector< std::vector< Real > > dphidxi_map
Definition: fe_map.h:759
std::vector< std::vector< Real > > d2etadxyz2_map
Definition: fe_map.h:742
bool calculate_dxyz
Definition: fe_map.h:884
VectorValue< Real > RealVectorValue
std::vector< RealGradient > d2xyzdzeta2_map
Definition: fe_map.h:671
std::vector< Real > dzetady_map
Definition: fe_map.h:723
std::vector< std::vector< Real > > d2xidxyz2_map
Definition: fe_map.h:736
Real dxdxi_map(const unsigned int p) const
Definition: fe_map.h:546
std::vector< std::vector< Real > > d2phideta2_map
Definition: fe_map.h:791
std::vector< std::vector< Real > > d2phidxidzeta_map
Definition: fe_map.h:786
std::vector< Real > dxidz_map
Definition: fe_map.h:691
std::vector< std::vector< Real > > phi_map
Definition: fe_map.h:754
Real dydxi_map(const unsigned int p) const
Definition: fe_map.h:554
std::vector< std::vector< Real > > d2phidetadzeta_map
Definition: fe_map.h:796
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
libmesh_assert(j)
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
std::vector< Real > dzetadx_map
Definition: fe_map.h:717
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
Real dydzeta_map(const unsigned int p) const
Definition: fe_map.h:602
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
std::vector< Real > dzetadz_map
Definition: fe_map.h:729
bool calculate_d2xyz
Definition: fe_map.h:889
std::vector< std::vector< Real > > d2zetadxyz2_map
Definition: fe_map.h:748
Real dzdxi_map(const unsigned int p) const
Definition: fe_map.h:562
bool calculations_started
Definition: fe_map.h:874
std::vector< RealGradient > d2xyzdetadzeta_map
Definition: fe_map.h:665
std::vector< Real > dxidx_map
Definition: fe_map.h:679
bool calculate_xyz
Definition: fe_map.h:879
std::vector< Real > dxidy_map
Definition: fe_map.h:685
OStreamProxy err(std::cerr)
Real dxdeta_map(const unsigned int p) const
Definition: fe_map.h:570
Real dydeta_map(const unsigned int p) const
Definition: fe_map.h:578
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector< std::vector< Real > > d2phidxideta_map
Definition: fe_map.h:781
std::vector< Point > xyz
Definition: fe_map.h:615
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< const Node * > elem_nodes
Definition: fe_map.h:908
static PetscErrorCode Mat * A
std::vector< Real > detady_map
Definition: fe_map.h:704
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
std::vector< Real > jac
Definition: fe_map.h:863
std::vector< Real > detadz_map
Definition: fe_map.h:710
std::vector< std::vector< Real > > d2phidzeta2_map
Definition: fe_map.h:801
std::vector< std::vector< Real > > d2phidxi2_map
Definition: fe_map.h:776
std::vector< Real > detadx_map
Definition: fe_map.h:698
std::vector< std::vector< Real > > dphideta_map
Definition: fe_map.h:764
Real dzdeta_map(const unsigned int p) const
Definition: fe_map.h:586
std::vector< RealGradient > d2xyzdxidzeta_map
Definition: fe_map.h:659
Real dzdzeta_map(const unsigned int p) const
Definition: fe_map.h:610
void libMesh::FEMap::determine_calculations ( )
inlineprotected

Determine which values are to be calculated

Definition at line 524 of file fe_map.h.

References calculate_d2xyz, calculate_dxyz, calculations_started, and resize_quadrature_map_vectors().

Referenced by compute_edge_map(), compute_face_map(), init_edge_shape_functions(), init_face_shape_functions(), init_reference_to_physical_map(), and resize_quadrature_map_vectors().

524  {
525  calculations_started = true;
526 
527 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
528  // Second derivative calculations currently have first derivative
529  // calculations as a prerequisite
530  if (calculate_d2xyz)
531  calculate_dxyz = true;
532 #endif
533  }
bool calculate_dxyz
Definition: fe_map.h:884
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
Real libMesh::FEMap::dxdeta_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The x value of the pth entry of the dxzydeta_map.

Definition at line 570 of file fe_map.h.

References dxyzdeta_map.

Referenced by libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and compute_single_point_map().

570 { return dxyzdeta_map[p](0); }
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
Real libMesh::FEMap::dxdxi_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The x value of the pth entry of the dxzydxi_map.

Definition at line 546 of file fe_map.h.

References dxyzdxi_map.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and compute_single_point_map().

546 { return dxyzdxi_map[p](0); }
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
Real libMesh::FEMap::dxdzeta_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The x value of the pth entry of the dxzydzeta_map.

Definition at line 594 of file fe_map.h.

References dxyzdzeta_map.

Referenced by compute_single_point_map().

594 { return dxyzdzeta_map[p](0); }
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
Real libMesh::FEMap::dydeta_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The y value of the pth entry of the dxzydeta_map.

Definition at line 578 of file fe_map.h.

References dxyzdeta_map.

Referenced by libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and compute_single_point_map().

578 { return dxyzdeta_map[p](1); }
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
Real libMesh::FEMap::dydxi_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The y value of the pth entry of the dxzydxi_map.

Definition at line 554 of file fe_map.h.

References dxyzdxi_map.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and compute_single_point_map().

554 { return dxyzdxi_map[p](1); }
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
Real libMesh::FEMap::dydzeta_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The y value of the pth entry of the dxzydzeta_map.

Definition at line 602 of file fe_map.h.

References dxyzdzeta_map.

Referenced by compute_single_point_map().

602 { return dxyzdzeta_map[p](1); }
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
Real libMesh::FEMap::dzdeta_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The z value of the pth entry of the dxzydeta_map.

Definition at line 586 of file fe_map.h.

References dxyzdeta_map.

Referenced by libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and compute_single_point_map().

586 { return dxyzdeta_map[p](2); }
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
Real libMesh::FEMap::dzdxi_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The z value of the pth entry of the dxzydxi_map.

Definition at line 562 of file fe_map.h.

References dxyzdxi_map.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and compute_single_point_map().

562 { return dxyzdxi_map[p](2); }
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
Real libMesh::FEMap::dzdzeta_map ( const unsigned int  p) const
inlineprotected

Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected.

Returns
The z value of the pth entry of the dxzydzeta_map.

Definition at line 610 of file fe_map.h.

References dxyzdzeta_map.

Referenced by compute_single_point_map().

610 { return dxyzdzeta_map[p](2); }
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
const std::vector<Real>& libMesh::FEMap::get_curvatures ( ) const
inline
Returns
The curvatures for use in face integration.

Definition at line 377 of file fe_map.h.

References calculate_d2xyz, calculations_started, curvatures, libMesh::libmesh_assert(), print_JxW(), and print_xyz().

379  calculate_d2xyz = true; return curvatures;}
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< Real > curvatures
Definition: fe_map.h:858
const std::vector<std::vector<Real> >& libMesh::FEMap::get_d2etadxyz2 ( ) const
inline

Second derivatives of "eta" reference coordinate wrt physical coordinates.

Definition at line 314 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2etadxyz2_map, and libMesh::libmesh_assert().

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

316  calculate_d2xyz = true; return d2etadxyz2_map; }
std::vector< std::vector< Real > > d2etadxyz2_map
Definition: fe_map.h:742
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phideta2_map ( )
inline
Returns
The reference to physical map 2nd derivative

Definition at line 490 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2phideta2_map, and libMesh::libmesh_assert().

492  calculate_d2xyz = true; return d2phideta2_map; }
std::vector< std::vector< Real > > d2phideta2_map
Definition: fe_map.h:791
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidetadzeta_map ( )
inline
Returns
The reference to physical map 2nd derivative

Definition at line 497 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2phidetadzeta_map, and libMesh::libmesh_assert().

499  calculate_d2xyz = true; return d2phidetadzeta_map; }
std::vector< std::vector< Real > > d2phidetadzeta_map
Definition: fe_map.h:796
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidxi2_map ( )
inline
Returns
The reference to physical map 2nd derivative

Definition at line 469 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2phidxi2_map, and libMesh::libmesh_assert().

471  calculate_d2xyz = true; return d2phidxi2_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > d2phidxi2_map
Definition: fe_map.h:776
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidxideta_map ( )
inline
Returns
The reference to physical map 2nd derivative

Definition at line 476 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2phidxideta_map, and libMesh::libmesh_assert().

478  calculate_d2xyz = true; return d2phidxideta_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > d2phidxideta_map
Definition: fe_map.h:781
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidxidzeta_map ( )
inline
Returns
The reference to physical map 2nd derivative

Definition at line 483 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2phidxidzeta_map, and libMesh::libmesh_assert().

485  calculate_d2xyz = true; return d2phidxidzeta_map; }
std::vector< std::vector< Real > > d2phidxidzeta_map
Definition: fe_map.h:786
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidzeta2_map ( )
inline
Returns
The reference to physical map 2nd derivative

Definition at line 504 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2phidzeta2_map, and libMesh::libmesh_assert().

506  calculate_d2xyz = true; return d2phidzeta2_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > d2phidzeta2_map
Definition: fe_map.h:801
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2psideta2 ( )
inline
Returns
The reference to physical map 2nd derivative for the side/edge

Definition at line 433 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2psideta2_map, and libMesh::libmesh_assert().

435  calculate_d2xyz = true; return d2psideta2_map; }
std::vector< std::vector< Real > > d2psideta2_map
Definition: fe_map.h:841
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2psidxi2 ( )
inline
Returns
The reference to physical map 2nd derivative for the side/edge

Definition at line 419 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2psidxi2_map, and libMesh::libmesh_assert().

421  calculate_d2xyz = true; return d2psidxi2_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > d2psidxi2_map
Definition: fe_map.h:827
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2psidxideta ( )
inline
Returns
The reference to physical map 2nd derivative for the side/edge

Definition at line 426 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2psidxideta_map, and libMesh::libmesh_assert().

428  calculate_d2xyz = true; return d2psidxideta_map; }
std::vector< std::vector< Real > > d2psidxideta_map
Definition: fe_map.h:834
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
const std::vector<std::vector<Real> >& libMesh::FEMap::get_d2xidxyz2 ( ) const
inline

Second derivatives of "xi" reference coordinate wrt physical coordinates.

Definition at line 307 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xidxyz2_map, and libMesh::libmesh_assert().

Referenced by libMesh::HCurlFETransformation< OutputShape >::init_map_d2phi(), libMesh::H1FETransformation< OutputShape >::init_map_d2phi(), and libMesh::H1FETransformation< OutputShape >::map_d2phi().

309  calculate_d2xyz = true; return d2xidxyz2_map; }
std::vector< std::vector< Real > > d2xidxyz2_map
Definition: fe_map.h:736
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdeta2 ( ) const
inline
Returns
The second partial derivatives in eta.

Definition at line 191 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xyzdeta2_map, and libMesh::libmesh_assert().

193  calculate_d2xyz = true; return d2xyzdeta2_map; }
libmesh_assert(j)
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdetadzeta ( ) const
inline
Returns
The second partial derivatives in eta-zeta.

Definition at line 225 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xyzdetadzeta_map, and libMesh::libmesh_assert().

227  calculate_d2xyz = true; return d2xyzdetadzeta_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< RealGradient > d2xyzdetadzeta_map
Definition: fe_map.h:665
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdxi2 ( ) const
inline
Returns
The second partial derivatives in xi.

Definition at line 184 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xyzdxi2_map, and libMesh::libmesh_assert().

186  calculate_d2xyz = true; return d2xyzdxi2_map; }
libmesh_assert(j)
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdxideta ( ) const
inline
Returns
The second partial derivatives in xi-eta.

Definition at line 209 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xyzdxideta_map, and libMesh::libmesh_assert().

211  calculate_d2xyz = true; return d2xyzdxideta_map; }
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdxidzeta ( ) const
inline
Returns
The second partial derivatives in xi-zeta.

Definition at line 218 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xyzdxidzeta_map, and libMesh::libmesh_assert().

220  calculate_d2xyz = true; return d2xyzdxidzeta_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
std::vector< RealGradient > d2xyzdxidzeta_map
Definition: fe_map.h:659
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdzeta2 ( ) const
inline
Returns
The second partial derivatives in zeta.

Definition at line 200 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2xyzdzeta2_map, and libMesh::libmesh_assert().

202  calculate_d2xyz = true; return d2xyzdzeta2_map; }
std::vector< RealGradient > d2xyzdzeta2_map
Definition: fe_map.h:671
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculations_started
Definition: fe_map.h:874
const std::vector<std::vector<Real> >& libMesh::FEMap::get_d2zetadxyz2 ( ) const
inline

Second derivatives of "zeta" reference coordinate wrt physical coordinates.

Definition at line 321 of file fe_map.h.

References calculate_d2xyz, calculations_started, d2zetadxyz2_map, and libMesh::libmesh_assert().

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

323  calculate_d2xyz = true; return d2zetadxyz2_map; }
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
std::vector< std::vector< Real > > d2zetadxyz2_map
Definition: fe_map.h:748
bool calculations_started
Definition: fe_map.h:874
const std::vector<Real>& libMesh::FEMap::get_detadx ( ) const
inline
const std::vector<Real>& libMesh::FEMap::get_detady ( ) const
inline
const std::vector<Real>& libMesh::FEMap::get_detadz ( ) const
inline
const std::vector<std::vector<Real> >& libMesh::FEMap::get_dphideta_map ( ) const
inline
Returns
The reference to physical map derivative

Definition at line 349 of file fe_map.h.

References calculate_dxyz, calculations_started, dphideta_map, and libMesh::libmesh_assert().

351  calculate_dxyz = true; return dphideta_map; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > dphideta_map
Definition: fe_map.h:764
std::vector<std::vector<Real> >& libMesh::FEMap::get_dphideta_map ( )
inline
Returns
The reference to physical map derivative

Definition at line 454 of file fe_map.h.

References calculate_dxyz, calculations_started, dphideta_map, and libMesh::libmesh_assert().

456  calculate_dxyz = true; return dphideta_map; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > dphideta_map
Definition: fe_map.h:764
const std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidxi_map ( ) const
inline
Returns
The reference to physical map derivative

Definition at line 342 of file fe_map.h.

References calculate_dxyz, calculations_started, dphidxi_map, and libMesh::libmesh_assert().

344  calculate_dxyz = true; return dphidxi_map; }
std::vector< std::vector< Real > > dphidxi_map
Definition: fe_map.h:759
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidxi_map ( )
inline
Returns
The reference to physical map derivative

Definition at line 447 of file fe_map.h.

References calculate_dxyz, calculations_started, dphidxi_map, and libMesh::libmesh_assert().

449  calculate_dxyz = true; return dphidxi_map; }
std::vector< std::vector< Real > > dphidxi_map
Definition: fe_map.h:759
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
const std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidzeta_map ( ) const
inline
Returns
The reference to physical map derivative

Definition at line 356 of file fe_map.h.

References calculate_dxyz, calculations_started, dphidzeta_map, and libMesh::libmesh_assert().

358  calculate_dxyz = true; return dphidzeta_map; }
std::vector< std::vector< Real > > dphidzeta_map
Definition: fe_map.h:769
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidzeta_map ( )
inline
Returns
The reference to physical map derivative

Definition at line 461 of file fe_map.h.

References calculate_dxyz, calculations_started, dphidzeta_map, and libMesh::libmesh_assert().

463  calculate_dxyz = true; return dphidzeta_map; }
std::vector< std::vector< Real > > dphidzeta_map
Definition: fe_map.h:769
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_dpsideta ( )
inline
Returns
The reference to physical map derivative for the side/edge

Definition at line 412 of file fe_map.h.

References calculate_dxyz, calculations_started, dpsideta_map, and libMesh::libmesh_assert().

414  calculate_dxyz = true; return dpsideta_map; }
std::vector< std::vector< Real > > dpsideta_map
Definition: fe_map.h:820
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector<std::vector<Real> >& libMesh::FEMap::get_dpsidxi ( )
inline
Returns
The reference to physical map derivative for the side/edge

Definition at line 405 of file fe_map.h.

References calculate_dxyz, calculations_started, dpsidxi_map, and libMesh::libmesh_assert().

407  calculate_dxyz = true; return dpsidxi_map; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Real > > dpsidxi_map
Definition: fe_map.h:814
const std::vector<Real>& libMesh::FEMap::get_dxidy ( ) const
inline
const std::vector<Real>& libMesh::FEMap::get_dxidz ( ) const
inline
const std::vector<RealGradient>& libMesh::FEMap::get_dxyzdeta ( ) const
inline
Returns
The element tangents in eta-direction at the quadrature points.

Definition at line 169 of file fe_map.h.

References calculate_dxyz, calculations_started, dxyzdeta_map, and libMesh::libmesh_assert().

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

171  calculate_dxyz = true; return dxyzdeta_map; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
const std::vector<RealGradient>& libMesh::FEMap::get_dxyzdxi ( ) const
inline
Returns
The element tangents in xi-direction at the quadrature points.

Definition at line 161 of file fe_map.h.

References calculate_dxyz, calculations_started, dxyzdxi_map, and libMesh::libmesh_assert().

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

163  calculate_dxyz = true; return dxyzdxi_map; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
bool calculations_started
Definition: fe_map.h:874
const std::vector<RealGradient>& libMesh::FEMap::get_dxyzdzeta ( ) const
inline
Returns
The element tangents in zeta-direction at the quadrature points.

Definition at line 177 of file fe_map.h.

References calculate_dxyz, calculations_started, dxyzdzeta_map, and libMesh::libmesh_assert().

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

179  calculate_dxyz = true; return dxyzdzeta_map; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
bool calculations_started
Definition: fe_map.h:874
const std::vector<Real>& libMesh::FEMap::get_dzetadx ( ) const
inline
const std::vector<Real>& libMesh::FEMap::get_dzetady ( ) const
inline
const std::vector<Real>& libMesh::FEMap::get_dzetadz ( ) const
inline
const std::vector<Real>& libMesh::FEMap::get_jacobian ( ) const
inline
Returns
The element Jacobian for each quadrature point.

Definition at line 145 of file fe_map.h.

References calculate_dxyz, calculations_started, jac, and libMesh::libmesh_assert().

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

147  calculate_dxyz = true; return jac; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< Real > jac
Definition: fe_map.h:863
const std::vector<Real>& libMesh::FEMap::get_JxW ( ) const
inline
Returns
The element Jacobian times the quadrature weight for each quadrature point.

Definition at line 153 of file fe_map.h.

References calculate_dxyz, calculations_started, JxW, and libMesh::libmesh_assert().

155  calculate_dxyz = true; return JxW; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector<Real>& libMesh::FEMap::get_JxW ( )
inline
Returns
Writable reference to the element Jacobian times the quadrature weight for each quadrature point.

Definition at line 515 of file fe_map.h.

References calculate_dxyz, calculations_started, JxW, and libMesh::libmesh_assert().

517  calculate_dxyz = true; return JxW; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< Real > JxW
Definition: fe_map.h:868
const std::vector<Point>& libMesh::FEMap::get_normals ( ) const
inline
Returns
The outward pointing normal vectors for face integration.

Definition at line 370 of file fe_map.h.

References calculate_dxyz, calculations_started, libMesh::libmesh_assert(), and normals.

372  calculate_dxyz = true; return normals; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< Point > normals
Definition: fe_map.h:851
const std::vector<std::vector<Real> >& libMesh::FEMap::get_phi_map ( ) const
inline
Returns
The reference to physical map for the element

Definition at line 335 of file fe_map.h.

References calculate_xyz, calculations_started, libMesh::libmesh_assert(), and phi_map.

337  calculate_xyz = true; return phi_map; }
std::vector< std::vector< Real > > phi_map
Definition: fe_map.h:754
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
bool calculate_xyz
Definition: fe_map.h:879
std::vector<std::vector<Real> >& libMesh::FEMap::get_phi_map ( )
inline
Returns
The reference to physical map for the element

Definition at line 440 of file fe_map.h.

References calculate_xyz, calculations_started, libMesh::libmesh_assert(), and phi_map.

442  calculate_xyz = true; return phi_map; }
std::vector< std::vector< Real > > phi_map
Definition: fe_map.h:754
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
bool calculate_xyz
Definition: fe_map.h:879
const std::vector<std::vector<Real> >& libMesh::FEMap::get_psi ( ) const
inline
Returns
The reference to physical map for the side/edge

Definition at line 329 of file fe_map.h.

References psi_map.

330  { return psi_map; }
std::vector< std::vector< Real > > psi_map
Definition: fe_map.h:808
std::vector<std::vector<Real> >& libMesh::FEMap::get_psi ( )
inline
Returns
The reference to physical map for the side/edge

Definition at line 399 of file fe_map.h.

References psi_map.

400  { return psi_map; }
std::vector< std::vector< Real > > psi_map
Definition: fe_map.h:808
const std::vector<std::vector<Point> >& libMesh::FEMap::get_tangents ( ) const
inline
Returns
The tangent vectors for face integration.

Definition at line 363 of file fe_map.h.

References calculate_dxyz, calculations_started, libMesh::libmesh_assert(), and tangents.

365  calculate_dxyz = true; return tangents; }
bool calculate_dxyz
Definition: fe_map.h:884
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
std::vector< std::vector< Point > > tangents
Definition: fe_map.h:846
const std::vector<Point>& libMesh::FEMap::get_xyz ( ) const
inline
Returns
The xyz spatial locations of the quadrature points on the element.

Definition at line 138 of file fe_map.h.

References calculate_xyz, calculations_started, libMesh::libmesh_assert(), and xyz.

140  calculate_xyz = true; return xyz; }
libmesh_assert(j)
bool calculations_started
Definition: fe_map.h:874
bool calculate_xyz
Definition: fe_map.h:879
std::vector< Point > xyz
Definition: fe_map.h:615
template<unsigned int Dim>
void libMesh::FEMap::init_edge_shape_functions ( const std::vector< Point > &  qp,
const Elem edge 
)

Same as before, but for an edge. This is used for some projection operators.

Definition at line 513 of file fe_boundary.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, d2psidxi2_map, libMesh::Elem::default_order(), determine_calculations(), dpsidxi_map, libMesh::libmesh_assert(), libMesh::FE< Dim, T >::n_shape_functions(), psi_map, libMesh::FE< Dim, T >::shape(), libMesh::FE< Dim, T >::shape_deriv(), libMesh::FE< Dim, T >::shape_second_deriv(), and libMesh::Elem::type().

Referenced by compute_edge_map(), and ~FEMap().

515 {
516  // Start logging the shape function initialization
517  LOG_SCOPE("init_edge_shape_functions()", "FEMap");
518 
519  libmesh_assert(edge);
520 
521  // We're calculating now!
522  this->determine_calculations();
523 
524  // The element type and order to use in
525  // the map
526  const Order mapping_order (edge->default_order());
527  const ElemType mapping_elem_type (edge->type());
528 
529  // The number of quadrature points.
530  const unsigned int n_qp = cast_int<unsigned int>(qp.size());
531 
532  const unsigned int n_mapping_shape_functions =
533  FE<Dim,LAGRANGE>::n_shape_functions (mapping_elem_type,
534  mapping_order);
535 
536  // resize the vectors to hold current data
537  // Psi are the shape functions used for the FE mapping
538  if (calculate_xyz)
539  this->psi_map.resize (n_mapping_shape_functions);
540  if (calculate_dxyz)
541  this->dpsidxi_map.resize (n_mapping_shape_functions);
542  if (calculate_d2xyz)
543  this->d2psidxi2_map.resize (n_mapping_shape_functions);
544 
545  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
546  {
547  // Allocate space to store the values of the shape functions
548  // and their first and second derivatives at the quadrature points.
549  if (calculate_xyz)
550  this->psi_map[i].resize (n_qp);
551  if (calculate_dxyz)
552  this->dpsidxi_map[i].resize (n_qp);
553  if (calculate_d2xyz)
554  this->d2psidxi2_map[i].resize (n_qp);
555 
556  // Compute the value of shape function i, and its first and
557  // second derivatives at quadrature point p
558  // (Lagrange shape functions are used for the mapping)
559  for (unsigned int p=0; p<n_qp; p++)
560  {
561  if (calculate_xyz)
562  this->psi_map[i][p] = FE<1,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
563  if (calculate_dxyz)
564  this->dpsidxi_map[i][p] = FE<1,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
565  if (calculate_d2xyz)
566  this->d2psidxi2_map[i][p] = FE<1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 0, qp[p]);
567  }
568  }
569 }
std::vector< std::vector< Real > > psi_map
Definition: fe_map.h:808
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
bool calculate_dxyz
Definition: fe_map.h:884
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculate_xyz
Definition: fe_map.h:879
std::vector< std::vector< Real > > d2psidxi2_map
Definition: fe_map.h:827
std::vector< std::vector< Real > > dpsidxi_map
Definition: fe_map.h:814
virtual unsigned int n_shape_functions() const libmesh_override
Definition: fe.C:36
void determine_calculations()
Definition: fe_map.h:524
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
template<unsigned int Dim>
void libMesh::FEMap::init_face_shape_functions ( const std::vector< Point > &  qp,
const Elem side 
)

Initalizes the reference to physical element map for a side. This is used for boundary integration.

Definition at line 408 of file fe_boundary.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, d2psideta2_map, d2psidxi2_map, d2psidxideta_map, libMesh::Elem::default_order(), determine_calculations(), dpsideta_map, dpsidxi_map, libMesh::libmesh_assert(), libMesh::FE< Dim, T >::n_shape_functions(), psi_map, libMesh::FE< Dim, T >::shape(), libMesh::FE< Dim, T >::shape_deriv(), libMesh::FE< Dim, T >::shape_second_deriv(), and libMesh::Elem::type().

Referenced by compute_edge_map(), and ~FEMap().

410 {
411  // Start logging the shape function initialization
412  LOG_SCOPE("init_face_shape_functions()", "FEMap");
413 
415 
416  // We're calculating now!
417  this->determine_calculations();
418 
419  // The element type and order to use in
420  // the map
421  const Order mapping_order (side->default_order());
422  const ElemType mapping_elem_type (side->type());
423 
424  // The number of quadrature points.
425  const unsigned int n_qp = cast_int<unsigned int>(qp.size());
426 
427  const unsigned int n_mapping_shape_functions =
428  FE<Dim,LAGRANGE>::n_shape_functions (mapping_elem_type,
429  mapping_order);
430 
431  // resize the vectors to hold current data
432  // Psi are the shape functions used for the FE mapping
433  if (calculate_xyz)
434  this->psi_map.resize (n_mapping_shape_functions);
435 
436  if (Dim > 1)
437  {
438  if (calculate_dxyz)
439  this->dpsidxi_map.resize (n_mapping_shape_functions);
440  if (calculate_d2xyz)
441  this->d2psidxi2_map.resize (n_mapping_shape_functions);
442  }
443 
444  if (Dim == 3)
445  {
446  if (calculate_dxyz)
447  this->dpsideta_map.resize (n_mapping_shape_functions);
448  if (calculate_d2xyz)
449  {
450  this->d2psidxideta_map.resize (n_mapping_shape_functions);
451  this->d2psideta2_map.resize (n_mapping_shape_functions);
452  }
453  }
454 
455  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
456  {
457  // Allocate space to store the values of the shape functions
458  // and their first and second derivatives at the quadrature points.
459  if (calculate_xyz)
460  this->psi_map[i].resize (n_qp);
461  if (Dim > 1)
462  {
463  if (calculate_dxyz)
464  this->dpsidxi_map[i].resize (n_qp);
465  if (calculate_d2xyz)
466  this->d2psidxi2_map[i].resize (n_qp);
467  }
468  if (Dim == 3)
469  {
470  if (calculate_dxyz)
471  this->dpsideta_map[i].resize (n_qp);
472  if (calculate_d2xyz)
473  {
474  this->d2psidxideta_map[i].resize (n_qp);
475  this->d2psideta2_map[i].resize (n_qp);
476  }
477  }
478 
479  // Compute the value of shape function i, and its first and
480  // second derivatives at quadrature point p
481  // (Lagrange shape functions are used for the mapping)
482  for (unsigned int p=0; p<n_qp; p++)
483  {
484  if (calculate_xyz)
485  this->psi_map[i][p] = FE<Dim-1,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
486  if (Dim > 1)
487  {
488  if (calculate_dxyz)
489  this->dpsidxi_map[i][p] = FE<Dim-1,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
490  if (calculate_d2xyz)
491  this->d2psidxi2_map[i][p] = FE<Dim-1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 0, qp[p]);
492  }
493  // libMesh::out << "this->d2psidxi2_map["<<i<<"][p]=" << d2psidxi2_map[i][p] << std::endl;
494 
495  // If we are in 3D, then our sides are 2D faces.
496  // For the second derivatives, we must also compute the cross
497  // derivative d^2() / dxi deta
498  if (Dim == 3)
499  {
500  if (calculate_dxyz)
501  this->dpsideta_map[i][p] = FE<Dim-1,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
502  if (calculate_d2xyz)
503  {
504  this->d2psidxideta_map[i][p] = FE<Dim-1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 1, qp[p]);
505  this->d2psideta2_map[i][p] = FE<Dim-1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 2, qp[p]);
506  }
507  }
508  }
509  }
510 }
std::vector< std::vector< Real > > dpsideta_map
Definition: fe_map.h:820
std::vector< std::vector< Real > > psi_map
Definition: fe_map.h:808
std::vector< std::vector< Real > > d2psideta2_map
Definition: fe_map.h:841
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
bool calculate_dxyz
Definition: fe_map.h:884
unsigned short int side
Definition: xdr_io.C:49
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
std::vector< std::vector< Real > > d2psidxideta_map
Definition: fe_map.h:834
libmesh_assert(j)
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculate_xyz
Definition: fe_map.h:879
std::vector< std::vector< Real > > d2psidxi2_map
Definition: fe_map.h:827
std::vector< std::vector< Real > > dpsidxi_map
Definition: fe_map.h:814
virtual unsigned int n_shape_functions() const libmesh_override
Definition: fe.C:36
void determine_calculations()
Definition: fe_map.h:524
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
template<unsigned int Dim>
template void libMesh::FEMap::init_reference_to_physical_map< 3 > ( const std::vector< Point > &  qp,
const Elem elem 
)

Definition at line 68 of file fe_map.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, d2phideta2_map, d2phidetadzeta_map, d2phidxi2_map, d2phidxideta_map, d2phidxidzeta_map, d2phidzeta2_map, libMesh::Elem::default_order(), determine_calculations(), dphideta_map, dphidxi_map, dphidzeta_map, libMesh::Elem::is_linear(), libMesh::FE< Dim, T >::n_shape_functions(), phi_map, libMesh::FE< Dim, T >::shape(), libMesh::FE< Dim, T >::shape_deriv(), libMesh::FE< Dim, T >::shape_second_deriv(), and libMesh::Elem::type().

Referenced by ~FEMap().

70 {
71  // Start logging the reference->physical map initialization
72  LOG_SCOPE("init_reference_to_physical_map()", "FEMap");
73 
74  // We're calculating now!
75  this->determine_calculations();
76 
77  // The number of quadrature points.
78  const std::size_t n_qp = qp.size();
79 
80  // The element type and order to use in
81  // the map
82  const Order mapping_order (elem->default_order());
83  const ElemType mapping_elem_type (elem->type());
84 
85  // Number of shape functions used to construt the map
86  // (Lagrange shape functions are used for mapping)
87  const unsigned int n_mapping_shape_functions =
88  FE<Dim,LAGRANGE>::n_shape_functions (mapping_elem_type,
89  mapping_order);
90 
91  if (calculate_xyz)
92  this->phi_map.resize (n_mapping_shape_functions);
93  if (Dim > 0)
94  {
95  if (calculate_dxyz)
96  this->dphidxi_map.resize (n_mapping_shape_functions);
97 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
98  if (calculate_d2xyz)
99  this->d2phidxi2_map.resize (n_mapping_shape_functions);
100 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
101  }
102 
103  if (Dim > 1)
104  {
105  if (calculate_dxyz)
106  this->dphideta_map.resize (n_mapping_shape_functions);
107 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
108  if (calculate_d2xyz)
109  {
110  this->d2phidxideta_map.resize (n_mapping_shape_functions);
111  this->d2phideta2_map.resize (n_mapping_shape_functions);
112  }
113 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
114  }
115 
116  if (Dim > 2)
117  {
118  if (calculate_dxyz)
119  this->dphidzeta_map.resize (n_mapping_shape_functions);
120 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
121  if (calculate_d2xyz)
122  {
123  this->d2phidxidzeta_map.resize (n_mapping_shape_functions);
124  this->d2phidetadzeta_map.resize (n_mapping_shape_functions);
125  this->d2phidzeta2_map.resize (n_mapping_shape_functions);
126  }
127 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
128  }
129 
130 
131  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
132  {
133  if (calculate_xyz)
134  this->phi_map[i].resize (n_qp);
135  if (Dim > 0)
136  {
137  if (calculate_dxyz)
138  this->dphidxi_map[i].resize (n_qp);
139 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
140  if (calculate_d2xyz)
141  {
142  this->d2phidxi2_map[i].resize (n_qp);
143  if (Dim > 1)
144  {
145  this->d2phidxideta_map[i].resize (n_qp);
146  this->d2phideta2_map[i].resize (n_qp);
147  }
148  if (Dim > 2)
149  {
150  this->d2phidxidzeta_map[i].resize (n_qp);
151  this->d2phidetadzeta_map[i].resize (n_qp);
152  this->d2phidzeta2_map[i].resize (n_qp);
153  }
154  }
155 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
156 
157  if (Dim > 1 && calculate_dxyz)
158  this->dphideta_map[i].resize (n_qp);
159 
160  if (Dim > 2 && calculate_dxyz)
161  this->dphidzeta_map[i].resize (n_qp);
162  }
163  }
164 
165  // Optimize for the *linear* geometric elements case:
166  bool is_linear = elem->is_linear();
167 
168  switch (Dim)
169  {
170 
171  //------------------------------------------------------------
172  // 0D
173  case 0:
174  {
175  if (calculate_xyz)
176  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
177  for (std::size_t p=0; p<n_qp; p++)
178  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
179 
180  break;
181  }
182 
183  //------------------------------------------------------------
184  // 1D
185  case 1:
186  {
187  // Compute the value of the mapping shape function i at quadrature point p
188  // (Lagrange shape functions are used for mapping)
189  if (is_linear)
190  {
191  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
192  {
193  if (calculate_xyz)
194  this->phi_map[i][0] =
195  FE<Dim,LAGRANGE>::shape(mapping_elem_type,
196  mapping_order,
197  i,
198  qp[0]);
199 
200  if (calculate_dxyz)
201  this->dphidxi_map[i][0] =
202  FE<Dim,LAGRANGE>::shape_deriv(mapping_elem_type,
203  mapping_order,
204  i,
205  0,
206  qp[0]);
207 
208 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
209  if (calculate_d2xyz)
210  this->d2phidxi2_map[i][0] =
211  FE<Dim,LAGRANGE>::shape_second_deriv(mapping_elem_type,
212  mapping_order,
213  i,
214  0,
215  qp[0]);
216 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
217  for (std::size_t p=1; p<n_qp; p++)
218  {
219  if (calculate_xyz)
220  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
221  if (calculate_dxyz)
222  this->dphidxi_map[i][p] = this->dphidxi_map[i][0];
223 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
224  if (calculate_d2xyz)
225  this->d2phidxi2_map[i][p] = this->d2phidxi2_map[i][0];
226 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
227  }
228  }
229  }
230  else
231  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
232  for (std::size_t p=0; p<n_qp; p++)
233  {
234  if (calculate_xyz)
235  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
236  if (calculate_dxyz)
237  this->dphidxi_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
238 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
239  if (calculate_d2xyz)
240  this->d2phidxi2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
241 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
242  }
243 
244  break;
245  }
246  //------------------------------------------------------------
247  // 2D
248  case 2:
249  {
250  // Compute the value of the mapping shape function i at quadrature point p
251  // (Lagrange shape functions are used for mapping)
252  if (is_linear)
253  {
254  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
255  {
256  if (calculate_xyz)
257  this->phi_map[i][0] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[0]);
258  if (calculate_dxyz)
259  {
260  this->dphidxi_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
261  this->dphideta_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
262  }
263 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
264  if (calculate_d2xyz)
265  {
266  this->d2phidxi2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
267  this->d2phidxideta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
268  this->d2phideta2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[0]);
269  }
270 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
271  for (std::size_t p=1; p<n_qp; p++)
272  {
273  if (calculate_xyz)
274  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
275  if (calculate_dxyz)
276  {
277  this->dphidxi_map[i][p] = this->dphidxi_map[i][0];
278  this->dphideta_map[i][p] = this->dphideta_map[i][0];
279  }
280 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
281  if (calculate_d2xyz)
282  {
283  this->d2phidxi2_map[i][p] = this->d2phidxi2_map[i][0];
284  this->d2phidxideta_map[i][p] = this->d2phidxideta_map[i][0];
285  this->d2phideta2_map[i][p] = this->d2phideta2_map[i][0];
286  }
287 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
288  }
289  }
290  }
291  else
292  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
293  for (std::size_t p=0; p<n_qp; p++)
294  {
295  if (calculate_xyz)
296  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
297  if (calculate_dxyz)
298  {
299  this->dphidxi_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
300  this->dphideta_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
301  }
302 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
303  if (calculate_d2xyz)
304  {
305  this->d2phidxi2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
306  this->d2phidxideta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
307  this->d2phideta2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[p]);
308  }
309 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
310  }
311 
312  break;
313  }
314 
315  //------------------------------------------------------------
316  // 3D
317  case 3:
318  {
319  // Compute the value of the mapping shape function i at quadrature point p
320  // (Lagrange shape functions are used for mapping)
321  if (is_linear)
322  {
323  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
324  {
325  if (calculate_xyz)
326  this->phi_map[i][0] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[0]);
327  if (calculate_dxyz)
328  {
329  this->dphidxi_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
330  this->dphideta_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
331  this->dphidzeta_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 2, qp[0]);
332  }
333 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
334  if (calculate_d2xyz)
335  {
336  this->d2phidxi2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
337  this->d2phidxideta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
338  this->d2phideta2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[0]);
339  this->d2phidxidzeta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 3, qp[0]);
340  this->d2phidetadzeta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 4, qp[0]);
341  this->d2phidzeta2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 5, qp[0]);
342  }
343 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
344  for (std::size_t p=1; p<n_qp; p++)
345  {
346  if (calculate_xyz)
347  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
348  if (calculate_dxyz)
349  {
350  this->dphidxi_map[i][p] = this->dphidxi_map[i][0];
351  this->dphideta_map[i][p] = this->dphideta_map[i][0];
352  this->dphidzeta_map[i][p] = this->dphidzeta_map[i][0];
353  }
354 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
355  if (calculate_d2xyz)
356  {
357  this->d2phidxi2_map[i][p] = this->d2phidxi2_map[i][0];
358  this->d2phidxideta_map[i][p] = this->d2phidxideta_map[i][0];
359  this->d2phideta2_map[i][p] = this->d2phideta2_map[i][0];
360  this->d2phidxidzeta_map[i][p] = this->d2phidxidzeta_map[i][0];
361  this->d2phidetadzeta_map[i][p] = this->d2phidetadzeta_map[i][0];
362  this->d2phidzeta2_map[i][p] = this->d2phidzeta2_map[i][0];
363  }
364 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
365  }
366  }
367  }
368  else
369  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
370  for (std::size_t p=0; p<n_qp; p++)
371  {
372  if (calculate_xyz)
373  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
374  if (calculate_dxyz)
375  {
376  this->dphidxi_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
377  this->dphideta_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
378  this->dphidzeta_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 2, qp[p]);
379  }
380 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
381  if (calculate_d2xyz)
382  {
383  this->d2phidxi2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
384  this->d2phidxideta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
385  this->d2phideta2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[p]);
386  this->d2phidxidzeta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 3, qp[p]);
387  this->d2phidetadzeta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 4, qp[p]);
388  this->d2phidzeta2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 5, qp[p]);
389  }
390 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
391  }
392 
393  break;
394  }
395 
396  default:
397  libmesh_error_msg("Invalid Dim = " << Dim);
398  }
399 }
std::vector< std::vector< Real > > dphidzeta_map
Definition: fe_map.h:769
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
std::vector< std::vector< Real > > dphidxi_map
Definition: fe_map.h:759
bool calculate_dxyz
Definition: fe_map.h:884
std::vector< std::vector< Real > > d2phideta2_map
Definition: fe_map.h:791
std::vector< std::vector< Real > > d2phidxidzeta_map
Definition: fe_map.h:786
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
std::vector< std::vector< Real > > phi_map
Definition: fe_map.h:754
std::vector< std::vector< Real > > d2phidetadzeta_map
Definition: fe_map.h:796
bool calculate_d2xyz
Definition: fe_map.h:889
bool calculate_xyz
Definition: fe_map.h:879
std::vector< std::vector< Real > > d2phidxideta_map
Definition: fe_map.h:781
virtual unsigned int n_shape_functions() const libmesh_override
Definition: fe.C:36
void determine_calculations()
Definition: fe_map.h:524
std::vector< std::vector< Real > > d2phidzeta2_map
Definition: fe_map.h:801
std::vector< std::vector< Real > > d2phidxi2_map
Definition: fe_map.h:776
std::vector< std::vector< Real > > dphideta_map
Definition: fe_map.h:764
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
void libMesh::FEMap::print_JxW ( std::ostream &  os) const

Prints the Jacobian times the weight for each quadrature point.

Definition at line 1430 of file fe_map.C.

References JxW.

Referenced by get_curvatures().

1431 {
1432  for (std::size_t i=0; i<JxW.size(); ++i)
1433  os << " [" << i << "]: " << JxW[i] << std::endl;
1434 }
std::vector< Real > JxW
Definition: fe_map.h:868
void libMesh::FEMap::print_xyz ( std::ostream &  os) const

Prints the spatial location of each quadrature point (on the physical element).

Definition at line 1438 of file fe_map.C.

References xyz.

Referenced by get_curvatures().

1439 {
1440  for (std::size_t i=0; i<xyz.size(); ++i)
1441  os << " [" << i << "]: " << xyz[i];
1442 }
std::vector< Point > xyz
Definition: fe_map.h:615
void libMesh::FEMap::resize_quadrature_map_vectors ( const unsigned int  dim,
unsigned int  n_qp 
)
protected

A utility function for use by compute_*_map

Definition at line 1142 of file fe_map.C.

References calculate_d2xyz, calculate_dxyz, calculate_xyz, d2etadxyz2_map, d2xidxyz2_map, d2xyzdeta2_map, d2xyzdetadzeta_map, d2xyzdxi2_map, d2xyzdxideta_map, d2xyzdxidzeta_map, d2xyzdzeta2_map, d2zetadxyz2_map, detadx_map, detady_map, detadz_map, determine_calculations(), dxidx_map, dxidy_map, dxidz_map, dxyzdeta_map, dxyzdxi_map, dxyzdzeta_map, dzetadx_map, dzetady_map, dzetadz_map, jac, JxW, and xyz.

Referenced by compute_affine_map(), compute_map(), compute_null_map(), and determine_calculations().

1143 {
1144  // We're calculating now!
1145  this->determine_calculations();
1146 
1147  // Resize the vectors to hold data at the quadrature points
1148  if (calculate_xyz)
1149  xyz.resize(n_qp);
1150  if (calculate_dxyz)
1151  {
1152  dxyzdxi_map.resize(n_qp);
1153  dxidx_map.resize(n_qp);
1154  dxidy_map.resize(n_qp); // 1D element may live in 2D ...
1155  dxidz_map.resize(n_qp); // ... or 3D
1156  }
1157 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1158  if (calculate_d2xyz)
1159  {
1160  d2xyzdxi2_map.resize(n_qp);
1161 
1162  // Inverse map second derivatives
1163  d2xidxyz2_map.resize(n_qp);
1164  for (std::size_t i=0; i<d2xidxyz2_map.size(); ++i)
1165  d2xidxyz2_map[i].assign(6, 0.);
1166  }
1167 #endif
1168  if (dim > 1)
1169  {
1170  if (calculate_dxyz)
1171  {
1172  dxyzdeta_map.resize(n_qp);
1173  detadx_map.resize(n_qp);
1174  detady_map.resize(n_qp);
1175  detadz_map.resize(n_qp);
1176  }
1177 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1178  if (calculate_d2xyz)
1179  {
1180  d2xyzdxideta_map.resize(n_qp);
1181  d2xyzdeta2_map.resize(n_qp);
1182 
1183  // Inverse map second derivatives
1184  d2etadxyz2_map.resize(n_qp);
1185  for (std::size_t i=0; i<d2etadxyz2_map.size(); ++i)
1186  d2etadxyz2_map[i].assign(6, 0.);
1187  }
1188 #endif
1189  if (dim > 2)
1190  {
1191  if (calculate_dxyz)
1192  {
1193  dxyzdzeta_map.resize (n_qp);
1194  dzetadx_map.resize (n_qp);
1195  dzetady_map.resize (n_qp);
1196  dzetadz_map.resize (n_qp);
1197  }
1198 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1199  if (calculate_d2xyz)
1200  {
1201  d2xyzdxidzeta_map.resize(n_qp);
1202  d2xyzdetadzeta_map.resize(n_qp);
1203  d2xyzdzeta2_map.resize(n_qp);
1204 
1205  // Inverse map second derivatives
1206  d2zetadxyz2_map.resize(n_qp);
1207  for (std::size_t i=0; i<d2zetadxyz2_map.size(); ++i)
1208  d2zetadxyz2_map[i].assign(6, 0.);
1209  }
1210 #endif
1211  }
1212  }
1213 
1214  if (calculate_dxyz)
1215  {
1216  jac.resize(n_qp);
1217  JxW.resize(n_qp);
1218  }
1219 }
std::vector< std::vector< Real > > d2etadxyz2_map
Definition: fe_map.h:742
bool calculate_dxyz
Definition: fe_map.h:884
std::vector< RealGradient > d2xyzdzeta2_map
Definition: fe_map.h:671
std::vector< Real > dzetady_map
Definition: fe_map.h:723
std::vector< std::vector< Real > > d2xidxyz2_map
Definition: fe_map.h:736
std::vector< Real > dxidz_map
Definition: fe_map.h:691
std::vector< RealGradient > d2xyzdxideta_map
Definition: fe_map.h:645
std::vector< RealGradient > dxyzdzeta_map
Definition: fe_map.h:633
std::vector< Real > dzetadx_map
Definition: fe_map.h:717
std::vector< RealGradient > dxyzdxi_map
Definition: fe_map.h:621
std::vector< RealGradient > d2xyzdeta2_map
Definition: fe_map.h:651
std::vector< RealGradient > d2xyzdxi2_map
Definition: fe_map.h:639
std::vector< Real > dzetadz_map
Definition: fe_map.h:729
bool calculate_d2xyz
Definition: fe_map.h:889
std::vector< std::vector< Real > > d2zetadxyz2_map
Definition: fe_map.h:748
std::vector< RealGradient > d2xyzdetadzeta_map
Definition: fe_map.h:665
std::vector< Real > dxidx_map
Definition: fe_map.h:679
bool calculate_xyz
Definition: fe_map.h:879
std::vector< Real > dxidy_map
Definition: fe_map.h:685
std::vector< Real > JxW
Definition: fe_map.h:868
std::vector< Point > xyz
Definition: fe_map.h:615
std::vector< Real > detady_map
Definition: fe_map.h:704
std::vector< RealGradient > dxyzdeta_map
Definition: fe_map.h:627
void determine_calculations()
Definition: fe_map.h:524
std::vector< Real > jac
Definition: fe_map.h:863
std::vector< Real > detadz_map
Definition: fe_map.h:710
std::vector< Real > detadx_map
Definition: fe_map.h:698
std::vector< RealGradient > d2xyzdxidzeta_map
Definition: fe_map.h:659

Friends And Related Function Documentation

template<unsigned int Dim, FEFamily T>
friend class FE
friend

FE classes should be able to reset calculations_started in a few special cases.

Definition at line 896 of file fe_map.h.

Member Data Documentation

std::vector<Real> libMesh::FEMap::curvatures
protected

The mean curvature (= one half the sum of the principal curvatures) on the boundary at the quadrature points. The mean curvature is a scalar value.

Definition at line 858 of file fe_map.h.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and get_curvatures().

std::vector<std::vector<Real> > libMesh::FEMap::d2etadxyz2_map
protected

Second derivatives of "eta" reference coordinate wrt physical coordinates. At each qp: (eta_{xx}, eta_{xy}, eta_{xz}, eta_{yy}, eta_{yz}, eta_{zz})

Definition at line 742 of file fe_map.h.

Referenced by compute_inverse_map_second_derivs(), compute_single_point_map(), get_d2etadxyz2(), and resize_quadrature_map_vectors().

std::vector<std::vector<Real> > libMesh::FEMap::d2phideta2_map
protected

Map for the second derivative, d^2(phi)/d(eta)^2.

Definition at line 791 of file fe_map.h.

Referenced by compute_single_point_map(), get_d2phideta2_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::d2phidetadzeta_map
protected

Map for the second derivative, d^2(phi)/d(eta)d(zeta).

Definition at line 796 of file fe_map.h.

Referenced by compute_single_point_map(), get_d2phidetadzeta_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::d2phidxi2_map
protected

Map for the second derivative, d^2(phi)/d(xi)^2.

Definition at line 776 of file fe_map.h.

Referenced by compute_single_point_map(), get_d2phidxi2_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::d2phidxideta_map
protected

Map for the second derivative, d^2(phi)/d(xi)d(eta).

Definition at line 781 of file fe_map.h.

Referenced by compute_single_point_map(), get_d2phidxideta_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::d2phidxidzeta_map
protected

Map for the second derivative, d^2(phi)/d(xi)d(zeta).

Definition at line 786 of file fe_map.h.

Referenced by compute_single_point_map(), get_d2phidxidzeta_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::d2phidzeta2_map
protected

Map for the second derivative, d^2(phi)/d(zeta)^2.

Definition at line 801 of file fe_map.h.

Referenced by compute_single_point_map(), get_d2phidzeta2_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::d2psideta2_map
protected

Map for the second derivatives (in eta) of the side shape functions. Useful for computing the curvature at the quadrature points.

Definition at line 841 of file fe_map.h.

Referenced by libMesh::FEXYZMap::compute_face_map(), compute_face_map(), get_d2psideta2(), and init_face_shape_functions().

std::vector<std::vector<Real> > libMesh::FEMap::d2psidxi2_map
protected

Map for the second derivatives (in xi) of the side shape functions. Useful for computing the curvature at the quadrature points.

Definition at line 827 of file fe_map.h.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), get_d2psidxi2(), init_edge_shape_functions(), and init_face_shape_functions().

std::vector<std::vector<Real> > libMesh::FEMap::d2psidxideta_map
protected

Map for the second (cross) derivatives in xi, eta of the side shape functions. Useful for computing the curvature at the quadrature points.

Definition at line 834 of file fe_map.h.

Referenced by libMesh::FEXYZMap::compute_face_map(), compute_face_map(), get_d2psidxideta(), and init_face_shape_functions().

std::vector<std::vector<Real> > libMesh::FEMap::d2xidxyz2_map
protected

Second derivatives of "xi" reference coordinate wrt physical coordinates. At each qp: (xi_{xx}, xi_{xy}, xi_{xz}, xi_{yy}, xi_{yz}, xi_{zz})

Definition at line 736 of file fe_map.h.

Referenced by compute_inverse_map_second_derivs(), compute_single_point_map(), get_d2xidxyz2(), and resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::d2xyzdeta2_map
protected
std::vector<RealGradient> libMesh::FEMap::d2xyzdetadzeta_map
protected

Vector of mixed second partial derivatives in eta-zeta: d^2(x)/d(eta)d(zeta) d^2(y)/d(eta)d(zeta) d^2(z)/d(eta)d(zeta)

Definition at line 665 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_d2xyzdetadzeta(), and resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::d2xyzdxi2_map
protected

Vector of second partial derivatives in xi: d^2(x)/d(xi)^2, d^2(y)/d(xi)^2, d^2(z)/d(xi)^2

Definition at line 639 of file fe_map.h.

Referenced by compute_affine_map(), compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_d2xyzdxi2(), and resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::d2xyzdxideta_map
protected

Vector of mixed second partial derivatives in xi-eta: d^2(x)/d(xi)d(eta) d^2(y)/d(xi)d(eta) d^2(z)/d(xi)d(eta)

Definition at line 645 of file fe_map.h.

Referenced by compute_affine_map(), compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_d2xyzdxideta(), and resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::d2xyzdxidzeta_map
protected

Vector of second partial derivatives in xi-zeta: d^2(x)/d(xi)d(zeta), d^2(y)/d(xi)d(zeta), d^2(z)/d(xi)d(zeta)

Definition at line 659 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_d2xyzdxidzeta(), and resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::d2xyzdzeta2_map
protected

Vector of second partial derivatives in zeta: d^2(x)/d(zeta)^2

Definition at line 671 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_d2xyzdzeta2(), and resize_quadrature_map_vectors().

std::vector<std::vector<Real> > libMesh::FEMap::d2zetadxyz2_map
protected

Second derivatives of "zeta" reference coordinate wrt physical coordinates. At each qp: (zeta_{xx}, zeta_{xy}, zeta_{xz}, zeta_{yy}, zeta_{yz}, zeta_{zz})

Definition at line 748 of file fe_map.h.

Referenced by compute_inverse_map_second_derivs(), compute_single_point_map(), get_d2zetadxyz2(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::detadx_map
protected

Map for partial derivatives: d(eta)/d(x). Needed for the Jacobian.

Definition at line 698 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_detadx(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::detady_map
protected

Map for partial derivatives: d(eta)/d(y). Needed for the Jacobian.

Definition at line 704 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_detady(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::detadz_map
protected

Map for partial derivatives: d(eta)/d(z). Needed for the Jacobian.

Definition at line 710 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_detadz(), and resize_quadrature_map_vectors().

std::vector<std::vector<Real> > libMesh::FEMap::dphideta_map
protected

Map for the derivative, d(phi)/d(eta).

Definition at line 764 of file fe_map.h.

Referenced by compute_single_point_map(), get_dphideta_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::dphidxi_map
protected

Map for the derivative, d(phi)/d(xi).

Definition at line 759 of file fe_map.h.

Referenced by compute_single_point_map(), get_dphidxi_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::dphidzeta_map
protected

Map for the derivative, d(phi)/d(zeta).

Definition at line 769 of file fe_map.h.

Referenced by compute_single_point_map(), get_dphidzeta_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::dpsideta_map
protected

Map for the derivative of the side function, d(psi)/d(eta).

Definition at line 820 of file fe_map.h.

Referenced by libMesh::FEXYZMap::compute_face_map(), compute_face_map(), get_dpsideta(), and init_face_shape_functions().

std::vector<std::vector<Real> > libMesh::FEMap::dpsidxi_map
protected

Map for the derivative of the side functions, d(psi)/d(xi).

Definition at line 814 of file fe_map.h.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), get_dpsidxi(), init_edge_shape_functions(), and init_face_shape_functions().

std::vector<Real> libMesh::FEMap::dxidx_map
protected

Map for partial derivatives: d(xi)/d(x). Needed for the Jacobian.

Definition at line 679 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_dxidx(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::dxidy_map
protected

Map for partial derivatives: d(xi)/d(y). Needed for the Jacobian.

Definition at line 685 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_dxidy(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::dxidz_map
protected

Map for partial derivatives: d(xi)/d(z). Needed for the Jacobian.

Definition at line 691 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_dxidz(), and resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::dxyzdeta_map
protected
std::vector<RealGradient> libMesh::FEMap::dxyzdxi_map
protected
std::vector<RealGradient> libMesh::FEMap::dxyzdzeta_map
protected

Vector of parital derivatives: d(x)/d(zeta), d(y)/d(zeta), d(z)/d(zeta)

Definition at line 633 of file fe_map.h.

Referenced by compute_affine_map(), compute_null_map(), compute_single_point_map(), dxdzeta_map(), dydzeta_map(), dzdzeta_map(), get_dxyzdzeta(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::dzetadx_map
protected

Map for partial derivatives: d(zeta)/d(x). Needed for the Jacobian.

Definition at line 717 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_dzetadx(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::dzetady_map
protected

Map for partial derivatives: d(zeta)/d(y). Needed for the Jacobian.

Definition at line 723 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_dzetady(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::dzetadz_map
protected

Map for partial derivatives: d(zeta)/d(z). Needed for the Jacobian.

Definition at line 729 of file fe_map.h.

Referenced by compute_affine_map(), compute_inverse_map_second_derivs(), compute_null_map(), compute_single_point_map(), get_dzetadz(), and resize_quadrature_map_vectors().

std::vector<const Node *> libMesh::FEMap::elem_nodes
private

Work vector for compute_affine_map()

Definition at line 908 of file fe_map.h.

Referenced by compute_affine_map(), compute_map(), and ~FEMap().

std::vector<Real> libMesh::FEMap::jac
protected

Jacobian values at quadrature points

Definition at line 863 of file fe_map.h.

Referenced by compute_affine_map(), compute_null_map(), compute_single_point_map(), get_jacobian(), and resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::JxW
protected
std::vector<Point> libMesh::FEMap::normals
protected

Normal vectors on boundary at quadrature points

Definition at line 851 of file fe_map.h.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and get_normals().

std::vector<std::vector<Real> > libMesh::FEMap::phi_map
protected

Map for the shape function phi.

Definition at line 754 of file fe_map.h.

Referenced by compute_affine_map(), compute_single_point_map(), get_phi_map(), and init_reference_to_physical_map().

std::vector<std::vector<Real> > libMesh::FEMap::psi_map
protected

Map for the side shape functions, psi.

Definition at line 808 of file fe_map.h.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), get_psi(), init_edge_shape_functions(), and init_face_shape_functions().

std::vector<std::vector<Point> > libMesh::FEMap::tangents
protected

Tangent vectors on boundary at quadrature points.

Definition at line 846 of file fe_map.h.

Referenced by compute_edge_map(), libMesh::FEXYZMap::compute_face_map(), compute_face_map(), and get_tangents().

std::vector<Point> libMesh::FEMap::xyz
protected

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