libMesh
Public Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
libMesh::VariationalMeshSmoother Class Reference

This is an implementation of Larisa Branets' smoothing algorithms. More...

#include <mesh_smoother_vsmoother.h>

Inheritance diagram for libMesh::VariationalMeshSmoother:
[legend]

Public Member Functions

 VariationalMeshSmoother (UnstructuredMesh &mesh, Real dilation_weight=0.5, const bool preserve_subdomain_boundaries=true, const double relative_residual_tolerance=TOLERANCE *TOLERANCE, const double absolute_residual_tolerance=TOLERANCE *TOLERANCE, const bool solver_quiet=true, const bool solver_verbose=false)
 Constructor. More...
 
virtual ~VariationalMeshSmoother ()=default
 Destructor. More...
 
virtual void setup ()
 Setup method that creates equation systems, system, and constraints, to be called just prior to smoothing. More...
 
virtual void smooth () override
 Redefinition of the smooth function from the base class. More...
 
void smooth (unsigned int n_iterations)
 The actual smoothing function, gets called whenever the user specifies an actual number of smoothing iterations. More...
 
const MeshQualityInfoget_mesh_info () const
 Getter for the _system's _mesh_info attribute. More...
 

Protected Attributes

UnstructuredMesh_mesh
 

Private Member Functions

VariationalSmootherSystemsystem () const
 Getter for _system to protect against dangling pointers. More...
 

Private Attributes

const Real _dilation_weight
 Smoother control variables. More...
 
const bool _preserve_subdomain_boundaries
 Whether subdomain boundaries are subject to change via smoothing. More...
 
std::unique_ptr< DistributedMesh_mesh_copy
 Mesh copy to avoid multiple EquationSystems. More...
 
std::unique_ptr< EquationSystems_equation_systems
 EquationsSystems object associated with the smoother. More...
 
VariationalSmootherSystem_system
 
std::unique_ptr< VariationalSmootherConstraint_constraint
 Constraints imposed on the smoothing process. More...
 
bool _setup_called
 Attribute the keep track of whether the setup method has been called. More...
 
double _relative_residual_tolerance
 Solver relative residual tolerance. More...
 
double _absolute_residual_tolerance
 Solver absolute residual tolerance. More...
 
bool _solver_quiet
 Solver quiet setting. More...
 
bool _solver_verbose
 Solver verbose setting. More...
 

Detailed Description

This is an implementation of Larisa Branets' smoothing algorithms.

The initial implementation was done by her, the adaptation to libmesh was completed by Derek Gaston. The code was heavily refactored into something more closely resembling C++ by John Peterson in 2014. The code eventually fell out of use and stopped functioning. Patrick Behne reimplemented the smoother in 2025.

Here are the relevant publications: 1) L. Branets, G. Carey, "Extension of a mesh quality metric for elements with a curved boundary edge or surface", Journal of Computing and Information Science in Engineering, vol. 5(4), pp.302-308, 2005.

2) L. Branets, G. Carey, "A local cell quality metric and variational grid smoothing algorithm", Engineering with Computers, vol. 21, pp.19-28, 2005.

3) L. Branets, "A variational grid optimization algorithm based on a local cell quality metric", Ph.D. thesis, The University of Texas at Austin, 2005.

Notes:

1) The smoother supports tangled meshes. However, not all untangling solves converge. In other words, we are able to untangle SOME tangled meshes, but not ANY tangled mesh.

Author
Derek R. Gaston
Date
2006

Definition at line 75 of file mesh_smoother_vsmoother.h.

Constructor & Destructor Documentation

◆ VariationalMeshSmoother()

libMesh::VariationalMeshSmoother::VariationalMeshSmoother ( UnstructuredMesh mesh,
Real  dilation_weight = 0.5,
const bool  preserve_subdomain_boundaries = true,
const double  relative_residual_tolerance = TOLERANCE * TOLERANCE,
const double  absolute_residual_tolerance = TOLERANCE * TOLERANCE,
const bool  solver_quiet = true,
const bool  solver_verbose = false 
)

Constructor.

Parameters
meshMesh to smooth by minimizing the distortion-dilation metric.
dilation_wiehgtWeight to give the dilation metric. The distortion metric is given 1 - dilation_weight.
preserve_subdomain_boundariesWhether the smoother is to preserve or modify sobdomain boundaries.
relative_residual_toleranceSolver setting for the relative residual tolerance.
absolute_residual_toleranceSolver setting for the absolute residual tolerance.
solver_quietWhether to make the solver quiet.
solver_verboseWhether to make the solver verbose.

