]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add all_zero() to VectorSpaceVector.
authorBruno Turcksin <bruno.turcksin@gmail.com>
Tue, 9 May 2017 01:21:43 +0000 (21:21 -0400)
committerBruno Turcksin <bruno.turcksin@gmail.com>
Tue, 9 May 2017 20:51:08 +0000 (16:51 -0400)
include/deal.II/lac/cuda_vector.h
include/deal.II/lac/la_parallel_block_vector.h
include/deal.II/lac/la_parallel_vector.h
include/deal.II/lac/la_vector.h
include/deal.II/lac/la_vector.templates.h
include/deal.II/lac/trilinos_epetra_vector.h
include/deal.II/lac/vector_space_vector.h
source/lac/cuda_vector.cu
source/lac/trilinos_epetra_vector.cc

index c10ac4129ff06954710fc90d032871ccf6c99e56..4bf4c99711d1385ae42165040ea99a7cbc393c91 100644 (file)
@@ -169,6 +169,11 @@ namespace LinearAlgebra
        */
       virtual void equ(const Number a, const VectorSpaceVector<Number> &V) override;
 
+      /**
+       * Return wether the vector contains only elements with value zero.
+       */
+      virtual bool all_zero() const override;
+
       /**
        * Return the mean value of all the entries of this vector.
        */
index 7dbccdf7bd392765030f3b498fb9605db5713e97..a81fc411a3e79a01e9e5a2b1d43dc61be01a8e0b 100644 (file)
@@ -367,7 +367,7 @@ namespace LinearAlgebra
        * This function is mainly for internal consistency checks and should
        * seldom be used when not in debug mode since it uses quite some time.
        */
-      bool all_zero () const;
+      virtual bool all_zero () const override;
 
       /**
        * Compute the mean value of all the entries in the vector.
index 6a53c1e0b0abde39118df50ca1621b93d82c39c9..f09a287d6c087a3081bef815370723ed2e308654 100644 (file)
@@ -897,7 +897,7 @@ namespace LinearAlgebra
        * This is a collective operation. This function is expensive, because
        * potentially all elements have to be checked.
        */
-      bool all_zero () const;
+      virtual bool all_zero () const override;
 
       /**
        * Compute the mean value of all the entries in the vector.
index 380954c4344205d35cd1bc0b86c24300dd3fe8a5..bcb2d3b9c49fdbf5565655723c5a3592704b3bd3 100644 (file)
@@ -225,6 +225,11 @@ namespace LinearAlgebra
      */
     virtual void equ(const Number a, const VectorSpaceVector<Number> &V);
 
+    /**
+     * Return wether the vector contains only elements with value zero.
+     */
+    virtual bool all_zero() const override;
+
     /**
      * Return the mean value of all the entries of this vector.
      */
index c823b74777bb374d8556e4589a8c79f24965a138..f07f4e8237e52df2d12ab0c9b44a1403313a11a6 100644 (file)
@@ -327,6 +327,21 @@ namespace LinearAlgebra
 
 
 
+  template <typename Number>
+  bool Vector<Number>::all_zero() const
+  {
+    Assert(this->size(), ExcEmptyObject());
+
+    const size_type size = this->size();
+    for (size_type i=0; i<size; ++i)
+      if (this->val[i] != Number())
+        return false;
+
+    return true;
+  }
+
+
+
   template <typename Number>
   typename Vector<Number>::value_type Vector<Number>::mean_value() const
   {
index 2be2a661b0d9e18f815e45f4e55dde0c4974a42b..7c81db4d9dd899817e49c18d8316920ea1cb3f0a 100644 (file)
@@ -179,6 +179,11 @@ namespace LinearAlgebra
        */
       virtual void equ(const double a, const VectorSpaceVector<double> &V);
 
+      /**
+       * Return wether the vector contains only elements with value zero.
+       */
+      virtual bool all_zero() const override;
+
       /**
        * Return the mean value of the element of this vector.
        */
index 8a18ff9f657e379c219f6bc03b10634030939b9d..a70b4b52970926c87e1011845c0102047600a4b5 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2015 - 2016 by the deal.II authors
+// Copyright (C) 2015 - 2017 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -141,6 +141,11 @@ namespace LinearAlgebra
      */
     virtual void equ(const Number a, const VectorSpaceVector<Number> &V) = 0;
 
+    /**
+     * Return whether the vector contains only elements with value zero.
+     */
+    virtual bool all_zero() const = 0;
+
     /**
      * Return the mean value of all the entries of this vector.
      */
index 4d0921f2f0d6f6c0ec35c4cf081c5b5b1e8eca5c..42504b7fda2e508a71d78d463aa3708ee83a00c0 100644 (file)
@@ -915,6 +915,14 @@ namespace LinearAlgebra
 
 
 
+    template <typename Number>
+    bool Vector<Number>::all_zero() const
+    {
+      return (linfty_norm() == 0) ? true : false;
+    }
+
+
+
     template <typename Number>
     typename Vector<Number>::value_type Vector<Number>::mean_value() const
     {
index c44ff85cad135ae1c25c127af55a690a0eb3b68f..53476c7bf668afbb2d527fe84d9ad4eb35cb8893 100644 (file)
@@ -389,6 +389,35 @@ namespace LinearAlgebra
 
 
 
+    bool Vector::all_zero() const
+    {
+      // get a representation of the vector and
+      // loop over all the elements
+      double *start_ptr = (*vector)[0];
+      const double *ptr  = start_ptr,
+                    *eptr = start_ptr + vector->MyLength();
+      unsigned int flag = 0;
+      while (ptr != eptr)
+        {
+          if (*ptr != 0)
+            {
+              flag = 1;
+              break;
+            }
+          ++ptr;
+        }
+
+      // Check that the vector is zero on _all_ processors.
+      const Epetra_MpiComm *mpi_comm
+        = dynamic_cast<const Epetra_MpiComm *>(&vector->Map().Comm());
+      Assert(mpi_comm != nullptr, ExcInternalError());
+      unsigned int num_nonzero = Utilities::MPI::sum(flag, mpi_comm->Comm());
+
+      return num_nonzero == 0;
+    }
+
+
+
     double Vector::mean_value() const
     {
       double mean_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.