Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_BlockMultiVector_decl.hpp
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 // ************************************************************************
38 // @HEADER
39 
40 #ifndef TPETRA_BLOCKMULTIVECTOR_DECL_HPP
41 #define TPETRA_BLOCKMULTIVECTOR_DECL_HPP
42 
45 #include "Tpetra_MultiVector.hpp"
46 #include <memory>
47 #include <utility>
48 
49 namespace Tpetra {
50 
135 template<class Scalar,
136  class LO,
137  class GO,
138  class Node>
140  public Tpetra::DistObject<Scalar, LO, GO, Node>
141 {
142 private:
144  using STS = Teuchos::ScalarTraits<Scalar>;
145  using packet_type = typename dist_object_type::packet_type;
146 
147 public:
149 
150 
155 
157  using scalar_type = Scalar;
165  using local_ordinal_type = LO;
170 
172  using node_type = Node;
173 
176 
192  typedef Kokkos::View<impl_scalar_type *, device_type> little_vec_type;
193  typedef typename little_vec_type::HostMirror little_host_vec_type;
194 
200  typedef Kokkos::View<const impl_scalar_type *, device_type>
202 
204  host_device_type;
205  typedef Kokkos::View<const impl_scalar_type*, host_device_type>
206  const_little_host_vec_type;
207 
209 
211 
216  BlockMultiVector ();
217 
220 
223 
227 
231 
234  const Teuchos::DataAccess copyOrView);
235 
266  BlockMultiVector (const map_type& meshMap,
267  const LO blockSize,
268  const LO numVecs);
269 
274  BlockMultiVector (const map_type& meshMap,
275  const map_type& pointMap,
276  const LO blockSize,
277  const LO numVecs);
278 
291  BlockMultiVector (const mv_type& X_mv,
292  const map_type& meshMap,
293  const LO blockSize);
294 
300  const map_type& newMeshMap,
301  const map_type& newPointMap,
302  const size_t offset = 0);
303 
309  const map_type& newMeshMap,
310  const size_t offset = 0);
311 
313 
315 
322  static map_type
323  makePointMap (const map_type& meshMap, const LO blockSize);
324 
330  return pointMap_;
331  }
332 
334  LO getBlockSize () const {
335  return blockSize_;
336  }
337 
339  LO getNumVectors () const {
340  return static_cast<LO> (mv_.getNumVectors ());
341  }
342 
348  mv_type getMultiVectorView () const;
349 
351 
353 
355  void putScalar (const Scalar& val);
356 
358  void scale (const Scalar& val);
359 
366  void
367  update (const Scalar& alpha,
369  const Scalar& beta);
370 
392  void
393  blockWiseMultiply (const Scalar& alpha,
394  const Kokkos::View<const impl_scalar_type***,
395  device_type, Kokkos::MemoryUnmanaged>& D,
397 
428  void
429  blockJacobiUpdate (const Scalar& alpha,
430  const Kokkos::View<const impl_scalar_type***,
431  device_type, Kokkos::MemoryUnmanaged>& D,
434  const Scalar& beta);
435 
436 
438  template<class TargetMemorySpace>
439  bool need_sync () const {
440  return mv_.template need_sync<typename TargetMemorySpace::memory_space> ();
441  }
442 
444  bool need_sync_host() const {
445  return mv_.need_sync_host();
446  }
447 
449  bool need_sync_device() const {
450  return mv_.need_sync_device();
451  }
452 
454 
456 
474  bool replaceLocalValues (const LO localRowIndex, const LO colIndex, const Scalar vals[]);
475 
486  bool replaceGlobalValues (const GO globalRowIndex, const LO colIndex, const Scalar vals[]);
487 
498  bool sumIntoLocalValues (const LO localRowIndex, const LO colIndex, const Scalar vals[]);
499 
510  bool sumIntoGlobalValues (const GO globalRowIndex, const LO colIndex, const Scalar vals[]);
511 
512 
513  const_little_host_vec_type getLocalBlockHost(
514  const LO localRowIndex,
515  const LO colIndex,
516  const Access::ReadOnlyStruct) const;
517 
518  little_host_vec_type getLocalBlockHost(
519  const LO localRowIndex,
520  const LO colIndex,
521  const Access::ReadWriteStruct);
522 
526  little_host_vec_type getLocalBlockHost(
527  const LO localRowIndex,
528  const LO colIndex,
529  const Access::OverwriteAllStruct);
531 
532 protected:
538 
539 
540  virtual bool checkSizes (const Tpetra::SrcDistObject& source);
541 
542  virtual void
543  copyAndPermute
544  (const SrcDistObject& source,
545  const size_t numSameIDs,
546  const Kokkos::DualView<const local_ordinal_type*,
547  buffer_device_type>& permuteToLIDs,
548  const Kokkos::DualView<const local_ordinal_type*,
549  buffer_device_type>& permuteFromLIDs,
550  const CombineMode CM);
551 
552  virtual void
553  packAndPrepare
554  (const SrcDistObject& source,
555  const Kokkos::DualView<const local_ordinal_type*,
556  buffer_device_type>& exportLIDs,
557  Kokkos::DualView<packet_type*,
558  buffer_device_type>& exports,
559  Kokkos::DualView<size_t*,
560  buffer_device_type> numPacketsPerLID,
561  size_t& constantNumPackets);
562 
563  virtual void
564  unpackAndCombine
565  (const Kokkos::DualView<const local_ordinal_type*,
566  buffer_device_type>& importLIDs,
567  Kokkos::DualView<packet_type*,
568  buffer_device_type> imports,
569  Kokkos::DualView<size_t*,
570  buffer_device_type> numPacketsPerLID,
571  const size_t constantNumPackets,
572  const CombineMode combineMode);
573 
575 
576 protected:
577 
579  size_t getStrideX () const {
580  return static_cast<size_t> (1);
581  }
582 
584  size_t getStrideY () const {
585  return mv_.getStride ();
586  }
587 
590  bool isValidLocalMeshIndex (const LO meshLocalIndex) const;
591 
598 
599 private:
601  map_type pointMap_;
602 
603 protected:
606 
607 private:
608 
610  LO blockSize_;
611 
613  void
614  replaceLocalValuesImpl (const LO localRowIndex,
615  const LO colIndex,
616  const Scalar vals[]);
618  void
619  sumIntoLocalValuesImpl (const LO localRowIndex,
620  const LO colIndex,
621  const Scalar vals[]);
622 
623  static Teuchos::RCP<const mv_type>
624  getMultiVectorFromSrcDistObject (const Tpetra::SrcDistObject&);
625 
626  static Teuchos::RCP<const BlockMultiVector<Scalar, LO, GO, Node> >
627  getBlockMultiVectorFromSrcDistObject (const Tpetra::SrcDistObject& src);
628 };
629 
630 } // namespace Tpetra
631 
632 #endif // TPETRA_BLOCKMULTIVECTOR_DECL_HPP
Node node_type
The Kokkos Node type; a legacy thing that will go away at some point.
BlockMultiVector< Scalar, LO, GO, Node > & operator=(const BlockMultiVector< Scalar, LO, GO, Node > &)=default
Copy assigment (shallow copy).
void update(const Scalar &alpha, const BlockMultiVector< Scalar, LO, GO, Node > &X, const Scalar &beta)
Update: this = beta*this + alpha*X.
Kokkos::View< const impl_scalar_type *, device_type > const_little_vec_type
&quot;Const block view&quot; of all degrees of freedom at a mesh point, for a single column of the MultiVector...
void putScalar(const Scalar &val)
Fill all entries with the given value val.
bool need_sync_device() const
Whether this object needs synchronization to the device.
LO local_ordinal_type
The type of local indices.
size_t getNumVectors() const
Number of columns in the multivector.
typename map_type::device_type device_type
This class&#39; preferred Kokkos device type.
Kokkos::View< impl_scalar_type *, device_type > little_vec_type
&quot;Block view&quot; of all degrees of freedom at a mesh point, for a single column of the MultiVector...
void blockJacobiUpdate(const Scalar &alpha, const Kokkos::View< const impl_scalar_type ***, device_type, Kokkos::MemoryUnmanaged > &D, const BlockMultiVector< Scalar, LO, GO, Node > &X, BlockMultiVector< Scalar, LO, GO, Node > &Z, const Scalar &beta)
Block Jacobi update .
void blockWiseMultiply(const Scalar &alpha, const Kokkos::View< const impl_scalar_type ***, device_type, Kokkos::MemoryUnmanaged > &D, const BlockMultiVector< Scalar, LO, GO, Node > &X)
*this := alpha * D * X, where D is a block diagonal matrix.
size_t getStrideX() const
Stride between consecutive local entries in the same column.
Forward declaration of Tpetra::BlockCrsMatrix.
GO global_ordinal_type
The type of global indices.
virtual bool checkSizes(const Tpetra::SrcDistObject &source)
Compare the source and target (this) objects for compatibility.
typename Kokkos::Details::ArithTraits< Scalar >::val_type impl_scalar_type
The type used internally in place of Scalar.
bool replaceGlobalValues(const GO globalRowIndex, const LO colIndex, const Scalar vals[])
Replace all values at the given mesh point, using a global index.
size_t getStrideY() const
Stride between consecutive local entries in the same row.
bool replaceLocalValues(const LO localRowIndex, const LO colIndex, const Scalar vals[])
Replace all values at the given mesh point, using local row and column indices.
static map_type makePointMap(const map_type &meshMap, const LO blockSize)
Create and return the point Map corresponding to the given mesh Map and block size.
bool sumIntoLocalValues(const LO localRowIndex, const LO colIndex, const Scalar vals[])
Sum into all values at the given mesh point, using a local index.
bool need_sync_device() const
Whether this MultiVector needs synchronization to the device.
MultiVector for multiple degrees of freedom per mesh point.
bool need_sync_host() const
Whether this object needs synchronization to the host.
LO getBlockSize() const
Get the number of degrees of freedom per mesh point.
typename::Kokkos::Details::ArithTraits< Scalar >::val_type packet_type
The type of each datum being sent or received in an Import or Export.
bool isValidLocalMeshIndex(const LO meshLocalIndex) const
True if and only if meshLocalIndex is a valid local index in the mesh Map.
typename mv_type::device_type device_type
The Kokkos Device type.
CombineMode
Rule for combining data in an Import or Export.
map_type getPointMap() const
Get this BlockMultiVector&#39;s (previously computed) point Map.
Tpetra::MultiVector< Scalar, LO, GO, Node > mv_type
The specialization of Tpetra::MultiVector that this class uses.
bool sumIntoGlobalValues(const GO globalRowIndex, const LO colIndex, const Scalar vals[])
Sum into all values at the given mesh point, using a global index.
bool need_sync_host() const
Whether this MultiVector needs synchronization to the host.
mv_type mv_
The Tpetra::MultiVector used to represent the data.
Abstract base class for objects that can be the source of an Import or Export operation.
Forward declaration of Tpetra::BlockMultiVector.
Scalar scalar_type
The type of entries in the object.
bool need_sync() const
Whether this object needs synchronization to the given memory space.
LO getNumVectors() const
Get the number of columns (vectors) in the BlockMultiVector.
mv_type getMultiVectorView() const
Get a Tpetra::MultiVector that views this BlockMultiVector&#39;s data.
typename dist_object_type::buffer_device_type buffer_device_type
Kokkos::Device specialization used for communication buffers.
size_t getStride() const
Stride between columns in the multivector.
Kokkos::Device< typename device_type::execution_space, buffer_memory_space > buffer_device_type
void scale(const Scalar &val)
Multiply all entries in place by the given value val.
Base class for distributed Tpetra objects that support data redistribution.
typename mv_type::impl_scalar_type impl_scalar_type
The implementation type of entries in the object.
map_type meshMap_
Mesh Map given to constructor.