libMesh::FEXYZMap Class Reference

#include <fe_xyz_map.h>

Inheritance diagram for libMesh::FEXYZMap:

Public Member Functions

 FEXYZMap ()
 
virtual ~FEXYZMap ()
 
virtual void compute_face_map (int dim, const std::vector< Real > &qw, const Elem *side)
 
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)
 
virtual void compute_affine_map (const unsigned int dim, const std::vector< Real > &qw, const Elem *elem)
 
virtual void compute_map (const unsigned int dim, const std::vector< Real > &qw, const Elem *elem)
 
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
 
std::vector< Real > & get_JxW ()
 
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_psi () const
 
std::vector< std::vector< Real > > & get_psi ()
 
const std::vector< std::vector
< Real > > & 
get_phi_map () const
 
std::vector< std::vector< Real > > & get_phi_map ()
 
const std::vector< std::vector
< Real > > & 
get_dphidxi_map () const
 
std::vector< std::vector< Real > > & get_dphidxi_map ()
 
const std::vector< std::vector
< Real > > & 
get_dphideta_map () const
 
std::vector< std::vector< Real > > & get_dphideta_map ()
 
const std::vector< std::vector
< Real > > & 
get_dphidzeta_map () const
 
std::vector< std::vector< Real > > & get_dphidzeta_map ()
 
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_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_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 ()
 

Static Public Member Functions

static AutoPtr< FEMapbuild (FEType fe_type)
 

Protected Member Functions

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 > > 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
 

Detailed Description

Definition at line 30 of file fe_xyz_map.h.

Constructor & Destructor Documentation

libMesh::FEXYZMap::FEXYZMap ( )
inline

Definition at line 34 of file fe_xyz_map.h.

35  : FEMap(){}
virtual libMesh::FEXYZMap::~FEXYZMap ( )
inlinevirtual

Definition at line 37 of file fe_xyz_map.h.

37 {}

Member Function Documentation

AutoPtr< FEMap > libMesh::FEMap::build ( FEType  fe_type)
staticinherited

Definition at line 41 of file fe_map.C.

References libMesh::FEType::family, and libMeshEnums::XYZ.

42 {
43  switch( fe_type.family )
44  {
45  case XYZ:
46  {
47  AutoPtr<FEMap> ap( new FEXYZMap );
48  return ap;
49  }
50  default:
51  {
52  AutoPtr<FEMap> ap( new FEMap );
53  return ap;
54  }
55  }
56 
57  //Shouldn't ever get here
58  libmesh_error();
59  return AutoPtr<FEMap>();
60 }
void libMesh::FEMap::compute_affine_map ( const unsigned int  dim,
const std::vector< Real > &  qw,
const Elem elem 
)
virtualinherited

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 716 of file fe_map.C.

References libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::d2xyzdeta2_map, libMesh::FEMap::d2xyzdetadzeta_map, libMesh::FEMap::d2xyzdxi2_map, libMesh::FEMap::d2xyzdxideta_map, libMesh::FEMap::d2xyzdxidzeta_map, libMesh::FEMap::d2xyzdzeta2_map, libMesh::FEMap::detadx_map, libMesh::FEMap::detady_map, libMesh::FEMap::detadz_map, libMesh::FEMap::dxidx_map, libMesh::FEMap::dxidy_map, libMesh::FEMap::dxidz_map, libMesh::FEMap::dxyzdeta_map, libMesh::FEMap::dxyzdxi_map, libMesh::FEMap::dxyzdzeta_map, libMesh::FEMap::dzetadx_map, libMesh::FEMap::dzetady_map, libMesh::FEMap::dzetadz_map, libMesh::FEMap::jac, libMesh::FEMap::JxW, libMesh::libmesh_assert(), libMesh::FEMap::phi_map, libMesh::Elem::point(), libMesh::FEMap::resize_quadrature_map_vectors(), libMesh::START_LOG(), libMesh::STOP_LOG(), and libMesh::FEMap::xyz.

Referenced by libMesh::FEMap::compute_map().

719 {
720  // Start logging the map computation.
721  START_LOG("compute_affine_map()", "FEMap");
722 
723  libmesh_assert(elem);
724 
725  const unsigned int n_qp = libmesh_cast_int<unsigned int>(qw.size());
726 
727  // Resize the vectors to hold data at the quadrature points
728  this->resize_quadrature_map_vectors(dim, n_qp);
729 
730  // Compute map at quadrature point 0
731  this->compute_single_point_map(dim, qw, elem, 0);
732 
733  // Compute xyz at all other quadrature points
734  for (unsigned int p=1; p<n_qp; p++)
735  {
736  xyz[p].zero();
737  for (unsigned int i=0; i<phi_map.size(); i++) // sum over the nodes
738  xyz[p].add_scaled (elem->point(i), phi_map[i][p] );
739  }
740 
741  // Copy other map data from quadrature point 0
742  for (unsigned int p=1; p<n_qp; p++) // for each extra quadrature point
743  {
744  dxyzdxi_map[p] = dxyzdxi_map[0];
745  dxidx_map[p] = dxidx_map[0];
746  dxidy_map[p] = dxidy_map[0];
747  dxidz_map[p] = dxidz_map[0];
748 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
749  // The map should be affine, so second derivatives are zero
750  d2xyzdxi2_map[p] = 0.;
751 #endif
752  if (dim > 1)
753  {
754  dxyzdeta_map[p] = dxyzdeta_map[0];
755  detadx_map[p] = detadx_map[0];
756  detady_map[p] = detady_map[0];
757  detadz_map[p] = detadz_map[0];
758 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
759  d2xyzdxideta_map[p] = 0.;
760  d2xyzdeta2_map[p] = 0.;
761 #endif
762  if (dim > 2)
763  {
765  dzetadx_map[p] = dzetadx_map[0];
766  dzetady_map[p] = dzetady_map[0];
767  dzetadz_map[p] = dzetadz_map[0];
768 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
769  d2xyzdxidzeta_map[p] = 0.;
770  d2xyzdetadzeta_map[p] = 0.;
771  d2xyzdzeta2_map[p] = 0.;
772 #endif
773  }
774  }
775  jac[p] = jac[0];
776  JxW[p] = JxW[0] / qw[0] * qw[p];
777  }
778 
779  STOP_LOG("compute_affine_map()", "FEMap");
780 }
void libMesh::FEMap::compute_edge_map ( int  dim,
const std::vector< Real > &  qw,
const Elem side 
)
inherited

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

Definition at line 807 of file fe_boundary.C.

References libMesh::FEMap::compute_face_map(), libMesh::FEMap::curvatures, libMesh::FEMap::d2psidxi2_map, libMesh::FEMap::d2xyzdeta2_map, libMesh::FEMap::d2xyzdxi2_map, libMesh::FEMap::d2xyzdxideta_map, libMesh::FEMap::dpsidxi_map, libMesh::FEMap::dxdxi_map(), libMesh::FEMap::dxyzdeta_map, libMesh::FEMap::dxyzdxi_map, libMesh::FEMap::dydxi_map(), libMesh::FEMap::dzdxi_map(), libMesh::FEMap::JxW, libMesh::libmesh_assert(), libMesh::libmesh_assert_greater(), libMesh::FEMap::normals, libMesh::Elem::point(), libMesh::FEMap::psi_map, libMesh::Real, libMesh::START_LOG(), libMesh::STOP_LOG(), libMesh::FEMap::tangents, and libMesh::FEMap::xyz.

810 {
811  libmesh_assert(edge);
812 
813  if (dim == 2)
814  {
815  // A 2D finite element living in either 2D or 3D space.
816  // The edges here are the sides of the element, so the
817  // (misnamed) compute_face_map function does what we want
818  this->compute_face_map(dim, qw, edge);
819  return;
820  }
821 
822  libmesh_assert_equal_to (dim, 3); // 1D is unnecessary and currently unsupported
823 
824  START_LOG("compute_edge_map()", "FEMap");
825 
826  // The number of quadrature points.
827  const unsigned int n_qp = libmesh_cast_int<unsigned int>(qw.size());
828 
829  // Resize the vectors to hold data at the quadrature points
830  this->xyz.resize(n_qp);
831  this->dxyzdxi_map.resize(n_qp);
832  this->dxyzdeta_map.resize(n_qp);
833  this->d2xyzdxi2_map.resize(n_qp);
834  this->d2xyzdxideta_map.resize(n_qp);
835  this->d2xyzdeta2_map.resize(n_qp);
836  this->tangents.resize(n_qp);
837  this->normals.resize(n_qp);
838  this->curvatures.resize(n_qp);
839 
840  this->JxW.resize(n_qp);
841 
842  // Clear the entities that will be summed
843  for (unsigned int p=0; p<n_qp; p++)
844  {
845  this->tangents[p].resize(1);
846  this->xyz[p].zero();
847  this->dxyzdxi_map[p].zero();
848  this->dxyzdeta_map[p].zero();
849  this->d2xyzdxi2_map[p].zero();
850  this->d2xyzdxideta_map[p].zero();
851  this->d2xyzdeta2_map[p].zero();
852  }
853 
854  // compute x, dxdxi at the quadrature points
855  for (unsigned int i=0; i<this->psi_map.size(); i++) // sum over the nodes
856  {
857  const Point& edge_point = edge->point(i);
858 
859  for (unsigned int p=0; p<n_qp; p++) // for each quadrature point...
860  {
861  this->xyz[p].add_scaled (edge_point, this->psi_map[i][p]);
862  this->dxyzdxi_map[p].add_scaled (edge_point, this->dpsidxi_map[i][p]);
863  this->d2xyzdxi2_map[p].add_scaled (edge_point, this->d2psidxi2_map[i][p]);
864  }
865  }
866 
867  // Compute the tangents at the quadrature point
868  // FIXME: normals (plural!) and curvatures are uncalculated
869  for (unsigned int p=0; p<n_qp; p++)
870  {
871  const Point n = this->dxyzdxi_map[p].cross(this->dxyzdeta_map[p]);
872  this->tangents[p][0] = this->dxyzdxi_map[p].unit();
873 
874  // compute the jacobian at the quadrature points
875  const Real the_jac = std::sqrt(this->dxdxi_map(p)*this->dxdxi_map(p) +
876  this->dydxi_map(p)*this->dydxi_map(p) +
877  this->dzdxi_map(p)*this->dzdxi_map(p));
878 
879  libmesh_assert_greater (the_jac, 0.);
880 
881  this->JxW[p] = the_jac*qw[p];
882  }
883 
884  STOP_LOG("compute_edge_map()", "FEMap");
885 }
void libMesh::FEXYZMap::compute_face_map ( int  dim,
const std::vector< Real > &  qw,
const Elem side 
)
virtual

