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, No Fortran Support

 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 `MatCoarsenSetType(agg, "my_agg")` or at runtime via the option `-mat_coarsen_type my_agg`

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

 37: /*@
 38:   MatCoarsenGetType - Gets the Coarsen method type and name (as a string)
 39:   from the coarsen context.

 41:   Not Collective

 43:   Input Parameter:
 44: . coarsen - the coarsen context

 46:   Output Parameter:
 47: . type - coarsener type

 49:   Level: advanced

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

 62: /*@
 63:   MatCoarsenApply - Gets a coarsen for a matrix.

 65:   Collective

 67:   Input Parameter:
 68: . coarser - the coarsen

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

 74:   Level: advanced

 76:   Notes:
 77:   When the coarsening is used inside `PCGAMG` then the options database keys are prefixed with `-pc_gamg_`

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

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

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

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

103:   Collective

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

109:   Level: advanced

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

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

125:   Logically Collective

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

131:   Level: advanced

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

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

146:   Collective

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

151:   Level: advanced

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

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

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

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

178:   Collective

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

185:   Options Database Key:
186: . -name [viewertype][:...] - option name and values. See `PetscObjectViewFromOptions()` for the possible arguments

188:   Level: intermediate

190: .seealso: `MatCoarsen`, `MatCoarsenView()`, `PetscObjectViewFromOptions()`, `MatCoarsenCreate()`
191: @*/
192: PetscErrorCode MatCoarsenViewFromOptions(MatCoarsen A, PetscObject obj, const char name[])
193: {
194:   PetscFunctionBegin;
196:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
197:   PetscFunctionReturn(PETSC_SUCCESS);
198: }

200: /*@
201:   MatCoarsenView - Prints the coarsen data structure.

203:   Collective

205:   Input Parameters:
206: + agg    - the coarsen context
207: - viewer - optional visualization context

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

211:   Level: advanced

213: .seealso: `MatCoarsen`, `PetscViewer`, `PetscViewerASCIIOpen()`, `MatCoarsenViewFromOptions`
214: @*/
215: PetscErrorCode MatCoarsenView(MatCoarsen agg, PetscViewer viewer)
216: {
217:   PetscBool isascii;

219:   PetscFunctionBegin;
221:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)agg), &viewer));
223:   PetscCheckSameComm(agg, 1, viewer, 2);

225:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
226:   PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)agg, viewer));
227:   if (agg->ops->view) {
228:     PetscCall(PetscViewerASCIIPushTab(viewer));
229:     PetscUseTypeMethod(agg, view, viewer);
230:     PetscCall(PetscViewerASCIIPopTab(viewer));
231:   }
232:   if (agg->strength_index_size > 0) PetscCall(PetscViewerASCIIPrintf(viewer, " Using scalar strength-of-connection index[%" PetscInt_FMT "] = {%" PetscInt_FMT ", ..}\n", agg->strength_index_size, agg->strength_index[0]));
233:   PetscFunctionReturn(PETSC_SUCCESS);
234: }

