]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Apply patch by Andrew Baker.
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 14 Dec 2013 16:10:22 +0000 (16:10 +0000)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 14 Dec 2013 16:10:22 +0000 (16:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@32009 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/lac/trilinos_precondition.h
deal.II/source/lac/trilinos_precondition.cc
tests/trilinos/precondition_amg_smoother.cc [new file with mode: 0644]
tests/trilinos/precondition_amg_smoother.output [new file with mode: 0644]

index 1cdd7325ca4588fd70b1108dd40c1d063887ee52..afa759f9093cf46d78d5b596e017a649152e47be 100644 (file)
@@ -59,6 +59,12 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> New: It is now possible to select between different smoothers and coarse
+  solvers in the Trilinos AMG preconditioners by a string to the smoother's name.
+  <br>
+  (Andrew Baker, 2013/12/14)
+  </li>
+
 </ol>
 
 
index 0f295a99e12df14a8c5f10127e5d6555608c442e..eb8747a34de3c12e9c91b8b7033d710a50d0b96d 100644 (file)
@@ -1293,10 +1293,8 @@ namespace TrilinosWrappers
     struct AdditionalData
     {
       /**
-       * Constructor. By default, we
-       * pretend to work on elliptic
-       * problems with linear finite
-       * elements on a scalar equation.
+       * Constructor. By default, we pretend to work on elliptic problems with
+       * linear finite elements on a scalar equation.
        */
       AdditionalData (const bool                             elliptic = true,
                       const bool                             higher_order_elements = false,
@@ -1306,153 +1304,151 @@ namespace TrilinosWrappers
                       const std::vector<std::vector<bool> > &constant_modes = std::vector<std::vector<bool> > (1),
                       const unsigned int                     smoother_sweeps = 2,
                       const unsigned int                     smoother_overlap = 0,
-                      const bool                             output_details = false);
+                      const bool                             output_details = false,
+                      const char*                            smoother_type = "Chebyshev",
+                      const char*                            coarse_type = "Amesos-KLU");
 
       /**
-       * Determines whether the AMG
-       * preconditioner should be optimized
-       * for elliptic problems (ML option
-       * smoothed aggregation SA, using a
-       * Chebyshev smoother) or for
-       * non-elliptic problems (ML option
-       * non-symmetric smoothed aggregation
-       * NSSA, smoother is SSOR with
+       * Determines whether the AMG preconditioner should be optimized for
+       * elliptic problems (ML option smoothed aggregation SA, using a
+       * Chebyshev smoother) or for non-elliptic problems (ML option
+       * non-symmetric smoothed aggregation NSSA, smoother is SSOR with
        * underrelaxation).
        */
       bool elliptic;
 
       /**
-       * Determines whether the matrix that
-       * the preconditioner is built upon
-       * is generated from linear or
-       * higher-order elements.
+       * Determines whether the matrix that the preconditioner is built upon
+       * is generated from linear or higher-order elements.
        */
       bool higher_order_elements;
 
       /**
-       * Defines how many multigrid cycles
-       * should be performed by the
+       * Defines how many multigrid cycles should be performed by the
        * preconditioner.
        */
       unsigned int n_cycles;
 
       /**
-       * Defines whether a w-cycle should be
-       * used instead of the standard setting
-       * of a v-cycle.
+       * Defines whether a w-cycle should be used instead of the standard
+       * setting of a v-cycle.
        */
       bool w_cycle;
 
       /**
-       * This threshold tells the AMG setup
-       * how the coarsening should be
-       * performed. In the AMG used by ML,
-       * all points that strongly couple
-       * with the tentative coarse-level
-       * point form one aggregate. The term
-       * <em>strong coupling</em> is
-       * controlled by the variable
-       * <tt>aggregation_threshold</tt>,
-       * meaning that all elements that are
-       * not smaller than
-       * <tt>aggregation_threshold</tt>
-       * times the diagonal element do
-       * couple strongly.
+       * This threshold tells the AMG setup how the coarsening should be
+       * performed. In the AMG used by ML, all points that strongly couple
+       * with the tentative coarse-level point form one aggregate. The term
+       * <em>strong coupling</em> is controlled by the variable
+       * <tt>aggregation_threshold</tt>, meaning that all elements that are
+       * not smaller than <tt>aggregation_threshold</tt> times the diagonal
+       * element do couple strongly.
        */
       double aggregation_threshold;
 
       /**
-       * Specifies the constant modes (near
-       * null space) of the matrix. This
-       * parameter tells AMG whether we
-       * work on a scalar equation (where
-       * the near null space only consists
-       * of ones) or on a vector-valued
+       * Specifies the constant modes (near null space) of the matrix. This
+       * parameter tells AMG whether we work on a scalar equation (where the
+       * near null space only consists of ones) or on a vector-valued
        * equation.
        */
       std::vector<std::vector<bool> > constant_modes;
 
       /**
-       * Determines how many sweeps of the
-       * smoother should be performed. When
-       * the flag <tt>elliptic</tt> is set
-       * to <tt>true</tt>, i.e., for
-       * elliptic or almost elliptic
-       * problems, the polynomial degree of
-       * the Chebyshev smoother is set to
-       * <tt>smoother_sweeps</tt>. The term
-       * sweeps refers to the number of
-       * matrix-vector products performed
-       * in the Chebyshev case. In the
-       * non-elliptic case,
-       * <tt>smoother_sweeps</tt> sets the
-       * number of SSOR relaxation sweeps
-       * for post-smoothing to be
-       * performed.
+       * Determines how many sweeps of the smoother should be performed. When
+       * the flag <tt>elliptic</tt> is set to <tt>true</tt>, i.e., for
+       * elliptic or almost elliptic problems, the polynomial degree of the
+       * Chebyshev smoother is set to <tt>smoother_sweeps</tt>. The term
+       * sweeps refers to the number of matrix-vector products performed in
+       * the Chebyshev case. In the non-elliptic case,
+       * <tt>smoother_sweeps</tt> sets the number of SSOR relaxation sweeps
+       * for post-smoothing to be performed.
        */
       unsigned int smoother_sweeps;
 
       /**
-       * Determines the overlap in the
-       * SSOR/Chebyshev error smoother when
-       * run in parallel.
+       * Determines the overlap in the SSOR/Chebyshev error smoother when run
+       * in parallel.
        */
       unsigned int smoother_overlap;
 
       /**
-       * If this flag is set to
-       * <tt>true</tt>, then internal
-       * information from the ML
-       * preconditioner is printed to
-       * screen. This can be useful when
+       * If this flag is set to <tt>true</tt>, then internal information from
+       * the ML preconditioner is printed to screen. This can be useful when
        * debugging the preconditioner.
        */
       bool output_details;
+
+      /**
+       * Determines which smoother to use for the AMG cycle. Possibilities
+       * for smoother_type are the following:
+       *     "Aztec"
+       *     "IFPACK"
+       *     "Jacobi"
+       *     "ML symmetric Gauss-Seidel"
+       *     "symmetric Gauss-Seidel"
+       *     "ML Gauss-Seidel"
+       *     "Gauss-Seidel"
+       *     "block Gauss-Seidel"
+       *     "symmetric block Gauss-Seidel"
+       *     "Chebyshev"
+       *     "MLS"
+       *     "Hiptmair"
+       *     "Amesos-KLU"
+       *     "Amesos-Superlu"
+       *     "Amesos-UMFPACK"
+       *     "Amesos-Superludist"
+       *     "Amesos-MUMPS"
+       *     "user-defined"
+       *     "SuperLU"
+       *     "IFPACK-Chebyshev"
+       *     "self"
+       *     "do-nothing"
+       *     "IC"
+       *     "ICT"
+       *     "ILU"
+       *     "ILUT"
+       *     "Block Chebyshev"
+       *     "IFPACK-Block Chebyshev"
+       */
+      const char* smoother_type;
+
+      /**
+       * Determines which solver to use on the coarsest level. The same
+       * settings as for the smoother type are possible.
+       */
+      const char* coarse_type;
     };
 
