Actual source code: petscmatcoarsen.h
1: #pragma once
3: #include <petscmat.h>
5: /* SUBMANSEC = Mat */
7: PETSC_EXTERN PetscFunctionList MatCoarsenList;
9: /*S
10: MatCoarsen - Object for managing the coarsening of a graph (symmetric matrix)
12: Level: advanced
14: Note:
15: This is used by the `PCGAMG` to generate coarser representations of an algebraic problem
17: .seealso: [](ch_matrices), [](sec_graph), `Mat`, `MatCoarsenCreate()`, `MatCoarsenType`, `MatColoringType`, `MatPartitioningType`, `MatOrderingType`
18: `MatColoring`, `MatPartitioning`
19: S*/
20: typedef struct _p_MatCoarsen *MatCoarsen;
22: /*J
23: MatCoarsenType - String with the name of a PETSc matrix coarsening algorithm
25: Level: beginner
27: .seealso: [](ch_matrices), [](sec_graph), `Mat`, `MatCoarsenCreate()`, `MatCoarsen`, `MatColoringType`, `MatPartitioningType`, `MatOrderingType`
28: J*/
29: typedef const char *MatCoarsenType;
30: #define MATCOARSENMIS "mis"
31: #define MATCOARSENHEM "hem"
32: #define MATCOARSENMISK "misk"
34: /* linked list for aggregates */
35: typedef struct _PetscCDIntNd {
36: struct _PetscCDIntNd *next;
37: PetscInt gid;
38: } PetscCDIntNd;
40: /* only used by node pool */
41: typedef struct _PetscCDArrNd {
42: struct _PetscCDArrNd *next;
43: struct _PetscCDIntNd *array;
44: } PetscCDArrNd;
46: /* linked list data structure that encodes aggregates and C-F points with array[idx] == NULL for F point and array of indices in an aggregate or C point (first index is always global index my0 + idx */
47: typedef struct _PetscCoarsenData {
48: PetscCDArrNd pool_list; /* node pool */
49: PetscCDIntNd *new_node;
50: PetscInt new_left;
51: PetscInt chk_sz; /* chunk size */
52: PetscCDIntNd *extra_nodes;
53: PetscCDIntNd **array; /* Array of lists */
54: PetscInt size; /* size of 'array' */
55: Mat mat; /* cache a Mat for communication data */
56: } PetscCoarsenData;
58: PETSC_EXTERN PetscErrorCode MatCoarsenCreate(MPI_Comm, MatCoarsen *);
59: PETSC_EXTERN PetscErrorCode MatCoarsenSetType(MatCoarsen, MatCoarsenType);
60: PETSC_EXTERN PetscErrorCode MatCoarsenSetAdjacency(MatCoarsen, Mat);
61: PETSC_EXTERN PetscErrorCode MatCoarsenSetGreedyOrdering(MatCoarsen, const IS);
62: PETSC_EXTERN PetscErrorCode MatCoarsenSetStrictAggs(MatCoarsen, PetscBool);
63: PETSC_EXTERN PetscErrorCode MatCoarsenGetData(MatCoarsen, PetscCoarsenData **);
64: PETSC_EXTERN PetscErrorCode MatCoarsenApply(MatCoarsen);
65: PETSC_EXTERN PetscErrorCode MatCoarsenDestroy(MatCoarsen *);
66: PETSC_EXTERN PetscErrorCode MatCoarsenRegister(const char[], PetscErrorCode (*)(MatCoarsen));
67: PETSC_EXTERN PetscErrorCode MatCoarsenView(MatCoarsen, PetscViewer);
68: PETSC_EXTERN PetscErrorCode MatCoarsenSetFromOptions(MatCoarsen);
69: PETSC_EXTERN PetscErrorCode MatCoarsenGetType(MatCoarsen, MatCoarsenType *);
70: PETSC_EXTERN PetscErrorCode MatCoarsenViewFromOptions(MatCoarsen, PetscObject, const char[]);
72: PETSC_EXTERN PetscErrorCode MatCoarsenMISKSetDistance(MatCoarsen, PetscInt);
73: PETSC_EXTERN PetscErrorCode MatCoarsenMISKGetDistance(MatCoarsen, PetscInt *);
74: PETSC_EXTERN PetscErrorCode MatCoarsenSetMaximumIterations(MatCoarsen, PetscInt);
75: PETSC_EXTERN PetscErrorCode MatCoarsenSetThreshold(MatCoarsen, PetscReal);
76: PETSC_EXTERN PetscErrorCode MatCoarsenSetStrengthIndex(MatCoarsen, PetscInt, PetscInt[]);