Definition at line 51 of file mesh_smoother_vsmoother.C.

58  : MeshSmoother(mesh),
59  _dilation_weight(dilation_weight),
60  _preserve_subdomain_boundaries(preserve_subdomain_boundaries),
61  _setup_called(false),
62  _relative_residual_tolerance(relative_residual_tolerance),
63  _absolute_residual_tolerance(absolute_residual_tolerance),
64  _solver_quiet(solver_quiet),
65  _solver_verbose(solver_verbose)
66 {}
const Real _dilation_weight
Smoother control variables.
bool _solver_verbose
Solver verbose setting.
MeshBase & mesh
double _absolute_residual_tolerance
Solver absolute residual tolerance.
bool _setup_called
Attribute the keep track of whether the setup method has been called.
MeshSmoother(UnstructuredMesh &mesh)
Constructor.
Definition: mesh_smoother.h:46
bool _solver_quiet
Solver quiet setting.
double _relative_residual_tolerance
Solver relative residual tolerance.
const bool _preserve_subdomain_boundaries
Whether subdomain boundaries are subject to change via smoothing.

◆ ~VariationalMeshSmoother()

virtual libMesh::VariationalMeshSmoother::~VariationalMeshSmoother ( )
virtualdefault

Destructor.

Member Function Documentation

◆ get_mesh_info()

const MeshQualityInfo & libMesh::VariationalMeshSmoother::get_mesh_info ( ) const

Getter for the _system's _mesh_info attribute.

Definition at line 166 of file mesh_smoother_vsmoother.C.

References _setup_called, _system, and libMesh::VariationalSmootherSystem::get_mesh_info().

Referenced by MeshSmootherTest::testVariationalSmoother().

167 {
168  libmesh_error_msg_if(!_setup_called, "Need to first call the setup() method of "
169  << "this VariationalMeshSmoother object, then call get_mesh_info() prior "
170  << "to calling smooth().");
171  return _system->get_mesh_info();
172 }
VariationalSmootherSystem * _system
bool _setup_called
Attribute the keep track of whether the setup method has been called.
const MeshQualityInfo & get_mesh_info()
Getter for the _mesh_info attribute.

◆ setup()

void libMesh::VariationalMeshSmoother::setup ( )
virtual

Setup method that creates equation systems, system, and constraints, to be called just prior to smoothing.

Definition at line 68 of file mesh_smoother_vsmoother.C.

References _absolute_residual_tolerance, _constraint, _dilation_weight, _equation_systems, libMesh::MeshSmoother::_mesh, _mesh_copy, _preserve_subdomain_boundaries, _relative_residual_tolerance, _setup_called, _solver_quiet, _solver_verbose, _system, libMesh::System::attach_constraint_object(), libMesh::MeshBase::elem_dimensions(), libMesh::VariationalSmootherSystem::get_dilation_weight(), libMesh::DiffSolver::quiet, system(), libMesh::DifferentiableSystem::time_solver, and libMesh::DiffSolver::verbose.

Referenced by smooth(), and MeshSmootherTest::testVariationalSmoother().

