Actual source code: coarsen.c

  1: #include <petsc/private/matimpl.h>

  3: /* Logging support */
  4: PetscClassId MAT_COARSEN_CLASSID;

  6: PetscFunctionList MatCoarsenList              = NULL;
  7: PetscBool         MatCoarsenRegisterAllCalled = PETSC_FALSE;

  9: /*@C
 10:   MatCoarsenRegister - Adds a new sparse matrix coarsening algorithm to the matrix package.

 12:   Logically Collective

 14:   Input Parameters:
 15: + sname    - name of coarsen (for example `MATCOARSENMIS`)
 16: - function - function pointer that creates the coarsen type

 18:   Level: developer

 20:   Example Usage:
 21: .vb
 22:    MatCoarsenRegister("my_agg", MyAggCreate);
 23: .ve

 25:   Then, your aggregator can be chosen with the procedural interface via
 26: $     MatCoarsenSetType(agg, "my_agg")
 27:   or at runtime via the option
 28: $     -mat_coarsen_type my_agg

 30: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenSetType()`, `MatCoarsenCreate()`, `MatCoarsenRegisterDestroy()`, `MatCoarsenRegisterAll()`
 31: @*/
 32: PetscErrorCode MatCoarsenRegister(const char sname[], PetscErrorCode (*function)(MatCoarsen))
 33: {
 34:   PetscFunctionBegin;
 35:   PetscCall(MatInitializePackage());
 36:   PetscCall(PetscFunctionListAdd(&MatCoarsenList, sname, function));
 37:   PetscFunctionReturn(PETSC_SUCCESS);
 38: }

 40: /*@C
 41:   MatCoarsenGetType - Gets the Coarsen method type and name (as a string)
 42:   from the coarsen context.

 44:   Not Collective

 46:   Input Parameter:
 47: . coarsen - the coarsen context

 49:   Output Parameter:
 50: . type - coarsener type

 52:   Level: advanced

 54: .seealso: `MatCoarsen`, `MatCoarsenCreate()`, `MatCoarsenType`, `MatCoarsenSetType()`, `MatCoarsenRegister()`
 55: @*/
 56: PetscErrorCode MatCoarsenGetType(MatCoarsen coarsen, MatCoarsenType *type)
 57: {
 58:   PetscFunctionBegin;
 60:   PetscAssertPointer(type, 2);
 61:   *type = ((PetscObject)coarsen)->type_name;
 62:   PetscFunctionReturn(PETSC_SUCCESS);
 63: }

 65: /*@
 66:   MatCoarsenApply - Gets a coarsen for a matrix.

 68:   Collective

 70:   Input Parameter:
 71: . coarser - the coarsen

 73:   Options Database Keys:
 74: + -mat_coarsen_type mis|hem|misk - mis: maximal independent set based; misk: distance k MIS; hem: heavy edge matching
 75: - -mat_coarsen_view              - view the coarsening object

 77:   Level: advanced

 79:   Notes:
 80:   Use `MatCoarsenGetData()` to access the results of the coarsening

 82:   The user can define additional coarsens; see `MatCoarsenRegister()`.

 84: .seealso: `MatCoarsen`, `MatCoarseSetFromOptions()`, `MatCoarsenSetType()`, `MatCoarsenRegister()`, `MatCoarsenCreate()`,
 85:           `MatCoarsenDestroy()`, `MatCoarsenSetAdjacency()`
 86:           `MatCoarsenGetData()`
 87: @*/
 88: PetscErrorCode MatCoarsenApply(MatCoarsen coarser)
 89: {
 90:   PetscFunctionBegin;
 92:   PetscAssertPointer(coarser, 1);
 93:   PetscCheck(coarser->graph->assembled, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_WRONGSTATE, "Not for unassembled matrix");
 94:   PetscCheck(!coarser->graph->factortype, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_WRONGSTATE, "Not for factored matrix");
 95:   PetscCall(PetscLogEventBegin(MAT_Coarsen, coarser, 0, 0, 0));
 96:   PetscUseTypeMethod(coarser, apply);
 97:   PetscCall(PetscLogEventEnd(MAT_Coarsen, coarser, 0, 0, 0));
 98:   PetscFunctionReturn(PETSC_SUCCESS);
 99: }

