22 #ifndef _aspect_boundary_velocity_interface_h 23 #define _aspect_boundary_velocity_interface_h 30 #include <deal.II/base/point.h> 31 #include <deal.II/base/parameter_handler.h> 33 #include <boost/core/demangle.hpp> 45 namespace BoundaryVelocity
62 virtual ~Interface() =
default;
102 boundary_velocity (
const types::boundary_id boundary_indicator,
103 const Point<dim> &position)
const = 0;
123 parse_parameters (ParameterHandler &prm);
139 ~Manager ()
override;
160 boundary_velocity (
const types::boundary_id boundary_indicator,
161 const Point<dim> &position)
const;
182 register_boundary_velocity (
const std::string &name,
183 const std::string &description,
184 void (*declare_parameters_function) (ParameterHandler &),
185 std::unique_ptr<Interface<dim>> (*factory_function) ());
200 const std::map<types::boundary_id, std::pair<std::string,std::vector<std::string>>> &
201 get_active_boundary_velocity_names ()
const;
212 const std::map<types::boundary_id,std::vector<std::unique_ptr<BoundaryVelocity::Interface<dim>>>> &
213 get_active_boundary_velocity_conditions ()
const;
219 const std::set<types::boundary_id> &
220 get_zero_boundary_velocity_indicators ()
const;
226 const std::set<types::boundary_id> &
227 get_tangential_boundary_velocity_indicators ()
const;
243 parse_parameters (ParameterHandler &prm);
254 template <
typename BoundaryVelocityType,
255 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,BoundaryVelocityType>::value>>
257 has_matching_boundary_velocity_model ()
const;
270 template <
typename BoundaryVelocityType,
271 typename =
typename std::enable_if_t<std::is_base_of<Interface<dim>,BoundaryVelocityType>::value>>
272 const BoundaryVelocityType &
273 get_matching_boundary_velocity_model ()
const;
294 <<
"Could not find entry <" 296 <<
"> among the names of registered boundary velocity objects.");
330 template <
typename BoundaryVelocityType,
typename>
335 for (
const auto &boundary : boundary_velocity_objects)
336 for (
const auto &p : boundary.second)
337 if (Plugins::plugin_type_matches<BoundaryVelocityType>(*p))
344 template <
typename BoundaryVelocityType,
typename>
346 const BoundaryVelocityType &
349 AssertThrow(has_matching_boundary_velocity_model<BoundaryVelocityType> (),
350 ExcMessage(
"You asked BoundaryVelocity::Manager::get_boundary_velocity_model() for a " 351 "boundary velocity model of type <" + boost::core::demangle(
typeid(BoundaryVelocityType).name()) +
"> " 352 "that could not be found in the current model. Activate this " 353 "boundary velocity model in the input file."));
355 typename std::map<types::boundary_id,std::vector<std::unique_ptr<BoundaryVelocity::Interface<dim>>>>::const_iterator boundary_velocity_model;
356 for (
const auto &boundary : boundary_velocity_objects)
357 for (
const auto &p : boundary)
358 if (Plugins::plugin_type_matches<BoundaryVelocityType>(*p))
359 return Plugins::get_plugin_as_type<BoundaryVelocityType>(*p);
362 return Plugins::get_plugin_as_type<BoundaryVelocityType>(*(*boundary_velocity_model));
384 #define ASPECT_REGISTER_BOUNDARY_VELOCITY_MODEL(classname, name, description) \ 385 template class classname<2>; \ 386 template class classname<3>; \ 387 namespace ASPECT_REGISTER_BOUNDARY_VELOCITY_MODEL_ ## classname \ 389 aspect::internal::Plugins::RegisterHelper<aspect::BoundaryVelocity::Interface<2>,classname<2>> \ 390 dummy_ ## classname ## _2d (&aspect::BoundaryVelocity::Manager<2>::register_boundary_velocity, \ 391 name, description); \ 392 aspect::internal::Plugins::RegisterHelper<aspect::BoundaryVelocity::Interface<3>,classname<3>> \ 393 dummy_ ## classname ## _3d (&aspect::BoundaryVelocity::Manager<3>::register_boundary_velocity, \ 394 name, description); \
void write_plugin_graph(std::ostream &output_stream)
std::set< types::boundary_id > tangential_velocity_boundary_indicators
std::map< types::boundary_id, std::pair< std::string, std::vector< std::string > > > boundary_velocity_indicators
std::string get_valid_model_names_pattern()
std::map< types::boundary_id, std::vector< std::unique_ptr< BoundaryVelocity::Interface< dim > > > > boundary_velocity_objects
bool has_matching_boundary_velocity_model() const
const BoundaryVelocityType & get_matching_boundary_velocity_model() const
void declare_parameters(ParameterHandler &prm)
std::set< types::boundary_id > zero_velocity_boundary_indicators
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.")