Special implementation for XYZ finite elements

Reimplemented from libMesh::FEMap.

Definition at line 22 of file fe_xyz_map.C.

References libMesh::TypeVector< T >::cross(), libMesh::FEMap::curvatures, libMesh::FEMap::d2psideta2_map, libMesh::FEMap::d2psidxi2_map, libMesh::FEMap::d2psidxideta_map, libMesh::FEMap::d2xyzdeta2_map, libMesh::FEMap::d2xyzdxi2_map, libMesh::FEMap::d2xyzdxideta_map, libMesh::FEMap::dpsideta_map, libMesh::FEMap::dpsidxi_map, libMesh::FEMap::dxdeta_map(), libMesh::FEMap::dxdxi_map(), libMesh::FEMap::dxyzdeta_map, libMesh::FEMap::dxyzdxi_map, libMesh::FEMap::dydeta_map(), libMesh::FEMap::dydxi_map(), libMesh::FEMap::dzdeta_map(), libMesh::FEMap::dzdxi_map(), libMesh::FEMap::JxW, libMesh::libmesh_assert(), libMesh::libmesh_assert_greater(), libMesh::FEMap::normals, libMesh::Elem::point(), libMesh::FEMap::psi_map, libMesh::Real, libMesh::START_LOG(), libMesh::STOP_LOG(), libMesh::FEMap::tangents, libMesh::TypeVector< T >::unit(), and libMesh::FEMap::xyz.

23 {
25 
26  START_LOG("compute_face_map()", "FEXYZMap");
27 
28  // The number of quadrature points.
29  const unsigned int n_qp = libmesh_cast_int<unsigned int>(qw.size());
30 
31  switch(dim)
32  {
33  case 2:
34  {
35 
36  // Resize the vectors to hold data at the quadrature points
37  {
38  this->xyz.resize(n_qp);
39  this->dxyzdxi_map.resize(n_qp);
40  this->d2xyzdxi2_map.resize(n_qp);
41  this->tangents.resize(n_qp);
42  this->normals.resize(n_qp);
43  this->curvatures.resize(n_qp);
44 
45  this->JxW.resize(n_qp);
46  }
47 
48  // Clear the entities that will be summed
49  // Compute the tangent & normal at the quadrature point
50  for (unsigned int p=0; p<n_qp; p++)
51  {
52  this->tangents[p].resize(LIBMESH_DIM-1); // 1 Tangent in 2D, 2 in 3D
53  this->xyz[p].zero();
54  this->dxyzdxi_map[p].zero();
55  this->d2xyzdxi2_map[p].zero();
56  }
57 
58  // compute x, dxdxi at the quadrature points
59  for (unsigned int i=0; i<this->psi_map.size(); i++) // sum over the nodes
60  {
61  const Point& side_point = side->point(i);
62 
63  for (unsigned int p=0; p<n_qp; p++) // for each quadrature point...
64  {
65  this->xyz[p].add_scaled (side_point, this->psi_map[i][p]);
66  this->dxyzdxi_map[p].add_scaled (side_point, this->dpsidxi_map[i][p]);
67  this->d2xyzdxi2_map[p].add_scaled(side_point, this->d2psidxi2_map[i][p]);
68  }
69  }
70 
71  // Compute the tangent & normal at the quadrature point
72  for (unsigned int p=0; p<n_qp; p++)
73  {
74  const Point n(this->dxyzdxi_map[p](1), -this->dxyzdxi_map[p](0), 0.);
75 
76  this->normals[p] = n.unit();
77  this->tangents[p][0] = this->dxyzdxi_map[p].unit();
78 #if LIBMESH_DIM == 3 // Only good in 3D space
79  this->tangents[p][1] = this->dxyzdxi_map[p].cross(n).unit();
80 #endif
81  // The curvature is computed via the familiar Frenet formula:
82  // curvature = [d^2(x) / d (xi)^2] dot [normal]
83  // For a reference, see:
84  // F.S. Merritt, Mathematics Manual, 1962, McGraw-Hill, p. 310
85  //
86  // Note: The sign convention here is different from the
87  // 3D case. Concave-upward curves (smiles) have a positive
88  // curvature. Concave-downward curves (frowns) have a
89  // negative curvature. Be sure to take that into account!
90  const Real numerator = this->d2xyzdxi2_map[p] * this->normals[p];
91  const Real denominator = this->dxyzdxi_map[p].size_sq();
92  libmesh_assert_not_equal_to (denominator, 0);
93  this->curvatures[p] = numerator / denominator;
94  }
95 
96  // compute the jacobian at the quadrature points
97  for (unsigned int p=0; p<n_qp; p++)
98  {
99  const Real the_jac = std::sqrt(this->dxdxi_map(p)*this->dxdxi_map(p) +
100  this->dydxi_map(p)*this->dydxi_map(p));
101 
102  libmesh_assert_greater (the_jac, 0.);
103 
104  this->JxW[p] = the_jac*qw[p];
105  }
106 
107  break;
108  }
109 
110  case 3:
111  {
112  // Resize the vectors to hold data at the quadrature points
113  {
114  this->xyz.resize(n_qp);
115  this->dxyzdxi_map.resize(n_qp);
116  this->dxyzdeta_map.resize(n_qp);
117  this->d2xyzdxi2_map.resize(n_qp);
118  this->d2xyzdxideta_map.resize(n_qp);
119  this->d2xyzdeta2_map.resize(n_qp);
120  this->tangents.resize(n_qp);
121  this->normals.resize(n_qp);
122  this->curvatures.resize(n_qp);
123 
124  this->JxW.resize(n_qp);
125  }
126 
127  // Clear the entities that will be summed
128  for (unsigned int p=0; p<n_qp; p++)
129  {
130  this->tangents[p].resize(LIBMESH_DIM-1); // 1 Tangent in 2D, 2 in 3D
131  this->xyz[p].zero();
132  this->dxyzdxi_map[p].zero();
133  this->dxyzdeta_map[p].zero();
134  this->d2xyzdxi2_map[p].zero();
135  this->d2xyzdxideta_map[p].zero();
136  this->d2xyzdeta2_map[p].zero();
137  }
138 
139  // compute x, dxdxi at the quadrature points
140  for (unsigned int i=0; i<this->psi_map.size(); i++) // sum over the nodes
141  {
142  const Point& side_point = side->point(i);
143 
144  for (unsigned int p=0; p<n_qp; p++) // for each quadrature point...
145  {
146  this->xyz[p].add_scaled (side_point, this->psi_map[i][p]);
147  this->dxyzdxi_map[p].add_scaled (side_point, this->dpsidxi_map[i][p]);
148  this->dxyzdeta_map[p].add_scaled (side_point, this->dpsideta_map[i][p]);
149  this->d2xyzdxi2_map[p].add_scaled (side_point, this->d2psidxi2_map[i][p]);
150  this->d2xyzdxideta_map[p].add_scaled(side_point, this->d2psidxideta_map[i][p]);
151  this->d2xyzdeta2_map[p].add_scaled (side_point, this->d2psideta2_map[i][p]);
152  }
153  }
154 
155  // Compute the tangents, normal, and curvature at the quadrature point
156  for (unsigned int p=0; p<n_qp; p++)
157  {
158  const Point n = this->dxyzdxi_map[p].cross(this->dxyzdeta_map[p]);
159  this->normals[p] = n.unit();
160  this->tangents[p][0] = this->dxyzdxi_map[p].unit();
161  this->tangents[p][1] = n.cross(this->dxyzdxi_map[p]).unit();
162 
163  // Compute curvature using the typical nomenclature
164  // of the first and second fundamental forms.
165  // For reference, see:
166  // 1) http://mathworld.wolfram.com/MeanCurvature.html
167  // (note -- they are using inward normal)
168  // 2) F.S. Merritt, Mathematics Manual, 1962, McGraw-Hill
169  const Real L = -this->d2xyzdxi2_map[p] * this->normals[p];
170  const Real M = -this->d2xyzdxideta_map[p] * this->normals[p];
171  const Real N = -this->d2xyzdeta2_map[p] * this->normals[p];
172  const Real E = this->dxyzdxi_map[p].size_sq();
173  const Real F = this->dxyzdxi_map[p] * this->dxyzdeta_map[p];
174  const Real G = this->dxyzdeta_map[p].size_sq();
175 
176  const Real numerator = E*N -2.*F*M + G*L;
177  const Real denominator = E*G - F*F;
178  libmesh_assert_not_equal_to (denominator, 0.);
179  this->curvatures[p] = 0.5*numerator/denominator;
180  }
181 
182  // compute the jacobian at the quadrature points, see
183  // http://sp81.msi.umn.edu:999/fluent/fidap/help/theory/thtoc.htm
184  for (unsigned int p=0; p<n_qp; p++)
185  {
186  const Real g11 = (this->dxdxi_map(p)*this->dxdxi_map(p) +
187  this->dydxi_map(p)*this->dydxi_map(p) +
188  this->dzdxi_map(p)*this->dzdxi_map(p));
189 
190  const Real g12 = (this->dxdxi_map(p)*this->dxdeta_map(p) +
191  this->dydxi_map(p)*this->dydeta_map(p) +
192  this->dzdxi_map(p)*this->dzdeta_map(p));
193 
194  const Real g21 = g12;
195 
196  const Real g22 = (this->dxdeta_map(p)*this->dxdeta_map(p) +
197  this->dydeta_map(p)*this->dydeta_map(p) +
198  this->dzdeta_map(p)*this->dzdeta_map(p));
199 
200 
201  const Real the_jac = std::sqrt(g11*g22 - g12*g21);
202 
203  libmesh_assert_greater (the_jac, 0.);
204 
205  this->JxW[p] = the_jac*qw[p];
206  }
207 
208  break;
209  }
210  default:
211  libmesh_error();
212 
213  } // switch(dim)
214 
215  STOP_LOG("compute_face_map()", "FEXYZMap");
216 
217  return;
218 }
void libMesh::FEMap::compute_map ( const unsigned int  dim,
const std::vector< Real > &  qw,
const Elem elem 
)
virtualinherited

