From 82132f67995ee51aaf34186cdda6d61cf9e01577 Mon Sep 17 00:00:00 2001 From: Daniel Garcia-Sanchez Date: Fri, 2 Aug 2019 10:39:51 +0200 Subject: [PATCH] Add std::complex make_periodicity_constraints() --- include/deal.II/dofs/dof_tools.h | 16 ++++----- source/dofs/dof_tools_constraints.cc | 22 ++++++------ source/dofs/dof_tools_constraints.inst.in | 44 ++++++++++++++++++++++- 3 files changed, 62 insertions(+), 20 deletions(-) diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index 28bde4e4ab..1045467385 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -1153,12 +1153,12 @@ namespace DoFTools * * @author Matthias Maier, 2012 - 2015 */ - template + template void make_periodicity_constraints( const FaceIterator & face_1, const typename identity::type &face_2, - AffineConstraints & constraints, + AffineConstraints & constraints, const ComponentMask & component_mask = ComponentMask(), const bool face_orientation = true, const bool face_flip = false, @@ -1190,13 +1190,13 @@ namespace DoFTools * * @author Daniel Arndt, Matthias Maier, 2013 - 2015 */ - template + template void make_periodicity_constraints( const std::vector< GridTools::PeriodicFacePair> & periodic_faces, - AffineConstraints & constraints, + AffineConstraints & constraints, const ComponentMask & component_mask = ComponentMask(), const std::vector &first_vector_components = std::vector()); @@ -1233,14 +1233,14 @@ namespace DoFTools * * @author Matthias Maier, 2012 */ - template + template void make_periodicity_constraints( const DoFHandlerType & dof_handler, const types::boundary_id b_id1, const types::boundary_id b_id2, const int direction, - AffineConstraints &constraints, + AffineConstraints &constraints, const ComponentMask & component_mask = ComponentMask()); @@ -1269,13 +1269,13 @@ namespace DoFTools * @ref GlossPeriodicConstraints "Glossary entry on periodic boundary conditions" * for further information. */ - template + template void make_periodicity_constraints( const DoFHandlerType & dof_handler, const types::boundary_id b_id, const int direction, - AffineConstraints &constraints, + AffineConstraints &constraints, const ComponentMask & component_mask = ComponentMask()); /** diff --git a/source/dofs/dof_tools_constraints.cc b/source/dofs/dof_tools_constraints.cc index 44a85b5545..3a9e3f65b4 100644 --- a/source/dofs/dof_tools_constraints.cc +++ b/source/dofs/dof_tools_constraints.cc @@ -1817,13 +1817,13 @@ namespace DoFTools * such "identity constraints" if the opposite constraint already * exists. */ - template + template void set_periodicity_constraints( const FaceIterator & face_1, const typename identity::type &face_2, const FullMatrix & transformation, - AffineConstraints & affine_constraints, + AffineConstraints & affine_constraints, const ComponentMask & component_mask, const bool face_orientation, const bool face_flip, @@ -2073,7 +2073,7 @@ namespace DoFTools // a dependency cycle: bool constraints_are_cyclic = true; - double cycle_factor = factor; + number cycle_factor = factor; for (auto test_dof = dof_right; test_dof != dof_left;) { @@ -2220,12 +2220,12 @@ namespace DoFTools // Low level interface: - template + template void make_periodicity_constraints( const FaceIterator & face_1, const typename identity::type &face_2, - AffineConstraints & affine_constraints, + AffineConstraints & affine_constraints, const ComponentMask & component_mask, const bool face_orientation, const bool face_flip, @@ -2449,13 +2449,13 @@ namespace DoFTools - template + template void make_periodicity_constraints( const std::vector< GridTools::PeriodicFacePair> & periodic_faces, - AffineConstraints & constraints, + AffineConstraints & constraints, const ComponentMask & component_mask, const std::vector &first_vector_components) { @@ -2489,13 +2489,13 @@ namespace DoFTools // High level interface variants: - template + template void make_periodicity_constraints(const DoFHandlerType & dof_handler, const types::boundary_id b_id1, const types::boundary_id b_id2, const int direction, - dealii::AffineConstraints &constraints, + dealii::AffineConstraints &constraints, const ComponentMask &component_mask) { static const int space_dim = DoFHandlerType::space_dimension; @@ -2522,12 +2522,12 @@ namespace DoFTools - template + template void make_periodicity_constraints(const DoFHandlerType & dof_handler, const types::boundary_id b_id, const int direction, - AffineConstraints &constraints, + AffineConstraints &constraints, const ComponentMask & component_mask) { static const int dim = DoFHandlerType::dimension; diff --git a/source/dofs/dof_tools_constraints.inst.in b/source/dofs/dof_tools_constraints.inst.in index 931ef21f26..82c18c6f16 100644 --- a/source/dofs/dof_tools_constraints.inst.in +++ b/source/dofs/dof_tools_constraints.inst.in @@ -60,13 +60,36 @@ for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS) const FullMatrix &, const std::vector &); - template void DoFTools::make_periodicity_constraints>( +#ifdef DEAL_II_WITH_COMPLEX_VALUES + template void DoFTools::make_periodicity_constraints( + const DH::face_iterator &, + const DH::face_iterator &, + AffineConstraints> &, + const ComponentMask &, + bool, + bool, + bool, + const FullMatrix &, + const std::vector &); +#endif + + template void + DoFTools::make_periodicity_constraints, double>( const std::vector< GridTools::PeriodicFacePair::cell_iterator>> &, AffineConstraints &, const ComponentMask &, const std::vector &); +#ifdef DEAL_II_WITH_COMPLEX_VALUES + template void DoFTools::make_periodicity_constraints, + std::complex>( + const std::vector< + GridTools::PeriodicFacePair::cell_iterator>> &, + AffineConstraints> &, + const ComponentMask &, + const std::vector &); +#endif template void DoFTools::make_periodicity_constraints( const DH &, @@ -76,12 +99,31 @@ for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS) AffineConstraints &, const ComponentMask &); +#ifdef DEAL_II_WITH_COMPLEX_VALUES + template void DoFTools::make_periodicity_constraints( + const DH &, + const types::boundary_id, + const types::boundary_id, + const int, + AffineConstraints> &, + const ComponentMask &); +#endif + template void DoFTools::make_periodicity_constraints( const DH &, const types::boundary_id, const int, AffineConstraints &, const ComponentMask &); + +#ifdef DEAL_II_WITH_COMPLEX_VALUES + template void DoFTools::make_periodicity_constraints( + const DH &, + const types::boundary_id, + const int, + AffineConstraints> &, + const ComponentMask &); +#endif } for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS; -- 2.39.5