Actual source code: pcmat.c


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

  4: static PetscErrorCode PCApply_Mat(PC pc, Vec x, Vec y)
  5: {
  6:   PetscFunctionBegin;
  7:   PetscCall(MatMult(pc->pmat, x, y));
  8:   PetscFunctionReturn(PETSC_SUCCESS);
  9: }

 11: static PetscErrorCode PCMatApply_Mat(PC pc, Mat X, Mat Y)
 12: {
 13:   PetscFunctionBegin;
 14:   PetscCall(MatMatMult(pc->pmat, X, MAT_REUSE_MATRIX, PETSC_DEFAULT, &Y));
 15:   PetscFunctionReturn(PETSC_SUCCESS);
 16: }

 18: static PetscErrorCode PCApplyTranspose_Mat(PC pc, Vec x, Vec y)
 19: {
 20:   PetscFunctionBegin;
 21:   PetscCall(MatMultTranspose(pc->pmat, x, y));
 22:   PetscFunctionReturn(PETSC_SUCCESS);
 23: }

 25: static PetscErrorCode PCDestroy_Mat(PC pc)
 26: {
 27:   PetscFunctionBegin;
 28:   PetscFunctionReturn(PETSC_SUCCESS);
 29: }

 31: /*MC
 32:      PCMAT - A preconditioner obtained by multiplying by the preconditioner matrix supplied
 33:              in `PCSetOperators()` or `KSPSetOperators()`

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

 39:    Level: intermediate

 41: .seealso: `PCCreate()`, `PCSetType()`, `PCType`, `PC`,
 42:           `PCSHELL`
 43: M*/

 45: PETSC_EXTERN PetscErrorCode PCCreate_Mat(PC pc)
 46: {
 47:   PetscFunctionBegin;
 48:   pc->ops->apply               = PCApply_Mat;
 49:   pc->ops->matapply            = PCMatApply_Mat;
 50:   pc->ops->applytranspose      = PCApplyTranspose_Mat;
 51:   pc->ops->setup               = NULL;
 52:   pc->ops->destroy             = PCDestroy_Mat;
 53:   pc->ops->setfromoptions      = NULL;
 54:   pc->ops->view                = NULL;
 55:   pc->ops->applyrichardson     = NULL;
 56:   pc->ops->applysymmetricleft  = NULL;
 57:   pc->ops->applysymmetricright = NULL;
 58:   PetscFunctionReturn(PETSC_SUCCESS);
 59: }