Zoltan2
Zoltan2_TpetraRowMatrixAdapter.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Zoltan2: A package of combinatorial algorithms for scientific computing
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Karen Devine (kddevin@sandia.gov)
39 // Erik Boman (egboman@sandia.gov)
40 // Siva Rajamanickam (srajama@sandia.gov)
41 //
42 // ***********************************************************************
43 //
44 // @HEADER
45 
50 #ifndef _ZOLTAN2_TPETRAROWMATRIXADAPTER_HPP_
51 #define _ZOLTAN2_TPETRAROWMATRIXADAPTER_HPP_
52 
54 #include <Zoltan2_StridedData.hpp>
56 
57 #include <Tpetra_RowMatrix.hpp>
58 
59 namespace Zoltan2 {
60 
62 
75 template <typename User, typename UserCoord=User>
76  class TpetraRowMatrixAdapter : public MatrixAdapter<User,UserCoord> {
77 public:
78 
79 #ifndef DOXYGEN_SHOULD_SKIP_THIS
80  typedef typename InputTraits<User>::scalar_t scalar_t;
81  typedef typename InputTraits<User>::lno_t lno_t;
82  typedef typename InputTraits<User>::gno_t gno_t;
83  typedef typename InputTraits<User>::part_t part_t;
84  typedef typename InputTraits<User>::node_t node_t;
85  typedef User user_t;
86  typedef UserCoord userCoord_t;
87 #endif
88 
92 
98  TpetraRowMatrixAdapter(const RCP<const User> &inmatrix,
99  int nWeightsPerRow=0);
100 
113  void setWeights(const scalar_t *weightVal, int stride, int idx = 0);
114 
130  void setRowWeights(const scalar_t *weightVal, int stride, int idx = 0);
131 
137  void setWeightIsDegree(int idx);
138 
144  void setRowWeightIsNumberOfNonZeros(int idx);
145 
147  // The MatrixAdapter interface.
149 
150  size_t getLocalNumRows() const {
151  return matrix_->getNodeNumRows();
152  }
153 
154  size_t getLocalNumColumns() const {
155  return matrix_->getNodeNumCols();
156  }
157 
158  size_t getLocalNumEntries() const {
159  return matrix_->getNodeNumEntries();
160  }
161 
162  bool CRSViewAvailable() const { return true; }
163 
164  void getRowIDsView(const gno_t *&rowIds) const
165  {
166  ArrayView<const gno_t> rowView = rowMap_->getNodeElementList();
167  rowIds = rowView.getRawPtr();
168  }
169 
170  void getCRSView(const lno_t *&offsets, const gno_t *&colIds) const
171  {
172  offsets = offset_.getRawPtr();
173  colIds = columnIds_.getRawPtr();
174  }
175 
176  void getCRSView(const lno_t *&offsets, const gno_t *&colIds,
177  const scalar_t *&values) const
178  {
179  offsets = offset_.getRawPtr();
180  colIds = columnIds_.getRawPtr();
181  values = values_.getRawPtr();
182  }
183 
184 
185  int getNumWeightsPerRow() const { return nWeightsPerRow_; }
186 
187  void getRowWeightsView(const scalar_t *&weights, int &stride,
188  int idx = 0) const
189  {
190  env_->localInputAssertion(__FILE__, __LINE__,
191  "invalid weight index",
192  idx >= 0 && idx < nWeightsPerRow_, BASIC_ASSERTION);
193  size_t length;
194  rowWeights_[idx].getStridedList(length, weights, stride);
195  }
196 
197  bool useNumNonzerosAsRowWeight(int idx) const { return numNzWeight_[idx];}
198 
199  template <typename Adapter>
200  void applyPartitioningSolution(const User &in, User *&out,
201  const PartitioningSolution<Adapter> &solution) const;
202 
203  template <typename Adapter>
204  void applyPartitioningSolution(const User &in, RCP<User> &out,
205  const PartitioningSolution<Adapter> &solution) const;
206 
207 private:
208 
209  RCP<Environment> env_; // for error messages, etc.
210 
211  RCP<const User> matrix_;
212  RCP<const Tpetra::Map<lno_t, gno_t, node_t> > rowMap_;
213  RCP<const Tpetra::Map<lno_t, gno_t, node_t> > colMap_;
214  lno_t base_;
215  ArrayRCP<lno_t> offset_;
216  ArrayRCP<gno_t> columnIds_; // TODO: KDD Is it necessary to copy and store
217  ArrayRCP<scalar_t> values_; // TODO: the matrix here? Would prefer views.
218 
219  int nWeightsPerRow_;
220  ArrayRCP<StridedData<lno_t, scalar_t> > rowWeights_;
221  ArrayRCP<bool> numNzWeight_;
222 
223  bool mayHaveDiagonalEntries;
224 
225  RCP<User> doMigration(const User &from, size_t numLocalRows,
226  const gno_t *myNewRows) const;
227 };
228 
230 // Definitions
232 
233 template <typename User, typename UserCoord>
235  const RCP<const User> &inmatrix, int nWeightsPerRow):
236  env_(rcp(new Environment)),
237  matrix_(inmatrix), rowMap_(), colMap_(), base_(),
238  offset_(), columnIds_(),
239  nWeightsPerRow_(nWeightsPerRow), rowWeights_(), numNzWeight_(),
240  mayHaveDiagonalEntries(true)
241 {
242  typedef StridedData<lno_t,scalar_t> input_t;
243 
244  rowMap_ = matrix_->getRowMap();
245  colMap_ = matrix_->getColMap();
246  base_ = rowMap_->getIndexBase();
247 
248  size_t nrows = matrix_->getNodeNumRows();
249  size_t nnz = matrix_->getNodeNumEntries();
250  size_t maxnumentries =
251  matrix_->getNodeMaxNumRowEntries(); // Diff from CrsMatrix
252 
253  offset_.resize(nrows+1, 0);
254  columnIds_.resize(nnz);
255  values_.resize(nnz);
256  ArrayRCP<lno_t> indices(maxnumentries); // Diff from CrsMatrix
257  ArrayRCP<scalar_t> nzs(maxnumentries); // Diff from CrsMatrix
258  lno_t next = 0;
259  for (size_t i=0; i < nrows; i++){
260  lno_t row = i + base_;
261  matrix_->getLocalRowCopy(row, indices(), nzs(), nnz); // Diff from CrsMatrix
262  for (size_t j=0; j < nnz; j++){
263  values_[next] = nzs[j];
264  // TODO - this will be slow
265  // Is it possible that global columns ids might be stored in order?
266  columnIds_[next++] = colMap_->getGlobalElement(indices[j]);
267  }
268  offset_[i+1] = offset_[i] + nnz;
269  }
270 
271  if (nWeightsPerRow_ > 0){
272  rowWeights_ = arcp(new input_t [nWeightsPerRow_], 0, nWeightsPerRow_, true);
273  numNzWeight_ = arcp(new bool [nWeightsPerRow_], 0, nWeightsPerRow_, true);
274  for (int i=0; i < nWeightsPerRow_; i++)
275  numNzWeight_[i] = false;
276  }
277 }
278 
280 template <typename User, typename UserCoord>
282  const scalar_t *weightVal, int stride, int idx)
283 {
284  if (this->getPrimaryEntityType() == MATRIX_ROW)
285  setRowWeights(weightVal, stride, idx);
286  else {
287  // TODO: Need to allow weights for columns and/or nonzeros
288  std::ostringstream emsg;
289  emsg << __FILE__ << "," << __LINE__
290  << " error: setWeights not yet supported for"
291  << " columns or nonzeros."
292  << std::endl;
293  throw std::runtime_error(emsg.str());
294  }
295 }
296 
298 template <typename User, typename UserCoord>
300  const scalar_t *weightVal, int stride, int idx)
301 {
302  typedef StridedData<lno_t,scalar_t> input_t;
303  env_->localInputAssertion(__FILE__, __LINE__,
304  "invalid row weight index",
305  idx >= 0 && idx < nWeightsPerRow_, BASIC_ASSERTION);
306  size_t nvtx = getLocalNumRows();
307  ArrayRCP<const scalar_t> weightV(weightVal, 0, nvtx*stride, false);
308  rowWeights_[idx] = input_t(weightV, stride);
309 }
310 
312 template <typename User, typename UserCoord>
314  int idx)
315 {
316  if (this->getPrimaryEntityType() == MATRIX_ROW)
317  setRowWeightIsNumberOfNonZeros(idx);
318  else {
319  // TODO: Need to allow weights for columns and/or nonzeros
320  std::ostringstream emsg;
321  emsg << __FILE__ << "," << __LINE__
322  << " error: setWeightIsNumberOfNonZeros not yet supported for"
323  << " columns" << std::endl;
324  throw std::runtime_error(emsg.str());
325  }
326 }
327 
329 template <typename User, typename UserCoord>
331  int idx)
332 {
333  env_->localInputAssertion(__FILE__, __LINE__,
334  "invalid row weight index",
335  idx >= 0 && idx < nWeightsPerRow_, BASIC_ASSERTION);
336 
337  numNzWeight_[idx] = true;
338 }
339 
341 template <typename User, typename UserCoord>
342  template <typename Adapter>
344  const User &in, User *&out,
345  const PartitioningSolution<Adapter> &solution) const
346 {
347  // Get an import list (rows to be received)
348  size_t numNewRows;
349  ArrayRCP<gno_t> importList;
350  try{
351  numNewRows = Zoltan2::getImportList<Adapter,
353  (solution, this, importList);
354  }
356 
357  // Move the rows, creating a new matrix.
358  RCP<User> outPtr = doMigration(in, numNewRows, importList.getRawPtr());
359  out = outPtr.get();
360  outPtr.release();
361 }
362 
364 template <typename User, typename UserCoord>
365  template <typename Adapter>
367  const User &in, RCP<User> &out,
368  const PartitioningSolution<Adapter> &solution) const
369 {
370  // Get an import list (rows to be received)
371  size_t numNewRows;
372  ArrayRCP<gno_t> importList;
373  try{
374  numNewRows = Zoltan2::getImportList<Adapter,
376  (solution, this, importList);
377  }
379 
380  // Move the rows, creating a new matrix.
381  out = doMigration(in, numNewRows, importList.getRawPtr());
382 }
383 
384 
386 template < typename User, typename UserCoord>
388  const User &from,
389  size_t numLocalRows,
390  const gno_t *myNewRows
391 ) const
392 {
393  typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
394  typedef Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> tcrsmatrix_t;
395 
396  // We cannot create a Tpetra::RowMatrix, unless the underlying type is
397  // something we know (like Tpetra::CrsMatrix).
398  // If the underlying type is something different, the user probably doesn't
399  // want a Tpetra::CrsMatrix back, so we throw an error.
400 
401  // Try to cast "from" matrix to a TPetra::CrsMatrix
402  // If that fails we throw an error.
403  // We could cast as a ref which will throw std::bad_cast but with ptr
404  // approach it might be clearer what's going on here
405  const tcrsmatrix_t *pCrsMatrix = dynamic_cast<const tcrsmatrix_t *>(&from);
406 
407  if(!pCrsMatrix) {
408  throw std::logic_error("TpetraRowMatrixAdapter cannot migrate data for "
409  "your RowMatrix; it can migrate data only for "
410  "Tpetra::CrsMatrix. "
411  "You can inherit from TpetraRowMatrixAdapter and "
412  "implement migration for your RowMatrix.");
413  }
414 
415  lno_t base = 0;
416 
417  // source map
418  const RCP<const map_t> &smap = from.getRowMap();
419  gno_t numGlobalRows = smap->getGlobalNumElements();
420 
421  // target map
422  ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
423  const RCP<const Teuchos::Comm<int> > &comm = from.getComm();
424  RCP<const map_t> tmap = rcp(new map_t(numGlobalRows, rowList, base, comm));
425 
426  // importer
427  Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
428 
429  // target matrix
430  // Chris Siefert proposed using the following to make migration
431  // more efficient.
432  // By default, the Domain and Range maps are the same as in "from".
433  // As in the original code, we instead set them both to tmap.
434  // The assumption is a square matrix.
435  // TODO: what about rectangular matrices?
436  // TODO: Should choice of domain/range maps be an option to this function?
437 
438  // KDD 3/7/16: disabling Chris' new code to avoid dashboard failures;
439  // KDD 3/7/16: can re-enable when issue #114 is fixed.
440  // KDD 3/7/16: when re-enable CSIEFERT code, can comment out
441  // KDD 3/7/16: "Original way" code.
442  // CSIEFERT RCP<tcrsmatrix_t> M;
443  // CSIEFERT from.importAndFillComplete(M, importer, tmap, tmap);
444 
445  // Original way we did it:
446  //
447  int oldNumElts = smap->getNodeNumElements();
448  int newNumElts = numLocalRows;
449 
450  // number of non zeros in my new rows
451  typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> vector_t;
452  vector_t numOld(smap); // TODO These vectors should have scalar=size_t,
453  vector_t numNew(tmap); // but ETI does not yet support that.
454  for (int lid=0; lid < oldNumElts; lid++){
455  numOld.replaceGlobalValue(smap->getGlobalElement(lid),
456  scalar_t(from.getNumEntriesInLocalRow(lid)));
457  }
458  numNew.doImport(numOld, importer, Tpetra::INSERT);
459 
460  // TODO Could skip this copy if could declare vector with scalar=size_t.
461  ArrayRCP<size_t> nnz(newNumElts);
462  if (newNumElts > 0){
463  ArrayRCP<scalar_t> ptr = numNew.getDataNonConst(0);
464  for (int lid=0; lid < newNumElts; lid++){
465  nnz[lid] = static_cast<size_t>(ptr[lid]);
466  }
467  }
468 
469  RCP<tcrsmatrix_t> M = rcp(new tcrsmatrix_t(tmap, nnz,
470  Tpetra::StaticProfile));
471  M->doImport(from, importer, Tpetra::INSERT);
472  M->fillComplete();
473 
474  // End of original way we did it.
475  return Teuchos::rcp_dynamic_cast<User>(M);
476 }
477 
478 } //namespace Zoltan2
479 
480 #endif
InputTraits< User >::scalar_t scalar_t
size_t getLocalNumColumns() const
Returns the number of columns on this process.
Helper functions for Partitioning Problems.
fast typical checks for valid arguments
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
MatrixAdapter defines the adapter interface for matrices.
InputTraits< User >::gno_t gno_t
bool useNumNonzerosAsRowWeight(int idx) const
Indicate whether row weight with index idx should be the global number of nonzeros in the row...
default_part_t part_t
The data type to represent part numbers.
InputTraits< User >::lno_t lno_t
int getNumWeightsPerRow() const
Returns the number of weights per row (0 or greater). Row weights may be used when partitioning matri...
static ArrayRCP< ArrayRCP< zscalar_t > > weights
size_t getImportList(const PartitioningSolution< SolutionAdapter > &solution, const DataAdapter *const data, ArrayRCP< typename DataAdapter::gno_t > &imports)
From a PartitioningSolution, get a list of IDs to be imported. Assumes part numbers in PartitioningSo...
Provides access for Zoltan2 to Tpetra::RowMatrix data.
size_t getLocalNumRows() const
Returns the number of rows on this process.
void getCRSView(const lno_t *&offsets, const gno_t *&colIds, const scalar_t *&values) const
A PartitioningSolution is a solution to a partitioning problem.
size_t getLocalNumEntries() const
Returns the number of nonzeros on this process.
default_lno_t lno_t
The ordinal type (e.g., int, long, int64_t) that represents local counts and local indices...
The StridedData class manages lists of weights or coordinates.
void getRowWeightsView(const scalar_t *&weights, int &stride, int idx=0) const
InputTraits< User >::part_t part_t
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
void setRowWeights(const scalar_t *weightVal, int stride, int idx=0)
Specify a weight for each row.
default_gno_t gno_t
The ordinal type (e.g., int, long, int64_t) that can represent global counts and identifiers.
default_node_t node_t
The Kokkos node type. This is only meaningful for users of Tpetra objects.
void getCRSView(const lno_t *&offsets, const gno_t *&colIds) const
Sets pointers to this process&#39; matrix entries using compressed sparse row (CRS) format. All matrix adapters must implement either getCRSView or getCCSView, but implementation of both is not required.
bool CRSViewAvailable() const
Indicates whether the MatrixAdapter implements a view of the matrix in compressed sparse row (CRS) fo...
void setWeightIsDegree(int idx)
Specify an index for which the weight should be the degree of the entity.
Defines the MatrixAdapter interface.
void setWeights(const scalar_t *weightVal, int stride, int idx=0)
Specify a weight for each entity of the primaryEntityType.
void applyPartitioningSolution(const User &in, User *&out, const PartitioningSolution< Adapter > &solution) const
void setRowWeightIsNumberOfNonZeros(int idx)
Specify an index for which the row weight should be the global number of nonzeros in the row...
default_scalar_t scalar_t
The data type for weights and coordinates.
void getRowIDsView(const gno_t *&rowIds) const
Sets pointer to this process&#39; rows&#39; global IDs.
TpetraRowMatrixAdapter(const RCP< const User > &inmatrix, int nWeightsPerRow=0)
Constructor
This file defines the StridedData class.