]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Enable deal.II solver use with Trilinos LinearOperator
authorJean-Paul Pelteret <jppelteret@gmail.com>
Wed, 8 Feb 2017 08:04:22 +0000 (09:04 +0100)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Thu, 9 Feb 2017 09:23:10 +0000 (10:23 +0100)
By making a distinction between an `inverse_payload` created for
Trilinos solvers/preconditioners and all other solver/preconditioner
types, we allow a `TrilinosPayload` to be created when the
solver/preconditioner is incompatible with `Epetra_Multivector`. In
turn, it is now possible to use the deal.II solvers (e.g.
`dealii::SolverCG<TrilinosWrappers::MPI::Vector>`) with `Trilinos`
`LinearOperator`s as such a call no longer results in a compiler error.

include/deal.II/lac/trilinos_sparse_matrix.h
tests/lac/linear_operator_10a.cc [new file with mode: 0644]
tests/lac/linear_operator_10a.with_cxx11=on.with_trilinos=on.output [new file with mode: 0644]
tests/lac/linear_operator_12a.cc [new file with mode: 0644]
tests/lac/linear_operator_12a.with_cxx11=on.with_trilinos=on.with_mpi=true.mpirun=1.output [new file with mode: 0644]
tests/lac/linear_operator_12a.with_cxx11=on.with_trilinos=on.with_mpi=true.mpirun=2.output [new file with mode: 0644]

index 7d328c3812c5d5f22aac2d3e078e86f3deed7b1d..d7c55438cc83b3ae96887ec7e624cb654d813269 100644 (file)
@@ -31,6 +31,7 @@
 
 #ifdef DEAL_II_WITH_CXX11
 #include <deal.II/lac/vector_memory.h>
+#include <type_traits>
 #endif // DEAL_II_WITH_CXX11
 
 #  include <vector>
@@ -2049,6 +2050,7 @@ namespace TrilinosWrappers
 #ifdef DEAL_II_WITH_CXX11
 
   // forwards declarations
+  class SolverBase;
   class PreconditionBase;
 
   namespace internal
@@ -2216,17 +2218,50 @@ namespace TrilinosWrappers
         TrilinosPayload transpose_payload () const;
 
         /**
-        * Returns a payload configured for inverse operations
-        *
-        * Invoking this constructor will configure two additional functions,
-        * namely <tt>inv_vmult</tt> and <tt>inv_Tvmult</tt>, both of which wrap
-        * inverse operations.
-        * The <tt>vmult</tt> and <tt>Tvmult</tt> operations retain the standard
-        * definitions inherited from @p op.
-        */
+         * Returns a payload configured for inverse operations
+         *
+         * Invoking this factory function will configure two additional functions,
+         * namely <tt>inv_vmult</tt> and <tt>inv_Tvmult</tt>, both of which wrap
+         * inverse operations.
+         * The <tt>vmult</tt> and <tt>Tvmult</tt> operations retain the standard
+         * definitions inherited from @p op.
+         *
+         * @note This function is enabled only if the solver and preconditioner
+         * derive from the respective TrilinosWrappers base classes.
+         * The C++ compiler will therefore only consider this function if the
+         * following criterion are satisfied:
+         * 1. the @p Solver derives from TrilinosWrappers::SolverBase, and
+         * 2. the @p Preconditioner derives from TrilinosWrappers::PreconditionBase.
+         */
         template <typename Solver, typename Preconditioner>
-        TrilinosPayload inverse_payload (Solver &, const Preconditioner &) const;
+        typename std::enable_if<
+        std::is_base_of<TrilinosWrappers::SolverBase,Solver>::value &&
+        std::is_base_of<TrilinosWrappers::PreconditionBase,Preconditioner>::value,
+            TrilinosPayload>::type
+            inverse_payload (Solver &, const Preconditioner &) const;
 
