1 #define PETSCKSP_DLL 2 3 /* 4 5 */ 6 #include "private/pcimpl.h" /*I "petscpc.h" I*/ 7 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 }; 17 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; 28 29 #undef __FUNCT__ 30 #define __FUNCT__ "PCView_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; 38 39 PetscFunctionBegin; 40 ierr = PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&iascii);CHKERRQ(ierr); 41 if (iascii) { 42 ierr = PetscViewerASCIIPrintf(viewer," FieldSplit with %s composition: total splits = %D, blocksize = %D\n",PCCompositeTypes[jac->type],jac->nsplits,jac->bs);CHKERRQ(ierr); 43 ierr = PetscViewerASCIIPrintf(viewer," Solver info for each split is in the following KSP objects:\n");CHKERRQ(ierr); 44 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 45 for (i=0; i<jac->nsplits; i++) { 46 ierr = PetscViewerASCIIPrintf(viewer,"Split number %D Fields ",i);CHKERRQ(ierr); 47 ierr = PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);CHKERRQ(ierr); 48 for (j=0; j<ilink->nfields; j++) { 49 if (j > 0) { 50 ierr = PetscViewerASCIIPrintf(viewer,",");CHKERRQ(ierr); 51 } 52 ierr = PetscViewerASCIIPrintf(viewer," %D",ilink->fields[j]);CHKERRQ(ierr); 53 } 54 ierr = PetscViewerASCIIPrintf(viewer,"\n");CHKERRQ(ierr); 55 ierr = PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);CHKERRQ(ierr); 56 ierr = KSPView(ilink->ksp,viewer);CHKERRQ(ierr); 57 ilink = ilink->next; 58 } 59 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 60 } else { 61 SETERRQ1(PETSC_ERR_SUP,"Viewer type %s not supported for PCFieldSplit",((PetscObject)viewer)->type_name); 62 } 63 PetscFunctionReturn(0); 64 } 65 66 #undef __FUNCT__ 67 #define __FUNCT__ "PCFieldSplitSetDefaults" 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; 75 76 PetscFunctionBegin; 77 ierr = PetscOptionsGetTruth(pc->prefix,"-pc_fieldsplit_default",&flg,PETSC_NULL);CHKERRQ(ierr); 78 if (!ilink || flg) { 79 ierr = PetscInfo(pc,"Using default splitting of fields\n");CHKERRQ(ierr); 80 if (jac->bs <= 0) { 81 ierr = MatGetBlockSize(pc->pmat,&jac->bs);CHKERRQ(ierr); 82 } 83 ierr = PetscMalloc(jac->bs*sizeof(PetscTruth),&fields);CHKERRQ(ierr); 84 ierr = PetscMemzero(fields,jac->bs*sizeof(PetscTruth));CHKERRQ(ierr); 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 ierr = PCFieldSplitSetFields(pc,1,&i);CHKERRQ(ierr); 95 } else { 96 jac->defaultsplit = PETSC_FALSE; 97 } 98 } 99 } 100 PetscFunctionReturn(0); 101 } 102 103 104 #undef __FUNCT__ 105 #define __FUNCT__ "PCSetUp_FieldSplit" 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 114 PetscFunctionBegin; 115 ierr = PCFieldSplitSetDefaults(pc);CHKERRQ(ierr); 116 nsplit = jac->nsplits; 117 ilink = jac->head; 118 119 /* get the matrices for each split */ 120 if (!jac->is) { 121 PetscInt rstart,rend,nslots,bs; 122 123 bs = jac->bs; 124 ierr = MatGetOwnershipRange(pc->pmat,&rstart,&rend);CHKERRQ(ierr); 125 ierr = MatGetLocalSize(pc->pmat,PETSC_NULL,&ccsize);CHKERRQ(ierr); 126 nslots = (rend - rstart)/bs; 127 ierr = PetscMalloc(nsplit*sizeof(IS),&jac->is);CHKERRQ(ierr); 128 ierr = PetscMalloc(nsplit*sizeof(IS),&jac->cis);CHKERRQ(ierr); 129 ierr = PetscMalloc(nsplit*sizeof(PetscInt),&jac->csize);CHKERRQ(ierr); 130 for (i=0; i<nsplit; i++) { 131 if (jac->defaultsplit) { 132 ierr = ISCreateStride(pc->comm,nslots,rstart+i,nsplit,&jac->is[i]);CHKERRQ(ierr); 133 jac->csize[i] = ccsize/nsplit; 134 } else { 135 PetscInt *ii,j,k,nfields = ilink->nfields,*fields = ilink->fields; 136 PetscTruth sorted; 137 ierr = PetscMalloc(ilink->nfields*nslots*sizeof(PetscInt),&ii);CHKERRQ(ierr); 138 for (j=0; j<nslots; j++) { 139 for (k=0; k<nfields; k++) { 140 ii[nfields*j + k] = rstart + bs*j + fields[k]; 141 } 142 } 143 ierr = ISCreateGeneral(pc->comm,nslots*nfields,ii,&jac->is[i]);CHKERRQ(ierr); 144 jac->csize[i] = (ccsize/bs)*ilink->nfields; 145 ierr = ISSorted(jac->is[i],&sorted);CHKERRQ(ierr); 146 if (!sorted) SETERRQ(PETSC_ERR_USER,"Fields must be sorted when creating split"); 147 ierr = PetscFree(ii);CHKERRQ(ierr); 148 ilink = ilink->next; 149 } 150 ierr = ISAllGather(jac->is[i],&jac->cis[i]);CHKERRQ(ierr); 151 } 152 } 153 154 if (!jac->pmat) { 155 ierr = PetscMalloc(nsplit*sizeof(Mat),&jac->pmat);CHKERRQ(ierr); 156 for (i=0; i<nsplit; i++) { 157 ierr = MatGetSubMatrix(pc->pmat,jac->is[i],jac->cis[i],jac->csize[i],MAT_INITIAL_MATRIX,&jac->pmat[i]);CHKERRQ(ierr); 158 } 159 } else { 160 for (i=0; i<nsplit; i++) { 161 ierr = MatGetSubMatrix(pc->pmat,jac->is[i],jac->cis[i],jac->csize[i],MAT_REUSE_MATRIX,&jac->pmat[i]);CHKERRQ(ierr); 162 } 163 } 164 165 /* set up the individual PCs */ 166 i = 0; 167 ilink = jac->head; 168 while (ilink) { 169 ierr = KSPSetOperators(ilink->ksp,jac->pmat[i],jac->pmat[i],flag);CHKERRQ(ierr); 170 ierr = KSPSetFromOptions(ilink->ksp);CHKERRQ(ierr); 171 ierr = KSPSetUp(ilink->ksp);CHKERRQ(ierr); 172 i++; 173 ilink = ilink->next; 174 } 175 176 /* create work vectors for each split */ 177 if (!jac->x) { 178 Vec xtmp; 179 ierr = PetscMalloc2(nsplit,Vec,&jac->x,nsplit,Vec,&jac->y);CHKERRQ(ierr); 180 ilink = jac->head; 181 for (i=0; i<nsplit; i++) { 182 Vec *vl,*vr; 183 184 ierr = KSPGetVecs(ilink->ksp,1,&vr,1,&vl);CHKERRQ(ierr); 185 ilink->x = *vr; 186 ilink->y = *vl; 187 ierr = PetscFree(vr);CHKERRQ(ierr); 188 ierr = PetscFree(vl);CHKERRQ(ierr); 189 jac->x[i] = ilink->x; 190 jac->y[i] = ilink->y; 191 ilink = ilink->next; 192 } 193 /* compute scatter contexts needed by multiplicative versions and non-default splits */ 194 195 ilink = jac->head; 196 ierr = MatGetVecs(pc->pmat,&xtmp,PETSC_NULL);CHKERRQ(ierr); 197 for (i=0; i<nsplit; i++) { 198 ierr = VecScatterCreate(xtmp,jac->is[i],jac->x[i],PETSC_NULL,&ilink->sctx);CHKERRQ(ierr); 199 ilink = ilink->next; 200 } 201 ierr = VecDestroy(xtmp);CHKERRQ(ierr); 202 } 203 PetscFunctionReturn(0); 204 } 205 206 #define FieldSplitSplitSolveAdd(ilink,xx,yy) \ 207 (VecScatterBegin(xx,ilink->x,INSERT_VALUES,SCATTER_FORWARD,ilink->sctx) || \ 208 VecScatterEnd(xx,ilink->x,INSERT_VALUES,SCATTER_FORWARD,ilink->sctx) || \ 209 KSPSolve(ilink->ksp,ilink->x,ilink->y) || \ 210 VecScatterBegin(ilink->y,yy,ADD_VALUES,SCATTER_REVERSE,ilink->sctx) || \ 211 VecScatterEnd(ilink->y,yy,ADD_VALUES,SCATTER_REVERSE,ilink->sctx)) 212 213 #undef __FUNCT__ 214 #define __FUNCT__ "PCApply_FieldSplit" 215 static PetscErrorCode PCApply_FieldSplit(PC pc,Vec x,Vec y) 216 { 217 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 218 PetscErrorCode ierr; 219 PC_FieldSplitLink ilink = jac->head; 220 PetscInt bs; 221 222 PetscFunctionBegin; 223 CHKMEMQ; 224 ierr = VecGetBlockSize(x,&bs);CHKERRQ(ierr); 225 ierr = VecSetBlockSize(x,jac->bs);CHKERRQ(ierr); 226 ierr = VecSetBlockSize(y,jac->bs);CHKERRQ(ierr); 227 228 if (jac->type == PC_COMPOSITE_ADDITIVE) { 229 if (jac->defaultsplit) { 230 ierr = VecStrideGatherAll(x,jac->x,INSERT_VALUES);CHKERRQ(ierr); 231 while (ilink) { 232 ierr = KSPSolve(ilink->ksp,ilink->x,ilink->y);CHKERRQ(ierr); 233 ilink = ilink->next; 234 } 235 ierr = VecStrideScatterAll(jac->y,y,INSERT_VALUES);CHKERRQ(ierr); 236 } else { 237 PetscInt i = 0; 238 239 ierr = VecSet(y,0.0);CHKERRQ(ierr); 240 while (ilink) { 241 ierr = FieldSplitSplitSolveAdd(ilink,x,y);CHKERRQ(ierr); 242 ilink = ilink->next; 243 i++; 244 } 245 } 246 } else { 247 if (!jac->w1) { 248 ierr = VecDuplicate(x,&jac->w1);CHKERRQ(ierr); 249 ierr = VecDuplicate(x,&jac->w2);CHKERRQ(ierr); 250 } 251 ierr = VecSet(y,0.0);CHKERRQ(ierr); 252 ierr = FieldSplitSplitSolveAdd(ilink,x,y);CHKERRQ(ierr); 253 while (ilink->next) { 254 ilink = ilink->next; 255 ierr = MatMult(pc->pmat,y,jac->w1);CHKERRQ(ierr); 256 ierr = VecWAXPY(jac->w2,-1.0,jac->w1,x);CHKERRQ(ierr); 257 ierr = FieldSplitSplitSolveAdd(ilink,jac->w2,y);CHKERRQ(ierr); 258 } 259 if (jac->type == PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE) { 260 while (ilink->previous) { 261 ilink = ilink->previous; 262 ierr = MatMult(pc->pmat,y,jac->w1);CHKERRQ(ierr); 263 ierr = VecWAXPY(jac->w2,-1.0,jac->w1,x);CHKERRQ(ierr); 264 ierr = FieldSplitSplitSolveAdd(ilink,jac->w2,y);CHKERRQ(ierr); 265 } 266 } 267 } 268 CHKMEMQ; 269 PetscFunctionReturn(0); 270 } 271 272 #undef __FUNCT__ 273 #define __FUNCT__ "PCDestroy_FieldSplit" 274 static PetscErrorCode PCDestroy_FieldSplit(PC pc) 275 { 276 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 277 PetscErrorCode ierr; 278 PC_FieldSplitLink ilink = jac->head,next; 279 280 PetscFunctionBegin; 281 while (ilink) { 282 ierr = KSPDestroy(ilink->ksp);CHKERRQ(ierr); 283 if (ilink->x) {ierr = VecDestroy(ilink->x);CHKERRQ(ierr);} 284 if (ilink->y) {ierr = VecDestroy(ilink->y);CHKERRQ(ierr);} 285 if (ilink->sctx) {ierr = VecScatterDestroy(ilink->sctx);CHKERRQ(ierr);} 286 next = ilink->next; 287 ierr = PetscFree2(ilink,ilink->fields);CHKERRQ(ierr); 288 ilink = next; 289 } 290 ierr = PetscFree2(jac->x,jac->y);CHKERRQ(ierr); 291 if (jac->pmat) {ierr = MatDestroyMatrices(jac->nsplits,&jac->pmat);CHKERRQ(ierr);} 292 if (jac->is) { 293 PetscInt i; 294 for (i=0; i<jac->nsplits; i++) {ierr = ISDestroy(jac->is[i]);CHKERRQ(ierr);} 295 ierr = PetscFree(jac->is);CHKERRQ(ierr); 296 } 297 if (jac->cis) { 298 PetscInt i; 299 for (i=0; i<jac->nsplits; i++) {ierr = ISDestroy(jac->cis[i]);CHKERRQ(ierr);} 300 ierr = PetscFree(jac->cis);CHKERRQ(ierr); 301 } 302 if (jac->w1) {ierr = VecDestroy(jac->w1);CHKERRQ(ierr);} 303 if (jac->w2) {ierr = VecDestroy(jac->w2);CHKERRQ(ierr);} 304 ierr = PetscFree(jac->csize);CHKERRQ(ierr); 305 ierr = PetscFree(jac);CHKERRQ(ierr); 306 PetscFunctionReturn(0); 307 } 308 309 #undef __FUNCT__ 310 #define __FUNCT__ "PCSetFromOptions_FieldSplit" 311 static PetscErrorCode PCSetFromOptions_FieldSplit(PC pc) 312 { 313 PetscErrorCode ierr; 314 PetscInt i = 0,nfields,*fields,bs; 315 PetscTruth flg; 316 char optionname[128]; 317 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 318 319 PetscFunctionBegin; 320 ierr = PetscOptionsHead("FieldSplit options");CHKERRQ(ierr); 321 ierr = PetscOptionsInt("-pc_fieldsplit_block_size","Blocksize that defines number of fields","PCFieldSplitSetBlockSize",jac->bs,&bs,&flg);CHKERRQ(ierr); 322 if (flg) { 323 ierr = PCFieldSplitSetBlockSize(pc,bs);CHKERRQ(ierr); 324 } 325 if (jac->bs <= 0) { 326 ierr = PCFieldSplitSetBlockSize(pc,1);CHKERRQ(ierr); 327 } 328 ierr = PetscOptionsEnum("-pc_fieldsplit_type","Type of composition","PCFieldSplitSetType",PCCompositeTypes,(PetscEnum)jac->type,(PetscEnum*)&jac->type,&flg);CHKERRQ(ierr); 329 ierr = PetscMalloc(jac->bs*sizeof(PetscInt),&fields);CHKERRQ(ierr); 330 while (PETSC_TRUE) { 331 sprintf(optionname,"-pc_fieldsplit_%d_fields",(int)i++); 332 nfields = jac->bs; 333 ierr = PetscOptionsIntArray(optionname,"Fields in this split","PCFieldSplitSetFields",fields,&nfields,&flg);CHKERRQ(ierr); 334 if (!flg) break; 335 if (!nfields) SETERRQ(PETSC_ERR_USER,"Cannot list zero fields"); 336 ierr = PCFieldSplitSetFields(pc,nfields,fields);CHKERRQ(ierr); 337 } 338 ierr = PetscFree(fields);CHKERRQ(ierr); 339 ierr = PetscOptionsTail();CHKERRQ(ierr); 340 PetscFunctionReturn(0); 341 } 342 343 /*------------------------------------------------------------------------------------*/ 344 345 EXTERN_C_BEGIN 346 #undef __FUNCT__ 347 #define __FUNCT__ "PCFieldSplitSetFields_FieldSplit" 348 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitSetFields_FieldSplit(PC pc,PetscInt n,PetscInt *fields) 349 { 350 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 351 PetscErrorCode ierr; 352 PC_FieldSplitLink ilink,next = jac->head; 353 char prefix[128]; 354 PetscInt i; 355 356 PetscFunctionBegin; 357 if (n <= 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Negative number of fields requested"); 358 for (i=0; i<n; i++) { 359 if (fields[i] >= jac->bs) SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE,"Field %D requested but only %D exist",fields[i],jac->bs); 360 if (fields[i] < 0) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Negative field %D requested",fields[i]); 361 } 362 ierr = PetscMalloc2(1,struct _PC_FieldSplitLink,&ilink,n,PetscInt,&ilink->fields);CHKERRQ(ierr); 363 ierr = PetscMemcpy(ilink->fields,fields,n*sizeof(PetscInt));CHKERRQ(ierr); 364 ilink->nfields = n; 365 ilink->next = PETSC_NULL; 366 ierr = KSPCreate(pc->comm,&ilink->ksp);CHKERRQ(ierr); 367 ierr = KSPSetType(ilink->ksp,KSPPREONLY);CHKERRQ(ierr); 368 369 if (pc->prefix) { 370 sprintf(prefix,"%sfieldsplit_%d_",pc->prefix,(int)jac->nsplits); 371 } else { 372 sprintf(prefix,"fieldsplit_%d_",(int)jac->nsplits); 373 } 374 ierr = KSPSetOptionsPrefix(ilink->ksp,prefix);CHKERRQ(ierr); 375 376 if (!next) { 377 jac->head = ilink; 378 ilink->previous = PETSC_NULL; 379 } else { 380 while (next->next) { 381 next = next->next; 382 } 383 next->next = ilink; 384 ilink->previous = next; 385 } 386 jac->nsplits++; 387 PetscFunctionReturn(0); 388 } 389 EXTERN_C_END 390 391 392 EXTERN_C_BEGIN 393 #undef __FUNCT__ 394 #define __FUNCT__ "PCFieldSplitGetSubKSP_FieldSplit" 395 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitGetSubKSP_FieldSplit(PC pc,PetscInt *n,KSP **subksp) 396 { 397 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 398 PetscErrorCode ierr; 399 PetscInt cnt = 0; 400 PC_FieldSplitLink ilink = jac->head; 401 402 PetscFunctionBegin; 403 ierr = PetscMalloc(jac->nsplits*sizeof(KSP*),subksp);CHKERRQ(ierr); 404 while (ilink) { 405 (*subksp)[cnt++] = ilink->ksp; 406 ilink = ilink->next; 407 } 408 if (cnt != jac->nsplits) SETERRQ2(PETSC_ERR_PLIB,"Corrupt PCFIELDSPLIT object: number splits in linked list %D in object %D",cnt,jac->nsplits); 409 *n = jac->nsplits; 410 PetscFunctionReturn(0); 411 } 412 EXTERN_C_END 413 414 #undef __FUNCT__ 415 #define __FUNCT__ "PCFieldSplitSetFields" 416 /*@ 417 PCFieldSplitSetFields - Sets the fields for one particular split in the field split preconditioner 418 419 Collective on PC 420 421 Input Parameters: 422 + pc - the preconditioner context 423 . n - the number of fields in this split 424 . fields - the fields in this split 425 426 Level: intermediate 427 428 .seealso: PCFieldSplitGetSubKSP(), PCFIELDSPLIT, PCFieldSplitSetBlockSize() 429 430 @*/ 431 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitSetFields(PC pc,PetscInt n, PetscInt *fields) 432 { 433 PetscErrorCode ierr,(*f)(PC,PetscInt,PetscInt *); 434 435 PetscFunctionBegin; 436 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 437 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetFields_C",(void (**)(void))&f);CHKERRQ(ierr); 438 if (f) { 439 ierr = (*f)(pc,n,fields);CHKERRQ(ierr); 440 } 441 PetscFunctionReturn(0); 442 } 443 444 #undef __FUNCT__ 445 #define __FUNCT__ "PCFieldSplitSetBlockSize" 446 /*@ 447 PCFieldSplitSetBlockSize - Sets the block size for defining where fields start in the 448 fieldsplit preconditioner. If not set the matrix block size is used. 449 450 Collective on PC 451 452 Input Parameters: 453 + pc - the preconditioner context 454 - bs - the block size 455 456 Level: intermediate 457 458 .seealso: PCFieldSplitGetSubKSP(), PCFIELDSPLIT, PCFieldSplitSetFields() 459 460 @*/ 461 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitSetBlockSize(PC pc,PetscInt bs) 462 { 463 PetscErrorCode ierr,(*f)(PC,PetscInt); 464 465 PetscFunctionBegin; 466 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 467 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetBlockSize_C",(void (**)(void))&f);CHKERRQ(ierr); 468 if (f) { 469 ierr = (*f)(pc,bs);CHKERRQ(ierr); 470 } 471 PetscFunctionReturn(0); 472 } 473 474 #undef __FUNCT__ 475 #define __FUNCT__ "PCFieldSplitGetSubKSP" 476 /*@C 477 PCFieldSplitGetSubKSP - Gets the KSP contexts for all splits 478 479 Collective on KSP 480 481 Input Parameter: 482 . pc - the preconditioner context 483 484 Output Parameters: 485 + n - the number of split 486 - pc - the array of KSP contexts 487 488 Note: 489 After PCFieldSplitGetSubKSP() the array of KSPs IS to be freed 490 491 You must call KSPSetUp() before calling PCFieldSplitGetSubKSP(). 492 493 Level: advanced 494 495 .seealso: PCFIELDSPLIT 496 @*/ 497 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitGetSubKSP(PC pc,PetscInt *n,KSP *subksp[]) 498 { 499 PetscErrorCode ierr,(*f)(PC,PetscInt*,KSP **); 500 501 PetscFunctionBegin; 502 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 503 PetscValidIntPointer(n,2); 504 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitGetSubKSP_C",(void (**)(void))&f);CHKERRQ(ierr); 505 if (f) { 506 ierr = (*f)(pc,n,subksp);CHKERRQ(ierr); 507 } else { 508 SETERRQ(PETSC_ERR_ARG_WRONG,"Cannot get subksp for this type of PC"); 509 } 510 PetscFunctionReturn(0); 511 } 512 513 EXTERN_C_BEGIN 514 #undef __FUNCT__ 515 #define __FUNCT__ "PCFieldSplitSetType_FieldSplit" 516 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitSetType_FieldSplit(PC pc,PCCompositeType type) 517 { 518 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 519 520 PetscFunctionBegin; 521 jac->type = type; 522 PetscFunctionReturn(0); 523 } 524 EXTERN_C_END 525 526 EXTERN_C_BEGIN 527 #undef __FUNCT__ 528 #define __FUNCT__ "PCFieldSplitSetBlockSize_FieldSplit" 529 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitSetBlockSize_FieldSplit(PC pc,PetscInt bs) 530 { 531 PC_FieldSplit *jac = (PC_FieldSplit*)pc->data; 532 533 PetscFunctionBegin; 534 if (bs < 1) SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE,"Blocksize must be positive, you gave %D",bs); 535 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); 536 jac->bs = bs; 537 PetscFunctionReturn(0); 538 } 539 EXTERN_C_END 540 541 #undef __FUNCT__ 542 #define __FUNCT__ "PCFieldSplitSetType" 543 /*@C 544 PCFieldSplitSetType - Sets the type of fieldsplit preconditioner. 545 546 Collective on PC 547 548 Input Parameter: 549 . pc - the preconditioner context 550 . type - PC_COMPOSITE_ADDITIVE (default), PC_COMPOSITE_MULTIPLICATIVE 551 552 Options Database Key: 553 . -pc_fieldsplit_type <type: one of multiplicative, additive, special> - Sets fieldsplit preconditioner type 554 555 Level: Developer 556 557 .keywords: PC, set, type, composite preconditioner, additive, multiplicative 558 559 .seealso: PCCompositeSetType() 560 561 @*/ 562 PetscErrorCode PETSCKSP_DLLEXPORT PCFieldSplitSetType(PC pc,PCCompositeType type) 563 { 564 PetscErrorCode ierr,(*f)(PC,PCCompositeType); 565 566 PetscFunctionBegin; 567 PetscValidHeaderSpecific(pc,PC_COOKIE,1); 568 ierr = PetscObjectQueryFunction((PetscObject)pc,"PCFieldSplitSetType_C",(void (**)(void))&f);CHKERRQ(ierr); 569 if (f) { 570 ierr = (*f)(pc,type);CHKERRQ(ierr); 571 } 572 PetscFunctionReturn(0); 573 } 574 575 /* -------------------------------------------------------------------------------------*/ 576 /*MC 577 PCFIELDSPLIT - Preconditioner created by combining separate preconditioners for individual 578 fields or groups of fields 579 580 581 To set options on the solvers for each block append -sub_ to all the PC 582 options database keys. For example, -sub_pc_type ilu -sub_pc_factor_levels 1 583 584 To set the options on the solvers separate for each block call PCFieldSplitGetSubKSP() 585 and set the options directly on the resulting KSP object 586 587 Level: intermediate 588 589 Options Database Keys: 590 + -pc_splitfield_%d_fields <a,b,..> - indicates the fields to be used in the %d'th split 591 . -pc_splitfield_default - automatically add any fields to additional splits that have not 592 been supplied explicitly by -pc_splitfield_%d_fields 593 . -pc_splitfield_block_size <bs> - size of block that defines fields (i.e. there are bs fields) 594 - -pc_splitfield_type <additive,multiplicative> 595 596 Concepts: physics based preconditioners 597 598 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 599 PCFieldSplitGetSubKSP(), PCFieldSplitSetFields(),PCFieldSplitSetType() 600 M*/ 601 602 EXTERN_C_BEGIN 603 #undef __FUNCT__ 604 #define __FUNCT__ "PCCreate_FieldSplit" 605 PetscErrorCode PETSCKSP_DLLEXPORT PCCreate_FieldSplit(PC pc) 606 { 607 PetscErrorCode ierr; 608 PC_FieldSplit *jac; 609 610 PetscFunctionBegin; 611 ierr = PetscNew(PC_FieldSplit,&jac);CHKERRQ(ierr); 612 ierr = PetscLogObjectMemory(pc,sizeof(PC_FieldSplit));CHKERRQ(ierr); 613 jac->bs = -1; 614 jac->nsplits = 0; 615 jac->type = PC_COMPOSITE_ADDITIVE; 616 617 pc->data = (void*)jac; 618 619 pc->ops->apply = PCApply_FieldSplit; 620 pc->ops->setup = PCSetUp_FieldSplit; 621 pc->ops->destroy = PCDestroy_FieldSplit; 622 pc->ops->setfromoptions = PCSetFromOptions_FieldSplit; 623 pc->ops->view = PCView_FieldSplit; 624 pc->ops->applyrichardson = 0; 625 626 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitGetSubKSP_C","PCFieldSplitGetSubKSP_FieldSplit", 627 PCFieldSplitGetSubKSP_FieldSplit);CHKERRQ(ierr); 628 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetFields_C","PCFieldSplitSetFields_FieldSplit", 629 PCFieldSplitSetFields_FieldSplit);CHKERRQ(ierr); 630 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetType_C","PCFieldSplitSetType_FieldSplit", 631 PCFieldSplitSetType_FieldSplit);CHKERRQ(ierr); 632 ierr = PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCFieldSplitSetBlockSize_C","PCFieldSplitSetBlockSize_FieldSplit", 633 PCFieldSplitSetBlockSize_FieldSplit);CHKERRQ(ierr); 634 PetscFunctionReturn(0); 635 } 636 EXTERN_C_END 637 638 639