libMesh::VariationalSmootherConstraint Class Reference

Constraint class for the VariationalMeshSmoother. More...

#include <variational_smoother_constraint.h>

Inheritance diagram for libMesh::VariationalSmootherConstraint:

Public Member Functions
	VariationalSmootherConstraint (System &sys, const bool &preserve_subdomain_boundaries, const unsigned int verbosity)
	Constructor.
virtual	~VariationalSmootherConstraint () override
virtual void	constrain () override
	Constraint function.

Private Member Functions
void	fix_node (const Node &node)
	Constrain (i.e., fix) a node to not move during mesh smoothing.
void	constrain_node_to_plane (const Node &node, const Point &ref_normal_vec)
	Constrain a node to remain in the given plane during mesh smoothing.
void	constrain_node_to_line (const Node &node, const Point &line_vec)
	Constrain a node to remain on the given line during mesh smoothing.
void	impose_constraint (const Node &node, const ConstraintVariant &constraint)
	Applies a given constraint to a node (e.g., fixing it, restricting it to a line or plane).

Static Private Member Functions
static bool	nodes_share_boundary_id (const Node &boundary_node, const Node &neighbor_node, const Elem &containing_elem, const BoundaryInfo &boundary_info)
	Determines whether two neighboring nodes share a common boundary id.
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)
	Get the relevant nodal neighbors for a subdomain constraint.
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)
	Get the relevant nodal neighbors for an external boundary constraint.
static ConstraintVariant	determine_constraint (const Node &node, const unsigned int dim, const std::set< std::set< const Node * > > &side_grouped_boundary_neighbors)
	Determines the appropriate constraint (PointConstraint, LineConstraint, or PlaneConstraint) for a node based on its neighbors.

Private Attributes
const unsigned int	_verbosity
	Verbosity setting.
System &	_sys
const bool	_preserve_subdomain_boundaries
	Whether nodes on subdomain boundaries are subject to change via smoothing.

Detailed Description

Constraint class for the VariationalMeshSmoother.

Currently, all mesh boundary nodes are constrained to not move during smoothing. If requested (preserve_subdomain_boundaries = true), nodes on subdomain boundaries are also constrained to not move.

Definition at line 389 of file variational_smoother_constraint.h.

Constructor & Destructor Documentation

VariationalSmootherConstraint()

libMesh::VariationalSmootherConstraint::VariationalSmootherConstraint	(	System &	sys,
		const bool &	preserve_subdomain_boundaries,
		const unsigned int	verbosity
	)

Constructor.

Parameters

sys	System to constrain.
preserve_subdomain_boundaries	Whether to constrain nodes on subdomain boundaries to not move.

Definition at line 334 of file variational_smoother_constraint.C.

  : Constraint(), _verbosity(verbosity), _sys(sys), _preserve_subdomain_boundaries(preserve_subdomain_boundaries) {
}

~VariationalSmootherConstraint()

libMesh::VariationalSmootherConstraint::~VariationalSmootherConstraint ( )

overridevirtualdefault

Member Function Documentation

constrain()

void libMesh::VariationalSmootherConstraint::constrain ( )

overridevirtual

Constraint function.

This function will be called to constrain the system prior to a solve and must be provided by the user in a derived class.

Implements libMesh::System::Constraint.

Definition at line 341 of file variational_smoother_constraint.C.