Compute the jacobian and some other additional data fields. Takes the integration weights as input, along with a pointer to the element.

Definition at line 784 of file fe_map.C.

References libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::Elem::has_affine_map(), libMesh::libmesh_assert(), libMesh::FEMap::resize_quadrature_map_vectors(), libMesh::START_LOG(), and libMesh::STOP_LOG().

787 {
788  if (elem->has_affine_map())
789  {
790  compute_affine_map(dim, qw, elem);
791  return;
792  }
793 
794  // Start logging the map computation.
795  START_LOG("compute_map()", "FEMap");
796 
797  libmesh_assert(elem);
798 
799  const unsigned int n_qp = libmesh_cast_int<unsigned int>(qw.size());
800 
801  // Resize the vectors to hold data at the quadrature points
802  this->resize_quadrature_map_vectors(dim, n_qp);
803 
804  // Compute map at all quadrature points
805  for (unsigned int p=0; p!=n_qp; p++)
806  this->compute_single_point_map(dim, qw, elem, p);
807 
808  // Stop logging the map computation.
809  STOP_LOG("compute_map()", "FEMap");
810 }
void libMesh::FEMap::compute_single_point_map ( const unsigned int  dim,
const std::vector< Real > &  qw,
const Elem elem,
unsigned int  p 
)
inherited

Compute the jacobian and some other additional data fields at the single point with index p.

Definition at line 317 of file fe_map.C.

References libMesh::FEMap::d2phideta2_map, libMesh::FEMap::d2phidetadzeta_map, libMesh::FEMap::d2phidxi2_map, libMesh::FEMap::d2phidxideta_map, libMesh::FEMap::d2phidxidzeta_map, libMesh::FEMap::d2phidzeta2_map, libMesh::FEMap::d2xyzdeta2_map, libMesh::FEMap::d2xyzdetadzeta_map, libMesh::FEMap::d2xyzdxi2_map, libMesh::FEMap::d2xyzdxideta_map, libMesh::FEMap::d2xyzdxidzeta_map, libMesh::FEMap::d2xyzdzeta2_map, libMesh::FEMap::detadx_map, libMesh::FEMap::detady_map, libMesh::FEMap::detadz_map, libMesh::FEMap::dphideta_map, libMesh::FEMap::dphidxi_map, libMesh::FEMap::dphidzeta_map, libMesh::FEMap::dxdeta_map(), libMesh::FEMap::dxdxi_map(), libMesh::FEMap::dxdzeta_map(), libMesh::FEMap::dxidx_map, libMesh::FEMap::dxidy_map, libMesh::FEMap::dxidz_map, libMesh::FEMap::dxyzdeta_map, libMesh::FEMap::dxyzdxi_map, libMesh::FEMap::dxyzdzeta_map, libMesh::FEMap::dydeta_map(), libMesh::FEMap::dydxi_map(), libMesh::FEMap::dydzeta_map(), libMesh::FEMap::dzdeta_map(), libMesh::FEMap::dzdxi_map(), libMesh::FEMap::dzdzeta_map(), libMesh::FEMap::dzetadx_map, libMesh::FEMap::dzetady_map, libMesh::FEMap::dzetadz_map, libMesh::err, libMesh::DofObject::id(), libMesh::FEMap::jac, libMesh::FEMap::JxW, libMesh::libmesh_assert(), libMesh::FEMap::phi_map, libMesh::Elem::point(), libMesh::Real, and libMesh::FEMap::xyz.

Referenced by libMesh::FEMap::compute_affine_map(), and libMesh::FEMap::compute_map().

