Reference documentation for deal.II version GIT 69313620a0 2022-11-28 13:00: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_solver.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 
16 
17 #include <deal.II/base/logstream.h>
18 
20 
21 #ifdef DEAL_II_WITH_PETSC
22 
23 # include <deal.II/lac/exceptions.h>
28 
29 # include <petscversion.h>
30 
31 # include <cmath>
32 # include <memory>
33 
35 
36 namespace PETScWrappers
37 {
39  {
41  }
42 
43 
44 
46  : solver_control(cn)
47  {}
48 
49 
50 
51  void
53  VectorBase & x,
54  const VectorBase & b,
55  const PreconditionBase &preconditioner)
56  {
57  /*
58  TODO: PETSc duplicates communicators, so this does not work (you put
59  MPI_COMM_SELF in, but get something other out when you ask PETSc for the
60  communicator. This mainly fails due to the MatrixFree classes, that can not
61  ask PETSc for a communicator. //Timo Heister
62  Assert(A.get_mpi_communicator()==mpi_communicator, ExcMessage("PETSc Solver
63  and Matrix need to use the same MPI_Comm."));
64  Assert(x.get_mpi_communicator()==mpi_communicator, ExcMessage("PETSc Solver
65  and Vector need to use the same MPI_Comm."));
66  Assert(b.get_mpi_communicator()==mpi_communicator, ExcMessage("PETSc Solver
67  and Vector need to use the same MPI_Comm."));
68  */
69 
70  // first create a solver object if this
71  // is necessary
72  if (solver_data.get() == nullptr)
73  {
74  solver_data = std::make_unique<SolverData>();
75 
76  PetscErrorCode ierr =
77  KSPCreate(A.get_mpi_communicator(), &solver_data->ksp);
78  AssertThrow(ierr == 0, ExcPETScError(ierr));
79 
80  // let derived classes set the solver
81  // type, and the preconditioning
82  // object set the type of
83  // preconditioner
85 
86  ierr = KSPSetPC(solver_data->ksp, preconditioner.get_pc());
87  AssertThrow(ierr == 0, ExcPETScError(ierr));
88 
89  // make sure the preconditioner has an associated matrix set
90  const Mat B = preconditioner;
91  AssertThrow(B != nullptr,
92  ExcMessage("PETSc preconditioner should have an "
93  "associated matrix set to be used in solver."));
94 
95  // setting the preconditioner overwrites the used matrices.
96  // hence, we need to set the matrices after the preconditioner.
97  ierr = KSPSetOperators(solver_data->ksp, A, preconditioner);
98  AssertThrow(ierr == 0, ExcPETScError(ierr));
99 
100  // then a convergence monitor
101  // function. that function simply
102  // checks with the solver_control
103  // object we have in this object for
104  // convergence
105  ierr = KSPSetConvergenceTest(solver_data->ksp,
107  reinterpret_cast<void *>(&solver_control),
108  nullptr);
109  AssertThrow(ierr == 0, ExcPETScError(ierr));
110  }
111 
112  // set the command line option prefix name
113  PetscErrorCode ierr =
114  KSPSetOptionsPrefix(solver_data->ksp, prefix_name.c_str());
115  AssertThrow(ierr == 0, ExcPETScError(ierr));
116 
117  // set the command line options provided
118  // by the user to override the defaults
119  ierr = KSPSetFromOptions(solver_data->ksp);
120  AssertThrow(ierr == 0, ExcPETScError(ierr));
121 
122  // then do the real work: set up solver
123  // internal data and solve the
124  // system.
125  ierr = KSPSetUp(solver_data->ksp);
126  AssertThrow(ierr == 0, ExcPETScError(ierr));
127 
128  ierr = KSPSolve(solver_data->ksp, b, x);
129  AssertThrow(ierr == 0, ExcPETScError(ierr));
130 
131  // do not destroy solver object
132  // solver_data.reset ();
133 
134  // in case of failure: throw
135  // exception
137  AssertThrow(false,
140  // otherwise exit as normal
141  }
142 
143 
144  void
145  SolverBase::set_prefix(const std::string &prefix)
146  {
147  prefix_name = prefix;
148  }
149 
150 
151  void
153  {
154  solver_data.reset();
155  }
156 
157 
158  SolverControl &
160  {
161  return solver_control;
162  }
163 
164 
165  int
167  const PetscInt iteration,
168  const PetscReal residual_norm,
169  KSPConvergedReason *reason,
170  void * solver_control_x)
171  {
173  *reinterpret_cast<SolverControl *>(solver_control_x);
174 
175  const SolverControl::State state =
176  solver_control.check(iteration, residual_norm);
177 
178  switch (state)
179  {
180  case ::SolverControl::iterate:
181  *reason = KSP_CONVERGED_ITERATING;
182  break;
183 
184  case ::SolverControl::success:
185  *reason = static_cast<KSPConvergedReason>(1);
186  break;
187 
188  case ::SolverControl::failure:
190  *reason = KSP_DIVERGED_ITS;
191  else
192  *reason = KSP_DIVERGED_DTOL;
193  break;
194 
195  default:
196  Assert(false, ExcNotImplemented());
197  }
198 
199  // return without failure
200  return 0;
201  }
202 
203 
204 
205  void
207  {
208  PetscErrorCode ierr;
209 
210  solver_data = std::make_unique<SolverData>();
211 
212  ierr = KSPCreate(preconditioner.get_mpi_communicator(), &solver_data->ksp);
213  AssertThrow(ierr == 0, ExcPETScError(ierr));
214 
215  // let derived classes set the solver
216  // type, and the preconditioning
217  // object set the type of
218  // preconditioner
220 
221  ierr = KSPSetPC(solver_data->ksp, preconditioner.get_pc());
222  AssertThrow(ierr == 0, ExcPETScError(ierr));
223 
224  // then a convergence monitor
225  // function. that function simply
226  // checks with the solver_control
227  // object we have in this object for
228  // convergence
229  ierr = KSPSetConvergenceTest(solver_data->ksp,
231  reinterpret_cast<void *>(&solver_control),
232  nullptr);
233  AssertThrow(ierr == 0, ExcPETScError(ierr));
234 
235  // set the command line options provided
236  // by the user to override the defaults
237  ierr = KSPSetFromOptions(solver_data->ksp);
238  AssertThrow(ierr == 0, ExcPETScError(ierr));
239  }
240 
241 
242 
243  /* ---------------------- SolverRichardson ------------------------ */
244 
246  : omega(omega)
247  {}
248 
249 
250 
252  const AdditionalData &data)
253  : SolverBase(cn)
254  , additional_data(data)
255  {}
256 
257 
258 
260  const MPI_Comm &,
261  const AdditionalData &data)
262  : SolverRichardson(cn, data)
263  {}
264 
265 
266  void
268  {
269  PetscErrorCode ierr = KSPSetType(ksp, KSPRICHARDSON);
270  AssertThrow(ierr == 0, ExcPETScError(ierr));
271 
272  // set the damping factor from the data
273  ierr = KSPRichardsonSetScale(ksp, additional_data.omega);
274  AssertThrow(ierr == 0, ExcPETScError(ierr));
275 
276  // in the deal.II solvers, we always
277  // honor the initial guess in the
278  // solution vector. do so here as well:
279  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
280  AssertThrow(ierr == 0, ExcPETScError(ierr));
281 
282  // Hand over the absolute
283  // tolerance and the maximum
284  // iteration number to the PETSc
285  // convergence criterion. The
286  // custom deal.II SolverControl
287  // object is ignored by the PETSc
288  // Richardson method (when no
289  // PETSc monitoring is present),
290  // since in this case PETSc
291  // uses a faster version of
292  // the Richardson iteration,
293  // where no residual is
294  // available.
295  ierr = KSPSetTolerances(ksp,
296  PETSC_DEFAULT,
297  this->solver_control.tolerance(),
298  PETSC_DEFAULT,
299  this->solver_control.max_steps() + 1);
300  AssertThrow(ierr == 0, ExcPETScError(ierr));
301  }
302 
303 
304  /* ---------------------- SolverChebychev ------------------------ */
305 
307  const AdditionalData &data)
308  : SolverBase(cn)
309  , additional_data(data)
310  {}
311 
312 
314  const MPI_Comm &,
315  const AdditionalData &data)
316  : SolverChebychev(cn, data)
317  {}
318 
319 
320  void
322  {
323  PetscErrorCode ierr = KSPSetType(ksp, KSPCHEBYSHEV);
324  AssertThrow(ierr == 0, ExcPETScError(ierr));
325 
326  // in the deal.II solvers, we always
327  // honor the initial guess in the
328  // solution vector. do so here as well:
329  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
330  AssertThrow(ierr == 0, ExcPETScError(ierr));
331  }
332 
333 
334  /* ---------------------- SolverCG ------------------------ */
335 
337  : SolverBase(cn)
338  , additional_data(data)
339  {}
340 
341 
343  const MPI_Comm &,
344  const AdditionalData &data)
345  : SolverCG(cn, data)
346  {}
347 
348 
349  void
351  {
352  PetscErrorCode ierr = KSPSetType(ksp, KSPCG);
353  AssertThrow(ierr == 0, ExcPETScError(ierr));
354 
355  // in the deal.II solvers, we always
356  // honor the initial guess in the
357  // solution vector. do so here as well:
358  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
359  AssertThrow(ierr == 0, ExcPETScError(ierr));
360  }
361 
362 
363  /* ---------------------- SolverBiCG ------------------------ */
364 
366  : SolverBase(cn)
367  , additional_data(data)
368  {}
369 
370 
372  const MPI_Comm &,
373  const AdditionalData &data)
374  : SolverBiCG(cn, data)
375  {}
376 
377 
378  void
380  {
381  PetscErrorCode ierr = KSPSetType(ksp, KSPBICG);
382  AssertThrow(ierr == 0, ExcPETScError(ierr));
383 
384  // in the deal.II solvers, we always
385  // honor the initial guess in the
386  // solution vector. do so here as well:
387  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
388  AssertThrow(ierr == 0, ExcPETScError(ierr));
389  }
390 
391 
392  /* ---------------------- SolverGMRES ------------------------ */
393 
395  const unsigned int restart_parameter,
396  const bool right_preconditioning)
397  : restart_parameter(restart_parameter)
398  , right_preconditioning(right_preconditioning)
399  {}
400 
401 
402 
404  : SolverBase(cn)
405  , additional_data(data)
406  {}
407 
408 
410  const MPI_Comm &,
411  const AdditionalData &data)
412  : SolverGMRES(cn, data)
413  {}
414 
415 
416  void
418  {
419  PetscErrorCode ierr = KSPSetType(ksp, KSPGMRES);
420  AssertThrow(ierr == 0, ExcPETScError(ierr));
421 
422  ierr = KSPGMRESSetRestart(ksp, additional_data.restart_parameter);
423  AssertThrow(ierr == 0, ExcPETScError(ierr));
424 
425  // Set preconditioning side to right
427  {
428  ierr = KSPSetPCSide(ksp, PC_RIGHT);
429  AssertThrow(ierr == 0, ExcPETScError(ierr));
430  }
431 
432  // in the deal.II solvers, we always
433  // honor the initial guess in the
434  // solution vector. do so here as well:
435  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
436  AssertThrow(ierr == 0, ExcPETScError(ierr));
437  }
438 
439 
440  /* ---------------------- SolverBicgstab ------------------------ */
441 
443  : SolverBase(cn)
444  , additional_data(data)
445  {}
446 
447 
449  const MPI_Comm &,
450  const AdditionalData &data)
451  : SolverBicgstab(cn, data)
452  {}
453 
454 
455  void
457  {
458  PetscErrorCode ierr = KSPSetType(ksp, KSPBCGS);
459  AssertThrow(ierr == 0, ExcPETScError(ierr));
460 
461  // in the deal.II solvers, we always
462  // honor the initial guess in the
463  // solution vector. do so here as well:
464  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
465  AssertThrow(ierr == 0, ExcPETScError(ierr));
466  }
467 
468 
469  /* ---------------------- SolverCGS ------------------------ */
470 
472  : SolverBase(cn)
473  , additional_data(data)
474  {}
475 
476 
478  const MPI_Comm &,
479  const AdditionalData &data)
480  : SolverCGS(cn, data)
481  {}
482 
483 
484  void
486  {
487  PetscErrorCode ierr = KSPSetType(ksp, KSPCGS);
488  AssertThrow(ierr == 0, ExcPETScError(ierr));
489 
490  // in the deal.II solvers, we always
491  // honor the initial guess in the
492  // solution vector. do so here as well:
493  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
494  AssertThrow(ierr == 0, ExcPETScError(ierr));
495  }
496 
497 
498  /* ---------------------- SolverTFQMR ------------------------ */
499 
501  : SolverBase(cn)
502  , additional_data(data)
503  {}
504 
505 
507  const MPI_Comm &,
508  const AdditionalData &data)
509  : SolverTFQMR(cn, data)
510  {}
511 
512 
513  void
515  {
516  PetscErrorCode ierr = KSPSetType(ksp, KSPTFQMR);
517  AssertThrow(ierr == 0, ExcPETScError(ierr));
518 
519  // in the deal.II solvers, we always
520  // honor the initial guess in the
521  // solution vector. do so here as well:
522  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
523  AssertThrow(ierr == 0, ExcPETScError(ierr));
524  }
525 
526 
527  /* ---------------------- SolverTCQMR ------------------------ */
528 
530  : SolverBase(cn)
531  , additional_data(data)
532  {}
533 
534 
536  const MPI_Comm &,
537  const AdditionalData &data)
538  : SolverTCQMR(cn, data)
539  {}
540 
541 
542  void
544  {
545  PetscErrorCode ierr = KSPSetType(ksp, KSPTCQMR);
546  AssertThrow(ierr == 0, ExcPETScError(ierr));
547 
548  // in the deal.II solvers, we always
549  // honor the initial guess in the
550  // solution vector. do so here as well:
551  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
552  AssertThrow(ierr == 0, ExcPETScError(ierr));
553  }
554 
555 
556  /* ---------------------- SolverCR ------------------------ */
557 
559  : SolverBase(cn)
560  , additional_data(data)
561  {}
562 
563 
565  const MPI_Comm &,
566  const AdditionalData &data)
567  : SolverCR(cn, data)
568  {}
569 
570 
571  void
573  {
574  PetscErrorCode ierr = KSPSetType(ksp, KSPCR);
575  AssertThrow(ierr == 0, ExcPETScError(ierr));
576 
577  // in the deal.II solvers, we always
578  // honor the initial guess in the
579  // solution vector. do so here as well:
580  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
581  AssertThrow(ierr == 0, ExcPETScError(ierr));
582  }
583 
584 
585  /* ---------------------- SolverLSQR ------------------------ */
586 
588  : SolverBase(cn)
589  , additional_data(data)
590  {}
591 
592 
593 
595  const MPI_Comm &,
596  const AdditionalData &data)
597  : SolverLSQR(cn, data)
598  {}
599 
600 
601 
602  void
604  {
605  PetscErrorCode ierr = KSPSetType(ksp, KSPLSQR);
606  AssertThrow(ierr == 0, ExcPETScError(ierr));
607 
608  // in the deal.II solvers, we always
609  // honor the initial guess in the
610  // solution vector. do so here as well:
611  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
612  AssertThrow(ierr == 0, ExcPETScError(ierr));
613  }
614 
615 
616  /* ---------------------- SolverPreOnly ------------------------ */
617 
619  : SolverBase(cn)
620  , additional_data(data)
621  {}
622 
623 
625  const MPI_Comm &,
626  const AdditionalData &data)
627  : SolverPreOnly(cn, data)
628  {}
629 
630 
631  void
633  {
634  PetscErrorCode ierr = KSPSetType(ksp, KSPPREONLY);
635  AssertThrow(ierr == 0, ExcPETScError(ierr));
636 
637  // The KSPPREONLY solver of
638  // PETSc never calls the convergence
639  // monitor, which leads to failure
640  // even when everything was ok.
641  // Therefore the SolverControl status
642  // is set to some nice values, which
643  // guarantee a nice result at the end
644  // of the solution process.
645  solver_control.check(1, 0.0);
646 
647  // Using the PREONLY solver with
648  // a nonzero initial guess leads
649  // PETSc to produce some error messages.
650  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_FALSE);
651  AssertThrow(ierr == 0, ExcPETScError(ierr));
652  }
653 
654 
655  /* ---------------------- SparseDirectMUMPS------------------------ */
656 
658  const AdditionalData &data)
659  : SolverBase(cn)
660  , additional_data(data)
661  , symmetric_mode(false)
662  {}
663 
664 
665 
667  const MPI_Comm &,
668  const AdditionalData &data)
669  : SparseDirectMUMPS(cn, data)
670  {}
671 
672 
673 
675  {
677  // the 'pc' object is owned by the 'ksp' object, and consequently
678  // does not have to be destroyed explicitly here
679  }
680 
681 
682  void
684  {
685  /*
686  * KSPPREONLY implements a stub method that applies only the
687  * preconditioner. Its use is due to SparseDirectMUMPS being a direct
688  * (rather than iterative) solver
689  */
690  PetscErrorCode ierr = KSPSetType(ksp, KSPPREONLY);
691  AssertThrow(ierr == 0, ExcPETScError(ierr));
692 
693  /*
694  * The KSPPREONLY solver of PETSc never calls the convergence monitor,
695  * which leads to failure even when everything was ok. Therefore, the
696  * SolverControl status is set to some nice values, which guarantee a
697  * nice result at the end of the solution process.
698  */
699  solver_control.check(1, 0.0);
700 
701  /*
702  * Using a PREONLY solver with a nonzero initial guess leads PETSc to
703  * produce some error messages.
704  */
705  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_FALSE);
706  AssertThrow(ierr == 0, ExcPETScError(ierr));
707  }
708 
709  void
711  VectorBase & x,
712  const VectorBase &b)
713  {
714 # ifdef DEAL_II_PETSC_WITH_MUMPS
715  /*
716  * factorization matrix to be obtained from MUMPS
717  */
718  Mat F;
719 
720  /*
721  * setting MUMPS integer control parameters ICNTL to be passed to
722  * MUMPS. Setting entry 7 of MUMPS ICNTL array (of size 40) to a value
723  * of 2. This sets use of Approximate Minimum Fill (AMF)
724  */
725  PetscInt ival = 2, icntl = 7;
726  /*
727  * number of iterations to solution (should be 1) for a direct solver
728  */
729  PetscInt its;
730  /*
731  * norm of residual
732  */
733  PetscReal rnorm;
734 
735  /*
736  * creating a solver object if this is necessary
737  */
738  if (solver_data == nullptr)
739  {
740  solver_data = std::make_unique<SolverDataMUMPS>();
741 
742  /*
743  * creates the default KSP context and puts it in the location
744  * solver_data->ksp
745  */
746  PetscErrorCode ierr =
747  KSPCreate(A.get_mpi_communicator(), &solver_data->ksp);
748  AssertThrow(ierr == 0, ExcPETScError(ierr));
749 
750  /*
751  * set the matrices involved. the last argument is irrelevant here,
752  * since we use the solver only once anyway
753  */
754  ierr = KSPSetOperators(solver_data->ksp, A, A);
755  AssertThrow(ierr == 0, ExcPETScError(ierr));
756 
757  /*
758  * setting the solver type
759  */
761 
762  /*
763  * getting the associated preconditioner context
764  */
765  ierr = KSPGetPC(solver_data->ksp, &solver_data->pc);
766  AssertThrow(ierr == 0, ExcPETScError(ierr));
767 
768  /*
769  * build PETSc PC for particular PCLU or PCCHOLESKY preconditioner
770  * depending on whether the symmetric mode has been set
771  */
772  if (symmetric_mode)
773  ierr = PCSetType(solver_data->pc, PCCHOLESKY);
774  else
775  ierr = PCSetType(solver_data->pc, PCLU);
776  AssertThrow(ierr == 0, ExcPETScError(ierr));
777 
778  /*
779  * convergence monitor function that checks with the solver_control
780  * object for convergence
781  */
782  ierr = KSPSetConvergenceTest(solver_data->ksp,
784  reinterpret_cast<void *>(&solver_control),
785  PETSC_NULL);
786  AssertThrow(ierr == 0, ExcPETScError(ierr));
787 
788  /*
789  * set the software that is to be used to perform the lu
790  * factorization here we start to see differences with the base
791  * class solve function
792  */
793 # if DEAL_II_PETSC_VERSION_LT(3, 9, 0)
794  ierr = PCFactorSetMatSolverPackage(solver_data->pc, MATSOLVERMUMPS);
795 # else
796  ierr = PCFactorSetMatSolverType(solver_data->pc, MATSOLVERMUMPS);
797 # endif
798  AssertThrow(ierr == 0, ExcPETScError(ierr));
799 
800  /*
801  * set up the package to call for the factorization
802  */
803 # if DEAL_II_PETSC_VERSION_LT(3, 9, 0)
804  ierr = PCFactorSetUpMatSolverPackage(solver_data->pc);
805 # else
806  ierr = PCFactorSetUpMatSolverType(solver_data->pc);
807 # endif
808  AssertThrow(ierr == 0, ExcPETScError(ierr));
809 
810  /*
811  * get the factored matrix F from the preconditioner context. This
812  * routine is valid only for LU, ILU, Cholesky, and incomplete
813  * Cholesky
814  */
815  ierr = PCFactorGetMatrix(solver_data->pc, &F);
816  AssertThrow(ierr == 0, ExcPETScError(ierr));
817 
818  /*
819  * Passing the control parameters to MUMPS
820  */
821  ierr = MatMumpsSetIcntl(F, icntl, ival);
822  AssertThrow(ierr == 0, ExcPETScError(ierr));
823 
824  /*
825  * set the command line option prefix name
826  */
827  ierr = KSPSetOptionsPrefix(solver_data->ksp, prefix_name.c_str());
828  AssertThrow(ierr == 0, ExcPETScError(ierr));
829 
830  /*
831  * set the command line options provided by the user to override
832  * the defaults
833  */
834  ierr = KSPSetFromOptions(solver_data->ksp);
835  AssertThrow(ierr == 0, ExcPETScError(ierr));
836  }
837 
838  /*
839  * solve the linear system
840  */
841  PetscErrorCode ierr = KSPSolve(solver_data->ksp, b, x);
842  AssertThrow(ierr == 0, ExcPETScError(ierr));
843 
844  /*
845  * in case of failure throw exception
846  */
848  {
849  AssertThrow(false,
852  }
853  else
854  {
855  /*
856  * obtain convergence information. obtain the number of iterations
857  * and residual norm
858  */
859  ierr = KSPGetIterationNumber(solver_data->ksp, &its);
860  AssertThrow(ierr == 0, ExcPETScError(ierr));
861  ierr = KSPGetResidualNorm(solver_data->ksp, &rnorm);
862  AssertThrow(ierr == 0, ExcPETScError(ierr));
863  }
864 
865 # else // DEAL_II_PETSC_WITH_MUMPS
866  Assert(
867  false,
868  ExcMessage(
869  "Your PETSc installation does not include a copy of "
870  "the MUMPS package necessary for this solver. You will need to configure "
871  "PETSc so that it includes MUMPS, recompile it, and then re-configure "
872  "and recompile deal.II as well."));
873 
874  // Cast to void to silence compiler warnings
875  (void)A;
876  (void)x;
877  (void)b;
878 # endif
879  }
880 
881 
882 
883  PetscErrorCode
885  const PetscInt iteration,
886  const PetscReal residual_norm,
887  KSPConvergedReason *reason,
888  void * solver_control_x)
889  {
891  *reinterpret_cast<SolverControl *>(solver_control_x);
892 
893  const SolverControl::State state =
894  solver_control.check(iteration, residual_norm);
895 
896  switch (state)
897  {
898  case ::SolverControl::iterate:
899  *reason = KSP_CONVERGED_ITERATING;
900  break;
901 
902  case ::SolverControl::success:
903  *reason = static_cast<KSPConvergedReason>(1);
904  break;
905 
906  case ::SolverControl::failure:
908  *reason = KSP_DIVERGED_ITS;
909  else
910  *reason = KSP_DIVERGED_DTOL;
911  break;
912 
913  default:
914  Assert(false, ExcNotImplemented());
915  }
916 
917  return 0;
918  }
919 
920 
921 
922  void
924  {
925  symmetric_mode = flag;
926  }
927 
928 } // namespace PETScWrappers
929 
931 
932 #endif // DEAL_II_WITH_PETSC
virtual void set_solver_type(KSP &ksp) const =0
SolverControl & control() const
SolverControl & solver_control
Definition: petsc_solver.h:151
void initialize(const PreconditionBase &preconditioner)
SolverBase(SolverControl &cn)
Definition: petsc_solver.cc:45
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
void set_prefix(const std::string &prefix)
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionBase &preconditioner)
Definition: petsc_solver.cc:52
std::unique_ptr< SolverData > solver_data
Definition: petsc_solver.h:216
SolverBiCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
SolverBicgstab(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverCGS(SolverControl &cn, const AdditionalData &data=AdditionalData())
SolverCG(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
SolverCR(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverChebychev(SolverControl &cn, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:498
virtual void set_solver_type(KSP &ksp) const override
SolverGMRES(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverLSQR(SolverControl &cn, const AdditionalData &data=AdditionalData())
SolverPreOnly(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:277
SolverRichardson(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverTCQMR(SolverControl &cn, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverTFQMR(SolverControl &cn, const AdditionalData &data=AdditionalData())
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b)
virtual void set_solver_type(KSP &ksp) const override
void set_symmetric_mode(const bool flag)
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
SparseDirectMUMPS(SolverControl &cn, const AdditionalData &data=AdditionalData())
std::unique_ptr< SolverDataMUMPS > solver_data
Definition: petsc_solver.h:978
State last_check() const
unsigned int last_step() const
double last_value() const
virtual State check(const unsigned int step, const double check_value)
unsigned int max_steps() const
double tolerance() const
@ success
Stop iteration, goal reached.
#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 & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1611
static const char A
PetscErrorCode destroy_krylov_solver(KSP &krylov_solver)
Tensor< 2, dim, Number > F(const Tensor< 2, dim, Number > &Grad_u)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
AdditionalData(const unsigned int restart_parameter=30, const bool right_preconditioning=false)