Reference documentation for deal.II version GIT 44af803bb8 2022-12-03 23:15:02+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 -------------------- */
142 
144  : PreconditionBase(MPI_COMM_NULL)
145  {}
146 
147 
148 
150  const AdditionalData &additional_data_)
152  {
153  additional_data = additional_data_;
154 
155  PetscErrorCode ierr = PCCreate(comm, &pc);
156  AssertThrow(ierr == 0, ExcPETScError(ierr));
157 
158  initialize();
159  }
160 
161 
162 
164  const AdditionalData &additional_data)
165  : PreconditionBase(matrix.get_mpi_communicator())
166  {
168  }
169 
170 
171 
172  void
174  {
176 
177  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCJACOBI));
178  AssertThrow(ierr == 0, ExcPETScError(ierr));
179 
180  ierr = PCSetFromOptions(pc);
181  AssertThrow(ierr == 0, ExcPETScError(ierr));
182  }
183 
184 
185 
186  void
188  const AdditionalData &additional_data_)
189  {
190  clear();
191 
193 
194  matrix = static_cast<Mat>(matrix_);
195  additional_data = additional_data_;
196 
197  create_pc();
198  initialize();
199 
200  PetscErrorCode ierr = PCSetUp(pc);
201  AssertThrow(ierr == 0, ExcPETScError(ierr));
202  }
203 
204 
205  /* ----------------- PreconditionBlockJacobi -------------------- */
206 
208  : PreconditionBase(MPI_COMM_NULL)
209  {}
210 
212  const MPI_Comm & comm,
213  const AdditionalData &additional_data_)
215  {
216  additional_data = additional_data_;
217 
218  PetscErrorCode ierr = PCCreate(comm, &pc);
219  AssertThrow(ierr == 0, ExcPETScError(ierr));
220 
221  initialize();
222  }
223 
224 
225 
227  const MatrixBase & matrix,
228  const AdditionalData &additional_data)
229  : PreconditionBase(matrix.get_mpi_communicator())
230  {
232  }
233 
234 
235 
236  void
238  {
239  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCBJACOBI));
240  AssertThrow(ierr == 0, ExcPETScError(ierr));
241 
242  ierr = PCSetFromOptions(pc);
243  AssertThrow(ierr == 0, ExcPETScError(ierr));
244  }
245 
246 
247 
248  void
250  const AdditionalData &additional_data_)
251  {
252  clear();
253 
255 
256  matrix = static_cast<Mat>(matrix_);
257  additional_data = additional_data_;
258 
259  create_pc();
260  initialize();
261 
262  PetscErrorCode ierr = PCSetUp(pc);
263  AssertThrow(ierr == 0, ExcPETScError(ierr));
264  }
265 
266 
267  /* ----------------- PreconditionSOR -------------------- */
268 
270  : PreconditionBase(MPI_COMM_NULL)
271  {}
272 
273 
274 
276  : omega(omega)
277  {}
278 
279 
280 
284  {
286  }
287 
288 
289  void
291  const AdditionalData &additional_data_)
292  {
293  clear();
294 
296 
297  matrix = static_cast<Mat>(matrix_);
298  additional_data = additional_data_;
299 
300  create_pc();
301 
302  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCSOR));
303  AssertThrow(ierr == 0, ExcPETScError(ierr));
304 
305  // then set flags as given
306  ierr = PCSORSetOmega(pc, additional_data.omega);
307  AssertThrow(ierr == 0, ExcPETScError(ierr));
308 
309  ierr = PCSetFromOptions(pc);
310  AssertThrow(ierr == 0, ExcPETScError(ierr));
311 
312  ierr = PCSetUp(pc);
313  AssertThrow(ierr == 0, ExcPETScError(ierr));
314  }
315 
316 
317  /* ----------------- PreconditionSSOR -------------------- */
318 
320  : PreconditionBase(MPI_COMM_NULL)
321  {}
322 
323 
324 
326  : omega(omega)
327  {}
328 
329 
330 
334  {
336  }
337 
338 
339  void
341  const AdditionalData &additional_data_)
342  {
343  clear();
344 
346 
347  matrix = static_cast<Mat>(matrix_);
348  additional_data = additional_data_;
349 
350  create_pc();
351 
352  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCSOR));
353  AssertThrow(ierr == 0, ExcPETScError(ierr));
354 
355  // then set flags as given
356  ierr = PCSORSetOmega(pc, additional_data.omega);
357  AssertThrow(ierr == 0, ExcPETScError(ierr));
358 
359  // convert SOR to SSOR
360  ierr = PCSORSetSymmetric(pc, SOR_SYMMETRIC_SWEEP);
361  AssertThrow(ierr == 0, ExcPETScError(ierr));
362 
363  ierr = PCSetFromOptions(pc);
364  AssertThrow(ierr == 0, ExcPETScError(ierr));
365 
366  ierr = PCSetUp(pc);
367  AssertThrow(ierr == 0, ExcPETScError(ierr));
368  }
369 
370 
371  /* ----------------- PreconditionICC -------------------- */
372 
374  : PreconditionBase(MPI_COMM_NULL)
375  {}
376 
377 
378 
380  : levels(levels)
381  {}
382 
383 
384 
388  {
390  }
391 
392 
393  void
395  const AdditionalData &additional_data_)
396  {
397  clear();
398 
400 
401  matrix = static_cast<Mat>(matrix_);
402  additional_data = additional_data_;
403 
404  create_pc();
405 
406  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCICC));
407  AssertThrow(ierr == 0, ExcPETScError(ierr));
408 
409  // then set flags
410  ierr = PCFactorSetLevels(pc, additional_data.levels);
411  AssertThrow(ierr == 0, ExcPETScError(ierr));
412 
413  ierr = PCSetFromOptions(pc);
414  AssertThrow(ierr == 0, ExcPETScError(ierr));
415 
416  ierr = PCSetUp(pc);
417  AssertThrow(ierr == 0, ExcPETScError(ierr));
418  }
419 
420 
421  /* ----------------- PreconditionILU -------------------- */
422 
424  : PreconditionBase(MPI_COMM_NULL)
425  {}
426 
427 
428 
430  : levels(levels)
431  {}
432 
433 
434 
438  {
440  }
441 
442 
443  void
445  const AdditionalData &additional_data_)
446  {
447  clear();
448 
450 
451  matrix = static_cast<Mat>(matrix_);
452  additional_data = additional_data_;
453 
454  create_pc();
455 
456  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCILU));
457  AssertThrow(ierr == 0, ExcPETScError(ierr));
458 
459  // then set flags
460  ierr = PCFactorSetLevels(pc, additional_data.levels);
461  AssertThrow(ierr == 0, ExcPETScError(ierr));
462 
463  ierr = PCSetFromOptions(pc);
464  AssertThrow(ierr == 0, ExcPETScError(ierr));
465 
466  ierr = PCSetUp(pc);
467  AssertThrow(ierr == 0, ExcPETScError(ierr));
468  }
469 
470 
471  /* ----------------- PreconditionBoomerAMG -------------------- */
472 
474  const bool symmetric_operator,
475  const double strong_threshold,
476  const double max_row_sum,
477  const unsigned int aggressive_coarsening_num_levels,
478  const bool output_details,
479  const RelaxationType relaxation_type_up,
480  const RelaxationType relaxation_type_down,
481  const RelaxationType relaxation_type_coarse,
482  const unsigned int n_sweeps_coarse,
483  const double tol,
484  const unsigned int max_iter,
485  const bool w_cycle)
486  : symmetric_operator(symmetric_operator)
487  , strong_threshold(strong_threshold)
488  , max_row_sum(max_row_sum)
489  , aggressive_coarsening_num_levels(aggressive_coarsening_num_levels)
490  , output_details(output_details)
491  , relaxation_type_up(relaxation_type_up)
492  , relaxation_type_down(relaxation_type_down)
493  , relaxation_type_coarse(relaxation_type_coarse)
494  , n_sweeps_coarse(n_sweeps_coarse)
495  , tol(tol)
496  , max_iter(max_iter)
497  , w_cycle(w_cycle)
498  {}
499 
500 
501 
502 # ifdef DEAL_II_PETSC_WITH_HYPRE
503  namespace
504  {
509  std::string
510  to_string(
512  {
513  std::string string_type;
514 
515  switch (relaxation_type)
516  {
518  string_type = "Jacobi";
519  break;
522  string_type = "sequential-Gauss-Seidel";
523  break;
526  string_type = "seqboundary-Gauss-Seidel";
527  break;
529  string_type = "SOR/Jacobi";
530  break;
533  string_type = "backward-SOR/Jacobi";
534  break;
537  string_type = "symmetric-SOR/Jacobi";
538  break;
541  string_type = " l1scaled-SOR/Jacobi";
542  break;
545  string_type = "Gaussian-elimination";
546  break;
549  string_type = "l1-Gauss-Seidel";
550  break;
553  string_type = "backward-l1-Gauss-Seidel";
554  break;
556  string_type = "CG";
557  break;
559  string_type = "Chebyshev";
560  break;
562  string_type = "FCF-Jacobi";
563  break;
566  string_type = "l1scaled-Jacobi";
567  break;
569  string_type = "None";
570  break;
571  default:
572  Assert(false, ExcNotImplemented());
573  }
574  return string_type;
575  }
576  } // namespace
577 # endif
578 
579 
580 
582  : PreconditionBase(MPI_COMM_NULL)
583  {}
584 
585 
586 
588  const MPI_Comm & comm,
589  const AdditionalData &additional_data_)
591  {
592  additional_data = additional_data_;
593 
594  PetscErrorCode ierr = PCCreate(comm, &pc);
595  AssertThrow(ierr == 0, ExcPETScError(ierr));
596 
597 # ifdef DEAL_II_PETSC_WITH_HYPRE
598  initialize();
599 # else // DEAL_II_PETSC_WITH_HYPRE
600  (void)pc;
601  Assert(false,
602  ExcMessage("Your PETSc installation does not include a copy of "
603  "the hypre package necessary for this preconditioner."));
604 # endif
605  }
606 
607 
608 
610  const MatrixBase & matrix,
611  const AdditionalData &additional_data)
612  : PreconditionBase(matrix.get_mpi_communicator())
613  {
615  }
616 
617 
618 
619  void
621  {
622 # ifdef DEAL_II_PETSC_WITH_HYPRE
623  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCHYPRE));
624  AssertThrow(ierr == 0, ExcPETScError(ierr));
625 
626  ierr = PCHYPRESetType(pc, "boomeramg");
627  AssertThrow(ierr == 0, ExcPETScError(ierr));
628 
630  {
631  set_option_value("-pc_hypre_boomeramg_print_statistics", "1");
632  }
633 
634  set_option_value("-pc_hypre_boomeramg_agg_nl",
637 
638  set_option_value("-pc_hypre_boomeramg_max_row_sum",
640 
641  set_option_value("-pc_hypre_boomeramg_strong_threshold",
643 
644  // change to symmetric SOR/Jacobi when using a symmetric operator for
645  // backward compatibility
649  {
652  }
653 
657  {
660  }
661 
665  {
668  }
669 
670  auto relaxation_type_is_symmetric =
671  [](AdditionalData::RelaxationType relaxation_type) {
672  return relaxation_type == AdditionalData::RelaxationType::Jacobi ||
673  relaxation_type ==
675  relaxation_type ==
677  relaxation_type == AdditionalData::RelaxationType::None ||
678  relaxation_type ==
680  relaxation_type == AdditionalData::RelaxationType::CG ||
682  };
683 
685  !relaxation_type_is_symmetric(additional_data.relaxation_type_up))
686  Assert(false,
687  ExcMessage("Use a symmetric smoother for relaxation_type_up"));
688 
690  !relaxation_type_is_symmetric(additional_data.relaxation_type_down))
691  Assert(false,
692  ExcMessage("Use a symmetric smoother for relaxation_type_down"));
693 
695  !relaxation_type_is_symmetric(additional_data.relaxation_type_coarse))
696  Assert(false,
697  ExcMessage("Use a symmetric smoother for relaxation_type_coarse"));
698 
699  set_option_value("-pc_hypre_boomeramg_relax_type_up",
701  set_option_value("-pc_hypre_boomeramg_relax_type_down",
703  set_option_value("-pc_hypre_boomeramg_relax_type_coarse",
705  set_option_value("-pc_hypre_boomeramg_grid_sweeps_coarse",
707 
708  set_option_value("-pc_hypre_boomeramg_tol",
710  set_option_value("-pc_hypre_boomeramg_max_iter",
712 
714  {
715  set_option_value("-pc_hypre_boomeramg_cycle_type", "W");
716  }
717 
718  ierr = PCSetFromOptions(pc);
719  AssertThrow(ierr == 0, ExcPETScError(ierr));
720 # else
721  Assert(false,
722  ExcMessage("Your PETSc installation does not include a copy of "
723  "the hypre package necessary for this preconditioner."));
724 # endif
725  }
726 
727 
728 
729  void
731  const AdditionalData &additional_data_)
732  {
733 # ifdef DEAL_II_PETSC_WITH_HYPRE
734  clear();
735 
737 
738  matrix = static_cast<Mat>(matrix_);
739  additional_data = additional_data_;
740 
741  create_pc();
742  initialize();
743 
744  PetscErrorCode ierr = PCSetUp(pc);
745  AssertThrow(ierr == 0, ExcPETScError(ierr));
746 
747 # else // DEAL_II_PETSC_WITH_HYPRE
748  (void)matrix_;
749  (void)additional_data_;
750  Assert(false,
751  ExcMessage("Your PETSc installation does not include a copy of "
752  "the hypre package necessary for this preconditioner."));
753 # endif
754  }
755 
756 
757  /* ----------------- PreconditionParaSails -------------------- */
758 
760  const unsigned int symmetric,
761  const unsigned int n_levels,
762  const double threshold,
763  const double filter,
764  const bool output_details)
766  , n_levels(n_levels)
767  , threshold(threshold)
768  , filter(filter)
769  , output_details(output_details)
770  {}
771 
772 
773 
775  : PreconditionBase(MPI_COMM_NULL)
776  {}
777 
778 
779 
781  const MatrixBase & matrix,
782  const AdditionalData &additional_data)
783  : PreconditionBase(matrix.get_mpi_communicator())
784  {
786  }
787 
788 
789  void
791  const AdditionalData &additional_data_)
792  {
793  clear();
794 
796 
797  matrix = static_cast<Mat>(matrix_);
798  additional_data = additional_data_;
799 
800 # ifdef DEAL_II_PETSC_WITH_HYPRE
801  create_pc();
802 
803  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCHYPRE));
804  AssertThrow(ierr == 0, ExcPETScError(ierr));
805 
806  ierr = PCHYPRESetType(pc, "parasails");
807  AssertThrow(ierr == 0, ExcPETScError(ierr));
808 
810  {
811  set_option_value("-pc_hypre_parasails_logging", "1");
812  }
813 
816  ExcMessage(
817  "ParaSails parameter symmetric can only be equal to 0, 1, 2!"));
818 
819  std::stringstream ssStream;
820 
821  switch (additional_data.symmetric)
822  {
823  case 0:
824  {
825  ssStream << "nonsymmetric";
826  break;
827  }
828 
829  case 1:
830  {
831  ssStream << "SPD";
832  break;
833  }
834 
835  case 2:
836  {
837  ssStream << "nonsymmetric,SPD";
838  break;
839  }
840 
841  default:
842  Assert(
843  false,
844  ExcMessage(
845  "ParaSails parameter symmetric can only be equal to 0, 1, 2!"));
846  }
847 
848  set_option_value("-pc_hypre_parasails_sym", ssStream.str());
849 
850  set_option_value("-pc_hypre_parasails_nlevels",
852 
853  ssStream.str(""); // empty the stringstream
854  ssStream << additional_data.threshold;
855  set_option_value("-pc_hypre_parasails_thresh", ssStream.str());
856 
857  ssStream.str(""); // empty the stringstream
858  ssStream << additional_data.filter;
859  set_option_value("-pc_hypre_parasails_filter", ssStream.str());
860 
861  ierr = PCSetFromOptions(pc);
862  AssertThrow(ierr == 0, ExcPETScError(ierr));
863 
864  ierr = PCSetUp(pc);
865  AssertThrow(ierr == 0, ExcPETScError(ierr));
866 
867 # else // DEAL_II_PETSC_WITH_HYPRE
868  (void)pc;
869  Assert(false,
870  ExcMessage("Your PETSc installation does not include a copy of "
871  "the hypre package necessary for this preconditioner."));
872 # endif
873  }
874 
875 
876  /* ----------------- PreconditionNone ------------------------- */
877 
879  : PreconditionBase(MPI_COMM_NULL)
880  {}
881 
882 
883 
885  const AdditionalData &additional_data)
886  : PreconditionBase(matrix.get_mpi_communicator())
887  {
889  }
890 
891 
892  void
894  const AdditionalData &additional_data_)
895  {
896  clear();
897 
899 
900  matrix = static_cast<Mat>(matrix_);
901  additional_data = additional_data_;
902 
903  create_pc();
904 
905  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCNONE));
906  AssertThrow(ierr == 0, ExcPETScError(ierr));
907 
908  ierr = PCSetFromOptions(pc);
909  AssertThrow(ierr == 0, ExcPETScError(ierr));
910 
911  ierr = PCSetUp(pc);
912  AssertThrow(ierr == 0, ExcPETScError(ierr));
913  }
914 
915 
916  /* ----------------- PreconditionLU -------------------- */
917 
919  const double zero_pivot,
920  const double damping)
921  : pivoting(pivoting)
922  , zero_pivot(zero_pivot)
923  , damping(damping)
924  {}
925 
926 
927 
929  : PreconditionBase(MPI_COMM_NULL)
930  {}
931 
932 
933 
935  const AdditionalData &additional_data)
936  : PreconditionBase(matrix.get_mpi_communicator())
937  {
939  }
940 
941 
942  void
944  const AdditionalData &additional_data_)
945  {
946  clear();
947 
949 
950  matrix = static_cast<Mat>(matrix_);
951  additional_data = additional_data_;
952 
953  create_pc();
954 
955  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCLU));
956  AssertThrow(ierr == 0, ExcPETScError(ierr));
957 
958  // set flags as given
959  ierr = PCFactorSetColumnPivot(pc, additional_data.pivoting);
960  AssertThrow(ierr == 0, ExcPETScError(ierr));
961 
962  ierr = PCFactorSetZeroPivot(pc, additional_data.zero_pivot);
963  AssertThrow(ierr == 0, ExcPETScError(ierr));
964 
965  ierr = PCFactorSetShiftAmount(pc, additional_data.damping);
966  AssertThrow(ierr == 0, ExcPETScError(ierr));
967 
968  ierr = PCSetFromOptions(pc);
969  AssertThrow(ierr == 0, ExcPETScError(ierr));
970 
971  ierr = PCSetUp(pc);
972  AssertThrow(ierr == 0, ExcPETScError(ierr));
973  }
974 
975  /* ----------------- PreconditionBDDC -------------------- */
976 
977  template <int dim>
979  const bool use_vertices,
980  const bool use_edges,
981  const bool use_faces,
982  const bool symmetric,
983  const std::vector<Point<dim>> coords)
984  : use_vertices(use_vertices)
985  , use_edges(use_edges)
986  , use_faces(use_faces)
988  , coords(coords)
989  {}
990 
991 
992 
993  template <int dim>
995  : PreconditionBase(MPI_COMM_NULL)
996  {}
997 
998 
999 
1000  template <int dim>
1002  const MPI_Comm comm,
1003  const AdditionalData &additional_data_)
1005  {
1006  additional_data = additional_data_;
1007 
1008  PetscErrorCode ierr = PCCreate(comm, &pc);
1009  AssertThrow(ierr == 0, ExcPETScError(ierr));
1010 
1011  initialize();
1012  }
1013 
1014 
1015 
1016  template <int dim>
1018  const AdditionalData &additional_data)
1019  : PreconditionBase(matrix.get_mpi_communicator())
1020  {
1022  }
1023 
1024 
1025 
1026  template <int dim>
1027  void
1029  {
1030 # if DEAL_II_PETSC_VERSION_GTE(3, 10, 0)
1031  PetscErrorCode ierr = PCSetType(pc, const_cast<char *>(PCBDDC));
1032  AssertThrow(ierr == 0, ExcPETScError(ierr));
1033 
1034  // The matrix must be of IS type. We check for this to avoid the PETSc error
1035  // in order to suggest the correct matrix reinit method.
1036  {
1037  MatType current_type;
1038  ierr = MatGetType(matrix, &current_type);
1039  AssertThrow(ierr == 0, ExcPETScError(ierr));
1040  AssertThrow(
1041  strcmp(current_type, MATIS) == 0,
1042  ExcMessage(
1043  "Matrix must be of IS type. For this, the variant of reinit that includes the active dofs must be used."));
1044  }
1045 
1046 
1047  std::stringstream ssStream;
1048 
1049  if (additional_data.use_vertices)
1050  set_option_value("-pc_bddc_use_vertices", "true");
1051  else
1052  set_option_value("-pc_bddc_use_vertices", "false");
1053  if (additional_data.use_edges)
1054  set_option_value("-pc_bddc_use_edges", "true");
1055  else
1056  set_option_value("-pc_bddc_use_edges", "false");
1057  if (additional_data.use_faces)
1058  set_option_value("-pc_bddc_use_faces", "true");
1059  else
1060  set_option_value("-pc_bddc_use_faces", "false");
1061  if (additional_data.symmetric)
1062  set_option_value("-pc_bddc_symmetric", "true");
1063  else
1064  set_option_value("-pc_bddc_symmetric", "false");
1065  if (additional_data.coords.size() > 0)
1066  {
1067  set_option_value("-pc_bddc_corner_selection", "true");
1068  // Convert coords vector to PETSc data array
1069  std::vector<PetscReal> coords_petsc(additional_data.coords.size() *
1070  dim);
1071  for (unsigned int i = 0, j = 0; i < additional_data.coords.size(); ++i)
1072  {
1073  for (j = 0; j < dim; ++j)
1074  coords_petsc[dim * i + j] = additional_data.coords[i][j];
1075  }
1076 
1077  ierr = PCSetCoordinates(pc,
1078  dim,
1079  additional_data.coords.size(),
1080  coords_petsc.data());
1081  AssertThrow(ierr == 0, ExcPETScError(ierr));
1082  }
1083  else
1084  {
1085  set_option_value("-pc_bddc_corner_selection", "false");
1086  ierr = PCSetCoordinates(pc, 0, 0, NULL);
1087  AssertThrow(ierr == 0, ExcPETScError(ierr));
1088  }
1089 
1090 
1091  ierr = PCSetFromOptions(pc);
1092  AssertThrow(ierr == 0, ExcPETScError(ierr));
1093 # else
1094  AssertThrow(
1095  false, ExcMessage("BDDC preconditioner requires PETSc 3.10.0 or newer"));
1096 # endif
1097  }
1098 
1099 
1100 
1101  template <int dim>
1102  void
1104  const AdditionalData &additional_data_)
1105  {
1106  clear();
1107 
1108  mpi_communicator = matrix_.get_mpi_communicator();
1109 
1110  matrix = static_cast<Mat>(matrix_);
1111  additional_data = additional_data_;
1112 
1113  create_pc();
1114  initialize();
1115 
1116  PetscErrorCode ierr = PCSetUp(pc);
1117  AssertThrow(ierr == 0, ExcPETScError(ierr));
1118  }
1119 
1120 
1121 } // namespace PETScWrappers
1122 
1123 template class PETScWrappers::PreconditionBDDC<2>;
1124 template class PETScWrappers::PreconditionBDDC<3>;
1125 
1127 
1128 #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())
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