libMesh
variational_smoother_system.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 #ifndef LIBMESH_VARIATIONAL_SMOOTHER_SYSTEM_H
19 #define LIBMESH_VARIATIONAL_SMOOTHER_SYSTEM_H
20 
21 // libMesh includes
22 #include "libmesh/enum_fe_family.h"
23 #include "libmesh/fem_function_base.h"
24 #include "libmesh/fem_system.h"
25 #include "libmesh/libmesh_common.h"
26 
27 // C++ includes
28 #include <map>
29 #include <memory>
30 
31 namespace libMesh
32 {
33 
38 {
39  // dof_id_types will hold the id of the Elem where the metric occurs
40  // IMPORTANT: the Real should be the first entry of the pair so that taking
41  // the min/max accross processors will compare the numeric values instead of
42  // the element ids.
43 
44  std::pair<Real, dof_id_type> max_elem_distortion{std::numeric_limits<Real>::lowest(),
46  std::pair<Real, dof_id_type> min_elem_distortion{std::numeric_limits<Real>::max(),
49 
50  std::pair<Real, dof_id_type> max_elem_dilation{std::numeric_limits<Real>::lowest(),
52  std::pair<Real, dof_id_type> min_elem_dilation{std::numeric_limits<Real>::max(),
55 
56  std::pair<Real, dof_id_type> max_elem_combined{std::numeric_limits<Real>::lowest(),
58  std::pair<Real, dof_id_type> min_elem_combined{std::numeric_limits<Real>::max(),
61 
62  std::pair<Real, dof_id_type> max_elem_det_S{std::numeric_limits<Real>::lowest(),
64  std::pair<Real, dof_id_type> min_elem_det_S{std::numeric_limits<Real>::max(),
67  Real max_qp_det_S = std::numeric_limits<Real>::lowest();
68  Real min_qp_det_S = std::numeric_limits<Real>::max();
69 
70  bool mesh_is_tangled = false;
71  bool initialized = false;
72 };
73 
74 // FEMSystem, TimeSolver and NewtonSolver will handle most tasks,
75 // but we must specify element residuals
77 {
87 public:
89  const std::string & name,
90  const unsigned int number)
91  : libMesh::FEMSystem(es, name, number),
94  _ref_vol(0.),
95  _dilation_weight(0.5),
96  _untangling_solve(false)
97  {}
98 
99  // Default destructor
100  ~VariationalSmootherSystem() override;
101 
105  virtual void assembly (bool get_residual,
106  bool get_jacobian,
107  bool apply_heterogeneous_constraints = false,
108  bool apply_no_constraints = false) override;
109 
111 
117  virtual void solve() override;
118 
126  static std::pair<std::unique_ptr<Elem>, std::vector<std::unique_ptr<Node>>>
127  get_target_elem(const ElemType & type);
128 
137  static void get_target_to_reference_jacobian(const Elem * const target_elem,
138  const FEMContext & femcontext,
139  std::vector<RealTensor> & jacobians,
140  std::vector<Real> & jacobian_dets);
141 
146  const MeshQualityInfo & get_mesh_info();
147 
148  /*
149  * Computes information about the mesh quality and sets the _mesh_info attribute.
150  */
152 
153 protected:
154 
155  // System initialization
156  virtual void init_data () override;
157 
158  // Context initialization
159  virtual void init_context (libMesh::DiffContext & context) override;
160 
161  // Element residual and jacobian calculations
162  // Time dependent parts
163  virtual bool element_time_derivative (bool request_jacobian,
164  libMesh::DiffContext & context) override;
165 
166  /* Computes the element reference volume used in the dilation metric
167  * The reference value is set to the averaged value of all elements' average
168  * |J|. Also computes any applicable target element inverse Jacobians. Target
169  * elements are relavant when the reference element does not minimize the
170  * distortion metric.
171  */
172  void prepare_for_smoothing();
173 
178 
184 
189 
194 
195  /* Map to hold target qp-dependent element target-to-reference mapping
196  * Jacobians, if any
197  */
198  std::map<ElemType, std::vector<RealTensor>> _target_jacobians;
199 
200  /*
201  * Map to hold the determinants of _target_jacobians.
202  */
203  std::map<ElemType, std::vector<Real>> _target_jacobian_dets;
204 
209 
214 };
215 
216 } // namespace libMesh
217 
218 #endif // LIBMESH_VARIATIONAL_SMOOTHER_SYSTEM_H
virtual void init_data() override
Initializes the member data fields associated with the system, so that, e.g., assemble() may be used...
std::pair< Real, dof_id_type > min_elem_dilation
virtual void solve() override
Solves the system to smooth the mesh.
ElemType
Defines an enum for geometric element types.
This is the EquationSystems class.
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:55
std::pair< Real, dof_id_type > max_elem_distortion
static constexpr Real TOLERANCE
std::pair< Real, dof_id_type > max_elem_det_S
Real _epsilon_squared_assembly
Epsilon squared value determined at runtime during each assembly.
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
virtual void assembly(bool get_residual, bool get_jacobian, bool apply_heterogeneous_constraints=false, bool apply_no_constraints=false) override
Assembly method to update the mesh based on the smoother solve.
static std::pair< std::unique_ptr< Elem >, std::vector< std::unique_ptr< Node > > > get_target_elem(const ElemType &type)
Get the target element for a given element type.
static void get_target_to_reference_jacobian(const Elem *const target_elem, const FEMContext &femcontext, std::vector< RealTensor > &jacobians, std::vector< Real > &jacobian_dets)
Get the jacobians (and determinants) of the target-to-reference element mapping.
The libMesh namespace provides an interface to certain functionality in the library.
std::pair< Real, dof_id_type > max_elem_combined
This class provides a specific system class.
Definition: fem_system.h:53
Real _dilation_weight
The relative weight to give the dilation metric.
std::pair< Real, dof_id_type > min_elem_distortion
Struct to hold smoother-relevant information about the mesh quality.
unsigned int number() const
Definition: system.h:2350
const MeshQualityInfo & get_mesh_info()
Getter for the _mesh_info attribute.
static const dof_id_type invalid_id
An invalid id to distinguish an uninitialized DofObject.
Definition: dof_object.h:482
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:62
bool _untangling_solve
Flag to indicate if the current solve is to untangle or smooth.
virtual bool element_time_derivative(bool request_jacobian, libMesh::DiffContext &context) override
Adds the time derivative contribution on elem to elem_residual.
std::pair< Real, dof_id_type > min_elem_combined
Real _ref_vol
The reference volume for each element.
const Real _epsilon_squared
The small nonzero constant to prevent zero denominators (degenerate meshes only)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::map< ElemType, std::vector< Real > > _target_jacobian_dets
std::pair< Real, dof_id_type > min_elem_det_S
VariationalSmootherSystem(libMesh::EquationSystems &es, const std::string &name, const unsigned int number)
This is an FEMSystem to solve the optimization probelem posed by the VariationalMeshSmoother class...
virtual void init_context(libMesh::DiffContext &context) override
const std::string & name() const
Definition: system.h:2342
std::map< ElemType, std::vector< RealTensor > > _target_jacobians
MeshQualityInfo _mesh_info
Information about the mesh quality.
std::pair< Real, dof_id_type > max_elem_dilation