21 #include "libmesh/equation_systems.h" 22 #include "libmesh/kelly_error_estimator.h" 23 #include "libmesh/patch_recovery_error_estimator.h" 24 #include "libmesh/fourth_error_estimators.h" 25 #include "libmesh/parallel.h" 26 #include "libmesh/error_vector.h" 27 #include "libmesh/distributed_mesh.h" 31 #ifdef LIBMESH_ENABLE_AMR 37 _fe_problem(fe_problem),
38 _mesh(_fe_problem.
mesh()),
39 _mesh_refinement_on(false),
43 _print_mesh_changed(false),
44 _t(_fe_problem.time()),
45 _step(_fe_problem.timeStep()),
47 _start_time(-
std::numeric_limits<
Real>::
max()),
50 _use_new_system(false),
52 _recompute_markers_during_cycles(false)
60 const unsigned int initial_steps,
61 const bool p_refinement)
69 _error = std::make_unique<ErrorVector>();
112 if (error_estimator_name ==
"KellyErrorEstimator")
114 else if (error_estimator_name ==
"LaplacianErrorEstimator")
116 else if (error_estimator_name ==
"PatchRecoveryErrorEstimator")
119 mooseError(std::string(
"Unknown error_estimator selection: ") +
120 std::string(error_estimator_name));
126 mooseAssert(
_error_estimator,
"error_estimator not initialized. Did you call init_adaptivity()?");
133 TIME_SECTION(
"adaptMesh", 3,
"Adapting Mesh");
136 if (marker_name.empty())
139 bool mesh_changed =
false;
142 bool distributed_adaptivity =
false;
146 if (!marker_name.empty())
150 std::vector<Number> serialized_solution;
155 std::unique_ptr<ConstElemRange> all_elems;
160 if (distributed_mesh && !distributed_mesh->is_serial_on_zero())
164 distributed_adaptivity =
true;
176 all_elems = std::make_unique<ConstElemRange>(
184 distributed_adaptivity =
false;
199 Threads::parallel_reduce(*all_elems, fet);
229 if (distributed_adaptivity)
242 if (distributed_adaptivity)
249 bool displaced_mesh_changed =
254 mooseAssert(displaced_mesh_changed,
"Undisplaced mesh changed, but displaced mesh did not!");
276 mooseAssert(
mesh,
"Mesh pointer must not be NULL");
282 level =
mesh->uniformRefineLevel();
290 if (
mesh->skipDeletionRepartitionAfterRefine())
294 mesh->needsRemoteElemDeletion(
false);
303 TIME_SECTION(
"uniformRefineWithProjection", 2,
"Uniformly Refining and Reprojecting");
312 for (
unsigned int i = 0; i < level; i++)
321 displaced_mesh_refinement.uniformly_refine(1);
332 mooseError(
"Mesh adaptivity system not available");
368 return *insert_pair.first->second;
374 TIME_SECTION(
"updateErrorVectors", 5,
"Updating Error Vectors");
398 #endif // LIBMESH_ENABLE_AMR bool adaptMesh(std::string marker_name=std::string())
Adapts the mesh based on the error estimator used.
bool initialAdaptMesh()
Used during initial adaptivity.
virtual void meshChanged(bool intermediate_change, bool contract_mesh, bool clean_refinement_flags)
Update data after a mesh change.
std::unique_ptr< libMesh::ErrorVector > _error
Error vector for use with the error estimator.
libMesh::ConstElemRange * getActiveLocalElementRange()
Return pointers to range objects for various types of ranges (local nodes, boundary elems...
std::shared_ptr< DisplacedProblem > _displaced_problem
Real _stop_time
When adaptivity stops.
std::string _marker_variable_name
Name of the marker variable if using the new adaptivity system.
const unsigned int invalid_uint
bool _initialized
on/off flag reporting if the adaptivity system has been initialized
NumericVector< Number > & solution()
static void uniformRefine(MooseMesh *mesh, unsigned int level=libMesh::invalid_uint)
Performs uniform refinement of the passed Mesh object.
void skip_partitioning(bool skip)
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
std::unique_ptr< libMesh::ErrorEstimator > _error_estimator
Error estimator to be used by the apps.
void updateErrorVectors()
Update the ErrorVectors that have been requested through calls to getErrorVector().
std::unique_ptr< libMesh::MeshRefinement > _mesh_refinement
A mesh refinement object to be used either with initial refinement or with Adaptivity.
void setTimeActive(Real start_time, Real stop_time)
Sets the time when the adaptivity is active.
void init(const unsigned int steps, const unsigned int initial_steps, const bool p_refinement)
Initialize and turn on adaptivity for the simulation.
const Parallel::Communicator & comm() const
void setUseNewSystem()
Tells this object we're using the "new" adaptivity system.
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
std::pair< typename M::iterator, bool > moose_try_emplace(M &m, const typename M::key_type &k, Args &&... args)
Function to mirror the behavior of the C++17 std::map::try_emplace() method (no hint).
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
std::string _initial_marker_variable_name
Name of the initial marker variable if using the new adaptivity system.
auto max(const L &left, const R &right)
void update()
Update the system (doing libMesh magic)
bool _print_mesh_changed
True if we want to print out info when mesh has changed.
virtual libMesh::DofMap & dofMap()
Gets writeable reference to the dof map.
void allow_remote_element_removal(bool allow)
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
PeriodicBoundaries * get_periodic_boundaries()
bool _mesh_refinement_on
on/off flag reporting if the adaptivity is being used
Adaptivity(FEProblemBase &fe_problem)
virtual libMesh::EquationSystems & es() override
An inteface for the _console for outputting to the Console object.
virtual dof_id_type max_elem_id() const=0
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
void setAdaptivityOn(bool state)
Allow adaptivity to be toggled programatically.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
std::unique_ptr< libMesh::MeshRefinement > _displaced_mesh_refinement
A mesh refinement object for displaced mesh.
NonlinearSystemBase & getNonlinearSystemBase(const unsigned int sys_num)
unsigned int uniformRefineLevel() const
Returns the level of uniform refinement requested (zero if AMR is disabled).
void setErrorEstimator(const MooseEnum &error_estimator_name)
Set the error estimator.
AuxiliarySystem & getAuxiliarySystem()
unsigned int _steps
steps of adaptivity to perform
Interface for objects interacting with the PerfGraph.
void uniformRefineWithProjection()
Performs uniform refinement on the meshes in the current object.
void printInfo(std::ostream &os=libMesh::out, const unsigned int verbosity=0) const
Calls print_info() on the underlying Mesh.
unsigned int _initial_steps
the number of adaptivity steps to do at the beginning of simulation
bool _use_new_system
Whether or not to use the "new" adaptivity system.
std::map< std::string, std::unique_ptr< libMesh::ErrorVector > > _indicator_field_to_error_vector
Stores pointers to ErrorVectors associated with indicator field names.
unsigned int _interval
intreval between adaptivity runs
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual std::shared_ptr< const DisplacedProblem > getDisplacedProblem() const
libMesh::ErrorVector & getErrorVector(const std::string &indicator_field)
Get an ErrorVector that will be filled up with values corresponding to the indicator field name passe...
T & set(const std::string &)
virtual MooseMesh & mesh() override
Real _start_time
When adaptivity start.
bool isAdaptivityDue()
Query if an adaptivity step should be performed at the current time / time step.
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
void setMarkerVariableName(std::string marker_field)
Sets the name of the field variable to actually use to flag elements for refinement / coarsening...
void doingPRefinement(bool doing_p_refinement)
Indicate whether the kind of adaptivity we're doing is p-refinement.
unsigned int _max_h_level
The maximum number of refinement levels.
void setErrorNorm(libMesh::SystemNorm &sys_norm)
Set the error norm (FIXME: improve description)
void setInitialMarkerVariableName(std::string marker_field)
Sets the name of the field variable to actually use to flag elements for initial refinement / coarsen...
FEProblemBase & _fe_problem
void uniformly_refine(unsigned int n=1)
virtual void localize(std::vector< T > &v_local) const=0
virtual libMesh::System & system() override
Get the reference to the libMesh system.