Tpetra parallel linear algebra  Version of the Day
Tpetra_Vector_def.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Tpetra: Templated Linear Algebra Services Package
5 // Copyright (2008) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #ifndef TPETRA_VECTOR_DEF_HPP
43 #define TPETRA_VECTOR_DEF_HPP
44 
52 
53 #include "Tpetra_MultiVector.hpp"
54 #include "Tpetra_Details_gathervPrint.hpp"
55 #include "KokkosCompat_View.hpp"
56 #include "Kokkos_Blas1_MV.hpp"
57 
58 namespace Tpetra {
59 
60  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
62  Vector (const Teuchos::RCP<const map_type>& map,
63  const bool zeroOut)
64  : base_type (map, 1, zeroOut)
65  {}
66 
67  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
70  : base_type (source)
71  {}
72 
73  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
76  const Teuchos::DataAccess copyOrView)
77  : base_type (source, copyOrView)
78  {}
79 
80  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
82  Vector (const Teuchos::RCP<const map_type>& map,
83  const Teuchos::ArrayView<const Scalar>& values)
84  : base_type (map, values, values.size (), 1)
85  {}
86 
87  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
89  Vector (const Teuchos::RCP<const map_type>& map,
90  const dual_view_type& view)
91  : base_type (map, view)
92  {}
93 
94  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
96  Vector (const Teuchos::RCP<const map_type>& map,
97  const dual_view_type& view,
98  const dual_view_type& origView)
99  : base_type (map, view, origView)
100  {}
101 
102  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
105  {}
106 
107  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
108  void
110  replaceGlobalValue (const GlobalOrdinal globalRow, const Scalar& value) const {
111  this->base_type::replaceGlobalValue (globalRow, 0, value);
112  }
113 
114  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
115  void
117  sumIntoGlobalValue (const GlobalOrdinal globalRow,
118  const Scalar& value,
119  const bool atomic) const
120  {
121  this->base_type::sumIntoGlobalValue (globalRow, 0, value, atomic);
122  }
123 
124  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
125  void
127  replaceLocalValue (const LocalOrdinal myRow, const Scalar& value) const {
128  this->base_type::replaceLocalValue (myRow, 0, value);
129  }
130 
131  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
132  void
134  sumIntoLocalValue (const LocalOrdinal globalRow,
135  const Scalar& value,
136  const bool atomic) const
137  {
138  this->base_type::sumIntoLocalValue (globalRow, 0, value, atomic);
139  }
140 
141  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
142  void
144  get1dCopy (const Teuchos::ArrayView<Scalar>& A) const {
145  const size_t lda = this->getLocalLength ();
146  this->get1dCopy (A, lda);
147  }
148 
149  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
153  {
154  dot_type result;
155  this->dot (y, Teuchos::arrayView (&result, 1));
156  return result;
157  }
158 
159  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
160  Scalar
162  meanValue () const
163  {
164  Scalar mean;
165  this->meanValue (Teuchos::arrayView (&mean, 1));
166  return mean;
167  }
168 
169  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
172  norm1 () const
173  {
174  mag_type norm;
175  this->norm1 (Teuchos::arrayView (&norm, 1));
176  return norm;
177  }
178 
179  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
182  norm2 () const
183  {
184  mag_type norm;
185  this->norm2 (Teuchos::arrayView (&norm, 1));
186  return norm;
187  }
188 
189  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
192  normInf () const
193  {
194  mag_type norm;
195  this->normInf (Teuchos::arrayView (&norm, 1));
196  return norm;
197  }
198 
199  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
201  TPETRA_DEPRECATED
204  {
205  using Kokkos::ALL;
206  using Kokkos::subview;
207  using Teuchos::Comm;
208  using Teuchos::RCP;
209  using Teuchos::reduceAll;
210  using Teuchos::REDUCE_SUM;
211  typedef Kokkos::Details::ArithTraits<impl_scalar_type> ATS;
212  typedef Kokkos::Details::ArithTraits<mag_type> ATM;
213  typedef Kokkos::View<mag_type, device_type> norm_view_type; // just one
214  const char tfecfFuncName[] = "normWeighted: ";
215 
216 #ifdef HAVE_TPETRA_DEBUG
217  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
218  ! this->getMap ()->isCompatible (*weights.getMap ()), std::runtime_error,
219  "Vectors do not have compatible Maps:" << std::endl
220  << "this->getMap(): " << std::endl << *this->getMap()
221  << "weights.getMap(): " << std::endl << *weights.getMap() << std::endl);
222 #else
223  const size_t lclNumRows = this->getLocalLength ();
224  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
225  lclNumRows != weights.getLocalLength (), std::runtime_error,
226  "Vectors do not have the same local length.");
227 #endif // HAVE_TPETRA_DEBUG
228 
229  norm_view_type lclNrm ("lclNrm"); // local norm
230  mag_type gblNrm = ATM::zero (); // return value
231 
232  auto X_lcl = this->template getLocalView<device_type> ();
233  auto W_lcl = this->template getLocalView<device_type> ();
234  KokkosBlas::nrm2w_squared (lclNrm,
235  subview (X_lcl, ALL (), 0),
236  subview (W_lcl, ALL (), 0));
237  const mag_type OneOverN =
238  ATM::one () / static_cast<mag_type> (this->getGlobalLength ());
239  RCP<const Comm<int> > comm = this->getMap ().is_null () ?
240  Teuchos::null : this->getMap ()->getComm ();
241 
242  if (! comm.is_null () && this->isDistributed ()) {
243  // Assume that MPI can access device memory.
244  reduceAll<int, mag_type> (*comm, REDUCE_SUM, 1, lclNrm.ptr_on_device (),
245  &gblNrm);
246  gblNrm = ATM::sqrt (gblNrm * OneOverN);
247  }
248  else {
249  auto lclNrm_h = Kokkos::create_mirror_view (lclNrm);
250  Kokkos::deep_copy (lclNrm_h, lclNrm);
251  gblNrm = ATM::sqrt (ATS::magnitude (lclNrm_h()) * OneOverN);
252  }
253 
254  return gblNrm;
255  }
256 
257  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
259  description () const
260  {
261  return this->descriptionImpl ("Tpetra::Vector");
262  }
263 
264  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
266  describe (Teuchos::FancyOStream& out,
267  const Teuchos::EVerbosityLevel verbLevel) const
268  {
269  this->describeImpl (out, "Tpetra::Vector", verbLevel);
270  }
271 
272  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
275  {
276  // The 2-argument copy constructor with second argument =
277  // Teuchos::Copy does a deep copy of its input.
279 
280  // The Kokkos refactor version of Vector has view semantics, so
281  // returning the Vector directly, rather than through RCP, only
282  // does a shallow copy.
283  return dst;
284  }
285 
286  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
287  Teuchos::RCP<const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node, classic> >
288  Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node, classic>::
289  offsetView (const Teuchos::RCP<const map_type>& subMap,
290  const size_t offset) const
291  {
292  using Kokkos::ALL;
293  using Kokkos::subview;
294  using Teuchos::rcp;
295  typedef Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node, classic> V;
296 
297  const size_t newNumRows = subMap->getNodeNumElements ();
298  const bool tooManyElts = newNumRows + offset > this->getOrigNumLocalRows ();
299  if (tooManyElts) {
300  const int myRank = this->getMap ()->getComm ()->getRank ();
301  TEUCHOS_TEST_FOR_EXCEPTION(
302  newNumRows + offset > this->getLocalLength (), std::runtime_error,
303  "Tpetra::Vector::offsetView(NonConst): Invalid input Map. The input "
304  "Map owns " << newNumRows << " entries on process " << myRank << ". "
305  "offset = " << offset << ". Yet, the Vector contains only "
306  << this->getOrigNumLocalRows () << " rows on this process.");
307  }
308 
309  const std::pair<size_t, size_t> offsetPair (offset, offset + newNumRows);
310  // Need 'this->' to get view_ and origView_ from parent class.
311  return rcp (new V (subMap,
312  subview (this->view_, offsetPair, ALL ()),
313  this->origView_));
314  }
315 
316  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
317  Teuchos::RCP<Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node, classic> >
318  Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node, classic>::
319  offsetViewNonConst (const Teuchos::RCP<const map_type>& subMap,
320  const size_t offset)
321  {
322  typedef Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node, classic> V;
323  return Teuchos::rcp_const_cast<V> (this->offsetView (subMap, offset));
324  }
325 
326 } // namespace Tpetra
327 
336 #define TPETRA_VECTOR_INSTANT(SCALAR,LO,GO,NODE) \
337  template class Vector< SCALAR , LO , GO , NODE >; \
338  template Vector< SCALAR , LO , GO , NODE > createCopy (const Vector< SCALAR , LO , GO , NODE >& src);
339 
340 #endif // TPETRA_VECTOR_DEF_HPP
Vector(const Teuchos::RCP< const map_type > &map, const bool zeroOut=true)
Basic constructor.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
MultiVector< ST, LO, GO, NT, classic > createCopy(const MultiVector< ST, LO, GO, NT, classic > &src)
Return a deep copy of the given MultiVector.
virtual std::string description() const
Return a one-line description of this object.
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Describe this object in a human-readable way to the given output stream.
base_type::dot_type dot_type
Type of an inner ("dot") product result.
void deep_copy(MultiVector< DS, DL, DG, DN, dstClassic > &dst, const MultiVector< SS, SL, SG, SN, srcClassic > &src)
Copy the contents of the MultiVector src into dst.
Kokkos::Details::InnerProductSpaceTraits< impl_scalar_type >::dot_type dot_type
Type of an inner ("dot") product result.
base_type::mag_type mag_type
Type of a norm result.
mag_type normInf() const
Compute Inf-norm of this Vector.
mag_type norm1() const
Return 1-norm of this Vector.
size_t getLocalLength() const
Local number of rows on the calling process.
virtual ~Vector()
Destructor.
mag_type norm2() const
Compute 2-norm of this Vector.
dot_type dot(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node, classic > &y) const
Computes dot product of this Vector against input Vector x.
virtual Teuchos::RCP< const map_type > getMap() const
The Map describing the parallel distribution of this object.
void sumIntoLocalValue(const LocalOrdinal myRow, const Scalar &value, const bool atomic=base_type::useAtomicUpdatesByDefault) const
Add value to existing value, using local (row) index.
void get1dCopy(const Teuchos::ArrayView< Scalar > &A) const
Return multi-vector values in user-provided two-dimensional array (using Teuchos memory management cl...
void replaceGlobalValue(const GlobalOrdinal globalRow, const Scalar &value) const
Replace current value at the specified location with specified value.
Kokkos::Details::ArithTraits< impl_scalar_type >::mag_type mag_type
Type of a norm result.
void replaceLocalValue(const LocalOrdinal myRow, const Scalar &value) const
Replace current value at the specified location with specified values.
Kokkos::DualView< impl_scalar_type **, Kokkos::LayoutLeft, typename execution_space::execution_space > dual_view_type
Kokkos::DualView specialization used by this class.
void sumIntoGlobalValue(const GlobalOrdinal globalRow, const Scalar &value, const bool atomic=base_type::useAtomicUpdatesByDefault) const
Add value to existing value, using global (row) index.
A distributed dense vector.
mag_type TPETRA_DEPRECATED normWeighted(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node, classic > &weights) const
Compute Weighted 2-norm (RMS Norm) of this Vector.
Scalar meanValue() const
Compute mean (average) value of this Vector.