Tpetra parallel linear algebra  Version of the Day
Tpetra_CrsGraph_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_CRSGRAPH_DEF_HPP
43 #define TPETRA_CRSGRAPH_DEF_HPP
44 
52 
53 #include "Tpetra_Distributor.hpp"
54 #include "Teuchos_Assert.hpp"
55 #include "Teuchos_NullIteratorTraits.hpp"
56 #include "Teuchos_as.hpp"
57 #include "Teuchos_SerialDenseMatrix.hpp"
58 #include "KokkosCompat_ClassicNodeAPI_Wrapper.hpp"
59 
60 #include <algorithm>
61 #include <string>
62 #include <utility>
63 
64 
65 namespace Tpetra {
66 
67  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
68  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
69  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
70  size_t maxNumEntriesPerRow,
71  ProfileType pftype,
72  const Teuchos::RCP<Teuchos::ParameterList>& params) :
73  dist_object_type (rowMap)
74  , rowMap_ (rowMap)
75  , nodeNumEntries_ (0)
76  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
77  , pftype_ (pftype)
78  , numAllocForAllRows_ (maxNumEntriesPerRow)
79  , storageStatus_ (pftype == StaticProfile ?
80  Details::STORAGE_1D_UNPACKED :
81  Details::STORAGE_2D)
82  , indicesAreAllocated_ (false)
83  , indicesAreLocal_ (false)
84  , indicesAreGlobal_ (false)
85  , fillComplete_ (false)
86  , indicesAreSorted_ (true)
87  , noRedundancies_ (true)
88  , haveLocalConstants_ (false)
89  , haveGlobalConstants_ (false)
90  , sortGhostsAssociatedWithEachProcessor_ (true)
91  {
92  const char tfecfFuncName[] = "CrsGraph(rowMap,maxNumEntriesPerRow,"
93  "pftype,params): ";
94  staticAssertions ();
95  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
96  maxNumEntriesPerRow == Teuchos::OrdinalTraits<size_t>::invalid (),
97  std::invalid_argument, "The allocation hint maxNumEntriesPerRow must be "
98  "a valid size_t value, which in this case means it must not be "
99  "Teuchos::OrdinalTraits<size_t>::invalid().");
100  resumeFill (params);
101  checkInternalState ();
102  }
103 
104  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
105  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
106  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
107  const Teuchos::RCP<const map_type>& colMap,
108  const size_t maxNumEntriesPerRow,
109  const ProfileType pftype,
110  const Teuchos::RCP<Teuchos::ParameterList>& params) :
111  dist_object_type (rowMap)
112  , rowMap_ (rowMap)
113  , colMap_ (colMap)
114  , nodeNumEntries_ (0)
115  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
116  , pftype_ (pftype)
117  , numAllocForAllRows_ (maxNumEntriesPerRow)
118  , storageStatus_ (pftype == StaticProfile ?
119  Details::STORAGE_1D_UNPACKED :
120  Details::STORAGE_2D)
121  , indicesAreAllocated_ (false)
122  , indicesAreLocal_ (false)
123  , indicesAreGlobal_ (false)
124  , fillComplete_ (false)
125  , indicesAreSorted_ (true)
126  , noRedundancies_ (true)
127  , haveLocalConstants_ (false)
128  , haveGlobalConstants_ (false)
129  , sortGhostsAssociatedWithEachProcessor_ (true)
130  {
131  const char tfecfFuncName[] = "CrsGraph(rowMap,colMap,maxNumEntriesPerRow,"
132  "pftype,params): ";
133  staticAssertions ();
134  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
135  maxNumEntriesPerRow == Teuchos::OrdinalTraits<size_t>::invalid (),
136  std::invalid_argument, "The allocation hint maxNumEntriesPerRow must be "
137  "a valid size_t value, which in this case means it must not be "
138  "Teuchos::OrdinalTraits<size_t>::invalid().");
139  resumeFill (params);
140  checkInternalState ();
141  }
142 
143  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
144  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
145  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
146  const Teuchos::ArrayRCP<const size_t>& numEntPerRow,
147  const ProfileType pftype,
148  const Teuchos::RCP<Teuchos::ParameterList>& params) :
149  dist_object_type (rowMap)
150  , rowMap_ (rowMap)
151  , nodeNumEntries_ (0)
152  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
153  , pftype_ (pftype)
154  , numAllocForAllRows_ (0)
155  , storageStatus_ (pftype == StaticProfile ?
156  Details::STORAGE_1D_UNPACKED :
157  Details::STORAGE_2D)
158  , indicesAreAllocated_ (false)
159  , indicesAreLocal_ (false)
160  , indicesAreGlobal_ (false)
161  , fillComplete_ (false)
162  , indicesAreSorted_ (true)
163  , noRedundancies_ (true)
164  , haveLocalConstants_ (false)
165  , haveGlobalConstants_ (false)
166  , sortGhostsAssociatedWithEachProcessor_ (true)
167  {
168  const char tfecfFuncName[] = "CrsGraph(rowMap,numEntPerRow,pftype,params): ";
169  staticAssertions ();
170 
171  const size_t lclNumRows = rowMap.is_null () ?
172  static_cast<size_t> (0) : rowMap->getNodeNumElements ();
173  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
174  static_cast<size_t> (numEntPerRow.size ()) != lclNumRows,
175  std::invalid_argument, "numEntPerRow has length " << numEntPerRow.size ()
176  << " != the local number of rows " << lclNumRows << " as specified by "
177  "the input row Map.");
178 
179 #ifdef HAVE_TPETRA_DEBUG
180  for (size_t r = 0; r < lclNumRows; ++r) {
181  const size_t curRowCount = numEntPerRow[r];
182  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
183  curRowCount == Teuchos::OrdinalTraits<size_t>::invalid (),
184  std::invalid_argument, "numEntPerRow(" << r << ") specifies an invalid "
185  "number of entries (Teuchos::OrdinalTraits<size_t>::invalid()).");
186  }
187 #endif // HAVE_TPETRA_DEBUG
188 
189  // Allocate k_numAllocPerRow_ (upper bound on number of entries
190  // per row). We get a host view of the data, as an ArrayRCP.
191  // Create a nonowning Kokkos::View of it, copy into
192  // k_numAllocPerRow, and sync. Then assign to k_numAllocPerRow_
193  // (which is a const view, so we can't copy into it directly).
194  typedef Kokkos::DualView<size_t*, Kokkos::LayoutLeft, execution_space>
195  dual_view_type;
196  typedef typename dual_view_type::host_mirror_space host_type;
197  typedef Kokkos::View<const size_t*, Kokkos::LayoutLeft, host_type,
198  Kokkos::MemoryUnmanaged> in_view_type;
199  in_view_type numAllocPerRowIn (numEntPerRow.getRawPtr (), lclNumRows);
200  dual_view_type k_numAllocPerRow ("Tpetra::CrsGraph::numAllocPerRow",
201  lclNumRows);
202  k_numAllocPerRow.template modify<host_type> ();
203  Kokkos::deep_copy (k_numAllocPerRow.h_view, numAllocPerRowIn);
204  k_numAllocPerRow.template sync<execution_space> ();
205  k_numAllocPerRow_ = k_numAllocPerRow;
206 
207  resumeFill (params);
208  checkInternalState ();
209  }
210 
211  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
212  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
213  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
214  const Kokkos::DualView<const size_t*, execution_space>& numEntPerRow,
215  const ProfileType pftype,
216  const Teuchos::RCP<Teuchos::ParameterList>& params) :
217  dist_object_type (rowMap)
218  , rowMap_ (rowMap)
219  , nodeNumEntries_ (0)
220  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
221  , pftype_ (pftype)
222  , k_numAllocPerRow_ (numEntPerRow)
223  , numAllocForAllRows_ (0)
224  , storageStatus_ (pftype == StaticProfile ?
225  Details::STORAGE_1D_UNPACKED :
226  Details::STORAGE_2D)
227  , indicesAreAllocated_ (false)
228  , indicesAreLocal_ (false)
229  , indicesAreGlobal_ (false)
230  , fillComplete_ (false)
231  , indicesAreSorted_ (true)
232  , noRedundancies_ (true)
233  , haveLocalConstants_ (false)
234  , haveGlobalConstants_ (false)
235  , sortGhostsAssociatedWithEachProcessor_ (true)
236  {
237  const char tfecfFuncName[] = "CrsGraph(rowMap,numEntPerRow,pftype,params): ";
238  staticAssertions ();
239 
240  const size_t lclNumRows = rowMap.is_null () ?
241  static_cast<size_t> (0) : rowMap->getNodeNumElements ();
242  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
243  static_cast<size_t> (numEntPerRow.dimension_0 ()) != lclNumRows,
244  std::invalid_argument, "numEntPerRow has length " <<
245  numEntPerRow.dimension_0 () << " != the local number of rows " <<
246  lclNumRows << " as specified by " "the input row Map.");
247 
248 #ifdef HAVE_TPETRA_DEBUG
249  for (size_t r = 0; r < lclNumRows; ++r) {
250  const size_t curRowCount = numEntPerRow.h_view(r);
251  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
252  curRowCount == Teuchos::OrdinalTraits<size_t>::invalid (),
253  std::invalid_argument, "numEntPerRow(" << r << ") specifies an invalid "
254  "number of entries (Teuchos::OrdinalTraits<size_t>::invalid()).");
255  }
256 #endif // HAVE_TPETRA_DEBUG
257 
258  resumeFill (params);
259  checkInternalState ();
260  }
261 
262 
263  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
264  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
265  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
266  const Teuchos::RCP<const map_type>& colMap,
267  const Kokkos::DualView<const size_t*, execution_space>& numEntPerRow,
268  const ProfileType pftype,
269  const Teuchos::RCP<Teuchos::ParameterList>& params) :
270  dist_object_type (rowMap)
271  , rowMap_ (rowMap)
272  , colMap_ (colMap)
273  , nodeNumEntries_ (0)
274  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
275  , pftype_ (pftype)
276  , k_numAllocPerRow_ (numEntPerRow)
277  , numAllocForAllRows_ (0)
278  , storageStatus_ (pftype == StaticProfile ?
279  Details::STORAGE_1D_UNPACKED :
280  Details::STORAGE_2D)
281  , indicesAreAllocated_ (false)
282  , indicesAreLocal_ (false)
283  , indicesAreGlobal_ (false)
284  , fillComplete_ (false)
285  , indicesAreSorted_ (true)
286  , noRedundancies_ (true)
287  , haveLocalConstants_ (false)
288  , haveGlobalConstants_ (false)
289  , sortGhostsAssociatedWithEachProcessor_ (true)
290  {
291  const char tfecfFuncName[] = "CrsGraph(rowMap,colMap,numEntPerRow,pftype,params): ";
292  staticAssertions ();
293 
294  const size_t lclNumRows = rowMap.is_null () ?
295  static_cast<size_t> (0) : rowMap->getNodeNumElements ();
296  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
297  static_cast<size_t> (numEntPerRow.dimension_0 ()) != lclNumRows,
298  std::invalid_argument, "numEntPerRow has length " <<
299  numEntPerRow.dimension_0 () << " != the local number of rows " <<
300  lclNumRows << " as specified by " "the input row Map.");
301 
302 #ifdef HAVE_TPETRA_DEBUG
303  for (size_t r = 0; r < lclNumRows; ++r) {
304  const size_t curRowCount = numEntPerRow.h_view(r);
305  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
306  curRowCount == Teuchos::OrdinalTraits<size_t>::invalid (),
307  std::invalid_argument, "numEntPerRow(" << r << ") specifies an invalid "
308  "number of entries (Teuchos::OrdinalTraits<size_t>::invalid()).");
309  }
310 #endif // HAVE_TPETRA_DEBUG
311 
312  resumeFill (params);
313  checkInternalState ();
314  }
315 
316 
317  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
318  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
319  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
320  const Teuchos::RCP<const map_type>& colMap,
321  const Teuchos::ArrayRCP<const size_t>& numEntPerRow,
322  ProfileType pftype,
323  const Teuchos::RCP<Teuchos::ParameterList>& params) :
324  dist_object_type (rowMap)
325  , rowMap_ (rowMap)
326  , colMap_ (colMap)
327  , nodeNumEntries_ (0)
328  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
329  , pftype_ (pftype)
330  , numAllocForAllRows_ (0)
331  , storageStatus_ (pftype == StaticProfile ?
332  Details::STORAGE_1D_UNPACKED :
333  Details::STORAGE_2D)
334  , indicesAreAllocated_ (false)
335  , indicesAreLocal_ (false)
336  , indicesAreGlobal_ (false)
337  , fillComplete_ (false)
338  , indicesAreSorted_ (true)
339  , noRedundancies_ (true)
340  , haveLocalConstants_ (false)
341  , haveGlobalConstants_ (false)
342  , sortGhostsAssociatedWithEachProcessor_ (true)
343  {
344  const char tfecfFuncName[] = "CrsGraph(rowMap,colMap,numEntPerRow,pftype,"
345  "params): ";
346  staticAssertions ();
347 
348  const size_t lclNumRows = rowMap.is_null () ?
349  static_cast<size_t> (0) : rowMap->getNodeNumElements ();
350  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
351  static_cast<size_t> (numEntPerRow.size ()) != lclNumRows,
352  std::invalid_argument, "numEntPerRow has length " << numEntPerRow.size ()
353  << " != the local number of rows " << lclNumRows << " as specified by "
354  "the input row Map.");
355 
356 #ifdef HAVE_TPETRA_DEBUG
357  for (size_t r = 0; r < lclNumRows; ++r) {
358  const size_t curRowCount = numEntPerRow[r];
359  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
360  curRowCount == Teuchos::OrdinalTraits<size_t>::invalid (),
361  std::invalid_argument, "numEntPerRow(" << r << ") specifies an invalid "
362  "number of entries (Teuchos::OrdinalTraits<size_t>::invalid()).");
363  }
364 #endif // HAVE_TPETRA_DEBUG
365 
366  // Allocate k_numAllocPerRow_ (upper bound on number of entries
367  // per row). We get a host view of the data, as an ArrayRCP.
368  // Create a nonowning Kokkos::View of it, copy into
369  // k_numAllocPerRow, and sync. Then assign to k_numAllocPerRow_
370  // (which is a const view, so we can't copy into it directly).
371  typedef Kokkos::DualView<size_t*, Kokkos::LayoutLeft, execution_space>
372  dual_view_type;
373  typedef typename dual_view_type::host_mirror_space host_type;
374  typedef Kokkos::View<const size_t*, Kokkos::LayoutLeft, host_type,
375  Kokkos::MemoryUnmanaged> in_view_type;
376  in_view_type numAllocPerRowIn (numEntPerRow.getRawPtr (), lclNumRows);
377  dual_view_type k_numAllocPerRow ("Tpetra::CrsGraph::numAllocPerRow",
378  lclNumRows);
379  k_numAllocPerRow.template modify<host_type> ();
380  Kokkos::deep_copy (k_numAllocPerRow.h_view, numAllocPerRowIn);
381  k_numAllocPerRow.template sync<execution_space> ();
382  k_numAllocPerRow_ = k_numAllocPerRow;
383 
384  resumeFill (params);
385  checkInternalState ();
386  }
387 
388 
389  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
390  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
391  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
392  const Teuchos::RCP<const map_type>& colMap,
393  const typename local_graph_type::row_map_type& rowPointers,
394  const typename local_graph_type::entries_type::non_const_type& columnIndices,
395  const Teuchos::RCP<Teuchos::ParameterList>& params) :
396  dist_object_type (rowMap)
397  , rowMap_(rowMap)
398  , colMap_(colMap)
399  , globalNumEntries_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
400  , globalNumDiags_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
401  , globalMaxNumRowEntries_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
402  , nodeNumEntries_(0)
403  , nodeNumAllocated_(OrdinalTraits<size_t>::invalid())
404  , pftype_(StaticProfile)
405  , numAllocForAllRows_(0)
406  , storageStatus_ (Details::STORAGE_1D_PACKED)
407  , indicesAreAllocated_(true)
408  , indicesAreLocal_(true)
409  , indicesAreGlobal_(false)
410  , fillComplete_(false)
411  , indicesAreSorted_(true)
412  , noRedundancies_(true)
413  , haveLocalConstants_ (false)
414  , haveGlobalConstants_ (false)
415  , sortGhostsAssociatedWithEachProcessor_(true)
416  {
417  staticAssertions ();
418  setAllIndices (rowPointers, columnIndices);
419  checkInternalState ();
420  }
421 
422 
423  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
424  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
425  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
426  const Teuchos::RCP<const map_type>& colMap,
427  const Teuchos::ArrayRCP<size_t>& rowPointers,
428  const Teuchos::ArrayRCP<LocalOrdinal> & columnIndices,
429  const Teuchos::RCP<Teuchos::ParameterList>& params) :
430  dist_object_type (rowMap)
431  , rowMap_ (rowMap)
432  , colMap_ (colMap)
433  , globalNumEntries_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
434  , globalNumDiags_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
435  , globalMaxNumRowEntries_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
436  , nodeNumEntries_ (0)
437  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
438  , pftype_ (StaticProfile)
439  , numAllocForAllRows_ (0)
440  , storageStatus_ (Details::STORAGE_1D_PACKED)
441  , indicesAreAllocated_ (true)
442  , indicesAreLocal_ (true)
443  , indicesAreGlobal_ (false)
444  , fillComplete_ (false)
445  , indicesAreSorted_ (true)
446  , noRedundancies_ (true)
447  , haveLocalConstants_ (false)
448  , haveGlobalConstants_ (false)
449  , sortGhostsAssociatedWithEachProcessor_ (true)
450  {
451  staticAssertions ();
452  setAllIndices (rowPointers, columnIndices);
453  checkInternalState ();
454  }
455 
456 
457  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
458  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
459  CrsGraph (const Teuchos::RCP<const map_type>& rowMap,
460  const Teuchos::RCP<const map_type>& colMap,
461  const local_graph_type& k_local_graph_,
462  const Teuchos::RCP<Teuchos::ParameterList>& params)
463  : DistObject<GlobalOrdinal, LocalOrdinal, GlobalOrdinal, node_type> (rowMap)
464  , rowMap_ (rowMap)
465  , colMap_ (colMap)
466  , lclGraph_ (k_local_graph_)
467  , globalNumEntries_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
468  , globalNumDiags_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
469  , globalMaxNumRowEntries_ (Teuchos::OrdinalTraits<global_size_t>::invalid ())
470  , nodeNumEntries_ (0) // FIXME (mfh 17 Mar 2014) should get from lclGraph_ right now
471  , nodeNumAllocated_ (Teuchos::OrdinalTraits<size_t>::invalid ())
472  , pftype_ (StaticProfile)
473  , numAllocForAllRows_ (0)
474  , storageStatus_ (Details::STORAGE_1D_PACKED)
475  , indicesAreAllocated_ (true)
476  , indicesAreLocal_ (true)
477  , indicesAreGlobal_ (false)
478  , fillComplete_ (false)
479  , indicesAreSorted_ (true)
480  , noRedundancies_ (true)
481  , haveLocalConstants_ (false)
482  , haveGlobalConstants_ (false)
483  , sortGhostsAssociatedWithEachProcessor_(true)
484  {
485  using Teuchos::arcp;
486  using Teuchos::ArrayRCP;
487  using Teuchos::as;
488  using Teuchos::ParameterList;
489  using Teuchos::parameterList;
490  using Teuchos::rcp;
491  typedef GlobalOrdinal GO;
492  typedef LocalOrdinal LO;
493 
494  staticAssertions();
495  const char tfecfFuncName[] = "CrsGraph(Map,Map,Kokkos::LocalStaticCrsGraph)";
496 
497  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
498  colMap.is_null (), std::runtime_error,
499  ": The input column Map must be nonnull.");
500  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
501  k_local_graph_.numRows () != rowMap->getNodeNumElements (),
502  std::runtime_error,
503  ": The input row Map and the input local graph need to have the same "
504  "number of rows. The row Map claims " << rowMap->getNodeNumElements ()
505  << " row(s), but the local graph claims " << k_local_graph_.numRows ()
506  << " row(s).");
507  // NOTE (mfh 17 Mar 2014) getNodeNumRows() returns
508  // rowMap_->getNodeNumElements(), but it doesn't have to.
509  // TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
510  // k_local_graph_.numRows () != getNodeNumRows (), std::runtime_error,
511  // ": The input row Map and the input local graph need to have the same "
512  // "number of rows. The row Map claims " << getNodeNumRows () << " row(s), "
513  // "but the local graph claims " << k_local_graph_.numRows () << " row(s).");
514  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
515  k_lclInds1D_.dimension_0 () != 0 || k_gblInds1D_.dimension_0 () != 0, std::logic_error,
516  ": cannot have 1D data structures allocated.");
517  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
518  ! lclInds2D_.is_null () || ! gblInds2D_.is_null (), std::logic_error,
519  ": cannot have 2D data structures allocated.");
520 
521  nodeNumAllocated_ = k_local_graph_.row_map (getNodeNumRows ());
522  nodeNumEntries_ = k_local_graph_.row_map (getNodeNumRows ());
523 
524  // NOTE (mfh 17 Mar 2014) We also need a version of this CrsGraph
525  // constructor that takes a domain and range Map, as well as a row
526  // and column Map. In that case, we must pass the domain and
527  // range Map into the following method.
528  setDomainRangeMaps (rowMap_, rowMap_);
529  makeImportExport ();
530 
531  k_lclInds1D_ = lclGraph_.entries;
532  k_rowPtrs_ = lclGraph_.row_map;
533 
534  typename local_graph_type::row_map_type d_ptrs = lclGraph_.row_map;
535  typename local_graph_type::entries_type d_inds = lclGraph_.entries;
536 
537  // Reset local properties
538  upperTriangular_ = true;
539  lowerTriangular_ = true;
540  nodeMaxNumRowEntries_ = 0;
541  nodeNumDiags_ = 0;
542 
543  // Compute triangular properties
544  const size_t numLocalRows = getNodeNumRows ();
545  for (size_t localRow = 0; localRow < numLocalRows; ++localRow) {
546  const GO globalRow = rowMap_->getGlobalElement (localRow);
547  const LO rlcid = colMap_->getLocalElement (globalRow);
548 
549  // It's entirely possible that the local matrix has no entries
550  // in the column corresponding to the current row. In that
551  // case, the column Map may not necessarily contain that GID.
552  // This is why we check whether rlcid is "invalid" (which means
553  // that globalRow is not a GID in the column Map).
554  if (rlcid != Teuchos::OrdinalTraits<LO>::invalid ()) {
555  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
556  rlcid + 1 >= static_cast<LO> (d_ptrs.dimension_0 ()),
557  std::runtime_error, ": The given row Map and/or column Map is/are "
558  "not compatible with the provided local graph.");
559  if (d_ptrs(rlcid) != d_ptrs(rlcid + 1)) {
560  const size_t smallestCol =
561  static_cast<size_t> (d_inds(d_ptrs(rlcid)));
562  const size_t largestCol =
563  static_cast<size_t> (d_inds(d_ptrs(rlcid + 1)-1));
564  if (smallestCol < localRow) {
565  upperTriangular_ = false;
566  }
567  if (localRow < largestCol) {
568  lowerTriangular_ = false;
569  }
570  for (size_t i = d_ptrs(rlcid); i < d_ptrs(rlcid + 1); ++i) {
571  if (d_inds(i) == rlcid) {
572  ++nodeNumDiags_;
573  }
574  }
575  }
576  nodeMaxNumRowEntries_ =
577  std::max (static_cast<size_t> (d_ptrs(rlcid + 1) - d_ptrs(rlcid)),
578  nodeMaxNumRowEntries_);
579  }
580  }
581 
582  haveLocalConstants_ = true;
583  computeGlobalConstants ();
584 
585  fillComplete_ = true;
586  checkInternalState ();
587  }
588 
589 
590  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
591  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
592  ~CrsGraph ()
593  {}
594 
595 
596  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
597  Teuchos::RCP<const Teuchos::ParameterList>
598  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
599  getValidParameters () const
600  {
601  using Teuchos::RCP;
602  using Teuchos::ParameterList;
603  using Teuchos::parameterList;
604 
605  RCP<ParameterList> params = parameterList ("Tpetra::CrsGraph");
606 
607  // Make a sublist for the Import.
608  RCP<ParameterList> importSublist = parameterList ("Import");
609 
610  // FIXME (mfh 02 Apr 2012) We should really have the Import and
611  // Export objects fill in these lists. However, we don't want to
612  // create an Import or Export unless we need them. For now, we
613  // know that the Import and Export just pass the list directly to
614  // their Distributor, so we can create a Distributor here
615  // (Distributor's constructor is a lightweight operation) and have
616  // it fill in the list.
617 
618  // Fill in Distributor default parameters by creating a
619  // Distributor and asking it to do the work.
620  Distributor distributor (rowMap_->getComm (), importSublist);
621  params->set ("Import", *importSublist, "How the Import performs communication.");
622 
623  // Make a sublist for the Export. For now, it's a clone of the
624  // Import sublist. It's not a shallow copy, though, since we
625  // might like the Import to do communication differently than the
626  // Export.
627  params->set ("Export", *importSublist, "How the Export performs communication.");
628 
629  return params;
630  }
631 
632 
633  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
634  void
635  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
636  setParameterList (const Teuchos::RCP<Teuchos::ParameterList>& params)
637  {
638  Teuchos::RCP<const Teuchos::ParameterList> validParams =
639  getValidParameters ();
640  params->validateParametersAndSetDefaults (*validParams);
641  this->setMyParamList (params);
642  }
643 
644 
645  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
646  global_size_t
647  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
648  getGlobalNumRows () const
649  {
650  return rowMap_->getGlobalNumElements ();
651  }
652 
653 
654  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
655  global_size_t
656  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
657  getGlobalNumCols () const
658  {
659  const char tfecfFuncName[] = "getGlobalNumCols: ";
660  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
661  ! isFillComplete () || getDomainMap ().is_null (), std::runtime_error,
662  "The graph does not have a domain Map. You may not call this method in "
663  "that case.");
664  return getDomainMap ()->getGlobalNumElements ();
665  }
666 
667 
668  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
669  size_t
670  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
671  getNodeNumRows () const
672  {
673  return rowMap_.is_null () ? static_cast<size_t> (0) :
674  rowMap_->getNodeNumElements ();
675  }
676 
677 
678  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
679  size_t
680  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
681  getNodeNumCols () const
682  {
683  const char tfecfFuncName[] = "getNodeNumCols: ";
684  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
685  ! hasColMap (), std::runtime_error,
686  "The graph does not have a column Map. You may not call this method "
687  "unless the graph has a column Map. This requires either that a custom "
688  "column Map was given to the constructor, or that fillComplete() has "
689  "been called.");
690  return colMap_.is_null () ? static_cast<size_t> (0) :
691  colMap_->getNodeNumElements ();
692  }
693 
694 
695  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
696  size_t
697  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
698  getNodeNumDiags () const
699  {
700  return nodeNumDiags_;
701  }
702 
703 
704  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
705  global_size_t
706  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
707  getGlobalNumDiags () const
708  {
709  return globalNumDiags_;
710  }
711 
712 
713  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
714  Teuchos::RCP<Node>
715  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
716  getNode () const
717  {
718  return rowMap_.is_null () ? Teuchos::null : rowMap_->getNode ();
719  }
720 
721 
722  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
723  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
724  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
725  getRowMap () const
726  {
727  return rowMap_;
728  }
729 
730 
731  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
732  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
733  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
734  getColMap () const
735  {
736  return colMap_;
737  }
738 
739 
740  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
741  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
742  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
743  getDomainMap () const
744  {
745  return domainMap_;
746  }
747 
748 
749  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
750  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
751  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
752  getRangeMap () const
753  {
754  return rangeMap_;
755  }
756 
757  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
758  Teuchos::RCP<const Import<LocalOrdinal, GlobalOrdinal, Node> >
759  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
760  getImporter () const
761  {
762  return importer_;
763  }
764 
765 
766  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
767  Teuchos::RCP<const Export<LocalOrdinal, GlobalOrdinal, Node> >
768  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
769  getExporter () const
770  {
771  return exporter_;
772  }
773 
774 
775  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
776  bool
777  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
778  hasColMap () const
779  {
780  return ! colMap_.is_null ();
781  }
782 
783 
784  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
785  bool
786  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
787  isStorageOptimized () const
788  {
789  // FIXME (mfh 07 Aug 2014) Why wouldn't storage be optimized if
790  // getNodeNumRows() is zero?
791 
792  const bool isOpt = indicesAreAllocated_ &&
793  k_numRowEntries_.dimension_0 () == 0 &&
794  getNodeNumRows () > 0;
795 
796 #ifdef HAVE_TPETRA_DEBUG
797  const char tfecfFuncName[] = "isStorageOptimized";
798  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
799  isOpt && getProfileType () == DynamicProfile, std::logic_error,
800  ": The matrix claims to have optimized storage, but getProfileType() "
801  "returns DynamicProfile. This should never happen. Please report this "
802  "bug to the Tpetra developers.");
803 #endif // HAVE_TPETRA_DEBUG
804 
805  return isOpt;
806  }
807 
808 
809  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
810  ProfileType
811  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
812  getProfileType () const
813  {
814  return pftype_;
815  }
816 
817 
818  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
819  global_size_t
820  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
821  getGlobalNumEntries () const
822  {
823  return globalNumEntries_;
824  }
825 
826 
827  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
828  size_t
829  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
830  getNodeNumEntries () const
831  {
832  return nodeNumEntries_;
833  }
834 
835 
836  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
837  global_size_t
838  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
839  getGlobalMaxNumRowEntries () const
840  {
841  return globalMaxNumRowEntries_;
842  }
843 
844 
845  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
846  size_t
847  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
848  getNodeMaxNumRowEntries () const
849  {
850  return nodeMaxNumRowEntries_;
851  }
852 
853 
854  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
855  bool
856  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
857  isFillComplete () const
858  {
859  return fillComplete_;
860  }
861 
862 
863  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
864  bool
865  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
866  isFillActive () const
867  {
868  return ! fillComplete_;
869  }
870 
871 
872  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
873  bool
874  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
875  isUpperTriangular () const
876  {
877  return upperTriangular_;
878  }
879 
880 
881  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
882  bool
883  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
884  isLowerTriangular () const
885  {
886  return lowerTriangular_;
887  }
888 
889 
890  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
891  bool
892  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
893  isLocallyIndexed () const
894  {
895  return indicesAreLocal_;
896  }
897 
898 
899  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
900  bool
901  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
902  isGloballyIndexed () const
903  {
904  return indicesAreGlobal_;
905  }
906 
907 
908  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
909  size_t
910  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
911  getNodeAllocationSize () const
912  {
913  return nodeNumAllocated_;
914  }
915 
916 
917  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
918  Teuchos::RCP<const Teuchos::Comm<int> >
919  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
920  getComm () const
921  {
922  return rowMap_.is_null () ? Teuchos::null : rowMap_->getComm ();
923  }
924 
925 
926  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
927  GlobalOrdinal
928  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
929  getIndexBase () const
930  {
931  return rowMap_->getIndexBase ();
932  }
933 
934 
935  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
936  bool
937  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
938  indicesAreAllocated () const
939  {
940  return indicesAreAllocated_;
941  }
942 
943 
944  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
945  bool
946  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
947  isSorted () const
948  {
949  return indicesAreSorted_;
950  }
951 
952 
953  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
954  bool
955  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
956  isMerged () const
957  {
958  return noRedundancies_;
959  }
960 
961 
962  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
963  void
964  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
965  setLocallyModified ()
966  {
967  // FIXME (mfh 07 May 2013) How do we know that the change
968  // introduced a redundancy, or even that it invalidated the sorted
969  // order of indices? CrsGraph has always made this conservative
970  // guess. It could be a bit costly to check at insertion time,
971  // though.
972  indicesAreSorted_ = false;
973  noRedundancies_ = false;
974 
975  // We've modified the graph, so we'll have to recompute local
976  // constants like the number of diagonal entries on this process.
977  haveLocalConstants_ = false;
978  }
979 
980 
981  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
982  void
983  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
984  allocateIndices (const ELocalGlobal lg)
985  {
986  using Teuchos::arcp;
987  using Teuchos::Array;
988  using Teuchos::ArrayRCP;
989  typedef Teuchos::ArrayRCP<size_t>::size_type size_type;
990  typedef typename local_graph_type::row_map_type::non_const_type
991  non_const_row_map_type;
992  typedef typename local_graph_type::entries_type::non_const_type
993  lcl_col_inds_type;
994  typedef Kokkos::View<GlobalOrdinal*,
995  typename lcl_col_inds_type::array_layout,
996  device_type> gbl_col_inds_type;
997  const char tfecfFuncName[] = "allocateIndices: ";
998 
999  // This is a protected function, only callable by us. If it was
1000  // called incorrectly, it is our fault. That's why the tests
1001  // below throw std::logic_error instead of std::invalid_argument.
1002  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1003  isLocallyIndexed () && lg == GlobalIndices, std::logic_error,
1004  "The graph is locally indexed, but Tpetra code is calling this method "
1005  "with lg=GlobalIndices. Please report this bug to the Tpetra developers.");
1006  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1007  isGloballyIndexed () && lg == LocalIndices, std::logic_error,
1008  "The graph is globally indexed, but Tpetra code is calling this method "
1009  "with lg=LocalIndices. Please report this bug to the Tpetra developers.");
1010  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1011  indicesAreAllocated (), std::logic_error, "The graph's indices are "
1012  "already allocated, but Tpetra code is calling allocateIndices) again. "
1013  "Please report this bug to the Tpetra developers.");
1014 
1015  const size_t numRows = getNodeNumRows ();
1016 
1017  if (getProfileType () == StaticProfile) {
1018  //
1019  // STATIC ALLOCATION PROFILE
1020  //
1021  non_const_row_map_type k_rowPtrs ("Tpetra::CrsGraph::ptr", numRows + 1);
1022 
1023  if (k_numAllocPerRow_.dimension_0 () != 0) {
1024  // It's OK to throw std::invalid_argument here, because we
1025  // haven't incurred any side effects yet. Throwing that
1026  // exception (and not, say, std::logic_error) implies that the
1027  // instance can recover.
1028  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1029  k_numAllocPerRow_.dimension_0 () != numRows, std::invalid_argument,
1030  "k_numAllocPerRow_ is allocated (has length != 0), but its length = "
1031  << k_numAllocPerRow_.dimension_0 () << " != numRows = " << numRows
1032  << ".");
1033  // FIXME hack until we get parallel_scan in kokkos
1034  //
1035  // FIXME (mfh 11 Aug 2014) This assumes UVM, since k_rowPtrs_
1036  // is currently a device View. Should instead use a DualView.
1037  typename Kokkos::DualView<const size_t*, execution_space>::t_host h_numAllocPerRow =
1038  k_numAllocPerRow_.h_view; // use a host view for now, since we compute on host
1039  bool anyInvalidAllocSizes = false;
1040  for (size_t i = 0; i < numRows; ++i) {
1041  size_t allocSize = h_numAllocPerRow(i);
1042  if (allocSize == Teuchos::OrdinalTraits<size_t>::invalid ()) {
1043  anyInvalidAllocSizes = true;
1044  allocSize = 0;
1045  }
1046  k_rowPtrs(i+1) = k_rowPtrs(i) + allocSize;
1047  }
1048  // It's OK to throw std::invalid_argument here, because we
1049  // haven't incurred any side effects yet. Throwing that
1050  // exception (and not, say, std::logic_error) implies that the
1051  // instance can recover.
1052  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1053  anyInvalidAllocSizes, std::invalid_argument, "The input array of "
1054  "allocation sizes per row had at least one invalid (== "
1055  "Teuchos::OrdinalTraits<size_t>::invalid()) entry.");
1056  }
1057  else {
1058  // It's OK to throw std::invalid_argument here, because we
1059  // haven't incurred any side effects yet. Throwing that
1060  // exception (and not, say, std::logic_error) implies that the
1061  // instance can recover.
1062  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1063  numAllocForAllRows_ == Teuchos::OrdinalTraits<size_t>::invalid (),
1064  std::invalid_argument, "numAllocForAllRows_ has an invalid value, "
1065  "namely Teuchos::OrdinalTraits<size_t>::invalid() = " <<
1066  Teuchos::OrdinalTraits<size_t>::invalid () << ".");
1067  // FIXME hack until we get parallel_scan in kokkos
1068  //
1069  // FIXME (mfh 11 Aug 2014) This assumes UVM, since k_rowPtrs_
1070  // is currently a device View. Should instead use a DualView.
1071  for (size_t i = 0; i < numRows; ++i) {
1072  k_rowPtrs(i+1) = k_rowPtrs(i) + numAllocForAllRows_;
1073  }
1074  }
1075 
1076  // "Commit" the resulting row offsets.
1077  k_rowPtrs_ = k_rowPtrs;
1078 
1079  // FIXME (mfh 05,11 Aug 2014) This assumes UVM, since k_rowPtrs_
1080  // is currently a device View. Should instead use a DualView.
1081  const size_type numInds = static_cast<size_type> (k_rowPtrs_(numRows));
1082  if (lg == LocalIndices) {
1083  k_lclInds1D_ = lcl_col_inds_type ("Tpetra::CrsGraph::ind", numInds);
1084  }
1085  else {
1086  k_gblInds1D_ = gbl_col_inds_type ("Tpetra::CrsGraph::ind", numInds);
1087  }
1088  nodeNumAllocated_ = numInds;
1089  storageStatus_ = Details::STORAGE_1D_UNPACKED;
1090  }
1091  else {
1092  //
1093  // DYNAMIC ALLOCATION PROFILE
1094  //
1095 
1096  // Use the host view of k_numAllocPerRow_, since we have to
1097  // allocate 2-D storage on the host.
1098  typename Kokkos::DualView<const size_t*, execution_space>::t_host h_numAllocPerRow =
1099  k_numAllocPerRow_.h_view;
1100  const bool useNumAllocPerRow = (k_numAllocPerRow_.dimension_0 () != 0);
1101 
1102  if (lg == LocalIndices) {
1103  lclInds2D_ = arcp<Array<LocalOrdinal> > (numRows);
1104  nodeNumAllocated_ = 0;
1105  for (size_t i = 0; i < numRows; ++i) {
1106  const size_t howMany = useNumAllocPerRow ?
1107  h_numAllocPerRow(i) : numAllocForAllRows_;
1108  nodeNumAllocated_ += howMany;
1109  if (howMany > 0) {
1110  lclInds2D_[i].resize (howMany);
1111  }
1112  }
1113  }
1114  else { // allocate global indices
1115  gblInds2D_ = arcp<Array<GlobalOrdinal> > (numRows);
1116  nodeNumAllocated_ = 0;
1117  for (size_t i = 0; i < numRows; ++i) {
1118  const size_t howMany = useNumAllocPerRow ?
1119  h_numAllocPerRow(i) : numAllocForAllRows_;
1120  nodeNumAllocated_ += howMany;
1121  if (howMany > 0) {
1122  gblInds2D_[i].resize (howMany);
1123  }
1124  }
1125  }
1126  storageStatus_ = Details::STORAGE_2D;
1127  }
1128 
1129  indicesAreLocal_ = (lg == LocalIndices);
1130  indicesAreGlobal_ = (lg == GlobalIndices);
1131 
1132  if (numRows > 0) {
1133  k_numRowEntries_ =
1134  t_numRowEntries_ ("Tpetra::CrsGraph::numRowEntries", numRows);
1135 
1136  // Fill with zeros on the host. k_numRowEntries_ is a DualView.
1137  //
1138  // TODO (mfh 05 Aug 2014) Write a device kernel for this.
1139  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
1140  size_t* const hostPtr = k_numRowEntries_.h_view.ptr_on_device ();
1141  std::fill (hostPtr, hostPtr + numRows, static_cast<size_t> (0));
1142  k_numRowEntries_.template sync<execution_space> ();
1143  }
1144 
1145  // done with these
1146  numAllocForAllRows_ = 0;
1147  k_numAllocPerRow_ = Kokkos::DualView<const size_t*, execution_space> ();
1148  indicesAreAllocated_ = true;
1149  checkInternalState ();
1150  }
1151 
1152 
1153  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1154  template <class T>
1155  Teuchos::ArrayRCP<Teuchos::Array<T> >
1156  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1157  allocateValues2D () const
1158  {
1159  using Teuchos::arcp;
1160  using Teuchos::Array;
1161  using Teuchos::ArrayRCP;
1162  const char tfecfFuncName[] = "allocateValues2D: ";
1163  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1164  ! indicesAreAllocated (), std::runtime_error,
1165  "Graph indices must be allocated before values.");
1166  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1167  getProfileType () != DynamicProfile, std::runtime_error,
1168  "Graph indices must be allocated in a dynamic profile.");
1169 
1170  ArrayRCP<Array<T> > values2D;
1171  values2D = arcp<Array<T> > (getNodeNumRows ());
1172  if (lclInds2D_ != null) {
1173  const size_t numRows = lclInds2D_.size ();
1174  for (size_t r = 0; r < numRows; ++r) {
1175  values2D[r].resize (lclInds2D_[r].size ());
1176  }
1177  }
1178  else if (gblInds2D_ != null) {
1179  const size_t numRows = gblInds2D_.size ();
1180  for (size_t r = 0; r < numRows; ++r) {
1181  values2D[r].resize (gblInds2D_[r].size ());
1182  }
1183  }
1184  return values2D;
1185  }
1186 
1187 
1188  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1189  Teuchos::ArrayView<const LocalOrdinal>
1190  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1191  getLocalView (const RowInfo rowinfo) const
1192  {
1193  using Kokkos::subview;
1194  using Kokkos::View;
1195  typedef LocalOrdinal LO;
1196  typedef View<const LO*, execution_space, Kokkos::MemoryUnmanaged> row_view_type;
1197 
1198  if (rowinfo.allocSize == 0) {
1199  return Teuchos::ArrayView<const LO> ();
1200  }
1201  else { // nothing in the row to view
1202  if (k_lclInds1D_.dimension_0 () != 0) { // 1-D storage
1203  const size_t start = rowinfo.offset1D;
1204  const size_t len = rowinfo.allocSize;
1205  const std::pair<size_t, size_t> rng (start, start + len);
1206  row_view_type rowView = subview (k_lclInds1D_, rng);
1207 
1208  const LO* const rowViewRaw = (len == 0) ? NULL : rowView.ptr_on_device ();
1209  return Teuchos::ArrayView<const LO> (rowViewRaw, len, Teuchos::RCP_DISABLE_NODE_LOOKUP);
1210  }
1211  else if (! lclInds2D_[rowinfo.localRow].empty ()) { // 2-D storage
1212  return lclInds2D_[rowinfo.localRow] ();
1213  }
1214  else {
1215  return Teuchos::ArrayView<const LO> (); // nothing in the row to view
1216  }
1217  }
1218  }
1219 
1220 
1221  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1222  Teuchos::ArrayView<LocalOrdinal>
1223  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1224  getLocalViewNonConst (const RowInfo rowinfo)
1225  {
1226  using Kokkos::subview;
1227  using Kokkos::View;
1228  typedef LocalOrdinal LO;
1229  typedef View<LO*, execution_space, Kokkos::MemoryUnmanaged> row_view_type;
1230 
1231  if (rowinfo.allocSize == 0) { // nothing in the row to view
1232  return Teuchos::ArrayView<LO> ();
1233  }
1234  else {
1235  if (k_lclInds1D_.dimension_0 () != 0) { // 1-D storage
1236  const size_t start = rowinfo.offset1D;
1237  const size_t len = rowinfo.allocSize;
1238  const std::pair<size_t, size_t> rng (start, start + len);
1239  row_view_type rowView = subview (k_lclInds1D_, rng);
1240 
1241  LO* const rowViewRaw = (len == 0) ? NULL : rowView.ptr_on_device ();
1242  return Teuchos::ArrayView<LO> (rowViewRaw, len, Teuchos::RCP_DISABLE_NODE_LOOKUP);
1243  }
1244  else if (! lclInds2D_[rowinfo.localRow].empty ()) { // 2-D storage
1245  return lclInds2D_[rowinfo.localRow] ();
1246  }
1247  else {
1248  return Teuchos::ArrayView<LO> (); // nothing in the row to view
1249  }
1250  }
1251  }
1252 
1253 
1254  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1255  Teuchos::ArrayView<const GlobalOrdinal>
1256  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1257  getGlobalView (const RowInfo rowinfo) const
1258  {
1259  Teuchos::ArrayView<const GlobalOrdinal> view;
1260  if (rowinfo.allocSize > 0) {
1261  if (k_gblInds1D_.dimension_0 () != 0) {
1262  auto rng = std::make_pair (rowinfo.offset1D, rowinfo.offset1D + rowinfo.allocSize);
1263  view = Kokkos::Compat::getConstArrayView (Kokkos::subview (k_gblInds1D_, rng));
1264  }
1265  else if (! gblInds2D_[rowinfo.localRow].empty()) {
1266  view = gblInds2D_[rowinfo.localRow] ();
1267  }
1268  }
1269  return view;
1270  }
1271 
1272 
1273  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1274  Teuchos::ArrayView<GlobalOrdinal>
1275  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1276  getGlobalViewNonConst (const RowInfo rowinfo)
1277  {
1278  Teuchos::ArrayView<GlobalOrdinal> view;
1279  if (rowinfo.allocSize > 0) {
1280  if (k_gblInds1D_.dimension_0 () != 0) {
1281  auto rng = std::make_pair (rowinfo.offset1D, rowinfo.offset1D + rowinfo.allocSize);
1282  view = Kokkos::Compat::getArrayView (Kokkos::subview (k_gblInds1D_, rng));
1283  }
1284  else if (! gblInds2D_[rowinfo.localRow].empty()) {
1285  view = gblInds2D_[rowinfo.localRow] ();
1286  }
1287  }
1288  return view;
1289  }
1290 
1291 
1292  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1293  RowInfo
1294  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1295  getRowInfo (const size_t myRow) const
1296  {
1297 #ifdef HAVE_TPETRA_DEBUG
1298  const char tfecfFuncName[] = "getRowInfo";
1299  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1300  ! rowMap_->isNodeLocalElement (myRow), std::logic_error,
1301  ": The given (local) row index myRow = " << myRow
1302  << " does not belong to the graph's row Map. "
1303  "This probably indicates a bug in Tpetra::CrsGraph or Tpetra::CrsMatrix. "
1304  "Please report this bug to the Tpetra developers.");
1305  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1306  ! hasRowInfo (), std::logic_error,
1307  ": Late catch! Graph does not have row info anymore. "
1308  "Error should have been caught earlier. "
1309  "Please report this bug to the Tpetra developers.");
1310 #endif // HAVE_TPETRA_DEBUG
1311 
1312  const size_t STINV = Teuchos::OrdinalTraits<size_t>::invalid ();
1313  RowInfo ret;
1314  ret.localRow = myRow;
1315  if (nodeNumAllocated_ != 0 && nodeNumAllocated_ != STINV) {
1316  // graph data structures have the info that we need
1317  //
1318  // if static graph, offsets tell us the allocation size
1319  if (getProfileType() == StaticProfile) {
1320  ret.offset1D = k_rowPtrs_(myRow);
1321  ret.allocSize = k_rowPtrs_(myRow+1) - k_rowPtrs_(myRow);
1322  if (k_numRowEntries_.dimension_0 () == 0) {
1323  ret.numEntries = ret.allocSize;
1324  } else {
1325  ret.numEntries = k_numRowEntries_.h_view(myRow);
1326  }
1327  }
1328  else {
1329  ret.offset1D = STINV;
1330  if (isLocallyIndexed ()) {
1331  ret.allocSize = lclInds2D_[myRow].size ();
1332  }
1333  else {
1334  ret.allocSize = gblInds2D_[myRow].size ();
1335  }
1336  ret.numEntries = k_numRowEntries_.h_view(myRow);
1337  }
1338  }
1339  else if (nodeNumAllocated_ == 0) {
1340  // have performed allocation, but the graph has no allocation or entries
1341  ret.allocSize = 0;
1342  ret.numEntries = 0;
1343  ret.offset1D = STINV;
1344  }
1345  else if (! indicesAreAllocated ()) {
1346  // haven't performed allocation yet; probably won't hit this code
1347  //
1348  // FIXME (mfh 07 Aug 2014) We want graph's constructors to
1349  // allocate, rather than doing lazy allocation at first insert.
1350  // This will make k_numAllocPerRow_ obsolete.
1351  const bool useNumAllocPerRow = (k_numAllocPerRow_.dimension_0 () != 0);
1352  if (useNumAllocPerRow) {
1353  ret.allocSize = k_numAllocPerRow_.h_view(myRow);
1354  } else {
1355  ret.allocSize = numAllocForAllRows_;
1356  }
1357  ret.numEntries = 0;
1358  ret.offset1D = STINV;
1359  }
1360  else {
1361  // don't know how we ended up here...
1362  TEUCHOS_TEST_FOR_EXCEPT(true);
1363  }
1364  return ret;
1365  }
1366 
1367 
1368  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1369  void
1370  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1371  staticAssertions () const
1372  {
1373  using Teuchos::OrdinalTraits;
1374  typedef LocalOrdinal LO;
1375  typedef GlobalOrdinal GO;
1376  typedef global_size_t GST;
1377 
1378  // Assumption: sizeof(GlobalOrdinal) >= sizeof(LocalOrdinal):
1379  // This is so that we can store local indices in the memory
1380  // formerly occupied by global indices.
1381  static_assert (sizeof (GlobalOrdinal) >= sizeof (LocalOrdinal),
1382  "Tpetra::CrsGraph: sizeof(GlobalOrdinal) must be >= sizeof(LocalOrdinal).");
1383  // Assumption: max(size_t) >= max(LocalOrdinal)
1384  // This is so that we can represent any LocalOrdinal as a size_t.
1385  static_assert (sizeof (size_t) >= sizeof (LocalOrdinal),
1386  "Tpetra::CrsGraph: sizeof(size_t) must be >= sizeof(LocalOrdinal).");
1387  static_assert (sizeof(GST) >= sizeof(size_t),
1388  "Tpetra::CrsGraph: sizeof(Tpetra::global_size_t) must be >= sizeof(size_t).");
1389 
1390  // FIXME (mfh 30 Sep 2015) We're not using
1391  // Teuchos::CompileTimeAssert any more. Can we do these checks
1392  // with static_assert?
1393 
1394  // can't call max() with CompileTimeAssert, because it isn't a
1395  // constant expression; will need to make this a runtime check
1396  const char msg[] = "Tpetra::CrsGraph: Object cannot be created with the "
1397  "given template arguments: size assumptions are not valid.";
1398  TEUCHOS_TEST_FOR_EXCEPTION(
1399  static_cast<size_t> (OrdinalTraits<LO>::max ()) > OrdinalTraits<size_t>::max (),
1400  std::runtime_error, msg);
1401  TEUCHOS_TEST_FOR_EXCEPTION(
1402  static_cast<GST> (OrdinalTraits<LO>::max ()) > static_cast<GST> (OrdinalTraits<GO>::max ()),
1403  std::runtime_error, msg);
1404  TEUCHOS_TEST_FOR_EXCEPTION(
1405  static_cast<size_t> (OrdinalTraits<GO>::max ()) > OrdinalTraits<GST>::max(),
1406  std::runtime_error, msg);
1407  TEUCHOS_TEST_FOR_EXCEPTION(
1408  OrdinalTraits<size_t>::max () > OrdinalTraits<GST>::max (),
1409  std::runtime_error, msg);
1410  }
1411 
1412 
1413  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1414  size_t
1415  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1416  insertIndices (const RowInfo& rowinfo,
1417  const SLocalGlobalViews &newInds,
1418  const ELocalGlobal lg,
1419  const ELocalGlobal I)
1420  {
1421  using Teuchos::ArrayView;
1422 #ifdef HAVE_TPETRA_DEBUG
1423  TEUCHOS_TEST_FOR_EXCEPTION(
1424  lg != GlobalIndices && lg != LocalIndices, std::invalid_argument,
1425  "Tpetra::CrsGraph::insertIndices: lg must be either GlobalIndices or "
1426  "LocalIndices.");
1427 #endif // HAVE_TPETRA_DEBUG
1428  size_t numNewInds = 0;
1429  if (lg == GlobalIndices) { // input indices are global
1430  ArrayView<const GlobalOrdinal> new_ginds = newInds.ginds;
1431  numNewInds = new_ginds.size();
1432  if (I == GlobalIndices) { // store global indices
1433  ArrayView<GlobalOrdinal> gind_view = getGlobalViewNonConst(rowinfo);
1434  std::copy(new_ginds.begin(), new_ginds.end(), gind_view.begin()+rowinfo.numEntries);
1435  }
1436  else if (I == LocalIndices) { // store local indices
1437  ArrayView<LocalOrdinal> lind_view = getLocalViewNonConst(rowinfo);
1438  typename ArrayView<const GlobalOrdinal>::iterator in = new_ginds.begin();
1439  const typename ArrayView<const GlobalOrdinal>::iterator stop = new_ginds.end();
1440  typename ArrayView<LocalOrdinal>::iterator out = lind_view.begin()+rowinfo.numEntries;
1441  while (in != stop) {
1442  *out++ = colMap_->getLocalElement (*in++);
1443  }
1444  }
1445  }
1446  else if (lg == LocalIndices) { // input indices are local
1447  ArrayView<const LocalOrdinal> new_linds = newInds.linds;
1448  numNewInds = new_linds.size();
1449  if (I == LocalIndices) { // store local indices
1450  ArrayView<LocalOrdinal> lind_view = getLocalViewNonConst(rowinfo);
1451  std::copy(new_linds.begin(), new_linds.end(), lind_view.begin()+rowinfo.numEntries);
1452  }
1453  else if (I == GlobalIndices) {
1454  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Tpetra::CrsGraph::"
1455  "insertIndices: the case where the input indices are local and the "
1456  "indices to write are global (lg=LocalIndices, I=GlobalIndices) is "
1457  "not implemented, because it does not make sense." << std::endl <<
1458  "If you have correct local column indices, that means the graph has "
1459  "a column Map. In that case, you should be storing local indices.");
1460  }
1461  }
1462 
1463  // FIXME (mfh 07 Aug 2014) We should just sync at fillComplete,
1464  // but for now, for correctness, do the modify-sync cycle here.
1465  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
1466  k_numRowEntries_.h_view(rowinfo.localRow) += numNewInds;
1467  k_numRowEntries_.template sync<execution_space> ();
1468 
1469  nodeNumEntries_ += numNewInds;
1470  setLocallyModified ();
1471  return numNewInds;
1472  }
1473 
1474 
1475  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1476  void
1477  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1478  insertGlobalIndicesImpl (const LocalOrdinal myRow,
1479  const Teuchos::ArrayView<const GlobalOrdinal>& indices)
1480  {
1481  const char tfecfFuncName[] = "insertGlobalIndicesImpl";
1482 
1483  RowInfo rowInfo = getRowInfo(myRow);
1484  const size_t numNewInds = indices.size();
1485  const size_t newNumEntries = rowInfo.numEntries + numNewInds;
1486  if (newNumEntries > rowInfo.allocSize) {
1487  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1488  getProfileType() == StaticProfile, std::runtime_error,
1489  ": new indices exceed statically allocated graph structure.");
1490 
1491  // update allocation, doubling size to reduce # reallocations
1492  size_t newAllocSize = 2*rowInfo.allocSize;
1493  if (newAllocSize < newNumEntries) {
1494  newAllocSize = newNumEntries;
1495  }
1496  gblInds2D_[myRow].resize(newAllocSize);
1497  nodeNumAllocated_ += (newAllocSize - rowInfo.allocSize);
1498  }
1499 
1500  // Copy new indices at end of global index array
1501  if (k_gblInds1D_.dimension_0 () != 0) {
1502  const size_t numIndsToCopy = static_cast<size_t> (indices.size ());
1503  const size_t offset = rowInfo.offset1D + rowInfo.numEntries;
1504  for (size_t k = 0; k < numIndsToCopy; ++k) {
1505  k_gblInds1D_[offset + k] = indices[k];
1506  }
1507  }
1508  else {
1509  std::copy(indices.begin(), indices.end(),
1510  gblInds2D_[myRow].begin()+rowInfo.numEntries);
1511  }
1512 
1513  // FIXME (mfh 07 Aug 2014) We should just sync at fillComplete,
1514  // but for now, for correctness, do the modify-sync cycle here.
1515  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
1516  k_numRowEntries_.h_view(myRow) += numNewInds;
1517  k_numRowEntries_.template sync<execution_space> ();
1518 
1519  nodeNumEntries_ += numNewInds;
1520  setLocallyModified ();
1521 
1522 #ifdef HAVE_TPETRA_DEBUG
1523  {
1524  const size_t chkNewNumEntries = getNumEntriesInLocalRow (myRow);
1525  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1526  chkNewNumEntries != newNumEntries, std::logic_error,
1527  ": Internal logic error. Please contact Tpetra team.");
1528  }
1529 #endif
1530  }
1531 
1532 
1533  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1534  void
1535  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1536  insertLocalIndicesImpl (const LocalOrdinal myRow,
1537  const Teuchos::ArrayView<const LocalOrdinal>& indices)
1538  {
1539  using Kokkos::MemoryUnmanaged;
1540  using Kokkos::subview;
1541  using Kokkos::View;
1542  typedef LocalOrdinal LO;
1543  const char* tfecfFuncName ("insertLocallIndicesImpl");
1544 
1545  RowInfo rowInfo = getRowInfo(myRow);
1546  const size_t numNewInds = indices.size();
1547  const size_t newNumEntries = rowInfo.numEntries + numNewInds;
1548  if (newNumEntries > rowInfo.allocSize) {
1549  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1550  getProfileType() == StaticProfile, std::runtime_error,
1551  ": new indices exceed statically allocated graph structure.");
1552 
1553  // update allocation, doubling size to reduce number of reallocations
1554  size_t newAllocSize = 2*rowInfo.allocSize;
1555  if (newAllocSize < newNumEntries)
1556  newAllocSize = newNumEntries;
1557  lclInds2D_[myRow].resize(newAllocSize);
1558  nodeNumAllocated_ += (newAllocSize - rowInfo.allocSize);
1559  }
1560 
1561  // Store the new indices at the end of row myRow.
1562  if (k_lclInds1D_.dimension_0 () != 0) {
1563  typedef View<const LO*, execution_space, MemoryUnmanaged> input_view_type;
1564  typedef View<LO*, execution_space> row_view_type;
1565 
1566  input_view_type inputInds (indices.getRawPtr (), indices.size ());
1567  const size_t start = rowInfo.offset1D + rowInfo.numEntries; // end of row
1568  const std::pair<size_t, size_t> rng (start, start + newNumEntries);
1569  row_view_type myInds = subview (k_lclInds1D_, rng);
1570  Kokkos::deep_copy (myInds, inputInds);
1571  }
1572  else {
1573  std::copy (indices.begin (), indices.end (),
1574  lclInds2D_[myRow].begin () + rowInfo.numEntries);
1575  }
1576 
1577  // FIXME (mfh 07 Aug 2014) We should just sync at fillComplete,
1578  // but for now, for correctness, do the modify-sync cycle here.
1579  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
1580  k_numRowEntries_.h_view(myRow) += numNewInds;
1581  k_numRowEntries_.template sync<execution_space> ();
1582 
1583  nodeNumEntries_ += numNewInds;
1584  setLocallyModified ();
1585 #ifdef HAVE_TPETRA_DEBUG
1586  {
1587  const size_t chkNewNumEntries = getNumEntriesInLocalRow (myRow);
1588  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1589  chkNewNumEntries != newNumEntries, std::logic_error,
1590  ": Internal logic error. Please contact Tpetra team.");
1591  }
1592 #endif
1593  }
1594 
1595 
1596  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1597  template <class Scalar>
1598  void
1599  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1600  insertIndicesAndValues (const RowInfo& rowInfo,
1601  const SLocalGlobalViews& newInds,
1602  const Teuchos::ArrayView<Scalar>& oldRowVals,
1603  const Teuchos::ArrayView<const Scalar>& newRowVals,
1604  const ELocalGlobal lg,
1605  const ELocalGlobal I)
1606  {
1607 #ifdef HAVE_TPETRA_DEBUG
1608  size_t numNewInds = 0;
1609  try {
1610  numNewInds = insertIndices (rowInfo, newInds, lg, I);
1611  } catch (std::exception& e) {
1612  TEUCHOS_TEST_FOR_EXCEPTION(
1613  true, std::runtime_error, "Tpetra::CrsGraph::insertIndicesAndValues: "
1614  "insertIndices threw an exception: " << e.what ());
1615  }
1616  TEUCHOS_TEST_FOR_EXCEPTION(
1617  numNewInds > static_cast<size_t> (oldRowVals.size ()), std::runtime_error,
1618  "Tpetra::CrsGraph::insertIndicesAndValues: numNewInds (" << numNewInds
1619  << ") > oldRowVals.size() (" << oldRowVals.size () << ".");
1620 #else
1621  const size_t numNewInds = insertIndices (rowInfo, newInds, lg, I);
1622 #endif // HAVE_TPETRA_DEBUG
1623 
1624  typedef typename Teuchos::ArrayView<Scalar>::size_type size_type;
1625 
1626 #ifdef HAVE_TPETRA_DEBUG
1627  TEUCHOS_TEST_FOR_EXCEPTION(
1628  rowInfo.numEntries + numNewInds > static_cast<size_t> (oldRowVals.size ()),
1629  std::runtime_error, "Tpetra::CrsGraph::insertIndicesAndValues: rowInfo."
1630  "numEntries (" << rowInfo.numEntries << ") + numNewInds (" << numNewInds
1631  << ") > oldRowVals.size() (" << oldRowVals.size () << ").");
1632  TEUCHOS_TEST_FOR_EXCEPTION(
1633  static_cast<size_type> (numNewInds) > newRowVals.size (),
1634  std::runtime_error, "Tpetra::CrsGraph::insertIndicesAndValues: "
1635  "numNewInds (" << numNewInds << ") > newRowVals.size() ("
1636  << newRowVals.size () << ").");
1637 #endif // HAVE_TPETRA_DEBUG
1638 
1639  size_type oldInd = static_cast<size_type> (rowInfo.numEntries);
1640 
1641 #ifdef HAVE_TPETRA_DEBUG
1642  try {
1643 #endif // HAVE_TPETRA_DEBUG
1644  //NOTE: The code in the else branch fails on GCC 4.9 and newer in the assignement oldRowVals[oldInd] = newRowVals[newInd];
1645  //We supply a workaround n as well as other code variants which produce or not produce the error
1646  #if defined(__GNUC__) && defined(__GNUC_MINOR__) && defined(__GNUC_PATCHLEVEL__)
1647  #define GCC_VERSION __GNUC__*100+__GNUC_MINOR__*10+__GNUC_PATCHLEVEL__
1648  #if GCC_VERSION >= 490
1649  #define GCC_WORKAROUND
1650  #endif
1651  #endif
1652  #ifdef GCC_WORKAROUND
1653  size_type nNI = static_cast<size_type>(numNewInds);
1654  if (nNI > 0)
1655  memcpy(&oldRowVals[oldInd], &newRowVals[0], nNI*sizeof(Scalar));
1656  /*
1657  //Original Code Fails
1658  for (size_type newInd = 0; newInd < static_cast<size_type> (numNewInds);
1659  ++newInd, ++oldInd) {
1660  oldRowVals[oldInd] = newRowVals[newInd];
1661  }
1662 
1663 
1664  //char cast variant fails
1665  char* oldRowValPtr = (char*)&oldRowVals[oldInd];
1666  const char* newRowValPtr = (const char*) &newRowVals[0];
1667 
1668  for(size_type newInd = 0; newInd < (nNI * sizeof(Scalar)); newInd++) {
1669  oldRowValPtr[newInd] = newRowValPtr[newInd];
1670  }
1671 
1672  //Raw ptr variant fails
1673  Scalar* oldRowValPtr = &oldRowVals[oldInd];
1674  Scalar* newRowValPtr = const_cast<Scalar*>(&newRowVals[0]);
1675 
1676  for(size_type newInd = 0; newInd < nNI; newInd++) {
1677  oldRowValPtr[newInd] = newRowValPtr[newInd];
1678  }
1679 
1680  //memcpy works
1681  for (size_type newInd = 0; newInd < nNI; newInd++) {
1682  memcpy( &oldRowVals[oldInd+newInd], &newRowVals[newInd], sizeof(Scalar));
1683  }
1684 
1685  //just one loop index fails
1686  for (size_type newInd = 0; newInd < nNI; newInd++) {
1687  oldRowVals[oldInd+newInd] = newRowVals[newInd];
1688  }
1689 
1690  //inline increment fails
1691  for (size_type newInd = 0; newInd < numNewInds;) {
1692  oldRowVals[oldInd++] = newRowVals[newInd++];
1693  }
1694 
1695  */
1696 
1697  #else // GCC Workaround above
1698  for (size_type newInd = 0; newInd < static_cast<size_type> (numNewInds);
1699  ++newInd, ++oldInd) {
1700  oldRowVals[oldInd] = newRowVals[newInd];
1701  }
1702  #endif // GCC Workaround
1703 #ifdef HAVE_TPETRA_DEBUG
1704  }
1705  catch (std::exception& e) {
1706  TEUCHOS_TEST_FOR_EXCEPTION(
1707  true, std::runtime_error, "Tpetra::CrsGraph::insertIndicesAndValues: "
1708  "for loop for copying values threw an exception: " << e.what ());
1709  }
1710 #endif // HAVE_TPETRA_DEBUG
1711  }
1712 
1713 
1714  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1715  void
1716  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1717  sortRowIndices (const RowInfo rowinfo)
1718  {
1719  if (rowinfo.numEntries > 0) {
1720  Teuchos::ArrayView<LocalOrdinal> inds_view =
1721  this->getLocalViewNonConst (rowinfo);
1722  std::sort (inds_view.begin (), inds_view.begin () + rowinfo.numEntries);
1723  }
1724  }
1725 
1726 
1727  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1728  template <class Scalar>
1729  void
1730  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1731  sortRowIndicesAndValues (const RowInfo rowinfo,
1732  const Teuchos::ArrayView<Scalar>& values)
1733  {
1734  if (rowinfo.numEntries > 0) {
1735  Teuchos::ArrayView<LocalOrdinal> inds_view =
1736  this->getLocalViewNonConst (rowinfo);
1737  sort2 (inds_view.begin (), inds_view.begin () + rowinfo.numEntries,
1738  values.begin ());
1739  }
1740  }
1741 
1742 
1743  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1744  void
1745  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1746  mergeRowIndices (RowInfo rowinfo)
1747  {
1748  using Teuchos::ArrayView;
1749  const char tfecfFuncName[] = "mergeRowIndices: ";
1750  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1751  isStorageOptimized (), std::logic_error, "The graph is already storage "
1752  "optimized, so we shouldn't be merging any indices. "
1753  "Please report this bug to the Tpetra developers.");
1754 
1755  ArrayView<LocalOrdinal> inds_view = this->getLocalViewNonConst (rowinfo);
1756  typename ArrayView<LocalOrdinal>::iterator beg, end, newend;
1757  beg = inds_view.begin();
1758  end = inds_view.begin() + rowinfo.numEntries;
1759  newend = std::unique(beg,end);
1760  const size_t mergedEntries = newend - beg;
1761 #ifdef HAVE_TPETRA_DEBUG
1762  // merge should not have eliminated any entries; if so, the
1763  // assignment below will destroy the packed structure
1764  TEUCHOS_TEST_FOR_EXCEPT( isStorageOptimized () && mergedEntries != rowinfo.numEntries );
1765 #endif // HAVE_TPETRA_DEBUG
1766 
1767  // FIXME (mfh 07 Aug 2014) We should just sync at fillComplete,
1768  // but for now, for correctness, do the modify-sync cycle here.
1769  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
1770  k_numRowEntries_.h_view(rowinfo.localRow) = mergedEntries;
1771  k_numRowEntries_.template sync<execution_space> ();
1772 
1773  nodeNumEntries_ -= (rowinfo.numEntries - mergedEntries);
1774  }
1775 
1776 
1777  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1778  template<class Scalar>
1779  void
1780  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1781  mergeRowIndicesAndValues (RowInfo rowinfo,
1782  const Teuchos::ArrayView<Scalar>& rowValues)
1783  {
1784  using Teuchos::ArrayView;
1785  const char tfecfFuncName[] = "mergeRowIndicesAndValues: ";
1786  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1787  isStorageOptimized(), std::logic_error, "It is invalid to call this "
1788  "method if the graph's storage has already been optimized. Please "
1789  "report this bug to the Tpetra developers.");
1790 
1791  typedef typename ArrayView<Scalar>::iterator Iter;
1792  Iter rowValueIter = rowValues.begin ();
1793  ArrayView<LocalOrdinal> inds_view = getLocalViewNonConst (rowinfo);
1794  typename ArrayView<LocalOrdinal>::iterator beg, end, newend;
1795 
1796  // beg,end define a half-exclusive interval over which to iterate.
1797  beg = inds_view.begin();
1798  end = inds_view.begin() + rowinfo.numEntries;
1799  newend = beg;
1800  if (beg != end) {
1801  typename ArrayView<LocalOrdinal>::iterator cur = beg + 1;
1802  Iter vcur = rowValueIter + 1;
1803  Iter vend = rowValueIter;
1804  cur = beg+1;
1805  while (cur != end) {
1806  if (*cur != *newend) {
1807  // new entry; save it
1808  ++newend;
1809  ++vend;
1810  (*newend) = (*cur);
1811  (*vend) = (*vcur);
1812  }
1813  else {
1814  // old entry; merge it
1815  //(*vend) = f (*vend, *vcur);
1816  (*vend) += *vcur;
1817  }
1818  ++cur;
1819  ++vcur;
1820  }
1821  ++newend; // one past the last entry, per typical [beg,end) semantics
1822  }
1823  const size_t mergedEntries = newend - beg;
1824 #ifdef HAVE_TPETRA_DEBUG
1825  // merge should not have eliminated any entries; if so, the
1826  // assignment below will destroy the packed structure
1827  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
1828  isStorageOptimized() && mergedEntries != rowinfo.numEntries,
1829  std::logic_error,
1830  ": Merge was incorrect; it eliminated entries from the graph. "
1831  << std::endl << "Please report this bug to the Tpetra developers.");
1832 #endif // HAVE_TPETRA_DEBUG
1833 
1834  // FIXME (mfh 07 Aug 2014) We should just sync at fillComplete,
1835  // but for now, for correctness, do the modify-sync cycle here.
1836  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
1837  k_numRowEntries_.h_view(rowinfo.localRow) = mergedEntries;
1838  k_numRowEntries_.template sync<execution_space> ();
1839 
1840  nodeNumEntries_ -= (rowinfo.numEntries - mergedEntries);
1841  }
1842 
1843 
1844  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1845  void
1846  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1847  setDomainRangeMaps (const Teuchos::RCP<const map_type>& domainMap,
1848  const Teuchos::RCP<const map_type>& rangeMap)
1849  {
1850  // simple pointer comparison for equality
1851  if (domainMap_ != domainMap) {
1852  domainMap_ = domainMap;
1853  importer_ = null;
1854  }
1855  if (rangeMap_ != rangeMap) {
1856  rangeMap_ = rangeMap;
1857  exporter_ = null;
1858  }
1859  }
1860 
1861 
1862  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1863  size_t
1864  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1865  findLocalIndex (RowInfo rowinfo, LocalOrdinal ind, size_t hint) const
1866  {
1867  using Teuchos::ArrayView;
1868  ArrayView<const LocalOrdinal> colInds = this->getLocalView (rowinfo);
1869  return this->findLocalIndex (rowinfo, ind, colInds, hint);
1870  }
1871 
1872 
1873  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1874  size_t
1875  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1876  findLocalIndex (RowInfo rowinfo,
1877  LocalOrdinal ind,
1878  Teuchos::ArrayView<const LocalOrdinal> colInds,
1879  size_t hint) const
1880  {
1881  typedef typename Teuchos::ArrayView<const LocalOrdinal>::iterator IT;
1882 
1883  // If the hint was correct, then the hint is the offset to return.
1884  if (hint < rowinfo.numEntries && colInds[hint] == ind) {
1885  return hint;
1886  }
1887 
1888  // The hint was wrong, so we must search for the given column
1889  // index in the column indices for the given row. How we do the
1890  // search depends on whether the graph's column indices are
1891  // sorted.
1892  IT beg = colInds.begin ();
1893  IT end = beg + rowinfo.numEntries;
1894  IT ptr = beg + rowinfo.numEntries; // "null"
1895  bool found = true;
1896 
1897  if (isSorted ()) {
1898  std::pair<IT,IT> p = std::equal_range (beg, end, ind); // binary search
1899  if (p.first == p.second) {
1900  found = false;
1901  } else {
1902  ptr = p.first;
1903  }
1904  }
1905  else {
1906  ptr = std::find (beg, end, ind); // direct search
1907  if (ptr == end) {
1908  found = false;
1909  }
1910  }
1911 
1912  if (found) {
1913  return static_cast<size_t> (ptr - beg);
1914  }
1915  else {
1916  return Teuchos::OrdinalTraits<size_t>::invalid ();
1917  }
1918  }
1919 
1920 
1921  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1922  size_t
1923  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1924  findGlobalIndex (RowInfo rowinfo, GlobalOrdinal ind, size_t hint) const
1925  {
1926  using Teuchos::ArrayView;
1927  typedef typename ArrayView<const GlobalOrdinal>::iterator IT;
1928 
1929  // Don't let an invalid global column index through.
1930  if (ind == Teuchos::OrdinalTraits<GlobalOrdinal>::invalid ()) {
1931  return Teuchos::OrdinalTraits<size_t>::invalid ();
1932  }
1933 
1934  ArrayView<const GlobalOrdinal> indices = getGlobalView (rowinfo);
1935 
1936  // We don't actually require that the hint be a valid index.
1937  // If it is not in range, we just ignore it.
1938  if (hint < rowinfo.numEntries && indices[hint] == ind) {
1939  return hint;
1940  }
1941 
1942  IT beg = indices.begin ();
1943  IT end = indices.begin () + rowinfo.numEntries; // not indices.end()
1944  if (isSorted ()) { // use binary search
1945  const std::pair<IT,IT> p = std::equal_range (beg, end, ind);
1946  if (p.first == p.second) { // range of matching entries is empty
1947  return Teuchos::OrdinalTraits<size_t>::invalid ();
1948  } else {
1949  return p.first - beg;
1950  }
1951  }
1952  else { // not sorted; must use linear search
1953  const IT loc = std::find (beg, end, ind);
1954  if (loc == end) {
1955  return Teuchos::OrdinalTraits<size_t>::invalid ();
1956  } else {
1957  return loc - beg;
1958  }
1959  }
1960  }
1961 
1962 
1963  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1964  void
1965  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1966  clearGlobalConstants ()
1967  {
1968  globalNumEntries_ = Teuchos::OrdinalTraits<global_size_t>::invalid ();
1969  globalNumDiags_ = Teuchos::OrdinalTraits<global_size_t>::invalid ();
1970  globalMaxNumRowEntries_ = Teuchos::OrdinalTraits<global_size_t>::invalid ();
1971  haveGlobalConstants_ = false;
1972  }
1973 
1974 
1975  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
1976  void
1977  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
1978  checkInternalState () const
1979  {
1980 #ifdef HAVE_TPETRA_DEBUG
1981  const global_size_t GSTI = Teuchos::OrdinalTraits<global_size_t>::invalid ();
1982  const size_t STI = Teuchos::OrdinalTraits<size_t>::invalid ();
1983  const char err[] = "Tpetra::CrsGraph::checkInternalState: Likely internal "
1984  "logic error. Please contact Tpetra team.";
1985  // check the internal state of this data structure
1986  // this is called by numerous state-changing methods, in a debug build, to ensure that the object
1987  // always remains in a valid state
1988  // the graph should have been allocated with a row map
1989  TEUCHOS_TEST_FOR_EXCEPTION( rowMap_ == null, std::logic_error, err );
1990  // am either complete or active
1991  TEUCHOS_TEST_FOR_EXCEPTION( isFillActive() == isFillComplete(), std::logic_error, err );
1992  // if the graph has been fill completed, then all maps should be present
1993  TEUCHOS_TEST_FOR_EXCEPTION( isFillComplete() == true && (colMap_ == null || rangeMap_ == null || domainMap_ == null), std::logic_error, err );
1994  // if storage has been optimized, then indices should have been allocated (even if trivially so)
1995  TEUCHOS_TEST_FOR_EXCEPTION( isStorageOptimized() == true && indicesAreAllocated() == false, std::logic_error, err );
1996  // if storage has been optimized, then number of allocated is now the number of entries
1997  TEUCHOS_TEST_FOR_EXCEPTION( isStorageOptimized() == true && getNodeAllocationSize() != getNodeNumEntries(), std::logic_error, err );
1998  // if graph doesn't have the global constants, then they should all be marked as invalid
1999  TEUCHOS_TEST_FOR_EXCEPTION( haveGlobalConstants_ == false && ( globalNumEntries_ != GSTI || globalNumDiags_ != GSTI || globalMaxNumRowEntries_ != GSTI ), std::logic_error, err );
2000  // if the graph has global cosntants, then they should be valid.
2001  TEUCHOS_TEST_FOR_EXCEPTION( haveGlobalConstants_ == true && ( globalNumEntries_ == GSTI || globalNumDiags_ == GSTI || globalMaxNumRowEntries_ == GSTI ), std::logic_error, err );
2002  TEUCHOS_TEST_FOR_EXCEPTION( haveGlobalConstants_ == true && ( globalNumEntries_ < nodeNumEntries_ || globalNumDiags_ < nodeNumDiags_ || globalMaxNumRowEntries_ < nodeMaxNumRowEntries_ ),
2003  std::logic_error, err );
2004  // if indices are allocated, then the allocation specifications should have been released
2005  TEUCHOS_TEST_FOR_EXCEPTION(
2006  indicesAreAllocated () && (numAllocForAllRows_ != 0 ||
2007  k_numAllocPerRow_.dimension_0 () != 0),
2008  std::logic_error, err );
2009  // if indices are not allocated, then information dictating allocation quantities should be present
2010  TEUCHOS_TEST_FOR_EXCEPTION(
2011  ! indicesAreAllocated () && (nodeNumAllocated_ != STI ||
2012  nodeNumEntries_ != 0),
2013  std::logic_error, err );
2014  // if storage is optimized, then profile should be static
2015  TEUCHOS_TEST_FOR_EXCEPTION( isStorageOptimized() && pftype_ != StaticProfile, std::logic_error, err );
2016 
2017  // If k_rowPtrs_ exists (has nonzero size), it must have N+1 rows,
2018  // and k_rowPtrs_(N) must equal k_gblInds1D_.dimension_0() (if
2019  // globally indexed) or k_lclInds1D_.dimension_0() (if locally
2020  // indexed).
2021 
2022  TEUCHOS_TEST_FOR_EXCEPTION(
2023  isGloballyIndexed () && k_rowPtrs_.dimension_0 () != 0 &&
2024  (static_cast<size_t> (k_rowPtrs_.dimension_0 ()) != getNodeNumRows () + 1 ||
2025  k_rowPtrs_(getNodeNumRows ()) != static_cast<size_t> (k_gblInds1D_.dimension_0 ())),
2026  std::logic_error, err );
2027 
2028  TEUCHOS_TEST_FOR_EXCEPTION(
2029  isLocallyIndexed () && k_rowPtrs_.dimension_0 () != 0 &&
2030  (static_cast<size_t> (k_rowPtrs_.dimension_0 ()) != getNodeNumRows () + 1 ||
2031  k_rowPtrs_(getNodeNumRows ()) != static_cast<size_t> (k_lclInds1D_.dimension_0 ())),
2032  std::logic_error, err );
2033 
2034  // If profile is dynamic and indices are allocated, then 2-D
2035  // allocations of column index storage (either local or global)
2036  // must be present.
2037  TEUCHOS_TEST_FOR_EXCEPTION(
2038  pftype_ == DynamicProfile &&
2039  indicesAreAllocated () &&
2040  getNodeNumRows () > 0 &&
2041  lclInds2D_.is_null () && gblInds2D_.is_null (),
2042  std::logic_error, err );
2043 
2044  // If profile is dynamic and the calling process owns nonzero
2045  // rows, then k_numRowEntries_ and 2-D storage of column indices
2046  // (whether local or global) must be present.
2047  TEUCHOS_TEST_FOR_EXCEPTION(
2048  pftype_ == DynamicProfile &&
2049  indicesAreAllocated () &&
2050  getNodeNumRows () > 0 &&
2051  (k_numRowEntries_.dimension_0 () == 0 || (lclInds2D_.is_null () && gblInds2D_.is_null ())),
2052  std::logic_error, err );
2053 
2054  // if profile is dynamic, then 1D allocations should not be present
2055  TEUCHOS_TEST_FOR_EXCEPTION(
2056  pftype_ == DynamicProfile &&
2057  (k_lclInds1D_.dimension_0 () != 0 || k_gblInds1D_.dimension_0 () != 0),
2058  std::logic_error, err );
2059  // if profile is dynamic, then row offsets should not be present
2060  TEUCHOS_TEST_FOR_EXCEPTION(
2061  pftype_ == DynamicProfile && k_rowPtrs_.dimension_0 () != 0,
2062  std::logic_error, err );
2063  // if profile is static and we have allocated non-trivially, then
2064  // 1D allocations should be present
2065  TEUCHOS_TEST_FOR_EXCEPTION(
2066  pftype_ == StaticProfile && indicesAreAllocated () &&
2067  getNodeAllocationSize () > 0 && k_lclInds1D_.dimension_0 () == 0 &&
2068  k_gblInds1D_.dimension_0 () == 0,
2069  std::logic_error, err);
2070  // if profile is static, then 2D allocations should not be present
2071  TEUCHOS_TEST_FOR_EXCEPTION(
2072  pftype_ == StaticProfile && (lclInds2D_ != null || gblInds2D_ != null),
2073  std::logic_error, err );
2074 
2075  // if indices are not allocated, then none of the buffers should be.
2076  TEUCHOS_TEST_FOR_EXCEPTION(
2077  ! indicesAreAllocated () &&
2078  ((k_rowPtrs_.dimension_0 () != 0 || k_numRowEntries_.dimension_0 () != 0) ||
2079  k_lclInds1D_.dimension_0 () != 0 || lclInds2D_ != null ||
2080  k_gblInds1D_.dimension_0 () != 0 || gblInds2D_ != null),
2081  std::logic_error, err );
2082 
2083  // indices may be local or global only if they are allocated
2084  // (numAllocated is redundant; could simply be indicesAreLocal_ ||
2085  // indicesAreGlobal_)
2086  TEUCHOS_TEST_FOR_EXCEPTION( (indicesAreLocal_ || indicesAreGlobal_) && ! indicesAreAllocated_, std::logic_error, err );
2087  // indices may be local or global, but not both
2088  TEUCHOS_TEST_FOR_EXCEPTION( indicesAreLocal_ == true && indicesAreGlobal_ == true, std::logic_error, err );
2089  // if indices are local, then global allocations should not be present
2090  TEUCHOS_TEST_FOR_EXCEPTION(
2091  indicesAreLocal_ && (k_gblInds1D_.dimension_0 () != 0 || gblInds2D_ != null),
2092  std::logic_error, err );
2093  // if indices are global, then local allocations should not be present
2094  TEUCHOS_TEST_FOR_EXCEPTION(
2095  indicesAreGlobal_ && (k_lclInds1D_.dimension_0 () != 0 || lclInds2D_ != null),
2096  std::logic_error, err );
2097  // if indices are local, then local allocations should be present
2098  TEUCHOS_TEST_FOR_EXCEPTION(
2099  indicesAreLocal_ && getNodeAllocationSize () > 0 &&
2100  k_lclInds1D_.dimension_0 () == 0 && getNodeNumRows () > 0 &&
2101  lclInds2D_.is_null (),
2102  std::logic_error, err);
2103  // if indices are global, then global allocations should be present
2104  TEUCHOS_TEST_FOR_EXCEPTION(
2105  indicesAreGlobal_ && getNodeAllocationSize () > 0 &&
2106  k_gblInds1D_.dimension_0 () == 0 && getNodeNumRows () > 0 &&
2107  gblInds2D_.is_null (),
2108  std::logic_error, err);
2109  // if indices are allocated, then we should have recorded how many were allocated
2110  TEUCHOS_TEST_FOR_EXCEPTION( indicesAreAllocated() == true && nodeNumAllocated_ == STI, std::logic_error, err );
2111  // if indices are not allocated, then the allocation size should be marked invalid
2112  TEUCHOS_TEST_FOR_EXCEPTION( indicesAreAllocated() == false && nodeNumAllocated_ != STI, std::logic_error, err );
2113  // check the actual allocations
2114  if (indicesAreAllocated()) {
2115  size_t actualNumAllocated = 0;
2116  if (pftype_ == DynamicProfile) {
2117  if (isGloballyIndexed() && gblInds2D_ != null) {
2118  for (size_t r = 0; r < getNodeNumRows(); ++r) {
2119  actualNumAllocated += gblInds2D_[r].size();
2120  }
2121  }
2122  else if (isLocallyIndexed() && lclInds2D_ != null) {
2123  for (size_t r = 0; r < getNodeNumRows(); ++r) {
2124  actualNumAllocated += lclInds2D_[r].size();
2125  }
2126  }
2127  TEUCHOS_TEST_FOR_EXCEPTION(actualNumAllocated != nodeNumAllocated_, std::logic_error, err );
2128  }
2129  else if (k_rowPtrs_.dimension_0 () != 0) { // pftype_ == StaticProfile
2130  TEUCHOS_TEST_FOR_EXCEPTION(
2131  static_cast<size_t> (k_rowPtrs_.dimension_0 ()) != getNodeNumRows () + 1,
2132  std::logic_error, err);
2133 
2134  actualNumAllocated = k_rowPtrs_(getNodeNumRows ());
2135  TEUCHOS_TEST_FOR_EXCEPTION(
2136  isLocallyIndexed () &&
2137  static_cast<size_t> (k_lclInds1D_.dimension_0 ()) != actualNumAllocated,
2138  std::logic_error, err );
2139  TEUCHOS_TEST_FOR_EXCEPTION(
2140  isGloballyIndexed () &&
2141  static_cast<size_t> (k_gblInds1D_.dimension_0 ()) != actualNumAllocated,
2142  std::logic_error, err );
2143  TEUCHOS_TEST_FOR_EXCEPTION(actualNumAllocated != nodeNumAllocated_, std::logic_error, err );
2144  }
2145  }
2146 #endif // HAVE_TPETRA_DEBUG
2147  }
2148 
2149 
2150  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2151  size_t
2152  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2153  getNumEntriesInGlobalRow (GlobalOrdinal globalRow) const
2154  {
2155  using Teuchos::OrdinalTraits;
2156  const LocalOrdinal lrow = rowMap_->getLocalElement (globalRow);
2157  if (hasRowInfo () && lrow != OrdinalTraits<LocalOrdinal>::invalid ()) {
2158  const RowInfo rowinfo = getRowInfo (lrow);
2159  return rowinfo.numEntries;
2160  } else {
2161  return OrdinalTraits<size_t>::invalid ();
2162  }
2163  }
2164 
2165 
2166  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2167  size_t
2168  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2169  getNumEntriesInLocalRow (LocalOrdinal localRow) const
2170  {
2171  if (hasRowInfo () && rowMap_->isNodeLocalElement (localRow)) {
2172  const RowInfo rowinfo = getRowInfo (localRow);
2173  return rowinfo.numEntries;
2174  } else {
2175  return Teuchos::OrdinalTraits<size_t>::invalid ();
2176  }
2177  }
2178 
2179 
2180  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2181  size_t
2182  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2183  getNumAllocatedEntriesInGlobalRow (GlobalOrdinal globalRow) const
2184  {
2185  const LocalOrdinal lrow = rowMap_->getLocalElement (globalRow);
2186  if (hasRowInfo () && lrow != Teuchos::OrdinalTraits<LocalOrdinal>::invalid ()) {
2187  const RowInfo rowinfo = getRowInfo (lrow);
2188  return rowinfo.allocSize;
2189  } else {
2190  return Teuchos::OrdinalTraits<size_t>::invalid ();
2191  }
2192  }
2193 
2194 
2195  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2196  size_t
2197  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2198  getNumAllocatedEntriesInLocalRow (LocalOrdinal localRow) const
2199  {
2200  if (hasRowInfo () && rowMap_->isNodeLocalElement (localRow)) {
2201  const RowInfo rowinfo = getRowInfo (localRow);
2202  return rowinfo.allocSize;
2203  } else {
2204  return Teuchos::OrdinalTraits<size_t>::invalid();
2205  }
2206  }
2207 
2208 
2209  namespace { // (anonymous)
2210  template<class OutputViewType, class InputViewType>
2211  class CopyOffsets {
2212  public:
2213  typedef typename OutputViewType::execution_space execution_space;
2214 
2215  CopyOffsets (const OutputViewType& ptr_out,
2216  const InputViewType& ptr_in) :
2217  ptr_out_ (ptr_out),
2218  ptr_in_ (ptr_in)
2219  {}
2220 
2221  KOKKOS_INLINE_FUNCTION void operator () (const ptrdiff_t& i) const {
2222  typedef typename OutputViewType::non_const_value_type value_type;
2223  // FIXME (mfh 22 Mar 2015) Change this into parallel_reduce
2224  // and check for overflow.
2225  ptr_out_(i) = static_cast<value_type> (ptr_in_(i));
2226  }
2227 
2228  private:
2229  OutputViewType ptr_out_;
2230  InputViewType ptr_in_;
2231  };
2232  } // namespace (anonymous)
2233 
2234 
2235  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2236  Teuchos::ArrayRCP<const size_t>
2237  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2238  getNodeRowPtrs () const
2239  {
2240  using Kokkos::ViewAllocateWithoutInitializing;
2241  using Kokkos::create_mirror_view;
2242  using Teuchos::ArrayRCP;
2243  typedef typename local_graph_type::row_map_type row_map_type;
2244  typedef typename row_map_type::non_const_value_type row_offset_type;
2245 #ifdef HAVE_TPETRA_DEBUG
2246  const char prefix[] = "Tpetra::CrsGraph::getNodeRowPtrs: ";
2247  const char suffix[] = " Please report this bug to the Tpetra developers.";
2248 #endif // HAVE_TPETRA_DEBUG
2249  const size_t size = k_rowPtrs_.dimension_0 ();
2250  const bool same = Kokkos::Impl::is_same<size_t, row_offset_type>::value;
2251 
2252  if (size == 0) {
2253  return ArrayRCP<const size_t> ();
2254  }
2255 
2256  ArrayRCP<const row_offset_type> ptr_rot;
2257  ArrayRCP<const size_t> ptr_st;
2258  if (same) { // size_t == row_offset_type
2259  // NOTE (mfh 22 Mar 2015) In a debug build of Kokkos, the result
2260  // of create_mirror_view might actually be a new allocation.
2261  // This helps with debugging when there are two memory spaces.
2262  typename row_map_type::HostMirror ptr_h = create_mirror_view (k_rowPtrs_);
2263  Kokkos::deep_copy (ptr_h, k_rowPtrs_);
2264 #ifdef HAVE_TPETRA_DEBUG
2265  TEUCHOS_TEST_FOR_EXCEPTION(
2266  ptr_h.dimension_0 () != k_rowPtrs_.dimension_0 (), std::logic_error,
2267  prefix << "size_t == row_offset_type, but ptr_h.dimension_0() = "
2268  << ptr_h.dimension_0 () << " != k_rowPtrs_.dimension_0() = "
2269  << k_rowPtrs_.dimension_0 () << ".");
2270  TEUCHOS_TEST_FOR_EXCEPTION(
2271  same && size != 0 && k_rowPtrs_.ptr_on_device () == NULL, std::logic_error,
2272  prefix << "size_t == row_offset_type and k_rowPtrs_.dimension_0() = "
2273  << size << " != 0, but k_rowPtrs_.ptr_on_device() == NULL." << suffix);
2274  TEUCHOS_TEST_FOR_EXCEPTION(
2275  same && size != 0 && ptr_h.ptr_on_device () == NULL, std::logic_error,
2276  prefix << "size_t == row_offset_type and k_rowPtrs_.dimension_0() = "
2277  << size << " != 0, but create_mirror_view(k_rowPtrs_).ptr_on_device() "
2278  "== NULL." << suffix);
2279 #endif // HAVE_TPETRA_DEBUG
2280  ptr_rot = Kokkos::Compat::persistingView (ptr_h);
2281  }
2282  else { // size_t != row_offset_type
2283  typedef Kokkos::View<size_t*, device_type> ret_view_type;
2284  ret_view_type ptr_d (ViewAllocateWithoutInitializing ("ptr"), size);
2285  CopyOffsets<ret_view_type, row_map_type> functor (ptr_d, k_rowPtrs_);
2286  Kokkos::parallel_for (size, functor);
2287  typename ret_view_type::HostMirror ptr_h = create_mirror_view (ptr_d);
2288  Kokkos::deep_copy (ptr_h, ptr_d);
2289  ptr_st = Kokkos::Compat::persistingView (ptr_h);
2290  }
2291 #ifdef HAVE_TPETRA_DEBUG
2292  TEUCHOS_TEST_FOR_EXCEPTION(
2293  same && size != 0 && ptr_rot.is_null (), std::logic_error,
2294  prefix << "size_t == row_offset_type and size = " << size
2295  << " != 0, but ptr_rot is null." << suffix);
2296  TEUCHOS_TEST_FOR_EXCEPTION(
2297  ! same && size != 0 && ptr_st.is_null (), std::logic_error,
2298  prefix << "size_t != row_offset_type and size = " << size
2299  << " != 0, but ptr_st is null." << suffix);
2300 #endif // HAVE_TPETRA_DEBUG
2301 
2302  // If size_t == row_offset_type, return a persisting host view of
2303  // k_rowPtrs_. Otherwise, return a size_t host copy of k_rowPtrs_.
2304 #ifdef HAVE_TPETRA_DEBUG
2305  ArrayRCP<const size_t> retval =
2306  Kokkos::Impl::if_c<same,
2307  ArrayRCP<const row_offset_type>,
2308  ArrayRCP<const size_t> >::select (ptr_rot, ptr_st);
2309  TEUCHOS_TEST_FOR_EXCEPTION(
2310  size != 0 && retval.is_null (), std::logic_error,
2311  prefix << "size = " << size << " != 0, but retval is null." << suffix);
2312  return retval;
2313 #else
2314  return Kokkos::Impl::if_c<same,
2315  ArrayRCP<const row_offset_type>,
2316  ArrayRCP<const size_t> >::select (ptr_rot, ptr_st);
2317 #endif // HAVE_TPETRA_DEBUG
2318  }
2319 
2320 
2321  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2322  Teuchos::ArrayRCP<const LocalOrdinal>
2323  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2324  getNodePackedIndices () const
2325  {
2326  return Kokkos::Compat::persistingView (k_lclInds1D_);
2327  }
2328 
2329 
2330  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2331  void
2332  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2333  getLocalRowCopy (LocalOrdinal localRow,
2334  const Teuchos::ArrayView<LocalOrdinal>&indices,
2335  size_t& numEntries) const
2336  {
2337  using Teuchos::ArrayView;
2338  typedef LocalOrdinal LO;
2339  typedef GlobalOrdinal GO;
2340 
2341  TEUCHOS_TEST_FOR_EXCEPTION(
2342  isGloballyIndexed () && ! hasColMap (), std::runtime_error,
2343  "Tpetra::CrsGraph::getLocalRowCopy: The graph is globally indexed and "
2344  "does not have a column Map yet. That means we don't have local indices "
2345  "for columns yet, so it doesn't make sense to call this method. If the "
2346  "graph doesn't have a column Map yet, you should call fillComplete on "
2347  "it first.");
2348  TEUCHOS_TEST_FOR_EXCEPTION(
2349  ! hasRowInfo(), std::runtime_error,
2350  "Tpetra::CrsGraph::getLocalRowCopy: graph row information was deleted "
2351  "at fillComplete().");
2352 
2353  if (! getRowMap ()->isNodeLocalElement (localRow)) {
2354  numEntries = 0;
2355  return;
2356  }
2357 
2358  const RowInfo rowinfo = getRowInfo (localRow);
2359  const size_t theNumEntries = rowinfo.numEntries;
2360 
2361  TEUCHOS_TEST_FOR_EXCEPTION(
2362  static_cast<size_t> (indices.size ()) < theNumEntries,
2363  std::runtime_error,
2364  "Tpetra::CrsGraph::getLocalRowCopy: The given row " << localRow << " has "
2365  << theNumEntries << " entries, but indices.size() = " << indices.size ()
2366  << ", which does not suffice to store the row's indices.");
2367 
2368  numEntries = theNumEntries;
2369 
2370  if (isLocallyIndexed ()) {
2371  ArrayView<const LO> lview = getLocalView (rowinfo);
2372  std::copy (lview.begin (), lview.begin () + numEntries, indices.begin ());
2373  }
2374  else if (isGloballyIndexed ()) {
2375  ArrayView<const GO> gview = getGlobalView (rowinfo);
2376  const map_type& colMap = * (this->getColMap ());
2377  for (size_t j = 0; j < numEntries; ++j) {
2378  indices[j] = colMap.getLocalElement (gview[j]);
2379  }
2380  }
2381  else {
2382  // If the graph on the calling process is neither locally nor
2383  // globally indexed, that means it owns no column indices.
2384  //
2385  // FIXME (mfh 21 Oct 2013) It's not entirely clear to me whether
2386  // we can reach this branch, given the checks above. However,
2387  // if that is the case, it should still be correct to call this
2388  // function if the calling process owns no column indices.
2389  numEntries = 0;
2390  }
2391  }
2392 
2393 
2394  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2395  void
2396  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2397  getGlobalRowCopy (GlobalOrdinal globalRow,
2398  const Teuchos::ArrayView<GlobalOrdinal>& indices,
2399  size_t& NumIndices) const
2400  {
2401  using Teuchos::ArrayView;
2402  const char tfecfFuncName[] = "getGlobalRowCopy: ";
2403  // we either currently store global indices, or we have a column
2404  // map with which to transcribe our local indices for the user
2405  const LocalOrdinal lrow = rowMap_->getLocalElement (globalRow);
2406 
2407  // FIXME (mfh 22 Aug 2014) Instead of throwing an exception,
2408  // should just set NumIndices=0 and return. In that case, the
2409  // calling process owns no entries in that row, so the right thing
2410  // to do is to return an empty copy.
2411  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2412  lrow == Teuchos::OrdinalTraits<LocalOrdinal>::invalid (),
2413  std::runtime_error,
2414  "GlobalRow (== " << globalRow << ") does not belong to this process.");
2415  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2416  ! hasRowInfo (), std::runtime_error,
2417  "Graph row information was deleted at fillComplete().");
2418  const RowInfo rowinfo = this->getRowInfo (static_cast<size_t> (lrow));
2419  NumIndices = rowinfo.numEntries;
2420  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2421  static_cast<size_t> (indices.size ()) < NumIndices, std::runtime_error,
2422  "Specified storage (size==" << indices.size () << ") does not suffice "
2423  "to hold all entries for this row (NumIndices == " << NumIndices << ").");
2424  if (isLocallyIndexed ()) {
2425  ArrayView<const LocalOrdinal> lview = this->getLocalView (rowinfo);
2426  for (size_t j = 0; j < NumIndices; ++j) {
2427  indices[j] = colMap_->getGlobalElement (lview[j]);
2428  }
2429  }
2430  else if (isGloballyIndexed ()) {
2431  ArrayView<const GlobalOrdinal> gview = this->getGlobalView (rowinfo);
2432  std::copy (gview.begin (), gview.begin () + NumIndices, indices.begin ());
2433  }
2434  }
2435 
2436 
2437  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2438  void
2439  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2440  getLocalRowView (LocalOrdinal localRow,
2441  Teuchos::ArrayView<const LocalOrdinal>& indices) const
2442  {
2443  const char tfecfFuncName[] = "getLocalRowView: ";
2444  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2445  isGloballyIndexed (), std::runtime_error, "The graph's indices are "
2446  "currently stored as global indices, so we cannot return a view with "
2447  "local column indices, whether or not the graph has a column Map. If "
2448  "the graph _does_ have a column Map, use getLocalRowCopy() instead.");
2449  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2450  ! hasRowInfo (), std::runtime_error, "Graph row information was "
2451  "deleted at fillComplete().");
2452  indices = Teuchos::null;
2453  if (rowMap_->isNodeLocalElement (localRow)) {
2454  const RowInfo rowinfo = this->getRowInfo (localRow);
2455  if (rowinfo.numEntries > 0) {
2456  indices = this->getLocalView (rowinfo);
2457  // getLocalView returns a view of the _entire_ row, including
2458  // any extra space at the end (which 1-D unpacked storage
2459  // might have, for example). That's why we have to take a
2460  // subview of the returned view.
2461  indices = indices (0, rowinfo.numEntries);
2462  }
2463  }
2464 #ifdef HAVE_TPETRA_DEBUG
2465  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2466  static_cast<size_t> (indices.size ()) != this->getNumEntriesInLocalRow (localRow),
2467  std::logic_error, ": Violated stated post-conditions. Please contact Tpetra team.");
2468 #endif // HAVE_TPETRA_DEBUG
2469  }
2470 
2471 
2472  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2473  void
2474  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2475  getGlobalRowView (GlobalOrdinal globalRow,
2476  Teuchos::ArrayView<const GlobalOrdinal>& indices) const
2477  {
2478  const char tfecfFuncName[] = "getGlobalRowView: ";
2479  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2480  isLocallyIndexed (), std::runtime_error, "The graph's indices are "
2481  "currently stored as local indices, so we cannot return a view with "
2482  "global column indices. Use getGlobalRowCopy() instead.");
2483  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2484  ! hasRowInfo (), std::runtime_error,
2485  "Graph row information was deleted at fillComplete().");
2486 
2487  // isNodeGlobalElement() requires a global to local lookup anyway,
2488  // and getLocalElement() returns invalid() if the element wasn't found.
2489  const LocalOrdinal localRow = rowMap_->getLocalElement (globalRow);
2490  indices = Teuchos::null;
2491  if (localRow != Teuchos::OrdinalTraits<LocalOrdinal>::invalid ()) {
2492  const RowInfo rowInfo = getRowInfo (static_cast<size_t> (localRow));
2493  if (rowInfo.numEntries > 0) {
2494  indices = (this->getGlobalView (rowInfo)) (0, rowInfo.numEntries);
2495  }
2496  }
2497 #ifdef HAVE_TPETRA_DEBUG
2498  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2499  static_cast<size_t> (indices.size ()) != this->getNumEntriesInGlobalRow (globalRow),
2500  std::logic_error,
2501  "Violated stated postconditions: indices.size() = " << indices.size ()
2502  << " != getNumEntriesInGlobalRow(globalRow=" << globalRow
2503  << ") = " << this->getNumEntriesInGlobalRow (globalRow)
2504  << ". Please report this bug to the Tpetra developers.");
2505 #endif // HAVE_TPETRA_DEBUG
2506  }
2507 
2508 
2509  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2510  void
2511  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2512  insertLocalIndices (const LocalOrdinal localRow,
2513  const Teuchos::ArrayView<const LocalOrdinal>& indices)
2514  {
2515  const char tfecfFuncName[] = "insertLocalIndices";
2516 
2517  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2518  ! isFillActive (), std::runtime_error,
2519  ": requires that fill is active.");
2520  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2521  isGloballyIndexed (), std::runtime_error,
2522  ": graph indices are global; use insertGlobalIndices().");
2523  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2524  ! hasColMap (), std::runtime_error,
2525  ": cannot insert local indices without a column map.");
2526  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2527  ! rowMap_->isNodeLocalElement (localRow), std::runtime_error,
2528  ": row does not belong to this node.");
2529  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2530  ! hasRowInfo (), std::runtime_error,
2531  ": graph row information was deleted at fillComplete().");
2532  if (! indicesAreAllocated ()) {
2533  allocateIndices (LocalIndices);
2534  }
2535 
2536 #ifdef HAVE_TPETRA_DEBUG
2537  // In a debug build, if the graph has a column Map, test whether
2538  // any of the given column indices are not in the column Map.
2539  // Keep track of the invalid column indices so we can tell the
2540  // user about them.
2541  if (hasColMap ()) {
2542  using Teuchos::Array;
2543  using Teuchos::toString;
2544  using std::endl;
2545  typedef typename Teuchos::ArrayView<const LocalOrdinal>::size_type size_type;
2546 
2547  const map_type& colMap = * (getColMap ());
2548  Array<LocalOrdinal> badColInds;
2549  bool allInColMap = true;
2550  for (size_type k = 0; k < indices.size (); ++k) {
2551  if (! colMap.isNodeLocalElement (indices[k])) {
2552  allInColMap = false;
2553  badColInds.push_back (indices[k]);
2554  }
2555  }
2556  if (! allInColMap) {
2557  std::ostringstream os;
2558  os << "Tpetra::CrsMatrix::insertLocalIndices: You attempted to insert "
2559  "entries in owned row " << localRow << ", at the following column "
2560  "indices: " << toString (indices) << "." << endl;
2561  os << "Of those, the following indices are not in the column Map on "
2562  "this process: " << toString (badColInds) << "." << endl << "Since "
2563  "the graph has a column Map already, it is invalid to insert entries "
2564  "at those locations.";
2565  TEUCHOS_TEST_FOR_EXCEPTION(! allInColMap, std::invalid_argument, os.str ());
2566  }
2567  }
2568 #endif // HAVE_TPETRA_DEBUG
2569 
2570  insertLocalIndicesImpl (localRow, indices);
2571 
2572 #ifdef HAVE_TPETRA_DEBUG
2573  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2574  indicesAreAllocated() == false || isLocallyIndexed() == false,
2575  std::logic_error,
2576  ": Violated stated post-conditions. Please contact Tpetra team.");
2577 #endif
2578  }
2579 
2580 
2581  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2582  void
2583  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2584  insertLocalIndicesFiltered (const LocalOrdinal localRow,
2585  const Teuchos::ArrayView<const LocalOrdinal>& indices)
2586  {
2587  typedef LocalOrdinal LO;
2588  const char tfecfFuncName[] = "insertLocalIndicesFiltered";
2589 
2590  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2591  isFillActive() == false, std::runtime_error,
2592  ": requires that fill is active.");
2593  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2594  isGloballyIndexed() == true, std::runtime_error,
2595  ": graph indices are global; use insertGlobalIndices().");
2596  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2597  hasColMap() == false, std::runtime_error,
2598  ": cannot insert local indices without a column map.");
2599  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2600  rowMap_->isNodeLocalElement(localRow) == false, std::runtime_error,
2601  ": row does not belong to this node.");
2602  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2603  ! hasRowInfo (), std::runtime_error,
2604  ": graph row information was deleted at fillComplete().");
2605  if (! indicesAreAllocated ()) {
2606  allocateIndices (LocalIndices);
2607  }
2608 
2609  // If we have a column map, use it to filter the entries.
2610  if (hasColMap ()) {
2611  Teuchos::Array<LO> filtered_indices (indices);
2612  SLocalGlobalViews inds_view;
2613  SLocalGlobalNCViews inds_ncview;
2614  inds_ncview.linds = filtered_indices();
2615  const size_t numFilteredEntries =
2616  filterIndices<LocalIndices>(inds_ncview);
2617  inds_view.linds = filtered_indices (0, numFilteredEntries);
2618  insertLocalIndicesImpl(localRow, inds_view.linds);
2619  }
2620  else {
2621  insertLocalIndicesImpl(localRow, indices);
2622  }
2623 #ifdef HAVE_TPETRA_DEBUG
2624  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2625  indicesAreAllocated() == false || isLocallyIndexed() == false,
2626  std::logic_error,
2627  ": Violated stated post-conditions. Please contact Tpetra team.");
2628 #endif
2629  }
2630 
2631 
2632  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2633  void
2634  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2635  insertGlobalIndices (const GlobalOrdinal grow,
2636  const Teuchos::ArrayView<const GlobalOrdinal>& indices)
2637  {
2638  using Teuchos::ArrayView;
2639  typedef LocalOrdinal LO;
2640  typedef GlobalOrdinal GO;
2641  typedef typename ArrayView<const GO>::size_type size_type;
2642  const char tfecfFuncName[] = "insertGlobalIndices";
2643 
2644  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2645  isLocallyIndexed() == true, std::runtime_error,
2646  ": graph indices are local; use insertLocalIndices().");
2647  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2648  ! hasRowInfo (), std::runtime_error,
2649  ": graph row information was deleted at fillComplete().");
2650  // This can't really be satisfied for now, because if we are
2651  // fillComplete(), then we are local. In the future, this may
2652  // change. However, the rule that modification require active
2653  // fill will not change.
2654  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2655  isFillActive() == false, std::runtime_error,
2656  ": You are not allowed to call this method if fill is not active. "
2657  "If fillComplete has been called, you must first call resumeFill "
2658  "before you may insert indices.");
2659  if (! indicesAreAllocated ()) {
2660  allocateIndices (GlobalIndices);
2661  }
2662  const LO myRow = rowMap_->getLocalElement (grow);
2663  if (myRow != Teuchos::OrdinalTraits<LO>::invalid ()) {
2664 #ifdef HAVE_TPETRA_DEBUG
2665  if (hasColMap ()) {
2666  using std::endl;
2667  const map_type& colMap = * (getColMap ());
2668  // In a debug build, keep track of the nonowned ("bad") column
2669  // indices, so that we can display them in the exception
2670  // message. In a release build, just ditch the loop early if
2671  // we encounter a nonowned column index.
2672  Array<GO> badColInds;
2673  bool allInColMap = true;
2674  for (size_type k = 0; k < indices.size (); ++k) {
2675  if (! colMap.isNodeGlobalElement (indices[k])) {
2676  allInColMap = false;
2677  badColInds.push_back (indices[k]);
2678  }
2679  }
2680  if (! allInColMap) {
2681  std::ostringstream os;
2682  os << "Tpetra::CrsGraph::insertGlobalIndices: You attempted to insert "
2683  "entries in owned row " << grow << ", at the following column "
2684  "indices: " << toString (indices) << "." << endl;
2685  os << "Of those, the following indices are not in the column Map on "
2686  "this process: " << toString (badColInds) << "." << endl << "Since "
2687  "the matrix has a column Map already, it is invalid to insert "
2688  "entries at those locations.";
2689  TEUCHOS_TEST_FOR_EXCEPTION(! allInColMap, std::invalid_argument, os.str ());
2690  }
2691  }
2692 #endif // HAVE_TPETRA_DEBUG
2693  insertGlobalIndicesImpl (myRow, indices);
2694  }
2695  else { // a nonlocal row
2696  const size_type numIndices = indices.size ();
2697  // This creates the Array if it doesn't exist yet. std::map's
2698  // operator[] does a lookup each time, so it's better to pull
2699  // nonlocals_[grow] out of the loop.
2700  std::vector<GO>& nonlocalRow = nonlocals_[grow];
2701  for (size_type k = 0; k < numIndices; ++k) {
2702  nonlocalRow.push_back (indices[k]);
2703  }
2704  }
2705 #ifdef HAVE_TPETRA_DEBUG
2706  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2707  indicesAreAllocated() == false || isGloballyIndexed() == false,
2708  std::logic_error,
2709  ": Violated stated post-conditions. Please contact Tpetra team.");
2710 #endif
2711  }
2712 
2713 
2714  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2715  void
2716  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2717  insertGlobalIndicesFiltered (const GlobalOrdinal grow,
2718  const Teuchos::ArrayView<const GlobalOrdinal>& indices)
2719  {
2720  using Teuchos::Array;
2721  using Teuchos::ArrayView;
2722  typedef LocalOrdinal LO;
2723  typedef GlobalOrdinal GO;
2724  const char tfecfFuncName[] = "insertGlobalIndicesFiltered";
2725 
2726  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2727  isLocallyIndexed() == true, std::runtime_error,
2728  ": graph indices are local; use insertLocalIndices().");
2729  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2730  ! hasRowInfo (), std::runtime_error,
2731  ": graph row information was deleted at fillComplete().");
2732  // This can't really be satisfied for now, because if we are
2733  // fillComplete(), then we are local. In the future, this may
2734  // change. However, the rule that modification require active
2735  // fill will not change.
2736  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2737  isFillActive() == false, std::runtime_error,
2738  ": You are not allowed to call this method if fill is not active. "
2739  "If fillComplete has been called, you must first call resumeFill "
2740  "before you may insert indices.");
2741  if (! indicesAreAllocated ()) {
2742  allocateIndices (GlobalIndices);
2743  }
2744  const LO myRow = rowMap_->getLocalElement (grow);
2745  if (myRow != Teuchos::OrdinalTraits<LO>::invalid ()) {
2746  // If we have a column map, use it to filter the entries.
2747  if (hasColMap ()) {
2748  Array<GO> filtered_indices(indices);
2749  SLocalGlobalViews inds_view;
2750  SLocalGlobalNCViews inds_ncview;
2751  inds_ncview.ginds = filtered_indices();
2752  const size_t numFilteredEntries =
2753  filterIndices<GlobalIndices> (inds_ncview);
2754  inds_view.ginds = filtered_indices (0, numFilteredEntries);
2755  insertGlobalIndicesImpl(myRow, inds_view.ginds);
2756  }
2757  else {
2758  insertGlobalIndicesImpl(myRow, indices);
2759  }
2760  }
2761  else { // a nonlocal row
2762  typedef typename ArrayView<const GO>::size_type size_type;
2763  const size_type numIndices = indices.size ();
2764  for (size_type k = 0; k < numIndices; ++k) {
2765  nonlocals_[grow].push_back (indices[k]);
2766  }
2767  }
2768 #ifdef HAVE_TPETRA_DEBUG
2769  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2770  indicesAreAllocated() == false || isGloballyIndexed() == false,
2771  std::logic_error,
2772  ": Violated stated post-conditions. Please contact Tpetra team.");
2773 #endif
2774  }
2775 
2776 
2777  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2778  void
2779  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2780  removeLocalIndices (LocalOrdinal lrow)
2781  {
2782  const char tfecfFuncName[] = "removeLocalIndices: ";
2783  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2784  ! isFillActive (), std::runtime_error, "requires that fill is active.");
2785  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2786  isStorageOptimized (), std::runtime_error,
2787  "cannot remove indices after optimizeStorage() has been called.");
2788  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2789  isGloballyIndexed (), std::runtime_error, "graph indices are global.");
2790  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2791  ! rowMap_->isNodeLocalElement (lrow), std::runtime_error,
2792  "Local row " << lrow << " is not in the row Map on the calling process.");
2793  if (! indicesAreAllocated ()) {
2794  allocateIndices (LocalIndices);
2795  }
2796 
2797  // FIXME (mfh 13 Aug 2014) What if they haven't been cleared on
2798  // all processes?
2799  clearGlobalConstants ();
2800 
2801  if (k_numRowEntries_.dimension_0 () != 0) {
2802  const size_t oldNumEntries = k_numRowEntries_.h_view (lrow);
2803  nodeNumEntries_ -= oldNumEntries;
2804 
2805  // FIXME (mfh 07 Aug 2014) We should just sync at fillComplete,
2806  // but for now, for correctness, do the modify-sync cycle here.
2807  k_numRowEntries_.template modify<typename t_numRowEntries_::host_mirror_space> ();
2808  k_numRowEntries_.h_view(lrow) = 0;
2809  k_numRowEntries_.template sync<execution_space> ();
2810  }
2811 #ifdef HAVE_TPETRA_DEBUG
2812  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2813  getNumEntriesInLocalRow (lrow) != 0 ||
2814  ! indicesAreAllocated () ||
2815  ! isLocallyIndexed (), std::logic_error,
2816  ": Violated stated post-conditions. Please contact Tpetra team.");
2817 #endif // HAVE_TPETRA_DEBUG
2818  }
2819 
2820 
2821  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2822  void
2823  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2824  setAllIndices (const typename local_graph_type::row_map_type& rowPointers,
2825  const typename local_graph_type::entries_type::non_const_type& columnIndices)
2826  {
2827  const char tfecfFuncName[] = "setAllIndices: ";
2828  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2829  ! hasColMap () || getColMap ().is_null (), std::runtime_error,
2830  "The graph must have a column Map before you may call this method.");
2831  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2832  static_cast<size_t> (rowPointers.size ()) != this->getNodeNumRows () + 1,
2833  std::runtime_error, "rowPointers.size() = " << rowPointers.size () <<
2834  " != this->getNodeNumRows()+1 = " << (this->getNodeNumRows () + 1) <<
2835  ".");
2836 
2837  // FIXME (mfh 07 Aug 2014) We need to relax this restriction,
2838  // since the future model will be allocation at construction, not
2839  // lazy allocation on first insert.
2840  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2841  k_lclInds1D_.dimension_0 () != 0 || k_gblInds1D_.dimension_0 () != 0,
2842  std::runtime_error, "You may not call this method if 1-D data structures "
2843  "are already allocated.");
2844 
2845  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
2846  lclInds2D_ != Teuchos::null || gblInds2D_ != Teuchos::null,
2847  std::runtime_error, "You may not call this method if 2-D data structures "
2848  "are already allocated.");
2849 
2850  const size_t localNumEntries = rowPointers(getNodeNumRows ());
2851 
2852  indicesAreAllocated_ = true;
2853  indicesAreLocal_ = true;
2854  pftype_ = StaticProfile; // if the profile wasn't static before, it sure is now.
2855  k_lclInds1D_ = columnIndices;
2856  k_rowPtrs_ = rowPointers;
2857  nodeNumAllocated_ = localNumEntries;
2858  nodeNumEntries_ = localNumEntries;
2859 
2860  // These normally get cleared out at the end of allocateIndices.
2861  // It makes sense to clear them out here, because at the end of
2862  // this method, the graph is allocated on the calling process.
2863  numAllocForAllRows_ = 0;
2864  k_numAllocPerRow_ = Kokkos::DualView<const size_t*, execution_space> ();
2865 
2866  checkInternalState ();
2867  }
2868 
2869 
2870  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2871  void
2872  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2873  setAllIndices (const Teuchos::ArrayRCP<size_t>& rowPointers,
2874  const Teuchos::ArrayRCP<LocalOrdinal>& columnIndices)
2875  {
2876  using Kokkos::View;
2877  typedef typename local_graph_type::row_map_type row_map_type;
2878  typedef typename row_map_type::array_layout layout_type;
2879  typedef typename row_map_type::non_const_value_type row_offset_type;
2880  typedef View<size_t*, layout_type , Kokkos::HostSpace,
2881  Kokkos::MemoryUnmanaged> input_view_type;
2882  typedef typename row_map_type::non_const_type nc_row_map_type;
2883 
2884  const size_t size = static_cast<size_t> (rowPointers.size ());
2885  const bool same = Kokkos::Impl::is_same<size_t, row_offset_type>::value;
2886  input_view_type ptr_in (rowPointers.getRawPtr (), size);
2887 
2888  nc_row_map_type ptr_rot ("Tpetra::CrsGraph::ptr", size);
2889 
2890  if (same) { // size_t == row_offset_type
2891  // This compile-time logic ensures that the compiler never sees
2892  // an assignment of View<row_offset_type*, ...> to View<size_t*,
2893  // ...> unless size_t == row_offset_type.
2894  input_view_type ptr_decoy (rowPointers.getRawPtr (), size); // never used
2895  Kokkos::deep_copy (Kokkos::Impl::if_c<same,
2896  nc_row_map_type,
2897  input_view_type>::select (ptr_rot, ptr_decoy),
2898  ptr_in);
2899  }
2900  else { // size_t != row_offset_type
2901  // CudaUvmSpace != HostSpace, so this will be false in that case.
2902  const bool inHostMemory =
2903  Kokkos::Impl::is_same<typename row_map_type::memory_space,
2904  Kokkos::HostSpace>::value;
2905  if (inHostMemory) {
2906  // Copy (with cast from size_t to row_offset_type) to ptr_rot.
2907  typedef CopyOffsets<nc_row_map_type, input_view_type> functor_type;
2908  functor_type functor (ptr_rot, ptr_in);
2909  Kokkos::parallel_for (size, functor);
2910  }
2911  else { // Copy input row offsets to device first.
2912  //
2913  // FIXME (mfh 24 Mar 2015) If CUDA UVM, running in the host's
2914  // execution space would avoid the double copy.
2915  //
2916  View<size_t*, layout_type ,execution_space > ptr_st ("Tpetra::CrsGraph::ptr", size);
2917  Kokkos::deep_copy (ptr_st, ptr_in);
2918  // Copy on device (casting from size_t to row_offset_type) to
2919  // ptr_rot. This executes in the output View's execution
2920  // space, which is the same as execution_space.
2921  typedef CopyOffsets<nc_row_map_type,
2922  View<size_t*, layout_type, execution_space> > functor_type;
2923  functor_type functor (ptr_rot, ptr_st);
2924  Kokkos::parallel_for (size, functor);
2925  }
2926  }
2927 
2928  Kokkos::View<LocalOrdinal*, layout_type , execution_space > k_ind =
2929  Kokkos::Compat::getKokkosViewDeepCopy<device_type> (columnIndices ());
2930  setAllIndices (ptr_rot, k_ind);
2931  }
2932 
2933 
2934  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
2935  void
2936  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
2937  getNumEntriesPerLocalRowUpperBound (Teuchos::ArrayRCP<const size_t>& boundPerLocalRow,
2938  size_t& boundForAllLocalRows,
2939  bool& boundSameForAllLocalRows) const
2940  {
2941  // The three output arguments. We assign them to the actual
2942  // output arguments at the end, in order to implement
2943  // transactional semantics.
2944  Teuchos::ArrayRCP<const size_t> numEntriesPerRow;
2945  size_t numEntriesForAll = 0;
2946  bool allRowsSame = true;
2947 
2948  const ptrdiff_t numRows = static_cast<ptrdiff_t> (this->getNodeNumRows ());
2949 
2950  if (! this->indicesAreAllocated ()) {
2951  if (k_numAllocPerRow_.dimension_0 () != 0) {
2952  numEntriesPerRow = Kokkos::Compat::persistingView (k_numAllocPerRow_.h_view);
2953  allRowsSame = false; // conservatively; we don't check the array
2954  }
2955  else {
2956  numEntriesForAll = numAllocForAllRows_;
2957  allRowsSame = true;
2958  }
2959  }
2960  else if (k_numRowEntries_.dimension_0 () != 0) {
2961  numEntriesPerRow = Kokkos::Compat::persistingView (k_numRowEntries_.h_view);
2962  allRowsSame = false; // conservatively; we don't check the array
2963  }
2964  else if (this->nodeNumAllocated_ == 0) {
2965  numEntriesForAll = 0;
2966  allRowsSame = true;
2967  }
2968  else {
2969  // left with the case that we have optimized storage. in this
2970  // case, we have to construct a list of row sizes.
2971  TEUCHOS_TEST_FOR_EXCEPTION(
2972  this->getProfileType () != StaticProfile, std::logic_error,
2973  "Tpetra::CrsGraph::getNumEntriesPerRowUpperBound: "
2974  "The graph is not StaticProfile, but storage appears to be optimized. "
2975  "Please report this bug to the Tpetra developers.");
2976  TEUCHOS_TEST_FOR_EXCEPTION(
2977  numRows != 0 && k_rowPtrs_.dimension_0 () == 0, std::logic_error,
2978  "Tpetra::CrsGraph::getNumEntriesPerRowUpperBound: "
2979  "The graph has " << numRows << " (> 0) row" << (numRows != 1 ? "s" : "")
2980  << " on the calling process, but the k_rowPtrs_ array has zero entries. "
2981  "Please report this bug to the Tpetra developers.");
2982 
2983  Teuchos::ArrayRCP<size_t> numEnt;
2984  if (numRows != 0) {
2985  numEnt = Teuchos::arcp<size_t> (numRows);
2986  }
2987 
2988  // We have to iterate through the row offsets anyway, so we
2989  // might as well check whether all rows' bounds are the same.
2990  bool allRowsReallySame = false;
2991  for (ptrdiff_t i = 0; i < numRows; ++i) {
2992  numEnt[i] = k_rowPtrs_(i+1) - k_rowPtrs_(i);
2993  if (i != 0 && numEnt[i] != numEnt[i-1]) {
2994  allRowsReallySame = false;
2995  }
2996  }
2997  if (allRowsReallySame) {
2998  if (numRows == 0) {
2999  numEntriesForAll = 0;
3000  } else {
3001  numEntriesForAll = numEnt[1] - numEnt[0];
3002  }
3003  allRowsSame = true;
3004  }
3005  else {
3006  numEntriesPerRow = numEnt; // Teuchos::arcp_const_cast<const size_t> (numEnt);
3007  allRowsSame = false; // conservatively; we don't check the array
3008  }
3009  }
3010 
3011  TEUCHOS_TEST_FOR_EXCEPTION(
3012  numEntriesForAll != 0 && numEntriesPerRow.size () != 0, std::logic_error,
3013  "Tpetra::CrsGraph::getNumEntriesPerLocalRowUpperBound: "
3014  "numEntriesForAll and numEntriesPerRow are not consistent. The former "
3015  "is nonzero (" << numEntriesForAll << "), but the latter has nonzero "
3016  "size " << numEntriesPerRow.size () << ". "
3017  "Please report this bug to the Tpetra developers.");
3018  TEUCHOS_TEST_FOR_EXCEPTION(
3019  numEntriesForAll != 0 && ! allRowsSame, std::logic_error,
3020  "Tpetra::CrsGraph::getNumEntriesPerLocalRowUpperBound: "
3021  "numEntriesForAll and allRowsSame are not consistent. The former "
3022  "is nonzero (" << numEntriesForAll << "), but the latter is false. "
3023  "Please report this bug to the Tpetra developers.");
3024  TEUCHOS_TEST_FOR_EXCEPTION(
3025  numEntriesPerRow.size () != 0 && allRowsSame, std::logic_error,
3026  "Tpetra::CrsGraph::getNumEntriesPerLocalRowUpperBound: "
3027  "numEntriesPerRow and allRowsSame are not consistent. The former has "
3028  "nonzero length " << numEntriesForAll << ", but the latter is true. "
3029  "Please report this bug to the Tpetra developers.");
3030 
3031  boundPerLocalRow = numEntriesPerRow;
3032  boundForAllLocalRows = numEntriesForAll;
3033  boundSameForAllLocalRows = allRowsSame;
3034  }
3035 
3036 
3037  // TODO: in the future, globalAssemble() should use import/export functionality
3038  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3039  void
3040  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3041  globalAssemble ()
3042  {
3043  using Teuchos::Array;
3044  using Teuchos::as;
3045  using Teuchos::Comm;
3046  using Teuchos::gatherAll;
3047  using Teuchos::ireceive;
3048  using Teuchos::isend;
3049  using Teuchos::outArg;
3050  using Teuchos::REDUCE_MAX;
3051  using Teuchos::reduceAll;
3052  using Teuchos::toString;
3053  using Teuchos::waitAll;
3054  using std::endl;
3055  using std::make_pair;
3056  using std::pair;
3057  typedef GlobalOrdinal GO;
3058  typedef typename std::map<GO, std::vector<GO> >::const_iterator NLITER;
3059  typedef typename Array<GO>::size_type size_type;
3060 
3061  const char tfecfFuncName[] = "globalAssemble"; // for exception macro
3062  RCP<const Comm<int> > comm = getComm();
3063 
3064  const int numImages = comm->getSize();
3065  const int myImageID = comm->getRank();
3066 #ifdef HAVE_TPETRA_DEBUG
3067  Teuchos::barrier (*comm);
3068 #endif // HAVE_TPETRA_DEBUG
3069 
3070  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC( ! isFillActive(), std::runtime_error,
3071  ": requires that fill is active.");
3072  // Determine if any nodes have global entries to share
3073  {
3074  size_t MyNonlocals = nonlocals_.size(), MaxGlobalNonlocals;
3075  reduceAll (*comm, REDUCE_MAX, MyNonlocals, outArg (MaxGlobalNonlocals));
3076  if (MaxGlobalNonlocals == 0) {
3077  return; // no entries to share
3078  }
3079  }
3080 
3081  // compute a list of NLRs from nonlocals_ and use it to compute:
3082  // IdsAndRows: a vector of (id,row) pairs
3083  // NLR2Id: a map from NLR to the Id that owns it
3084  // globalNeighbors: a global graph of connectivity between images:
3085  // globalNeighbors(i,j) indicates that j sends to i
3086  // sendIDs: a list of all images I send to
3087  // recvIDs: a list of all images I receive from (constructed later)
3088  Array<pair<int, GO> > IdsAndRows;
3089  std::map<GO, int> NLR2Id;
3090  Teuchos::SerialDenseMatrix<int, char> globalNeighbors;
3091  Array<int> sendIDs, recvIDs;
3092  {
3093  // nonlocals_ contains the entries we are holding for all
3094  // nonowned rows. Compute list of rows for which we have data.
3095  Array<GO> NLRs;
3096  std::set<GO> setOfRows;
3097  for (NLITER iter = nonlocals_.begin (); iter != nonlocals_.end (); ++iter) {
3098  setOfRows.insert (iter->first);
3099  }
3100  // copy the elements in the set into an Array
3101  NLRs.resize (setOfRows.size ());
3102  std::copy (setOfRows.begin (), setOfRows.end (), NLRs.begin ());
3103 
3104  // get a list of ImageIDs for the non-local rows (NLRs)
3105  Array<int> NLRIds(NLRs.size());
3106  {
3107  const LookupStatus stat =
3108  rowMap_->getRemoteIndexList (NLRs (), NLRIds ());
3109  int lclerror = ( stat == IDNotPresent ? 1 : 0 );
3110  int gblerror;
3111  reduceAll<int, int> (*comm, REDUCE_MAX, lclerror, outArg (gblerror));
3112  if (gblerror != 0) {
3113  const int myRank = comm->getRank ();
3114  std::ostringstream os;
3115  os << "On one or more processes in the communicator, "
3116  << "there were insertions into rows of the graph that do not "
3117  << "exist in the row Map on any process in the communicator."
3118  << endl << "This process " << myRank << " is "
3119  << (lclerror == 0 ? "not " : "") << "one of those offending "
3120  << "processes." << endl;
3121  if (lclerror != 0) {
3122  // If NLRIds[k] is -1, then NLRs[k] is a row index not in
3123  // the row Map. Collect this list of invalid row indices
3124  // for display in the exception message.
3125  Array<GO> invalidNonlocalRows;
3126  for (size_type k = 0; k < NLRs.size (); ++k) {
3127  if (NLRIds[k] == -1) {
3128  invalidNonlocalRows.push_back (NLRs[k]);
3129  }
3130  }
3131  const size_type numInvalid = invalidNonlocalRows.size ();
3132  os << "On this process, " << numInvalid << " nonlocal row"
3133  << (numInvalid != 1 ? "s " : " ") << " were inserted that are "
3134  << "not in the row Map on any process." << endl;
3135  // Don't print _too_ many nonlocal rows.
3136  if (numInvalid <= 100) {
3137  os << "Offending row indices: "
3138  << toString (invalidNonlocalRows ()) << endl;
3139  }
3140  }
3141  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3142  gblerror != 0, std::runtime_error,
3143  ": nonlocal entries correspond to invalid rows."
3144  << endl << os.str ());
3145  }
3146  }
3147 
3148  // build up a list of neighbors, as well as a map between NLRs and Ids
3149  // localNeighbors[i] != 0 iff I have data to send to image i
3150  // put NLRs,Ids into an array of pairs
3151  IdsAndRows.reserve(NLRs.size());
3152  Array<char> localNeighbors(numImages, 0);
3153  typename Array<GO>::const_iterator nlr;
3154  typename Array<int>::const_iterator id;
3155  for (nlr = NLRs.begin(), id = NLRIds.begin();
3156  nlr != NLRs.end(); ++nlr, ++id) {
3157  NLR2Id[*nlr] = *id;
3158  localNeighbors[*id] = 1;
3159  // IdsAndRows.push_back(make_pair<int,GlobalOrdinal>(*id,*nlr));
3160  IdsAndRows.push_back(make_pair(*id,*nlr));
3161  }
3162  for (int j=0; j<numImages; ++j) {
3163  if (localNeighbors[j]) {
3164  sendIDs.push_back(j);
3165  }
3166  }
3167  // sort IdsAndRows, by Ids first, then rows
3168  std::sort(IdsAndRows.begin(),IdsAndRows.end());
3169  // gather from other nodes to form the full graph
3170  globalNeighbors.shapeUninitialized(numImages,numImages);
3171  gatherAll (*getComm(), numImages, localNeighbors.getRawPtr(),
3172  numImages * numImages, globalNeighbors.values());
3173  // globalNeighbors at this point contains (on all images) the
3174  // connectivity between the images.
3175  // globalNeighbors(i,j) != 0 means that j sends to i/that i receives from j
3176  }
3177 
3179  // FIGURE OUT WHO IS SENDING TO WHOM AND HOW MUCH
3180  // DO THIS IN THE PROCESS OF PACKING ALL OUTGOING DATA ACCORDING TO DESTINATION ID
3182 
3183  // loop over all columns to know from which images I can expect to receive something
3184  for (int j=0; j<numImages; ++j) {
3185  if (globalNeighbors(myImageID,j)) {
3186  recvIDs.push_back(j);
3187  }
3188  }
3189  const size_t numRecvs = recvIDs.size();
3190 
3191  // we know how many we're sending to already
3192  // form a contiguous list of all data to be sent
3193  // track the number of entries for each ID
3194  Array<pair<GO, GO> > IJSendBuffer;
3195  Array<size_t> sendSizes(sendIDs.size(), 0);
3196  size_t numSends = 0;
3197  for (typename Array<pair<int, GO> >::const_iterator IdAndRow = IdsAndRows.begin();
3198  IdAndRow != IdsAndRows.end(); ++IdAndRow) {
3199  int id = IdAndRow->first;
3200  GO row = IdAndRow->second;
3201  // have we advanced to a new send?
3202  if (sendIDs[numSends] != id) {
3203  numSends++;
3204  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(sendIDs[numSends] != id,
3205  std::logic_error, ": internal logic error. Contact Tpetra team.");
3206  }
3207  // copy data for row into contiguous storage
3208  for (typename std::vector<GO>::const_iterator j = nonlocals_[row].begin ();
3209  j != nonlocals_[row].end (); ++j) {
3210  IJSendBuffer.push_back (pair<GlobalOrdinal, GlobalOrdinal> (row, *j));
3211  sendSizes[numSends]++;
3212  }
3213  }
3214  if (IdsAndRows.size() > 0) {
3215  numSends++; // one last increment, to make it a count instead of an index
3216  }
3217  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3218  as<typename Array<int>::size_type> (numSends) != sendIDs.size (),
3219  std::logic_error, ": internal logic error. Contact Tpetra team.");
3220 
3221  // don't need this data anymore
3222  nonlocals_.clear();
3223 
3225  // TRANSMIT SIZE INFO BETWEEN SENDERS AND RECEIVERS
3227 
3228  // Array of pending nonblocking communication requests. It's OK
3229  // to mix nonblocking send and receive requests in the same
3230  // waitAll() call.
3231  Array<RCP<Teuchos::CommRequest<int> > > requests;
3232 
3233  // perform non-blocking sends: send sizes to our recipients
3234  for (size_t s = 0; s < numSends ; ++s) {
3235  // We're using a nonowning RCP because all communication
3236  // will be local to this method and the scope of our data
3237  requests.push_back (isend<int, size_t> (*comm,
3238  rcp (&sendSizes[s], false),
3239  sendIDs[s]));
3240  }
3241  // perform non-blocking receives: receive sizes from our senders
3242  Array<size_t> recvSizes (numRecvs);
3243  for (size_t r = 0; r < numRecvs; ++r) {
3244  // We're using a nonowning RCP because all communication
3245  // will be local to this method and the scope of our data
3246  requests.push_back (ireceive (*comm, rcp (&recvSizes[r], false), recvIDs[r]));
3247  }
3248  // Wait on all the nonblocking sends and receives.
3249  if (! requests.empty()) {
3250  waitAll (*comm, requests());
3251  }
3252 #ifdef HAVE_TPETRA_DEBUG
3253  Teuchos::barrier (*comm);
3254 #endif // HAVE_TPETRA_DEBUG
3255 
3256  // This doesn't necessarily deallocate the array.
3257  requests.resize (0);
3258 
3260  // NOW SEND/RECEIVE ALL ROW DATA
3262  // from the size info, build the ArrayViews into IJSendBuffer
3263  Array<ArrayView<pair<GO,GO> > > sendBuffers(numSends,null);
3264  {
3265  size_t cur = 0;
3266  for (size_t s=0; s<numSends; ++s) {
3267  sendBuffers[s] = IJSendBuffer(cur,sendSizes[s]);
3268  cur += sendSizes[s];
3269  }
3270  }
3271  // perform non-blocking sends
3272  for (size_t s=0; s < numSends ; ++s)
3273  {
3274  // We're using a nonowning RCP because all communication
3275  // will be local to this method and the scope of our data
3276  ArrayRCP<pair<GO,GO> > tmpSendBuf =
3277  arcp (sendBuffers[s].getRawPtr(), 0, sendBuffers[s].size(), false);
3278  requests.push_back (isend<int, pair<GO,GO> > (*comm, tmpSendBuf, sendIDs[s]));
3279  }
3280  // calculate amount of storage needed for receives
3281  // setup pointers for the receives as well
3282  size_t totalRecvSize = std::accumulate (recvSizes.begin(), recvSizes.end(), 0);
3283  Array<pair<GO,GO> > IJRecvBuffer (totalRecvSize);
3284  // from the size info, build the ArrayViews into IJRecvBuffer
3285  Array<ArrayView<pair<GO,GO> > > recvBuffers (numRecvs, null);
3286  {
3287  size_t cur = 0;
3288  for (size_t r=0; r<numRecvs; ++r) {
3289  recvBuffers[r] = IJRecvBuffer(cur,recvSizes[r]);
3290  cur += recvSizes[r];
3291  }
3292  }
3293  // perform non-blocking recvs
3294  for (size_t r = 0; r < numRecvs; ++r) {
3295  // We're using a nonowning RCP because all communication
3296  // will be local to this method and the scope of our data
3297  ArrayRCP<pair<GO,GO> > tmpRecvBuf =
3298  arcp (recvBuffers[r].getRawPtr(), 0, recvBuffers[r].size(), false);
3299  requests.push_back (ireceive (*comm, tmpRecvBuf, recvIDs[r]));
3300  }
3301  // perform waits
3302  if (! requests.empty()) {
3303  waitAll (*comm, requests());
3304  }
3305 #ifdef HAVE_TPETRA_DEBUG
3306  Teuchos::barrier (*comm);
3307 #endif // HAVE_TPETRA_DEBUG
3308 
3310  // NOW PROCESS THE RECEIVED ROW DATA
3312  // TODO: instead of adding one entry at a time, add one row at a time.
3313  // this requires resorting; they arrived sorted by sending node,
3314  // so that entries could be non-contiguous if we received
3315  // multiple entries for a particular row from different processors.
3316  // it also requires restoring the data, which may make it not worth the trouble.
3317  for (typename Array<pair<GO,GO> >::const_iterator ij = IJRecvBuffer.begin();
3318  ij != IJRecvBuffer.end(); ++ij)
3319  {
3320  insertGlobalIndicesFiltered (ij->first, tuple<GO> (ij->second));
3321  }
3322  checkInternalState();
3323  }
3324 
3325 
3326  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3327  void
3328  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3329  resumeFill (const Teuchos::RCP<Teuchos::ParameterList>& params)
3330  {
3331  const char tfecfFuncName[] = "resumeFill";
3332  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(! hasRowInfo(), std::runtime_error,
3333  ": Sorry, you cannot resume fill of the CrsGraph, since the graph's row "
3334  "information was deleted in fillComplete().");
3335 
3336 #ifdef HAVE_TPETRA_DEBUG
3337  Teuchos::barrier( *rowMap_->getComm() );
3338 #endif // HAVE_TPETRA_DEBUG
3339  clearGlobalConstants();
3340  if (params != null) this->setParameterList (params);
3341  lowerTriangular_ = false;
3342  upperTriangular_ = false;
3343  // either still sorted/merged or initially sorted/merged
3344  indicesAreSorted_ = true;
3345  noRedundancies_ = true;
3346  fillComplete_ = false;
3347 #ifdef HAVE_TPETRA_DEBUG
3348  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3349  ! isFillActive() || isFillComplete(), std::logic_error,
3350  "::resumeFill(): At end of method, either fill is not active or fill is "
3351  "complete. This violates stated post-conditions. Please report this bug "
3352  "to the Tpetra developers.");
3353 #endif // HAVE_TPETRA_DEBUG
3354  }
3355 
3356 
3357  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3358  void
3359  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3360  fillComplete (const Teuchos::RCP<Teuchos::ParameterList>& params)
3361  {
3362  // If the graph already has domain and range Maps, don't clobber
3363  // them. If it doesn't, use the current row Map for both the
3364  // domain and range Maps.
3365  //
3366  // NOTE (mfh 28 Sep 2014): If the graph was constructed without a
3367  // column Map, and column indices are inserted which are not in
3368  // the row Map on any process, this will cause troubles. However,
3369  // that is not a common case for most applications that we
3370  // encounter, and checking for it might require more
3371  // communication.
3372  Teuchos::RCP<const map_type> domMap = this->getDomainMap ();
3373  if (domMap.is_null ()) {
3374  domMap = this->getRowMap ();
3375  }
3376  Teuchos::RCP<const map_type> ranMap = this->getRangeMap ();
3377  if (ranMap.is_null ()) {
3378  ranMap = this->getRowMap ();
3379  }
3380  this->fillComplete (domMap, ranMap, params);
3381  }
3382 
3383 
3384  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3385  void
3386  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3387  fillComplete (const Teuchos::RCP<const map_type>& domainMap,
3388  const Teuchos::RCP<const map_type>& rangeMap,
3389  const Teuchos::RCP<Teuchos::ParameterList>& params)
3390  {
3391  const char tfecfFuncName[] = "fillComplete";
3392 
3393 #ifdef HAVE_TPETRA_DEBUG
3394  rowMap_->getComm ()->barrier ();
3395 #endif // HAVE_TPETRA_DEBUG
3396 
3397  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC( ! isFillActive() || isFillComplete(),
3398  std::runtime_error, ": Graph fill state must be active (isFillActive() "
3399  "must be true) before calling fillComplete().");
3400 
3401  const int numProcs = getComm ()->getSize ();
3402 
3403  //
3404  // Read and set parameters
3405  //
3406 
3407  // Does the caller want to sort remote GIDs (within those owned by
3408  // the same process) in makeColMap()?
3409  if (! params.is_null ()) {
3410  if (params->isParameter ("sort column map ghost gids")) {
3411  sortGhostsAssociatedWithEachProcessor_ =
3412  params->get<bool> ("sort column map ghost gids",
3413  sortGhostsAssociatedWithEachProcessor_);
3414  }
3415  else if (params->isParameter ("Sort column Map ghost GIDs")) {
3416  sortGhostsAssociatedWithEachProcessor_ =
3417  params->get<bool> ("Sort column Map ghost GIDs",
3418  sortGhostsAssociatedWithEachProcessor_);
3419  }
3420  }
3421 
3422  // If true, the caller promises that no process did nonlocal
3423  // changes since the last call to fillComplete.
3424  bool assertNoNonlocalInserts = false;
3425  if (! params.is_null ()) {
3426  assertNoNonlocalInserts =
3427  params->get<bool> ("No Nonlocal Changes", assertNoNonlocalInserts);
3428  }
3429 
3430  //
3431  // Allocate indices, if they haven't already been allocated
3432  //
3433  if (! indicesAreAllocated ()) {
3434  if (hasColMap ()) {
3435  // We have a column Map, so use local indices.
3436  allocateIndices (LocalIndices);
3437  } else {
3438  // We don't have a column Map, so use global indices.
3439  allocateIndices (GlobalIndices);
3440  }
3441  }
3442 
3443  //
3444  // Do global assembly, if requested and if the communicator
3445  // contains more than one process.
3446  //
3447  const bool mayNeedGlobalAssemble = ! assertNoNonlocalInserts && numProcs > 1;
3448  if (mayNeedGlobalAssemble) {
3449  // This first checks if we need to do global assembly.
3450  // The check costs a single all-reduce.
3451  globalAssemble ();
3452  }
3453  else {
3454  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3455  numProcs > 1 && nonlocals_.size() > 0, std::runtime_error,
3456  ":" << std::endl << "The graph's communicator contains only one "
3457  "process, but there are nonlocal entries. " << std::endl <<
3458  "This probably means that invalid entries were added to the graph.");
3459  }
3460 
3461  // Set domain and range Map. This may clear the Import / Export
3462  // objects if the new Maps differ from any old ones.
3463  setDomainRangeMaps (domainMap, rangeMap);
3464 
3465  // If the graph does not already have a column Map (either from
3466  // the user constructor calling the version of the constructor
3467  // that takes a column Map, or from a previous fillComplete call),
3468  // then create it.
3469  if (! hasColMap ()) {
3470  makeColMap ();
3471  }
3472 
3473  // Make indices local, if they aren't already.
3474  // The method doesn't do any work if the indices are already local.
3475  makeIndicesLocal ();
3476 
3477  if (! isSorted ()) {
3478  // If this process has no indices, then CrsGraph considers it
3479  // already trivially sorted. Thus, this method need not be
3480  // called on all processes in the row Map's communicator.
3481  sortAllIndices ();
3482  }
3483 
3484  if (! isMerged()) {
3485  mergeAllIndices ();
3486  }
3487  makeImportExport (); // Make Import and Export objects, if necessary
3488  computeGlobalConstants ();
3489  fillLocalGraph (params);
3490  fillComplete_ = true;
3491 
3492 #ifdef HAVE_TPETRA_DEBUG
3493  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3494  isFillActive() == true || isFillComplete() == false, std::logic_error,
3495  ": Violated stated post-conditions. Please contact Tpetra team.");
3496 #endif
3497 
3498  checkInternalState ();
3499  }
3500 
3501 
3502  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3503  void
3504  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3505  expertStaticFillComplete (const Teuchos::RCP<const map_type>& domainMap,
3506  const Teuchos::RCP<const map_type>& rangeMap,
3507  const Teuchos::RCP<const import_type>& importer,
3508  const Teuchos::RCP<const export_type>& exporter,
3509  const Teuchos::RCP<Teuchos::ParameterList>& params)
3510  {
3511  const char tfecfFuncName[] = "expertStaticFillComplete: ";
3512 #ifdef HAVE_TPETRA_MMM_TIMINGS
3513  std::string label;
3514  if(!params.is_null())
3515  label = params->get("Timer Label",label);
3516  std::string prefix = std::string("Tpetra ")+ label + std::string(": ");
3517  using Teuchos::TimeMonitor;
3518  Teuchos::RCP<Teuchos::TimeMonitor> MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-Setup"))));
3519 #endif
3520 
3521 
3522  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3523  domainMap.is_null () || rangeMap.is_null (),
3524  std::runtime_error, "The input domain Map and range Map must be nonnull.");
3525  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3526  pftype_ != StaticProfile, std::runtime_error, "You may not call this "
3527  "method unless the graph is StaticProfile.");
3528  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3529  isFillComplete () || ! hasColMap (), std::runtime_error, "You may not "
3530  "call this method unless the graph has a column Map.");
3531  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3532  getNodeNumRows () > 0 && k_rowPtrs_.dimension_0 () == 0,
3533  std::runtime_error, "The calling process has getNodeNumRows() = "
3534  << getNodeNumRows () << " > 0 rows, but the row offsets array has not "
3535  "been set.");
3536  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3537  static_cast<size_t> (k_rowPtrs_.dimension_0 ()) != getNodeNumRows () + 1,
3538  std::runtime_error, "The row offsets array has length " <<
3539  k_rowPtrs_.dimension_0 () << " != getNodeNumRows()+1 = " <<
3540  (getNodeNumRows () + 1) << ".");
3541 
3542  // Note: We don't need to do the following things which are normally done in fillComplete:
3543  // allocateIndices, globalAssemble, makeColMap, makeIndicesLocal, sortAllIndices, mergeAllIndices
3544 
3545  // Note: Need to do this so computeGlobalConstants & fillLocalGraph work
3546  //
3547  // The first assignment is always true if the graph has 1-D
3548  // storage (StaticProfile). The second assignment is only true if
3549  // storage is packed.
3550  nodeNumAllocated_ = k_rowPtrs_(getNodeNumRows ());
3551  nodeNumEntries_ = nodeNumAllocated_;
3552 
3553  // Constants from allocateIndices
3554  //
3555  // mfh 08 Aug 2014: numAllocForAllRows_ and k_numAllocPerRow_ go
3556  // away once the graph is allocated. expertStaticFillComplete
3557  // either presumes that the graph is allocated, or "allocates" it.
3558  //
3559  // FIXME (mfh 08 Aug 2014) The goal for the Kokkos refactor
3560  // version of CrsGraph is to allocate in the constructor, not
3561  // lazily on first insert. That will make both
3562  // numAllocForAllRows_ and k_numAllocPerRow_ obsolete.
3563  numAllocForAllRows_ = 0;
3564  k_numAllocPerRow_ = Kokkos::DualView<const size_t*, execution_space> ();
3565  indicesAreAllocated_ = true;
3566 
3567  // Constants from makeIndicesLocal
3568  //
3569  // The graph has a column Map, so its indices had better be local.
3570  indicesAreLocal_ = true;
3571  indicesAreGlobal_ = false;
3572 
3573  // set domain/range map: may clear the import/export objects
3574 #ifdef HAVE_TPETRA_MMM_TIMINGS
3575  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-Maps"))));
3576 #endif
3577  setDomainRangeMaps (domainMap, rangeMap);
3578 
3579  // Presume the user sorted and merged the arrays first
3580  indicesAreSorted_ = true;
3581  noRedundancies_ = true;
3582 
3583  // makeImportExport won't create a new importer/exporter if I set one here first.
3584 #ifdef HAVE_TPETRA_MMM_TIMINGS
3585  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-mIXcheckI"))));
3586 #endif
3587 
3588  importer_ = Teuchos::null;
3589  exporter_ = Teuchos::null;
3590  if (importer != Teuchos::null) {
3591  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3592  ! importer->getSourceMap ()->isSameAs (*getDomainMap ()) ||
3593  ! importer->getTargetMap ()->isSameAs (*getColMap ()),
3594  std::invalid_argument,": importer does not match matrix maps.");
3595  importer_ = importer;
3596 
3597  }
3598 
3599 #ifdef HAVE_TPETRA_MMM_TIMINGS
3600  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-mIXcheckE"))));
3601 #endif
3602 
3603  if (exporter != Teuchos::null) {
3604  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3605  ! exporter->getSourceMap ()->isSameAs (*getRowMap ()) ||
3606  ! exporter->getTargetMap ()->isSameAs (*getRangeMap ()),
3607  std::invalid_argument,": exporter does not match matrix maps.");
3608  exporter_ = exporter;
3609  }
3610 
3611 #ifdef HAVE_TPETRA_MMM_TIMINGS
3612  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-mIXmake"))));
3613 #endif
3614 
3615  makeImportExport ();
3616 
3617  // Compute the constants
3618 #ifdef HAVE_TPETRA_MMM_TIMINGS
3619  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-cGC"))));
3620 #endif
3621  computeGlobalConstants ();
3622 
3623  // Since we have a StaticProfile, fillLocalGraph will do the right thing...
3624 #ifdef HAVE_TPETRA_MMM_TIMINGS
3625  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-fLG"))));
3626 #endif
3627  fillLocalGraph (params);
3628  fillComplete_ = true;
3629 
3630 #ifdef HAVE_TPETRA_MMM_TIMINGS
3631  MM = Teuchos::rcp(new TimeMonitor(*TimeMonitor::getNewTimer(prefix + std::string("ESFC-G-cIS"))));
3632 #endif
3633  checkInternalState ();
3634  }
3635 
3636 
3637  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3638  void
3639  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3640  fillLocalGraph (const Teuchos::RCP<Teuchos::ParameterList>& params)
3641  {
3642  using Kokkos::create_mirror_view;
3643  using Teuchos::ArrayRCP;
3644  using Teuchos::null;
3645  using Teuchos::RCP;
3646  using Teuchos::rcp;
3647  typedef ArrayRCP<size_t>::size_type size_type;
3648  typedef t_numRowEntries_ row_entries_type;
3649  typedef typename local_graph_type::row_map_type row_map_type;
3650  typedef typename row_map_type::non_const_type non_const_row_map_type;
3651  typedef typename local_graph_type::entries_type::non_const_type lclinds_1d_type;
3652 
3653  const size_t lclNumRows = this->getNodeNumRows ();
3654 
3655  // This method's goal is to fill in the two arrays (compressed
3656  // sparse row format) that define the sparse graph's structure.
3657  //
3658  // Use the nonconst version of row_map_type for k_ptrs, because
3659  // the latter is const and we need to modify k_ptrs here.
3660  non_const_row_map_type k_ptrs;
3661  row_map_type k_ptrs_const;
3662  lclinds_1d_type k_inds;
3663 
3664  // The number of entries in each locally owned row. This is a
3665  // DualView. 2-D storage lives on host and is currently not
3666  // thread-safe for parallel kernels even on host, so we have to
3667  // work sequentially with host storage in that case.
3668  row_entries_type k_numRowEnt = k_numRowEntries_;
3669  typename row_entries_type::t_host h_numRowEnt = k_numRowEnt.h_view;
3670 
3671  bool requestOptimizedStorage = true;
3672  if (! params.is_null () && ! params->get ("Optimize Storage", true)) {
3673  requestOptimizedStorage = false;
3674  }
3675  if (getProfileType () == DynamicProfile) {
3676  // Pack 2-D storage (DynamicProfile) into 1-D packed storage.
3677  //
3678  // DynamicProfile means that the graph's column indices are
3679  // currently stored in a 2-D "unpacked" format, in the
3680  // arrays-of-arrays lclInds2D_. We allocate 1-D storage
3681  // (k_inds) and then copy from 2-D storage (lclInds2D_) into 1-D
3682  // storage (k_inds).
3683  TEUCHOS_TEST_FOR_EXCEPTION(
3684  static_cast<size_t> (k_numRowEnt.dimension_0 ()) != lclNumRows,
3685  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph (called from "
3686  "fillComplete or expertStaticFillComplete): For the DynamicProfile "
3687  "branch, k_numRowEnt has the wrong length. "
3688  "k_numRowEnt.dimension_0() = " << k_numRowEnt.dimension_0 ()
3689  << " != getNodeNumRows() = " << lclNumRows << "");
3690 
3691  // Pack the row offsets into k_ptrs, by doing a sum-scan of
3692  // the array of valid entry counts per row (h_numRowEnt).
3693  //
3694  // Total number of entries in the matrix on the calling
3695  // process. We will compute this in the loop below. It's
3696  // cheap to compute and useful as a sanity check.
3697  size_t lclTotalNumEntries = 0;
3698  // This will be a host view of packed row offsets.
3699  typename non_const_row_map_type::HostMirror h_ptrs;
3700  {
3701  // Allocate the packed row offsets array.
3702  k_ptrs = non_const_row_map_type ("Tpetra::CrsGraph::ptr", lclNumRows+1);
3703  k_ptrs_const = k_ptrs;
3704  //
3705  // FIXME hack until we get parallel_scan in kokkos
3706  //
3707  h_ptrs = create_mirror_view (k_ptrs);
3708  h_ptrs(0) = 0;
3709  for (size_type i = 0; i < static_cast<size_type> (lclNumRows); ++i) {
3710  const size_t numEnt = h_numRowEnt(i);
3711  lclTotalNumEntries += numEnt;
3712  h_ptrs(i+1) = h_ptrs(i) + numEnt;
3713  }
3714  Kokkos::deep_copy (k_ptrs, h_ptrs);
3715  }
3716 
3717  TEUCHOS_TEST_FOR_EXCEPTION(
3718  static_cast<size_t> (k_ptrs.dimension_0 ()) != lclNumRows + 1,
3719  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: In DynamicProfile "
3720  "branch, after packing k_ptrs, k_ptrs.dimension_0() = "
3721  << k_ptrs.dimension_0 () << " != (lclNumRows+1) = "
3722  << (lclNumRows+1) << ".");
3723  TEUCHOS_TEST_FOR_EXCEPTION(
3724  static_cast<size_t> (h_ptrs.dimension_0 ()) != lclNumRows + 1,
3725  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: In DynamicProfile "
3726  "branch, after packing h_ptrs, h_ptrs.dimension_0() = "
3727  << h_ptrs.dimension_0 () << " != (lclNumRows+1) = "
3728  << (lclNumRows+1) << ".");
3729  // FIXME (mfh 08 Aug 2014) This assumes UVM.
3730  TEUCHOS_TEST_FOR_EXCEPTION(
3731  k_ptrs(lclNumRows) != lclTotalNumEntries, std::logic_error,
3732  "Tpetra::CrsGraph::fillLocalGraph: In DynamicProfile branch, after "
3733  "packing k_ptrs, k_ptrs(lclNumRows = " << lclNumRows << ") = " <<
3734  k_ptrs(lclNumRows) << " != total number of entries on the calling "
3735  "process = " << lclTotalNumEntries << ".");
3736 
3737  // Allocate the array of packed column indices.
3738  k_inds = lclinds_1d_type ("Tpetra::CrsGraph::ind", lclTotalNumEntries);
3739  // We need a host view of the above, since 2-D storage lives on host.
3740  typename lclinds_1d_type::HostMirror h_inds = create_mirror_view (k_inds);
3741  // Pack the column indices.
3742  for (size_t row = 0; row < lclNumRows; ++row) {
3743  const size_t numEnt = h_numRowEnt(row);
3744  std::copy (lclInds2D_[row].begin (),
3745  lclInds2D_[row].begin () + numEnt,
3746  h_inds.ptr_on_device () + h_ptrs(row));
3747  }
3748  Kokkos::deep_copy (k_inds, h_inds);
3749 
3750  // Sanity check of packed row offsets.
3751  if (k_ptrs.dimension_0 () != 0) {
3752  const size_t numOffsets = static_cast<size_t> (k_ptrs.dimension_0 ());
3753  TEUCHOS_TEST_FOR_EXCEPTION(
3754  static_cast<size_t> (k_ptrs(numOffsets-1)) != k_inds.dimension_0 (),
3755  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: "
3756  "In DynamicProfile branch, after packing, k_ptrs(" << (numOffsets-1)
3757  << ") = " << k_ptrs(numOffsets-1) << " != k_inds.dimension_0() = "
3758  << k_inds.dimension_0 () << ".");
3759  }
3760  }
3761  else if (getProfileType () == StaticProfile) {
3762  // StaticProfile means that the graph's column indices are
3763  // currently stored in a 1-D format, with row offsets in
3764  // k_rowPtrs_ and local column indices in k_lclInds1D_.
3765 
3766  // StaticProfile also means that the graph's array of row
3767  // offsets must already be allocated.
3768  TEUCHOS_TEST_FOR_EXCEPTION(
3769  k_rowPtrs_.dimension_0 () == 0, std::logic_error,
3770  "k_rowPtrs_ has size zero, but shouldn't");
3771  TEUCHOS_TEST_FOR_EXCEPTION(
3772  k_rowPtrs_.dimension_0 () != lclNumRows + 1, std::logic_error,
3773  "Tpetra::CrsGraph::fillLocalGraph: k_rowPtrs_ has size "
3774  << k_rowPtrs_.dimension_0 () << " != (lclNumRows + 1) = "
3775  << (lclNumRows + 1) << ".")
3776  {
3777  const size_t numOffsets = k_rowPtrs_.dimension_0 ();
3778  // FIXME (mfh 08 Aug 2014) This relies on UVM.
3779  TEUCHOS_TEST_FOR_EXCEPTION(
3780  numOffsets != 0 &&
3781  k_lclInds1D_.dimension_0 () != k_rowPtrs_(numOffsets-1),
3782  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: "
3783  "numOffsets = " << numOffsets << " != 0 and "
3784  "k_lclInds1D_.dimension_0() = " << k_lclInds1D_.dimension_0 ()
3785  << " != k_rowPtrs_(" << numOffsets << ") = "
3786  << k_rowPtrs_(numOffsets-1) << ".");
3787  }
3788 
3789  if (nodeNumEntries_ != nodeNumAllocated_) {
3790  // The graph's current 1-D storage is "unpacked." This means
3791  // the row offsets may differ from what the final row offsets
3792  // should be. This could happen, for example, if the user
3793  // specified StaticProfile in the constructor and set an upper
3794  // bound on the number of entries in each row, but didn't fill
3795  // all those entries.
3796  TEUCHOS_TEST_FOR_EXCEPTION(
3797  static_cast<size_t> (k_numRowEnt.dimension_0 ()) != lclNumRows,
3798  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph (called from "
3799  "fillComplete or expertStaticFillComplete): In StaticProfile "
3800  "unpacked branch, k_numRowEnt has the wrong length. "
3801  "k_numRowEnt.dimension_0() = " << k_numRowEnt.dimension_0 ()
3802  << " != getNodeNumRows() = " << lclNumRows << "");
3803 
3804  if (k_rowPtrs_.dimension_0 () != 0) {
3805  const size_t numOffsets =
3806  static_cast<size_t> (k_rowPtrs_.dimension_0 ());
3807  TEUCHOS_TEST_FOR_EXCEPTION(
3808  k_rowPtrs_(numOffsets-1) != static_cast<size_t> (k_lclInds1D_.dimension_0 ()),
3809  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: "
3810  "In StaticProfile branch, before allocating or packing, "
3811  "k_rowPtrs_(" << (numOffsets-1) << ") = "
3812  << k_rowPtrs_(numOffsets-1) << " != k_lclInds1D_.dimension_0() = "
3813  << k_lclInds1D_.dimension_0 () << ".");
3814  }
3815 
3816  // Pack the row offsets into k_ptrs, by doing a sum-scan of
3817  // the array of valid entry counts per row (h_numRowEnt).
3818 
3819  // Total number of entries in the matrix on the calling
3820  // process. We will compute this in the loop below. It's
3821  // cheap to compute and useful as a sanity check.
3822  size_t lclTotalNumEntries = 0;
3823  {
3824  // Allocate the packed row offsets array.
3825  k_ptrs = non_const_row_map_type ("Tpetra::CrsGraph::ptr", lclNumRows + 1);
3826  k_ptrs_const = k_ptrs;
3827  //
3828  // FIXME hack until we get parallel_scan in kokkos
3829  //
3830  // Unlike in the 2-D storage case above, we don't need the
3831  // host view of the packed row offsets array after packing
3832  // the row offsets.
3833  typename non_const_row_map_type::HostMirror h_k_ptrs =
3834  create_mirror_view (k_ptrs);
3835  h_k_ptrs(0) = 0;
3836  for (size_t i = 0; i < lclNumRows; ++i) {
3837  const size_t numEnt = h_numRowEnt(i);
3838  lclTotalNumEntries += numEnt;
3839  h_k_ptrs(i+1) = h_k_ptrs(i) + numEnt;
3840  }
3841  Kokkos::deep_copy (k_ptrs, h_k_ptrs);
3842  }
3843 
3844  TEUCHOS_TEST_FOR_EXCEPTION(
3845  static_cast<size_t> (k_ptrs.dimension_0 ()) != lclNumRows + 1,
3846  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: In "
3847  "StaticProfile unpacked branch, after allocating k_ptrs, "
3848  "k_ptrs.dimension_0() = " << k_ptrs.dimension_0 () << " != "
3849  "lclNumRows+1 = " << (lclNumRows+1) << ".");
3850  // FIXME (mfh 08 Aug 2014) This assumes UVM.
3851  TEUCHOS_TEST_FOR_EXCEPTION(
3852  k_ptrs(lclNumRows) != lclTotalNumEntries, std::logic_error,
3853  "Tpetra::CrsGraph::fillLocalGraph: In StaticProfile unpacked "
3854  "branch, after filling k_ptrs, k_ptrs(lclNumRows=" << lclNumRows
3855  << ") = " << k_ptrs(lclNumRows) << " != total number of entries on "
3856  "the calling process = " << lclTotalNumEntries << ".");
3857 
3858  // Allocate the array of packed column indices.
3859  k_inds = lclinds_1d_type ("Tpetra::CrsGraph::ind", lclTotalNumEntries);
3860 
3861  // k_rowPtrs_ and k_lclInds1D_ are currently unpacked. Pack
3862  // them, using the packed row offsets array k_ptrs that we
3863  // created above.
3864  //
3865  // FIXME (mfh 08 Aug 2014) If "Optimize Storage" is false (in
3866  // CrsMatrix?), we need to keep around the unpacked row
3867  // offsets and column indices.
3868 
3869  // Pack the column indices from unpacked k_lclInds1D_ into
3870  // packed k_inds. We will replace k_lclInds1D_ below.
3871  typedef pack_functor<typename local_graph_type::entries_type::non_const_type, row_map_type> inds_packer_type;
3872  inds_packer_type f (k_inds, k_lclInds1D_, k_ptrs, k_rowPtrs_);
3873  Kokkos::parallel_for (lclNumRows, f);
3874 
3875  TEUCHOS_TEST_FOR_EXCEPTION(
3876  k_ptrs.dimension_0 () == 0, std::logic_error, "Tpetra::CrsGraph::"
3877  "fillLocalGraph: In StaticProfile \"Optimize Storage\" = true branch,"
3878  " after packing, k_ptrs.dimension_0() = 0. This probably means that "
3879  "k_rowPtrs_ was never allocated.");
3880  if (k_ptrs.dimension_0 () != 0) {
3881  const size_t numOffsets = static_cast<size_t> (k_ptrs.dimension_0 ());
3882  TEUCHOS_TEST_FOR_EXCEPTION(
3883  static_cast<size_t> (k_ptrs(numOffsets - 1)) != k_inds.dimension_0 (),
3884  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: "
3885  "In StaticProfile \"Optimize Storage\"=true branch, after packing, "
3886  "k_ptrs(" << (numOffsets-1) << ") = " << k_ptrs(numOffsets-1) <<
3887  " != k_inds.dimension_0() = " << k_inds.dimension_0 () << ".");
3888  }
3889  }
3890  else { // We don't have to pack, so just set the pointers.
3891  k_ptrs_const = k_rowPtrs_;
3892  k_inds = k_lclInds1D_;
3893 
3894  TEUCHOS_TEST_FOR_EXCEPTION(
3895  k_ptrs_const.dimension_0 () == 0, std::logic_error, "Tpetra::CrsGraph::"
3896  "fillLocalGraph: In StaticProfile \"Optimize Storage\" = "
3897  "false branch, k_ptrs_const.dimension_0() = 0. This probably means that "
3898  "k_rowPtrs_ was never allocated.");
3899  if (k_ptrs_const.dimension_0 () != 0) {
3900  const size_t numOffsets = static_cast<size_t> (k_ptrs_const.dimension_0 ());
3901  TEUCHOS_TEST_FOR_EXCEPTION(
3902  static_cast<size_t> (k_ptrs_const(numOffsets - 1)) != k_inds.dimension_0 (),
3903  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: "
3904  "In StaticProfile \"Optimize Storage\" = false branch, "
3905  "k_ptrs_const(" << (numOffsets-1) << ") = " << k_ptrs_const(numOffsets - 1)
3906  << " != k_inds.dimension_0() = " << k_inds.dimension_0 () << ".");
3907  }
3908  }
3909  }
3910 
3911  // Extra sanity checks.
3912  TEUCHOS_TEST_FOR_EXCEPTION(
3913  static_cast<size_t> (k_ptrs_const.dimension_0 ()) != lclNumRows + 1,
3914  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: After packing, "
3915  "k_ptrs_const.dimension_0() = " << k_ptrs_const.dimension_0 ()
3916  << " != lclNumRows+1 = " << (lclNumRows+1) << ".");
3917  if (k_ptrs_const.dimension_0 () != 0) {
3918  const size_t numOffsets = static_cast<size_t> (k_ptrs_const.dimension_0 ());
3919  TEUCHOS_TEST_FOR_EXCEPTION(
3920  static_cast<size_t> (k_ptrs_const(numOffsets - 1)) != k_inds.dimension_0 (),
3921  std::logic_error, "Tpetra::CrsGraph::fillLocalGraph: After packing, "
3922  "k_ptrs_const(" << (numOffsets-1) << ") = " << k_ptrs_const(numOffsets-1)
3923  << " != k_inds.dimension_0() = " << k_inds.dimension_0 () << ".");
3924  }
3925 
3926  if (requestOptimizedStorage) {
3927  // With optimized storage, we don't need to store the 2-D column
3928  // indices array-of-arrays, or the array of row entry counts.
3929 
3930  // Free graph data structures that are only needed for 2-D or
3931  // unpacked 1-D storage.
3932  lclInds2D_ = null;
3933  k_numRowEntries_ = row_entries_type ();
3934 
3935  // Keep the new 1-D packed allocations.
3936  k_rowPtrs_ = k_ptrs_const;
3937  k_lclInds1D_ = k_inds;
3938 
3939  // Storage is packed now, so the number of allocated entries is
3940  // the same as the actual number of entries.
3941  nodeNumAllocated_ = nodeNumEntries_;
3942  // The graph is definitely StaticProfile now, whether or not it
3943  // was before.
3944  pftype_ = StaticProfile;
3945  }
3946 
3947  // FIXME (mfh 28 Aug 2014) "Local Graph" sublist no longer used.
3948 
3949  // Build the local graph.
3950  lclGraph_ = local_graph_type (k_inds, k_ptrs_const);
3951 
3952  // TODO (mfh 13 Mar 2014) getNodeNumDiags(), isUpperTriangular(),
3953  // and isLowerTriangular() depend on computeGlobalConstants(), in
3954  // particular the part where it looks at the local matrix. You
3955  // have to use global indices to determine which entries are
3956  // diagonal, or above or below the diagonal. However, lower or
3957  // upper triangularness is a local property.
3958  }
3959 
3960  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3961  void
3962  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3963  replaceColMap (const Teuchos::RCP<const map_type>& newColMap)
3964  {
3965  // NOTE: This safety check matches the code, but not the documentation of Crsgraph
3966  //
3967  // FIXME (mfh 18 Aug 2014) This will break if the calling process
3968  // has no entries, because in that case, currently it is neither
3969  // locally nor globally indexed. This will change once we get rid
3970  // of lazy allocation (so that the constructor allocates indices
3971  // and therefore commits to local vs. global).
3972  const char tfecfFuncName[] = "replaceColMap: ";
3973  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3974  isLocallyIndexed () || isGloballyIndexed (), std::runtime_error,
3975  "Requires matching maps and non-static graph.");
3976  colMap_ = newColMap;
3977  }
3978 
3979  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
3980  void
3981  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
3982  reindexColumns (const Teuchos::RCP<const map_type>& newColMap,
3983  const Teuchos::RCP<const import_type>& newImport,
3984  const bool sortIndicesInEachRow)
3985  {
3986  using Teuchos::REDUCE_MIN;
3987  using Teuchos::reduceAll;
3988  using Teuchos::RCP;
3989  typedef GlobalOrdinal GO;
3990  typedef LocalOrdinal LO;
3991  typedef typename local_graph_type::entries_type::non_const_type col_inds_type;
3992  const char tfecfFuncName[] = "reindexColumns: ";
3993 
3994  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
3995  isFillComplete (), std::runtime_error, "The graph is fill complete "
3996  "(isFillComplete() returns true). You must call resumeFill() before "
3997  "you may call this method.");
3998 
3999  // mfh 19 Aug 2014: This method does NOT redistribute data; it
4000  // doesn't claim to do the work of an Import or Export. This
4001  // means that for all processes, the calling process MUST own all
4002  // column indices, in both the old column Map (if it exists) and
4003  // the new column Map. We check this via an all-reduce.
4004  //
4005  // Some processes may be globally indexed, others may be locally
4006  // indexed, and others (that have no graph entries) may be
4007  // neither. This method will NOT change the graph's current
4008  // state. If it's locally indexed, it will stay that way, and
4009  // vice versa. It would easy to add an option to convert indices
4010  // from global to local, so as to save a global-to-local
4011  // conversion pass. However, we don't do this here. The intended
4012  // typical use case is that the graph already has a column Map and
4013  // is locally indexed, and this is the case for which we optimize.
4014 
4015  const size_t lclNumRows = getNodeNumRows ();
4016 
4017  // Attempt to convert indices to the new column Map's version of
4018  // local. This will fail if on the calling process, the graph has
4019  // indices that are not on that process in the new column Map.
4020  // After the local conversion attempt, we will do an all-reduce to
4021  // see if any processes failed.
4022 
4023  // If this is false, then either the graph contains a column index
4024  // which is invalid in the CURRENT column Map, or the graph is
4025  // locally indexed but currently has no column Map. In either
4026  // case, there is no way to convert the current local indices into
4027  // global indices, so that we can convert them into the new column
4028  // Map's local indices. It's possible for this to be true on some
4029  // processes but not others, due to replaceColMap.
4030  bool allCurColIndsValid = true;
4031  // On the calling process, are all valid current column indices
4032  // also in the new column Map on the calling process? In other
4033  // words, does local reindexing suffice, or should the user have
4034  // done an Import or Export instead?
4035  bool localSuffices = true;
4036 
4037  // Final arrays for the local indices. We will allocate exactly
4038  // one of these ONLY if the graph is locally indexed on the
4039  // calling process, and ONLY if the graph has one or more entries
4040  // (is not empty) on the calling process. In that case, we
4041  // allocate the first (1-D storage) if the graph has a static
4042  // profile, else we allocate the second (2-D storage).
4043  typename local_graph_type::entries_type::non_const_type newLclInds1D;
4044  Teuchos::ArrayRCP<Teuchos::Array<LO> > newLclInds2D;
4045 
4046  // If indices aren't allocated, that means the calling process
4047  // owns no entries in the graph. Thus, there is nothing to
4048  // convert, and it trivially succeeds locally.
4049  if (indicesAreAllocated ()) {
4050  if (isLocallyIndexed ()) {
4051  if (hasColMap ()) { // locally indexed, and currently has a column Map
4052  const map_type& oldColMap = * (getColMap ());
4053  if (pftype_ == StaticProfile) {
4054  // Allocate storage for the new local indices.
4055  RCP<node_type> node = getRowMap ()->getNode ();
4056  newLclInds1D =
4057  col_inds_type ("Tpetra::CrsGraph::ind", nodeNumAllocated_);
4058  // Attempt to convert the new indices locally.
4059  for (size_t lclRow = 0; lclRow < lclNumRows; ++lclRow) {
4060  const RowInfo rowInfo = getRowInfo (lclRow);
4061  const size_t beg = rowInfo.offset1D;
4062  const size_t end = beg + rowInfo.numEntries;
4063  for (size_t k = beg; k < end; ++k) {
4064  // FIXME (mfh 21 Aug 2014) This assumes UVM. Should
4065  // use a DualView instead.
4066  const LO oldLclCol = k_lclInds1D_(k);
4067  if (oldLclCol == Teuchos::OrdinalTraits<LO>::invalid ()) {
4068  allCurColIndsValid = false;
4069  break; // Stop at the first invalid index
4070  }
4071  const GO gblCol = oldColMap.getGlobalElement (oldLclCol);
4072 
4073  // The above conversion MUST succeed. Otherwise, the
4074  // current local index is invalid, which means that
4075  // the graph was constructed incorrectly.
4076  if (gblCol == Teuchos::OrdinalTraits<GO>::invalid ()) {
4077  allCurColIndsValid = false;
4078  break; // Stop at the first invalid index
4079  }
4080  else {
4081  const LO newLclCol = newColMap->getLocalElement (gblCol);
4082  if (newLclCol == Teuchos::OrdinalTraits<LO>::invalid ()) {
4083  localSuffices = false;
4084  break; // Stop at the first invalid index
4085  }
4086  // FIXME (mfh 21 Aug 2014) This assumes UVM. Should
4087  // use a DualView instead.
4088  newLclInds1D(k) = newLclCol;
4089  }
4090  } // for each entry in the current row
4091  } // for each locally owned row
4092  }
4093  else { // pftype_ == DynamicProfile
4094  // Allocate storage for the new local indices. We only
4095  // allocate the outer array here; we will allocate the
4096  // inner arrays below.
4097  newLclInds2D = Teuchos::arcp<Teuchos::Array<LO> > (lclNumRows);
4098 
4099  // Attempt to convert the new indices locally.
4100  for (size_t lclRow = 0; lclRow < lclNumRows; ++lclRow) {
4101  const RowInfo rowInfo = getRowInfo (lclRow);
4102  newLclInds2D.resize (rowInfo.allocSize);
4103 
4104  Teuchos::ArrayView<const LO> oldLclRowView = getLocalView (rowInfo);
4105  Teuchos::ArrayView<LO> newLclRowView = (newLclInds2D[lclRow]) ();
4106 
4107  for (size_t k = 0; k < rowInfo.numEntries; ++k) {
4108  const LO oldLclCol = oldLclRowView[k];
4109  if (oldLclCol == Teuchos::OrdinalTraits<LO>::invalid ()) {
4110  allCurColIndsValid = false;
4111  break; // Stop at the first invalid index
4112  }
4113  const GO gblCol = oldColMap.getGlobalElement (oldLclCol);
4114 
4115  // The above conversion MUST succeed. Otherwise, the
4116  // local index is invalid and the graph is wrong.
4117  if (gblCol == Teuchos::OrdinalTraits<GO>::invalid ()) {
4118  allCurColIndsValid = false;
4119  break; // Stop at the first invalid index
4120  }
4121  else {
4122  const LO newLclCol = newColMap->getLocalElement (gblCol);
4123  if (newLclCol == Teuchos::OrdinalTraits<LO>::invalid ()) {
4124  localSuffices = false;
4125  break; // Stop at the first invalid index.
4126  }
4127  newLclRowView[k] = newLclCol;
4128  }
4129  } // for each entry in the current row
4130  } // for each locally owned row
4131  } // pftype_
4132  }
4133  else { // locally indexed, but no column Map
4134  // This case is only possible if replaceColMap() was called
4135  // with a null argument on the calling process. It's
4136  // possible, but it means that this method can't possibly
4137  // succeed, since we have no way of knowing how to convert
4138  // the current local indices to global indices.
4139  allCurColIndsValid = false;
4140  }
4141  }
4142  else { // globally indexed
4143  // If the graph is globally indexed, we don't need to save
4144  // local indices, but we _do_ need to know whether the current
4145  // global indices are valid in the new column Map. We may
4146  // need to do a getRemoteIndexList call to find this out.
4147  //
4148  // In this case, it doesn't matter whether the graph currently
4149  // has a column Map. We don't need the old column Map to
4150  // convert from global indices to the _new_ column Map's local
4151  // indices. Furthermore, we can use the same code, whether
4152  // the graph is static or dynamic profile.
4153 
4154  // Test whether the current global indices are in the new
4155  // column Map on the calling process.
4156  for (size_t lclRow = 0; lclRow < lclNumRows; ++lclRow) {
4157  const RowInfo rowInfo = getRowInfo (lclRow);
4158  Teuchos::ArrayView<const GO> oldGblRowView = getGlobalView (rowInfo);
4159  for (size_t k = 0; k < rowInfo.numEntries; ++k) {
4160  const GO gblCol = oldGblRowView[k];
4161  if (! newColMap->isNodeGlobalElement (gblCol)) {
4162  localSuffices = false;
4163  break; // Stop at the first invalid index
4164  }
4165  } // for each entry in the current row
4166  } // for each locally owned row
4167  } // locally or globally indexed
4168  } // whether indices are allocated
4169 
4170  // Do an all-reduce to check both possible error conditions.
4171  int lclSuccess[2];
4172  lclSuccess[0] = allCurColIndsValid ? 1 : 0;
4173  lclSuccess[1] = localSuffices ? 1 : 0;
4174  int gblSuccess[2];
4175  gblSuccess[0] = 0;
4176  gblSuccess[1] = 0;
4177  RCP<const Teuchos::Comm<int> > comm =
4178  getRowMap ().is_null () ? Teuchos::null : getRowMap ()->getComm ();
4179  if (! comm.is_null ()) {
4180  reduceAll<int, int> (*comm, REDUCE_MIN, 2, lclSuccess, gblSuccess);
4181  }
4182 
4183  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4184  gblSuccess[0] == 0, std::runtime_error, "It is not possible to continue."
4185  " The most likely reason is that the graph is locally indexed, but the "
4186  "column Map is missing (null) on some processes, due to a previous call "
4187  "to replaceColMap().");
4188 
4189  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4190  gblSuccess[1] == 0, std::runtime_error, "On some process, the graph "
4191  "contains column indices that are in the old column Map, but not in the "
4192  "new column Map (on that process). This method does NOT redistribute "
4193  "data; it does not claim to do the work of an Import or Export operation."
4194  " This means that for all processess, the calling process MUST own all "
4195  "column indices, in both the old column Map and the new column Map. In "
4196  "this case, you will need to do an Import or Export operation to "
4197  "redistribute data.");
4198 
4199  // Commit the results.
4200  if (isLocallyIndexed ()) {
4201  if (pftype_ == StaticProfile) {
4202  k_lclInds1D_ = newLclInds1D;
4203  } else { // dynamic profile
4204  lclInds2D_ = newLclInds2D;
4205  }
4206  // We've reindexed, so we don't know if the indices are sorted.
4207  //
4208  // FIXME (mfh 17 Sep 2014) It could make sense to check this,
4209  // since we're already going through all the indices above. We
4210  // could also sort each row in place; that way, we would only
4211  // have to make one pass over the rows.
4212  indicesAreSorted_ = false;
4213  if (sortIndicesInEachRow) {
4214  // NOTE (mfh 17 Sep 2014) The graph must be locally indexed in
4215  // order to call this method.
4216  //
4217  // FIXME (mfh 17 Sep 2014) This violates the strong exception
4218  // guarantee. It would be better to sort the new index arrays
4219  // before committing them.
4220  sortAllIndices ();
4221  }
4222  }
4223  colMap_ = newColMap;
4224 
4225  if (newImport.is_null ()) {
4226  // FIXME (mfh 19 Aug 2014) Should use the above all-reduce to
4227  // check whether the input Import is null on any process.
4228  //
4229  // If the domain Map hasn't been set yet, we can't compute a new
4230  // Import object. Leave it what it is; it should be null, but
4231  // it doesn't matter. If the domain Map _has_ been set, then
4232  // compute a new Import object if necessary.
4233  if (! domainMap_.is_null ()) {
4234  if (! domainMap_->isSameAs (* newColMap)) {
4235  importer_ = Teuchos::rcp (new import_type (domainMap_, newColMap));
4236  } else {
4237  importer_ = Teuchos::null; // don't need an Import
4238  }
4239  }
4240  } else {
4241  // The caller gave us an Import object. Assume that it's valid.
4242  importer_ = newImport;
4243  }
4244  }
4245 
4246 
4247  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4248  void
4249  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4250  replaceDomainMapAndImporter (const Teuchos::RCP<const map_type>& newDomainMap,
4251  const Teuchos::RCP<const import_type>& newImporter)
4252  {
4253  const char prefix[] = "Tpetra::CrsGraph::replaceDomainMapAndImporter: ";
4254  TEUCHOS_TEST_FOR_EXCEPTION(
4255  colMap_.is_null (), std::invalid_argument, prefix << "You may not call "
4256  "this method unless the graph already has a column Map.");
4257  TEUCHOS_TEST_FOR_EXCEPTION(
4258  newDomainMap.is_null (), std::invalid_argument,
4259  prefix << "The new domain Map must be nonnull.");
4260 
4261 #ifdef HAVE_TPETRA_DEBUG
4262  if (newImporter.is_null ()) {
4263  // It's not a good idea to put expensive operations in a macro
4264  // clause, even if they are side effect - free, because macros
4265  // don't promise that they won't evaluate their arguments more
4266  // than once. It's polite for them to do so, but not required.
4267  const bool colSameAsDom = colMap_->isSameAs (*newDomainMap);
4268  TEUCHOS_TEST_FOR_EXCEPTION(
4269  colSameAsDom, std::invalid_argument, "If the new Import is null, "
4270  "then the new domain Map must be the same as the current column Map.");
4271  }
4272  else {
4273  const bool colSameAsTgt =
4274  colMap_->isSameAs (* (newImporter->getTargetMap ()));
4275  const bool newDomSameAsSrc =
4276  newDomainMap->isSameAs (* (newImporter->getSourceMap ()));
4277  TEUCHOS_TEST_FOR_EXCEPTION(
4278  ! colSameAsTgt || ! newDomSameAsSrc, std::invalid_argument, "If the "
4279  "new Import is nonnull, then the current column Map must be the same "
4280  "as the new Import's target Map, and the new domain Map must be the "
4281  "same as the new Import's source Map.");
4282  }
4283 #endif // HAVE_TPETRA_DEBUG
4284 
4285  domainMap_ = newDomainMap;
4286  importer_ = Teuchos::rcp_const_cast<import_type> (newImporter);
4287  }
4288 
4289  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4290  typename CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::local_graph_type
4291  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4292  getLocalGraph () const
4293  {
4294  return lclGraph_;
4295  }
4296 
4297  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4298  typename CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::local_graph_type
4299  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4300  getLocalGraph_Kokkos () const
4301  {
4302  return lclGraph_;
4303  }
4304 
4305  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4306  void
4307  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4308  computeGlobalConstants ()
4309  {
4310  using Teuchos::as;
4311  using Teuchos::outArg;
4312  using Teuchos::reduceAll;
4313  typedef LocalOrdinal LO;
4314  typedef GlobalOrdinal GO;
4315  typedef global_size_t GST;
4316 
4317 #ifdef HAVE_TPETRA_DEBUG
4318  TEUCHOS_TEST_FOR_EXCEPTION(! hasColMap(), std::logic_error, "Tpetra::"
4319  "CrsGraph::computeGlobalConstants: At this point, the graph should have "
4320  "a column Map, but it does not. Please report this bug to the Tpetra "
4321  "developers.");
4322 #endif // HAVE_TPETRA_DEBUG
4323 
4324  // If necessary, (re)compute the local constants: nodeNumDiags_,
4325  // lowerTriangular_, upperTriangular_, and nodeMaxNumRowEntries_.
4326  if (! haveLocalConstants_) {
4327  // We have actually already computed nodeNumEntries_.
4328  // nodeNumEntries_ gets updated by methods that insert or remove
4329  // indices (including setAllIndices and
4330  // expertStaticFillComplete). Before fillComplete, its count
4331  // may include duplicate column indices in the same row.
4332  // However, mergeRowIndices and mergeRowIndicesAndValues both
4333  // subtract off merged indices in each row from the total count.
4334  // Thus, nodeNumEntries_ _should_ be accurate at this point,
4335  // meaning that we don't have to re-count it here.
4336 
4337  // Reset local properties
4338  upperTriangular_ = true;
4339  lowerTriangular_ = true;
4340  nodeMaxNumRowEntries_ = 0;
4341  nodeNumDiags_ = 0;
4342 
4343  // At this point, we know that we have both a row Map and a column Map.
4344  const map_type& rowMap = *rowMap_;
4345  const map_type& colMap = *colMap_;
4346 
4347  // Go through all the entries of the graph. Count the number of
4348  // diagonal elements we encounter, and figure out whether the
4349  // graph is lower or upper triangular. Diagonal elements are
4350  // determined using global indices, with respect to the whole
4351  // graph. However, lower or upper triangularity is a local
4352  // property, and is determined using local indices.
4353  //
4354  // At this point, indices have already been sorted in each row.
4355  // That makes finding out whether the graph is lower / upper
4356  // triangular easier.
4357  if (indicesAreAllocated () && nodeNumAllocated_ > 0) {
4358  const size_t numLocalRows = getNodeNumRows ();
4359  for (size_t localRow = 0; localRow < numLocalRows; ++localRow) {
4360  const GO globalRow = rowMap.getGlobalElement (localRow);
4361  // Find the local (column) index for the diagonal entry.
4362  // This process might not necessarily own _any_ entries in
4363  // the current row. If it doesn't, skip this row. It won't
4364  // affect any of the attributes (nodeNumDiagons_,
4365  // upperTriangular_, lowerTriangular_, or
4366  // nodeMaxNumRowEntries_) which this loop sets.
4367  const LO rlcid = colMap.getLocalElement (globalRow);
4368  // This process owns one or more entries in the current row.
4369  RowInfo rowInfo = getRowInfo (localRow);
4370  ArrayView<const LO> rview = getLocalView (rowInfo);
4371  typename ArrayView<const LO>::iterator beg, end, cur;
4372  beg = rview.begin();
4373  end = beg + rowInfo.numEntries;
4374  if (beg != end) {
4375  for (cur = beg; cur != end; ++cur) {
4376  // is this the diagonal?
4377  if (rlcid == *cur) ++nodeNumDiags_;
4378  }
4379  // Local column indices are sorted in each row. That means
4380  // the smallest column index in this row (on this process)
4381  // is *beg, and the largest column index in this row (on
4382  // this process) is *(end - 1). We know that end - 1 is
4383  // valid because beg != end.
4384  const size_t smallestCol = static_cast<size_t> (*beg);
4385  const size_t largestCol = static_cast<size_t> (*(end - 1));
4386 
4387  if (smallestCol < localRow) {
4388  upperTriangular_ = false;
4389  }
4390  if (localRow < largestCol) {
4391  lowerTriangular_ = false;
4392  }
4393  }
4394  // Update the max number of entries over all rows.
4395  nodeMaxNumRowEntries_ = std::max (nodeMaxNumRowEntries_, rowInfo.numEntries);
4396  }
4397  }
4398  haveLocalConstants_ = true;
4399  } // if my process doesn't have local constants
4400 
4401  // Compute global constants from local constants. Processes that
4402  // already have local constants still participate in the
4403  // all-reduces, using their previously computed values.
4404  if (haveGlobalConstants_ == false) {
4405  // Promote all the nodeNum* and nodeMaxNum* quantities from
4406  // size_t to global_size_t, when doing the all-reduces for
4407  // globalNum* / globalMaxNum* results.
4408  //
4409  // FIXME (mfh 07 May 2013) Unfortunately, we either have to do
4410  // this in two all-reduces (one for the sum and the other for
4411  // the max), or use a custom MPI_Op that combines the sum and
4412  // the max. The latter might even be slower than two
4413  // all-reduces on modern network hardware. It would also be a
4414  // good idea to use nonblocking all-reduces (MPI 3), so that we
4415  // don't have to wait around for the first one to finish before
4416  // starting the second one.
4417  GST lcl[2], gbl[2];
4418  lcl[0] = static_cast<GST> (nodeNumEntries_);
4419  lcl[1] = static_cast<GST> (nodeNumDiags_);
4420  reduceAll<int,GST> (*getComm (), Teuchos::REDUCE_SUM,
4421  2, lcl, gbl);
4422  globalNumEntries_ = gbl[0];
4423  globalNumDiags_ = gbl[1];
4424  reduceAll<int,GST> (*getComm (), Teuchos::REDUCE_MAX,
4425  static_cast<GST> (nodeMaxNumRowEntries_),
4426  outArg (globalMaxNumRowEntries_));
4427  haveGlobalConstants_ = true;
4428  }
4429  }
4430 
4431 
4432  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4433  void
4434  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4435  makeIndicesLocal ()
4436  {
4437  using Teuchos::arcp;
4438  using Teuchos::Array;
4439  typedef LocalOrdinal LO;
4440  typedef GlobalOrdinal GO;
4441  typedef typename local_graph_type::entries_type::non_const_type
4442  lcl_col_inds_type;
4443  typedef Kokkos::View<GlobalOrdinal*,
4444  typename lcl_col_inds_type::array_layout,
4445  device_type> gbl_col_inds_type;
4446  const char tfecfFuncName[] = "makeIndicesLocal";
4447 
4448  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4449  ! this->hasColMap (), std::logic_error, ": The graph does not have a "
4450  "column Map yet. This method should never be called in that case. "
4451  "Please report this bug to the Tpetra developers.");
4452  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4453  this->getColMap ().is_null (), std::logic_error, ": The graph claims "
4454  "that it has a column Map, because hasColMap() returns true. However, "
4455  "the result of getColMap() is null. This should never happen. Please "
4456  "report this bug to the Tpetra developers.");
4457 
4458  const size_t lclNumRows = this->getNodeNumRows ();
4459  const map_type& colMap = * (this->getColMap ());
4460 
4461  if (isGloballyIndexed () && lclNumRows > 0) {
4462  typename t_numRowEntries_::t_host h_numRowEnt = k_numRowEntries_.h_view;
4463 
4464  // allocate data for local indices
4465  if (getProfileType () == StaticProfile) {
4466  // If GO and LO are the same size, we can reuse the existing
4467  // array of 1-D index storage to convert column indices from
4468  // GO to LO. Otherwise, we'll just allocate a new buffer.
4469  if (nodeNumAllocated_ && Kokkos::Impl::is_same<LO,GO>::value) {
4470  k_lclInds1D_ = Kokkos::Impl::if_c<Kokkos::Impl::is_same<LO,GO>::value,
4471  t_GlobalOrdinal_1D,
4472  lcl_col_inds_type >::select (k_gblInds1D_, k_lclInds1D_);
4473  }
4474  else {
4475  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4476  k_rowPtrs_.dimension_0 () == 0, std::logic_error, ": This should "
4477  "never happen at this point. Please report this bug to the Tpetra "
4478  "developers.");
4479  const size_t numEnt = k_rowPtrs_[lclNumRows];
4480 
4481  k_lclInds1D_ = lcl_col_inds_type ("Tpetra::CrsGraph::lclind", numEnt);
4482  }
4483 
4484  for (size_t r = 0; r < lclNumRows; ++r) {
4485  const size_t offset = k_rowPtrs_(r);
4486  const size_t numentry = h_numRowEnt(r);
4487  for (size_t j = 0; j < numentry; ++j) {
4488  const GO gid = k_gblInds1D_(offset + j);
4489  const LO lid = colMap.getLocalElement (gid);
4490  k_lclInds1D_(offset + j) = lid;
4491 #ifdef HAVE_TPETRA_DEBUG
4492  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4493  k_lclInds1D_(offset + j) == Teuchos::OrdinalTraits<LO>::invalid(),
4494  std::logic_error,
4495  ": In local row r=" << r << ", global column " << gid << " is "
4496  "not in the column Map. This should never happen. Please "
4497  "report this bug to the Tpetra developers.");
4498 #endif // HAVE_TPETRA_DEBUG
4499  }
4500  }
4501  // We've converted column indices from global to local, so we
4502  // can deallocate the global column indices (which we know are
4503  // in 1-D storage, because the graph has static profile).
4504  k_gblInds1D_ = gbl_col_inds_type ();
4505  }
4506  else { // the graph has dynamic profile (2-D index storage)
4507  lclInds2D_ = arcp<Array<LO> > (lclNumRows);
4508  for (size_t r = 0; r < lclNumRows; ++r) {
4509  if (! gblInds2D_[r].empty ()) {
4510  const GO* const ginds = gblInds2D_[r].getRawPtr ();
4511  const size_t rna = gblInds2D_[r].size ();
4512  const size_t numentry = h_numRowEnt(r);
4513  lclInds2D_[r].resize (rna);
4514  LO* const linds = lclInds2D_[r].getRawPtr ();
4515  for (size_t j = 0; j < numentry; ++j) {
4516  const GO gid = ginds[j];
4517  const LO lid = colMap.getLocalElement (gid);
4518  linds[j] = lid;
4519 #ifdef HAVE_TPETRA_DEBUG
4520  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4521  linds[j] == Teuchos::OrdinalTraits<LO>::invalid(),
4522  std::logic_error,
4523  ": Global column ginds[j=" << j << "]=" << ginds[j]
4524  << " of local row r=" << r << " is not in the column Map. "
4525  "This should never happen. Please report this bug to the "
4526  "Tpetra developers.");
4527 #endif // HAVE_TPETRA_DEBUG
4528  }
4529  }
4530  }
4531  gblInds2D_ = Teuchos::null;
4532  }
4533  } // globallyIndexed() && lclNumRows > 0
4534 
4535  lclGraph_ = local_graph_type (k_lclInds1D_, k_rowPtrs_);
4536  indicesAreLocal_ = true;
4537  indicesAreGlobal_ = false;
4538  checkInternalState ();
4539  }
4540 
4541 
4542  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4543  void
4544  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4545  sortAllIndices ()
4546  {
4547  // this should be called only after makeIndicesLocal()
4548  TEUCHOS_TEST_FOR_EXCEPT( isGloballyIndexed () );
4549  if (isSorted () == false) {
4550  // FIXME (mfh 06 Mar 2014) This would be a good place for a
4551  // thread-parallel kernel.
4552  for (size_t row = 0; row < getNodeNumRows (); ++row) {
4553  RowInfo rowInfo = getRowInfo (row);
4554  sortRowIndices (rowInfo);
4555  }
4556  }
4557  indicesAreSorted_ = true; // we just sorted every row
4558  }
4559 
4560 
4561  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4562  void
4563  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4564  makeColMap ()
4565  {
4566  using Teuchos::Array;
4567  using Teuchos::ArrayView;
4568  using Teuchos::rcp;
4569  using Teuchos::REDUCE_MAX;
4570  using Teuchos::reduceAll;
4571  using std::endl;
4572  typedef LocalOrdinal LO;
4573  typedef GlobalOrdinal GO;
4574  const char tfecfFuncName[] = "makeColMap";
4575 
4576  if (hasColMap ()) { // The graph already has a column Map.
4577  // FIXME (mfh 26 Feb 2013): This currently prevents structure
4578  // changes that affect the column Map.
4579  return;
4580  }
4581 
4582  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4583  isLocallyIndexed (), std::runtime_error,
4584  ": The graph is locally indexed. Calling makeColMap() to make the "
4585  "column Map requires that the graph be globally indexed.");
4586 
4587  // After the calling process is done going through all of the rows
4588  // it owns, myColumns will contain the list of indices owned by
4589  // this process in the column Map.
4590  Array<GO> myColumns;
4591 
4592  // If we reach this point, we don't have a column Map yet, so the
4593  // graph can't be locally indexed. Thus, isGloballyIndexed() ==
4594  // false means that the graph is empty on this process, so
4595  // myColumns will be left empty.
4596  if (isGloballyIndexed ()) {
4597  // Go through all the rows, finding the populated column indices.
4598  //
4599  // Our final list of indices for the column Map constructor will
4600  // have the following properties (all of which are with respect
4601  // to the calling process):
4602  //
4603  // 1. Indices in the domain Map go first.
4604  // 2. Indices not in the domain Map follow, ordered first
4605  // contiguously by their owning process rank (in the domain
4606  // Map), then in increasing order within that.
4607  // 3. No duplicate indices.
4608  //
4609  // This imitates the ordering used by Aztec(OO) and Epetra.
4610  // Storing indices owned by the same process (in the domain Map)
4611  // contiguously permits the use of contiguous send and receive
4612  // buffers.
4613  //
4614  // We begin by partitioning the column indices into "local" GIDs
4615  // (owned by the domain Map) and "remote" GIDs (not owned by the
4616  // domain Map). We use the same order for local GIDs as the
4617  // domain Map, so we track them in place in their array. We use
4618  // an std::set (RemoteGIDSet) to keep track of remote GIDs, so
4619  // that we don't have to merge duplicates later.
4620  const LO LINV = Teuchos::OrdinalTraits<LO>::invalid ();
4621  size_t numLocalColGIDs = 0, numRemoteColGIDs = 0;
4622 
4623  // GIDisLocal[lid] == 0 if and only if local index lid in the
4624  // domain Map is remote (not local).
4625  Array<char> GIDisLocal (domainMap_->getNodeNumElements (), 0);
4626  std::set<GO> RemoteGIDSet;
4627  // This preserves the not-sorted Epetra order of GIDs.
4628  std::vector<GO> RemoteGIDUnorderedVector;
4629  const size_t myNumRows = getNodeNumRows ();
4630  for (size_t r = 0; r < myNumRows; ++r) {
4631  RowInfo rowinfo = getRowInfo (r);
4632  if (rowinfo.numEntries > 0) {
4633  // NOTE (mfh 02 Sep 2014) getGlobalView() returns a view of
4634  // all the space in the row, not just the occupied entries.
4635  // (This matters for the case of unpacked 1-D storage. We
4636  // might not have packed it yet.) That's why we need to
4637  // take a subview.
4638  ArrayView<const GO> rowGids = getGlobalView (rowinfo);
4639  rowGids = rowGids (0, rowinfo.numEntries);
4640 
4641  for (size_t k = 0; k < rowinfo.numEntries; ++k) {
4642  const GO gid = rowGids[k];
4643  const LO lid = domainMap_->getLocalElement (gid);
4644  if (lid != LINV) {
4645  const char alreadyFound = GIDisLocal[lid];
4646  if (alreadyFound == 0) {
4647  GIDisLocal[lid] = static_cast<char> (1);
4648  ++numLocalColGIDs;
4649  }
4650  }
4651  else {
4652  const bool notAlreadyFound = RemoteGIDSet.insert (gid).second;
4653  if (notAlreadyFound) { // gid did not exist in the set before
4654  if (! sortGhostsAssociatedWithEachProcessor_) {
4655  // The user doesn't want to sort remote GIDs (for
4656  // each remote process); they want us to keep remote
4657  // GIDs in their original order. We do this by
4658  // stuffing each remote GID into an array as we
4659  // encounter it for the first time. The std::set
4660  // helpfully tracks first encounters.
4661  RemoteGIDUnorderedVector.push_back (gid);
4662  }
4663  ++numRemoteColGIDs;
4664  }
4665  }
4666  } // for each entry k in row r
4667  } // if row r contains > 0 entries
4668  } // for each locally owned row r
4669 
4670  // Possible short-circuit for serial scenario:
4671  //
4672  // If all domain GIDs are present as column indices, then set
4673  // ColMap=DomainMap. By construction, LocalGIDs is a subset of
4674  // DomainGIDs.
4675  //
4676  // If we have
4677  // * Number of remote GIDs is 0, so that ColGIDs == LocalGIDs,
4678  // and
4679  // * Number of local GIDs is number of domain GIDs
4680  // then
4681  // * LocalGIDs \subset DomainGIDs && size(LocalGIDs) ==
4682  // size(DomainGIDs) => DomainGIDs == LocalGIDs == ColGIDs
4683  // on the calling process.
4684  //
4685  // We will concern ourselves only with the special case of a
4686  // serial DomainMap, obviating the need for communication.
4687  //
4688  // If
4689  // * DomainMap has a serial communicator
4690  // then we can set the column Map as the domain Map
4691  // return. Benefit: this graph won't need an Import object
4692  // later.
4693  //
4694  // Note, for a serial domain map, there can be no RemoteGIDs,
4695  // because there are no remote processes. Likely explanations
4696  // for this are:
4697  // * user submitted erroneous column indices
4698  // * user submitted erroneous domain Map
4699  if (domainMap_->getComm ()->getSize () == 1) {
4700  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4701  numRemoteColGIDs != 0, std::runtime_error,
4702  ": " << numRemoteColGIDs << " column "
4703  << (numRemoteColGIDs != 1 ? "indices are" : "index is")
4704  << " not in the domain Map." << endl
4705  << "Either these indices are invalid or the domain Map is invalid."
4706  << endl << "Remember that nonsquare matrices, or matrices where the "
4707  "row and range Maps are different, require calling the version of "
4708  "fillComplete that takes the domain and range Maps as input.");
4709  if (numLocalColGIDs == domainMap_->getNodeNumElements()) {
4710  colMap_ = domainMap_;
4711  checkInternalState ();
4712  return;
4713  }
4714  }
4715 
4716  // Populate myColumns with a list of all column GIDs. Put
4717  // locally owned (in the domain Map) GIDs at the front: they
4718  // correspond to "same" and "permuted" entries between the
4719  // column Map and the domain Map. Put remote GIDs at the back.
4720  myColumns.resize (numLocalColGIDs + numRemoteColGIDs);
4721  // get pointers into myColumns for each part
4722  ArrayView<GO> LocalColGIDs = myColumns (0, numLocalColGIDs);
4723  ArrayView<GO> RemoteColGIDs = myColumns (numLocalColGIDs, numRemoteColGIDs);
4724 
4725  // Copy the remote GIDs into myColumns
4726  if (sortGhostsAssociatedWithEachProcessor_) {
4727  // The std::set puts GIDs in increasing order.
4728  std::copy (RemoteGIDSet.begin(), RemoteGIDSet.end(),
4729  RemoteColGIDs.begin());
4730  } else {
4731  // Respect the originally encountered order.
4732  std::copy (RemoteGIDUnorderedVector.begin(),
4733  RemoteGIDUnorderedVector.end(), RemoteColGIDs.begin());
4734  }
4735 
4736  // Make a list of process ranks corresponding to the remote GIDs.
4737  Array<int> RemoteImageIDs (numRemoteColGIDs);
4738  // Look up the remote process' ranks in the domain Map.
4739  {
4740  const LookupStatus stat =
4741  domainMap_->getRemoteIndexList (RemoteColGIDs, RemoteImageIDs ());
4742 #ifdef HAVE_TPETRA_DEBUG
4743  // If any process returns IDNotPresent, then at least one of
4744  // the remote indices was not present in the domain Map. This
4745  // means that the Import object cannot be constructed, because
4746  // of incongruity between the column Map and domain Map.
4747  // This has two likely causes:
4748  // - The user has made a mistake in the column indices
4749  // - The user has made a mistake with respect to the domain Map
4750  const int missingID_lcl = (stat == IDNotPresent ? 1 : 0);
4751  int missingID_gbl = 0;
4752  reduceAll<int, int> (*getComm (), REDUCE_MAX, missingID_lcl,
4753  outArg (missingID_gbl));
4754  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4755  missingID_gbl == 1, std::runtime_error,
4756  ": Some column indices are not in the domain Map." << endl
4757  << "Either these column indices are invalid or the domain Map is "
4758  "invalid." << endl << "Likely cause: For a nonsquare matrix, you "
4759  "must give the domain and range Maps as input to fillComplete.");
4760 #else
4761  (void) stat; // forestall compiler warning for unused variable
4762 #endif // HAVE_TPETRA_DEBUG
4763  }
4764  // Sort incoming remote column indices by their owning process
4765  // rank, so that all columns coming from a given remote process
4766  // are contiguous. This means the Import's Distributor doesn't
4767  // need to reorder data.
4768  //
4769  // NOTE (mfh 02 Sep 2014) This needs to be a stable sort, so
4770  // that it respects either of the possible orderings of GIDs
4771  // (sorted, or original order) specified above.
4772  sort2 (RemoteImageIDs.begin(), RemoteImageIDs.end(), RemoteColGIDs.begin());
4773 
4774  // Copy the local GIDs into myColumns. Two cases:
4775  // 1. If the number of Local column GIDs is the same as the
4776  // number of Local domain GIDs, we can simply read the domain
4777  // GIDs into the front part of ColIndices (see logic above
4778  // from the serial short circuit case)
4779  // 2. We step through the GIDs of the DomainMap, checking to see
4780  // if each domain GID is a column GID. We want to do this to
4781  // maintain a consistent ordering of GIDs between the columns
4782  // and the domain.
4783 
4784  const size_t numDomainElts = domainMap_->getNodeNumElements ();
4785  if (numLocalColGIDs == numDomainElts) {
4786  // If the number of locally owned GIDs are the same as the
4787  // number of local domain Map elements, then the local domain
4788  // Map elements are the same as the locally owned GIDs.
4789  if (domainMap_->isContiguous ()) {
4790  // NOTE (mfh 03 Mar 2013, 02 Sep 2014) In the common case
4791  // that the domain Map is contiguous, it's more efficient to
4792  // avoid calling getNodeElementList(), since that
4793  // permanently constructs and caches the GID list in the
4794  // contiguous Map.
4795  GO curColMapGid = domainMap_->getMinGlobalIndex ();
4796  for (size_t k = 0; k < numLocalColGIDs; ++k, ++curColMapGid) {
4797  LocalColGIDs[k] = curColMapGid;
4798  }
4799  }
4800  else {
4801  ArrayView<const GO> domainElts = domainMap_->getNodeElementList ();
4802  std::copy (domainElts.begin(), domainElts.end(), LocalColGIDs.begin());
4803  }
4804  }
4805  else {
4806  // Count the number of locally owned GIDs, both to keep track
4807  // of the current array index, and as a sanity check.
4808  size_t numLocalCount = 0;
4809  if (domainMap_->isContiguous ()) {
4810  // NOTE (mfh 03 Mar 2013, 02 Sep 2014) In the common case
4811  // that the domain Map is contiguous, it's more efficient to
4812  // avoid calling getNodeElementList(), since that
4813  // permanently constructs and caches the GID list in the
4814  // contiguous Map.
4815  GO curColMapGid = domainMap_->getMinGlobalIndex ();
4816  for (size_t i = 0; i < numDomainElts; ++i, ++curColMapGid) {
4817  if (GIDisLocal[i]) {
4818  LocalColGIDs[numLocalCount++] = curColMapGid;
4819  }
4820  }
4821  }
4822  else {
4823  ArrayView<const GO> domainElts = domainMap_->getNodeElementList ();
4824  for (size_t i = 0; i < numDomainElts; ++i) {
4825  if (GIDisLocal[i]) {
4826  LocalColGIDs[numLocalCount++] = domainElts[i];
4827  }
4828  }
4829  }
4830  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
4831  numLocalCount != numLocalColGIDs, std::logic_error,
4832  ": numLocalCount = " << numLocalCount << " != numLocalColGIDs = "
4833  << numLocalColGIDs << ". This should never happen. Please report "
4834  "this bug to the Tpetra developers.");
4835  }
4836 
4837  // FIXME (mfh 03 Apr 2013) Now would be a good time to use the
4838  // information we collected above to construct the Import. In
4839  // particular, building an Import requires:
4840  //
4841  // 1. numSameIDs (length of initial contiguous sequence of GIDs
4842  // on this process that are the same in both Maps; this
4843  // equals the number of domain Map elements on this process)
4844  //
4845  // 2. permuteToLIDs and permuteFromLIDs (both empty in this
4846  // case, since there's no permutation going on; the column
4847  // Map starts with the domain Map's GIDs, and immediately
4848  // after them come the remote GIDs)
4849  //
4850  // 3. remoteGIDs (exactly those GIDs that we found out above
4851  // were not in the domain Map) and remoteLIDs (which we could
4852  // have gotten above by using the three-argument version of
4853  // getRemoteIndexList() that computes local indices as well
4854  // as process ranks, instead of the two-argument version that
4855  // was used above)
4856  //
4857  // 4. remotePIDs (which we have from the getRemoteIndexList()
4858  // call above)
4859  //
4860  // 5. Sorting remotePIDs, and applying that permutation to
4861  // remoteGIDs and remoteLIDs (by calling sort3 above instead
4862  // of sort2)
4863  //
4864  // 6. Everything after the sort3 call in Import::setupExport():
4865  // a. Create the Distributor via createFromRecvs(), which
4866  // computes exportGIDs and exportPIDs
4867  // b. Compute exportLIDs from exportGIDs (by asking the
4868  // source Map, in this case the domain Map, to convert
4869  // global to local)
4870  //
4871  // Steps 1-5 come for free, since we must do that work anyway in
4872  // order to compute the column Map. In particular, Step 3 is
4873  // even more expensive than Step 6a, since it involves both
4874  // creating and using a new Distributor object.
4875 
4876  } // if the graph is globally indexed
4877 
4878  const global_size_t gstInv =
4879  Teuchos::OrdinalTraits<global_size_t>::invalid ();
4880  // FIXME (mfh 05 Mar 2014) Doesn't the index base of a Map have to
4881  // be the same as the Map's min GID? If the first column is empty
4882  // (contains no entries), then the column Map's min GID won't
4883  // necessarily be the same as the domain Map's index base.
4884  const GO indexBase = domainMap_->getIndexBase ();
4885  colMap_ = rcp (new map_type (gstInv, myColumns, indexBase,
4886  domainMap_->getComm (),
4887  domainMap_->getNode ()));
4888 
4889  checkInternalState ();
4890  }
4891 
4892 
4893  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4894  void
4895  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4896  mergeAllIndices ()
4897  {
4898  TEUCHOS_TEST_FOR_EXCEPT( isGloballyIndexed() ); // call only after makeIndicesLocal()
4899  TEUCHOS_TEST_FOR_EXCEPT( ! isSorted() ); // call only after sortIndices()
4900  if (! isMerged ()) {
4901  for (size_t row=0; row < getNodeNumRows(); ++row) {
4902  RowInfo rowInfo = getRowInfo(row);
4903  mergeRowIndices(rowInfo);
4904  }
4905  // we just merged every row
4906  noRedundancies_ = true;
4907  }
4908  }
4909 
4910 
4911  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4912  void
4913  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4914  makeImportExport ()
4915  {
4916  using Teuchos::ParameterList;
4917  using Teuchos::RCP;
4918  using Teuchos::rcp;
4919 
4920  TEUCHOS_TEST_FOR_EXCEPTION(! hasColMap (), std::logic_error, "Tpetra::"
4921  "CrsGraph::makeImportExport: This method may not be called unless the "
4922  "graph has a column Map.");
4923  RCP<ParameterList> params = this->getNonconstParameterList (); // could be null
4924 
4925  // Don't do any checks to see if we need to create the Import, if
4926  // it exists already.
4927  //
4928  // FIXME (mfh 25 Mar 2013) This will become incorrect if we
4929  // change CrsGraph in the future to allow changing the column
4930  // Map after fillComplete. For now, the column Map is fixed
4931  // after the first fillComplete call.
4932  if (importer_.is_null ()) {
4933  // Create the Import instance if necessary.
4934  if (domainMap_ != colMap_ && (! domainMap_->isSameAs (*colMap_))) {
4935  if (params.is_null () || ! params->isSublist ("Import")) {
4936  importer_ = rcp (new import_type (domainMap_, colMap_));
4937  } else {
4938  RCP<ParameterList> importSublist = sublist (params, "Import", true);
4939  importer_ = rcp (new import_type (domainMap_, colMap_, importSublist));
4940  }
4941  }
4942  }
4943 
4944  // Don't do any checks to see if we need to create the Export, if
4945  // it exists already.
4946  if (exporter_.is_null ()) {
4947  // Create the Export instance if necessary.
4948  if (rangeMap_ != rowMap_ && ! rangeMap_->isSameAs (*rowMap_)) {
4949  if (params.is_null () || ! params->isSublist ("Export")) {
4950  exporter_ = rcp (new export_type (rowMap_, rangeMap_));
4951  }
4952  else {
4953  RCP<ParameterList> exportSublist = sublist (params, "Export", true);
4954  exporter_ = rcp (new export_type (rowMap_, rangeMap_, exportSublist));
4955  }
4956  }
4957  }
4958  }
4959 
4960 
4961  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4962  std::string
4963  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4964  description () const
4965  {
4966  std::ostringstream oss;
4967  oss << dist_object_type::description ();
4968  if (isFillComplete ()) {
4969  oss << "{status = fill complete"
4970  << ", global rows = " << getGlobalNumRows()
4971  << ", global cols = " << getGlobalNumCols()
4972  << ", global num entries = " << getGlobalNumEntries()
4973  << "}";
4974  }
4975  else {
4976  oss << "{status = fill not complete"
4977  << ", global rows = " << getGlobalNumRows()
4978  << "}";
4979  }
4980  return oss.str();
4981  }
4982 
4983 
4984  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
4985  void
4986  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
4987  describe (Teuchos::FancyOStream &out,
4988  const Teuchos::EVerbosityLevel verbLevel) const
4989  {
4990  const char tfecfFuncName[] = "describe()";
4991  using std::endl;
4992  using std::setw;
4993  using Teuchos::VERB_DEFAULT;
4994  using Teuchos::VERB_NONE;
4995  using Teuchos::VERB_LOW;
4996  using Teuchos::VERB_MEDIUM;
4997  using Teuchos::VERB_HIGH;
4998  using Teuchos::VERB_EXTREME;
4999  Teuchos::EVerbosityLevel vl = verbLevel;
5000  if (vl == VERB_DEFAULT) vl = VERB_LOW;
5001  RCP<const Comm<int> > comm = this->getComm();
5002  const int myImageID = comm->getRank(),
5003  numImages = comm->getSize();
5004  size_t width = 1;
5005  for (size_t dec=10; dec<getGlobalNumRows(); dec *= 10) {
5006  ++width;
5007  }
5008  width = std::max<size_t> (width, static_cast<size_t> (11)) + 2;
5009  Teuchos::OSTab tab (out);
5010  // none: print nothing
5011  // low: print O(1) info from node 0
5012  // medium: print O(P) info, num entries per node
5013  // high: print O(N) info, num entries per row
5014  // extreme: print O(NNZ) info: print graph indices
5015  //
5016  // for medium and higher, print constituent objects at specified verbLevel
5017  if (vl != VERB_NONE) {
5018  if (myImageID == 0) out << this->description() << std::endl;
5019  // O(1) globals, minus what was already printed by description()
5020  if (isFillComplete() && myImageID == 0) {
5021  out << "Global number of diagonals = " << globalNumDiags_ << std::endl;
5022  out << "Global max number of row entries = " << globalMaxNumRowEntries_ << std::endl;
5023  }
5024  // constituent objects
5025  if (vl == VERB_MEDIUM || vl == VERB_HIGH || vl == VERB_EXTREME) {
5026  if (myImageID == 0) out << "\nRow map: " << std::endl;
5027  rowMap_->describe(out,vl);
5028  if (colMap_ != null) {
5029  if (myImageID == 0) out << "\nColumn map: " << std::endl;
5030  colMap_->describe(out,vl);
5031  }
5032  if (domainMap_ != null) {
5033  if (myImageID == 0) out << "\nDomain map: " << std::endl;
5034  domainMap_->describe(out,vl);
5035  }
5036  if (rangeMap_ != null) {
5037  if (myImageID == 0) out << "\nRange map: " << std::endl;
5038  rangeMap_->describe(out,vl);
5039  }
5040  }
5041  // O(P) data
5042  if (vl == VERB_MEDIUM || vl == VERB_HIGH || vl == VERB_EXTREME) {
5043  for (int imageCtr = 0; imageCtr < numImages; ++imageCtr) {
5044  if (myImageID == imageCtr) {
5045  out << "Node ID = " << imageCtr << std::endl
5046  << "Node number of entries = " << nodeNumEntries_ << std::endl
5047  << "Node number of diagonals = " << nodeNumDiags_ << std::endl
5048  << "Node max number of entries = " << nodeMaxNumRowEntries_ << std::endl;
5049  if (indicesAreAllocated()) {
5050  out << "Node number of allocated entries = " << nodeNumAllocated_ << std::endl;
5051  }
5052  else {
5053  out << "Indices are not allocated." << std::endl;
5054  }
5055  }
5056  comm->barrier();
5057  comm->barrier();
5058  comm->barrier();
5059  }
5060  }
5061  // O(N) and O(NNZ) data
5062  if (vl == VERB_HIGH || vl == VERB_EXTREME) {
5063  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5064  ! hasRowInfo (), std::runtime_error, ": reduce verbosity level; "
5065  "graph row information was deleted at fillComplete().");
5066  for (int imageCtr = 0; imageCtr < numImages; ++imageCtr) {
5067  if (myImageID == imageCtr) {
5068  out << std::setw(width) << "Node ID"
5069  << std::setw(width) << "Global Row"
5070  << std::setw(width) << "Num Entries";
5071  if (vl == VERB_EXTREME) {
5072  out << " Entries";
5073  }
5074  out << std::endl;
5075  for (size_t r=0; r < getNodeNumRows(); ++r) {
5076  RowInfo rowinfo = getRowInfo(r);
5077  GlobalOrdinal gid = rowMap_->getGlobalElement(r);
5078  out << std::setw(width) << myImageID
5079  << std::setw(width) << gid
5080  << std::setw(width) << rowinfo.numEntries;
5081  if (vl == VERB_EXTREME) {
5082  out << " ";
5083  if (isGloballyIndexed()) {
5084  ArrayView<const GlobalOrdinal> rowview = getGlobalView(rowinfo);
5085  for (size_t j=0; j < rowinfo.numEntries; ++j) out << rowview[j] << " ";
5086  }
5087  else if (isLocallyIndexed()) {
5088  ArrayView<const LocalOrdinal> rowview = getLocalView(rowinfo);
5089  for (size_t j=0; j < rowinfo.numEntries; ++j) out << colMap_->getGlobalElement(rowview[j]) << " ";
5090  }
5091  }
5092  out << std::endl;
5093  }
5094  }
5095  comm->barrier();
5096  comm->barrier();
5097  comm->barrier();
5098  }
5099  }
5100  }
5101  }
5102 
5103 
5104  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5105  bool
5106  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5107  checkSizes (const SrcDistObject& source)
5108  {
5109  (void) source; // forestall "unused variable" compiler warnings
5110 
5111  // It's not clear what kind of compatibility checks on sizes can
5112  // be performed here. Epetra_CrsGraph doesn't check any sizes for
5113  // compatibility.
5114  return true;
5115  }
5116 
5117 
5118  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5119  void
5120  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5121  copyAndPermute (const SrcDistObject& source,
5122  size_t numSameIDs,
5123  const Teuchos::ArrayView<const LocalOrdinal> &permuteToLIDs,
5124  const Teuchos::ArrayView<const LocalOrdinal> &permuteFromLIDs)
5125  {
5126  using Teuchos::Array;
5127  using Teuchos::ArrayView;
5128  typedef LocalOrdinal LO;
5129  typedef GlobalOrdinal GO;
5130  const char tfecfFuncName[] = "copyAndPermute";
5131  typedef CrsGraph<LO, GO, node_type> this_type;
5132  typedef RowGraph<LO, GO, node_type> row_graph_type;
5133 
5134  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5135  permuteToLIDs.size() != permuteFromLIDs.size(), std::runtime_error,
5136  ": permuteToLIDs and permuteFromLIDs must have the same size.");
5137  // Make sure that the source object has the right type. We only
5138  // actually need it to be a RowGraph, with matching first three
5139  // template parameters. If it's a CrsGraph, we can use view mode
5140  // instead of copy mode to get each row's data.
5141  //
5142  // FIXME (mfh 07 Jul 2013) It should not be necessary for any of
5143  // the template parameters but GO to match. GO has to match
5144  // because the graph has to send indices as global ordinals, if
5145  // the source and target graphs do not have the same column Map.
5146  // If LO doesn't match, the graphs could communicate using global
5147  // indices. It could be possible that Node affects the graph's
5148  // storage format, but packAndPrepare should assume a common
5149  // communication format in any case.
5150  const row_graph_type* srcRowGraph = dynamic_cast<const row_graph_type*> (&source);
5151  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5152  srcRowGraph == NULL, std::invalid_argument,
5153  ": The source object must be a RowGraph with matching first three "
5154  "template parameters.");
5155 
5156  // If the source object is actually a CrsGraph, we can use view
5157  // mode instead of copy mode to access the entries in each row,
5158  // if the graph is not fill complete.
5159  const this_type* srcCrsGraph = dynamic_cast<const this_type*> (&source);
5160 
5161  const map_type& srcRowMap = * (srcRowGraph->getRowMap ());
5162  const map_type& tgtRowMap = * (this->getRowMap ());
5163  const bool src_filled = srcRowGraph->isFillComplete ();
5164  Array<GO> row_copy;
5165  LO myid = 0;
5166 
5167  //
5168  // "Copy" part of "copy and permute."
5169  //
5170  if (src_filled || srcCrsGraph == NULL) {
5171  // If the source graph is fill complete, we can't use view mode,
5172  // because the data might be stored in a different format not
5173  // compatible with the expectations of view mode. Also, if the
5174  // source graph is not a CrsGraph, we can't use view mode,
5175  // because RowGraph only provides copy mode access to the data.
5176  for (size_t i = 0; i < numSameIDs; ++i, ++myid) {
5177  const GO gid = srcRowMap.getGlobalElement (myid);
5178  size_t row_length = srcRowGraph->getNumEntriesInGlobalRow (gid);
5179  row_copy.resize (row_length);
5180  size_t check_row_length = 0;
5181  srcRowGraph->getGlobalRowCopy (gid, row_copy (), check_row_length);
5182  this->insertGlobalIndices (gid, row_copy ());
5183  }
5184  } else {
5185  for (size_t i = 0; i < numSameIDs; ++i, ++myid) {
5186  const GO gid = srcRowMap.getGlobalElement (myid);
5187  ArrayView<const GO> row;
5188  srcCrsGraph->getGlobalRowView (gid, row);
5189  this->insertGlobalIndices (gid, row);
5190  }
5191  }
5192 
5193  //
5194  // "Permute" part of "copy and permute."
5195  //
5196  if (src_filled || srcCrsGraph == NULL) {
5197  for (LO i = 0; i < permuteToLIDs.size (); ++i) {
5198  const GO mygid = tgtRowMap.getGlobalElement (permuteToLIDs[i]);
5199  const GO srcgid = srcRowMap.getGlobalElement (permuteFromLIDs[i]);
5200  size_t row_length = srcRowGraph->getNumEntriesInGlobalRow (srcgid);
5201  row_copy.resize (row_length);
5202  size_t check_row_length = 0;
5203  srcRowGraph->getGlobalRowCopy (srcgid, row_copy (), check_row_length);
5204  this->insertGlobalIndices (mygid, row_copy ());
5205  }
5206  } else {
5207  for (LO i = 0; i < permuteToLIDs.size (); ++i) {
5208  const GO mygid = tgtRowMap.getGlobalElement (permuteToLIDs[i]);
5209  const GO srcgid = srcRowMap.getGlobalElement (permuteFromLIDs[i]);
5210  ArrayView<const GO> row;
5211  srcCrsGraph->getGlobalRowView (srcgid, row);
5212  this->insertGlobalIndices (mygid, row);
5213  }
5214  }
5215  }
5216 
5217 
5218  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5219  void
5220  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5221  packAndPrepare (const SrcDistObject& source,
5222  const Teuchos::ArrayView<const LocalOrdinal> &exportLIDs,
5223  Teuchos::Array<GlobalOrdinal> &exports,
5224  const Teuchos::ArrayView<size_t> & numPacketsPerLID,
5225  size_t& constantNumPackets,
5226  Distributor &distor)
5227  {
5228  using Teuchos::Array;
5229  const char tfecfFuncName[] = "packAndPrepare";
5230  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5231  exportLIDs.size() != numPacketsPerLID.size(), std::runtime_error,
5232  ": exportLIDs and numPacketsPerLID must have the same size.");
5233  typedef RowGraph<LocalOrdinal, GlobalOrdinal, node_type> row_graph_type;
5234  const row_graph_type& srcGraph = dynamic_cast<const row_graph_type&> (source);
5235 
5236  // We don't check whether src_graph has had fillComplete called,
5237  // because it doesn't matter whether the *source* graph has been
5238  // fillComplete'd. The target graph can not be fillComplete'd yet.
5239  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5240  this->isFillComplete (), std::runtime_error,
5241  ": The target graph of an Import or Export must not be fill complete.");
5242 
5243  srcGraph.pack (exportLIDs, exports, numPacketsPerLID, constantNumPackets, distor);
5244  }
5245 
5246 
5247  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5248  void
5249  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5250  pack (const Teuchos::ArrayView<const LocalOrdinal>& exportLIDs,
5251  Teuchos::Array<GlobalOrdinal>& exports,
5252  const Teuchos::ArrayView<size_t>& numPacketsPerLID,
5253  size_t& constantNumPackets,
5254  Distributor& distor) const
5255  {
5256  using Teuchos::Array;
5257  typedef LocalOrdinal LO;
5258  typedef GlobalOrdinal GO;
5259  const char tfecfFuncName[] = "pack";
5260  (void) distor; // forestall "unused argument" compiler warning
5261 
5262  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5263  exportLIDs.size() != numPacketsPerLID.size(), std::runtime_error,
5264  ": exportLIDs and numPacketsPerLID must have the same size.");
5265 
5266  const map_type& srcMap = * (this->getMap ());
5267  constantNumPackets = 0;
5268 
5269  // Set numPacketsPerLID[i] to the number of entries owned by the
5270  // calling process in (local) row exportLIDs[i] of the graph, that
5271  // the caller wants us to send out. Compute the total number of
5272  // packets (that is, entries) owned by this process in all the
5273  // rows that the caller wants us to send out.
5274  size_t totalNumPackets = 0;
5275  Array<GO> row;
5276  for (LO i = 0; i < exportLIDs.size (); ++i) {
5277  const GO GID = srcMap.getGlobalElement (exportLIDs[i]);
5278  size_t row_length = this->getNumEntriesInGlobalRow (GID);
5279  numPacketsPerLID[i] = row_length;
5280  totalNumPackets += row_length;
5281  }
5282 
5283  exports.resize (totalNumPackets);
5284 
5285  // Loop again over the rows to export, and pack rows of indices
5286  // into the output buffer.
5287  size_t exportsOffset = 0;
5288  for (LO i = 0; i < exportLIDs.size (); ++i) {
5289  const GO GID = srcMap.getGlobalElement (exportLIDs[i]);
5290  size_t row_length = this->getNumEntriesInGlobalRow (GID);
5291  row.resize (row_length);
5292  size_t check_row_length = 0;
5293  this->getGlobalRowCopy (GID, row (), check_row_length);
5294  typename Array<GO>::const_iterator row_iter = row.begin();
5295  typename Array<GO>::const_iterator row_end = row.end();
5296  size_t j = 0;
5297  for (; row_iter != row_end; ++row_iter, ++j) {
5298  exports[exportsOffset+j] = *row_iter;
5299  }
5300  exportsOffset += row.size();
5301  }
5302  }
5303 
5304 
5305  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5306  void
5307  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5308  unpackAndCombine (const Teuchos::ArrayView<const LocalOrdinal> &importLIDs,
5309  const Teuchos::ArrayView<const GlobalOrdinal> &imports,
5310  const Teuchos::ArrayView<size_t> &numPacketsPerLID,
5311  size_t constantNumPackets,
5312  Distributor& /* distor */,
5313  CombineMode /* CM */)
5314  {
5315  using Teuchos::ArrayView;
5316  typedef LocalOrdinal LO;
5317  typedef GlobalOrdinal GO;
5318 
5319  // FIXME (mfh 02 Apr 2012) REPLACE combine mode has a perfectly
5320  // reasonable meaning, whether or not the matrix is fill complete.
5321  // It's just more work to implement.
5322 
5323  // We are not checking the value of the CombineMode input
5324  // argument. For CrsGraph, we only support import/export
5325  // operations if fillComplete has not yet been called. Any
5326  // incoming column-indices are inserted into the target graph. In
5327  // this context, CombineMode values of ADD vs INSERT are
5328  // equivalent. What is the meaning of REPLACE for CrsGraph? If a
5329  // duplicate column-index is inserted, it will be compressed out
5330  // when fillComplete is called.
5331  //
5332  // Note: I think REPLACE means that an existing row is replaced by
5333  // the imported row, i.e., the existing indices are cleared. CGB,
5334  // 6/17/2010
5335 
5336  const char tfecfFuncName[] = "unpackAndCombine: ";
5337  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5338  importLIDs.size() != numPacketsPerLID.size(), std::runtime_error,
5339  "importLIDs and numPacketsPerLID must have the same size.");
5340  TEUCHOS_TEST_FOR_EXCEPTION_CLASS_FUNC(
5341  isFillComplete (), std::runtime_error,
5342  "Import or Export operations are not allowed on the destination "
5343  "CrsGraph if it is fill complete.");
5344  size_t importsOffset = 0;
5345 
5346  typedef typename ArrayView<const LO>::const_iterator iter_type;
5347  iter_type impLIDiter = importLIDs.begin();
5348  iter_type impLIDend = importLIDs.end();
5349 
5350  for (size_t i = 0; impLIDiter != impLIDend; ++impLIDiter, ++i) {
5351  LO row_length = numPacketsPerLID[i];
5352 
5353  const GO* const row_raw = (row_length == 0) ? NULL : &imports[importsOffset];
5354  ArrayView<const GlobalOrdinal> row (row_raw, row_length);
5355  insertGlobalIndicesFiltered (this->getMap ()->getGlobalElement (*impLIDiter), row);
5356  importsOffset += row_length;
5357  }
5358  }
5359 
5360 
5361  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5362  void
5363  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5364  removeEmptyProcessesInPlace (const Teuchos::RCP<const map_type>& newMap)
5365  {
5366  using Teuchos::Comm;
5367  using Teuchos::null;
5368  using Teuchos::ParameterList;
5369  using Teuchos::RCP;
5370 
5371  // We'll set all the state "transactionally," so that this method
5372  // satisfies the strong exception guarantee. This object's state
5373  // won't be modified until the end of this method.
5374  RCP<const map_type> rowMap, domainMap, rangeMap, colMap;
5375  RCP<import_type> importer;
5376  RCP<export_type> exporter;
5377 
5378  rowMap = newMap;
5379  RCP<const Comm<int> > newComm =
5380  (newMap.is_null ()) ? null : newMap->getComm ();
5381 
5382  if (! domainMap_.is_null ()) {
5383  if (domainMap_.getRawPtr () == rowMap_.getRawPtr ()) {
5384  // Common case: original domain and row Maps are identical.
5385  // In that case, we need only replace the original domain Map
5386  // with the new Map. This ensures that the new domain and row
5387  // Maps _stay_ identical.
5388  domainMap = newMap;
5389  } else {
5390  domainMap = domainMap_->replaceCommWithSubset (newComm);
5391  }
5392  }
5393  if (! rangeMap_.is_null ()) {
5394  if (rangeMap_.getRawPtr () == rowMap_.getRawPtr ()) {
5395  // Common case: original range and row Maps are identical. In
5396  // that case, we need only replace the original range Map with
5397  // the new Map. This ensures that the new range and row Maps
5398  // _stay_ identical.
5399  rangeMap = newMap;
5400  } else {
5401  rangeMap = rangeMap_->replaceCommWithSubset (newComm);
5402  }
5403  }
5404  if (! colMap.is_null ()) {
5405  colMap = colMap_->replaceCommWithSubset (newComm);
5406  }
5407 
5408  // (Re)create the Export and / or Import if necessary.
5409  if (! newComm.is_null ()) {
5410  RCP<ParameterList> params = this->getNonconstParameterList (); // could be null
5411  //
5412  // The operations below are collective on the new communicator.
5413  //
5414  // (Re)create the Export object if necessary. If I haven't
5415  // called fillComplete yet, I don't have a rangeMap, so I must
5416  // first check if the _original_ rangeMap is not null. Ditto
5417  // for the Import object and the domain Map.
5418  if (! rangeMap_.is_null () &&
5419  rangeMap != rowMap &&
5420  ! rangeMap->isSameAs (*rowMap)) {
5421  if (params.is_null () || ! params->isSublist ("Export")) {
5422  exporter = rcp (new export_type (rowMap, rangeMap));
5423  }
5424  else {
5425  RCP<ParameterList> exportSublist = sublist (params, "Export", true);
5426  exporter = rcp (new export_type (rowMap, rangeMap, exportSublist));
5427  }
5428  }
5429  // (Re)create the Import object if necessary.
5430  if (! domainMap_.is_null () &&
5431  domainMap != colMap &&
5432  ! domainMap->isSameAs (*colMap)) {
5433  if (params.is_null () || ! params->isSublist ("Import")) {
5434  importer = rcp (new import_type (domainMap, colMap));
5435  } else {
5436  RCP<ParameterList> importSublist = sublist (params, "Import", true);
5437  importer = rcp (new import_type (domainMap, colMap, importSublist));
5438  }
5439  }
5440  } // if newComm is not null
5441 
5442  // Defer side effects until the end. If no destructors throw
5443  // exceptions (they shouldn't anyway), then this method satisfies
5444  // the strong exception guarantee.
5445  exporter_ = exporter;
5446  importer_ = importer;
5447  rowMap_ = rowMap;
5448  // mfh 31 Mar 2013: DistObject's map_ is the row Map of a CrsGraph
5449  // or CrsMatrix. CrsGraph keeps a redundant pointer (rowMap_) to
5450  // the same object. We might want to get rid of this redundant
5451  // pointer sometime, but for now, we'll leave it alone and just
5452  // set map_ to the same object.
5453  this->map_ = rowMap;
5454  domainMap_ = domainMap;
5455  rangeMap_ = rangeMap;
5456  colMap_ = colMap;
5457  }
5458 
5459  template <class LocalOrdinal, class GlobalOrdinal, class Node, const bool classic>
5460  bool
5461  CrsGraph<LocalOrdinal, GlobalOrdinal, Node, classic>::
5462  hasRowInfo () const
5463  {
5464  if (indicesAreAllocated () &&
5465  getProfileType () == StaticProfile &&
5466  k_rowPtrs_.dimension_0 () == 0) {
5467  return false;
5468  } else {
5469  return true;
5470  }
5471  }
5472 
5473 } // namespace Tpetra
5474 
5475 //
5476 // Explicit instantiation macro
5477 //
5478 // Must be expanded from within the Tpetra namespace!
5479 //
5480 #define TPETRA_CRSGRAPH_GRAPH_INSTANT(LO,GO,NODE) template class CrsGraph< LO , GO , NODE >;
5481 #define TPETRA_CRSGRAPH_SORTROWINDICESANDVALUES_INSTANT(S,LO,GO,NODE) template void CrsGraph< LO , GO , NODE >::sortRowIndicesAndValues< S >(const RowInfo, const Teuchos::ArrayView< S >& );
5482 #define TPETRA_CRSGRAPH_MERGEROWINDICESANDVALUES_INSTANT(S,LO,GO,NODE) template void CrsGraph< LO , GO , NODE >::mergeRowIndicesAndValues< S >(RowInfo, const ArrayView< S >& );
5483 #define TPETRA_CRSGRAPH_ALLOCATEVALUES1D_INSTANT(S,LO,GO,NODE)
5484 #define TPETRA_CRSGRAPH_ALLOCATEVALUES2D_INSTANT(S,LO,GO,NODE) template ArrayRCP< Array< S > > CrsGraph< LO , GO , NODE >::allocateValues2D< S >() const;
5485 
5486 #define TPETRA_CRSGRAPH_INSTANT(S,LO,GO,NODE) \
5487  TPETRA_CRSGRAPH_SORTROWINDICESANDVALUES_INSTANT(S,LO,GO,NODE) \
5488  TPETRA_CRSGRAPH_MERGEROWINDICESANDVALUES_INSTANT(S,LO,GO,NODE) \
5489  TPETRA_CRSGRAPH_ALLOCATEVALUES1D_INSTANT(S,LO,GO,NODE) \
5490  TPETRA_CRSGRAPH_ALLOCATEVALUES2D_INSTANT(S,LO,GO,NODE) \
5491  template void CrsGraph< LO , GO , NODE >::insertIndicesAndValues< S > (const RowInfo&, const SLocalGlobalViews&, const Teuchos::ArrayView< S >&, const Teuchos::ArrayView<const S >&, const ELocalGlobal, const ELocalGlobal);
5492 
5493 #endif // TPETRA_CRSGRAPH_DEF_HPP
Tpetra::CrsGraph::getNodePackedIndices
Teuchos::ArrayRCP< const LocalOrdinal > getNodePackedIndices() const
Get an Teuchos::ArrayRCP of the packed column-indices.
Definition: Tpetra_CrsGraph_def.hpp:2324
Tpetra::Import
Communication plan for data redistribution from a uniquely-owned to a (possibly) multiply-owned distr...
Definition: Tpetra_Import_decl.hpp:124
Tpetra::CrsGraph::isGloballyIndexed
bool isGloballyIndexed() const
If graph indices are in the global range, this function returns true. Otherwise, this function return...
Definition: Tpetra_CrsGraph_def.hpp:902
Tpetra::CrsGraph::getComm
Teuchos::RCP< const Comm< int > > getComm() const
Returns the communicator.
Definition: Tpetra_CrsGraph_def.hpp:920
Tpetra
Namespace Tpetra contains the class and methods constituting the Tpetra library.
Tpetra::CrsGraph::haveGlobalConstants_
bool haveGlobalConstants_
Whether all processes have computed global constants.
Definition: Tpetra_CrsGraph_decl.hpp:1951
Tpetra::CrsGraph::copyAndPermute
virtual void copyAndPermute(const SrcDistObject &source, size_t numSameIDs, const Teuchos::ArrayView< const LocalOrdinal > &permuteToLIDs, const Teuchos::ArrayView< const LocalOrdinal > &permuteFromLIDs)
Perform copies and permutations that are local to this process.
Definition: Tpetra_CrsGraph_def.hpp:5121
Tpetra::DistObject< GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node >::description
virtual std::string description() const
One-line descriptiion of this object.
Tpetra::CrsGraph::description
std::string description() const
Return a simple one-line description of this object.
Definition: Tpetra_CrsGraph_def.hpp:4964
Tpetra::CrsGraph::expertStaticFillComplete
void expertStaticFillComplete(const Teuchos::RCP< const map_type > &domainMap, const Teuchos::RCP< const map_type > &rangeMap, const Teuchos::RCP< const import_type > &importer=Teuchos::null, const Teuchos::RCP< const export_type > &exporter=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Perform a fillComplete on a graph that already has data, via setAllIndices().
Definition: Tpetra_CrsGraph_def.hpp:3505
Tpetra::CrsGraph::getRowMap
Teuchos::RCP< const map_type > getRowMap() const
Returns the Map that describes the row distribution in this graph.
Definition: Tpetra_CrsGraph_def.hpp:725
Tpetra::CrsGraph::insertLocalIndicesFiltered
void insertLocalIndicesFiltered(const LocalOrdinal localRow, const Teuchos::ArrayView< const LocalOrdinal > &indices)
Like insertLocalIndices(), but with column Map filtering.
Definition: Tpetra_CrsGraph_def.hpp:2584
Tpetra::CrsGraph::isFillActive
bool isFillActive() const
Returns true if resumeFill() has been called and the graph is in edit mode.
Definition: Tpetra_CrsGraph_def.hpp:866
Tpetra::CrsGraph::getLocalViewNonConst
Teuchos::ArrayView< LocalOrdinal > getLocalViewNonConst(const RowInfo rowinfo)
Get a nonconst, nonowned, locally indexed view of the locally owned row myRow, such that rowinfo = ge...
Definition: Tpetra_CrsGraph_def.hpp:1224
Tpetra::CrsGraph::removeEmptyProcessesInPlace
virtual void removeEmptyProcessesInPlace(const Teuchos::RCP< const map_type > &newMap)
Remove processes owning zero rows from the Maps and their communicator.
Definition: Tpetra_CrsGraph_def.hpp:5364
Tpetra::CrsGraph::insertGlobalIndices
void insertGlobalIndices(GlobalOrdinal globalRow, const Teuchos::ArrayView< const GlobalOrdinal > &indices)
Insert global indices into the graph.
Definition: Tpetra_CrsGraph_def.hpp:2635
Tpetra::CrsGraph::getRangeMap
Teuchos::RCP< const map_type > getRangeMap() const
Returns the Map associated with the domain of this graph.
Definition: Tpetra_CrsGraph_def.hpp:752
Tpetra::Map::isNodeGlobalElement
bool isNodeGlobalElement(GlobalOrdinal globalIndex) const
Whether the given global index is owned by this Map on the calling process.
Definition: Tpetra_Map_def.hpp:719
Tpetra::CrsGraph::indicesAreSorted_
bool indicesAreSorted_
Whether the graph&#39;s indices are sorted in each row, on this process.
Definition: Tpetra_CrsGraph_decl.hpp:1944
Tpetra::CrsGraph::getNodeNumCols
size_t getNodeNumCols() const
Returns the number of columns connected to the locally owned rows of this graph.
Definition: Tpetra_CrsGraph_def.hpp:681
Tpetra::CrsGraph::hasRowInfo
bool hasRowInfo() const
Whether it is correct to call getRowInfo().
Definition: Tpetra_CrsGraph_def.hpp:5462
Tpetra::CrsGraph::removeLocalIndices
void removeLocalIndices(LocalOrdinal localRow)
Remove all graph indices from the specified local row.
Definition: Tpetra_CrsGraph_def.hpp:2780
Tpetra::CrsGraph::noRedundancies_
bool noRedundancies_
Whether the graph&#39;s indices are non-redundant (merged) in each row, on this process.
Definition: Tpetra_CrsGraph_decl.hpp:1947
Tpetra::CrsGraph::getLocalGraph
local_graph_type getLocalGraph() const
Get the local graph.
Definition: Tpetra_CrsGraph_def.hpp:4292
Tpetra::CrsGraph::sortGhostsAssociatedWithEachProcessor_
bool sortGhostsAssociatedWithEachProcessor_
Whether to require makeColMap() (and therefore fillComplete()) to order column Map GIDs associated wi...
Definition: Tpetra_CrsGraph_decl.hpp:1970
Tpetra::CrsGraph::device_type
Node::device_type device_type
This class&#39; Kokkos device type.
Definition: Tpetra_CrsGraph_decl.hpp:259
Tpetra::CrsGraph::hasColMap
bool hasColMap() const
Whether the graph has a column Map.
Definition: Tpetra_CrsGraph_def.hpp:778
Tpetra::CrsGraph::rangeMap_
Teuchos::RCP< const map_type > rangeMap_
The Map describing the range of the (matrix corresponding to the) graph.
Definition: Tpetra_CrsGraph_decl.hpp:1757
Tpetra::CrsGraph::pftype_
ProfileType pftype_
Whether the graph was allocated with static or dynamic profile.
Definition: Tpetra_CrsGraph_decl.hpp:1786
Tpetra::CrsGraph::getGlobalRowCopy
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const Teuchos::ArrayView< GlobalOrdinal > &Indices, size_t &NumIndices) const
Get a copy of the given row, using global indices.
Definition: Tpetra_CrsGraph_def.hpp:2397
Tpetra::CrsGraph::mergeRowIndices
void mergeRowIndices(RowInfo rowinfo)
Merge duplicate row indices in the given row.
Definition: Tpetra_CrsGraph_def.hpp:1746
Tpetra::CrsGraph::getGlobalMaxNumRowEntries
size_t getGlobalMaxNumRowEntries() const
Maximum number of entries in all rows over all processes.
Definition: Tpetra_CrsGraph_def.hpp:839
Tpetra::LookupStatus
LookupStatus
Return status of Map remote index lookup (getRemoteIndexList()).
Definition: Tpetra_ConfigDefs.hpp:167
Tpetra::CrsGraph::getGlobalNumRows
global_size_t getGlobalNumRows() const
Returns the number of global rows in the graph.
Definition: Tpetra_CrsGraph_def.hpp:648
Tpetra::CrsGraph::isSorted
bool isSorted() const
Whether graph indices in all rows are known to be sorted.
Definition: Tpetra_CrsGraph_def.hpp:947
Tpetra::CrsGraph::k_numRowEntries_
t_numRowEntries_ k_numRowEntries_
The number of local entries in each locally owned row.
Definition: Tpetra_CrsGraph_decl.hpp:1919
Tpetra::CrsGraph::setAllIndices
void setAllIndices(const typename local_graph_type::row_map_type &rowPointers, const typename local_graph_type::entries_type::non_const_type &columnIndices)
Set the graph&#39;s data directly, using 1-D storage.
Definition: Tpetra_CrsGraph_def.hpp:2824
Tpetra::CrsGraph::getProfileType
ProfileType getProfileType() const
Returns true if the graph was allocated with static data structures.
Definition: Tpetra_CrsGraph_def.hpp:812
Tpetra::CrsGraph::k_lclInds1D_
local_graph_type::entries_type::non_const_type k_lclInds1D_
Local column indices for all rows.
Definition: Tpetra_CrsGraph_decl.hpp:1828
Tpetra::CrsGraph::mergeRowIndicesAndValues
void mergeRowIndicesAndValues(RowInfo rowinfo, const Teuchos::ArrayView< Scalar > &rowValues)
Merge duplicate row indices in the given row, along with their corresponding values.
Definition: Tpetra_CrsGraph_def.hpp:1781
Tpetra::CrsGraph::gblInds2D_
Teuchos::ArrayRCP< Teuchos::Array< GlobalOrdinal > > gblInds2D_
Global column indices for all rows.
Definition: Tpetra_CrsGraph_decl.hpp:1909
Tpetra::deep_copy
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.
Definition: Tpetra_MultiVector_decl.hpp:2120
Tpetra::CrsGraph::resumeFill
void resumeFill(const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Definition: Tpetra_CrsGraph_def.hpp:3329
Tpetra::CrsGraph::getImporter
Teuchos::RCP< const import_type > getImporter() const
Returns the importer associated with this graph.
Definition: Tpetra_CrsGraph_def.hpp:760
Tpetra::CrsGraph::checkInternalState
void checkInternalState() const
Throw an exception if the internal state is not consistent.
Definition: Tpetra_CrsGraph_def.hpp:1978
Tpetra::CrsGraph::isStorageOptimized
bool isStorageOptimized() const
Returns true if storage has been optimized.
Definition: Tpetra_CrsGraph_def.hpp:787
Tpetra::CrsGraph::fillComplete
void fillComplete(const Teuchos::RCP< const map_type > &domainMap, const Teuchos::RCP< const map_type > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Signal that data entry is complete, specifying domain and range maps.
Definition: Tpetra_CrsGraph_def.hpp:3387
Tpetra::CrsGraph::upperTriangular_
bool upperTriangular_
Whether the graph is locally upper triangular.
Definition: Tpetra_CrsGraph_decl.hpp:1942
Tpetra::CrsGraph::setParameterList
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &params)
Set the given list of parameters (must be nonnull).
Definition: Tpetra_CrsGraph_def.hpp:636
Tpetra::CrsGraph::insertGlobalIndicesFiltered
void insertGlobalIndicesFiltered(const GlobalOrdinal localRow, const Teuchos::ArrayView< const GlobalOrdinal > &indices)
Like insertGlobalIndices(), but with column Map filtering.
Definition: Tpetra_CrsGraph_def.hpp:2717
Tpetra::CrsGraph::domainMap_
Teuchos::RCP< const map_type > domainMap_
The Map describing the domain of the (matrix corresponding to the) graph.
Definition: Tpetra_CrsGraph_decl.hpp:1759
Tpetra::RowInfo
Allocation information for a locally owned row in a CrsGraph or CrsMatrix.
Definition: Tpetra_CrsGraph_decl.hpp:108
Tpetra::CrsGraph::getNodeMaxNumRowEntries
size_t getNodeMaxNumRowEntries() const
Maximum number of entries in all rows owned by the calling process.
Definition: Tpetra_CrsGraph_def.hpp:848
Tpetra::Map::isNodeLocalElement
bool isNodeLocalElement(LocalOrdinal localIndex) const
Whether the given local index is valid for this Map on the calling process.
Definition: Tpetra_Map_def.hpp:707
Tpetra::CrsGraph::isFillComplete
bool isFillComplete() const
Returns true if fillComplete() has been called and the graph is in compute mode.
Definition: Tpetra_CrsGraph_def.hpp:857
Tpetra::CrsGraph::replaceColMap
void replaceColMap(const Teuchos::RCP< const map_type > &newColMap)
Replace the graph&#39;s current column Map with the given Map.
Definition: Tpetra_CrsGraph_def.hpp:3963
Tpetra::CrsGraph::colMap_
Teuchos::RCP< const map_type > colMap_
The Map describing the distribution of columns of the graph.
Definition: Tpetra_CrsGraph_decl.hpp:1755
Tpetra::CrsGraph::getNodeAllocationSize
size_t getNodeAllocationSize() const
Returns the total number of indices allocated for the graph, across all rows on this node...
Definition: Tpetra_CrsGraph_def.hpp:911
Details
Implementation details of Tpetra.
Tpetra::CrsGraph::getGlobalNumEntries
global_size_t getGlobalNumEntries() const
Returns the global number of entries in the graph.
Definition: Tpetra_CrsGraph_def.hpp:821
Tpetra::CrsGraph::numAllocForAllRows_
size_t numAllocForAllRows_
The maximum number of entries to allow in each locally owned row.
Definition: Tpetra_CrsGraph_decl.hpp:1816
Tpetra::global_size_t
size_t global_size_t
Global size_t object.
Definition: Tpetra_ConfigDefs.hpp:152
Tpetra::CrsGraph::local_graph_type
Kokkos::StaticCrsGraph< LocalOrdinal, Kokkos::LayoutLeft, execution_space > local_graph_type
The type of the part of the sparse graph on each MPI process.
Definition: Tpetra_CrsGraph_decl.hpp:266
Tpetra::CrsGraph::isLocallyIndexed
bool isLocallyIndexed() const
If graph indices are in the local range, this function returns true. Otherwise, this function returns...
Definition: Tpetra_CrsGraph_def.hpp:893
Tpetra::CrsGraph::getDomainMap
Teuchos::RCP< const map_type > getDomainMap() const
Returns the Map associated with the domain of this graph.
Definition: Tpetra_CrsGraph_def.hpp:743
Tpetra::CrsGraph::haveLocalConstants_
bool haveLocalConstants_
Whether this process has computed local constants.
Definition: Tpetra_CrsGraph_decl.hpp:1949
Tpetra::CrsGraph::isMerged
bool isMerged() const
Whether duplicate column indices in each row have been merged.
Definition: Tpetra_CrsGraph_def.hpp:956
Tpetra::CrsGraph::k_gblInds1D_
t_GlobalOrdinal_1D k_gblInds1D_
Global column indices for all rows.
Definition: Tpetra_CrsGraph_decl.hpp:1840
Tpetra::CrsGraph::getNode
Teuchos::RCP< node_type > getNode() const
Returns the underlying node.
Definition: Tpetra_CrsGraph_def.hpp:716
Tpetra::CrsGraph::getIndexBase
GlobalOrdinal getIndexBase() const
Returns the index base for global indices for this graph.
Definition: Tpetra_CrsGraph_def.hpp:929
Tpetra::CrsGraph::findGlobalIndex
size_t findGlobalIndex(RowInfo rowinfo, GlobalOrdinal ind, size_t hint=0) const
Find the column offset corresponding to the given (global) column index.
Definition: Tpetra_CrsGraph_def.hpp:1924
Tpetra::CrsGraph::lclGraph_
local_graph_type lclGraph_
Local graph; only initialized after first fillComplete() call.
Definition: Tpetra_CrsGraph_decl.hpp:1777
Tpetra::CrsGraph::globalAssemble
void globalAssemble()
Communicate non-local contributions to other processes.
Definition: Tpetra_CrsGraph_def.hpp:3041
Tpetra::Export
Communication plan for data redistribution from a (possibly) multiply-owned to a uniquely-owned distr...
Definition: Tpetra_Export_decl.hpp:124
Tpetra::CrsGraph::sortRowIndicesAndValues
void sortRowIndicesAndValues(const RowInfo rowinfo, const Teuchos::ArrayView< Scalar > &values)
Sort the column indices and their values in the given row.
Definition: Tpetra_CrsGraph_def.hpp:1731
Tpetra::CrsGraph::insertIndicesAndValues
void insertIndicesAndValues(const RowInfo &rowInfo, const SLocalGlobalViews &newInds, const Teuchos::ArrayView< Scalar > &oldRowVals, const Teuchos::ArrayView< const Scalar > &newRowVals, const ELocalGlobal lg, const ELocalGlobal I)
Insert indices and their values into the given row.
Definition: Tpetra_CrsGraph_def.hpp:1600
Tpetra::CrsGraph::getLocalGraph_Kokkos
TPETRA_DEPRECATED local_graph_type getLocalGraph_Kokkos() const
Get the local graph (DEPRECATED: call getLocalGraph() instead).
Definition: Tpetra_CrsGraph_def.hpp:4300
Tpetra::CrsGraph::isUpperTriangular
bool isUpperTriangular() const
Whether the graph is locally upper triangular.
Definition: Tpetra_CrsGraph_def.hpp:875
Tpetra::Distributor
Sets up and executes a communication plan for a Tpetra DistObject.
Definition: Tpetra_Distributor.hpp:186
Tpetra::Map::getGlobalElement
GlobalOrdinal getGlobalElement(LocalOrdinal localIndex) const
The global index corresponding to the given local index.
Definition: Tpetra_Map_def.hpp:689
Tpetra::CombineMode
CombineMode
Rule for combining data in an Import or Export.
Definition: Tpetra_CombineMode.hpp:94
Tpetra::DistObject< GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node >::map_
Teuchos::RCP< const map_type > map_
The Map over which this object is distributed.
Definition: Tpetra_DistObject_decl.hpp:695
Tpetra::CrsGraph::makeColMap
void makeColMap()
Make the graph&#39;s column Map, if it does not already have one.
Definition: Tpetra_CrsGraph_def.hpp:4564
Tpetra::CrsGraph::getLocalRowCopy
void getLocalRowCopy(LocalOrdinal LocalRow, const Teuchos::ArrayView< LocalOrdinal > &indices, size_t &NumIndices) const
Get a copy of the given row, using local indices.
Definition: Tpetra_CrsGraph_def.hpp:2333
Tpetra::CrsGraph::k_numAllocPerRow_
Kokkos::DualView< const size_t *, Kokkos::LayoutLeft, execution_space > k_numAllocPerRow_
The maximum number of entries to allow in each locally owned row, per row.
Definition: Tpetra_CrsGraph_decl.hpp:1805
Tpetra::CrsGraph::getGlobalNumCols
global_size_t getGlobalNumCols() const
Returns the number of global columns in the graph.
Definition: Tpetra_CrsGraph_def.hpp:657
Tpetra::SrcDistObject
Abstract base class for objects that can be the source of an Import or Export operation.
Definition: Tpetra_SrcDistObject.hpp:89
Tpetra::CrsGraph::importer_
Teuchos::RCP< const import_type > importer_
The Import from the domain Map to the column Map.
Definition: Tpetra_CrsGraph_decl.hpp:1767
Tpetra::CrsGraph::node_type
Node node_type
This class&#39; Kokkos Node type.
Definition: Tpetra_CrsGraph_decl.hpp:256
Tpetra::CrsGraph::reindexColumns
void reindexColumns(const Teuchos::RCP< const map_type > &newColMap, const Teuchos::RCP< const import_type > &newImport=Teuchos::null, const bool sortIndicesInEachRow=true)
Reindex the column indices in place, and replace the column Map. Optionally, replace the Import objec...
Definition: Tpetra_CrsGraph_def.hpp:3982
Tpetra::CrsGraph::getColMap
Teuchos::RCP< const map_type > getColMap() const
Returns the Map that describes the column distribution in this graph.
Definition: Tpetra_CrsGraph_def.hpp:734
Tpetra::CrsGraph::insertIndices
size_t insertIndices(const RowInfo &rowInfo, const SLocalGlobalViews &newInds, const ELocalGlobal lg, const ELocalGlobal I)
Insert indices into the given row.
Definition: Tpetra_CrsGraph_def.hpp:1416
Tpetra::CrsGraph::getGlobalViewNonConst
Teuchos::ArrayView< GlobalOrdinal > getGlobalViewNonConst(const RowInfo rowinfo)
Get a nonconst, nonowned, globally indexed view of the locally owned row myRow, such that rowinfo = g...
Definition: Tpetra_CrsGraph_def.hpp:1276
Tpetra::CrsGraph::rowMap_
Teuchos::RCP< const map_type > rowMap_
The Map describing the distribution of rows of the graph.
Definition: Tpetra_CrsGraph_decl.hpp:1753
Tpetra::CrsGraph::getNodeNumEntries
size_t getNodeNumEntries() const
Returns the local number of entries in the graph.
Definition: Tpetra_CrsGraph_def.hpp:830
Tpetra::CrsGraph::pack
virtual void pack(const Teuchos::ArrayView< const LocalOrdinal > &exportLIDs, Teuchos::Array< GlobalOrdinal > &exports, const Teuchos::ArrayView< size_t > &numPacketsPerLID, size_t &constantNumPackets, Distributor &distor) const
Pack this object&#39;s data for Import or Export.
Definition: Tpetra_CrsGraph_def.hpp:5250
Tpetra::CrsGraph::lclInds2D_
Teuchos::ArrayRCP< Teuchos::Array< LocalOrdinal > > lclInds2D_
Local column indices for all rows.
Definition: Tpetra_CrsGraph_decl.hpp:1896
Tpetra::CrsGraph::getGlobalNumDiags
global_size_t getGlobalNumDiags() const
Returns the number of global diagonal entries, based on global row/column index comparisons.
Definition: Tpetra_CrsGraph_def.hpp:707
Tpetra::CrsGraph::findLocalIndex
size_t findLocalIndex(RowInfo rowinfo, LocalOrdinal ind, size_t hint=0) const
Find the column offset corresponding to the given (local) column index.
Definition: Tpetra_CrsGraph_def.hpp:1865
Tpetra::Map::getLocalElement
LocalOrdinal getLocalElement(GlobalOrdinal globalIndex) const
The local index corresponding to the given global index.
Definition: Tpetra_Map_def.hpp:666
Tpetra::sort2
void sort2(const IT1 &first1, const IT1 &last1, const IT2 &first2)
Sort the first array, and apply the resulting permutation to the second array.
Definition: Tpetra_Util.hpp:545
Tpetra::CrsGraph
A distributed graph accessed by rows (adjacency lists) and stored sparsely.
Definition: Tpetra_CrsGraph_decl.hpp:230
Tpetra::CrsGraph::getNumAllocatedEntriesInGlobalRow
size_t getNumAllocatedEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of allocated entries for this node in the specified global row ...
Definition: Tpetra_CrsGraph_def.hpp:2183
Tpetra::Map
Describes a parallel distribution of objects over processes.
Definition: Tpetra_DirectoryImpl_decl.hpp:69
Tpetra::CrsGraph::setDomainRangeMaps
void setDomainRangeMaps(const Teuchos::RCP< const map_type > &domainMap, const Teuchos::RCP< const map_type > &rangeMap)
Definition: Tpetra_CrsGraph_def.hpp:1847
Tpetra::CrsGraph::getNumAllocatedEntriesInLocalRow
size_t getNumAllocatedEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of allocated entries on this node in the specified local row...
Definition: Tpetra_CrsGraph_def.hpp:2198
Tpetra::CrsGraph::getNumEntriesInLocalRow
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
Definition: Tpetra_CrsGraph_def.hpp:2169
Tpetra::CrsGraph::storageStatus_
Details::EStorageStatus storageStatus_
Status of the graph&#39;s storage, when not in a fill-complete state.
Definition: Tpetra_CrsGraph_decl.hpp:1932
Tpetra::CrsGraph::getNodeNumDiags
size_t getNodeNumDiags() const
Returns the number of local diagonal entries, based on global row/column index comparisons.
Definition: Tpetra_CrsGraph_def.hpp:698
Tpetra::CrsGraph::execution_space
device_type::execution_space execution_space
This class&#39; Kokkos execution space.
Definition: Tpetra_CrsGraph_decl.hpp:261
Tpetra::CrsGraph::k_rowPtrs_
local_graph_type::row_map_type::const_type k_rowPtrs_
Row offsets for "1-D" storage.
Definition: Tpetra_CrsGraph_decl.hpp:1869
Tpetra::CrsGraph::getNodeNumRows
size_t getNodeNumRows() const
Returns the number of graph rows owned on the calling node.
Definition: Tpetra_CrsGraph_def.hpp:671
Tpetra::CrsGraph::getGlobalRowView
void getGlobalRowView(GlobalOrdinal GlobalRow, Teuchos::ArrayView< const GlobalOrdinal > &Indices) const
Extract a const, non-persisting view of global indices in a specified row of the graph.
Definition: Tpetra_CrsGraph_def.hpp:2475
Tpetra::CrsGraph::nonlocals_
std::map< GlobalOrdinal, std::vector< GlobalOrdinal > > nonlocals_
Nonlocal data given to insertGlobalValues or sumIntoGlobalValues.
Definition: Tpetra_CrsGraph_decl.hpp:1954
Tpetra::DistObject< GlobalOrdinal, LocalOrdinal, GlobalOrdinal, Node >::getMap
virtual Teuchos::RCP< const map_type > getMap() const
The Map describing the parallel distribution of this object.
Definition: Tpetra_DistObject_decl.hpp:347
Tpetra::CrsGraph::replaceDomainMapAndImporter
void replaceDomainMapAndImporter(const Teuchos::RCP< const map_type > &newDomainMap, const Teuchos::RCP< const import_type > &newImporter)
Replace the current domain Map and Import with the given parameters.
Definition: Tpetra_CrsGraph_def.hpp:4250
Tpetra::CrsGraph::mergeAllIndices
void mergeAllIndices()
Merge duplicate row indices in all of the rows.
Definition: Tpetra_CrsGraph_def.hpp:4896
Tpetra::CrsGraph::getLocalView
Teuchos::ArrayView< const LocalOrdinal > getLocalView(const RowInfo rowinfo) const
Get a const, nonowned, locally indexed view of the locally owned row myRow, such that rowinfo = getRo...
Definition: Tpetra_CrsGraph_def.hpp:1191
Tpetra::CrsGraph::setLocallyModified
void setLocallyModified()
Report that we made a local modification to its structure.
Definition: Tpetra_CrsGraph_def.hpp:965
Tpetra::CrsGraph::~CrsGraph
virtual ~CrsGraph()
Destructor.
Definition: Tpetra_CrsGraph_def.hpp:592
Tpetra::CrsGraph::t_GlobalOrdinal_1D
Kokkos::View< GlobalOrdinal *, execution_space > t_GlobalOrdinal_1D
Type of the k_gblInds1D_ array of global column indices.
Definition: Tpetra_CrsGraph_decl.hpp:1831
Tpetra::CrsGraph::getValidParameters
Teuchos::RCP< const ParameterList > getValidParameters() const
Default parameter list suitable for validation.
Definition: Tpetra_CrsGraph_def.hpp:599
Tpetra::CrsGraph::sortAllIndices
void sortAllIndices()
Sort the column indices in all the rows.
Definition: Tpetra_CrsGraph_def.hpp:4545
Tpetra::CrsGraph::isLowerTriangular
bool isLowerTriangular() const
Whether the graph is locally lower triangular.
Definition: Tpetra_CrsGraph_def.hpp:884
Tpetra::CrsGraph::sortRowIndices
void sortRowIndices(const RowInfo rowinfo)
Sort the column indices in the given row.
Definition: Tpetra_CrsGraph_def.hpp:1717
Tpetra::CrsGraph::insertLocalIndices
void insertLocalIndices(const LocalOrdinal localRow, const Teuchos::ArrayView< const LocalOrdinal > &indices)
Insert local indices into the graph.
Definition: Tpetra_CrsGraph_def.hpp:2512
Tpetra::CrsGraph::exporter_
Teuchos::RCP< const export_type > exporter_
The Export from the row Map to the range Map.
Definition: Tpetra_CrsGraph_decl.hpp:1774
Tpetra::CrsGraph::describe
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object to the given output stream with given verbosity level.
Definition: Tpetra_CrsGraph_def.hpp:4987
Tpetra::CrsGraph::getGlobalView
Teuchos::ArrayView< const GlobalOrdinal > getGlobalView(const RowInfo rowinfo) const
Get a const, nonowned, globally indexed view of the locally owned row myRow, such that rowinfo = getR...
Definition: Tpetra_CrsGraph_def.hpp:1257
Tpetra::CrsGraph::lowerTriangular_
bool lowerTriangular_
Whether the graph is locally lower triangular.
Definition: Tpetra_CrsGraph_decl.hpp:1940
Tpetra::CrsGraph::getExporter
Teuchos::RCP< const export_type > getExporter() const
Returns the exporter associated with this graph.
Definition: Tpetra_CrsGraph_def.hpp:769
Tpetra::CrsGraph::getNumEntriesPerLocalRowUpperBound
void getNumEntriesPerLocalRowUpperBound(Teuchos::ArrayRCP< const size_t > &boundPerLocalRow, size_t &boundForAllLocalRows, bool &boundSameForAllLocalRows) const
Get an upper bound on the number of entries that can be stored in each row.
Definition: Tpetra_CrsGraph_def.hpp:2937
Tpetra::CrsGraph::checkSizes
virtual bool checkSizes(const SrcDistObject &source)
Compare the source and target (this) objects for compatibility.
Definition: Tpetra_CrsGraph_def.hpp:5107
Tpetra::CrsGraph::getNodeRowPtrs
Teuchos::ArrayRCP< const size_t > getNodeRowPtrs() const
Get a host view of the row offsets.
Definition: Tpetra_CrsGraph_def.hpp:2238
Tpetra::CrsGraph::getRowInfo
RowInfo getRowInfo(const size_t myRow) const
Get information about the locally owned row with local index myRow.
Definition: Tpetra_CrsGraph_def.hpp:1295
Tpetra::CrsGraph::getLocalRowView
void getLocalRowView(LocalOrdinal LocalRow, Teuchos::ArrayView< const LocalOrdinal > &indices) const
Extract a const, non-persisting view of local indices in a specified row of the graph.
Definition: Tpetra_CrsGraph_def.hpp:2440
Tpetra::CrsGraph::getNumEntriesInGlobalRow
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries on this node in the specified global row.
Definition: Tpetra_CrsGraph_def.hpp:2153
Generated by   doxygen 1.8.11