From: Daniel Arndt Date: Thu, 16 Apr 2020 18:20:26 +0000 (-0400) Subject: Fix petsc/copy_to_dealvec_block X-Git-Tag: v9.2.0-rc1~215^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F9910%2Fhead;p=dealii.git Fix petsc/copy_to_dealvec_block --- diff --git a/include/deal.II/lac/la_parallel_block_vector.templates.h b/include/deal.II/lac/la_parallel_block_vector.templates.h index 59935321b0..7cfff39de5 100644 --- a/include/deal.II/lac/la_parallel_block_vector.templates.h +++ b/include/deal.II/lac/la_parallel_block_vector.templates.h @@ -213,6 +213,38 @@ namespace LinearAlgebra #ifdef DEAL_II_WITH_PETSC + namespace petsc_helpers + { + template + void + copy_petsc_vector(const PETSC_Number *petsc_start_ptr, + const PETSC_Number *petsc_end_ptr, + Number * ptr) + { + std::copy(petsc_start_ptr, petsc_end_ptr, ptr); + } + + template + void + copy_petsc_vector(const std::complex *petsc_start_ptr, + const std::complex *petsc_end_ptr, + std::complex * ptr) + { + std::copy(petsc_start_ptr, petsc_end_ptr, ptr); + } + + template + void + copy_petsc_vector(const std::complex * /*petsc_start_ptr*/, + const std::complex * /*petsc_end_ptr*/, + Number * /*ptr*/) + { + AssertThrow(false, ExcMessage("Tried to copy complex -> real")); + } + } // namespace petsc_helpers + + + template BlockVector & BlockVector:: @@ -221,14 +253,34 @@ namespace LinearAlgebra AssertDimension(this->n_blocks(), petsc_vec.n_blocks()); for (unsigned int i = 0; i < this->n_blocks(); ++i) { - const auto &partitioner = this->block(i).get_partitioner(); - IndexSet combined_set = partitioner->locally_owned_range(); - combined_set.add_indices(partitioner->ghost_indices()); - ReadWriteVector rw_vector(combined_set); - rw_vector.import(petsc_vec.block(i), VectorOperation::insert); - this->block(i).import(rw_vector, VectorOperation::insert); - - if (this->block(i).has_ghost_elements() || + // We would like to use the same compact infrastructure as for the + // Trilinos vector below, but the interface through ReadWriteVector + // does not support overlapping (ghosted) PETSc vectors, which we need + // for backward compatibility. + + Assert(petsc_vec.block(i).locally_owned_elements() == + this->block(i).locally_owned_elements(), + StandardExceptions::ExcInvalidState()); + + // get a representation of the vector and copy it + PetscScalar * start_ptr; + PetscErrorCode ierr = + VecGetArray(static_cast(petsc_vec.block(i)), + &start_ptr); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + + const size_type vec_size = this->block(i).local_size(); + petsc_helpers::copy_petsc_vector(start_ptr, + start_ptr + vec_size, + this->block(i).begin()); + + // restore the representation of the vector + ierr = VecRestoreArray(static_cast(petsc_vec.block(i)), + &start_ptr); + AssertThrow(ierr == 0, ExcPETScError(ierr)); + + // spread ghost values between processes? + if (this->block(i).vector_is_ghosted || petsc_vec.block(i).has_ghost_elements()) this->block(i).update_ghost_values(); } diff --git a/include/deal.II/lac/la_parallel_vector.templates.h b/include/deal.II/lac/la_parallel_vector.templates.h index 0e33633825..ad6d4f41f0 100644 --- a/include/deal.II/lac/la_parallel_vector.templates.h +++ b/include/deal.II/lac/la_parallel_vector.templates.h @@ -749,42 +749,6 @@ namespace LinearAlgebra -#ifdef DEAL_II_WITH_PETSC - - namespace petsc_helpers - { - template - void - copy_petsc_vector(const PETSC_Number *petsc_start_ptr, - const PETSC_Number *petsc_end_ptr, - Number * ptr) - { - std::copy(petsc_start_ptr, petsc_end_ptr, ptr); - } - - template - void - copy_petsc_vector(const std::complex *petsc_start_ptr, - const std::complex *petsc_end_ptr, - std::complex * ptr) - { - std::copy(petsc_start_ptr, petsc_end_ptr, ptr); - } - - template - void - copy_petsc_vector(const std::complex * /*petsc_start_ptr*/, - const std::complex * /*petsc_end_ptr*/, - Number * /*ptr*/) - { - AssertThrow(false, ExcMessage("Tried to copy complex -> real")); - } - } // namespace petsc_helpers - -#endif - - - template void Vector::compress(