]> https://gitweb.dealii.org/ - dealii.git/commitdiff
test for PETSc Vector reinit helper
authorESeNonFossiIo <esenonfossiio@gmail.com>
Tue, 12 Apr 2016 20:11:33 +0000 (22:11 +0200)
committerESeNonFossiIo <esenonfossiio@gmail.com>
Thu, 14 Apr 2016 08:41:23 +0000 (10:41 +0200)
include/deal.II/lac/petsc_parallel_sparse_matrix.h
include/deal.II/lac/petsc_sparse_matrix.h
source/lac/petsc_parallel_sparse_matrix.cc
source/lac/petsc_sparse_matrix.cc
tests/lac/linear_operator_09.cc [new file with mode: 0644]
tests/lac/linear_operator_09.with_cxx11=on.with_petsc=on.output [new file with mode: 0644]

index fd5905527fc13ce8f97b46fd297be4fe666daf0d..6d17876da623328b7ee22e005ffe237d9e085995 100644 (file)
@@ -379,19 +379,19 @@ namespace PETScWrappers
       PetscScalar matrix_scalar_product (const Vector &u,
                                          const Vector &v) const;
 
-       /**
-        * Return the partitioning of the domain space of this matrix, i.e., the
-        * partitioning of the vectors this matrix has to be multiplied with.
-        */
-       IndexSet locally_owned_domain_indices() const;
-
-       /**
-        * Return the partitioning of the range space of this matrix, i.e., the
-        * partitioning of the vectors that are result from matrix-vector
-        * products.
-        */
-       IndexSet locally_owned_range_indices() const;
-                                         
+      /**
+       * Return the partitioning of the domain space of this matrix, i.e., the
+       * partitioning of the vectors this matrix has to be multiplied with.
+       */
+      IndexSet locally_owned_domain_indices() const;
+
+      /**
+       * Return the partitioning of the range space of this matrix, i.e., the
+       * partitioning of the vectors that are result from matrix-vector
+       * products.
+       */
+      IndexSet locally_owned_range_indices() const;
+
     private:
 
       /**
index 4c6dbc037cf62459b43b2f285f98b61031044f82..41a1481f71c10a7b538a0b1d605b3b4886eb9aa3 100644 (file)
@@ -225,15 +225,15 @@ namespace PETScWrappers
     PetscScalar matrix_scalar_product (const VectorBase &u,
                                        const VectorBase &v) const;
 
-     /**
-      * Return the number of rows of this matrix.
-      */
-     size_t m() const;
-
-     /**
-      * Return the number of coloumns of this matrix.
-      */
-     size_t n() const;
+    /**
+     * Return the number of rows of this matrix.
+     */
+    size_t m() const;
+
+    /**
+     * Return the number of coloumns of this matrix.
+     */
+    size_t n() const;
 
   private:
 
index c2837beb54270bb16d329b7eea149da4f65d63b4..043a8f5e0073d39416146a095b4476151688c808 100644 (file)
@@ -895,16 +895,16 @@ namespace PETScWrappers
     {
       PetscInt n_rows, n_cols, min, max;
       PetscErrorCode ierr;
-      ISrows = nullptr;
-      IScols = nullptr;
-      
+      IS *rows = nullptr;
+      IS *cols = nullptr;
+
       ierr = MatGetSize (matrix, &n_rows, &n_cols);
       ierr = MatGetOwnershipIS(matrix, rows, cols);
       ierr = ISGetMinMax(*rows, &min, &max);
 
       IndexSet locally_owned_domain_indices(n_rows);
       locally_owned_domain_indices.add_range(min, max);
-      
+
       return locally_owned_domain_indices;
     }
 
@@ -913,19 +913,19 @@ namespace PETScWrappers
     {
       PetscInt n_rows, n_cols, min, max;
       PetscErrorCode ierr;
-      ISrows = nullptr;
-      IScols = nullptr;
-      
+      IS *rows = nullptr;
+      IS *cols = nullptr;
+
       ierr = MatGetSize (matrix, &n_rows, &n_cols);
       ierr = MatGetOwnershipIS(matrix, rows, cols);
       ierr = ISGetMinMax(*cols, &min, &max);
-      
+
       IndexSet locally_owned_range_indices(n_cols);
       locally_owned_range_indices.add_range(min, max);
-      
+
       return locally_owned_range_indices;
     }
-    
+
   }
 }
 
index d0b42978248c0370ba7157b570783cd9d1642ca2..ec0eaa8471fc028b2a03de9204836806989927f8 100644 (file)
@@ -312,8 +312,8 @@ namespace PETScWrappers
   SparseMatrix::m () const
   {
     PetscInt m,n;
-    PetscErrorCode ierr = MatGetSize(matrix, &m, &n); 
-    
+    PetscErrorCode ierr = MatGetSize(matrix, &m, &n);
+
     return m;
   }
 
@@ -321,11 +321,11 @@ namespace PETScWrappers
   SparseMatrix::n () const
   {
     PetscInt m,n;
-    PetscErrorCode ierr = MatGetSize(matrix, &m, &n); 
-    
+    PetscErrorCode ierr = MatGetSize(matrix, &m, &n);
+
     return n;
   }
-  
+
   // Explicit instantiations
   //
   template
diff --git a/tests/lac/linear_operator_09.cc b/tests/lac/linear_operator_09.cc
new file mode 100644 (file)
index 0000000..41c48cc
--- /dev/null
@@ -0,0 +1,56 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 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 that it is possible to instantiate a LinearOperator object for all
+// different kinds of Trilinos matrices and vectors
+// TODO: A bit more tests...
+
+#include "../tests.h"
+
+
+#include <deal.II/lac/linear_operator.h>
+
+#include <deal.II/lac/petsc_vector.h>
+#include <deal.II/lac/petsc_sparse_matrix.h>
+
+#include <deal.II/lac/petsc_parallel_vector.h>
+#include <deal.II/lac/petsc_parallel_sparse_matrix.h>
+
+using namespace dealii;
+
+int main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1);
+
+  initlog();
+  deallog << std::setprecision(10);
+
+  {
+    PETScWrappers::SparseMatrix a;
+    auto op_a  = linear_operator<PETScWrappers::Vector>(a);
+  }
+
+  {
+    PETScWrappers::MPI::SparseMatrix a;
+    auto op_a  = linear_operator<PETScWrappers::MPI::Vector>(a);
+  }
+
+
+  deallog << "OK" << std::endl;
+
+  return 0;
+}
+
+
diff --git a/tests/lac/linear_operator_09.with_cxx11=on.with_petsc=on.output b/tests/lac/linear_operator_09.with_cxx11=on.with_petsc=on.output
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::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.