Actual source code: pcmat.c


  2: #include <petsc/private/pcimpl.h>

  4: static PetscErrorCode PCApply_Mat(PC pc, Vec x, Vec y)
  5: {
  6:   MatMult(pc->pmat, x, y);
  7:   return 0;
  8: }

 10: static PetscErrorCode PCMatApply_Mat(PC pc, Mat X, Mat Y)
 11: {
 12:   MatMatMult(pc->pmat, X, MAT_REUSE_MATRIX, PETSC_DEFAULT, &Y);
 13:   return 0;
 14: }

 16: static PetscErrorCode PCApplyTranspose_Mat(PC pc, Vec x, Vec y)
 17: {
 18:   MatMultTranspose(pc->pmat, x, y);
 19:   return 0;
 20: }

 22: static PetscErrorCode PCDestroy_Mat(PC pc)
 23: {
 24:   return 0;
 25: }

 27: /*MC
 28:      PCMAT - A preconditioner obtained by multiplying by the preconditioner matrix supplied
 29:              in `PCSetOperators()` or `KSPSetOperators()`

 31:    Note:
 32:     This one is a little strange. One rarely has an explicit matrix that approximates the
 33:          inverse of the matrix they wish to solve for.

 35:    Level: intermediate

 37: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
 38:           `PCSHELL`
 39: M*/

 41: PETSC_EXTERN PetscErrorCode PCCreate_Mat(PC pc)
 42: {
 43:   pc->ops->apply               = PCApply_Mat;
 44:   pc->ops->matapply            = PCMatApply_Mat;
 45:   pc->ops->applytranspose      = PCApplyTranspose_Mat;
 46:   pc->ops->setup               = NULL;
 47:   pc->ops->destroy             = PCDestroy_Mat;
 48:   pc->ops->setfromoptions      = NULL;
 49:   pc->ops->view                = NULL;
 50:   pc->ops->applyrichardson     = NULL;
 51:   pc->ops->applysymmetricleft  = NULL;
 52:   pc->ops->applysymmetricright = NULL;
 53:   return 0;
 54: }