{
  const auto & mesh = _sys.get_mesh();
  const auto dim = mesh.mesh_dimension();
  const auto & proc_id = mesh.processor_id();
 
  // Only compute the node to elem map once
  std::unordered_map<dof_id_type, std::vector<const Elem *>> nodes_to_elem_map;
  MeshTools::build_nodes_to_elem_map(mesh, nodes_to_elem_map);
 
  const auto & boundary_info = mesh.get_boundary_info();
 
  const auto boundary_node_ids = MeshTools::find_boundary_nodes(mesh);
 
  // Identify/constrain subdomain boundary nodes, if requested
  std::unordered_map<dof_id_type, ConstraintVariant> subdomain_boundary_map;
  if (_preserve_subdomain_boundaries)
    {
      for (const auto * elem : mesh.active_element_ptr_range())
        {
          const auto & sub_id1 = elem->subdomain_id();
          for (const auto side : elem->side_index_range())
            {
              const auto * neighbor = elem->neighbor_ptr(side);
              if (neighbor == nullptr)
                continue;
 
              const auto & sub_id2 = neighbor->subdomain_id();
              if (sub_id1 == sub_id2)
                continue;
 
              // elem and neighbor are in different subdomains, and share nodes
              // that need to be constrained
              for (const auto local_node_id : elem->nodes_on_side(side))
                {
                  const auto & node = mesh.node_ref(elem->node_id(local_node_id));
                  // Make sure we haven't already processed this node.
                  // We leave the responsibility of computing the constraint to
                  // the proc that owns the node.
                  if (subdomain_boundary_map.count(node.id()) ||
                      node.processor_id() != proc_id)
                    continue;
 
                  // Get the relevant nodal neighbors for the subdomain constraint
                  const auto side_grouped_boundary_neighbors =
                      get_neighbors_for_subdomain_constraint(
                          mesh, node, sub_id1, nodes_to_elem_map);
 
                  // Determine which constraint should be imposed
                  const auto subdomain_constraint =
                      determine_constraint(node, dim, side_grouped_boundary_neighbors);
 
                  // This subdomain boundary node does not lie on an external boundary,
                  // go ahead and impose constraint
                  if (boundary_node_ids.find(node.id()) == boundary_node_ids.end())
                    {
                      if (_verbosity > 20)
                        libMesh::out << "Imposing subdomain constraint on "
                                     << std::endl << "  " << node << std::endl
                                     << "    " << subdomain_constraint << std::endl;
 
                      this->impose_constraint(node, subdomain_constraint);
                    }
 
                  // This subdomain boundary node could lie on an external boundary, save it
                  // for later to combine with the external boundary constraint.
                  // We also save constraints for non-boundary nodes so we don't try to
                  // re-constrain the node when accessed from the neighboring elem.
                  // See subdomain_boundary_map.count call above.
                  subdomain_boundary_map[node.id()] = subdomain_constraint;
 
                } // for local_node_id
 
            } // for side
        } // for elem
    }
 
  // Loop through boundary nodes and impose constraints
  for (const auto & bid : boundary_node_ids)
    {
      const auto & node = mesh.node_ref(bid);
      // We leave the responsibility of computing the constraint to the proc
      // that owns the node.
      if (node.processor_id() != proc_id)
        continue;
 
      // Get the relevant nodal neighbors for the boundary constraint
      const auto side_grouped_boundary_neighbors = get_neighbors_for_boundary_constraint(
          mesh, node, boundary_node_ids, boundary_info, nodes_to_elem_map);
 
      // Determine which constraint should be imposed
      const auto boundary_constraint =
          determine_constraint(node, dim, side_grouped_boundary_neighbors);
 
      // Check for the case where this boundary node is also part of a subdomain id boundary
      if (const auto it = subdomain_boundary_map.find(bid); it != subdomain_boundary_map.end())
        {
          const auto & subdomain_constraint = it->second;
          // Combine current boundary constraint with previously determined
          // subdomain_constraint
          auto combined_constraint = intersect_constraints(subdomain_constraint, boundary_constraint);
 
          // This will catch cases where constraints have no intersection
          // Fall back to fixed node constraint
          if (std::holds_alternative<InvalidConstraint>(combined_constraint))
            combined_constraint = PointConstraint(node);
 
          if (_verbosity > 20)
            libMesh::out << "Imposing boundary/subdomain constraint on "
                         << std::endl << "  " << node << std::endl
                         << "    " << combined_constraint << std::endl;
 
          this->impose_constraint(node, combined_constraint);
        }
 
      else
        {
          if (_verbosity > 20)
            libMesh::out << "Imposing boundary constraint on "
                         << std::endl << node << std::endl
                         << "    " << boundary_constraint << std::endl;
 
          this->impose_constraint(node, boundary_constraint);
        }
 
    } // end bid
}

