18 #ifdef DEAL_II_WITH_PETSC
35 VectorReference::operator PetscScalar()
const
49 VecGetOwnershipRange(vector.vector, &
begin, &
end);
52 Vec locally_stored_elements =
nullptr;
53 ierr = VecGhostGetLocalForm(vector.vector, &locally_stored_elements);
57 ierr = VecGetSize(locally_stored_elements, &lsize);
60 const PetscScalar *ptr;
61 ierr = VecGetArrayRead(locally_stored_elements, &ptr);
70 value = *(ptr + index -
begin);
75 Assert(vector.ghost_indices.is_element(index),
77 "You are trying to access an element of a vector "
78 "that is neither a locally owned element nor a "
79 "ghost element of the vector."));
81 vector.ghost_indices.index_within_set(index);
84 value = *(ptr + ghostidx +
end -
begin);
87 ierr = VecRestoreArrayRead(locally_stored_elements, &ptr);
91 VecGhostRestoreLocalForm(vector.vector, &locally_stored_elements);
103 PetscErrorCode ierr = VecGetOwnershipRange(vector.vector, &
begin, &
end);
108 ExcAccessToNonlocalElement(index,
begin,
end - 1));
110 const PetscScalar *ptr;
112 ierr = VecGetArrayRead(vector.vector, &ptr);
114 value = *(ptr + index -
begin);
115 ierr = VecRestoreArrayRead(vector.vector, &ptr);
134 , ghost_indices(v.ghost_indices)
153 const PetscErrorCode ierr =
154 PetscObjectReference(
reinterpret_cast<PetscObject
>(
vector));
163 const PetscErrorCode ierr = VecDestroy(&
vector);
176 PetscErrorCode ierr =
177 PetscObjectReference(
reinterpret_cast<PetscObject
>(v));
179 ierr = VecDestroy(&
vector);
187 const PetscErrorCode ierr = VecDestroy(&
vector);
202 PetscErrorCode ierr = VecCopy(v,
vector);
217 PetscErrorCode ierr = VecSet(
vector, s);
223 ierr = VecGhostGetLocalForm(
vector, &ghost);
226 ierr = VecSet(ghost, s);
229 ierr = VecGhostRestoreLocalForm(
vector, &ghost);
244 const PetscErrorCode ierr = VecEqual(
vector, v.
vector, &flag);
247 return (flag == PETSC_TRUE);
258 const PetscErrorCode ierr = VecEqual(
vector, v.
vector, &flag);
261 return (flag == PETSC_FALSE);
270 const PetscErrorCode ierr = VecGetSize(
vector, &sz);
282 const PetscErrorCode ierr = VecGetLocalSize(
vector, &sz);
294 const PetscErrorCode ierr = VecGetLocalSize(
vector, &sz);
302 std::pair<VectorBase::size_type, VectorBase::size_type>
306 const PetscErrorCode ierr =
307 VecGetOwnershipRange(
static_cast<const Vec &
>(
vector), &
begin, &
end);
317 const std::vector<PetscScalar> &
values)
320 ExcMessage(
"Function called with arguments of different sizes"));
328 const std::vector<PetscScalar> &
values)
331 ExcMessage(
"Function called with arguments of different sizes"));
339 const ::Vector<PetscScalar> &
values)
342 ExcMessage(
"Function called with arguments of different sizes"));
351 const PetscScalar *
values)
372 const PetscErrorCode ierr = VecDot(vec.
vector,
vector, &result);
396 # ifdef DEAL_II_WITH_MPI
400 int all_int_last_action;
402 const int ierr = MPI_Allreduce(&my_int_last_action,
403 &all_int_last_action,
412 ExcMessage(
"Error: not all processors agree on the last "
413 "VectorOperation before this compress() call."));
421 "Missing compress() or calling with wrong VectorOperation argument."));
437 PetscErrorCode ierr = VecAssemblyBegin(
vector);
439 ierr = VecAssemblyEnd(
vector);
467 const PetscErrorCode ierr = VecSum(
vector, &
sum);
469 return sum /
static_cast<PetscReal
>(
size());
474 const PetscScalar *start_ptr;
475 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
478 PetscScalar
mean = 0;
480 PetscScalar sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
485 const PetscScalar *ptr = start_ptr;
486 const PetscScalar *eptr = ptr + (
size() / 4) * 4;
495 while (ptr != start_ptr +
size())
498 mean = (sum0 + sum1 + sum2 + sum3) /
static_cast<PetscReal
>(
size());
503 ierr = VecRestoreArrayRead(
vector, &start_ptr);
515 const PetscErrorCode ierr = VecNorm(
vector, NORM_1, &
d);
528 const PetscErrorCode ierr = VecNorm(
vector, NORM_2, &
d);
541 const PetscScalar *start_ptr;
542 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
547 real_type sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
552 const PetscScalar *ptr = start_ptr;
553 const PetscScalar *eptr = ptr + (
size() / 4) * 4;
562 while (ptr != start_ptr +
size())
565 norm = std::pow(sum0 + sum1 + sum2 + sum3, 1. / p);
570 ierr = VecRestoreArrayRead(
vector, &start_ptr);
583 const PetscErrorCode ierr = VecNorm(
vector, NORM_INFINITY, &
d);
597 const PetscErrorCode ierr = VecMin(
vector, &p, &
d);
610 const PetscErrorCode ierr = VecMax(
vector, &p, &
d);
623 const PetscScalar *start_ptr;
624 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
627 const PetscScalar *ptr = start_ptr,
642 ierr = VecRestoreArrayRead(
vector, &start_ptr);
651 template <
typename T>
660 template <
typename T>
666 "whether it is non-negative."))
return true;
677 const PetscScalar *start_ptr;
678 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
681 const PetscScalar *ptr = start_ptr,
696 ierr = VecRestoreArrayRead(
vector, &start_ptr);
710 const PetscErrorCode ierr = VecScale(
vector, a);
724 const PetscScalar factor = 1. / a;
727 const PetscErrorCode ierr = VecScale(
vector, factor);
739 const PetscErrorCode ierr = VecAXPY(
vector, 1, v);
751 const PetscErrorCode ierr = VecAXPY(
vector, -1, v);
765 const PetscErrorCode ierr = VecShift(
vector, s);
777 const PetscErrorCode ierr = VecAXPY(
vector, a, v);
793 const PetscScalar weights[2] = {a,
b};
794 Vec addends[2] = {v.
vector,
w.vector};
796 const PetscErrorCode ierr = VecMAXPY(
vector, 2, weights, addends);
808 const PetscErrorCode ierr = VecAYPX(
vector, s, v);
836 const PetscErrorCode ierr = VecPointwiseMult(
vector, factors,
vector);
865 MPI_Comm
comm = PetscObjectComm((PetscObject)
vector);
868 PetscErrorCode ierr =
869 PetscViewerSetFormat(PETSC_VIEWER_STDOUT_(
comm), format);
873 ierr = VecView(
vector, PETSC_VIEWER_STDOUT_(
comm));
881 const unsigned int precision,
882 const bool scientific,
883 const bool across)
const
889 const PetscScalar *val;
890 PetscErrorCode ierr = VecGetArrayRead(
vector, &val);
895 const std::ios::fmtflags old_flags = out.flags();
896 const unsigned int old_precision = out.precision(precision);
898 out.precision(precision);
900 out.setf(std::ios::scientific, std::ios::floatfield);
902 out.setf(std::ios::fixed, std::ios::floatfield);
906 out << val[i] <<
' ';
909 out << val[i] << std::endl;
913 out.flags(old_flags);
914 out.precision(old_precision);
918 ierr = VecRestoreArrayRead(
vector, &val);
934 VectorBase::operator
const Vec &()
const
950 std::size_t mem =
sizeof(Vec) +
sizeof(
last_action) +
972 const PetscScalar *
values,
973 const bool add_values)
978 internal::VectorReference::ExcWrongMode(action,
last_action));
988 const PetscInt *petsc_indices =
989 reinterpret_cast<const PetscInt *
>(indices);
991 const InsertMode mode = (add_values ? ADD_VALUES : INSERT_VALUES);
992 const PetscErrorCode ierr =
993 VecSetValues(
vector, n_elements, petsc_indices,
values, mode);
size_type n_elements() const
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
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)
const MPI_Comm & get_mpi_communicator() const
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)
VectorBase & operator=(const VectorBase &)
real_type l2_norm() const
types::global_dof_index size_type
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
#define AssertIsFinite(number)
#define AssertThrowMPI(error_code)
#define AssertNothrow(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcGhostsPresent()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
types::global_dof_index size_type
double norm(const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim >>> &Du)
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
bool is_non_negative(const T &t)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
T sum(const T &t, const MPI_Comm &mpi_communicator)