321 {
322  libmesh_assert(elem);
323 
324  switch (dim)
325  {
326  //--------------------------------------------------------------------
327  // 0D
328  case 0:
329  {
330  xyz[p] = elem->point(0);
331  jac[p] = 1.0;
332  JxW[p] = qw[p];
333  break;
334  }
335 
336  //--------------------------------------------------------------------
337  // 1D
338  case 1:
339  {
340  // Clear the entities that will be summed
341  xyz[p].zero();
342  dxyzdxi_map[p].zero();
343 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
344  d2xyzdxi2_map[p].zero();
345 #endif
346 
347  // compute x, dx, d2x at the quadrature point
348  for (unsigned int i=0; i<phi_map.size(); i++) // sum over the nodes
349  {
350  // Reference to the point, helps eliminate
351  // exessive temporaries in the inner loop
352  const Point& elem_point = elem->point(i);
353 
354  xyz[p].add_scaled (elem_point, phi_map[i][p] );
355  dxyzdxi_map[p].add_scaled (elem_point, dphidxi_map[i][p]);
356 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
357  d2xyzdxi2_map[p].add_scaled(elem_point, d2phidxi2_map[i][p]);
358 #endif
359  }
360 
361  // Compute the jacobian
362  //
363  // 1D elements can live in 2D or 3D space.
364  // The transformation matrix from local->global
365  // coordinates is
366  //
367  // T = | dx/dxi |
368  // | dy/dxi |
369  // | dz/dxi |
370  //
371  // The generalized determinant of T (from the
372  // so-called "normal" eqns.) is
373  // jac = "det(T)" = sqrt(det(T'T))
374  //
375  // where T'= transpose of T, so
376  //
377  // jac = sqrt( (dx/dxi)^2 + (dy/dxi)^2 + (dz/dxi)^2 )
378  jac[p] = dxyzdxi_map[p].size();
379 
380  if (jac[p] <= 0.)
381  {
382  libMesh::err << "ERROR: negative Jacobian: "
383  << jac[p]
384  << " in element "
385  << elem->id()
386  << std::endl;
387  libmesh_error();
388  }
389 
390  // The inverse Jacobian entries also come from the
391  // generalized inverse of T (see also the 2D element
392  // living in 3D code).
393  const Real jacm2 = 1./jac[p]/jac[p];
394  dxidx_map[p] = jacm2*dxdxi_map(p);
395 #if LIBMESH_DIM > 1
396  dxidy_map[p] = jacm2*dydxi_map(p);
397 #endif
398 #if LIBMESH_DIM > 2
399  dxidz_map[p] = jacm2*dzdxi_map(p);
400 #endif
401 
402  JxW[p] = jac[p]*qw[p];
403 
404  // done computing the map
405  break;
406  }
407 
408 
409  //--------------------------------------------------------------------
410  // 2D
411  case 2:
412  {
413  //------------------------------------------------------------------
414  // Compute the (x,y) values at the quadrature points,
415  // the Jacobian at the quadrature points
416 
417  xyz[p].zero();
418 
419  dxyzdxi_map[p].zero();
420  dxyzdeta_map[p].zero();
421 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
422  d2xyzdxi2_map[p].zero();
423  d2xyzdxideta_map[p].zero();
424  d2xyzdeta2_map[p].zero();
425 #endif
426 
427 
428  // compute (x,y) at the quadrature points, derivatives once
429  for (unsigned int i=0; i<phi_map.size(); i++) // sum over the nodes
430  {
431  // Reference to the point, helps eliminate
432  // exessive temporaries in the inner loop
433  const Point& elem_point = elem->point(i);
434 
435  xyz[p].add_scaled (elem_point, phi_map[i][p] );
436 
437  dxyzdxi_map[p].add_scaled (elem_point, dphidxi_map[i][p] );
438  dxyzdeta_map[p].add_scaled (elem_point, dphideta_map[i][p]);
439 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
440  d2xyzdxi2_map[p].add_scaled (elem_point, d2phidxi2_map[i][p]);
441  d2xyzdxideta_map[p].add_scaled (elem_point, d2phidxideta_map[i][p]);
442  d2xyzdeta2_map[p].add_scaled (elem_point, d2phideta2_map[i][p]);
443 #endif
444  }
445 
446  // compute the jacobian once
447  const Real dx_dxi = dxdxi_map(p),
448  dx_deta = dxdeta_map(p),
449  dy_dxi = dydxi_map(p),
450  dy_deta = dydeta_map(p);
451 
452 #if LIBMESH_DIM == 2
453  // Compute the Jacobian. This assumes the 2D face
454  // lives in 2D space
455  //
456  // Symbolically, the matrix determinant is
457  //
458  // | dx/dxi dx/deta |
459  // jac = | dy/dxi dy/deta |
460  //
461  // jac = dx/dxi*dy/deta - dx/deta*dy/dxi
462  jac[p] = (dx_dxi*dy_deta - dx_deta*dy_dxi);
463 
464  if (jac[p] <= 0.)
465  {
466  libMesh::err << "ERROR: negative Jacobian: "
467  << jac[p]
468  << " in element "
469  << elem->id()
470  << std::endl;
471  libmesh_error();
472  }
473 
474  JxW[p] = jac[p]*qw[p];
475 
476  // Compute the shape function derivatives wrt x,y at the
477  // quadrature points
478  const Real inv_jac = 1./jac[p];
479 
480  dxidx_map[p] = dy_deta*inv_jac; //dxi/dx = (1/J)*dy/deta
481  dxidy_map[p] = -dx_deta*inv_jac; //dxi/dy = -(1/J)*dx/deta
482  detadx_map[p] = -dy_dxi* inv_jac; //deta/dx = -(1/J)*dy/dxi
483  detady_map[p] = dx_dxi* inv_jac; //deta/dy = (1/J)*dx/dxi
484 
485  dxidz_map[p] = detadz_map[p] = 0.;
486 #else
487 
488  const Real dz_dxi = dzdxi_map(p),
489  dz_deta = dzdeta_map(p);
490 
491  // Compute the Jacobian. This assumes a 2D face in
492  // 3D space.
493  //
494  // The transformation matrix T from local to global
495  // coordinates is
496  //
497  // | dx/dxi dx/deta |
498  // T = | dy/dxi dy/deta |
499  // | dz/dxi dz/deta |
500  // note det(T' T) = det(T')det(T) = det(T)det(T)
501  // so det(T) = std::sqrt(det(T' T))
502  //
503  //----------------------------------------------
504  // Notes:
505  //
506  // dX = R dXi -> R'dX = R'R dXi
507  // (R^-1)dX = dXi [(R'R)^-1 R']dX = dXi
508  //
509  // so R^-1 = (R'R)^-1 R'
510  //
511  // and R^-1 R = (R'R)^-1 R'R = I.
512  //
513  const Real g11 = (dx_dxi*dx_dxi +
514  dy_dxi*dy_dxi +
515  dz_dxi*dz_dxi);
516 
517  const Real g12 = (dx_dxi*dx_deta +
518  dy_dxi*dy_deta +
519  dz_dxi*dz_deta);
520 
521  const Real g21 = g12;
522 
523  const Real g22 = (dx_deta*dx_deta +
524  dy_deta*dy_deta +
525  dz_deta*dz_deta);
526 
527  const Real det = (g11*g22 - g12*g21);
528 
529  if (det <= 0.)
530  {
531  libMesh::err << "ERROR: negative Jacobian! "
532  << " in element "
533  << elem->id()
534  << std::endl;
535  libmesh_error();
536  }
537 
538  const Real inv_det = 1./det;
539  jac[p] = std::sqrt(det);
540 
541  JxW[p] = jac[p]*qw[p];
542 
543  const Real g11inv = g22*inv_det;
544  const Real g12inv = -g12*inv_det;
545  const Real g21inv = -g21*inv_det;
546  const Real g22inv = g11*inv_det;
547 
548  dxidx_map[p] = g11inv*dx_dxi + g12inv*dx_deta;
549  dxidy_map[p] = g11inv*dy_dxi + g12inv*dy_deta;
550  dxidz_map[p] = g11inv*dz_dxi + g12inv*dz_deta;
551 
552  detadx_map[p] = g21inv*dx_dxi + g22inv*dx_deta;
553  detady_map[p] = g21inv*dy_dxi + g22inv*dy_deta;
554  detadz_map[p] = g21inv*dz_dxi + g22inv*dz_deta;
555 
556 #endif
557  // done computing the map
558  break;
559  }
560 
561 
562 
563  //--------------------------------------------------------------------
564  // 3D
565  case 3:
566  {
567  //------------------------------------------------------------------
568  // Compute the (x,y,z) values at the quadrature points,
569  // the Jacobian at the quadrature point
570 
571  // Clear the entities that will be summed
572  xyz[p].zero ();
573  dxyzdxi_map[p].zero ();
574  dxyzdeta_map[p].zero ();
575  dxyzdzeta_map[p].zero ();
576 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
577  d2xyzdxi2_map[p].zero();
578  d2xyzdxideta_map[p].zero();
579  d2xyzdxidzeta_map[p].zero();
580  d2xyzdeta2_map[p].zero();
581  d2xyzdetadzeta_map[p].zero();
582  d2xyzdzeta2_map[p].zero();
583 #endif
584 
585 
586  // compute (x,y,z) at the quadrature points,
587  // dxdxi, dydxi, dzdxi,
588  // dxdeta, dydeta, dzdeta,
589  // dxdzeta, dydzeta, dzdzeta all once
590  for (unsigned int i=0; i<phi_map.size(); i++) // sum over the nodes
591  {
592  // Reference to the point, helps eliminate
593  // exessive temporaries in the inner loop
594  const Point& elem_point = elem->point(i);
595 
596  xyz[p].add_scaled (elem_point, phi_map[i][p] );
597  dxyzdxi_map[p].add_scaled (elem_point, dphidxi_map[i][p] );
598  dxyzdeta_map[p].add_scaled (elem_point, dphideta_map[i][p] );
599  dxyzdzeta_map[p].add_scaled (elem_point, dphidzeta_map[i][p]);
600 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
601  d2xyzdxi2_map[p].add_scaled (elem_point,
602  d2phidxi2_map[i][p]);
603  d2xyzdxideta_map[p].add_scaled (elem_point,
604  d2phidxideta_map[i][p]);
605  d2xyzdxidzeta_map[p].add_scaled (elem_point,
606  d2phidxidzeta_map[i][p]);
607  d2xyzdeta2_map[p].add_scaled (elem_point,
608  d2phideta2_map[i][p]);
609  d2xyzdetadzeta_map[p].add_scaled (elem_point,
610  d2phidetadzeta_map[i][p]);
611  d2xyzdzeta2_map[p].add_scaled (elem_point,
612  d2phidzeta2_map[i][p]);
613 #endif
614  }
615 
616  // compute the jacobian
617  const Real
618  dx_dxi = dxdxi_map(p), dy_dxi = dydxi_map(p), dz_dxi = dzdxi_map(p),
619  dx_deta = dxdeta_map(p), dy_deta = dydeta_map(p), dz_deta = dzdeta_map(p),
620  dx_dzeta = dxdzeta_map(p), dy_dzeta = dydzeta_map(p), dz_dzeta = dzdzeta_map(p);
621 
622  // Symbolically, the matrix determinant is
623  //
624  // | dx/dxi dy/dxi dz/dxi |
625  // jac = | dx/deta dy/deta dz/deta |
626  // | dx/dzeta dy/dzeta dz/dzeta |
627  //
628  // jac = dx/dxi*(dy/deta*dz/dzeta - dz/deta*dy/dzeta) +
629  // dy/dxi*(dz/deta*dx/dzeta - dx/deta*dz/dzeta) +
630  // dz/dxi*(dx/deta*dy/dzeta - dy/deta*dx/dzeta)
631 
632  jac[p] = (dx_dxi*(dy_deta*dz_dzeta - dz_deta*dy_dzeta) +
633  dy_dxi*(dz_deta*dx_dzeta - dx_deta*dz_dzeta) +
634  dz_dxi*(dx_deta*dy_dzeta - dy_deta*dx_dzeta));
635 
636  if (jac[p] <= 0.)
637  {
638  libMesh::err << "ERROR: negative Jacobian: "
639  << jac[p]
640  << " in element "
641  << elem->id()
642  << std::endl;
643  libmesh_error();
644  }
645 
646  JxW[p] = jac[p]*qw[p];
647 
648  // Compute the shape function derivatives wrt x,y at the
649  // quadrature points
650  const Real inv_jac = 1./jac[p];
651 
652  dxidx_map[p] = (dy_deta*dz_dzeta - dz_deta*dy_dzeta)*inv_jac;
653  dxidy_map[p] = (dz_deta*dx_dzeta - dx_deta*dz_dzeta)*inv_jac;
654  dxidz_map[p] = (dx_deta*dy_dzeta - dy_deta*dx_dzeta)*inv_jac;
655 
656  detadx_map[p] = (dz_dxi*dy_dzeta - dy_dxi*dz_dzeta )*inv_jac;
657  detady_map[p] = (dx_dxi*dz_dzeta - dz_dxi*dx_dzeta )*inv_jac;
658  detadz_map[p] = (dy_dxi*dx_dzeta - dx_dxi*dy_dzeta )*inv_jac;
659 
660  dzetadx_map[p] = (dy_dxi*dz_deta - dz_dxi*dy_deta )*inv_jac;
661  dzetady_map[p] = (dz_dxi*dx_deta - dx_dxi*dz_deta )*inv_jac;
662  dzetadz_map[p] = (dx_dxi*dy_deta - dy_dxi*dx_deta )*inv_jac;
663 
664  // done computing the map
665  break;
666  }
667 
668  default:
669  libmesh_error();
670  }
671 }
Real libMesh::FEMap::dxdeta_map ( const unsigned int  p) const
inlineprotectedinherited

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 458 of file fe_map.h.

