max_element (const TrilinosWrappers::Vector &criteria)
{
TrilinosScalar m = 0;
- criteria.vector->MaxValue(&m);
+ criteria.trilinos_vector().MaxValue(&m);
return m;
}
min_element (const TrilinosWrappers::Vector &criteria)
{
TrilinosScalar m = 0;
- criteria.vector->MinValue(&m);
+ criteria.trilinos_vector().MinValue(&m);
return m;
}
#endif
stokes_matrix (&S),
m_inverse (&Mpinv),
a_preconditioner (Apreconditioner),
- tmp (stokes_matrix->block(1,1).matrix->RowMap())
+ tmp (stokes_matrix->block(1,1).row_partitioner())
{}
std::vector<double> rhs_values(n_q_points);
- TrilinosWrappers::MPI::Vector rhs (temperature_mass_matrix.matrix->RowMap()),
- sol (temperature_mass_matrix.matrix->RowMap());
+ TrilinosWrappers::MPI::Vector rhs (temperature_mass_matrix.row_partitioner()),
+ sol (temperature_mass_matrix.row_partitioner());
for (; cell!=endc; ++cell)
if (cell->subdomain_id() ==
// deal.II wrapper class itself. So
// we access to the actual Trilinos
// matrix in the Trilinos wrapper
- // class, and create a plain pointer
- // when passing it in.
+ // class by the command
+ // trilinos_matrix(). Trilinos wants
+ // the matrix to be non-constant, so
+ // we have to manually remove the
+ // constantness using a const_cast.
case Parameters::Solver::gmres:
{
- Epetra_Vector x(View, system_matrix.matrix->RowMap(),
+ Epetra_Vector x(View, system_matrix.domain_partitioner(),
newton_update.begin());
- Epetra_Vector b(View, system_matrix.matrix->RowMap(),
+ Epetra_Vector b(View, system_matrix.range_partitioner(),
right_hand_side.begin());
AztecOO solver;
solver.SetAztecParam(AZ_athresh, parameters.ilut_atol);
solver.SetAztecParam(AZ_rthresh, parameters.ilut_rtol);
- solver.SetUserMatrix(&*system_matrix.matrix);
+ solver.SetUserMatrix(const_cast<Epetra_CrsMatrix*>
+ (&system_matrix.trilinos_matrix()));
solver.Iterate(parameters.max_iterations, parameters.linear_residual);
TrilinosScalar frobenius_norm () const;
//@}
+/**
+ * @name Access to underlying Trilinos data
+ */
+//@{
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos
+ * Epetra_CrsMatrix data.
+ */
+ const Epetra_CrsMatrix & trilinos_matrix () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos
+ * Epetra_CrsGraph data that stores
+ * the sparsity pattern of the
+ * matrix.
+ */
+ const Epetra_CrsGraph & trilinos_sparsity_pattern () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * domain space of this matrix, i.e.,
+ * the partitioning of the vectors
+ * this matrix has to be multiplied
+ * with.
+ */
+ const Epetra_Map & domain_partitioner () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * range space of this matrix, i.e.,
+ * the partitioning of the vectors
+ * that are result from matrix-vector
+ * products.
+ */
+ const Epetra_Map & range_partitioner () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * matrix rows. Equal to the
+ * partitioning of the range.
+ */
+ const Epetra_Map & row_partitioner () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * matrix columns. This is in general
+ * not equal to the partitioner
+ * Epetra_Map for the domain because
+ * of overlap in the matrix.
+ */
+ const Epetra_Map & col_partitioner () const;
+//@}
/**
* @name Iterators
*/
*/
std::vector<TrilinosScalar> column_values;
- public:
/**
* A sparse matrix object in
* Trilinos to be used for
* elements. The actual type,
* a sparse matrix, is set in
* the constructor.
- *
- * TODO: This object should
- * finally become private.
*/
std::auto_ptr<Epetra_FECrsMatrix> matrix;
}
+
+ inline
+ const Epetra_CrsMatrix &
+ SparseMatrix::trilinos_matrix () const
+ {
+ return static_cast<Epetra_CrsMatrix&>(*matrix);
+ }
+
+
+
+ inline
+ const Epetra_CrsGraph &
+ SparseMatrix::trilinos_sparsity_pattern () const
+ {
+ return matrix->Graph();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparseMatrix::domain_partitioner () const
+ {
+ return matrix->DomainMap();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparseMatrix::range_partitioner () const
+ {
+ return matrix->RangeMap();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparseMatrix::row_partitioner () const
+ {
+ return matrix->RowMap();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparseMatrix::col_partitioner () const
+ {
+ return matrix->ColMap();
+ }
+
#endif // DOXYGEN
}
const unsigned int n_cols,
const unsigned int *col_indices);
//@}
+/**
+ * @name Access of underlying Trilinos data
+ */
+//@{
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos
+ * Epetra_CrsGraph data that stores
+ * the sparsity pattern.
+ */
+ const Epetra_CrsGraph & trilinos_sparsity_pattern () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the parallel
+ * partitioning of the domain space
+ * of this sparsity pattern, i.e.,
+ * the partitioning of the vectors
+ * matrices based on this sparsity
+ * pattern are multiplied with.
+ */
+ const Epetra_Map & domain_partitioner () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * range space of this sparsity
+ * pattern, i.e., the partitioning of
+ * the vectors that are result from
+ * matrix-vector products.
+ */
+ const Epetra_Map & range_partitioner () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * sparsity pattern rows. Equal to
+ * the partitioning of the range.
+ */
+ const Epetra_Map & row_partitioner () const;
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the partitioning of the
+ * sparsity pattern columns. This is
+ * in general not equal to the
+ * partitioner Epetra_Map for the
+ * domain because of overlap in the
+ * matrix.
+ */
+ const Epetra_Map & col_partitioner () const;
+//@}
/**
* @name Iterators
*/
*/
unsigned int row_in_cache;
- friend class SparseMatrix;
friend class SparsityPatternIterators::const_iterator;
friend class SparsityPatternIterators::const_iterator::Accessor;
};
if (n_cols == 0)
return;
- int * col_index_ptr = (int*)col_indices;
+ int * col_index_ptr = (int*)(col_indices);
compressed = false;
int ierr;
// Epetra_CrsGraph input function, which
// is much faster than the
// Epetra_FECrsGraph function.
- if (row_map.MyGID(row) == true)
- ierr = graph->Epetra_CrsGraph::InsertGlobalIndices(row,
- n_cols,
- col_index_ptr);
- else
+ //if (row_map.MyGID(row) == true)
+ //ierr = graph->Epetra_CrsGraph::InsertGlobalIndices(row,
+ // n_cols,
+ // col_index_ptr);
+ //else
{
// When we're at off-processor data, we
// have to stick with the standard
- // SumIntoGlobalValues
- // function. Nevertheless, the way we
- // call it is the fastest one (any other
- // will lead to repeated allocation and
- // deallocation of memory in order to
- // call the function we already use,
- // which is very unefficient if writing
- // one element at a time).
-
- ierr = graph->InsertGlobalIndices (1, (int*)&row, n_cols,
- col_index_ptr);
+ // SumIntoGlobalValues method.
+
+ ierr = graph->InsertGlobalIndices (1, (int*)&row, n_cols, col_index_ptr);
}
//Assert (ierr <= 0, ExcAccessToNonPresentElement(row, col_index_ptr[0]));
}
+
+ inline
+ const Epetra_CrsGraph &
+ SparsityPattern::trilinos_sparsity_pattern () const
+ {
+ return static_cast<Epetra_CrsGraph&>(*graph);
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparsityPattern::domain_partitioner () const
+ {
+ return (Epetra_Map &) graph->DomainMap();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparsityPattern::range_partitioner () const
+ {
+ return (Epetra_Map &) graph->RangeMap();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparsityPattern::row_partitioner () const
+ {
+ return (Epetra_Map &) graph->RowMap();
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ SparsityPattern::col_partitioner () const
+ {
+ return (Epetra_Map &) graph->ColMap();
+ }
+
+
#endif // DOXYGEN
}
*/
//@{
+ /**
+ * Return a const reference to the
+ * underlying Trilinos
+ * Epetra_MultiVector class.
+ */
+ const Epetra_MultiVector & trilinos_vector () const;
+
+ /**
+ * Return a (modifyable) reference to
+ * the underlying Trilinos
+ * Epetra_FEVector class.
+ */
+ Epetra_FEVector & trilinos_vector ();
+
+ /**
+ * Return a const reference to the
+ * underlying Trilinos Epetra_Map
+ * that sets the parallel
+ * partitioning of the vector.
+ */
+ const Epetra_Map & vector_partitioner () const;
+
/**
* Output of vector in
* user-defined format in analogy
*/
bool compressed;
- public:
/**
* An Epetra distibuted vector
* type. Requires an existing
* Epetra_Map for storing data.
- * TODO: Should become private
- * at some point.
*/
std::auto_ptr<Epetra_FEVector> vector;
}
+
+ inline
+ const Epetra_MultiVector &
+ VectorBase::trilinos_vector () const
+ {
+ return static_cast<Epetra_MultiVector &>(*vector);
+ }
+
+
+
+ inline
+ Epetra_FEVector &
+ VectorBase::trilinos_vector ()
+ {
+ return *vector;
+ }
+
+
+
+ inline
+ const Epetra_Map &
+ VectorBase::vector_partitioner () const
+ {
+ return (Epetra_Map &)vector->Map();
+ }
+
#endif // DOXYGEN
}
// and copy it
TrilinosScalar **start_ptr;
- int ierr = localized_vector.vector->ExtractView (&start_ptr);
+ int ierr = localized_vector.trilinos_vector().ExtractView (&start_ptr);
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
std::copy (start_ptr[0], start_ptr[0]+vec_size, begin());
// and copy it
TrilinosScalar **start_ptr;
- int ierr = v.vector->ExtractView (&start_ptr);
+ int ierr = v.trilinos_vector().ExtractView (&start_ptr);
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
std::copy (start_ptr[0], start_ptr[0]+vec_size, begin());
// get a representation of the vector
// and copy it
TrilinosScalar **start_ptr;
- int ierr = v.vector->ExtractView (&start_ptr);
+ int ierr = v.trilinos_vector().ExtractView (&start_ptr);
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
std::copy (start_ptr[0], start_ptr[0]+vec_size, begin());
PreconditionBase::vmult (VectorBase &dst,
const VectorBase &src) const
{
- Assert (dst.vector->Map().SameAs(preconditioner->OperatorRangeMap()),
+ Assert (dst.vector_partitioner().SameAs(preconditioner->OperatorRangeMap()),
ExcNonMatchingMaps("dst"));
- Assert (src.vector->Map().SameAs(preconditioner->OperatorDomainMap()),
+ Assert (src.vector_partitioner().SameAs(preconditioner->OperatorDomainMap()),
ExcNonMatchingMaps("src"));
- const int ierr = preconditioner->ApplyInverse (*src.vector, *dst.vector);
+ const int ierr = preconditioner->ApplyInverse (src.trilinos_vector(),
+ dst.trilinos_vector());
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
}
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("point relaxation",
- &*matrix.matrix, 0));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("point relaxation",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ 0));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("point relaxation",
- &*matrix.matrix,
- additional_data.overlap));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("point relaxation",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ additional_data.overlap));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("point relaxation",
- &*matrix.matrix,
- additional_data.overlap));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("point relaxation",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ additional_data.overlap));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("IC", &*matrix.matrix,
- additional_data.overlap));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("IC",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ additional_data.overlap));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("ILU", &*matrix.matrix,
- additional_data.overlap));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("ILU",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ additional_data.overlap));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("ILUT", &*matrix.matrix,
- additional_data.overlap));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("ILUT",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ additional_data.overlap));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (Ifpack().Create ("Amesos", &*matrix.matrix,
- additional_data.overlap));
+ ifpack = Teuchos::rcp (Ifpack().Create
+ ("Amesos",
+ const_cast<Epetra_CrsMatrix*>(&matrix.trilinos_matrix()),
+ additional_data.overlap));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
preconditioner.release();
ifpack.release();
- ifpack = Teuchos::rcp (new Ifpack_Chebyshev(&*matrix.matrix));
+ ifpack = Teuchos::rcp (new Ifpack_Chebyshev (&matrix.trilinos_matrix()));
Assert (&*ifpack != 0, ExcMessage ("Trilinos could not create this "
"preconditioner"));
else
parameter_list.set("ML output", 0);
- const Epetra_Map & domain_map = matrix.matrix->DomainMap();
+ const Epetra_Map & domain_map = matrix.domain_partitioner();
Epetra_MultiVector distributed_constant_modes (domain_map,
constant_modes_dimension);
}
multilevel_operator = Teuchos::rcp (new ML_Epetra::MultiLevelPreconditioner(
- *matrix.matrix, parameter_list, true));
+ matrix.trilinos_matrix(), parameter_list,
+ true));
if (additional_data.output_details)
multilevel_operator->PrintUnused(0);
// (Thyra wrapper around the
// matrix).
Teuchos::RCP<const Thyra::LinearOpBase<double> > tmpM =
- Thyra::epetraLinearOp(Teuchos::rcp(&*input_M.matrix, false));
+ Thyra::epetraLinearOp(Teuchos::rcp(&input_M.trilinos_matrix(), false));
Teuchos::RCP<const Thyra::LinearOpSourceBase<double> > M
= Teuchos::rcp(new Thyra::DefaultLinearOpSource<double>(tmpM), true);
else
{
Teuchos::RCP<const Thyra::LinearOpBase<double> >
- tmpFp = Thyra::epetraLinearOp(Teuchos::rcp(&*Fp_matrix.matrix,false));
+ tmpFp = Thyra::epetraLinearOp(Teuchos::rcp(&Fp_matrix.trilinos_matrix(),
+ false));
Fp_op = tmpFp;
}
// gradient and divergence
// operators in Thyra format.
Teuchos::RCP<const Thyra::LinearOpBase<double> > tmpS01 =
- Thyra::epetraLinearOp(Teuchos::rcp(&*system_matrix.block(0,1).matrix,
+ Thyra::epetraLinearOp(Teuchos::rcp(&system_matrix.block(0,1).trilinos_matrix(),
false));
Thyra::ConstLinearOperator<double> S01 = tmpS01;
Teuchos::RCP<const Thyra::LinearOpBase<double> > tmpS10 =
- Thyra::epetraLinearOp(Teuchos::rcp(&*system_matrix.block(1,0).matrix,
+ Thyra::epetraLinearOp(Teuchos::rcp(&system_matrix.block(1,0).trilinos_matrix(),
false));
Thyra::ConstLinearOperator<double> S10 = tmpS10;
// to let the AztecOO solver
// know about the matrix and
// vectors.
- linear_problem = std::auto_ptr<Epetra_LinearProblem> (
- new Epetra_LinearProblem(&*(A.matrix), &*x.vector,
- &*b.vector));
+ linear_problem = std::auto_ptr<Epetra_LinearProblem>
+ (new Epetra_LinearProblem(const_cast<Epetra_CrsMatrix*>(&A.trilinos_matrix()),
+ &x.trilinos_vector(),
+ const_cast<Epetra_MultiVector*>(&b.trilinos_vector())));
// Next we can allocate the
// AztecOO solver...
ExcDimensionMismatch(b.size(), A.m()));
Assert (A.local_range ().second == A.m(),
ExcMessage ("Can only work in serial when using deal.II vectors."));
- Assert (A.matrix->Filled(),
+ Assert (A.trilinos_matrix().Filled(),
ExcMessage ("Matrix is not compressed. Call compress() method."));
- Epetra_Vector ep_x (View, A.matrix->DomainMap(), x.begin());
- Epetra_Vector ep_b (View, A.matrix->RangeMap(), const_cast<double*>(b.begin()));
+ Epetra_Vector ep_x (View, A.domain_partitioner(), x.begin());
+ Epetra_Vector ep_b (View, A.range_partitioner(), const_cast<double*>(b.begin()));
// We need an
// Epetra_LinearProblem object
// know about the matrix and
// vectors.
linear_problem = std::auto_ptr<Epetra_LinearProblem>
- (new Epetra_LinearProblem(&*(A.matrix), &ep_x, &ep_b));
+ (new Epetra_LinearProblem
+ (const_cast<Epetra_CrsMatrix*>(&A.trilinos_matrix()), &ep_x, &ep_b));
// Next we can allocate the
// AztecOO solver...
// to let the AztecOO solver
// know about the matrix and
// vectors.
- linear_problem = std::auto_ptr<Epetra_LinearProblem> (
- new Epetra_LinearProblem(&*(A.matrix), &*x.vector,
- &*b.vector));
+ linear_problem = std::auto_ptr<Epetra_LinearProblem>
+ (new Epetra_LinearProblem(const_cast<Epetra_CrsMatrix*>(&A.trilinos_matrix()),
+ &x.trilinos_vector(),
+ const_cast<Epetra_MultiVector*>(&b.trilinos_vector())));
// Next we can allocate the
// AztecOO solver...
ExcDimensionMismatch(b.size(), A.m()));
Assert (A.local_range ().second == A.m(),
ExcMessage ("Can only work in serial when using deal.II vectors."));
- Epetra_Vector ep_x (View, A.matrix->DomainMap(), x.begin());
- Epetra_Vector ep_b (View, A.matrix->RangeMap(), const_cast<double*>(b.begin()));
+ Epetra_Vector ep_x (View, A.domain_partitioner(), x.begin());
+ Epetra_Vector ep_b (View, A.range_partitioner(), const_cast<double*>(b.begin()));
// We need an
// Epetra_LinearProblem object
// know about the matrix and
// vectors.
linear_problem = std::auto_ptr<Epetra_LinearProblem>
- (new Epetra_LinearProblem(&*(A.matrix), &ep_x, &ep_b));
+ (new Epetra_LinearProblem
+ (const_cast<Epetra_CrsMatrix*>(&A.trilinos_matrix()), &ep_x, &ep_b));
// Next we can allocate the
// AztecOO solver...
MPI::BlockVector tmpb (b.n_blocks());
for (unsigned int i=0; i<x.n_blocks(); ++i)
{
- tmpx.block(i).reinit (A.block(i,i).matrix->DomainMap(), false);
+ tmpx.block(i).reinit (A.block(i,i).domain_partitioner(), false);
tmpx.block(i) = x.block(i);
- tmpb.block(i).reinit (A.block(i,i).matrix->RangeMap(), false);
+ tmpb.block(i).reinit (A.block(i,i).range_partitioner(), false);
tmpb.block(i) = b.block(i);
}
tmpx.collect_sizes();
for (unsigned int i=0; i<n_rows; ++i)
{
- AssertThrow (input_x.block(i).vector->Map().UniqueGIDs() == true,
+ AssertThrow (input_x.block(i).vector_partitioner().UniqueGIDs() == true,
ExcOverlappingMaps("vector", "x"));
- AssertThrow (input_b.block(i).vector->Map().UniqueGIDs() == true,
+ AssertThrow (input_b.block(i).vector_partitioner().UniqueGIDs() == true,
ExcOverlappingMaps("vector", "b"));
for (unsigned int j=0; j<n_cols; ++j)
{
i_str << i;
std::ostringstream error_component;
error_component << "x.block(" << i_str.str() << ")";
- AssertThrow (input_x.block(j).vector->Map().SameAs(
- input_A.block(i,j).matrix->DomainMap()) == true,
+ AssertThrow (input_x.block(j).vector_partitioner().SameAs(
+ input_A.block(i,j).domain_partitioner()) == true,
ExcNonMatchingMaps (error_component.str()));
}
{
i_str << j;
std::ostringstream error_component;
error_component << "b.block(" << i_str.str() << ")";
- AssertThrow (input_b.block(i).vector->Map().SameAs(
- input_A.block(i,j).matrix->RangeMap()) == true,
+ AssertThrow (input_b.block(i).vector_partitioner().SameAs(
+ input_A.block(i,j).range_partitioner()) == true,
ExcNonMatchingMaps (error_component.str()));
}
}
{
Teuchos::RCP<const Thyra::VectorSpaceBase<double> > tmp_space
= Thyra::create_VectorSpace(
- Teuchos::rcp(&input_A.block(i,i).matrix->DomainMap(),
- false));
+ Teuchos::rcp(&input_A.block(i,i).domain_partitioner(),
+ false));
epetra_vector_spaces.push_back(tmp_space);
}
for (unsigned int j=0; j<n_cols; ++j)
{
Teuchos::RCP<const Thyra::LinearOpBase<double> >
- A_ij = Thyra::epetraLinearOp(Teuchos::rcp(&*input_A.block(i,j).matrix,
+ A_ij = Thyra::epetraLinearOp(Teuchos::rcp(&input_A.block(i,j).trilinos_matrix(),
false));
tmpA->setBlock(i, j, A_ij);
}
Thyra::Vector<double> sol = A.domain().createMember();
for (unsigned int i=0; i<n_rows; ++i)
{
- Epetra_Vector *block_ptr = (*input_b.block(i).vector)(0);
+ Epetra_Vector *block_ptr = const_cast<Epetra_Vector*>
+ ((input_b.block(i).trilinos_vector())(0));
Teuchos::RCP<Thyra::VectorBase<double> > tmp_rhs_i
= Thyra::create_Vector(Teuchos::rcp(block_ptr, false),
epetra_vector_spaces[i]);
const Thyra::Vector<double> rhs_i = tmp_rhs_i;
rhs.setBlock(i, rhs_i);
- block_ptr = (*input_x.block(i).vector)(0);
+ block_ptr = (input_x.block(i).trilinos_vector())(0);
Teuchos::RCP<Thyra::VectorBase<double> > tmp_sol_i
= Thyra::create_Vector(Teuchos::rcp(block_ptr, false),
epetra_vector_spaces[i]);
TrilinosScalar *values = new TrilinosScalar(colnums);
int ierr;
- ierr = matrix->matrix->ExtractGlobalRowCopy((int)this->a_row, colnums,
- ncols, &(values[0]));
+ ierr = matrix->trilinos_matrix().ExtractGlobalRowCopy((int)this->a_row,
+ colnums,
+ ncols, &(values[0]));
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
// copy it into our caches if the
SparseMatrix::SparseMatrix (const SparsityPattern &InputSP)
:
Subscriptor(),
- row_map (InputSP.row_map),
- col_map (InputSP.col_map),
+ row_map (InputSP.range_partitioner()),
+ col_map (InputSP.domain_partitioner()),
last_action (Zero),
compressed (true),
matrix (std::auto_ptr<Epetra_FECrsMatrix>
- (new Epetra_FECrsMatrix(Copy, *InputSP.graph, false)))
+ (new Epetra_FECrsMatrix(Copy,
+ InputSP.trilinos_sparsity_pattern(),
+ false)))
{
- Assert(InputSP.graph->Filled() == true,
+ Assert(InputSP.trilinos_sparsity_pattern().Filled() == true,
ExcMessage("The Trilinos sparsity pattern has not been compressed."));
compress();
}
{
matrix.reset();
- row_map = sparsity_pattern.row_map;
- col_map = sparsity_pattern.col_map;
+ row_map = sparsity_pattern.range_partitioner();
+ col_map = sparsity_pattern.domain_partitioner();
- Assert (sparsity_pattern.graph->Filled() == true,
+ Assert (sparsity_pattern.trilinos_sparsity_pattern().Filled() == true,
ExcMessage("The Trilinos sparsity pattern has not been compressed"));
matrix = std::auto_ptr<Epetra_FECrsMatrix>
- (new Epetra_FECrsMatrix(Copy, *sparsity_pattern.graph, false));
+ (new Epetra_FECrsMatrix(Copy, sparsity_pattern.trilinos_sparsity_pattern(),
+ false));
compress();
}
row_map = sparse_matrix.row_map;
col_map = sparse_matrix.col_map;
matrix = std::auto_ptr<Epetra_FECrsMatrix>
- (new Epetra_FECrsMatrix(Copy, sparse_matrix.matrix->Graph(), false));
+ (new Epetra_FECrsMatrix(Copy, sparse_matrix.trilinos_sparsity_pattern(),
+ false));
compress();
}
if (matrix->Filled() == false)
matrix->FillComplete(col_map, row_map, true);
- Assert (src.vector->Map().SameAs(matrix->DomainMap()) == true,
+ Assert (src.vector_partitioner().SameAs(matrix->DomainMap()) == true,
ExcMessage ("Column map of matrix does not fit with vector map!"));
- Assert (dst.vector->Map().SameAs(matrix->RangeMap()) == true,
+ Assert (dst.vector_partitioner().SameAs(matrix->RangeMap()) == true,
ExcMessage ("Row map of matrix does not fit with vector map!"));
- const int ierr = matrix->Multiply (false, *(src.vector), *(dst.vector));
+ const int ierr = matrix->Multiply (false, src.trilinos_vector(),
+ dst.trilinos_vector());
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
}
if (matrix->Filled() == false)
matrix->FillComplete(col_map, row_map, true);
- Assert (src.vector->Map().SameAs(matrix->RangeMap()) == true,
+ Assert (src.vector_partitioner().SameAs(matrix->RangeMap()) == true,
ExcMessage ("Column map of matrix does not fit with vector map!"));
- Assert (dst.vector->Map().SameAs(matrix->DomainMap()) == true,
+ Assert (dst.vector_partitioner().SameAs(matrix->DomainMap()) == true,
ExcMessage ("Row map of matrix does not fit with vector map!"));
- const int ierr = matrix->Multiply (true, *(src.vector), *(dst.vector));
+ const int ierr = matrix->Multiply (true, src.trilinos_vector(),
+ dst.trilinos_vector());
AssertThrow (ierr == 0, ExcTrilinosError(ierr));
}
Vector::import_nonlocal_data_for_fe (const TrilinosWrappers::SparseMatrix &m,
const Vector &v)
{
- Assert (m.matrix->Filled() == true,
+ Assert (m.trilinos_matrix().Filled() == true,
ExcMessage ("Matrix is not compressed. "
"Cannot find exchange information!"));
Assert (v.vector->Map().UniqueGIDs() == true,
ExcMessage ("The input vector has overlapping data, "
"which is not allowed."));
- if (vector->Map().SameAs(m.matrix->ColMap()) == false)
+ if (vector->Map().SameAs(m.col_partitioner()) == false)
{
- map = m.matrix->ColMap();
+ map = m.col_partitioner();
vector = std::auto_ptr<Epetra_FEVector> (new Epetra_FEVector(map));
}
projects / dealii-svn.git / commitdiff