|
libMesh
|
This is an implementation of Larisa Branets' smoothing algorithms. More...
#include <mesh_smoother_vsmoother.h>
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 MeshQualityInfo & | get_mesh_info () const |
| Getter for the _system's _mesh_info attribute. More... | |
Protected Attributes | |
| UnstructuredMesh & | _mesh |
Private Member Functions | |
| VariationalSmootherSystem * | system () 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... | |
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.
Definition at line 75 of file mesh_smoother_vsmoother.h.
| 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.
| mesh | Mesh to smooth by minimizing the distortion-dilation metric. |
| dilation_wiehgt | Weight to give the dilation metric. The distortion metric is given 1 - dilation_weight. |
| preserve_subdomain_boundaries | Whether the smoother is to preserve or modify sobdomain boundaries. |
| relative_residual_tolerance | Solver setting for the relative residual tolerance. |
| absolute_residual_tolerance | Solver setting for the absolute residual tolerance. |
| solver_quiet | Whether to make the solver quiet. |
| solver_verbose | Whether to make the solver verbose. |
Definition at line 51 of file mesh_smoother_vsmoother.C.
|
virtualdefault |
Destructor.
| 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().
|
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().
|
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().
| 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().
|
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().
|
private |
Solver absolute residual tolerance.
Definition at line 193 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
private |
Constraints imposed on the smoothing process.
Definition at line 178 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
private |
Smoother control variables.
Definition at line 134 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
private |
EquationsSystems object associated with the smoother.
Definition at line 154 of file mesh_smoother_vsmoother.h.
|
protectedinherited |
Definition at line 61 of file mesh_smoother.h.
Referenced by libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::LaplaceMeshSmoother::init(), setup(), libMesh::LaplaceMeshSmoother::smooth(), and smooth().
|
private |
Mesh copy to avoid multiple EquationSystems.
Definition at line 148 of file mesh_smoother_vsmoother.h.
|
private |
Whether subdomain boundaries are subject to change via smoothing.
Definition at line 139 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
private |
Solver relative residual tolerance.
Definition at line 188 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
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().
|
private |
Solver quiet setting.
Definition at line 198 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
private |
Solver verbose setting.
Definition at line 203 of file mesh_smoother_vsmoother.h.
Referenced by setup().
|
private |
Definition at line 163 of file mesh_smoother_vsmoother.h.
Referenced by get_mesh_info(), setup(), and system().
1.8.14