References libMesh::FEMap::dxyzdeta_map.

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

458 { return dxyzdeta_map[p](0); }
Real libMesh::FEMap::dxdxi_map ( const unsigned int  p) const
inlineprotectedinherited

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 437 of file fe_map.h.

References libMesh::FEMap::dxyzdxi_map.

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

437 { return dxyzdxi_map[p](0); }
Real libMesh::FEMap::dxdzeta_map ( const unsigned int  p) const
inlineprotectedinherited

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 479 of file fe_map.h.

References libMesh::FEMap::dxyzdzeta_map.

Referenced by libMesh::FEMap::compute_single_point_map().

479 { return dxyzdzeta_map[p](0); }
Real libMesh::FEMap::dydeta_map ( const unsigned int  p) const
inlineprotectedinherited

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 465 of file fe_map.h.

References libMesh::FEMap::dxyzdeta_map.

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

465 { return dxyzdeta_map[p](1); }
Real libMesh::FEMap::dydxi_map ( const unsigned int  p) const
inlineprotectedinherited

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 444 of file fe_map.h.

References libMesh::FEMap::dxyzdxi_map.

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

444 { return dxyzdxi_map[p](1); }
Real libMesh::FEMap::dydzeta_map ( const unsigned int  p) const
inlineprotectedinherited

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 486 of file fe_map.h.

References libMesh::FEMap::dxyzdzeta_map.

Referenced by libMesh::FEMap::compute_single_point_map().

486 { return dxyzdzeta_map[p](1); }
Real libMesh::FEMap::dzdeta_map ( const unsigned int  p) const
inlineprotectedinherited

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 472 of file fe_map.h.

References libMesh::FEMap::dxyzdeta_map.

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

472 { return dxyzdeta_map[p](2); }
Real libMesh::FEMap::dzdxi_map ( const unsigned int  p) const
inlineprotectedinherited

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 451 of file fe_map.h.

References libMesh::FEMap::dxyzdxi_map.

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

451 { return dxyzdxi_map[p](2); }
Real libMesh::FEMap::dzdzeta_map ( const unsigned int  p) const
inlineprotectedinherited

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 493 of file fe_map.h.

References libMesh::FEMap::dxyzdzeta_map.

Referenced by libMesh::FEMap::compute_single_point_map().

493 { return dxyzdzeta_map[p](2); }
const std::vector<Real>& libMesh::FEMap::get_curvatures ( ) const
inlineinherited
Returns
the curvatures for use in face integration.

Definition at line 300 of file fe_map.h.

References libMesh::FEMap::curvatures.

301  { return curvatures;}
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phideta2_map ( )
inlineinherited
Returns
the reference to physical map 2nd derivative

Definition at line 400 of file fe_map.h.

References libMesh::FEMap::d2phideta2_map.

