Reference documentation for deal.II version GIT e22cb6a53e 2023-09-21 14: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\}}\)
trilinos_precondition.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 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_TRILINOS
19 
22 # include <deal.II/lac/vector.h>
23 
24 # include <Epetra_MultiVector.h>
25 # include <Ifpack.h>
26 # include <Ifpack_Chebyshev.h>
27 # include <Teuchos_ParameterList.hpp>
28 # include <Teuchos_RCP.hpp>
29 
31 
32 namespace TrilinosWrappers
33 {
35  : communicator(MPI_COMM_SELF)
36  {}
37 
38 
39  void
41  {
42  preconditioner.reset();
43  communicator = MPI_COMM_SELF;
44  }
45 
46 
47  MPI_Comm
49  {
50  return communicator.Comm();
51  }
52 
53 
56  {
57  AssertThrow(!preconditioner.is_null(),
58  ExcMessage("Trying to dereference a null pointer."));
59  return (*preconditioner);
60  }
61 
62 
63  IndexSet
65  {
66  return IndexSet(preconditioner->OperatorDomainMap());
67  }
68 
69 
70  IndexSet
72  {
73  return IndexSet(preconditioner->OperatorRangeMap());
74  }
75 
76  /* -------------------------- PreconditionJacobi -------------------------- */
77 
79  const double omega,
80  const double min_diagonal,
81  const unsigned int n_sweeps)
82  : omega(omega)
83  , min_diagonal(min_diagonal)
84  , n_sweeps(n_sweeps)
85  {}
86 
87 
88 
89  void
91  const AdditionalData &additional_data)
92  {
93  // release memory before reallocation
94  preconditioner.reset();
95  preconditioner.reset(
96  Ifpack().Create("point relaxation",
97  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
98  0));
99 
100  Ifpack_Preconditioner *ifpack =
101  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
102  Assert(ifpack != nullptr,
103  ExcMessage("Trilinos could not create this "
104  "preconditioner"));
105 
106  int ierr;
107 
108  Teuchos::ParameterList parameter_list;
109  parameter_list.set("relaxation: sweeps",
110  static_cast<int>(additional_data.n_sweeps));
111  parameter_list.set("relaxation: type", "Jacobi");
112  parameter_list.set("relaxation: damping factor", additional_data.omega);
113  parameter_list.set("relaxation: min diagonal value",
114  additional_data.min_diagonal);
115 
116  ierr = ifpack->SetParameters(parameter_list);
117  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
118 
119  ierr = ifpack->Initialize();
120  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
121 
122  ierr = ifpack->Compute();
123  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
124  }
125 
126 
127 
128  /* -------------------------- PreconditionSSOR -------------------------- */
129 
131  const double min_diagonal,
132  const unsigned int overlap,
133  const unsigned int n_sweeps)
134  : omega(omega)
135  , min_diagonal(min_diagonal)
136  , overlap(overlap)
137  , n_sweeps(n_sweeps)
138  {}
139 
140 
141 
142  void
144  const AdditionalData &additional_data)
145  {
146  preconditioner.reset();
147  preconditioner.reset(
148  Ifpack().Create("point relaxation",
149  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
150  additional_data.overlap));
151 
152  Ifpack_Preconditioner *ifpack =
153  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
154  Assert(ifpack != nullptr,
155  ExcMessage("Trilinos could not create this "
156  "preconditioner"));
157 
158  int ierr;
159 
160  Teuchos::ParameterList parameter_list;
161  parameter_list.set("relaxation: sweeps",
162  static_cast<int>(additional_data.n_sweeps));
163  parameter_list.set("relaxation: type", "symmetric Gauss-Seidel");
164  parameter_list.set("relaxation: damping factor", additional_data.omega);
165  parameter_list.set("relaxation: min diagonal value",
166  additional_data.min_diagonal);
167  parameter_list.set("schwarz: combine mode", "Add");
168 
169  ierr = ifpack->SetParameters(parameter_list);
170  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
171 
172  ierr = ifpack->Initialize();
173  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
174 
175  ierr = ifpack->Compute();
176  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
177  }
178 
179 
180 
181  /* -------------------------- PreconditionSOR -------------------------- */
182 
184  const double min_diagonal,
185  const unsigned int overlap,
186  const unsigned int n_sweeps)
187  : omega(omega)
188  , min_diagonal(min_diagonal)
189  , overlap(overlap)
190  , n_sweeps(n_sweeps)
191  {}
192 
193 
194 
195  void
197  const AdditionalData &additional_data)
198  {
199  preconditioner.reset();
200  preconditioner.reset(
201  Ifpack().Create("point relaxation",
202  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
203  additional_data.overlap));
204 
205  Ifpack_Preconditioner *ifpack =
206  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
207  Assert(ifpack != nullptr,
208  ExcMessage("Trilinos could not create this "
209  "preconditioner"));
210 
211  int ierr;
212 
213  Teuchos::ParameterList parameter_list;
214  parameter_list.set("relaxation: sweeps",
215  static_cast<int>(additional_data.n_sweeps));
216  parameter_list.set("relaxation: type", "Gauss-Seidel");
217  parameter_list.set("relaxation: damping factor", additional_data.omega);
218  parameter_list.set("relaxation: min diagonal value",
219  additional_data.min_diagonal);
220  parameter_list.set("schwarz: combine mode", "Add");
221 
222  ierr = ifpack->SetParameters(parameter_list);
223  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
224 
225  ierr = ifpack->Initialize();
226  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
227 
228  ierr = ifpack->Compute();
229  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
230  }
231 
232 
233 
234  /* ----------------------- PreconditionBlockJacobi ---------------------- */
235 
237  const unsigned int block_size,
238  const std::string &block_creation_type,
239  const double omega,
240  const double min_diagonal,
241  const unsigned int n_sweeps)
243  , block_creation_type(block_creation_type)
244  , omega(omega)
245  , min_diagonal(min_diagonal)
246  , n_sweeps(n_sweeps)
247  {}
248 
249 
250 
251  void
253  const AdditionalData &additional_data)
254  {
255  // release memory before reallocation
256  preconditioner.reset();
257 
258  // Block relaxation setup fails if we have no locally owned rows. As a
259  // work-around we just pretend to use point relaxation on those processors:
260  preconditioner.reset(
261  Ifpack().Create((matrix.trilinos_matrix().NumMyRows() == 0) ?
262  "point relaxation" :
263  "block relaxation",
264  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
265  0));
266 
267  Ifpack_Preconditioner *ifpack =
268  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
269  Assert(ifpack != nullptr,
270  ExcMessage("Trilinos could not create this "
271  "preconditioner"));
272 
273  int ierr;
274 
275  Teuchos::ParameterList parameter_list;
276  parameter_list.set("relaxation: sweeps",
277  static_cast<int>(additional_data.n_sweeps));
278  parameter_list.set("relaxation: type", "Jacobi");
279  parameter_list.set("relaxation: damping factor", additional_data.omega);
280  parameter_list.set("relaxation: min diagonal value",
281  additional_data.min_diagonal);
282  parameter_list.set("partitioner: type",
283  additional_data.block_creation_type);
284  int n_local_parts =
285  (matrix.trilinos_matrix().NumMyRows() + additional_data.block_size - 1) /
286  additional_data.block_size;
287  parameter_list.set("partitioner: local parts", n_local_parts);
288 
289  ierr = ifpack->SetParameters(parameter_list);
290  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
291 
292  ierr = ifpack->Initialize();
293  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
294 
295  ierr = ifpack->Compute();
296  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
297  }
298 
299 
300 
301  /* ----------------------- PreconditionBlockSSOR ------------------------ */
302 
304  const unsigned int block_size,
305  const std::string &block_creation_type,
306  const double omega,
307  const double min_diagonal,
308  const unsigned int overlap,
309  const unsigned int n_sweeps)
311  , block_creation_type(block_creation_type)
312  , omega(omega)
313  , min_diagonal(min_diagonal)
314  , overlap(overlap)
315  , n_sweeps(n_sweeps)
316  {}
317 
318 
319 
320  void
322  const AdditionalData &additional_data)
323  {
324  preconditioner.reset();
325 
326  // Block relaxation setup fails if we have no locally owned rows. As a
327  // work-around we just pretend to use point relaxation on those processors:
328  preconditioner.reset(
329  Ifpack().Create((matrix.trilinos_matrix().NumMyRows() == 0) ?
330  "point relaxation" :
331  "block relaxation",
332  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
333  additional_data.overlap));
334 
335  Ifpack_Preconditioner *ifpack =
336  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
337  Assert(ifpack != nullptr,
338  ExcMessage("Trilinos could not create this "
339  "preconditioner"));
340 
341  int ierr;
342 
343  Teuchos::ParameterList parameter_list;
344  parameter_list.set("relaxation: sweeps",
345  static_cast<int>(additional_data.n_sweeps));
346  parameter_list.set("relaxation: type", "symmetric Gauss-Seidel");
347  parameter_list.set("relaxation: damping factor", additional_data.omega);
348  parameter_list.set("relaxation: min diagonal value",
349  additional_data.min_diagonal);
350  parameter_list.set("schwarz: combine mode", "Add");
351  parameter_list.set("partitioner: type",
352  additional_data.block_creation_type);
353  int n_local_parts =
354  (matrix.trilinos_matrix().NumMyRows() + additional_data.block_size - 1) /
355  additional_data.block_size;
356  parameter_list.set("partitioner: local parts", n_local_parts);
357 
358  ierr = ifpack->SetParameters(parameter_list);
359  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
360 
361  ierr = ifpack->Initialize();
362  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
363 
364  ierr = ifpack->Compute();
365  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
366  }
367 
368 
369 
370  /* ------------------------ PreconditionBlockSOR ------------------------ */
371 
373  const unsigned int block_size,
374  const std::string &block_creation_type,
375  const double omega,
376  const double min_diagonal,
377  const unsigned int overlap,
378  const unsigned int n_sweeps)
380  , block_creation_type(block_creation_type)
381  , omega(omega)
382  , min_diagonal(min_diagonal)
383  , overlap(overlap)
384  , n_sweeps(n_sweeps)
385  {}
386 
387 
388 
389  void
391  const AdditionalData &additional_data)
392  {
393  preconditioner.reset();
394 
395  // Block relaxation setup fails if we have no locally owned rows. As a
396  // work-around we just pretend to use point relaxation on those processors:
397  preconditioner.reset(
398  Ifpack().Create((matrix.trilinos_matrix().NumMyRows() == 0) ?
399  "point relaxation" :
400  "block relaxation",
401  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
402  additional_data.overlap));
403 
404  Ifpack_Preconditioner *ifpack =
405  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
406  Assert(ifpack != nullptr,
407  ExcMessage("Trilinos could not create this "
408  "preconditioner"));
409 
410  int ierr;
411 
412  Teuchos::ParameterList parameter_list;
413  parameter_list.set("relaxation: sweeps",
414  static_cast<int>(additional_data.n_sweeps));
415  parameter_list.set("relaxation: type", "Gauss-Seidel");
416  parameter_list.set("relaxation: damping factor", additional_data.omega);
417  parameter_list.set("relaxation: min diagonal value",
418  additional_data.min_diagonal);
419  parameter_list.set("schwarz: combine mode", "Add");
420  parameter_list.set("partitioner: type",
421  additional_data.block_creation_type);
422  int n_local_parts =
423  (matrix.trilinos_matrix().NumMyRows() + additional_data.block_size - 1) /
424  additional_data.block_size;
425  parameter_list.set("partitioner: local parts", n_local_parts);
426 
427  ierr = ifpack->SetParameters(parameter_list);
428  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
429 
430  ierr = ifpack->Initialize();
431  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
432 
433  ierr = ifpack->Compute();
434  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
435  }
436 
437 
438 
439  /* -------------------------- PreconditionIC -------------------------- */
440 
442  const double ic_atol,
443  const double ic_rtol,
444  const unsigned int overlap)
445  : ic_fill(ic_fill)
446  , ic_atol(ic_atol)
447  , ic_rtol(ic_rtol)
448  , overlap(overlap)
449  {}
450 
451 
452 
453  void
455  const AdditionalData &additional_data)
456  {
457  preconditioner.reset();
458  preconditioner.reset(
459  Ifpack().Create("IC",
460  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
461  additional_data.overlap));
462 
463  Ifpack_Preconditioner *ifpack =
464  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
465  Assert(ifpack != nullptr,
466  ExcMessage("Trilinos could not create this "
467  "preconditioner"));
468 
469  int ierr;
470 
471  Teuchos::ParameterList parameter_list;
472  parameter_list.set("fact: level-of-fill", additional_data.ic_fill);
473  parameter_list.set("fact: absolute threshold", additional_data.ic_atol);
474  parameter_list.set("fact: relative threshold", additional_data.ic_rtol);
475  parameter_list.set("schwarz: combine mode", "Add");
476 
477  ierr = ifpack->SetParameters(parameter_list);
478  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
479 
480  ierr = ifpack->Initialize();
481  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
482 
483  ierr = ifpack->Compute();
484  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
485  }
486 
487 
488 
489  /* -------------------------- PreconditionILU -------------------------- */
490 
492  const double ilu_atol,
493  const double ilu_rtol,
494  const unsigned int overlap)
495  : ilu_fill(ilu_fill)
496  , ilu_atol(ilu_atol)
497  , ilu_rtol(ilu_rtol)
498  , overlap(overlap)
499  {}
500 
501 
502 
503  void
505  const AdditionalData &additional_data)
506  {
507  preconditioner.reset();
508  preconditioner.reset(
509  Ifpack().Create("ILU",
510  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
511  additional_data.overlap));
512 
513  Ifpack_Preconditioner *ifpack =
514  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
515  Assert(ifpack != nullptr,
516  ExcMessage("Trilinos could not create this "
517  "preconditioner"));
518 
519  int ierr;
520 
521  Teuchos::ParameterList parameter_list;
522  parameter_list.set("fact: level-of-fill",
523  static_cast<int>(additional_data.ilu_fill));
524  parameter_list.set("fact: absolute threshold", additional_data.ilu_atol);
525  parameter_list.set("fact: relative threshold", additional_data.ilu_rtol);
526  parameter_list.set("schwarz: combine mode", "Add");
527 
528  ierr = ifpack->SetParameters(parameter_list);
529  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
530 
531  ierr = ifpack->Initialize();
532  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
533 
534  ierr = ifpack->Compute();
535  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
536  }
537 
538 
539 
540  /* -------------------------- PreconditionILUT -------------------------- */
541 
543  const unsigned int ilut_fill,
544  const double ilut_atol,
545  const double ilut_rtol,
546  const unsigned int overlap)
547  : ilut_drop(ilut_drop)
548  , ilut_fill(ilut_fill)
549  , ilut_atol(ilut_atol)
550  , ilut_rtol(ilut_rtol)
551  , overlap(overlap)
552  {}
553 
554 
555 
556  void
558  const AdditionalData &additional_data)
559  {
560  preconditioner.reset();
561  preconditioner.reset(
562  Ifpack().Create("ILUT",
563  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
564  additional_data.overlap));
565 
566  Ifpack_Preconditioner *ifpack =
567  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
568  Assert(ifpack != nullptr,
569  ExcMessage("Trilinos could not create this "
570  "preconditioner"));
571 
572  int ierr;
573 
574  Teuchos::ParameterList parameter_list;
575  parameter_list.set("fact: drop value", additional_data.ilut_drop);
576  parameter_list.set("fact: level-of-fill",
577  static_cast<int>(additional_data.ilut_fill));
578  parameter_list.set("fact: absolute threshold", additional_data.ilut_atol);
579  parameter_list.set("fact: relative threshold", additional_data.ilut_rtol);
580  parameter_list.set("schwarz: combine mode", "Add");
581 
582  ierr = ifpack->SetParameters(parameter_list);
583  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
584 
585  ierr = ifpack->Initialize();
586  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
587 
588  ierr = ifpack->Compute();
589  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
590  }
591 
592 
593 
594  /* ---------------------- PreconditionBlockDirect --------------------- */
595 
597  const unsigned int overlap)
598  : overlap(overlap)
599  {}
600 
601 
602 
603  void
605  const AdditionalData &additional_data)
606  {
607  preconditioner.reset();
608  preconditioner.reset(
609  Ifpack().Create("Amesos",
610  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
611  additional_data.overlap));
612 
613  Ifpack_Preconditioner *ifpack =
614  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
615  Assert(ifpack != nullptr,
616  ExcMessage("Trilinos could not create this "
617  "preconditioner"));
618 
619  int ierr;
620 
621  Teuchos::ParameterList parameter_list;
622  parameter_list.set("schwarz: combine mode", "Add");
623 
624  ierr = ifpack->SetParameters(parameter_list);
625  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
626 
627  ierr = ifpack->Initialize();
628  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
629 
630  ierr = ifpack->Compute();
631  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
632  }
633 
634 
635 
636  /* ---------------------- PreconditionBlockDirect --------------------- */
637 
639  const unsigned int degree,
640  const double max_eigenvalue,
641  const double eigenvalue_ratio,
642  const double min_eigenvalue,
643  const double min_diagonal,
644  const bool nonzero_starting)
645  : degree(degree)
646  , max_eigenvalue(max_eigenvalue)
647  , eigenvalue_ratio(eigenvalue_ratio)
648  , min_eigenvalue(min_eigenvalue)
649  , min_diagonal(min_diagonal)
650  , nonzero_starting(nonzero_starting)
651  {}
652 
653 
654 
655  void
657  const AdditionalData &additional_data)
658  {
660  Teuchos::rcp(new Ifpack_Chebyshev(&matrix.trilinos_matrix()));
661 
662  Ifpack_Chebyshev *ifpack =
663  dynamic_cast<Ifpack_Chebyshev *>(preconditioner.get());
664  Assert(ifpack != nullptr,
665  ExcMessage("Trilinos could not create this "
666  "preconditioner"));
667 
668  int ierr;
669 
670  Teuchos::ParameterList parameter_list;
671  parameter_list.set("chebyshev: ratio eigenvalue",
672  additional_data.eigenvalue_ratio);
673  parameter_list.set("chebyshev: min eigenvalue",
674  additional_data.min_eigenvalue);
675  parameter_list.set("chebyshev: max eigenvalue",
676  additional_data.max_eigenvalue);
677  parameter_list.set("chebyshev: degree",
678  static_cast<int>(additional_data.degree));
679  parameter_list.set("chebyshev: min diagonal value",
680  additional_data.min_diagonal);
681  parameter_list.set("chebyshev: zero starting solution",
682  !additional_data.nonzero_starting);
683 
684  ierr = ifpack->SetParameters(parameter_list);
685  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
686 
687  ierr = ifpack->Initialize();
688  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
689 
690  ierr = ifpack->Compute();
691  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
692  }
693 
694 
695 
696  /* -------------------------- PreconditionIdentity --------------------- */
697 
698  void
700  const AdditionalData &)
701  {
702  // What follows just configures a dummy preconditioner that
703  // sets up the domain and range maps, as well as the communicator.
704  // It is never used as the vmult, Tvmult operations are
705  // given a custom definition.
706  // Note: This is only required in order to wrap this
707  // preconditioner in a LinearOperator without an exemplar
708  // matrix.
709 
710  // From PreconditionJacobi:
711  // release memory before reallocation
712  preconditioner.reset();
713  preconditioner.reset(
714  Ifpack().Create("point relaxation",
715  const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()),
716  0));
717 
718  Ifpack_Preconditioner *ifpack =
719  dynamic_cast<Ifpack_Preconditioner *>(preconditioner.get());
720  Assert(ifpack != nullptr,
721  ExcMessage("Trilinos could not create this "
722  "preconditioner"));
723 
724  int ierr;
725 
726  Teuchos::ParameterList parameter_list;
727  parameter_list.set("relaxation: sweeps", 1);
728  parameter_list.set("relaxation: type", "Jacobi");
729  parameter_list.set("relaxation: damping factor", 1.0);
730  parameter_list.set("relaxation: min diagonal value", 0.0);
731 
732  ierr = ifpack->SetParameters(parameter_list);
733  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
734 
735  ierr = ifpack->Initialize();
736  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
737 
738  ierr = ifpack->Compute();
739  AssertThrow(ierr == 0, ExcTrilinosError(ierr));
740  }
741 
742  void
744  {
745  dst = src;
746  }
747 
748  void
750  {
751  dst = src;
752  }
753 
754  void
756  const ::Vector<double> &src) const
757  {
758  dst = src;
759  }
760 
761  void
763  const ::Vector<double> &src) const
764  {
765  dst = src;
766  }
767 
768 # ifndef DOXYGEN
769  void
773  {
774  dst = src;
775  }
776 
777  void
781  {
782  dst = src;
783  }
784 # endif // DOXYGEN
785 } // namespace TrilinosWrappers
786 
788 
789 #endif // DEAL_II_WITH_TRILINOS
Epetra_Operator & trilinos_operator() const
Teuchos::RCP< Epetra_Operator > preconditioner
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void vmult(MPI::Vector &dst, const MPI::Vector &src) const override
void Tvmult(MPI::Vector &dst, const MPI::Vector &src) const override
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
void initialize(const SparseMatrix &matrix, const AdditionalData &additional_data=AdditionalData())
Definition: vector.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcTrilinosError(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1705
constexpr int block_size
Definition: cuda_size.h:29
@ matrix
Contents is actually a matrix.
AdditionalData(const unsigned int block_size=1, const std::string &block_creation_type="linear", const double omega=1, const double min_diagonal=0, const unsigned int n_sweeps=1)
AdditionalData(const unsigned int block_size=1, const std::string &block_creation_type="linear", const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)
AdditionalData(const unsigned int block_size=1, const std::string &block_creation_type="linear", const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)
AdditionalData(const unsigned int degree=1, const double max_eigenvalue=10., const double eigenvalue_ratio=30., const double min_eigenvalue=1., const double min_diagonal=1e-12, const bool nonzero_starting=false)
AdditionalData(const unsigned int ic_fill=0, const double ic_atol=0., const double ic_rtol=1., const unsigned int overlap=0)
AdditionalData(const double ilut_drop=0., const unsigned int ilut_fill=0, const double ilut_atol=0., const double ilut_rtol=1., const unsigned int overlap=0)
AdditionalData(const unsigned int ilu_fill=0, const double ilu_atol=0., const double ilu_rtol=1., const unsigned int overlap=0)
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int n_sweeps=1)
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)
AdditionalData(const double omega=1, const double min_diagonal=0, const unsigned int overlap=0, const unsigned int n_sweeps=1)