]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add test
authorTimo Heister <timo.heister@gmail.com>
Sun, 16 Jun 2013 17:14:22 +0000 (17:14 +0000)
committerTimo Heister <timo.heister@gmail.com>
Sun, 16 Jun 2013 17:14:22 +0000 (17:14 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_unify_linear_algebra@29822 0785d39b-7218-0410-832d-ea1e28bc413d

tests/gla/mat_04.cc [new file with mode: 0644]
tests/gla/mat_04/ncpu_2/cmp/generic [new file with mode: 0644]

diff --git a/tests/gla/mat_04.cc b/tests/gla/mat_04.cc
new file mode 100644 (file)
index 0000000..c7c20e2
--- /dev/null
@@ -0,0 +1,170 @@
+//---------------------------------------------------------------------------
+//    $Id: ghost_01.cc 28440 2013-02-17 14:35:08Z heister $
+//    Version: $Name$ 
+//
+//    Copyright (C) 2004, 2005, 2010 by the deal.II authors
+//
+//    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.
+//
+//---------------------------------------------------------------------------
+
+
+// LA::MPI::SparseMatrix, make sure it has the same
+// entries in Trilinos and PETSc (also when using
+// a TrilinosWrappers::SparsityPattern instead)
+
+#include "../tests.h"
+#include <deal.II/lac/generic_linear_algebra.h>
+#include <deal.II/base/index_set.h>
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <vector>
+
+#include <deal.II/lac/full_matrix.h>
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/lac/sparsity_tools.h>
+
+#include "gla.h"
+
+template <class LA, int dim>
+void test ()
+{
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  unsigned int numproc = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+  
+  if (myid==0)
+    deallog << "numproc=" << numproc << std::endl;
+
+  ConstraintMatrix cm;
+  cm.close();
+
+  parallel::distributed::Triangulation<dim>
+    triangulation (MPI_COMM_WORLD,
+                     typename Triangulation<dim>::MeshSmoothing
+                     (Triangulation<dim>::smoothing_on_refinement |
+                      Triangulation<dim>::smoothing_on_coarsening));
+  const double R0      = 6371000.-2890000;
+  const double R1      = 6371000.-  35000.;
+  GridGenerator::hyper_shell (triangulation,
+                                  Point<dim>(),
+                                  R0,
+                                  R1,
+                                  (dim==3) ? 96 : 12,
+                                  true);
+
+  FE_Q<dim>                                 temperature_fe(1);
+
+  DoFHandler<dim> dof_handler(triangulation);
+  dof_handler.distribute_dofs (temperature_fe);
+
+  IndexSet owned = dof_handler.locally_owned_dofs();
+  IndexSet relevant;
+  DoFTools::extract_locally_relevant_dofs (dof_handler, relevant);
+
+  CompressedSimpleSparsityPattern sp (relevant);
+  typename LA::MPI::SparseMatrix matrix;
+  DoFTools::make_sparsity_pattern (dof_handler, sp,
+                                       cm, false,
+                                       Utilities::MPI::
+                                       this_mpi_process(MPI_COMM_WORLD));
+  SparsityTools::distribute_sparsity_pattern(sp,
+      dof_handler.n_locally_owned_dofs_per_processor(),
+      MPI_COMM_WORLD, relevant);
+  sp.compress();
+  matrix.reinit (owned, owned, sp, MPI_COMM_WORLD);
+
+
+  matrix.print(deallog.get_file_stream());
+  
+
+                                  // done
+  if (myid==0)
+    deallog << "OK" << std::endl;
+}
+
+
+
+template <int dim>
+void test_trilinos_alternative ()
+{
+  typedef LA_Trilinos LA;
+   
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+  unsigned int numproc = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD);
+  
+  if (myid==0)
+    deallog << "numproc=" << numproc << std::endl;
+
+  ConstraintMatrix cm;
+  cm.close();
+
+  parallel::distributed::Triangulation<dim>
+    triangulation (MPI_COMM_WORLD,
+                     typename Triangulation<dim>::MeshSmoothing
+                     (Triangulation<dim>::smoothing_on_refinement |
+                      Triangulation<dim>::smoothing_on_coarsening));
+  const double R0      = 6371000.-2890000;
+  const double R1      = 6371000.-  35000.;
+  GridGenerator::hyper_shell (triangulation,
+                                  Point<dim>(),
+                                  R0,
+                                  R1,
+                                  (dim==3) ? 96 : 12,
+                                  true);
+
+  FE_Q<dim>                                 temperature_fe(1);
+
+  DoFHandler<dim> dof_handler(triangulation);
+  dof_handler.distribute_dofs (temperature_fe);
+
+  IndexSet owned = dof_handler.locally_owned_dofs();
+  IndexSet relevant;
+  DoFTools::extract_locally_relevant_dofs (dof_handler, relevant);
+
+  TrilinosWrappers::SparsityPattern sp (owned, MPI_COMM_WORLD);
+  typename LA::MPI::SparseMatrix matrix;
+  DoFTools::make_sparsity_pattern (dof_handler, sp,
+                                       cm, false,
+                                       Utilities::MPI::
+                                       this_mpi_process(MPI_COMM_WORLD));
+  sp.compress();
+  matrix.reinit (sp);
+
+  matrix.print(deallog.get_file_stream());
+
+                                  // done
+  if (myid==0)
+    deallog << "OK" << std::endl;
+}
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll log(__FILE__);
+
+  {    
+    deallog.push("PETSc");
+    test<LA_PETSc,2>();
+    deallog.pop();     
+    deallog.push("Trilinos");
+    test<LA_Trilinos,2>();
+    deallog.pop();     
+    deallog.push("Trilinos_alt");
+    test_trilinos_alternative<2>();
+    deallog.pop();  }
+  
+  // compile, don't run
+  //if (myid==9999)
+    //  test<LA_Dummy>();
+  
+
+}
diff --git a/tests/gla/mat_04/ncpu_2/cmp/generic b/tests/gla/mat_04/ncpu_2/cmp/generic
new file mode 100644 (file)
index 0000000..5272187
--- /dev/null
@@ -0,0 +1,441 @@
+
+DEAL:0:PETSc::numproc=2
+(0,0) 0.00000
+(0,1) 0.00000
+(0,2) 0.00000
+(0,3) 0.00000
+(0,11) 0.00000
+(0,13) 0.00000
+(1,0) 0.00000
+(1,1) 0.00000
+(1,2) 0.00000
+(1,3) 0.00000
+(1,4) 0.00000
+(1,5) 0.00000
+(2,0) 0.00000
+(2,1) 0.00000
+(2,2) 0.00000
+(2,3) 0.00000
+(2,11) 0.00000
+(2,13) 0.00000
+(3,0) 0.00000
+(3,1) 0.00000
+(3,2) 0.00000
+(3,3) 0.00000
+(3,4) 0.00000
+(3,5) 0.00000
+(4,1) 0.00000
+(4,3) 0.00000
+(4,4) 0.00000
+(4,5) 0.00000
+(4,6) 0.00000
+(4,7) 0.00000
+(5,1) 0.00000
+(5,3) 0.00000
+(5,4) 0.00000
+(5,5) 0.00000
+(5,6) 0.00000
+(5,7) 0.00000
+(6,4) 0.00000
+(6,5) 0.00000
+(6,6) 0.00000
+(6,7) 0.00000
+(6,8) 0.00000
+(6,9) 0.00000
+(7,4) 0.00000
+(7,5) 0.00000
+(7,6) 0.00000
+(7,7) 0.00000
+(7,8) 0.00000
+(7,9) 0.00000
+(8,6) 0.00000
+(8,7) 0.00000
+(8,8) 0.00000
+(8,9) 0.00000
+(8,14) 0.00000
+(8,15) 0.00000
+(9,6) 0.00000
+(9,7) 0.00000
+(9,8) 0.00000
+(9,9) 0.00000
+(9,14) 0.00000
+(9,15) 0.00000
+(10,10) 0.00000
+(10,11) 0.00000
+(10,12) 0.00000
+(10,13) 0.00000
+(10,22) 0.00000
+(10,23) 0.00000
+(11,0) 0.00000
+(11,2) 0.00000
+(11,10) 0.00000
+(11,11) 0.00000
+(11,12) 0.00000
+(11,13) 0.00000
+(12,10) 0.00000
+(12,11) 0.00000
+(12,12) 0.00000
+(12,13) 0.00000
+(12,22) 0.00000
+(12,23) 0.00000
+(13,0) 0.00000
+(13,2) 0.00000
+(13,10) 0.00000
+(13,11) 0.00000
+(13,12) 0.00000
+(13,13) 0.00000
+DEAL:0:PETSc::OK
+DEAL:0:Trilinos::numproc=2
+(0,0) 0.00000
+(0,1) 0.00000
+(0,2) 0.00000
+(0,3) 0.00000
+(0,11) 0.00000
+(0,13) 0.00000
+(1,0) 0.00000
+(1,1) 0.00000
+(1,2) 0.00000
+(1,3) 0.00000
+(1,4) 0.00000
+(1,5) 0.00000
+(2,0) 0.00000
+(2,1) 0.00000
+(2,2) 0.00000
+(2,3) 0.00000
+(2,11) 0.00000
+(2,13) 0.00000
+(3,0) 0.00000
+(3,1) 0.00000
+(3,2) 0.00000
+(3,3) 0.00000
+(3,4) 0.00000
+(3,5) 0.00000
+(4,1) 0.00000
+(4,3) 0.00000
+(4,4) 0.00000
+(4,5) 0.00000
+(4,6) 0.00000
+(4,7) 0.00000
+(5,1) 0.00000
+(5,3) 0.00000
+(5,4) 0.00000
+(5,5) 0.00000
+(5,6) 0.00000
+(5,7) 0.00000
+(6,4) 0.00000
+(6,5) 0.00000
+(6,6) 0.00000
+(6,7) 0.00000
+(6,8) 0.00000
+(6,9) 0.00000
+(7,4) 0.00000
+(7,5) 0.00000
+(7,6) 0.00000
+(7,7) 0.00000
+(7,8) 0.00000
+(7,9) 0.00000
+(8,6) 0.00000
+(8,7) 0.00000
+(8,8) 0.00000
+(8,9) 0.00000
+(8,14) 0.00000
+(8,15) 0.00000
+(9,6) 0.00000
+(9,7) 0.00000
+(9,8) 0.00000
+(9,9) 0.00000
+(9,14) 0.00000
+(9,15) 0.00000
+(10,10) 0.00000
+(10,11) 0.00000
+(10,12) 0.00000
+(10,13) 0.00000
+(10,22) 0.00000
+(10,23) 0.00000
+(11,0) 0.00000
+(11,2) 0.00000
+(11,10) 0.00000
+(11,11) 0.00000
+(11,12) 0.00000
+(11,13) 0.00000
+(12,10) 0.00000
+(12,11) 0.00000
+(12,12) 0.00000
+(12,13) 0.00000
+(12,22) 0.00000
+(12,23) 0.00000
+(13,0) 0.00000
+(13,2) 0.00000
+(13,10) 0.00000
+(13,11) 0.00000
+(13,12) 0.00000
+(13,13) 0.00000
+DEAL:0:Trilinos::OK
+DEAL:0:Trilinos_alt::numproc=2
+(0,0) 0.00000
+(0,1) 0.00000
+(0,2) 0.00000
+(0,3) 0.00000
+(0,11) 0.00000
+(0,13) 0.00000
+(1,0) 0.00000
+(1,1) 0.00000
+(1,2) 0.00000
+(1,3) 0.00000
+(1,4) 0.00000
+(1,5) 0.00000
+(2,0) 0.00000
+(2,1) 0.00000
+(2,2) 0.00000
+(2,3) 0.00000
+(2,11) 0.00000
+(2,13) 0.00000
+(3,0) 0.00000
+(3,1) 0.00000
+(3,2) 0.00000
+(3,3) 0.00000
+(3,4) 0.00000
+(3,5) 0.00000
+(4,1) 0.00000
+(4,3) 0.00000
+(4,4) 0.00000
+(4,5) 0.00000
+(4,6) 0.00000
+(4,7) 0.00000
+(5,1) 0.00000
+(5,3) 0.00000
+(5,4) 0.00000
+(5,5) 0.00000
+(5,6) 0.00000
+(5,7) 0.00000
+(6,4) 0.00000
+(6,5) 0.00000
+(6,6) 0.00000
+(6,7) 0.00000
+(6,8) 0.00000
+(6,9) 0.00000
+(7,4) 0.00000
+(7,5) 0.00000
+(7,6) 0.00000
+(7,7) 0.00000
+(7,8) 0.00000
+(7,9) 0.00000
+(8,6) 0.00000
+(8,7) 0.00000
+(8,8) 0.00000
+(8,9) 0.00000
+(8,14) 0.00000
+(8,15) 0.00000
+(9,6) 0.00000
+(9,7) 0.00000
+(9,8) 0.00000
+(9,9) 0.00000
+(9,14) 0.00000
+(9,15) 0.00000
+(10,10) 0.00000
+(10,11) 0.00000
+(10,12) 0.00000
+(10,13) 0.00000
+(10,22) 0.00000
+(10,23) 0.00000
+(11,0) 0.00000
+(11,2) 0.00000
+(11,10) 0.00000
+(11,11) 0.00000
+(11,12) 0.00000
+(11,13) 0.00000
+(12,10) 0.00000
+(12,11) 0.00000
+(12,12) 0.00000
+(12,13) 0.00000
+(12,22) 0.00000
+(12,23) 0.00000
+(13,0) 0.00000
+(13,2) 0.00000
+(13,10) 0.00000
+(13,11) 0.00000
+(13,12) 0.00000
+(13,13) 0.00000
+DEAL:0:Trilinos_alt::OK
+
+(14,8) 0.00000
+(14,9) 0.00000
+(14,14) 0.00000
+(14,15) 0.00000
+(14,16) 0.00000
+(14,17) 0.00000
+(15,8) 0.00000
+(15,9) 0.00000
+(15,14) 0.00000
+(15,15) 0.00000
+(15,16) 0.00000
+(15,17) 0.00000
+(16,14) 0.00000
+(16,15) 0.00000
+(16,16) 0.00000
+(16,17) 0.00000
+(16,18) 0.00000
+(16,19) 0.00000
+(17,14) 0.00000
+(17,15) 0.00000
+(17,16) 0.00000
+(17,17) 0.00000
+(17,18) 0.00000
+(17,19) 0.00000
+(18,16) 0.00000
+(18,17) 0.00000
+(18,18) 0.00000
+(18,19) 0.00000
+(18,20) 0.00000
+(18,21) 0.00000
+(19,16) 0.00000
+(19,17) 0.00000
+(19,18) 0.00000
+(19,19) 0.00000
+(19,20) 0.00000
+(19,21) 0.00000
+(20,18) 0.00000
+(20,19) 0.00000
+(20,20) 0.00000
+(20,21) 0.00000
+(20,22) 0.00000
+(20,23) 0.00000
+(21,18) 0.00000
+(21,19) 0.00000
+(21,20) 0.00000
+(21,21) 0.00000
+(21,22) 0.00000
+(21,23) 0.00000
+(22,10) 0.00000
+(22,12) 0.00000
+(22,20) 0.00000
+(22,21) 0.00000
+(22,22) 0.00000
+(22,23) 0.00000
+(23,10) 0.00000
+(23,12) 0.00000
+(23,20) 0.00000
+(23,21) 0.00000
+(23,22) 0.00000
+(23,23) 0.00000
+(14,14) 0.00000
+(14,15) 0.00000
+(14,16) 0.00000
+(14,17) 0.00000
+(14,8) 0.00000
+(14,9) 0.00000
+(15,14) 0.00000
+(15,15) 0.00000
+(15,16) 0.00000
+(15,17) 0.00000
+(15,8) 0.00000
+(15,9) 0.00000
+(16,14) 0.00000
+(16,15) 0.00000
+(16,16) 0.00000
+(16,17) 0.00000
+(16,18) 0.00000
+(16,19) 0.00000
+(17,14) 0.00000
+(17,15) 0.00000
+(17,16) 0.00000
+(17,17) 0.00000
+(17,18) 0.00000
+(17,19) 0.00000
+(18,16) 0.00000
+(18,17) 0.00000
+(18,18) 0.00000
+(18,19) 0.00000
+(18,20) 0.00000
+(18,21) 0.00000
+(19,16) 0.00000
+(19,17) 0.00000
+(19,18) 0.00000
+(19,19) 0.00000
+(19,20) 0.00000
+(19,21) 0.00000
+(20,18) 0.00000
+(20,19) 0.00000
+(20,20) 0.00000
+(20,21) 0.00000
+(20,22) 0.00000
+(20,23) 0.00000
+(21,18) 0.00000
+(21,19) 0.00000
+(21,20) 0.00000
+(21,21) 0.00000
+(21,22) 0.00000
+(21,23) 0.00000
+(22,20) 0.00000
+(22,21) 0.00000
+(22,22) 0.00000
+(22,23) 0.00000
+(22,10) 0.00000
+(22,12) 0.00000
+(23,20) 0.00000
+(23,21) 0.00000
+(23,22) 0.00000
+(23,23) 0.00000
+(23,10) 0.00000
+(23,12) 0.00000
+(14,14) 0.00000
+(14,15) 0.00000
+(14,16) 0.00000
+(14,17) 0.00000
+(14,8) 0.00000
+(14,9) 0.00000
+(15,14) 0.00000
+(15,15) 0.00000
+(15,16) 0.00000
+(15,17) 0.00000
+(15,8) 0.00000
+(15,9) 0.00000
+(16,14) 0.00000
+(16,15) 0.00000
+(16,16) 0.00000
+(16,17) 0.00000
+(16,18) 0.00000
+(16,19) 0.00000
+(17,14) 0.00000
+(17,15) 0.00000
+(17,16) 0.00000
+(17,17) 0.00000
+(17,18) 0.00000
+(17,19) 0.00000
+(18,16) 0.00000
+(18,17) 0.00000
+(18,18) 0.00000
+(18,19) 0.00000
+(18,20) 0.00000
+(18,21) 0.00000
+(19,16) 0.00000
+(19,17) 0.00000
+(19,18) 0.00000
+(19,19) 0.00000
+(19,20) 0.00000
+(19,21) 0.00000
+(20,18) 0.00000
+(20,19) 0.00000
+(20,20) 0.00000
+(20,21) 0.00000
+(20,22) 0.00000
+(20,23) 0.00000
+(21,18) 0.00000
+(21,19) 0.00000
+(21,20) 0.00000
+(21,21) 0.00000
+(21,22) 0.00000
+(21,23) 0.00000
+(22,20) 0.00000
+(22,21) 0.00000
+(22,22) 0.00000
+(22,23) 0.00000
+(22,10) 0.00000
+(22,12) 0.00000
+(23,20) 0.00000
+(23,21) 0.00000
+(23,22) 0.00000
+(23,23) 0.00000
+(23,10) 0.00000
+(23,12) 0.00000
+

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.