401  { return d2phideta2_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidetadzeta_map ( )
inlineinherited
Returns
the reference to physical map 2nd derivative

Definition at line 406 of file fe_map.h.

References libMesh::FEMap::d2phidetadzeta_map.

407  { return d2phidetadzeta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidxi2_map ( )
inlineinherited
Returns
the reference to physical map 2nd derivative

Definition at line 382 of file fe_map.h.

References libMesh::FEMap::d2phidxi2_map.

383  { return d2phidxi2_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidxideta_map ( )
inlineinherited
Returns
the reference to physical map 2nd derivative

Definition at line 388 of file fe_map.h.

References libMesh::FEMap::d2phidxideta_map.

389  { return d2phidxideta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidxidzeta_map ( )
inlineinherited
Returns
the reference to physical map 2nd derivative

Definition at line 394 of file fe_map.h.

References libMesh::FEMap::d2phidxidzeta_map.

395  { return d2phidxidzeta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2phidzeta2_map ( )
inlineinherited
Returns
the reference to physical map 2nd derivative

Definition at line 412 of file fe_map.h.

References libMesh::FEMap::d2phidzeta2_map.

413  { return d2phidzeta2_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2psideta2 ( )
inlineinherited
Returns
the reference to physical map 2nd derivative for the side/edge

Definition at line 351 of file fe_map.h.

References libMesh::FEMap::d2psideta2_map.

352  { return d2psideta2_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2psidxi2 ( )
inlineinherited
Returns
the reference to physical map 2nd derivative for the side/edge

Definition at line 339 of file fe_map.h.

References libMesh::FEMap::d2psidxi2_map.

340  { return d2psidxi2_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_d2psidxideta ( )
inlineinherited
Returns
the reference to physical map 2nd derivative for the side/edge

Definition at line 345 of file fe_map.h.

References libMesh::FEMap::d2psidxideta_map.

346  { return d2psidxideta_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdeta2 ( ) const
inlineinherited
Returns
the second partial derivatives in eta.

Definition at line 157 of file fe_map.h.

References libMesh::FEMap::d2xyzdeta2_map.

158  { return d2xyzdeta2_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdetadzeta ( ) const
inlineinherited
Returns
the second partial derivatives in eta-zeta.

Definition at line 187 of file fe_map.h.

References libMesh::FEMap::d2xyzdetadzeta_map.

188  { return d2xyzdetadzeta_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdxi2 ( ) const
inlineinherited
Returns
the second partial derivatives in xi.

Definition at line 151 of file fe_map.h.

References libMesh::FEMap::d2xyzdxi2_map.

152  { return d2xyzdxi2_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdxideta ( ) const
inlineinherited
Returns
the second partial derivatives in xi-eta.

Definition at line 173 of file fe_map.h.

References libMesh::FEMap::d2xyzdxideta_map.

174  { return d2xyzdxideta_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdxidzeta ( ) const
inlineinherited
Returns
the second partial derivatives in xi-zeta.

Definition at line 181 of file fe_map.h.

References libMesh::FEMap::d2xyzdxidzeta_map.

182  { return d2xyzdxidzeta_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_d2xyzdzeta2 ( ) const
inlineinherited
Returns
the second partial derivatives in zeta.

Definition at line 165 of file fe_map.h.

References libMesh::FEMap::d2xyzdzeta2_map.

166  { return d2xyzdzeta2_map; }
const std::vector<Real>& libMesh::FEMap::get_detadx ( ) const
inlineinherited
Returns
the deta/dx entry in the transformation matrix from physical to local coordinates.

Definition at line 217 of file fe_map.h.

References libMesh::FEMap::detadx_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

218  { return detadx_map; }
const std::vector<Real>& libMesh::FEMap::get_detady ( ) const
inlineinherited
Returns
the deta/dy entry in the transformation matrix from physical to local coordinates.

Definition at line 224 of file fe_map.h.

References libMesh::FEMap::detady_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

225  { return detady_map; }
const std::vector<Real>& libMesh::FEMap::get_detadz ( ) const
inlineinherited
Returns
the deta/dz entry in the transformation matrix from physical to local coordinates.

Definition at line 231 of file fe_map.h.

References libMesh::FEMap::detadz_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

232  { return detadz_map; }
const std::vector<std::vector<Real> >& libMesh::FEMap::get_dphideta_map ( ) const
inlineinherited
Returns
the reference to physical map derivative

Definition at line 276 of file fe_map.h.

References libMesh::FEMap::dphideta_map.

277  { return dphideta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_dphideta_map ( )
inlineinherited
Returns
the reference to physical map derivative

Definition at line 369 of file fe_map.h.

References libMesh::FEMap::dphideta_map.

370  { return dphideta_map; }
const std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidxi_map ( ) const
inlineinherited
Returns
the reference to physical map derivative

Definition at line 270 of file fe_map.h.

References libMesh::FEMap::dphidxi_map.

271  { return dphidxi_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidxi_map ( )
inlineinherited
Returns
the reference to physical map derivative

Definition at line 363 of file fe_map.h.

References libMesh::FEMap::dphidxi_map.

364  { return dphidxi_map; }
const std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidzeta_map ( ) const
inlineinherited
Returns
the reference to physical map derivative

Definition at line 282 of file fe_map.h.

References libMesh::FEMap::dphidzeta_map.

283  { return dphidzeta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_dphidzeta_map ( )
inlineinherited
Returns
the reference to physical map derivative

Definition at line 375 of file fe_map.h.

References libMesh::FEMap::dphidzeta_map.

376  { return dphidzeta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_dpsideta ( )
inlineinherited
Returns
the reference to physical map derivative for the side/edge

Definition at line 333 of file fe_map.h.

References libMesh::FEMap::dpsideta_map.

334  { return dpsideta_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_dpsidxi ( )
inlineinherited
Returns
the reference to physical map derivative for the side/edge

Definition at line 327 of file fe_map.h.

References libMesh::FEMap::dpsidxi_map.

328  { return dpsidxi_map; }
const std::vector<Real>& libMesh::FEMap::get_dxidx ( ) const
inlineinherited
Returns
the dxi/dx entry in the transformation matrix from physical to local coordinates.

Definition at line 196 of file fe_map.h.

References libMesh::FEMap::dxidx_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

197  { return dxidx_map; }
const std::vector<Real>& libMesh::FEMap::get_dxidy ( ) const
inlineinherited
Returns
the dxi/dy entry in the transformation matrix from physical to local coordinates.

Definition at line 203 of file fe_map.h.

References libMesh::FEMap::dxidy_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

204  { return dxidy_map; }
const std::vector<Real>& libMesh::FEMap::get_dxidz ( ) const
inlineinherited
Returns
the dxi/dz entry in the transformation matrix from physical to local coordinates.

Definition at line 210 of file fe_map.h.

References libMesh::FEMap::dxidz_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

211  { return dxidz_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_dxyzdeta ( ) const
inlineinherited
Returns
the element tangents in eta-direction at the quadrature points.

Definition at line 138 of file fe_map.h.

References libMesh::FEMap::dxyzdeta_map.

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

139  { return dxyzdeta_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_dxyzdxi ( ) const
inlineinherited
Returns
the element tangents in xi-direction at the quadrature points.

Definition at line 131 of file fe_map.h.

References libMesh::FEMap::dxyzdxi_map.

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

132  { return dxyzdxi_map; }
const std::vector<RealGradient>& libMesh::FEMap::get_dxyzdzeta ( ) const
inlineinherited
Returns
the element tangents in zeta-direction at the quadrature points.

Definition at line 145 of file fe_map.h.

References libMesh::FEMap::dxyzdzeta_map.

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

146  { return dxyzdzeta_map; }
const std::vector<Real>& libMesh::FEMap::get_dzetadx ( ) const
inlineinherited
Returns
the dzeta/dx entry in the transformation matrix from physical to local coordinates.

Definition at line 238 of file fe_map.h.

References libMesh::FEMap::dzetadx_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

239  { return dzetadx_map; }
const std::vector<Real>& libMesh::FEMap::get_dzetady ( ) const
inlineinherited
Returns
the dzeta/dy entry in the transformation matrix from physical to local coordinates.

Definition at line 245 of file fe_map.h.

References libMesh::FEMap::dzetady_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

246  { return dzetady_map; }
const std::vector<Real>& libMesh::FEMap::get_dzetadz ( ) const
inlineinherited
Returns
the dzeta/dz entry in the transformation matrix from physical to local coordinates.

Definition at line 252 of file fe_map.h.

References libMesh::FEMap::dzetadz_map.

Referenced by libMesh::H1FETransformation< T >::map_curl(), libMesh::H1FETransformation< T >::map_d2phi(), libMesh::H1FETransformation< T >::map_div(), libMesh::H1FETransformation< T >::map_dphi(), and libMesh::HCurlFETransformation< T >::map_phi().

253  { return dzetadz_map; }
const std::vector<Real>& libMesh::FEMap::get_jacobian ( ) const
inlineinherited
Returns
the element Jacobian for each quadrature point.

Definition at line 117 of file fe_map.h.

References libMesh::FEMap::jac.

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

118  { return jac; }
const std::vector<Real>& libMesh::FEMap::get_JxW ( ) const
inlineinherited
Returns
the element Jacobian times the quadrature weight for each quadrature point.

Definition at line 124 of file fe_map.h.

References libMesh::FEMap::JxW.

125  { return JxW; }
std::vector<Real>& libMesh::FEMap::get_JxW ( )
inlineinherited
Returns
writable reference to the element Jacobian times the quadrature weight for each quadrature point.

Definition at line 422 of file fe_map.h.

References libMesh::FEMap::JxW.

423  { return JxW; }
const std::vector<Point>& libMesh::FEMap::get_normals ( ) const
inlineinherited
Returns
the normal vectors for face integration.

Definition at line 294 of file fe_map.h.

References libMesh::FEMap::normals.

295  { return normals; }
const std::vector<std::vector<Real> >& libMesh::FEMap::get_phi_map ( ) const
inlineinherited
Returns
the reference to physical map for the element

Definition at line 264 of file fe_map.h.

References libMesh::FEMap::phi_map.

265  { return phi_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_phi_map ( )
inlineinherited
Returns
the reference to physical map for the element

Definition at line 357 of file fe_map.h.

References libMesh::FEMap::phi_map.

358  { return phi_map; }
const std::vector<std::vector<Real> >& libMesh::FEMap::get_psi ( ) const
inlineinherited
Returns
the reference to physical map for the side/edge

Definition at line 258 of file fe_map.h.

References libMesh::FEMap::psi_map.

259  { return psi_map; }
std::vector<std::vector<Real> >& libMesh::FEMap::get_psi ( )
inlineinherited
Returns
the reference to physical map for the side/edge

Definition at line 321 of file fe_map.h.

References libMesh::FEMap::psi_map.

322  { return psi_map; }
const std::vector<std::vector<Point> >& libMesh::FEMap::get_tangents ( ) const
inlineinherited
Returns
the tangent vectors for face integration.

Definition at line 288 of file fe_map.h.

References libMesh::FEMap::tangents.

289  { return tangents; }
const std::vector<Point>& libMesh::FEMap::get_xyz ( ) const
inlineinherited
Returns
the xyz spatial locations of the quadrature points on the element.

Definition at line 111 of file fe_map.h.

References libMesh::FEMap::xyz.

112  { return xyz; }
template<unsigned int Dim>
void libMesh::FEMap::init_edge_shape_functions ( const std::vector< Point > &  qp,
const Elem edge 
)
inherited

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

Start logging the shape function initialization

Stop logging the shape function initialization

Definition at line 475 of file fe_boundary.C.

References libMesh::FEMap::d2psidxi2_map, libMesh::Elem::default_order(), libMesh::FEMap::dpsidxi_map, libMesh::libmesh_assert(), libMesh::FE< Dim, T >::n_shape_functions(), libMesh::FEMap::psi_map, libMesh::FE< Dim, T >::shape(), libMesh::FE< Dim, T >::shape_deriv(), libMesh::FE< Dim, T >::shape_second_deriv(), libMesh::START_LOG(), libMesh::STOP_LOG(), and libMesh::Elem::type().

477 {
478  libmesh_assert(edge);
479 
483  START_LOG("init_edge_shape_functions()", "FEMap");
484 
485  // The element type and order to use in
486  // the map
487  const Order mapping_order (edge->default_order());
488  const ElemType mapping_elem_type (edge->type());
489 
490  // The number of quadrature points.
491  const unsigned int n_qp = libmesh_cast_int<unsigned int>(qp.size());
492 
493  const unsigned int n_mapping_shape_functions =
494  FE<Dim,LAGRANGE>::n_shape_functions (mapping_elem_type,
495  mapping_order);
496 
497  // resize the vectors to hold current data
498  // Psi are the shape functions used for the FE mapping
499  this->psi_map.resize (n_mapping_shape_functions);
500  this->dpsidxi_map.resize (n_mapping_shape_functions);
501  this->d2psidxi2_map.resize (n_mapping_shape_functions);
502 
503  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
504  {
505  // Allocate space to store the values of the shape functions
506  // and their first and second derivatives at the quadrature points.
507  this->psi_map[i].resize (n_qp);
508  this->dpsidxi_map[i].resize (n_qp);
509  this->d2psidxi2_map[i].resize (n_qp);
510 
511  // Compute the value of shape function i, and its first and
512  // second derivatives at quadrature point p
513  // (Lagrange shape functions are used for the mapping)
514  for (unsigned int p=0; p<n_qp; p++)
515  {
516  this->psi_map[i][p] = FE<1,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
517  this->dpsidxi_map[i][p] = FE<1,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
518  this->d2psidxi2_map[i][p] = FE<1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 0, qp[p]);
519  }
520  }
521 
525  STOP_LOG("init_edge_shape_functions()", "FEMap");
526 }
template<unsigned int Dim>
void libMesh::FEMap::init_face_shape_functions ( const std::vector< Point > &  qp,
const Elem side 
)
inherited

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

Start logging the shape function initialization

Stop logging the shape function initialization

Definition at line 386 of file fe_boundary.C.

References libMesh::FEMap::d2psideta2_map, libMesh::FEMap::d2psidxi2_map, libMesh::FEMap::d2psidxideta_map, libMesh::Elem::default_order(), libMesh::FEMap::dpsideta_map, libMesh::FEMap::dpsidxi_map, libMesh::libmesh_assert(), libMesh::FE< Dim, T >::n_shape_functions(), libMesh::FEMap::psi_map, libMesh::FE< Dim, T >::shape(), libMesh::FE< Dim, T >::shape_deriv(), libMesh::FE< Dim, T >::shape_second_deriv(), libMesh::START_LOG(), libMesh::STOP_LOG(), and libMesh::Elem::type().

388 {
390 
394  START_LOG("init_face_shape_functions()", "FEMap");
395 
396  // The element type and order to use in
397  // the map
398  const Order mapping_order (side->default_order());
399  const ElemType mapping_elem_type (side->type());
400 
401  // The number of quadrature points.
402  const unsigned int n_qp = libmesh_cast_int<unsigned int>(qp.size());
403 
404  const unsigned int n_mapping_shape_functions =
405  FE<Dim,LAGRANGE>::n_shape_functions (mapping_elem_type,
406  mapping_order);
407 
408  // resize the vectors to hold current data
409  // Psi are the shape functions used for the FE mapping
410  this->psi_map.resize (n_mapping_shape_functions);
411 
412  if (Dim > 1)
413  {
414  this->dpsidxi_map.resize (n_mapping_shape_functions);
415  this->d2psidxi2_map.resize (n_mapping_shape_functions);
416  }
417 
418  if (Dim == 3)
419  {
420  this->dpsideta_map.resize (n_mapping_shape_functions);
421  this->d2psidxideta_map.resize (n_mapping_shape_functions);
422  this->d2psideta2_map.resize (n_mapping_shape_functions);
423  }
424 
425  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
426  {
427  // Allocate space to store the values of the shape functions
428  // and their first and second derivatives at the quadrature points.
429  this->psi_map[i].resize (n_qp);
430  if (Dim > 1)
431  {
432  this->dpsidxi_map[i].resize (n_qp);
433  this->d2psidxi2_map[i].resize (n_qp);
434  }
435  if (Dim == 3)
436  {
437  this->dpsideta_map[i].resize (n_qp);
438  this->d2psidxideta_map[i].resize (n_qp);
439  this->d2psideta2_map[i].resize (n_qp);
440  }
441 
442  // Compute the value of shape function i, and its first and
443  // second derivatives at quadrature point p
444  // (Lagrange shape functions are used for the mapping)
445  for (unsigned int p=0; p<n_qp; p++)
446  {
447  this->psi_map[i][p] = FE<Dim-1,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
448  if (Dim > 1)
449  {
450  this->dpsidxi_map[i][p] = FE<Dim-1,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
451  this->d2psidxi2_map[i][p] = FE<Dim-1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 0, qp[p]);
452  }
453  // libMesh::out << "this->d2psidxi2_map["<<i<<"][p]=" << d2psidxi2_map[i][p] << std::endl;
454 
455  // If we are in 3D, then our sides are 2D faces.
456  // For the second derivatives, we must also compute the cross
457  // derivative d^2() / dxi deta
458  if (Dim == 3)
459  {
460  this->dpsideta_map[i][p] = FE<Dim-1,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
461  this->d2psidxideta_map[i][p] = FE<Dim-1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 1, qp[p]);
462  this->d2psideta2_map[i][p] = FE<Dim-1,LAGRANGE>::shape_second_deriv(mapping_elem_type, mapping_order, i, 2, qp[p]);
463  }
464  }
465  }
466 
467 
471  STOP_LOG("init_face_shape_functions()", "FEMap");
472 }
template<unsigned int Dim>
template void libMesh::FEMap::init_reference_to_physical_map< 3 > ( const std::vector< Point > &  qp,
const Elem elem 
)
inherited

Definition at line 65 of file fe_map.C.

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

67 {
68  // Start logging the reference->physical map initialization
69  START_LOG("init_reference_to_physical_map()", "FEMap");
70 
71  // The number of quadrature points.
72  const std::size_t n_qp = qp.size();
73 
74  // The element type and order to use in
75  // the map
76  const Order mapping_order (elem->default_order());
77  const ElemType mapping_elem_type (elem->type());
78 
79  // Number of shape functions used to construt the map
80  // (Lagrange shape functions are used for mapping)
81  const unsigned int n_mapping_shape_functions =
82  FE<Dim,LAGRANGE>::n_shape_functions (mapping_elem_type,
83  mapping_order);
84 
85  this->phi_map.resize (n_mapping_shape_functions);
86  if (Dim > 0)
87  {
88  this->dphidxi_map.resize (n_mapping_shape_functions);
89 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
90  this->d2phidxi2_map.resize (n_mapping_shape_functions);
91 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
92  }
93 
94  if (Dim > 1)
95  {
96  this->dphideta_map.resize (n_mapping_shape_functions);
97 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
98  this->d2phidxideta_map.resize (n_mapping_shape_functions);
99  this->d2phideta2_map.resize (n_mapping_shape_functions);
100 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
101  }
102 
103  if (Dim > 2)
104  {
105  this->dphidzeta_map.resize (n_mapping_shape_functions);
106 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
107  this->d2phidxidzeta_map.resize (n_mapping_shape_functions);
108  this->d2phidetadzeta_map.resize (n_mapping_shape_functions);
109  this->d2phidzeta2_map.resize (n_mapping_shape_functions);
110 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
111  }
112 
113 
114  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
115  {
116  this->phi_map[i].resize (n_qp);
117  if (Dim > 0)
118  {
119  this->dphidxi_map[i].resize (n_qp);
120 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
121  this->d2phidxi2_map[i].resize (n_qp);
122  if (Dim > 1)
123  {
124  this->d2phidxideta_map[i].resize (n_qp);
125  this->d2phideta2_map[i].resize (n_qp);
126  }
127  if (Dim > 2)
128  {
129  this->d2phidxidzeta_map[i].resize (n_qp);
130  this->d2phidetadzeta_map[i].resize (n_qp);
131  this->d2phidzeta2_map[i].resize (n_qp);
132  }
133 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
134 
135  if (Dim > 1)
136  this->dphideta_map[i].resize (n_qp);
137 
138  if (Dim > 2)
139  this->dphidzeta_map[i].resize (n_qp);
140  }
141  }
142 
143  // Optimize for the *linear* geometric elements case:
144  bool is_linear = elem->is_linear();
145 
146  switch (Dim)
147  {
148 
149  //------------------------------------------------------------
150  // 0D
151  case 0:
152  {
153  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
154  for (std::size_t p=0; p<n_qp; p++)
155  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
156 
157  break;
158  }
159 
160  //------------------------------------------------------------
161  // 1D
162  case 1:
163  {
164  // Compute the value of the mapping shape function i at quadrature point p
165  // (Lagrange shape functions are used for mapping)
166  if (is_linear)
167  {
168  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
169  {
170  this->phi_map[i][0] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[0]);
171  this->dphidxi_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
172 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
173  this->d2phidxi2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
174 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
175  for (std::size_t p=1; p<n_qp; p++)
176  {
177  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
178  this->dphidxi_map[i][p] = this->dphidxi_map[i][0];
179 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
180  this->d2phidxi2_map[i][p] = this->d2phidxi2_map[i][0];
181 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
182  }
183  }
184  }
185  else
186  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
187  for (std::size_t p=0; p<n_qp; p++)
188  {
189  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
190  this->dphidxi_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
191 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
192  this->d2phidxi2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
193 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
194  }
195 
196  break;
197  }
198  //------------------------------------------------------------
199  // 2D
200  case 2:
201  {
202  // Compute the value of the mapping shape function i at quadrature point p
203  // (Lagrange shape functions are used for mapping)
204  if (is_linear)
205  {
206  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
207  {
208  this->phi_map[i][0] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[0]);
209  this->dphidxi_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
210  this->dphideta_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
211 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
212  this->d2phidxi2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
213  this->d2phidxideta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
214  this->d2phideta2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[0]);
215 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
216  for (std::size_t p=1; p<n_qp; p++)
217  {
218  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
219  this->dphidxi_map[i][p] = this->dphidxi_map[i][0];
220  this->dphideta_map[i][p] = this->dphideta_map[i][0];
221 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
222  this->d2phidxi2_map[i][p] = this->d2phidxi2_map[i][0];
223  this->d2phidxideta_map[i][p] = this->d2phidxideta_map[i][0];
224  this->d2phideta2_map[i][p] = this->d2phideta2_map[i][0];
225 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
226  }
227  }
228  }
229  else
230  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
231  for (std::size_t p=0; p<n_qp; p++)
232  {
233  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
234  this->dphidxi_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
235  this->dphideta_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
236 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
237  this->d2phidxi2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
238  this->d2phidxideta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
239  this->d2phideta2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[p]);
240 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
241  }
242 
243  break;
244  }
245 
246  //------------------------------------------------------------
247  // 3D
248  case 3:
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  this->phi_map[i][0] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[0]);
257  this->dphidxi_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
258  this->dphideta_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
259  this->dphidzeta_map[i][0] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 2, qp[0]);
260 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
261  this->d2phidxi2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[0]);
262  this->d2phidxideta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[0]);
263  this->d2phideta2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[0]);
264  this->d2phidxidzeta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 3, qp[0]);
265  this->d2phidetadzeta_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 4, qp[0]);
266  this->d2phidzeta2_map[i][0] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 5, qp[0]);
267 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
268  for (std::size_t p=1; p<n_qp; p++)
269  {
270  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
271  this->dphidxi_map[i][p] = this->dphidxi_map[i][0];
272  this->dphideta_map[i][p] = this->dphideta_map[i][0];
273  this->dphidzeta_map[i][p] = this->dphidzeta_map[i][0];
274 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
275  this->d2phidxi2_map[i][p] = this->d2phidxi2_map[i][0];
276  this->d2phidxideta_map[i][p] = this->d2phidxideta_map[i][0];
277  this->d2phideta2_map[i][p] = this->d2phideta2_map[i][0];
278  this->d2phidxidzeta_map[i][p] = this->d2phidxidzeta_map[i][0];
279  this->d2phidetadzeta_map[i][p] = this->d2phidetadzeta_map[i][0];
280  this->d2phidzeta2_map[i][p] = this->d2phidzeta2_map[i][0];
281 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
282  }
283  }
284  }
285  else
286  for (unsigned int i=0; i<n_mapping_shape_functions; i++)
287  for (std::size_t p=0; p<n_qp; p++)
288  {
289  this->phi_map[i][p] = FE<Dim,LAGRANGE>::shape (mapping_elem_type, mapping_order, i, qp[p]);
290  this->dphidxi_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
291  this->dphideta_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
292  this->dphidzeta_map[i][p] = FE<Dim,LAGRANGE>::shape_deriv (mapping_elem_type, mapping_order, i, 2, qp[p]);
293 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
294  this->d2phidxi2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 0, qp[p]);
295  this->d2phidxideta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 1, qp[p]);
296  this->d2phideta2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 2, qp[p]);
297  this->d2phidxidzeta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 3, qp[p]);
298  this->d2phidetadzeta_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 4, qp[p]);
299  this->d2phidzeta2_map[i][p] = FE<Dim,LAGRANGE>::shape_second_deriv (mapping_elem_type, mapping_order, i, 5, qp[p]);
300 #endif // ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
301  }
302 
303  break;
304  }
305 
306  default:
307  libmesh_error();
308  }
309 
310  // Stop logging the reference->physical map initialization
311  STOP_LOG("init_reference_to_physical_map()", "FEMap");
312  return;
313 }
void libMesh::FEMap::print_JxW ( std::ostream &  os) const
inherited