References _preserve_subdomain_boundaries, _sys, _verbosity, libMesh::MeshTools::build_nodes_to_elem_map(), determine_constraint(), dim, libMesh::MeshTools::find_boundary_nodes(), libMesh::System::get_mesh(), get_neighbors_for_boundary_constraint(), get_neighbors_for_subdomain_constraint(), impose_constraint(), libMesh::intersect_constraints(), mesh, libMesh::MeshBase::mesh_dimension(), and libMesh::out.

constrain_node_to_line()

void libMesh::VariationalSmootherConstraint::constrain_node_to_line ( const Node &  node, const Point &  line_vec  )

private

Constrain a node to remain on the given line during mesh smoothing.

Parameters

node	Node to constrain
line_vec	vector parallel to the constraining line. This, along with the coordinates of node, are used to define the constraining line.

Definition at line 533 of file variational_smoother_constraint.C.

{
  const auto dim = _sys.get_mesh().mesh_dimension();
 
  // We will free the dimension most parallel to line_vec to keep the
  // constraint coefficients small
  const std::vector<Real> line_vec_coefs{line_vec(0), line_vec(1), line_vec(2)};
  auto it = std::max_element(line_vec_coefs.begin(), line_vec_coefs.end(), [](double a, double b) {
    return std::abs(a) < std::abs(b);
  });
  const unsigned int free_dim = std::distance(line_vec_coefs.begin(), it);
  const auto free_dof_index = node.dof_number(_sys.number(), free_dim, 0);
 
  // A line is parameterized as r(t) = node + t * line_vec, so
  // x(t) = node(x) + t * line_vec(x)
  // y(t) = node(y) + t * line_vec(y)
  // z(t) = node(z) + t * line_vec(z)
  // Let's say we leave x free. Then t = (x(t) - node(x)) / line_vec(x)
  // Then y and z can be constrained as
  // y = node(y) + line_vec_y * (x(t) - node(x)) / line_vec(x)
  //   = x(t) * line_vec(y) / line_vec(x) + (node(y) - node(x) * line_vec(y) / line_vec(x))
  // z = x(t) * line_vec(z) / line_vec(x) + (node(z) - node(x) * line_vec(z) / line_vec(x))
 
  libmesh_assert(!relative_fuzzy_equals(line_vec(free_dim), 0.));
  for (const auto constrained_dim : make_range(dim))
    {
      if (constrained_dim == free_dim)
        continue;
 
      DofConstraintRow constraint_row;
      constraint_row[free_dof_index] = line_vec(constrained_dim) / line_vec(free_dim);
      const auto inhomogeneous_part =
          node(constrained_dim) - node(free_dim) * line_vec(constrained_dim) / line_vec(free_dim);
      const auto constrained_dof_index = node.dof_number(_sys.number(), constrained_dim, 0);
      _sys.get_dof_map().add_constraint_row(
          constrained_dof_index, constraint_row, inhomogeneous_part, true);
    }
}