101: /*@
102:   MatCoarsenSetAdjacency - Sets the adjacency graph (matrix) of the thing to be coarsened.

104:   Collective

106:   Input Parameters:
107: + agg - the coarsen context
108: - adj - the adjacency matrix

110:   Level: advanced

112: .seealso: `MatCoarsen`, `MatCoarsenSetFromOptions()`, `Mat`, `MatCoarsenCreate()`, `MatCoarsenApply()`
113: @*/
114: PetscErrorCode MatCoarsenSetAdjacency(MatCoarsen agg, Mat adj)
115: {
116:   PetscFunctionBegin;
119:   agg->graph = adj;
120:   PetscFunctionReturn(PETSC_SUCCESS);
121: }

123: /*@
124:   MatCoarsenSetStrictAggs - Set whether to keep strict (non overlapping) aggregates in the linked list of aggregates for a coarsen context

126:   Logically Collective

128:   Input Parameters:
129: + agg - the coarsen context
130: - str - `PETSC_TRUE` keep strict aggregates, `PETSC_FALSE` allow overlap

132:   Level: advanced

134: .seealso: `MatCoarsen`, `MatCoarsenCreate()`, `MatCoarsenSetFromOptions()`
135: @*/
136: PetscErrorCode MatCoarsenSetStrictAggs(MatCoarsen agg, PetscBool str)
137: {
138:   PetscFunctionBegin;
140:   agg->strict_aggs = str;
141:   PetscFunctionReturn(PETSC_SUCCESS);
142: }

144: /*@
145:   MatCoarsenDestroy - Destroys the coarsen context.

147:   Collective

149:   Input Parameter:
150: . agg - the coarsen context

152:   Level: advanced

154: .seealso: `MatCoarsen`, `MatCoarsenCreate()`
155: @*/
156: PetscErrorCode MatCoarsenDestroy(MatCoarsen *agg)
157: {
158:   PetscFunctionBegin;
159:   if (!*agg) PetscFunctionReturn(PETSC_SUCCESS);
161:   if (--((PetscObject)*agg)->refct > 0) {
162:     *agg = NULL;
163:     PetscFunctionReturn(PETSC_SUCCESS);
164:   }

166:   PetscTryTypeMethod(*agg, destroy);
167:   if ((*agg)->agg_lists) PetscCall(PetscCDDestroy((*agg)->agg_lists));
168:   PetscCall(PetscObjectComposeFunction((PetscObject)*agg, "MatCoarsenSetMaximumIterations_C", NULL));
169:   PetscCall(PetscObjectComposeFunction((PetscObject)*agg, "MatCoarsenSetThreshold_C", NULL));
170:   PetscCall(PetscObjectComposeFunction((PetscObject)*agg, "MatCoarsenSetStrengthIndex_C", NULL));

172:   PetscCall(PetscHeaderDestroy(agg));
173:   PetscFunctionReturn(PETSC_SUCCESS);
174: }

176: /*@C
177:   MatCoarsenViewFromOptions - View the coarsener from the options database

179:   Collective

181:   Input Parameters:
182: + A    - the coarsen context
183: . obj  - Optional object that provides the prefix for the option name
184: - name - command line option (usually `-mat_coarsen_view`)

186:   Options Database Key:
187: . -mat_coarsen_view [viewertype]:... - the viewer and its options

189:   Note:
190: .vb
191:     If no value is provided ascii:stdout is used
192:        ascii[:[filename][:[format][:append]]]    defaults to stdout - format can be one of ascii_info, ascii_info_detail, or ascii_matlab,
193:                                                   for example ascii::ascii_info prints just the information about the object not all details
194:                                                   unless :append is given filename opens in write mode, overwriting what was already there
195:        binary[:[filename][:[format][:append]]]   defaults to the file binaryoutput
196:        draw[:drawtype[:filename]]                for example, draw:tikz, draw:tikz:figure.tex  or draw:x
197:        socket[:port]                             defaults to the standard output port
198:        saws[:communicatorname]                    publishes object to the Scientific Application Webserver (SAWs)
199: .ve

201:   Level: intermediate

203: .seealso: `MatCoarsen`, `MatCoarsenView`, `PetscObjectViewFromOptions()`, `MatCoarsenCreate()`
204: @*/
205: PetscErrorCode MatCoarsenViewFromOptions(MatCoarsen A, PetscObject obj, const char name[])
206: {
207:   PetscFunctionBegin;
209:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
210:   PetscFunctionReturn(PETSC_SUCCESS);
211: }

