From: bangerth Date: Sat, 12 Jul 2014 12:57:51 +0000 (+0000) Subject: Work around a bug in PETSc whereby PETSc does not zero-initialize ghost elements... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b496a19fdde8945c1a6f4d0fb5f3c80dac9ef6b5;p=dealii-svn.git Work around a bug in PETSc whereby PETSc does not zero-initialize ghost elements of a vector. git-svn-id: https://svn.dealii.org/trunk@33142 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 635384ed35..5b1d5b9809 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -148,6 +148,13 @@ inconvenience this causes.

Specific improvements

    +
  1. Fixed: PETSc up to at least version 3.5 has a bug where it does + not zero-initialize the ghost elements of a newly created ghosted + parallel vector. This is now worked around inside deal.II. +
    + (Wolfgang Bangerth, Michal Wichrowski, 2014/07/12) +
  2. +
  3. Improved: The Trilinos direct solver, TrilinosWrappers::SolverDirect, now takes a string to select among the available Amesos solvers. Moreover, the solver now also supports deal.II's distributed vectors. diff --git a/deal.II/source/lac/petsc_parallel_vector.cc b/deal.II/source/lac/petsc_parallel_vector.cc index ce866b663d..749e464410 100644 --- a/deal.II/source/lac/petsc_parallel_vector.cc +++ b/deal.II/source/lac/petsc_parallel_vector.cc @@ -1,7 +1,7 @@ // --------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2004 - 2013 by the deal.II authors +// Copyright (C) 2004 - 2014 by the deal.II authors // // This file is part of the deal.II library. // @@ -332,29 +332,42 @@ namespace PETScWrappers ExcDimensionMismatch (size(), n)); #if DEBUG - // test ghost allocation in debug mode - PetscInt begin, end; + { + // test ghost allocation in debug mode + PetscInt begin, end; - ierr = VecGetOwnershipRange (vector, &begin, &end); + ierr = VecGetOwnershipRange (vector, &begin, &end); - Assert(local_size==(size_type)(end-begin), ExcInternalError()); + Assert(local_size==(size_type)(end-begin), ExcInternalError()); - Vec l; - ierr = VecGhostGetLocalForm(vector, &l); - AssertThrow (ierr == 0, ExcPETScError(ierr)); - - PetscInt lsize; - ierr = VecGetSize(l, &lsize); - AssertThrow (ierr == 0, ExcPETScError(ierr)); + Vec l; + ierr = VecGhostGetLocalForm(vector, &l); + AssertThrow (ierr == 0, ExcPETScError(ierr)); - ierr = VecGhostRestoreLocalForm(vector, &l); - AssertThrow (ierr == 0, ExcPETScError(ierr)); + PetscInt lsize; + ierr = VecGetSize(l, &lsize); + AssertThrow (ierr == 0, ExcPETScError(ierr)); - Assert( lsize==end-begin+(PetscInt)ghost_indices.n_elements() ,ExcInternalError()); + ierr = VecGhostRestoreLocalForm(vector, &l); + AssertThrow (ierr == 0, ExcPETScError(ierr)); + Assert (lsize==end-begin+(PetscInt)ghost_indices.n_elements(), + ExcInternalError()); + } #endif + // in PETSc versions up to 3.5, VecCreateGhost zeroed out the locally + // owned vector elements but forgot about the ghost elements. we need to + // do this ourselves + // + // see https://code.google.com/p/dealii/issues/detail?id=233 +#if DEAL_II_PETSC_VERSION_LT(3,6,0) + PETScWrappers::MPI::Vector zero; + zero.reinit (communicator, this->size(), local_size); + *this = zero; +#endif + }