57 #include <Teuchos_CommandLineProcessor.hpp>
66 using Teuchos::ArrayView;
67 using Teuchos::ArrayRCP;
68 using Teuchos::CommandLineProcessor;
70 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>
tMVector_t;
71 typedef Tpetra::Map<zlno_t, zgno_t, znode_t>
tMap_t;
81 const RCP<
const Teuchos::Comm<int> > & comm,
88 ArrayRCP<zscalar_t>
weights(wgts, 0, len,
true);
92 for (
zlno_t i=0; i < len; i++)
96 for (
int i=0; i < 10; i++){
98 for (
zlno_t j=i; j < len; j += 10)
104 for (
zlno_t i=0; i < len; i++)
116 const RCP<
const Teuchos::Comm<int> > & comm,
119 int rank = comm->getRank();
120 int nprocs = comm->getSize();
122 double k = log(numGlobalCoords) / 3;
123 double xdimf = exp(k) + 0.5;
124 ssize_t xdim =
static_cast<ssize_t
>(floor(xdimf));
126 ssize_t zdim = numGlobalCoords / (xdim*ydim);
127 ssize_t num=xdim*ydim*zdim;
128 ssize_t diff = numGlobalCoords - num;
132 if (zdim > xdim && zdim > ydim){
134 newdiff = diff - (xdim*ydim);
139 else if (ydim > xdim && ydim > zdim){
141 newdiff = diff - (xdim*zdim);
148 newdiff = diff - (ydim*zdim);
158 diff = numGlobalCoords - num;
160 diff /= -numGlobalCoords;
162 diff /= numGlobalCoords;
166 std::cout <<
"Warning: Difference " << diff*100 <<
" percent" << std::endl;
167 std::cout <<
"Mesh size: " << xdim <<
"x" << ydim <<
"x" <<
168 zdim <<
", " << num <<
" vertices." << std::endl;
173 ssize_t numLocalCoords = num / nprocs;
174 ssize_t leftOver = num % nprocs;
177 if (rank <= leftOver)
178 gid0 =
zgno_t(rank) * (numLocalCoords+1);
180 gid0 = (leftOver * (numLocalCoords+1)) +
181 ((
zgno_t(rank) - leftOver) * numLocalCoords);
186 ssize_t gid1 = gid0 + numLocalCoords;
191 ArrayRCP<zgno_t> idArray(ids, 0, numLocalCoords,
true);
193 for (ssize_t i=gid0; i < gid1; i++)
196 RCP<const tMap_t> idMap = rcp(
197 new tMap_t(num, idArray.view(0, numLocalCoords), 0, comm));
204 ArrayRCP<zscalar_t> coordArray(x, 0, numLocalCoords*3,
true);
211 zgno_t xyPlane = xdim*ydim;
212 zgno_t zStart = gid0 / xyPlane;
213 zgno_t rem = gid0 % xyPlane;
220 for (
zscalar_t zval=zStart; next < numLocalCoords && zval < zdim; zval++){
221 for (
zscalar_t yval=yStart; next < numLocalCoords && yval < ydim; yval++){
222 for (
zscalar_t xval=xStart; next < numLocalCoords && xval < xdim; xval++){
233 ArrayView<const zscalar_t> xArray(x, numLocalCoords);
234 ArrayView<const zscalar_t> yArray(y, numLocalCoords);
235 ArrayView<const zscalar_t> zArray(z, numLocalCoords);
236 ArrayRCP<ArrayView<const zscalar_t> >
coordinates =
237 arcp(
new ArrayView<const zscalar_t> [3], 0, 3);
238 coordinates[0] = xArray;
239 coordinates[1] = yArray;
240 coordinates[2] = zArray;
242 ArrayRCP<const ArrayView<const zscalar_t> > constCoords =
243 coordinates.getConst();
245 RCP<tMVector_t> meshCoords = rcp(
new tMVector_t(
246 idMap, constCoords.view(0,3), 3));
252 int main(
int narg,
char *arg[])
256 Tpetra::ScopeGuard tscope(&narg, &arg);
257 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
259 int rank = comm->getRank();
260 int nprocs = comm->getSize();
265 std::cout <<
"Number of processes: " << nprocs << std::endl;
268 double numGlobalCoords = 1000;
271 string timingType(
"no_timers");
272 string debugLevel(
"basic_status");
273 string memoryOn(
"memoryOn");
274 string memoryOff(
"memoryOff");
275 string memoryProcs(
"0");
277 int numGlobalParts = nprocs;
279 CommandLineProcessor commandLine(
false,
true);
280 commandLine.setOption(
"size", &numGlobalCoords,
281 "Approximate number of global coordinates.");
282 commandLine.setOption(
"testCuts", &numTestCuts,
283 "Number of test cuts to make when looking for bisector.");
284 commandLine.setOption(
"numParts", &numGlobalParts,
285 "Number of parts (default is one per proc).");
286 commandLine.setOption(
"nWeights", &nWeights,
287 "Number of weights per coordinate, zero implies uniform weights.");
289 string balanceCount(
"balance_object_count");
290 string balanceWeight(
"balance_object_weight");
291 string mcnorm1(
"multicriteria_minimize_total_weight");
292 string mcnorm2(
"multicriteria_balance_total_maximum");
293 string mcnorm3(
"multicriteria_minimize_maximum_weight");
295 string objective(balanceWeight);
297 string doc(balanceCount); doc.append(
": ignore weights\n");
299 doc.append(balanceWeight); doc.append(
": balance on first weight\n");
302 doc.append(
": given multiple weights, balance their total.\n");
305 doc.append(
": given multiple weights, balance the maximum for each coordinate.\n");
308 doc.append(
": given multiple weights, balance the L2 norm of the weights.\n");
310 commandLine.setOption(
"objective", &objective, doc.c_str());
312 commandLine.setOption(
"timers", &timingType,
313 "no_timers, micro_timers, macro_timers, both_timers, test_timers");
315 commandLine.setOption(
"debug", &debugLevel,
316 "no_status, basic_status, detailed_status, verbose_detailed_status");
318 commandLine.setOption(memoryOn.c_str(), memoryOff.c_str(), &doMemory,
319 "do memory profiling");
321 commandLine.setOption(
"memoryProcs", &memoryProcs,
322 "list of processes that output memory usage");
324 CommandLineProcessor::EParseCommandLineReturn rc =
325 commandLine.parse(narg, arg);
327 if (rc != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL){
328 if (rc == Teuchos::CommandLineProcessor::PARSE_HELP_PRINTED){
330 std::cout <<
"PASS" << std::endl;
335 std::cout <<
"FAIL" << std::endl;
342 zgno_t globalSize =
static_cast<zgno_t>(numGlobalCoords);
345 size_t numLocalCoords = coordinates->getLocalLength();
349 for (
int p=0; p < nprocs; p++){
351 std::cout <<
"Rank " << rank <<
", " << numLocalCoords <<
"coords" << std::endl;
352 const zscalar_t *x = coordinates->getData(0).getRawPtr();
353 const zscalar_t *y = coordinates->getData(1).getRawPtr();
354 const zscalar_t *z = coordinates->getData(2).getRawPtr();
355 for (
zlno_t i=0; i < numLocalCoords; i++)
356 std::cout <<
" " << x[i] <<
" " << y[i] <<
" " << z[i] << std::endl;
363 Array<ArrayRCP<zscalar_t> >
weights(nWeights);
368 for (
int i=0; i < nWeights; i++){
378 MEMORY_CHECK(doMemory && rank==0,
"After creating input");
381 const RCP<const tMap_t> &coordmap = coordinates->getMap();
382 ArrayView<const zgno_t> ids = coordmap->getLocalElementList();
383 const zgno_t *globalIds = ids.getRawPtr();
385 size_t localCount = coordinates->getLocalLength();
387 RCP<inputAdapter_t> ia;
390 ia = rcp(
new inputAdapter_t (localCount, globalIds,
391 coordinates->getData(0).getRawPtr(), coordinates->getData(1).getRawPtr(),
392 coordinates->getData(2).getRawPtr(), 1,1,1));
395 vector<const zscalar_t *> values(3);
396 for (
int i=0; i < 3; i++)
397 values[i] = coordinates->getData(i).getRawPtr();
398 vector<int> valueStrides(0);
399 vector<const zscalar_t *> weightPtrs(nWeights);
400 for (
int i=0; i < nWeights; i++)
401 weightPtrs[i] = weights[i].getRawPtr();
402 vector<int> weightStrides(0);
404 ia = rcp(
new inputAdapter_t (localCount, globalIds,
405 values, valueStrides, weightPtrs, weightStrides));
408 MEMORY_CHECK(doMemory && rank==0,
"After creating input adapter");
412 Teuchos::ParameterList params;
414 if (timingType !=
"no_timers"){
415 params.set(
"timer_output_stream" ,
"std::cout");
416 params.set(
"timer_type" , timingType);
420 params.set(
"memory_output_stream" ,
"std::cout");
421 params.set(
"memory_procs" , memoryProcs);
424 params.set(
"debug_output_stream" ,
"std::cerr");
425 params.set(
"debug_procs" ,
"0");
427 if (debugLevel !=
"basic_status"){
428 params.set(
"debug_level" , debugLevel);
431 params.set(
"algorithm",
"rcb");
432 params.set(
"partitioning_objective", objective);
433 double tolerance = 1.1;
434 params.set(
"imbalance_tolerance", tolerance );
436 if (numGlobalParts != nprocs)
437 params.set(
"num_global_parts" , numGlobalParts);
440 std::cout <<
"Number of parts: " << numGlobalParts << std::endl;
456 std::cout <<
"PASS" << std::endl;
int main(int narg, char **arg)
Defines the PartitioningSolution class.
common code used by tests
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
void printTimers() const
Return the communicator passed to the problem.
Tpetra::Map::local_ordinal_type zlno_t
PartitioningProblem sets up partitioning problems for the user.
Defines the PartitioningProblem class.
Defines the BasicVectorAdapter class.
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
Tpetra::Map::global_ordinal_type zgno_t
void solve(bool updateInputData=true)
Direct the problem to create a solution.
#define MEMORY_CHECK(iPrint, msg)