Zoltan2
Loading...
Searching...
No Matches
XpetraCrsMatrixInput.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::XpetraCrsMatrixAdapter
47
53#include <string>
54
58
59#include <Teuchos_DefaultComm.hpp>
60#include <Teuchos_RCP.hpp>
61#include <Teuchos_Comm.hpp>
62#include <Teuchos_CommHelpers.hpp>
63
64using Teuchos::RCP;
65using Teuchos::rcp;
66using Teuchos::rcp_const_cast;
67using Teuchos::Comm;
68
69typedef Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> tmatrix_t;
70typedef Xpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> xmatrix_t;
71
72template<typename offset_t>
73void printMatrix(RCP<const Comm<int> > &comm, zlno_t nrows,
74 const zgno_t *rowIds, const offset_t *offsets, const zgno_t *colIds)
75{
76 int rank = comm->getRank();
77 int nprocs = comm->getSize();
78 comm->barrier();
79 for (int p=0; p < nprocs; p++){
80 if (p == rank){
81 std::cout << rank << ":" << std::endl;
82 for (zlno_t i=0; i < nrows; i++){
83 std::cout << " row " << rowIds[i] << ": ";
84 for (offset_t j=offsets[i]; j < offsets[i+1]; j++){
85 std::cout << colIds[j] << " ";
86 }
87 std::cout << std::endl;
88 }
89 std::cout.flush();
90 }
91 comm->barrier();
92 }
93 comm->barrier();
94}
95
96template <typename User>
99{
100 typedef typename Zoltan2::InputTraits<User>::offset_t offset_t;
101
102 RCP<const Comm<int> > comm = M.getComm();
103 int fail = 0, gfail=0;
104
105 if (!fail && ia.getLocalNumRows() != M.getLocalNumRows())
106 fail = 4;
107
108 if (M.getLocalNumRows()){
109 if (!fail && ia.getLocalNumColumns() != M.getLocalNumCols())
110 fail = 6;
111 }
112
113 gfail = globalFail(*comm, fail);
114
115 const zgno_t *rowIds=NULL;
116 ArrayRCP<const zgno_t> colIds;
117 ArrayRCP<const offset_t> offsets;
118 size_t nrows=0;
119
120 if (!gfail){
121
122 nrows = ia.getLocalNumRows();
123 ia.getRowIDsView(rowIds);
124 ia.getCRSView(offsets, colIds);
125
126 if (nrows != M.getLocalNumRows())
127 fail = 8;
128
129 gfail = globalFail(*comm, fail);
130
131 if (gfail == 0){
132 printMatrix<offset_t>(comm, nrows, rowIds, offsets.getRawPtr(), colIds.getRawPtr());
133 }
134 else{
135 if (!fail) fail = 10;
136 }
137 }
138 return fail;
139}
140
141int main(int narg, char *arg[])
142{
143 Tpetra::ScopeGuard tscope(&narg, &arg);
144 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
145
146 int rank = comm->getRank();
147 int fail = 0, gfail=0;
148 bool aok = true;
149
150 // Create object that can give us test Tpetra, Xpetra
151 // and Epetra matrices for testing.
152
153 RCP<UserInputForTests> uinput;
154 Teuchos::ParameterList params;
155 params.set("input file", "simple");
156 params.set("file type", "Chaco");
157
158 try{
159 uinput = rcp(new UserInputForTests(params, comm));
160 }
161 catch(std::exception &e){
162 aok = false;
163 std::cout << e.what() << std::endl;
164 }
165 TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
166
167 RCP<tmatrix_t> tM; // original matrix (for checking)
168 RCP<tmatrix_t> newM; // migrated matrix
169
170 tM = uinput->getUITpetraCrsMatrix();
171 size_t nrows = tM->getLocalNumRows();
172
173 // To test migration in the input adapter we need a Solution object.
174
175 RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment(comm));
176
177 int nWeights = 1;
178
181 typedef adapter_t::part_t part_t;
182
183 part_t *p = new part_t [nrows];
184 memset(p, 0, sizeof(part_t) * nrows);
185 ArrayRCP<part_t> solnParts(p, 0, nrows, true);
186
187 soln_t solution(env, comm, nWeights);
188 solution.setParts(solnParts);
189
191 // User object is Tpetra::CrsMatrix
192 if (!gfail){
193 if (rank==0)
194 std::cout << "Input adapter for Tpetra::CrsMatrix" << std::endl;
195
196 RCP<const tmatrix_t> ctM = rcp_const_cast<const tmatrix_t>(tM);
197 RCP<Zoltan2::XpetraCrsMatrixAdapter<tmatrix_t> > tMInput;
198
199 try {
200 tMInput =
202 }
203 catch (std::exception &e){
204 aok = false;
205 std::cout << e.what() << std::endl;
206 }
207 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter ", 1);
208
209 fail = verifyInputAdapter<tmatrix_t>(*tMInput, *tM);
210
211 gfail = globalFail(*comm, fail);
212
213 if (!gfail){
214 tmatrix_t *mMigrate = NULL;
215 try{
216 tMInput->applyPartitioningSolution(*tM, mMigrate, solution);
217 newM = rcp(mMigrate);
218 }
219 catch (std::exception &e){
220 fail = 11;
221 std::cout << "Error caught: " << e.what() << std::endl;
222 }
223
224 gfail = globalFail(*comm, fail);
225
226 if (!gfail){
227 RCP<const tmatrix_t> cnewM = rcp_const_cast<const tmatrix_t>(newM);
228 RCP<Zoltan2::XpetraCrsMatrixAdapter<tmatrix_t> > newInput;
229 try{
230 newInput = rcp(new Zoltan2::XpetraCrsMatrixAdapter<tmatrix_t>(cnewM));
231 }
232 catch (std::exception &e){
233 aok = false;
234 std::cout << e.what() << std::endl;
235 }
236 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 2 ", 1);
237
238 if (rank==0){
239 std::cout <<
240 "Input adapter for Tpetra::CrsMatrix migrated to proc 0" <<
241 std::endl;
242 }
243 fail = verifyInputAdapter<tmatrix_t>(*newInput, *newM);
244 if (fail) fail += 100;
245 gfail = globalFail(*comm, fail);
246 }
247 }
248 if (gfail){
249 printFailureCode(*comm, fail);
250 }
251 }
252
254 // User object is Xpetra::CrsMatrix
255 if (!gfail){
256 if (rank==0)
257 std::cout << "Input adapter for Xpetra::CrsMatrix" << std::endl;
258
259 RCP<xmatrix_t> xM = uinput->getUIXpetraCrsMatrix();
260 RCP<const xmatrix_t> cxM = rcp_const_cast<const xmatrix_t>(xM);
261 RCP<Zoltan2::XpetraCrsMatrixAdapter<xmatrix_t> > xMInput;
262
263 try {
264 xMInput =
266 }
267 catch (std::exception &e){
268 aok = false;
269 std::cout << e.what() << std::endl;
270 }
271 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 3 ", 1);
272
273 fail = verifyInputAdapter<xmatrix_t>(*xMInput, *tM);
274
275 gfail = globalFail(*comm, fail);
276
277 if (!gfail){
278 xmatrix_t *mMigrate =NULL;
279 try{
280 xMInput->applyPartitioningSolution(*xM, mMigrate, solution);
281 }
282 catch (std::exception &e){
283 std::cout << "Error caught: " << e.what() << std::endl;
284 fail = 11;
285 }
286
287 gfail = globalFail(*comm, fail);
288
289 if (!gfail){
290 RCP<const xmatrix_t> cnewM(mMigrate);
291 RCP<Zoltan2::XpetraCrsMatrixAdapter<xmatrix_t> > newInput;
292 try{
293 newInput =
295 }
296 catch (std::exception &e){
297 aok = false;
298 std::cout << e.what() << std::endl;
299 }
300 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 4 ", 1);
301
302 if (rank==0){
303 std::cout <<
304 "Input adapter for Xpetra::CrsMatrix migrated to proc 0" <<
305 std::endl;
306 }
307 fail = verifyInputAdapter<xmatrix_t>(*newInput, *newM);
308 if (fail) fail += 100;
309 gfail = globalFail(*comm, fail);
310 }
311 }
312 if (gfail){
313 printFailureCode(*comm, fail);
314 }
315 }
316
317#ifdef HAVE_EPETRA_DATA_TYPES
319 // User object is Epetra_CrsMatrix
320 typedef Epetra_CrsMatrix ematrix_t;
321 if (!gfail){
322 if (rank==0)
323 std::cout << "Input adapter for Epetra_CrsMatrix" << std::endl;
324
325 RCP<ematrix_t> eM = uinput->getUIEpetraCrsMatrix();
326 RCP<const ematrix_t> ceM = rcp_const_cast<const ematrix_t>(eM);
327 RCP<Zoltan2::XpetraCrsMatrixAdapter<ematrix_t> > eMInput;
328
329 bool goodAdapter = true;
330 try {
331 eMInput =
333 }
334 catch (std::exception &e){
335 if (std::is_same<znode_t, Xpetra::EpetraNode>::value) {
336 aok = false;
337 goodAdapter = false;
338 std::cout << e.what() << std::endl;
339 }
340 else {
341 // We expect an error from Xpetra when znode_t != Xpetra::EpetraNode
342 // Ignore it, but skip tests using this matrix adapter.
343 std::cout << "Node type is not supported by Xpetra's Epetra interface;"
344 << " Skipping this test." << std::endl;
345 std::cout << "FYI: Here's the exception message: " << std::endl
346 << e.what() << std::endl;
347 goodAdapter = false;
348 }
349 }
350 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 5 ", 1);
351
352 if (goodAdapter) {
353 fail = verifyInputAdapter<ematrix_t>(*eMInput, *tM);
354
355 gfail = globalFail(*comm, fail);
356
357 if (!gfail){
358 ematrix_t *mMigrate =NULL;
359 try{
360 eMInput->applyPartitioningSolution(*eM, mMigrate, solution);
361 }
362 catch (std::exception &e){
363 std::cout << "Error caught: " << e.what() << std::endl;
364 fail = 11;
365 }
366
367 gfail = globalFail(*comm, fail);
368
369 if (!gfail){
370 RCP<const ematrix_t> cnewM(mMigrate, true);
371 RCP<Zoltan2::XpetraCrsMatrixAdapter<ematrix_t> > newInput;
372 try{
373 newInput =
375 }
376 catch (std::exception &e){
377 aok = false;
378 std::cout << e.what() << std::endl;
379 }
380 TEST_FAIL_AND_EXIT(*comm, aok, "XpetraCrsMatrixAdapter 6 ", 1);
381
382 if (rank==0){
383 std::cout <<
384 "Input adapter for Epetra_CrsMatrix migrated to proc 0" <<
385 std::endl;
386 }
387 fail = verifyInputAdapter<ematrix_t>(*newInput, *newM);
388 if (fail) fail += 100;
389 gfail = globalFail(*comm, fail);
390 }
391 }
392 if (gfail){
393 printFailureCode(*comm, fail);
394 }
395 }
396 }
397#endif
398
400 // DONE
401
402 if (rank==0)
403 std::cout << "PASS" << std::endl;
404}
int globalFail(const Comm< int > &comm, int fail)
void printFailureCode(const Comm< int > &comm, int fail)
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
Xpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > xmatrix_t
int verifyInputAdapter(Zoltan2::XpetraCrsMatrixAdapter< User > &ia, tmatrix_t &M)
void printMatrix(RCP< const Comm< int > > &comm, zlno_t nrows, const zgno_t *rowIds, const offset_t *offsets, const zgno_t *colIds)
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tmatrix_t
Traits for application input objects.
common code used by tests
Tpetra::Map ::local_ordinal_type zlno_t
Tpetra::Map ::global_ordinal_type zgno_t
Defines the XpetraCrsMatrixAdapter class.
int main()
The user parameters, debug, timing and memory profiling output objects, and error checking methods.
A PartitioningSolution is a solution to a partitioning problem.
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
static const std::string fail
SparseMatrixAdapter_t::part_t part_t
default_offset_t offset_t
The data type to represent offsets.