From f53ef9ce461049ae4ca02ee43eb06a04304af70c Mon Sep 17 00:00:00 2001 From: kronbichler Date: Wed, 2 Mar 2011 17:10:07 +0000 Subject: [PATCH] Cleanup: replace Teuchos::rcp by shared_ptr in one place to be consistent with the rest of Trilinos part. Move one Trilinos header file to .cc file. git-svn-id: https://svn.dealii.org/trunk@23457 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/base/tensor.h | 54 ++++++++++++------------ deal.II/include/deal.II/base/utilities.h | 1 - deal.II/source/base/utilities.cc | 1 + deal.II/source/lac/trilinos_vector.cc | 17 +++++--- 4 files changed, 38 insertions(+), 35 deletions(-) diff --git a/deal.II/include/deal.II/base/tensor.h b/deal.II/include/deal.II/base/tensor.h index 448d7c823c..16fd1ac05c 100644 --- a/deal.II/include/deal.II/base/tensor.h +++ b/deal.II/include/deal.II/base/tensor.h @@ -65,7 +65,7 @@ class Tensor * the outside world. */ static const unsigned int rank = rank_; - + /** * Type of stored objects. This * is a tensor of lower rank. @@ -85,13 +85,13 @@ class Tensor * to zero. */ Tensor (); - + /** * Copy constructor, where the data is * copied from a C-style array. */ Tensor (const array_type &initializer); - + /** * Read-Write access operator. */ @@ -135,7 +135,7 @@ class Tensor * Add another tensor. */ Tensor & operator += (const Tensor &); - + /** * Subtract another tensor. */ @@ -196,7 +196,7 @@ class Tensor * contexts. */ double norm_square () const; - + /** * Fill a vector with all tensor elements. * @@ -244,9 +244,9 @@ class Tensor << "Invalid tensor index " << arg1); /** - * Read or write the data of this object to or + * Read or write the data of this object to or * from a stream for the purpose of serialization - */ + */ template void serialize(Archive & ar, const unsigned int version); @@ -306,7 +306,7 @@ typename Tensor::value_type& Tensor::operator[] (const unsigned int i) { Assert (i::value_type& Tensor::operator[] (const unsigned int i) const { Assert (i Tensor::operator + (const Tensor &t) const { Tensor tmp(*this); - + for (unsigned int i=0; i Tensor::operator - (const Tensor &t) const { Tensor tmp(*this); - + for (unsigned int i=0; i Tensor::operator - () const { Tensor tmp; - + for (unsigned int i=0; i::memory_consumption () template template inline -void +void Tensor::serialize(Archive & ar, const unsigned int) { ar & subtensor; -} +} #endif // DOXYGEN /* ----------------- Non-member functions operating on tensors. ------------ */ @@ -924,7 +924,7 @@ void contract (Tensor<3,dim> &dest, Assert (false, (typename Tensor<2,dim>::ExcInvalidTensorIndex (index2))); } - + break; case 2: switch (index2) @@ -947,7 +947,7 @@ void contract (Tensor<3,dim> &dest, Assert (false, (typename Tensor<2,dim>::ExcInvalidTensorIndex (index2))); } - + break; case 3: switch (index2) @@ -970,7 +970,7 @@ void contract (Tensor<3,dim> &dest, Assert (false, (typename Tensor<2,dim>::ExcInvalidTensorIndex (index2))); } - + break; default: Assert (false, @@ -978,7 +978,7 @@ void contract (Tensor<3,dim> &dest, } } - + /** * Multiplication operator performing a contraction of the last index * of the first argument and the first index of the second @@ -1322,7 +1322,7 @@ cross_product (Tensor<1,dim> &dst, const Tensor<1,dim> &src) { Assert (dim==2, ExcInternalError()); - + dst[0] = src[1]; dst[1] = -src[0]; } @@ -1346,7 +1346,7 @@ cross_product (Tensor<1,dim> &dst, const Tensor<1,dim> &src2) { Assert (dim==3, ExcInternalError()); - + dst[0] = src1[1]*src2[2] - src1[2]*src2[1]; dst[1] = src1[2]*src2[0] - src1[0]*src2[2]; dst[2] = src1[0]*src2[1] - src1[1]*src2[0]; @@ -1483,7 +1483,7 @@ double determinant (const Tensor<2,dim> &t) // for some algorithmic simplicity, we use // the expansion along the last row double det = 0; - + for (unsigned int k=0; k minor; @@ -1496,7 +1496,7 @@ double determinant (const Tensor<2,dim> &t) det += t[dim-1][k] * cofactor; } - + return std::pow (-1., static_cast(dim)) * det; } @@ -1535,7 +1535,7 @@ Tensor<2,dim> invert (const Tensor<2,dim> &t) { Tensor<2,dim> return_tensor; - switch (dim) + switch (dim) { case 1: return_tensor[0][0] = 1.0/t[0][0]; @@ -1552,7 +1552,7 @@ invert (const Tensor<2,dim> &t) return_tensor[1][1] = t[0][0]*t4; break; } - + case 3: { const double t4 = t[0][0]*t[1][1], @@ -1572,17 +1572,17 @@ invert (const Tensor<2,dim> &t) return_tensor[2][0] = -(-t[1][0]*t[2][1]+t[1][1]*t[2][0])*t07; return_tensor[2][1] = -(t[0][0]*t[2][1]-t01)*t07; return_tensor[2][2] = (t4-t8)*t07; - + break; } // if desired, take over the // inversion of a 4x4 tensor // from the FullMatrix - + default: AssertThrow (false, ExcNotImplemented()); - } + } return return_tensor; } diff --git a/deal.II/include/deal.II/base/utilities.h b/deal.II/include/deal.II/base/utilities.h index 6d8c067a0d..f99177a858 100644 --- a/deal.II/include/deal.II/base/utilities.h +++ b/deal.II/include/deal.II/base/utilities.h @@ -27,7 +27,6 @@ typedef int MPI_Comm; #endif #ifdef DEAL_II_USE_TRILINOS -# include # include # include # ifdef DEAL_II_COMPILER_SUPPORTS_MPI diff --git a/deal.II/source/base/utilities.cc b/deal.II/source/base/utilities.cc index e96654b336..d899364ce8 100644 --- a/deal.II/source/base/utilities.cc +++ b/deal.II/source/base/utilities.cc @@ -39,6 +39,7 @@ # include # include # endif +# include "Teuchos_RCP.hpp" # include "Epetra_SerialComm.h" #endif diff --git a/deal.II/source/lac/trilinos_vector.cc b/deal.II/source/lac/trilinos_vector.cc index a2ecd7ec61..a165a83680 100644 --- a/deal.II/source/lac/trilinos_vector.cc +++ b/deal.II/source/lac/trilinos_vector.cc @@ -18,10 +18,11 @@ # include # include -# include # include # include +# include + DEAL_II_NAMESPACE_OPEN @@ -240,12 +241,14 @@ namespace TrilinosWrappers Epetra_Map new_map (v.size(), n_elements, &global_ids[0], 0, v.block(0).vector_partitioner().Comm()); - if (import_data == false) - vector.reset (new Epetra_FEVector (new_map)); - - Teuchos::RCP actual_vec = (import_data == true) ? - Teuchos::rcp (new Epetra_FEVector (new_map), true) : - Teuchos::rcp (vector.get(), false); + std_cxx1x::shared_ptr actual_vec; + if ( import_data == true ) + actual_vec.reset (new Epetra_FEVector (new_map)); + else + { + vector.reset (new Epetra_FEVector (new_map)); + actual_vec = vector; + } TrilinosScalar* entries = (*actual_vec)[0]; block_offset = 0; -- 2.39.5