libMesh
fe_map.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_FE_MAP_H
21 #define LIBMESH_FE_MAP_H
22 
23 #include "libmesh/reference_counted_object.h"
24 #include "libmesh/point.h"
25 #include "libmesh/vector_value.h"
26 #include "libmesh/enum_elem_type.h"
27 #include "libmesh/fe_type.h"
28 #include "libmesh/auto_ptr.h"
29 
30 namespace libMesh
31 {
32 
33 // forward declarations
34 class Elem;
35 class Node;
36 
45 class FEMap
46 {
47 public:
48 
49  FEMap();
50  virtual ~FEMap(){}
51 
52  static UniquePtr<FEMap> build(FEType fe_type);
53 
54  template<unsigned int Dim>
55  void init_reference_to_physical_map(const std::vector<Point> & qp,
56  const Elem * elem);
57 
65  void compute_single_point_map(const unsigned int dim,
66  const std::vector<Real> & qw,
67  const Elem * elem,
68  unsigned int p,
69  const std::vector<const Node *> & elem_nodes,
70  bool compute_second_derivatives);
71 
78  virtual void compute_affine_map(const unsigned int dim,
79  const std::vector<Real> & qw,
80  const Elem * elem);
81 
87  virtual void compute_null_map(const unsigned int dim,
88  const std::vector<Real> & qw);
89 
90 
99  virtual void compute_map(const unsigned int dim,
100  const std::vector<Real> & qw,
101  const Elem * elem,
102  bool calculate_d2phi);
103 
107  virtual void compute_face_map(int dim,
108  const std::vector<Real> & qw,
109  const Elem * side);
110 
114  void compute_edge_map(int dim,
115  const std::vector<Real> & qw,
116  const Elem * side);
117 
122  template<unsigned int Dim>
123  void init_face_shape_functions(const std::vector<Point> & qp,
124  const Elem * side);
125 
130  template<unsigned int Dim>
131  void init_edge_shape_functions(const std::vector<Point> & qp,
132  const Elem * edge);
133 
138  const std::vector<Point> & get_xyz() const
140  calculate_xyz = true; return xyz; }
141 
145  const std::vector<Real> & get_jacobian() const
147  calculate_dxyz = true; return jac; }
148 
153  const std::vector<Real> & get_JxW() const
155  calculate_dxyz = true; return JxW; }
156 
161  const std::vector<RealGradient> & get_dxyzdxi() const
163  calculate_dxyz = true; return dxyzdxi_map; }
164 
169  const std::vector<RealGradient> & get_dxyzdeta() const
171  calculate_dxyz = true; return dxyzdeta_map; }
172 
177  const std::vector<RealGradient> & get_dxyzdzeta() const
179  calculate_dxyz = true; return dxyzdzeta_map; }
180 
184  const std::vector<RealGradient> & get_d2xyzdxi2() const
186  calculate_d2xyz = true; return d2xyzdxi2_map; }
187 
191  const std::vector<RealGradient> & get_d2xyzdeta2() const
193  calculate_d2xyz = true; return d2xyzdeta2_map; }
194 
195 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
196 
200  const std::vector<RealGradient> & get_d2xyzdzeta2() const
202  calculate_d2xyz = true; return d2xyzdzeta2_map; }
203 
204 #endif
205 
209  const std::vector<RealGradient> & get_d2xyzdxideta() const
211  calculate_d2xyz = true; return d2xyzdxideta_map; }
212 
213 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
214 
218  const std::vector<RealGradient> & get_d2xyzdxidzeta() const
220  calculate_d2xyz = true; return d2xyzdxidzeta_map; }
221 
225  const std::vector<RealGradient> & get_d2xyzdetadzeta() const
227  calculate_d2xyz = true; return d2xyzdetadzeta_map; }
228 
229 #endif
230 
235  const std::vector<Real> & get_dxidx() const
237  calculate_dxyz = true; return dxidx_map; }
238 
243  const std::vector<Real> & get_dxidy() const
245  calculate_dxyz = true; return dxidy_map; }
246 
251  const std::vector<Real> & get_dxidz() const
253  calculate_dxyz = true; return dxidz_map; }
254 
259  const std::vector<Real> & get_detadx() const
261  calculate_dxyz = true; return detadx_map; }
262 
267  const std::vector<Real> & get_detady() const
269  calculate_dxyz = true; return detady_map; }
270 
275  const std::vector<Real> & get_detadz() const
277  calculate_dxyz = true; return detadz_map; }
278 
283  const std::vector<Real> & get_dzetadx() const
285  calculate_dxyz = true; return dzetadx_map; }
286 
291  const std::vector<Real> & get_dzetady() const
293  calculate_dxyz = true; return dzetady_map; }
294 
299  const std::vector<Real> & get_dzetadz() const
301  calculate_dxyz = true; return dzetadz_map; }
302 
303 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
304 
307  const std::vector<std::vector<Real>> & get_d2xidxyz2() const
309  calculate_d2xyz = true; return d2xidxyz2_map; }
310 
314  const std::vector<std::vector<Real>> & get_d2etadxyz2() const
316  calculate_d2xyz = true; return d2etadxyz2_map; }
317 
321  const std::vector<std::vector<Real>> & get_d2zetadxyz2() const
323  calculate_d2xyz = true; return d2zetadxyz2_map; }
324 #endif
325 
329  const std::vector<std::vector<Real>> & get_psi() const
330  { return psi_map; }
331 
335  const std::vector<std::vector<Real>> & get_phi_map() const
337  calculate_xyz = true; return phi_map; }
338 
342  const std::vector<std::vector<Real>> & get_dphidxi_map() const
344  calculate_dxyz = true; return dphidxi_map; }
345 
349  const std::vector<std::vector<Real>> & get_dphideta_map() const
351  calculate_dxyz = true; return dphideta_map; }
352 
356  const std::vector<std::vector<Real>> & get_dphidzeta_map() const
358  calculate_dxyz = true; return dphidzeta_map; }
359 
363  const std::vector<std::vector<Point>> & get_tangents() const
365  calculate_dxyz = true; return tangents; }
366 
370  const std::vector<Point> & get_normals() const
372  calculate_dxyz = true; return normals; }
373 
377  const std::vector<Real> & get_curvatures() const
379  calculate_d2xyz = true; return curvatures;}
380 
384  void print_JxW(std::ostream & os) const;
385 
390  void print_xyz(std::ostream & os) const;
391 
392  /* FIXME: PB: The public functions below break encapsulation! I'm not happy
393  about it, but I don't have time to redo the infinite element routines.
394  These are needed in inf_fe_boundary.C and inf_fe.C. A proper implementation
395  would subclass this class and hide all the radial function business. */
399  std::vector<std::vector<Real>> & get_psi()
400  { return psi_map; }
401 
405  std::vector<std::vector<Real>> & get_dpsidxi()
407  calculate_dxyz = true; return dpsidxi_map; }
408 
412  std::vector<std::vector<Real>> & get_dpsideta()
414  calculate_dxyz = true; return dpsideta_map; }
415 
419  std::vector<std::vector<Real>> & get_d2psidxi2()
421  calculate_d2xyz = true; return d2psidxi2_map; }
422 
426  std::vector<std::vector<Real>> & get_d2psidxideta()
428  calculate_d2xyz = true; return d2psidxideta_map; }
429 
433  std::vector<std::vector<Real>> & get_d2psideta2()
435  calculate_d2xyz = true; return d2psideta2_map; }
436 
440  std::vector<std::vector<Real>> & get_phi_map()
442  calculate_xyz = true; return phi_map; }
443 
447  std::vector<std::vector<Real>> & get_dphidxi_map()
449  calculate_dxyz = true; return dphidxi_map; }
450 
454  std::vector<std::vector<Real>> & get_dphideta_map()
456  calculate_dxyz = true; return dphideta_map; }
457 
461  std::vector<std::vector<Real>> & get_dphidzeta_map()
463  calculate_dxyz = true; return dphidzeta_map; }
464 
465 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
466 
469  std::vector<std::vector<Real>> & get_d2phidxi2_map()
471  calculate_d2xyz = true; return d2phidxi2_map; }
472 
476  std::vector<std::vector<Real>> & get_d2phidxideta_map()
478  calculate_d2xyz = true; return d2phidxideta_map; }
479 
483  std::vector<std::vector<Real>> & get_d2phidxidzeta_map()
485  calculate_d2xyz = true; return d2phidxidzeta_map; }
486 
490  std::vector<std::vector<Real>> & get_d2phideta2_map()
492  calculate_d2xyz = true; return d2phideta2_map; }
493 
497  std::vector<std::vector<Real>> & get_d2phidetadzeta_map()
499  calculate_d2xyz = true; return d2phidetadzeta_map; }
500 
504  std::vector<std::vector<Real>> & get_d2phidzeta2_map()
506  calculate_d2xyz = true; return d2phidzeta2_map; }
507 #endif
508 
509  /* FIXME: PB: This function breaks encapsulation! Needed in FE<>::reinit and
510  InfFE<>::reinit. Not sure yet if the algorithm can be redone to avoid this. */
515  std::vector<Real> & get_JxW()
517  calculate_dxyz = true; return JxW; }
518 
519 protected:
520 
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  }
534 
538  void resize_quadrature_map_vectors(const unsigned int dim, unsigned int n_qp);
539 
546  Real dxdxi_map(const unsigned int p) const { return dxyzdxi_map[p](0); }
547 
554  Real dydxi_map(const unsigned int p) const { return dxyzdxi_map[p](1); }
555 
562  Real dzdxi_map(const unsigned int p) const { return dxyzdxi_map[p](2); }
563 
570  Real dxdeta_map(const unsigned int p) const { return dxyzdeta_map[p](0); }
571 
578  Real dydeta_map(const unsigned int p) const { return dxyzdeta_map[p](1); }
579 
586  Real dzdeta_map(const unsigned int p) const { return dxyzdeta_map[p](2); }
587 
594  Real dxdzeta_map(const unsigned int p) const { return dxyzdzeta_map[p](0); }
595 
602  Real dydzeta_map(const unsigned int p) const { return dxyzdzeta_map[p](1); }
603 
610  Real dzdzeta_map(const unsigned int p) const { return dxyzdzeta_map[p](2); }
611 
615  std::vector<Point> xyz;
616 
621  std::vector<RealGradient> dxyzdxi_map;
622 
627  std::vector<RealGradient> dxyzdeta_map;
628 
633  std::vector<RealGradient> dxyzdzeta_map;
634 
639  std::vector<RealGradient> d2xyzdxi2_map;
640 
645  std::vector<RealGradient> d2xyzdxideta_map;
646 
651  std::vector<RealGradient> d2xyzdeta2_map;
652 
653 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
654 
659  std::vector<RealGradient> d2xyzdxidzeta_map;
660 
665  std::vector<RealGradient> d2xyzdetadzeta_map;
666 
671  std::vector<RealGradient> d2xyzdzeta2_map;
672 
673 #endif
674 
679  std::vector<Real> dxidx_map;
680 
685  std::vector<Real> dxidy_map;
686 
691  std::vector<Real> dxidz_map;
692 
693 
698  std::vector<Real> detadx_map;
699 
704  std::vector<Real> detady_map;
705 
710  std::vector<Real> detadz_map;
711 
712 
717  std::vector<Real> dzetadx_map;
718 
723  std::vector<Real> dzetady_map;
724 
729  std::vector<Real> dzetadz_map;
730 
731 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
732 
736  std::vector<std::vector<Real>> d2xidxyz2_map;
737 
742  std::vector<std::vector<Real>> d2etadxyz2_map;
743 
748  std::vector<std::vector<Real>> d2zetadxyz2_map;
749 #endif
750 
754  std::vector<std::vector<Real>> phi_map;
755 
759  std::vector<std::vector<Real>> dphidxi_map;
760 
764  std::vector<std::vector<Real>> dphideta_map;
765 
769  std::vector<std::vector<Real>> dphidzeta_map;
770 
771 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
772 
776  std::vector<std::vector<Real>> d2phidxi2_map;
777 
781  std::vector<std::vector<Real>> d2phidxideta_map;
782 
786  std::vector<std::vector<Real>> d2phidxidzeta_map;
787 
791  std::vector<std::vector<Real>> d2phideta2_map;
792 
796  std::vector<std::vector<Real>> d2phidetadzeta_map;
797 
801  std::vector<std::vector<Real>> d2phidzeta2_map;
802 
803 #endif
804 
808  std::vector<std::vector<Real>> psi_map;
809 
814  std::vector<std::vector<Real>> dpsidxi_map;
815 
820  std::vector<std::vector<Real>> dpsideta_map;
821 
827  std::vector<std::vector<Real>> d2psidxi2_map;
828 
834  std::vector<std::vector<Real>> d2psidxideta_map;
835 
841  std::vector<std::vector<Real>> d2psideta2_map;
842 
846  std::vector<std::vector<Point>> tangents;
847 
851  std::vector<Point> normals;
852 
858  std::vector<Real> curvatures;
859 
863  std::vector<Real> jac;
864 
868  std::vector<Real> JxW;
869 
874  mutable bool calculations_started;
875 
879  mutable bool calculate_xyz;
880 
884  mutable bool calculate_dxyz;
885 
889  mutable bool calculate_d2xyz;
890 
895  template <unsigned int Dim, FEFamily T>
896  friend class FE;
897 
898 private:
903  void compute_inverse_map_second_derivs(unsigned p);
904 
908  std::vector<const Node *> elem_nodes;
909 };
910 
911 }
912 
913 #endif // LIBMESH_FE_MAP_H
const std::vector< Real > & get_JxW() const
Definition: fe_map.h:153
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
virtual void compute_affine_map(const unsigned int dim, const std::vector< Real > &qw, const Elem *elem)
Compute the jacobian and some other additional data fields.
Definition: fe_map.C:1223
std::vector< std::vector< Real > > & get_dphidzeta_map()
Definition: fe_map.h:461
void init_edge_shape_functions(const std::vector< Point > &qp, const Elem *edge)
Same as before, but for an edge.
Definition: fe_boundary.C:513
const std::vector< std::vector< Real > > & get_phi_map() const
Definition: fe_map.h:335
std::vector< std::vector< Real > > dpsideta_map
Map for the derivative of the side function, d(psi)/d(eta).
Definition: fe_map.h:820
std::vector< std::vector< Real > > psi_map
Map for the side shape functions, psi.
Definition: fe_map.h:808
const std::vector< std::vector< Real > > & get_d2xidxyz2() const
Second derivatives of "xi" reference coordinate wrt physical coordinates.
Definition: fe_map.h:307
std::vector< std::vector< Real > > d2psideta2_map
Map for the second derivatives (in eta) of the side shape functions.
Definition: fe_map.h:841
std::vector< std::vector< Real > > & get_d2phidxideta_map()
Definition: fe_map.h:476
std::vector< std::vector< Real > > & get_d2phideta2_map()
Definition: fe_map.h:490
std::vector< std::vector< Real > > dphidzeta_map
Map for the derivative, d(phi)/d(zeta).
Definition: fe_map.h:769
void compute_inverse_map_second_derivs(unsigned p)
A helper function used by FEMap::compute_single_point_map() to compute second derivatives of the inve...
Definition: fe_map.C:1446
virtual void compute_map(const unsigned int dim, const std::vector< Real > &qw, const Elem *elem, bool calculate_d2phi)
Compute the jacobian and some other additional data fields.
Definition: fe_map.C:1379
Real dxdzeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:594
std::vector< std::vector< Real > > dphidxi_map
Map for the derivative, d(phi)/d(xi).
Definition: fe_map.h:759
std::vector< std::vector< Real > > d2etadxyz2_map
Second derivatives of "eta" reference coordinate wrt physical coordinates.
Definition: fe_map.h:742
const std::vector< std::vector< Real > > & get_d2zetadxyz2() const
Second derivatives of "zeta" reference coordinate wrt physical coordinates.
Definition: fe_map.h:321
bool calculate_dxyz
Should we calculate mapping gradients?
Definition: fe_map.h:884
void compute_edge_map(int dim, const std::vector< Real > &qw, const Elem *side)
Same as before, but for an edge.
Definition: fe_boundary.C:921
std::vector< std::vector< Real > > & get_d2phidxi2_map()
Definition: fe_map.h:469
unsigned int dim
std::vector< RealGradient > d2xyzdzeta2_map
Vector of second partial derivatives in zeta: d^2(x)/d(zeta)^2.
Definition: fe_map.h:671
std::vector< Real > dzetady_map
Map for partial derivatives: d(zeta)/d(y).
Definition: fe_map.h:723
std::vector< std::vector< Real > > d2xidxyz2_map
Second derivatives of "xi" reference coordinate wrt physical coordinates.
Definition: fe_map.h:736
const std::vector< RealGradient > & get_d2xyzdxi2() const
Definition: fe_map.h:184
void init_reference_to_physical_map(const std::vector< Point > &qp, const Elem *elem)
Definition: fe_map.C:68
Real dxdxi_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:546
unsigned short int side
Definition: xdr_io.C:49
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
std::vector< std::vector< Real > > d2phideta2_map
Map for the second derivative, d^2(phi)/d(eta)^2.
Definition: fe_map.h:791
std::vector< std::vector< Real > > d2phidxidzeta_map
Map for the second derivative, d^2(phi)/d(xi)d(zeta).
Definition: fe_map.h:786
const std::vector< Real > & get_dxidz() const
Definition: fe_map.h:251
std::vector< Real > dxidz_map
Map for partial derivatives: d(xi)/d(z).
Definition: fe_map.h:691
std::vector< Real > & get_JxW()
Definition: fe_map.h:515
std::vector< std::vector< Real > > phi_map
Map for the shape function phi.
Definition: fe_map.h:754
The libMesh namespace provides an interface to certain functionality in the library.
Real dydxi_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:554
std::vector< std::vector< Real > > d2phidetadzeta_map
Map for the second derivative, d^2(phi)/d(eta)d(zeta).
Definition: fe_map.h:796
std::vector< std::vector< Real > > d2psidxideta_map
Map for the second (cross) derivatives in xi, eta of the side shape functions.
Definition: fe_map.h:834
std::vector< RealGradient > d2xyzdxideta_map
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(...
Definition: fe_map.h:645
libmesh_assert(j)
void print_JxW(std::ostream &os) const
Prints the Jacobian times the weight for each quadrature point.
Definition: fe_map.C:1430
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const std::vector< Real > & get_dxidy() const
Definition: fe_map.h:243
std::vector< RealGradient > dxyzdzeta_map
Vector of partial derivatives: d(x)/d(zeta), d(y)/d(zeta), d(z)/d(zeta)
Definition: fe_map.h:633
std::vector< std::vector< Real > > & get_dpsidxi()
Definition: fe_map.h:405
const std::vector< RealGradient > & get_dxyzdxi() const
Definition: fe_map.h:161
const std::vector< RealGradient > & get_d2xyzdetadzeta() const
Definition: fe_map.h:225
std::vector< std::vector< Real > > & get_d2phidxidzeta_map()
Definition: fe_map.h:483
std::vector< Real > dzetadx_map
Map for partial derivatives: d(zeta)/d(x).
Definition: fe_map.h:717
std::vector< RealGradient > dxyzdxi_map
Vector of partial derivatives: d(x)/d(xi), d(y)/d(xi), d(z)/d(xi)
Definition: fe_map.h:621
Real dydzeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:602
const std::vector< RealGradient > & get_d2xyzdeta2() const
Definition: fe_map.h:191
A specific instantiation of the FEBase class.
Definition: fe.h:89
std::vector< std::vector< Real > > & get_d2phidetadzeta_map()
Definition: fe_map.h:497
const std::vector< Real > & get_jacobian() const
Definition: fe_map.h:145
std::vector< std::vector< Real > > & get_d2psideta2()
Definition: fe_map.h:433
std::vector< RealGradient > d2xyzdeta2_map
Vector of second partial derivatives in eta: d^2(x)/d(eta)^2.
Definition: fe_map.h:651
const std::vector< std::vector< Real > > & get_d2etadxyz2() const
Second derivatives of "eta" reference coordinate wrt physical coordinates.
Definition: fe_map.h:314
std::vector< std::vector< Real > > & get_d2psidxideta()
Definition: fe_map.h:426
std::vector< RealGradient > d2xyzdxi2_map
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: fe_map.h:639
std::vector< std::vector< Real > > & get_dphidxi_map()
Definition: fe_map.h:447
std::vector< Real > dzetadz_map
Map for partial derivatives: d(zeta)/d(z).
Definition: fe_map.h:729
bool calculate_d2xyz
Should we calculate mapping hessians?
Definition: fe_map.h:889
std::vector< std::vector< Real > > d2zetadxyz2_map
Second derivatives of "zeta" reference coordinate wrt physical coordinates.
Definition: fe_map.h:748
Real dzdxi_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:562
bool calculations_started
Have calculations with this object already been started? Then all get_* functions should already have...
Definition: fe_map.h:874
std::vector< RealGradient > d2xyzdetadzeta_map
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...
Definition: fe_map.h:665
std::vector< Real > dxidx_map
Map for partial derivatives: d(xi)/d(x).
Definition: fe_map.h:679
bool calculate_xyz
Should we calculate physical point locations?
Definition: fe_map.h:879
std::vector< Real > dxidy_map
Map for partial derivatives: d(xi)/d(y).
Definition: fe_map.h:685
const std::vector< RealGradient > & get_dxyzdzeta() const
Definition: fe_map.h:177
Real dxdeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:570
const std::vector< RealGradient > & get_d2xyzdxideta() const
Definition: fe_map.h:209
Real dydeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:578
std::vector< Real > curvatures
The mean curvature (= one half the sum of the principal curvatures) on the boundary at the quadrature...
Definition: fe_map.h:858
const std::vector< std::vector< Real > > & get_dphideta_map() const
Definition: fe_map.h:349
std::vector< std::vector< Real > > d2psidxi2_map
Map for the second derivatives (in xi) of the side shape functions.
Definition: fe_map.h:827
const std::vector< std::vector< Real > > & get_psi() const
Definition: fe_map.h:329
const std::vector< Real > & get_dzetadz() const
Definition: fe_map.h:299
std::vector< std::vector< Real > > & get_d2psidxi2()
Definition: fe_map.h:419
const std::vector< Real > & get_curvatures() const
Definition: fe_map.h:377
std::vector< Real > JxW
Jacobian*Weight values at quadrature points.
Definition: fe_map.h:868
std::vector< std::vector< Real > > d2phidxideta_map
Map for the second derivative, d^2(phi)/d(xi)d(eta).
Definition: fe_map.h:781
std::vector< Point > normals
Normal vectors on boundary at quadrature points.
Definition: fe_map.h:851
const std::vector< std::vector< Real > > & get_dphidzeta_map() const
Definition: fe_map.h:356
std::vector< Point > xyz
The spatial locations of the quadrature points.
Definition: fe_map.h:615
const std::vector< Real > & get_detadz() const
Definition: fe_map.h:275
void print_xyz(std::ostream &os) const
Prints the spatial location of each quadrature point (on the physical element).
Definition: fe_map.C:1438
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::vector< Real > & get_dzetady() const
Definition: fe_map.h:291
std::vector< const Node * > elem_nodes
Work vector for compute_affine_map()
Definition: fe_map.h:908
std::vector< std::vector< Real > > & get_d2phidzeta2_map()
Definition: fe_map.h:504
virtual void compute_null_map(const unsigned int dim, const std::vector< Real > &qw)
Assign a fake jacobian and some other additional data fields.
Definition: fe_map.C:1304
static UniquePtr< FEMap > build(FEType fe_type)
Definition: fe_map.C:50
std::vector< std::vector< Real > > dpsidxi_map
Map for the derivative of the side functions, d(psi)/d(xi).
Definition: fe_map.h:814
std::vector< Real > detady_map
Map for partial derivatives: d(eta)/d(y).
Definition: fe_map.h:704
const std::vector< std::vector< Real > > & get_dphidxi_map() const
Definition: fe_map.h:342
const std::vector< RealGradient > & get_dxyzdeta() const
Definition: fe_map.h:169
std::vector< std::vector< Point > > tangents
Tangent vectors on boundary at quadrature points.
Definition: fe_map.h:846
const std::vector< Real > & get_detady() const
Definition: fe_map.h:267
std::vector< std::vector< Real > > & get_dpsideta()
Definition: fe_map.h:412
std::vector< RealGradient > dxyzdeta_map
Vector of partial derivatives: d(x)/d(eta), d(y)/d(eta), d(z)/d(eta)
Definition: fe_map.h:627
void resize_quadrature_map_vectors(const unsigned int dim, unsigned int n_qp)
A utility function for use by compute_*_map.
Definition: fe_map.C:1142
std::vector< std::vector< Real > > & get_psi()
Definition: fe_map.h:399
const std::vector< Point > & get_normals() const
Definition: fe_map.h:370
void determine_calculations()
Determine which values are to be calculated.
Definition: fe_map.h:524
const std::vector< std::vector< Point > > & get_tangents() const
Definition: fe_map.h:363
const std::vector< Real > & get_dxidx() const
Definition: fe_map.h:235
virtual void compute_face_map(int dim, const std::vector< Real > &qw, const Elem *side)
Same as compute_map, but for a side.
Definition: fe_boundary.C:573
std::vector< Real > jac
Jacobian values at quadrature points.
Definition: fe_map.h:863
std::vector< Real > detadz_map
Map for partial derivatives: d(eta)/d(z).
Definition: fe_map.h:710
std::vector< std::vector< Real > > d2phidzeta2_map
Map for the second derivative, d^2(phi)/d(zeta)^2.
Definition: fe_map.h:801
std::vector< std::vector< Real > > & get_phi_map()
Definition: fe_map.h:440
const std::vector< RealGradient > & get_d2xyzdxidzeta() const
Definition: fe_map.h:218
const std::vector< RealGradient > & get_d2xyzdzeta2() const
Definition: fe_map.h:200
const std::vector< Real > & get_detadx() const
Definition: fe_map.h:259
std::vector< std::vector< Real > > d2phidxi2_map
Map for the second derivative, d^2(phi)/d(xi)^2.
Definition: fe_map.h:776
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)
Compute the jacobian and some other additional data fields at the single point with index p...
Definition: fe_map.C:403
const std::vector< Real > & get_dzetadx() const
Definition: fe_map.h:283
Class contained in FE that encapsulates mapping (i.e.
Definition: fe_map.h:45
const std::vector< Point > & get_xyz() const
Definition: fe_map.h:138
std::vector< Real > detadx_map
Map for partial derivatives: d(eta)/d(x).
Definition: fe_map.h:698
std::vector< std::vector< Real > > dphideta_map
Map for the derivative, d(phi)/d(eta).
Definition: fe_map.h:764
std::vector< std::vector< Real > > & get_dphideta_map()
Definition: fe_map.h:454
void init_face_shape_functions(const std::vector< Point > &qp, const Elem *side)
Initializes the reference to physical element map for a side.
Definition: fe_boundary.C:408
Real dzdeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:586
virtual ~FEMap()
Definition: fe_map.h:50
std::vector< RealGradient > d2xyzdxidzeta_map
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: fe_map.h:659
Real dzdzeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Definition: fe_map.h:610