Reference documentation for deal.II version GIT d9d8a449a2 2022-08-17 08:45: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)
48  {}
49 
50 
51 
52  void
54  VectorBase & x,
55  const VectorBase & b,
56  const PreconditionBase &preconditioner)
57  {
58  /*
59  TODO: PETSc duplicates communicators, so this does not work (you put
60  MPI_COMM_SELF in, but get something other out when you ask PETSc for the
61  communicator. This mainly fails due to the MatrixFree classes, that can not
62  ask PETSc for a communicator. //Timo Heister
63  Assert(A.get_mpi_communicator()==mpi_communicator, ExcMessage("PETSc Solver
64  and Matrix need to use the same MPI_Comm."));
65  Assert(x.get_mpi_communicator()==mpi_communicator, ExcMessage("PETSc Solver
66  and Vector need to use the same MPI_Comm."));
67  Assert(b.get_mpi_communicator()==mpi_communicator, ExcMessage("PETSc Solver
68  and Vector need to use the same MPI_Comm."));
69  */
70 
71  // first create a solver object if this
72  // is necessary
73  if (solver_data.get() == nullptr)
74  {
75  solver_data = std::make_unique<SolverData>();
76 
77  PetscErrorCode ierr = KSPCreate(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  void
205  {
206  PetscErrorCode ierr;
207 
208  solver_data = std::make_unique<SolverData>();
209 
210  ierr = KSPCreate(mpi_communicator, &solver_data->ksp);
211  AssertThrow(ierr == 0, ExcPETScError(ierr));
212 
213  // let derived classes set the solver
214  // type, and the preconditioning
215  // object set the type of
216  // preconditioner
218 
219  ierr = KSPSetPC(solver_data->ksp, preconditioner.get_pc());
220  AssertThrow(ierr == 0, ExcPETScError(ierr));
221 
222  // then a convergence monitor
223  // function. that function simply
224  // checks with the solver_control
225  // object we have in this object for
226  // convergence
227  ierr = KSPSetConvergenceTest(solver_data->ksp,
229  reinterpret_cast<void *>(&solver_control),
230  nullptr);
231  AssertThrow(ierr == 0, ExcPETScError(ierr));
232 
233  // set the command line options provided
234  // by the user to override the defaults
235  ierr = KSPSetFromOptions(solver_data->ksp);
236  AssertThrow(ierr == 0, ExcPETScError(ierr));
237  }
238 
239 
240 
241  /* ---------------------- SolverRichardson ------------------------ */
242 
244  : omega(omega)
245  {}
246 
247 
248 
250  const MPI_Comm & mpi_communicator,
251  const AdditionalData &data)
253  , additional_data(data)
254  {}
255 
256 
257  void
259  {
260  PetscErrorCode ierr = KSPSetType(ksp, KSPRICHARDSON);
261  AssertThrow(ierr == 0, ExcPETScError(ierr));
262 
263  // set the damping factor from the data
264  ierr = KSPRichardsonSetScale(ksp, additional_data.omega);
265  AssertThrow(ierr == 0, ExcPETScError(ierr));
266 
267  // in the deal.II solvers, we always
268  // honor the initial guess in the
269  // solution vector. do so here as well:
270  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
271  AssertThrow(ierr == 0, ExcPETScError(ierr));
272 
273  // Hand over the absolute
274  // tolerance and the maximum
275  // iteration number to the PETSc
276  // convergence criterion. The
277  // custom deal.II SolverControl
278  // object is ignored by the PETSc
279  // Richardson method (when no
280  // PETSc monitoring is present),
281  // since in this case PETSc
282  // uses a faster version of
283  // the Richardson iteration,
284  // where no residual is
285  // available.
286  ierr = KSPSetTolerances(ksp,
287  PETSC_DEFAULT,
288  this->solver_control.tolerance(),
289  PETSC_DEFAULT,
290  this->solver_control.max_steps() + 1);
291  AssertThrow(ierr == 0, ExcPETScError(ierr));
292  }
293 
294 
295  /* ---------------------- SolverChebychev ------------------------ */
296 
298  const MPI_Comm & mpi_communicator,
299  const AdditionalData &data)
300  : SolverBase(cn, mpi_communicator)
301  , additional_data(data)
302  {}
303 
304 
305  void
307  {
308  PetscErrorCode ierr = KSPSetType(ksp, KSPCHEBYSHEV);
309  AssertThrow(ierr == 0, ExcPETScError(ierr));
310 
311  // in the deal.II solvers, we always
312  // honor the initial guess in the
313  // solution vector. do so here as well:
314  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
315  AssertThrow(ierr == 0, ExcPETScError(ierr));
316  }
317 
318 
319  /* ---------------------- SolverCG ------------------------ */
320 
322  const MPI_Comm & mpi_communicator,
323  const AdditionalData &data)
324  : SolverBase(cn, mpi_communicator)
325  , additional_data(data)
326  {}
327 
328 
329  void
331  {
332  PetscErrorCode ierr = KSPSetType(ksp, KSPCG);
333  AssertThrow(ierr == 0, ExcPETScError(ierr));
334 
335  // in the deal.II solvers, we always
336  // honor the initial guess in the
337  // solution vector. do so here as well:
338  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
339  AssertThrow(ierr == 0, ExcPETScError(ierr));
340  }
341 
342 
343  /* ---------------------- SolverBiCG ------------------------ */
344 
346  const MPI_Comm & mpi_communicator,
347  const AdditionalData &data)
348  : SolverBase(cn, mpi_communicator)
349  , additional_data(data)
350  {}
351 
352 
353  void
355  {
356  PetscErrorCode ierr = KSPSetType(ksp, KSPBICG);
357  AssertThrow(ierr == 0, ExcPETScError(ierr));
358 
359  // in the deal.II solvers, we always
360  // honor the initial guess in the
361  // solution vector. do so here as well:
362  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
363  AssertThrow(ierr == 0, ExcPETScError(ierr));
364  }
365 
366 
367  /* ---------------------- SolverGMRES ------------------------ */
368 
370  const unsigned int restart_parameter,
371  const bool right_preconditioning)
372  : restart_parameter(restart_parameter)
373  , right_preconditioning(right_preconditioning)
374  {}
375 
376 
377 
379  const MPI_Comm & mpi_communicator,
380  const AdditionalData &data)
382  , additional_data(data)
383  {}
384 
385 
386  void
388  {
389  PetscErrorCode ierr = KSPSetType(ksp, KSPGMRES);
390  AssertThrow(ierr == 0, ExcPETScError(ierr));
391 
392  ierr = KSPGMRESSetRestart(ksp, additional_data.restart_parameter);
393  AssertThrow(ierr == 0, ExcPETScError(ierr));
394 
395  // Set preconditioning side to right
397  {
398  ierr = KSPSetPCSide(ksp, PC_RIGHT);
399  AssertThrow(ierr == 0, ExcPETScError(ierr));
400  }
401 
402  // in the deal.II solvers, we always
403  // honor the initial guess in the
404  // solution vector. do so here as well:
405  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
406  AssertThrow(ierr == 0, ExcPETScError(ierr));
407  }
408 
409 
410  /* ---------------------- SolverBicgstab ------------------------ */
411 
413  const MPI_Comm & mpi_communicator,
414  const AdditionalData &data)
415  : SolverBase(cn, mpi_communicator)
416  , additional_data(data)
417  {}
418 
419 
420  void
422  {
423  PetscErrorCode ierr = KSPSetType(ksp, KSPBCGS);
424  AssertThrow(ierr == 0, ExcPETScError(ierr));
425 
426  // in the deal.II solvers, we always
427  // honor the initial guess in the
428  // solution vector. do so here as well:
429  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
430  AssertThrow(ierr == 0, ExcPETScError(ierr));
431  }
432 
433 
434  /* ---------------------- SolverCGS ------------------------ */
435 
437  const MPI_Comm & mpi_communicator,
438  const AdditionalData &data)
439  : SolverBase(cn, mpi_communicator)
440  , additional_data(data)
441  {}
442 
443 
444  void
446  {
447  PetscErrorCode ierr = KSPSetType(ksp, KSPCGS);
448  AssertThrow(ierr == 0, ExcPETScError(ierr));
449 
450  // in the deal.II solvers, we always
451  // honor the initial guess in the
452  // solution vector. do so here as well:
453  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
454  AssertThrow(ierr == 0, ExcPETScError(ierr));
455  }
456 
457 
458  /* ---------------------- SolverTFQMR ------------------------ */
459 
461  const MPI_Comm & mpi_communicator,
462  const AdditionalData &data)
463  : SolverBase(cn, mpi_communicator)
464  , additional_data(data)
465  {}
466 
467 
468  void
470  {
471  PetscErrorCode ierr = KSPSetType(ksp, KSPTFQMR);
472  AssertThrow(ierr == 0, ExcPETScError(ierr));
473 
474  // in the deal.II solvers, we always
475  // honor the initial guess in the
476  // solution vector. do so here as well:
477  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
478  AssertThrow(ierr == 0, ExcPETScError(ierr));
479  }
480 
481 
482  /* ---------------------- SolverTCQMR ------------------------ */
483 
485  const MPI_Comm & mpi_communicator,
486  const AdditionalData &data)
487  : SolverBase(cn, mpi_communicator)
488  , additional_data(data)
489  {}
490 
491 
492  void
494  {
495  PetscErrorCode ierr = KSPSetType(ksp, KSPTCQMR);
496  AssertThrow(ierr == 0, ExcPETScError(ierr));
497 
498  // in the deal.II solvers, we always
499  // honor the initial guess in the
500  // solution vector. do so here as well:
501  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
502  AssertThrow(ierr == 0, ExcPETScError(ierr));
503  }
504 
505 
506  /* ---------------------- SolverCR ------------------------ */
507 
509  const MPI_Comm & mpi_communicator,
510  const AdditionalData &data)
511  : SolverBase(cn, mpi_communicator)
512  , additional_data(data)
513  {}
514 
515 
516  void
518  {
519  PetscErrorCode ierr = KSPSetType(ksp, KSPCR);
520  AssertThrow(ierr == 0, ExcPETScError(ierr));
521 
522  // in the deal.II solvers, we always
523  // honor the initial guess in the
524  // solution vector. do so here as well:
525  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
526  AssertThrow(ierr == 0, ExcPETScError(ierr));
527  }
528 
529 
530  /* ---------------------- SolverLSQR ------------------------ */
531 
533  const MPI_Comm & mpi_communicator,
534  const AdditionalData &data)
535  : SolverBase(cn, mpi_communicator)
536  , additional_data(data)
537  {}
538 
539 
540  void
542  {
543  PetscErrorCode ierr = KSPSetType(ksp, KSPLSQR);
544  AssertThrow(ierr == 0, ExcPETScError(ierr));
545 
546  // in the deal.II solvers, we always
547  // honor the initial guess in the
548  // solution vector. do so here as well:
549  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_TRUE);
550  AssertThrow(ierr == 0, ExcPETScError(ierr));
551  }
552 
553 
554  /* ---------------------- SolverPreOnly ------------------------ */
555 
557  const MPI_Comm & mpi_communicator,
558  const AdditionalData &data)
559  : SolverBase(cn, mpi_communicator)
560  , additional_data(data)
561  {}
562 
563 
564  void
566  {
567  PetscErrorCode ierr = KSPSetType(ksp, KSPPREONLY);
568  AssertThrow(ierr == 0, ExcPETScError(ierr));
569 
570  // The KSPPREONLY solver of
571  // PETSc never calls the convergence
572  // monitor, which leads to failure
573  // even when everything was ok.
574  // Therefore the SolverControl status
575  // is set to some nice values, which
576  // guarantee a nice result at the end
577  // of the solution process.
578  solver_control.check(1, 0.0);
579 
580  // Using the PREONLY solver with
581  // a nonzero initial guess leads
582  // PETSc to produce some error messages.
583  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_FALSE);
584  AssertThrow(ierr == 0, ExcPETScError(ierr));
585  }
586 
587 
588  /* ---------------------- SparseDirectMUMPS------------------------ */
589 
591  {
593  // the 'pc' object is owned by the 'ksp' object, and consequently
594  // does not have to be destroyed explicitly here
595  }
596 
597 
599  const MPI_Comm & mpi_communicator,
600  const AdditionalData &data)
602  , additional_data(data)
603  , symmetric_mode(false)
604  {}
605 
606 
607  void
609  {
610  /*
611  * KSPPREONLY implements a stub method that applies only the
612  * preconditioner. Its use is due to SparseDirectMUMPS being a direct
613  * (rather than iterative) solver
614  */
615  PetscErrorCode ierr = KSPSetType(ksp, KSPPREONLY);
616  AssertThrow(ierr == 0, ExcPETScError(ierr));
617 
618  /*
619  * The KSPPREONLY solver of PETSc never calls the convergence monitor,
620  * which leads to failure even when everything was ok. Therefore, the
621  * SolverControl status is set to some nice values, which guarantee a
622  * nice result at the end of the solution process.
623  */
624  solver_control.check(1, 0.0);
625 
626  /*
627  * Using a PREONLY solver with a nonzero initial guess leads PETSc to
628  * produce some error messages.
629  */
630  ierr = KSPSetInitialGuessNonzero(ksp, PETSC_FALSE);
631  AssertThrow(ierr == 0, ExcPETScError(ierr));
632  }
633 
634  void
636  VectorBase & x,
637  const VectorBase &b)
638  {
639 # ifdef DEAL_II_PETSC_WITH_MUMPS
640  /*
641  * factorization matrix to be obtained from MUMPS
642  */
643  Mat F;
644 
645  /*
646  * setting MUMPS integer control parameters ICNTL to be passed to
647  * MUMPS. Setting entry 7 of MUMPS ICNTL array (of size 40) to a value
648  * of 2. This sets use of Approximate Minimum Fill (AMF)
649  */
650  PetscInt ival = 2, icntl = 7;
651  /*
652  * number of iterations to solution (should be 1) for a direct solver
653  */
654  PetscInt its;
655  /*
656  * norm of residual
657  */
658  PetscReal rnorm;
659 
660  /*
661  * creating a solver object if this is necessary
662  */
663  if (solver_data == nullptr)
664  {
665  solver_data = std::make_unique<SolverDataMUMPS>();
666 
667  /*
668  * creates the default KSP context and puts it in the location
669  * solver_data->ksp
670  */
671  PetscErrorCode ierr = KSPCreate(mpi_communicator, &solver_data->ksp);
672  AssertThrow(ierr == 0, ExcPETScError(ierr));
673 
674  /*
675  * set the matrices involved. the last argument is irrelevant here,
676  * since we use the solver only once anyway
677  */
678  ierr = KSPSetOperators(solver_data->ksp, A, A);
679  AssertThrow(ierr == 0, ExcPETScError(ierr));
680 
681  /*
682  * setting the solver type
683  */
685 
686  /*
687  * getting the associated preconditioner context
688  */
689  ierr = KSPGetPC(solver_data->ksp, &solver_data->pc);
690  AssertThrow(ierr == 0, ExcPETScError(ierr));
691 
692  /*
693  * build PETSc PC for particular PCLU or PCCHOLESKY preconditioner
694  * depending on whether the symmetric mode has been set
695  */
696  if (symmetric_mode)
697  ierr = PCSetType(solver_data->pc, PCCHOLESKY);
698  else
699  ierr = PCSetType(solver_data->pc, PCLU);
700  AssertThrow(ierr == 0, ExcPETScError(ierr));
701 
702  /*
703  * convergence monitor function that checks with the solver_control
704  * object for convergence
705  */
706  ierr = KSPSetConvergenceTest(solver_data->ksp,
708  reinterpret_cast<void *>(&solver_control),
709  PETSC_NULL);
710  AssertThrow(ierr == 0, ExcPETScError(ierr));
711 
712  /*
713  * set the software that is to be used to perform the lu
714  * factorization here we start to see differences with the base
715  * class solve function
716  */
717 # if DEAL_II_PETSC_VERSION_LT(3, 9, 0)
718  ierr = PCFactorSetMatSolverPackage(solver_data->pc, MATSOLVERMUMPS);
719 # else
720  ierr = PCFactorSetMatSolverType(solver_data->pc, MATSOLVERMUMPS);
721 # endif
722  AssertThrow(ierr == 0, ExcPETScError(ierr));
723 
724  /*
725  * set up the package to call for the factorization
726  */
727 # if DEAL_II_PETSC_VERSION_LT(3, 9, 0)
728  ierr = PCFactorSetUpMatSolverPackage(solver_data->pc);
729 # else
730  ierr = PCFactorSetUpMatSolverType(solver_data->pc);
731 # endif
732  AssertThrow(ierr == 0, ExcPETScError(ierr));
733 
734  /*
735  * get the factored matrix F from the preconditioner context. This
736  * routine is valid only for LU, ILU, Cholesky, and incomplete
737  * Cholesky
738  */
739  ierr = PCFactorGetMatrix(solver_data->pc, &F);
740  AssertThrow(ierr == 0, ExcPETScError(ierr));
741 
742  /*
743  * Passing the control parameters to MUMPS
744  */
745  ierr = MatMumpsSetIcntl(F, icntl, ival);
746  AssertThrow(ierr == 0, ExcPETScError(ierr));
747 
748  /*
749  * set the command line option prefix name
750  */
751  ierr = KSPSetOptionsPrefix(solver_data->ksp, prefix_name.c_str());
752  AssertThrow(ierr == 0, ExcPETScError(ierr));
753 
754  /*
755  * set the command line options provided by the user to override
756  * the defaults
757  */
758  ierr = KSPSetFromOptions(solver_data->ksp);
759  AssertThrow(ierr == 0, ExcPETScError(ierr));
760  }
761 
762  /*
763  * solve the linear system
764  */
765  PetscErrorCode ierr = KSPSolve(solver_data->ksp, b, x);
766  AssertThrow(ierr == 0, ExcPETScError(ierr));
767 
768  /*
769  * in case of failure throw exception
770  */
772  {
773  AssertThrow(false,
776  }
777  else
778  {
779  /*
780  * obtain convergence information. obtain the number of iterations
781  * and residual norm
782  */
783  ierr = KSPGetIterationNumber(solver_data->ksp, &its);
784  AssertThrow(ierr == 0, ExcPETScError(ierr));
785  ierr = KSPGetResidualNorm(solver_data->ksp, &rnorm);
786  AssertThrow(ierr == 0, ExcPETScError(ierr));
787  }
788 
789 # else // DEAL_II_PETSC_WITH_MUMPS
790  Assert(
791  false,
792  ExcMessage(
793  "Your PETSc installation does not include a copy of "
794  "the MUMPS package necessary for this solver. You will need to configure "
795  "PETSc so that it includes MUMPS, recompile it, and then re-configure "
796  "and recompile deal.II as well."));
797 
798  // Cast to void to silence compiler warnings
799  (void)A;
800  (void)x;
801  (void)b;
802 # endif
803  }
804 
805  PetscErrorCode
807  const PetscInt iteration,
808  const PetscReal residual_norm,
809  KSPConvergedReason *reason,
810  void * solver_control_x)
811  {
813  *reinterpret_cast<SolverControl *>(solver_control_x);
814 
815  const SolverControl::State state =
816  solver_control.check(iteration, residual_norm);
817 
818  switch (state)
819  {
820  case ::SolverControl::iterate:
821  *reason = KSP_CONVERGED_ITERATING;
822  break;
823 
824  case ::SolverControl::success:
825  *reason = static_cast<KSPConvergedReason>(1);
826  break;
827 
828  case ::SolverControl::failure:
830  *reason = KSP_DIVERGED_ITS;
831  else
832  *reason = KSP_DIVERGED_DTOL;
833  break;
834 
835  default:
836  Assert(false, ExcNotImplemented());
837  }
838 
839  return 0;
840  }
841 
842  void
844  {
845  symmetric_mode = flag;
846  }
847 
848 } // namespace PETScWrappers
849 
851 
852 #endif // DEAL_II_WITH_PETSC
virtual void set_solver_type(KSP &ksp) const =0
SolverControl & control() const
SolverControl & solver_control
Definition: petsc_solver.h:181
void initialize(const PreconditionBase &preconditioner)
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
const MPI_Comm mpi_communicator
Definition: petsc_solver.h:186
void set_prefix(const std::string &prefix)
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionBase &preconditioner)
Definition: petsc_solver.cc:53
std::unique_ptr< SolverData > solver_data
Definition: petsc_solver.h:251
SolverBase(SolverControl &cn, const MPI_Comm &mpi_communicator)
Definition: petsc_solver.cc:45
SolverBiCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverBicgstab(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverCGS(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, 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 MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverChebychev(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
Definition: petsc_solver.h:529
SolverGMRES(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
virtual void set_solver_type(KSP &ksp) const override
SolverLSQR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverPreOnly(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
SolverRichardson(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:311
virtual void set_solver_type(KSP &ksp) const override
SolverTCQMR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverTFQMR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const override
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b)
virtual void set_solver_type(KSP &ksp) const override
const AdditionalData additional_data
Definition: petsc_solver.h:956
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 MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
std::unique_ptr< SolverDataMUMPS > solver_data
Definition: petsc_solver.h:991
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:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
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)