Add tests for the HSL MAXX solvers.

git-svn-id: https://svn.dealii.org/trunk@9552 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/tests/bits/Makefile b/tests/bits/Makefile
index b22a1499769008ad2b1b467f92e18660280a4291..37f3067f49606e93ae8234d174f464174675230c 100644 (file)
--- a/tests/bits/Makefile
+++ b/tests/bits/Makefile
@@ -70,10 +70,17 @@ tests_x = anna_? \
 ifeq ($(USE_CONTRIB_PETSC),yes)
   tests_x += petsc_*
 endif
+
 ifeq ($(USE_CONTRIB_METIS),yes)
   tests_x += metis_*
 endif
 
+ifeq ($(USE_CONTRIB_HSL),yes)
+  tests_x += hsl_*
+endif
+
+
+
 # from above list of regular expressions, generate the real set of
 # tests
 expand = $(shell echo $(addsuffix .cc,$(1)) \

diff --git a/tests/bits/hsl_ma27_01.cc b/tests/bits/hsl_ma27_01.cc
new file mode 100644 (file)
index 0000000..3db27e9
--- /dev/null
+++ b/tests/bits/hsl_ma27_01.cc
@@ -0,0 +1,111 @@
+//----------------------------  hsl_ma27_01.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  hsl_ma27_01.cc  ---------------------------
+
+// test the hsl_ma27 sparse direct solver on a mass matrix
+
+#include "../tests.h"
+#include <iostream>
+#include <fstream>
+#include <cstdlib>
+
+#include <base/quadrature_lib.h>
+#include <base/function.h>
+
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+#include <dofs/dof_tools.h>
+
+#include <lac/sparse_matrix.h>
+#include <lac/sparsity_pattern.h>
+#include <lac/vector.h>
+#include <lac/sparse_direct.h>
+
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+
+#include <numerics/vectors.h>
+#include <numerics/matrices.h>
+
+
+
+template <int dim>
+void test ()
+{
+  deallog << dim << 'd' << std::endl;
+  
+  Triangulation<dim> tria;  
+  GridGenerator::hyper_cube (tria,0,1);
+  tria.refine_global (1);
+
+                                   // destroy the uniformity of the matrix by
+                                   // refining one cell
+  tria.begin_active()->set_refine_flag ();
+  tria.execute_coarsening_and_refinement ();
+
+  FE_Q<dim> fe (1);
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+    
+  SparsityPattern sparsity_pattern;
+  sparsity_pattern.reinit (dof_handler.n_dofs(),
+                           dof_handler.n_dofs(),
+                           dof_handler.max_couplings_between_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern);
+  sparsity_pattern.compress ();
+  
+  SparseMatrix<double> B;
+  B.reinit (sparsity_pattern);  
+  
+  QGauss<dim> qr (2);
+  MatrixTools::create_mass_matrix (dof_handler, qr, B);
+
+                                   // for a number of different solution
+                                   // vectors, make up a matching rhs vector
+                                   // and check what the MA27 solver finds
+  for (unsigned int i=0; i<3; ++i)
+    {
+      Vector<double> solution (dof_handler.n_dofs());
+      Vector<double> x (dof_handler.n_dofs());
+      Vector<double> b (dof_handler.n_dofs());
+
+      for (unsigned int j=0; j<dof_handler.n_dofs(); ++j)
+        solution(j) = j+j*(i+1)*(i+1);
+
+      B.vmult (b, solution);
+      x = b;
+      SparseDirectMA27().solve (B, x);
+
+      x -= solution;
+      deallog << "relative norm distance = "
+              << x.l2_norm() / solution.l2_norm()
+              << std::endl;
+      deallog << "absolute norms = "
+              << x.l2_norm() << ' ' << solution.l2_norm()
+              << std::endl;
+      Assert (x.l2_norm() / solution.l2_norm() < 1e-8,
+              ExcInternalError());
+    }
+}
+
+
+int main ()
+{
+  std::ofstream logfile("hsl_ma27_01.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+}

diff --git a/tests/bits/hsl_ma27_02.cc b/tests/bits/hsl_ma27_02.cc
new file mode 100644 (file)
index 0000000..41f85a9
--- /dev/null
+++ b/tests/bits/hsl_ma27_02.cc
@@ -0,0 +1,117 @@
+//----------------------------  hsl_ma27_02.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  hsl_ma27_02.cc  ---------------------------
+
+// test the hsl_ma27 sparse direct solver on a mass matrix; test that the
+// functions that allow for repeated solves yield the correct result even
+// though the matrix is decomposed only once
+
+#include "../tests.h"
+#include <iostream>
+#include <fstream>
+#include <cstdlib>
+
+#include <base/quadrature_lib.h>
+#include <base/function.h>
+
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+#include <dofs/dof_tools.h>
+
+#include <lac/sparse_matrix.h>
+#include <lac/sparsity_pattern.h>
+#include <lac/vector.h>
+#include <lac/sparse_direct.h>
+
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+
+#include <numerics/vectors.h>
+#include <numerics/matrices.h>
+
+
+
+template <int dim>
+void test ()
+{
+  deallog << dim << 'd' << std::endl;
+  
+  Triangulation<dim> tria;  
+  GridGenerator::hyper_cube (tria,0,1);
+  tria.refine_global (1);
+
+                                   // destroy the uniformity of the matrix by
+                                   // refining one cell
+  tria.begin_active()->set_refine_flag ();
+  tria.execute_coarsening_and_refinement ();
+
+  FE_Q<dim> fe (1);
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+    
+  SparsityPattern sparsity_pattern;
+  sparsity_pattern.reinit (dof_handler.n_dofs(),
+                           dof_handler.n_dofs(),
+                           dof_handler.max_couplings_between_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern);
+  sparsity_pattern.compress ();
+  
+  SparseMatrix<double> B;
+  B.reinit (sparsity_pattern);  
+  
+  QGauss<dim> qr (2);
+  MatrixTools::create_mass_matrix (dof_handler, qr, B);
+
+                                   // compute a decomposition of the matrix
+  SparseDirectMA27 Binv;
+  Binv.factorize (B);
+      
+                                   // for a number of different solution
+                                   // vectors, make up a matching rhs vector
+                                   // and check what the MA27 solver finds
+  for (unsigned int i=0; i<3; ++i)
+    {
+      Vector<double> solution (dof_handler.n_dofs());
+      Vector<double> x (dof_handler.n_dofs());
+      Vector<double> b (dof_handler.n_dofs());
+
+      for (unsigned int j=0; j<dof_handler.n_dofs(); ++j)
+        solution(j) = j+j*(i+1)*(i+1);
+
+      B.vmult (b, solution);
+      x = b;
+      Binv.solve (x);
+
+      x -= solution;
+      deallog << "relative norm distance = "
+              << x.l2_norm() / solution.l2_norm()
+              << std::endl;
+      deallog << "absolute norms = "
+              << x.l2_norm() << ' ' << solution.l2_norm()
+              << std::endl;
+      Assert (x.l2_norm() / solution.l2_norm() < 1e-8,
+              ExcInternalError());
+    }
+}
+
+
+int main ()
+{
+  std::ofstream logfile("hsl_ma27_02.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+}

diff --git a/tests/bits/hsl_ma47_01.cc b/tests/bits/hsl_ma47_01.cc
new file mode 100644 (file)
index 0000000..6d84927
--- /dev/null
+++ b/tests/bits/hsl_ma47_01.cc
@@ -0,0 +1,111 @@
+//----------------------------  hsl_ma47_01.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  hsl_ma47_01.cc  ---------------------------
+
+// test the hsl_ma47 sparse direct solver on a mass matrix
+
+#include "../tests.h"
+#include <iostream>
+#include <fstream>
+#include <cstdlib>
+
+#include <base/quadrature_lib.h>
+#include <base/function.h>
+
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+#include <dofs/dof_tools.h>
+
+#include <lac/sparse_matrix.h>
+#include <lac/sparsity_pattern.h>
+#include <lac/vector.h>
+#include <lac/sparse_direct.h>
+
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+
+#include <numerics/vectors.h>
+#include <numerics/matrices.h>
+
+
+
+template <int dim>
+void test ()
+{
+  deallog << dim << 'd' << std::endl;
+  
+  Triangulation<dim> tria;  
+  GridGenerator::hyper_cube (tria,0,1);
+  tria.refine_global (1);
+
+                                   // destroy the uniformity of the matrix by
+                                   // refining one cell
+  tria.begin_active()->set_refine_flag ();
+  tria.execute_coarsening_and_refinement ();
+
+  FE_Q<dim> fe (1);
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+    
+  SparsityPattern sparsity_pattern;
+  sparsity_pattern.reinit (dof_handler.n_dofs(),
+                           dof_handler.n_dofs(),
+                           dof_handler.max_couplings_between_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern);
+  sparsity_pattern.compress ();
+  
+  SparseMatrix<double> B;
+  B.reinit (sparsity_pattern);  
+  
+  QGauss<dim> qr (2);
+  MatrixTools::create_mass_matrix (dof_handler, qr, B);
+
+                                   // for a number of different solution
+                                   // vectors, make up a matching rhs vector
+                                   // and check what the MA47 solver finds
+  for (unsigned int i=0; i<3; ++i)
+    {
+      Vector<double> solution (dof_handler.n_dofs());
+      Vector<double> x (dof_handler.n_dofs());
+      Vector<double> b (dof_handler.n_dofs());
+
+      for (unsigned int j=0; j<dof_handler.n_dofs(); ++j)
+        solution(j) = j+j*(i+1)*(i+1);
+
+      B.vmult (b, solution);
+      x = b;
+      SparseDirectMA47().solve (B, x);
+
+      x -= solution;
+      deallog << "relative norm distance = "
+              << x.l2_norm() / solution.l2_norm()
+              << std::endl;
+      deallog << "absolute norms = "
+              << x.l2_norm() << ' ' << solution.l2_norm()
+              << std::endl;
+      Assert (x.l2_norm() / solution.l2_norm() < 1e-8,
+              ExcInternalError());
+    }
+}
+
+
+int main ()
+{
+  std::ofstream logfile("hsl_ma47_01.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+}

diff --git a/tests/bits/hsl_ma47_02.cc b/tests/bits/hsl_ma47_02.cc
new file mode 100644 (file)
index 0000000..6e005c9
--- /dev/null
+++ b/tests/bits/hsl_ma47_02.cc
@@ -0,0 +1,117 @@
+//----------------------------  hsl_ma47_02.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  hsl_ma47_02.cc  ---------------------------
+
+// test the hsl_ma47 sparse direct solver on a mass matrix; test that the
+// functions that allow for repeated solves yield the correct result even
+// though the matrix is decomposed only once
+
+#include "../tests.h"
+#include <iostream>
+#include <fstream>
+#include <cstdlib>
+
+#include <base/quadrature_lib.h>
+#include <base/function.h>
+
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+
+#include <dofs/dof_tools.h>
+
+#include <lac/sparse_matrix.h>
+#include <lac/sparsity_pattern.h>
+#include <lac/vector.h>
+#include <lac/sparse_direct.h>
+
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+
+#include <numerics/vectors.h>
+#include <numerics/matrices.h>
+
+
+
+template <int dim>
+void test ()
+{
+  deallog << dim << 'd' << std::endl;
+  
+  Triangulation<dim> tria;  
+  GridGenerator::hyper_cube (tria,0,1);
+  tria.refine_global (1);
+
+                                   // destroy the uniformity of the matrix by
+                                   // refining one cell
+  tria.begin_active()->set_refine_flag ();
+  tria.execute_coarsening_and_refinement ();
+
+  FE_Q<dim> fe (1);
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+    
+  SparsityPattern sparsity_pattern;
+  sparsity_pattern.reinit (dof_handler.n_dofs(),
+                           dof_handler.n_dofs(),
+                           dof_handler.max_couplings_between_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern);
+  sparsity_pattern.compress ();
+  
+  SparseMatrix<double> B;
+  B.reinit (sparsity_pattern);  
+  
+  QGauss<dim> qr (2);
+  MatrixTools::create_mass_matrix (dof_handler, qr, B);
+
+                                   // compute a decomposition of the matrix
+  SparseDirectMA47 Binv;
+  Binv.factorize (B);
+      
+                                   // for a number of different solution
+                                   // vectors, make up a matching rhs vector
+                                   // and check what the MA47 solver finds
+  for (unsigned int i=0; i<3; ++i)
+    {
+      Vector<double> solution (dof_handler.n_dofs());
+      Vector<double> x (dof_handler.n_dofs());
+      Vector<double> b (dof_handler.n_dofs());
+
+      for (unsigned int j=0; j<dof_handler.n_dofs(); ++j)
+        solution(j) = j+j*(i+1)*(i+1);
+
+      B.vmult (b, solution);
+      x = b;
+      Binv.solve (x);
+
+      x -= solution;
+      deallog << "relative norm distance = "
+              << x.l2_norm() / solution.l2_norm()
+              << std::endl;
+      deallog << "absolute norms = "
+              << x.l2_norm() << ' ' << solution.l2_norm()
+              << std::endl;
+      Assert (x.l2_norm() / solution.l2_norm() < 1e-8,
+              ExcInternalError());
+    }
+}
+
+
+int main ()
+{
+  std::ofstream logfile("hsl_ma47_02.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  test<1> ();
+  test<2> ();
+  test<3> ();
+}

diff --git a/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma27_01.output b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma27_01.output
new file mode 100644 (file)
index 0000000..655f55b
--- /dev/null
+++ b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma27_01.output
@@ -0,0 +1,22 @@
+
+DEAL::1d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 7.48331
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 18.7083
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 37.4166
+DEAL::2d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 57.2364
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 143.091
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 286.182
+DEAL::3d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 354.373
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 885.932
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 1771.86

diff --git a/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma27_02.output b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma27_02.output
new file mode 100644 (file)
index 0000000..655f55b
--- /dev/null
+++ b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma27_02.output
@@ -0,0 +1,22 @@
+
+DEAL::1d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 7.48331
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 18.7083
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 37.4166
+DEAL::2d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 57.2364
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 143.091
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 286.182
+DEAL::3d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 354.373
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 885.932
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 1771.86

diff --git a/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma47_01.output b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma47_01.output
new file mode 100644 (file)
index 0000000..655f55b
--- /dev/null
+++ b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma47_01.output
@@ -0,0 +1,22 @@
+
+DEAL::1d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 7.48331
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 18.7083
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 37.4166
+DEAL::2d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 57.2364
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 143.091
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 286.182
+DEAL::3d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 354.373
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 885.932
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 1771.86

diff --git a/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma47_02.output b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma47_02.output
new file mode 100644 (file)
index 0000000..655f55b
--- /dev/null
+++ b/tests/results/i686-pc-linux-gnu+gcc3.2/bits/hsl_ma47_02.output
@@ -0,0 +1,22 @@
+
+DEAL::1d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 7.48331
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 18.7083
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 37.4166
+DEAL::2d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 57.2364
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 143.091
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 286.182
+DEAL::3d
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 354.373
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 885.932
+DEAL::relative norm distance = 0.00000
+DEAL::absolute norms = 0.00000 1771.86