Actual source code: bicgstabcusp.cu
2: /* -------------------------------------------------------------------- */
4: /*
5: Include files needed for the CUSP BiCGSTAB preconditioner:
6: pcimpl.h - private include file intended for use by all preconditioners
7: */
9: #include <private/pcimpl.h> /*I "petscpc.h" I*/
10: #include <../src/mat/impls/aij/seq/aij.h>
11: #include <cusp/monitor.h>
12: #undef VecType
13: #include <cusp/krylov/bicgstab.h>
14: #define VecType char*
15: #include <../src/vec/vec/impls/dvecimpl.h>
16: #include <../src/mat/impls/aij/seq/seqcusp/cuspmatimpl.h>
19: /*
20: Private context (data structure) for the CUSP BiCGStab preconditioner.
21: */
22: typedef struct {
23: PetscInt maxits;
24: PetscReal rtol;
25: PetscBool monitorverbose;
26: CUSPMATRIX* mat;
27: } PC_BiCGStabCUSP;
31: static PetscErrorCode PCBiCGStabCUSPSetTolerance_BiCGStabCUSP(PC pc,PetscReal rtol)
32: {
33: PC_BiCGStabCUSP *bicg = (PC_BiCGStabCUSP*)pc->data;
36: bicg->rtol = rtol;
37: return(0);
38: }
42: static PetscErrorCode PCBiCGStabCUSPSetUseVerboseMonitor_BiCGStabCUSP(PC pc, PetscBool useverbose)
43: {
44: PC_BiCGStabCUSP *bicg = (PC_BiCGStabCUSP*)pc->data;
47: bicg->monitorverbose = useverbose;
48: return(0);
49: }
53: PetscErrorCode PCBiCGStabCUSPSetUseVerboseMonitor(PC pc, PetscBool useverbose)
54: {
59: PetscTryMethod(pc, "PCBiCGStabCUSPSetUseVerboseMonitors_C",(PC,PetscBool),(pc,useverbose));
60: return(0);
61: }
65: static PetscErrorCode PCBiCGStabCUSPSetIterations_BiCGStabCUSP(PC pc, PetscInt its)
66: {
67: PC_BiCGStabCUSP *bicg = (PC_BiCGStabCUSP*)pc->data;
70: bicg->maxits = its;
71: return(0);
72: }
76: PetscErrorCode PCBiCGStabCUSPSetITerations(PC pc, PetscInt its)
77: {
82: PetscTryMethod(pc, "PCBiCGStabCUSPSetIterations_C",(PC,PetscInt),(pc,its));
83: return(0);
84: }
88: PetscErrorCode PCBiCGStabCUSPSetTolerance(PC pc, PetscReal rtol)
89: {
94: PetscTryMethod(pc, "PCBiCGStabCUSPSetTolerance_C",(PC,PetscReal),(pc,rtol));
95: return(0);
96: }
98: /* -------------------------------------------------------------------------- */
99: /*
100: PCSetUp_BiCGStabCUSP - Prepares for the use of the CUSP BiCGStab preconditioner
101: by setting data structures and options.
103: Input Parameter:
104: . pc - the preconditioner context
106: Application Interface Routine: PCSetUp()
108: Notes:
109: The interface routine PCSetUp() is not usually called directly by
110: the user, but instead is called by PCApply() if necessary.
111: */
114: static PetscErrorCode PCSetUp_BiCGStabCUSP(PC pc)
115: {
116: PC_BiCGStabCUSP *bicg = (PC_BiCGStabCUSP*)pc->data;
117: PetscBool flg = PETSC_FALSE;
118: Mat_SeqAIJCUSP *gpustruct;
119: PetscErrorCode ierr;
122: PetscTypeCompare((PetscObject)pc->pmat,MATSEQAIJCUSP,&flg);
123: if (!flg) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_SUP,"Currently only handles CUSP matrices");
124: try{
125: MatCUSPCopyToGPU(pc->pmat);CHKERRCUSP(ierr);
126: gpustruct = (Mat_SeqAIJCUSP *)(pc->pmat->spptr);
127: bicg->mat = (CUSPMATRIX*)gpustruct->mat;
128: } catch(char* ex) {
129: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"CUSP error: %s",ex);
130: }
131: return(0);
132: }
134: /* -------------------------------------------------------------------------- */
135: /*
136: PCApply_BiCGStabCUSP - Applies the BiCGStabCUSP preconditioner to a vector.
138: Input Parameters:
139: . pc - the preconditioner context
140: . x - input vector
142: Output Parameter:
143: . y - output vector
145: Application Interface Routine: PCApply()
146: */
149: static PetscErrorCode PCApply_BiCGStabCUSP(PC pc,Vec x,Vec y)
150: {
151: PC_BiCGStabCUSP *bicg = (PC_BiCGStabCUSP*)pc->data;
152: PetscErrorCode ierr;
153: PetscBool flg1,flg2;
154: CUSPARRAY *xarray,*yarray;
157: PetscTypeCompare((PetscObject)x,VECSEQCUSP,&flg1);
158: PetscTypeCompare((PetscObject)y,VECSEQCUSP,&flg2);
159: if (!(flg1 && flg2)) SETERRQ(((PetscObject)pc)->comm,PETSC_ERR_SUP, "Currently only handles CUSP vectors");
160: if (!bicg->mat) {
161: PCSetUp_BiCGStabCUSP(pc);
162: }
163: VecSet(y,0.0);
164: VecCUSPGetArrayRead(x,&xarray);
165: VecCUSPGetArrayWrite(y,&yarray);
166: try {
167: cusp::default_monitor<PetscScalar> monitor(*xarray,bicg->maxits,bicg->rtol);
168: if (bicg->monitorverbose){
169: cusp::verbose_monitor<PetscScalar> verbosemonitor(*xarray,bicg->maxits,bicg->rtol);
170: cusp::krylov::bicgstab(*bicg->mat,*yarray,*xarray,verbosemonitor);
171: } else {
172: cusp::krylov::bicgstab(*bicg->mat,*yarray,*xarray,monitor);
173: }
174: } catch(char* ex) {
175: SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"CUSP error: %s", ex);
176: }
177: VecCUSPRestoreArrayRead(x,&xarray);
178: VecCUSPRestoreArrayWrite(y,&yarray);
179: PetscObjectStateIncrease((PetscObject)y);
180: return(0);
181: }
182: /* -------------------------------------------------------------------------- */
183: /*
184: PCDestroy_BiCGStabCUSP - Destroys the private context for the BiCGStabCUSP preconditioner
185: that was created with PCCreate_BiCGStabCUSP().
187: Input Parameter:
188: . pc - the preconditioner context
190: Application Interface Routine: PCDestroy()
191: */
194: static PetscErrorCode PCDestroy_BiCGStabCUSP(PC pc)
195: {
196: PetscErrorCode ierr;
199: /*
200: Free the private data structure that was hanging off the PC
201: */
202: PetscFree(pc->data);
203: return(0);
204: }
208: static PetscErrorCode PCSetFromOptions_BiCGStabCUSP(PC pc)
209: {
210: PC_BiCGStabCUSP *bicg = (PC_BiCGStabCUSP*)pc->data;
211: PetscErrorCode ierr;
214: PetscOptionsHead("BiCGStabCUSP options");
215: PetscOptionsReal("-pc_bicgstabcusp_rtol","relative tolerance for BiCGStabCUSP preconditioner","PCBiCGStabCUSPSetTolerance",bicg->rtol,&bicg->rtol,0);
216: PetscOptionsInt("-pc_bicgstabcusp_max_it","maximum iterations for BiCGStabCUSP preconditioner","PCBiCGStabCUSPSetIterations",bicg->maxits,&bicg->maxits,0);
217: PetscOptionsBool("-pc_bicgstabcusp_monitor_verbose","Print information about GPU BiCGStabCUSP iterations","PCBiCGStabCUSPSetUseVerboseMonitor",bicg->monitorverbose,&bicg->monitorverbose,0);
218: PetscOptionsTail();
219: return(0);
220: }
222: /* -------------------------------------------------------------------------- */
228: PetscErrorCode PCCreate_BiCGStabCUSP(PC pc)
229: {
230: PC_BiCGStabCUSP *bicg;
231: PetscErrorCode ierr;
234: /*
235: Creates the private data structure for this preconditioner and
236: attach it to the PC object.
237: */
238: PetscNewLog(pc,PC_BiCGStabCUSP,&bicg);
239: /*
240: Set default values. We don't actually want to set max iterations as far as I know, but the Cusp monitor requires them so we use a large number.
241: */
242: bicg->maxits = 1000;
243: bicg->rtol = 1.e-1;
244: bicg->monitorverbose = PETSC_FALSE;
245: pc->data = (void*)bicg;
246: /*
247: Set the pointers for the functions that are provided above.
248: Now when the user-level routines (such as PCApply(), PCDestroy(), etc.)
249: are called, they will automatically call these functions. Note we
250: choose not to provide a couple of these functions since they are
251: not needed.
252: */
253: pc->ops->apply = PCApply_BiCGStabCUSP;
254: pc->ops->applytranspose = 0;
255: pc->ops->setup = PCSetUp_BiCGStabCUSP;
256: pc->ops->destroy = PCDestroy_BiCGStabCUSP;
257: pc->ops->setfromoptions = PCSetFromOptions_BiCGStabCUSP;
258: pc->ops->view = 0;
259: pc->ops->applyrichardson = 0;
260: pc->ops->applysymmetricleft = 0;
261: pc->ops->applysymmetricright = 0;
262: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCBiCGStabCUSPSetTolerance_C","PCBiCGStabCUSPSetTolerance_BiCGStabCUSP",PCBiCGStabCUSPSetTolerance_BiCGStabCUSP);
263: PetscObjectComposeFunctionDynamic((PetscObject)pc, "PCBiCGStabCUSPSetIterations_C","PCBiCGStabCUSPSetIterations_BiCGStabCUSP", PCBiCGStabCUSPSetIterations_BiCGStabCUSP);
264: PetscObjectComposeFunctionDynamic((PetscObject)pc, "PCBiCGStabCUSPSetUseVerboseMonitor_C", "PCBiCGStabCUSPSetUseVerboseMonitor_BiCGStabCUSP", PCBiCGStabCUSPSetUseVerboseMonitor_BiCGStabCUSP);
265: return(0);
266: }