MFEM  v4.6.0
Finite element discretization library
eltrans.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_ELEMENTTRANSFORM
13 #define MFEM_ELEMENTTRANSFORM
14 
15 #include "../config/config.hpp"
16 #include "../linalg/linalg.hpp"
17 #include "intrules.hpp"
18 #include "fe.hpp"
19 
20 namespace mfem
21 {
22 
24 {
25 protected:
29  double Wght;
30  int EvalState;
32  {
39  };
41 
42  /** @brief Evaluate the Jacobian of the transformation at the IntPoint and
43  store it in dFdx. */
44  virtual const DenseMatrix &EvalJacobian() = 0;
45 
46  /** @brief Evaluate the Hessian of the transformation at the IntPoint and
47  store it in d2Fdx2. */
48  virtual const DenseMatrix &EvalHessian() = 0;
49 
50  double EvalWeight();
51  const DenseMatrix &EvalAdjugateJ();
53  const DenseMatrix &EvalInverseJ();
54 
55 public:
56 
57  /** This enumeration declares the values stored in
58  ElementTransformation::ElementType and indicates which group of objects
59  the index stored in ElementTransformation::ElementNo refers:
60 
61  | ElementType | Range of ElementNo
62  +-------------+-------------------------
63  | ELEMENT | [0, Mesh::GetNE() )
64  | BDR_ELEMENT | [0, Mesh::GetNBE() )
65  | EDGE | [0, Mesh::GetNEdges() )
66  | FACE | [0, Mesh::GetNFaces() )
67  | BDR_FACE | [0, Mesh::GetNBE() )
68  */
69  enum
70  {
71  ELEMENT = 1,
73  EDGE = 3,
74  FACE = 4,
76  };
77 
79 
80  /// The Mesh object containing the element.
81  /** If the element transformation belongs to a mesh, this will point to the
82  containing Mesh object. ElementNo will be the number of the element in
83  this Mesh. This will be NULL if the element does not belong to a mesh. */
84  class Mesh *mesh;
85 
87 
88  /** @brief Force the reevaluation of the Jacobian in the next call. */
89  void Reset() { EvalState = 0; }
90 
91  /** @brief Set the integration point @a ip that weights and Jacobians will
92  be evaluated at. */
93  void SetIntPoint(const IntegrationPoint *ip)
94  { IntPoint = ip; EvalState = 0; }
95 
96  /** @brief Get a const reference to the currently set integration point. This
97  will return NULL if no integration point is set. */
98  const IntegrationPoint &GetIntPoint() { return *IntPoint; }
99 
100  /** @brief Transform integration point from reference coordinates to
101  physical coordinates and store them in the vector. */
102  virtual void Transform(const IntegrationPoint &, Vector &) = 0;
103 
104  /** @brief Transform all the integration points from the integration rule
105  from reference coordinates to physical
106  coordinates and store them as column vectors in the matrix. */
107  virtual void Transform(const IntegrationRule &, DenseMatrix &) = 0;
108 
109  /** @brief Transform all the integration points from the column vectors
110  of @a matrix from reference coordinates to physical
111  coordinates and store them as column vectors in @a result. */
112  virtual void Transform(const DenseMatrix &matrix, DenseMatrix &result) = 0;
113 
114  /** @brief Return the Jacobian matrix of the transformation at the currently
115  set IntegrationPoint, using the method SetIntPoint(). */
116  /** The dimensions of the Jacobian matrix are physical-space-dim by
117  reference-space-dim. The first column contains the x derivatives of the
118  transformation, the second -- the y derivatives, etc. */
120  { return (EvalState & JACOBIAN_MASK) ? dFdx : EvalJacobian(); }
121 
122 
123  /** @brief Return the Hessian matrix of the transformation at the currently
124  set IntegrationPoint, using the method SetIntPoint(). */
126  { return (EvalState & HESSIAN_MASK) ? d2Fdx2 : EvalHessian(); }
127 
128  /** @brief Return the weight of the Jacobian matrix of the transformation
129  at the currently set IntegrationPoint.
130  The Weight evaluates to \f$ \sqrt{\lvert J^T J \rvert} \f$. */
131  double Weight() { return (EvalState & WEIGHT_MASK) ? Wght : EvalWeight(); }
132 
133  /** @brief Return the adjugate of the Jacobian matrix of the transformation
134  at the currently set IntegrationPoint. */
136  { return (EvalState & ADJUGATE_MASK) ? adjJ : EvalAdjugateJ(); }
137 
138  /** @brief Return the transpose of the adjugate of the Jacobian matrix of
139  the transformation at the currently set IntegrationPoint. */
142 
143  /** @brief Return the inverse of the Jacobian matrix of the transformation
144  at the currently set IntegrationPoint. */
146  { return (EvalState & INVERSE_MASK) ? invJ : EvalInverseJ(); }
147 
148  /// Return the order of the current element we are using for the transformation.
149  virtual int Order() const = 0;
150 
151  /// Return the order of the elements of the Jacobian of the transformation.
152  virtual int OrderJ() const = 0;
153 
154  /** @brief Return the order of the determinant of the Jacobian (weight)
155  of the transformation. */
156  virtual int OrderW() const = 0;
157 
158  /// Return the order of \f$ adj(J)^T \nabla fi \f$
159  virtual int OrderGrad(const FiniteElement *fe) const = 0;
160 
161  /// Return the Geometry::Type of the reference element.
162  Geometry::Type GetGeometryType() const { return geom; }
163 
164  /// Return the topological dimension of the reference element.
165  int GetDimension() const { return Geometry::Dimension[geom]; }
166 
167  /// Get the dimension of the target (physical) space.
168  /** We support 2D meshes embedded in 3D; in this case the function will
169  return "3". */
170  virtual int GetSpaceDim() const = 0;
171 
172  /** @brief Transform a point @a pt from physical space to a point @a ip in
173  reference space. */
174  /** Attempt to find the IntegrationPoint that is transformed into the given
175  point in physical space. If the inversion fails a non-zero value is
176  returned. This method is not 100 percent reliable for non-linear
177  transformations. */
178  virtual int TransformBack(const Vector &pt, IntegrationPoint &ip) = 0;
179 
181 };
182 
183 
184 /// The inverse transformation of a given ElementTransformation.
186 {
187 public:
188  /// Algorithms for selecting an initial guess.
190  {
191  Center = 0, ///< Use the center of the reference element.
192  ClosestPhysNode = 1, /**<
193  Use the point returned by FindClosestPhysPoint() from a reference-space
194  grid of type and size controlled by SetInitGuessPointsType() and
195  SetInitGuessRelOrder(), respectively. */
196  ClosestRefNode = 2, /**<
197  Use the point returned by FindClosestRefPoint() from a reference-space
198  grid of type and size controlled by SetInitGuessPointsType() and
199  SetInitGuessRelOrder(), respectively. */
200  GivenPoint = 3 ///< Use a specific point, set with SetInitialGuess().
201  };
202 
203  /// Solution strategy.
205  {
206  Newton = 0, /**<
207  Use Newton's algorithm, without restricting the reference-space points
208  (iterates) to the reference element. */
210  Use Newton's algorithm, restricting the reference-space points to the
211  reference element by scaling back the Newton increments, i.e.
212  projecting new iterates, x_new, lying outside the element, to the
213  intersection of the line segment [x_old, x_new] with the boundary. */
215  Use Newton's algorithm, restricting the reference-space points to the
216  reference element by projecting new iterates, x_new, lying outside the
217  element, to the point on the boundary closest (in reference-space) to
218  x_new. */
219  };
220 
221  /// Values returned by Transform().
223  {
224  Inside = 0, ///< The point is inside the element
225  Outside = 1, ///< The point is _probably_ outside the element
226  Unknown = 2 ///< The algorithm failed to determine where the point is
227  };
228 
229 protected:
230  // Pointer to the forward transformation. Not owned.
232 
233  // Parameters of the inversion algorithms:
235  int init_guess_type; // algorithm to use
236  int qpts_type; // Quadrature1D type for the initial guess type
237  int rel_qpts_order; // num_1D_qpts = max(trans_order+rel_qpts_order,0)+1
238  int solver_type; // solution strategy to use
239  int max_iter; // max. number of Newton iterations
240  double ref_tol; // reference space tolerance
241  double phys_rtol; // physical space tolerance (relative)
242  double ip_tol; // tolerance for checking if a point is inside the ref. elem.
244 
245  void NewtonPrint(int mode, double val);
246  void NewtonPrintPoint(const char *prefix, const Vector &pt,
247  const char *suffix);
248  int NewtonSolve(const Vector &pt, IntegrationPoint &ip);
249 
250 public:
251  /// Construct the InverseElementTransformation with default parameters.
252  /** Some practical considerations regarding the choice of initial guess type
253  and solver type:
254  1. The combination of #Center and #NewtonSegmentProject provides the
255  fastest way to estimate if a point lies inside an element, assuming
256  that most queried elements are not very deformed, e.g. if most
257  elements are convex.
258  2. [Default] The combination of #Center and #NewtonElementProject
259  provides a somewhat slower alternative to 1 with the benefit of being
260  more reliable in the case when the query point is inside the element
261  but potentially slower in the case when the query point is outside the
262  element.
263  3. The combination of #ClosestPhysNode and #NewtonElementProject is
264  slower than 1 and 2 but much more reliable, especially in the case of
265  highly distorted elements which do not have very high aspect ratios.
266  4. The combination of #ClosestRefNode and #NewtonElementProject should
267  generally be the most reliable, coming at a bit higher computational
268  cost than 3 while performing better on distorted meshes with elements
269  having high aspect ratios.
270  @note None of these choices provide a guarantee that if a query point is
271  inside the element then it will be found. The only guarantee is that if
272  the Transform() method returns #Inside then the point lies inside the
273  element up to one of the specified physical- or reference-space
274  tolerances. */
276  : T(Trans),
277  ip0(NULL),
279  qpts_type(Quadrature1D::OpenHalfUniform),
280  rel_qpts_order(-1),
282  max_iter(16),
283  ref_tol(1e-15),
284  phys_rtol(1e-15),
285  ip_tol(1e-8),
286  print_level(-1)
287  { }
288 
290 
291  /// Set a new forward ElementTransformation, @a Trans.
292  void SetTransformation(ElementTransformation &Trans) { T = &Trans; }
293 
294  /** @brief Choose how the initial guesses for subsequent calls to Transform()
295  will be selected. */
297 
298  /** @brief Set the initial guess for subsequent calls to Transform(),
299  switching to the #GivenPoint #InitGuessType at the same time. */
300  void SetInitialGuess(const IntegrationPoint &init_ip)
301  { ip0 = &init_ip; SetInitialGuessType(GivenPoint); }
302 
303  /// Set the Quadrature1D type used for the `Closest*` initial guess types.
304  void SetInitGuessPointsType(int q_type) { qpts_type = q_type; }
305 
306  /// Set the relative order used for the `Closest*` initial guess types.
307  /** The number of points in each spatial direction is given by the formula
308  max(trans_order+order,0)+1, where trans_order is the order of the current
309  ElementTransformation. */
310  void SetInitGuessRelOrder(int order) { rel_qpts_order = order; }
311 
312  /** @brief Specify which algorithm to use for solving the transformation
313  equation, i.e. when calling the Transform() method. */
314  void SetSolverType(SolverType stype) { solver_type = stype; }
315 
316  /// Set the maximum number of iterations when solving for a reference point.
317  void SetMaxIter(int max_it) { max_iter = max_it; }
318 
319  /// Set the reference-space convergence tolerance.
320  void SetReferenceTol(double ref_sp_tol) { ref_tol = ref_sp_tol; }
321 
322  /// Set the relative physical-space convergence tolerance.
323  void SetPhysicalRelTol(double phys_rel_tol) { phys_rtol = phys_rel_tol; }
324 
325  /** @brief Set the tolerance used to determine if a point lies inside or
326  outside of the reference element. */
327  /** This tolerance is used only with the pure #Newton solver. */
328  void SetElementTol(double el_tol) { ip_tol = el_tol; }
329 
330  /// Set the desired print level, useful for debugging.
331  /** The valid options are: -1 - never print (default); 0 - print only errors;
332  1 - print the first and last iterations; 2 - print every iteration;
333  and 3 - print every iteration including point coordinates. */
334  void SetPrintLevel(int pr_level) { print_level = pr_level; }
335 
336  /** @brief Find the IntegrationPoint mapped closest to @a pt. */
337  /** This function uses the given IntegrationRule, @a ir, maps its points to
338  physical space and finds the one that is closest to the point @a pt.
339 
340  @param pt The query point.
341  @param ir The IntegrationRule, i.e. the set of reference points to map
342  to physical space and check.
343  @return The index of the IntegrationPoint in @a ir whose mapped point is
344  closest to @a pt.
345  @see FindClosestRefPoint(). */
346  int FindClosestPhysPoint(const Vector& pt, const IntegrationRule &ir);
347 
348  /** @brief Find the IntegrationPoint mapped closest to @a pt, using a norm
349  that approximates the (unknown) distance in reference coordinates. */
350  /** @see FindClosestPhysPoint(). */
351  int FindClosestRefPoint(const Vector& pt, const IntegrationRule &ir);
352 
353  /** @brief Given a point, @a pt, in physical space, find its reference
354  coordinates, @a ip.
355 
356  @returns A value of type #TransformResult. */
357  virtual int Transform(const Vector &pt, IntegrationPoint &ip);
358 };
359 
360 /// A standard isoparametric element transformation
362 {
363 private:
364  DenseMatrix dshape,d2shape;
365  Vector shape;
366 
367  const FiniteElement *FElem;
368  DenseMatrix PointMat; // dim x dof
369 
370  /** @brief Evaluate the Jacobian of the transformation at the IntPoint and
371  store it in dFdx. */
372  virtual const DenseMatrix &EvalJacobian();
373  // Evaluate the Hessian of the transformation at the IntPoint and store it
374  // in d2Fdx2.
375  virtual const DenseMatrix &EvalHessian();
376 
377 public:
378  IsoparametricTransformation() : FElem(NULL) {}
379 
380  /// Set the element that will be used to compute the transformations
381  void SetFE(const FiniteElement *FE)
382  {
383  MFEM_ASSERT(FE != NULL, "Must provide a valid FiniteElement object!");
384  EvalState = (FE != FElem) ? 0 : EvalState;
385  FElem = FE; geom = FE->GetGeomType();
386  }
387 
388  /// Get the current element used to compute the transformations
389  const FiniteElement* GetFE() const { return FElem; }
390 
391  /// @brief Set the underlying point matrix describing the transformation.
392  /** The dimensions of the matrix are space-dim x dof. The transformation is
393  defined as
394  \f$ x = F( \hat x ) = P \phi( \hat x ) \f$
395 
396  where \f$ \hat x \f$ is the reference point, @a x is the corresponding
397  physical point, @a P is the point matrix, and \f$ \phi( \hat x ) \f$ is
398  the column-vector of all basis functions evaluated at \f$ \hat x \f$ .
399  The columns of @a P represent the control points in physical space
400  defining the transformation. */
401  void SetPointMat(const DenseMatrix &pm) { PointMat = pm; EvalState = 0; }
402 
403  /// Return the stored point matrix.
404  const DenseMatrix &GetPointMat() const { return PointMat; }
405 
406  /// @brief Write access to the stored point matrix. Use with caution.
407  /** If the point matrix is altered using this member function the Reset
408  function should also be called to force the reevaluation of the
409  Jacobian, etc.. */
410  DenseMatrix &GetPointMat() { return PointMat; }
411 
412  /// Set the FiniteElement Geometry for the reference elements being used.
414 
415  /** @brief Transform integration point from reference coordinates to
416  physical coordinates and store them in the vector. */
417  virtual void Transform(const IntegrationPoint &, Vector &);
418 
419  /** @brief Transform all the integration points from the integration rule
420  from reference coordinates to physical
421  coordinates and store them as column vectors in the matrix. */
422  virtual void Transform(const IntegrationRule &, DenseMatrix &);
423 
424  /** @brief Transform all the integration points from the column vectors
425  of @a matrix from reference coordinates to physical
426  coordinates and store them as column vectors in @a result. */
427  virtual void Transform(const DenseMatrix &matrix, DenseMatrix &result);
428 
429  /// Return the order of the current element we are using for the transformation.
430  virtual int Order() const { return FElem->GetOrder(); }
431 
432  /// Return the order of the elements of the Jacobian of the transformation.
433  virtual int OrderJ() const;
434 
435  /** @brief Return the order of the determinant of the Jacobian (weight)
436  of the transformation. */
437  virtual int OrderW() const;
438 
439  /// Return the order of \f$ adj(J)^T \nabla fi \f$
440  virtual int OrderGrad(const FiniteElement *fe) const;
441 
442  virtual int GetSpaceDim() const { return PointMat.Height(); }
443 
444  /** @brief Transform a point @a pt from physical space to a point @a ip in
445  reference space. */
446  /** Attempt to find the IntegrationPoint that is transformed into the given
447  point in physical space. If the inversion fails a non-zero value is
448  returned. This method is not 100 percent reliable for non-linear
449  transformations. */
450  virtual int TransformBack(const Vector & v, IntegrationPoint & ip)
451  {
452  InverseElementTransformation inv_tr(this);
453  return inv_tr.Transform(v, ip);
454  }
455 
457 
458  MFEM_DEPRECATED void FinalizeTransformation() {}
459 };
460 
462 {
463 public:
465  void Transform (const IntegrationPoint &, IntegrationPoint &);
466  void Transform (const IntegrationRule &, IntegrationRule &);
467 };
468 
469 /** @brief A specialized ElementTransformation class representing a face and
470  its two neighboring elements.
471 
472  This class can be used as a container for the element transformation data
473  needed for integrating discontinuous fields on element interfaces in a
474  Discontinuous Galerkin (DG) context.
475 
476  The secondary purpose of this class is to enable the
477  GridFunction::GetValue function, and various related functions, to properly
478  evaluate fields with limited continuity on boundary elements.
479 */
481 {
482 private:
483 
484  // Bitwise OR of ConfigMasks
485  int mask;
486 
487  IntegrationPoint eip1, eip2;
488 
489 protected: // interface for Mesh to be able to configure this object.
490 
491  friend class Mesh;
492 #ifdef MFEM_USE_MPI
493  friend class ParMesh;
494 #endif
495 
496  /// Set the mask indicating which portions of the object have been setup
497  /** The argument @a m is a bitmask used in
498  Mesh::GetFaceElementTransformations to indicate which portions of the
499  FaceElementTransformations object have been configured.
500 
501  mask & 1: Elem1 is configured
502  mask & 2: Elem2 is configured
503  mask & 4: Loc1 is configured
504  mask & 8: Loc2 is configured
505  mask & 16: The Face transformation itself is configured
506  */
507  void SetConfigurationMask(int m) { mask = m; }
508 
509 public:
510 
512  {
513  HAVE_ELEM1 = 1, ///< Element on side 1 is configured
514  HAVE_ELEM2 = 2, ///< Element on side 2 is configured
515  HAVE_LOC1 = 4, ///< Point transformation for side 1 is configured
516  HAVE_LOC2 = 8, ///< Point transformation for side 2 is configured
517  HAVE_FACE = 16 ///< Face transformation is configured
518  };
519 
521  Geometry::Type &FaceGeom; ///< @deprecated Use GetGeometryType instead
523  ElementTransformation *Face; ///< @deprecated No longer necessary
525 
527 
528  /** @brief Method to set the geometry type of the face.
529 
530  @note This method is designed to be used when
531  [Par]Mesh::GetFaceTransformation will not be called i.e. when the face
532  transformation will not be needed but the neighboring element
533  transformations will be. Using this method to override the GeometryType
534  should only be done with great care.
535  */
537 
538  /** @brief Return the mask defining the configuration state.
539 
540  The mask value indicates which portions of FaceElementTransformations
541  object have been configured.
542 
543  mask & 1: Elem1 is configured
544  mask & 2: Elem2 is configured
545  mask & 4: Loc1 is configured
546  mask & 8: Loc2 is configured
547  mask & 16: The Face transformation itself is configured
548  */
549  int GetConfigurationMask() const { return mask; }
550 
551  /** @brief Set the integration point in the Face and the two neighboring
552  elements, if present.
553 
554  The point @a face_ip must be in the reference coordinate system of the
555  face.
556  */
557  void SetIntPoint(const IntegrationPoint *face_ip);
558 
559  /** @brief Set the integration point in the Face and the two neighboring
560  elements, if present.
561 
562  This is a more expressive member function name than SetIntPoint, which
563  in this special case, does the same thing. This function can be used for
564  greater code clarity.
565  */
566  inline void SetAllIntPoints(const IntegrationPoint *face_ip)
568 
569  /** @brief Get a const reference to the integration point in neighboring
570  element 1 corresponding to the currently set integration point on the
571  face.
572 
573  This IntegrationPoint object will only contain up-to-date data if
574  SetIntPoint or SetAllIntPoints has been called with the latest
575  integration point for the face and the appropriate point transformation
576  has been configured. */
577  const IntegrationPoint &GetElement1IntPoint() { return eip1; }
578 
579  /** @brief Get a const reference to the integration point in neighboring
580  element 2 corresponding to the currently set integration point on the
581  face.
582 
583  This IntegrationPoint object will only contain up-to-date data if
584  SetIntPoint or SetAllIntPoints has been called with the latest
585  integration point for the face and the appropriate point transformation
586  has been configured. */
587  const IntegrationPoint &GetElement2IntPoint() { return eip2; }
588 
589  virtual void Transform(const IntegrationPoint &, Vector &);
590  virtual void Transform(const IntegrationRule &, DenseMatrix &);
591  virtual void Transform(const DenseMatrix &matrix, DenseMatrix &result);
592 
597 
598  /** @brief Check for self-consistency: compares the result of mapping the
599  reference face vertices to physical coordinates using the three
600  transformations: face, element 1, and element 2.
601 
602  @param[in] print_level If set to a positive number, print the physical
603  coordinates of the face vertices computed through
604  all available transformations: face, element 1,
605  and/or element 2.
606  @param[in,out] out The output stream to use for printing.
607 
608  @returns A maximal distance between physical coordinates of face vertices
609  that should coincide. A successful check should return a small
610  number relative to the mesh extents. If less than 2 of the three
611  transformations are set, returns 0.
612 
613  @warning This check will generally fail on periodic boundary faces.
614  */
615  double CheckConsistency(int print_level = 0,
616  std::ostream &out = mfem::out);
617 };
618 
619 /** Elem1(Loc1(x)) = Face(x) = Elem2(Loc2(x))
620 
621 
622  Physical Space
623 
624  *--------* ^ *--------*
625  Elem1No / / \ / \ / \ \ Elem2No
626  / / \ / \ / \ \
627  / / n \ / \ / \ \
628  *--------* ==> * ( ) * *--------*
629  \ \ / \ / \ / /
630  \ \ / \ / \ / /
631  \ \ / \ / \ / /
632  *--------* v *--------*
633 
634  ^ ^
635  | ^ |
636  Elem1 | | | Elem2
637  | | Face |
638  |
639  *--------* *--------*
640  / /| / /|
641  1 *--------* | 1 *--------* 1 *--------* |
642  | | | Loc1 | | Loc2 | | |
643  | | * <----- | x | -----> | | *
644  | |/ | | | |/
645  *--------* *--------* *--------*
646  0 1 0 1 0 1
647 
648  Reference Space
649 */
650 
651 }
652 
653 #endif
Abstract class for all finite elements.
Definition: fe_base.hpp:233
virtual ~ElementTransformation()
Definition: eltrans.hpp:180
IntegrationPointTransformation & GetIntPoint1Transformation()
Definition: eltrans.cpp:607
void SetPrintLevel(int pr_level)
Set the desired print level, useful for debugging.
Definition: eltrans.hpp:334
int NewtonSolve(const Vector &pt, IntegrationPoint &ip)
Definition: eltrans.cpp:172
const DenseMatrix & AdjugateJacobian()
Return the adjugate of the Jacobian matrix of the transformation at the currently set IntegrationPoin...
Definition: eltrans.hpp:135
ElementTransformation & GetElement2Transformation()
Definition: eltrans.cpp:599
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:96
virtual int OrderGrad(const FiniteElement *fe) const =0
Return the order of .
DenseMatrix & GetPointMat()
Write access to the stored point matrix. Use with caution.
Definition: eltrans.hpp:410
void SetIntPoint(const IntegrationPoint *face_ip)
Set the integration point in the Face and the two neighboring elements, if present.
Definition: eltrans.cpp:568
Element on side 1 is configured.
Definition: eltrans.hpp:513
Use a specific point, set with SetInitialGuess().
Definition: eltrans.hpp:200
virtual int OrderW() const =0
Return the order of the determinant of the Jacobian (weight) of the transformation.
void SetPhysicalRelTol(double phys_rel_tol)
Set the relative physical-space convergence tolerance.
Definition: eltrans.hpp:323
void SetInitialGuess(const IntegrationPoint &init_ip)
Set the initial guess for subsequent calls to Transform(), switching to the GivenPoint InitGuessType ...
Definition: eltrans.hpp:300
const IntegrationPoint * IntPoint
Definition: eltrans.hpp:26
ElementTransformation * T
Definition: eltrans.hpp:231
virtual void Transform(const IntegrationPoint &, Vector &)
Transform integration point from reference coordinates to physical coordinates and store them in the ...
Definition: eltrans.cpp:622
Use the center of the reference element.
Definition: eltrans.hpp:191
virtual int Transform(const Vector &pt, IntegrationPoint &ip)
Given a point, pt, in physical space, find its reference coordinates, ip.
Definition: eltrans.cpp:336
virtual int TransformBack(const Vector &v, IntegrationPoint &ip)
Transform a point pt from physical space to a point ip in reference space.
Definition: eltrans.hpp:450
IntegrationPointTransformation & GetIntPoint2Transformation()
Definition: eltrans.cpp:615
void SetIntPoint(const IntegrationPoint *ip)
Set the integration point ip that weights and Jacobians will be evaluated at.
Definition: eltrans.hpp:93
A specialized ElementTransformation class representing a face and its two neighboring elements...
Definition: eltrans.hpp:480
void NewtonPrint(int mode, double val)
Definition: eltrans.cpp:129
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
void SetInitialGuessType(InitGuessType itype)
Choose how the initial guesses for subsequent calls to Transform() will be selected.
Definition: eltrans.hpp:296
Point transformation for side 2 is configured.
Definition: eltrans.hpp:516
IntegrationPointTransformation Loc2
Definition: eltrans.hpp:524
const DenseMatrix & EvalAdjugateJ()
Definition: eltrans.cpp:36
const DenseMatrix & InverseJacobian()
Return the inverse of the Jacobian matrix of the transformation at the currently set IntegrationPoint...
Definition: eltrans.hpp:145
const IntegrationPoint & GetElement2IntPoint()
Get a const reference to the integration point in neighboring element 2 corresponding to the currentl...
Definition: eltrans.hpp:587
A class container for 1D quadrature type constants.
Definition: intrules.hpp:390
virtual void Transform(const IntegrationPoint &, Vector &)
Transform integration point from reference coordinates to physical coordinates and store them in the ...
Definition: eltrans.cpp:492
const IntegrationPoint & GetIntPoint()
Get a const reference to the currently set integration point. This will return NULL if no integration...
Definition: eltrans.hpp:98
Element on side 2 is configured.
Definition: eltrans.hpp:514
The inverse transformation of a given ElementTransformation.
Definition: eltrans.hpp:185
ElementTransformation * Elem2
Definition: eltrans.hpp:522
void SetSolverType(SolverType stype)
Specify which algorithm to use for solving the transformation equation, i.e. when calling the Transfo...
Definition: eltrans.hpp:314
virtual int OrderJ() const =0
Return the order of the elements of the Jacobian of the transformation.
void NewtonPrintPoint(const char *prefix, const Vector &pt, const char *suffix)
Definition: eltrans.cpp:158
int GetConfigurationMask() const
Return the mask defining the configuration state.
Definition: eltrans.hpp:549
void SetTransformation(ElementTransformation &Trans)
Set a new forward ElementTransformation, Trans.
Definition: eltrans.hpp:292
Geometry::Type GetGeomType() const
Returns the Geometry::Type of the reference element.
Definition: fe_base.hpp:320
int FindClosestRefPoint(const Vector &pt, const IntegrationRule &ir)
Find the IntegrationPoint mapped closest to pt, using a norm that approximates the (unknown) distance...
Definition: eltrans.cpp:97
void SetPointMat(const DenseMatrix &pm)
Set the underlying point matrix describing the transformation.
Definition: eltrans.hpp:401
The point is inside the element.
Definition: eltrans.hpp:224
void SetGeometryType(Geometry::Type g)
Method to set the geometry type of the face.
Definition: eltrans.hpp:536
static const int Dimension[NumGeom]
Definition: geom.hpp:47
Geometry::Type GetGeometryType() const
Return the Geometry::Type of the reference element.
Definition: eltrans.hpp:162
IntegrationPointTransformation Loc1
Definition: eltrans.hpp:524
The point is probably outside the element.
Definition: eltrans.hpp:225
const DenseMatrix & TransposeAdjugateJacobian()
Return the transpose of the adjugate of the Jacobian matrix of the transformation at the currently se...
Definition: eltrans.hpp:140
virtual int Order() const
Return the order of the current element we are using for the transformation.
Definition: eltrans.hpp:430
void SetAllIntPoints(const IntegrationPoint *face_ip)
Set the integration point in the Face and the two neighboring elements, if present.
Definition: eltrans.hpp:566
const DenseMatrix & Jacobian()
Return the Jacobian matrix of the transformation at the currently set IntegrationPoint, using the method SetIntPoint().
Definition: eltrans.hpp:119
void Reset()
Force the reevaluation of the Jacobian in the next call.
Definition: eltrans.hpp:89
void SetReferenceTol(double ref_sp_tol)
Set the reference-space convergence tolerance.
Definition: eltrans.hpp:320
void SetInitGuessPointsType(int q_type)
Set the Quadrature1D type used for the Closest* initial guess types.
Definition: eltrans.hpp:304
InverseElementTransformation(ElementTransformation *Trans=NULL)
Construct the InverseElementTransformation with default parameters.
Definition: eltrans.hpp:275
double Weight()
Return the weight of the Jacobian matrix of the transformation at the currently set IntegrationPoint...
Definition: eltrans.hpp:131
virtual const DenseMatrix & EvalHessian()=0
Evaluate the Hessian of the transformation at the IntPoint and store it in d2Fdx2.
void SetElementTol(double el_tol)
Set the tolerance used to determine if a point lies inside or outside of the reference element...
Definition: eltrans.hpp:328
int GetDimension() const
Return the topological dimension of the reference element.
Definition: eltrans.hpp:165
Point transformation for side 1 is configured.
Definition: eltrans.hpp:515
const DenseMatrix & EvalInverseJ()
Definition: eltrans.cpp:57
void SetConfigurationMask(int m)
Set the mask indicating which portions of the object have been setup.
Definition: eltrans.hpp:507
virtual int OrderW() const
Return the order of the determinant of the Jacobian (weight) of the transformation.
Definition: eltrans.cpp:457
void SetIdentityTransformation(Geometry::Type GeomType)
Set the FiniteElement Geometry for the reference elements being used.
Definition: eltrans.cpp:383
double CheckConsistency(int print_level=0, std::ostream &out=mfem::out)
Check for self-consistency: compares the result of mapping the reference face vertices to physical co...
Definition: eltrans.cpp:646
virtual int Order() const =0
Return the order of the current element we are using for the transformation.
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
void Transform(const IntegrationPoint &, IntegrationPoint &)
Definition: eltrans.cpp:546
const IntegrationPoint & GetElement1IntPoint()
Get a const reference to the integration point in neighboring element 1 corresponding to the currentl...
Definition: eltrans.hpp:577
int Height() const
Get the height (size of output) of the Operator. Synonym with NumRows().
Definition: operator.hpp:66
virtual int OrderGrad(const FiniteElement *fe) const
Return the order of .
Definition: eltrans.cpp:471
MFEM_DEPRECATED void FinalizeTransformation()
Definition: eltrans.hpp:458
virtual int GetSpaceDim() const
Get the dimension of the target (physical) space.
Definition: eltrans.hpp:442
virtual int GetSpaceDim() const =0
Get the dimension of the target (physical) space.
Face transformation is configured.
Definition: eltrans.hpp:517
Class for integration point with weight.
Definition: intrules.hpp:31
The algorithm failed to determine where the point is.
Definition: eltrans.hpp:226
SolverType
Solution strategy.
Definition: eltrans.hpp:204
A standard isoparametric element transformation.
Definition: eltrans.hpp:361
TransformResult
Values returned by Transform().
Definition: eltrans.hpp:222
int FindClosestPhysPoint(const Vector &pt, const IntegrationRule &ir)
Find the IntegrationPoint mapped closest to pt.
Definition: eltrans.cpp:70
void SetFE(const FiniteElement *FE)
Set the element that will be used to compute the transformations.
Definition: eltrans.hpp:381
ElementTransformation & GetElement1Transformation()
Definition: eltrans.cpp:591
IsoparametricTransformation Transf
Definition: eltrans.hpp:464
ElementTransformation * Elem1
Definition: eltrans.hpp:522
ElementTransformation * Face
Definition: eltrans.hpp:523
virtual const DenseMatrix & EvalJacobian()=0
Evaluate the Jacobian of the transformation at the IntPoint and store it in dFdx. ...
InitGuessType
Algorithms for selecting an initial guess.
Definition: eltrans.hpp:189
virtual int OrderJ() const
Return the order of the elements of the Jacobian of the transformation.
Definition: eltrans.cpp:443
Vector data type.
Definition: vector.hpp:58
const DenseMatrix & EvalTransAdjugateJ()
Definition: eltrans.cpp:46
virtual void Transform(const IntegrationPoint &, Vector &)=0
Transform integration point from reference coordinates to physical coordinates and store them in the ...
const DenseMatrix & GetPointMat() const
Return the stored point matrix.
Definition: eltrans.hpp:404
const IntegrationPoint * ip0
Definition: eltrans.hpp:234
const DenseMatrix & Hessian()
Return the Hessian matrix of the transformation at the currently set IntegrationPoint, using the method SetIntPoint().
Definition: eltrans.hpp:125
const FiniteElement * GetFE() const
Get the current element used to compute the transformations.
Definition: eltrans.hpp:389
void SetMaxIter(int max_it)
Set the maximum number of iterations when solving for a reference point.
Definition: eltrans.hpp:317
Class for parallel meshes.
Definition: pmesh.hpp:32
void SetInitGuessRelOrder(int order)
Set the relative order used for the Closest* initial guess types.
Definition: eltrans.hpp:310
int GetOrder() const
Returns the order of the finite element. In the case of anisotropic orders, returns the maximum order...
Definition: fe_base.hpp:327
class Mesh * mesh
The Mesh object containing the element.
Definition: eltrans.hpp:84
virtual int TransformBack(const Vector &pt, IntegrationPoint &ip)=0
Transform a point pt from physical space to a point ip in reference space.