From befcdc95471907a06d0c24180b3b5c488347fdce Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 6 Dec 2017 13:00:15 -0700 Subject: [PATCH] Provide more instantiations for std::complex data types. --- include/deal.II/lac/constraint_matrix.h | 8 +- .../deal.II/lac/constraint_matrix.templates.h | 37 +-- .../numerics/matrix_creator.templates.h | 110 ++++----- include/deal.II/numerics/matrix_tools.h | 32 +-- .../deal.II/numerics/vector_tools.templates.h | 38 ++- source/lac/constraint_matrix.cc | 63 ++--- source/lac/constraint_matrix.inst.in | 66 +++-- source/lac/sparse_matrix.inst.in | 13 + source/numerics/matrix_creator.inst.in | 230 ++++-------------- source/numerics/vector_tools_boundary.inst.in | 41 ++-- 10 files changed, 289 insertions(+), 349 deletions(-) diff --git a/include/deal.II/lac/constraint_matrix.h b/include/deal.II/lac/constraint_matrix.h index a58fe540fc..7c12c0d0c5 100644 --- a/include/deal.II/lac/constraint_matrix.h +++ b/include/deal.II/lac/constraint_matrix.h @@ -1464,13 +1464,13 @@ private: /** * Internal helper function for distribute_local_to_global function. */ - template - LocalType + template + typename ProductType::type resolve_vector_entry (const size_type i, const internals::GlobalRowsFromLocal &global_rows, - const Vector &local_vector, + const Vector &local_vector, const std::vector &local_dof_indices, - const FullMatrix &local_matrix) const; + const FullMatrix &local_matrix) const; }; diff --git a/include/deal.II/lac/constraint_matrix.templates.h b/include/deal.II/lac/constraint_matrix.templates.h index 6c397f2e24..4c732d2e2f 100644 --- a/include/deal.II/lac/constraint_matrix.templates.h +++ b/include/deal.II/lac/constraint_matrix.templates.h @@ -1399,7 +1399,7 @@ namespace internals * Access to the scratch data is only through the accessor class which * handles the access as well as marking the data as used. */ - template + template class ConstraintMatrixData { public: @@ -1434,7 +1434,7 @@ namespace internals /** * Temporary array for column values */ - std::vector values; + std::vector values; /** * Temporary array for block start indices @@ -1449,7 +1449,7 @@ namespace internals /** * Temporary array for vector values */ - std::vector vector_values; + std::vector vector_values; /** * Data array for reorder row/column indices. @@ -2187,19 +2187,19 @@ make_sorted_row_list (const std::vector &local_dof_indices, // Resolve the constraints from the vector and apply inhomogeneities. -template < typename LocalType> +template inline -LocalType +typename ProductType::type ConstraintMatrix:: resolve_vector_entry (const size_type i, const internals::GlobalRowsFromLocal &global_rows, - const Vector &local_vector, + const Vector &local_vector, const std::vector &local_dof_indices, - const FullMatrix &local_matrix) const + const FullMatrix &local_matrix) const { const size_type loc_row = global_rows.local_row(i); const size_type n_inhomogeneous_rows = global_rows.n_inhomogeneities(); - LocalType val = 0; + typename ProductType::type val = 0; // has a direct contribution from some local entry. If we have inhomogeneous // constraints, compute the contribution of the inhomogeneity in the current // row. @@ -2217,7 +2217,7 @@ resolve_vector_entry (const size_type i, for (size_type q=0; q::type add_this = local_vector (loc_row_q); for (size_type k=0; k::ScratchDataAccessor + typename internals::ConstraintMatrixData::ScratchDataAccessor scratch_data; internals::GlobalRowsFromLocal &global_rows = scratch_data->global_rows; @@ -2282,7 +2282,7 @@ ConstraintMatrix::distribute_local_to_global ( std::vector &vals = scratch_data->values; // create arrays for writing into the vector as well std::vector &vector_indices = scratch_data->vector_indices; - std::vector &vector_values = scratch_data->vector_values; + std::vector &vector_values = scratch_data->vector_values; vector_indices.resize(n_actual_dofs); vector_values.resize(n_actual_dofs); SparseMatrix *sparse_matrix @@ -2328,13 +2328,14 @@ ConstraintMatrix::distribute_local_to_global ( // hand side. if (use_vectors == true) { - const number val = resolve_vector_entry (i, global_rows, - local_vector, - local_dof_indices, - local_matrix); + const typename VectorType::value_type + val = resolve_vector_entry (i, global_rows, + local_vector, + local_dof_indices, + local_matrix); AssertIsFinite(val); - if (val != number ()) + if (val != typename VectorType::value_type ()) { vector_indices[local_row_n] = row; vector_values[local_row_n] = val; @@ -2406,7 +2407,7 @@ distribute_local_to_global ( } Assert (sorted == true, ExcMatrixNotClosed()); - typename internals::ConstraintMatrixData::ScratchDataAccessor + typename internals::ConstraintMatrixData::ScratchDataAccessor scratch_data; const size_type n_local_dofs = local_dof_indices.size(); @@ -2528,7 +2529,7 @@ ConstraintMatrix::distribute_local_to_global ( const size_type n_local_row_dofs = row_indices.size(); const size_type n_local_col_dofs = col_indices.size(); - typename internals::ConstraintMatrixData::ScratchDataAccessor + typename internals::ConstraintMatrixData::ScratchDataAccessor scratch_data; internals::GlobalRowsFromLocal &global_rows = scratch_data->global_rows; global_rows.reinit(n_local_row_dofs); diff --git a/include/deal.II/numerics/matrix_creator.templates.h b/include/deal.II/numerics/matrix_creator.templates.h index 578a8c054b..39b48200ba 100644 --- a/include/deal.II/numerics/matrix_creator.templates.h +++ b/include/deal.II/numerics/matrix_creator.templates.h @@ -142,7 +142,7 @@ namespace MatrixCreator { Scratch (const ::dealii::hp::FECollection &fe, const UpdateFlags update_flags, - const Function *coefficient, + const Function::real_type> *coefficient, const Function *rhs_function, const ::dealii::hp::QCollection &quadrature, const ::dealii::hp::MappingCollection &mapping) @@ -156,7 +156,7 @@ namespace MatrixCreator update_flags), coefficient_values(quadrature_collection.max_n_quadrature_points()), coefficient_vector_values (quadrature_collection.max_n_quadrature_points(), - dealii::Vector (fe_collection.n_components())), + dealii::Vector::real_type> (fe_collection.n_components())), rhs_values(quadrature_collection.max_n_quadrature_points()), rhs_vector_values(quadrature_collection.max_n_quadrature_points(), dealii::Vector (fe_collection.n_components())), @@ -196,13 +196,13 @@ namespace MatrixCreator ::dealii::hp::FEValues x_fe_values; - std::vector coefficient_values; - std::vector > coefficient_vector_values; - std::vector rhs_values; - std::vector > rhs_vector_values; + std::vector::real_type> coefficient_values; + std::vector::real_type> > coefficient_vector_values; + std::vector rhs_values; + std::vector > rhs_vector_values; - const Function *coefficient; - const Function *rhs_function; + const Function::real_type> *coefficient; + const Function *rhs_function; const UpdateFlags update_flags; }; @@ -211,10 +211,10 @@ namespace MatrixCreator template struct CopyData { - std::vector dof_indices; - FullMatrix cell_matrix; - dealii::Vector cell_rhs; - const ConstraintMatrix *constraints; + std::vector dof_indices; + FullMatrix::real_type> cell_matrix; + dealii::Vector cell_rhs; + const ConstraintMatrix *constraints; }; } @@ -225,7 +225,7 @@ namespace MatrixCreator typename number> void mass_assembler (const CellIterator &cell, MatrixCreator::internal::AssemblerData::Scratch &data, - MatrixCreator::internal::AssemblerData::CopyData ©_data) + MatrixCreator::internal::AssemblerData::CopyData ©_data) { data.x_fe_values.reinit (cell); const FEValues &fe_values = data.x_fe_values.get_present_fe_values (); @@ -280,14 +280,13 @@ namespace MatrixCreator else { data.coefficient_vector_values.resize (n_q_points, - dealii::Vector(n_components)); + dealii::Vector::real_type>(n_components)); data.coefficient->vector_value_list (fe_values.get_quadrature_points(), data.coefficient_vector_values); } } - number add_data; const std::vector &JxW = fe_values.get_JxW_values(); for (unsigned int i=0; i::real_type add_data = 0; if (use_coefficient) { if (data.coefficient->n_components==1) @@ -331,7 +329,7 @@ namespace MatrixCreator if (use_rhs_function) { - add_data = 0; + number add_data = 0; if (data.rhs_function->n_components==1) for (unsigned int point=0; point::real_type add_data = 0; for (unsigned int comp_i = 0; comp_i < n_components; ++comp_i) if (fe.get_nonzero_components(i)[comp_i] && fe.get_nonzero_components(j)[comp_i]) @@ -381,7 +379,7 @@ namespace MatrixCreator if (use_rhs_function) { - add_data = 0; + number add_data = 0; for (unsigned int comp_i = 0; comp_i < n_components; ++comp_i) if (fe.get_nonzero_components(i)[comp_i]) { @@ -408,7 +406,7 @@ namespace MatrixCreator typename CellIterator> void laplace_assembler (const CellIterator &cell, MatrixCreator::internal::AssemblerData::Scratch &data, - MatrixCreator::internal::AssemblerData::CopyData ©_data) + MatrixCreator::internal::AssemblerData::CopyData ©_data) { data.x_fe_values.reinit (cell); const FEValues &fe_values = data.x_fe_values.get_present_fe_values (); @@ -630,12 +628,12 @@ namespace MatrixCreator CopyData(CopyData const &data); - unsigned int dofs_per_cell; - std::vector dofs; - std::vector > dof_is_on_face; - typename DoFHandlerType::active_cell_iterator cell; - std::vector > cell_matrix; - std::vector > cell_vector; + unsigned int dofs_per_cell; + std::vector dofs; + std::vector > dof_is_on_face; + typename DoFHandlerType::active_cell_iterator cell; + std::vector::real_type> > cell_matrix; + std::vector > cell_vector; }; @@ -723,10 +721,10 @@ namespace MatrixCreator void create_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const coefficient, + const Function::real_type> *const coefficient, const ConstraintMatrix &constraints) { Assert (matrix.m() == dof.n_dofs(), @@ -754,7 +752,7 @@ namespace MatrixCreator static_cast::active_cell_iterator>(dof.end()), &MatrixCreator::internal::mass_assembler::active_cell_iterator,number>, std::bind(&MatrixCreator::internal:: - copy_local_to_global, Vector >, + copy_local_to_global::real_type>, Vector >, std::placeholders::_1, &matrix, &rhs_vector), assembler_data, copy_data); @@ -765,10 +763,10 @@ namespace MatrixCreator template void create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const coefficient, + const Function::real_type> *const coefficient, const ConstraintMatrix &constraints) { create_mass_matrix(StaticMappingQ1::mapping, @@ -833,10 +831,10 @@ namespace MatrixCreator void create_mass_matrix (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const coefficient, + const Function::real_type> *const coefficient, const ConstraintMatrix &constraints) { Assert (matrix.m() == dof.n_dofs(), @@ -861,7 +859,7 @@ namespace MatrixCreator static_cast::active_cell_iterator>(dof.end()), &MatrixCreator::internal::mass_assembler::active_cell_iterator,number>, std::bind (&MatrixCreator::internal:: - copy_local_to_global, Vector >, + copy_local_to_global::real_type>, Vector >, std::placeholders::_1, &matrix, &rhs_vector), assembler_data, copy_data); @@ -872,10 +870,10 @@ namespace MatrixCreator template void create_mass_matrix (const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const coefficient, + const Function::real_type> *const coefficient, const ConstraintMatrix &constraints) { create_mass_matrix(hp::StaticMappingQ1::mapping_collection, dof, q, @@ -892,12 +890,12 @@ namespace MatrixCreator create_boundary_mass_matrix_1 (typename DoFHandler::active_cell_iterator const &cell, MatrixCreator::internal::AssemblerBoundary::Scratch const &, MatrixCreator::internal::AssemblerBoundary::CopyData,number > ©_data, + spacedim>,number> ©_data, Mapping const &mapping, FiniteElement const &fe, Quadrature const &q, std::map*> const &boundary_functions, - Function const *const coefficient, + Function::real_type> const *const coefficient, std::vector const &component_mapping) { @@ -921,9 +919,9 @@ namespace MatrixCreator // two variables for the coefficient, one for the two cases // indicated in the name - std::vector coefficient_values (fe_values.n_quadrature_points, 1.); - std::vector > coefficient_vector_values (fe_values.n_quadrature_points, - Vector(n_components)); + std::vector::real_type> coefficient_values (fe_values.n_quadrature_points, 1.); + std::vector::real_type> > coefficient_vector_values (fe_values.n_quadrature_points, + Vector::real_type>(n_components)); const bool coefficient_is_vector = (coefficient != nullptr && coefficient->n_components != 1); std::vector rhs_values_scalar (fe_values.n_quadrature_points); @@ -1096,7 +1094,7 @@ namespace MatrixCreator spacedim>,number > const ©_data, std::map*> const &boundary_functions, std::vector const &dof_to_boundary_mapping, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, Vector &rhs_vector) { // now transfer cell matrix and vector to the whole boundary matrix @@ -1204,11 +1202,11 @@ namespace MatrixCreator create_boundary_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &boundary_functions, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const coefficient, + const Function::real_type> *const coefficient, std::vector component_mapping) { // what would that be in 1d? the @@ -1283,7 +1281,7 @@ namespace MatrixCreator hp::FECollection const &fe_collection, hp::QCollection const &q, const std::map*> &boundary_functions, - Function const *const coefficient, + Function::real_type> const *const coefficient, std::vector const &component_mapping) { const unsigned int n_components = fe_collection.n_components(); @@ -1308,8 +1306,8 @@ namespace MatrixCreator // two variables for the coefficient, // one for the two cases indicated in // the name - std::vector coefficient_values; - std::vector > coefficient_vector_values; + std::vector::real_type> coefficient_values; + std::vector::real_type> > coefficient_vector_values; const bool coefficient_is_vector = (coefficient != nullptr && coefficient->n_components != 1); @@ -1347,7 +1345,7 @@ namespace MatrixCreator // FE has several components { coefficient_vector_values.resize (fe_values.n_quadrature_points, - Vector(n_components)); + Vector::real_type>(n_components)); coefficient_values.resize (fe_values.n_quadrature_points, 1.); rhs_values_system.resize (fe_values.n_quadrature_points, @@ -1498,7 +1496,7 @@ namespace MatrixCreator ::CopyData,number > const ©_data, std::map*> const &boundary_functions, std::vector const &dof_to_boundary_mapping, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, Vector &rhs_vector) { // now transfer cell matrix and vector to the whole boundary matrix @@ -1591,11 +1589,11 @@ namespace MatrixCreator template void create_boundary_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &rhs, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const a, + const Function::real_type> *const a, std::vector component_mapping) { create_boundary_mass_matrix(StaticMappingQ1::mapping, dof, q, @@ -1609,11 +1607,11 @@ namespace MatrixCreator create_boundary_mass_matrix (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &boundary_functions, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const coefficient, + const Function::real_type> *const coefficient, std::vector component_mapping) { // what would that be in 1d? the @@ -1678,11 +1676,11 @@ namespace MatrixCreator template void create_boundary_mass_matrix (const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &rhs, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const a, + const Function::real_type> *const a, std::vector component_mapping) { create_boundary_mass_matrix(hp::StaticMappingQ1::mapping_collection, dof, q, diff --git a/include/deal.II/numerics/matrix_tools.h b/include/deal.II/numerics/matrix_tools.h index 06d6ff3bc0..b3fe2c22e9 100644 --- a/include/deal.II/numerics/matrix_tools.h +++ b/include/deal.II/numerics/matrix_tools.h @@ -277,10 +277,10 @@ namespace MatrixCreator void create_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, const ConstraintMatrix &constraints = ConstraintMatrix()); /** @@ -290,10 +290,10 @@ namespace MatrixCreator template void create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, const ConstraintMatrix &constraints = ConstraintMatrix()); /** @@ -324,10 +324,10 @@ namespace MatrixCreator void create_mass_matrix (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, const ConstraintMatrix &constraints = ConstraintMatrix()); /** @@ -336,10 +336,10 @@ namespace MatrixCreator template void create_mass_matrix (const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, const ConstraintMatrix &constraints = ConstraintMatrix()); @@ -372,11 +372,11 @@ namespace MatrixCreator void create_boundary_mass_matrix (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &boundary_functions, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const weight = 0, + const Function::real_type> *const weight = 0, std::vector component_mapping = std::vector()); @@ -387,11 +387,11 @@ namespace MatrixCreator template void create_boundary_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &boundary_functions, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, std::vector component_mapping = std::vector()); /** @@ -401,11 +401,11 @@ namespace MatrixCreator void create_boundary_mass_matrix (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &boundary_functions, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, std::vector component_mapping = std::vector()); /** @@ -414,11 +414,11 @@ namespace MatrixCreator template void create_boundary_mass_matrix (const hp::DoFHandler &dof, const hp::QCollection &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &boundary_functions, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function *const a = nullptr, + const Function::real_type> *const a = nullptr, std::vector component_mapping = std::vector()); /** diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index 65b671e924..a134b95429 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -22,6 +22,7 @@ #include #include #include +#include #include #include #include @@ -37,6 +38,7 @@ #include #include #include +#include #include #include #include @@ -922,7 +924,7 @@ namespace VectorTools } template - void invert_mass_matrix(const SparseMatrix > &/*mass_matrix*/, + void invert_mass_matrix(const SparseMatrix &/*mass_matrix*/, const Vector > &/*rhs*/, Vector > &/*solution*/) { @@ -2854,9 +2856,25 @@ namespace VectorTools return a > b; } + template + bool real_part_bigger_than(const number1 a, + const std::complex b) + { + Assert(std::abs(b.imag()) <= 1e-15*std::abs(b), ExcInternalError()); + return a > b.real(); + } + + template + bool real_part_bigger_than(const std::complex a, + const number2 b) + { + Assert(std::abs(a.imag()) <= 1e-15*std::abs(a), ExcInternalError()); + return a.real() > b; + } + template bool real_part_bigger_than(const std::complex a, - const std::complex &b) + const std::complex b) { Assert(std::abs(a.imag()) <= 1e-15*std::abs(a), ExcInternalError()); Assert(std::abs(b.imag()) <= 1e-15*std::abs(b), ExcInternalError()); @@ -2906,13 +2924,13 @@ namespace VectorTools template void invert_mass_matrix - (const SparseMatrix > &/*mass_matrix*/, - const FilteredMatrix > > &/*filtered_mass_matrix*/, - FilteredMatrix > > &/*filtered_preconditioner*/, - const Vector > &/*rhs*/, - Vector > &/*boundary_projection*/) + (const SparseMatrix &mass_matrix, + const FilteredMatrix > > &filtered_mass_matrix, + FilteredMatrix > > &filtered_preconditioner, + const Vector > &rhs, + Vector > &boundary_projection) { - Assert(false, ExcNotImplemented()); + Assert (false, ExcInternalError()); } @@ -3026,14 +3044,14 @@ namespace VectorTools // make mass matrix and right hand side - SparseMatrix mass_matrix(sparsity); + SparseMatrix::real_type> mass_matrix(sparsity); Vector rhs(sparsity.n_rows()); MatrixCreator::create_boundary_mass_matrix (mapping, dof, q, mass_matrix, boundary_functions, rhs, dof_to_boundary_mapping, - (const Function *) nullptr, + (const Function::real_type> *) nullptr, component_mapping); // For certain weird elements, diff --git a/source/lac/constraint_matrix.cc b/source/lac/constraint_matrix.cc index 6b159a292f..5939293efa 100644 --- a/source/lac/constraint_matrix.cc +++ b/source/lac/constraint_matrix.cc @@ -1324,13 +1324,13 @@ PARALLEL_VECTOR_FUNCTIONS(TrilinosWrappers::MPI::BlockVector); VectorType &, \ bool , \ std::integral_constant) const -#define MATRIX_FUNCTIONS(MatrixType) \ +#define MATRIX_FUNCTIONS(MatrixType,VectorScalar) \ template void ConstraintMatrix:: \ - distribute_local_to_global > (const FullMatrix &, \ - const Vector &, \ + distribute_local_to_global > (const FullMatrix &, \ + const Vector &, \ const std::vector &, \ MatrixType &, \ - Vector &, \ + Vector &, \ bool , \ std::integral_constant) const #define BLOCK_MATRIX_VECTOR_FUNCTIONS(MatrixType, VectorType) \ @@ -1352,37 +1352,42 @@ PARALLEL_VECTOR_FUNCTIONS(TrilinosWrappers::MPI::BlockVector); bool , \ std::integral_constant) const -MATRIX_FUNCTIONS(SparseMatrix); -MATRIX_FUNCTIONS(SparseMatrix); -MATRIX_FUNCTIONS(FullMatrix); -MATRIX_FUNCTIONS(FullMatrix); -MATRIX_FUNCTIONS(FullMatrix >); -MATRIX_FUNCTIONS(SparseMatrix >); -MATRIX_FUNCTIONS(SparseMatrix >); +MATRIX_FUNCTIONS(FullMatrix,double); +MATRIX_FUNCTIONS(FullMatrix,float); +MATRIX_FUNCTIONS(FullMatrix,std::complex); +MATRIX_FUNCTIONS(FullMatrix >,std::complex); + +MATRIX_FUNCTIONS(SparseMatrix,double); +MATRIX_FUNCTIONS(SparseMatrix,float); +MATRIX_FUNCTIONS(SparseMatrix,std::complex); +MATRIX_FUNCTIONS(SparseMatrix,std::complex); +MATRIX_FUNCTIONS(SparseMatrix >,std::complex); +MATRIX_FUNCTIONS(SparseMatrix >,std::complex); + +MATRIX_FUNCTIONS(SparseMatrixEZ,double); +MATRIX_FUNCTIONS(SparseMatrixEZ,float); +MATRIX_FUNCTIONS(ChunkSparseMatrix,double); +MATRIX_FUNCTIONS(ChunkSparseMatrix,float); + BLOCK_MATRIX_FUNCTIONS(BlockSparseMatrix); BLOCK_MATRIX_FUNCTIONS(BlockSparseMatrix); BLOCK_MATRIX_VECTOR_FUNCTIONS(BlockSparseMatrix, BlockVector); BLOCK_MATRIX_VECTOR_FUNCTIONS(BlockSparseMatrix, BlockVector); -MATRIX_FUNCTIONS(SparseMatrixEZ); -MATRIX_FUNCTIONS(SparseMatrixEZ); -MATRIX_FUNCTIONS(ChunkSparseMatrix); -MATRIX_FUNCTIONS(ChunkSparseMatrix); - // BLOCK_MATRIX_FUNCTIONS(BlockSparseMatrixEZ); // BLOCK_MATRIX_VECTOR_FUNCTIONS(BlockSparseMatrixEZ, Vector); #ifdef DEAL_II_WITH_PETSC -MATRIX_FUNCTIONS(PETScWrappers::SparseMatrix); -MATRIX_FUNCTIONS(PETScWrappers::MPI::SparseMatrix); +MATRIX_FUNCTIONS(PETScWrappers::SparseMatrix,PetscScalar); +MATRIX_FUNCTIONS(PETScWrappers::MPI::SparseMatrix,PetscScalar); BLOCK_MATRIX_FUNCTIONS(PETScWrappers::MPI::BlockSparseMatrix); MATRIX_VECTOR_FUNCTIONS(PETScWrappers::MPI::SparseMatrix, PETScWrappers::MPI::Vector); BLOCK_MATRIX_VECTOR_FUNCTIONS(PETScWrappers::MPI::BlockSparseMatrix,PETScWrappers::MPI::BlockVector); #endif #ifdef DEAL_II_WITH_TRILINOS -MATRIX_FUNCTIONS(TrilinosWrappers::SparseMatrix); +MATRIX_FUNCTIONS(TrilinosWrappers::SparseMatrix,double); BLOCK_MATRIX_FUNCTIONS(TrilinosWrappers::BlockSparseMatrix); MATRIX_VECTOR_FUNCTIONS(TrilinosWrappers::SparseMatrix, TrilinosWrappers::MPI::Vector); BLOCK_MATRIX_VECTOR_FUNCTIONS(TrilinosWrappers::BlockSparseMatrix, TrilinosWrappers::MPI::BlockVector); @@ -1468,15 +1473,17 @@ ONLY_MATRIX_FUNCTIONS(PETScWrappers::MPI::BlockSparseMatrix); // constructor of scratch_data (it won't allow one to be constructed in place). namespace internals { -#define SCRATCH_INITIALIZER(Number,Name) \ - ConstraintMatrixData::ScratchData scratch_data_initializer_##Name; \ - template <> Threads::ThreadLocalStorage::ScratchData> \ - ConstraintMatrixData::scratch_data(scratch_data_initializer_##Name) - - SCRATCH_INITIALIZER(double,double); - SCRATCH_INITIALIZER(float,float); - SCRATCH_INITIALIZER(std::complex,cdouble); - SCRATCH_INITIALIZER(std::complex,cfloat); +#define SCRATCH_INITIALIZER(MatrixScalar,VectorScalar,Name) \ + ConstraintMatrixData::ScratchData scratch_data_initializer_##Name; \ + template <> Threads::ThreadLocalStorage::ScratchData> \ + ConstraintMatrixData::scratch_data(scratch_data_initializer_##Name) + + SCRATCH_INITIALIZER(double,double,dd); + SCRATCH_INITIALIZER(float,float,ff); + SCRATCH_INITIALIZER(std::complex,std::complex,zz); + SCRATCH_INITIALIZER(std::complex,std::complex,cc); + SCRATCH_INITIALIZER(double,std::complex,dz); + SCRATCH_INITIALIZER(float,std::complex,fc); #undef SCRATCH_INITIALIZER } diff --git a/source/lac/constraint_matrix.inst.in b/source/lac/constraint_matrix.inst.in index 9244eed2fb..b3c08c4829 100644 --- a/source/lac/constraint_matrix.inst.in +++ b/source/lac/constraint_matrix.inst.in @@ -29,21 +29,57 @@ for (S: REAL_SCALARS; T : DEAL_II_VEC_TEMPLATES) { template void ConstraintMatrix::condense >(const LinearAlgebra::distributed::T &, LinearAlgebra::distributed::T &) const; template void ConstraintMatrix::condense >(LinearAlgebra::distributed::T &vec) const; - template void ConstraintMatrix::distribute_local_to_global > ( - const Vector&, const std::vector &, LinearAlgebra::distributed::T &, const FullMatrix&) const; - template void ConstraintMatrix::distribute_local_to_global > ( - const Vector&, const std::vector &, const std::vector &, LinearAlgebra::distributed::T &, const FullMatrix&, bool) const; - template void ConstraintMatrix::set_zero >(LinearAlgebra::distributed::T &) const; - template void ConstraintMatrix::distribute_local_to_global > > ( - const FullMatrix &, const std::vector< size_type > &, DiagonalMatrix > &) const; - template void ConstraintMatrix::distribute_local_to_global >, LinearAlgebra::distributed::T > ( - const FullMatrix &, const Vector&, const std::vector< size_type > &, - DiagonalMatrix > &, LinearAlgebra::distributed::T&, - bool, std::integral_constant) const; - template void ConstraintMatrix::distribute_local_to_global >, T > ( - const FullMatrix &, const Vector&, const std::vector< size_type > &, - DiagonalMatrix > &, T&, - bool, std::integral_constant) const; + + template + void + ConstraintMatrix::distribute_local_to_global > + (const Vector&, + const std::vector &, + LinearAlgebra::distributed::T &, + const FullMatrix&) const; + + template + void + ConstraintMatrix::distribute_local_to_global > + (const Vector&, + const std::vector &, + const std::vector &, + LinearAlgebra::distributed::T &, + const FullMatrix&, + bool) const; + + template + void + ConstraintMatrix::distribute_local_to_global > > + (const FullMatrix &, + const std::vector< size_type > &, + DiagonalMatrix > &) const; + + template + void + ConstraintMatrix::distribute_local_to_global >, LinearAlgebra::distributed::T > + (const FullMatrix &, + const Vector&, + const std::vector< size_type > &, + DiagonalMatrix > &, + LinearAlgebra::distributed::T&, + bool, + std::integral_constant) const; + + template + void + ConstraintMatrix::distribute_local_to_global >, T > + (const FullMatrix &, + const Vector&, + const std::vector< size_type > &, + DiagonalMatrix > &, + T&, + bool, + std::integral_constant) const; + + template + void + ConstraintMatrix::set_zero >(LinearAlgebra::distributed::T &) const; } diff --git a/source/lac/sparse_matrix.inst.in b/source/lac/sparse_matrix.inst.in index 1a54ee2587..bdb65baec8 100644 --- a/source/lac/sparse_matrix.inst.in +++ b/source/lac/sparse_matrix.inst.in @@ -137,6 +137,19 @@ for (S1, S2, S3 : REAL_SCALARS; Tvmult_add (V1 &, const V2 &) const; } +for (S1 : REAL_SCALARS; S2, S3 : COMPLEX_SCALARS; + V1, V2 : DEAL_II_VEC_TEMPLATES) +{ + template void SparseMatrix:: + vmult (V1 &, const V2 &) const; + template void SparseMatrix:: + Tvmult (V1 &, const V2 &) const; + template void SparseMatrix:: + vmult_add (V1 &, const V2 &) const; + template void SparseMatrix:: + Tvmult_add (V1 &, const V2 &) const; +} + for (S1 : REAL_SCALARS) { template void SparseMatrix:: diff --git a/source/numerics/matrix_creator.inst.in b/source/numerics/matrix_creator.inst.in index 2ddbd7e66d..e5cb415035 100644 --- a/source/numerics/matrix_creator.inst.in +++ b/source/numerics/matrix_creator.inst.in @@ -14,6 +14,8 @@ // --------------------------------------------------------------------- + +// MatrixCreator::create_mass_matrix for the matrix only --> only real-valued scalars for (scalar: REAL_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension <= deal_II_space_dimension @@ -33,29 +35,30 @@ for (scalar: REAL_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimensi SparseMatrix &matrix, const Function * const coefficient, const ConstraintMatrix &constraints); + + // hp versions of functions template void MatrixCreator::create_mass_matrix - (const Mapping &mapping, - const DoFHandler &dof, - const Quadrature &q, + (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const hp::QCollection &q, SparseMatrix &matrix, - const Function &rhs, - Vector &rhs_vector, const Function * const coefficient, const ConstraintMatrix &constraints); + template - void MatrixCreator::create_mass_matrix - (const DoFHandler &dof, - const Quadrature &q, + void MatrixCreator::create_mass_matrix + (const hp::DoFHandler &dof, + const hp::QCollection &q, SparseMatrix &matrix, - const Function &rhs, - Vector &rhs_vector, const Function * const coefficient, const ConstraintMatrix &constraints); #endif } -for (scalar: COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) + +// MatrixCreator::create_mass_matrix for matrix + rhs --> real- and complex-valued scalars +for (scalar: REAL_AND_COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension <= deal_II_space_dimension // non-hp version of create_mass_matrix @@ -64,103 +67,58 @@ for (scalar: COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dime (const Mapping &mapping, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, - const Function * const coefficient, + SparseMatrix::real_type> &matrix, + const Function &rhs, + Vector &rhs_vector, + const Function::real_type> * const coefficient, const ConstraintMatrix &constraints); template void MatrixCreator::create_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, - const Function * const coefficient, + SparseMatrix::real_type> &matrix, + const Function &rhs, + Vector &rhs_vector, + const Function::real_type> * const coefficient, const ConstraintMatrix &constraints); + + // hp version template void MatrixCreator::create_mass_matrix - (const Mapping &mapping, - const DoFHandler &dof, - const Quadrature &q, - SparseMatrix &matrix, + (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const hp::QCollection &q, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function * const coefficient, + const Function::real_type> * const coefficient, const ConstraintMatrix &constraints); + template - void MatrixCreator::create_mass_matrix - (const DoFHandler &dof, - const Quadrature &q, - SparseMatrix &matrix, + void MatrixCreator::create_mass_matrix + (const hp::DoFHandler &dof, + const hp::QCollection &q, + SparseMatrix::real_type> &matrix, const Function &rhs, Vector &rhs_vector, - const Function * const coefficient, + const Function::real_type> * const coefficient, const ConstraintMatrix &constraints); #endif } -for (scalar: REAL_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) -{ -#if deal_II_dimension <= deal_II_space_dimension - template - void MatrixCreator::create_boundary_mass_matrix - (const DoFHandler &dof, - const Quadrature &q, - SparseMatrix &matrix, - const std::map*> &rhs, - Vector &rhs_vector, - std::vector &dof_to_boundary_mapping, - const Function * const a, - std::vector); - - template - void MatrixCreator::create_boundary_mass_matrix - (const Mapping &, - const DoFHandler &dof, - const Quadrature &q, - SparseMatrix &matrix, - const std::map*> &rhs, - Vector &rhs_vector, - std::vector &dof_to_boundary_mapping, - const Function * const a, - std::vector); - - template - void - MatrixCreator::create_boundary_mass_matrix - (const hp::MappingCollection&, - const hp::DoFHandler&, - const hp::QCollection&, - SparseMatrix&, - const std::map*>&, - Vector&, - std::vector&, - const Function * const, - std::vector); - - template - void MatrixCreator::create_boundary_mass_matrix - (const hp::DoFHandler&, - const hp::QCollection&, - SparseMatrix&, - const std::map*>&, - Vector&, - std::vector&, - const Function * const, - std::vector); -#endif -} - -for (scalar: COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +for (scalar: REAL_AND_COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { #if deal_II_dimension <= deal_II_space_dimension template void MatrixCreator::create_boundary_mass_matrix (const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &rhs, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function * const a, + const Function::real_type> * const a, std::vector); template @@ -168,11 +126,11 @@ for (scalar: COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dime (const Mapping &, const DoFHandler &dof, const Quadrature &q, - SparseMatrix &matrix, + SparseMatrix::real_type> &matrix, const std::map*> &rhs, Vector &rhs_vector, std::vector &dof_to_boundary_mapping, - const Function * const a, + const Function::real_type> * const a, std::vector); template @@ -181,127 +139,27 @@ for (scalar: COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dime (const hp::MappingCollection&, const hp::DoFHandler&, const hp::QCollection&, - SparseMatrix&, + SparseMatrix::real_type>&, const std::map*>&, Vector&, std::vector&, - const Function * const, + const Function::real_type> * const, std::vector); template void MatrixCreator::create_boundary_mass_matrix (const hp::DoFHandler&, const hp::QCollection&, - SparseMatrix&, + SparseMatrix::real_type>&, const std::map*>&, Vector&, std::vector&, - const Function * const, + const Function::real_type> * const, std::vector); #endif } -//TODO[SP]: replace by -// where applicable and move to codimension cases above also when applicable -for (scalar : REAL_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) -{ -// hp versions of functions -#if deal_II_dimension <= deal_II_space_dimension - template - void MatrixCreator::create_mass_matrix - (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function * const coefficient, - const ConstraintMatrix &constraints); - - template - void MatrixCreator::create_mass_matrix - (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function &rhs, - Vector &rhs_vector, - const Function * const coefficient, - const ConstraintMatrix &constraints); - -#endif - -#if deal_II_dimension == deal_II_space_dimension - - - template - void MatrixCreator::create_mass_matrix - (const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function * const coefficient, - const ConstraintMatrix &constraints); - - template - void MatrixCreator::create_mass_matrix - (const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function &rhs, - Vector &rhs_vector, - const Function * const coefficient, - const ConstraintMatrix &constraints); - -#endif -} - -for (scalar : COMPLEX_SCALARS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) -{ -// hp versions of functions -#if deal_II_dimension <= deal_II_space_dimension - template - void MatrixCreator::create_mass_matrix - (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function * const coefficient, - const ConstraintMatrix &constraints); - - template - void MatrixCreator::create_mass_matrix - (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function &rhs, - Vector &rhs_vector, - const Function * const coefficient, - const ConstraintMatrix &constraints); - -#endif -#if deal_II_dimension == deal_II_space_dimension - - - template - void MatrixCreator::create_mass_matrix - (const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function * const coefficient, - const ConstraintMatrix &constraints); - - template - void MatrixCreator::create_mass_matrix - (const hp::DoFHandler &dof, - const hp::QCollection &q, - SparseMatrix &matrix, - const Function &rhs, - Vector &rhs_vector, - const Function * const coefficient, - const ConstraintMatrix &constraints); - -#endif -} for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { diff --git a/source/numerics/vector_tools_boundary.inst.in b/source/numerics/vector_tools_boundary.inst.in index dc07a99915..27488e0150 100644 --- a/source/numerics/vector_tools_boundary.inst.in +++ b/source/numerics/vector_tools_boundary.inst.in @@ -109,27 +109,46 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) #endif } -for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; number : REAL_AND_COMPLEX_SCALARS) { - namespace VectorTools \{ #if deal_II_dimension == deal_II_space_dimension + namespace VectorTools \{ template void project_boundary_values (const Mapping &, const DoFHandler &, - const std::map*> &, + const std::map*> &, const Quadrature&, - std::map&, std::vector); + std::map&, + std::vector); template void project_boundary_values (const DoFHandler &, - const std::map*> &, + const std::map*> &, const Quadrature&, - std::map&, + std::map&, std::vector); + template + void project_boundary_values + (const hp::DoFHandler &, + const std::map*> &, + const hp::QCollection &, + std::map &, + std::vector); + + \} +#endif +} + +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) +{ + namespace VectorTools \{ +#if deal_II_dimension == deal_II_space_dimension + + template void project_boundary_values (const Mapping &, @@ -146,16 +165,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) ConstraintMatrix&, std::vector); - template - void project_boundary_values - (const hp::DoFHandler &, - const std::map*> &, - const hp::QCollection &, - std::map &, - std::vector); - - - #if deal_II_dimension != 1 template void project_boundary_values_curl_conforming -- 2.39.5