18#ifdef DEAL_II_WITH_PETSC
36 VectorReference::operator PetscScalar()
const
50 VecGetOwnershipRange(vector.vector, &begin, &end);
53 Vec locally_stored_elements = PETSC_NULL;
54 ierr = VecGhostGetLocalForm(vector.vector, &locally_stored_elements);
58 ierr = VecGetSize(locally_stored_elements, &lsize);
62 ierr = VecGetArray(locally_stored_elements, &ptr);
67 if (index >=
static_cast<size_type
>(begin) &&
68 index <
static_cast<size_type
>(end))
71 value = *(ptr + index - begin);
76 Assert(vector.ghost_indices.is_element(index),
78 "You are trying to access an element of a vector "
79 "that is neither a locally owned element nor a "
80 "ghost element of the vector."));
81 const size_type ghostidx =
82 vector.ghost_indices.index_within_set(index);
85 value = *(ptr + ghostidx + end - begin);
88 ierr = VecRestoreArray(locally_stored_elements, &ptr);
92 VecGhostRestoreLocalForm(vector.vector, &locally_stored_elements);
104 PetscErrorCode ierr = VecGetOwnershipRange(vector.vector, &begin, &end);
107 AssertThrow((index >=
static_cast<size_type
>(begin)) &&
108 (index <
static_cast<size_type
>(end)),
109 ExcAccessToNonlocalElement(index, begin, end - 1));
111 PetscInt idx = index;
113 ierr = VecGetValues(vector.vector, 1, &idx, &value);
125 , obtained_ownership(true)
128 ExcMessage(
"PETSc does not support multi-threaded access, set "
129 "the thread limit to 1 in MPI_InitFinalize()."));
137 , ghost_indices(v.ghost_indices)
139 , obtained_ownership(true)
142 ExcMessage(
"PETSc does not support multi-threaded access, set "
143 "the thread limit to 1 in MPI_InitFinalize()."));
159 , obtained_ownership(false)
162 ExcMessage(
"PETSc does not support multi-threaded access, set "
163 "the thread limit to 1 in MPI_InitFinalize()."));
172 const PetscErrorCode ierr = VecDestroy(&
vector);
185 const PetscErrorCode ierr = VecDestroy(&
vector);
204 PetscErrorCode ierr = VecSet(
vector, s);
209 Vec ghost = PETSC_NULL;
210 ierr = VecGhostGetLocalForm(
vector, &ghost);
213 ierr = VecSet(ghost, s);
216 ierr = VecGhostRestoreLocalForm(
vector, &ghost);
231 const PetscErrorCode ierr = VecEqual(
vector, v.
vector, &flag);
234 return (flag == PETSC_TRUE);
245 const PetscErrorCode ierr = VecEqual(
vector, v.
vector, &flag);
248 return (flag == PETSC_FALSE);
257 const PetscErrorCode ierr = VecGetSize(
vector, &sz);
269 const PetscErrorCode ierr = VecGetLocalSize(
vector, &sz);
281 const PetscErrorCode ierr = VecGetLocalSize(
vector, &sz);
289 std::pair<VectorBase::size_type, VectorBase::size_type>
293 const PetscErrorCode ierr =
294 VecGetOwnershipRange(
static_cast<const Vec &
>(
vector), &begin, &end);
297 return std::make_pair(begin, end);
304 const std::vector<PetscScalar> &values)
306 Assert(indices.size() == values.size(),
307 ExcMessage(
"Function called with arguments of different sizes"));
315 const std::vector<PetscScalar> &values)
317 Assert(indices.size() == values.size(),
318 ExcMessage(
"Function called with arguments of different sizes"));
326 const ::Vector<PetscScalar> &values)
328 Assert(indices.size() == values.size(),
329 ExcMessage(
"Function called with arguments of different sizes"));
338 const PetscScalar *values)
359 const PetscErrorCode ierr = VecDot(vec.
vector,
vector, &result);
383# ifdef DEAL_II_WITH_MPI
387 int all_int_last_action;
389 const int ierr = MPI_Allreduce(&my_int_last_action,
390 &all_int_last_action,
399 ExcMessage(
"Error: not all processors agree on the last "
400 "VectorOperation before this compress() call."));
409 "Missing compress() or calling with wrong VectorOperation argument."));
425 PetscErrorCode ierr = VecAssemblyBegin(
vector);
427 ierr = VecAssemblyEnd(
vector);
455 const PetscErrorCode ierr = VecSum(
vector, &sum);
457 return sum /
static_cast<PetscReal
>(
size());
462 PetscScalar * start_ptr;
463 PetscErrorCode ierr = VecGetArray(
vector, &start_ptr);
466 PetscScalar mean = 0;
468 PetscScalar sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
473 const PetscScalar *ptr = start_ptr;
474 const PetscScalar *eptr = ptr + (
size() / 4) * 4;
483 while (ptr != start_ptr +
size())
486 mean = (sum0 + sum1 + sum2 + sum3) /
static_cast<PetscReal
>(
size());
491 ierr = VecRestoreArray(
vector, &start_ptr);
503 const PetscErrorCode ierr = VecNorm(
vector, NORM_1, &d);
516 const PetscErrorCode ierr = VecNorm(
vector, NORM_2, &d);
529 PetscScalar * start_ptr;
530 PetscErrorCode ierr = VecGetArray(
vector, &start_ptr);
535 real_type sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
540 const PetscScalar *ptr = start_ptr;
541 const PetscScalar *eptr = ptr + (
size() / 4) * 4;
550 while (ptr != start_ptr +
size())
553 norm =
std::pow(sum0 + sum1 + sum2 + sum3, 1. / p);
558 ierr = VecRestoreArray(
vector, &start_ptr);
571 const PetscErrorCode ierr = VecNorm(
vector, NORM_INFINITY, &d);
585 const PetscErrorCode ierr = VecMin(
vector, &p, &d);
598 const PetscErrorCode ierr = VecMax(
vector, &p, &d);
611 PetscScalar * start_ptr;
612 PetscErrorCode ierr = VecGetArray(
vector, &start_ptr);
615 const PetscScalar *ptr = start_ptr,
630 ierr = VecRestoreArray(
vector, &start_ptr);
639 template <
typename T>
648 template <
typename T>
654 "whether it is non-negative."))
return true;
665 PetscScalar * start_ptr;
666 PetscErrorCode ierr = VecGetArray(
vector, &start_ptr);
669 const PetscScalar *ptr = start_ptr,
684 ierr = VecRestoreArray(
vector, &start_ptr);
698 const PetscErrorCode ierr = VecScale(
vector, a);
712 const PetscScalar factor = 1. / a;
715 const PetscErrorCode ierr = VecScale(
vector, factor);
727 const PetscErrorCode ierr = VecAXPY(
vector, 1, v);
739 const PetscErrorCode ierr = VecAXPY(
vector, -1, v);
753 const PetscErrorCode ierr = VecShift(
vector, s);
765 const PetscErrorCode ierr = VecAXPY(
vector, a, v);
781 const PetscScalar weights[2] = {a, b};
782 Vec addends[2] = {v.
vector, w.vector};
784 const PetscErrorCode ierr = VecMAXPY(
vector, 2, weights, addends);
796 const PetscErrorCode ierr = VecAYPX(
vector, s, v);
824 const PetscErrorCode ierr = VecPointwiseMult(
vector, factors,
vector);
855 PetscErrorCode ierr =
856 PetscViewerSetFormat(PETSC_VIEWER_STDOUT_WORLD, format);
860 ierr = VecView(
vector, PETSC_VIEWER_STDOUT_WORLD);
868 const unsigned int precision,
869 const bool scientific,
870 const bool across)
const
877 PetscErrorCode ierr = VecGetArray(
vector, &val);
882 const std::ios::fmtflags old_flags = out.flags();
883 const unsigned int old_precision = out.precision(precision);
885 out.precision(precision);
887 out.setf(std::ios::scientific, std::ios::floatfield);
889 out.setf(std::ios::fixed, std::ios::floatfield);
893 out << val[i] <<
' ';
896 out << val[i] << std::endl;
900 out.flags(old_flags);
901 out.precision(old_precision);
905 ierr = VecRestoreArray(
vector, &val);
922 VectorBase::operator
const Vec &()
const
931 std::size_t mem =
sizeof(Vec) +
sizeof(
last_action) +
953 const PetscScalar *values,
954 const bool add_values)
961 internal::VectorReference::ExcWrongMode(action,
last_action));
971 const PetscInt *petsc_indices =
972 reinterpret_cast<const PetscInt *
>(indices);
974 const InsertMode mode = (add_values ? ADD_VALUES : INSERT_VALUES);
975 const PetscErrorCode ierr =
976 VecSetValues(
vector, n_elements, petsc_indices, values, mode);
size_type n_elements() const
static bool is_running_single_threaded()
real_type lp_norm(const real_type p) const
real_type l1_norm() const
VectorBase & operator+=(const VectorBase &V)
PetscScalar mean_value() const
size_type local_size() const
VectorOperation::values last_action
PetscScalar operator*(const VectorBase &vec) const
bool operator==(const VectorBase &v) const
VectorBase & operator*=(const PetscScalar factor)
VectorBase & operator-=(const VectorBase &V)
bool operator!=(const VectorBase &v) const
std::size_t memory_consumption() const
VectorBase & operator/=(const PetscScalar factor)
bool is_non_negative() const
virtual const MPI_Comm & get_mpi_communicator() const
void scale(const VectorBase &scaling_factors)
std::pair< size_type, size_type > local_range() const
void sadd(const PetscScalar s, const VectorBase &V)
real_type norm_sqr() const
void compress(const VectorOperation::values operation)
void set(const std::vector< size_type > &indices, const std::vector< PetscScalar > &values)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
PetscScalar add_and_dot(const PetscScalar a, const VectorBase &V, const VectorBase &W)
virtual ~VectorBase() override
real_type linfty_norm() const
void add(const std::vector< size_type > &indices, const std::vector< PetscScalar > &values)
void write_ascii(const PetscViewerFormat format=PETSC_VIEWER_DEFAULT)
bool has_ghost_elements() const
size_type locally_owned_size() const
void equ(const PetscScalar a, const VectorBase &V)
void do_set_add_operation(const size_type n_elements, const size_type *indices, const PetscScalar *values, const bool add_values)
real_type l2_norm() const
VectorBase & operator=(const VectorBase &)=delete
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcGhostsPresent()
#define Assert(cond, exc)
#define AssertIsFinite(number)
#define AssertThrowMPI(error_code)
#define AssertNothrow(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
bool is_non_negative(const T &t)
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)