16 #ifndef dealii__petsc_vector_base_h 17 #define dealii__petsc_vector_base_h 20 #include <deal.II/base/config.h> 22 #ifdef DEAL_II_WITH_PETSC 24 # include <deal.II/base/subscriptor.h> 25 # include <deal.II/lac/exceptions.h> 26 # include <deal.II/lac/vector.h> 31 # include <petscvec.h> 32 # include <deal.II/base/index_set.h> 34 DEAL_II_NAMESPACE_OPEN
37 template <
typename number>
class Vector;
88 VectorReference (
const VectorBase &vector,
89 const size_type index);
105 const VectorReference &operator = (
const VectorReference &r)
const;
112 VectorReference &operator = (
const VectorReference &r);
117 const VectorReference &operator = (
const PetscScalar &s)
const;
122 const VectorReference &operator += (
const PetscScalar &s)
const;
127 const VectorReference &operator -= (
const PetscScalar &s)
const;
132 const VectorReference &operator *= (
const PetscScalar &s)
const;
137 const VectorReference &operator /= (
const PetscScalar &s)
const;
142 PetscReal real ()
const;
150 PetscReal imag ()
const;
156 operator PetscScalar ()
const;
162 <<
"You tried to access element " << arg1
163 <<
" of a distributed vector, but only elements " 164 << arg2 <<
" through " << arg3
165 <<
" are stored locally and can be accessed.");
171 <<
"You tried to do a " 176 <<
" operation but the vector is currently in " 181 <<
" mode. You first have to call 'compress()'.");
187 const VectorBase &vector;
192 const size_type index;
198 friend class ::PETScWrappers::VectorBase;
239 typedef PetscReal real_type;
241 typedef internal::VectorReference reference;
242 typedef const internal::VectorReference const_reference;
272 virtual void clear ();
318 size_type size ()
const;
328 size_type local_size ()
const;
338 std::pair<size_type, size_type>
339 local_range ()
const;
345 bool in_local_range (
const size_type index)
const;
360 IndexSet locally_owned_elements ()
const;
368 bool has_ghost_elements()
const;
374 operator () (
const size_type index);
380 operator () (
const size_type index)
const;
388 operator [] (
const size_type index);
396 operator [] (
const size_type index)
const;
404 void set (
const std::vector<size_type> &indices,
405 const std::vector<PetscScalar> &values);
413 void extract_subvector_to (
const std::vector<size_type> &indices,
414 std::vector<PetscScalar> &values)
const;
420 template <
typename ForwardIterator,
typename OutputIterator>
421 void extract_subvector_to (
const ForwardIterator indices_begin,
422 const ForwardIterator indices_end,
423 OutputIterator values_begin)
const;
429 void add (
const std::vector<size_type> &indices,
430 const std::vector<PetscScalar> &values);
436 void add (
const std::vector<size_type> &indices,
437 const ::Vector<PetscScalar> &values);
444 void add (
const size_type n_elements,
445 const size_type *indices,
446 const PetscScalar *values);
454 PetscScalar operator * (
const VectorBase &vec)
const;
459 real_type norm_sqr ()
const;
464 PetscScalar mean_value ()
const;
473 real_type l1_norm ()
const;
479 real_type l2_norm ()
const;
485 real_type lp_norm (
const real_type p)
const;
491 real_type linfty_norm ()
const;
508 PetscScalar add_and_dot (
const PetscScalar a,
518 real_type normalize () const DEAL_II_DEPRECATED;
523 real_type min () const;
528 real_type max () const;
575 bool all_zero () const;
582 bool is_non_negative () const;
587 VectorBase &operator *= (const PetscScalar factor);
592 VectorBase &operator /= (const PetscScalar factor);
608 void add (const PetscScalar s);
615 void add (const
VectorBase &V) DEAL_II_DEPRECATED;
620 void add (const PetscScalar a, const
VectorBase &V);
625 void add (const PetscScalar a, const
VectorBase &V,
631 void sadd (const PetscScalar s,
637 void sadd (const PetscScalar s,
646 void sadd (const PetscScalar s,
658 void sadd (const PetscScalar s,
671 void scale (const
VectorBase &scaling_factors);
676 void equ (const PetscScalar a, const
VectorBase &V);
683 void equ (const PetscScalar a, const
VectorBase &V,
684 const PetscScalar b, const
VectorBase &W) DEAL_II_DEPRECATED;
706 void write_ascii (const PetscViewerFormat format = PETSC_VIEWER_DEFAULT) ;
715 void print (
std::ostream &out,
716 const
unsigned int precision = 3,
717 const
bool scientific = true,
718 const
bool across = true) const;
741 operator const Vec &() const;
746 std::
size_t memory_consumption () const;
752 virtual const MPI_Comm &get_mpi_communicator () const;
792 bool attained_ownership;
799 void do_set_add_operation (const size_type n_elements,
800 const size_type *indices,
801 const PetscScalar *values,
802 const
bool add_values);
826 void swap (VectorBase &u, VectorBase &v)
835 VectorReference::VectorReference (
const VectorBase &vector,
836 const size_type index)
844 const VectorReference &
845 VectorReference::operator = (
const VectorReference &r)
const 851 *
this =
static_cast<PetscScalar
> (r);
860 VectorReference::operator = (
const VectorReference &r)
866 *
this =
static_cast<PetscScalar
> (r);
874 const VectorReference &
875 VectorReference::operator = (
const PetscScalar &value)
const 885 const PetscInt petsc_i = index;
887 const PetscErrorCode ierr = VecSetValues (vector, 1, &petsc_i, &value,
899 const VectorReference &
900 VectorReference::operator += (
const PetscScalar &value)
const 919 if (value == PetscScalar())
923 const PetscInt petsc_i = index;
924 const PetscErrorCode ierr = VecSetValues (vector, 1, &petsc_i, &value,
935 const VectorReference &
936 VectorReference::operator -= (
const PetscScalar &value)
const 955 if (value == PetscScalar())
960 const PetscInt petsc_i = index;
961 const PetscScalar subtractand = -value;
962 const PetscErrorCode ierr = VecSetValues (vector, 1, &petsc_i, &subtractand,
972 const VectorReference &
973 VectorReference::operator *= (
const PetscScalar &value)
const 995 const PetscInt petsc_i = index;
996 const PetscScalar new_value
997 =
static_cast<PetscScalar
>(*this) * value;
999 const PetscErrorCode ierr = VecSetValues (vector, 1, &petsc_i, &new_value,
1009 const VectorReference &
1010 VectorReference::operator /= (
const PetscScalar &value)
const 1032 const PetscInt petsc_i = index;
1033 const PetscScalar new_value
1034 =
static_cast<PetscScalar
>(*this) / value;
1036 const PetscErrorCode ierr = VecSetValues (vector, 1, &petsc_i, &new_value,
1047 VectorReference::real ()
const 1049 #ifndef PETSC_USE_COMPLEX 1050 return static_cast<PetscScalar
>(*this);
1052 return PetscRealPart (static_cast<PetscScalar>(*
this));
1060 VectorReference::imag ()
const 1062 #ifndef PETSC_USE_COMPLEX 1063 return PetscReal (0);
1065 return PetscImaginaryPart (static_cast<PetscScalar>(*
this));
1075 PetscInt begin, end;
1076 const PetscErrorCode ierr = VecGetOwnershipRange (static_cast<const Vec &>(vector),
1080 return ((index >= static_cast<size_type>(begin)) &&
1081 (index < static_cast<size_type>(end)));
1092 const std::pair<size_type, size_type> x = local_range();
1109 internal::VectorReference
1112 return internal::VectorReference (*
this, index);
1121 return static_cast<PetscScalar
>(internal::VectorReference (*
this, index));
1127 internal::VectorReference
1130 return operator()(index);
1139 return operator()(index);
1146 static MPI_Comm comm;
1147 PetscObjectGetComm((PetscObject)vector, &comm);
1153 std::vector<PetscScalar> &values)
const 1155 extract_subvector_to(&(indices[0]), &(indices[0]) + indices.size(), &(values[0]));
1158 template <
typename ForwardIterator,
typename OutputIterator>
1161 const ForwardIterator indices_end,
1162 OutputIterator values_begin)
const 1164 const PetscInt n_idx =
static_cast<PetscInt
>(indices_end - indices_begin);
1189 PetscInt begin, end;
1190 PetscErrorCode ierr = VecGetOwnershipRange (vector, &begin, &end);
1193 Vec locally_stored_elements = PETSC_NULL;
1194 ierr = VecGhostGetLocalForm(vector, &locally_stored_elements);
1198 ierr = VecGetSize(locally_stored_elements, &lsize);
1202 ierr = VecGetArray(locally_stored_elements, &ptr);
1205 for (PetscInt i=0; i<n_idx; ++i)
1207 const unsigned int index = *(indices_begin+i);
1208 if ( index>=static_cast<unsigned int>(begin)
1209 && index<static_cast<unsigned int>(end) )
1212 *(values_begin+i) = *(ptr+index-begin);
1217 const unsigned int ghostidx
1218 = ghost_indices.index_within_set(index);
1221 *(values_begin+i) = *(ptr+ghostidx+end-begin);
1225 ierr = VecRestoreArray(locally_stored_elements, &ptr);
1228 ierr = VecGhostRestoreLocalForm(vector, &locally_stored_elements);
1237 PetscInt begin, end;
1238 PetscErrorCode ierr = VecGetOwnershipRange (vector, &begin, &end);
1242 ierr = VecGetArray(vector, &ptr);
1245 for (PetscInt i=0; i<n_idx; ++i)
1247 const unsigned int index = *(indices_begin+i);
1249 Assert(index>=static_cast<unsigned int>(begin)
1252 *(values_begin+i) = *(ptr+index-begin);
1255 ierr = VecRestoreArray(vector, &ptr);
1264 DEAL_II_NAMESPACE_CLOSE
1266 #endif // DEAL_II_WITH_PETSC reference operator[](const size_type index)
types::global_dof_index size_type
#define DeclException2(Exception2, type1, type2, outsequence)
bool has_ghost_elements() const
#define AssertThrow(cond, exc)
virtual const MPI_Comm & get_mpi_communicator() const
unsigned int global_dof_index
VectorOperation::values last_action
#define Assert(cond, exc)
void extract_subvector_to(const std::vector< size_type > &indices, std::vector< PetscScalar > &values) const
static::ExceptionBase & ExcGhostsPresent()
void add_range(const size_type begin, const size_type end)
reference operator()(const size_type index)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
IndexSet locally_owned_elements() const
static::ExceptionBase & ExcInternalError()
bool in_local_range(const size_type index) const