Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.4.1
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
petsc_vector_base.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 2022 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
17
18#ifdef DEAL_II_WITH_PETSC
19
22
26
27# include <cmath>
28
30
31namespace PETScWrappers
32{
33 namespace internal
34 {
35# ifndef DOXYGEN
36 VectorReference::operator PetscScalar() const
37 {
38 AssertIndexRange(index, vector.size());
39
40 // The vector may have ghost entries. In that case, we first need to
41 // figure out which elements we own locally, then get a pointer to the
42 // elements that are stored here (both the ones we own as well as the
43 // ghost elements). In this array, the locally owned elements come first
44 // followed by the ghost elements whose position we can get from an
45 // index set.
46 if (vector.ghosted)
47 {
48 PetscInt begin, end;
49 PetscErrorCode ierr =
50 VecGetOwnershipRange(vector.vector, &begin, &end);
51 AssertThrow(ierr == 0, ExcPETScError(ierr));
52
53 Vec locally_stored_elements = PETSC_NULL;
54 ierr = VecGhostGetLocalForm(vector.vector, &locally_stored_elements);
55 AssertThrow(ierr == 0, ExcPETScError(ierr));
56
57 PetscInt lsize;
58 ierr = VecGetSize(locally_stored_elements, &lsize);
59 AssertThrow(ierr == 0, ExcPETScError(ierr));
60
61 PetscScalar *ptr;
62 ierr = VecGetArray(locally_stored_elements, &ptr);
63 AssertThrow(ierr == 0, ExcPETScError(ierr));
64
65 PetscScalar value;
66
67 if (index >= static_cast<size_type>(begin) &&
68 index < static_cast<size_type>(end))
69 {
70 // local entry
71 value = *(ptr + index - begin);
72 }
73 else
74 {
75 // ghost entry
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);
83
84 AssertIndexRange(ghostidx + end - begin, lsize);
85 value = *(ptr + ghostidx + end - begin);
86 }
87
88 ierr = VecRestoreArray(locally_stored_elements, &ptr);
89 AssertThrow(ierr == 0, ExcPETScError(ierr));
90
91 ierr =
92 VecGhostRestoreLocalForm(vector.vector, &locally_stored_elements);
93 AssertThrow(ierr == 0, ExcPETScError(ierr));
94
95 return value;
96 }
97
98
99 // first verify that the requested
100 // element is actually locally
101 // available
102 PetscInt begin, end;
103
104 PetscErrorCode ierr = VecGetOwnershipRange(vector.vector, &begin, &end);
105 AssertThrow(ierr == 0, ExcPETScError(ierr));
106
107 AssertThrow((index >= static_cast<size_type>(begin)) &&
108 (index < static_cast<size_type>(end)),
109 ExcAccessToNonlocalElement(index, begin, end - 1));
110
111 PetscInt idx = index;
112 PetscScalar value;
113 ierr = VecGetValues(vector.vector, 1, &idx, &value);
114 AssertThrow(ierr == 0, ExcPETScError(ierr));
115
116 return value;
117 }
118# endif
119 } // namespace internal
120
122 : vector(nullptr)
123 , ghosted(false)
124 , last_action(::VectorOperation::unknown)
125 , obtained_ownership(true)
126 {
128 ExcMessage("PETSc does not support multi-threaded access, set "
129 "the thread limit to 1 in MPI_InitFinalize()."));
130 }
131
132
133
135 : Subscriptor()
136 , ghosted(v.ghosted)
137 , ghost_indices(v.ghost_indices)
138 , last_action(::VectorOperation::unknown)
139 , obtained_ownership(true)
140 {
142 ExcMessage("PETSc does not support multi-threaded access, set "
143 "the thread limit to 1 in MPI_InitFinalize()."));
144
145 PetscErrorCode ierr = VecDuplicate(v.vector, &vector);
146 AssertThrow(ierr == 0, ExcPETScError(ierr));
147
148 ierr = VecCopy(v.vector, vector);
149 AssertThrow(ierr == 0, ExcPETScError(ierr));
150 }
151
152
153
155 : Subscriptor()
156 , vector(v)
157 , ghosted(false)
158 , last_action(::VectorOperation::unknown)
159 , obtained_ownership(false)
160 {
162 ExcMessage("PETSc does not support multi-threaded access, set "
163 "the thread limit to 1 in MPI_InitFinalize()."));
164 }
165
166
167
169 {
171 {
172 const PetscErrorCode ierr = VecDestroy(&vector);
173 AssertNothrow(ierr == 0, ExcPETScError(ierr));
174 (void)ierr;
175 }
176 }
177
178
179
180 void
182 {
184 {
185 const PetscErrorCode ierr = VecDestroy(&vector);
186 AssertThrow(ierr == 0, ExcPETScError(ierr));
187 }
188
189 ghosted = false;
192 obtained_ownership = true;
193 }
194
195
196
197 VectorBase &
198 VectorBase::operator=(const PetscScalar s)
199 {
201
202 // TODO[TH]: assert(is_compressed())
203
204 PetscErrorCode ierr = VecSet(vector, s);
205 AssertThrow(ierr == 0, ExcPETScError(ierr));
206
207 if (has_ghost_elements())
208 {
209 Vec ghost = PETSC_NULL;
210 ierr = VecGhostGetLocalForm(vector, &ghost);
211 AssertThrow(ierr == 0, ExcPETScError(ierr));
212
213 ierr = VecSet(ghost, s);
214 AssertThrow(ierr == 0, ExcPETScError(ierr));
215
216 ierr = VecGhostRestoreLocalForm(vector, &ghost);
217 AssertThrow(ierr == 0, ExcPETScError(ierr));
218 }
219
220 return *this;
221 }
222
223
224
225 bool
227 {
228 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
229
230 PetscBool flag;
231 const PetscErrorCode ierr = VecEqual(vector, v.vector, &flag);
232 AssertThrow(ierr == 0, ExcPETScError(ierr));
233
234 return (flag == PETSC_TRUE);
235 }
236
237
238
239 bool
241 {
242 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
243
244 PetscBool flag;
245 const PetscErrorCode ierr = VecEqual(vector, v.vector, &flag);
246 AssertThrow(ierr == 0, ExcPETScError(ierr));
247
248 return (flag == PETSC_FALSE);
249 }
250
251
252
255 {
256 PetscInt sz;
257 const PetscErrorCode ierr = VecGetSize(vector, &sz);
258 AssertThrow(ierr == 0, ExcPETScError(ierr));
259
260 return sz;
261 }
262
263
264
267 {
268 PetscInt sz;
269 const PetscErrorCode ierr = VecGetLocalSize(vector, &sz);
270 AssertThrow(ierr == 0, ExcPETScError(ierr));
271
272 return sz;
273 }
274
275
276
279 {
280 PetscInt sz;
281 const PetscErrorCode ierr = VecGetLocalSize(vector, &sz);
282 AssertThrow(ierr == 0, ExcPETScError(ierr));
283
284 return sz;
285 }
286
287
288
289 std::pair<VectorBase::size_type, VectorBase::size_type>
291 {
292 PetscInt begin, end;
293 const PetscErrorCode ierr =
294 VecGetOwnershipRange(static_cast<const Vec &>(vector), &begin, &end);
295 AssertThrow(ierr == 0, ExcPETScError(ierr));
296
297 return std::make_pair(begin, end);
298 }
299
300
301
302 void
303 VectorBase::set(const std::vector<size_type> & indices,
304 const std::vector<PetscScalar> &values)
305 {
306 Assert(indices.size() == values.size(),
307 ExcMessage("Function called with arguments of different sizes"));
308 do_set_add_operation(indices.size(), indices.data(), values.data(), false);
309 }
310
311
312
313 void
314 VectorBase::add(const std::vector<size_type> & indices,
315 const std::vector<PetscScalar> &values)
316 {
317 Assert(indices.size() == values.size(),
318 ExcMessage("Function called with arguments of different sizes"));
319 do_set_add_operation(indices.size(), indices.data(), values.data(), true);
320 }
321
322
323
324 void
325 VectorBase::add(const std::vector<size_type> & indices,
326 const ::Vector<PetscScalar> &values)
327 {
328 Assert(indices.size() == values.size(),
329 ExcMessage("Function called with arguments of different sizes"));
330 do_set_add_operation(indices.size(), indices.data(), values.begin(), true);
331 }
332
333
334
335 void
336 VectorBase::add(const size_type n_elements,
337 const size_type * indices,
338 const PetscScalar *values)
339 {
340 do_set_add_operation(n_elements, indices, values, true);
341 }
342
343
344
345 PetscScalar
347 {
348 Assert(size() == vec.size(), ExcDimensionMismatch(size(), vec.size()));
349
350 PetscScalar result;
351
352 // For complex vectors, VecDot() computes
353 // val = (x,y) = y^H x,
354 // where y^H denotes the conjugate transpose of y.
355 // Note that this corresponds to the usual "mathematicians'"
356 // complex inner product where the SECOND argument gets the
357 // complex conjugate, which is also how we document this
358 // operation.
359 const PetscErrorCode ierr = VecDot(vec.vector, vector, &result);
360 AssertThrow(ierr == 0, ExcPETScError(ierr));
361
362 return result;
363 }
364
365
366
367 PetscScalar
368 VectorBase::add_and_dot(const PetscScalar a,
369 const VectorBase &V,
370 const VectorBase &W)
371 {
372 this->add(a, V);
373 return *this * W;
374 }
375
376
377
378 void
380 {
381 {
382# ifdef DEBUG
383# ifdef DEAL_II_WITH_MPI
384 // Check that all processors agree that last_action is the same (or none!)
385
386 int my_int_last_action = last_action;
387 int all_int_last_action;
388
389 const int ierr = MPI_Allreduce(&my_int_last_action,
390 &all_int_last_action,
391 1,
392 MPI_INT,
393 MPI_BOR,
395 AssertThrowMPI(ierr);
396
397 AssertThrow(all_int_last_action != (::VectorOperation::add |
399 ExcMessage("Error: not all processors agree on the last "
400 "VectorOperation before this compress() call."));
401# endif
402# endif
403 }
404
407 last_action == operation,
409 "Missing compress() or calling with wrong VectorOperation argument."));
410
411 // note that one may think that
412 // we only need to do something
413 // if in fact the state is
414 // anything but
415 // last_action::unknown. but
416 // that's not true: one
417 // frequently gets into
418 // situations where only one
419 // processor (or a subset of
420 // processors) actually writes
421 // something into a vector, but
422 // we still need to call
423 // VecAssemblyBegin/End on all
424 // processors.
425 PetscErrorCode ierr = VecAssemblyBegin(vector);
426 AssertThrow(ierr == 0, ExcPETScError(ierr));
427 ierr = VecAssemblyEnd(vector);
428 AssertThrow(ierr == 0, ExcPETScError(ierr));
429
430 // reset the last action field to
431 // indicate that we're back to a
432 // pristine state
434 }
435
436
437
440 {
441 const real_type d = l2_norm();
442 return d * d;
443 }
444
445
446
447 PetscScalar
449 {
450 // We can only use our more efficient
451 // routine in the serial case.
452 if (dynamic_cast<const PETScWrappers::MPI::Vector *>(this) != nullptr)
453 {
454 PetscScalar sum;
455 const PetscErrorCode ierr = VecSum(vector, &sum);
456 AssertThrow(ierr == 0, ExcPETScError(ierr));
457 return sum / static_cast<PetscReal>(size());
458 }
459
460 // get a representation of the vector and
461 // loop over all the elements
462 PetscScalar * start_ptr;
463 PetscErrorCode ierr = VecGetArray(vector, &start_ptr);
464 AssertThrow(ierr == 0, ExcPETScError(ierr));
465
466 PetscScalar mean = 0;
467 {
468 PetscScalar sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
469
470 // use modern processors better by
471 // allowing pipelined commands to be
472 // executed in parallel
473 const PetscScalar *ptr = start_ptr;
474 const PetscScalar *eptr = ptr + (size() / 4) * 4;
475 while (ptr != eptr)
476 {
477 sum0 += *ptr++;
478 sum1 += *ptr++;
479 sum2 += *ptr++;
480 sum3 += *ptr++;
481 }
482 // add up remaining elements
483 while (ptr != start_ptr + size())
484 sum0 += *ptr++;
485
486 mean = (sum0 + sum1 + sum2 + sum3) / static_cast<PetscReal>(size());
487 }
488
489 // restore the representation of the
490 // vector
491 ierr = VecRestoreArray(vector, &start_ptr);
492 AssertThrow(ierr == 0, ExcPETScError(ierr));
493
494 return mean;
495 }
496
497
500 {
501 real_type d;
502
503 const PetscErrorCode ierr = VecNorm(vector, NORM_1, &d);
504 AssertThrow(ierr == 0, ExcPETScError(ierr));
505
506 return d;
507 }
508
509
510
513 {
514 real_type d;
515
516 const PetscErrorCode ierr = VecNorm(vector, NORM_2, &d);
517 AssertThrow(ierr == 0, ExcPETScError(ierr));
518
519 return d;
520 }
521
522
523
526 {
527 // get a representation of the vector and
528 // loop over all the elements
529 PetscScalar * start_ptr;
530 PetscErrorCode ierr = VecGetArray(vector, &start_ptr);
531 AssertThrow(ierr == 0, ExcPETScError(ierr));
532
533 real_type norm = 0;
534 {
535 real_type sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
536
537 // use modern processors better by
538 // allowing pipelined commands to be
539 // executed in parallel
540 const PetscScalar *ptr = start_ptr;
541 const PetscScalar *eptr = ptr + (size() / 4) * 4;
542 while (ptr != eptr)
543 {
548 }
549 // add up remaining elements
550 while (ptr != start_ptr + size())
552
553 norm = std::pow(sum0 + sum1 + sum2 + sum3, 1. / p);
554 }
555
556 // restore the representation of the
557 // vector
558 ierr = VecRestoreArray(vector, &start_ptr);
559 AssertThrow(ierr == 0, ExcPETScError(ierr));
560
561 return norm;
562 }
563
564
565
568 {
569 real_type d;
570
571 const PetscErrorCode ierr = VecNorm(vector, NORM_INFINITY, &d);
572 AssertThrow(ierr == 0, ExcPETScError(ierr));
573
574 return d;
575 }
576
577
578
581 {
582 PetscInt p;
583 real_type d;
584
585 const PetscErrorCode ierr = VecMin(vector, &p, &d);
586 AssertThrow(ierr == 0, ExcPETScError(ierr));
587
588 return d;
589 }
590
591
594 {
595 PetscInt p;
596 real_type d;
597
598 const PetscErrorCode ierr = VecMax(vector, &p, &d);
599 AssertThrow(ierr == 0, ExcPETScError(ierr));
600
601 return d;
602 }
603
604
605
606 bool
608 {
609 // get a representation of the vector and
610 // loop over all the elements
611 PetscScalar * start_ptr;
612 PetscErrorCode ierr = VecGetArray(vector, &start_ptr);
613 AssertThrow(ierr == 0, ExcPETScError(ierr));
614
615 const PetscScalar *ptr = start_ptr,
616 *eptr = start_ptr + locally_owned_size();
617 bool flag = true;
618 while (ptr != eptr)
619 {
620 if (*ptr != value_type())
621 {
622 flag = false;
623 break;
624 }
625 ++ptr;
626 }
627
628 // restore the representation of the
629 // vector
630 ierr = VecRestoreArray(vector, &start_ptr);
631 AssertThrow(ierr == 0, ExcPETScError(ierr));
632
633 return flag;
634 }
635
636
637 namespace internal
638 {
639 template <typename T>
640 bool
641 is_non_negative(const T &t)
642 {
643 return t >= 0;
644 }
645
646
647
648 template <typename T>
649 bool
650 is_non_negative(const std::complex<T> &)
651 {
652 Assert(false,
653 ExcMessage("You can't ask a complex value "
654 "whether it is non-negative.")) return true;
655 }
656 } // namespace internal
657
658
659
660 bool
662 {
663 // get a representation of the vector and
664 // loop over all the elements
665 PetscScalar * start_ptr;
666 PetscErrorCode ierr = VecGetArray(vector, &start_ptr);
667 AssertThrow(ierr == 0, ExcPETScError(ierr));
668
669 const PetscScalar *ptr = start_ptr,
670 *eptr = start_ptr + locally_owned_size();
671 bool flag = true;
672 while (ptr != eptr)
673 {
674 if (!internal::is_non_negative(*ptr))
675 {
676 flag = false;
677 break;
678 }
679 ++ptr;
680 }
681
682 // restore the representation of the
683 // vector
684 ierr = VecRestoreArray(vector, &start_ptr);
685 AssertThrow(ierr == 0, ExcPETScError(ierr));
686
687 return flag;
688 }
689
690
691
692 VectorBase &
693 VectorBase::operator*=(const PetscScalar a)
694 {
697
698 const PetscErrorCode ierr = VecScale(vector, a);
699 AssertThrow(ierr == 0, ExcPETScError(ierr));
700
701 return *this;
702 }
703
704
705
706 VectorBase &
707 VectorBase::operator/=(const PetscScalar a)
708 {
711
712 const PetscScalar factor = 1. / a;
713 AssertIsFinite(factor);
714
715 const PetscErrorCode ierr = VecScale(vector, factor);
716 AssertThrow(ierr == 0, ExcPETScError(ierr));
717
718 return *this;
719 }
720
721
722
723 VectorBase &
725 {
727 const PetscErrorCode ierr = VecAXPY(vector, 1, v);
728 AssertThrow(ierr == 0, ExcPETScError(ierr));
729
730 return *this;
731 }
732
733
734
735 VectorBase &
737 {
739 const PetscErrorCode ierr = VecAXPY(vector, -1, v);
740 AssertThrow(ierr == 0, ExcPETScError(ierr));
741
742 return *this;
743 }
744
745
746
747 void
748 VectorBase::add(const PetscScalar s)
749 {
752
753 const PetscErrorCode ierr = VecShift(vector, s);
754 AssertThrow(ierr == 0, ExcPETScError(ierr));
755 }
756
757
758
759 void
760 VectorBase::add(const PetscScalar a, const VectorBase &v)
761 {
764
765 const PetscErrorCode ierr = VecAXPY(vector, a, v);
766 AssertThrow(ierr == 0, ExcPETScError(ierr));
767 }
768
769
770
771 void
772 VectorBase::add(const PetscScalar a,
773 const VectorBase &v,
774 const PetscScalar b,
775 const VectorBase &w)
776 {
780
781 const PetscScalar weights[2] = {a, b};
782 Vec addends[2] = {v.vector, w.vector};
783
784 const PetscErrorCode ierr = VecMAXPY(vector, 2, weights, addends);
785 AssertThrow(ierr == 0, ExcPETScError(ierr));
786 }
787
788
789
790 void
791 VectorBase::sadd(const PetscScalar s, const VectorBase &v)
792 {
795
796 const PetscErrorCode ierr = VecAYPX(vector, s, v);
797 AssertThrow(ierr == 0, ExcPETScError(ierr));
798 }
799
800
801
802 void
803 VectorBase::sadd(const PetscScalar s,
804 const PetscScalar a,
805 const VectorBase &v)
806 {
810
811 // there is nothing like a AXPAY
812 // operation in Petsc, so do it in two
813 // steps
814 *this *= s;
815 add(a, v);
816 }
817
818
819
820 void
822 {
824 const PetscErrorCode ierr = VecPointwiseMult(vector, factors, vector);
825 AssertThrow(ierr == 0, ExcPETScError(ierr));
826 }
827
828
829
830 void
831 VectorBase::equ(const PetscScalar a, const VectorBase &v)
832 {
835
836 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
837
838 // there is no simple operation for this
839 // in PETSc. there are multiple ways to
840 // emulate it, we choose this one:
841 const PetscErrorCode ierr = VecCopy(v.vector, vector);
842 AssertThrow(ierr == 0, ExcPETScError(ierr));
843
844 *this *= a;
845 }
846
847
848
849 void
850 VectorBase::write_ascii(const PetscViewerFormat format)
851 {
852 // TODO[TH]:assert(is_compressed())
853
854 // Set options
855 PetscErrorCode ierr =
856 PetscViewerSetFormat(PETSC_VIEWER_STDOUT_WORLD, format);
857 AssertThrow(ierr == 0, ExcPETScError(ierr));
858
859 // Write to screen
860 ierr = VecView(vector, PETSC_VIEWER_STDOUT_WORLD);
861 AssertThrow(ierr == 0, ExcPETScError(ierr));
862 }
863
864
865
866 void
867 VectorBase::print(std::ostream & out,
868 const unsigned int precision,
869 const bool scientific,
870 const bool across) const
871 {
872 AssertThrow(out.fail() == false, ExcIO());
873
874 // get a representation of the vector and
875 // loop over all the elements
876 PetscScalar * val;
877 PetscErrorCode ierr = VecGetArray(vector, &val);
878
879 AssertThrow(ierr == 0, ExcPETScError(ierr));
880
881 // save the state of out stream
882 const std::ios::fmtflags old_flags = out.flags();
883 const unsigned int old_precision = out.precision(precision);
884
885 out.precision(precision);
886 if (scientific)
887 out.setf(std::ios::scientific, std::ios::floatfield);
888 else
889 out.setf(std::ios::fixed, std::ios::floatfield);
890
891 if (across)
892 for (size_type i = 0; i < locally_owned_size(); ++i)
893 out << val[i] << ' ';
894 else
895 for (size_type i = 0; i < locally_owned_size(); ++i)
896 out << val[i] << std::endl;
897 out << std::endl;
898
899 // reset output format
900 out.flags(old_flags);
901 out.precision(old_precision);
902
903 // restore the representation of the
904 // vector
905 ierr = VecRestoreArray(vector, &val);
906 AssertThrow(ierr == 0, ExcPETScError(ierr));
907
908 AssertThrow(out.fail() == false, ExcIO());
909 }
910
911
912
913 void
915 {
916 const PetscErrorCode ierr = VecSwap(vector, v.vector);
917 AssertThrow(ierr == 0, ExcPETScError(ierr));
918 }
919
920
921
922 VectorBase::operator const Vec &() const
923 {
924 return vector;
925 }
926
927
928 std::size_t
930 {
931 std::size_t mem = sizeof(Vec) + sizeof(last_action) +
934
935 // TH: I am relatively sure that PETSc is
936 // storing the local data in a contiguous
937 // block without indices:
938 mem += locally_owned_size() * sizeof(PetscScalar);
939 // assume that PETSc is storing one index
940 // and one double per ghost element
941 if (ghosted)
942 mem += ghost_indices.n_elements() * (sizeof(PetscScalar) + sizeof(int));
943
944 // TODO[TH]: size of constant memory for PETSc?
945 return mem;
946 }
947
948
949
950 void
952 const size_type * indices,
953 const PetscScalar *values,
954 const bool add_values)
955 {
957 (add_values ? ::VectorOperation::add :
959 Assert((last_action == action) ||
961 internal::VectorReference::ExcWrongMode(action, last_action));
963 // VecSetValues complains if we
964 // come with an empty
965 // vector. however, it is not a
966 // collective operation, so we
967 // can skip the call if necessary
968 // (unlike the above calls)
969 if (n_elements != 0)
970 {
971 const PetscInt *petsc_indices =
972 reinterpret_cast<const PetscInt *>(indices);
973
974 const InsertMode mode = (add_values ? ADD_VALUES : INSERT_VALUES);
975 const PetscErrorCode ierr =
976 VecSetValues(vector, n_elements, petsc_indices, values, mode);
977 AssertThrow(ierr == 0, ExcPETScError(ierr));
978 }
979
980 // set the mode here, independent of whether we have actually
981 // written elements or whether the list was empty
982 last_action = action;
983 }
984
985} // namespace PETScWrappers
986
988
989#endif // DEAL_II_WITH_PETSC
size_type n_elements() const
Definition: index_set.h:1834
void clear()
Definition: index_set.h:1612
static bool is_running_single_threaded()
real_type lp_norm(const real_type p) const
VectorBase & operator+=(const VectorBase &V)
PetscScalar mean_value() 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)
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)
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
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 &)=delete
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcGhostsPresent()
#define Assert(cond, exc)
Definition: exceptions.h:1473
#define AssertIsFinite(number)
Definition: exceptions.h:1758
#define AssertThrowMPI(error_code)
Definition: exceptions.h:1790
#define AssertNothrow(cond, exc)
Definition: exceptions.h:1536
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
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)