46 #ifndef XPETRA_TPETRAVECTOR_DEF_HPP
47 #define XPETRA_TPETRAVECTOR_DEF_HPP
55 #include "Xpetra_TpetraMultiVector.hpp"
57 #include "Tpetra_MultiVector.hpp"
58 #include "Tpetra_Vector.hpp"
62 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
67 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
69 :
TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> (map,A,map->getNodeNumElements(),1)
72 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
77 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
79 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::replaceGlobalValue"); getTpetra_Vector()->replaceGlobalValue(globalRow, value); }
81 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
83 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::sumIntoGlobalValue"); getTpetra_Vector()->sumIntoGlobalValue(globalRow, value); }
85 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
87 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::replaceLocalValue"); getTpetra_Vector()->replaceLocalValue(myRow, value); }
89 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
91 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::sumIntoLocalValue"); getTpetra_Vector()->sumIntoLocalValue(myRow, value); }
93 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
95 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::norm1");
return getTpetra_Vector()->norm1(); }
97 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
99 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::norm2");
return getTpetra_Vector()->norm2(); }
101 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
103 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::normInf");
return getTpetra_Vector()->normInf(); }
105 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
107 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::meanValue");
return getTpetra_Vector()->meanValue(); }
109 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
111 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::description");
return getTpetra_Vector()->description(); }
113 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
115 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::describe"); getTpetra_Vector()->describe(out, verbLevel); }
117 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
119 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::dot");
return getTpetra_Vector()->dot(*
toTpetra(a)); }
121 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
123 {
XPETRA_MONITOR(
"TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::normWeighted");
return getTpetra_Vector()->normWeighted(*
toTpetra(weights)); }
125 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
129 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
130 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
134 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
136 return this->TpetraMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getHostLocalView();
139 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
141 return this->TpetraMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getDeviceLocalView();
145 template<
class TargetDeviceType>
146 typename Kokkos::Impl::if_c<
147 Kokkos::Impl::is_same<
148 typename dual_view_type::t_dev_um::execution_space::memory_space,
149 typename TargetDeviceType::memory_space>::value,
150 typename dual_view_type::t_dev_um,
151 typename dual_view_type::t_host_um>::type
152 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
153 TpetraVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getLocalView ()
const {
154 return this->TpetraMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::template getLocalView<TargetDeviceType>();
161 #ifdef HAVE_XPETRA_EPETRA
163 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
164 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
167 template <
class Scalar>
169 :
public virtual Vector<Scalar,int,int,EpetraNode>,
176 #undef XPETRA_TPETRAMULTIVECTOR_SHORT
274 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
278 typename dual_view_type::t_host_um getHostLocalView ()
const {
279 typename dual_view_type::t_host_um ret;
283 typename dual_view_type::t_dev_um getDeviceLocalView()
const {
284 typename dual_view_type::t_dev_um ret;
298 template<
class TargetDeviceType>
299 typename Kokkos::Impl::if_c<
300 Kokkos::Impl::is_same<
301 typename dual_view_type::t_dev_um::execution_space::memory_space,
302 typename TargetDeviceType::memory_space>::value,
303 typename dual_view_type::t_dev_um,
304 typename dual_view_type::t_host_um>::type
305 getLocalView ()
const {
306 typename Kokkos::Impl::if_c<
307 Kokkos::Impl::is_same<
308 typename dual_view_type::t_dev_um::execution_space::memory_space,
309 typename TargetDeviceType::memory_space>::value,
310 typename dual_view_type::t_dev_um,
311 typename dual_view_type::t_host_um>::type ret;
321 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))) || \
322 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))))
325 template <
class Scalar>
327 :
public virtual Vector<Scalar,int,long long,EpetraNode>,
334 #undef XPETRA_TPETRAMULTIVECTOR_SHORT
432 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
436 typename dual_view_type::t_host_um getHostLocalView ()
const {
437 typename dual_view_type::t_host_um ret;
441 typename dual_view_type::t_dev_um getDeviceLocalView()
const {
442 typename dual_view_type::t_dev_um ret;
456 template<
class TargetDeviceType>
457 typename Kokkos::Impl::if_c<
458 Kokkos::Impl::is_same<
459 typename dual_view_type::t_dev_um::execution_space::memory_space,
460 typename TargetDeviceType::memory_space>::value,
461 typename dual_view_type::t_dev_um,
462 typename dual_view_type::t_host_um>::type
463 getLocalView ()
const {
464 typename Kokkos::Impl::if_c<
465 Kokkos::Impl::is_same<
466 typename dual_view_type::t_dev_um::execution_space::memory_space,
467 typename TargetDeviceType::memory_space>::value,
468 typename dual_view_type::t_dev_um,
469 typename dual_view_type::t_host_um>::type ret;
478 #endif // HAVE_XPETRA_EPETRA
Teuchos::ScalarTraits< Scalar >::magnitudeType norm1() const
Return 1-norm of this Vector.
std::string description() const
Return a simple one-line description of this object.
virtual ~TpetraVector()
Destructor.
void sumIntoLocalValue(LocalOrdinal myRow, const Scalar &value)
Adds specified value to existing value at the specified location.
TpetraVector(const Teuchos::RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &vec)
TpetraMultiVector constructor to wrap a Tpetra::MultiVector object.
virtual ~TpetraVector()
Destructor.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
void sumIntoGlobalValue(GlobalOrdinal globalRow, const Scalar &value)
Adds specified value to existing value at the specified location.
Scalar meanValue() const
Compute mean (average) value of this Vector.
Teuchos::ScalarTraits< Scalar >::magnitudeType normInf() const
Compute Inf-norm of this Vector.
Teuchos::ScalarTraits< Scalar >::magnitudeType norm2() const
Compute 2-norm of this Vector.
#define XPETRA_TPETRA_ETI_EXCEPTION(cl, obj, go, node)
void replaceGlobalValue(GlobalOrdinal globalRow, const Scalar &value)
Replace current value at the specified location with specified value.
TpetraVector(const Teuchos::RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &vec)
TpetraMultiVector constructor to wrap a Tpetra::MultiVector object.
void sumIntoGlobalValue(GlobalOrdinal globalRow, const Scalar &value)
Adds specified value to existing value at the specified location.
TpetraVector(const Teuchos::RCP< const Map > &map, const Teuchos::ArrayView< const Scalar > &A)
Set multi-vector values from an array using Teuchos memory management classes. (copy) ...
void replaceLocalValue(LocalOrdinal myRow, const Scalar &value)
Replace current value at the specified location with specified values.
Teuchos::ScalarTraits< Scalar >::magnitudeType normWeighted(const Vector &weights) const
Compute Weighted 2-norm (RMS Norm) of this Vector.
Scalar dot(const Vector &a) const
Computes dot product of this Vector against input Vector x.
void replaceGlobalValue(GlobalOrdinal globalRow, const Scalar &value)
Replace current value at the specified location with specified value.
Teuchos::ScalarTraits< Scalar >::magnitudeType norm1() const
Return 1-norm of this Vector.
void sumIntoLocalValue(LocalOrdinal myRow, const Scalar &value)
Adds specified value to existing value at the specified location.
Teuchos::ScalarTraits< Scalar >::magnitudeType norm2() const
Compute 2-norm of this Vector.
Teuchos::ScalarTraits< Scalar >::magnitudeType norm2() const
Compute 2-norm of this Vector.
std::string description() const
Return a simple one-line description of this object.
Teuchos::ScalarTraits< Scalar >::magnitudeType norm1() const
Return 1-norm of this Vector.
Scalar meanValue() const
Compute mean (average) value of this Vector.
static const EVerbosityLevel verbLevel_default
void sumIntoGlobalValue(GlobalOrdinal globalRow, const Scalar &value)
Adds specified value to existing value at the specified location.
void replaceLocalValue(LocalOrdinal myRow, const Scalar &value)
Replace current value at the specified location with specified values.
static magnitudeType magnitude(T a)
RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_Vector() const
Get the underlying Tpetra multivector.
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_Vector() const
Get the underlying Tpetra multivector.
RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_Vector() const
Get the underlying Tpetra multivector.
Teuchos::ScalarTraits< Scalar >::magnitudeType normWeighted(const Vector &weights) const
Compute Weighted 2-norm (RMS Norm) of this Vector.
Teuchos::ScalarTraits< Scalar >::magnitudeType normInf() const
Compute Inf-norm of this Vector.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
Teuchos::ScalarTraits< Scalar >::magnitudeType normInf() const
Compute Inf-norm of this Vector.
void sumIntoLocalValue(LocalOrdinal myRow, const Scalar &value)
Adds specified value to existing value at the specified location.
void replaceLocalValue(LocalOrdinal myRow, const Scalar &value)
Replace current value at the specified location with specified values.
std::string description() const
Return a simple one-line description of this object.
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_MultiVector() const
Get the underlying Tpetra multivector.
#define XPETRA_MONITOR(funcName)
TpetraVector(const Teuchos::RCP< const Map > &map, const Teuchos::ArrayView< const Scalar > &A)
Set multi-vector values from an array using Teuchos memory management classes. (copy) ...
TpetraVector(const Teuchos::RCP< const Map > &map, bool zeroOut=true)
Sets all vector entries to zero.
Scalar dot(const Vector &a) const
Computes dot product of this Vector against input Vector x.
void replaceGlobalValue(GlobalOrdinal globalRow, const Scalar &value)
Replace current value at the specified location with specified value.
Scalar dot(const Vector &a) const
Computes dot product of this Vector against input Vector x.
Scalar meanValue() const
Compute mean (average) value of this Vector.
TpetraVector(const Teuchos::RCP< const Map > &map, bool zeroOut=true)
Sets all vector entries to zero.
virtual ~TpetraVector()
Destructor.
Teuchos::ScalarTraits< Scalar >::magnitudeType normWeighted(const Vector &weights) const
Compute Weighted 2-norm (RMS Norm) of this Vector.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
TpetraVector(const Teuchos::RCP< const Map > &map, bool zeroOut=true)
Sets all vector entries to zero.