Tpetra parallel linear algebra  Version of the Day
Tpetra_Experimental_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 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #ifndef TPETRA_EXPERIMENTAL_BLOCKMULTIVECTOR_DECL_HPP
43 #define TPETRA_EXPERIMENTAL_BLOCKMULTIVECTOR_DECL_HPP
44 
45 #include <Tpetra_MultiVector.hpp>
47 
48 namespace Tpetra {
49 namespace Experimental {
50 
51 #ifndef DOXYGEN_SHOULD_SKIP_THIS
52 // Forward declaration of BlockCrsMatrix, needed for debugging.
53 template<class S, class LO, class GO, class N> class BlockCrsMatrix;
54 #endif // DOXYGEN_SHOULD_SKIP_THIS
55 
144 template<class Scalar = Details::DefaultTypes::scalar_type,
146  class GO = Details::DefaultTypes::global_ordinal_type,
149  public Tpetra::DistObject<Scalar, LO, GO, Node>
150 {
151 private:
153  typedef Teuchos::ScalarTraits<Scalar> STS;
154 
155 protected:
156  typedef Scalar packet_type;
157 
158 public:
160 
161 
166 
168  typedef Scalar scalar_type;
177  typedef LO local_ordinal_type;
181  typedef Node node_type;
182 
199 
206 
208 
210 
241  BlockMultiVector (const map_type& meshMap,
242  const LO blockSize,
243  const LO numVecs);
244 
249  BlockMultiVector (const map_type& meshMap,
250  const map_type& pointMap,
251  const LO blockSize,
252  const LO numVecs);
253 
266  BlockMultiVector (const mv_type& X_mv,
267  const map_type& meshMap,
268  const LO blockSize);
269 
275  const map_type& newMeshMap,
276  const map_type& newPointMap,
277  const size_t offset = 0);
278 
284  const map_type& newMeshMap,
285  const size_t offset = 0);
286 
291  BlockMultiVector ();
292 
294 
296 
303  static map_type
304  makePointMap (const map_type& meshMap, const LO blockSize);
305 
310  map_type getPointMap () const {
311  return pointMap_;
312  }
313 
315  LO getBlockSize () const {
316  return blockSize_;
317  }
318 
320  LO getNumVectors () const {
321  return static_cast<LO> (mv_.getNumVectors ());
322  }
323 
329  mv_type getMultiVectorView ();
330 
332 
334 
336  void putScalar (const Scalar& val);
337 
339  void scale (const Scalar& val);
340 
342 
344 
362  bool replaceLocalValues (const LO localRowIndex, const LO colIndex, const Scalar vals[]) const;
363 
374  bool replaceGlobalValues (const GO globalRowIndex, const LO colIndex, const Scalar vals[]) const;
375 
386  bool sumIntoLocalValues (const LO localRowIndex, const LO colIndex, const Scalar vals[]) const;
387 
398  bool sumIntoGlobalValues (const GO globalRowIndex, const LO colIndex, const Scalar vals[]) const;
399 
410  bool getLocalRowView (const LO localRowIndex, const LO colIndex, Scalar*& vals) const;
411 
422  bool getGlobalRowView (const GO globalRowIndex, const LO colIndex, Scalar*& vals) const;
423 
435  little_vec_type
436  getLocalBlock (const LO localRowIndex, const LO colIndex) const;
438 
439 protected:
446 
447 
448  virtual bool checkSizes (const Tpetra::SrcDistObject& source);
449 
450  virtual void
451  copyAndPermute (const Tpetra::SrcDistObject& source,
452  size_t numSameIDs,
453  const Teuchos::ArrayView<const LO>& permuteToLIDs,
454  const Teuchos::ArrayView<const LO>& permuteFromLIDs);
455 
456  virtual void
457  packAndPrepare (const Tpetra::SrcDistObject& source,
458  const Teuchos::ArrayView<const LO>& exportLIDs,
459  Teuchos::Array<impl_scalar_type>& exports,
460  const Teuchos::ArrayView<size_t>& numPacketsPerLID,
461  size_t& constantNumPackets,
462  Tpetra::Distributor& distor);
463 
464  virtual void
465  unpackAndCombine (const Teuchos::ArrayView<const LO> &importLIDs,
466  const Teuchos::ArrayView<const impl_scalar_type> &imports,
467  const Teuchos::ArrayView<size_t> &numPacketsPerLID,
468  size_t constantNumPackets,
469  Tpetra::Distributor& distor,
472 
473 protected:
475  impl_scalar_type* getRawPtr () const {
476  return mvData_;
477  }
478 
480  size_t getStrideX () const {
481  return static_cast<size_t> (1);
482  }
483 
485  size_t getStrideY () const {
486  return mv_.getStride ();
487  }
488 
491  bool isValidLocalMeshIndex (const LO meshLocalIndex) const {
492  return meshLocalIndex != Teuchos::OrdinalTraits<LO>::invalid () &&
493  meshMap_.isNodeLocalElement (meshLocalIndex);
494  }
495 
496 private:
497  // mfh 20 May 2014: I'm only using this for debugging.
498  template<class Scalar2, class LO2, class GO2, class Node2>
499  friend class BlockCrsMatrix;
500 
506  map_type meshMap_;
507 
509  map_type pointMap_;
510 
511 protected:
513  mv_type mv_;
514 
515 private:
524  impl_scalar_type* mvData_;
525 
527  LO blockSize_;
528 
530  void
531  replaceLocalValuesImpl (const LO localRowIndex,
532  const LO colIndex,
533  const Scalar vals[]) const;
535  void
536  sumIntoLocalValuesImpl (const LO localRowIndex,
537  const LO colIndex,
538  const Scalar vals[]) const;
539 
540  static Teuchos::RCP<const mv_type>
541  getMultiVectorFromSrcDistObject (const Tpetra::SrcDistObject&);
542 
543  static Teuchos::RCP<const BlockMultiVector<Scalar, LO, GO, Node> >
544  getBlockMultiVectorFromSrcDistObject (const Tpetra::SrcDistObject& src);
545 };
546 
547 } // namespace Experimental
548 } // namespace Tpetra
549 
550 #endif // TPETRA_EXPERIMENTAL_BLOCKMULTIVECTOR_DECL_HPP
Namespace Tpetra contains the class and methods constituting the Tpetra library.
mv_type mv_
The Tpetra::MultiVector used to represent the data.
KokkosClassic::DefaultNode::DefaultNodeType node_type
Default value of Node template parameter.
void putScalar(const Scalar &val)
Fill all entries with the given value val.
mv_type getMultiVectorView()
Get a Tpetra::MultiVector that views this BlockMultiVector&#39;s data.
Tpetra::Map< LO, GO, Node > map_type
The specialization of Tpetra::Map that this class uses.
bool replaceGlobalValues(const GO globalRowIndex, const LO colIndex, const Scalar vals[]) const
Replace all values at the given mesh point, using a global index.
bool sumIntoGlobalValues(const GO globalRowIndex, const LO colIndex, const Scalar vals[]) const
Sum into all values at the given mesh point, using a global index.
MultiVector for multiple degrees of freedom per mesh point.
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.
virtual void unpackAndCombine(const Teuchos::ArrayView< const LO > &importLIDs, const Teuchos::ArrayView< const impl_scalar_type > &imports, const Teuchos::ArrayView< size_t > &numPacketsPerLID, size_t constantNumPackets, Tpetra::Distributor &distor, Tpetra::CombineMode CM)
Perform any unpacking and combining after communication.
int local_ordinal_type
Default value of LocalOrdinal template parameter.
void scale(const Scalar &val)
Multiply all entries in place by the given value val.
mv_type::impl_scalar_type impl_scalar_type
The implementation type of entries in the matrix.
little_vec_type getLocalBlock(const LO localRowIndex, const LO colIndex) const
Get a view of the degrees of freedom at the given mesh point.
bool isNodeLocalElement(LocalOrdinal localIndex) const
Whether the given local index is valid for this Map on the calling process.
LO getBlockSize() const
Get the number of degrees of freedom per mesh point.
impl_scalar_type * getRawPtr() const
Raw pointer to the MultiVector&#39;s data.
LittleBlock, LittleVector, and kernels.
Nonowning view of a set of degrees of freedom corresponding to a mesh point in a block vector or mult...
LO getNumVectors() const
Get the number of columns (vectors) in the BlockMultiVector.
Sets up and executes a communication plan for a Tpetra DistObject.
CombineMode
Rule for combining data in an Import or Export.
bool getLocalRowView(const LO localRowIndex, const LO colIndex, Scalar *&vals) const
Get a writeable view of the entries at the given mesh point, using a local index. ...
Abstract base class for objects that can be the source of an Import or Export operation.
double scalar_type
Default value of Scalar template parameter.
bool isValidLocalMeshIndex(const LO meshLocalIndex) const
True if and only if meshLocalIndex is a valid local index in the mesh Map.
size_t getStride() const
Stride between columns in the multivector.
Kokkos::Details::ArithTraits< Scalar >::val_type impl_scalar_type
The type used internally in place of Scalar.
map_type getPointMap() const
Get this BlockMultiVector&#39;s (previously computed) point Map.
virtual bool checkSizes(const Tpetra::SrcDistObject &source)
Compare the source and target (this) objects for compatibility.
LittleVector< impl_scalar_type, LO > little_vec_type
"Block view" of all degrees of freedom at a mesh point, for a single column of the MultiVector...
size_t getStrideY() const
Stride between consecutive local entries in the same row.
LittleVector< const impl_scalar_type, LO > const_little_vec_type
"Const block view" of all degrees of freedom at a mesh point, for a single column of the MultiVector...
virtual void copyAndPermute(const Tpetra::SrcDistObject &source, size_t numSameIDs, const Teuchos::ArrayView< const LO > &permuteToLIDs, const Teuchos::ArrayView< const LO > &permuteFromLIDs)
Perform copies and permutations that are local to this process.
size_t getNumVectors() const
Number of columns in the multivector.
bool replaceLocalValues(const LO localRowIndex, const LO colIndex, const Scalar vals[]) const
Replace all values at the given mesh point, using local row and column indices.
size_t getStrideX() const
Stride between consecutive local entries in the same column.
virtual void packAndPrepare(const Tpetra::SrcDistObject &source, const Teuchos::ArrayView< const LO > &exportLIDs, Teuchos::Array< impl_scalar_type > &exports, const Teuchos::ArrayView< size_t > &numPacketsPerLID, size_t &constantNumPackets, Tpetra::Distributor &distor)
Perform any packing or preparation required for communication.
Base class for distributed Tpetra objects that support data redistribution.
bool sumIntoLocalValues(const LO localRowIndex, const LO colIndex, const Scalar vals[]) const
Sum into all values at the given mesh point, using a local index.
Tpetra::MultiVector< Scalar, LO, GO, Node > mv_type
The specialization of Tpetra::MultiVector that this class uses.
bool getGlobalRowView(const GO globalRowIndex, const LO colIndex, Scalar *&vals) const
Get a writeable view of the entries at the given mesh point, using a global index.