]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
update/add tests
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 18 Apr 2013 23:18:23 +0000 (23:18 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 18 Apr 2013 23:18:23 +0000 (23:18 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_unify_linear_algebra@29338 0785d39b-7218-0410-832d-ea1e28bc413d

tests/gla/gla.h
tests/gla/mat_01.cc [new file with mode: 0644]
tests/gla/mat_01/ncpu_10/cmp/generic [new file with mode: 0644]
tests/gla/mat_01/ncpu_4/cmp/generic [new file with mode: 0644]

index f8199275e4578148139cf2696379abe5222e510b..74c9f99a9b09bb8076e023fe8a3e21cafd4225f0 100644 (file)
@@ -8,6 +8,7 @@ class LA_PETSc
     {
       public:
        typedef LinearAlgebraPETSc::MPI::Vector Vector;
+       typedef LinearAlgebraPETSc::MPI::SparseMatrix SparseMatrix;
     };
 };
 
@@ -17,7 +18,8 @@ class LA_Trilinos
     class MPI
     {
       public:
-       typedef LinearAlgebraPETSc::MPI::Vector Vector;
+       typedef LinearAlgebraTrilinos::MPI::Vector Vector;
+       typedef LinearAlgebraTrilinos::MPI::SparseMatrix SparseMatrix;
     };
 };
 
@@ -79,6 +81,26 @@ class LA_Dummy
            
            
        };
+
+       class SparseMatrix
+       {
+         public:
+           template <typename SP>
+           SparseMatrix(const IndexSet & local,
+                        const IndexSet &,
+                        SP & sp,
+                        const MPI_Comm &comm=MPI_COMM_WORLD)
+             {}
+           
+           void set(unsigned int, unsigned int, double);
+           
+           const double & operator()(unsigned int, unsigned int) const
+             {
+               static double d;
+               return d;
+             }
+           
+       };
        
        
     };
diff --git a/tests/gla/mat_01.cc b/tests/gla/mat_01.cc
new file mode 100644 (file)
index 0000000..685d1a5
--- /dev/null
@@ -0,0 +1,125 @@
+//---------------------------------------------------------------------------
+//    $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.
+//
+//---------------------------------------------------------------------------
+
+
+// create, size, and reinit of LA::MPI::SparseMatrix
+
+#include "../tests.h"
+#include <deal.II/lac/abstract_linear_algebra.h>
+#include <deal.II/base/index_set.h>
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <vector>
+
+#include "gla.h"
+
+template <class LA> 
+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;
+
+                                  // each processor owns 2 indices and all
+                                  // are ghosting Element 1 (the second)
+
+  IndexSet local_active(numproc*2);
+  local_active.add_range(myid*2,myid*2+2);
+  IndexSet local_relevant= local_active;
+  local_relevant.add_range(1,2);
+
+  CompressedSimpleSparsityPattern csp (local_relevant);
+  
+  for (unsigned int i=0;i<2*numproc;++i)
+    if (local_relevant.is_element(i))
+      csp.add(i,i);
+
+  if (myid==0)
+    csp.add(0,1);
+
+  typename LA::MPI::SparseMatrix mat;
+  mat.reinit (local_active, local_active, csp, MPI_COMM_WORLD);
+
+  Assert(mat.n()==numproc*2, ExcInternalError());
+  Assert(mat.m()==numproc*2, ExcInternalError());
+
+                                  // set local values
+  mat.set(myid*2,myid*2, myid*2.0);
+  mat.set(myid*2+1,myid*2+1, myid*2.0+1.0);
+
+  mat.compress(VectorOperation::insert);
+
+  mat.add(0,1,1.0);
+  
+  mat.compress(VectorOperation::add);
+  
+                                  // check local values
+  if (myid==0)
+    {
+      deallog << myid*2 << ": " << mat(myid*2,myid*2) << std::endl;
+      deallog << myid*2+1 << ": " << mat(myid*2+1,myid*2+1) << std::endl;
+      deallog << "0,1 : " << mat(0,1) << std::endl;
+    }  
+
+                                  // done
+  if (myid==0)
+    deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv);
+  unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD);
+
+  deallog.push(Utilities::int_to_string(myid));
+
+  if (myid == 0)
+    {
+      std::ofstream logfile(output_file_for_mpi("mat_01").c_str());
+      deallog.attach(logfile);
+      deallog << std::setprecision(4);
+      deallog.depth_console(0);
+      deallog.threshold_double(1.e-10);
+
+      {        
+       deallog.push("PETSc");
+       test<LA_PETSc>();
+       deallog.pop();  
+       deallog.push("Trilinos");
+       test<LA_Trilinos>();
+       deallog.pop();  
+      }
+      
+    }
+  else
+      {        
+       deallog.push("PETSc");
+       test<LA_PETSc>();
+       deallog.pop();  
+       deallog.push("Trilinos");
+       test<LA_Trilinos>();
+       deallog.pop();  
+      }
+
+  // compile, don't run
+  //if (myid==9999)
+    //  test<LA_Dummy>();
+  
+
+}
diff --git a/tests/gla/mat_01/ncpu_10/cmp/generic b/tests/gla/mat_01/ncpu_10/cmp/generic
new file mode 100644 (file)
index 0000000..a748da5
--- /dev/null
@@ -0,0 +1,11 @@
+
+DEAL:0:PETSc::numproc=10
+DEAL:0:PETSc::0: 0
+DEAL:0:PETSc::1: 1.000
+DEAL:0:PETSc::0,1 : 10.00
+DEAL:0:PETSc::OK
+DEAL:0:Trilinos::numproc=10
+DEAL:0:Trilinos::0: 0
+DEAL:0:Trilinos::1: 1.000
+DEAL:0:Trilinos::0,1 : 10.00
+DEAL:0:Trilinos::OK
diff --git a/tests/gla/mat_01/ncpu_4/cmp/generic b/tests/gla/mat_01/ncpu_4/cmp/generic
new file mode 100644 (file)
index 0000000..88d363c
--- /dev/null
@@ -0,0 +1,11 @@
+
+DEAL:0:PETSc::numproc=4
+DEAL:0:PETSc::0: 0
+DEAL:0:PETSc::1: 1.000
+DEAL:0:PETSc::0,1 : 4.000
+DEAL:0:PETSc::OK
+DEAL:0:Trilinos::numproc=4
+DEAL:0:Trilinos::0: 0
+DEAL:0:Trilinos::1: 1.000
+DEAL:0:Trilinos::0,1 : 4.000
+DEAL:0:Trilinos::OK

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.