1 2 #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 3 4 static PetscErrorCode PCApply_Mat(PC pc,Vec x,Vec y) 5 { 6 PetscFunctionBegin; 7 PetscCall(MatMult(pc->pmat,x,y)); 8 PetscFunctionReturn(0); 9 } 10 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(0); 16 } 17 18 static PetscErrorCode PCApplyTranspose_Mat(PC pc,Vec x,Vec y) 19 { 20 PetscFunctionBegin; 21 PetscCall(MatMultTranspose(pc->pmat,x,y)); 22 PetscFunctionReturn(0); 23 } 24 25 static PetscErrorCode PCDestroy_Mat(PC pc) 26 { 27 PetscFunctionBegin; 28 PetscFunctionReturn(0); 29 } 30 31 /*MC 32 PCMAT - A preconditioner obtained by multiplying by the preconditioner matrix supplied 33 in PCSetOperators() or KSPSetOperators() 34 35 Notes: 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. 38 39 Level: intermediate 40 41 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 42 PCSHELL 43 44 M*/ 45 46 PETSC_EXTERN PetscErrorCode PCCreate_Mat(PC pc) 47 { 48 PetscFunctionBegin; 49 pc->ops->apply = PCApply_Mat; 50 pc->ops->matapply = PCMatApply_Mat; 51 pc->ops->applytranspose = PCApplyTranspose_Mat; 52 pc->ops->setup = NULL; 53 pc->ops->destroy = PCDestroy_Mat; 54 pc->ops->setfromoptions = NULL; 55 pc->ops->view = NULL; 56 pc->ops->applyrichardson = NULL; 57 pc->ops->applysymmetricleft = NULL; 58 pc->ops->applysymmetricright = NULL; 59 PetscFunctionReturn(0); 60 } 61