+
     /**
-     * Let Trilinos compute a multilevel
-     * hierarchy for the solution of a
-     * linear system with the given
-     * matrix. The function uses the
-     * matrix format specified in
-     * TrilinosWrappers::SparseMatrix.
+     * Let Trilinos compute a multilevel hierarchy for the solution of a
+     * linear system with the given matrix. The function uses the matrix
+     * format specified in TrilinosWrappers::SparseMatrix.
      */
     void initialize (const SparseMatrix                    &matrix,
                      const AdditionalData &additional_data = AdditionalData());
 
     /**
-     * Let Trilinos compute a multilevel
-     * hierarchy for the solution of a
-     * linear system with the given
-     * matrix. The function uses the
-     * matrix format specified in
-     * TrilinosWrappers::SparseMatrix.
+     * Let Trilinos compute a multilevel hierarchy for the solution of a
+     * linear system with the given matrix. The function uses the matrix
+     * format specified in TrilinosWrappers::SparseMatrix.
      *
-     * This function is similar to the one
-     * above, but allows the user to set
-     * all the options of the Trilinos ML
-     * preconditioner. In order to find out
-     * about all the options for ML, we
-     * refer to the <a
-     * href=http://trilinos.sandia.gov/packages/ml/mlguide5.pdf>ML
-     * user's guide</a>. In particular,
-     * users need to follow the ML
-     * instructions in case a vector-valued
-     * problem ought to be solved.
+     * This function is similar to the one above, but allows the user to set
+     * all the options of the Trilinos ML preconditioner. In order to find out
+     * about all the options for ML, we refer to the <a
+     * href=http://trilinos.sandia.gov/packages/ml/mlguide5.pdf>ML user's
+     * guide</a>. In particular, users need to follow the ML instructions in
+     * case a vector-valued problem ought to be solved.
      */
     void initialize (const SparseMatrix           &matrix,
                      const Teuchos::ParameterList &ml_parameters);
 
     /**
-     * Let Trilinos compute a multilevel
-     * hierarchy for the solution of a
-     * linear system with the given
-     * matrix. This function takes a
-     * deal.ii matrix and copies the
-     * content into a Trilinos matrix, so
-     * the function can be considered
-     * rather inefficient.
+     * Let Trilinos compute a multilevel hierarchy for the solution of a
+     * linear system with the given matrix. This function takes a deal.ii
+     * matrix and copies the content into a Trilinos matrix, so the function
+     * can be considered rather inefficient.
      */
     template <typename number>
     void initialize (const ::dealii::SparseMatrix<number> &deal_ii_sparse_matrix,
@@ -1461,47 +1457,33 @@ namespace TrilinosWrappers
                      const ::dealii::SparsityPattern      *use_this_sparsity = 0);
 
     /**
-     * This function can be used for a
-     * faster recalculation of the
-     * preconditioner construction when
-     * the matrix entries underlying the
-     * preconditioner have changed, but
-     * the matrix sparsity pattern has
-     * remained the same. What this
-     * function does is taking the
-     * already generated coarsening
-     * structure, computing the AMG
-     * prolongation and restriction
-     * according to a smoothed
-     * aggregation strategy and then
-     * building the whole multilevel
-     * hiearchy. This function can be
-     * considerably faster than the
-     * initialize function, since the
-     * coarsening pattern is usually the
-     * most difficult thing to do when
-     * setting up the AMG ML
-     * preconditioner.
+     * This function can be used for a faster recalculation of the
+     * preconditioner construction when the matrix entries underlying the
+     * preconditioner have changed, but the matrix sparsity pattern has
+     * remained the same. What this function does is taking the already
+     * generated coarsening structure, computing the AMG prolongation and
+     * restriction according to a smoothed aggregation strategy and then
+     * building the whole multilevel hiearchy. This function can be
+     * considerably faster than the initialize function, since the coarsening
+     * pattern is usually the most difficult thing to do when setting up the
+     * AMG ML preconditioner.
      */
     void reinit ();
 
     /**
-     * Destroys the preconditioner, leaving
-     * an object like just after having
+     * Destroys the preconditioner, leaving an object like just after having
      * called the constructor.
      */
     void clear ();
 
     /**
-     * Prints an estimate of the memory
-     * consumption of this class.
+     * Prints an estimate of the memory consumption of this class.
      */
     size_type memory_consumption () const;
 
   private:
     /**
-     * A copy of the deal.II matrix into
-     * Trilinos format.
+     * A copy of the deal.II matrix into Trilinos format.
      */
     std_cxx1x::shared_ptr<SparseMatrix> trilinos_matrix;
   };