References _sys, libMesh::DofMap::add_constraint_row(), b, dim, libMesh::DofObject::dof_number(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::libmesh_assert(), libMesh::make_range(), libMesh::MeshBase::mesh_dimension(), libMesh::System::number(), and libMesh::relative_fuzzy_equals().

Referenced by impose_constraint().

constrain_node_to_plane()

void libMesh::VariationalSmootherConstraint::constrain_node_to_plane ( const Node &  node, const Point &  ref_normal_vec  )

private

Constrain a node to remain in the given plane during mesh smoothing.

Parameters

node	Node to constrain
ref_normal_vec	Reference normal vector to the constraining plane. This, along with the coordinates of node, are used to define the constraining plane.

Definition at line 483 of file variational_smoother_constraint.C.

{
  const auto dim = _sys.get_mesh().mesh_dimension();
  // determine equation of plane: c_x * x + c_y * y + c_z * z + c = 0
  std::vector<Real> xyz_coefs; // vector to hold c_x, c_y, c_z
  Real c = 0.;
 
  // We choose to constrain the dimension with the largest magnitude coefficient
  // This approach ensures the coefficients added to the constraint_row
  // (i.e., -c_xyz / c_max) have as small magnitude as possible
 
  // We initialize this to avoid maybe-uninitialized compiler error
  unsigned int constrained_dim = 0;
 
  // Let's assert that we have a nonzero normal to ensure that constrained_dim
  // is always set
  libmesh_assert(ref_normal_vec.norm() > TOLERANCE);
 
  Real max_abs_coef = 0.;
  for (const auto d : make_range(dim))
    {
      const auto coef = ref_normal_vec(d);
      xyz_coefs.push_back(coef);
      c -= coef * node(d);
 
      const auto coef_abs = std::abs(coef);
      if (coef_abs > max_abs_coef)
        {
          max_abs_coef = coef_abs;
          constrained_dim = d;
        }
    }
 
  DofConstraintRow constraint_row;
  for (const auto free_dim : make_range(dim))
    {
      if (free_dim == constrained_dim)
        continue;
 
      const auto free_dof_index = node.dof_number(_sys.number(), free_dim, 0);
      constraint_row[free_dof_index] = -xyz_coefs[free_dim] / xyz_coefs[constrained_dim];
    }
 
  const auto inhomogeneous_part = -c / xyz_coefs[constrained_dim];
  const auto constrained_dof_index = node.dof_number(_sys.number(), constrained_dim, 0);
  _sys.get_dof_map().add_constraint_row(
      constrained_dof_index, constraint_row, inhomogeneous_part, true);
}

References _sys, libMesh::DofMap::add_constraint_row(), dim, libMesh::DofObject::dof_number(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::libmesh_assert(), libMesh::make_range(), libMesh::MeshBase::mesh_dimension(), libMesh::TypeVector< T >::norm(), libMesh::System::number(), libMesh::Real, and libMesh::TOLERANCE.

Referenced by impose_constraint().

determine_constraint()

ConstraintVariant libMesh::VariationalSmootherConstraint::determine_constraint ( const Node &  node, const unsigned int  dim, const std::set< std::set< const Node * > > &  side_grouped_boundary_neighbors  )

staticprivate

Determines the appropriate constraint (PointConstraint, LineConstraint, or PlaneConstraint) for a node based on its neighbors.

Parameters

node	The node to constrain.
dim	The mesh dimension.

Returns

The best-fit constraint for the given geometry.

Definition at line 798 of file variational_smoother_constraint.C.

{
  // Determines the appropriate geometric constraint for a node based on its
  // neighbors.
 
  // Extract neighbors in flat vector
  std::vector<const Node *> neighbors;
  for (const auto & side : side_grouped_boundary_neighbors)
    neighbors.insert(neighbors.end(), side.begin(), side.end());
 
  // Constrain the node to it's current location
  // Note that neighbors.size() <= 1 is used here instead of == 1 so we don't
  // crash when trying to access *neighbors[0] a few lines below in the
  // dim == 2 case
  if (dim == 1 || neighbors.size() <= 1)
    return PointConstraint{node};
 
  if (dim == 2 || neighbors.size() == 2)
    {
      // Determine whether node + all neighbors are colinear
      bool all_colinear = true;
      const Point ref_dir = (*neighbors[0] - node).unit();
      for (const auto i : make_range(std::size_t(1), neighbors.size()))
        {
          const Point delta = *(neighbors[i]) - node;
          libmesh_assert(delta.norm() >= TOLERANCE);
          const Point dir = delta.unit();
          if (!dir.relative_fuzzy_equals(ref_dir) && !dir.relative_fuzzy_equals(-ref_dir))
          {
            all_colinear = false;
            break;
          }
        }
 
      if (all_colinear)
        return LineConstraint{node, ref_dir};
 
      return PointConstraint{node};
    }
 
  // dim == 3, neighbors.size() >= 3
  std::set<PlaneConstraint> valid_planes;
  for (const auto & side_nodes : side_grouped_boundary_neighbors)
    {
      std::vector<const Node *> side_nodes_vec(side_nodes.begin(), side_nodes.end());
      for (const auto i : index_range(side_nodes_vec))
        {
          const Point vec_i = (*side_nodes_vec[i] - node);
          for (const auto j : make_range(i))
            {
              const Point vec_j = (*side_nodes_vec[j] - node);
              Point candidate_normal = vec_i.cross(vec_j);
              if (candidate_normal.norm() <= TOLERANCE)
                continue;
 
              const PlaneConstraint candidate_plane{node, candidate_normal};
              valid_planes.emplace(candidate_plane);
            }
        }
    }
 
  // Fall back to point constraint
  if (valid_planes.empty())
    return PointConstraint(node);
 
  // Combine all the planes together to get a common constraint
  auto it = valid_planes.begin();
  ConstraintVariant current = *it++;
  for (; it != valid_planes.end(); ++it)
    {
      current = intersect_constraints(current, *it);
 
      // This will catch cases where constraints have no intersection
      // (i.e., the element surface is non-planar)
      // Fall back to fixed node constraint
      if (std::holds_alternative<InvalidConstraint>(current))
        {
          current = PointConstraint(node);
          break;
        }
    }
 
  return current;
}

References libMesh::TypeVector< T >::cross(), dim, libMesh::index_range(), libMesh::intersect_constraints(), libMesh::libmesh_assert(), libMesh::make_range(), libMesh::TypeVector< T >::norm(), libMesh::TypeVector< T >::relative_fuzzy_equals(), libMesh::TOLERANCE, and libMesh::TypeVector< T >::unit().

Referenced by constrain().

fix_node()

void libMesh::VariationalSmootherConstraint::fix_node ( const Node &  node )

private

Constrain (i.e., fix) a node to not move during mesh smoothing.

Parameters

node	Node to fix.

Definition at line 469 of file variational_smoother_constraint.C.

{
  for (const auto d : make_range(_sys.get_mesh().mesh_dimension()))
    {
      const auto constrained_dof_index = node.dof_number(_sys.number(), d, 0);
      DofConstraintRow constraint_row;
      // Leave the constraint row as all zeros so this dof is independent from other dofs
      const auto constrained_value = node(d);
      // Simply constrain this dof to retain it's current value
      _sys.get_dof_map().add_constraint_row(
          constrained_dof_index, constraint_row, constrained_value, true);
    }
}

References _sys, libMesh::DofMap::add_constraint_row(), libMesh::DofObject::dof_number(), libMesh::System::get_dof_map(), libMesh::System::get_mesh(), libMesh::make_range(), libMesh::MeshBase::mesh_dimension(), and libMesh::System::number().

Referenced by impose_constraint().

get_neighbors_for_boundary_constraint()

std::set< std::set< const Node * > > libMesh::VariationalSmootherConstraint::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  )

staticprivate

Get the relevant nodal neighbors for an external boundary constraint.

Parameters

mesh	The mesh being smoothed.
node	The node (on the external boundary) being constrained.
boundary_node_ids	The set of mesh's external boundary node ids.
boundary_info	The mesh's BoundaryInfo.
nodes_to_elem_map	A mapping from node id to containing element ids.

Returns

A set of node pointer sets containing nodal neighbors to 'node' on the external boundary. The subsets are grouped by element faces that form the external boundary. Note that 'node' itself does not appear in this set.

Definition at line 717 of file variational_smoother_constraint.C.

{
 
  // Find all the nodal neighbors... that is the nodes directly connected
  // to this node through one edge, or if none exists, use other nodes on the
  // containing face
  std::vector<const Node *> neighbors;
  MeshTools::find_nodal_or_face_neighbors(
      mesh, node, nodes_to_elem_map, neighbors);
 
  // Each constituent set corresponds to neighbors sharing a face on the
  // boundary
  std::set<std::set<const Node *>> side_grouped_boundary_neighbors;
 
  for (const auto * neigh : neighbors)
    {
      const bool is_neighbor_boundary_node =
          boundary_node_ids.find(neigh->id()) != boundary_node_ids.end();
      if (!is_neighbor_boundary_node)
        continue;
 
      // Determine whether nodes share a common boundary id
      // First, find the common element that both node and neigh belong to
      const auto & elems_containing_node = libmesh_map_find(nodes_to_elem_map, node.id());
      const auto & elems_containing_neigh = libmesh_map_find(nodes_to_elem_map, neigh->id());
      const Elem * common_elem = nullptr;
      for (const auto * neigh_elem : elems_containing_neigh)
        {
          const bool is_neigh_common =
              std::find(elems_containing_node.begin(), elems_containing_node.end(), neigh_elem) !=
              elems_containing_node.end();
          if (!is_neigh_common)
            continue;
          common_elem = neigh_elem;
          // Keep this in the neigh_elem loop because there can be multiple common
          // elements Now, determine whether node and neigh share a common boundary
          // id
          const bool nodes_have_common_bid = VariationalSmootherConstraint::nodes_share_boundary_id(
              node, *neigh, *common_elem, boundary_info);
          if (nodes_have_common_bid)
            {
              // Find the side coinciding with the shared boundary
              for (const auto side : common_elem->side_index_range())
                {
                  // We only care about external boundaries here, make sure side doesn't
                  // have a neighbor
                  if (common_elem->neighbor_ptr(side))
                    continue;
 
                  bool node_found_on_side = false;
                  bool neigh_found_on_side = false;
                  const auto nodes_on_side = common_elem->nodes_on_side(side);
                  for (const auto local_node_id : nodes_on_side)
                    {
                      if (common_elem->node_id(local_node_id) == node.id())
                        node_found_on_side = true;
                      else if (common_elem->node_id(local_node_id) == neigh->id())
                        neigh_found_on_side = true;
                    }
                  if (!(node_found_on_side && neigh_found_on_side))
                    continue;
 
                  std::set<const Node *> node_ptrs_on_side;
                  for (const auto local_node_id : nodes_on_side)
                    node_ptrs_on_side.insert(common_elem->node_ptr(local_node_id));
                  node_ptrs_on_side.erase(node_ptrs_on_side.find(&node));
                  side_grouped_boundary_neighbors.insert(node_ptrs_on_side);
                }
              continue;
            }
        }
    }
 
  return side_grouped_boundary_neighbors;
}

References libMesh::MeshTools::find_nodal_or_face_neighbors(), libMesh::DofObject::id(), mesh, libMesh::Elem::neighbor_ptr(), libMesh::Elem::node_id(), libMesh::Elem::node_ptr(), libMesh::Elem::nodes_on_side(), nodes_share_boundary_id(), and libMesh::Elem::side_index_range().

Referenced by constrain().

get_neighbors_for_subdomain_constraint()

std::set< std::set< const Node * > > libMesh::VariationalSmootherConstraint::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  )