69 {
70  // Check for multiple dimensions
71  if (_mesh.elem_dimensions().size() > 1)
72  libmesh_not_implemented_msg("Meshes containing elements of differing dimension are not yet supported.");
73 
74  /*
75  Ideally we'd want to update _mesh directly even if it already has an
76  EquationSystems attached and doing work on it ... and that's thwarted by the
77  problems:
78 
79  1. We can't easily tell if there's already an EquationSystems attached.
80  2. We can't attach a second EquationSystems safely (if it already has a system 0)
81  because the DoF indexing will need to be overwritten by our system.
82  3. The destructor of es won't even clean up after itself (we generally expect
83  a mesh to go unused after its EquationSystems is destroyed), much less know
84  how to restore anything from a previous EquationSystems.
85 
86  To avoid these issues, we'll just construct a new DistributedMesh _mesh_copy
87  from _mesh, then do the solve on _mesh_copy, then copy its node locations back
88  to _mesh after the solve is done. That'll be slightly less memory-efficient
89  and somewhat more CPU-efficient in the case where _mesh is serial, though
90  it'll be significantly less memory-efficient when _mesh is already distributed,
91  but either way the robustness is probably worth it.
92  */
93 
94  // Create a new mesh, EquationSystems, and System
95  _mesh_copy = std::make_unique<DistributedMesh>(_mesh);
96  _equation_systems = std::make_unique<EquationSystems>(*_mesh_copy);
97  _system =
98  &(_equation_systems->add_system<VariationalSmootherSystem>("variational_smoother_system"));
99 
100  // Set this to something > 0 to add more quadrature points than the default
101  // rule that integrates order 2 * fe_order + 1 polynomials exactly.
102  // Using higher quadrature orders has not had a significant effect on observed solutions.
103  //system()->extra_quadrature_order = 0;
104 
105  // Uncomment these to debug
106  //system()->print_element_solutions=true;
107  //system()->print_element_residuals=true;
108  //system()->print_element_jacobians=true;
109 
110  // Add boundary node and hanging node constraints
111  _constraint =
112  std::make_unique<VariationalSmootherConstraint>(*_system, _preserve_subdomain_boundaries);
114 
115  // Set system parameters
117 
118  // Set up solver
119  system()->time_solver = std::make_unique<SteadySolver>(*_system);
120 
121  // Uncomment this line and use -snes_test_jacobian and -snes_test_jacobian_view
122  // flags to compare the hand-coded jacobian in VariationalSmootherSystem
123  // to finite difference jacobians.
124  //system()->time_solver->diff_solver() = std::make_unique<PetscDiffSolver>(*_system);
125 
126  _equation_systems->init();
127 
128  // Solver verbosity
129  DiffSolver & solver = *(system()->time_solver->diff_solver().get());
130  solver.quiet = _solver_quiet;
131  solver.verbose = _solver_verbose;
132 
133  // Solver convergence tolerances
134  system()->time_solver->diff_solver()->relative_residual_tolerance = _relative_residual_tolerance;
135  system()->time_solver->diff_solver()->absolute_residual_tolerance = _absolute_residual_tolerance;
136 
137  _setup_called = true;
138 }
const Real _dilation_weight
Smoother control variables.
bool _solver_verbose
Solver verbose setting.
std::unique_ptr< EquationSystems > _equation_systems
EquationsSystems object associated with the smoother.
std::unique_ptr< TimeSolver > time_solver
A pointer to the solver object we&#39;re going to use.
Definition: diff_system.h:253
VariationalSmootherSystem * _system
UnstructuredMesh & _mesh
Definition: mesh_smoother.h:61
VariationalSmootherSystem * system() const
Getter for _system to protect against dangling pointers.
std::unique_ptr< VariationalSmootherConstraint > _constraint
Constraints imposed on the smoothing process.
void attach_constraint_object(Constraint &constrain)
Register a user object for imposing constraints.
Definition: system.C:1987
std::unique_ptr< DistributedMesh > _mesh_copy
Mesh copy to avoid multiple EquationSystems.
double _absolute_residual_tolerance
Solver absolute residual tolerance.
bool _setup_called
Attribute the keep track of whether the setup method has been called.
const std::set< unsigned char > & elem_dimensions() const
Definition: mesh_base.h:282
bool _solver_quiet
Solver quiet setting.
double _relative_residual_tolerance
Solver relative residual tolerance.
const bool _preserve_subdomain_boundaries
Whether subdomain boundaries are subject to change via smoothing.

◆ smooth() [1/2]

virtual void libMesh::VariationalMeshSmoother::smooth ( )
inlineoverridevirtual

Redefinition of the smooth function from the base class.

All this does is call the smooth function in this class which takes an int, using a default value of 1.

Implements libMesh::MeshSmoother.

Definition at line 115 of file mesh_smoother_vsmoother.h.

References smooth().

Referenced by main(), smooth(), MeshSmootherTest::testVariationalSmoother(), and MeshSmootherTest::testVariationalSmootherRegression().

115 {this->smooth(1); }
virtual void smooth() override
Redefinition of the smooth function from the base class.

◆ smooth() [2/2]

void libMesh::VariationalMeshSmoother::smooth ( unsigned int  n_iterations)

The actual smoothing function, gets called whenever the user specifies an actual number of smoothing iterations.

Definition at line 140 of file mesh_smoother_vsmoother.C.

References _equation_systems, libMesh::MeshSmoother::_mesh, _mesh_copy, _setup_called, libMesh::ParallelObject::comm(), libMesh::make_range(), libMesh::MeshBase::node_ref(), setup(), libMesh::VariationalSmootherSystem::solve(), libMesh::Parallel::sync_dofobject_data_by_id(), and system().

