Assert (&src != &dst, ExcSourceEqualsDestination());
Assert (matrix->Filled(), ExcMatrixNotCompressed());
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ const int n_global_elements = matrix->DomainMap().NumGlobalElements();
+#else
+ const long long int n_global_elements = matrix->DomainMap().NumGlobalElements64();
+#endif
+
AssertDimension (static_cast<size_type>(matrix->DomainMap().NumMyElements()),
- static_cast<size_type>(matrix->DomainMap().NumGlobalElements()));
+ static_cast<size_type>(n_global_elements));
AssertDimension (dst.size(), static_cast<size_type>(matrix->RangeMap().NumMyElements()));
AssertDimension (src.size(), static_cast<size_type>(matrix->DomainMap().NumMyElements()));
Assert (&src != &dst, ExcSourceEqualsDestination());
Assert (matrix->Filled(), ExcMatrixNotCompressed());
+#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
+ const int n_global_elements = matrix->DomainMap().NumGlobalElements();
+#else
+ const long long int n_global_elements = matrix->DomainMap().NumGlobalElements64();
+#endif
+
AssertDimension (static_cast<size_type>(matrix->DomainMap().NumMyElements()),
- static_cast<size_type>(matrix->DomainMap().NumGlobalElements()));
+ static_cast<size_type>(n_global_elements));
AssertDimension (dst.size(), static_cast<size_type>(matrix->DomainMap().NumMyElements()));
AssertDimension (src.size(), static_cast<size_type>(matrix->RangeMap().NumMyElements()));
// $Id$
// Version: $Name$
//
-// Copyright (C) 2008, 2009, 2011, 2012 by the deal.II authors
+// Copyright (C) 2008, 2009, 2011, 2012, 2013 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
namespace TrilinosWrappers
{
+ namespace
+ {
+ // 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();
+ }
+#else
+ long long int n_global_elements (const Epetra_BlockMap &map)
+ {
+ return map.NumGlobalElements64();
+ }
+#endif
+ }
+
+
namespace MPI
{
BlockVector &
}
-
+
BlockVector &
BlockVector::operator = (const BlockVector &v)
{
for (size_type i=0; i<no_blocks; ++i)
{
- block_sizes[i] = input_maps[i].NumGlobalElements();
+ block_sizes[i] = n_global_elements(input_maps[i]);
}
this->block_indices.reinit (block_sizes);
std::vector<size_type> block_sizes (no_blocks);
for (size_type i=0; i<no_blocks; ++i)
- block_sizes[i] = input_maps[i].NumGlobalElements();
+ block_sizes[i] = n_global_elements(input_maps[i]);
this->block_indices.reinit (block_sizes);
namespace TrilinosWrappers
{
+ namespace
+ {
+ // 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();
+ }
+#else
+ long long int n_global_elements (const Epetra_BlockMap &map)
+ {
+ return map.NumGlobalElements64();
+ }
+#endif
+ }
+
+
namespace SparseMatrixIterators
{
void
if (input_row_map.Comm().MyPID() == 0)
{
AssertDimension (sparsity_pattern.n_rows(),
- static_cast<size_type>(input_row_map.NumGlobalElements()));
+ static_cast<size_type>(n_global_elements(input_row_map)));
AssertDimension (sparsity_pattern.n_cols(),
- static_cast<size_type>(input_col_map.NumGlobalElements()));
+ static_cast<size_type>(n_global_elements(input_col_map)));
}
column_space_map.reset (new Epetra_Map (input_col_map));
const size_type n_rows = dealii_sparse_matrix.m();
-#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
- Assert (static_cast<size_type>(input_row_map.NumGlobalElements()) == n_rows,
- ExcDimensionMismatch (input_row_map.NumGlobalElements(),
+ Assert (static_cast<size_type>(n_global_elements(input_row_map)) == n_rows,
+ ExcDimensionMismatch (n_global_elements(input_row_map),
n_rows));
-#else
- Assert (static_cast<size_type>(input_row_map.NumGlobalElements64()) == n_rows,
- ExcDimensionMismatch (input_row_map.NumGlobalElements64(),
- n_rows));
-#endif
- Assert (input_row_map.NumGlobalElements() == (TrilinosWrappers::types::int_type)n_rows,
- ExcDimensionMismatch (input_row_map.NumGlobalElements(),
+ Assert (n_global_elements(input_row_map) == (TrilinosWrappers::types::int_type)n_rows,
+ ExcDimensionMismatch (n_global_elements(input_row_map),
n_rows));
- Assert (input_col_map.NumGlobalElements() == (TrilinosWrappers::types::int_type)dealii_sparse_matrix.n(),
- ExcDimensionMismatch (input_col_map.NumGlobalElements(),
+ Assert (n_global_elements(input_col_map) == (TrilinosWrappers::types::int_type)dealii_sparse_matrix.n(),
+ ExcDimensionMismatch (n_global_elements(input_col_map),
dealii_sparse_matrix.n()));
const ::dealii::SparsityPattern &sparsity_pattern =
namespace TrilinosWrappers
{
+ namespace
+ {
+ // 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();
+ }
+#else
+ long long int n_global_elements (const Epetra_BlockMap &map)
+ {
+ return map.NumGlobalElements64();
+ }
+#endif
+ }
+
+
namespace SparsityPatternIterators
{
void
// release memory before reallocation
graph.reset ();
AssertDimension (n_entries_per_row.size(),
- static_cast<size_type>(input_row_map.NumGlobalElements()));
+ 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.NumGlobalElements()));
+ static_cast<size_type>(n_global_elements(input_row_map)));
AssertDimension (sp.n_cols(),
- static_cast<size_type>(input_col_map.NumGlobalElements()));
+ static_cast<size_type>(n_global_elements(input_col_map)));
column_space_map.reset (new Epetra_Map (input_col_map));
compressed = false;
// Generate the view and make
// sure that we have not generated
// an error.
- int ierr = graph->ExtractMyRowView(static_cast<int>(trilinos_i),
+ int ierr = graph->ExtractMyRowView(static_cast<int>(trilinos_i),
nnz_extracted, col_indices);
Assert (ierr==0, ExcTrilinosError(ierr));
ExcDimensionMismatch(nnz_present, nnz_extracted));
// Search the index
- int *el_find = std::find(col_indices, col_indices + nnz_present,
+ int *el_find = std::find(col_indices, col_indices + nnz_present,
static_cast<int>(trilinos_j));
int local_col_index = (int)(el_find - col_indices);
if (graph->Filled() == true)
n_cols = graph -> NumGlobalCols();
else
- n_cols = column_space_map->NumGlobalElements();
+ n_cols = n_global_elements(*column_space_map);
return n_cols;
}
namespace TrilinosWrappers
{
+ namespace
+ {
+ // 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();
+ }
+#else
+ long long int n_global_elements (const Epetra_BlockMap &map)
+ {
+ return map.NumGlobalElements64();
+ }
+#endif
+ }
+
namespace MPI
{
:
VectorBase()
{
- AssertThrow (input_map.NumGlobalElements() == v.vector->Map().NumGlobalElements(),
- ExcDimensionMismatch (input_map.NumGlobalElements(),
- v.vector->Map().NumGlobalElements()));
+ 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 ((int)parallel_partitioner.size() == v.vector->Map().NumGlobalElements(),
+ AssertThrow ((size_type)parallel_partitioner.size() == n_global_elements(v.vector->Map()),
ExcDimensionMismatch (parallel_partitioner.size(),
- v.vector->Map().NumGlobalElements()));
+ n_global_elements(v.vector->Map())));
last_action = Zero;
Vector::Vector (const Epetra_Map &input_map)
{
last_action = Zero;
- Epetra_LocalMap map (input_map.NumGlobalElements(),
+ 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().NumGlobalElements(),
+ 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().NumGlobalElements() != input_map.NumGlobalElements())
+ if (n_global_elements(vector->Map()) != n_global_elements(input_map))
{
- Epetra_LocalMap map (input_map.NumGlobalElements(),
+ 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().NumGlobalElements() !=
+ 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 (size() != v.size())
{
- Epetra_LocalMap map (v.vector->Map().NumGlobalElements(),
+ 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())
{
-#ifndef DEAL_II_USE_LARGE_INDEX_TYPE
- const int n_global_elements
- = v.vector->Map().NumGlobalElements();
-#else
- const long long int n_global_elements
- = v.vector->Map().NumGlobalElements64();
-#endif
- Epetra_LocalMap map (n_global_elements,
+ Epetra_LocalMap map (n_global_elements(v.vector->Map()),
v.vector->Map().IndexBase(),
v.vector->Comm());
vector.reset (new Epetra_FEVector(map));