Reference documentation for deal.II version 9.3.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
petsc_vector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_petsc_vector_h
17 # define dealii_petsc_vector_h
18 
19 
20 # include <deal.II/base/config.h>
21 
22 # ifdef DEAL_II_WITH_PETSC
23 
24 # include <deal.II/base/index_set.h>
26 
27 # include <deal.II/lac/exceptions.h>
29 # include <deal.II/lac/vector.h>
32 
34 
35 
39 namespace PETScWrappers
40 {
47  namespace MPI
48  {
156  class Vector : public VectorBase
157  {
158  public:
163 
167  Vector();
168 
185  explicit Vector(const MPI_Comm &communicator,
186  const size_type n,
188 
199  template <typename Number>
200  explicit Vector(const MPI_Comm & communicator,
201  const ::Vector<Number> &v,
203 
204 
218  explicit Vector(const MPI_Comm & communicator,
219  const VectorBase &v,
220  const size_type local_size);
221 
245  Vector(const IndexSet &local,
246  const IndexSet &ghost,
247  const MPI_Comm &communicator);
248 
260  explicit Vector(const IndexSet &local, const MPI_Comm &communicator);
261 
265  Vector(const Vector &v);
266 
271  virtual void
272  clear() override;
273 
278  Vector &
279  operator=(const Vector &v);
280 
287  Vector &
288  operator=(const PetscScalar s);
289 
299  template <typename number>
300  Vector &
301  operator=(const ::Vector<number> &v);
302 
319  void
321  const size_type N,
323  const bool omit_zeroing_entries = false);
324 
334  void
335  reinit(const Vector &v, const bool omit_zeroing_entries = false);
336 
344  void
345  reinit(const IndexSet &local,
346  const IndexSet &ghost,
347  const MPI_Comm &communicator);
348 
356  void
357  reinit(const IndexSet &local, const MPI_Comm &communicator);
358 
363  const MPI_Comm &
364  get_mpi_communicator() const override;
365 
377  void
378  print(std::ostream & out,
379  const unsigned int precision = 3,
380  const bool scientific = true,
381  const bool across = true) const;
382 
389  bool
390  all_zero() const;
391 
392  protected:
399  virtual void
401 
402 
403 
409  virtual void
410  create_vector(const size_type n,
412  const IndexSet &ghostnodes);
413 
414 
415  private:
420  };
421 
422 
423  // ------------------ template and inline functions -------------
424 
425 
433  inline void
435  {
436  u.swap(v);
437  }
438 
439 
440 # ifndef DOXYGEN
441 
442  template <typename number>
444  const ::Vector<number> &v,
446  : communicator(communicator)
447  {
449 
450  *this = v;
451  }
452 
453 
454 
455  inline Vector &
456  Vector::operator=(const PetscScalar s)
457  {
459 
460  return *this;
461  }
462 
463 
464 
465  template <typename number>
466  inline Vector &
467  Vector::operator=(const ::Vector<number> &v)
468  {
469  Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
470 
471  // FIXME: the following isn't necessarily fast, but this is due to
472  // the fact that PETSc doesn't offer an inlined access operator.
473  //
474  // if someone wants to contribute some code: to make this code
475  // faster, one could either first convert all values to PetscScalar,
476  // and then set them all at once using VecSetValues. This has the
477  // drawback that it could take quite some memory, if the vector is
478  // large, and it would in addition allocate memory on the heap, which
479  // is expensive. an alternative would be to split the vector into
480  // chunks of, say, 128 elements, convert a chunk at a time and set it
481  // in the output vector using VecSetValues. since 128 elements is
482  // small enough, this could easily be allocated on the stack (as a
483  // local variable) which would make the whole thing much more
484  // efficient.
485  //
486  // a second way to make things faster is for the special case that
487  // number==PetscScalar. we could then declare a specialization of
488  // this template, and omit the conversion. the problem with this is
489  // that the best we can do is to use VecSetValues, but this isn't
490  // very efficient either: it wants to see an array of indices, which
491  // in this case a) again takes up a whole lot of memory on the heap,
492  // and b) is totally dumb since its content would simply be the
493  // sequence 0,1,2,3,...,n. the best of all worlds would probably be a
494  // function in Petsc that would take a pointer to an array of
495  // PetscScalar values and simply copy n elements verbatim into the
496  // vector...
497  for (size_type i = 0; i < v.size(); ++i)
498  (*this)(i) = v(i);
499 
501 
502  return *this;
503  }
504 
505 
506 
507  inline const MPI_Comm &
509  {
510  return communicator;
511  }
512 
513 # endif // DOXYGEN
514  } // namespace MPI
515 } // namespace PETScWrappers
516 
517 namespace internal
518 {
519  namespace LinearOperatorImplementation
520  {
521  template <typename>
522  class ReinitHelper;
523 
528  template <>
530  {
531  public:
532  template <typename Matrix>
533  static void
534  reinit_range_vector(const Matrix & matrix,
536  bool /*omit_zeroing_entries*/)
537  {
538  v.reinit(matrix.locally_owned_range_indices(),
539  matrix.get_mpi_communicator());
540  }
541 
542  template <typename Matrix>
543  static void
544  reinit_domain_vector(const Matrix & matrix,
546  bool /*omit_zeroing_entries*/)
547  {
548  v.reinit(matrix.locally_owned_domain_indices(),
549  matrix.get_mpi_communicator());
550  }
551  };
552 
553  } // namespace LinearOperatorImplementation
554 } /* namespace internal */
555 
562 template <>
563 struct is_serial_vector<PETScWrappers::MPI::Vector> : std::false_type
564 {};
565 
566 
568 
569 # endif // DEAL_II_WITH_PETSC
570 
571 #endif
572 /*------------------------- petsc_vector.h -------------------------*/
virtual void create_vector(const size_type n, const size_type locally_owned_size)
Contents is actually a matrix.
size_type locally_owned_size() const
static void reinit_range_vector(const Matrix &matrix, PETScWrappers::MPI::Vector &v, bool)
Definition: petsc_vector.h:534
virtual void clear() override
VectorBase & operator=(const VectorBase &)=delete
void reinit(const MPI_Comm &communicator, const size_type N, const size_type locally_owned_size, const bool omit_zeroing_entries=false)
void compress(const VectorOperation::values operation)
static void reinit_domain_vector(const Matrix &matrix, PETScWrappers::MPI::Vector &v, bool)
Definition: petsc_vector.h:544
#define Assert(cond, exc)
Definition: exceptions.h:1465
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
Vector & operator=(const Vector &v)
unsigned int global_dof_index
Definition: types.h:76
const MPI_Comm & get_mpi_communicator() const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
static const char N
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
#define DEAL_II_DEPRECATED
Definition: config.h:157
void swap(Vector &u, Vector &v)
Definition: petsc_vector.h:434