236: /*@
237:   MatCoarsenSetType - Sets the type of aggregator to use

239:   Collective

241:   Input Parameters:
242: + coarser - the coarsen context.
243: - type    - a known coarsening method

245:   Options Database Key:
246: . -mat_coarsen_type  type - maximal independent set based; distance k MIS; heavy edge matching

248:   Level: advanced

250: .seealso: `MatCoarsen`, `MatCoarsenCreate()`, `MatCoarsenApply()`, `MatCoarsenType`, `MatCoarsenGetType()`
251: @*/
252: PetscErrorCode MatCoarsenSetType(MatCoarsen coarser, MatCoarsenType type)
253: {
254:   PetscBool match;
255:   PetscErrorCode (*r)(MatCoarsen);

257:   PetscFunctionBegin;
259:   PetscAssertPointer(type, 2);

261:   PetscCall(PetscObjectTypeCompare((PetscObject)coarser, type, &match));
262:   if (match) PetscFunctionReturn(PETSC_SUCCESS);

264:   PetscTryTypeMethod(coarser, destroy);
265:   coarser->ops->destroy = NULL;
266:   PetscCall(PetscMemzero(coarser->ops, sizeof(struct _MatCoarsenOps)));

268:   PetscCall(PetscFunctionListFind(MatCoarsenList, type, &r));
269:   PetscCheck(r, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown coarsen type %s", type);
270:   PetscCall((*r)(coarser));

272:   PetscCall(PetscFree(((PetscObject)coarser)->type_name));
273:   PetscCall(PetscStrallocpy(type, &((PetscObject)coarser)->type_name));
274:   PetscFunctionReturn(PETSC_SUCCESS);
275: }

277: /*@
278:   MatCoarsenSetGreedyOrdering - Sets the ordering of the vertices to use with a greedy coarsening method

280:   Logically Collective

282:   Input Parameters:
283: + coarser - the coarsen context
284: - perm    - vertex ordering of (greedy) algorithm

286:   Level: advanced

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

291: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
292: @*/
293: PetscErrorCode MatCoarsenSetGreedyOrdering(MatCoarsen coarser, const IS perm)
294: {
295:   PetscFunctionBegin;
297:   coarser->perm = perm;
298:   PetscFunctionReturn(PETSC_SUCCESS);
299: }

301: /*@C
302:   MatCoarsenGetData - Gets the weights for vertices for a coarsener.

304:   Logically Collective, No Fortran Support

306:   Input Parameter:
307: . coarser - the coarsen context

309:   Output Parameter:
310: . llist - linked list of aggregates

312:   Level: advanced

314:   Note:
315:   This passes ownership to the caller and nullifies the value of weights (`PetscCoarsenData`) within the `MatCoarsen`

317: .seealso: `MatCoarsen`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`, `PetscCoarsenData`
318: @*/
319: PetscErrorCode MatCoarsenGetData(MatCoarsen coarser, PetscCoarsenData **llist)
320: {
321:   PetscFunctionBegin;
323:   PetscCheck(coarser->agg_lists, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_WRONGSTATE, "No linked list - generate it or call ApplyCoarsen");
324:   *llist             = coarser->agg_lists;
325:   coarser->agg_lists = NULL; /* giving up ownership */
326:   PetscFunctionReturn(PETSC_SUCCESS);
327: }

329: /*@
330:   MatCoarsenSetFromOptions - Sets various coarsen options from the options database.

332:   Collective

334:   Input Parameter:
335: . coarser - the coarsen context.

337:   Options Database Key:
338: + -mat_coarsen_type  type - see `MatCoarsenType`
339: - -mat_coarsen_max_it its - number of iterations to use in the coarsening process, see `MatCoarsenSetMaximumIterations()`

341:   Level: advanced

343:   Notes:
344:   When the coarsening is used inside `PCGAMG` then the options database keys are prefixed with `-pc_gamg_`

346:   Sets the `MatCoarsenType` to `MATCOARSENMISK` if has not been set previously

348: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`,
349:           `MatCoarsenSetMaximumIterations()`
350: @*/
351: PetscErrorCode MatCoarsenSetFromOptions(MatCoarsen coarser)
352: {
353:   PetscBool   flag;
354:   char        type[256];
355:   const char *def;

357:   PetscFunctionBegin;
358:   PetscObjectOptionsBegin((PetscObject)coarser);
359:   if (!((PetscObject)coarser)->type_name) {
360:     def = MATCOARSENMISK;
361:   } else {
362:     def = ((PetscObject)coarser)->type_name;
363:   }
364:   PetscCall(PetscOptionsFList("-mat_coarsen_type", "Type of aggregator", "MatCoarsenSetType", MatCoarsenList, def, type, 256, &flag));
365:   if (flag) PetscCall(MatCoarsenSetType(coarser, type));

367:   PetscCall(PetscOptionsInt("-mat_coarsen_max_it", "Number of iterations (for HEM)", "MatCoarsenSetMaximumIterations", coarser->max_it, &coarser->max_it, NULL));
368:   PetscCall(PetscOptionsInt("-mat_coarsen_threshold", "Threshold (for HEM)", "MatCoarsenSetThreshold", coarser->max_it, &coarser->max_it, NULL));
369:   coarser->strength_index_size = MAT_COARSEN_STRENGTH_INDEX_SIZE;
370:   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));
371:   /*
372:    Set the type if it was never set.
373:    */
374:   if (!((PetscObject)coarser)->type_name) PetscCall(MatCoarsenSetType(coarser, def));

376:   PetscTryTypeMethod(coarser, setfromoptions, PetscOptionsObject);
377:   PetscOptionsEnd();
378:   PetscFunctionReturn(PETSC_SUCCESS);
379: }

