Actual source code: coarsen.c


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

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

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

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

 13:    Logically Collective

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

 19:    Level: developer

 21:    Sample usage:
 22: .vb
 23:    MatCoarsenRegister("my_agg",MyAggCreate);
 24: .ve

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

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

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

 45:    Not Collective

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

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

 53:    Level: advanced

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

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

 69:    Collective

 71:    Input Parameter:
 72: .   coarsen - the coarsen

 74:    Options Database Keys:
 75:    To specify the coarsen through the options database, use one of
 76:    the following
 77: $    -mat_coarsen_type mis|hem|misk
 78:    To see the coarsen result
 79: $    -mat_coarsen_view

 81:    Level: advanced

 83:    Notes:
 84:     Use `MatCoarsenGetData()` to access the results of the coarsening

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

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

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

108:    Collective

110:    Input Parameters:
111: +  agg - the coarsen context
112: -  adj - the adjacency matrix

114:    Level: advanced

116: .seealso: `MatCoarsen`, `MatCoarsenSetFromOptions()`, `Mat`, `MatCoarsenCreate()`, `MatCoarsenApply()`
117: @*/
118: PetscErrorCode MatCoarsenSetAdjacency(MatCoarsen agg, Mat adj)
119: {
120:   PetscFunctionBegin;
123:   agg->graph = adj;
124:   PetscFunctionReturn(PETSC_SUCCESS);
125: }

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

130:    Logically Collective

132:    Input Parameters:
133: +  agg - the coarsen context
134: -  str - `PETSC_TRUE` keep strict aggregates, `PETSC_FALSE` allow overlap
135:    Level: advanced

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

147: /*@
148:    MatCoarsenDestroy - Destroys the coarsen context.

150:    Collective

152:    Input Parameters:
153: .  agg - the coarsen context

155:    Level: advanced

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

169:   if ((*agg)->ops->destroy) PetscCall((*(*agg)->ops->destroy)((*agg)));

171:   if ((*agg)->agg_lists) PetscCall(PetscCDDestroy((*agg)->agg_lists));

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

177: /*@
178:    MatCoarsenCreate - Creates a coarsen context.

180:    Collective

182:    Input Parameter:
183: .   comm - MPI communicator

185:    Output Parameter:
186: .  newcrs - location to put the context

188:    Level: advanced

190: .seealso: `MatCoarsen`, `MatCoarsenSetType()`, `MatCoarsenApply()`, `MatCoarsenDestroy()`,
191:           `MatCoarsenSetAdjacency()`, `MatCoarsenGetData()`

193: @*/
194: PetscErrorCode MatCoarsenCreate(MPI_Comm comm, MatCoarsen *newcrs)
195: {
196:   MatCoarsen agg;

198:   PetscFunctionBegin;
199:   *newcrs = NULL;

201:   PetscCall(MatInitializePackage());
202:   PetscCall(PetscHeaderCreate(agg, MAT_COARSEN_CLASSID, "MatCoarsen", "Matrix/graph coarsen", "MatCoarsen", comm, MatCoarsenDestroy, MatCoarsenView));

204:   *newcrs = agg;
205:   PetscFunctionReturn(PETSC_SUCCESS);
206: }

208: /*@C
209:    MatCoarsenViewFromOptions - View the coarsener from the options database

211:    Collective

213:    Input Parameters:
214: +  A - the coarsen context
215: .  obj - Optional object that provides the prefix for the option name
216: -  name - command line option (usually `-mat_coarsen_view`)

218:   Options Database:
219: .  -mat_coarsen_view [viewertype]:... - the viewer and its options

221:   Note:
222: .vb
223:     If no value is provided ascii:stdout is used
224:        ascii[:[filename][:[format][:append]]]    defaults to stdout - format can be one of ascii_info, ascii_info_detail, or ascii_matlab,
225:                                                   for example ascii::ascii_info prints just the information about the object not all details
226:                                                   unless :append is given filename opens in write mode, overwriting what was already there
227:        binary[:[filename][:[format][:append]]]   defaults to the file binaryoutput
228:        draw[:drawtype[:filename]]                for example, draw:tikz, draw:tikz:figure.tex  or draw:x
229:        socket[:port]                             defaults to the standard output port
230:        saws[:communicatorname]                    publishes object to the Scientific Application Webserver (SAWs)
231: .ve

233:    Level: intermediate

235: .seealso: `MatCoarsen`, `MatCoarsenView`, `PetscObjectViewFromOptions()`, `MatCoarsenCreate()`
236: @*/
237: PetscErrorCode MatCoarsenViewFromOptions(MatCoarsen A, PetscObject obj, const char name[])
238: {
239:   PetscFunctionBegin;
241:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
242:   PetscFunctionReturn(PETSC_SUCCESS);
243: }

245: /*@C
246:    MatCoarsenView - Prints the coarsen data structure.

248:    Collective

250:    Input Parameters:
251: +  agg - the coarsen context
252: -  viewer - optional visualization context

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

256:    Level: advanced

258: .seealso: `MatCoarsen`, `PetscViewer`, `PetscViewerASCIIOpen()`, `MatCoarsenViewFromOptions`
259: @*/
260: PetscErrorCode MatCoarsenView(MatCoarsen agg, PetscViewer viewer)
261: {
262:   PetscBool iascii;

264:   PetscFunctionBegin;
266:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)agg), &viewer));
268:   PetscCheckSameComm(agg, 1, viewer, 2);

270:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
271:   PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)agg, viewer));
272:   if (agg->ops->view) {
273:     PetscCall(PetscViewerASCIIPushTab(viewer));
274:     PetscUseTypeMethod(agg, view, viewer);
275:     PetscCall(PetscViewerASCIIPopTab(viewer));
276:   }
277:   PetscFunctionReturn(PETSC_SUCCESS);
278: }

280: /*@C
281:    MatCoarsenSetType - Sets the type of aggregator to use

283:    Collective

285:    Input Parameters:
286: +  coarser - the coarsen context.
287: -  type - a known coarsening method

289:    Options Database Key:
290: .  -mat_coarsen_type  <type> - (for instance, misk), use -help for a list of available methods

292:    Level: advanced

294: .seealso: `MatCoarsen`, `MatCoarsenCreate()`, `MatCoarsenApply()`, `MatCoarsenType`, `MatCoarsenGetType()`
295: @*/
296: PetscErrorCode MatCoarsenSetType(MatCoarsen coarser, MatCoarsenType type)
297: {
298:   PetscBool match;
299:   PetscErrorCode (*r)(MatCoarsen);

301:   PetscFunctionBegin;

305:   PetscCall(PetscObjectTypeCompare((PetscObject)coarser, type, &match));
306:   if (match) PetscFunctionReturn(PETSC_SUCCESS);

308:   PetscTryTypeMethod(coarser, destroy);
309:   coarser->ops->destroy = NULL;
310:   PetscCall(PetscMemzero(coarser->ops, sizeof(struct _MatCoarsenOps)));

312:   PetscCall(PetscFunctionListFind(MatCoarsenList, type, &r));
313:   PetscCheck(r, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown coarsen type %s", type);
314:   PetscCall((*r)(coarser));

316:   PetscCall(PetscFree(((PetscObject)coarser)->type_name));
317:   PetscCall(PetscStrallocpy(type, &((PetscObject)coarser)->type_name));
318:   PetscFunctionReturn(PETSC_SUCCESS);
319: }

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

324:    Logically Collective

326:    Input Parameters:
327: +  coarser - the coarsen context
328: -  perm - vertex ordering of (greedy) algorithm

330:    Level: advanced

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

335: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
336: @*/
337: PetscErrorCode MatCoarsenSetGreedyOrdering(MatCoarsen coarser, const IS perm)
338: {
339:   PetscFunctionBegin;
341:   coarser->perm = perm;
342:   PetscFunctionReturn(PETSC_SUCCESS);
343: }

345: /*@C
346:    MatCoarsenGetData - Gets the weights for vertices for a coarsener.

348:    Logically Collective

350:    Input Parameter:
351: .  coarser - the coarsen context

353:    Output Parameter:
354: .  llist - linked list of aggregates

356:    Level: advanced

358: .seealso: `MatCoarsen`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
359: @*/
360: PetscErrorCode MatCoarsenGetData(MatCoarsen coarser, PetscCoarsenData **llist)
361: {
362:   PetscFunctionBegin;
364:   PetscCheck(coarser->agg_lists, PetscObjectComm((PetscObject)coarser), PETSC_ERR_ARG_WRONGSTATE, "No linked list - generate it or call ApplyCoarsen");
365:   *llist             = coarser->agg_lists;
366:   coarser->agg_lists = NULL; /* giving up ownership */
367:   PetscFunctionReturn(PETSC_SUCCESS);
368: }

370: /*@
371:    MatCoarsenSetFromOptions - Sets various coarsen options from the options database.

373:    Collective

375:    Input Parameter:
376: .  coarser - the coarsen context.

378:    Options Database Key:
379: .  -mat_coarsen_type  <type> - (for instance, mis), use -help for a list of available methods

381:    Level: advanced

383:    Note:
384:    Sets the `MatCoarsenType` to `MATCOARSENMISK` if has not been set previously

386: .seealso: `MatCoarsen`, `MatCoarsenType`, `MatCoarsenApply()`, `MatCoarsenCreate()`, `MatCoarsenSetType()`
387: @*/
388: PetscErrorCode MatCoarsenSetFromOptions(MatCoarsen coarser)
389: {
390:   PetscBool   flag;
391:   char        type[256];
392:   const char *def;

394:   PetscFunctionBegin;
395:   PetscObjectOptionsBegin((PetscObject)coarser);
396:   if (!((PetscObject)coarser)->type_name) {
397:     def = MATCOARSENMISK;
398:   } else {
399:     def = ((PetscObject)coarser)->type_name;
400:   }

402:   PetscCall(PetscOptionsFList("-mat_coarsen_type", "Type of aggregator", "MatCoarsenSetType", MatCoarsenList, def, type, 256, &flag));
403:   if (flag) PetscCall(MatCoarsenSetType(coarser, type));
404:   /*
405:    Set the type if it was never set.
406:    */
407:   if (!((PetscObject)coarser)->type_name) PetscCall(MatCoarsenSetType(coarser, def));

409:   PetscTryTypeMethod(coarser, setfromoptions, PetscOptionsObject);
410:   PetscOptionsEnd();

412:   PetscFunctionReturn(PETSC_SUCCESS);
413: }