]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Move some template functions into the .cc file. There is a problem that in static...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 26 Oct 2010 19:43:59 +0000 (19:43 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 26 Oct 2010 19:43:59 +0000 (19:43 +0000)
git-svn-id: https://svn.dealii.org/trunk@22506 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/lac/vector_memory.h
deal.II/source/lac/vector_memory.cc

index b30a447ab846de2dd2571faab9c4582860f00f29..0e3619d90623f763186c63f7c0d9361b3bab0701 100644 (file)
@@ -248,19 +248,6 @@ class PrimitiveVectorMemory : public VectorMemory<VECTOR>
  * since they are reused, this should be of no concern. Additionally,
  * the destructor of the Pool warns about memory leaks.
  *
- * @note Instantiations for this class are provided for the types
- * Vector and BlockVector with number types <tt>float</tt>,
- * <tt>double</tt>, <tt>long double</tt>, <tt>@<std::complex@<float@>@></tt>,
- * <tt>@<std::complex@<double@>@></tt> and <tt>@<std::complex@<long
- * double@>@></tt>. In order to generate additional instances, it is
- * sufficient to define the Pool variable somewhere by
- * @code
- * #include <lac/vector_memory.h>
- *
- * template GrowingVectorMemory<MyVector>::Pool
- *   GrowingVectorMemory<MyVector>::pool;
- * @endcode
- *
  * @author Guido Kanschat, 1999, 2007
  */
 template<class VECTOR = dealii::Vector<double> >
@@ -354,13 +341,27 @@ class GrowingVectorMemory : public VectorMemory<VECTOR>
                                      */
     struct Pool
     {
-/// Standard constructor creating an empty pool
+                                        /**
+                                         * Standard constructor
+                                         * creating an empty pool
+                                         */
        Pool();
-/// Destructor. Frees memory and warns about memory leaks
+                                        /**
+                                         * Destructor. Frees memory
+                                         * and warns about memory
+                                         * leaks
+                                         */
        ~Pool();
-/// Create data vector; does nothing after first initialization
+                                        /**
+                                         * Create data vector; does
+                                         * nothing after first
+                                         * initialization
+                                         */
        void initialize(const unsigned int size);
-/// Pointer to the storage object
+                                        /**
+                                         * Pointer to the storage
+                                         * object
+                                         */
        std::vector<entry_type>* data;
     };
 
@@ -447,185 +448,6 @@ VECTOR * VectorMemory<VECTOR>::Pointer::operator -> () const
 }
 
 
-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
-  unsigned int n=0;
-  for (typename std::vector<entry_type>::iterator i=data->begin();
-       i != data->end();
-       ++i)
-    {
-      if (i->first == true)
-       ++n;
-      delete i->second;
-    }
-  delete data;
-}
-
-
-template <typename VECTOR>
-inline
-void
-GrowingVectorMemory<VECTOR>::Pool::initialize(const unsigned int 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 unsigned int initial_size,
-                                                 const bool log_statistics)
-
-               :
-               total_alloc(0),
-               current_alloc(0),
-               log_statistics(log_statistics)
-{
-  Threads::ThreadMutex::ScopedLock lock(mutex);
-  pool.initialize(initial_size);
-}
-
-
-template<typename VECTOR>
-inline
-GrowingVectorMemory<VECTOR>::~GrowingVectorMemory()
-{
-  AssertThrow(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::ThreadMutex::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::ThreadMutex::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::ThreadMutex::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
-unsigned int
-GrowingVectorMemory<VECTOR>::memory_consumption () const
-{
-  Threads::ThreadMutex::ScopedLock lock(mutex);
-
-  unsigned int 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;
-}
-
 #endif // DOXYGEN
 
 DEAL_II_NAMESPACE_CLOSE
index 0fbbb45406b4c020ef1dd84aeac1dcfa194657ca..267dca6fb8b1e26c673573a20379c7d24d5b97be 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2007, 2008 by the deal.II authors
+//    Copyright (C) 2007, 2008, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -31,10 +31,190 @@ typename GrowingVectorMemory<VECTOR>::Pool GrowingVectorMemory<VECTOR>::pool;
 template <typename VECTOR>
 Threads::ThreadMutex 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
+  unsigned int n=0;
+  for (typename std::vector<entry_type>::iterator i=data->begin();
+       i != data->end();
+       ++i)
+    {
+      if (i->first == true)
+       ++n;
+      delete i->second;
+    }
+  delete data;
+}
+
+
+template <typename VECTOR>
+inline
+void
+GrowingVectorMemory<VECTOR>::Pool::initialize(const unsigned int 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 unsigned int initial_size,
+                                                 const bool log_statistics)
+
+               :
+               total_alloc(0),
+               current_alloc(0),
+               log_statistics(log_statistics)
+{
+  Threads::ThreadMutex::ScopedLock lock(mutex);
+  pool.initialize(initial_size);
+}
+
+
+template<typename VECTOR>
+inline
+GrowingVectorMemory<VECTOR>::~GrowingVectorMemory()
+{
+  AssertThrow(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::ThreadMutex::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::ThreadMutex::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::ThreadMutex::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
+unsigned int
+GrowingVectorMemory<VECTOR>::memory_consumption () const
+{
+  Threads::ThreadMutex::ScopedLock lock(mutex);
+
+  unsigned int 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
 
+
 #include "vector_memory.inst"
 
 //TODO: Fold this into the list of vectors to be instantiated

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.