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 PetscErrorCode ierr; 7 8 PetscFunctionBegin; 9 ierr = MatMult(pc->pmat,x,y);CHKERRQ(ierr); 10 PetscFunctionReturn(0); 11 } 12 13 static PetscErrorCode PCMatApply_Mat(PC pc,Mat X,Mat Y) 14 { 15 PetscErrorCode ierr; 16 17 PetscFunctionBegin; 18 ierr = MatMatMult(pc->pmat,X,MAT_REUSE_MATRIX,PETSC_DEFAULT,&Y);CHKERRQ(ierr); 19 PetscFunctionReturn(0); 20 } 21 22 static PetscErrorCode PCApplyTranspose_Mat(PC pc,Vec x,Vec y) 23 { 24 PetscErrorCode ierr; 25 26 PetscFunctionBegin; 27 ierr = MatMultTranspose(pc->pmat,x,y);CHKERRQ(ierr); 28 PetscFunctionReturn(0); 29 } 30 31 static PetscErrorCode PCDestroy_Mat(PC pc) 32 { 33 PetscFunctionBegin; 34 PetscFunctionReturn(0); 35 } 36 37 /*MC 38 PCMAT - A preconditioner obtained by multiplying by the preconditioner matrix supplied 39 in PCSetOperators() or KSPSetOperators() 40 41 Notes: 42 This one is a little strange. One rarely has an explicit matrix that approximates the 43 inverse of the matrix they wish to solve for. 44 45 Level: intermediate 46 47 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC, 48 PCSHELL 49 50 M*/ 51 52 PETSC_EXTERN PetscErrorCode PCCreate_Mat(PC pc) 53 { 54 PetscFunctionBegin; 55 pc->ops->apply = PCApply_Mat; 56 pc->ops->matapply = PCMatApply_Mat; 57 pc->ops->applytranspose = PCApplyTranspose_Mat; 58 pc->ops->setup = NULL; 59 pc->ops->destroy = PCDestroy_Mat; 60 pc->ops->setfromoptions = NULL; 61 pc->ops->view = NULL; 62 pc->ops->applyrichardson = NULL; 63 pc->ops->applysymmetricleft = NULL; 64 pc->ops->applysymmetricright = NULL; 65 PetscFunctionReturn(0); 66 } 67 68