213: /*@C
214:   MatCoarsenView - Prints the coarsen data structure.

216:   Collective

218:   Input Parameters:
219: + agg    - the coarsen context
220: - viewer - optional visualization context

222:    For viewing the options database see `MatCoarsenViewFromOptions()`

224:   Level: advanced

226: .seealso: `MatCoarsen`, `PetscViewer`, `PetscViewerASCIIOpen()`, `MatCoarsenViewFromOptions`
227: @*/
228: PetscErrorCode MatCoarsenView(MatCoarsen agg, PetscViewer viewer)
229: {
230:   PetscBool iascii;

232:   PetscFunctionBegin;
234:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)agg), &viewer));
236:   PetscCheckSameComm(agg, 1, viewer, 2);

238:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
239:   PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)agg, viewer));
240:   if (agg->ops->view) {
241:     PetscCall(PetscViewerASCIIPushTab(viewer));
242:     PetscUseTypeMethod(agg, view, viewer);
243:     PetscCall(PetscViewerASCIIPopTab(viewer));
244:   }
245:   if (agg->strength_index_size > 0) PetscCall(PetscViewerASCIIPrintf(viewer, " Using scalar strength-of-connection index index[%d] = {%d, ..}\n", (int)agg->strength_index_size, (int)agg->strength_index[0]));
246:   PetscFunctionReturn(PETSC_SUCCESS);
247: }

