ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2024 by the authors of the ASPECT code.
3 
4  This file is part of ASPECT.
5 
6  ASPECT is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2, or (at your option)
9  any later version.
10 
11  ASPECT is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with ASPECT; see the file LICENSE. If not see
18  <http://www.gnu.org/licenses/>.
19 */
20 
21 
22 #ifndef _aspect_mesh_refinement_interface_h
23 #define _aspect_mesh_refinement_interface_h
24 
25 #include <aspect/global.h>
26 #include <aspect/plugins.h>
28 
29 #include <memory>
30 #include <deal.II/base/table_handler.h>
31 #include <deal.II/base/parameter_handler.h>
32 #include <deal.II/distributed/tria.h>
33 
34 
35 namespace aspect
36 {
37  using namespace dealii;
38 
39  template <int dim> class Simulator;
40  template <int dim> class SimulatorAccess;
41 
42 
49  namespace MeshRefinement
50  {
51 
73  template <int dim>
75  {
76  public:
86  virtual
87  void
88  execute (Vector<float> &error_indicators) const;
89 
100  virtual
101  void
102  tag_additional_cells () const;
103  };
104 
105 
106 
107 
108 
109 
116  template <int dim>
117  class Manager : public Plugins::ManagerBase<Interface<dim>>, public SimulatorAccess<dim>
118  {
119  public:
126  virtual
127  void
128  execute (Vector<float> &error_indicators) const;
129 
135  virtual
136  void
137  tag_additional_cells () const;
138 
143  static
144  void
145  declare_parameters (ParameterHandler &prm);
146 
152  void
153  parse_parameters (ParameterHandler &prm) override;
154 
169  template <typename MeshRefinementType,
170  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,MeshRefinementType>::value>>
171  DEAL_II_DEPRECATED
172  bool
173  has_matching_mesh_refinement_strategy () const;
174 
191  template <typename MeshRefinementType,
192  typename = typename std::enable_if_t<std::is_base_of<Interface<dim>,MeshRefinementType>::value>>
193  DEAL_II_DEPRECATED
194  const MeshRefinementType &
195  get_matching_mesh_refinement_strategy () const;
196 
214  static
215  void
216  register_mesh_refinement_criterion (const std::string &name,
217  const std::string &description,
218  void (*declare_parameters_function) (ParameterHandler &),
219  std::unique_ptr<Interface<dim>> (*factory_function) ());
220 
230  static
231  void
232  write_plugin_graph (std::ostream &output_stream);
233 
237  DeclException1 (ExcMeshRefinementNameNotFound,
238  std::string,
239  << "Could not find entry <"
240  << arg1
241  << "> among the names of registered mesh refinement objects.");
242  private:
248  { plus, max };
249 
254 
260 
265  std::vector<double> scaling_factors;
266  };
267 
268 
269 
270  template <int dim>
271  template <typename MeshRefinementType, typename>
272  inline
273  bool
275  {
276  return this->template has_matching_active_plugin<MeshRefinementType>();
277  }
278 
279 
280 
281  template <int dim>
282  template <typename MeshRefinementType, typename>
283  inline
284  const MeshRefinementType &
286  {
287  return this->template get_matching_active_plugin<MeshRefinementType>();
288  }
289 
290 
291 
299 #define ASPECT_REGISTER_MESH_REFINEMENT_CRITERION(classname,name,description) \
300  template class classname<2>; \
301  template class classname<3>; \
302  namespace ASPECT_REGISTER_MESH_REFINEMENT_CRITERION_ ## classname \
303  { \
304  aspect::internal::Plugins::RegisterHelper<aspect::MeshRefinement::Interface<2>,classname<2>> \
305  dummy_ ## classname ## _2d (&aspect::MeshRefinement::Manager<2>::register_mesh_refinement_criterion, \
306  name, description); \
307  aspect::internal::Plugins::RegisterHelper<aspect::MeshRefinement::Interface<3>,classname<3>> \
308  dummy_ ## classname ## _3d (&aspect::MeshRefinement::Manager<3>::register_mesh_refinement_criterion, \
309  name, description); \
310  }
311  }
312 }
313 
314 
315 #endif
MergeOperation merge_operation
Definition: interface.h:253
void write_plugin_graph(std::ostream &output_stream)
std::vector< double > scaling_factors
Definition: interface.h:265
Definition: compat.h:59
void declare_parameters(ParameterHandler &prm)
Definition: compat.h:42
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.")