381: /*@
382:   MatCoarsenSetMaximumIterations - Maximum `MATCOARSENHEM` iterations to use

384:   Logically Collective

386:   Input Parameters:
387: + coarse - the coarsen context
388: - n      - number of HEM iterations

390:   Options Database Key:
391: . -mat_coarsen_max_it n - Maximum `MATCOARSENHEM` iterations to use

393:   Level: intermediate

395:   Note:
396:   When the coarsening is used inside `PCGAMG` then the options database keys are prefixed with `-pc_gamg_`

398: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
399: @*/
400: PetscErrorCode MatCoarsenSetMaximumIterations(MatCoarsen coarse, PetscInt n)
401: {
402:   PetscFunctionBegin;
405:   PetscTryMethod(coarse, "MatCoarsenSetMaximumIterations_C", (MatCoarsen, PetscInt), (coarse, n));
406:   PetscFunctionReturn(PETSC_SUCCESS);
407: }

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

416: /*@
417:   MatCoarsenSetStrengthIndex -  Index array to use for index to use for strength of connection

419:   Logically Collective

421:   Input Parameters:
422: + coarse - the coarsen context
423: . n      - number of indices
424: - idx    - array of indices

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

429:   Level: intermediate

431:   Note:
432:   When the coarsening is used inside `PCGAMG` then the options database keys are prefixed with `-pc_gamg_`

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 - Set the threshold for HEM

456:   Logically Collective

458:   Input Parameters:
459: + coarse - the coarsen context
460: - b      - threshold value, default is -1

462:   Options Database Key:
463: . -mat_coarsen_threshold b - threshold

465:   Level: intermediate

467:   Note:
468:   When the coarsening is used inside `PCGAMG` then the options database keys are prefixed with `-pc_gamg_`

470:   Developer Note:
471:   It is not documented how this threshold is used

473: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
474: @*/
475: PetscErrorCode MatCoarsenSetThreshold(MatCoarsen coarse, PetscReal b)
476: {
477:   PetscFunctionBegin;
480:   PetscTryMethod(coarse, "MatCoarsenSetThreshold_C", (MatCoarsen, PetscReal), (coarse, b));
481:   PetscFunctionReturn(PETSC_SUCCESS);
482: }

484: static PetscErrorCode MatCoarsenSetThreshold_MATCOARSEN(MatCoarsen coarse, PetscReal b)
485: {
486:   PetscFunctionBegin;
487:   coarse->threshold = b;
488:   PetscFunctionReturn(PETSC_SUCCESS);
489: }

491: /*@
492:   MatCoarsenCreate - Creates a coarsen context.

494:   Collective

496:   Input Parameter:
497: . comm - MPI communicator

499:   Output Parameter:
500: . newcrs - location to put the context

502:   Level: advanced

504: .seealso: `MatCoarsen`, `MatCoarsenSetType()`, `MatCoarsenApply()`, `MatCoarsenDestroy()`,
505:           `MatCoarsenSetAdjacency()`, `MatCoarsenGetData()`
506: @*/
507: PetscErrorCode MatCoarsenCreate(MPI_Comm comm, MatCoarsen *newcrs)
508: {
509:   MatCoarsen agg;

511:   PetscFunctionBegin;
512:   PetscAssertPointer(newcrs, 2);
513:   PetscCall(MatInitializePackage());

515:   PetscCall(PetscHeaderCreate(agg, MAT_COARSEN_CLASSID, "MatCoarsen", "Matrix/graph coarsen", "MatCoarsen", comm, MatCoarsenDestroy, MatCoarsenView));
516:   PetscCall(PetscObjectComposeFunction((PetscObject)agg, "MatCoarsenSetMaximumIterations_C", MatCoarsenSetMaximumIterations_MATCOARSEN));
517:   PetscCall(PetscObjectComposeFunction((PetscObject)agg, "MatCoarsenSetThreshold_C", MatCoarsenSetThreshold_MATCOARSEN));
518:   PetscCall(PetscObjectComposeFunction((PetscObject)agg, "MatCoarsenSetStrengthIndex_C", MatCoarsenSetStrengthIndex_MATCOARSEN));
519:   agg->strength_index_size = 0;
520:   *newcrs                  = agg;
521:   PetscFunctionReturn(PETSC_SUCCESS);
522: }