From 5086b968edaa8f7ce65753d35e89644f9bf2e214 Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 9 Apr 2014 09:44:44 +0000 Subject: [PATCH] Fix bug in memory consumption of AlignedVector and in some other places. Mention change in Table class and use the table in step-37 git-svn-id: https://svn.dealii.org/trunk@32748 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 7 ++++ deal.II/examples/step-37/step-37.cc | 32 +++++++++---------- deal.II/include/deal.II/base/aligned_vector.h | 11 ++++--- .../include/deal.II/base/memory_consumption.h | 24 +++++++++++++- deal.II/include/deal.II/base/table.h | 2 +- .../deal.II/matrix_free/dof_info.templates.h | 11 ++++--- .../matrix_free/mapping_info.templates.h | 2 +- .../matrix_free/matrix_free.templates.h | 5 +-- 8 files changed, 64 insertions(+), 30 deletions(-) diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 365644606b..5c2f730fe7 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -148,6 +148,13 @@ inconvenience this causes.

Specific improvements

    +
  1. Changed: TableBase now uses AlignedVector for storing data + instead of std::vector, which allows its use for VectorizedArray + data fields which require more alignment. +
    + (Martin Kronbichler, 2014/04/09) +
  2. +
  3. Improved: Piola transformation for FE_BDM is now active.
    (Guido Kanschat, 2014/04/09) diff --git a/deal.II/examples/step-37/step-37.cc b/deal.II/examples/step-37/step-37.cc index 9cf1b00a84..2f07151750 100644 --- a/deal.II/examples/step-37/step-37.cc +++ b/deal.II/examples/step-37/step-37.cc @@ -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 AlignedVector > - // structure. One would think that one can use - // std::vector > 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 + // VectorizedArray 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. + // std::vector > 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 class LaplaceOperator : public Subscriptor { @@ -271,7 +272,7 @@ namespace Step37 void evaluate_coefficient(const Coefficient &function); MatrixFree data; - AlignedVector > coefficient; + Table<2, VectorizedArray > coefficient; Vector diagonal_values; bool diagonal_is_available; @@ -389,12 +390,12 @@ namespace Step37 { const unsigned int n_cells = data.n_macro_cells(); FEEvaluation phi (data); - coefficient.resize (n_cells * phi.n_q_points); + coefficient.reinit (n_cells, phi.n_q_points); for (unsigned int cell=0; cell &cell_range) const { FEEvaluation phi (data); - AssertDimension (coefficient.size(), data.n_macro_cells() * phi.n_q_points); for (unsigned int cell=cell_range.first; cell::value == true && types_are_equal::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::size_type AlignedVector::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; } diff --git a/deal.II/include/deal.II/base/memory_consumption.h b/deal.II/include/deal.II/base/memory_consumption.h index a42f4ca401..b88f670ca2 100644 --- a/deal.II/include/deal.II/base/memory_consumption.h +++ b/deal.II/include/deal.II/base/memory_consumption.h @@ -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. // @@ -29,6 +29,10 @@ DEAL_II_NAMESPACE_OPEN +// forward declaration +template 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 &); + /** + * Determine the amount of memory in bytes consumed by a + * VectorizedArray variable. + */ + template + inline + std::size_t memory_consumption (const VectorizedArray &); + /** * Determine an estimate of the * amount of memory in bytes @@ -522,6 +534,7 @@ namespace MemoryConsumption } + template inline std::size_t memory_consumption (const std::complex &) @@ -531,6 +544,15 @@ namespace MemoryConsumption + template + inline + std::size_t memory_consumption (const VectorizedArray &) + { + return sizeof(VectorizedArray); + } + + + inline std::size_t memory_consumption (const std::string &s) { diff --git a/deal.II/include/deal.II/base/table.h b/deal.II/include/deal.II/base/table.h index af4aa8d9f4..a861fa0089 100644 --- a/deal.II/include/deal.II/base/table.h +++ b/deal.II/include/deal.II/base/table.h @@ -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. // diff --git a/deal.II/include/deal.II/matrix_free/dof_info.templates.h b/deal.II/include/deal.II/matrix_free/dof_info.templates.h index 3431748971..9dcd2d8aad 100644 --- a/deal.II/include/deal.II/matrix_free/dof_info.templates.h +++ b/deal.II/include/deal.II/matrix_free/dof_info.templates.h @@ -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- diff --git a/deal.II/include/deal.II/matrix_free/mapping_info.templates.h b/deal.II/include/deal.II/matrix_free/mapping_info.templates.h index 74051a1c85..f93060e5b7 100644 --- a/deal.II/include/deal.II/matrix_free/mapping_info.templates.h +++ b/deal.II/include/deal.II/matrix_free/mapping_info.templates.h @@ -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; } diff --git a/deal.II/include/deal.II/matrix_free/matrix_free.templates.h b/deal.II/include/deal.II/matrix_free/matrix_free.templates.h index ad08fa64aa..c01b4bc4d1 100644 --- a/deal.II/include/deal.II/matrix_free/matrix_free.templates.h +++ b/deal.II/include/deal.II/matrix_free/matrix_free.templates.h @@ -762,7 +762,7 @@ std::size_t MatrixFree::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) + -- 2.39.5