Actual source code: fieldsplit.c
1: #define PETSCKSP_DLL
3: /*
5: */
6: #include private/pcimpl.h
8: typedef struct _PC_FieldSplitLink *PC_FieldSplitLink;
9: struct _PC_FieldSplitLink {
10: KSP ksp;
11: Vec x,y;
12: PetscInt nfields;
13: PetscInt *fields;
14: VecScatter sctx;
15: PC_FieldSplitLink next,previous;
16: };
18: typedef struct {
19: PCCompositeType type; /* additive or multiplicative */
20: PetscTruth defaultsplit;
21: PetscInt bs;
22: PetscInt nsplits,*csize;
23: Vec *x,*y,w1,w2;
24: Mat *pmat;
25: IS *is,*cis;
26: PC_FieldSplitLink head;
27: } PC_FieldSplit;
31: static PetscErrorCode PCView_FieldSplit(PC pc,PetscViewer viewer)
32: {
33: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
34: PetscErrorCode ierr;
35: PetscTruth iascii;
36: PetscInt i,j;
37: PC_FieldSplitLink ilink = jac->head;
40: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);
41: if (iascii) {
42: PetscViewerASCIIPrintf(viewer," FieldSplit with %s composition: total splits = %D, blocksize = %D\n",PCCompositeTypes[jac->type],jac->nsplits,jac->bs);
43: PetscViewerASCIIPrintf(viewer," Solver info for each split is in the following KSP objects:\n");
44: PetscViewerASCIIPushTab(viewer);
45: for (i=0; i<jac->nsplits; i++) {
46: PetscViewerASCIIPrintf(viewer,"Split number %D Fields ",i);
47: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
48: for (j=0; j<ilink->nfields; j++) {
49: if (j > 0) {
50: PetscViewerASCIIPrintf(viewer,",");
51: }
52: PetscViewerASCIIPrintf(viewer," %D",ilink->fields[j]);
53: }
54: PetscViewerASCIIPrintf(viewer,"\n");
55: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
56: KSPView(ilink->ksp,viewer);
57: ilink = ilink->next;
58: }
59: PetscViewerASCIIPopTab(viewer);
60: } else {
61: SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCFieldSplit",((PetscObject)viewer)->type_name);
62: }
63: return(0);
64: }
68: static PetscErrorCode PCFieldSplitSetDefaults(PC pc)
69: {
70: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
71: PetscErrorCode ierr;
72: PC_FieldSplitLink ilink = jac->head;
73: PetscInt i;
74: PetscTruth flg = PETSC_FALSE,*fields;
77: PetscOptionsGetTruth(pc->prefix,"-pc_fieldsplit_default",&flg,PETSC_NULL);
78: if (!ilink || flg) {
79: PetscInfo(pc,"Using default splitting of fields\n");
80: if (jac->bs <= 0) {
81: MatGetBlockSize(pc->pmat,&jac->bs);
82: }
83: PetscMalloc(jac->bs*sizeof(PetscTruth),&fields);
84: PetscMemzero(fields,jac->bs*sizeof(PetscTruth));
85: while (ilink) {
86: for (i=0; i<ilink->nfields; i++) {
87: fields[ilink->fields[i]] = PETSC_TRUE;
88: }
89: ilink = ilink->next;
90: }
91: jac->defaultsplit = PETSC_TRUE;
92: for (i=0; i<jac->bs; i++) {
93: if (!fields[i]) {
94: PCFieldSplitSetFields(pc,1,&i);
95: } else {
96: jac->defaultsplit = PETSC_FALSE;
97: }
98: }
99: }
100: return(0);
101: }
106: static PetscErrorCode PCSetUp_FieldSplit(PC pc)
107: {
108: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
109: PetscErrorCode ierr;
110: PC_FieldSplitLink ilink;
111: PetscInt i,nsplit,ccsize;
112: MatStructure flag = pc->flag;
113: PetscTruth sorted;
116: PCFieldSplitSetDefaults(pc);
117: nsplit = jac->nsplits;
118: ilink = jac->head;
120: /* get the matrices for each split */
121: if (!jac->is) {
122: PetscInt rstart,rend,nslots,bs;
124: bs = jac->bs;
125: MatGetOwnershipRange(pc->pmat,&rstart,&rend);
126: MatGetLocalSize(pc->pmat,PETSC_NULL,&ccsize);
127: nslots = (rend - rstart)/bs;
128: PetscMalloc(nsplit*sizeof(IS),&jac->is);
129: PetscMalloc(nsplit*sizeof(IS),&jac->cis);
130: PetscMalloc(nsplit*sizeof(PetscInt),&jac->csize);
131: for (i=0; i<nsplit; i++) {
132: if (jac->defaultsplit) {
133: ISCreateStride(pc->comm,nslots,rstart+i,nsplit,&jac->is[i]);
134: jac->csize[i] = ccsize/nsplit;
135: } else {
136: if (ilink->nfields > 1) {
137: PetscInt *ii,j,k,nfields = ilink->nfields,*fields = ilink->fields;
138: PetscMalloc(ilink->nfields*nslots*sizeof(PetscInt),&ii);
139: for (j=0; j<nslots; j++) {
140: for (k=0; k<nfields; k++) {
141: ii[nfields*j + k] = rstart + bs*j + fields[k];
142: }
143: }
144: ISCreateGeneral(pc->comm,nslots*nfields,ii,&jac->is[i]);
145: PetscFree(ii);
146: } else {
147: ISCreateStride(pc->comm,nslots,ilink->fields[0],bs,&jac->is[i]);
148: }
149: jac->csize[i] = (ccsize/bs)*ilink->nfields;
150: ISSorted(jac->is[i],&sorted);
151: if (!sorted) SETERRQ(PETSC_ERR_USER,"Fields must be sorted when creating split");
152: ilink = ilink->next;
153: }
154: ISAllGather(jac->is[i],&jac->cis[i]);
155: }
156: }
157:
158: if (!jac->pmat) {
159: PetscMalloc(nsplit*sizeof(Mat),&jac->pmat);
160: for (i=0; i<nsplit; i++) {
161: MatGetSubMatrix(pc->pmat,jac->is[i],jac->cis[i],jac->csize[i],MAT_INITIAL_MATRIX,&jac->pmat[i]);
162: }
163: } else {
164: for (i=0; i<nsplit; i++) {
165: MatGetSubMatrix(pc->pmat,jac->is[i],jac->cis[i],jac->csize[i],MAT_REUSE_MATRIX,&jac->pmat[i]);
166: }
167: }
169: /* set up the individual PCs */
170: i = 0;
171: ilink = jac->head;
172: while (ilink) {
173: KSPSetOperators(ilink->ksp,jac->pmat[i],jac->pmat[i],flag);
174: KSPSetFromOptions(ilink->ksp);
175: KSPSetUp(ilink->ksp);
176: i++;
177: ilink = ilink->next;
178: }
179:
180: /* create work vectors for each split */
181: if (!jac->x) {
182: Vec xtmp;
183: PetscMalloc2(nsplit,Vec,&jac->x,nsplit,Vec,&jac->y);
184: ilink = jac->head;
185: for (i=0; i<nsplit; i++) {
186: Vec *vl,*vr;
188: KSPGetVecs(ilink->ksp,1,&vr,1,&vl);
189: ilink->x = *vr;
190: ilink->y = *vl;
191: PetscFree(vr);
192: PetscFree(vl);
193: jac->x[i] = ilink->x;
194: jac->y[i] = ilink->y;
195: ilink = ilink->next;
196: }
197: /* compute scatter contexts needed by multiplicative versions and non-default splits */
198:
199: ilink = jac->head;
200: MatGetVecs(pc->pmat,&xtmp,PETSC_NULL);
201: for (i=0; i<nsplit; i++) {
202: VecScatterCreate(xtmp,jac->is[i],jac->x[i],PETSC_NULL,&ilink->sctx);
203: ilink = ilink->next;
204: }
205: VecDestroy(xtmp);
206: }
207: return(0);
208: }
210: #define FieldSplitSplitSolveAdd(ilink,xx,yy) \
211: (VecScatterBegin(xx,ilink->x,INSERT_VALUES,SCATTER_FORWARD,ilink->sctx) || \
212: VecScatterEnd(xx,ilink->x,INSERT_VALUES,SCATTER_FORWARD,ilink->sctx) || \
213: KSPSolve(ilink->ksp,ilink->x,ilink->y) || \
214: VecScatterBegin(ilink->y,yy,ADD_VALUES,SCATTER_REVERSE,ilink->sctx) || \
215: VecScatterEnd(ilink->y,yy,ADD_VALUES,SCATTER_REVERSE,ilink->sctx))
219: static PetscErrorCode PCApply_FieldSplit(PC pc,Vec x,Vec y)
220: {
221: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
222: PetscErrorCode ierr;
223: PC_FieldSplitLink ilink = jac->head;
224: PetscInt bs;
227: CHKMEMQ;
228: VecGetBlockSize(x,&bs);
229: VecSetBlockSize(x,jac->bs);
230: VecSetBlockSize(y,jac->bs);
232: if (jac->type == PC_COMPOSITE_ADDITIVE) {
233: if (jac->defaultsplit) {
234: VecStrideGatherAll(x,jac->x,INSERT_VALUES);
235: while (ilink) {
236: KSPSolve(ilink->ksp,ilink->x,ilink->y);
237: ilink = ilink->next;
238: }
239: VecStrideScatterAll(jac->y,y,INSERT_VALUES);
240: } else {
241: VecSet(y,0.0);
242: while (ilink) {
243: FieldSplitSplitSolveAdd(ilink,x,y);
244: ilink = ilink->next;
245: }
246: }
247: } else {
248: if (!jac->w1) {
249: VecDuplicate(x,&jac->w1);
250: VecDuplicate(x,&jac->w2);
251: }
252: VecSet(y,0.0);
253: FieldSplitSplitSolveAdd(ilink,x,y);
254: while (ilink->next) {
255: ilink = ilink->next;
256: MatMult(pc->pmat,y,jac->w1);
257: VecWAXPY(jac->w2,-1.0,jac->w1,x);
258: FieldSplitSplitSolveAdd(ilink,jac->w2,y);
259: }
260: if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
261: while (ilink->previous) {
262: ilink = ilink->previous;
263: MatMult(pc->pmat,y,jac->w1);
264: VecWAXPY(jac->w2,-1.0,jac->w1,x);
265: FieldSplitSplitSolveAdd(ilink,jac->w2,y);
266: }
267: }
268: }
269: CHKMEMQ;
270: return(0);
271: }
273: #define FieldSplitSplitSolveAddTranspose(ilink,xx,yy) \
274: (VecScatterBegin(xx,ilink->y,INSERT_VALUES,SCATTER_FORWARD,ilink->sctx) || \
275: VecScatterEnd(xx,ilink->y,INSERT_VALUES,SCATTER_FORWARD,ilink->sctx) || \
276: KSPSolveTranspose(ilink->ksp,ilink->y,ilink->x) || \
277: VecScatterBegin(ilink->x,yy,ADD_VALUES,SCATTER_REVERSE,ilink->sctx) || \
278: VecScatterEnd(ilink->x,yy,ADD_VALUES,SCATTER_REVERSE,ilink->sctx))
282: static PetscErrorCode PCApplyTranspose_FieldSplit(PC pc,Vec x,Vec y)
283: {
284: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
285: PetscErrorCode ierr;
286: PC_FieldSplitLink ilink = jac->head;
287: PetscInt bs;
290: CHKMEMQ;
291: VecGetBlockSize(x,&bs);
292: VecSetBlockSize(x,jac->bs);
293: VecSetBlockSize(y,jac->bs);
295: if (jac->type == PC_COMPOSITE_ADDITIVE) {
296: if (jac->defaultsplit) {
297: VecStrideGatherAll(x,jac->x,INSERT_VALUES);
298: while (ilink) {
299: KSPSolveTranspose(ilink->ksp,ilink->x,ilink->y);
300: ilink = ilink->next;
301: }
302: VecStrideScatterAll(jac->y,y,INSERT_VALUES);
303: } else {
304: VecSet(y,0.0);
305: while (ilink) {
306: FieldSplitSplitSolveAddTranspose(ilink,x,y);
307: ilink = ilink->next;
308: }
309: }
310: } else {
311: if (!jac->w1) {
312: VecDuplicate(x,&jac->w1);
313: VecDuplicate(x,&jac->w2);
314: }
315: VecSet(y,0.0);
316: if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) {
317: FieldSplitSplitSolveAddTranspose(ilink,x,y);
318: while (ilink->next) {
319: ilink = ilink->next;
320: MatMultTranspose(pc->pmat,y,jac->w1);
321: VecWAXPY(jac->w2,-1.0,jac->w1,x);
322: FieldSplitSplitSolveAddTranspose(ilink,jac->w2,y);
323: }
324: while (ilink->previous) {
325: ilink = ilink->previous;
326: MatMultTranspose(pc->pmat,y,jac->w1);
327: VecWAXPY(jac->w2,-1.0,jac->w1,x);
328: FieldSplitSplitSolveAddTranspose(ilink,jac->w2,y);
329: }
330: } else {
331: while (ilink->next) { /* get to last entry in linked list */
332: ilink = ilink->next;
333: }
334: FieldSplitSplitSolveAddTranspose(ilink,x,y);
335: while (ilink->previous) {
336: ilink = ilink->previous;
337: MatMultTranspose(pc->pmat,y,jac->w1);
338: VecWAXPY(jac->w2,-1.0,jac->w1,x);
339: FieldSplitSplitSolveAddTranspose(ilink,jac->w2,y);
340: }
341: }
342: }
343: CHKMEMQ;
344: return(0);
345: }
349: static PetscErrorCode PCDestroy_FieldSplit(PC pc)
350: {
351: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
352: PetscErrorCode ierr;
353: PC_FieldSplitLink ilink = jac->head,next;
356: while (ilink) {
357: KSPDestroy(ilink->ksp);
358: if (ilink->x) {VecDestroy(ilink->x);}
359: if (ilink->y) {VecDestroy(ilink->y);}
360: if (ilink->sctx) {VecScatterDestroy(ilink->sctx);}
361: next = ilink->next;
362: PetscFree2(ilink,ilink->fields);
363: ilink = next;
364: }
365: PetscFree2(jac->x,jac->y);
366: if (jac->pmat) {MatDestroyMatrices(jac->nsplits,&jac->pmat);}
367: if (jac->is) {
368: PetscInt i;
369: for (i=0; i<jac->nsplits; i++) {ISDestroy(jac->is[i]);}
370: PetscFree(jac->is);
371: }
372: if (jac->cis) {
373: PetscInt i;
374: for (i=0; i<jac->nsplits; i++) {ISDestroy(jac->cis[i]);}
375: PetscFree(jac->cis);
376: }
377: if (jac->w1) {VecDestroy(jac->w1);}
378: if (jac->w2) {VecDestroy(jac->w2);}
379: PetscFree(jac->csize);
380: PetscFree(jac);
381: return(0);
382: }
386: static PetscErrorCode PCSetFromOptions_FieldSplit(PC pc)
387: {
389: PetscInt i = 0,nfields,*fields,bs;
390: PetscTruth flg;
391: char optionname[128];
392: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
395: PetscOptionsHead("FieldSplit options");
396: PetscOptionsInt("-pc_fieldsplit_block_size","Blocksize that defines number of fields","PCFieldSplitSetBlockSize",jac->bs,&bs,&flg);
397: if (flg) {
398: PCFieldSplitSetBlockSize(pc,bs);
399: }
400: if (jac->bs <= 0) {
401: PCFieldSplitSetBlockSize(pc,1);
402: }
403: PetscOptionsEnum("-pc_fieldsplit_type","Type of composition","PCFieldSplitSetType",PCCompositeTypes,(PetscEnum)jac->type,(PetscEnum*)&jac->type,&flg);
404: PetscMalloc(jac->bs*sizeof(PetscInt),&fields);
405: while (PETSC_TRUE) {
406: sprintf(optionname,"-pc_fieldsplit_%d_fields",(int)i++);
407: nfields = jac->bs;
408: PetscOptionsIntArray(optionname,"Fields in this split","PCFieldSplitSetFields",fields,&nfields,&flg);
409: if (!flg) break;
410: if (!nfields) SETERRQ(PETSC_ERR_USER,"Cannot list zero fields");
411: PCFieldSplitSetFields(pc,nfields,fields);
412: }
413: PetscFree(fields);
414: PetscOptionsTail();
415: return(0);
416: }
418: /*------------------------------------------------------------------------------------*/
423: PetscErrorCode PCFieldSplitSetFields_FieldSplit(PC pc,PetscInt n,PetscInt *fields)
424: {
425: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
426: PetscErrorCode ierr;
427: PC_FieldSplitLink ilink,next = jac->head;
428: char prefix[128];
429: PetscInt i;
432: if (n <= 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Negative number of fields requested");
433: for (i=0; i<n; i++) {
434: if (fields[i] >= jac->bs) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Field %D requested but only %D exist",fields[i],jac->bs);
435: if (fields[i] < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Negative field %D requested",fields[i]);
436: }
437: PetscMalloc2(1,struct _PC_FieldSplitLink,&ilink,n,PetscInt,&ilink->fields);
438: PetscMemcpy(ilink->fields,fields,n*sizeof(PetscInt));
439: ilink->nfields = n;
440: ilink->next = PETSC_NULL;
441: KSPCreate(pc->comm,&ilink->ksp);
442: KSPSetType(ilink->ksp,KSPPREONLY);
444: if (pc->prefix) {
445: sprintf(prefix,"%sfieldsplit_%d_",pc->prefix,(int)jac->nsplits);
446: } else {
447: sprintf(prefix,"fieldsplit_%d_",(int)jac->nsplits);
448: }
449: KSPSetOptionsPrefix(ilink->ksp,prefix);
451: if (!next) {
452: jac->head = ilink;
453: ilink->previous = PETSC_NULL;
454: } else {
455: while (next->next) {
456: next = next->next;
457: }
458: next->next = ilink;
459: ilink->previous = next;
460: }
461: jac->nsplits++;
462: return(0);
463: }
470: PetscErrorCode PCFieldSplitGetSubKSP_FieldSplit(PC pc,PetscInt *n,KSP **subksp)
471: {
472: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
473: PetscErrorCode ierr;
474: PetscInt cnt = 0;
475: PC_FieldSplitLink ilink = jac->head;
478: PetscMalloc(jac->nsplits*sizeof(KSP*),subksp);
479: while (ilink) {
480: (*subksp)[cnt++] = ilink->ksp;
481: ilink = ilink->next;
482: }
483: if (cnt != jac->nsplits) SETERRQ2(PETSC_ERR_PLIB,"Corrupt PCFIELDSPLIT object: number splits in linked list %D in object %D",cnt,jac->nsplits);
484: *n = jac->nsplits;
485: return(0);
486: }
491: /*@
492: PCFieldSplitSetFields - Sets the fields for one particular split in the field split preconditioner
494: Collective on PC
496: Input Parameters:
497: + pc - the preconditioner context
498: . n - the number of fields in this split
499: . fields - the fields in this split
501: Level: intermediate
503: .seealso: PCFieldSplitGetSubKSP(), PCFIELDSPLIT, PCFieldSplitSetBlockSize()
505: @*/
506: PetscErrorCode PCFieldSplitSetFields(PC pc,PetscInt n, PetscInt *fields)
507: {
508: PetscErrorCode ierr,(*f)(PC,PetscInt,PetscInt *);
512: PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetFields_C",(void (**)(void))&f);
513: if (f) {
514: (*f)(pc,n,fields);
515: }
516: return(0);
517: }
521: /*@
522: PCFieldSplitSetBlockSize - Sets the block size for defining where fields start in the
523: fieldsplit preconditioner. If not set the matrix block size is used.
525: Collective on PC
527: Input Parameters:
528: + pc - the preconditioner context
529: - bs - the block size
531: Level: intermediate
533: .seealso: PCFieldSplitGetSubKSP(), PCFIELDSPLIT, PCFieldSplitSetFields()
535: @*/
536: PetscErrorCode PCFieldSplitSetBlockSize(PC pc,PetscInt bs)
537: {
538: PetscErrorCode ierr,(*f)(PC,PetscInt);
542: PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetBlockSize_C",(void (**)(void))&f);
543: if (f) {
544: (*f)(pc,bs);
545: }
546: return(0);
547: }
551: /*@C
552: PCFieldSplitGetSubKSP - Gets the KSP contexts for all splits
553:
554: Collective on KSP
556: Input Parameter:
557: . pc - the preconditioner context
559: Output Parameters:
560: + n - the number of split
561: - pc - the array of KSP contexts
563: Note:
564: After PCFieldSplitGetSubKSP() the array of KSPs IS to be freed
566: You must call KSPSetUp() before calling PCFieldSplitGetSubKSP().
568: Level: advanced
570: .seealso: PCFIELDSPLIT
571: @*/
572: PetscErrorCode PCFieldSplitGetSubKSP(PC pc,PetscInt *n,KSP *subksp[])
573: {
574: PetscErrorCode ierr,(*f)(PC,PetscInt*,KSP **);
579: PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitGetSubKSP_C",(void (**)(void))&f);
580: if (f) {
581: (*f)(pc,n,subksp);
582: } else {
583: SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot get subksp for this type of PC");
584: }
585: return(0);
586: }
591: PetscErrorCode PCFieldSplitSetType_FieldSplit(PC pc,PCCompositeType type)
592: {
593: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
596: jac->type = type;
597: return(0);
598: }
604: PetscErrorCode PCFieldSplitSetBlockSize_FieldSplit(PC pc,PetscInt bs)
605: {
606: PC_FieldSplit *jac = (PC_FieldSplit*)pc->data;
609: if (bs < 1) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Blocksize must be positive, you gave %D",bs);
610: if (jac->bs > 0 && jac->bs != bs) SETERRQ2(PETSC_ERR_ARG_WRONGSTATE,"Cannot change fieldsplit blocksize from %D to %D after it has been set",jac->bs,bs);
611: jac->bs = bs;
612: return(0);
613: }
618: /*@C
619: PCFieldSplitSetType - Sets the type of fieldsplit preconditioner.
620:
621: Collective on PC
623: Input Parameter:
624: . pc - the preconditioner context
625: . type - PC_COMPOSITE_ADDITIVE (default), PC_COMPOSITE_MULTIPLICATIVE
627: Options Database Key:
628: . -pc_fieldsplit_type <type: one of multiplicative, additive, special> - Sets fieldsplit preconditioner type
630: Level: Developer
632: .keywords: PC, set, type, composite preconditioner, additive, multiplicative
634: .seealso: PCCompositeSetType()
636: @*/
637: PetscErrorCode PCFieldSplitSetType(PC pc,PCCompositeType type)
638: {
639: PetscErrorCode ierr,(*f)(PC,PCCompositeType);
643: PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetType_C",(void (**)(void))&f);
644: if (f) {
645: (*f)(pc,type);
646: }
647: return(0);
648: }
650: /* -------------------------------------------------------------------------------------*/
651: /*MC
652: PCFIELDSPLIT - Preconditioner created by combining separate preconditioners for individual
653: fields or groups of fields
656: To set options on the solvers for each block append -sub_ to all the PC
657: options database keys. For example, -sub_pc_type ilu -sub_pc_factor_levels 1
658:
659: To set the options on the solvers separate for each block call PCFieldSplitGetSubKSP()
660: and set the options directly on the resulting KSP object
662: Level: intermediate
664: Options Database Keys:
665: + -pc_splitfield_%d_fields <a,b,..> - indicates the fields to be used in the %d'th split
666: . -pc_splitfield_default - automatically add any fields to additional splits that have not
667: been supplied explicitly by -pc_splitfield_%d_fields
668: . -pc_splitfield_block_size <bs> - size of block that defines fields (i.e. there are bs fields)
669: - -pc_splitfield_type <additive,multiplicative>
671: Concepts: physics based preconditioners
673: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
674: PCFieldSplitGetSubKSP(), PCFieldSplitSetFields(),PCFieldSplitSetType()
675: M*/
680: PetscErrorCode PCCreate_FieldSplit(PC pc)
681: {
683: PC_FieldSplit *jac;
686: PetscNew(PC_FieldSplit,&jac);
687: PetscLogObjectMemory(pc,sizeof(PC_FieldSplit));
688: jac->bs = -1;
689: jac->nsplits = 0;
690: jac->type = PC_COMPOSITE_ADDITIVE;
692: pc->data = (void*)jac;
694: pc->ops->apply = PCApply_FieldSplit;
695: pc->ops->applytranspose = PCApplyTranspose_FieldSplit;
696: pc->ops->setup = PCSetUp_FieldSplit;
697: pc->ops->destroy = PCDestroy_FieldSplit;
698: pc->ops->setfromoptions = PCSetFromOptions_FieldSplit;
699: pc->ops->view = PCView_FieldSplit;
700: pc->ops->applyrichardson = 0;
702: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitGetSubKSP_C","PCFieldSplitGetSubKSP_FieldSplit",
703: PCFieldSplitGetSubKSP_FieldSplit);
704: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetFields_C","PCFieldSplitSetFields_FieldSplit",
705: PCFieldSplitSetFields_FieldSplit);
706: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetType_C","PCFieldSplitSetType_FieldSplit",
707: PCFieldSplitSetType_FieldSplit);
708: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetBlockSize_C","PCFieldSplitSetBlockSize_FieldSplit",
709: PCFieldSplitSetBlockSize_FieldSplit);
710: return(0);
711: }