#include <deal.II/dofs/dof_handler.h>
#include <deal.II/matrix_free/helper_functions.h>
+#include <deal.II/base/std_cxx1x/array.h>
+
#include <memory>
*/
struct DoFInfo
{
-
- /**
- * size_type of the dof_indicies object.
- */
-
- typedef std::vector<types::global_dof_index>::size_type size_dof;
-
- /**
- * size_type of the constraint_indicators object.
- */
-
- typedef std::vector<std::pair<unsigned short, unsigned short> >::size_type size_constraint;
-
/**
* Default empty constructor.
*/
/**
* 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.
* 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
*/
void compute_renumber_serial (const std::vector<unsigned int> &boundary_cells,
const SizeInfo &size_info,
- std::vector<types::global_dof_index> &renumbering);
+ std::vector<unsigned int> &renumbering);
/**
* Reorganizes cells in the hp case without parallelism such that all
* DoFHandlers.
*/
void compute_renumber_hp_serial (SizeInfo &size_info,
- std::vector<types::global_dof_index> &renumbering,
+ std::vector<unsigned int> &renumbering,
std::vector<unsigned int> &irregular_cells);
/**
*/
void compute_renumber_parallel (const std::vector<unsigned int> &boundary_cells,
SizeInfo &size_info,
- std::vector<types::global_dof_index> &renumbering);
+ std::vector<unsigned int> &renumbering);
/**
* This method reorders the way cells are gone through based on a given
* vectorization.
*/
void reorder_cells (const SizeInfo &size_info,
- const std::vector<types::global_dof_index> &renumbering,
+ const std::vector<unsigned int> &renumbering,
const std::vector<unsigned int> &constraint_pool_row_index,
const std::vector<unsigned int> &irregular_cells,
const unsigned int vectorization_length);
void
make_thread_graph_partition_color (SizeInfo &size_info,
TaskInfo &task_info,
- std::vector<types::global_dof_index> &renumbering,
+ std::vector<unsigned int> &renumbering,
std::vector<unsigned int> &irregular_cells,
const bool hp_bool);
void
make_thread_graph_partition_partition (SizeInfo &size_info,
TaskInfo &task_info,
- std::vector<types::global_dof_index> &renumbering,
+ std::vector<unsigned int> &renumbering,
std::vector<unsigned int> &irregular_cells,
const bool hp_bool);
void
make_connectivity_graph (const SizeInfo &size_info,
const TaskInfo &task_info,
- const std::vector<types::global_dof_index> &renumbering,
+ const std::vector<unsigned int> &renumbering,
const std::vector<unsigned int> &irregular_cells,
const bool do_blocking,
CompressedSimpleSparsityPattern &connectivity) const;
* certain structure in the indices, like indices for vector-valued
* problems or for cells where not all vector components are filled.
*/
- std::vector<std_cxx1x::tuple<size_dof,
- size_constraint,
- unsigned int> > row_starts;
+ std::vector<std_cxx1x::array<unsigned int, 3> > 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 <tt>0</tt> and
+ * <tt>n_locally_owned_dofs-1</tt> and ghost indices in the range
+ * <tt>n_locally_owned_dofs</tt> to
+ * <tt>n_locally_owned_dofs+n_ghost_dofs</tt>. 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<types::global_dof_index> dof_indices;
+ std::vector<unsigned int> dof_indices;
/**
* This variable describes the position of constraints in terms of the
* This stores a (sorted) list of all locally owned degrees of freedom
* that are constrained.
*/
- std::vector<types::global_dof_index> constrained_dofs;
+ std::vector<unsigned int> constrained_dofs;
/**
* Stores the rowstart indices of the compressed row storage in the @p
* plain_dof_indices fields.
*/
- std::vector<types::global_dof_index> row_starts_plain_indices;
+ std::vector<unsigned int> row_starts_plain_indices;
/**
* Stores the indices of the degrees of freedom for each cell. This
* contiguous memory region and store the rowstart in the variable @p
* row_starts_plain_indices.
*/
- std::vector<types::global_dof_index> plain_dof_indices;
+ std::vector<unsigned int> plain_dof_indices;
/**
* Stores the dimension of the underlying DoFHandler. Since the indices
#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;
}
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;
}
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]);
}
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;
}
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;
}
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);
}
inline
- const types::global_dof_index *
+ const unsigned int *
DoFInfo::end_indices_plain (const unsigned int row) const
{
return begin_indices_plain(row) +
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<unsigned int>::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<unsigned short,unsigned short> constraint_iterator (0,0);
unsigned int dofs_this_cell = (cell_active_fe_index.empty()) ?
(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
// not ghost, so transform to the local index space
// directly
dof_indices.push_back
- (constraint_indices[j] - first_owned);
+ (static_cast<unsigned int>(constraint_indices[j] -
+ first_owned));
}
}
}
else
current_dof -= first_owned;
- dof_indices.push_back (current_dof);
+ dof_indices.push_back (static_cast<unsigned int>(current_dof));
// make sure constraint_iterator.first is always within the
// bounds of unsigned short
constraint_iterator.first++;
}
}
- row_starts[cell_number+1] = std_cxx1x::tuple<size_dof,
- size_constraint,
- unsigned int>
- (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
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<dofs_this_cell; ++i)
}
else
current_dof -= first_owned;
- plain_dof_indices.push_back (current_dof);
+ plain_dof_indices.push_back (static_cast<unsigned int>
+ (current_dof));
}
}
}
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<types::global_dof_index>::iterator dof = dof_indices.begin();
+ for (std::vector<unsigned int>::iterator dof = dof_indices.begin();
dof!=dof_indices.end(); ++dof)
AssertIndexRange (*dof, n_owned+n_ghosts);
#endif
// replace the temporary numbering of ghosts by the real number in
// the index set, we need to store these values
std::vector<std::pair<types::global_dof_index,unsigned int> > ghost_origin(n_ghosts);
- for (size_dof i=0; i<n_ghosts; ++i)
+ for (std::size_t i=0; i<n_ghosts; ++i)
{
ghost_origin[i].first = ghost_dofs[i];
ghost_origin[i].second = i;
types::global_dof_index last_contiguous_start = ghost_origin[0].first;
ghost_numbering[ghost_origin[0].second] = 0;
- for (size_dof i=1; i<n_ghosts; i++)
+ for (std::size_t i=1; i<n_ghosts; i++)
{
- if (ghost_origin[i].first>ghost_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);
// dofs. the ghost index set should store the same number
{
AssertDimension (n_unique_ghosts, ghost_indices.n_elements());
- for (size_dof i=0; i<n_ghosts; ++i)
+ for (std::size_t i=0; i<n_ghosts; ++i)
Assert (ghost_numbering[i] ==
ghost_indices.index_within_set(ghost_dofs[i]),
ExcInternalError());
const unsigned int n_boundary_cells = boundary_cells.size();
for (unsigned int i=0; i<n_boundary_cells; ++i)
{
- types::global_dof_index *data_ptr = const_cast<types::global_dof_index *>
- (begin_indices(boundary_cells[i]));
+ unsigned int *data_ptr = const_cast<unsigned int *> (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)
?
{
if (row_length_indicators(boundary_cells[i]) > 0)
{
- types::global_dof_index *data_ptr = const_cast<types::global_dof_index *>
- (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<unsigned int *> (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)
?
}
}
- std::vector<types::global_dof_index> new_ghosts;
- ghost_dofs.swap(new_ghosts);
+ std::vector<types::global_dof_index> 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
void
DoFInfo::compute_renumber_serial (const std::vector<unsigned int> &boundary_cells,
const SizeInfo &size_info,
- std::vector<types::global_dof_index> &renumbering)
+ std::vector<unsigned int> &renumbering)
{
- std::vector<types::global_dof_index> reverse_numbering (size_info.n_active_cells,
- numbers::invalid_unsigned_int);
+ std::vector<unsigned int> 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<n_boundary_cells; ++j)
reverse_numbering[boundary_cells[j]] =
void
DoFInfo::compute_renumber_hp_serial (SizeInfo &size_info,
- std::vector<types::global_dof_index> &renumbering,
+ std::vector<unsigned int> &renumbering,
std::vector<unsigned int> &irregular_cells)
{
if (max_fe_index < 2)
void
DoFInfo::compute_renumber_parallel (const std::vector<unsigned int> &boundary_cells,
SizeInfo &size_info,
- std::vector<types::global_dof_index> &renumbering)
+ std::vector<unsigned int> &renumbering)
{
- std::vector<types::global_dof_index> reverse_numbering (size_info.n_active_cells,
- numbers::invalid_unsigned_int);
+ std::vector<unsigned int> 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<n_boundary_cells; ++j)
reverse_numbering[boundary_cells[j]] = j;
void
DoFInfo::reorder_cells (const SizeInfo &size_info,
- const std::vector<types::global_dof_index> &renumbering,
+ const std::vector<unsigned int> &renumbering,
const std::vector<unsigned int> &constraint_pool_row_index,
const std::vector<unsigned int> &irregular_cells,
const unsigned int vectorization_length)
std::swap (new_active_fe_index, cell_active_fe_index);
}
- std::vector<std_cxx1x::tuple<size_dof, size_constraint,
- unsigned int> > new_row_starts;
- std::vector<types::global_dof_index> new_dof_indices;
+ std::vector<std_cxx1x::array<unsigned int, 3> > new_row_starts;
+ std::vector<unsigned int> new_dof_indices;
std::vector<std::pair<unsigned short,unsigned short> >
new_constraint_indicator;
- std::vector<types::global_dof_index> new_plain_indices, new_rowstart_plain;
+ std::vector<unsigned int> 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());
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());
}
// 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<const types::global_dof_index *> glob_indices (vectorization_length);
- std::vector<const types::global_dof_index *> plain_glob_indices (vectorization_length);
+ std::vector<const unsigned int *> glob_indices (vectorization_length);
+ std::vector<const unsigned int *> plain_glob_indices (vectorization_length);
std::vector<const std::pair<unsigned short,unsigned short>*>
constr_ind(vectorization_length), constr_end(vectorization_length);
std::vector<unsigned int> index(vectorization_length);
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<size_dof, size_constraint, unsigned int>
- (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);
}
AssertDimension (position_cell+1, row_starts.size());
- new_row_starts[size_info.n_macro_cells] =
- std_cxx1x::tuple<size_dof, size_constraint, unsigned int>
- (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(),
#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; i<dof_indices.size(); ++i)
+ const unsigned int index_range = (vector_partitioner->local_range().second-
+ vector_partitioner->local_range().first)
+ + vector_partitioner->ghost_indices().n_elements();
+ for (std::size_t i=0; i<dof_indices.size(); ++i)
AssertIndexRange (dof_indices[i], index_range);
// sanity check 2: for the constraint indicators, the first index should
// sanity check 3: all non-boundary cells should have indices that only
// refer to the locally owned range
- const types::global_dof_index local_size = (vector_partitioner->local_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<size_info.boundary_cells_start; ++row)
{
- const types::global_dof_index *ptr = begin_indices(row);
- const types::global_dof_index *end_ptr = end_indices (row);
+ const unsigned int *ptr = begin_indices(row);
+ const unsigned int *end_ptr = end_indices (row);
for ( ; ptr != end_ptr; ++ptr)
AssertIndexRange (*ptr, local_size);
}
for (unsigned int row=size_info.boundary_cells_end;
row<size_info.n_macro_cells; ++row)
{
- const types::global_dof_index *ptr = begin_indices(row);
- const types::global_dof_index *end_ptr = end_indices (row);
+ const unsigned int *ptr = begin_indices(row);
+ const unsigned int *end_ptr = end_indices (row);
for ( ; ptr != end_ptr; ++ptr)
AssertIndexRange (*ptr, local_size);
}
void DoFInfo::make_thread_graph_partition_color
- (SizeInfo &size_info,
- TaskInfo &task_info,
- std::vector<types::global_dof_index> &renumbering,
- std::vector<unsigned int> &irregular_cells,
- const bool hp_bool)
+ (SizeInfo &size_info,
+ TaskInfo &task_info,
+ std::vector<unsigned int> &renumbering,
+ std::vector<unsigned int> &irregular_cells,
+ const bool hp_bool)
{
if (size_info.n_macro_cells == 0)
return;
// yet assigned a partition.
std::vector<unsigned int> cell_partition(task_info.n_blocks,
size_info.n_macro_cells);
- std::vector<types::global_dof_index> neighbor_list;
- std::vector<types::global_dof_index> neighbor_neighbor_list;
+ std::vector<unsigned int> neighbor_list;
+ std::vector<unsigned int> 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<types::global_dof_index> partition_list (task_info.n_blocks,0);
+ std::vector<unsigned int> partition_list (task_info.n_blocks,0);
std::vector<unsigned int> partition_color_list(task_info.n_blocks,0);
// This vector points to the start of each partition.
// check that the renumbering is one-to-one
#ifdef DEBUG
{
- std::vector<types::global_dof_index> sorted_renumbering (renumbering);
+ std::vector<unsigned int> sorted_renumbering (renumbering);
std::sort(sorted_renumbering.begin(), sorted_renumbering.end());
for (unsigned int i=0; i<sorted_renumbering.size(); ++i)
Assert(sorted_renumbering[i] == i, ExcInternalError());
void
DoFInfo::make_thread_graph_partition_partition
- (SizeInfo &size_info,
- TaskInfo &task_info,
- std::vector<types::global_dof_index> &renumbering,
- std::vector<unsigned int> &irregular_cells,
- const bool hp_bool)
+ (SizeInfo &size_info,
+ TaskInfo &task_info,
+ std::vector<unsigned int> &renumbering,
+ std::vector<unsigned int> &irregular_cells,
+ const bool hp_bool)
{
if (size_info.n_macro_cells == 0)
return;
// 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<unsigned int> partition_list(size_info.n_active_cells,0);
- std::vector<types::global_dof_index> partition_partition_list(size_info.n_active_cells,0);
+ std::vector<unsigned int> partition_partition_list(size_info.n_active_cells,0);
// This vector points to the start of each partition.
std::vector<unsigned int> partition_size(2,0);
missing_macros = 0;
std::vector<unsigned int> remaining_per_macro_cell
(max_fe_index);
- std::vector<std::vector<types::global_dof_index> >
+ std::vector<std::vector<unsigned int> >
renumbering_fe_index;
unsigned int cell;
bool filled = true;
reserve (2000);
}
- void reserve (const types::global_dof_index size)
+ void reserve (const std::size_t size)
{
if (size > 0)
this->std::vector<types::global_dof_index>::reserve (size);
// 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<types::global_dof_index>::iterator &dat)
{
- AssertIndexRange (static_cast<types::global_dof_index>(dat - begin()), size()+1);
- AssertIndexRange (static_cast<types::global_dof_index>(end() - dat), size()+1);
+ AssertIndexRange (static_cast<std::size_t>(dat - begin()), size()+1);
+ AssertIndexRange (static_cast<std::size_t>(end() - dat), size()+1);
AssertIndexRange (size(), capacity());
while (dat != end() && *dat < entry)
++dat;
DoFInfo::make_connectivity_graph
(const SizeInfo &size_info,
const TaskInfo &task_info,
- const std::vector<types::global_dof_index> &renumbering,
+ const std::vector<unsigned int> &renumbering,
const std::vector<unsigned int> &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<unsigned int> row_lengths(n_rows);
unsigned int cell_start = 0, mcell_start = 0;
- std::vector<types::global_dof_index> scratch;
+ std::vector<unsigned int> scratch;
for (unsigned int block = 0; block < n_blocks; ++block)
{
// if we have the blocking variant (used in the coloring scheme), we
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<n_unique; ++i)
+ for (unsigned int i=0; i<n_unique; ++i)
row_lengths[scratch[i]]++;
}
}
// disregard dofs that only sit on one cell
- for (types::global_dof_index row=0; row<n_rows; ++row)
+ for (unsigned int row=0; row<n_rows; ++row)
if (row_lengths[row] == 1)
row_lengths[row] = 0;
for (unsigned int cell = cell_start; cell < cell_start+n_comp;
++cell)
{
- const types::global_dof_index
+ const unsigned int
*it = begin_indices (renumbering[cell]),
*end_cell = end_indices (renumbering[cell]);
for ( ; it != end_cell; ++it)
}
else
{
- const types::global_dof_index
+ const unsigned int
*it = begin_indices (block),
*end_cell = end_indices (block);
for ( ; it != end_cell; ++it)
++cell)
{
// apply renumbering when we do blocking
- const types::global_dof_index
+ const unsigned int
*it = begin_indices (renumbering[cell]),
*end_cell = end_indices (renumbering[cell]);
for ( ; it != end_cell; ++it)
}
else
{
- const types::global_dof_index *it = begin_indices (block),
- * end_cell = end_indices (block);
+ const unsigned int *it = begin_indices (block),
+ * end_cell = end_indices (block);
for ( ; it != end_cell; ++it)
if (row_lengths[*it] > 0)
{
renumbering.resize (local_size, numbers::invalid_dof_index);
types::global_dof_index counter = 0;
- std::vector<types::global_dof_index>::iterator dof_ind = dof_indices.begin(),
- end_ind = dof_indices.end();
+ std::vector<unsigned int>::iterator dof_ind = dof_indices.begin(),
+ end_ind = dof_indices.end();
for ( ; dof_ind != end_ind; ++dof_ind)
{
if (*dof_ind < local_size)
}
AssertIndexRange (counter, local_size+1);
- for (size_dof i=0; i<renumbering.size(); ++i)
+ for (std::size_t i=0; i<renumbering.size(); ++i)
if (renumbering[i] == numbers::invalid_dof_index)
renumbering[i] = counter++;
// adjust the constrained DoFs
- std::vector<types::global_dof_index> new_constrained_dofs (constrained_dofs.size());
- for (size_dof i=0; i<constrained_dofs.size(); ++i)
+ std::vector<unsigned int> new_constrained_dofs (constrained_dofs.size());
+ for (std::size_t i=0; i<constrained_dofs.size(); ++i)
new_constrained_dofs[i] = renumbering[constrained_dofs[i]];
// the new constrained DoFs should be sorted already as they are not
// contained in dof_indices and then get contiguous numbers
#ifdef DEBUG
- for (size_dof i=1; i<new_constrained_dofs.size(); ++i)
+ for (std::size_t i=1; i<new_constrained_dofs.size(); ++i)
Assert (new_constrained_dofs[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; i<renumbering.size(); ++i)
+ renumbering[i] = vector_partitioner->local_to_global(renumbering[i]);
+
AssertDimension (counter, renumbering.size());
}
DoFInfo::memory_consumption () const
{
std::size_t memory = sizeof(*this);
- memory += (row_starts.capacity()*sizeof(std_cxx1x::tuple<size_dof,
- size_constraint, unsigned int>));
+ memory += (row_starts.capacity()*sizeof(std_cxx1x::array<unsigned int,3>));
memory += MemoryConsumption::memory_consumption (dof_indices);
memory += MemoryConsumption::memory_consumption (row_starts_plain_indices);
memory += MemoryConsumption::memory_consumption (plain_dof_indices);
{
out << " Memory row starts indices: ";
size_info.print_memory_statistics
- (out, (row_starts.capacity()*sizeof(std_cxx1x::tuple<size_dof, size_constraint, unsigned int>)));
+ (out, (row_starts.capacity()*sizeof(std_cxx1x::array<unsigned int, 3>)));
out << " Memory dof indices: ";
size_info.print_memory_statistics
(out, MemoryConsumption::memory_consumption (dof_indices));
for (unsigned int row=0 ; row<n_rows ; ++row)
{
out << "Entries row " << row << ": ";
- const types::global_dof_index
- *glob_indices = begin_indices(row),
- *end_row = end_indices(row);
+ const unsigned int *glob_indices = begin_indices(row),
+ *end_row = end_indices(row);
unsigned int index = 0;
const std::pair<unsigned short,unsigned short>
*con_it = begin_indicators(row),