Prints the Jacobian times the weight for each quadrature point.

Definition at line 814 of file fe_map.C.

References libMesh::FEMap::JxW.

815 {
816  for (unsigned int i=0; i<JxW.size(); ++i)
817  os << " [" << i << "]: " << JxW[i] << std::endl;
818 }
void libMesh::FEMap::print_xyz ( std::ostream &  os) const
inherited

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

Definition at line 822 of file fe_map.C.

References libMesh::FEMap::xyz.

823 {
824  for (unsigned int i=0; i<xyz.size(); ++i)
825  os << " [" << i << "]: " << xyz[i];
826 }
void libMesh::FEMap::resize_quadrature_map_vectors ( const unsigned int  dim,
unsigned int  n_qp 
)
protectedinherited

A utility function for use by compute_*_map

Definition at line 675 of file fe_map.C.

References libMesh::FEMap::d2xyzdeta2_map, libMesh::FEMap::d2xyzdetadzeta_map, libMesh::FEMap::d2xyzdxi2_map, libMesh::FEMap::d2xyzdxideta_map, libMesh::FEMap::d2xyzdxidzeta_map, libMesh::FEMap::d2xyzdzeta2_map, libMesh::FEMap::detadx_map, libMesh::FEMap::detady_map, libMesh::FEMap::detadz_map, libMesh::FEMap::dxidx_map, libMesh::FEMap::dxidy_map, libMesh::FEMap::dxidz_map, libMesh::FEMap::dxyzdeta_map, libMesh::FEMap::dxyzdxi_map, libMesh::FEMap::dxyzdzeta_map, libMesh::FEMap::dzetadx_map, libMesh::FEMap::dzetady_map, libMesh::FEMap::dzetadz_map, libMesh::FEMap::jac, libMesh::FEMap::JxW, and libMesh::FEMap::xyz.

