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