Tpetra parallel linear algebra Version of the Day
TpetraExt_MMHelpers_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_MMHELPERS_DEF_HPP
43#define TPETRA_MMHELPERS_DEF_HPP
44
46#include "Teuchos_VerboseObject.hpp"
47
52namespace Tpetra {
53
54template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
55CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsMatrixStruct()
56{
57}
58
59template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
60CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>::~CrsMatrixStruct()
61{
62 deleteContents();
63}
64
65template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
66void CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
67deleteContents ()
68{
69 importMatrix.reset();
70 origMatrix = Teuchos::null;
71}
72
73template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
74int dumpCrsMatrixStruct (const CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& M)
75{
76 std::cout << "proc " << M.rowMap->Comm().MyPID()<<std::endl;
77 std::cout << "numRows: " << M.numRows<<std::endl;
78 for(LocalOrdinal i=0; i<M.numRows; ++i) {
79 for(LocalOrdinal j=0; j<M.numEntriesPerRow[i]; ++j) {
80 std::cout << " "<<M.rowMap->GID(i)<<" "
81 <<M.colMap->GID(M.indices[i][j])<<" "<<M.values[i][j]<<std::endl;
82 }
83 }
84
85 return 0;
86}
87
88template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
89CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
90CrsWrapper_CrsMatrix (CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& crsmatrix)
91 : crsmat_ (crsmatrix)
92{
93}
94
95template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
96CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::~CrsWrapper_CrsMatrix()
97{
98}
99
100template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
101Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
102CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getRowMap() const
103{
104 return crsmat_.getRowMap();
105}
106
107template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
108bool CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
109isFillComplete ()
110{
111 return crsmat_.isFillComplete ();
112}
113
114template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
115void
116CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
117insertGlobalValues (GlobalOrdinal globalRow,
118 const Teuchos::ArrayView<const GlobalOrdinal> &indices,
119 const Teuchos::ArrayView<const Scalar> &values)
120{
121 crsmat_.insertGlobalValues (globalRow, indices, values);
122}
123
124template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
125void
126CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
127sumIntoGlobalValues (GlobalOrdinal globalRow,
128 const Teuchos::ArrayView<const GlobalOrdinal> &indices,
129 const Teuchos::ArrayView<const Scalar> &values)
130{
131 crsmat_.sumIntoGlobalValues (globalRow, indices, values);
132}
133
134
135
136template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
137CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
138CrsWrapper_GraphBuilder (const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >& map)
139 : graph_(),
140 rowmap_(map),
141 max_row_length_(0)
142{
143 Teuchos::ArrayView<const GlobalOrdinal> rows = map->getNodeElementList ();
144 const LocalOrdinal numRows = static_cast<LocalOrdinal> (rows.size ());
145 for (LocalOrdinal i = 0; i < numRows; ++i) {
146 graph_[rows[i]] = new std::set<GlobalOrdinal>;
147 }
148}
149
150template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
151CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
152~CrsWrapper_GraphBuilder ()
153{
154 typename std::map<GlobalOrdinal,std::set<GlobalOrdinal>*>::iterator
155 iter = graph_.begin(), iter_end = graph_.end();
156 for (; iter != iter_end; ++iter) {
157 delete iter->second;
158 }
159 graph_.clear ();
160}
161
162template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
163bool CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::isFillComplete()
164{
165 return false;
166}
167
168template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
169void
170CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
171insertGlobalValues (GlobalOrdinal globalRow,
172 const Teuchos::ArrayView<const GlobalOrdinal> &indices,
173 const Teuchos::ArrayView<const Scalar> &/* values */)
174{
175 typename std::map<GlobalOrdinal,std::set<GlobalOrdinal>*>::iterator
176 iter = graph_.find (globalRow);
177
178 TEUCHOS_TEST_FOR_EXCEPTION(
179 iter == graph_.end(), std::runtime_error,
180 "Tpetra::CrsWrapper_GraphBuilder::insertGlobalValues could not find row "
181 << globalRow << " in the graph. Super bummer man. Hope you figure it out.");
182
183 std::set<GlobalOrdinal>& cols = * (iter->second);
184
185 for (typename Teuchos::ArrayView<const GlobalOrdinal>::size_type i = 0;
186 i < indices.size (); ++i) {
187 cols.insert (indices[i]);
188 }
189
190 const global_size_t row_length = static_cast<global_size_t> (cols.size ());
191 if (row_length > max_row_length_) {
192 max_row_length_ = row_length;
193 }
194}
195
196template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
197void
198CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
199sumIntoGlobalValues (GlobalOrdinal globalRow,
200 const Teuchos::ArrayView<const GlobalOrdinal> &indices,
201 const Teuchos::ArrayView<const Scalar> &values)
202{
203 insertGlobalValues (globalRow, indices, values);
204}
205
206template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
207std::map<GlobalOrdinal,std::set<GlobalOrdinal>*>&
208CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::get_graph ()
209{
210 return graph_;
211}
212
213template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
214void
215insert_matrix_locations (CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>& graphbuilder,
216 CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C)
217{
218 global_size_t max_row_length = graphbuilder.get_max_row_length();
219 if (max_row_length < 1) return;
220
221 Teuchos::Array<GlobalOrdinal> indices(max_row_length);
222 Teuchos::Array<Scalar> zeros(max_row_length, Teuchos::ScalarTraits<Scalar>::zero());
223
224 typedef std::map<GlobalOrdinal,std::set<GlobalOrdinal>*> Graph;
225 typedef typename Graph::iterator GraphIter;
226 Graph& graph = graphbuilder.get_graph ();
227
228 const GraphIter iter_end = graph.end ();
229 for (GraphIter iter = graph.begin (); iter != iter_end; ++iter) {
230 const GlobalOrdinal row = iter->first;
231 const std::set<GlobalOrdinal>& cols = * (iter->second);
232 // "copy" entries out of set into contiguous array storage
233 const size_t num_entries = std::copy (cols.begin (), cols.end (), indices.begin ()) - indices.begin ();
234 // insert zeros into the result matrix at the appropriate locations
235 C.insertGlobalValues (row, indices (0, num_entries), zeros (0, num_entries));
236 }
237}
238
239} // namespace Tpetra
240
241//
242// Explicit instantiation macro
243//
244// Must be expanded from within the Tpetra namespace!
245//
246
247#define TPETRA_CRSMATRIXSTRUCT_INSTANT(SCALAR,LO,GO,NODE) \
248 \
249 template class CrsMatrixStruct< SCALAR , LO , GO , NODE >;
250
251#define TPETRA_CRSWRAPPER_INSTANT(SCALAR,LO,GO,NODE) \
252 \
253 template class CrsWrapper< SCALAR , LO , GO , NODE >;
254
255#define TPETRA_CRSWRAPPER_CRSMATRIX_INSTANT(SCALAR,LO,GO,NODE) \
256 \
257 template class CrsWrapper_CrsMatrix< SCALAR , LO , GO , NODE >;
258
259#define TPETRA_CRSWRAPPER_GRAPHBUILDER_INSTANT(SCALAR,LO,GO,NODE) \
260 \
261 template class CrsWrapper_GraphBuilder< SCALAR , LO , GO , NODE >;
262
263#endif // TPETRA_MMHELPERS_DEF_HPP
Declaration of Tpetra::MMMultiMultiply and nonmember constructors.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
size_t global_size_t
Global size_t object.