22 #ifndef _aspect_heating_model_interface_h 23 #define _aspect_heating_model_interface_h 30 #include <deal.II/base/point.h> 31 #include <deal.II/base/parameter_handler.h> 33 #include <boost/core/demangle.hpp> 39 template <
int dim>
class SimulatorAccess;
46 namespace HeatingModel
69 const unsigned int n_comp);
191 parse_parameters (ParameterHandler &prm);
237 adiabatic_heating_enabled()
const;
244 shear_heating_enabled()
const;
261 parse_parameters (ParameterHandler &prm);
313 register_heating_model (
const std::string &name,
314 const std::string &description,
315 void (*declare_parameters_function) (ParameterHandler &),
323 const std::vector<std::string> &
324 get_active_heating_model_names ()
const;
330 const std::list<std::unique_ptr<Interface<dim>>> &
331 get_active_heating_models ()
const;
342 template <
typename HeatingModelType,
343 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,HeatingModelType>::value>>
345 has_matching_heating_model ()
const;
358 template <
typename HeatingModelType,
359 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,HeatingModelType>::value>>
360 const HeatingModelType &
361 get_matching_heating_model ()
const;
382 <<
"Could not find entry <" 384 <<
"> among the names of registered heating model objects.");
402 template <
typename HeatingModelType,
typename>
407 for (
const auto &p : heating_model_objects)
408 if (Plugins::plugin_type_matches<HeatingModelType>(*p))
415 template <
typename HeatingModelType,
typename>
417 const HeatingModelType &
420 AssertThrow(has_matching_heating_model<HeatingModelType> (),
421 ExcMessage(
"You asked HeatingModel::Manager::get_heating_model() for a " 422 "heating model of type <" + boost::core::demangle(
typeid(HeatingModelType).name()) +
"> " 423 "that could not be found in the current model. Activate this " 424 "heating model in the input file."));
426 typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator heating_model;
427 for (
typename std::list<std::unique_ptr<
Interface<dim>>>::const_iterator
428 p = heating_model_objects.begin();
429 p != heating_model_objects.end(); ++p)
430 if (Plugins::plugin_type_matches<HeatingModelType>(*(*p)))
431 return Plugins::get_plugin_as_type<HeatingModelType>(*(*p));
434 return Plugins::get_plugin_as_type<HeatingModelType>(*(*heating_model));
456 #define ASPECT_REGISTER_HEATING_MODEL(classname,name,description) \ 457 template class classname<2>; \ 458 template class classname<3>; \ 459 namespace ASPECT_REGISTER_HEATING_MODEL_ ## classname \ 461 aspect::internal::Plugins::RegisterHelper<aspect::HeatingModel::Interface<2>,classname<2>> \ 462 dummy_ ## classname ## _2d (&aspect::HeatingModel::Manager<2>::register_heating_model, \ 463 name, description); \ 464 aspect::internal::Plugins::RegisterHelper<aspect::HeatingModel::Interface<3>,classname<3>> \ 465 dummy_ ## classname ## _3d (&aspect::HeatingModel::Manager<3>::register_heating_model, \ 466 name, description); \ std::list< std::unique_ptr< Interface< dim > > > heating_model_objects
std::vector< double > heating_source_terms
std::vector< double > rates_of_temperature_change
bool has_matching_heating_model() const
void write_plugin_graph(std::ostream &output_stream)
const HeatingModelType & get_matching_heating_model() const
std::string get_valid_model_names_pattern()
std::vector< double > lhs_latent_heat_terms
void declare_parameters(ParameterHandler &prm)
std::vector< std::string > model_names
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.")