SparseMatrix::in_local_range (const size_type index) const
{
TrilinosWrappers::types::int_type begin, end;
- begin = matrix->RowMap().MinMyGID64();
- end = matrix->RowMap().MaxMyGID64()+1;
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ begin = matrix->RowMap().MinMyGID();
+ end = matrix->RowMap().MaxMyGID()+1;
+#else
+ begin = matrix->RowMap().MinMyGID();
+ end = matrix->RowMap().MaxMyGID()+1;
+#endif
return ((index >= static_cast<size_type>(begin)) &&
(index < static_cast<size_type>(end)));
SparseMatrix::size_type
SparseMatrix::m () const
{
- return matrix -> NumGlobalRows64();
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ return matrix->NumGlobalRows();
+#else
+ return matrix->NumGlobalRows64();
+#endif
}
SparseMatrix::size_type
SparseMatrix::n () const
{
- return matrix -> NumGlobalCols64();
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ return matrix->NumGlobalCols();
+#else
+ return matrix->NumGlobalCols64();
+#endif
}
SparseMatrix::local_range () const
{
size_type begin, end;
- begin = matrix -> RowMap().MinMyGID64();
- end = matrix -> RowMap().MaxMyGID64()+1;
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ begin = matrix->RowMap().MinMyGID();
+ end = matrix->RowMap().MaxMyGID()+1;
+#else
+ begin = matrix->RowMap().MinMyGID64();
+ end = matrix->RowMap().MaxMyGID64()+1;
+#endif
return std::make_pair (begin, end);
}
SparseMatrix::size_type
SparseMatrix::n_nonzero_elements () const
{
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ return matrix->NumGlobalNonzeros();
+#else
return matrix->NumGlobalNonzeros64();
+#endif
}
SparsityPattern::in_local_range (const size_type index) const
{
TrilinosWrappers::types::int_type begin, end;
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ begin = graph->RowMap().MinMyGID();
+ end = graph->RowMap().MaxMyGID()+1;
+#else
begin = graph->RowMap().MinMyGID64();
end = graph->RowMap().MaxMyGID64()+1;
+#endif
return ((index >= static_cast<size_type>(begin)) &&
(index < static_cast<size_type>(end)));
namespace TrilinosWrappers
{
class SparseMatrix;
+
+ namespace
+ {
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ // define a helper function that queries the global ID of local ID of
+ // an Epetra_BlockMap object by calling either the 32- or 64-bit
+ // function necessary.
+ int gid(const Epetra_BlockMap &map, int i)
+ {
+ return map.GID(i);
+ }
+#else
+ // define a helper function that queries the global ID of local ID of
+ // an Epetra_BlockMap object by calling either the 32- or 64-bit
+ // function necessary.
+ long long int gid(const Epetra_BlockMap &map, int i)
+ {
+ return map.GID64(i);
+ }
+#endif
+ }
+
/**
* Namespace for Trilinos vector classes that work in parallel over
* MPI. This namespace is restricted to vectors only, whereas matrices
// deal.II might not use doubles, so
// that a direct access is not possible.
for (int i=0; i<size; ++i)
- (*vector)[0][i] = v(parallel_partitioner.GID64(i));
+ (*vector)[0][i] = v(gid(parallel_partitioner,i));
}
// forward declaration
class VectorBase;
-
/**
* @cond internal
*/
VectorBase::size_type
VectorBase::size () const
{
- return (size_type) (vector->Map().MaxAllGID64() + 1 -
- vector->Map().MinAllGID64());
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ return (size_type) (vector->Map().MaxAllGID() + 1 - vector->Map().MinAllGID());
+#else
+ return (size_type) (vector->Map().MaxAllGID64() + 1 - vector->Map().MinAllGID64());
+#endif
}
VectorBase::local_range () const
{
TrilinosWrappers::types::int_type begin, end;
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ begin = vector->Map().MinMyGID();
+ end = vector->Map().MaxMyGID()+1;
+#else
begin = vector->Map().MinMyGID64();
end = vector->Map().MaxMyGID64()+1;
+#endif
return std::make_pair (begin, end);
}
namespace TrilinosWrappers
{
+ namespace
+ {
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ int global_length (const Epetra_MultiVector &vector)
+ {
+ return vector.GlobalLength();
+ }
+#else
+ long long int global_length (const Epetra_MultiVector &vector)
+ {
+ return vector.GlobalLength64();
+ }
+#endif
+ }
PreconditionBase::PreconditionBase()
#ifdef DEAL_II_WITH_MPI
Assert (n_relevant_rows == my_size,
ExcDimensionMismatch(n_relevant_rows, my_size));
Assert (n_rows ==
- static_cast<size_type>(distributed_constant_modes.GlobalLength64()),
+ static_cast<size_type>(global_length(distributed_constant_modes)),
ExcDimensionMismatch(n_rows,
- distributed_constant_modes.GlobalLength64()));
+ global_length(distributed_constant_modes)));
// Reshape null space as a
// contiguous vector of
{
namespace
{
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
// define a helper function that queries the size of an Epetra_Map object
// by calling either the 32- or 64-bit function necessary, and returns the
// result in the correct data type so that we can use it in calling other
// Epetra member functions that are overloaded by index type
-#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
int n_global_elements (const Epetra_BlockMap &map)
{
return map.NumGlobalElements();
{
return map.MaxMyGID();
}
+
+ int n_global_cols(const Epetra_CrsGraph &graph)
+ {
+ return graph.NumGlobalCols();
+ }
+
+ int global_column_index(const Epetra_CrsMatrix &matrix, int i)
+ {
+ return matrix.GCID(i);
+ }
+
+ int global_row_index(const Epetra_CrsMatrix &matrix, int i)
+ {
+ return matrix.GRID(i);
+ }
#else
+ // define a helper function that queries the size of an Epetra_Map object
+ // by calling either the 32- or 64-bit function necessary, and returns the
+ // result in the correct data type so that we can use it in calling other
+ // Epetra member functions that are overloaded by index type
long long int n_global_elements (const Epetra_BlockMap &map)
{
return map.NumGlobalElements64();
{
return map.MaxMyGID64();
}
+
+ long long int n_global_cols(const Epetra_CrsGraph &graph)
+ {
+ return graph.NumGlobalCols64();
+ }
+
+ long long int global_column_index(const Epetra_CrsMatrix &matrix, int i)
+ {
+ return matrix.GCID64(i);
+ }
+
+ long long int global_row_index(const Epetra_CrsMatrix &matrix, int i)
+ {
+ return matrix.GRID64(i);
+ }
#endif
}
graph->OptimizeStorage();
// check whether we got the number of columns right.
- AssertDimension (sparsity_pattern.n_cols(),static_cast<size_type>(graph->NumGlobalCols64()));
+ AssertDimension (sparsity_pattern.n_cols(),static_cast<size_type>(
+ n_global_cols(*graph)));
// And now finally generate the matrix.
matrix.reset (new Epetra_FECrsMatrix(Copy, *graph, false));
{
matrix->ExtractMyRowView (i, num_entries, values, indices);
for (TrilinosWrappers::types::int_type j=0; j<num_entries; ++j)
- out << "(" << matrix->GRID64(i) << "," << matrix->GCID64(indices[j]) << ") "
+ out << "(" << global_row_index(*matrix,i) << ","
+ << global_column_index(*matrix,indices[j]) << ") "
<< values[j] << std::endl;
}
}
return map.MaxMyGID();
}
- int n_global_rows(const Epetra_FECrsGraph&graph)
+ int n_global_rows(const Epetra_CrsGraph &graph)
{
return graph.NumGlobalRows();
}
- int n_global_cols(const Epetra_FECrsGraph&graph)
+ int n_global_cols(const Epetra_CrsGraph &graph)
{
return graph.NumGlobalCols();
}
- int n_global_entries(const Epetra_FECrsGraph&graph)
+ int n_global_entries(const Epetra_CrsGraph &graph)
{
return graph.NumGlobalEntries();
}
+
+ int global_row_index(const Epetra_CrsGraph &graph, int i)
+ {
+ return graph.GRID(i);
+ }
#else
long long int n_global_elements (const Epetra_BlockMap &map)
{
return map.MaxMyGID64();
}
- long long int n_global_rows(const Epetra_FECrsGraph&graph)
+ long long int n_global_rows(const Epetra_CrsGraph &graph)
{
return graph.NumGlobalRows64();
}
- long long int n_global_cols(const Epetra_FECrsGraph&graph)
+ long long int n_global_cols(const Epetra_CrsGraph &graph)
{
return graph.NumGlobalCols64();
}
- long long int n_global_entries(const Epetra_FECrsGraph&graph)
+ long long int n_global_entries(const Epetra_CrsGraph &graph)
{
- return graph.NumGlobalEntries();
+ return graph.NumGlobalEntries64();
+ }
+
+ long long int global_row_index(const Epetra_CrsGraph &graph, int i)
+ {
+ return graph.GRID64(i);
}
#endif
}
// release memory before reallocation
graph.reset ();
AssertDimension (n_entries_per_row.size(),
- static_cast<size_type>(input_row_map.NumGlobalElements64()));
+ static_cast<size_type>(n_global_elements(input_row_map)));
column_space_map.reset (new Epetra_Map (input_col_map));
compressed = false;
graph.reset ();
AssertDimension (sp.n_rows(),
- static_cast<size_type>(input_row_map.NumGlobalElements64()));
+ static_cast<size_type>(n_global_elements(input_row_map)));
AssertDimension (sp.n_cols(),
- static_cast<size_type>(input_col_map.NumGlobalElements64()));
+ static_cast<size_type>(n_global_elements(input_col_map)));
column_space_map.reset (new Epetra_Map (input_col_map));
compressed = false;
if (graph->Filled() == true)
n_cols = n_global_cols(*graph);
else
- n_cols = column_space_map->NumGlobalElements64();
+ n_cols = n_global_elements(*column_space_map);
return n_cols;
}
{
graph->ExtractMyRowView (i, num_entries, indices);
for (TrilinosWrappers::types::int_type j=0; j<num_entries; ++j)
- out << "(" << i << "," << indices[graph->GRID64(j)] << ") "
+ out << "(" << i << "," << indices[global_row_index(*graph,j)] << ") "
<< std::endl;
}
}
// j horizontal, gnuplot output is
// x-y, that is we have to exchange
// the order of output
- out << indices[graph->GRID64(static_cast<TrilinosWrappers::types::int_type>(j))]
+ out << indices[global_row_index(*graph,static_cast<int>(j))]
<< " " << -static_cast<TrilinosWrappers::types::int_type>(row)
<< std::endl;
}
namespace TrilinosWrappers
{
+ namespace
+ {
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ // define a helper function that queries the size of an Epetra_BlockMap object
+ // by calling either the 32- or 64-bit function necessary, and returns the
+ // result in the correct data type so that we can use it in calling other
+ // Epetra member functions that are overloaded by index type
+ int n_global_elements (const Epetra_BlockMap &map)
+ {
+ return map.NumGlobalElements();
+ }
+ // define a helper function that queries the pointer to internal array
+ // containing list of global IDs assigned to the calling processor
+ // by calling either the 32- or 64-bit function necessary, and returns the
+ // result in the correct data type so that we can use it in calling other
+ // Epetra member functions that are overloaded by index type
+ int* my_global_elements(const Epetra_BlockMap &map)
+ {
+ return map.MyGlobalElements();
+ }
+ // define a helper function that queries the global vector length of an
+ // Epetra_FEVector object by calling either the 32- or 64-bit
+ // function necessary.
+ int global_length(const Epetra_FEVector &vector)
+ {
+ return vector.GlobalLength();
+ }
+#else
+ // define a helper function that queries the size of an Epetra_BlockMap object
+ // by calling either the 32- or 64-bit function necessary, and returns the
+ // result in the correct data type so that we can use it in calling other
+ // Epetra member functions that are overloaded by index type
+ long long int n_global_elements (const Epetra_BlockMap &map)
+ {
+ return map.NumGlobalElements64();
+ }
+ // define a helper function that queries the pointer to internal array
+ // containing list of global IDs assigned to the calling processor
+ // by calling either the 32- or 64-bit function necessary, and returns the
+ // result in the correct data type so that we can use it in calling other
+ // Epetra member functions that are overloaded by index type
+ long long int* my_global_elements(const Epetra_BlockMap &map)
+ {
+ return map.MyGlobalElements64();
+ }
+ // define a helper function that queries the global vector length of an
+ // Epetra_FEVector object by calling either the 32- or 64-bit
+ // function necessary.
+ long long int global_length(const Epetra_FEVector &vector)
+ {
+ return vector.GlobalLength64();
+ }
+#endif
+ }
+
namespace MPI
{
:
VectorBase()
{
- AssertThrow (input_map.NumGlobalElements64() == v.vector->Map().NumGlobalElements64(),
- ExcDimensionMismatch (input_map.NumGlobalElements64(),
- v.vector->Map().NumGlobalElements64()));
+ AssertThrow (n_global_elements(input_map) == n_global_elements(v.vector->Map()),
+ ExcDimensionMismatch (n_global_elements(input_map),
+ n_global_elements(v.vector->Map())));
last_action = Zero;
VectorBase()
{
AssertThrow (parallel_partitioner.size() ==
- static_cast<size_type>(v.vector->Map().NumGlobalElements64()),
+ static_cast<size_type>(n_global_elements(v.vector->Map())),
ExcDimensionMismatch (parallel_partitioner.size(),
- v.vector->Map().NumGlobalElements64()));
+ n_global_elements(v.vector->Map())));
last_action = Zero;
for (size_type block=0; block<v.n_blocks(); ++block)
{
TrilinosWrappers::types::int_type *glob_elements =
- v.block(block).vector_partitioner().MyGlobalElements64();
+ my_global_elements(v.block(block).vector_partitioner());
for (size_type i=0; i<v.block(block).local_size(); ++i)
global_ids[added_elements++] = glob_elements[i] + block_offset;
block_offset += v.block(block).size();
if (import_data == true)
{
- AssertThrow (static_cast<size_type>(actual_vec->GlobalLength64())
+ AssertThrow (static_cast<size_type>(global_length(*actual_vec))
== v.size(),
- ExcDimensionMismatch (actual_vec->GlobalLength64(),
+ ExcDimensionMismatch (global_length(*actual_vec),
v.size()));
Epetra_Import data_exchange (vector->Map(), actual_vec->Map());
Vector::Vector (const Epetra_Map &input_map)
{
last_action = Zero;
- Epetra_LocalMap map (input_map.NumGlobalElements64(),
+ Epetra_LocalMap map (n_global_elements(input_map),
input_map.IndexBase(),
input_map.Comm());
vector.reset (new Epetra_FEVector(map));
Vector::Vector (const VectorBase &v)
{
last_action = Zero;
- Epetra_LocalMap map (v.vector->Map().NumGlobalElements64(),
+ Epetra_LocalMap map (n_global_elements(v.vector->Map()),
v.vector->Map().IndexBase(),
v.vector->Map().Comm());
vector.reset (new Epetra_FEVector(map));
Vector::reinit (const Epetra_Map &input_map,
const bool fast)
{
- if (vector->Map().NumGlobalElements64() != input_map.NumGlobalElements64())
+ if (n_global_elements(vector->Map()) != n_global_elements(input_map))
{
- Epetra_LocalMap map (input_map.NumGlobalElements64(),
+ Epetra_LocalMap map (n_global_elements(input_map),
input_map.IndexBase(),
input_map.Comm());
vector.reset (new Epetra_FEVector (map));
const MPI_Comm &communicator,
const bool fast)
{
- if (vector->Map().NumGlobalElements64() !=
+ if (n_global_elements(vector->Map()) !=
static_cast<TrilinosWrappers::types::int_type>(partitioning.size()))
{
Epetra_LocalMap map (static_cast<TrilinosWrappers::types::int_type>(partitioning.size()),
{
if (local_range() != v.local_range())
{
- Epetra_LocalMap map (v.vector->GlobalLength64(),
+ Epetra_LocalMap map (global_length(*(v.vector)),
v.vector->Map().IndexBase(),
v.vector->Comm());
vector.reset (new Epetra_FEVector(map));
{
if (size() != v.size())
{
- Epetra_LocalMap map (v.vector->Map().NumGlobalElements64(),
+ Epetra_LocalMap map (n_global_elements(v.vector->Map()),
v.vector->Map().IndexBase(),
v.vector->Comm());
vector.reset (new Epetra_FEVector(map));
{
if (size() != v.size())
{
- Epetra_LocalMap map (v.vector->Map().NumGlobalElements64(),
+ Epetra_LocalMap map (n_global_elements(v.vector->Map()),
v.vector->Map().IndexBase(),
v.vector->Comm());
vector.reset (new Epetra_FEVector(map));
namespace TrilinosWrappers
{
+ namespace
+ {
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ // define a helper function that queries the global vector length of an
+ // Epetra_FEVector object by calling either the 32- or 64-bit
+ // function necessary.
+ int global_length(const Epetra_FEVector &vector)
+ {
+ return vector.GlobalLength();
+ }
+#else
+ // define a helper function that queries the global vector length of an
+ // Epetra_FEVector object by calling either the 32- or 64-bit
+ // function necessary.
+ long long int global_length(const Epetra_FEVector &vector)
+ {
+ return vector.GlobalLength64();
+ }
+#endif
+ }
+
namespace internal
{
VectorReference::operator TrilinosScalar () const
void
VectorBase::print (const char *format) const
{
- Assert (vector->GlobalLength64()!=0, ExcEmptyObject());
+ Assert (global_length(*vector)!=0, ExcEmptyObject());
for (size_type j=0; j<size(); ++j)
{