staticprivate

Get the relevant nodal neighbors for a subdomain constraint.

Parameters

mesh	The mesh being smoothed.
node	The node (on the subdomain boundary) being constrained.
sub_id	The subdomain id of the block on one side of the subdomain boundary.
nodes_to_elem_map	A mapping from node id to containing element ids.

Returns

A set of node pointer sets containing nodal neighbors to 'node' on the sub_id1-sub_id2 boundary. The subsets are grouped by element faces that form the subdomain boundary. Note that 'node' itself does not appear in this set.

Definition at line 632 of file variational_smoother_constraint.C.

{
 
  // Find all the nodal neighbors... that is the nodes directly connected
  // to this node through one edge, or if none exists, use other nodes on the
  // containing face
  std::vector<const Node *> neighbors;
  MeshTools::find_nodal_or_face_neighbors(
      mesh, node, nodes_to_elem_map, neighbors);
 
  // Each constituent set corresponds to neighbors sharing a face on the
  // subdomain boundary
  std::set<std::set<const Node *>> side_grouped_boundary_neighbors;
 
  for (const auto * neigh : neighbors)
    {
      // Determine whether the neighbor is on the subdomain boundary
      // First, find the common elements that both node and neigh belong to
      const auto & elems_containing_node = libmesh_map_find(nodes_to_elem_map, node.id());
      const auto & elems_containing_neigh = libmesh_map_find(nodes_to_elem_map, neigh->id());
      const Elem * common_elem = nullptr;
      for (const auto * neigh_elem : elems_containing_neigh)
        {
          if ((std::find(elems_containing_node.begin(), elems_containing_node.end(), neigh_elem) !=
               elems_containing_node.end())
              // We should be able to find a common element on the sub_id boundary
              && (neigh_elem->subdomain_id() == sub_id))
            common_elem = neigh_elem;
          else
            continue;
 
          // Now, determine whether node and neigh are on a side coincident
          // with the subdomain boundary
          for (const auto common_side : common_elem->side_index_range())
            {
              bool node_found_on_side = false;
              bool neigh_found_on_side = false;
              for (const auto local_node_id : common_elem->nodes_on_side(common_side))
                {
                  if (common_elem->node_id(local_node_id) == node.id())
                    node_found_on_side = true;
                  else if (common_elem->node_id(local_node_id) == neigh->id())
                    neigh_found_on_side = true;
                }
 
              if (!(node_found_on_side && neigh_found_on_side &&
                    common_elem->neighbor_ptr(common_side)))
                continue;
 
              const auto matched_side = common_side;
              // There could be multiple matched sides, so keep this next part
              // inside the common_side loop
              //
              // Does matched_side, containing both node and neigh, lie on the
              // sub_id subdomain boundary?
              const auto matched_neighbor_sub_id =
                  common_elem->neighbor_ptr(matched_side)->subdomain_id();
              const bool is_matched_side_on_subdomain_boundary = matched_neighbor_sub_id != sub_id;
 
              if (is_matched_side_on_subdomain_boundary)
                {
                  // Store all nodes that live on this side
                  const auto nodes_on_side = common_elem->nodes_on_side(common_side);
                  std::set<const Node *> node_ptrs_on_side;
                  for (const auto local_node_id : nodes_on_side)
                    node_ptrs_on_side.insert(common_elem->node_ptr(local_node_id));
                  node_ptrs_on_side.erase(node_ptrs_on_side.find(&node));
                  side_grouped_boundary_neighbors.insert(node_ptrs_on_side);
 
                  continue;
                }
 
            } // for common_side
 
        } // for neigh_elem
    }
 
  return side_grouped_boundary_neighbors;
}

