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_boundary_temperature_interface_h
23 #define _aspect_boundary_temperature_interface_h
24 
25 #include <aspect/plugins.h>
27 #include <aspect/utilities.h>
29 
30 #include <deal.II/base/parameter_handler.h>
31 #include <deal.II/distributed/tria.h>
32 
33 #include <boost/core/demangle.hpp>
34 #include <typeinfo>
35 
36 
37 namespace aspect
38 {
39  template <int dim> class SimulatorAccess;
40 
47  namespace BoundaryTemperature
48  {
49  using namespace dealii;
50 
56  template <int dim>
57  class Interface
58  {
59  public:
64  virtual ~Interface() = default;
65 
71  virtual void initialize ();
72 
84  virtual
85  double boundary_temperature (const types::boundary_id boundary_indicator,
86  const Point<dim> &position) const;
87 
95  virtual
96  double minimal_temperature (const std::set<types::boundary_id> &fixed_boundary_ids
97  = std::set<types::boundary_id>()) const = 0;
98 
106  virtual
107  double maximal_temperature (const std::set<types::boundary_id> &fixed_boundary_ids
108  = std::set<types::boundary_id>()) const = 0;
109 
121  virtual
122  void
123  update ();
124 
131  static
132  void
133  declare_parameters (ParameterHandler &prm);
134 
141  virtual
142  void
143  parse_parameters (ParameterHandler &prm);
144  };
145 
146 
152  template <int dim>
153  class Manager : public ::aspect::SimulatorAccess<dim>
154  {
155  public:
160  ~Manager () override;
161 
171  virtual
172  void
173  update ();
174 
179  static
180  void
181  declare_parameters (ParameterHandler &prm);
182 
188  void
189  parse_parameters (ParameterHandler &prm);
190 
196  double
197  boundary_temperature (const types::boundary_id boundary_indicator,
198  const Point<dim> &position) const;
199 
204  double minimal_temperature (const std::set<types::boundary_id> &fixed_boundary_ids
205  = std::set<types::boundary_id>()) const;
206 
211  double maximal_temperature (const std::set<types::boundary_id> &fixed_boundary_ids
212  = std::set<types::boundary_id>()) const;
213 
231  static
232  void
233  register_boundary_temperature (const std::string &name,
234  const std::string &description,
235  void (*declare_parameters_function) (ParameterHandler &),
236  std::unique_ptr<Interface<dim>> (*factory_function) ());
237 
238 
243  const std::vector<std::string> &
244  get_active_boundary_temperature_names () const;
245 
250  const std::vector<std::unique_ptr<Interface<dim>>> &
251  get_active_boundary_temperature_conditions () const;
252 
260  template <typename BoundaryTemperatureType>
262  BoundaryTemperatureType *
263  find_boundary_temperature_model () const;
264 
271  template <typename BoundaryTemperatureType>
272  bool
273  has_matching_boundary_temperature_model () const;
274 
283  template <typename BoundaryTemperatureType>
284  const BoundaryTemperatureType &
285  get_matching_boundary_temperature_model () const;
286 
287  /*
288  * Return a set of boundary indicators for which boundary
289  * temperatures are prescribed.
290  */
291  const std::set<types::boundary_id> &
292  get_fixed_temperature_boundary_indicators() const;
293 
294  /*
295  * Return whether Dirichlet boundary conditions will be applied
296  * on parts of the boundaries where material flows out.
297  */
298  bool
299  allows_fixed_temperature_on_outflow_boundaries() const;
300 
310  static
311  void
312  write_plugin_graph (std::ostream &output_stream);
313 
314 
318  DeclException1 (ExcBoundaryTemperatureNameNotFound,
319  std::string,
320  << "Could not find entry <"
321  << arg1
322  << "> among the names of registered boundary temperature objects.");
323  private:
328  std::vector<std::unique_ptr<Interface<dim>>> boundary_temperature_objects;
329 
334  std::vector<std::string> model_names;
335 
342  std::vector<aspect::Utilities::Operator> model_operators;
343 
348  std::set<types::boundary_id> fixed_temperature_boundary_indicators;
349 
355  };
356 
357 
358 
359  template <int dim>
360  template <typename BoundaryTemperatureType>
361  inline
362  BoundaryTemperatureType *
364  {
365  for (const auto &p : boundary_temperature_objects)
366  if (BoundaryTemperatureType *x = dynamic_cast<BoundaryTemperatureType *> ( p.get()) )
367  return x;
368  return nullptr;
369  }
370 
371 
372  template <int dim>
373  template <typename BoundaryTemperatureType>
374  inline
375  bool
377  {
378  for (const auto &p : boundary_temperature_objects)
379  if (Plugins::plugin_type_matches<BoundaryTemperatureType>(*p))
380  return true;
381  return false;
382  }
383 
384 
385  template <int dim>
386  template <typename BoundaryTemperatureType>
387  inline
388  const BoundaryTemperatureType &
390  {
391  AssertThrow(has_matching_boundary_temperature_model<BoundaryTemperatureType> (),
392  ExcMessage("You asked BoundaryTemperature::Manager::get_boundary_temperature_model() for a "
393  "boundary temperature model of type <" + boost::core::demangle(typeid(BoundaryTemperatureType).name()) + "> "
394  "that could not be found in the current model. Activate this "
395  "boundary temperature model in the input file."));
396 
397  for (auto &p : boundary_temperature_objects)
398  if (Plugins::plugin_type_matches<BoundaryTemperatureType>(*p))
399  return Plugins::get_plugin_as_type<BoundaryTemperatureType>(*p);
400 
401  // We will never get here, because we had the Assert above. Just to avoid warnings.
402  typename std::vector<std::unique_ptr<Interface<dim>>>::const_iterator boundary_temperature_model;
403  return Plugins::get_plugin_as_type<BoundaryTemperatureType>(*(*boundary_temperature_model));
404  }
405 
406 
413  template <int dim>
414  std::string
416 
417 
425 #define ASPECT_REGISTER_BOUNDARY_TEMPERATURE_MODEL(classname, name, description) \
426  template class classname<2>; \
427  template class classname<3>; \
428  namespace ASPECT_REGISTER_BOUNDARY_TEMPERATURE_MODEL_ ## classname \
429  { \
430  aspect::internal::Plugins::RegisterHelper<aspect::BoundaryTemperature::Interface<2>,classname<2>> \
431  dummy_ ## classname ## _2d (&aspect::BoundaryTemperature::Manager<2>::register_boundary_temperature, \
432  name, description); \
433  aspect::internal::Plugins::RegisterHelper<aspect::BoundaryTemperature::Interface<3>,classname<3>> \
434  dummy_ ## classname ## _3d (&aspect::BoundaryTemperature::Manager<3>::register_boundary_temperature, \
435  name, description); \
436  }
437  }
438 }
439 
440 
441 #endif
const BoundaryTemperatureType & get_matching_boundary_temperature_model() const
Definition: interface.h:389
std::vector< std::unique_ptr< Interface< dim > > > boundary_temperature_objects
Definition: interface.h:328
StructuredDataLookup< dim > DEAL_II_DEPRECATED
std::string get_valid_model_names_pattern()
void write_plugin_graph(std::ostream &output_stream)
std::vector< std::string > model_names
Definition: interface.h:334
std::vector< aspect::Utilities::Operator > model_operators
Definition: interface.h:342
void declare_parameters(ParameterHandler &prm)
DEAL_II_DEPRECATED BoundaryTemperatureType * find_boundary_temperature_model() const
Definition: compat.h:88
std::set< types::boundary_id > fixed_temperature_boundary_indicators
Definition: interface.h:348
bool has_matching_boundary_temperature_model() const
Definition: interface.h:376
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.")