Referenced by libMesh::FEMap::compute_affine_map(), and libMesh::FEMap::compute_map().

676 {
677  // Resize the vectors to hold data at the quadrature points
678  xyz.resize(n_qp);
679  dxyzdxi_map.resize(n_qp);
680  dxidx_map.resize(n_qp);
681  dxidy_map.resize(n_qp); // 1D element may live in 2D ...
682  dxidz_map.resize(n_qp); // ... or 3D
683 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
684  d2xyzdxi2_map.resize(n_qp);
685 #endif
686  if (dim > 1)
687  {
688  dxyzdeta_map.resize(n_qp);
689  detadx_map.resize(n_qp);
690  detady_map.resize(n_qp);
691  detadz_map.resize(n_qp);
692 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
693  d2xyzdxideta_map.resize(n_qp);
694  d2xyzdeta2_map.resize(n_qp);
695 #endif
696  if (dim > 2)
697  {
698  dxyzdzeta_map.resize (n_qp);
699  dzetadx_map.resize (n_qp);
700  dzetady_map.resize (n_qp);
701  dzetadz_map.resize (n_qp);
702 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
703  d2xyzdxidzeta_map.resize(n_qp);
704  d2xyzdetadzeta_map.resize(n_qp);
705  d2xyzdzeta2_map.resize(n_qp);
706 #endif
707  }
708  }
709 
710  jac.resize(n_qp);
711  JxW.resize(n_qp);
712 }

Member Data Documentation

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

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 721 of file fe_map.h.

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

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

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

Definition at line 654 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2phideta2_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 659 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2phidetadzeta_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 639 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2phidxi2_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 644 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2phidxideta_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 649 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2phidxidzeta_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 664 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2phidzeta2_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 704 of file fe_map.h.

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

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

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

Definition at line 690 of file fe_map.h.

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

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

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 697 of file fe_map.h.

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

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

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 548 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2xyzdetadzeta(), and libMesh::FEMap::resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::d2xyzdxi2_map
protectedinherited
std::vector<RealGradient> libMesh::FEMap::d2xyzdxideta_map
protectedinherited

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 528 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_edge_map(), compute_face_map(), libMesh::FEMap::compute_face_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2xyzdxideta(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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 542 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2xyzdxidzeta(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 554 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_d2xyzdzeta2(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 581 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_detadx(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 587 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_detady(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 593 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_detadz(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 627 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dphideta_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 622 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dphidxi_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 632 of file fe_map.h.

Referenced by libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dphidzeta_map(), and libMesh::FEMap::init_reference_to_physical_map().

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

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

Definition at line 683 of file fe_map.h.

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

std::vector<std::vector<Real> > libMesh::FEMap::dpsidxi_map
protectedinherited
std::vector<Real> libMesh::FEMap::dxidx_map
protectedinherited

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

Definition at line 562 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dxidx(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 568 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dxidy(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 574 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dxidz(), and libMesh::FEMap::resize_quadrature_map_vectors().

std::vector<RealGradient> libMesh::FEMap::dxyzdzeta_map
protectedinherited
std::vector<Real> libMesh::FEMap::dzetadx_map
protectedinherited

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

Definition at line 600 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dzetadx(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 606 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dzetady(), and libMesh::FEMap::resize_quadrature_map_vectors().

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

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

Definition at line 612 of file fe_map.h.

Referenced by libMesh::FEMap::compute_affine_map(), libMesh::FEMap::compute_single_point_map(), libMesh::FEMap::get_dzetadz(), and libMesh::FEMap::resize_quadrature_map_vectors().

std::vector<Real> libMesh::FEMap::jac
protectedinherited
std::vector<Point> libMesh::FEMap::normals
protectedinherited

Normal vectors on boundary at quadrature points

Definition at line 714 of file fe_map.h.

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

std::vector<std::vector<Real> > libMesh::FEMap::phi_map
protectedinherited
std::vector<std::vector<Real> > libMesh::FEMap::psi_map
protectedinherited
std::vector<std::vector<Point> > libMesh::FEMap::tangents
protectedinherited

Tangent vectors on boundary at quadrature points.

Definition at line 709 of file fe_map.h.

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


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:58:01 UTC

Hosted By:
SourceForge.net Logo