Actual source code: supportgraph.c
2: /* --------------------------------------------------------------------
4: This file implements a SupportGraph preconditioner in PETSc as part of PC.
5: You can use this as a starting point for implementing your own
6: preconditioner that is not provided with PETSc. (You might also consider
7: just using PCSHELL)
9: The following basic routines are required for each preconditioner.
10: PCCreate_XXX() - Creates a preconditioner context
11: PCSetFromOptions_XXX() - Sets runtime options
12: PCApply_XXX() - Applies the preconditioner
13: PCDestroy_XXX() - Destroys the preconditioner context
14: where the suffix "_XXX" denotes a particular implementation, in
15: this case we use _SupportGraph (e.g., PCCreate_SupportGraph, PCApply_SupportGraph).
16: These routines are actually called via the common user interface
17: routines PCCreate(), PCSetFromOptions(), PCApply(), and PCDestroy(),
18: so the application code interface remains identical for all
19: preconditioners.
21: Another key routine is:
22: PCSetUp_XXX() - Prepares for the use of a preconditioner
23: by setting data structures and options. The interface routine PCSetUp()
24: is not usually called directly by the user, but instead is called by
25: PCApply() if necessary.
27: Additional basic routines are:
28: PCView_XXX() - Prints details of runtime options that
29: have actually been used.
30: These are called by application codes via the interface routines
31: PCView().
33: The various types of solvers (preconditioners, Krylov subspace methods,
34: nonlinear solvers, timesteppers) are all organized similarly, so the
35: above description applies to these categories also. One exception is
36: that the analogues of PCApply() for these components are KSPSolve(),
37: SNESSolve(), and TSSolve().
39: Additional optional functionality unique to preconditioners is left and
40: right symmetric preconditioner application via PCApplySymmetricLeft()
41: and PCApplySymmetricRight(). The SupportGraph implementation is
42: PCApplySymmetricLeftOrRight_SupportGraph().
44: -------------------------------------------------------------------- */
46: /*
47: Include files needed for the SupportGraph preconditioner:
48: pcimpl.h - private include file intended for use by all preconditioners
49: adjacency_list.hpp
50: */
52: #include <private/pcimpl.h> /*I "petscpc.h" I*/
54: /*
55: Private context (data structure) for the SupportGraph preconditioner.
56: */
57: typedef struct {
58: Mat pre; /* Cholesky factored preconditioner matrix */
59: PetscBool augment; /* whether to augment the spanning tree */
60: PetscReal maxCong; /* create subgraph with at most this much congestion (only used with augment) */
61: PetscReal tol; /* throw out entries smaller than this */
62: } PC_SupportGraph;
66: static PetscErrorCode PCView_SupportGraph(PC pc,PetscViewer viewer)
67: {
68: PC_SupportGraph *sg = (PC_SupportGraph*)pc->data;
69: PetscErrorCode ierr;
70: PetscBool iascii;
73: PetscTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
74: if (iascii) {
75: PetscViewerASCIIPrintf(viewer," SupportGraph: maxCong = %f\n",sg->maxCong);
76: PetscViewerASCIIPrintf(viewer," SupportGraph: tol = %f\n",sg->tol);
77: PetscViewerASCIIPrintf(viewer," Factored Matrix:\n");
78: PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_INFO);
79: PetscViewerASCIIPushTab(viewer);
80: PetscViewerASCIIPushTab(viewer);
81: MatView(sg->pre, viewer);
82: PetscViewerASCIIPopTab(viewer);
83: PetscViewerASCIIPopTab(viewer);
84: PetscViewerPopFormat(viewer);
85: } else {
86: SETERRQ1(((PetscObject)pc)->comm,PETSC_ERR_SUP,"Viewer type %s not supported for PCSupportGraph",((PetscObject)viewer)->type_name);
87: }
88: return(0);
89: }
93: /* -------------------------------------------------------------------------- */
94: /*
95: PCSetUp_SupportGraph - Prepares for the use of the SupportGraph preconditioner
96: by setting data structures and options.
98: Input Parameter:
99: . pc - the preconditioner context
101: Application Interface Routine: PCSetUp()
103: Notes:
104: The interface routine PCSetUp() is not usually called directly by
105: the user, but instead is called by PCApply() if necessary.
106: */
109: static PetscErrorCode PCSetUp_SupportGraph(PC pc)
110: {
111: PC_SupportGraph *sg = (PC_SupportGraph*)pc->data;
112: PetscErrorCode ierr;
113: /*
114: Vec diag;
115: PetscInt n,i;
116: PetscScalar *x;
117: PetscBool zeroflag = PETSC_FALSE;
118: */
121: if(!pc->setupcalled) {
122: if (!MatIsSymmetric(pc->pmat, 1.0e-9)) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_ARG_WRONG,"matrix must be symmetric");
123: /* note that maxCong is being updated */
124: AugmentedLowStretchSpanningTree(pc->pmat, &sg->pre, sg->augment, sg->tol, sg->maxCong);
125: }
126: return(0);
127: }
130: /* -------------------------------------------------------------------------- */
131: /*
132: PCApply_SupportGraph - Applies the SupportGraph preconditioner to a vector.
134: Input Parameters:
135: . pc - the preconditioner context
136: . x - input vector
138: Output Parameter:
139: . y - output vector
141: Application Interface Routine: PCApply()
142: */
145: static PetscErrorCode PCApply_SupportGraph(PC pc,Vec x,Vec y)
146: {
147: PC_SupportGraph *sg = (PC_SupportGraph*)pc->data;
151: MatSolve(sg->pre,x,y);
152: return(0);
153: }
154: /* -------------------------------------------------------------------------- */
155: /*
156: PCDestroy_SupportGraph - Destroys the private context for the SupportGraph preconditioner
157: that was created with PCCreate_SupportGraph().
159: Input Parameter:
160: . pc - the preconditioner context
162: Application Interface Routine: PCDestroy()
163: */
166: static PetscErrorCode PCDestroy_SupportGraph(PC pc)
167: {
168: PC_SupportGraph *sg = (PC_SupportGraph*)pc->data;
169: PetscErrorCode ierr;
172: MatDestroy(&sg->pre);
173: /*
174: Free the private data structure that was hanging off the PC
175: */
176: PetscFree(pc->data);
177: return(0);
178: }
182: static PetscErrorCode PCSetFromOptions_SupportGraph(PC pc)
183: {
184: PC_SupportGraph *sg = (PC_SupportGraph*)pc->data;
188: PetscOptionsHead("SupportGraph options");
189: PetscOptionsBool("-pc_sg_augment","Max congestion","",sg->augment,&sg->augment,0);
190: PetscOptionsReal("-pc_sg_cong","Max congestion","",sg->maxCong,&sg->maxCong,0);
191: PetscOptionsReal("-pc_sg_tol","Smallest usable value","",sg->tol,&sg->tol,0);
192: PetscOptionsTail();
193: return(0);
194: }
196: /* -------------------------------------------------------------------------- */
197: /*
198: PCCreate_SupportGraph - Creates a SupportGraph preconditioner context, PC_SupportGraph,
199: and sets this as the private data within the generic preconditioning
200: context, PC, that was created within PCCreate().
202: Input Parameter:
203: . pc - the preconditioner context
205: Application Interface Routine: PCCreate()
206: */
208: /*MC
209: PCSUPPORTGRAPH - SupportGraph (i.e. diagonal scaling preconditioning)
211: Options Database Key:
212: . -pc_supportgraph_augment - augment the spanning tree
214: Level: beginner
216: Concepts: SupportGraph, diagonal scaling, preconditioners
218: Notes: Zero entries along the diagonal are replaced with the value 1.0
220: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC
221: M*/
226: PetscErrorCode PCCreate_SupportGraph(PC pc)
227: {
228: PC_SupportGraph *sg;
232: /*
233: Creates the private data structure for this preconditioner and
234: attach it to the PC object.
235: */
236: PetscNewLog(pc,PC_SupportGraph,&sg);
237: pc->data = (void*)sg;
239: sg->pre = 0;
240: sg->augment = PETSC_TRUE;
241: sg->maxCong = 3.0;
242: sg->tol = 0;
243:
245: /*
246: Set the pointers for the functions that are provided above.
247: Now when the user-level routines (such as PCApply(), PCDestroy(), etc.)
248: are called, they will automatically call these functions. Note we
249: choose not to provide a couple of these functions since they are
250: not needed.
251: */
252: pc->ops->apply = PCApply_SupportGraph;
253: pc->ops->applytranspose = 0;
254: pc->ops->setup = PCSetUp_SupportGraph;
255: pc->ops->destroy = PCDestroy_SupportGraph;
256: pc->ops->setfromoptions = PCSetFromOptions_SupportGraph;
257: pc->ops->view = PCView_SupportGraph;
258: pc->ops->applyrichardson = 0;
259: pc->ops->applysymmetricleft = 0;
260: pc->ops->applysymmetricright = 0;
261: return(0);
262: }