54 #include <Tpetra_Map.hpp>
60 template <
typename User>
75 : nids(nids_), gids(gids_), dim(dim_), coords(coords_), weights(weights_)
85 { wgt = weights; stride = 1; }
91 coo = &(coords[idx*nids]);
105 template <
typename User>
118 : nids(nids_), gids(gids_), dim(dim_), coords(coords_), weights(weights_)
128 { wgt = weights; stride = 1; }
134 coo = &(coords[idx]);
150 template <
typename User>
165 : nids(nids_), dim(dim_)
169 typedef Kokkos::DualView<gno_t *, device_t> view_ids_t;
170 kokkos_gids = view_ids_t(
171 Kokkos::ViewAllocateWithoutInitializing(
"gids"), nids);
173 auto host_gids = kokkos_gids.h_view;
174 for(
size_t n = 0; n < nids; ++n) {
175 host_gids(n) = gids_[n];
178 kokkos_gids.template modify<typename view_ids_t::host_mirror_space>();
179 kokkos_gids.sync_host();
180 kokkos_gids.template sync<device_t>();
186 typedef Kokkos::DualView<scalar_t **, device_t> view_weights_t;
187 kokkos_weights = view_weights_t(
188 Kokkos::ViewAllocateWithoutInitializing(
"weights"), nids, 0);
189 auto host_kokkos_weights = kokkos_weights.h_view;
190 for(
size_t n = 0; n < nids; ++n) {
191 host_kokkos_weights(n,0) = weights_[n];
194 kokkos_weights.template modify<typename view_weights_t::host_mirror_space>();
195 kokkos_weights.sync_host();
196 kokkos_weights.template sync<device_t>();
201 typedef Kokkos::DualView<scalar_t **, Kokkos::LayoutLeft, device_t> kokkos_coords_t;
202 kokkos_coords = kokkos_coords_t(
203 Kokkos::ViewAllocateWithoutInitializing(
"coords"), nids, dim);
204 auto host_kokkos_coords = kokkos_coords.h_view;
206 for(
size_t n = 0; n < nids; ++n) {
207 for(
int idx = 0; idx < dim; ++idx) {
208 host_kokkos_coords(n,idx) = coords_[n+idx*nids];
212 kokkos_coords.template modify<typename kokkos_coords_t::host_mirror_space>();
213 kokkos_coords.sync_host();
214 kokkos_coords.template sync<device_t>();
221 auto kokkosIds = kokkos_gids.view_host();
222 ids = kokkosIds.data();
226 ids = kokkos_gids.template view<device_t>();
232 int idx = 0)
const override
234 auto h_wgts_2d = kokkos_weights.view_host();
236 wgt = Kokkos::subview(h_wgts_2d, Kokkos::ALL, idx).data();
241 wgt = kokkos_weights.template view<device_t>();
247 int &stride,
int idx = 0)
const override {
248 elements = kokkos_coords.view_host().data();
253 Kokkos::LayoutLeft,
device_t> & coo)
const {
254 coo = kokkos_coords.template view<device_t>();
259 Kokkos::DualView<gno_t *, device_t> kokkos_gids;
261 Kokkos::DualView<scalar_t **, Kokkos::LayoutLeft, device_t> kokkos_coords;
262 Kokkos::DualView<scalar_t **, device_t> kokkos_weights;
270 const Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
272 int rank = comm->getRank();
273 int nprocs = comm->getSize();
275 typedef Tpetra::Map<> Map_t;
276 typedef Map_t::local_ordinal_type localId_t;
277 typedef Map_t::global_ordinal_type globalId_t;
278 typedef double scalar_t;
289 size_t localCount = 40;
293 scalar_t *cStrided =
new scalar_t [dim * localCount];
295 for (
size_t i = 0; i < localCount; i++)
296 for (
int d = 0; d < dim; d++)
297 cStrided[cnt++] = d*1000 + rank*100 + i;
300 scalar_t *cContig =
new scalar_t [dim * localCount];
302 for (
int d = 0; d < dim; d++)
303 for (
size_t i = 0; i < localCount; i++)
304 cContig[cnt++] = d*1000 + rank*100 + i;
307 globalId_t *globalIds =
new globalId_t [localCount];
308 globalId_t offset = rank * localCount;
309 for (
size_t i=0; i < localCount; i++) globalIds[i] = offset++;
314 Teuchos::ParameterList params(
"test params");
315 params.set(
"debug_level",
"basic_status");
316 params.set(
"error_check_level",
"debug_mode_assertions");
318 params.set(
"algorithm", algorithm);
319 params.set(
"num_global_parts", nprocs+1);
325 stridedAdapter_t *ia1 =
326 new stridedAdapter_t(localCount,globalIds,dim,cStrided);
333 quality_t *metricObject1 =
new quality_t(ia1, ¶ms, comm,
341 std::cout <<
"no weights -- balance satisfied: " << imb << std::endl;
343 std::cout <<
"no weights -- balance failure: " << imb << std::endl;
346 std::cout << std::endl;
350 contigAdapter_t *ia2 =
new contigAdapter_t(localCount,globalIds,dim,cContig);
358 kokkosAdapter_t *ia3 =
new kokkosAdapter_t(localCount,globalIds,dim,cContig);
367 for (
size_t i = 0; i < localCount; i++) {
368 if((problem1->
getSolution().getPartListView()[i] !=
373 std::cout << rank <<
" Error: differing parts for index " << i
374 << problem1->
getSolution().getPartListView()[i] <<
" "
375 << problem2->
getSolution().getPartListView()[i] <<
" "
376 << problem3->
getSolution().getPartListView()[i] << std::endl;
381 if (ndiff > 0) nFail++;
382 else if (rank == 0) std::cout <<
"no weights -- comparisons OK " << std::endl;
384 delete metricObject1;
396 scalar_t *
weights =
new scalar_t [localCount];
397 for (
size_t i=0; i < localCount; i++) weights[i] = 1 + scalar_t(rank);
400 ia1 =
new stridedAdapter_t(localCount, globalIds, dim, cStrided, weights);
414 std::cout <<
"weighted -- balance satisfied " << imb << std::endl;
416 std::cout <<
"weighted -- balance failed " << imb << std::endl;
419 std::cout << std::endl;
423 ia2 =
new contigAdapter_t(localCount, globalIds, dim, cContig, weights);
431 for (
size_t i = 0; i < localCount; i++) {
432 if (problem1->
getSolution().getPartListView()[i] !=
434 std::cout << rank <<
" Error: differing parts for index " << i
435 << problem1->
getSolution().getPartListView()[i] <<
" "
436 << problem2->
getSolution().getPartListView()[i] << std::endl;
441 if (ndiff > 0) nFail++;
442 else if (rank == 0) std::cout <<
"weighted -- comparisons OK " << std::endl;
444 delete metricObject1;
451 if (weights)
delete []
weights;
452 if (cStrided)
delete [] cStrided;
453 if (cContig)
delete [] cContig;
454 if (globalIds)
delete [] globalIds;
459 Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM, 1, &nFail, &gnFail);
464 int main(
int narg,
char *arg[])
466 Tpetra::ScopeGuard scope(&narg, &arg);
467 const Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
474 if (comm->getRank() == 0) {
475 if (err == 0) std::cout <<
"PASS" << std::endl;
476 else std::cout <<
"FAIL: " << err <<
" tests failed" << std::endl;
OldSchoolVectorAdapterStrided(const size_t nids_, const gno_t *gids_, const int dim_, const scalar_t *coords_, const scalar_t *weights_=NULL)
void printMetrics(std::ostream &os) const
Print all metrics.
size_t getLocalNumIDs() const
Returns the number of objects on this process.
Zoltan2::InputTraits< User >::gno_t gno_t
void getEntriesView(const scalar_t *&coo, int &stride, int idx=0) const
int getNumEntriesPerID() const
Return the number of vectors.
map_t::global_ordinal_type gno_t
Zoltan2::InputTraits< User >::scalar_t scalar_t
A simple class that can be the User template argument for an InputAdapter.
Defines the VectorAdapter interface.
void getWeightsView(const scalar_t *&wgt, int &stride, int idx=0) const override
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t
int main(int narg, char **arg)
Defines the PartitioningSolution class.
virtual void getWeightsKokkosView(Kokkos::View< scalar_t **, device_t > &wgt) const
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
void getEntriesView(const scalar_t *&elements, int &stride, int idx=0) const override
Zoltan2::InputTraits< User >::node_t node_t
scalar_t getWeightImbalance(int weightIndex) const
Return the imbalance for the requested weight.
void getIDsView(const gno_t *&ids) const
void getWeightsView(const scalar_t *&wgt, int &stride, int idx=0) const
virtual void getEntriesKokkosView(Kokkos::View< scalar_t **, Kokkos::LayoutLeft, device_t > &coo) const
void getEntriesView(const scalar_t *&coo, int &stride, int idx=0) const
void getIDsView(const gno_t *&ids) const
virtual void getIDsKokkosView(Kokkos::View< const gno_t *, device_t > &ids) const
size_t getLocalNumIDs() const
Returns the number of objects on this process.
VectorAdapter defines the interface for vector input.
int run_test_strided_versus_contig(const std::string &algorithm)
int getNumEntriesPerID() const
Return the number of vectors.
Zoltan2::InputTraits< User >::gno_t gno_t
node_t::device_type device_t
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
void getIDsView(const gno_t *&ids) const override
Tpetra::Map::node_type node_t
void getWeightsView(const scalar_t *&wgt, int &stride, int idx=0) const
Defines the PartitioningProblem class.
Zoltan2::InputTraits< User >::gno_t gno_t
A class that computes and returns quality metrics.
size_t getLocalNumIDs() const
Returns the number of objects on this process.
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
void solve(bool updateInputData=true)
Direct the problem to create a solution.
node_t::device_type device_t
int getNumEntriesPerID() const
Return the number of vectors.
OldSchoolVectorAdapterContig(const size_t nids_, const gno_t *gids_, const int dim_, const scalar_t *coords_, const scalar_t *weights_=NULL)
scalar_t getObjectCountImbalance() const
Return the object count imbalance.
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
Zoltan2::InputTraits< User >::scalar_t scalar_t
KokkosVectorAdapter(const size_t nids_, const gno_t *gids_, const int dim_, const scalar_t *coords_, const scalar_t *weights_=NULL)
Zoltan2::InputTraits< User >::scalar_t scalar_t