249: /*@C
250:   MatCoarsenSetType - Sets the type of aggregator to use

252:   Collective

254:   Input Parameters:
255: + coarser - the coarsen context.
256: - type    - a known coarsening method

258:   Options Database Key:
259: . -mat_coarsen_type  <type> - maximal independent set based; distance k MIS; heavy edge matching

261:   Level: advanced

263: .seealso: `MatCoarsen`, `MatCoarsenCreate()`, `MatCoarsenApply()`, `MatCoarsenType`, `MatCoarsenGetType()`
264: @*/
265: PetscErrorCode MatCoarsenSetType(MatCoarsen coarser, MatCoarsenType type)
266: {
267:   PetscBool match;
268:   PetscErrorCode (*r)(MatCoarsen);

270:   PetscFunctionBegin;
272:   PetscAssertPointer(type, 2);

274:   PetscCall(PetscObjectTypeCompare((PetscObject)coarser, type, &match));
275:   if (match) PetscFunctionReturn(PETSC_SUCCESS);

277:   PetscTryTypeMethod(coarser, destroy);
278:   coarser->ops->destroy = NULL;
279:   PetscCall(PetscMemzero(coarser->ops, sizeof(struct _MatCoarsenOps)));

281:   PetscCall(PetscFunctionListFind(MatCoarsenList, type, &r));
282:   PetscCheck(r, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown coarsen type %s", type);
283:   PetscCall((*r)(coarser));

285:   PetscCall(PetscFree(((PetscObject)coarser)->type_name));
286:   PetscCall(PetscStrallocpy(type, &((PetscObject)coarser)->type_name));
287:   PetscFunctionReturn(PETSC_SUCCESS);
288: }

290: /*@C
291:   MatCoarsenSetGreedyOrdering - Sets the ordering of the vertices to use with a greedy coarsening method

293:   Logically Collective

295:   Input Parameters:
296: + coarser - the coarsen context
297: - perm    - vertex ordering of (greedy) algorithm

299:   Level: advanced

301:   Note:
302:   The `IS` weights is freed by PETSc, the user should not destroy it or change it after this call

304: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
305: @*/
306: PetscErrorCode MatCoarsenSetGreedyOrdering(MatCoarsen coarser, const IS perm)
307: {
308:   PetscFunctionBegin;
310:   coarser->perm = perm;
311:   PetscFunctionReturn(PETSC_SUCCESS);
312: }

314: /*@C
315:   MatCoarsenGetData - Gets the weights for vertices for a coarsener.

317:   Logically Collective

319:   Input Parameter:
320: . coarser - the coarsen context

322:   Output Parameter:
323: . llist - linked list of aggregates

325:   Level: advanced

327: .seealso: `MatCoarsen`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
328: @*/
329: PetscErrorCode MatCoarsenGetData(MatCoarsen coarser, PetscCoarsenData **llist)
330: {
331:   PetscFunctionBegin;
333:   PetscCheck(coarser->agg_lists, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_WRONGSTATE, "No linked list - generate it or call ApplyCoarsen");
334:   *llist             = coarser->agg_lists;
335:   coarser->agg_lists = NULL; /* giving up ownership */
336:   PetscFunctionReturn(PETSC_SUCCESS);
337: }

339: /*@
340:   MatCoarsenSetFromOptions - Sets various coarsen options from the options database.

342:   Collective

344:   Input Parameter:
345: . coarser - the coarsen context.

347:   Options Database Key:
348: . -mat_coarsen_type  <type> - mis: maximal independent set based; misk: distance k MIS; hem: heavy edge matching

350:   Level: advanced

352:   Note:
353:   Sets the `MatCoarsenType` to `MATCOARSENMISK` if has not been set previously

355: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
356: @*/
357: PetscErrorCode MatCoarsenSetFromOptions(MatCoarsen coarser)
358: {
359:   PetscBool   flag;
360:   char        type[256];
361:   const char *def;

363:   PetscFunctionBegin;
364:   PetscObjectOptionsBegin((PetscObject)coarser);
365:   if (!((PetscObject)coarser)->type_name) {
366:     def = MATCOARSENMISK;
367:   } else {
368:     def = ((PetscObject)coarser)->type_name;
369:   }
370:   PetscCall(PetscOptionsFList("-mat_coarsen_type", "Type of aggregator", "MatCoarsenSetType", MatCoarsenList, def, type, 256, &flag));
371:   if (flag) PetscCall(MatCoarsenSetType(coarser, type));

373:   PetscCall(PetscOptionsInt("-mat_coarsen_max_it", "Number of iterations (for HEM)", "MatCoarsenSetMaximumIterations", coarser->max_it, &coarser->max_it, NULL));
374:   PetscCall(PetscOptionsInt("-mat_coarsen_threshold", "Threshold (for HEM)", "MatCoarsenSetThreshold", coarser->max_it, &coarser->max_it, NULL));
375:   coarser->strength_index_size = MAT_COARSEN_STRENGTH_INDEX_SIZE;
376:   PetscCall(PetscOptionsIntArray("-mat_coarsen_strength_index", "Array of indices to use strength of connection measure (default is all indices)", "MatCoarsenSetStrengthIndex", coarser->strength_index, &coarser->strength_index_size, NULL));
377:   /*
378:    Set the type if it was never set.
379:    */
380:   if (!((PetscObject)coarser)->type_name) PetscCall(MatCoarsenSetType(coarser, def));

382:   PetscTryTypeMethod(coarser, setfromoptions, PetscOptionsObject);
383:   PetscOptionsEnd();
384:   PetscFunctionReturn(PETSC_SUCCESS);
385: }

387: /*@
388:   MatCoarsenSetMaximumIterations - Max HEM iterations

390:   Logically Collective

392:   Input Parameters:
393: + coarse - the coarsen context
394: - n      - number of HEM iterations

396:   Options Database Key:
397: . -mat_coarsen_max_it <default=4> - Max HEM iterations

399:   Level: intermediate

401: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
402: @*/
403: PetscErrorCode MatCoarsenSetMaximumIterations(MatCoarsen coarse, PetscInt n)
404: {
405:   PetscFunctionBegin;
408:   PetscTryMethod(coarse, "MatCoarsenSetMaximumIterations_C", (MatCoarsen, PetscInt), (coarse, n));
409:   PetscFunctionReturn(PETSC_SUCCESS);
410: }

412: static PetscErrorCode MatCoarsenSetMaximumIterations_MATCOARSEN(MatCoarsen coarse, PetscInt b)
413: {
414:   PetscFunctionBegin;
415:   coarse->max_it = b;
416:   PetscFunctionReturn(PETSC_SUCCESS);
417: }

419: /*@
420:   MatCoarsenSetStrengthIndex -  Index array to use for index to use for strength of connection

422:   Logically Collective

424:   Input Parameters:
425: + coarse - the coarsen context
426: . n      - number of indices
427: - idx    - array of indices

429:   Options Database Key:
430: . -mat_coarsen_strength_index - array of subset of variables per vertex to use for strength norm, -1 for using all (default)

432:   Level: intermediate

434: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
435: @*/
436: PetscErrorCode MatCoarsenSetStrengthIndex(MatCoarsen coarse, PetscInt n, PetscInt idx[])
437: {
438:   PetscFunctionBegin;
441:   PetscTryMethod(coarse, "MatCoarsenSetStrengthIndex_C", (MatCoarsen, PetscInt, PetscInt[]), (coarse, n, idx));
442:   PetscFunctionReturn(PETSC_SUCCESS);
443: }

445: static PetscErrorCode MatCoarsenSetStrengthIndex_MATCOARSEN(MatCoarsen coarse, PetscInt n, PetscInt idx[])
446: {
447:   PetscFunctionBegin;
448:   coarse->strength_index_size = n;
449:   for (int iii = 0; iii < n; iii++) coarse->strength_index[iii] = idx[iii];
450:   PetscFunctionReturn(PETSC_SUCCESS);
451: }

453: /*@
454:   MatCoarsenSetThreshold - Max HEM iterations

456:   Logically Collective

458:   Input Parameters:
459: + coarse - the coarsen context
460: - b      - threshold value

462:   Options Database Key:
463: . -mat_coarsen_threshold <-1> - Max HEM iterations

465:   Level: intermediate

467: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
468: @*/
469: PetscErrorCode MatCoarsenSetThreshold(MatCoarsen coarse, PetscReal b)
470: {
471:   PetscFunctionBegin;
474:   PetscTryMethod(coarse, "MatCoarsenSetThreshold_C", (MatCoarsen, PetscReal), (coarse, b));
475:   PetscFunctionReturn(PETSC_SUCCESS);
476: }

478: static PetscErrorCode MatCoarsenSetThreshold_MATCOARSEN(MatCoarsen coarse, PetscReal b)
479: {
480:   PetscFunctionBegin;
481:   coarse->threshold = b;
482:   PetscFunctionReturn(PETSC_SUCCESS);
483: }

485: /*@
486:   MatCoarsenCreate - Creates a coarsen context.

488:   Collective

490:   Input Parameter:
491: . comm - MPI communicator

493:   Output Parameter:
494: . newcrs - location to put the context

496:   Level: advanced

498: .seealso: `MatCoarsen`, `MatCoarsenSetType()`, `MatCoarsenApply()`, `MatCoarsenDestroy()`,
499:           `MatCoarsenSetAdjacency()`, `MatCoarsenGetData()`

501: @*/
502: PetscErrorCode MatCoarsenCreate(MPI_Comm comm, MatCoarsen *newcrs)
503: {
504:   MatCoarsen agg;

506:   PetscFunctionBegin;
507:   *newcrs = NULL;

509:   PetscCall(MatInitializePackage());
510:   PetscCall(PetscHeaderCreate(agg, MAT_COARSEN_CLASSID, "MatCoarsen", "Matrix/graph coarsen", "MatCoarsen", comm, MatCoarsenDestroy, MatCoarsenView));
511:   PetscCall(PetscObjectComposeFunction((PetscObject)agg, "MatCoarsenSetMaximumIterations_C", MatCoarsenSetMaximumIterations_MATCOARSEN));
512:   PetscCall(PetscObjectComposeFunction((PetscObject)agg, "MatCoarsenSetThreshold_C", MatCoarsenSetThreshold_MATCOARSEN));
513:   PetscCall(PetscObjectComposeFunction((PetscObject)agg, "MatCoarsenSetStrengthIndex_C", MatCoarsenSetStrengthIndex_MATCOARSEN));

515:   agg->strength_index_size = 0;

517:   *newcrs = agg;
518:   PetscFunctionReturn(PETSC_SUCCESS);
519: }