variational_smoother_constraint.h

Go to the documentation of this file.

// The libMesh Finite Element Library.
// Copyright (C) 2002-2026 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
 
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
 
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
 
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 
#ifndef LIBMESH_VARIATIONAL_SMOOTHER_CONSTRAINT_H
#define LIBMESH_VARIATIONAL_SMOOTHER_CONSTRAINT_H
 
// Local Includes
#include "libmesh/system.h"
#include "libmesh/dof_map.h"
 
// C++ includes
#include <variant>
 
namespace libMesh
{
 
// Forward declarations
class PointConstraint;
class LineConstraint;
class PlaneConstraint;
class InvalidConstraint;
 
using ConstraintVariant = std::variant<PointConstraint, LineConstraint,
                                       PlaneConstraint, InvalidConstraint>;
 
class PointConstraint
{
 
public:
  PointConstraint() = default;
 
  PointConstraint(const Point &point, const Real &tol = TOLERANCE * TOLERANCE);
 
  bool operator<(const PointConstraint &other) const;
 
  bool operator==(const PointConstraint &other) const;
 
  bool contains_point(const PointConstraint &p) const { return *this == p; }
 
  ConstraintVariant intersect(const ConstraintVariant &other) const;
 
  const Point &point() const { return _point; }
 
  const Real &tol() const { return _tol; }
 
private:
  // Life is easier if we don't make this const
  Point _point;
 
  Real _tol;
};
 
class LineConstraint
{
public:
  LineConstraint() = default;
 
  LineConstraint(const Point &point, const Point &direction,
                 const Real &tol = TOLERANCE * TOLERANCE);
 
  bool operator<(const LineConstraint &other) const;
 
  bool operator==(const LineConstraint &other) const;
 
  bool contains_point(const PointConstraint &p) const;
 
  bool is_parallel(const LineConstraint &l) const;
 
  bool is_parallel(const PlaneConstraint &p) const;
 
  ConstraintVariant intersect(const ConstraintVariant &other) const;
 
  const Point &point() const { return _point; }
 
  const Point &direction() const { return _direction; }
 
  const Real &tol() const { return _tol; }
 
private:
  // Life is easier if we don't make these const
  Point _point;
 
  Point _direction;
 
  Real _tol;
};
 
class PlaneConstraint
{
 
public:
  PlaneConstraint() = default;
 
  PlaneConstraint(const Point &point, const Point &normal,
                  const Real &tol = TOLERANCE * TOLERANCE);
 
  bool operator<(const PlaneConstraint &other) const;
 
  bool operator==(const PlaneConstraint &other) const;
 
  bool contains_point(const PointConstraint &p) const;
 
  bool contains_line(const LineConstraint &l) const;
 
  bool is_parallel(const PlaneConstraint &p) const;
 
  bool is_parallel(const LineConstraint &l) const;
 
  ConstraintVariant intersect(const ConstraintVariant &other) const;
 
  const Point &point() const { return _point; }
 
  const Point &normal() const { return _normal; }
 
  const Real &tol() const { return _tol; }
 
private:
  // Life is easier if we don't make these const
  Point _point;
 
  Point _normal;
 
  Real _tol;
};
 
class InvalidConstraint
{
 
public:
  InvalidConstraint()
    : _err_msg("We should never get here! The InvalidConstraint object should be "
               "detected and replaced with a valid ConstraintVariant prior to calling "
               "any class methods.")
  {
  }
 
  ConstraintVariant intersect(const ConstraintVariant &) const {
    libmesh_assert_msg(false, _err_msg);
    return *this;
  }
 
  bool contains_point(const PointConstraint &) const {
    libmesh_assert_msg(false, _err_msg);
    return false;
  }
 
private:
  std::string _err_msg;
};
 
inline ConstraintVariant intersect_constraints(const ConstraintVariant &a,
                                               const ConstraintVariant &b) {
  // std::visit applies the visitor v (a Callable that can be called with any
  // combination of types from Variants) to the active value inside a
  // std::Variant. This circumvents the issue that the literal ConstraintVariant
  // type does not have a method called 'intersect' (but the types defining
  // ConstraintVariant do)
  return std::visit(
      [](const auto &lhs, const auto &rhs) -> ConstraintVariant {
        return lhs.intersect(rhs);
      },
      a, b);
}
 
class VariationalSmootherConstraint : public System::Constraint
{
private:
 
  const unsigned int _verbosity;
 
  System & _sys;
 
  const bool _preserve_subdomain_boundaries;
 
  void fix_node(const Node & node);
 
  void constrain_node_to_plane(const Node & node, const Point & ref_normal_vec);
 
  void constrain_node_to_line(const Node & node, const Point & line_vec);
 
  static bool nodes_share_boundary_id(
      const Node & boundary_node,
      const Node & neighbor_node,
      const Elem & containing_elem,
      const BoundaryInfo & boundary_info);
 
  static std::set<std::set<const Node *>>
  get_neighbors_for_subdomain_constraint(
      const MeshBase &mesh, const Node &node, const subdomain_id_type sub_id,
      const std::unordered_map<dof_id_type, std::vector<const Elem *>>
          &nodes_to_elem_map);
 
  static std::set<std::set<const Node *>> get_neighbors_for_boundary_constraint(
      const MeshBase &mesh, const Node &node,
      const std::unordered_set<dof_id_type> &boundary_node_ids,
      const BoundaryInfo &boundary_info,
      const std::unordered_map<dof_id_type, std::vector<const Elem *>>
          &nodes_to_elem_map);
 
  static ConstraintVariant determine_constraint(
      const Node &node, const unsigned int dim,
      const std::set<std::set<const Node *>> &side_grouped_boundary_neighbors);
 
  void impose_constraint(const Node &node, const ConstraintVariant &constraint);
 
public:
  VariationalSmootherConstraint(System & sys,
                                const bool & preserve_subdomain_boundaries,
                                const unsigned int verbosity);
 
  virtual ~VariationalSmootherConstraint() override;
 
  virtual void constrain() override;
};
 
 
} // namespace libMesh
 
#endif // LIBMESH_VARIATIONAL_SMOOTHER_CONSTRAINT_H