+        /**
+         * Returns a payload configured for inverse operations
+         *
+         * Invoking this factory function will configure two additional functions,
+         * namely <tt>inv_vmult</tt> and <tt>inv_Tvmult</tt>, both of which
+         * are disabled because the @p Solver or @p Preconditioner are not
+         * compatible with Epetra_MultiVector.
+         * The <tt>vmult</tt> and <tt>Tvmult</tt> operations retain the standard
+         * definitions inherited from @p op.
+         *
+         * @note The C++ compiler will only consider this function if the
+         * following criterion are satisfied:
+         * 1. the @p Solver does not derive from TrilinosWrappers::SolverBase, and
+         * 2. the @p Preconditioner does not derive from
+         * TrilinosWrappers::PreconditionBase.
+         */
+        template <typename Solver, typename Preconditioner>
+        typename std::enable_if<
+        !(std::is_base_of<TrilinosWrappers::SolverBase,Solver>::value &&
+          std::is_base_of<TrilinosWrappers::PreconditionBase,Preconditioner>::value),
+              TrilinosPayload>::type
+              inverse_payload (Solver &, const Preconditioner &) const;
 
 //@}
 
@@ -3144,10 +3179,13 @@ namespace TrilinosWrappers
     namespace LinearOperator
     {
       template <typename Solver, typename Preconditioner>
-      TrilinosPayload
-      TrilinosPayload::inverse_payload (
-        Solver                &solver,
-        const Preconditioner  &preconditioner) const
+      typename std::enable_if<
+      std::is_base_of<TrilinosWrappers::SolverBase,Solver>::value &&
+      std::is_base_of<TrilinosWrappers::PreconditionBase,Preconditioner>::value,
+          TrilinosPayload>::type
+          TrilinosPayload::inverse_payload (
+            Solver                &solver,
+            const Preconditioner  &preconditioner) const
       {
         const auto &payload = *this;
 
@@ -3192,6 +3230,36 @@ namespace TrilinosWrappers
 
         return return_op;
       }
+
+      template <typename Solver, typename Preconditioner>
+      typename std::enable_if<
+      !(std::is_base_of<TrilinosWrappers::SolverBase,Solver>::value &&
+        std::is_base_of<TrilinosWrappers::PreconditionBase,Preconditioner>::value),
+            TrilinosPayload>::type
+            TrilinosPayload::inverse_payload (
+              Solver                &solver,
+              const Preconditioner  &preconditioner) const
+      {
+        TrilinosPayload return_op(*this);
+
+        return_op.inv_vmult = [](TrilinosPayload::Domain &,
+                                 const TrilinosPayload::Range &)
+        {
+          AssertThrow(false,
+                      ExcMessage("Payload inv_vmult disabled because of "
+                                 "incompatible solver/preconditioner choice."));
+        };
+
+        return_op.inv_Tvmult = [](TrilinosPayload::Range &,
+                                  const TrilinosPayload::Domain &)
+        {
+          AssertThrow(false,
+                      ExcMessage("Payload inv_vmult disabled because of "
+                                 "incompatible solver/preconditioner choice."));
+        };
+
+        return return_op;
+      }
     } // namespace LinearOperator
   } // namespace internal
 #endif // DEAL_II_WITH_CXX11
