Tpetra parallel linear algebra  Version of the Day
Tpetra_Experimental_BlockVector_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_BLOCKVECTOR_DECL_HPP
43 #define TPETRA_EXPERIMENTAL_BLOCKVECTOR_DECL_HPP
44 
45 #include <Tpetra_Experimental_BlockMultiVector.hpp>
46 #include <Tpetra_Vector.hpp>
47 
48 namespace Tpetra {
49 namespace Experimental {
50 
76 template<class Scalar = ::Tpetra::Details::DefaultTypes::scalar_type,
78  class GO = ::Tpetra::Details::DefaultTypes::global_ordinal_type,
80 class BlockVector : public BlockMultiVector<Scalar, LO, GO, Node> {
81 private:
83  typedef Teuchos::ScalarTraits<Scalar> STS;
84 
85 public:
87 
88 
98  typedef typename base_type::node_type node_type;
100  typedef typename Node::device_type device_type;
101 
108 
123  typedef Kokkos::View<impl_scalar_type*,
124  Kokkos::LayoutRight,
125  device_type,
126  Kokkos::MemoryTraits<Kokkos::Unmanaged> >
128 
133  typedef Kokkos::View<const impl_scalar_type*,
134  Kokkos::LayoutRight,
135  device_type,
136  Kokkos::MemoryTraits<Kokkos::Unmanaged> >
138 
140 
142 
172  BlockVector (const map_type& meshMap, const LO blockSize);
173 
178  BlockVector (const map_type& meshMap,
179  const map_type& pointMap,
180  const LO blockSize);
181 
195  BlockVector (const mv_type& X_mv,
196  const map_type& meshMap,
197  const LO blockSize);
198 
204  const map_type& newMeshMap,
205  const map_type& newPointMap,
206  const size_t offset = 0);
207 
213  const map_type& newMeshMap,
214  const size_t offset = 0);
215 
220  BlockVector ();
221 
223 
225 
231 
233 
235 
251  bool replaceLocalValues (const LO localRowIndex, const Scalar vals[]) const;
252 
262  bool replaceGlobalValues (const GO globalRowIndex, const Scalar vals[]) const;
263 
273  bool sumIntoLocalValues (const LO localRowIndex, const Scalar vals[]) const;
274 
284  bool sumIntoGlobalValues (const GO globalRowIndex, const Scalar vals[]) const;
285 
296  bool getLocalRowView (const LO localRowIndex, Scalar*& vals) const;
297 
308  bool getGlobalRowView (const GO globalRowIndex, Scalar*& vals) const;
309 
322  little_vec_type getLocalBlock (const LO localRowIndex) const;
324 };
325 
326 } // namespace Experimental
327 } // namespace Tpetra
328 
329 #endif // TPETRA_EXPERIMENTAL_BLOCKMULTIVECTOR_DECL_HPP
Namespace Tpetra contains the class and methods constituting the Tpetra library.
base_type::impl_scalar_type impl_scalar_type
The implementation type of entries in the vector.
bool sumIntoGlobalValues(const GO globalRowIndex, const Scalar vals[]) const
Sum into all values at the given mesh point, using a global index.
bool getGlobalRowView(const GO globalRowIndex, Scalar *&vals) const
Get a writeable view of the entries at the given mesh point, using a global index.
KokkosClassic::DefaultNode::DefaultNodeType node_type
Default value of Node template parameter.
MultiVector for multiple degrees of freedom per mesh point.
Vector for multiple degrees of freedom per mesh point.
Kokkos::View< impl_scalar_type *, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > little_vec_type
"Block view" of all degrees of freedom at a mesh point, for a single column of the MultiVector...
int local_ordinal_type
Default value of LocalOrdinal template parameter.
mv_type::impl_scalar_type impl_scalar_type
The implementation type of entries in the matrix.
Tpetra::Vector< Scalar, LO, GO, Node > vec_type
The specialization of Tpetra::Vector that this class uses.
base_type::node_type node_type
The Kokkos Node type.
Kokkos::View< impl_scalar_type *, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > little_vec_type
"Block view" of all degrees of freedom at a mesh point.
base_type::local_ordinal_type local_ordinal_type
The type of local indices.
bool sumIntoLocalValues(const LO localRowIndex, const Scalar vals[]) const
Sum into all values at the given mesh point, using a local index.
Tpetra::Map< LO, GO, Node > map_type
The specialization of Tpetra::Map that this class uses.
bool replaceLocalValues(const LO localRowIndex, const Scalar vals[]) const
Replace all values at the given mesh point, using a local index.
base_type::scalar_type scalar_type
The type of entries in the vector.
base_type::global_ordinal_type global_ordinal_type
The type of global indices.
vec_type getVectorView()
Get a Tpetra::Vector that views this BlockVector&#39;s data.
Kokkos::View< const impl_scalar_type *, Kokkos::LayoutRight, device_type, Kokkos::MemoryTraits< Kokkos::Unmanaged > > const_little_vec_type
"Const block view" of all degrees of freedom at a mesh point.
Node::device_type device_type
The Kokkos Device type.
double scalar_type
Default value of Scalar template parameter.
bool getLocalRowView(const LO localRowIndex, Scalar *&vals) const
Get a writeable view of the entries at the given mesh point, using a local index. ...
little_vec_type getLocalBlock(const LO localRowIndex) const
Get a view of the degrees of freedom at the given mesh point, using a local index.
bool replaceGlobalValues(const GO globalRowIndex, const Scalar vals[]) const
Replace all values at the given mesh point, using a global index.
A distributed dense vector.
Tpetra::MultiVector< Scalar, LO, GO, Node > mv_type
The specialization of Tpetra::MultiVector that this class uses.