index db67419c435f2843d03554be96e634ac5b5615c0..be35b106638acc08ece69e42870509160172f4fb 100644 (file)
@@ -516,7 +516,9 @@ namespace TrilinosWrappers
                   const std::vector<std::vector<bool> > &constant_modes,
                   const unsigned int                     smoother_sweeps,
                   const unsigned int                     smoother_overlap,
-                  const bool                             output_details)
+                  const bool                             output_details,
+                  const char*                            smoother_type,
+                  const char*                            coarse_type)
     :
     elliptic (elliptic),
     higher_order_elements (higher_order_elements),
@@ -526,7 +528,9 @@ namespace TrilinosWrappers
     constant_modes (constant_modes),
     smoother_sweeps (smoother_sweeps),
     smoother_overlap (smoother_overlap),
-    output_details (output_details)
+    output_details (output_details),
+    smoother_type (smoother_type),
+    coarse_type (coarse_type)
   {}
 
 
@@ -543,26 +547,17 @@ namespace TrilinosWrappers
     if (additional_data.elliptic == true)
       {
         ML_Epetra::SetDefaults("SA",parameter_list);
-        parameter_list.set("smoother: type", "Chebyshev");
-
-        // uncoupled mode can give a lot of
-        // warnings or even fail when there
-        // are too many entries per row and
-        // aggreggation gets complicated, but
-        // MIS does not work if too few
-        // elements are located on one
-        // processor. work around these
-        // warnings by choosing the different
-        // strategies in different
-        // situations: for low order, always
-        // use the standard choice
-        // uncoupled. if higher order, right
-        // now we also just use Uncoupled,
-        // but we should be aware that maybe
-        // MIS might be needed
+
+        // uncoupled mode can give a lot of warnings or even fail when there
+        // are too many entries per row and aggreggation gets complicated, but
+        // MIS does not work if too few elements are located on one
+        // processor. work around these warnings by choosing the different
+        // strategies in different situations: for low order, always use the
+        // standard choice uncoupled. if higher order, right now we also just
+        // use Uncoupled, but we should be aware that maybe MIS might be
+        // needed
         //
-        // TODO: Maybe there are any
-        // other/better options?
+        // TODO: Maybe there are any other/better options?
         if (additional_data.higher_order_elements)
           {
             //if (matrix.m()/matrix.matrix->Comm().NumProc() < 50000)
@@ -578,6 +573,9 @@ namespace TrilinosWrappers
         parameter_list.set("aggregation: block scaling", true);
       }
 
+    parameter_list.set("smoother: type", additional_data.smoother_type);
+    parameter_list.set("coarse: type", additional_data.coarse_type);
+
     parameter_list.set("smoother: sweeps",
                        static_cast<int>(additional_data.smoother_sweeps));
     parameter_list.set("cycle applications",
@@ -622,10 +620,8 @@ namespace TrilinosWrappers
                 ExcDimensionMismatch(n_rows,
                                      global_length(distributed_constant_modes)));
 
-        // Reshape null space as a
-        // contiguous vector of
-        // doubles so that Trilinos
-        // can read from it.
+        // Reshape null space as a contiguous vector of doubles so that
+        // Trilinos can read from it.
         for (size_type d=0; d<constant_modes_dimension; ++d)
           for (size_type row=0; row<my_size; ++row)
             {
@@ -641,10 +637,8 @@ namespace TrilinosWrappers
         if (my_size > 0)
           parameter_list.set("null space: vectors",
                              distributed_constant_modes.Values());
-        // We need to set a valid pointer to data even
-        // if there is no data on the current
-        // processor. Therefore, pass a dummy in that
-        // case
+        // We need to set a valid pointer to data even if there is no data on
+        // the current processor. Therefore, pass a dummy in that case
         else
           parameter_list.set("null space: vectors",
                              &dummy[0]);
diff --git a/tests/trilinos/precondition_amg_smoother.cc b/tests/trilinos/precondition_amg_smoother.cc
new file mode 100644 (file)
index 0000000..11c1754
--- /dev/null
@@ -0,0 +1,302 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 2013 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// solves a 2D Poisson equation for linear elements with AMG preconditioner
+// and various smoothers
+
+#include "../tests.h"
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/lac/compressed_sparsity_pattern.h>
+#include <deal.II/lac/solver_cg.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/lac/trilinos_precondition.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/base/function.h>
+#include <deal.II/grid/tria.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+template <int dim>
+class Step4
+{
+public:
+  Step4 ();
+  void run ();
+
+private:
+  void make_grid ();
+  void setup_system();
+  void assemble_system ();
+  void solve ();
+
+  Triangulation<dim>   triangulation;
+  FE_Q<dim>            fe;
+  DoFHandler<dim>      dof_handler;
+
+  ConstraintMatrix     constraints;
+
+  TrilinosWrappers::SparseMatrix system_matrix;
+
+  Vector<double>       solution;
+  Vector<double>       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 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 += 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 ()
+  :
+  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();
+
+  CompressedSparsityPattern c_sparsity(dof_handler.n_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, c_sparsity, constraints, false);
+  system_matrix.reinit (c_sparsity);
+
+  solution.reinit (dof_handler.n_dofs());
+  system_rhs.reinit (dof_handler.n_dofs());
+}
+
+
+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)
+    {
+      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);
+}
+
+
+
+template <int dim>
+void Step4<dim>::solve ()
+{
+
+  // variant 1: solve with AMG
+  deallog.push(Utilities::int_to_string(dof_handler.n_dofs(),5));
+  deallog.push("Chebyshev");
+  TrilinosWrappers::PreconditionAMG preconditioner;
+  TrilinosWrappers::PreconditionAMG::AdditionalData data;
+  data.coarse_type = "Amesos-KLU";
+  data.smoother_type = "Chebyshev";
+  data.output_details = true;
+  data.aggregation_threshold = 1e-3;
+  data.smoother_sweeps = 3;
+  {
+    solution = 0;
+    SolverControl           solver_control (1000, 1e-10);
+    SolverCG<>              solver (solver_control);
+    preconditioner.initialize(system_matrix, data);
+    solver.solve (system_matrix, solution, system_rhs,
+                  preconditioner);
+  }
+  deallog.pop();
+
+  deallog.push("SGS");
+  data.smoother_type = "symmetric Gauss-Seidel";
+  data.smoother_sweeps = 2;
+  {
+    solution = 0;
+    SolverControl           solver_control (1000, 1e-12);
+    SolverCG<>              solver (solver_control);
+    preconditioner.initialize(system_matrix, data);
+    solver.solve (system_matrix, solution, system_rhs,
+                  preconditioner);
+  }
+  deallog.pop();
+  deallog.pop();
+}
+
+
+
+template <int dim>
+void Step4<dim>::run()
+{
+  for (unsigned int cycle = 0; cycle < 2; ++cycle)
+    {
+      if (cycle == 0)
+        make_grid();
+      else
+        triangulation.refine_global(1);
+
+      setup_system();
+      assemble_system();
+      solve();
+    }
+}
+
+
+int main (int argc, char **argv)
+{
+  std::ofstream logfile("output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  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/trilinos/precondition_amg_smoother.output b/tests/trilinos/precondition_amg_smoother.output
new file mode 100644 (file)
index 0000000..ba3157c
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL:04225:Chebyshev:cg::Starting value 21.8299
+DEAL:04225:Chebyshev:cg::Convergence step 10 value 0
+DEAL:04225:SGS:cg::Starting value 21.8299
+DEAL:04225:SGS:cg::Convergence step 10 value 0
+DEAL:16641:Chebyshev:cg::Starting value 30.8770
+DEAL:16641:Chebyshev:cg::Convergence step 8 value 0
+DEAL:16641:SGS:cg::Starting value 30.8770
+DEAL:16641:SGS:cg::Convergence step 7 value 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.