diff --git a/tests/lac/linear_operator_10a.cc b/tests/lac/linear_operator_10a.cc
new file mode 100644 (file)
index 0000000..035bcd8
--- /dev/null
@@ -0,0 +1,385 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Test internal preconditioner and solver options
+
+#include "../tests.h"
+
+#include <deal.II/lac/trilinos_linear_operator.h>
+#include <deal.II/lac/packaged_operation.h>
+
+#include <deal.II/lac/solver_bicgstab.h>
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/lac/solver_gmres.h>
+#include <deal.II/lac/solver_minres.h>
+#include <deal.II/lac/solver_qmrs.h>
+#include <deal.II/lac/solver_richardson.h>
+
+#include <deal.II/lac/trilinos_block_sparse_matrix.h>
+#include <deal.II/lac/trilinos_block_vector.h>
+#include <deal.II/lac/trilinos_precondition.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_vector.h>
+
+
+using namespace dealii;
+
+template<typename VECTOR>
+void
+print (const VECTOR &vec)
+{
+  for (types::global_dof_index i=0; i < vec.size(); ++i)
+    {
+      deallog << vec(i) << " ";
+    }
+  deallog << std::endl;
+}
+
+
+template<class PRECONDITIONER, class MATRIX, class VECTOR,
+         class ADDITIONAL_DATA = typename PRECONDITIONER::AdditionalData>
+void
+test_preconditioner (const MATRIX &A,
+                     const VECTOR &b,
+                     const ADDITIONAL_DATA &data = ADDITIONAL_DATA())
+{
+  const auto lo_A = linear_operator<VECTOR>(A);
+  // Note: The above should be equivalent to the following:
+  //
+  //  typedef dealii::TrilinosWrappers::internal::LinearOperator::TrilinosPayload PAYLOAD;
+  //  const auto lo_A = linear_operator<VECTOR,VECTOR,PAYLOAD>(A);
+
+  PRECONDITIONER preconditioner;
+  preconditioner.initialize(A, data);
+
+  typedef SolverCG<VECTOR> SOLVER;
+  SolverControl solver_control (100, 1.0e-10, false,false);
+  SOLVER solver (solver_control);
+
+  // Exact inverse
+  const auto lo_A_inv = inverse_operator(lo_A, solver, preconditioner);
+  // Note: The above should be equivalent to the following:
+  //
+  //  const auto lo_A_inv = inverse_operator<PAYLOAD,SOLVER,
+  //             PRECONDITIONER,
+  //             VECTOR,VECTOR>(lo_A,
+  //                            solver,
+  //                            preconditioner);
+
+  // Singular operation
+  {
+    deallog.push("S_Op");
+    const VECTOR x = lo_A_inv*b;
+    print(x);
+    deallog.pop();
+  }
+
+  // Composite operation
+  {
+    deallog.push("C_Op");
+    const VECTOR x = (lo_A_inv*lo_A*lo_A_inv)*b;
+    print(x);
+    deallog.pop();
+  }
+
+  // Approximate inverse
+  deallog.push("Approx");
+  {
+    // Using exemplar matrix
+    deallog.push("Exemp");
+    const auto lo_A_inv_approx = linear_operator<VECTOR,VECTOR>(A, preconditioner);
+    // Note: The above should be equivalent to the following:
+    //
+    //    typedef dealii::TrilinosWrappers::internal::LinearOperator::TrilinosPayload PAYLOAD;
+    //    const auto lo_A_inv_approx = linear_operator<VECTOR,VECTOR,PAYLOAD>(A, preconditioner);
+
+    // Singular operation
+    {
+      deallog.push("S_Op");
+      const VECTOR x_approx = lo_A_inv_approx*b;
+      print(x_approx);
+      deallog.pop();
+    }
+
+    // Composite operation
+    {
+      deallog.push("C_Op");
+      const VECTOR x_approx = (lo_A_inv_approx*lo_A*lo_A_inv_approx)*b;
+      print(x_approx);
+      deallog.pop();
+    }
+
+    deallog.pop();
+  }
+  {
+    // Stand-alone
+    deallog.push("S.A.");
+    typedef dealii::TrilinosWrappers::internal::LinearOperator::TrilinosPayload PAYLOAD;
+    const auto lo_A_inv_approx = linear_operator<VECTOR,VECTOR,PAYLOAD>(preconditioner);
+
+    // Singular operation
+    {
+      deallog.push("S_Op");
+      const VECTOR x_approx = lo_A_inv_approx*b;
+      print(x_approx);
+      deallog.pop();
+    }
+
+    // Composite operation
+    {
+      deallog.push("C_Op");
+      const VECTOR x_approx = (lo_A_inv_approx*lo_A*lo_A_inv_approx)*b;
+      print(x_approx);
+      deallog.pop();
+    }
+
+    deallog.pop();
+  }
+  deallog.pop();
+}
+
+template<class SOLVER, class MATRIX, class VECTOR>
+void
+test_solver (const MATRIX &A,
+             const VECTOR &b)
+{
+  const auto lo_A = linear_operator<VECTOR>(A);
+  // Note: The above should be equivalent to the following:
+  //
+  //  typedef dealii::TrilinosWrappers::internal::LinearOperator::TrilinosPayload PAYLOAD;
+  //  const auto lo_A = linear_operator<VECTOR,VECTOR,PAYLOAD>(A);
+
+  SolverControl solver_control (100, 1.0e-10, false,false);
+  SOLVER solver (solver_control);
+
+  typedef TrilinosWrappers::PreconditionJacobi PRECONDITIONER;
+  PRECONDITIONER preconditioner;
+  preconditioner.initialize(A);
+
+  {
+    const auto lo_A_inv = inverse_operator(lo_A, solver, preconditioner);
+    // Note: The above should be equivalent to the following:
+    //
+    //  const auto lo_A_inv = inverse_operator<PAYLOAD,SOLVER,
+    //             PRECONDITIONER,
+    //             VECTOR,VECTOR>(lo_A,
+    //                            solver,
+    //                            preconditioner);
+
+    // Singular operation
+    {
+      deallog.push("S_Op");
+      const VECTOR x_approx = lo_A_inv*b;
+      print(x_approx);
+      deallog.pop();
+    }
+
+    // Composite operation
+    {
+      deallog.push("C_Op");
+      const VECTOR x_approx = (lo_A_inv*lo_A*lo_A_inv)*b;
+      print(x_approx);
+      deallog.pop();
+    }
+  }
+
+  // Composite operation 2
+  {
+    deallog.push("C_Op2");
+    SolverControl solver_control_1 (100, 1.0e-10, false,false);
+    SOLVER solver_1 (solver_control_1);
+    const auto lo_A_inv_1 = inverse_operator(lo_A, solver_1, preconditioner);
+    SolverControl solver_control_2 (100, 1.0e-10, false,false);
+    SOLVER solver_2 (solver_control_2);
+    const auto lo_A_inv_2 = inverse_operator(lo_A, solver_2, preconditioner);
+    const VECTOR x_approx = (lo_A_inv_2*lo_A*lo_A_inv_1)*b;
+    print(x_approx);
+    deallog.pop();
+  }
+}
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll all;
+
+  deallog.depth_console(0);
+  deallog << std::setprecision(10);
+
+  // TrilinosWrappers::SparseMatrix
+  {
+    const unsigned int rc=10;
+    TrilinosWrappers::SparsityPattern sparsity_pattern (rc, rc, /*n_entries_per_row =*/ 1);
+    for (unsigned int i=0; i < rc; ++i)
+      {
+        sparsity_pattern.add(i,i);
+      }
+    sparsity_pattern.compress();
+
+    TrilinosWrappers::SparseMatrix A (sparsity_pattern);
+    TrilinosWrappers::Vector b (A.domain_partitioner());
+    TrilinosWrappers::MPI::Vector c (A.domain_partitioner());
+    for (unsigned int i=0; i < rc; ++i)
+      {
+        A.set(i,i,2.0);
+        b(i) = i;
+      }
+
+    // === PRECONDITIONERS ===
+    deallog << "PRECONDITIONERS" << std::endl;
+    deallog.push("Preconditioners");
+
+    {
+      deallog.push("PreconditionAMG");
+      typedef TrilinosWrappers::PreconditionAMG PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionAMGMueLu");
+      typedef TrilinosWrappers::PreconditionAMGMueLu PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionChebyshev");
+      typedef TrilinosWrappers::PreconditionChebyshev PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionIC");
+      typedef TrilinosWrappers::PreconditionIC PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionIdentity");
+      typedef TrilinosWrappers::PreconditionIdentity PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionILU");
+      typedef TrilinosWrappers::PreconditionILU PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionILUT");
+      typedef TrilinosWrappers::PreconditionILUT PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionJacobi");
+      typedef TrilinosWrappers::PreconditionJacobi PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionSOR");
+      typedef TrilinosWrappers::PreconditionSOR PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("PreconditionSSOR");
+      typedef TrilinosWrappers::PreconditionSSOR PREC;
+      test_preconditioner<PREC>(A, b);
+      test_preconditioner<PREC>(A, c);
+      deallog.pop();
+    }
+
+    deallog.pop();
+
+    // === SOLVERS ===
+    deallog << std::endl;
+    deallog << "SOLVERS" << std::endl;
+    deallog.push("Solvers");
+
+    {
+      deallog.push("SolverBicgstab");
+      typedef SolverBicgstab<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("SolverCG");
+      typedef SolverCG<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("SolverGMRES");
+      typedef SolverGMRES<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("SolverFGMRES");
+      typedef SolverFGMRES<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("SolverMinRes");
+      typedef SolverMinRes<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("SolverQMRS");
+      typedef SolverQMRS<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    {
+      deallog.push("SolverRichardson");
+      typedef SolverRichardson<TrilinosWrappers::Vector> SLVR;
+      test_solver<SLVR> (A, b);
+      deallog.pop();
+    }
+
+    deallog.pop();
+    deallog << "TrilinosWrappers::SparseMatrix OK" << std::endl;
+  } // TrilinosWrappers::SparseMatrix
+
+}
diff --git a/tests/lac/linear_operator_10a.with_cxx11=on.with_trilinos=on.output b/tests/lac/linear_operator_10a.with_cxx11=on.with_trilinos=on.output
new file mode 100644 (file)
index 0000000..0e3e7af
--- /dev/null
@@ -0,0 +1,146 @@
+
+DEAL:0::PRECONDITIONERS
+DEAL:0:Preconditioners:PreconditionAMG:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMG:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMG:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMG:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMG:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionAMGMueLu:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:Exemp:S_Op::0.000000000 0.08823529412 0.1764705882 0.2647058824 0.3529411765 0.4411764706 0.5294117647 0.6176470588 0.7058823529 0.7941176471 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:Exemp:C_Op::0.000000000 0.01557093426 0.03114186851 0.04671280277 0.06228373702 0.07785467128 0.09342560554 0.1089965398 0.1245674740 0.1401384083 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:S.A.:S_Op::0.000000000 0.08823529412 0.1764705882 0.2647058824 0.3529411765 0.4411764706 0.5294117647 0.6176470588 0.7058823529 0.7941176471 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:S.A.:C_Op::0.000000000 0.01557093426 0.03114186851 0.04671280277 0.06228373702 0.07785467128 0.09342560554 0.1089965398 0.1245674740 0.1401384083 
+DEAL:0:Preconditioners:PreconditionChebyshev:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionChebyshev:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIC:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIC:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIC:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIC:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIC:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIdentity:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:Exemp:S_Op::0.000000000 1.000000000 2.000000000 3.000000000 4.000000000 5.000000000 6.000000000 7.000000000 8.000000000 9.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:Exemp:C_Op::0.000000000 2.000000000 4.000000000 6.000000000 8.000000000 10.00000000 12.00000000 14.00000000 16.00000000 18.00000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:S.A.:S_Op::0.000000000 1.000000000 2.000000000 3.000000000 4.000000000 5.000000000 6.000000000 7.000000000 8.000000000 9.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:S.A.:C_Op::0.000000000 2.000000000 4.000000000 6.000000000 8.000000000 10.00000000 12.00000000 14.00000000 16.00000000 18.00000000 
+DEAL:0:Preconditioners:PreconditionIdentity:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionIdentity:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILU:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILU:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILU:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILU:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILU:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILUT:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILUT:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionILUT:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILUT:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionILUT:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionJacobi:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionJacobi:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionJacobi:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionJacobi:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionJacobi:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSOR:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSOR:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSOR:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSOR:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSOR:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSSOR:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSSOR:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:Exemp:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:Exemp:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:S.A.:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:S.A.:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Preconditioners:PreconditionSSOR:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSSOR:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:Exemp:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:Exemp:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:S.A.:S_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0:Preconditioners:PreconditionSSOR:Approx:S.A.:C_Op::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 
+DEAL:0::
+DEAL:0::SOLVERS
+DEAL:0:Solvers:SolverBicgstab:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverBicgstab:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverBicgstab:C_Op2::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverCG:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverCG:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverCG:C_Op2::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverGMRES:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverGMRES:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverGMRES:C_Op2::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverFGMRES:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverFGMRES:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverFGMRES:C_Op2::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverMinRes:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverMinRes:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverMinRes:C_Op2::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverQMRS:S_Op::0.000000000 1.000000000 2.000000000 3.000000000 4.000000000 5.000000000 6.000000000 7.000000000 8.000000000 9.000000000 
+DEAL:0:Solvers:SolverQMRS:C_Op::0.000000000 2.000000000 4.000000000 6.000000000 8.000000000 10.00000000 12.00000000 14.00000000 16.00000000 18.00000000 
+DEAL:0:Solvers:SolverQMRS:C_Op2::0.000000000 2.000000000 4.000000000 6.000000000 8.000000000 10.00000000 12.00000000 14.00000000 16.00000000 18.00000000 
+DEAL:0:Solvers:SolverRichardson:S_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverRichardson:C_Op::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0:Solvers:SolverRichardson:C_Op2::0.000000000 0.5000000000 1.000000000 1.500000000 2.000000000 2.500000000 3.000000000 3.500000000 4.000000000 4.500000000 
+DEAL:0::TrilinosWrappers::SparseMatrix OK
diff --git a/tests/lac/linear_operator_12a.cc b/tests/lac/linear_operator_12a.cc
new file mode 100644 (file)
index 0000000..756900d
--- /dev/null
@@ -0,0 +1,389 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// tests Trilinos direct solvers on a 2D Poisson equation for linear elements
+// Note: This test is a modified version of tests/trilinos/direct_solver_2.cc
+
+#include "../tests.h"
+
+#include <deal.II/lac/trilinos_linear_operator.h>
+#include <deal.II/lac/packaged_operation.h>
+
+#include <deal.II/base/function.h>
+#include <deal.II/base/utilities.h>
+#include <deal.II/base/index_set.h>
+#include <deal.II/base/conditional_ostream.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_values.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/lac/sparsity_tools.h>
+
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/lac/trilinos_solver.h>
+#include <deal.II/lac/trilinos_precondition.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_vector.h>
+
+#include <deal.II/numerics/vector_tools.h>
+
+#include <fstream>
+#include <iomanip>
+
+using namespace dealii;
+
+template <int dim>
+class Step4
+{
+public:
+  Step4 ();
+  void run ();
+
+private:
+  void make_grid ();
+  void setup_system();
+  void assemble_system ();
+  void solve ();
+
+  parallel::distributed::Triangulation<dim>   triangulation;
+  FE_Q<dim>                                   fe;
+  DoFHandler<dim>                             dof_handler;
+
+  ConstraintMatrix                            constraints;
+  SparsityPattern                             sparsity_pattern;
+
+  TrilinosWrappers::SparseMatrix              system_matrix;
+
+  TrilinosWrappers::MPI::Vector               solution;
+  TrilinosWrappers::MPI::Vector               system_rhs;
+};
+
+
+template <int dim>
+class RightHandSide : public Function<dim>
+{
+public:
+  RightHandSide () : Function<dim>() {}
+
+  virtual double value (const Point<dim>   &p,
+                        const unsigned int  component = 0) const;
+};
+
+
+template <int dim>
+class RightHandSideTwo : public Function<dim>
+{
+public:
+  RightHandSideTwo () : Function<dim>() {}
+
+  virtual double value (const Point<dim>   &p,
+                        const unsigned int  component = 0) const;
+};
+
+
+
+
+template <int dim>
+class BoundaryValues : public Function<dim>
+{
+public:
+  BoundaryValues () : Function<dim>() {}
+
+  virtual double value (const Point<dim>   &p,
+                        const unsigned int  component = 0) const;
+};
+
+
+
+
+template <int dim>
+double RightHandSide<dim>::value (const Point<dim> &p,
+                                  const unsigned int /*component*/) const
+{
+  double return_value = 0;
+  for (unsigned int i=0; i<dim; ++i)
+    return_value += 2*std::pow(p(i), 2);
+
+  return return_value;
+}
+
+
+template <int dim>
+double RightHandSideTwo<dim>::value (const Point<dim> &p,
+                                     const unsigned int /*component*/) const
+{
+  double return_value = 0;
+  for (unsigned int i=0; i<dim; ++i)
+    return_value += 4*std::pow(p(i), 4);
+
+  return return_value;
+}
+
+
+template <int dim>
+double BoundaryValues<dim>::value (const Point<dim> &p,
+                                   const unsigned int /*component*/) const
+{
+  return p.square();
+}
+
+
+
+template <int dim>
+Step4<dim>::Step4 ()
+  :
+  triangulation(MPI_COMM_WORLD,
+                typename Triangulation<dim>::MeshSmoothing
+                (Triangulation<dim>::smoothing_on_refinement |
+                 Triangulation<dim>::smoothing_on_coarsening)),
+  fe (1),
+  dof_handler (triangulation)
+{}
+
+
+template <int dim>
+void Step4<dim>::make_grid ()
+{
+  GridGenerator::hyper_cube (triangulation, -1, 1);
+  triangulation.refine_global (6);
+}
+
+
+
+template <int dim>
+void Step4<dim>::setup_system ()
+{
+  dof_handler.distribute_dofs (fe);
+
+  constraints.clear();
+  std::map<unsigned int,double> boundary_values;
+  VectorTools::interpolate_boundary_values (dof_handler,
+                                            0,
+                                            BoundaryValues<dim>(),
+                                            constraints);
+  constraints.close();
+
+  IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs();
+  IndexSet locally_relevant_dofs;
+
+  DoFTools::extract_locally_relevant_dofs(dof_handler,
+                                          locally_relevant_dofs);
+
+
+  DynamicSparsityPattern dsp(dof_handler.n_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, dsp, constraints, false);
+  SparsityTools::distribute_sparsity_pattern(dsp,
+                                             dof_handler.n_locally_owned_dofs_per_processor(),
+                                             MPI_COMM_WORLD,
+                                             locally_relevant_dofs);
+
+  system_matrix.reinit (locally_owned_dofs,
+                        locally_owned_dofs,
+                        dsp,
+                        MPI_COMM_WORLD);
+
+  solution.reinit (locally_relevant_dofs,
+                   MPI_COMM_WORLD);
+
+  system_rhs.reinit (locally_owned_dofs,
+                     locally_relevant_dofs,
+                     MPI_COMM_WORLD,
+                     true);
+
+}
+
+
+template <int dim>
+void Step4<dim>::assemble_system ()
+{
+  QGauss<dim>  quadrature_formula(fe.degree+1);
+
+  const RightHandSide<dim> right_hand_side;
+
+  FEValues<dim> fe_values (fe, quadrature_formula,
+                           update_values   | update_gradients |
+                           update_quadrature_points | update_JxW_values);
+
+  const unsigned int   dofs_per_cell = fe.dofs_per_cell;
+  const unsigned int   n_q_points    = quadrature_formula.size();
+
+  FullMatrix<double>   cell_matrix (dofs_per_cell, dofs_per_cell);
+  Vector<double>       cell_rhs (dofs_per_cell);
+
+  std::vector<types::global_dof_index> local_dof_indices (dofs_per_cell);
+
+  typename DoFHandler<dim>::active_cell_iterator
+  cell = dof_handler.begin_active(),
+  endc = dof_handler.end();
+
+  for (; cell!=endc; ++cell)
+    {
+      if (cell->is_locally_owned())
+        {
+          fe_values.reinit (cell);
+          cell_matrix = 0;
+          cell_rhs = 0;
+
+          for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
+            for (unsigned int i=0; i<dofs_per_cell; ++i)
+              {
+                for (unsigned int j=0; j<dofs_per_cell; ++j)
+                  cell_matrix(i,j) += (fe_values.shape_grad (i, q_point) *
+                                       fe_values.shape_grad (j, q_point) *
+                                       fe_values.JxW (q_point));
+
+                cell_rhs(i) += (fe_values.shape_value (i, q_point) *
+                                right_hand_side.value (fe_values.quadrature_point (q_point)) *
+                                fe_values.JxW (q_point));
+
+              }
+
+          cell->get_dof_indices (local_dof_indices);
+          constraints.distribute_local_to_global(cell_matrix, cell_rhs,
+                                                 local_dof_indices,
+                                                 system_matrix, system_rhs);
+        }
+    }
+  system_matrix.compress(VectorOperation::add);
+  system_rhs.compress(VectorOperation::add);
+}
+
+
+
+template <int dim>
+void Step4<dim>::solve ()
+{
+  typedef TrilinosWrappers::MPI::Vector VectorType;
+
+  // Compute 'reference' solution with direct solver
+  VectorType temp_solution(system_rhs);
+  {
+    deallog.push("DirectKLU");
+    temp_solution = 0;
+    TrilinosWrappers::SolverDirect::AdditionalData data;
+    data.solver_type = "Amesos_Klu";
+    SolverControl           solver_control (1000, 1e-10);
+    TrilinosWrappers::SolverDirect solver(solver_control, data);
+    solver.solve (system_matrix, temp_solution, system_rhs);
+    constraints.distribute(temp_solution);
+    solution = temp_solution;
+    deallog.pop();
+  }
+
+  VectorType output(system_rhs);
+  {
+    deallog.push("deal_II_CG_SSOR");
+    output = 0;
+    SolverControl                           solver_control (1000, 1e-12);
+    SolverCG<TrilinosWrappers::MPI::Vector> solver (solver_control);
+    TrilinosWrappers::PreconditionSSOR preconditioner;
+    preconditioner.initialize(system_matrix);
+    solver.solve (system_matrix, output, system_rhs,
+                  preconditioner);
+    constraints.distribute(output);
+    output -= temp_solution;
+    const double local_error = output.l2_norm();
+    const double global_error = std::sqrt(Utilities::MPI::sum(
+                                            local_error * local_error,
+                                            MPI_COMM_WORLD));
+    deallog << "Norm of error in standard solve: " << global_error
+            << std::endl;
+    deallog.pop();
+  }
+
+  {
+    deallog.push("LinearOperator_deal_II_CG_SSOR");
+    output = 0;
+    SolverControl                           solver_control (1000, 1e-12);
+    SolverCG<TrilinosWrappers::MPI::Vector> solver (solver_control);
+    TrilinosWrappers::PreconditionSSOR preconditioner;
+    preconditioner.initialize(system_matrix);
+    const auto lo_A = linear_operator<VectorType>(system_matrix);
+    const auto lo_A_inv = inverse_operator(lo_A, solver, preconditioner);
+    output = lo_A_inv*system_rhs;
+    constraints.distribute(output);
+    output -= temp_solution;
+    const double local_error = output.l2_norm();
+    const double global_error = std::sqrt(Utilities::MPI::sum(
+                                            local_error * local_error,
+                                            MPI_COMM_WORLD));
+    deallog << "Norm of error in LinearOperator solve: " << global_error
+            << std::endl;
+    deallog.pop();
+  }
+
+}
+
+
+
+template <int dim>
+void Step4<dim>::run()
+{
+  make_grid();
+  setup_system();
+  assemble_system();
+  solve();
+}
+
+
+int main (int argc, char **argv)
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.threshold_double(1.e-10);
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  try
+    {
+      Step4<2> test;
+      test.run();
+    }
+  catch (std::exception &exc)
+    {
+      deallog << std::endl << std::endl
+              << "----------------------------------------------------"
+              << std::endl;
+      deallog << "Exception on processing: " << std::endl
+              << exc.what() << std::endl
+              << "Aborting!" << std::endl
+              << "----------------------------------------------------"
+              << std::endl;
+
+      return 1;
+    }
+  catch (...)
+    {
+      deallog << std::endl << std::endl
+              << "----------------------------------------------------"
+              << std::endl;
+      deallog << "Unknown exception!" << std::endl
+              << "Aborting!" << std::endl
+              << "----------------------------------------------------"
+              << std::endl;
+      return 1;
+    };
+}
diff --git a/tests/lac/linear_operator_12a.with_cxx11=on.with_trilinos=on.with_mpi=true.mpirun=1.output b/tests/lac/linear_operator_12a.with_cxx11=on.with_trilinos=on.with_mpi=true.mpirun=1.output
new file mode 100644 (file)
index 0000000..a53349a
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL:DirectKLU::Starting value 0
+DEAL:DirectKLU::Convergence step 0 value 0
+DEAL:deal_II_CG_SSOR:cg::Starting value 21.8027
+DEAL:deal_II_CG_SSOR:cg::Convergence step 86 value 0
+DEAL:deal_II_CG_SSOR::Norm of error in standard solve: 0
+DEAL:LinearOperator_deal_II_CG_SSOR:cg::Starting value 21.8027
+DEAL:LinearOperator_deal_II_CG_SSOR:cg::Convergence step 86 value 0
+DEAL:LinearOperator_deal_II_CG_SSOR::Norm of error in LinearOperator solve: 0
diff --git a/tests/lac/linear_operator_12a.with_cxx11=on.with_trilinos=on.with_mpi=true.mpirun=2.output b/tests/lac/linear_operator_12a.with_cxx11=on.with_trilinos=on.with_mpi=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..9ed8f08
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL:DirectKLU::Starting value 0
+DEAL:DirectKLU::Convergence step 0 value 0
+DEAL:deal_II_CG_SSOR:cg::Starting value 21.8027
+DEAL:deal_II_CG_SSOR:cg::Convergence step 105 value 0
+DEAL:deal_II_CG_SSOR::Norm of error in standard solve: 0
+DEAL:LinearOperator_deal_II_CG_SSOR:cg::Starting value 21.8027
+DEAL:LinearOperator_deal_II_CG_SSOR:cg::Convergence step 105 value 0
+DEAL:LinearOperator_deal_II_CG_SSOR::Norm of error in LinearOperator solve: 0

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.