libMesh::InfFEMap Class Reference

Class that encapsulates mapping (i.e. More...

#include <inf_fe_map.h>

Static Public Member Functions
static Point	map (const unsigned int dim, const Elem *inf_elem, const Point &reference_point)
static Point	inverse_map (const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
static void	inverse_map (const unsigned int dim, const Elem *elem, const std::vector< Point > &physical_points, std::vector< Point > &reference_points, const Real tolerance=TOLERANCE, const bool secure=true)
	Takes a number points in physical space (in the physical_points vector) and finds their location on the reference element for the input element elem. More...
static Real	eval (Real v, Order o, unsigned int i)
static Real	eval_deriv (Real v, Order o, unsigned int i)

Detailed Description

Class that encapsulates mapping (i.e.

from physical space to reference space and vice-versa) operations on infinite elements.

Definition at line 47 of file inf_fe_map.h.

Member Function Documentation

eval()

Real libMesh::InfFEMap::eval ( Real  v, Order  o, unsigned int  i  )

static

Returns

The value of the \( i^{th} \) polynomial evaluated at v. This method provides the approximation in radial direction for the overall shape functions, which is defined in InfFE::shape(). This method is allowed to be static, since it is independent of dimension and base_family. It is templated, though, w.r.t. to radial FEFamily.

The value of the \( i^{th} \) mapping shape function in radial direction evaluated at v when T_radial == INFINITE_MAP. Currently, only one specific mapping shape is used. Namely the one by Marques JMMC, Owen DRJ: Infinite elements in quasi-static materially nonlinear problems, Computers and Structures, 1984.

Definition at line 27 of file inf_fe_map_eval.C.

{
  libmesh_assert (-1.-1.e-5 <= v && v < 1.);
 
  switch (i)
    {
    case 0:
      return -2.*v/(1.-v);
    case 1:
      return (1.+v)/(1.-v);
 
    default:
      libmesh_error_msg("bad index i = " << i);
    }
}

References libMesh::libmesh_assert().

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::eval(), libMesh::InfFE< Dim, T_radial, T_map >::init_radial_shape_functions(), and map().

eval_deriv()

Real libMesh::InfFEMap::eval_deriv ( Real  v, Order  o, unsigned int  i  )

static

Returns

The value of the first derivative of the \( i^{th} \) polynomial at coordinate v. See eval for details.

Definition at line 45 of file inf_fe_map_eval.C.

{
  libmesh_assert (-1.-1.e-5 <= v && v < 1.);
 
  switch (i)
    {
    case 0:
      return -2./((1.-v)*(1.-v));
    case 1:
      return 2./((1.-v)*(1.-v));
 
    default:
      libmesh_error_msg("bad index i = " << i);
    }
}

References libMesh::libmesh_assert().

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::eval_deriv(), and libMesh::InfFE< Dim, T_radial, T_map >::init_radial_shape_functions().

inverse_map() [1/2]

Point libMesh::InfFEMap::inverse_map ( const unsigned int  dim, const Elem *  elem, const Point &  p, const Real  tolerance = TOLERANCE, const bool  secure = true  )

static

Returns

The location (on the reference element) of the point p located in physical space. First, the location in the base face is computed. This requires inverting the (possibly nonlinear) transformation map in the base, so it is not trivial. The optional parameter tolerance defines how close is "good enough." The map inversion iteration computes the sequence \( \{ p_n \} \), and the iteration is terminated when \( \|p - p_n\| < \mbox{\texttt{tolerance}} \). Once the base face point is determined, the radial local coordinate is directly evaluated. If interpolated is true, the interpolated distance from the base element to the infinite element origin is used for the map in radial direction.

Definition at line 88 of file inf_fe_map.C.

{
  libmesh_assert(inf_elem);
  libmesh_assert_greater_equal (tolerance, 0.);
  libmesh_assert(dim > 0);
 
  // Start logging the map inversion.
  LOG_SCOPE("inverse_map()", "InfFEMap");
 
  // The strategy is:
  // compute the intersection of the line
  // physical_point - origin with the base element,
  // find its internal coordinatels using FEMap::inverse_map():
  // The radial part can then be computed directly later on.
 
  // 1.)
  // build a base element to do the map inversion in the base face
  std::unique_ptr<Elem> base_elem (InfFEBase::build_elem (inf_elem));
 
  // The point on the reference element (which we are looking for).
  // start with an invalid guess:
  Point p;
  p(dim-1)=-2.;
 
  // 2.)
  // Now find the intersection of a plane represented by the base
  // element nodes and the line given by the origin of the infinite
  // element and the physical point.
  Point intersection;
 
  // the origin of the infinite element
  const Point o = inf_elem->origin();
 
  switch (dim)
    {
      // unnecessary for 1D
    case 1:
      break;
 
    case 2:
      libmesh_error_msg("ERROR: InfFE::inverse_map is not yet implemented in 2d");
 
    case 3:
      {
        // references to the nodal points of the base element
        const Point & p0 = base_elem->point(0);
        const Point & p1 = base_elem->point(1);
        const Point & p2 = base_elem->point(2);
 
        // a reference to the physical point
        const Point & fp = physical_point;
 
        // The intersection of the plane and the line is given by
        // can be computed solving a linear 3x3 system
        // a*({p1}-{p0})+b*({p2}-{p0})-c*({fp}-{o})={fp}-{p0}.
        const Real c_factor = -(p1(0)*fp(1)*p0(2)-p1(0)*fp(2)*p0(1)
                                +fp(0)*p1(2)*p0(1)-p0(0)*fp(1)*p1(2)
                                +p0(0)*fp(2)*p1(1)+p2(0)*fp(2)*p0(1)
                                -p2(0)*fp(1)*p0(2)-fp(0)*p2(2)*p0(1)
                                +fp(0)*p0(2)*p2(1)+p0(0)*fp(1)*p2(2)
                                -p0(0)*fp(2)*p2(1)-fp(0)*p0(2)*p1(1)
                                +p0(2)*p2(0)*p1(1)-p0(1)*p2(0)*p1(2)
                                -fp(0)*p1(2)*p2(1)+p2(1)*p0(0)*p1(2)
                                -p2(0)*fp(2)*p1(1)-p1(0)*fp(1)*p2(2)
                                +p2(2)*p1(0)*p0(1)+p1(0)*fp(2)*p2(1)
                                -p0(2)*p1(0)*p2(1)-p2(2)*p0(0)*p1(1)
                                +fp(0)*p2(2)*p1(1)+p2(0)*fp(1)*p1(2))/
          (fp(0)*p1(2)*p0(1)-p1(0)*fp(2)*p0(1)
           +p1(0)*fp(1)*p0(2)-p1(0)*o(1)*p0(2)
           +o(0)*p2(2)*p0(1)-p0(0)*fp(2)*p2(1)
           +p1(0)*o(1)*p2(2)+fp(0)*p0(2)*p2(1)
           -fp(0)*p1(2)*p2(1)-p0(0)*o(1)*p2(2)
           +p0(0)*fp(1)*p2(2)-o(0)*p0(2)*p2(1)
           +o(0)*p1(2)*p2(1)-p2(0)*fp(2)*p1(1)
           +fp(0)*p2(2)*p1(1)-p2(0)*fp(1)*p0(2)
           -o(2)*p0(0)*p1(1)-fp(0)*p0(2)*p1(1)
           +p0(0)*o(1)*p1(2)+p0(0)*fp(2)*p1(1)
           -p0(0)*fp(1)*p1(2)-o(0)*p1(2)*p0(1)
           -p2(0)*o(1)*p1(2)-o(2)*p2(0)*p0(1)
           -o(2)*p1(0)*p2(1)+o(2)*p0(0)*p2(1)
           -fp(0)*p2(2)*p0(1)+o(2)*p2(0)*p1(1)
           +p2(0)*o(1)*p0(2)+p2(0)*fp(1)*p1(2)
           +p2(0)*fp(2)*p0(1)-p1(0)*fp(1)*p2(2)
           +p1(0)*fp(2)*p2(1)-o(0)*p2(2)*p1(1)
           +o(2)*p1(0)*p0(1)+o(0)*p0(2)*p1(1));
 
 
        // Check whether the above system is ill-posed.  It should
        // only happen when \p physical_point is not in \p
        // inf_elem. In this case, any point that is not in the
        // element is a valid answer.
        if (libmesh_isinf(c_factor))
          return p;
 
        // The intersection should always be between the origin and the physical point.
        // It can become positive close to the border, but should
        // be only very small.
        //  So as in the case above, we can be sufficiently sure here
        // that \p fp is not in this element:
        if (c_factor > 0.01)
          return p;
 
        // Compute the intersection with
        // {intersection} = {fp} + c*({fp}-{o}).
        intersection.add_scaled(fp,1.+c_factor);
        intersection.add_scaled(o,-c_factor);
 
        break;
      }
 
    default:
      libmesh_error_msg("Invalid dim = " << dim);
    }
 
  // 3.)
  // Now we have the intersection-point (projection of physical point onto base-element).
  // Lets compute its internal coordinates (being p(0) and p(1)):
  p= FEMap::inverse_map(dim-1, base_elem.get(), intersection,
                        tolerance, secure);
 
  // 4.
  // Now that we have the local coordinates in the base,
  // we compute the radial distance with Newton iteration.
 
  // distance from the physical point to the ifem origin
  const Real fp_o_dist = Point(o-physical_point).norm();
 
  // the distance from the intersection on the
  // base to the origin
  const Real a_dist = Point(o-intersection).norm();
 
  // element coordinate in radial direction
 
  // fp_o_dist is at infinity.
  if (libmesh_isinf(fp_o_dist))
    {
      p(dim-1)=1;
      return p;
    }
 
  // when we are somewhere in this element:
  Real v = 0;
 
  // For now we're sticking with T_map == CARTESIAN
  // if (T_map == CARTESIAN)
  v = 1.-2.*a_dist/fp_o_dist;
  // else
  //   libmesh_not_implemented();
 
  // do not put the point back into the element, otherwise the contains_point-function
  // gives false positives!
  //if (v <= -1.-1e-5)
  //  v=-1.;
  //if (v >= 1.)
  //  v=1.-1e-5;
 
  p(dim-1)=v;
#ifdef DEBUG
  // first check whether we are in the reference-element:
  if (-1.-1.e-5 < v && v < 1.)
    {
      const Point check = map (dim, inf_elem, p);
      const Point diff  = physical_point - check;
 
      if (diff.norm() > tolerance)
        libmesh_warning("WARNING:  diff is " << diff.norm());
    }
#endif
 
  return p;
}

References libMesh::TypeVector< T >::add_scaled(), libMesh::InfFEBase::build_elem(), dim, libMesh::FEMap::inverse_map(), libMesh::libmesh_assert(), libMesh::libmesh_isinf(), map(), libMesh::TypeVector< T >::norm(), libMesh::Elem::origin(), and libMesh::Real.

Referenced by libMesh::InfQuad4::contains_point(), libMesh::InfPrism::contains_point(), libMesh::InfHex::contains_point(), inverse_map(), libMesh::FEMap::inverse_map(), and libMesh::InfFE< Dim, T_radial, T_map >::inverse_map().

inverse_map() [2/2]

void libMesh::InfFEMap::inverse_map ( const unsigned int  dim, const Elem *  elem, const std::vector< Point > &  physical_points, std::vector< Point > &  reference_points, const Real  tolerance = TOLERANCE, const bool  secure = true  )

static

Takes a number points in physical space (in the physical_points vector) and finds their location on the reference element for the input element elem.

The values on the reference element are returned in the vector reference_points

Definition at line 266 of file inf_fe_map.C.

{
  // The number of points to find the
  // inverse map of
  const std::size_t n_points = physical_points.size();
 
  // Resize the vector to hold the points
  // on the reference element
  reference_points.resize(n_points);
 
  // Find the coordinates on the reference
  // element of each point in physical space
  for (unsigned int p=0; p<n_points; p++)
    reference_points[p] =
      inverse_map (dim, elem, physical_points[p], tolerance, secure);
}

References dim, inverse_map(), and libMesh::TypeVector< T >::size().

map()

Point libMesh::InfFEMap::map ( const unsigned int  dim, const Elem *  inf_elem, const Point &  reference_point  )

static

Returns

The location (in physical space) of the point p located on the reference element.

Definition at line 40 of file inf_fe_map.C.

{
  libmesh_assert(inf_elem);
  libmesh_assert_not_equal_to (dim, 0);
 
  std::unique_ptr<Elem>      base_elem (InfFEBase::build_elem (inf_elem));
 
  const Order        radial_mapping_order (InfFERadial::mapping_order());
  const Real         v                    (reference_point(dim-1));
 
  // map in the base face
  Point base_point;
  switch (dim)
    {
    case 1:
      base_point = inf_elem->point(0);
      break;
    case 2:
    case 3:
      base_point = FEMap::map (dim-1, base_elem.get(), reference_point);
      break;
    default:
#ifdef DEBUG
      libmesh_error_msg("Unknown dim = " << dim);
#endif
      break;
    }
 
 
  // map in the outer node face not necessary. Simply
  // compute the outer_point = base_point + (base_point-origin)
  const Point outer_point (base_point * 2. - inf_elem->origin());
 
  Point p;
 
  // there are only two mapping shapes in radial direction
  p.add_scaled (base_point,  eval (v, radial_mapping_order, 0));
  p.add_scaled (outer_point, eval (v, radial_mapping_order, 1));
 
  return p;
}

References libMesh::TypeVector< T >::add_scaled(), libMesh::InfFEBase::build_elem(), dim, eval(), libMesh::libmesh_assert(), libMesh::FEMap::map(), libMesh::InfFERadial::mapping_order(), libMesh::Elem::origin(), libMesh::Elem::point(), and libMesh::Real.

Referenced by inverse_map(), libMesh::FEMap::map(), and libMesh::InfFE< Dim, T_radial, T_map >::map().

---

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

• include/fe/inf_fe_map.h
• src/fe/inf_fe_map.C
• src/fe/inf_fe_map_eval.C