]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Put implementation (Growing)VectorMemory into .templates.h file 566/head
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 18 Feb 2015 16:35:18 +0000 (17:35 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 18 Feb 2015 21:03:48 +0000 (22:03 +0100)
This enables users to provide instantiations of the class for non-deal.II
vector types and their use in Solver classes.

Reported by Karl Ljungkvist.

doc/news/changes.h
include/deal.II/lac/solver.h
include/deal.II/lac/vector_memory.templates.h [new file with mode: 0644]
source/lac/vector_memory.cc

index c1c72b0fe95c3cd054b41cdffcaac36976bafb7c..23518e1e709fadbe0af1d438b5101ec13d434a45 100644 (file)
@@ -350,6 +350,15 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> Fixed: The implementation of the class GrowingVectorMemory has been
+  moved from source/lac/vector_memory.cc to the new file
+  include/deal.II/lac/vector_memory.templates.h. This allows users to
+  create instantiations of GrowingVectorMemory for their own vector classes
+  in case they intend to use them for the deal.II solvers.
+  <br>
+  (Martin Kronbichler, 2015/02/18)
+  </li>
+
   <li> Changed: All members of MultithreadInfo are now static so it is no
   longer necessary to use the global instance multithread_info (now
   deprecated) or create your own instance (which does not work correctly
index 49e5f51481ff1e24de9803bb99acd9157a7058fe..2388cd2aad5eb4014f4d59fb3d5f9ffa3e974e78 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 1998 - 2014 by the deal.II authors
+// Copyright (C) 1998 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -120,6 +120,24 @@ template <typename number> class Vector;
  * <tt>swap(VECTOR &a, VECTOR &b)</tt> that exchanges the values of the two
  * vectors.
  *
+ * Finally, the solvers also expect an instantiation of
+ * GrowingVectorMemory<VECTOR>. These instantiations are provided by the
+ * deal.II library for the built-in vector types, but must be explicitly added
+ * for user-provided vector classes. Otherwise, the linker will complain that
+ * it cannout find the constructors and destructors of GrowingVectorMemory
+ * that happen in the @p Solver class below.
+ *
+ * @code
+ * #include <deal.II/lac/vector_memory.templates.h>
+ *
+ * // Definition and implementation of vector class
+ * class UserVector { ... };
+ *
+ * // Create explicit instantiation for the vector class
+ * template class VectorMemory<UserVector >;
+ * template class GrowingVectorMemory<UserVector >;
+ * @endcode
+ *
  * The preconditioners used must have the same interface as matrices, i.e. in
  * particular they have to provide a member function @p vmult which denotes
  * the application of the preconditioner.
diff --git a/include/deal.II/lac/vector_memory.templates.h b/include/deal.II/lac/vector_memory.templates.h
new file mode 100644 (file)
index 0000000..2bc6391
--- /dev/null
@@ -0,0 +1,210 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2007 - 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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef __deal2__vector_memory_templates_h
+#define __deal2__vector_memory_templates_h
+
+
+#include <deal.II/lac/vector_memory.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+
+template <typename VECTOR>
+typename GrowingVectorMemory<VECTOR>::Pool GrowingVectorMemory<VECTOR>::pool;
+
+template <typename VECTOR>
+Threads::Mutex GrowingVectorMemory<VECTOR>::mutex;
+
+template <typename VECTOR>
+inline
+GrowingVectorMemory<VECTOR>::Pool::Pool()
+  :
+  data(0)
+{}
+
+
+
+template <typename VECTOR>
+inline
+GrowingVectorMemory<VECTOR>::Pool::~Pool()
+{
+  // Nothing to do if memory was unused.
+  if (data == 0) return;
+
+  // First, delete all remaining
+  // vectors. Actually, there should
+  // be none, if there is no memory
+  // leak
+  for (typename std::vector<entry_type>::iterator i=data->begin();
+       i != data->end();
+       ++i)
+    {
+      delete i->second;
+    }
+  delete data;
+}
+
+
+template <typename VECTOR>
+inline
+void
+GrowingVectorMemory<VECTOR>::Pool::initialize(const size_type size)
+{
+  if (data == 0)
+    {
+      data = new std::vector<entry_type>(size);
+
+      for (typename std::vector<entry_type>::iterator i= data->begin();
+           i != data->end();
+           ++i)
+        {
+          i->first = false;
+          i->second = new VECTOR;
+        }
+    }
+}
+
+
+template <typename VECTOR>
+inline
+GrowingVectorMemory<VECTOR>::GrowingVectorMemory (const size_type initial_size,
+                                                  const bool log_statistics)
+
+  :
+  total_alloc(0),
+  current_alloc(0),
+  log_statistics(log_statistics)
+{
+  Threads::Mutex::ScopedLock lock(mutex);
+  pool.initialize(initial_size);
+}
+
+
+template<typename VECTOR>
+inline
+GrowingVectorMemory<VECTOR>::~GrowingVectorMemory()
+{
+  AssertNothrow(current_alloc == 0,
+                StandardExceptions::ExcMemoryLeak(current_alloc));
+  if (log_statistics)
+    {
+      deallog << "GrowingVectorMemory:Overall allocated vectors: "
+              << total_alloc << std::endl;
+      deallog << "GrowingVectorMemory:Maximum allocated vectors: "
+              << pool.data->size() << std::endl;
+    }
+}
+
+
+
+template<typename VECTOR>
+inline
+VECTOR *
+GrowingVectorMemory<VECTOR>::alloc ()
+{
+  Threads::Mutex::ScopedLock lock(mutex);
+  ++total_alloc;
+  ++current_alloc;
+  // see if there is a free vector
+  // available in our list
+  for (typename std::vector<entry_type>::iterator i=pool.data->begin();
+       i != pool.data->end(); ++i)
+    {
+      if (i->first == false)
+        {
+          i->first = true;
+          return (i->second);
+        }
+    }
+
+  // no free vector found, so let's
+  // just allocate a new one
+  const entry_type t (true, new VECTOR);
+  pool.data->push_back(t);
+
+  return t.second;
+}
+
+
+
+template<typename VECTOR>
+inline
+void
+GrowingVectorMemory<VECTOR>::free(const VECTOR *const v)
+{
+  Threads::Mutex::ScopedLock lock(mutex);
+  for (typename std::vector<entry_type>::iterator i=pool.data->begin();
+       i != pool.data->end(); ++i)
+    {
+      if (v == (i->second))
+        {
+          i->first = false;
+          --current_alloc;
+          return;
+        }
+    }
+  Assert(false, typename VectorMemory<VECTOR>::ExcNotAllocatedHere());
+}
+
+
+
+template<typename VECTOR>
+inline
+void
+GrowingVectorMemory<VECTOR>::release_unused_memory ()
+{
+  Threads::Mutex::ScopedLock lock(mutex);
+
+  std::vector<entry_type> new_data;
+
+  if (pool.data != 0)
+    {
+      const typename std::vector<entry_type>::const_iterator
+      end = pool.data->end();
+      for (typename std::vector<entry_type>::const_iterator
+           i = pool.data->begin(); i != end ; ++i)
+        if (i->first == false)
+          delete i->second;
+        else
+          new_data.push_back (*i);
+
+      *pool.data = new_data;
+    }
+}
+
+
+
+template<typename VECTOR>
+inline
+std::size_t
+GrowingVectorMemory<VECTOR>::memory_consumption () const
+{
+  Threads::Mutex::ScopedLock lock(mutex);
+
+  std::size_t result = sizeof (*this);
+  const typename std::vector<entry_type>::const_iterator
+  end = pool.data->end();
+  for (typename std::vector<entry_type>::const_iterator
+       i = pool.data->begin(); i != end ; ++i)
+    result += sizeof (*i) + i->second->memory_consumption();
+
+  return result;
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 96d0687b16a9537ed0edb364ed282179186cda15..9e38ad9b3971d38add8bf8e7426067c915a2b32b 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2007 - 2014 by the deal.II authors
+// Copyright (C) 2007 - 2015 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -13,7 +13,6 @@
 //
 // ---------------------------------------------------------------------
 
-#include <deal.II/lac/vector_memory.h>
 #include <deal.II/lac/vector.h>
 #include <deal.II/lac/block_vector.h>
 #include <deal.II/lac/parallel_vector.h>
 #include <deal.II/lac/trilinos_vector.h>
 #include <deal.II/lac/trilinos_block_vector.h>
 
-DEAL_II_NAMESPACE_OPEN
-
-
-template <typename VECTOR>
-typename GrowingVectorMemory<VECTOR>::Pool GrowingVectorMemory<VECTOR>::pool;
-
-template <typename VECTOR>
-Threads::Mutex GrowingVectorMemory<VECTOR>::mutex;
-
-template <typename VECTOR>
-inline
-GrowingVectorMemory<VECTOR>::Pool::Pool()
-  :
-  data(0)
-{}
-
-
-
-template <typename VECTOR>
-inline
-GrowingVectorMemory<VECTOR>::Pool::~Pool()
-{
-  // Nothing to do if memory was unused.
-  if (data == 0) return;
-
-  // First, delete all remaining
-  // vectors. Actually, there should
-  // be none, if there is no memory
-  // leak
-  for (typename std::vector<entry_type>::iterator i=data->begin();
-       i != data->end();
-       ++i)
-    {
-      delete i->second;
-    }
-  delete data;
-}
-
-
-template <typename VECTOR>
-inline
-void
-GrowingVectorMemory<VECTOR>::Pool::initialize(const size_type size)
-{
-  if (data == 0)
-    {
-      data = new std::vector<entry_type>(size);
-
-      for (typename std::vector<entry_type>::iterator i= data->begin();
-           i != data->end();
-           ++i)
-        {
-          i->first = false;
-          i->second = new VECTOR;
-        }
-    }
-}
-
-
-template <typename VECTOR>
-inline
-GrowingVectorMemory<VECTOR>::GrowingVectorMemory (const size_type initial_size,
-                                                  const bool log_statistics)
-
-  :
-  total_alloc(0),
-  current_alloc(0),
-  log_statistics(log_statistics)
-{
-  Threads::Mutex::ScopedLock lock(mutex);
-  pool.initialize(initial_size);
-}
-
-
-template<typename VECTOR>
-inline
-GrowingVectorMemory<VECTOR>::~GrowingVectorMemory()
-{
-  AssertNothrow(current_alloc == 0,
-                StandardExceptions::ExcMemoryLeak(current_alloc));
-  if (log_statistics)
-    {
-      deallog << "GrowingVectorMemory:Overall allocated vectors: "
-              << total_alloc << std::endl;
-      deallog << "GrowingVectorMemory:Maximum allocated vectors: "
-              << pool.data->size() << std::endl;
-    }
-}
-
-
+#include <deal.II/lac/vector_memory.templates.h>
 
-template<typename VECTOR>
-inline
-VECTOR *
-GrowingVectorMemory<VECTOR>::alloc ()
-{
-  Threads::Mutex::ScopedLock lock(mutex);
-  ++total_alloc;
-  ++current_alloc;
-  // see if there is a free vector
-  // available in our list
-  for (typename std::vector<entry_type>::iterator i=pool.data->begin();
-       i != pool.data->end(); ++i)
-    {
-      if (i->first == false)
-        {
-          i->first = true;
-          return (i->second);
-        }
-    }
-
-  // no free vector found, so let's
-  // just allocate a new one
-  const entry_type t (true, new VECTOR);
-  pool.data->push_back(t);
-
-  return t.second;
-}
-
-
-
-template<typename VECTOR>
-inline
-void
-GrowingVectorMemory<VECTOR>::free(const VECTOR *const v)
-{
-  Threads::Mutex::ScopedLock lock(mutex);
-  for (typename std::vector<entry_type>::iterator i=pool.data->begin();
-       i != pool.data->end(); ++i)
-    {
-      if (v == (i->second))
-        {
-          i->first = false;
-          --current_alloc;
-          return;
-        }
-    }
-  Assert(false, typename VectorMemory<VECTOR>::ExcNotAllocatedHere());
-}
-
-
-
-template<typename VECTOR>
-inline
-void
-GrowingVectorMemory<VECTOR>::release_unused_memory ()
-{
-  Threads::Mutex::ScopedLock lock(mutex);
-
-  std::vector<entry_type> new_data;
-
-  if (pool.data != 0)
-    {
-      const typename std::vector<entry_type>::const_iterator
-      end = pool.data->end();
-      for (typename std::vector<entry_type>::const_iterator
-           i = pool.data->begin(); i != end ; ++i)
-        if (i->first == false)
-          delete i->second;
-        else
-          new_data.push_back (*i);
-
-      *pool.data = new_data;
-    }
-}
-
-
-
-template<typename VECTOR>
-inline
-std::size_t
-GrowingVectorMemory<VECTOR>::memory_consumption () const
-{
-  Threads::Mutex::ScopedLock lock(mutex);
-
-  std::size_t result = sizeof (*this);
-  const typename std::vector<entry_type>::const_iterator
-  end = pool.data->end();
-  for (typename std::vector<entry_type>::const_iterator
-       i = pool.data->begin(); i != end ; ++i)
-    result += sizeof (*i) + i->second->memory_consumption();
-
-  return result;
-}
-
-
-// -------------------------------------------------------------
-// explicit instantiations
 
+DEAL_II_NAMESPACE_OPEN
 
 #include "vector_memory.inst"
 

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.