From 7485388545ff9459faa83b53d4ca3bf5576763fc Mon Sep 17 00:00:00 2001 From: kronbichler Date: Mon, 22 Jul 2013 09:43:59 +0000 Subject: [PATCH] Fix global_dof_index in matrix-free class: It should only be used when dealing with deal.II's global dof indices, not the matrix-free dof indices which are transformed into a process-local index space and start at zero on all processors. This limits the MPI-local problem size to 4 billion, but i) that limit is still far away and ii) there were many internal uses of unsigned int left that would have needed re-organization. This change saves memory as dof indices can be saved in 32-bit integer instead of 64-bit integers and still allows problems with >4g global dofs. git-svn-id: https://svn.dealii.org/trunk@30090 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-37/step-37.cc | 4 +- .../include/deal.II/matrix_free/dof_info.h | 85 ++++--- .../deal.II/matrix_free/dof_info.templates.h | 215 +++++++++--------- .../deal.II/matrix_free/fe_evaluation.h | 10 +- .../include/deal.II/matrix_free/matrix_free.h | 28 +-- .../matrix_free/matrix_free.templates.h | 20 +- tests/matrix_free/step-37.cc | 4 +- 7 files changed, 177 insertions(+), 189 deletions(-) diff --git a/deal.II/examples/step-37/step-37.cc b/deal.II/examples/step-37/step-37.cc index e4f4ddb568..62c51b30ac 100644 --- a/deal.II/examples/step-37/step-37.cc +++ b/deal.II/examples/step-37/step-37.cc @@ -596,7 +596,7 @@ namespace Step37 { data.cell_loop (&LaplaceOperator::local_apply, this, dst, src); - const std::vector & + const std::vector & constrained_dofs = data.get_constrained_dofs(); for (unsigned int i=0; i & + const std::vector & constrained_dofs = data.get_constrained_dofs(); for (unsigned int i=0; i #include +#include + #include @@ -55,19 +57,6 @@ namespace internal */ struct DoFInfo { - - /** - * size_type of the dof_indicies object. - */ - - typedef std::vector::size_type size_dof; - - /** - * size_type of the constraint_indicators object. - */ - - typedef std::vector >::size_type size_constraint; - /** * Default empty constructor. */ @@ -87,13 +76,13 @@ namespace internal /** * Returns a pointer to the first index in the DoF row @p row. */ - const types::global_dof_index *begin_indices (const unsigned int row) const; + const unsigned int *begin_indices (const unsigned int row) const; /** * Returns a pointer to the one past the last DoF index in the row @p * row. */ - const types::global_dof_index *end_indices (const unsigned int row) const; + const unsigned int *end_indices (const unsigned int row) const; /** * Returns the number of entries in the indices field for the given row. @@ -124,13 +113,13 @@ namespace internal * Returns a pointer to the first index in the DoF row @p row for plain * indices (i.e., the entries where constraints are not embedded). */ - const types::global_dof_index *begin_indices_plain (const unsigned int row) const; + const unsigned int *begin_indices_plain (const unsigned int row) const; /** * Returns a pointer to the one past the last DoF index in the row @p * row (i.e., the entries where constraints are not embedded). */ - const types::global_dof_index *end_indices_plain (const unsigned int row) const; + const unsigned int *end_indices_plain (const unsigned int row) const; /** * Returns the FE index for a given finite element degree. If not in hp @@ -182,7 +171,7 @@ namespace internal */ void compute_renumber_serial (const std::vector &boundary_cells, const SizeInfo &size_info, - std::vector &renumbering); + std::vector &renumbering); /** * Reorganizes cells in the hp case without parallelism such that all @@ -191,7 +180,7 @@ namespace internal * DoFHandlers. */ void compute_renumber_hp_serial (SizeInfo &size_info, - std::vector &renumbering, + std::vector &renumbering, std::vector &irregular_cells); /** @@ -201,7 +190,7 @@ namespace internal */ void compute_renumber_parallel (const std::vector &boundary_cells, SizeInfo &size_info, - std::vector &renumbering); + std::vector &renumbering); /** * This method reorders the way cells are gone through based on a given @@ -210,7 +199,7 @@ namespace internal * vectorization. */ void reorder_cells (const SizeInfo &size_info, - const std::vector &renumbering, + const std::vector &renumbering, const std::vector &constraint_pool_row_index, const std::vector &irregular_cells, const unsigned int vectorization_length); @@ -241,7 +230,7 @@ namespace internal void make_thread_graph_partition_color (SizeInfo &size_info, TaskInfo &task_info, - std::vector &renumbering, + std::vector &renumbering, std::vector &irregular_cells, const bool hp_bool); @@ -262,7 +251,7 @@ namespace internal void make_thread_graph_partition_partition (SizeInfo &size_info, TaskInfo &task_info, - std::vector &renumbering, + std::vector &renumbering, std::vector &irregular_cells, const bool hp_bool); @@ -275,7 +264,7 @@ namespace internal void make_connectivity_graph (const SizeInfo &size_info, const TaskInfo &task_info, - const std::vector &renumbering, + const std::vector &renumbering, const std::vector &irregular_cells, const bool do_blocking, CompressedSimpleSparsityPattern &connectivity) const; @@ -317,19 +306,25 @@ namespace internal * certain structure in the indices, like indices for vector-valued * problems or for cells where not all vector components are filled. */ - std::vector > row_starts; + std::vector > row_starts; /** - * Stores the (global) indices of the degrees of freedom for each cell. This - * array also includes the indirect contributions from constraints, + * Stores the indices of the degrees of freedom for each cell. These + * indices are computed in MPI-local index space, i.e., each processor + * stores the locally owned indices as numbers between 0 and + * n_locally_owned_dofs-1 and ghost indices in the range + * n_locally_owned_dofs to + * n_locally_owned_dofs+n_ghost_dofs. The translation between + * this MPI-local index space and the global numbering of degrees of + * freedom is stored in the @p vector_partitioner data structure. + + * This array also includes the indirect contributions from constraints, * which are described by the @p constraint_indicator field. Because of * variable lengths of rows, this would be a vector of a * vector. However, we use one contiguous memory region and store the * rowstart in the variable @p row_starts. */ - std::vector dof_indices; + std::vector dof_indices; /** * This variable describes the position of constraints in terms of the @@ -354,13 +349,13 @@ namespace internal * This stores a (sorted) list of all locally owned degrees of freedom * that are constrained. */ - std::vector constrained_dofs; + std::vector constrained_dofs; /** * Stores the rowstart indices of the compressed row storage in the @p * plain_dof_indices fields. */ - std::vector row_starts_plain_indices; + std::vector row_starts_plain_indices; /** * Stores the indices of the degrees of freedom for each cell. This @@ -370,7 +365,7 @@ namespace internal * contiguous memory region and store the rowstart in the variable @p * row_starts_plain_indices. */ - std::vector plain_dof_indices; + std::vector plain_dof_indices; /** * Stores the dimension of the underlying DoFHandler. Since the indices @@ -432,11 +427,11 @@ namespace internal #ifndef DOXYGEN inline - const types::global_dof_index * + const unsigned int * DoFInfo::begin_indices (const unsigned int row) const { AssertIndexRange (row, row_starts.size()-1); - const unsigned int index = std_cxx1x::get<0>(row_starts[row]); + const unsigned int index = row_starts[row][0]; AssertIndexRange(index, dof_indices.size()+1); return dof_indices.empty() ? 0 : &dof_indices[0] + index; } @@ -444,11 +439,11 @@ namespace internal inline - const types::global_dof_index * + const unsigned int * DoFInfo::end_indices (const unsigned int row) const { AssertIndexRange (row, row_starts.size()-1); - const unsigned int index = std_cxx1x::get<0>(row_starts[row+1]); + const unsigned int index = row_starts[row+1][0]; AssertIndexRange(index, dof_indices.size()+1); return dof_indices.empty() ? 0 : &dof_indices[0] + index; } @@ -460,8 +455,7 @@ namespace internal DoFInfo::row_length_indices (const unsigned int row) const { AssertIndexRange (row, row_starts.size()-1); - return (std_cxx1x::get<0>(row_starts[row+1]) - - std_cxx1x::get<0>(row_starts[row])); + return (row_starts[row+1][0] - row_starts[row][0]); } @@ -471,7 +465,7 @@ namespace internal DoFInfo::begin_indicators (const unsigned int row) const { AssertIndexRange (row, row_starts.size()-1); - const unsigned int index = std_cxx1x::get<1>(row_starts[row]); + const unsigned int index = row_starts[row][1]; AssertIndexRange (index, constraint_indicator.size()+1); return constraint_indicator.empty() ? 0 : &constraint_indicator[0] + index; } @@ -483,7 +477,7 @@ namespace internal DoFInfo::end_indicators (const unsigned int row) const { AssertIndexRange (row, row_starts.size()-1); - const unsigned int index = std_cxx1x::get<1>(row_starts[row+1]); + const unsigned int index = row_starts[row+1][1]; AssertIndexRange (index, constraint_indicator.size()+1); return constraint_indicator.empty() ? 0 : &constraint_indicator[0] + index; } @@ -495,20 +489,19 @@ namespace internal DoFInfo::row_length_indicators (const unsigned int row) const { AssertIndexRange (row, row_starts.size()-1); - return (std_cxx1x::get<1>(row_starts[row+1]) - - std_cxx1x::get<1>(row_starts[row])); + return (row_starts[row+1][1] - row_starts[row][1]); } inline - const types::global_dof_index * + const unsigned int * DoFInfo::begin_indices_plain (const unsigned int row) const { // if we have no constraints, should take the data from dof_indices if (row_length_indicators(row) == 0) { - Assert (row_starts_plain_indices[row] == numbers::invalid_dof_index, + Assert (row_starts_plain_indices[row]==numbers::invalid_unsigned_int, ExcInternalError()); return begin_indices(row); } @@ -524,7 +517,7 @@ namespace internal inline - const types::global_dof_index * + const unsigned int * DoFInfo::end_indices_plain (const unsigned int row) const { return begin_indices_plain(row) + 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 6b6975afe3..16c53c268d 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 @@ -184,7 +184,10 @@ namespace internal const unsigned int n_mpi_procs = vector_partitioner->n_mpi_processes(); const types::global_dof_index first_owned = vector_partitioner->local_range().first; const types::global_dof_index last_owned = vector_partitioner->local_range().second; - const types::global_dof_index n_owned = last_owned - first_owned; + Assert (last_owned-first_owned < std::numeric_limits::max(), + ExcMessage("The size local range of owned indices must not " + "exceed the size of unsigned int")); + const unsigned int n_owned = last_owned - first_owned; std::pair constraint_iterator (0,0); unsigned int dofs_this_cell = (cell_active_fe_index.empty()) ? @@ -242,8 +245,7 @@ namespace internal (constraint_indices[j] < first_owned || constraint_indices[j] >= last_owned)) { - dof_indices.push_back - (n_owned + ghost_dofs.size()); + dof_indices.push_back (n_owned + ghost_dofs.size()); // collect ghosts so that we can later construct an // IndexSet for them. also store whether the current @@ -255,7 +257,8 @@ namespace internal // not ghost, so transform to the local index space // directly dof_indices.push_back - (constraint_indices[j] - first_owned); + (static_cast(constraint_indices[j] - + first_owned)); } } } @@ -276,7 +279,7 @@ no_constraint: else current_dof -= first_owned; - dof_indices.push_back (current_dof); + dof_indices.push_back (static_cast(current_dof)); // make sure constraint_iterator.first is always within the // bounds of unsigned short @@ -286,10 +289,9 @@ no_constraint: constraint_iterator.first++; } } - row_starts[cell_number+1] = std_cxx1x::tuple - (dof_indices.size(), constraint_indicator.size(), 0); + row_starts[cell_number+1][0] = dof_indices.size(); + row_starts[cell_number+1][1] = constraint_indicator.size(); + row_starts[cell_number+1][2] = 0; // now to the plain indices: in case we have constraints on this cell, // store the indices without the constraints resolve once again @@ -298,9 +300,8 @@ no_constraint: if (cell_number == 0) row_starts_plain_indices.resize (row_starts.size()); row_starts_plain_indices[cell_number] = plain_dof_indices.size(); - bool cell_has_constraints = - std_cxx1x::get<1>(row_starts[cell_number+1]) > - std_cxx1x::get<1>(row_starts[cell_number]); + bool cell_has_constraints = (row_starts[cell_number+1][1] > + row_starts[cell_number][1]); if (cell_has_constraints == true) { for (unsigned int i=0; i + (current_dof)); } } } @@ -331,12 +333,12 @@ no_constraint: Assert (boundary_cells.size() < row_starts.size(), ExcInternalError()); // sort ghost dofs and compress out duplicates - const types::global_dof_index n_owned = (vector_partitioner->local_range().second- - vector_partitioner->local_range().first); - const size_dof n_ghosts = ghost_dofs.size(); - unsigned int n_unique_ghosts= 0; + const unsigned int n_owned = (vector_partitioner->local_range().second- + vector_partitioner->local_range().first); + const std::size_t n_ghosts = ghost_dofs.size(); + unsigned int n_unique_ghosts= 0; #ifdef DEBUG - for (std::vector::iterator dof = dof_indices.begin(); + for (std::vector::iterator dof = dof_indices.begin(); dof!=dof_indices.end(); ++dof) AssertIndexRange (*dof, n_owned+n_ghosts); #endif @@ -349,7 +351,7 @@ no_constraint: // replace the temporary numbering of ghosts by the real number in // the index set, we need to store these values std::vector > ghost_origin(n_ghosts); - for (size_dof i=0; ighost_origin[i-1].first+1) + if (ghost_origin[i].first > ghost_origin[i-1].first+1) { ghost_indices.add_range (last_contiguous_start, ghost_origin[i-1].first+1); @@ -379,7 +381,7 @@ no_constraint: // dofs. the ghost index set should store the same number { AssertDimension (n_unique_ghosts, ghost_indices.n_elements()); - for (size_dof i=0; i - (begin_indices(boundary_cells[i])); + unsigned int *data_ptr = const_cast (begin_indices(boundary_cells[i])); - const types::global_dof_index *row_end = end_indices(boundary_cells[i]); + const unsigned int *row_end = end_indices(boundary_cells[i]); for ( ; data_ptr != row_end; ++data_ptr) *data_ptr = ((*data_ptr < n_owned) ? @@ -409,10 +410,8 @@ no_constraint: { if (row_length_indicators(boundary_cells[i]) > 0) { - types::global_dof_index *data_ptr = const_cast - (begin_indices_plain(boundary_cells[i])); - const types::global_dof_index *row_end = - (end_indices_plain(boundary_cells[i])); + unsigned int *data_ptr = const_cast (begin_indices_plain(boundary_cells[i])); + const unsigned int *row_end = end_indices_plain(boundary_cells[i]); for ( ; data_ptr != row_end; ++data_ptr) *data_ptr = ((*data_ptr < n_owned) ? @@ -425,8 +424,8 @@ no_constraint: } } - std::vector new_ghosts; - ghost_dofs.swap(new_ghosts); + std::vector empty; + ghost_dofs.swap(empty); // set the ghost indices now. need to cast away constness here, but that // is uncritical since we reset the Partitioner in the same initialize @@ -441,10 +440,10 @@ no_constraint: void DoFInfo::compute_renumber_serial (const std::vector &boundary_cells, const SizeInfo &size_info, - std::vector &renumbering) + std::vector &renumbering) { - std::vector reverse_numbering (size_info.n_active_cells, - numbers::invalid_unsigned_int); + std::vector reverse_numbering (size_info.n_active_cells, + numbers::invalid_unsigned_int); const unsigned int n_boundary_cells = boundary_cells.size(); for (unsigned int j=0; j &renumbering, + std::vector &renumbering, std::vector &irregular_cells) { if (max_fe_index < 2) @@ -557,10 +556,10 @@ no_constraint: void DoFInfo::compute_renumber_parallel (const std::vector &boundary_cells, SizeInfo &size_info, - std::vector &renumbering) + std::vector &renumbering) { - std::vector reverse_numbering (size_info.n_active_cells, - numbers::invalid_unsigned_int); + std::vector reverse_numbering (size_info.n_active_cells, + numbers::invalid_unsigned_int); const unsigned int n_boundary_cells = boundary_cells.size(); for (unsigned int j=0; j &renumbering, + const std::vector &renumbering, const std::vector &constraint_pool_row_index, const std::vector &irregular_cells, const unsigned int vectorization_length) @@ -610,12 +609,11 @@ no_constraint: std::swap (new_active_fe_index, cell_active_fe_index); } - std::vector > new_row_starts; - std::vector new_dof_indices; + std::vector > new_row_starts; + std::vector new_dof_indices; std::vector > new_constraint_indicator; - std::vector new_plain_indices, new_rowstart_plain; + std::vector new_plain_indices, new_rowstart_plain; unsigned int position_cell = 0; new_row_starts.resize (size_info.n_macro_cells + 1); new_dof_indices.reserve (dof_indices.size()); @@ -623,7 +621,7 @@ no_constraint: if (store_plain_indices == true) { new_rowstart_plain.resize (size_info.n_macro_cells + 1, - numbers::invalid_dof_index); + numbers::invalid_unsigned_int); new_plain_indices.reserve (plain_dof_indices.size()); } @@ -632,8 +630,8 @@ no_constraint: // vectors are adjacent, i.e., first dof index 0 for all vectors, then // dof index 1 for all vectors, and so on. This involves some extra // resorting. - std::vector glob_indices (vectorization_length); - std::vector plain_glob_indices (vectorization_length); + std::vector glob_indices (vectorization_length); + std::vector plain_glob_indices (vectorization_length); std::vector*> constr_ind(vectorization_length), constr_end(vectorization_length); std::vector index(vectorization_length); @@ -642,10 +640,9 @@ no_constraint: const unsigned int dofs_mcell = dofs_per_cell[cell_active_fe_index.size() == 0 ? 0 : cell_active_fe_index[i]] * vectorization_length; - new_row_starts[i] = - std_cxx1x::tuple - (new_dof_indices.size(), new_constraint_indicator.size(), - irregular_cells[i]); + new_row_starts[i][0] = new_dof_indices.size(); + new_row_starts[i][1] = new_constraint_indicator.size(); + new_row_starts[i][2] = irregular_cells[i]; const unsigned int n_comp = (irregular_cells[i]>0 ? irregular_cells[i] : vectorization_length); @@ -726,9 +723,9 @@ no_constraint: } AssertDimension (position_cell+1, row_starts.size()); - new_row_starts[size_info.n_macro_cells] = - std_cxx1x::tuple - (new_dof_indices.size(), new_constraint_indicator.size(), 0); + new_row_starts[size_info.n_macro_cells][0] = new_dof_indices.size(); + new_row_starts[size_info.n_macro_cells][1] = new_constraint_indicator.size(); + new_row_starts[size_info.n_macro_cells][2] = 0; AssertDimension(dof_indices.size(), new_dof_indices.size()); AssertDimension(constraint_indicator.size(), @@ -743,10 +740,10 @@ no_constraint: #ifdef DEBUG // sanity check 1: all indices should be smaller than the number of dofs // locally owned plus the number of ghosts - const types::global_dof_index index_range = (vector_partitioner->local_range().second- - vector_partitioner->local_range().first) - + vector_partitioner->ghost_indices().n_elements(); - for (size_dof i=0; ilocal_range().second- + vector_partitioner->local_range().first) + + vector_partitioner->ghost_indices().n_elements(); + for (std::size_t i=0; ilocal_range().second- - vector_partitioner->local_range().first); + const unsigned int local_size = (vector_partitioner->local_range().second- + vector_partitioner->local_range().first); for (unsigned int row=0; row &renumbering, - std::vector &irregular_cells, - const bool hp_bool) + (SizeInfo &size_info, + TaskInfo &task_info, + std::vector &renumbering, + std::vector &irregular_cells, + const bool hp_bool) { if (size_info.n_macro_cells == 0) return; @@ -929,12 +926,12 @@ no_constraint: // yet assigned a partition. std::vector cell_partition(task_info.n_blocks, size_info.n_macro_cells); - std::vector neighbor_list; - std::vector neighbor_neighbor_list; + std::vector neighbor_list; + std::vector neighbor_neighbor_list; // In element j of this variable, one puts the old number of the block // that should be the jth block in the new numeration. - std::vector partition_list (task_info.n_blocks,0); + std::vector partition_list (task_info.n_blocks,0); std::vector partition_color_list(task_info.n_blocks,0); // This vector points to the start of each partition. @@ -1127,7 +1124,7 @@ no_constraint: // check that the renumbering is one-to-one #ifdef DEBUG { - std::vector sorted_renumbering (renumbering); + std::vector sorted_renumbering (renumbering); std::sort(sorted_renumbering.begin(), sorted_renumbering.end()); for (unsigned int i=0; i &renumbering, - std::vector &irregular_cells, - const bool hp_bool) + (SizeInfo &size_info, + TaskInfo &task_info, + std::vector &renumbering, + std::vector &irregular_cells, + const bool hp_bool) { if (size_info.n_macro_cells == 0) return; @@ -1187,7 +1184,7 @@ no_constraint: // In element j of this variable, one puts the old number of the block // that should be the jth block in the new numeration. std::vector partition_list(size_info.n_active_cells,0); - std::vector partition_partition_list(size_info.n_active_cells,0); + std::vector partition_partition_list(size_info.n_active_cells,0); // This vector points to the start of each partition. std::vector partition_size(2,0); @@ -1442,7 +1439,7 @@ not_connect: missing_macros = 0; std::vector remaining_per_macro_cell (max_fe_index); - std::vector > + std::vector > renumbering_fe_index; unsigned int cell; bool filled = true; @@ -1655,7 +1652,7 @@ not_connect: reserve (2000); } - void reserve (const types::global_dof_index size) + void reserve (const std::size_t size) { if (size > 0) this->std::vector::reserve (size); @@ -1664,11 +1661,11 @@ not_connect: // insert a given entry. dat is a pointer within this vector (the user // needs to make sure that it really stays there) - void insert (const types::global_dof_index entry, + void insert (const unsigned int entry, std::vector::iterator &dat) { - AssertIndexRange (static_cast(dat - begin()), size()+1); - AssertIndexRange (static_cast(end() - dat), size()+1); + AssertIndexRange (static_cast(dat - begin()), size()+1); + AssertIndexRange (static_cast(end() - dat), size()+1); AssertIndexRange (size(), capacity()); while (dat != end() && *dat < entry) ++dat; @@ -1694,23 +1691,23 @@ not_connect: DoFInfo::make_connectivity_graph (const SizeInfo &size_info, const TaskInfo &task_info, - const std::vector &renumbering, + const std::vector &renumbering, const std::vector &irregular_cells, const bool do_blocking, CompressedSimpleSparsityPattern &connectivity) const { AssertDimension (row_starts.size()-1, size_info.n_active_cells); - const types::global_dof_index n_rows = + const unsigned int n_rows = (vector_partitioner->local_range().second- vector_partitioner->local_range().first) + vector_partitioner->ghost_indices().n_elements(); - const types::global_dof_index n_blocks = (do_blocking == true) ? - task_info.n_blocks : size_info.n_active_cells; + const unsigned int n_blocks = (do_blocking == true) ? + task_info.n_blocks : size_info.n_active_cells; // first determine row lengths std::vector row_lengths(n_rows); unsigned int cell_start = 0, mcell_start = 0; - std::vector scratch; + std::vector scratch; for (unsigned int block = 0; block < n_blocks; ++block) { // if we have the blocking variant (used in the coloring scheme), we @@ -1747,15 +1744,15 @@ not_connect: scratch.insert(scratch.end(), begin_indices(block), end_indices(block)); std::sort(scratch.begin(), scratch.end()); - const types::global_dof_index n_unique = + const unsigned int n_unique = std::unique(scratch.begin(), scratch.end())-scratch.begin(); - for (types::global_dof_index i=0; i 0) { @@ -1883,8 +1880,8 @@ not_connect: renumbering.resize (local_size, numbers::invalid_dof_index); types::global_dof_index counter = 0; - std::vector::iterator dof_ind = dof_indices.begin(), - end_ind = dof_indices.end(); + std::vector::iterator dof_ind = dof_indices.begin(), + end_ind = dof_indices.end(); for ( ; dof_ind != end_ind; ++dof_ind) { if (*dof_ind < local_size) @@ -1896,23 +1893,27 @@ not_connect: } AssertIndexRange (counter, local_size+1); - for (size_dof i=0; i new_constrained_dofs (constrained_dofs.size()); - for (size_dof i=0; i new_constrained_dofs (constrained_dofs.size()); + for (std::size_t i=0; i new_constrained_dofs[i-1], ExcInternalError()); #endif std::swap (constrained_dofs, new_constrained_dofs); + // transform indices to global index space + for (std::size_t i=0; ilocal_to_global(renumbering[i]); + AssertDimension (counter, renumbering.size()); } @@ -1922,8 +1923,7 @@ not_connect: DoFInfo::memory_consumption () const { std::size_t memory = sizeof(*this); - memory += (row_starts.capacity()*sizeof(std_cxx1x::tuple)); + memory += (row_starts.capacity()*sizeof(std_cxx1x::array)); memory += MemoryConsumption::memory_consumption (dof_indices); memory += MemoryConsumption::memory_consumption (row_starts_plain_indices); memory += MemoryConsumption::memory_consumption (plain_dof_indices); @@ -1941,7 +1941,7 @@ not_connect: { out << " Memory row starts indices: "; size_info.print_memory_statistics - (out, (row_starts.capacity()*sizeof(std_cxx1x::tuple))); + (out, (row_starts.capacity()*sizeof(std_cxx1x::array))); out << " Memory dof indices: "; size_info.print_memory_statistics (out, MemoryConsumption::memory_consumption (dof_indices)); @@ -1969,9 +1969,8 @@ not_connect: for (unsigned int row=0 ; row *con_it = begin_indicators(row), diff --git a/deal.II/include/deal.II/matrix_free/fe_evaluation.h b/deal.II/include/deal.II/matrix_free/fe_evaluation.h index 34a849fc3a..0e46cdb59c 100644 --- a/deal.II/include/deal.II/matrix_free/fe_evaluation.h +++ b/deal.II/include/deal.II/matrix_free/fe_evaluation.h @@ -1821,15 +1821,14 @@ FEEvaluationBase // loop over all local dofs. ind_local holds local number on cell, index // iterates over the elements of index_local_to_global and dof_indices // points to the global indices stored in index_local_to_global - const types::global_dof_index *dof_indices = dof_info.begin_indices(cell); + const unsigned int *dof_indices = dof_info.begin_indices(cell); const std::pair *indicators = dof_info.begin_indicators(cell); const std::pair *indicators_end = dof_info.end_indicators(cell); unsigned int ind_local = 0; - const unsigned int n_irreg_components_filled = - std_cxx1x::get<2>(dof_info.row_starts[cell]); + const unsigned int n_irreg_components_filled = dof_info.row_starts[cell][2]; const bool at_irregular_cell = n_irreg_components_filled > 0; // scalar case (or case when all components have the same degrees of freedom @@ -2443,10 +2442,9 @@ FEEvaluationBase // loop over all local dofs. ind_local holds local number on cell, index // iterates over the elements of index_local_to_global and dof_indices // points to the global indices stored in index_local_to_global - const types::global_dof_index *dof_indices = dof_info.begin_indices_plain(cell); + const unsigned int *dof_indices = dof_info.begin_indices_plain(cell); - const unsigned int n_irreg_components_filled = - std_cxx1x::get<2>(dof_info.row_starts[cell]); + const unsigned int n_irreg_components_filled = dof_info.row_starts[cell][2]; const bool at_irregular_cell = n_irreg_components_filled > 0; // scalar case (or case when all components have the same degrees of freedom diff --git a/deal.II/include/deal.II/matrix_free/matrix_free.h b/deal.II/include/deal.II/matrix_free/matrix_free.h index c76a95dd32..c3fe2de0f1 100644 --- a/deal.II/include/deal.II/matrix_free/matrix_free.h +++ b/deal.II/include/deal.II/matrix_free/matrix_free.h @@ -558,15 +558,18 @@ public: /** * Returns a list of all degrees of freedom that are constrained. The list - * is returned in local index space for the locally owned range of the - * vector, not in global numbers. In addition, it only returns the indices - * for degrees of freedom that are owned locally, not for ghosts. - */ - const std::vector & + * is returned in MPI-local index space for the locally owned range of the + * vector, not in global MPI index space that spans all MPI processors. To + * get numbers in global index space, call + * get_vector_partitioner()->local_to_global on an entry of the + * vector. In addition, it only returns the indices for degrees of freedom + * that are owned locally, not for ghosts. + */ + const std::vector & get_constrained_dofs (const unsigned int fe_component = 0) const; /** - * Calls renumber_dofs function in dof info which renumbers the degrees + * Calls renumber_dofs function in dof_info which renumbers the degrees * of freedom according to the ordering for parallelization. */ void renumber_dofs (std::vector &renumbering, @@ -976,7 +979,7 @@ MatrixFree::get_vector_partitioner (const unsigned int comp) const template inline -const std::vector & +const std::vector & MatrixFree::get_constrained_dofs (const unsigned int comp) const { AssertIndexRange (comp, n_components()); @@ -1207,8 +1210,7 @@ MatrixFree::get_cell_iterator(const unsigned int macro_cell_number, AssertIndexRange (dof_index, dof_handlers.n_dof_handlers); AssertIndexRange (macro_cell_number, size_info.n_macro_cells); AssertIndexRange (vector_number, vectorization_length); - const unsigned int irreg_filled = - std_cxx1x::get<2>(dof_info[dof_index].row_starts[macro_cell_number]); + const unsigned int irreg_filled = dof_info[dof_index].row_starts[macro_cell_number][2]; if (irreg_filled > 0) AssertIndexRange (vector_number, irreg_filled); #endif @@ -1246,8 +1248,7 @@ MatrixFree::get_hp_cell_iterator(const unsigned int macro_cell_numbe AssertIndexRange (dof_index, dof_handlers.n_dof_handlers); AssertIndexRange (macro_cell_number, size_info.n_macro_cells); AssertIndexRange (vector_number, vectorization_length); - const unsigned int irreg_filled = - std_cxx1x::get<2>(dof_info[dof_index].row_starts[macro_cell_number]); + const unsigned int irreg_filled = dof_info[dof_index].row_starts[macro_cell_number][2]; if (irreg_filled > 0) AssertIndexRange (vector_number, irreg_filled); #endif @@ -1269,7 +1270,7 @@ bool MatrixFree::at_irregular_cell (const unsigned int macro_cell) const { AssertIndexRange (macro_cell, size_info.n_macro_cells); - return std_cxx1x::get<2>(dof_info[0].row_starts[macro_cell]) > 0; + return dof_info[0].row_starts[macro_cell][2] > 0; } @@ -1280,8 +1281,7 @@ unsigned int MatrixFree::n_components_filled (const unsigned int macro_cell) const { AssertIndexRange (macro_cell, size_info.n_macro_cells); - const unsigned int n_filled = - std_cxx1x::get<2>(dof_info[0].row_starts[macro_cell]); + const unsigned int n_filled = dof_info[0].row_starts[macro_cell][2]; if (n_filled == 0) return VectorizedArray::n_array_elements; else 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 f01f109faf..d2477b4298 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 @@ -178,7 +178,7 @@ internal_reinit(const Mapping &mapping, dof_info[i].n_components = dof_handler[i]->get_fe().element_multiplicity(0); dof_info[i].dofs_per_cell.push_back(dof_handler[i]->get_fe().dofs_per_cell); dof_info[i].row_starts.resize(size_info.n_macro_cells+1); - std_cxx1x::get<2>(dof_info[i].row_starts.back()) = + dof_info[i].row_starts.back()[2] = cell_level_index.size() % VectorizedArray::n_array_elements; // if indices are not initialized, the cell_level_index might not be @@ -300,7 +300,7 @@ internal_reinit(const Mapping &mapping, dof_info[i].n_components = dof_handler[i]->get_fe()[0].element_multiplicity(0); dof_info[i].dofs_per_cell.push_back(dof_handler[i]->get_fe()[0].dofs_per_cell); dof_info[i].row_starts.resize(size_info.n_macro_cells+1); - std_cxx1x::get<2>(dof_info[i].row_starts.back()) = + dof_info[i].row_starts.back()[2] = cell_level_index.size() % VectorizedArray::n_array_elements; // if indices are not initialized, the cell_level_index might not be @@ -470,7 +470,6 @@ void MatrixFree::initialize_indices AssertDimension (n_fe, constraint.size()); std::vector local_dof_indices; - std::vector > ghost_dofs(n_fe); std::vector > > lexicographic_inv(n_fe); internal::MatrixFreeFunctions::ConstraintValues constraint_values; @@ -581,8 +580,9 @@ void MatrixFree::initialize_indices // initialize the arrays for indices dof_info[no].row_starts.resize (n_active_cells+1); - dof_info[no].row_starts[0] = std_cxx1x::tuple (0,0,0); + dof_info[no].row_starts[0][0] = 0; + dof_info[no].row_starts[0][1] = 0; + dof_info[no].row_starts[0][2] = 0; dof_info[no].dof_indices.reserve ((n_active_cells*dof_info[no].dofs_per_cell[0]*3)/2); @@ -593,11 +593,9 @@ void MatrixFree::initialize_indices end_index = dof_info[no].vector_partitioner->local_range().second; for (types::global_dof_index i=start_index; iis_constrained(i)==true) - dof_info[no].constrained_dofs.push_back(i-start_index); + dof_info[no].constrained_dofs. + push_back(static_cast(i-start_index)); } - - if (n_mpi_procs > 1) - ghost_dofs[no].reserve (locally_owned_set[no].n_elements()/10+1); } // extract all the global indices associated with the computation, and form @@ -689,7 +687,7 @@ void MatrixFree::initialize_indices // reorganize the indices in order to overlap communication in MPI with // computations: Place all cells with ghost indices into one chunk. Also // reorder cells so that we can parallelize by threads - std::vector renumbering; + std::vector renumbering; if (task_info.use_multithreading == true) { dof_info[0].compute_renumber_parallel (boundary_cells, size_info, @@ -719,7 +717,7 @@ void MatrixFree::initialize_indices // operations will then be done simultaneously). #ifdef DEBUG { - std::vector sorted_renumbering (renumbering); + std::vector sorted_renumbering (renumbering); std::sort (sorted_renumbering.begin(), sorted_renumbering.end()); for (unsigned int i=0; i & + const std::vector & constrained_dofs = data.get_constrained_dofs(); for (unsigned int i=0; i & + const std::vector & constrained_dofs = data.get_constrained_dofs(); for (unsigned int i=0; i