From: Wolfgang Bangerth Date: Mon, 30 Oct 2017 23:06:02 +0000 (-0600) Subject: Instantiate functions and classes also for complex-valued vectors. X-Git-Tag: v9.0.0-rc1~821^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=48e9efd18180b74e40df0f55d455bcce0b37b56a;p=dealii.git Instantiate functions and classes also for complex-valued vectors. --- diff --git a/cmake/config/template-arguments.in b/cmake/config/template-arguments.in index a07bc908c2..657dde43a2 100644 --- a/cmake/config/template-arguments.in +++ b/cmake/config/template-arguments.in @@ -12,6 +12,12 @@ REAL_SCALARS := { double; float } COMPLEX_SCALARS := { std::complex; std::complex } +// real scalar floating point types +REAL_AND_COMPLEX_SCALARS := { double; + float; + std::complex; + std::complex } + // differentiable scalar types DIFFERENTIABLE_REAL_SCALARS := { @DEAL_II_EXPAND_TRILINOS_SACADO_TYPES@; @@ -32,29 +38,46 @@ MPI_SCALARS := { int; unsigned long long int; float; double; - long double } + long double; + std::complex; + std::complex } // template names for serial vectors that we can instantiate as T where // S=REAL_SCALARS for example DEAL_II_VEC_TEMPLATES := { Vector; BlockVector } -// Serial vector types +// Serial vector types, based on real or complex scalars // TODO: why are parallel vectors in here? SERIAL_VECTORS := { Vector; Vector ; + Vector >; + Vector >; + BlockVector; BlockVector; + BlockVector >; + BlockVector >; + LinearAlgebra::Vector; LinearAlgebra::Vector ; + LinearAlgebra::Vector >; + LinearAlgebra::Vector > ; + LinearAlgebra::distributed::Vector; LinearAlgebra::distributed::Vector ; + LinearAlgebra::distributed::Vector >; + LinearAlgebra::distributed::Vector > ; + LinearAlgebra::distributed::BlockVector; LinearAlgebra::distributed::BlockVector ; + LinearAlgebra::distributed::BlockVector >; + LinearAlgebra::distributed::BlockVector > ; + @DEAL_II_EXPAND_TRILINOS_MPI_VECTOR@; @DEAL_II_EXPAND_EPETRA_VECTOR@; @DEAL_II_EXPAND_PETSC_MPI_VECTOR@; @@ -63,7 +86,8 @@ SERIAL_VECTORS := { Vector; @DEAL_II_EXPAND_PETSC_MPI_BLOCKVECTOR@; } -// same as SERIAL_VECTORS but only with real-valued PETSc vectors +// same as SERIAL_VECTORS but only real-valued vectors (and only with PETSc +// vectors if the PETScScalar data type is real-valued) REAL_SERIAL_VECTORS := { Vector; Vector ; diff --git a/source/base/mpi.inst.in b/source/base/mpi.inst.in index bbbc0b8807..c96e796efd 100644 --- a/source/base/mpi.inst.in +++ b/source/base/mpi.inst.in @@ -13,14 +13,12 @@ // // --------------------------------------------------------------------- -for (S : REAL_SCALARS) + +for (S : MPI_SCALARS) { template void sum (const LAPACKFullMatrix &, const MPI_Comm &, LAPACKFullMatrix &); -} -for (S : MPI_SCALARS) -{ template void sum (const Vector &, const MPI_Comm &, Vector &); @@ -70,24 +68,6 @@ for (S : MPI_SCALARS) } -for (S : COMPLEX_SCALARS) -{ - template - void sum (const Vector &, const MPI_Comm &, Vector &); - - template - void sum (const FullMatrix &, const MPI_Comm &, FullMatrix &); - - template - S sum (const S &, const MPI_Comm &); - - template - void sum (const std::vector &, const MPI_Comm &, std::vector &); - - template - void sum (const ArrayView &, const MPI_Comm &, const ArrayView &); -} - for (S : REAL_SCALARS; rank: RANKS; dim : SPACE_DIMENSIONS) { template diff --git a/source/base/partitioner.inst.in b/source/base/partitioner.inst.in index 111d4ed477..3ebb55585c 100644 --- a/source/base/partitioner.inst.in +++ b/source/base/partitioner.inst.in @@ -37,27 +37,3 @@ for (SCALAR : MPI_SCALARS) std::vector &) const; #endif } - - -for (SCALAR : COMPLEX_SCALARS) -{ -#ifdef DEAL_II_WITH_MPI - template void Utilities::MPI::Partitioner::export_to_ghosted_array_start(const unsigned int , - const ArrayView &, - const ArrayView &, - const ArrayView &, - std::vector &) const; - template void Utilities::MPI::Partitioner::export_to_ghosted_array_finish(const ArrayView &, - std::vector &) const; - template void Utilities::MPI::Partitioner::import_from_ghosted_array_start(const VectorOperation::values , - const unsigned int , - const ArrayView &, - const ArrayView &, - std::vector &) const; - template void Utilities::MPI::Partitioner::import_from_ghosted_array_finish(const VectorOperation::values , - const ArrayView &, - const ArrayView &, - const ArrayView &, - std::vector &) const; -#endif -} diff --git a/source/lac/block_vector.inst.in b/source/lac/block_vector.inst.in index b6d45a48e1..ba0f278ffb 100644 --- a/source/lac/block_vector.inst.in +++ b/source/lac/block_vector.inst.in @@ -15,7 +15,7 @@ -for (S : REAL_SCALARS) +for (S : REAL_AND_COMPLEX_SCALARS) { template class BlockVector; } @@ -26,13 +26,14 @@ for (S1, S2 : REAL_SCALARS) const bool); } - -for (S : COMPLEX_SCALARS) +for (S1, S2 : COMPLEX_SCALARS) { - template class BlockVector; + template void BlockVector::reinit(const BlockVector&, + const bool); } -for (S1, S2 : COMPLEX_SCALARS) + +for (S1 : COMPLEX_SCALARS; S2 : REAL_SCALARS) { template void BlockVector::reinit(const BlockVector&, const bool); diff --git a/source/lac/constraint_matrix.inst.in b/source/lac/constraint_matrix.inst.in index 1e62dd37c0..b99e631eb9 100644 --- a/source/lac/constraint_matrix.inst.in +++ b/source/lac/constraint_matrix.inst.in @@ -76,8 +76,3 @@ for (Vec : SERIAL_VECTORS) { template void ConstraintMatrix::distribute(Vec &) const; } - -for (S: COMPLEX_SCALARS; T : DEAL_II_VEC_TEMPLATES) -{ - template void ConstraintMatrix::distribute >(T &) const; -} diff --git a/source/lac/full_matrix.cc b/source/lac/full_matrix.cc index 4d6b9a8a0d..7cecdfcbc9 100644 --- a/source/lac/full_matrix.cc +++ b/source/lac/full_matrix.cc @@ -32,9 +32,6 @@ template void FullMatrix::vmult(Vector &, template void FullMatrix::Tvmult(Vector &, const Vector &, bool) const; template void FullMatrix::add (const long double, const FullMatrix &); -// This is needed if PETSc was compiled with complex, though, it may -// be used elsewhere too. -template void dealii::FullMatrix::vmult >(dealii::Vector > &, dealii::Vector > const &, bool) const; // do a few functions that currently don't fit the scheme because they have // two template arguments that need to be different (the case of same diff --git a/source/lac/full_matrix.inst.in b/source/lac/full_matrix.inst.in index 68716b0d74..65201f6582 100644 --- a/source/lac/full_matrix.inst.in +++ b/source/lac/full_matrix.inst.in @@ -15,7 +15,7 @@ -for (S : REAL_SCALARS) +for (S : REAL_AND_COMPLEX_SCALARS) { template class FullMatrix; @@ -32,7 +32,11 @@ for (S : REAL_SCALARS) template void FullMatrix::copy_from<3>( const Tensor<2,3>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type); +} + +for (S : REAL_SCALARS) +{ template void FullMatrix::copy_to<1>( Tensor<2,1>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type) const; @@ -43,6 +47,7 @@ for (S : REAL_SCALARS) Tensor<2,3>&, const size_type, const size_type, const size_type, const size_type, const size_type, const size_type) const; } + for (S1, S2 : REAL_SCALARS) { template @@ -98,6 +103,19 @@ for (S1, S2 : REAL_SCALARS) const FullMatrix&, const std::vector&, const std::vector&); + + template + void FullMatrix::cholesky (const FullMatrix&); + + template + void FullMatrix::outer_product (const Vector&, + const Vector&); +} + + +// real matrices can be multiplied by real or complex vectors +for (S1 : REAL_SCALARS; S2 : REAL_AND_COMPLEX_SCALARS) +{ template void FullMatrix::vmult( Vector&, const Vector&, bool) const; @@ -120,15 +138,38 @@ for (S1, S2 : REAL_SCALARS) template void FullMatrix::precondition_Jacobi ( Vector &, const Vector &, const S1) const; +} + +// complex matrices can be multiplied only by complex vectors +for (S1 : COMPLEX_SCALARS; S2 : COMPLEX_SCALARS) +{ template - void FullMatrix::cholesky (const FullMatrix&); + void FullMatrix::vmult( + Vector&, const Vector&, bool) const; + template + void FullMatrix::Tvmult( + Vector&, const Vector&, bool) const; + template + S2 FullMatrix::matrix_norm_square ( + const Vector &) const; + template + S2 FullMatrix::matrix_scalar_product( + const Vector&, const Vector&) const; + template + void FullMatrix::forward( + Vector&, const Vector&) const; + template + void FullMatrix::backward( + Vector&, const Vector&) const; template - void FullMatrix::outer_product (const Vector&, - const Vector&); + void FullMatrix::precondition_Jacobi ( + Vector &, const Vector &, const S1) const; } + + for (S1, S2, S3 : REAL_SCALARS) { template @@ -140,20 +181,6 @@ for (S1, S2, S3 : REAL_SCALARS) - - -// same for complex scalars - -for (S : COMPLEX_SCALARS) -{ - template class FullMatrix; - - template void FullMatrix::print( - LogStream&, const unsigned int, const unsigned int) const; - template void FullMatrix::print( - std::ostream&, const unsigned int, const unsigned int) const; -} - for (S1, S2 : COMPLEX_SCALARS) { template @@ -203,28 +230,6 @@ for (S1, S2 : COMPLEX_SCALARS) const FullMatrix&, const std::vector&, const std::vector&); - template - void FullMatrix::vmult( - Vector&, const Vector&, bool) const; - template - void FullMatrix::Tvmult( - Vector&, const Vector&, bool) const; - template - S2 FullMatrix::matrix_norm_square ( - const Vector &) const; - template - S2 FullMatrix::matrix_scalar_product( - const Vector&, const Vector&) const; - template - void FullMatrix::forward( - Vector&, const Vector&) const; - template - void FullMatrix::backward( - Vector&, const Vector&) const; - - template - void FullMatrix::precondition_Jacobi ( - Vector &, const Vector &, const S1) const; } for (S1, S2, S3 : COMPLEX_SCALARS) diff --git a/source/lac/la_parallel_block_vector.inst.in b/source/lac/la_parallel_block_vector.inst.in index 41ecec5924..7d8f41cd7c 100644 --- a/source/lac/la_parallel_block_vector.inst.in +++ b/source/lac/la_parallel_block_vector.inst.in @@ -15,7 +15,7 @@ -for (SCALAR : REAL_SCALARS) +for (SCALAR : REAL_AND_COMPLEX_SCALARS) { namespace LinearAlgebra \{ @@ -32,7 +32,23 @@ for (SCALAR : REAL_SCALARS) \} } -for (S1, S2 : REAL_SCALARS) +for (S1 : REAL_AND_COMPLEX_SCALARS; S2 : REAL_SCALARS) +{ + namespace LinearAlgebra + \{ + namespace distributed + \{ + template void BlockVector::reinit (const BlockVector&, + const bool); + template void BlockVector::add (const std::vector &, + const ::dealii::Vector&); + \} + \} +} + + + +for (S1, S2 : COMPLEX_SCALARS) { namespace LinearAlgebra \{ diff --git a/source/lac/la_parallel_vector.inst.in b/source/lac/la_parallel_vector.inst.in index cdb4b357b0..fcd3d32e10 100644 --- a/source/lac/la_parallel_vector.inst.in +++ b/source/lac/la_parallel_vector.inst.in @@ -15,7 +15,7 @@ -for (SCALAR : REAL_SCALARS) +for (SCALAR : REAL_AND_COMPLEX_SCALARS) { namespace LinearAlgebra \{ @@ -26,7 +26,7 @@ for (SCALAR : REAL_SCALARS) \} } -for (S1, S2 : REAL_SCALARS) +for (S1 : REAL_AND_COMPLEX_SCALARS; S2 : REAL_SCALARS) { namespace LinearAlgebra \{ @@ -34,22 +34,12 @@ for (S1, S2 : REAL_SCALARS) \{ template void Vector::reinit (const Vector&, const bool); - template void Vector::copy_locally_owned_data_from (const Vector&); template S1 Vector::inner_product_local (const Vector&) const; + template void Vector::copy_locally_owned_data_from (const Vector&); \} \} } -for (SCALAR : COMPLEX_SCALARS) -{ - namespace LinearAlgebra - \{ - namespace distributed - \{ - template class Vector; - \} - \} -} for (S1, S2 : COMPLEX_SCALARS) { @@ -59,7 +49,9 @@ for (S1, S2 : COMPLEX_SCALARS) \{ template void Vector::reinit (const Vector&, const bool); + template S1 Vector::inner_product_local (const Vector&) const; template void Vector::copy_locally_owned_data_from (const Vector&); \} \} } + diff --git a/source/lac/lapack_full_matrix.inst.in b/source/lac/lapack_full_matrix.inst.in index 5277991229..fbc567da81 100644 --- a/source/lac/lapack_full_matrix.inst.in +++ b/source/lac/lapack_full_matrix.inst.in @@ -15,13 +15,13 @@ -for (S : REAL_SCALARS) +for (S : REAL_AND_COMPLEX_SCALARS) { template class LAPACKFullMatrix; template class PreconditionLU; } -for (S1, S2 : REAL_SCALARS) +for (S1 : REAL_AND_COMPLEX_SCALARS; S2 : REAL_SCALARS) { template LAPACKFullMatrix & LAPACKFullMatrix::operator = (const FullMatrix &M); diff --git a/source/lac/vector.inst.in b/source/lac/vector.inst.in index 25df5d3658..6588fbfb01 100644 --- a/source/lac/vector.inst.in +++ b/source/lac/vector.inst.in @@ -15,7 +15,7 @@ -for (SCALAR : REAL_SCALARS) +for (SCALAR : REAL_AND_COMPLEX_SCALARS) { template class Vector; } @@ -32,13 +32,6 @@ for (S1, S2 : REAL_SCALARS) void Vector::reinit(const Vector&, const bool); } - - -for (SCALAR : COMPLEX_SCALARS) -{ - template class Vector; -} - for (S1, S2 : COMPLEX_SCALARS) { template @@ -51,8 +44,10 @@ for (S1, S2 : COMPLEX_SCALARS) void Vector::reinit(const Vector&, const bool); } -for (S1: REAL_SCALARS; S2: COMPLEX_SCALARS) +for (S1: COMPLEX_SCALARS; S2: REAL_SCALARS) { template - Vector& Vector::operator= (const Vector &); + Vector& Vector::operator= (const Vector &); + template + void Vector::reinit(const Vector&, const bool); } diff --git a/source/lac/vector_memory.inst.in b/source/lac/vector_memory.inst.in index 75fd686f20..28f51735c9 100644 --- a/source/lac/vector_memory.inst.in +++ b/source/lac/vector_memory.inst.in @@ -20,12 +20,3 @@ for (VECTOR : SERIAL_VECTORS) template class VectorMemory; template class GrowingVectorMemory; } - -for (SCALAR : COMPLEX_SCALARS) -{ - template class VectorMemory >; - template class GrowingVectorMemory >; - - template class VectorMemory >; - template class GrowingVectorMemory >; -} diff --git a/source/lac/vector_memory_release.inst.in b/source/lac/vector_memory_release.inst.in index 0acf5031c8..cbf8653ae9 100644 --- a/source/lac/vector_memory_release.inst.in +++ b/source/lac/vector_memory_release.inst.in @@ -19,9 +19,3 @@ for (VECTOR : SERIAL_VECTORS) { dealii::GrowingVectorMemory::release_unused_memory(); } - -for (SCALAR : COMPLEX_SCALARS) -{ - dealii::GrowingVectorMemory >::release_unused_memory(); - dealii::GrowingVectorMemory >::release_unused_memory(); -} diff --git a/source/numerics/data_out_dof_data.inst.in b/source/numerics/data_out_dof_data.inst.in index 8ceb65aacb..be26c154ed 100644 --- a/source/numerics/data_out_dof_data.inst.in +++ b/source/numerics/data_out_dof_data.inst.in @@ -14,7 +14,7 @@ // --------------------------------------------------------------------- -for (VEC : REAL_SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS) +for (VEC : SERIAL_VECTORS; DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS) { template void DataOut_DoFData,deal_II_dimension,deal_II_dimension>::