Actual source code: precon.c

      2: /*
      3:     The PC (preconditioner) interface routines, callable by users.
      4: */
      5: #include <petsc/private/pcimpl.h>
      6: #include <petscdm.h>

      8: /* Logging support */
      9: PetscClassId  PC_CLASSID;
     10: PetscLogEvent PC_SetUp, PC_SetUpOnBlocks, PC_Apply, PC_MatApply, PC_ApplyCoarse, PC_ApplyMultiple, PC_ApplySymmetricLeft;
     11: PetscLogEvent PC_ApplySymmetricRight, PC_ModifySubMatrices, PC_ApplyOnBlocks, PC_ApplyTransposeOnBlocks;
     12: PetscInt      PetscMGLevelId;

     14: PetscErrorCode PCGetDefaultType_Private(PC pc,const char *type[])
     15: {
     17:   PetscMPIInt    size;
     18:   PetscBool      hasop,flg1,flg2,set,flg3;

     21:   MPI_Comm_size(PetscObjectComm((PetscObject)pc),&size);
     22:   if (pc->pmat) {
     23:     MatHasOperation(pc->pmat,MATOP_GET_DIAGONAL_BLOCK,&hasop);
     24:     if (size == 1) {
     25:       MatGetFactorAvailable(pc->pmat,"petsc",MAT_FACTOR_ICC,&flg1);
     26:       MatGetFactorAvailable(pc->pmat,"petsc",MAT_FACTOR_ILU,&flg2);
     27:       MatIsSymmetricKnown(pc->pmat,&set,&flg3);
     28:       if (flg1 && (!flg2 || (set && flg3))) {
     29:         *type = PCICC;
     30:       } else if (flg2) {
     31:         *type = PCILU;
     32:       } else if (hasop) { /* likely is a parallel matrix run on one processor */
     33:         *type = PCBJACOBI;
     34:       } else {
     35:         *type = PCNONE;
     36:       }
     37:     } else {
     38:        if (hasop) {
     39:         *type = PCBJACOBI;
     40:       } else {
     41:         *type = PCNONE;
     42:       }
     43:     }
     44:   } else {
     45:     if (size == 1) {
     46:       *type = PCILU;
     47:     } else {
     48:       *type = PCBJACOBI;
     49:     }
     50:   }
     51:   return(0);
     52: }

     54: /*@
     55:    PCReset - Resets a PC context to the pcsetupcalled = 0 state and removes any allocated Vecs and Mats

     57:    Collective on PC

     59:    Input Parameter:
     60: .  pc - the preconditioner context

     62:    Level: developer

     64:    Notes:
     65:     This allows a PC to be reused for a different sized linear system but using the same options that have been previously set in the PC

     67: .seealso: PCCreate(), PCSetUp()
     68: @*/
     69: PetscErrorCode  PCReset(PC pc)
     70: {

     75:   if (pc->ops->reset) {
     76:     (*pc->ops->reset)(pc);
     77:   }
     78:   VecDestroy(&pc->diagonalscaleright);
     79:   VecDestroy(&pc->diagonalscaleleft);
     80:   MatDestroy(&pc->pmat);
     81:   MatDestroy(&pc->mat);

     83:   pc->setupcalled = 0;
     84:   return(0);
     85: }

     87: /*@C
     88:    PCDestroy - Destroys PC context that was created with PCCreate().

     90:    Collective on PC

     92:    Input Parameter:
     93: .  pc - the preconditioner context

     95:    Level: developer

     97: .seealso: PCCreate(), PCSetUp()
     98: @*/
     99: PetscErrorCode  PCDestroy(PC *pc)
    100: {

    104:   if (!*pc) return(0);
    106:   if (--((PetscObject)(*pc))->refct > 0) {*pc = NULL; return(0);}

    108:   PCReset(*pc);

    110:   /* if memory was published with SAWs then destroy it */
    111:   PetscObjectSAWsViewOff((PetscObject)*pc);
    112:   if ((*pc)->ops->destroy) {(*(*pc)->ops->destroy)((*pc));}
    113:   DMDestroy(&(*pc)->dm);
    114:   PetscHeaderDestroy(pc);
    115:   return(0);
    116: }

    118: /*@C
    119:    PCGetDiagonalScale - Indicates if the preconditioner applies an additional left and right
    120:       scaling as needed by certain time-stepping codes.

    122:    Logically Collective on PC

    124:    Input Parameter:
    125: .  pc - the preconditioner context

    127:    Output Parameter:
    128: .  flag - PETSC_TRUE if it applies the scaling

    130:    Level: developer

    132:    Notes:
    133:     If this returns PETSC_TRUE then the system solved via the Krylov method is
    134: $           D M A D^{-1} y = D M b  for left preconditioning or
    135: $           D A M D^{-1} z = D b for right preconditioning

    137: .seealso: PCCreate(), PCSetUp(), PCDiagonalScaleLeft(), PCDiagonalScaleRight(), PCSetDiagonalScale()
    138: @*/
    139: PetscErrorCode  PCGetDiagonalScale(PC pc,PetscBool  *flag)
    140: {
    144:   *flag = pc->diagonalscale;
    145:   return(0);
    146: }

    148: /*@
    149:    PCSetDiagonalScale - Indicates the left scaling to use to apply an additional left and right
    150:       scaling as needed by certain time-stepping codes.

    152:    Logically Collective on PC

    154:    Input Parameters:
    155: +  pc - the preconditioner context
    156: -  s - scaling vector

    158:    Level: intermediate

    160:    Notes:
    161:     The system solved via the Krylov method is
    162: $           D M A D^{-1} y = D M b  for left preconditioning or
    163: $           D A M D^{-1} z = D b for right preconditioning

    165:    PCDiagonalScaleLeft() scales a vector by D. PCDiagonalScaleRight() scales a vector by D^{-1}.

    167: .seealso: PCCreate(), PCSetUp(), PCDiagonalScaleLeft(), PCDiagonalScaleRight(), PCGetDiagonalScale()
    168: @*/
    169: PetscErrorCode  PCSetDiagonalScale(PC pc,Vec s)
    170: {

    176:   pc->diagonalscale     = PETSC_TRUE;

    178:   PetscObjectReference((PetscObject)s);
    179:   VecDestroy(&pc->diagonalscaleleft);

    181:   pc->diagonalscaleleft = s;

    183:   VecDuplicate(s,&pc->diagonalscaleright);
    184:   VecCopy(s,pc->diagonalscaleright);
    185:   VecReciprocal(pc->diagonalscaleright);
    186:   return(0);
    187: }

    189: /*@
    190:    PCDiagonalScaleLeft - Scales a vector by the left scaling as needed by certain time-stepping codes.

    192:    Logically Collective on PC

    194:    Input Parameters:
    195: +  pc - the preconditioner context
    196: .  in - input vector
    197: -  out - scaled vector (maybe the same as in)

    199:    Level: intermediate

    201:    Notes:
    202:     The system solved via the Krylov method is
    203: $           D M A D^{-1} y = D M b  for left preconditioning or
    204: $           D A M D^{-1} z = D b for right preconditioning

    206:    PCDiagonalScaleLeft() scales a vector by D. PCDiagonalScaleRight() scales a vector by D^{-1}.

    208:    If diagonal scaling is turned off and in is not out then in is copied to out

    210: .seealso: PCCreate(), PCSetUp(), PCDiagonalScaleSet(), PCDiagonalScaleRight(), PCDiagonalScale()
    211: @*/
    212: PetscErrorCode  PCDiagonalScaleLeft(PC pc,Vec in,Vec out)
    213: {

    220:   if (pc->diagonalscale) {
    221:     VecPointwiseMult(out,pc->diagonalscaleleft,in);
    222:   } else if (in != out) {
    223:     VecCopy(in,out);
    224:   }
    225:   return(0);
    226: }

    228: /*@
    229:    PCDiagonalScaleRight - Scales a vector by the right scaling as needed by certain time-stepping codes.

    231:    Logically Collective on PC

    233:    Input Parameters:
    234: +  pc - the preconditioner context
    235: .  in - input vector
    236: -  out - scaled vector (maybe the same as in)

    238:    Level: intermediate

    240:    Notes:
    241:     The system solved via the Krylov method is
    242: $           D M A D^{-1} y = D M b  for left preconditioning or
    243: $           D A M D^{-1} z = D b for right preconditioning

    245:    PCDiagonalScaleLeft() scales a vector by D. PCDiagonalScaleRight() scales a vector by D^{-1}.

    247:    If diagonal scaling is turned off and in is not out then in is copied to out

    249: .seealso: PCCreate(), PCSetUp(), PCDiagonalScaleLeft(), PCDiagonalScaleSet(), PCDiagonalScale()
    250: @*/
    251: PetscErrorCode  PCDiagonalScaleRight(PC pc,Vec in,Vec out)
    252: {

    259:   if (pc->diagonalscale) {
    260:     VecPointwiseMult(out,pc->diagonalscaleright,in);
    261:   } else if (in != out) {
    262:     VecCopy(in,out);
    263:   }
    264:   return(0);
    265: }

    267: /*@
    268:    PCSetUseAmat - Sets a flag to indicate that when the preconditioner needs to apply (part of) the
    269:    operator during the preconditioning process it applies the Amat provided to TSSetRHSJacobian(),
    270:    TSSetIJacobian(), SNESSetJacobian(), KSPSetOperator() or PCSetOperator() not the Pmat.

    272:    Logically Collective on PC

    274:    Input Parameters:
    275: +  pc - the preconditioner context
    276: -  flg - PETSC_TRUE to use the Amat, PETSC_FALSE to use the Pmat (default is false)

    278:    Options Database Key:
    279: .  -pc_use_amat <true,false>

    281:    Notes:
    282:    For the common case in which the linear system matrix and the matrix used to construct the
    283:    preconditioner are identical, this routine is does nothing.

    285:    Level: intermediate

    287: .seealso: PCGetUseAmat(), PCBJACOBI, PGMG, PCFIELDSPLIT, PCCOMPOSITE
    288: @*/
    289: PetscErrorCode  PCSetUseAmat(PC pc,PetscBool flg)
    290: {
    293:   pc->useAmat = flg;
    294:   return(0);
    295: }

    297: /*@
    298:    PCSetErrorIfFailure - Causes PC to generate an error if a FPE, for example a zero pivot, is detected.

    300:    Logically Collective on PC

    302:    Input Parameters:
    303: +  pc - iterative context obtained from PCCreate()
    304: -  flg - PETSC_TRUE indicates you want the error generated

    306:    Level: advanced

    308:    Notes:
    309:     Normally PETSc continues if a linear solver fails due to a failed setup of a preconditioner, you can call KSPGetConvergedReason() after a KSPSolve()
    310:     to determine if it has converged or failed. Or use -ksp_error_if_not_converged to cause the program to terminate as soon as lack of convergence is
    311:     detected.

    313:     This is propagated into KSPs used by this PC, which then propagate it into PCs used by those KSPs

    315: .seealso: PCGetInitialGuessNonzero(), PCSetInitialGuessKnoll(), PCGetInitialGuessKnoll()
    316: @*/
    317: PetscErrorCode  PCSetErrorIfFailure(PC pc,PetscBool flg)
    318: {
    322:   pc->erroriffailure = flg;
    323:   return(0);
    324: }

    326: /*@
    327:    PCGetUseAmat - Gets a flag to indicate that when the preconditioner needs to apply (part of) the
    328:    operator during the preconditioning process it applies the Amat provided to TSSetRHSJacobian(),
    329:    TSSetIJacobian(), SNESSetJacobian(), KSPSetOperator() or PCSetOperator() not the Pmat.

    331:    Logically Collective on PC

    333:    Input Parameter:
    334: .  pc - the preconditioner context

    336:    Output Parameter:
    337: .  flg - PETSC_TRUE to use the Amat, PETSC_FALSE to use the Pmat (default is false)

    339:    Notes:
    340:    For the common case in which the linear system matrix and the matrix used to construct the
    341:    preconditioner are identical, this routine is does nothing.

    343:    Level: intermediate

    345: .seealso: PCSetUseAmat(), PCBJACOBI, PGMG, PCFIELDSPLIT, PCCOMPOSITE
    346: @*/
    347: PetscErrorCode  PCGetUseAmat(PC pc,PetscBool *flg)
    348: {
    351:   *flg = pc->useAmat;
    352:   return(0);
    353: }

    355: /*@
    356:    PCCreate - Creates a preconditioner context.

    358:    Collective

    360:    Input Parameter:
    361: .  comm - MPI communicator

    363:    Output Parameter:
    364: .  pc - location to put the preconditioner context

    366:    Notes:
    367:    The default preconditioner for sparse matrices is PCILU or PCICC with 0 fill on one process and block Jacobi with PCILU or PCICC
    368:    in parallel. For dense matrices it is always PCNONE.

    370:    Level: developer

    372: .seealso: PCSetUp(), PCApply(), PCDestroy()
    373: @*/
    374: PetscErrorCode  PCCreate(MPI_Comm comm,PC *newpc)
    375: {
    376:   PC             pc;

    381:   *newpc = NULL;
    382:   PCInitializePackage();

    384:   PetscHeaderCreate(pc,PC_CLASSID,"PC","Preconditioner","PC",comm,PCDestroy,PCView);

    386:   pc->mat                  = NULL;
    387:   pc->pmat                 = NULL;
    388:   pc->setupcalled          = 0;
    389:   pc->setfromoptionscalled = 0;
    390:   pc->data                 = NULL;
    391:   pc->diagonalscale        = PETSC_FALSE;
    392:   pc->diagonalscaleleft    = NULL;
    393:   pc->diagonalscaleright   = NULL;

    395:   pc->modifysubmatrices  = NULL;
    396:   pc->modifysubmatricesP = NULL;

    398:   *newpc = pc;
    399:   return(0);

    401: }

    403: /* -------------------------------------------------------------------------------*/

    405: /*@
    406:    PCApply - Applies the preconditioner to a vector.

    408:    Collective on PC

    410:    Input Parameters:
    411: +  pc - the preconditioner context
    412: -  x - input vector

    414:    Output Parameter:
    415: .  y - output vector

    417:    Level: developer

    419: .seealso: PCApplyTranspose(), PCApplyBAorAB()
    420: @*/
    421: PetscErrorCode  PCApply(PC pc,Vec x,Vec y)
    422: {
    424:   PetscInt       m,n,mv,nv;

    430:   if (x == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
    431:   if (pc->erroriffailure) {VecValidValues(x,2,PETSC_TRUE);}
    432:   /* use pmat to check vector sizes since for KSPLSQR the pmat may be of a different size than mat */
    433:   MatGetLocalSize(pc->pmat,&m,&n);
    434:   VecGetLocalSize(x,&mv);
    435:   VecGetLocalSize(y,&nv);
    436:   /* check pmat * y = x is feasible */
    437:   if (mv != m) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Preconditioner number of local rows %D does not equal input vector size %D",m,mv);
    438:   if (nv != n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Preconditioner number of local columns %D does not equal output vector size %D",n,nv);
    439:   VecSetErrorIfLocked(y,3);

    441:   PCSetUp(pc);
    442:   if (!pc->ops->apply) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PC does not have apply");
    443:   VecLockReadPush(x);
    444:   PetscLogEventBegin(PC_Apply,pc,x,y,0);
    445:   (*pc->ops->apply)(pc,x,y);
    446:   PetscLogEventEnd(PC_Apply,pc,x,y,0);
    447:   if (pc->erroriffailure) {VecValidValues(y,3,PETSC_FALSE);}
    448:   VecLockReadPop(x);
    449:   return(0);
    450: }

    452: /*@
    453:    PCMatApply - Applies the preconditioner to multiple vectors stored as a MATDENSE. Like PCApply(), Y and X must be different matrices.

    455:    Collective on PC

    457:    Input Parameters:
    458: +  pc - the preconditioner context
    459: -  X - block of input vectors

    461:    Output Parameter:
    462: .  Y - block of output vectors

    464:    Level: developer

    466: .seealso: PCApply(), KSPMatSolve()
    467: @*/
    468: PetscErrorCode  PCMatApply(PC pc,Mat X,Mat Y)
    469: {
    470:   Mat            A;
    471:   Vec            cy, cx;
    472:   PetscInt       m1, M1, m2, M2, n1, N1, n2, N2, m3, M3, n3, N3;
    473:   PetscBool      match;

    482:   if (Y == X) SETERRQ(PetscObjectComm((PetscObject)pc), PETSC_ERR_ARG_IDN, "Y and X must be different matrices");
    483:   PCGetOperators(pc, NULL, &A);
    484:   MatGetLocalSize(A, &m3, &n3);
    485:   MatGetLocalSize(X, &m2, &n2);
    486:   MatGetLocalSize(Y, &m1, &n1);
    487:   MatGetSize(A, &M3, &N3);
    488:   MatGetSize(X, &M2, &N2);
    489:   MatGetSize(Y, &M1, &N1);
    490:   if (n1 != n2 || N1 != N2) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Incompatible number of columns between block of input vectors (n,N) = (%D,%D) and block of output vectors (n,N) = (%D,%D)", n2, N2, n1, N1);
    491:   if (m2 != m3 || M2 != M3) SETERRQ6(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Incompatible layout between block of input vectors (m,M) = (%D,%D) and Pmat (m,M)x(n,N) = (%D,%D)x(%D,%D)", m2, M2, m3, M3, n3, N3);
    492:   if (m1 != n3 || M1 != N3) SETERRQ6(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Incompatible layout between block of output vectors (m,M) = (%D,%D) and Pmat (m,M)x(n,N) = (%D,%D)x(%D,%D)", m1, M1, m3, M3, n3, N3);
    493:   PetscObjectBaseTypeCompareAny((PetscObject)Y, &match, MATSEQDENSE, MATMPIDENSE, "");
    494:   if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Provided block of output vectors not stored in a dense Mat");
    495:   PetscObjectBaseTypeCompareAny((PetscObject)X, &match, MATSEQDENSE, MATMPIDENSE, "");
    496:   if (!match) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Provided block of input vectors not stored in a dense Mat");
    497:   PCSetUp(pc);
    498:   if (pc->ops->matapply) {
    499:     PetscLogEventBegin(PC_MatApply, pc, X, Y, 0);
    500:     (*pc->ops->matapply)(pc, X, Y);
    501:     PetscLogEventEnd(PC_MatApply, pc, X, Y, 0);
    502:   } else {
    503:     PetscInfo1(pc, "PC type %s applying column by column\n", ((PetscObject)pc)->type_name);
    504:     for (n1 = 0; n1 < N1; ++n1) {
    505:       MatDenseGetColumnVecRead(X, n1, &cx);
    506:       MatDenseGetColumnVecWrite(Y, n1, &cy);
    507:       PCApply(pc, cx, cy);
    508:       MatDenseRestoreColumnVecWrite(Y, n1, &cy);
    509:       MatDenseRestoreColumnVecRead(X, n1, &cx);
    510:     }
    511:   }
    512:   return(0);
    513: }

    515: /*@
    516:    PCApplySymmetricLeft - Applies the left part of a symmetric preconditioner to a vector.

    518:    Collective on PC

    520:    Input Parameters:
    521: +  pc - the preconditioner context
    522: -  x - input vector

    524:    Output Parameter:
    525: .  y - output vector

    527:    Notes:
    528:    Currently, this routine is implemented only for PCICC and PCJACOBI preconditioners.

    530:    Level: developer

    532: .seealso: PCApply(), PCApplySymmetricRight()
    533: @*/
    534: PetscErrorCode  PCApplySymmetricLeft(PC pc,Vec x,Vec y)
    535: {

    542:   if (x == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
    543:   if (pc->erroriffailure) {VecValidValues(x,2,PETSC_TRUE);}
    544:   PCSetUp(pc);
    545:   if (!pc->ops->applysymmetricleft) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PC does not have left symmetric apply");
    546:   VecLockReadPush(x);
    547:   PetscLogEventBegin(PC_ApplySymmetricLeft,pc,x,y,0);
    548:   (*pc->ops->applysymmetricleft)(pc,x,y);
    549:   PetscLogEventEnd(PC_ApplySymmetricLeft,pc,x,y,0);
    550:   VecLockReadPop(x);
    551:   if (pc->erroriffailure) {VecValidValues(y,3,PETSC_FALSE);}
    552:   return(0);
    553: }

    555: /*@
    556:    PCApplySymmetricRight - Applies the right part of a symmetric preconditioner to a vector.

    558:    Collective on PC

    560:    Input Parameters:
    561: +  pc - the preconditioner context
    562: -  x - input vector

    564:    Output Parameter:
    565: .  y - output vector

    567:    Level: developer

    569:    Notes:
    570:    Currently, this routine is implemented only for PCICC and PCJACOBI preconditioners.

    572: .seealso: PCApply(), PCApplySymmetricLeft()
    573: @*/
    574: PetscErrorCode  PCApplySymmetricRight(PC pc,Vec x,Vec y)
    575: {

    582:   if (x == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
    583:   if (pc->erroriffailure) {VecValidValues(x,2,PETSC_TRUE);}
    584:   PCSetUp(pc);
    585:   if (!pc->ops->applysymmetricright) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PC does not have left symmetric apply");
    586:   VecLockReadPush(x);
    587:   PetscLogEventBegin(PC_ApplySymmetricRight,pc,x,y,0);
    588:   (*pc->ops->applysymmetricright)(pc,x,y);
    589:   PetscLogEventEnd(PC_ApplySymmetricRight,pc,x,y,0);
    590:   VecLockReadPop(x);
    591:   if (pc->erroriffailure) {VecValidValues(y,3,PETSC_FALSE);}
    592:   return(0);
    593: }

    595: /*@
    596:    PCApplyTranspose - Applies the transpose of preconditioner to a vector.

    598:    Collective on PC

    600:    Input Parameters:
    601: +  pc - the preconditioner context
    602: -  x - input vector

    604:    Output Parameter:
    605: .  y - output vector

    607:    Notes:
    608:     For complex numbers this applies the non-Hermitian transpose.

    610:    Developer Notes:
    611:     We need to implement a PCApplyHermitianTranspose()

    613:    Level: developer

    615: .seealso: PCApply(), PCApplyBAorAB(), PCApplyBAorABTranspose(), PCApplyTransposeExists()
    616: @*/
    617: PetscErrorCode  PCApplyTranspose(PC pc,Vec x,Vec y)
    618: {

    625:   if (x == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
    626:   if (pc->erroriffailure) {VecValidValues(x,2,PETSC_TRUE);}
    627:   PCSetUp(pc);
    628:   if (!pc->ops->applytranspose) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PC does not have apply transpose");
    629:   VecLockReadPush(x);
    630:   PetscLogEventBegin(PC_Apply,pc,x,y,0);
    631:   (*pc->ops->applytranspose)(pc,x,y);
    632:   PetscLogEventEnd(PC_Apply,pc,x,y,0);
    633:   VecLockReadPop(x);
    634:   if (pc->erroriffailure) {VecValidValues(y,3,PETSC_FALSE);}
    635:   return(0);
    636: }

    638: /*@
    639:    PCApplyTransposeExists - Test whether the preconditioner has a transpose apply operation

    641:    Collective on PC

    643:    Input Parameters:
    644: .  pc - the preconditioner context

    646:    Output Parameter:
    647: .  flg - PETSC_TRUE if a transpose operation is defined

    649:    Level: developer

    651: .seealso: PCApplyTranspose()
    652: @*/
    653: PetscErrorCode  PCApplyTransposeExists(PC pc,PetscBool  *flg)
    654: {
    658:   if (pc->ops->applytranspose) *flg = PETSC_TRUE;
    659:   else *flg = PETSC_FALSE;
    660:   return(0);
    661: }

    663: /*@
    664:    PCApplyBAorAB - Applies the preconditioner and operator to a vector. y = B*A*x or y = A*B*x.

    666:    Collective on PC

    668:    Input Parameters:
    669: +  pc - the preconditioner context
    670: .  side - indicates the preconditioner side, one of PC_LEFT, PC_RIGHT, or PC_SYMMETRIC
    671: .  x - input vector
    672: -  work - work vector

    674:    Output Parameter:
    675: .  y - output vector

    677:    Level: developer

    679:    Notes:
    680:     If the PC has had PCSetDiagonalScale() set then D M A D^{-1} for left preconditioning or  D A M D^{-1} is actually applied. Note that the
    681:    specific KSPSolve() method must also be written to handle the post-solve "correction" for the diagonal scaling.

    683: .seealso: PCApply(), PCApplyTranspose(), PCApplyBAorABTranspose()
    684: @*/
    685: PetscErrorCode  PCApplyBAorAB(PC pc,PCSide side,Vec x,Vec y,Vec work)
    686: {

    698:   if (x == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
    699:   if (side != PC_LEFT && side != PC_SYMMETRIC && side != PC_RIGHT) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Side must be right, left, or symmetric");
    700:   if (pc->diagonalscale && side == PC_SYMMETRIC) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot include diagonal scaling with symmetric preconditioner application");
    701:   if (pc->erroriffailure) {VecValidValues(x,3,PETSC_TRUE);}

    703:   PCSetUp(pc);
    704:   if (pc->diagonalscale) {
    705:     if (pc->ops->applyBA) {
    706:       Vec work2; /* this is expensive, but to fix requires a second work vector argument to PCApplyBAorAB() */
    707:       VecDuplicate(x,&work2);
    708:       PCDiagonalScaleRight(pc,x,work2);
    709:       (*pc->ops->applyBA)(pc,side,work2,y,work);
    710:       PCDiagonalScaleLeft(pc,y,y);
    711:       VecDestroy(&work2);
    712:     } else if (side == PC_RIGHT) {
    713:       PCDiagonalScaleRight(pc,x,y);
    714:       PCApply(pc,y,work);
    715:       MatMult(pc->mat,work,y);
    716:       PCDiagonalScaleLeft(pc,y,y);
    717:     } else if (side == PC_LEFT) {
    718:       PCDiagonalScaleRight(pc,x,y);
    719:       MatMult(pc->mat,y,work);
    720:       PCApply(pc,work,y);
    721:       PCDiagonalScaleLeft(pc,y,y);
    722:     } else if (side == PC_SYMMETRIC) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot provide diagonal scaling with symmetric application of preconditioner");
    723:   } else {
    724:     if (pc->ops->applyBA) {
    725:       (*pc->ops->applyBA)(pc,side,x,y,work);
    726:     } else if (side == PC_RIGHT) {
    727:       PCApply(pc,x,work);
    728:       MatMult(pc->mat,work,y);
    729:     } else if (side == PC_LEFT) {
    730:       MatMult(pc->mat,x,work);
    731:       PCApply(pc,work,y);
    732:     } else if (side == PC_SYMMETRIC) {
    733:       /* There's an extra copy here; maybe should provide 2 work vectors instead? */
    734:       PCApplySymmetricRight(pc,x,work);
    735:       MatMult(pc->mat,work,y);
    736:       VecCopy(y,work);
    737:       PCApplySymmetricLeft(pc,work,y);
    738:     }
    739:   }
    740:   if (pc->erroriffailure) {VecValidValues(y,4,PETSC_FALSE);}
    741:   return(0);
    742: }

    744: /*@
    745:    PCApplyBAorABTranspose - Applies the transpose of the preconditioner
    746:    and operator to a vector. That is, applies tr(B) * tr(A) with left preconditioning,
    747:    NOT tr(B*A) = tr(A)*tr(B).

    749:    Collective on PC

    751:    Input Parameters:
    752: +  pc - the preconditioner context
    753: .  side - indicates the preconditioner side, one of PC_LEFT, PC_RIGHT, or PC_SYMMETRIC
    754: .  x - input vector
    755: -  work - work vector

    757:    Output Parameter:
    758: .  y - output vector

    761:    Notes:
    762:     this routine is used internally so that the same Krylov code can be used to solve A x = b and A' x = b, with a preconditioner
    763:       defined by B'. This is why this has the funny form that it computes tr(B) * tr(A)

    765:     Level: developer

    767: .seealso: PCApply(), PCApplyTranspose(), PCApplyBAorAB()
    768: @*/
    769: PetscErrorCode  PCApplyBAorABTranspose(PC pc,PCSide side,Vec x,Vec y,Vec work)
    770: {

    778:   if (x == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"x and y must be different vectors");
    779:   if (pc->erroriffailure) {VecValidValues(x,3,PETSC_TRUE);}
    780:   if (pc->ops->applyBAtranspose) {
    781:     (*pc->ops->applyBAtranspose)(pc,side,x,y,work);
    782:     if (pc->erroriffailure) {VecValidValues(y,4,PETSC_FALSE);}
    783:     return(0);
    784:   }
    785:   if (side != PC_LEFT && side != PC_RIGHT) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_OUTOFRANGE,"Side must be right or left");

    787:   PCSetUp(pc);
    788:   if (side == PC_RIGHT) {
    789:     PCApplyTranspose(pc,x,work);
    790:     MatMultTranspose(pc->mat,work,y);
    791:   } else if (side == PC_LEFT) {
    792:     MatMultTranspose(pc->mat,x,work);
    793:     PCApplyTranspose(pc,work,y);
    794:   }
    795:   /* add support for PC_SYMMETRIC */
    796:   if (pc->erroriffailure) {VecValidValues(y,4,PETSC_FALSE);}
    797:   return(0);
    798: }

    800: /* -------------------------------------------------------------------------------*/

    802: /*@
    803:    PCApplyRichardsonExists - Determines whether a particular preconditioner has a
    804:    built-in fast application of Richardson's method.

    806:    Not Collective

    808:    Input Parameter:
    809: .  pc - the preconditioner

    811:    Output Parameter:
    812: .  exists - PETSC_TRUE or PETSC_FALSE

    814:    Level: developer

    816: .seealso: PCApplyRichardson()
    817: @*/
    818: PetscErrorCode  PCApplyRichardsonExists(PC pc,PetscBool  *exists)
    819: {
    823:   if (pc->ops->applyrichardson) *exists = PETSC_TRUE;
    824:   else *exists = PETSC_FALSE;
    825:   return(0);
    826: }

    828: /*@
    829:    PCApplyRichardson - Applies several steps of Richardson iteration with
    830:    the particular preconditioner. This routine is usually used by the
    831:    Krylov solvers and not the application code directly.

    833:    Collective on PC

    835:    Input Parameters:
    836: +  pc  - the preconditioner context
    837: .  b   - the right hand side
    838: .  w   - one work vector
    839: .  rtol - relative decrease in residual norm convergence criteria
    840: .  abstol - absolute residual norm convergence criteria
    841: .  dtol - divergence residual norm increase criteria
    842: .  its - the number of iterations to apply.
    843: -  guesszero - if the input x contains nonzero initial guess

    845:    Output Parameter:
    846: +  outits - number of iterations actually used (for SOR this always equals its)
    847: .  reason - the reason the apply terminated
    848: -  y - the solution (also contains initial guess if guesszero is PETSC_FALSE

    850:    Notes:
    851:    Most preconditioners do not support this function. Use the command
    852:    PCApplyRichardsonExists() to determine if one does.

    854:    Except for the multigrid PC this routine ignores the convergence tolerances
    855:    and always runs for the number of iterations

    857:    Level: developer

    859: .seealso: PCApplyRichardsonExists()
    860: @*/
    861: PetscErrorCode  PCApplyRichardson(PC pc,Vec b,Vec y,Vec w,PetscReal rtol,PetscReal abstol, PetscReal dtol,PetscInt its,PetscBool guesszero,PetscInt *outits,PCRichardsonConvergedReason *reason)
    862: {

    870:   if (b == y) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_IDN,"b and y must be different vectors");
    871:   PCSetUp(pc);
    872:   if (!pc->ops->applyrichardson) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PC does not have apply richardson");
    873:   (*pc->ops->applyrichardson)(pc,b,y,w,rtol,abstol,dtol,its,guesszero,outits,reason);
    874:   return(0);
    875: }

    877: /*@
    878:    PCSetFailedReason - Sets the reason a PCSetUp() failed or PC_NOERROR if it did not fail

    880:    Logically Collective on PC

    882:    Input Parameter:
    883: +  pc - the preconditioner context
    884: -  reason - the reason it failedx

    886:    Level: advanced

    888: .seealso: PCCreate(), PCApply(), PCDestroy(), PCFailedReason
    889: @*/
    890: PetscErrorCode PCSetFailedReason(PC pc,PCFailedReason reason)
    891: {
    893:   pc->failedreason = reason;
    894:   return(0);
    895: }

    897: /*@
    898:    PCGetFailedReason - Gets the reason a PCSetUp() failed or PC_NOERROR if it did not fail

    900:    Logically Collective on PC

    902:    Input Parameter:
    903: .  pc - the preconditioner context

    905:    Output Parameter:
    906: .  reason - the reason it failed

    908:    Level: advanced

    910:    Notes: This is the maximum over reason over all ranks in the PC communicator. It is only valid after
    911:    a call KSPCheckDot() or  KSPCheckNorm() inside a KSPSolve(). It is not valid immediately after a PCSetUp()
    912:    or PCApply(), then use PCGetFailedReasonRank()

    914: .seealso: PCCreate(), PCApply(), PCDestroy(), PCGetFailedReasonRank(), PCSetFailedReason()
    915: @*/
    916: PetscErrorCode PCGetFailedReason(PC pc,PCFailedReason *reason)
    917: {
    919:   if (pc->setupcalled < 0) *reason = (PCFailedReason)pc->setupcalled;
    920:   else *reason = pc->failedreason;
    921:   return(0);
    922: }

    924: /*@
    925:    PCGetFailedReasonRank - Gets the reason a PCSetUp() failed or PC_NOERROR if it did not fail on this MPI rank

    927:   Not Collective on PC

    929:    Input Parameter:
    930: .  pc - the preconditioner context

    932:    Output Parameter:
    933: .  reason - the reason it failed

    935:    Notes:
    936:      Different ranks may have different reasons or no reason, see PCGetFailedReason()

    938:    Level: advanced

    940: .seealso: PCCreate(), PCApply(), PCDestroy(), PCGetFailedReason(), PCSetFailedReason()
    941: @*/
    942: PetscErrorCode PCGetFailedReasonRank(PC pc,PCFailedReason *reason)
    943: {
    945:   if (pc->setupcalled < 0) *reason = (PCFailedReason)pc->setupcalled;
    946:   else *reason = pc->failedreason;
    947:   return(0);
    948: }

    950: /*  Next line needed to deactivate KSP_Solve logging */
    951: #include <petsc/private/kspimpl.h>

    953: /*
    954:       a setupcall of 0 indicates never setup,
    955:                      1 indicates has been previously setup
    956:                     -1 indicates a PCSetUp() was attempted and failed
    957: */
    958: /*@
    959:    PCSetUp - Prepares for the use of a preconditioner.

    961:    Collective on PC

    963:    Input Parameter:
    964: .  pc - the preconditioner context

    966:    Level: developer

    968: .seealso: PCCreate(), PCApply(), PCDestroy()
    969: @*/
    970: PetscErrorCode  PCSetUp(PC pc)
    971: {
    972:   PetscErrorCode   ierr;
    973:   const char       *def;
    974:   PetscObjectState matstate, matnonzerostate;

    978:   if (!pc->mat) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_ARG_WRONGSTATE,"Matrix must be set first");

    980:   if (pc->setupcalled && pc->reusepreconditioner) {
    981:     PetscInfo(pc,"Leaving PC with identical preconditioner since reuse preconditioner is set\n");
    982:     return(0);
    983:   }

    985:   PetscObjectStateGet((PetscObject)pc->pmat,&matstate);
    986:   MatGetNonzeroState(pc->pmat,&matnonzerostate);
    987:   if (!pc->setupcalled) {
    988:     PetscInfo(pc,"Setting up PC for first time\n");
    989:     pc->flag = DIFFERENT_NONZERO_PATTERN;
    990:   } else if (matstate == pc->matstate) {
    991:     PetscInfo(pc,"Leaving PC with identical preconditioner since operator is unchanged\n");
    992:     return(0);
    993:   } else {
    994:     if (matnonzerostate > pc->matnonzerostate) {
    995:       PetscInfo(pc,"Setting up PC with different nonzero pattern\n");
    996:       pc->flag = DIFFERENT_NONZERO_PATTERN;
    997:     } else {
    998:       PetscInfo(pc,"Setting up PC with same nonzero pattern\n");
    999:       pc->flag = SAME_NONZERO_PATTERN;
    1000:     }
    1001:   }
    1002:   pc->matstate        = matstate;
    1003:   pc->matnonzerostate = matnonzerostate;

    1005:   if (!((PetscObject)pc)->type_name) {
    1006:     PCGetDefaultType_Private(pc,&def);
    1007:     PCSetType(pc,def);
    1008:   }

    1010:   MatSetErrorIfFailure(pc->pmat,pc->erroriffailure);
    1011:   MatSetErrorIfFailure(pc->mat,pc->erroriffailure);
    1012:   PetscLogEventBegin(PC_SetUp,pc,0,0,0);
    1013:   if (pc->ops->setup) {
    1014:     /* do not log solves and applications of preconditioners while constructing preconditioners; perhaps they should be logged separately from the regular solves */
    1015:     PetscLogEventDeactivatePush(KSP_Solve);
    1016:     PetscLogEventDeactivatePush(PC_Apply);
    1017:     (*pc->ops->setup)(pc);
    1018:     PetscLogEventDeactivatePop(KSP_Solve);
    1019:     PetscLogEventDeactivatePop(PC_Apply);
    1020:   }
    1021:   PetscLogEventEnd(PC_SetUp,pc,0,0,0);
    1022:   if (!pc->setupcalled) pc->setupcalled = 1;
    1023:   return(0);
    1024: }

    1026: /*@
    1027:    PCSetUpOnBlocks - Sets up the preconditioner for each block in
    1028:    the block Jacobi, block Gauss-Seidel, and overlapping Schwarz
    1029:    methods.

    1031:    Collective on PC

    1033:    Input Parameters:
    1034: .  pc - the preconditioner context

    1036:    Level: developer

    1038: .seealso: PCCreate(), PCApply(), PCDestroy(), PCSetUp()
    1039: @*/
    1040: PetscErrorCode  PCSetUpOnBlocks(PC pc)
    1041: {

    1046:   if (!pc->ops->setuponblocks) return(0);
    1047:   PetscLogEventBegin(PC_SetUpOnBlocks,pc,0,0,0);
    1048:   (*pc->ops->setuponblocks)(pc);
    1049:   PetscLogEventEnd(PC_SetUpOnBlocks,pc,0,0,0);
    1050:   return(0);
    1051: }

    1053: /*@C
    1054:    PCSetModifySubMatrices - Sets a user-defined routine for modifying the
    1055:    submatrices that arise within certain subdomain-based preconditioners.
    1056:    The basic submatrices are extracted from the preconditioner matrix as
    1057:    usual; the user can then alter these (for example, to set different boundary
    1058:    conditions for each submatrix) before they are used for the local solves.

    1060:    Logically Collective on PC

    1062:    Input Parameters:
    1063: +  pc - the preconditioner context
    1064: .  func - routine for modifying the submatrices
    1065: -  ctx - optional user-defined context (may be null)

    1067:    Calling sequence of func:
    1068: $     func (PC pc,PetscInt nsub,IS *row,IS *col,Mat *submat,void *ctx);

    1070: +  row - an array of index sets that contain the global row numbers
    1071:          that comprise each local submatrix
    1072: .  col - an array of index sets that contain the global column numbers
    1073:          that comprise each local submatrix
    1074: .  submat - array of local submatrices
    1075: -  ctx - optional user-defined context for private data for the
    1076:          user-defined func routine (may be null)

    1078:    Notes:
    1079:    PCSetModifySubMatrices() MUST be called before KSPSetUp() and
    1080:    KSPSolve().

    1082:    A routine set by PCSetModifySubMatrices() is currently called within
    1083:    the block Jacobi (PCBJACOBI) and additive Schwarz (PCASM)
    1084:    preconditioners.  All other preconditioners ignore this routine.

    1086:    Level: advanced

    1088: .seealso: PCModifySubMatrices()
    1089: @*/
    1090: PetscErrorCode  PCSetModifySubMatrices(PC pc,PetscErrorCode (*func)(PC,PetscInt,const IS[],const IS[],Mat[],void*),void *ctx)
    1091: {
    1094:   pc->modifysubmatrices  = func;
    1095:   pc->modifysubmatricesP = ctx;
    1096:   return(0);
    1097: }

    1099: /*@C
    1100:    PCModifySubMatrices - Calls an optional user-defined routine within
    1101:    certain preconditioners if one has been set with PCSetModifySubMatrices().

    1103:    Collective on PC

    1105:    Input Parameters:
    1106: +  pc - the preconditioner context
    1107: .  nsub - the number of local submatrices
    1108: .  row - an array of index sets that contain the global row numbers
    1109:          that comprise each local submatrix
    1110: .  col - an array of index sets that contain the global column numbers
    1111:          that comprise each local submatrix
    1112: .  submat - array of local submatrices
    1113: -  ctx - optional user-defined context for private data for the
    1114:          user-defined routine (may be null)

    1116:    Output Parameter:
    1117: .  submat - array of local submatrices (the entries of which may
    1118:             have been modified)

    1120:    Notes:
    1121:    The user should NOT generally call this routine, as it will
    1122:    automatically be called within certain preconditioners (currently
    1123:    block Jacobi, additive Schwarz) if set.

    1125:    The basic submatrices are extracted from the preconditioner matrix
    1126:    as usual; the user can then alter these (for example, to set different
    1127:    boundary conditions for each submatrix) before they are used for the
    1128:    local solves.

    1130:    Level: developer

    1132: .seealso: PCSetModifySubMatrices()
    1133: @*/
    1134: PetscErrorCode  PCModifySubMatrices(PC pc,PetscInt nsub,const IS row[],const IS col[],Mat submat[],void *ctx)
    1135: {

    1140:   if (!pc->modifysubmatrices) return(0);
    1141:   PetscLogEventBegin(PC_ModifySubMatrices,pc,0,0,0);
    1142:   (*pc->modifysubmatrices)(pc,nsub,row,col,submat,ctx);
    1143:   PetscLogEventEnd(PC_ModifySubMatrices,pc,0,0,0);
    1144:   return(0);
    1145: }

    1147: /*@
    1148:    PCSetOperators - Sets the matrix associated with the linear system and
    1149:    a (possibly) different one associated with the preconditioner.

    1151:    Logically Collective on PC

    1153:    Input Parameters:
    1154: +  pc - the preconditioner context
    1155: .  Amat - the matrix that defines the linear system
    1156: -  Pmat - the matrix to be used in constructing the preconditioner, usually the same as Amat.

    1158:    Notes:
    1159:     Passing a NULL for Amat or Pmat removes the matrix that is currently used.

    1161:     If you wish to replace either Amat or Pmat but leave the other one untouched then
    1162:     first call KSPGetOperators() to get the one you wish to keep, call PetscObjectReference()
    1163:     on it and then pass it back in in your call to KSPSetOperators().

    1165:    More Notes about Repeated Solution of Linear Systems:
    1166:    PETSc does NOT reset the matrix entries of either Amat or Pmat
    1167:    to zero after a linear solve; the user is completely responsible for
    1168:    matrix assembly.  See the routine MatZeroEntries() if desiring to
    1169:    zero all elements of a matrix.

    1171:    Level: intermediate

    1173: .seealso: PCGetOperators(), MatZeroEntries()
    1174:  @*/
    1175: PetscErrorCode  PCSetOperators(PC pc,Mat Amat,Mat Pmat)
    1176: {
    1177:   PetscErrorCode   ierr;
    1178:   PetscInt         m1,n1,m2,n2;

    1186:   if (pc->setupcalled && pc->mat && pc->pmat && Amat && Pmat) {
    1187:     MatGetLocalSize(Amat,&m1,&n1);
    1188:     MatGetLocalSize(pc->mat,&m2,&n2);
    1189:     if (m1 != m2 || n1 != n2) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Cannot change local size of Amat after use old sizes %D %D new sizes %D %D",m2,n2,m1,n1);
    1190:     MatGetLocalSize(Pmat,&m1,&n1);
    1191:     MatGetLocalSize(pc->pmat,&m2,&n2);
    1192:     if (m1 != m2 || n1 != n2) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Cannot change local size of Pmat after use old sizes %D %D new sizes %D %D",m2,n2,m1,n1);
    1193:   }

    1195:   if (Pmat != pc->pmat) {
    1196:     /* changing the operator that defines the preconditioner thus reneed to clear current states so new preconditioner is built */
    1197:     pc->matnonzerostate = -1;
    1198:     pc->matstate        = -1;
    1199:   }

    1201:   /* reference first in case the matrices are the same */
    1202:   if (Amat) {PetscObjectReference((PetscObject)Amat);}
    1203:   MatDestroy(&pc->mat);
    1204:   if (Pmat) {PetscObjectReference((PetscObject)Pmat);}
    1205:   MatDestroy(&pc->pmat);
    1206:   pc->mat  = Amat;
    1207:   pc->pmat = Pmat;
    1208:   return(0);
    1209: }

    1211: /*@
    1212:    PCSetReusePreconditioner - reuse the current preconditioner even if the operator in the preconditioner has changed.

    1214:    Logically Collective on PC

    1216:    Input Parameters:
    1217: +  pc - the preconditioner context
    1218: -  flag - PETSC_TRUE do not compute a new preconditioner, PETSC_FALSE do compute a new preconditioner

    1220:     Level: intermediate

    1222: .seealso: PCGetOperators(), MatZeroEntries(), PCGetReusePreconditioner(), KSPSetReusePreconditioner()
    1223:  @*/
    1224: PetscErrorCode  PCSetReusePreconditioner(PC pc,PetscBool flag)
    1225: {
    1229:   pc->reusepreconditioner = flag;
    1230:   return(0);
    1231: }

    1233: /*@
    1234:    PCGetReusePreconditioner - Determines if the PC reuses the current preconditioner even if the operator in the preconditioner has changed.

    1236:    Not Collective

    1238:    Input Parameter:
    1239: .  pc - the preconditioner context

    1241:    Output Parameter:
    1242: .  flag - PETSC_TRUE do not compute a new preconditioner, PETSC_FALSE do compute a new preconditioner

    1244:    Level: intermediate

    1246: .seealso: PCGetOperators(), MatZeroEntries(), PCSetReusePreconditioner()
    1247:  @*/
    1248: PetscErrorCode  PCGetReusePreconditioner(PC pc,PetscBool *flag)
    1249: {
    1253:   *flag = pc->reusepreconditioner;
    1254:   return(0);
    1255: }

    1257: /*@
    1258:    PCGetOperators - Gets the matrix associated with the linear system and
    1259:    possibly a different one associated with the preconditioner.

    1261:    Not collective, though parallel Mats are returned if the PC is parallel

    1263:    Input Parameter:
    1264: .  pc - the preconditioner context

    1266:    Output Parameters:
    1267: +  Amat - the matrix defining the linear system
    1268: -  Pmat - the matrix from which the preconditioner is constructed, usually the same as Amat.

    1270:    Level: intermediate

    1272:    Notes:
    1273:     Does not increase the reference count of the matrices, so you should not destroy them

    1275:    Alternative usage: If the operators have NOT been set with KSP/PCSetOperators() then the operators
    1276:       are created in PC and returned to the user. In this case, if both operators
    1277:       mat and pmat are requested, two DIFFERENT operators will be returned. If
    1278:       only one is requested both operators in the PC will be the same (i.e. as
    1279:       if one had called KSP/PCSetOperators() with the same argument for both Mats).
    1280:       The user must set the sizes of the returned matrices and their type etc just
    1281:       as if the user created them with MatCreate(). For example,

    1283: $         KSP/PCGetOperators(ksp/pc,&Amat,NULL); is equivalent to
    1284: $           set size, type, etc of Amat

    1286: $         MatCreate(comm,&mat);
    1287: $         KSP/PCSetOperators(ksp/pc,Amat,Amat);
    1288: $         PetscObjectDereference((PetscObject)mat);
    1289: $           set size, type, etc of Amat

    1291:      and

    1293: $         KSP/PCGetOperators(ksp/pc,&Amat,&Pmat); is equivalent to
    1294: $           set size, type, etc of Amat and Pmat

    1296: $         MatCreate(comm,&Amat);
    1297: $         MatCreate(comm,&Pmat);
    1298: $         KSP/PCSetOperators(ksp/pc,Amat,Pmat);
    1299: $         PetscObjectDereference((PetscObject)Amat);
    1300: $         PetscObjectDereference((PetscObject)Pmat);
    1301: $           set size, type, etc of Amat and Pmat

    1303:     The rational for this support is so that when creating a TS, SNES, or KSP the hierarchy
    1304:     of underlying objects (i.e. SNES, KSP, PC, Mat) and their livespans can be completely
    1305:     managed by the top most level object (i.e. the TS, SNES, or KSP). Another way to look
    1306:     at this is when you create a SNES you do not NEED to create a KSP and attach it to
    1307:     the SNES object (the SNES object manages it for you). Similarly when you create a KSP
    1308:     you do not need to attach a PC to it (the KSP object manages the PC object for you).
    1309:     Thus, why should YOU have to create the Mat and attach it to the SNES/KSP/PC, when
    1310:     it can be created for you?

    1313: .seealso: PCSetOperators(), KSPGetOperators(), KSPSetOperators(), PCGetOperatorsSet()
    1314: @*/
    1315: PetscErrorCode  PCGetOperators(PC pc,Mat *Amat,Mat *Pmat)
    1316: {

    1321:   if (Amat) {
    1322:     if (!pc->mat) {
    1323:       if (pc->pmat && !Pmat) {  /* Apmat has been set, but user did not request it, so use for Amat */
    1324:         pc->mat = pc->pmat;
    1325:         PetscObjectReference((PetscObject)pc->mat);
    1326:       } else {                  /* both Amat and Pmat are empty */
    1327:         MatCreate(PetscObjectComm((PetscObject)pc),&pc->mat);
    1328:         if (!Pmat) { /* user did NOT request Pmat, so make same as Amat */
    1329:           pc->pmat = pc->mat;
    1330:           PetscObjectReference((PetscObject)pc->pmat);
    1331:         }
    1332:       }
    1333:     }
    1334:     *Amat = pc->mat;
    1335:   }
    1336:   if (Pmat) {
    1337:     if (!pc->pmat) {
    1338:       if (pc->mat && !Amat) {    /* Amat has been set but was not requested, so use for pmat */
    1339:         pc->pmat = pc->mat;
    1340:         PetscObjectReference((PetscObject)pc->pmat);
    1341:       } else {
    1342:         MatCreate(PetscObjectComm((PetscObject)pc),&pc->pmat);
    1343:         if (!Amat) { /* user did NOT request Amat, so make same as Pmat */
    1344:           pc->mat = pc->pmat;
    1345:           PetscObjectReference((PetscObject)pc->mat);
    1346:         }
    1347:       }
    1348:     }
    1349:     *Pmat = pc->pmat;
    1350:   }
    1351:   return(0);
    1352: }

    1354: /*@C
    1355:    PCGetOperatorsSet - Determines if the matrix associated with the linear system and
    1356:    possibly a different one associated with the preconditioner have been set in the PC.

    1358:    Not collective, though the results on all processes should be the same

    1360:    Input Parameter:
    1361: .  pc - the preconditioner context

    1363:    Output Parameters:
    1364: +  mat - the matrix associated with the linear system was set
    1365: -  pmat - matrix associated with the preconditioner was set, usually the same

    1367:    Level: intermediate

    1369: .seealso: PCSetOperators(), KSPGetOperators(), KSPSetOperators(), PCGetOperators()
    1370: @*/
    1371: PetscErrorCode  PCGetOperatorsSet(PC pc,PetscBool  *mat,PetscBool  *pmat)
    1372: {
    1375:   if (mat) *mat = (pc->mat)  ? PETSC_TRUE : PETSC_FALSE;
    1376:   if (pmat) *pmat = (pc->pmat) ? PETSC_TRUE : PETSC_FALSE;
    1377:   return(0);
    1378: }

    1380: /*@
    1381:    PCFactorGetMatrix - Gets the factored matrix from the
    1382:    preconditioner context.  This routine is valid only for the LU,
    1383:    incomplete LU, Cholesky, and incomplete Cholesky methods.

    1385:    Not Collective on PC though Mat is parallel if PC is parallel

    1387:    Input Parameters:
    1388: .  pc - the preconditioner context

    1390:    Output parameters:
    1391: .  mat - the factored matrix

    1393:    Level: advanced

    1395:    Notes:
    1396:     Does not increase the reference count for the matrix so DO NOT destroy it

    1398: @*/
    1399: PetscErrorCode  PCFactorGetMatrix(PC pc,Mat *mat)
    1400: {

    1406:   if (pc->ops->getfactoredmatrix) {
    1407:     (*pc->ops->getfactoredmatrix)(pc,mat);
    1408:   } else SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"PC type does not support getting factor matrix");
    1409:   return(0);
    1410: }

    1412: /*@C
    1413:    PCSetOptionsPrefix - Sets the prefix used for searching for all
    1414:    PC options in the database.

    1416:    Logically Collective on PC

    1418:    Input Parameters:
    1419: +  pc - the preconditioner context
    1420: -  prefix - the prefix string to prepend to all PC option requests

    1422:    Notes:
    1423:    A hyphen (-) must NOT be given at the beginning of the prefix name.
    1424:    The first character of all runtime options is AUTOMATICALLY the
    1425:    hyphen.

    1427:    Level: advanced

    1429: .seealso: PCAppendOptionsPrefix(), PCGetOptionsPrefix()
    1430: @*/
    1431: PetscErrorCode  PCSetOptionsPrefix(PC pc,const char prefix[])
    1432: {

    1437:   PetscObjectSetOptionsPrefix((PetscObject)pc,prefix);
    1438:   return(0);
    1439: }

    1441: /*@C
    1442:    PCAppendOptionsPrefix - Appends to the prefix used for searching for all
    1443:    PC options in the database.

    1445:    Logically Collective on PC

    1447:    Input Parameters:
    1448: +  pc - the preconditioner context
    1449: -  prefix - the prefix string to prepend to all PC option requests

    1451:    Notes:
    1452:    A hyphen (-) must NOT be given at the beginning of the prefix name.
    1453:    The first character of all runtime options is AUTOMATICALLY the
    1454:    hyphen.

    1456:    Level: advanced

    1458: .seealso: PCSetOptionsPrefix(), PCGetOptionsPrefix()
    1459: @*/
    1460: PetscErrorCode  PCAppendOptionsPrefix(PC pc,const char prefix[])
    1461: {

    1466:   PetscObjectAppendOptionsPrefix((PetscObject)pc,prefix);
    1467:   return(0);
    1468: }

    1470: /*@C
    1471:    PCGetOptionsPrefix - Gets the prefix used for searching for all
    1472:    PC options in the database.

    1474:    Not Collective

    1476:    Input Parameters:
    1477: .  pc - the preconditioner context

    1479:    Output Parameters:
    1480: .  prefix - pointer to the prefix string used, is returned

    1482:    Notes:
    1483:     On the fortran side, the user should pass in a string 'prifix' of
    1484:    sufficient length to hold the prefix.

    1486:    Level: advanced

    1488: .seealso: PCSetOptionsPrefix(), PCAppendOptionsPrefix()
    1489: @*/
    1490: PetscErrorCode  PCGetOptionsPrefix(PC pc,const char *prefix[])
    1491: {

    1497:   PetscObjectGetOptionsPrefix((PetscObject)pc,prefix);
    1498:   return(0);
    1499: }

    1501: /*
    1502:    Indicates the right hand side will be changed by KSPSolve(), this occurs for a few
    1503:   preconditioners including BDDC and Eisentat that transform the equations before applying
    1504:   the Krylov methods
    1505: */
    1506: PETSC_INTERN PetscErrorCode  PCPreSolveChangeRHS(PC pc,PetscBool *change)
    1507: {

    1513:   *change = PETSC_FALSE;
    1514:   PetscTryMethod(pc,"PCPreSolveChangeRHS_C",(PC,PetscBool*),(pc,change));
    1515:   return(0);
    1516: }

    1518: /*@
    1519:    PCPreSolve - Optional pre-solve phase, intended for any
    1520:    preconditioner-specific actions that must be performed before
    1521:    the iterative solve itself.

    1523:    Collective on PC

    1525:    Input Parameters:
    1526: +  pc - the preconditioner context
    1527: -  ksp - the Krylov subspace context

    1529:    Level: developer

    1531:    Sample of Usage:
    1532: .vb
    1533:     PCPreSolve(pc,ksp);
    1534:     KSPSolve(ksp,b,x);
    1535:     PCPostSolve(pc,ksp);
    1536: .ve

    1538:    Notes:
    1539:    The pre-solve phase is distinct from the PCSetUp() phase.

    1541:    KSPSolve() calls this directly, so is rarely called by the user.

    1543: .seealso: PCPostSolve()
    1544: @*/
    1545: PetscErrorCode  PCPreSolve(PC pc,KSP ksp)
    1546: {
    1548:   Vec            x,rhs;

    1553:   pc->presolvedone++;
    1554:   if (pc->presolvedone > 2) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Cannot embed PCPreSolve() more than twice");
    1555:   KSPGetSolution(ksp,&x);
    1556:   KSPGetRhs(ksp,&rhs);

    1558:   if (pc->ops->presolve) {
    1559:     (*pc->ops->presolve)(pc,ksp,rhs,x);
    1560:   }
    1561:   return(0);
    1562: }

    1564: /*@
    1565:    PCPostSolve - Optional post-solve phase, intended for any
    1566:    preconditioner-specific actions that must be performed after
    1567:    the iterative solve itself.

    1569:    Collective on PC

    1571:    Input Parameters:
    1572: +  pc - the preconditioner context
    1573: -  ksp - the Krylov subspace context

    1575:    Sample of Usage:
    1576: .vb
    1577:     PCPreSolve(pc,ksp);
    1578:     KSPSolve(ksp,b,x);
    1579:     PCPostSolve(pc,ksp);
    1580: .ve

    1582:    Note:
    1583:    KSPSolve() calls this routine directly, so it is rarely called by the user.

    1585:    Level: developer

    1587: .seealso: PCPreSolve(), KSPSolve()
    1588: @*/
    1589: PetscErrorCode  PCPostSolve(PC pc,KSP ksp)
    1590: {
    1592:   Vec            x,rhs;

    1597:   pc->presolvedone--;
    1598:   KSPGetSolution(ksp,&x);
    1599:   KSPGetRhs(ksp,&rhs);
    1600:   if (pc->ops->postsolve) {
    1601:      (*pc->ops->postsolve)(pc,ksp,rhs,x);
    1602:   }
    1603:   return(0);
    1604: }

    1606: /*@C
    1607:   PCLoad - Loads a PC that has been stored in binary  with PCView().

    1609:   Collective on PetscViewer

    1611:   Input Parameters:
    1612: + newdm - the newly loaded PC, this needs to have been created with PCCreate() or
    1613:            some related function before a call to PCLoad().
    1614: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen()

    1616:    Level: intermediate

    1618:   Notes:
    1619:    The type is determined by the data in the file, any type set into the PC before this call is ignored.

    1621:   Notes for advanced users:
    1622:   Most users should not need to know the details of the binary storage
    1623:   format, since PCLoad() and PCView() completely hide these details.
    1624:   But for anyone who's interested, the standard binary matrix storage
    1625:   format is
    1626: .vb
    1627:      has not yet been determined
    1628: .ve

    1630: .seealso: PetscViewerBinaryOpen(), PCView(), MatLoad(), VecLoad()
    1631: @*/
    1632: PetscErrorCode  PCLoad(PC newdm, PetscViewer viewer)
    1633: {
    1635:   PetscBool      isbinary;
    1636:   PetscInt       classid;
    1637:   char           type[256];

    1642:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
    1643:   if (!isbinary) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Invalid viewer; open viewer with PetscViewerBinaryOpen()");

    1645:   PetscViewerBinaryRead(viewer,&classid,1,NULL,PETSC_INT);
    1646:   if (classid != PC_FILE_CLASSID) SETERRQ(PetscObjectComm((PetscObject)newdm),PETSC_ERR_ARG_WRONG,"Not PC next in file");
    1647:   PetscViewerBinaryRead(viewer,type,256,NULL,PETSC_CHAR);
    1648:   PCSetType(newdm, type);
    1649:   if (newdm->ops->load) {
    1650:     (*newdm->ops->load)(newdm,viewer);
    1651:   }
    1652:   return(0);
    1653: }

    1655: #include <petscdraw.h>
    1656: #if defined(PETSC_HAVE_SAWS)
    1657: #include <petscviewersaws.h>
    1658: #endif

    1660: /*@C
    1661:    PCViewFromOptions - View from Options

    1663:    Collective on PC

    1665:    Input Parameters:
    1666: +  A - the PC context
    1667: .  obj - Optional object
    1668: -  name - command line option

    1670:    Level: intermediate
    1671: .seealso:  PC, PCView, PetscObjectViewFromOptions(), PCCreate()
    1672: @*/
    1673: PetscErrorCode  PCViewFromOptions(PC A,PetscObject obj,const char name[])
    1674: {

    1679:   PetscObjectViewFromOptions((PetscObject)A,obj,name);
    1680:   return(0);
    1681: }

    1683: /*@C
    1684:    PCView - Prints the PC data structure.

    1686:    Collective on PC

    1688:    Input Parameters:
    1689: +  PC - the PC context
    1690: -  viewer - optional visualization context

    1692:    Note:
    1693:    The available visualization contexts include
    1694: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
    1695: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
    1696:          output where only the first processor opens
    1697:          the file.  All other processors send their
    1698:          data to the first processor to print.

    1700:    The user can open an alternative visualization contexts with
    1701:    PetscViewerASCIIOpen() (output to a specified file).

    1703:    Level: developer

    1705: .seealso: KSPView(), PetscViewerASCIIOpen()
    1706: @*/
    1707: PetscErrorCode  PCView(PC pc,PetscViewer viewer)
    1708: {
    1709:   PCType         cstr;
    1711:   PetscBool      iascii,isstring,isbinary,isdraw;
    1712: #if defined(PETSC_HAVE_SAWS)
    1713:   PetscBool      issaws;
    1714: #endif

    1718:   if (!viewer) {
    1719:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pc),&viewer);
    1720:   }

    1724:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
    1725:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
    1726:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
    1727:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
    1728: #if defined(PETSC_HAVE_SAWS)
    1729:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSAWS,&issaws);
    1730: #endif

    1732:   if (iascii) {
    1733:     PetscObjectPrintClassNamePrefixType((PetscObject)pc,viewer);
    1734:     if (!pc->setupcalled) {
    1735:       PetscViewerASCIIPrintf(viewer,"  PC has not been set up so information may be incomplete\n");
    1736:     }
    1737:     if (pc->ops->view) {
    1738:       PetscViewerASCIIPushTab(viewer);
    1739:       (*pc->ops->view)(pc,viewer);
    1740:       PetscViewerASCIIPopTab(viewer);
    1741:     }
    1742:     if (pc->mat) {
    1743:       PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
    1744:       if (pc->pmat == pc->mat) {
    1745:         PetscViewerASCIIPrintf(viewer,"  linear system matrix = precond matrix:\n");
    1746:         PetscViewerASCIIPushTab(viewer);
    1747:         MatView(pc->mat,viewer);
    1748:         PetscViewerASCIIPopTab(viewer);
    1749:       } else {
    1750:         if (pc->pmat) {
    1751:           PetscViewerASCIIPrintf(viewer,"  linear system matrix followed by preconditioner matrix:\n");
    1752:         } else {
    1753:           PetscViewerASCIIPrintf(viewer,"  linear system matrix:\n");
    1754:         }
    1755:         PetscViewerASCIIPushTab(viewer);
    1756:         MatView(pc->mat,viewer);
    1757:         if (pc->pmat) {MatView(pc->pmat,viewer);}
    1758:         PetscViewerASCIIPopTab(viewer);
    1759:       }
    1760:       PetscViewerPopFormat(viewer);
    1761:     }
    1762:   } else if (isstring) {
    1763:     PCGetType(pc,&cstr);
    1764:     PetscViewerStringSPrintf(viewer," PCType: %-7.7s",cstr);
    1765:     if (pc->ops->view) {(*pc->ops->view)(pc,viewer);}
    1766:     if (pc->mat) {MatView(pc->mat,viewer);}
    1767:     if (pc->pmat && pc->pmat != pc->mat) {MatView(pc->pmat,viewer);}
    1768:   } else if (isbinary) {
    1769:     PetscInt    classid = PC_FILE_CLASSID;
    1770:     MPI_Comm    comm;
    1771:     PetscMPIInt rank;
    1772:     char        type[256];

    1774:     PetscObjectGetComm((PetscObject)pc,&comm);
    1775:     MPI_Comm_rank(comm,&rank);
    1776:     if (!rank) {
    1777:       PetscViewerBinaryWrite(viewer,&classid,1,PETSC_INT);
    1778:       PetscStrncpy(type,((PetscObject)pc)->type_name,256);
    1779:       PetscViewerBinaryWrite(viewer,type,256,PETSC_CHAR);
    1780:     }
    1781:     if (pc->ops->view) {
    1782:       (*pc->ops->view)(pc,viewer);
    1783:     }
    1784:   } else if (isdraw) {
    1785:     PetscDraw draw;
    1786:     char      str[25];
    1787:     PetscReal x,y,bottom,h;
    1788:     PetscInt  n;

    1790:     PetscViewerDrawGetDraw(viewer,0,&draw);
    1791:     PetscDrawGetCurrentPoint(draw,&x,&y);
    1792:     if (pc->mat) {
    1793:       MatGetSize(pc->mat,&n,NULL);
    1794:       PetscSNPrintf(str,25,"PC: %s (%D)",((PetscObject)pc)->type_name,n);
    1795:     } else {
    1796:       PetscSNPrintf(str,25,"PC: %s",((PetscObject)pc)->type_name);
    1797:     }
    1798:     PetscDrawStringBoxed(draw,x,y,PETSC_DRAW_RED,PETSC_DRAW_BLACK,str,NULL,&h);
    1799:     bottom = y - h;
    1800:     PetscDrawPushCurrentPoint(draw,x,bottom);
    1801:     if (pc->ops->view) {
    1802:       (*pc->ops->view)(pc,viewer);
    1803:     }
    1804:     PetscDrawPopCurrentPoint(draw);
    1805: #if defined(PETSC_HAVE_SAWS)
    1806:   } else if (issaws) {
    1807:     PetscMPIInt rank;

    1809:     PetscObjectName((PetscObject)pc);
    1810:     MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
    1811:     if (!((PetscObject)pc)->amsmem && !rank) {
    1812:       PetscObjectViewSAWs((PetscObject)pc,viewer);
    1813:     }
    1814:     if (pc->mat) {MatView(pc->mat,viewer);}
    1815:     if (pc->pmat && pc->pmat != pc->mat) {MatView(pc->pmat,viewer);}
    1816: #endif
    1817:   }
    1818:   return(0);
    1819: }

    1821: /*@C
    1822:   PCRegister -  Adds a method to the preconditioner package.

    1824:    Not collective

    1826:    Input Parameters:
    1827: +  name_solver - name of a new user-defined solver
    1828: -  routine_create - routine to create method context

    1830:    Notes:
    1831:    PCRegister() may be called multiple times to add several user-defined preconditioners.

    1833:    Sample usage:
    1834: .vb
    1835:    PCRegister("my_solver", MySolverCreate);
    1836: .ve

    1838:    Then, your solver can be chosen with the procedural interface via
    1839: $     PCSetType(pc,"my_solver")
    1840:    or at runtime via the option
    1841: $     -pc_type my_solver

    1843:    Level: advanced

    1845: .seealso: PCRegisterAll()
    1846: @*/
    1847: PetscErrorCode  PCRegister(const char sname[],PetscErrorCode (*function)(PC))
    1848: {

    1852:   PCInitializePackage();
    1853:   PetscFunctionListAdd(&PCList,sname,function);
    1854:   return(0);
    1855: }

    1857: static PetscErrorCode MatMult_PC(Mat A,Vec X,Vec Y)
    1858: {
    1859:   PC             pc;

    1863:   MatShellGetContext(A,&pc);
    1864:   PCApply(pc,X,Y);
    1865:   return(0);
    1866: }

    1868: /*@
    1869:     PCComputeOperator - Computes the explicit preconditioned operator.

    1871:     Collective on PC

    1873:     Input Parameter:
    1874: +   pc - the preconditioner object
    1875: -   mattype - the matrix type to be used for the operator

    1877:     Output Parameter:
    1878: .   mat - the explict preconditioned operator

    1880:     Notes:
    1881:     This computation is done by applying the operators to columns of the identity matrix.
    1882:     This routine is costly in general, and is recommended for use only with relatively small systems.
    1883:     Currently, this routine uses a dense matrix format when mattype == NULL

    1885:     Level: advanced

    1887: .seealso: KSPComputeOperator(), MatType

    1889: @*/
    1890: PetscErrorCode  PCComputeOperator(PC pc,MatType mattype,Mat *mat)
    1891: {
    1893:   PetscInt       N,M,m,n;
    1894:   Mat            A,Apc;

    1899:   PCGetOperators(pc,&A,NULL);
    1900:   MatGetLocalSize(A,&m,&n);
    1901:   MatGetSize(A,&M,&N);
    1902:   MatCreateShell(PetscObjectComm((PetscObject)pc),m,n,M,N,pc,&Apc);
    1903:   MatShellSetOperation(Apc,MATOP_MULT,(void (*)(void))MatMult_PC);
    1904:   MatComputeOperator(Apc,mattype,mat);
    1905:   MatDestroy(&Apc);
    1906:   return(0);
    1907: }

    1909: /*@
    1910:    PCSetCoordinates - sets the coordinates of all the nodes on the local process

    1912:    Collective on PC

    1914:    Input Parameters:
    1915: +  pc - the solver context
    1916: .  dim - the dimension of the coordinates 1, 2, or 3
    1917: .  nloc - the blocked size of the coordinates array
    1918: -  coords - the coordinates array

    1920:    Level: intermediate

    1922:    Notes:
    1923:    coords is an array of the dim coordinates for the nodes on
    1924:    the local processor, of size dim*nloc.
    1925:    If there are 108 equation on a processor
    1926:    for a displacement finite element discretization of elasticity (so
    1927:    that there are nloc = 36 = 108/3 nodes) then the array must have 108
    1928:    double precision values (ie, 3 * 36).  These x y z coordinates
    1929:    should be ordered for nodes 0 to N-1 like so: [ 0.x, 0.y, 0.z, 1.x,
    1930:    ... , N-1.z ].

    1932: .seealso: MatSetNearNullSpace()
    1933: @*/
    1934: PetscErrorCode PCSetCoordinates(PC pc, PetscInt dim, PetscInt nloc, PetscReal coords[])
    1935: {

    1941:   PetscTryMethod(pc,"PCSetCoordinates_C",(PC,PetscInt,PetscInt,PetscReal*),(pc,dim,nloc,coords));
    1942:   return(0);
    1943: }

    1945: /*@
    1946:    PCGetInterpolations - Gets interpolation matrices for all levels (except level 0)

    1948:    Logically Collective on PC

    1950:    Input Parameters:
    1951: +  pc - the precondition context

    1953:    Output Parameter:
    1954: -  num_levels - the number of levels
    1955: .  interpolations - the interpolation matrices (size of num_levels-1)

    1957:    Level: advanced

    1959: .keywords: MG, GAMG, BoomerAMG, multigrid, interpolation, level

    1961: .seealso: PCMGGetRestriction(), PCMGSetInterpolation(), PCMGGetInterpolation(), PCGetCoarseOperators()
    1962: @*/
    1963: PetscErrorCode PCGetInterpolations(PC pc,PetscInt *num_levels,Mat *interpolations[])
    1964: {

    1971:   PetscUseMethod(pc,"PCGetInterpolations_C",(PC,PetscInt*,Mat*[]),(pc,num_levels,interpolations));
    1972:   return(0);
    1973: }

    1975: /*@
    1976:    PCGetCoarseOperators - Gets coarse operator matrices for all levels (except the finest level)

    1978:    Logically Collective on PC

    1980:    Input Parameters:
    1981: +  pc - the precondition context

    1983:    Output Parameter:
    1984: -  num_levels - the number of levels
    1985: .  coarseOperators - the coarse operator matrices (size of num_levels-1)

    1987:    Level: advanced

    1989: .keywords: MG, GAMG, BoomerAMG, get, multigrid, interpolation, level

    1991: .seealso: PCMGGetRestriction(), PCMGSetInterpolation(), PCMGGetRScale(), PCMGGetInterpolation(), PCGetInterpolations()
    1992: @*/
    1993: PetscErrorCode PCGetCoarseOperators(PC pc,PetscInt *num_levels,Mat *coarseOperators[])
    1994: {

    2001:   PetscUseMethod(pc,"PCGetCoarseOperators_C",(PC,PetscInt*,Mat*[]),(pc,num_levels,coarseOperators));
    2002:   return(0);
    2003: }