References libMesh::MeshTools::find_nodal_or_face_neighbors(), libMesh::DofObject::id(), mesh, libMesh::Elem::neighbor_ptr(), libMesh::Elem::node_id(), libMesh::Elem::node_ptr(), libMesh::Elem::nodes_on_side(), libMesh::Elem::side_index_range(), and libMesh::Elem::subdomain_id().

Referenced by constrain().

impose_constraint()

void libMesh::VariationalSmootherConstraint::impose_constraint ( const Node &  node, const ConstraintVariant &  constraint  )

private

Applies a given constraint to a node (e.g., fixing it, restricting it to a line or plane).

Parameters

node	The node to constrain.
constraint	The geometric constraint variant to apply.

Exceptions

libMesh::logicError

If the constraint cannot be imposed.

Definition at line 888 of file variational_smoother_constraint.C.

{
 
  libmesh_assert_msg(!std::holds_alternative<InvalidConstraint>(constraint),
                     "Cannot impose constraint using InvalidConstraint.");
 
  if (std::holds_alternative<PointConstraint>(constraint))
    fix_node(node);
  else if (std::holds_alternative<LineConstraint>(constraint))
    constrain_node_to_line(node, std::get<LineConstraint>(constraint).direction());
  else if (std::holds_alternative<PlaneConstraint>(constraint))
    constrain_node_to_plane(node, std::get<PlaneConstraint>(constraint).normal());
 
  else
    libmesh_assert_msg(false, "Unknown constraint type.");
}

