1 /* 2 3 */ 4 #include "src/ksp/pc/pcimpl.h" /*I "petscpc.h" I*/ 5 6 typedef struct _PC_FieldSplitLink *PC_FieldSplitLink; 7 struct _PC_FieldSplitLink { 8 PC pc; 9 Vec x,y; 10 PetscInt nfields; 11 PetscInt *fields; 12 PC_FieldSplitLink next; 13 }; 14 15 typedef struct { 16 PetscInt bs; 17 PetscInt nsplits; 18 PC_FieldSplitLink head; 19 Vec *x,*y; 20 } PC_FieldSplit; 21 22 #undef __FUNCT__ 23 #define __FUNCT__ "PCView_FieldSplit" 24 static PetscErrorCode PCView_FieldSplit(PC pc,PetscViewer viewer) 25 { 26 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 27 PetscErrorCode ierr; 28 PetscTruth iascii; 29 PetscInt i,j; 30 PC_FieldSplitLink link = jac->head; 31 32 PetscFunctionBegin; 33 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 34 if (iascii) { 35 ierr = PetscViewerASCIIPrintf(viewer," FieldSplit: total splits = %D",jac->nsplits);CHKERRQ(ierr); 36 ierr = PetscViewerASCIIPrintf(viewer," Solver info for each split is in the following PC objects:\n");CHKERRQ(ierr); 37 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 38 for (i=0; i<jac->nsplits; i++) { 39 ierr = PetscViewerASCIIPrintf(viewer,"Split number %D Fields ",i);CHKERRQ(ierr); 40 for (j=0; j<link->nfields; j++) { 41 ierr = PetscViewerASCIIPrintf(viewer,"%D \n",link->fields[j]);CHKERRQ(ierr); 42 } 43 ierr = PetscViewerASCIIPrintf(viewer,"\n");CHKERRQ(ierr); 44 ierr = PCView(link->pc,viewer);CHKERRQ(ierr); 45 link = link->next; 46 } 47 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 48 } else { 49 SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCFieldSplit",((PetscObject)viewer)->type_name); 50 } 51 PetscFunctionReturn(0); 52 } 53 54 #undef __FUNCT__ 55 #define __FUNCT__ "PCSetUp_FieldSplit" 56 static PetscErrorCode PCSetUp_FieldSplit(PC pc) 57 { 58 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 59 PetscErrorCode ierr; 60 PC_FieldSplitLink link = jac->head; 61 62 PetscFunctionBegin; 63 while (link) { 64 ierr = PCSetUp(link->pc);CHKERRQ(ierr); 65 link = link->next; 66 } 67 PetscFunctionReturn(0); 68 } 69 70 #undef __FUNCT__ 71 #define __FUNCT__ "PCApply_FieldSplit" 72 static PetscErrorCode PCApply_FieldSplit(PC pc,Vec x,Vec y) 73 { 74 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 75 PetscErrorCode ierr; 76 PC_FieldSplitLink link = jac->head; 77 78 PetscFunctionBegin; 79 ierr = VecStrideGatherAll(x,jac->x,INSERT_VALUES);CHKERRQ(ierr); 80 while (link) { 81 ierr = PCApply(link->pc,link->x,link->y);CHKERRQ(ierr); 82 link = link->next; 83 } 84 ierr = VecStrideScatterAll(jac->y,y,INSERT_VALUES);CHKERRQ(ierr); 85 PetscFunctionReturn(0); 86 } 87 88 #undef __FUNCT__ 89 #define __FUNCT__ "PCDestroy_FieldSplit" 90 static PetscErrorCode PCDestroy_FieldSplit(PC pc) 91 { 92 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 93 PetscErrorCode ierr; 94 PC_FieldSplitLink link = jac->head; 95 96 PetscFunctionBegin; 97 while (link) { 98 ierr = PCDestroy(link->pc);CHKERRQ(ierr); 99 ierr = PetscFree2(link,link->fields);CHKERRQ(ierr); 100 } 101 if (jac->x) {ierr = PetscFree2(jac->x,jac->y);CHKERRQ(ierr);} 102 ierr = PetscFree(jac);CHKERRQ(ierr); 103 PetscFunctionReturn(0); 104 } 105 106 #undef __FUNCT__ 107 #define __FUNCT__ "PCSetFromOptions_FieldSplit" 108 static PetscErrorCode PCSetFromOptions_FieldSplit(PC pc) 109 { 110 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 111 PetscErrorCode ierr; 112 PC_FieldSplitLink link = jac->head; 113 114 PetscFunctionBegin; 115 if (!link) { /* user never set fields, so set them ourselves */ 116 link = jac->head; 117 } 118 while (link) { 119 ierr = PCApply(link->pc,link->x,link->y);CHKERRQ(ierr); 120 link = link->next; 121 } 122 PetscFunctionReturn(0); 123 } 124 125 /*------------------------------------------------------------------------------------*/ 126 127 EXTERN_C_BEGIN 128 #undef __FUNCT__ 129 #define __FUNCT__ "PCFieldSplitSetFields_FieldSplit" 130 PetscErrorCode PCFieldSplitSetFields_FieldSplit(PC pc,PetscInt n,PetscInt *fields) 131 { 132 PC_FieldSplit *jac; 133 PetscErrorCode ierr; 134 PetscInt *myfields; 135 PC_FieldSplitLink link,next = jac->head; 136 137 PetscFunctionBegin; 138 if (n <= 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Negative number of fields requested"); 139 ierr = PetscMalloc2(1,struct _PC_FieldSplitLink,&link,n,PetscInt,&myfields);CHKERRQ(ierr); 140 ierr = PetscMemcpy(myfields,fields,n*sizeof(PetscInt));CHKERRQ(ierr); 141 link->fields = myfields; 142 link->nfields = n; 143 144 if (!next) { 145 jac->head = link; 146 } else { 147 while (next->next) { 148 next = next->next; 149 } 150 next->next = link; 151 } 152 jac->nsplits++; 153 PetscFunctionReturn(0); 154 } 155 EXTERN_C_END 156 157 158 EXTERN_C_BEGIN 159 #undef __FUNCT__ 160 #define __FUNCT__ "PCFieldSplitGetSubPC_FieldSplit" 161 PetscErrorCode PCFieldSplitGetSubPC_FieldSplit(PC pc,PetscInt *n,PC **subpc) 162 { 163 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 164 PetscErrorCode ierr; 165 PetscInt cnt = 0; 166 PC_FieldSplitLink link; 167 168 PetscFunctionBegin; 169 ierr = PetscMalloc(jac->nsplits*sizeof(PC*),subpc);CHKERRQ(ierr); 170 while (link) { 171 (*subpc)[cnt++] = link->pc; 172 link = link->next; 173 } 174 if (cnt != jac->nsplits) SETERRQ2(PETSC_ERR_PLIB,"Corrupt PCFIELDSPLIT object: number splits in linked list %D in object %D",cnt,jac->nsplits); 175 *n = jac->nsplits; 176 PetscFunctionReturn(0); 177 } 178 EXTERN_C_END 179 180 #undef __FUNCT__ 181 #define __FUNCT__ "PCFieldSplitSetFields" 182 /*@ 183 PCFieldSplitSetFields - Sets the fields for one particular split in the field split preconditioner 184 185 Collective on PC 186 187 Input Parameters: 188 + pc - the preconditioner context 189 . n - the number of fields in this split 190 . fields - the fields in this split 191 192 Level: intermediate 193 194 .seealso: PCFieldSplitGetSubKSP(), PCFIELDSPLIT 195 196 @*/ 197 PetscErrorCode PCFieldSplitSetFields(PC pc,PetscInt n, PetscInt *fields) 198 { 199 PetscErrorCode ierr,(*f)(PC,PetscInt,PetscInt *); 200 201 PetscFunctionBegin; 202 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 203 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetFields_C",(void (**)(void))&f);CHKERRQ(ierr); 204 if (f) { 205 ierr = (*f)(pc,n,fields);CHKERRQ(ierr); 206 } 207 PetscFunctionReturn(0); 208 } 209 210 #undef __FUNCT__ 211 #define __FUNCT__ "PCFieldSplitGetSubPC" 212 /*@C 213 PCFieldSplitGetSubPC - Gets the PC contexts for all splits 214 215 Collective on PC 216 217 Input Parameter: 218 . pc - the preconditioner context 219 220 Output Parameters: 221 + n - the number of split 222 - pc - the array of PC contexts 223 224 Note: 225 After PCFieldSplitGetSubPC() the array of PCs IS to be freed 226 227 You must call KSPSetUp() before calling PCFieldSplitGetSubPC(). 228 229 Level: advanced 230 231 .keywords: PC, 232 233 .seealso: PCFIELDSPLIT 234 @*/ 235 PetscErrorCode PCFieldSplitGetSubPC(PC pc,PetscInt *n,PC *subpc[]) 236 { 237 PetscErrorCode ierr,(*f)(PC,PetscInt*,PC **); 238 239 PetscFunctionBegin; 240 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 241 PetscValidIntPointer(n,2); 242 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitGetSubPC_C",(void (**)(void))&f);CHKERRQ(ierr); 243 if (f) { 244 ierr = (*f)(pc,n,subpc);CHKERRQ(ierr); 245 } else { 246 SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot get subpc for this type of PC"); 247 } 248 249 PetscFunctionReturn(0); 250 } 251 252 /* -------------------------------------------------------------------------------------*/ 253 /*MC 254 PCFieldSplit - Preconditioner created by combining seperate preconditioners for individual 255 fields or groups of fields 256 257 258 To set options on the solvers for each block append -sub_ to all the PC 259 options database keys. For example, -sub_pc_type ilu -sub_pc_ilu_levels 1 260 261 To set the options on the solvers seperate for each block call PCFieldSplitGetSubPC() 262 and set the options directly on the resulting PC object 263 264 Level: intermediate 265 266 Concepts: physics based preconditioners 267 268 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 269 PCFieldSplitGetSubPC(), PCFieldSplitSetFields() 270 M*/ 271 272 EXTERN_C_BEGIN 273 #undef __FUNCT__ 274 #define __FUNCT__ "PCCreate_FieldSplit" 275 PetscErrorCode PCCreate_FieldSplit(PC pc) 276 { 277 PetscErrorCode ierr; 278 PC_FieldSplit *jac; 279 280 PetscFunctionBegin; 281 ierr = PetscNew(PC_FieldSplit,&jac);CHKERRQ(ierr); 282 PetscLogObjectMemory(pc,sizeof(PC_FieldSplit)); 283 ierr = PetscMemzero(jac,sizeof(PC_FieldSplit));CHKERRQ(ierr); 284 jac->bs = -1; 285 jac->nsplits = 0; 286 pc->data = (void*)jac; 287 288 pc->ops->apply = PCApply_FieldSplit; 289 pc->ops->setup = PCSetUp_FieldSplit; 290 pc->ops->destroy = PCDestroy_FieldSplit; 291 pc->ops->setfromoptions = PCSetFromOptions_FieldSplit; 292 pc->ops->view = PCView_FieldSplit; 293 pc->ops->applyrichardson = 0; 294 295 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitGetSubPC_C","PCFieldSplitGetSubPC_FieldSplit", 296 PCFieldSplitGetSubPC_FieldSplit);CHKERRQ(ierr); 297 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetFields_C","PCFieldSplitSetFields_FieldSplit", 298 PCFieldSplitSetFields_FieldSplit);CHKERRQ(ierr); 299 PetscFunctionReturn(0); 300 } 301 EXTERN_C_END 302 303 304