22 #ifndef _aspect_initial_composition_interface_h 23 #define _aspect_initial_composition_interface_h 29 #include <deal.II/base/point.h> 30 #include <deal.II/base/parameter_handler.h> 32 #include <boost/core/demangle.hpp> 38 template <
int dim>
class SimulatorAccess;
46 namespace InitialComposition
63 virtual ~Interface() =
default;
78 double initial_composition (
const Point<dim> &position,
const unsigned int n_comp)
const = 0;
99 parse_parameters (ParameterHandler &prm);
118 ~Manager ()
override;
134 parse_parameters (ParameterHandler &prm);
143 initial_composition (
const Point<dim> &position,
144 const unsigned int n_comp)
const;
165 register_initial_composition (
const std::string &name,
166 const std::string &description,
167 void (*declare_parameters_function) (ParameterHandler &),
168 std::unique_ptr<Interface<dim>> (*factory_function) ());
175 const std::vector<std::string> &
176 get_active_initial_composition_names ()
const;
182 const std::list<std::unique_ptr<Interface<dim>>> &
183 get_active_initial_composition_conditions ()
const;
194 template <
typename InitialCompositionType,
195 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,InitialCompositionType>::value>>
197 has_matching_initial_composition_model ()
const;
210 template <
typename InitialCompositionType,
211 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,InitialCompositionType>::value>>
212 const InitialCompositionType &
213 get_matching_initial_composition_model ()
const;
233 <<
"Could not find entry <" 235 <<
"> among the names of registered initial composition objects.");
261 template <
typename InitialCompositionType,
typename>
266 for (
const auto &p : initial_composition_objects)
267 if (Plugins::plugin_type_matches<InitialCompositionType>(*p))
274 template <
typename InitialCompositionType,
typename>
276 const InitialCompositionType &
279 AssertThrow(has_matching_initial_composition_model<InitialCompositionType> (),
280 ExcMessage(
"You asked InitialComposition::Manager::get_initial_composition_model() for a " 281 "initial composition model of type <" + boost::core::demangle(
typeid(InitialCompositionType).name()) +
"> " 282 "that could not be found in the current model. Activate this " 283 "initial composition model in the input file."));
285 typename std::list<std::unique_ptr<Interface<dim>>>::const_iterator initial_composition_model;
286 for (
const auto &p : initial_composition_objects)
287 if (Plugins::plugin_type_matches<InitialCompositionType>(*p))
288 return Plugins::get_plugin_as_type<InitialCompositionType>(*p);
291 return Plugins::get_plugin_as_type<InitialCompositionType>(*(*initial_composition_model));
316 #define ASPECT_REGISTER_INITIAL_COMPOSITION_MODEL(classname,name,description) \ 317 template class classname<2>; \ 318 template class classname<3>; \ 319 namespace ASPECT_REGISTER_INITIAL_COMPOSITION_MODEL_ ## classname \ 321 aspect::internal::Plugins::RegisterHelper<aspect::InitialComposition::Interface<2>,classname<2>> \ 322 dummy_ ## classname ## _2d (&aspect::InitialComposition::Manager<2>::register_initial_composition, \ 323 name, description); \ 324 aspect::internal::Plugins::RegisterHelper<aspect::InitialComposition::Interface<3>,classname<3>> \ 325 dummy_ ## classname ## _3d (&aspect::InitialComposition::Manager<3>::register_initial_composition, \ 326 name, description); \
std::list< std::unique_ptr< Interface< dim > > > initial_composition_objects
void write_plugin_graph(std::ostream &output_stream)
const InitialCompositionType & get_matching_initial_composition_model() const
std::string get_valid_model_names_pattern()
void declare_parameters(ParameterHandler &prm)
std::vector< std::string > model_names
bool has_matching_initial_composition_model() const
std::vector< aspect::Utilities::Operator > model_operators
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.")