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: }