Class for handling inclusions in an immersed boundary method. More...
#include <inclusions.h>
Public Member Functions | |
Inclusions (const unsigned int n_vector_components=1) | |
Class for computing the inclusions of a given mesh. | |
types::global_dof_index | n_dofs () const |
types::global_dof_index | n_particles () const |
types::global_dof_index | n_inclusions () const |
Returns the number of inclusions in the mesh. | |
unsigned int | n_dofs_per_inclusion () const |
Number of degrees of freedom associated to each inclusion. | |
void | initialize () |
std::vector< types::global_dof_index > | get_dof_indices (const types::global_dof_index &quadrature_id) const |
Returns the degrees of freedom indices associated with a given quadrature point. | |
void | setup_inclusions_particles (const parallel::distributed::Triangulation< spacedim > &tria) |
Sets up the inclusions particles for the given triangulation. | |
types::global_dof_index | get_inclusion_id (const types::global_dof_index &quadrature_id) const |
unsigned int | get_component (const types::global_dof_index &dof_index) const |
Get the ith component for the given dof index. | |
unsigned int | get_fourier_component (const types::global_dof_index &dof_index) const |
Get the ith Fourier component for the given dof index. | |
double | get_inclusion_data (const types::global_dof_index &inclusion_id, const types::global_dof_index &dof_index) const |
Get the Fourier data for the given local dof index. | |
double | get_inclusion_data (const types::global_dof_index &inclusion_id, const types::global_dof_index &dof_index, const Point< spacedim > &point) const |
Get the Fourier data for the given local dof index. | |
const Tensor< 1, spacedim > & | get_normal (const types::global_dof_index &quadrature_id) const |
Get the normal. | |
double | get_JxW (const types::global_dof_index &particle_id) const |
return weight for integration normal direction of an inclusion at a quadrature point | |
const std::vector< double > & | get_fe_values (const types::global_dof_index particle_id) const |
Point< spacedim > | get_center (const types::global_dof_index &inclusion_id) const |
Get the center of the inclusion. | |
double | get_radius (const types::global_dof_index &inclusion_id) const |
Get the radius of the inclusion. | |
double | get_section_measure (const types::global_dof_index &inclusion_id) const |
Get the measure of the section of the inclusion. | |
Tensor< 1, spacedim > | get_direction (const types::global_dof_index &inclusion_id) const |
Get the direction of the inclusion. | |
Tensor< 2, spacedim > | get_rotation (const types::global_dof_index &inclusion_id) const |
Get the rotation of the inclusion. | |
const std::vector< Point< spacedim > > & | get_current_support_points (const types::global_dof_index &inclusion_id) const |
void | output_particles (const std::string &filename) const |
print the inclusions in parallel on a .vtu file | |
void | update_displacement_hdf5 () |
Update the displacement data after the initialization reading from a hdf5 file. | |
void | update_inclusions_data (std::vector< double > new_data) |
Update the displacement data after the initialization, 3D case only. | |
void | update_inclusions_data (std::vector< std::vector< double > > new_data) |
int | get_vesselID (const types::global_dof_index &inclusion_id) const |
3D return the vessel that a given inclusion belongs to, 2D return 0 | |
void | update_single_inclusion_data_along_normal (const types::global_dof_index &inclusion_id, const double nd) |
void | update_single_vessel_data (const types::global_dof_index &, const std::vector< double >) |
unsigned int | get_n_vessels () const |
unsigned int | get_n_coefficients () const |
unsigned int | get_offset_coefficients () const |
unsigned int | get_inclusions_in_vessel (unsigned int vessel_id) const |
void | compute_rotated_inclusion_data () |
std::vector< double > | get_rotated_inclusion_data (const types::global_dof_index &inclusion_id) const |
return the rotate the data of a given inclusion | |
void | set_n_q_points (unsigned int n_q) |
void | set_n_coefficients (unsigned int n_coefficients) |
![]() | |
ParameterAcceptor (const std::string §ion_name="") | |
unsigned int | get_acceptor_id () const |
virtual | ~ParameterAcceptor () override |
virtual void | declare_parameters (ParameterHandler &prm) |
virtual void | parse_parameters (ParameterHandler &prm) |
std::string | get_section_name () const |
std::vector< std::string > | get_section_path () const |
void | add_parameter (const std::string &entry, ParameterType ¶meter, const std::string &documentation="", ParameterHandler &prm_=prm, const Patterns::PatternBase &pattern=*Patterns::Tools::Convert< ParameterType >::to_pattern()) |
void | enter_subsection (const std::string &subsection) |
void | leave_subsection () |
void | enter_my_subsection (ParameterHandler &prm) |
void | leave_my_subsection (ParameterHandler &prm) |
void | serialize (Archive &ar, const unsigned int version) |
unsigned int | n_subscriptions () const |
void | list_subscribers (StreamType &stream) const |
void | list_subscribers () const |
unsigned int | n_subscriptions () const |
void | list_subscribers (StreamType &stream) const |
void | list_subscribers () const |
![]() | |
EnableObserverPointer () | |
EnableObserverPointer (const EnableObserverPointer &) | |
EnableObserverPointer (EnableObserverPointer &&) noexcept | |
virtual | ~EnableObserverPointer () |
EnableObserverPointer & | operator= (const EnableObserverPointer &) |
EnableObserverPointer & | operator= (EnableObserverPointer &&) noexcept |
void | serialize (Archive &ar, const unsigned int version) |
unsigned int | n_subscriptions () const |
void | list_subscribers (StreamType &stream) const |
void | list_subscribers () const |
unsigned int | n_subscriptions () const |
void | list_subscribers (StreamType &stream) const |
void | list_subscribers () const |
Public Attributes | |
ParameterAcceptorProxy< Functions::ParsedFunction< spacedim > > | inclusions_rhs |
Particles::ParticleHandler< spacedim > | inclusions_as_particles |
std::vector< std::vector< double > > | inclusions |
std::string | data_file = "" |
std::vector< std::vector< double > > | inclusions_data |
std::vector< double > | reference_inclusion_data |
std::vector< std::vector< double > > | rotated_inclusion_data |
std::map< unsigned int, std::vector< types::global_dof_index > > | map_vessel_inclusions |
![]() | |
boost::signals2::signal< void()> | declare_parameters_call_back |
boost::signals2::signal< void()> | parse_parameters_call_back |
Private Member Functions | |
void | check_vessels () |
Check that all vesselsID are present and create the map vessel_inclusions. | |
Private Attributes | |
unsigned int | n_q_points = 100 |
unsigned int | n_coefficients = 1 |
unsigned int | offset_coefficients = 0 |
std::vector< unsigned int > | selected_coefficients |
const unsigned int | n_vector_components |
MPI_Comm | mpi_communicator |
std::vector< Point< spacedim > > | support_points |
std::vector< double > | theta |
std::vector< Tensor< 1, spacedim > > | normals |
std::vector< Point< spacedim > > | current_support_points |
std::vector< double > | current_fe_values |
unsigned int | n_vessels = 1 |
std::string | inclusions_file = "" |
unsigned int | rtree_extraction_level = 1 |
Friends | |
template<int dim, typename number, int n_components> | |
class | CouplingOperator |
Additional Inherited Members | |
![]() | |
static void | initialize (const std::string &filename="", const std::string &output_filename="", const ParameterHandler::OutputStyle output_style_for_output_filename=ParameterHandler::Short, ParameterHandler &prm=ParameterAcceptor::prm, const ParameterHandler::OutputStyle output_style_for_filename=ParameterHandler::DefaultStyle) |
static void | initialize (std::istream &input_stream, ParameterHandler &prm=ParameterAcceptor::prm) |
static void | clear () |
static void | parse_all_parameters (ParameterHandler &prm=ParameterAcceptor::prm) |
static void | declare_all_parameters (ParameterHandler &prm=ParameterAcceptor::prm) |
static ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
![]() | |
static ::ExceptionBase & | ExcInUse (int arg1, std::string arg2, std::string arg3) |
static ::ExceptionBase & | ExcNoSubscriber (std::string arg1, std::string arg2) |
![]() | |
static ParameterHandler | prm |
![]() | |
const std::string | section_name |
std::vector< std::string > | subsections |
Class for handling inclusions in an immersed boundary method.
This class provides functionality for handling inclusions in an immersed boundary method. It stores the positions, radii, and Fourier coefficients of the inclusions, and provides methods for computing their normals, center, and direction. It also provides methods for initializing the inclusions from a file, setting up particles for the inclusions, and getting degrees of freedom associated with each inclusion.
Definition at line 64 of file inclusions.h.
|
inline |
Class for computing the inclusions of a given mesh.
n_vector_components | Number of vector components. |
Definition at line 75 of file inclusions.h.
References ParameterAcceptor::add_parameter(), data_file, inclusions, inclusions_file, inclusions_rhs, n_coefficients, n_q_points, n_vector_components, ParameterAcceptor::ParameterAcceptor(), ParameterAcceptor::prm, reference_inclusion_data, rtree_extraction_level, and selected_coefficients.
|
inlineprivate |
Check that all vesselsID are present and create the map vessel_inclusions.
Definition at line 1065 of file inclusions.h.
References get_vesselID(), inclusions, map_vessel_inclusions, and n_vessels.
Referenced by initialize().
|
inline |
Definition at line 941 of file inclusions.h.
References get_rotation(), inclusions_data, rotated_inclusion_data, and Tensor< int rank_, int dim, typename Number >::unroll().
|
inline |
Get the center of the inclusion.
inclusion_id |
Definition at line 551 of file inclusions.h.
References inclusions.
Referenced by get_current_support_points().
|
inline |
Get the ith component for the given dof index.
dof_index | A number in [0,n_dofs()) |
Definition at line 366 of file inclusions.h.
References n_dofs(), n_dofs_per_inclusion(), n_vector_components, and selected_coefficients.
|
inline |
Definition at line 683 of file inclusions.h.
References current_support_points, get_center(), get_radius(), get_rotation(), n_q_points, and support_points.
Referenced by setup_inclusions_particles().
|
inline |
Get the direction of the inclusion.
inclusion_id |
Definition at line 613 of file inclusions.h.
References inclusions, and Tensor< int rank_, int dim, typename Number >::norm().
Referenced by get_JxW(), get_rotation(), and get_section_measure().
|
inline |
Returns the degrees of freedom indices associated with a given quadrature point.
quadrature_id | The global index of the quadrature point. |
Definition at line 277 of file inclusions.h.
References n_dofs_per_inclusion(), n_particles(), and n_q_points.
|
inline |
Get a list of fe values
particle_id |
Definition at line 492 of file inclusions.h.
References std::cos(), current_fe_values, get_radius(), inclusions, n_coefficients, n_particles(), n_q_points, n_vector_components, selected_coefficients, std::sin(), std::sqrt(), and theta.
|
inline |
Get the ith Fourier component for the given dof index.
dof_index | A number in [0,n_dofs()) |
Definition at line 382 of file inclusions.h.
References n_dofs(), n_dofs_per_inclusion(), and selected_coefficients.
Referenced by get_inclusion_data(), and get_inclusion_data().
|
inline |
Get the Fourier data for the given local dof index.
inclusion_id | A number in [0,n_inclusions()) |
dof_index | A number in [0,n_dofs_per_inclusion()) |
Definition at line 398 of file inclusions.h.
References get_fourier_component(), get_rotated_inclusion_data(), n_dofs(), and n_inclusions().
|
inline |
Get the Fourier data for the given local dof index.
inclusion_id | A number in [0,n_inclusions()) |
dof_index | A number in [0,n_dofs_per_inclusion()) |
point | The inclusion quadrature point location |
Definition at line 417 of file inclusions.h.
References get_fourier_component(), get_rotated_inclusion_data(), inclusions_data, inclusions_rhs, n_dofs(), n_inclusions(), and n_vector_components.
|
inline |
Quadrature id to inclusion id.
quadrature_id |
Definition at line 353 of file inclusions.h.
References n_particles(), and n_q_points.
Referenced by get_normal().
|
inline |
Definition at line 933 of file inclusions.h.
References map_vessel_inclusions, and n_vessels.
|
inline |
return weight for integration normal direction of an inclusion at a quadrature point
quadrature_id |
Definition at line 475 of file inclusions.h.
References get_direction(), get_radius(), inclusions, n_particles(), n_q_points, Tensor< int rank_, int dim, typename Number >::norm(), and numbers::PI.
|
inline |
Definition at line 921 of file inclusions.h.
References n_coefficients.
|
inline |
Definition at line 915 of file inclusions.h.
References n_vessels.
|
inline |
Get the normal.
quadrature_id |
Definition at line 455 of file inclusions.h.
References get_inclusion_id(), get_rotation(), n_particles(), n_q_points, and normals.
|
inline |
Definition at line 927 of file inclusions.h.
References offset_coefficients.
|
inline |
Get the radius of the inclusion.
inclusion_id |
Definition at line 570 of file inclusions.h.
References inclusions.
Referenced by get_current_support_points(), get_fe_values(), get_JxW(), and get_section_measure().
|
inline |
return the rotate the data of a given inclusion
inclusion_id |
Definition at line 999 of file inclusions.h.
References inclusions, inclusions_data, and rotated_inclusion_data.
Referenced by get_inclusion_data(), and get_inclusion_data().
|
inline |
Get the rotation of the inclusion.
inclusion_id |
Definition at line 644 of file inclusions.h.
References std::abs(), and get_direction().
Referenced by compute_rotated_inclusion_data(), get_current_support_points(), and get_normal().
|
inline |
Get the measure of the section of the inclusion.
inclusion_id |
Definition at line 593 of file inclusions.h.
References get_direction(), get_radius(), Tensor< int rank_, int dim, typename Number >::norm(), and numbers::PI.
|
inline |
3D return the vessel that a given inclusion belongs to, 2D return 0
inclusion_id |
Definition at line 871 of file inclusions.h.
References inclusions.
Referenced by check_vessels().
|
inline |
Reinit all reference points and normals, and read inclusions file.
Definition at line 168 of file inclusions.h.
References check_vessels(), std::cos(), current_fe_values, current_support_points, data_file, inclusions, inclusions_data, inclusions_file, mpi_communicator, n_coefficients, n_dofs_per_inclusion(), n_inclusions(), n_q_points, n_vector_components, normals, numbers::PI, reference_inclusion_data, selected_coefficients, std::sin(), support_points, theta, and Utilities::MPI::this_mpi_process().
Referenced by setup_inclusions_particles().
|
inline |
Returns the number of degrees of freedom in the system.
Definition at line 123 of file inclusions.h.
References inclusions, and n_dofs_per_inclusion().
Referenced by get_component(), get_fourier_component(), get_inclusion_data(), get_inclusion_data(), and setup_inclusions_particles().
|
inline |
Number of degrees of freedom associated to each inclusion.
Definition at line 158 of file inclusions.h.
References selected_coefficients.
Referenced by get_component(), get_dof_indices(), get_fourier_component(), initialize(), and n_dofs().
|
inline |
Returns the number of inclusions in the mesh.
Definition at line 147 of file inclusions.h.
References inclusions.
Referenced by get_inclusion_data(), get_inclusion_data(), initialize(), setup_inclusions_particles(), and update_displacement_hdf5().
|
inline |
Returns the number of particles in the system.
Definition at line 135 of file inclusions.h.
References inclusions, and n_q_points.
Referenced by get_dof_indices(), get_fe_values(), get_inclusion_id(), get_JxW(), get_normal(), and setup_inclusions_particles().
|
inline |
print the inclusions in parallel on a .vtu file
filename |
Definition at line 702 of file inclusions.h.
References Particles::DataOut< int dim, int spacedim >::build_patches(), inclusions_as_particles, mpi_communicator, and Particles::DataOut< int dim, int spacedim >::write_vtu_in_parallel().
|
inline |
Definition at line 1018 of file inclusions.h.
References n_coefficients.
|
inline |
Definition at line 1011 of file inclusions.h.
References n_q_points.
|
inline |
Sets up the inclusions particles for the given triangulation.
tria | The triangulation to set up the inclusions particles for. |
Definition at line 293 of file inclusions.h.
References parallel::distributed::Triangulation< int dim, int spacedim >::active_cell_iterators(), Utilities::MPI::all_gather(), Utilities::MPI::create_evenly_distributed_partitioning(), parallel::distributed::Triangulation< int dim, int spacedim >::get_communicator(), get_current_support_points(), inclusions_as_particles, initialize(), StaticMappingQ1< int dim, int spacedim >::mapping, mpi_communicator, n_dofs(), n_inclusions(), parallel::distributed::Triangulation< int dim, int spacedim >::n_locally_owned_active_cells(), n_particles(), and rtree_extraction_level.
|
inline |
Update the displacement data after the initialization reading from a hdf5 file.
Definition at line 713 of file inclusions.h.
References Utilities::MPI::create_evenly_distributed_partitioning(), data_file, inclusions_data, mpi_communicator, and n_inclusions().
|
inline |
Update the displacement data after the initialization, 3D case only.
new_data | vector of lenght equal to the number of inclusions or to the number of vessels, elements of the vector are the new values to be assigned for normal expansion |
Definition at line 759 of file inclusions.h.
References inclusions, map_vessel_inclusions, n_vessels, and update_single_inclusion_data_along_normal().
|
inline |
Definition at line 795 of file inclusions.h.
References map_vessel_inclusions, n_vessels, and update_single_inclusion_data_along_normal().
|
inline |
Definition at line 887 of file inclusions.h.
References inclusions_data.
Referenced by update_inclusions_data(), and update_inclusions_data().
|
inline |
Definition at line 910 of file inclusions.h.
|
friend |
Definition at line 68 of file inclusions.h.
References CouplingOperator.
Referenced by CouplingOperator.
|
mutableprivate |
Definition at line 1053 of file inclusions.h.
Referenced by get_fe_values(), and initialize().
|
mutableprivate |
Definition at line 1052 of file inclusions.h.
Referenced by get_current_support_points(), and initialize().
std::string Inclusions< spacedim >::data_file = "" |
Definition at line 1028 of file inclusions.h.
Referenced by Inclusions(), initialize(), and update_displacement_hdf5().
std::vector<std::vector<double> > Inclusions< spacedim >::inclusions |
Definition at line 1026 of file inclusions.h.
Referenced by check_vessels(), get_center(), get_direction(), get_fe_values(), get_JxW(), get_radius(), get_rotated_inclusion_data(), get_vesselID(), Inclusions(), initialize(), n_dofs(), n_inclusions(), n_particles(), and update_inclusions_data().
Particles::ParticleHandler<spacedim> Inclusions< spacedim >::inclusions_as_particles |
Definition at line 1025 of file inclusions.h.
Referenced by output_particles(), and setup_inclusions_particles().
std::vector<std::vector<double> > Inclusions< spacedim >::inclusions_data |
Definition at line 1032 of file inclusions.h.
Referenced by compute_rotated_inclusion_data(), get_inclusion_data(), get_rotated_inclusion_data(), initialize(), update_displacement_hdf5(), and update_single_inclusion_data_along_normal().
|
private |
Definition at line 1057 of file inclusions.h.
Referenced by Inclusions(), and initialize().
ParameterAcceptorProxy<Functions::ParsedFunction<spacedim> > Inclusions< spacedim >::inclusions_rhs |
Definition at line 1024 of file inclusions.h.
Referenced by get_inclusion_data(), and Inclusions().
std::map<unsigned int, std::vector<types::global_dof_index> > Inclusions< spacedim >::map_vessel_inclusions |
Definition at line 1037 of file inclusions.h.
Referenced by check_vessels(), get_inclusions_in_vessel(), update_inclusions_data(), and update_inclusions_data().
|
private |
Definition at line 1046 of file inclusions.h.
Referenced by initialize(), output_particles(), setup_inclusions_particles(), and update_displacement_hdf5().
|
private |
Definition at line 1041 of file inclusions.h.
Referenced by get_fe_values(), get_n_coefficients(), Inclusions(), initialize(), and set_n_coefficients().
|
private |
Definition at line 1040 of file inclusions.h.
Referenced by get_current_support_points(), get_dof_indices(), get_fe_values(), get_inclusion_id(), get_JxW(), get_normal(), Inclusions(), initialize(), n_particles(), and set_n_q_points().
|
private |
Definition at line 1045 of file inclusions.h.
Referenced by get_component(), get_fe_values(), get_inclusion_data(), Inclusions(), and initialize().
|
private |
Definition at line 1055 of file inclusions.h.
Referenced by check_vessels(), get_inclusions_in_vessel(), get_n_vessels(), update_inclusions_data(), and update_inclusions_data().
|
mutableprivate |
Definition at line 1051 of file inclusions.h.
Referenced by get_normal(), and initialize().
|
private |
Definition at line 1042 of file inclusions.h.
Referenced by get_offset_coefficients().
std::vector<double> Inclusions< spacedim >::reference_inclusion_data |
Definition at line 1033 of file inclusions.h.
Referenced by Inclusions(), and initialize().
std::vector<std::vector<double> > Inclusions< spacedim >::rotated_inclusion_data |
Definition at line 1034 of file inclusions.h.
Referenced by compute_rotated_inclusion_data(), and get_rotated_inclusion_data().
|
private |
Definition at line 1058 of file inclusions.h.
Referenced by CouplingOperator< dim, number, n_components >::CouplingOperator(), Inclusions(), and setup_inclusions_particles().
|
private |
Definition at line 1043 of file inclusions.h.
Referenced by get_component(), get_fe_values(), get_fourier_component(), Inclusions(), initialize(), and n_dofs_per_inclusion().
|
private |
Definition at line 1047 of file inclusions.h.
Referenced by get_current_support_points(), and initialize().
|
private |
Definition at line 1048 of file inclusions.h.
Referenced by get_fe_values(), and initialize().