18 #ifdef DEAL_II_WITH_PETSC
36 VectorReference::operator PetscScalar()
const
50 VecGetOwnershipRange(vector.vector, &
begin, &
end);
53 Vec locally_stored_elements =
nullptr;
54 ierr = VecGhostGetLocalForm(vector.vector, &locally_stored_elements);
58 ierr = VecGetSize(locally_stored_elements, &lsize);
61 const PetscScalar *ptr;
62 ierr = VecGetArrayRead(locally_stored_elements, &ptr);
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."));
82 vector.ghost_indices.index_within_set(index);
85 value = *(ptr + ghostidx +
end -
begin);
88 ierr = VecRestoreArrayRead(locally_stored_elements, &ptr);
92 VecGhostRestoreLocalForm(vector.vector, &locally_stored_elements);
104 PetscErrorCode ierr = VecGetOwnershipRange(vector.vector, &
begin, &
end);
109 ExcAccessToNonlocalElement(index,
begin,
end - 1));
111 const PetscScalar *ptr;
113 ierr = VecGetArrayRead(vector.vector, &ptr);
115 value = *(ptr + index -
begin);
116 ierr = VecRestoreArrayRead(vector.vector, &ptr);
130 ExcMessage(
"PETSc does not support multi-threaded access, set "
131 "the thread limit to 1 in MPI_InitFinalize()."));
139 , ghost_indices(v.ghost_indices)
143 ExcMessage(
"PETSc does not support multi-threaded access, set "
144 "the thread limit to 1 in MPI_InitFinalize()."));
163 ExcMessage(
"PETSc does not support multi-threaded access, set "
164 "the thread limit to 1 in MPI_InitFinalize()."));
166 const PetscErrorCode ierr =
167 PetscObjectReference(
reinterpret_cast<PetscObject
>(
vector));
176 const PetscErrorCode ierr = VecDestroy(&
vector);
189 PetscErrorCode ierr =
190 PetscObjectReference(
reinterpret_cast<PetscObject
>(v));
192 ierr = VecDestroy(&
vector);
200 const PetscErrorCode ierr = VecDestroy(&
vector);
215 PetscErrorCode ierr = VecCopy(v,
vector);
230 PetscErrorCode ierr = VecSet(
vector, s);
236 ierr = VecGhostGetLocalForm(
vector, &ghost);
239 ierr = VecSet(ghost, s);
242 ierr = VecGhostRestoreLocalForm(
vector, &ghost);
257 const PetscErrorCode ierr = VecEqual(
vector, v.
vector, &flag);
260 return (flag == PETSC_TRUE);
271 const PetscErrorCode ierr = VecEqual(
vector, v.
vector, &flag);
274 return (flag == PETSC_FALSE);
283 const PetscErrorCode ierr = VecGetSize(
vector, &sz);
295 const PetscErrorCode ierr = VecGetLocalSize(
vector, &sz);
307 const PetscErrorCode ierr = VecGetLocalSize(
vector, &sz);
315 std::pair<VectorBase::size_type, VectorBase::size_type>
319 const PetscErrorCode ierr =
320 VecGetOwnershipRange(
static_cast<const Vec &
>(
vector), &
begin, &
end);
330 const std::vector<PetscScalar> &
values)
333 ExcMessage(
"Function called with arguments of different sizes"));
341 const std::vector<PetscScalar> &
values)
344 ExcMessage(
"Function called with arguments of different sizes"));
352 const ::Vector<PetscScalar> &
values)
355 ExcMessage(
"Function called with arguments of different sizes"));
364 const PetscScalar *
values)
385 const PetscErrorCode ierr = VecDot(vec.
vector,
vector, &result);
409 # ifdef DEAL_II_WITH_MPI
413 int all_int_last_action;
415 const int ierr = MPI_Allreduce(&my_int_last_action,
416 &all_int_last_action,
425 ExcMessage(
"Error: not all processors agree on the last "
426 "VectorOperation before this compress() call."));
435 "Missing compress() or calling with wrong VectorOperation argument."));
451 PetscErrorCode ierr = VecAssemblyBegin(
vector);
453 ierr = VecAssemblyEnd(
vector);
481 const PetscErrorCode ierr = VecSum(
vector, &
sum);
483 return sum /
static_cast<PetscReal
>(
size());
488 const PetscScalar *start_ptr;
489 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
492 PetscScalar
mean = 0;
494 PetscScalar sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
499 const PetscScalar *ptr = start_ptr;
500 const PetscScalar *eptr = ptr + (
size() / 4) * 4;
509 while (ptr != start_ptr +
size())
512 mean = (sum0 + sum1 + sum2 + sum3) /
static_cast<PetscReal
>(
size());
517 ierr = VecRestoreArrayRead(
vector, &start_ptr);
529 const PetscErrorCode ierr = VecNorm(
vector, NORM_1, &
d);
542 const PetscErrorCode ierr = VecNorm(
vector, NORM_2, &
d);
555 const PetscScalar *start_ptr;
556 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
561 real_type sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
566 const PetscScalar *ptr = start_ptr;
567 const PetscScalar *eptr = ptr + (
size() / 4) * 4;
576 while (ptr != start_ptr +
size())
579 norm = std::pow(sum0 + sum1 + sum2 + sum3, 1. / p);
584 ierr = VecRestoreArrayRead(
vector, &start_ptr);
597 const PetscErrorCode ierr = VecNorm(
vector, NORM_INFINITY, &
d);
611 const PetscErrorCode ierr = VecMin(
vector, &p, &
d);
624 const PetscErrorCode ierr = VecMax(
vector, &p, &
d);
637 const PetscScalar *start_ptr;
638 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
641 const PetscScalar *ptr = start_ptr,
656 ierr = VecRestoreArrayRead(
vector, &start_ptr);
665 template <
typename T>
674 template <
typename T>
680 "whether it is non-negative."))
return true;
691 const PetscScalar *start_ptr;
692 PetscErrorCode ierr = VecGetArrayRead(
vector, &start_ptr);
695 const PetscScalar *ptr = start_ptr,
710 ierr = VecRestoreArrayRead(
vector, &start_ptr);
724 const PetscErrorCode ierr = VecScale(
vector, a);
738 const PetscScalar factor = 1. / a;
741 const PetscErrorCode ierr = VecScale(
vector, factor);
753 const PetscErrorCode ierr = VecAXPY(
vector, 1, v);
765 const PetscErrorCode ierr = VecAXPY(
vector, -1, v);
779 const PetscErrorCode ierr = VecShift(
vector, s);
791 const PetscErrorCode ierr = VecAXPY(
vector, a, v);
807 const PetscScalar weights[2] = {a,
b};
808 Vec addends[2] = {v.
vector,
w.vector};
810 const PetscErrorCode ierr = VecMAXPY(
vector, 2, weights, addends);
822 const PetscErrorCode ierr = VecAYPX(
vector, s, v);
850 const PetscErrorCode ierr = VecPointwiseMult(
vector, factors,
vector);
879 MPI_Comm
comm = PetscObjectComm((PetscObject)
vector);
882 PetscErrorCode ierr =
883 PetscViewerSetFormat(PETSC_VIEWER_STDOUT_(
comm), format);
887 ierr = VecView(
vector, PETSC_VIEWER_STDOUT_(
comm));
895 const unsigned int precision,
896 const bool scientific,
897 const bool across)
const
903 const PetscScalar *val;
904 PetscErrorCode ierr = VecGetArrayRead(
vector, &val);
909 const std::ios::fmtflags old_flags = out.flags();
910 const unsigned int old_precision = out.precision(precision);
912 out.precision(precision);
914 out.setf(std::ios::scientific, std::ios::floatfield);
916 out.setf(std::ios::fixed, std::ios::floatfield);
920 out << val[i] <<
' ';
923 out << val[i] << std::endl;
927 out.flags(old_flags);
928 out.precision(old_precision);
932 ierr = VecRestoreArrayRead(
vector, &val);
948 VectorBase::operator
const Vec &()
const
964 std::size_t mem =
sizeof(Vec) +
sizeof(
last_action) +
986 const PetscScalar *
values,
987 const bool add_values)
994 internal::VectorReference::ExcWrongMode(action,
last_action));
1002 if (n_elements != 0)
1004 const PetscInt *petsc_indices =
1005 reinterpret_cast<const PetscInt *
>(indices);
1007 const InsertMode mode = (add_values ? ADD_VALUES : INSERT_VALUES);
1008 const PetscErrorCode ierr =
1009 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
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)