8 #include <Teuchos_RCP.hpp> 9 #include <Teuchos_Array.hpp> 10 #include <Teuchos_GlobalMPISession.hpp> 11 #include <Teuchos_ParameterList.hpp> 12 #include "Teuchos_XMLParameterListHelpers.hpp" 14 #include "Tpetra_MultiVector.hpp" 15 #include <Tpetra_CrsGraph.hpp> 16 #include <Tpetra_Map.hpp> 21 typedef Tpetra::CrsGraph<zlno_t, zgno_t, znode_t>
tcrsGraph_t;
22 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>
tMVector_t;
26 int main(
int argc,
char *argv[]){
30 Teuchos::GlobalMPISession session(&argc, &argv);
31 Teuchos::RCP<const Teuchos::Comm<int> > tcomm = Teuchos::DefaultComm<int>::getComm();
32 int nx = 2, ny = 2, nz = 2;
33 for (
int i = 1 ; i < argc ; ++i ) {
34 if ( 0 == strcasecmp( argv[i] ,
"NX" ) ) {
35 nx = atoi( argv[++i] );
37 else if ( 0 == strcasecmp( argv[i] ,
"NY" ) ) {
38 ny = atoi( argv[++i] );
40 else if ( 0 == strcasecmp( argv[i] ,
"NZ" ) ) {
41 nz = atoi( argv[++i] );
44 std::cerr <<
"Unrecognized command line argument #" << i <<
": " << argv[i] << std::endl ;
53 int rank = tcomm->getRank();
54 part_t numProcs = tcomm->getSize();
57 zgno_t numGlobalTasks = nx*ny*nz;
59 zgno_t myTasks = numGlobalTasks / numProcs;
60 zgno_t taskLeftOver = numGlobalTasks % numProcs;
61 if (rank < taskLeftOver ) ++myTasks;
63 zgno_t myTaskBegin = (numGlobalTasks / numProcs) * rank;
64 myTaskBegin += (taskLeftOver < rank ? taskLeftOver : rank);
65 zgno_t myTaskEnd = myTaskBegin + myTasks;
69 for(
int i = 0; i < coordDim; ++i){
75 zlno_t *task_communication_xadj_ =
new zlno_t [myTasks+1];
76 zgno_t *task_communication_adj_ =
new zgno_t [myTasks * 6];
79 task_communication_xadj_[0] = 0;
80 for (
zlno_t i = myTaskBegin; i < myTaskEnd; ++i) {
81 task_gnos[i - myTaskBegin] = i;
84 zlno_t y = (i / (nx)) % ny;
85 zlno_t z = (i / (nx)) / ny;
86 partCenters[0][i - myTaskBegin] = x;
87 partCenters[1][i - myTaskBegin] = y;
88 partCenters[2][i - myTaskBegin] = z;
93 task_communication_adj_[prevNCount++] = i - 1;
96 task_communication_adj_[prevNCount++] = i + 1;
99 task_communication_adj_[prevNCount++] = i - nx;
102 task_communication_adj_[prevNCount++] = i + nx;
105 task_communication_adj_[prevNCount++] = i - nx * ny;
108 task_communication_adj_[prevNCount++] = i + nx * ny;
110 task_communication_xadj_[i+1 - myTaskBegin] = prevNCount;
113 RCP<my_adapter_t> ia;
114 typedef Tpetra::Map<>::node_type mytest_znode_t;
115 typedef Tpetra::Map<zlno_t, zgno_t, mytest_znode_t> map_t;
116 RCP<const map_t> map = rcp (
new map_t (numGlobalTasks, myTasks, 0, tcomm));
118 RCP<tcrsGraph_t> TpetraCrsGraph(
new tcrsGraph_t (map, 0));
121 for (
zlno_t lclRow = 0; lclRow < myTasks; ++lclRow) {
122 const zgno_t gblRow = map->getGlobalElement (lclRow);
123 zgno_t begin = task_communication_xadj_[lclRow];
124 zgno_t end = task_communication_xadj_[lclRow + 1];
125 const ArrayView< const zgno_t > indices(task_communication_adj_+begin, end-begin);
126 TpetraCrsGraph->insertGlobalIndices(gblRow, indices);
128 TpetraCrsGraph->fillComplete ();
129 RCP<const tcrsGraph_t> const_data = rcp_const_cast<
const tcrsGraph_t>(TpetraCrsGraph);
133 const int coord_dim = 3;
134 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(coord_dim);
137 Teuchos::ArrayView<const zscalar_t> a(partCenters[0], myTasks);
139 Teuchos::ArrayView<const zscalar_t> b(partCenters[1], myTasks);
141 Teuchos::ArrayView<const zscalar_t> c(partCenters[2], myTasks);
145 Teuchos::ArrayView<const zscalar_t> a;
150 RCP<tMVector_t> coords(
new tMVector_t(map, coordView.view(0, coord_dim), coord_dim));
151 RCP<const tMVector_t> const_coords = rcp_const_cast<
const tMVector_t>(coords);
153 ia->setCoordinateInput(adapter.getRawPtr());
168 ParameterList zoltan2_parameters;
169 zoltan2_parameters.set(
"compute_metrics",
true);
170 zoltan2_parameters.set(
"imbalance_tolerance", 1.0);
171 zoltan2_parameters.set(
"num_global_parts", tcomm->getSize());
172 zoltan2_parameters.set(
"algorithm",
"multijagged");
173 zoltan2_parameters.set(
"mj_keep_part_boxes",
false);
174 zoltan2_parameters.set(
"mj_recursion_depth", 3);
175 RCP<xcrsGraph_problem_t> partition_problem;
176 partition_problem = RCP<xcrsGraph_problem_t> (
new xcrsGraph_problem_t(ia.getRawPtr(),&zoltan2_parameters,tcomm));
179 partition_problem->solve();
181 RCP<const Zoltan2::Environment> env = partition_problem->getEnvironment();
183 RCP<quality_t>metricObject =
184 rcp(
new quality_t(ia.getRawPtr(), &zoltan2_parameters, tcomm,
185 &partition_problem->getSolution()));
187 if (tcomm->getRank() == 0){
188 metricObject->printMetrics(std::cout);
190 partition_problem->printTimers();
192 part_t *proc_to_task_xadj_ = NULL;
193 part_t *proc_to_task_adj_ = NULL;
201 Teuchos::rcpFromRef(mach),
203 rcpFromRef(partition_problem->getSolution()),
207 ctm.
getProcTask(proc_to_task_xadj_, proc_to_task_adj_);
209 if (tcomm->getRank() == 0){
210 for (part_t i = 0; i < numProcs; ++i){
211 std::cout <<
"\nProc i:" << i <<
" ";
212 for (part_t j = proc_to_task_xadj_[i]; j < proc_to_task_xadj_[i+1]; ++j){
213 std::cout <<
" " << proc_to_task_adj_[j];
216 std::cout << std::endl;
223 zlno_t prevNCount_tmp = 0;
224 zgno_t *task_communication_adj_tmp =
new zgno_t [numGlobalTasks * 6];
225 zlno_t *task_communication_xadj_tmp =
new zlno_t [numGlobalTasks+1];
226 task_communication_xadj_tmp[0] = 0;
227 for (
zlno_t i = 0; i < numGlobalTasks; ++i) {
229 zlno_t y = (i / (nx)) % ny;
230 zlno_t z = (i / (nx)) / ny;
233 task_communication_adj_tmp[prevNCount_tmp++] = i - 1;
236 task_communication_adj_tmp[prevNCount_tmp++] = i + 1;
239 task_communication_adj_tmp[prevNCount_tmp++] = i - nx;
242 task_communication_adj_tmp[prevNCount_tmp++] = i + nx;
245 task_communication_adj_tmp[prevNCount_tmp++] = i - nx * ny;
248 task_communication_adj_tmp[prevNCount_tmp++] = i + nx * ny;
250 task_communication_xadj_tmp[i+1] = prevNCount_tmp;
254 std::vector <int> mach_extent(mach_coord_dim);
257 std::vector <part_t> all_parts(numGlobalTasks), copy(numGlobalTasks, 0);
259 const part_t *parts = partition_problem->getSolution().getPartListView();
264 for (part_t i = 0; i < myTasks; ++i){
265 zgno_t g = map->getGlobalElement(i);
269 reduceAll<int, part_t>(
281 int *machine_extent =
new int [mach_coord_dim];
282 bool *machine_extent_wrap_around =
new bool[mach_coord_dim];
287 for (
zlno_t i = 0; i < numGlobalTasks; ++i){
288 zlno_t b = task_communication_xadj_tmp[i];
289 zlno_t e = task_communication_xadj_tmp[i+1];
293 for (
zlno_t j = b; j < e; ++j){
294 zgno_t n = task_communication_adj_tmp[j];
300 for (
int k = 0 ; k < mach_coord_dim ; ++k){
301 part_t distance =
ZOLTAN2_ABS(proc_coords[k][procId1] - proc_coords[k][procId2]);
302 if (machine_extent_wrap_around[k]){
303 if (machine_extent[k] - distance < distance){
304 distance = machine_extent[k] - distance;
311 delete [] machine_extent_wrap_around;
312 delete [] machine_extent;
314 if (tcomm->getRank() == 0)
315 std::cout <<
"HOPS:" << hops <<
" HOPS2:" << hops2 << std::endl;
317 delete [] task_communication_xadj_tmp;
318 delete [] task_communication_adj_tmp;
346 if (tcomm->getRank() == 0){
347 cout <<
"PASS" << endl;
351 for (
int i = 0; i < coordDim; i++)
delete [] partCenters[i];
352 delete [] partCenters;
355 catch(std::string &s){
362 catch(
char const * s){
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
Zoltan2::XpetraCrsGraphAdapter< tcrsGraph_t, tMVector_t > my_adapter_t
bool getHopCount(int rank1, int rank2, pcoord_t &hops) const
Provides access for Zoltan2 to Xpetra::CrsGraph data.
common code used by tests
part_t getAssignedProcForTask(part_t taskId)
getAssignedProcForTask function, returns the assigned processor id for the given task ...
int main(int argc, char *argv[])
Defines the XpetraMultiVectorAdapter.
Defines XpetraCrsGraphAdapter class.
void getProcTask(part_t *&proc_to_task_xadj_, part_t *&proc_to_task_adj_)
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t
bool getMachineExtentWrapArounds(bool *wrap_around) const
if the machine has a wrap-around tourus link in each dimension. return true if the information is ava...
int getMachineDim() const
returns the dimension (number of coords per node) in the machine
InputTraits< User >::part_t part_t
bool getAllMachineCoordinatesView(pcoord_t **&allCoords) const
getProcDim function set the coordinates of all ranks allCoords[i][j], i=0,...,getMachineDim(), j=0,...,getNumRanks(), is the i-th dimensional coordinate for rank j. return true if coordinates are available for all ranks
An adapter for Xpetra::MultiVector.
MachineRepresentation Class Base class for representing machine coordinates, networks, etc.
PartitioningProblem sets up partitioning problems for the user.
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
bool getMachineExtent(int *nxyz) const
sets the number of unique coordinates in each machine dimension return true if coordinates are availa...
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t