References constrain_node_to_line(), constrain_node_to_plane(), and fix_node().

Referenced by constrain().

nodes_share_boundary_id()

bool libMesh::VariationalSmootherConstraint::nodes_share_boundary_id ( const Node &  boundary_node, const Node &  neighbor_node, const Elem &  containing_elem, const BoundaryInfo &  boundary_info  )

staticprivate

Determines whether two neighboring nodes share a common boundary id.

Parameters

boundary_node	The first of the two prospective nodes.
neighbor_node	The second of the two prospective nodes.
containing_elem	The element containing node1 and node2.
boundary_info	The mesh's BoundaryInfo.

Returns

nodes_share_bid Whether node1 and node2 share a common boundary id.

Definition at line 600 of file variational_smoother_constraint.C.

{
  bool nodes_share_bid = false;
 
  // Node ids local to containing_elem
  const auto node_id = containing_elem.get_node_index(&boundary_node);
  const auto neighbor_id = containing_elem.get_node_index(&neighbor_node);
 
  for (const auto side_id : containing_elem.side_index_range())
    {
      // We don't care about this side if it doesn't contain our boundary and neighbor nodes
      if (!(containing_elem.is_node_on_side(node_id, side_id) &&
            containing_elem.is_node_on_side(neighbor_id, side_id)))
        continue;
 
      // If the current side, containing boundary_node and neighbor_node, lies on a boundary,
      // we can say that boundary_node and neighbor_node have a common boundary id.
      std::vector<boundary_id_type> boundary_ids;
      boundary_info.boundary_ids(&containing_elem, side_id, boundary_ids);
      if (boundary_ids.size())
        {
          nodes_share_bid = true;
          break;
        }
    }
  return nodes_share_bid;
}

References libMesh::BoundaryInfo::boundary_ids(), libMesh::Elem::get_node_index(), libMesh::Elem::is_node_on_side(), side_id, and libMesh::Elem::side_index_range().

Referenced by get_neighbors_for_boundary_constraint().

Member Data Documentation

_preserve_subdomain_boundaries

const bool libMesh::VariationalSmootherConstraint::_preserve_subdomain_boundaries

private

Whether nodes on subdomain boundaries are subject to change via smoothing.

Definition at line 410 of file variational_smoother_constraint.h.

Referenced by constrain().

_sys

System& libMesh::VariationalSmootherConstraint::_sys

private

Definition at line 405 of file variational_smoother_constraint.h.

Referenced by constrain(), constrain_node_to_line(), constrain_node_to_plane(), and fix_node().

_verbosity

const unsigned int libMesh::VariationalSmootherConstraint::_verbosity

private

Verbosity setting.

The verbosity levels and the corresponding information output are as follows:

verbosity = 0: No information.

20 < verbosity: Prints:

• Constraint information on all boundary and subdomain nodes

Definition at line 403 of file variational_smoother_constraint.h.

Referenced by constrain().

---

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

• include/systems/variational_smoother_constraint.h
• src/systems/variational_smoother_constraint.C