ASPECT
interface.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2011 - 2022 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>
74  class Interface
75  {
76  public:
81  virtual ~Interface () = default;
82 
88  virtual void initialize ();
89 
99  virtual
100  void
101  update ();
102 
112  virtual
113  void
114  execute (Vector<float> &error_indicators) const;
115 
126  virtual
127  void
128  tag_additional_cells () const;
129 
143  static
144  void
145  declare_parameters (ParameterHandler &prm);
146 
153  virtual
154  void
155  parse_parameters (ParameterHandler &prm);
156  };
157 
158 
159 
160 
161 
162 
169  template <int dim>
170  class Manager : public ::aspect::SimulatorAccess<dim>
171  {
172  public:
177  ~Manager () override;
178 
185  virtual
186  void
187  update ();
188 
195  virtual
196  void
197  execute (Vector<float> &error_indicators) const;
198 
204  virtual
205  void
206  tag_additional_cells () const;
207 
212  static
213  void
214  declare_parameters (ParameterHandler &prm);
215 
221  void
222  parse_parameters (ParameterHandler &prm);
223 
230  template <typename MeshRefinementType>
231  bool
232  has_matching_mesh_refinement_strategy () const;
233 
242  template <typename MeshRefinementType>
243  const MeshRefinementType &
244  get_matching_mesh_refinement_strategy () const;
245 
263  static
264  void
265  register_mesh_refinement_criterion (const std::string &name,
266  const std::string &description,
267  void (*declare_parameters_function) (ParameterHandler &),
268  std::unique_ptr<Interface<dim>> (*factory_function) ());
269 
279  static
280  void
281  write_plugin_graph (std::ostream &output_stream);
282 
286  DeclException1 (ExcMeshRefinementNameNotFound,
287  std::string,
288  << "Could not find entry <"
289  << arg1
290  << "> among the names of registered mesh refinement objects.");
291  private:
297  { plus, max };
298 
303 
309 
314  std::vector<double> scaling_factors;
315 
320  std::list<std::unique_ptr<Interface<dim>>> mesh_refinement_objects;
321  };
322 
323 
324 
325  template <int dim>
326  template <typename MeshRefinementType>
327  inline
328  bool
330  {
331  for (typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator
332  p = mesh_refinement_objects.begin();
333  p != mesh_refinement_objects.end(); ++p)
334  if (Plugins::plugin_type_matches<MeshRefinementType>(*(*p)))
335  return true;
336 
337  return false;
338  }
339 
340 
341 
342  template <int dim>
343  template <typename MeshRefinementType>
344  inline
345  const MeshRefinementType &
347  {
348  AssertThrow(has_matching_mesh_refinement_strategy<MeshRefinementType> (),
349  ExcMessage("You asked MeshRefinement::Manager::get_matching_mesh_refinement_strategy() for a "
350  "mesh refinement strategy of type <" + boost::core::demangle(typeid(MeshRefinementType).name()) + "> "
351  "that could not be found in the current model. Activate this "
352  "mesh refinement strategy in the input file."));
353 
354  for (typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator
355  p = mesh_refinement_objects.begin();
356  p != mesh_refinement_objects.end(); ++p)
357  if (Plugins::plugin_type_matches<MeshRefinementType>(*(*p)))
358  return Plugins::get_plugin_as_type<MeshRefinementType>(*(*p));
359 
360  // We will never get here, because we had the Assert above. Just to avoid warnings.
361  typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator mesh_refinement_strategy;
362  return Plugins::get_plugin_as_type<MeshRefinementType>(*(*mesh_refinement_strategy));
363  }
364 
365 
366 
374 #define ASPECT_REGISTER_MESH_REFINEMENT_CRITERION(classname,name,description) \
375  template class classname<2>; \
376  template class classname<3>; \
377  namespace ASPECT_REGISTER_MESH_REFINEMENT_CRITERION_ ## classname \
378  { \
379  aspect::internal::Plugins::RegisterHelper<aspect::MeshRefinement::Interface<2>,classname<2>> \
380  dummy_ ## classname ## _2d (&aspect::MeshRefinement::Manager<2>::register_mesh_refinement_criterion, \
381  name, description); \
382  aspect::internal::Plugins::RegisterHelper<aspect::MeshRefinement::Interface<3>,classname<3>> \
383  dummy_ ## classname ## _3d (&aspect::MeshRefinement::Manager<3>::register_mesh_refinement_criterion, \
384  name, description); \
385  }
386  }
387 }
388 
389 
390 #endif
MergeOperation merge_operation
Definition: interface.h:302
void write_plugin_graph(std::ostream &output_stream)
std::vector< double > scaling_factors
Definition: interface.h:314
void declare_parameters(ParameterHandler &prm)
Definition: compat.h:88
std::list< std::unique_ptr< Interface< dim > > > mesh_refinement_objects
Definition: interface.h:320
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.")