]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Fix bug in memory consumption of AlignedVector and in some other places. Mention...
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 9 Apr 2014 09:44:44 +0000 (09:44 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 9 Apr 2014 09:44:44 +0000 (09:44 +0000)
git-svn-id: https://svn.dealii.org/trunk@32748 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/examples/step-37/step-37.cc
deal.II/include/deal.II/base/aligned_vector.h
deal.II/include/deal.II/base/memory_consumption.h
deal.II/include/deal.II/base/table.h
deal.II/include/deal.II/matrix_free/dof_info.templates.h
deal.II/include/deal.II/matrix_free/mapping_info.templates.h
deal.II/include/deal.II/matrix_free/matrix_free.templates.h

index 365644606b93f7a4f96e457bbd420235ae9d2744..5c2f730fe7506acae6721b2e5edcf4fef84a490a 100644 (file)
@@ -148,6 +148,13 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> Changed: TableBase<N,T> now uses AlignedVector for storing data
+  instead of std::vector, which allows its use for VectorizedArray<Number>
+  data fields which require more alignment.
+  <br>
+  (Martin Kronbichler, 2014/04/09)
+  </li>
+
   <li> Improved: Piola transformation for FE_BDM is now active.
   <br>
   (Guido Kanschat, 2014/04/09)
index 9cf1b00a84537b230e8394a9c74b3969437a4fbc..2f071517505a9600968066857ce08ca68b0adf21 100644 (file)
@@ -221,17 +221,18 @@ namespace Step37
   // twice. Rather, we would keep this object in the main class and simply
   // store a reference.
   //
-  // @note Observe how we store the values for the coefficient: We use a
-  // vector type <code>AlignedVector<VectorizedArray<number> ></code>
-  // structure. One would think that one can use
-  // <code>std::vector<VectorizedArray<number> ></code> as well, but there are
-  // some technicalities with vectorization: A certain alignment of the data
-  // with the memory address boundaries is required (essentially, a
-  // VectorizedArray of 16 bytes length as in SSE needs to start at a memory
-  // address that is divisible by 16). The chosen class makes sure that this
-  // alignment is respected, whereas std::vector can in general not, which may
-  // lead to segmentation faults at strange places for some systems or
-  // suboptimal performance for other systems.
+  // @note Note that storing values of type
+  // <code>VectorizedArray<number></code> requires care: Here, we use the
+  // deal.II table class which is prepared to hold the data with correct
+  // alignment. However, storing it in e.g.
+  // <code>std::vector<VectorizedArray<number> ></code> is not possible with
+  // vectorization: A certain alignment of the data with the memory address
+  // boundaries is required (essentially, a VectorizedArray of 16 bytes length
+  // as in SSE needs to start at a memory address that is divisible by
+  // 16). The table class (as well as the AlignedVector class it is based on)
+  // makes sure that this alignment is respected, whereas std::vector can in
+  // general not, which may lead to segmentation faults at strange places for
+  // some systems or suboptimal performance for other systems.
   template <int dim, int fe_degree, typename number>
   class LaplaceOperator : public Subscriptor
   {
@@ -271,7 +272,7 @@ namespace Step37
     void evaluate_coefficient(const Coefficient<dim> &function);
 
     MatrixFree<dim,number>      data;
-    AlignedVector<VectorizedArray<number> > coefficient;
+    Table<2, VectorizedArray<number> > coefficient;
 
     Vector<number>  diagonal_values;
     bool            diagonal_is_available;
@@ -389,12 +390,12 @@ namespace Step37
   {
     const unsigned int n_cells = data.n_macro_cells();
     FEEvaluation<dim,fe_degree,fe_degree+1,1,number> phi (data);
-    coefficient.resize (n_cells * phi.n_q_points);
+    coefficient.reinit (n_cells, phi.n_q_points);
     for (unsigned int cell=0; cell<n_cells; ++cell)
       {
         phi.reinit (cell);
         for (unsigned int q=0; q<phi.n_q_points; ++q)
-          coefficient[cell*phi.n_q_points+q] =
+          coefficient(cell,q) =
             coefficient_function.value(phi.quadrature_point(q));
       }
   }
@@ -499,7 +500,6 @@ namespace Step37
                const std::pair<unsigned int,unsigned int> &cell_range) const
   {
     FEEvaluation<dim,fe_degree,fe_degree+1,1,number> phi (data);
-    AssertDimension (coefficient.size(), data.n_macro_cells() * phi.n_q_points);
 
     for (unsigned int cell=cell_range.first; cell<cell_range.second; ++cell)
       {
@@ -507,7 +507,7 @@ namespace Step37
         phi.read_dof_values(src);
         phi.evaluate (false,true,false);
         for (unsigned int q=0; q<phi.n_q_points; ++q)
-          phi.submit_gradient (coefficient[cell*phi.n_q_points+q] *
+          phi.submit_gradient (coefficient(cell,q) *
                                phi.get_gradient(q), q);
         phi.integrate (false,true);
         phi.distribute_local_to_global (dst);
index 54540bf352f319aa25aa9794d10996b979339ea5..a9331d151d67b3144256a65144c70c686aa25409 100644 (file)
@@ -364,8 +364,9 @@ namespace internal
       if (size == 0)
         return;
 
-      // do not use memset for long double because on some systems it does not
-      // completely fill its memory
+      // do not use memcmp for long double because on some systems it does not
+      // completely fill its memory and may lead to false positives in
+      // e.g. valgrind
       if (std_cxx1x::is_trivial<T>::value == true &&
           types_are_equal<T,long double>::value == false)
         {
@@ -406,6 +407,7 @@ namespace internal
     mutable T *destination_;
     bool trivial_element;
   };
+
 } // end of namespace internal
 
 
@@ -793,8 +795,9 @@ inline
 typename AlignedVector<T>::size_type
 AlignedVector<T>::memory_consumption () const
 {
-  size_type memory = sizeof(this);
-  memory += sizeof(T) * capacity();
+  size_type memory = sizeof(*this);
+  for (const T* t = _data ; t != _end_data; ++t)
+    memory += dealii::MemoryConsumption::memory_consumption(*t);
   return memory;
 }
 
index a42f4ca40167dbbd833f7080028bc16880512f67..b88f670ca2118ab274abdb2d3ed4712f8579cc88 100644 (file)
@@ -1,7 +1,7 @@
 // ---------------------------------------------------------------------
 // $Id$
 //
-// Copyright (C) 2000 - 2013 by the deal.II authors
+// Copyright (C) 2000 - 2014 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
 DEAL_II_NAMESPACE_OPEN
 
 
+// forward declaration
+template <typename T> class VectorizedArray;
+
+
 /**
  * This namespace provides functions helping to determine the amount
  * of memory used by objects. The goal is not necessarily to give the
@@ -198,6 +202,14 @@ namespace MemoryConsumption
   inline
   std::size_t memory_consumption (const std::complex<T> &);
 
+  /**
+   * Determine the amount of memory in bytes consumed by a
+   * <tt>VectorizedArray</tt> variable.
+   */
+  template <typename T>
+  inline
+  std::size_t memory_consumption (const VectorizedArray<T> &);
+
   /**
    * Determine an estimate of the
    * amount of memory in bytes
@@ -522,6 +534,7 @@ namespace MemoryConsumption
   }
 
 
+
   template <typename T>
   inline
   std::size_t memory_consumption (const std::complex<T> &)
@@ -531,6 +544,15 @@ namespace MemoryConsumption
 
 
 
+  template <typename T>
+  inline
+  std::size_t memory_consumption (const VectorizedArray<T> &)
+  {
+    return sizeof(VectorizedArray<T>);
+  }
+
+
+
   inline
   std::size_t memory_consumption (const std::string &s)
   {
index af4aa8d9f4d4e1d61f29c166ca9d031d50a0b733..a861fa0089c23ab8725ad668a5c15912169a5ccf 100644 (file)
@@ -1,7 +1,7 @@
 // ---------------------------------------------------------------------
 // $Id$
 //
-// Copyright (C) 2002 - 2013 by the deal.II authors
+// Copyright (C) 2002 - 2014 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
index 34317489718d1d8edc2f0b00b92b6e2ef68e5139..9dcd2d8aad2732d5953e8d5352ae50421dcf6778 100644 (file)
@@ -830,9 +830,8 @@ no_constraint:
       // set up partitions. if we just use coloring without partitions, do
       // nothing here, assume all cells to belong to the zero partition (that
       // we otherwise use for MPI boundary cells)
-      unsigned int partition = 0, start_up = 0, counter = 0;
-      unsigned int start_nonboundary = numbers::invalid_unsigned_int;
-      bool work = true;
+      unsigned int start_up = 0,
+        start_nonboundary = numbers::invalid_unsigned_int;
       if (task_info.use_coloring_only == false)
         {
           start_nonboundary =
@@ -920,6 +919,8 @@ no_constraint:
                                true, connectivity);
 
       // Create cell-block  partitioning.
+      unsigned int partition = 0, counter = 0;
+      bool work = true;
 
       // For each block of cells, this variable saves to which partitions the
       // block belongs. Initialize all to n_macro_cells to mark them as not
@@ -1131,8 +1132,8 @@ no_constraint:
       }
 #endif
       AssertDimension(counter,size_info.n_active_cells);
-      task_info.evens = (partition+1)>>1;
-      task_info.odds  = (partition)>>1;
+      task_info.evens = (partition+1)/2;
+      task_info.odds  = (partition)/2;
       task_info.n_blocked_workers = task_info.odds-
                                     (task_info.odds+task_info.evens+1)%2;
       task_info.n_workers = task_info.partition_color_blocks_data.size()-1-
index 74051a1c85a9373102e52e4774bc1d3675af8d68..f93060e5b782e1e91fcde50aefc9e61e1a8d15ff 100644 (file)
@@ -850,7 +850,7 @@ namespace internal
       memory += MemoryConsumption::memory_consumption (affine_data);
       memory += MemoryConsumption::memory_consumption (cartesian_data);
       memory += MemoryConsumption::memory_consumption (cell_type);
-      memory += sizeof (this);
+      memory += sizeof (*this);
       return memory;
     }
 
index ad08fa64aaed6e5fd6a713a26fad90bd1830d09b..c01b4bc4d11e3a63aafa24f113fde289751f7395 100644 (file)
@@ -762,7 +762,7 @@ std::size_t MatrixFree<dim,Number>::memory_consumption () const
   memory += MemoryConsumption::memory_consumption (constraint_pool_data);
   memory += MemoryConsumption::memory_consumption (constraint_pool_row_index);
   memory += MemoryConsumption::memory_consumption (task_info);
-  memory += sizeof(this);
+  memory += sizeof(*this);
   memory += mapping_info.memory_consumption();
   return memory;
 }
@@ -853,7 +853,8 @@ namespace internal
     std::size_t
     TaskInfo::memory_consumption () const
     {
-      return (MemoryConsumption::memory_consumption (partition_color_blocks_row_index) +
+      return (sizeof(*this)+
+              MemoryConsumption::memory_consumption (partition_color_blocks_row_index) +
               MemoryConsumption::memory_consumption (partition_color_blocks_data)+
               MemoryConsumption::memory_consumption (partition_evens) +
               MemoryConsumption::memory_consumption (partition_odds) +

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.