Zoltan2
TpetraRowMatrixInput.cpp
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 //
46 // Basic testing of Zoltan2::TpetraRowMatrixAdapter
47 
53 #include <string>
54 
56 #include <Zoltan2_InputTraits.hpp>
57 #include <Zoltan2_TestHelpers.hpp>
58 
59 #include <Teuchos_GlobalMPISession.hpp>
60 #include <Teuchos_DefaultComm.hpp>
61 #include <Teuchos_RCP.hpp>
62 #include <Teuchos_Comm.hpp>
63 #include <Teuchos_CommHelpers.hpp>
64 
65 using namespace std;
66 using Teuchos::RCP;
67 using Teuchos::rcp;
68 using Teuchos::rcp_const_cast;
69 using Teuchos::rcp_dynamic_cast;
70 using Teuchos::Comm;
71 using Teuchos::DefaultComm;
72 
73 typedef Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> ztcrsmatrix_t;
74 typedef Tpetra::RowMatrix<zscalar_t, zlno_t, zgno_t, znode_t> ztrowmatrix_t;
75 
77 
78 void printMatrix(RCP<const Comm<int> > &comm, zlno_t nrows,
79  const zgno_t *rowIds, const zlno_t *offsets, const zgno_t *colIds)
80 {
81  int rank = comm->getRank();
82  int nprocs = comm->getSize();
83  comm->barrier();
84  for (int p=0; p < nprocs; p++){
85  if (p == rank){
86  std::cout << rank << ":" << std::endl;
87  for (zlno_t i=0; i < nrows; i++){
88  std::cout << " row " << rowIds[i] << ": ";
89  for (zlno_t j=offsets[i]; j < offsets[i+1]; j++){
90  std::cout << colIds[j] << " ";
91  }
92  std::cout << std::endl;
93  }
94  std::cout.flush();
95  }
96  comm->barrier();
97  }
98  comm->barrier();
99 }
100 
102 
103 template <typename User>
106 {
107  RCP<const Comm<int> > comm = M.getComm();
108  int fail = 0, gfail=0;
109 
110  if (!fail && ia.getLocalNumRows() != M.getNodeNumRows())
111  fail = 4;
112 
113  if (M.getNodeNumRows()){
114  if (!fail && ia.getLocalNumColumns() != M.getNodeNumCols())
115  fail = 6;
116  }
117 
118  gfail = globalFail(comm, fail);
119 
120  const zgno_t *rowIds=NULL, *colIds=NULL;
121  const zlno_t *offsets=NULL;
122  size_t nrows=0;
123 
124  if (!gfail){
125 
126  nrows = ia.getLocalNumRows();
127  ia.getRowIDsView(rowIds);
128  ia.getCRSView(offsets, colIds);
129 
130  if (nrows != M.getNodeNumRows())
131  fail = 8;
132 
133  gfail = globalFail(comm, fail);
134 
135  if (gfail == 0){
136  printMatrix(comm, nrows, rowIds, offsets, colIds);
137  }
138  else{
139  if (!fail) fail = 10;
140  }
141  }
142  return fail;
143 }
144 
145 
146 int main(int argc, char *argv[])
147 {
148  Teuchos::GlobalMPISession session(&argc, &argv);
149  RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
150  int rank = comm->getRank();
151  int fail = 0, gfail=0;
152  bool aok = true;
153 
154  // Create object that can give us Tpetra matrices for testing.
155 
156  RCP<UserInputForTests> uinput;
157 
158  try{
159  uinput =
160  rcp(new UserInputForTests(
161  testDataFilePath,std::string("simple"), comm, true));
162  }
163  catch(std::exception &e){
164  aok = false;
165  std::cout << e.what() << std::endl;
166  }
167  TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
168 
169  // Input matrix and row matrix cast from it.
170  RCP<ztcrsmatrix_t> tM = uinput->getUITpetraCrsMatrix();
171  RCP<ztrowmatrix_t> trM = rcp_dynamic_cast<ztrowmatrix_t>(tM);
172 
173  RCP<ztrowmatrix_t> newM; // migrated matrix
174 
175  size_t nrows = trM->getNodeNumRows();
176 
177  // To test migration in the input adapter we need a Solution object.
178 
179  RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
180 
181  int nWeights = 1;
182 
185  typedef adapter_t::part_t part_t;
186 
187  part_t *p = new part_t [nrows];
188  memset(p, 0, sizeof(part_t) * nrows);
189  ArrayRCP<part_t> solnParts(p, 0, nrows, true);
190 
191  soln_t solution(env, comm, nWeights);
192  solution.setParts(solnParts);
193 
195  // User object is Tpetra::RowMatrix
196  if (!gfail){
197  if (rank==0)
198  std::cout << "Input adapter for Tpetra::RowMatrix" << std::endl;
199 
200  RCP<const ztrowmatrix_t> ctrM = rcp_const_cast<const ztrowmatrix_t>(
201  rcp_dynamic_cast<ztrowmatrix_t>(tM));
202  RCP<adapter_t> trMInput;
203 
204  try {
205  trMInput = rcp(new adapter_t(ctrM));
206  }
207  catch (std::exception &e){
208  aok = false;
209  std::cout << e.what() << std::endl;
210  }
211  TEST_FAIL_AND_EXIT(*comm, aok, "TpetraRowMatrixAdapter ", 1);
212 
213  fail = verifyInputAdapter<ztrowmatrix_t>(*trMInput, *trM);
214 
215  gfail = globalFail(comm, fail);
216 
217  if (!gfail){
218  ztrowmatrix_t *mMigrate = NULL;
219  try{
220  trMInput->applyPartitioningSolution(*trM, mMigrate, solution);
221  newM = rcp(mMigrate);
222  }
223  catch (std::exception &e){
224  fail = 11;
225  cout << "Error caught: " << e.what() << endl;
226  }
227 
228  gfail = globalFail(comm, fail);
229 
230  if (!gfail){
231  RCP<const ztrowmatrix_t> cnewM =
232  rcp_const_cast<const ztrowmatrix_t>(newM);
233  RCP<adapter_t> newInput;
234  try{
235  newInput = rcp(new adapter_t(cnewM));
236  }
237  catch (std::exception &e){
238  aok = false;
239  std::cout << e.what() << std::endl;
240  }
241  TEST_FAIL_AND_EXIT(*comm, aok, "TpetraRowMatrixAdapter 2 ", 1);
242 
243  if (rank==0){
244  std::cout <<
245  "Input adapter for Tpetra::RowMatrix migrated to proc 0" <<
246  std::endl;
247  }
248  fail = verifyInputAdapter<ztrowmatrix_t>(*newInput, *newM);
249  if (fail) fail += 100;
250  gfail = globalFail(comm, fail);
251  }
252  }
253  if (gfail){
254  printFailureCode(comm, fail);
255  }
256  }
257 
259  // DONE
260 
261  if (rank==0)
262  std::cout << "PASS" << std::endl;
263 }
int verifyInputAdapter(Zoltan2::TpetraRowMatrixAdapter< User > &ia, ztrowmatrix_t &M)
size_t getLocalNumColumns() const
Returns the number of columns on this process.
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > ztcrsmatrix_t
int globalFail(const RCP< const Comm< int > > &comm, int fail)
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
int zlno_t
common code used by tests
Provides access for Zoltan2 to Tpetra::RowMatrix data.
size_t getLocalNumRows() const
Returns the number of rows on this process.
void printMatrix(RCP< const Comm< int > > &comm, zlno_t nrows, const zgno_t *rowIds, const zlno_t *offsets, const zgno_t *colIds)
A PartitioningSolution is a solution to a partitioning problem.
Traits for application input objects.
Defines the TpetraRowMatrixAdapter class.
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
static const std::string fail
int zgno_t
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.
void printFailureCode(const RCP< const Comm< int > > &comm, int fail)
Tpetra::RowMatrix< zscalar_t, zlno_t, zgno_t, znode_t > ztrowmatrix_t
int main(int argc, char *argv[])
void getRowIDsView(const gno_t *&rowIds) const
Sets pointer to this process&#39; rows&#39; global IDs.
std::string testDataFilePath(".")