Reference documentation for deal.II version GIT 85a786e179 2022-11-27 23:15:01+00:00
\(\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\}}\)
petsc_precondition.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2021 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 
20 # include <deal.II/base/utilities.h>
21 
22 # include <deal.II/lac/exceptions.h>
27 
28 # include <petscconf.h>
29 
30 # include <cmath>
31 
33 
34 namespace PETScWrappers
35 {
37  : mpi_communicator(comm)
38  , pc(nullptr)
39  , matrix(nullptr)
40  {}
41 
42 
43 
45  : PreconditionBase(MPI_COMM_NULL)
46  {}
47 
48 
49 
51  {
52  try
53  {
54  clear();
55  }
56  catch (...)
57  {}
58  }
59 
60 
61  void
63  {
64  mpi_communicator = MPI_COMM_NULL;
65 
66  matrix = nullptr;
67 
68  if (pc != nullptr)
69  {
70  PetscErrorCode ierr = PCDestroy(&pc);
71  pc = nullptr;
72  AssertThrow(ierr == 0, ExcPETScError(ierr));
73  }
74  }
75 
76 
77  void
79  {
81 
82  const PetscErrorCode ierr = PCApply(pc, src, dst);
83  AssertThrow(ierr == 0, ExcPETScError(ierr));
84  }
85 
86 
87 
88  void
90  {
92 
93  const PetscErrorCode ierr = PCApplyTranspose(pc, src, dst);
94  AssertThrow(ierr == 0, ExcPETScError(ierr));
95  }
96 
97 
98 
99  MPI_Comm
101  {
102  return mpi_communicator;
103  }
104 
105 
106 
107  void
109  {
110  // only allow the creation of the
111  // preconditioner once
113 
114  MPI_Comm comm = get_mpi_communicator();
115 
116  PetscErrorCode ierr = PCCreate(comm, &pc);
117  AssertThrow(ierr == 0, ExcPETScError(ierr));
118 
119 # if DEAL_II_PETSC_VERSION_LT(3, 5, 0)
120  ierr = PCSetOperators(pc, matrix, matrix, SAME_PRECONDITIONER);
121 # else
122  ierr = PCSetOperators(pc, matrix, matrix);
123 # endif
124  AssertThrow(ierr == 0, ExcPETScError(ierr));
125  }
126 
127 
128  const PC &
130  {
131  return pc;
132  }
133 
134 
135  PreconditionBase::operator Mat() const
136  {
137  return matrix;
138  }
139 
140 
141  /* ----------------- PreconditionJacobi -------------------- */
143  const AdditionalData &additional_data_)
145  {
146  additional_data = additional_data_;
147 
148  PetscErrorCode ierr = PCCreate(comm, &pc);
149  AssertThrow(ierr == 0, ExcPETScError(ierr));
150 
151  initialize();
152  }
153 
154 
155 
157  const AdditionalData &additional_data)
158  : PreconditionBase(matrix.get_mpi_communicator())
159  {
161  }
162 
163 
164 
165  void
167  {
169 
170  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCJACOBI));
171  AssertThrow(ierr == 0, ExcPETScError(ierr));
172 
173  ierr = PCSetFromOptions(pc);
174  AssertThrow(ierr == 0, ExcPETScError(ierr));
175  }
176 
177 
178 
179  void
181  const AdditionalData &additional_data_)
182  {
183  clear();
184 
186 
187  matrix = static_cast<Mat>(matrix_);
188  additional_data = additional_data_;
189 
190  create_pc();
191  initialize();
192 
193  PetscErrorCode ierr = PCSetUp(pc);
194  AssertThrow(ierr == 0, ExcPETScError(ierr));
195  }
196 
197 
198  /* ----------------- PreconditionBlockJacobi -------------------- */
199 
200 
202  const MPI_Comm & comm,
203  const AdditionalData &additional_data_)
205  {
206  additional_data = additional_data_;
207 
208  PetscErrorCode ierr = PCCreate(comm, &pc);
209  AssertThrow(ierr == 0, ExcPETScError(ierr));
210 
211  initialize();
212  }
213 
214 
215 
217  const MatrixBase & matrix,
218  const AdditionalData &additional_data)
219  : PreconditionBase(matrix.get_mpi_communicator())
220  {
222  }
223 
224 
225 
226  void
228  {
229  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCBJACOBI));
230  AssertThrow(ierr == 0, ExcPETScError(ierr));
231 
232  ierr = PCSetFromOptions(pc);
233  AssertThrow(ierr == 0, ExcPETScError(ierr));
234  }
235 
236 
237 
238  void
240  const AdditionalData &additional_data_)
241  {
242  clear();
243 
245 
246  matrix = static_cast<Mat>(matrix_);
247  additional_data = additional_data_;
248 
249  create_pc();
250  initialize();
251 
252  PetscErrorCode ierr = PCSetUp(pc);
253  AssertThrow(ierr == 0, ExcPETScError(ierr));
254  }
255 
256 
257  /* ----------------- PreconditionSOR -------------------- */
258 
260  : omega(omega)
261  {}
262 
263 
264 
268  {
270  }
271 
272 
273  void
275  const AdditionalData &additional_data_)
276  {
277  clear();
278 
280 
281  matrix = static_cast<Mat>(matrix_);
282  additional_data = additional_data_;
283 
284  create_pc();
285 
286  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCSOR));
287  AssertThrow(ierr == 0, ExcPETScError(ierr));
288 
289  // then set flags as given
290  ierr = PCSORSetOmega(pc, additional_data.omega);
291  AssertThrow(ierr == 0, ExcPETScError(ierr));
292 
293  ierr = PCSetFromOptions(pc);
294  AssertThrow(ierr == 0, ExcPETScError(ierr));
295 
296  ierr = PCSetUp(pc);
297  AssertThrow(ierr == 0, ExcPETScError(ierr));
298  }
299 
300 
301  /* ----------------- PreconditionSSOR -------------------- */
302 
304  : omega(omega)
305  {}
306 
307 
308 
312  {
314  }
315 
316 
317  void
319  const AdditionalData &additional_data_)
320  {
321  clear();
322 
324 
325  matrix = static_cast<Mat>(matrix_);
326  additional_data = additional_data_;
327 
328  create_pc();
329 
330  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCSOR));
331  AssertThrow(ierr == 0, ExcPETScError(ierr));
332 
333  // then set flags as given
334  ierr = PCSORSetOmega(pc, additional_data.omega);
335  AssertThrow(ierr == 0, ExcPETScError(ierr));
336 
337  // convert SOR to SSOR
338  ierr = PCSORSetSymmetric(pc, SOR_SYMMETRIC_SWEEP);
339  AssertThrow(ierr == 0, ExcPETScError(ierr));
340 
341  ierr = PCSetFromOptions(pc);
342  AssertThrow(ierr == 0, ExcPETScError(ierr));
343 
344  ierr = PCSetUp(pc);
345  AssertThrow(ierr == 0, ExcPETScError(ierr));
346  }
347 
348 
349  /* ----------------- PreconditionICC -------------------- */
350 
351 
353  : levels(levels)
354  {}
355 
356 
357 
361  {
363  }
364 
365 
366  void
368  const AdditionalData &additional_data_)
369  {
370  clear();
371 
373 
374  matrix = static_cast<Mat>(matrix_);
375  additional_data = additional_data_;
376 
377  create_pc();
378 
379  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCICC));
380  AssertThrow(ierr == 0, ExcPETScError(ierr));
381 
382  // then set flags
383  ierr = PCFactorSetLevels(pc, additional_data.levels);
384  AssertThrow(ierr == 0, ExcPETScError(ierr));
385 
386  ierr = PCSetFromOptions(pc);
387  AssertThrow(ierr == 0, ExcPETScError(ierr));
388 
389  ierr = PCSetUp(pc);
390  AssertThrow(ierr == 0, ExcPETScError(ierr));
391  }
392 
393 
394  /* ----------------- PreconditionILU -------------------- */
395 
397  : levels(levels)
398  {}
399 
400 
401 
405  {
407  }
408 
409 
410  void
412  const AdditionalData &additional_data_)
413  {
414  clear();
415 
417 
418  matrix = static_cast<Mat>(matrix_);
419  additional_data = additional_data_;
420 
421  create_pc();
422 
423  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCILU));
424  AssertThrow(ierr == 0, ExcPETScError(ierr));
425 
426  // then set flags
427  ierr = PCFactorSetLevels(pc, additional_data.levels);
428  AssertThrow(ierr == 0, ExcPETScError(ierr));
429 
430  ierr = PCSetFromOptions(pc);
431  AssertThrow(ierr == 0, ExcPETScError(ierr));
432 
433  ierr = PCSetUp(pc);
434  AssertThrow(ierr == 0, ExcPETScError(ierr));
435  }
436 
437 
438  /* ----------------- PreconditionBoomerAMG -------------------- */
439 
441  const bool symmetric_operator,
442  const double strong_threshold,
443  const double max_row_sum,
444  const unsigned int aggressive_coarsening_num_levels,
445  const bool output_details,
446  const RelaxationType relaxation_type_up,
447  const RelaxationType relaxation_type_down,
448  const RelaxationType relaxation_type_coarse,
449  const unsigned int n_sweeps_coarse,
450  const double tol,
451  const unsigned int max_iter,
452  const bool w_cycle)
453  : symmetric_operator(symmetric_operator)
454  , strong_threshold(strong_threshold)
455  , max_row_sum(max_row_sum)
456  , aggressive_coarsening_num_levels(aggressive_coarsening_num_levels)
457  , output_details(output_details)
458  , relaxation_type_up(relaxation_type_up)
459  , relaxation_type_down(relaxation_type_down)
460  , relaxation_type_coarse(relaxation_type_coarse)
461  , n_sweeps_coarse(n_sweeps_coarse)
462  , tol(tol)
463  , max_iter(max_iter)
464  , w_cycle(w_cycle)
465  {}
466 
467 
468 
469 # ifdef DEAL_II_PETSC_WITH_HYPRE
470  namespace
471  {
476  std::string
477  to_string(
479  {
480  std::string string_type;
481 
482  switch (relaxation_type)
483  {
485  string_type = "Jacobi";
486  break;
489  string_type = "sequential-Gauss-Seidel";
490  break;
493  string_type = "seqboundary-Gauss-Seidel";
494  break;
496  string_type = "SOR/Jacobi";
497  break;
500  string_type = "backward-SOR/Jacobi";
501  break;
504  string_type = "symmetric-SOR/Jacobi";
505  break;
508  string_type = " l1scaled-SOR/Jacobi";
509  break;
512  string_type = "Gaussian-elimination";
513  break;
516  string_type = "l1-Gauss-Seidel";
517  break;
520  string_type = "backward-l1-Gauss-Seidel";
521  break;
523  string_type = "CG";
524  break;
526  string_type = "Chebyshev";
527  break;
529  string_type = "FCF-Jacobi";
530  break;
533  string_type = "l1scaled-Jacobi";
534  break;
536  string_type = "None";
537  break;
538  default:
539  Assert(false, ExcNotImplemented());
540  }
541  return string_type;
542  }
543  } // namespace
544 # endif
545 
547  const MPI_Comm & comm,
548  const AdditionalData &additional_data_)
550  {
551  additional_data = additional_data_;
552 
553  PetscErrorCode ierr = PCCreate(comm, &pc);
554  AssertThrow(ierr == 0, ExcPETScError(ierr));
555 
556 # ifdef DEAL_II_PETSC_WITH_HYPRE
557  initialize();
558 # else // DEAL_II_PETSC_WITH_HYPRE
559  (void)pc;
560  Assert(false,
561  ExcMessage("Your PETSc installation does not include a copy of "
562  "the hypre package necessary for this preconditioner."));
563 # endif
564  }
565 
566 
567 
569  const MatrixBase & matrix,
570  const AdditionalData &additional_data)
571  : PreconditionBase(matrix.get_mpi_communicator())
572  {
574  }
575 
576 
577 
578  void
580  {
581 # ifdef DEAL_II_PETSC_WITH_HYPRE
582  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCHYPRE));
583  AssertThrow(ierr == 0, ExcPETScError(ierr));
584 
585  ierr = PCHYPRESetType(pc, "boomeramg");
586  AssertThrow(ierr == 0, ExcPETScError(ierr));
587 
589  {
590  set_option_value("-pc_hypre_boomeramg_print_statistics", "1");
591  }
592 
593  set_option_value("-pc_hypre_boomeramg_agg_nl",
596 
597  set_option_value("-pc_hypre_boomeramg_max_row_sum",
599 
600  set_option_value("-pc_hypre_boomeramg_strong_threshold",
602 
603  // change to symmetric SOR/Jacobi when using a symmetric operator for
604  // backward compatibility
608  {
611  }
612 
616  {
619  }
620 
624  {
627  }
628 
629  auto relaxation_type_is_symmetric =
630  [](AdditionalData::RelaxationType relaxation_type) {
631  return relaxation_type == AdditionalData::RelaxationType::Jacobi ||
632  relaxation_type ==
634  relaxation_type ==
636  relaxation_type == AdditionalData::RelaxationType::None ||
637  relaxation_type ==
639  relaxation_type == AdditionalData::RelaxationType::CG ||
641  };
642 
644  !relaxation_type_is_symmetric(additional_data.relaxation_type_up))
645  Assert(false,
646  ExcMessage("Use a symmetric smoother for relaxation_type_up"));
647 
649  !relaxation_type_is_symmetric(additional_data.relaxation_type_down))
650  Assert(false,
651  ExcMessage("Use a symmetric smoother for relaxation_type_down"));
652 
654  !relaxation_type_is_symmetric(additional_data.relaxation_type_coarse))
655  Assert(false,
656  ExcMessage("Use a symmetric smoother for relaxation_type_coarse"));
657 
658  set_option_value("-pc_hypre_boomeramg_relax_type_up",
660  set_option_value("-pc_hypre_boomeramg_relax_type_down",
662  set_option_value("-pc_hypre_boomeramg_relax_type_coarse",
664  set_option_value("-pc_hypre_boomeramg_grid_sweeps_coarse",
666 
667  set_option_value("-pc_hypre_boomeramg_tol",
669  set_option_value("-pc_hypre_boomeramg_max_iter",
671 
673  {
674  set_option_value("-pc_hypre_boomeramg_cycle_type", "W");
675  }
676 
677  ierr = PCSetFromOptions(pc);
678  AssertThrow(ierr == 0, ExcPETScError(ierr));
679 # else
680  Assert(false,
681  ExcMessage("Your PETSc installation does not include a copy of "
682  "the hypre package necessary for this preconditioner."));
683 # endif
684  }
685 
686 
687 
688  void
690  const AdditionalData &additional_data_)
691  {
692 # ifdef DEAL_II_PETSC_WITH_HYPRE
693  clear();
694 
696 
697  matrix = static_cast<Mat>(matrix_);
698  additional_data = additional_data_;
699 
700  create_pc();
701  initialize();
702 
703  PetscErrorCode ierr = PCSetUp(pc);
704  AssertThrow(ierr == 0, ExcPETScError(ierr));
705 
706 # else // DEAL_II_PETSC_WITH_HYPRE
707  (void)matrix_;
708  (void)additional_data_;
709  Assert(false,
710  ExcMessage("Your PETSc installation does not include a copy of "
711  "the hypre package necessary for this preconditioner."));
712 # endif
713  }
714 
715 
716  /* ----------------- PreconditionParaSails -------------------- */
717 
719  const unsigned int symmetric,
720  const unsigned int n_levels,
721  const double threshold,
722  const double filter,
723  const bool output_details)
725  , n_levels(n_levels)
726  , threshold(threshold)
727  , filter(filter)
728  , output_details(output_details)
729  {}
730 
731 
732 
734  const MatrixBase & matrix,
737  {
739  }
740 
741 
742  void
744  const AdditionalData &additional_data_)
745  {
746  clear();
747 
749 
750  matrix = static_cast<Mat>(matrix_);
751  additional_data = additional_data_;
752 
753 # ifdef DEAL_II_PETSC_WITH_HYPRE
754  create_pc();
755 
756  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCHYPRE));
757  AssertThrow(ierr == 0, ExcPETScError(ierr));
758 
759  ierr = PCHYPRESetType(pc, "parasails");
760  AssertThrow(ierr == 0, ExcPETScError(ierr));
761 
763  {
764  set_option_value("-pc_hypre_parasails_logging", "1");
765  }
766 
769  ExcMessage(
770  "ParaSails parameter symmetric can only be equal to 0, 1, 2!"));
771 
772  std::stringstream ssStream;
773 
774  switch (additional_data.symmetric)
775  {
776  case 0:
777  {
778  ssStream << "nonsymmetric";
779  break;
780  }
781 
782  case 1:
783  {
784  ssStream << "SPD";
785  break;
786  }
787 
788  case 2:
789  {
790  ssStream << "nonsymmetric,SPD";
791  break;
792  }
793 
794  default:
795  Assert(
796  false,
797  ExcMessage(
798  "ParaSails parameter symmetric can only be equal to 0, 1, 2!"));
799  }
800 
801  set_option_value("-pc_hypre_parasails_sym", ssStream.str());
802 
803  set_option_value("-pc_hypre_parasails_nlevels",
805 
806  ssStream.str(""); // empty the stringstream
807  ssStream << additional_data.threshold;
808  set_option_value("-pc_hypre_parasails_thresh", ssStream.str());
809 
810  ssStream.str(""); // empty the stringstream
811  ssStream << additional_data.filter;
812  set_option_value("-pc_hypre_parasails_filter", ssStream.str());
813 
814  ierr = PCSetFromOptions(pc);
815  AssertThrow(ierr == 0, ExcPETScError(ierr));
816 
817  ierr = PCSetUp(pc);
818  AssertThrow(ierr == 0, ExcPETScError(ierr));
819 
820 # else // DEAL_II_PETSC_WITH_HYPRE
821  (void)pc;
822  Assert(false,
823  ExcMessage("Your PETSc installation does not include a copy of "
824  "the hypre package necessary for this preconditioner."));
825 # endif
826  }
827 
828 
829  /* ----------------- PreconditionNone ------------------------- */
830 
832  const AdditionalData &additional_data)
833  : PreconditionBase(matrix.get_mpi_communicator())
834  {
836  }
837 
838 
839  void
841  const AdditionalData &additional_data_)
842  {
843  clear();
844 
846 
847  matrix = static_cast<Mat>(matrix_);
848  additional_data = additional_data_;
849 
850  create_pc();
851 
852  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCNONE));
853  AssertThrow(ierr == 0, ExcPETScError(ierr));
854 
855  ierr = PCSetFromOptions(pc);
856  AssertThrow(ierr == 0, ExcPETScError(ierr));
857 
858  ierr = PCSetUp(pc);
859  AssertThrow(ierr == 0, ExcPETScError(ierr));
860  }
861 
862 
863  /* ----------------- PreconditionLU -------------------- */
864 
866  const double zero_pivot,
867  const double damping)
868  : pivoting(pivoting)
869  , zero_pivot(zero_pivot)
870  , damping(damping)
871  {}
872 
873 
874 
878  {
880  }
881 
882 
883  void
885  const AdditionalData &additional_data_)
886  {
887  clear();
888 
890 
891  matrix = static_cast<Mat>(matrix_);
892  additional_data = additional_data_;
893 
894  create_pc();
895 
896  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCLU));
897  AssertThrow(ierr == 0, ExcPETScError(ierr));
898 
899  // set flags as given
900  ierr = PCFactorSetColumnPivot(pc, additional_data.pivoting);
901  AssertThrow(ierr == 0, ExcPETScError(ierr));
902 
903  ierr = PCFactorSetZeroPivot(pc, additional_data.zero_pivot);
904  AssertThrow(ierr == 0, ExcPETScError(ierr));
905 
906  ierr = PCFactorSetShiftAmount(pc, additional_data.damping);
907  AssertThrow(ierr == 0, ExcPETScError(ierr));
908 
909  ierr = PCSetFromOptions(pc);
910  AssertThrow(ierr == 0, ExcPETScError(ierr));
911 
912  ierr = PCSetUp(pc);
913  AssertThrow(ierr == 0, ExcPETScError(ierr));
914  }
915 
916  /* ----------------- PreconditionBDDC -------------------- */
917 
918  template <int dim>
920  const bool use_vertices,
921  const bool use_edges,
922  const bool use_faces,
923  const bool symmetric,
924  const std::vector<Point<dim>> coords)
925  : use_vertices(use_vertices)
926  , use_edges(use_edges)
927  , use_faces(use_faces)
929  , coords(coords)
930  {}
931 
932 
933 
934  template <int dim>
936  const MPI_Comm comm,
937  const AdditionalData &additional_data_)
939  {
940  additional_data = additional_data_;
941 
942  PetscErrorCode ierr = PCCreate(comm, &pc);
943  AssertThrow(ierr == 0, ExcPETScError(ierr));
944 
945  initialize();
946  }
947 
948 
949 
950  template <int dim>
952  const AdditionalData &additional_data)
953  : PreconditionBase(matrix.get_mpi_communicator())
954  {
956  }
957 
958 
959 
960  template <int dim>
961  void
963  {
964 # if DEAL_II_PETSC_VERSION_GTE(3, 10, 0)
965  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCBDDC));
966  AssertThrow(ierr == 0, ExcPETScError(ierr));
967 
968  // The matrix must be of IS type. We check for this to avoid the PETSc error
969  // in order to suggest the correct matrix reinit method.
970  {
971  MatType current_type;
972  ierr = MatGetType(matrix, &current_type);
973  AssertThrow(ierr == 0, ExcPETScError(ierr));
974  AssertThrow(
975  strcmp(current_type, MATIS) == 0,
976  ExcMessage(
977  "Matrix must be of IS type. For this, the variant of reinit that includes the active dofs must be used."));
978  }
979 
980 
981  std::stringstream ssStream;
982 
983  if (additional_data.use_vertices)
984  set_option_value("-pc_bddc_use_vertices", "true");
985  else
986  set_option_value("-pc_bddc_use_vertices", "false");
987  if (additional_data.use_edges)
988  set_option_value("-pc_bddc_use_edges", "true");
989  else
990  set_option_value("-pc_bddc_use_edges", "false");
991  if (additional_data.use_faces)
992  set_option_value("-pc_bddc_use_faces", "true");
993  else
994  set_option_value("-pc_bddc_use_faces", "false");
995  if (additional_data.symmetric)
996  set_option_value("-pc_bddc_symmetric", "true");
997  else
998  set_option_value("-pc_bddc_symmetric", "false");
999  if (additional_data.coords.size() > 0)
1000  {
1001  set_option_value("-pc_bddc_corner_selection", "true");
1002  // Convert coords vector to PETSc data array
1003  std::vector<PetscReal> coords_petsc(additional_data.coords.size() *
1004  dim);
1005  for (unsigned int i = 0, j = 0; i < additional_data.coords.size(); ++i)
1006  {
1007  for (j = 0; j < dim; ++j)
1008  coords_petsc[dim * i + j] = additional_data.coords[i][j];
1009  }
1010 
1011  ierr = PCSetCoordinates(pc,
1012  dim,
1013  additional_data.coords.size(),
1014  coords_petsc.data());
1015  AssertThrow(ierr == 0, ExcPETScError(ierr));
1016  }
1017  else
1018  {
1019  set_option_value("-pc_bddc_corner_selection", "false");
1020  ierr = PCSetCoordinates(pc, 0, 0, NULL);
1021  AssertThrow(ierr == 0, ExcPETScError(ierr));
1022  }
1023 
1024 
1025  ierr = PCSetFromOptions(pc);
1026  AssertThrow(ierr == 0, ExcPETScError(ierr));
1027 # else
1028  AssertThrow(
1029  false, ExcMessage("BDDC preconditioner requires PETSc 3.10.0 or newer"));
1030 # endif
1031  }
1032 
1033 
1034 
1035  template <int dim>
1036  void
1038  const AdditionalData &additional_data_)
1039  {
1040  clear();
1041 
1042  mpi_communicator = matrix_.get_mpi_communicator();
1043 
1044  matrix = static_cast<Mat>(matrix_);
1045  additional_data = additional_data_;
1046 
1047  create_pc();
1048  initialize();
1049 
1050  PetscErrorCode ierr = PCSetUp(pc);
1051  AssertThrow(ierr == 0, ExcPETScError(ierr));
1052  }
1053 
1054 
1055 } // namespace PETScWrappers
1056 
1057 template class PETScWrappers::PreconditionBDDC<2>;
1058 template class PETScWrappers::PreconditionBDDC<3>;
1059 
1061 
1062 #endif // DEAL_II_WITH_PETSC
virtual const MPI_Comm & get_mpi_communicator() const =0
void Tvmult(VectorBase &dst, const VectorBase &src) const
void vmult(VectorBase &dst, const VectorBase &src) const
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const MatrixBase &matrix, const AdditionalData &additional_data=AdditionalData())
Definition: point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
#define Assert(cond, exc)
Definition: exceptions.h:1501
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcInvalidState()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1611
@ matrix
Contents is actually a matrix.
@ symmetric
Matrix is symmetric.
void set_option_value(const std::string &name, const std::string &value)
std::string to_string(const T &t)
Definition: patterns.h:2394
AdditionalData(const bool use_vertices=true, const bool use_edges=false, const bool use_faces=false, const bool symmetric=false, const std::vector< Point< dim >> coords={})
AdditionalData(const bool symmetric_operator=false, const double strong_threshold=0.25, const double max_row_sum=0.9, const unsigned int aggressive_coarsening_num_levels=0, const bool output_details=false, const RelaxationType relaxation_type_up=RelaxationType::SORJacobi, const RelaxationType relaxation_type_down=RelaxationType::SORJacobi, const RelaxationType relaxation_type_coarse=RelaxationType::GaussianElimination, const unsigned int n_sweeps_coarse=1, const double tol=0.0, const unsigned int max_iter=1, const bool w_cycle=false)
AdditionalData(const double pivoting=1.e-6, const double zero_pivot=1.e-12, const double damping=0.0)
AdditionalData(const unsigned int symmetric=1, const unsigned int n_levels=1, const double threshold=0.1, const double filter=0.05, const bool output_details=false)
const MPI_Comm & comm