Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
MatrixAdapter.cpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Zoltan2: A package of combinatorial algorithms for scientific computing
4 //
5 // Copyright 2012 NTESS and the Zoltan2 contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 //
11 // Testing of GraphAdapter built from Xpetra matrix input adapters.
12 //
13 
18 #include "Kokkos_UnorderedMap.hpp"
20 #include <Zoltan2_InputTraits.hpp>
21 #include <Zoltan2_TestHelpers.hpp>
24 
25 #include <bitset>
26 #include <iostream>
27 #include <string>
28 
29 #include <Teuchos_ArrayView.hpp>
30 #include <Teuchos_Comm.hpp>
31 #include <Teuchos_DefaultComm.hpp>
32 #include <Teuchos_TestingHelpers.hpp>
33 
34 using Teuchos::ArrayView;
35 using Teuchos::Comm;
36 using Teuchos::RCP;
37 using Teuchos::rcp;
38 using Teuchos::rcp_const_cast;
39 using Teuchos::rcp_dynamic_cast;
40 
42 
43 using tcrsMatrix_t = Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t>;
44 using trowMatrix_t = Tpetra::RowMatrix<zscalar_t, zlno_t, zgno_t, znode_t>;
45 using tmap_t = Tpetra::Map<zlno_t, zgno_t, znode_t>;
46 
51 
53 int main(int narg, char *arg[]) {
54  Tpetra::ScopeGuard tscope(&narg, &arg);
55  Teuchos::RCP<const Teuchos::Comm<int>> comm = Tpetra::getDefaultComm();
56 
57  int rank = comm->getRank();
58 
59  int nVtxWeights = 2;
60  int nnzWgtIdx = -1;
61  std::string fname("simple");
62 
63  if (rank == 0)
64  std::cout << "TESTING base case (global)" << std::endl;
65 
66  // Input generator
67  UserInputForTests *uinput;
68 
69  uinput = new UserInputForTests(testDataFilePath, fname, comm, true);
70 
71  RCP<tcrsMatrix_t> M = uinput->getUITpetraCrsMatrix();
72  RCP<trowMatrix_t> trM = rcp_dynamic_cast<trowMatrix_t>(M);
73  RCP<const trowMatrix_t> ctrM = rcp_const_cast<const trowMatrix_t>(trM);
74  zlno_t nLocalRows = M->getLocalNumRows();
75  std::cout << "nLocalRows: " << nLocalRows << std::endl;
76 
77  // Weights:
78  zscalar_t **rowWeights = nullptr;
79  // create as many 1-D weights views as nVtxWeights
81  nLocalRows);
83  nLocalRows);
84 
85  auto wgts0Host = Kokkos::create_mirror_view(wgts0);
86  auto wgts1Host = Kokkos::create_mirror_view(wgts1);
87 
88  for (zlno_t i = 0; i < nLocalRows; i++) {
89  wgts0Host(i) = i;
90  wgts1Host(i) = 200000 + i;
91  }
92 
93  if (nVtxWeights > 0) {
94  rowWeights = new zscalar_t *[nVtxWeights];
95  for (int i = 0; i < nVtxWeights; i++) {
96  if (nnzWgtIdx == i) {
97  rowWeights[i] = nullptr;
98  } else {
99  rowWeights[i] = new zscalar_t[nLocalRows];
100  for (zlno_t j = 0; j < nLocalRows; j++) {
101  rowWeights[i][j] = 200000 * i + j;
102  }
103  }
104  }
105  }
106 
107  Kokkos::deep_copy(wgts0, wgts0Host);
108  Kokkos::deep_copy(wgts1, wgts1Host);
109 
110  tRowMAdapter_t tmi(ctrM, nVtxWeights);
111  for (int i = 0; i < nVtxWeights; i++) {
112  tmi.setWeights(rowWeights[i], 1, i);
113  }
114  tmi.setRowWeightsDevice(wgts0, 0);
115  tmi.setRowWeightsDevice(wgts1, 1);
116 
117  simpleVAdapter_t *via = nullptr;
118 
119  // Set up some fake input
120  zscalar_t **coords = nullptr;
121  int coordDim = 3;
122 
123  if (coordDim > 0) {
124  coords = new zscalar_t *[coordDim];
125  for (int i = 0; i < coordDim; i++) {
126  coords[i] = new zscalar_t[nLocalRows];
127  for (zlno_t j = 0; j < nLocalRows; j++) {
128  coords[i][j] = 100000 * i + j;
129  }
130  }
131  }
132 
133  zgno_t *gids = nullptr;
134  if (coordDim > 0) {
135  gids = new zgno_t[nLocalRows];
136  for (zlno_t i = 0; i < nLocalRows; i++)
137  gids[i] = M->getRowMap()->getGlobalElement(i);
138  via = new simpleVAdapter_t(nLocalRows, gids, coords[0],
139  (coordDim > 1 ? coords[1] : nullptr),
140  (coordDim > 2 ? coords[2] : nullptr), 1, 1, 1);
141  tmi.setCoordinateInput(via);
142  }
143 
144  // TEST of getIDsView, getIDsKokkosView, getRowIDsView, getRowIDsHostView and
145  // getRowIDsDeviceView
146 
147  const zgno_t *ids;
148  const zgno_t *rowIds;
152 
153  tmi.getIDsView(ids);
154  tmi.getRowIDsView(rowIds);
155  tmi.getIDsKokkosView(kIds);
156  auto kIdsHost = Kokkos::create_mirror_view(kIds);
157  Kokkos::deep_copy(kIdsHost, kIds);
158 
159  tmi.getRowIDsHostView(kHostIds);
160  tmi.getRowIDsDeviceView(kDeviceIds);
161 
162  auto kDeviceIdsHost = Kokkos::create_mirror_view(kDeviceIds);
163  Kokkos::deep_copy(kDeviceIdsHost, kDeviceIds);
164 
165  bool success = true;
166  for (size_t i = 0; i < tmi.getLocalNumIDs(); ++i) {
167  TEUCHOS_TEST_EQUALITY(ids[i], kIdsHost(i), std::cout, success);
168  TEUCHOS_TEST_EQUALITY(ids[i], rowIds[i], std::cout, success);
169  TEUCHOS_TEST_EQUALITY(kIdsHost(i), kHostIds(i), std::cout, success);
170  TEUCHOS_TEST_EQUALITY(kIdsHost(i), kDeviceIdsHost(i), std::cout, success);
171  }
172  TEST_FAIL_AND_EXIT(*comm, success, "ids != vertexIds != kIds", 1)
173 
174  // TEST of getCRSView, getCRSHostView and getCRSDeviceView
175  // const zoffset_t *offsets;
176  ArrayRCP<const zoffset_t> offsets;
177  ArrayRCP<const zgno_t> colIds;
182  tmi.getCRSView(offsets, colIds);
183  tmi.getCRSHostView(kHostOffsets, kHostColIds);
184  tmi.getCRSDeviceView(kDeviceOffsets, kDeviceColIds);
185 
186  auto kDeviceColIdsHost = Kokkos::create_mirror_view(kDeviceColIds);
187  Kokkos::deep_copy(kDeviceColIdsHost, kDeviceColIds);
188 
189  auto kDeviceOffsetsHost = Kokkos::create_mirror_view(kDeviceOffsets);
190  Kokkos::deep_copy(kDeviceOffsetsHost, kDeviceOffsets);
191 
192  for (int i = 0; success && i < colIds.size(); i++) {
193  TEUCHOS_TEST_EQUALITY(colIds[i], kHostColIds(i), std::cout, success);
194  TEUCHOS_TEST_EQUALITY(colIds[i], kDeviceColIdsHost(i), std::cout, success);
195  }
196  TEST_FAIL_AND_EXIT(*comm, success, "colIds != kHostColIds != kDeviceColIds",
197  1)
198 
199  for (int i = 0; success && i < offsets.size(); i++) {
200  TEUCHOS_TEST_EQUALITY(offsets[i], kHostOffsets(i), std::cout, success);
201  TEUCHOS_TEST_EQUALITY(offsets[i], kDeviceOffsetsHost(i), std::cout,
202  success);
203  }
204  TEST_FAIL_AND_EXIT(*comm, success,
205  "offsets != kHostOffsets != kDeviceOffsets", 1)
206 
207  ArrayRCP<const zscalar_t> values;
210 
211  tmi.getCRSView(offsets, colIds, values);
212  tmi.getCRSHostView(kHostOffsets, kHostColIds, kHostValues);
213  tmi.getCRSDeviceView(kDeviceOffsets, kDeviceColIds, kDeviceValues);
214  auto kDeviceValuesHost = Kokkos::create_mirror_view(kDeviceValues);
215  Kokkos::deep_copy(kDeviceValuesHost, kDeviceValues);
216 
217  for (int i = 0; success && i < colIds.size(); i++) {
218  TEUCHOS_TEST_EQUALITY(colIds[i], kHostColIds(i), std::cout, success);
219  TEUCHOS_TEST_EQUALITY(colIds[i], kDeviceColIdsHost(i), std::cout, success);
220  }
221  TEST_FAIL_AND_EXIT(*comm, success, "colIds != kHostColIds != kDeviceColIds",
222  1)
223 
224  for (int i = 0; success && i < offsets.size(); i++) {
225  TEUCHOS_TEST_EQUALITY(offsets[i], kHostOffsets(i), std::cout, success);
226  TEUCHOS_TEST_EQUALITY(offsets[i], kDeviceOffsetsHost(i), std::cout,
227  success);
228  }
229  TEST_FAIL_AND_EXIT(*comm, success,
230  "offsets != kHostOffsets != kDeviceOffsets", 1)
231 
232  for (int i = 0; success && i < values.size(); i++) {
233  TEUCHOS_TEST_EQUALITY(values[i], kHostValues(i), std::cout, success);
234  TEUCHOS_TEST_EQUALITY(values[i], kDeviceValuesHost(i), std::cout, success);
235  }
236  TEST_FAIL_AND_EXIT(*comm, success, "values != kHostValues != kDeviceValues",
237  1)
238 
239  // TEST of getRowWeightsView, getRowWeightsHost0View and
240  // getRowWeightsDeviceView
241 
246 
247  tmi.getRowWeightsHostView(weightsHost0, 0);
248  tmi.getRowWeightsDeviceView(weightsDevice0, 0);
249 
250  tmi.getRowWeightsHostView(weightsHost1, 1);
251  tmi.getRowWeightsDeviceView(weightsDevice1, 1);
252 
253  auto hostWeightsDevice0 = Kokkos::create_mirror_view(weightsDevice0);
254  Kokkos::deep_copy(hostWeightsDevice0, weightsDevice0);
255 
256  auto hostWeightsDevice1 = Kokkos::create_mirror_view(weightsDevice1);
257  Kokkos::deep_copy(hostWeightsDevice1, weightsDevice1);
258 
259  for (int w = 0; success && w < nVtxWeights; w++) {
260  const zscalar_t *wgts;
261 
262  Kokkos::View<zscalar_t *, typename znode_t::device_type> wkgts;
263  int stride;
264  tmi.getRowWeightsView(wgts, stride, w);
265 
266  if (w == 0) {
267 
268  for (zlno_t i = 0; success && i < nLocalRows; ++i) {
269  TEUCHOS_TEST_EQUALITY(wgts[stride * i], weightsHost0(i), std::cout,
270  success);
271  TEUCHOS_TEST_EQUALITY(wgts[stride * i], hostWeightsDevice0(i), std::cout,
272  success);
273  }
274  } else {
275  for (zlno_t i = 0; success && i < nLocalRows; ++i) {
276  TEUCHOS_TEST_EQUALITY(wgts[stride * i], weightsHost1(i), std::cout,
277  success);
278  TEUCHOS_TEST_EQUALITY(wgts[stride * i], hostWeightsDevice1(i), std::cout,
279  success);
280  }
281  }
282  }
283  TEST_FAIL_AND_EXIT(*comm, success, "wgts != vwgts != wkgts", 1)
284 
285  // TEST of useNumNonzerosAsRowWeight, useNumNonzerosAsRowWeightHost and
286  // useNumNonzerosAsRowWeightDevice
287  // for (int w = 0; success && w < nVtxWeights; w++) {
288  // TEUCHOS_TEST_EQUALITY(tmi.useNumNonzerosAsRowWeight(w),
289  // tmi.useNumNonzerosAsRowWeightHost(w), std::cout,
290  // success);
291  // TEUCHOS_TEST_EQUALITY(tmi.useNumNonzerosAsRowWeight(w),
292  // tmi.useNumNonzerosAsRowWeightDevice(w), std::cout,
293  // success);
294  // }
295  // TEST_FAIL_AND_EXIT(
296  // *comm, success,
297  // "useNumNonzerosAsRowWeight != useNumNonzerosAsRowWeightHost != "
298  // "useNumNonzerosAsRowWeightDevice",
299  // 1)
300 
301  // clean
302  if (nVtxWeights > 0) {
303  for (int i = 0; i < nVtxWeights; i++) {
304  if (rowWeights[i])
305  delete[] rowWeights[i];
306  }
307  delete[] rowWeights;
308  }
309 
310  if (coordDim > 0) {
311  delete via;
312  delete[] gids;
313  for (int i = 0; i < coordDim; i++) {
314  if (coords[i])
315  delete[] coords[i];
316  }
317  delete[] coords;
318  }
319 
320  delete uinput;
321 
322  if (rank == 0)
323  std::cout << "PASS" << std::endl;
324 
325  return 0;
326 }
typename baseMAdapter_t::offset_t zoffset_t
typename Zoltan2::TpetraRowMatrixAdapter< trowMatrix_t, simpleUser_t > tRowMAdapter_t
MatrixAdapter defines the adapter interface for matrices.
RCP< tcrsMatrix_t > getUITpetraCrsMatrix()
Kokkos::View< scalar_t *, device_t > WeightsDeviceView1D
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
int main(int narg, char **arg)
Definition: coloring1.cpp:164
common code used by tests
Kokkos::View< const gno_t *, device_t > ConstIdsDeviceView
typename ConstScalarsDeviceView::HostMirror ConstScalarsHostView
Provides access for Zoltan2 to Tpetra::RowMatrix data.
list fname
Begin.
Definition: validXML.py:19
Defines the XpetraCrsMatrixAdapter class.
Kokkos::View< const scalar_t *, device_t > ConstScalarsDeviceView
Zoltan2::BasicVectorAdapter< simpleUser_t > simpleVAdapter_t
typename ConstOffsetsDeviceView::HostMirror ConstOffsetsHostView
typename ConstIdsDeviceView::HostMirror ConstIdsHostView
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
Traits for application input objects.
Defines the TpetraRowMatrixAdapter class.
Tpetra::Map::local_ordinal_type zlno_t
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tcrsMatrix_t
Tpetra::Map< zlno_t, zgno_t, znode_t > tmap_t
typename InputTraits< User >::offset_t offset_t
typename WeightsDeviceView1D::HostMirror WeightsHostView1D
float zscalar_t
Defines the BasicVectorAdapter class.
Kokkos::View< const offset_t *, device_t > ConstOffsetsDeviceView
Tpetra::Map::global_ordinal_type zgno_t
Tpetra::RowMatrix< zscalar_t, zlno_t, zgno_t, znode_t > trowMatrix_t
std::string testDataFilePath(".")