141 {
142  if (!_setup_called)
143  setup();
144 
145  system()->solve();
146 
147  // Update _mesh from _mesh_copy
148  for (auto * node_copy : _mesh_copy->local_node_ptr_range())
149  {
150  auto & node = _mesh.node_ref(node_copy->id());
151  for (const auto d : make_range(_mesh_copy->mesh_dimension()))
152  node(d) = (*node_copy)(d);
153  }
154 
155  SyncNodalPositions sync_object(_mesh);
156  Parallel::sync_dofobject_data_by_id (_mesh.comm(), _mesh.nodes_begin(), _mesh.nodes_end(), sync_object);
157 
158  // Release memory occupied by _mesh_copy
159  // Destruct this before _mesh_copy because it references _mesh_copy
160  _equation_systems.reset();
161  _mesh_copy.reset();
162  // We'll need to call setup again since we'll have to reconstruct _mesh_copy
163  _setup_called = false;
164 }
virtual void solve() override
Solves the system to smooth the mesh.
std::unique_ptr< EquationSystems > _equation_systems
EquationsSystems object associated with the smoother.
const Parallel::Communicator & comm() const
UnstructuredMesh & _mesh
Definition: mesh_smoother.h:61
VariationalSmootherSystem * system() const
Getter for _system to protect against dangling pointers.
std::unique_ptr< DistributedMesh > _mesh_copy
Mesh copy to avoid multiple EquationSystems.
bool _setup_called
Attribute the keep track of whether the setup method has been called.
void sync_dofobject_data_by_id(const Communicator &comm, const Iterator &range_begin, const Iterator &range_end, SyncFunctor &sync)
Request data about a range of ghost dofobjects uniquely identified by their id.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
Definition: int_range.h:140
virtual const Node & node_ref(const dof_id_type i) const
Definition: mesh_base.h:596
virtual void setup()
Setup method that creates equation systems, system, and constraints, to be called just prior to smoot...

◆ system()

VariationalSmootherSystem* libMesh::VariationalMeshSmoother::system ( ) const
inlineprivate

Getter for _system to protect against dangling pointers.

Definition at line 168 of file mesh_smoother_vsmoother.h.

References _system, and libMesh::libmesh_assert().

Referenced by setup(), and smooth().

169  {
171  return _system;
172  }
VariationalSmootherSystem * _system
libmesh_assert(ctx)

Member Data Documentation

◆ _absolute_residual_tolerance

double libMesh::VariationalMeshSmoother::_absolute_residual_tolerance
private

Solver absolute residual tolerance.

Definition at line 193 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _constraint

std::unique_ptr<VariationalSmootherConstraint> libMesh::VariationalMeshSmoother::_constraint
private

Constraints imposed on the smoothing process.

Definition at line 178 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _dilation_weight

const Real libMesh::VariationalMeshSmoother::_dilation_weight
private

Smoother control variables.

Definition at line 134 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _equation_systems

std::unique_ptr<EquationSystems> libMesh::VariationalMeshSmoother::_equation_systems
private

EquationsSystems object associated with the smoother.

Definition at line 154 of file mesh_smoother_vsmoother.h.

Referenced by setup(), and smooth().

◆ _mesh

UnstructuredMesh& libMesh::MeshSmoother::_mesh
protectedinherited

◆ _mesh_copy

std::unique_ptr<DistributedMesh> libMesh::VariationalMeshSmoother::_mesh_copy
private

Mesh copy to avoid multiple EquationSystems.

Definition at line 148 of file mesh_smoother_vsmoother.h.

Referenced by setup(), and smooth().

◆ _preserve_subdomain_boundaries

const bool libMesh::VariationalMeshSmoother::_preserve_subdomain_boundaries
private

Whether subdomain boundaries are subject to change via smoothing.

Definition at line 139 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _relative_residual_tolerance

double libMesh::VariationalMeshSmoother::_relative_residual_tolerance
private

Solver relative residual tolerance.

Definition at line 188 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _setup_called

bool libMesh::VariationalMeshSmoother::_setup_called
private

Attribute the keep track of whether the setup method has been called.

Definition at line 183 of file mesh_smoother_vsmoother.h.

Referenced by get_mesh_info(), setup(), and smooth().

◆ _solver_quiet

bool libMesh::VariationalMeshSmoother::_solver_quiet
private

Solver quiet setting.

Definition at line 198 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _solver_verbose

bool libMesh::VariationalMeshSmoother::_solver_verbose
private

Solver verbose setting.

Definition at line 203 of file mesh_smoother_vsmoother.h.

Referenced by setup().

◆ _system

VariationalSmootherSystem* libMesh::VariationalMeshSmoother::_system
private

Definition at line 163 of file mesh_smoother_vsmoother.h.

Referenced by get_mesh_info(), setup(), and system().


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