]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update the GMRES implementations to use VectorMemory::Pointer.
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 24 Aug 2017 17:14:52 +0000 (11:14 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 24 Aug 2017 19:11:16 +0000 (13:11 -0600)
This facilitates automatic memory allocation cleanup.

include/deal.II/lac/solver_gmres.h

index 4143379a3d1aa5206d71e766fad5a8a785a58e7f..447c7089b58e066f34865a02c79817e999baef53 100644 (file)
@@ -28,6 +28,8 @@
 #include <deal.II/lac/lapack_full_matrix.h>
 #include <deal.II/lac/vector.h>
 
+#include <deal.II/base/std_cxx14/memory.h>
+
 #include <vector>
 #include <cmath>
 #include <algorithm>
@@ -64,9 +66,9 @@ namespace internal
                  VectorMemory<VectorType> &vmem);
 
       /**
-       * Delete all allocated vectors.
+       * Destructor. Delete all allocated vectors.
        */
-      ~TmpVectors();
+      ~TmpVectors() = default;
 
       /**
        * Get vector number @p i. If this vector was unused before, an error
@@ -99,7 +101,7 @@ namespace internal
       /**
        * Field for storing the vectors.
        */
-      std::vector<VectorType *> data;
+      std::vector<typename VectorMemory<VectorType>::Pointer> data;
 
       /**
        * Offset of the first vector. This is for later when vector rotation
@@ -517,21 +519,11 @@ namespace internal
                 VectorMemory<VectorType> &vmem)
       :
       mem(vmem),
-      data (max_size, nullptr),
+      data (max_size),
       offset(0)
     {}
 
 
-    template <class VectorType>
-    inline
-    TmpVectors<VectorType>::~TmpVectors ()
-    {
-      for (typename std::vector<VectorType *>::iterator v = data.begin();
-           v != data.end(); ++v)
-        if (*v != nullptr)
-          mem.free(*v);
-    }
-
 
     template <class VectorType>
     inline VectorType &
@@ -545,6 +537,7 @@ namespace internal
     }
 
 
+
     template <class VectorType>
     inline VectorType &
     TmpVectors<VectorType>::operator() (const unsigned int i,
@@ -554,13 +547,14 @@ namespace internal
               ExcIndexRange(i,-offset, data.size()-offset));
       if (data[i-offset] == nullptr)
         {
-          data[i-offset] = mem.alloc();
+          data[i-offset] = std::move(typename VectorMemory<VectorType>::Pointer(mem));
           data[i-offset]->reinit(temp);
         }
       return *data[i-offset];
     }
 
 
+
     template <class VectorType>
     unsigned int
     TmpVectors<VectorType>::size() const
@@ -569,6 +563,7 @@ namespace internal
     }
 
 
+
     // A comparator for better printing eigenvalues
     inline
     bool complex_less_pred(const std::complex<double> &x,
@@ -822,20 +817,19 @@ SolverGMRES<VectorType>::solve (const MatrixType         &A,
   VectorType &v = tmp_vectors(0, x);
   VectorType &p = tmp_vectors(n_tmp_vectors-1, x);
 
-  // Following vectors are needed
-  // when not the default residuals
-  // are used as stopping criterion
-  VectorType *r=nullptr;
-  VectorType *x_=nullptr;
-  dealii::Vector<double> *gamma_=nullptr;
+  // Following vectors are needed when we are not using the default residuals
+  // as stopping criterion
+  typename VectorMemory<VectorType>::Pointer r;
+  typename VectorMemory<VectorType>::Pointer x_;
+  std::unique_ptr<dealii::Vector<double> > gamma_;
   if (!use_default_residual)
     {
-      r=this->memory.alloc();
-      x_=this->memory.alloc();
+      r  = std::move(typename VectorMemory<VectorType>::Pointer(this->memory));
+      x_ = std::move(typename VectorMemory<VectorType>::Pointer(this->memory));
       r->reinit(x);
       x_->reinit(x);
 
-      gamma_ = new dealii::Vector<double> (gamma.size());
+      gamma_ = std_cxx14::make_unique<dealii::Vector<double> >(gamma.size());
     }
 
   bool re_orthogonalize = additional_data.force_re_orthogonalization;
@@ -1038,14 +1032,6 @@ SolverGMRES<VectorType>::solve (const MatrixType         &A,
   if (!krylov_space_signal.empty())
     krylov_space_signal(tmp_vectors);
 
-  if (!use_default_residual)
-    {
-      this->memory.free(r);
-      this->memory.free(x_);
-
-      delete gamma_;
-    }
-
   deallog.pop();
 
   // in case of failure: throw exception
@@ -1186,7 +1172,7 @@ SolverFGMRES<VectorType>::solve (const MatrixType         &A,
   // Iteration starts here
   double res = -std::numeric_limits<double>::max();
 
-  VectorType *aux = this->memory.alloc();
+  typename VectorMemory<VectorType>::Pointer aux (this->memory);
   aux->reinit(x);
   do
     {
@@ -1247,8 +1233,6 @@ SolverFGMRES<VectorType>::solve (const MatrixType         &A,
     }
   while (iteration_state == SolverControl::iterate);
 
-  this->memory.free(aux);
-
   deallog.pop();
   // in case of failure: throw exception
   if (iteration_state != SolverControl::success)

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.