22 #ifndef _aspect_mesh_refinement_interface_h 23 #define _aspect_mesh_refinement_interface_h 30 #include <deal.II/base/table_handler.h> 31 #include <deal.II/base/parameter_handler.h> 32 #include <deal.II/distributed/tria.h> 39 template <
int dim>
class Simulator;
40 template <
int dim>
class SimulatorAccess;
49 namespace MeshRefinement
114 execute (Vector<float> &error_indicators)
const;
128 tag_additional_cells ()
const;
155 parse_parameters (ParameterHandler &prm);
197 execute (Vector<float> &error_indicators)
const;
206 tag_additional_cells ()
const;
222 parse_parameters (ParameterHandler &prm);
233 template <
typename MeshRefinementType,
234 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,MeshRefinementType>::value>>
236 has_matching_mesh_refinement_strategy ()
const;
249 template <
typename MeshRefinementType,
250 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,MeshRefinementType>::value>>
251 const MeshRefinementType &
252 get_matching_mesh_refinement_strategy ()
const;
273 register_mesh_refinement_criterion (
const std::string &name,
274 const std::string &description,
275 void (*declare_parameters_function) (ParameterHandler &),
296 <<
"Could not find entry <" 298 <<
"> among the names of registered mesh refinement objects.");
334 template <
typename MeshRefinementType,
typename>
339 for (
typename std::list<std::unique_ptr<
Interface<dim>>>::const_iterator
340 p = mesh_refinement_objects.begin();
341 p != mesh_refinement_objects.end(); ++p)
342 if (Plugins::plugin_type_matches<MeshRefinementType>(*(*p)))
351 template <
typename MeshRefinementType,
typename>
353 const MeshRefinementType &
356 AssertThrow(has_matching_mesh_refinement_strategy<MeshRefinementType> (),
357 ExcMessage(
"You asked MeshRefinement::Manager::get_matching_mesh_refinement_strategy() for a " 358 "mesh refinement strategy of type <" + boost::core::demangle(
typeid(MeshRefinementType).name()) +
"> " 359 "that could not be found in the current model. Activate this " 360 "mesh refinement strategy in the input file."));
362 for (
typename std::list<std::unique_ptr<
Interface<dim>>>::const_iterator
363 p = mesh_refinement_objects.begin();
364 p != mesh_refinement_objects.end(); ++p)
365 if (Plugins::plugin_type_matches<MeshRefinementType>(*(*p)))
366 return Plugins::get_plugin_as_type<MeshRefinementType>(*(*p));
369 typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator mesh_refinement_strategy;
370 return Plugins::get_plugin_as_type<MeshRefinementType>(*(*mesh_refinement_strategy));
382 #define ASPECT_REGISTER_MESH_REFINEMENT_CRITERION(classname,name,description) \ 383 template class classname<2>; \ 384 template class classname<3>; \ 385 namespace ASPECT_REGISTER_MESH_REFINEMENT_CRITERION_ ## classname \ 387 aspect::internal::Plugins::RegisterHelper<aspect::MeshRefinement::Interface<2>,classname<2>> \ 388 dummy_ ## classname ## _2d (&aspect::MeshRefinement::Manager<2>::register_mesh_refinement_criterion, \ 389 name, description); \ 390 aspect::internal::Plugins::RegisterHelper<aspect::MeshRefinement::Interface<3>,classname<3>> \ 391 dummy_ ## classname ## _3d (&aspect::MeshRefinement::Manager<3>::register_mesh_refinement_criterion, \ 392 name, description); \ MergeOperation merge_operation
void write_plugin_graph(std::ostream &output_stream)
std::vector< double > scaling_factors
void declare_parameters(ParameterHandler &prm)
std::list< std::unique_ptr< Interface< dim > > > mesh_refinement_objects
DeclException1(ProbabilityFunctionNegative, Point< dim >,<< "Your probability density function in the particle generator " "returned a negative probability density for the following position: "<< arg1<< ". Please check your function expression.")