1 2 /* -------------------------------------------------------------------- */ 3 4 /* 5 Include files needed for the ViennaCL Smoothed Aggregation preconditioner: 6 pcimpl.h - private include file intended for use by all preconditioners 7 */ 8 #define PETSC_SKIP_SPINLOCK 9 #include <petsc/private/pcimpl.h> /*I "petscpc.h" I*/ 10 #include <../src/mat/impls/aij/seq/aij.h> 11 #include <../src/vec/vec/impls/dvecimpl.h> 12 #include <../src/mat/impls/aij/seq/seqviennacl/viennaclmatimpl.h> 13 #include <../src/vec/vec/impls/seq/seqviennacl/viennaclvecimpl.h> 14 #include <viennacl/linalg/amg.hpp> 15 16 /* 17 Private context (data structure) for the SAVIENNACL preconditioner. 18 */ 19 typedef struct { 20 viennacl::linalg::amg_precond<viennacl::compressed_matrix<PetscScalar> > *SAVIENNACL; 21 } PC_SAVIENNACL; 22 23 24 /* -------------------------------------------------------------------------- */ 25 /* 26 PCSetUp_SAVIENNACL - Prepares for the use of the SAVIENNACL preconditioner 27 by setting data structures and options. 28 29 Input Parameter: 30 . pc - the preconditioner context 31 32 Application Interface Routine: PCSetUp() 33 34 Notes: 35 The interface routine PCSetUp() is not usually called directly by 36 the user, but instead is called by PCApply() if necessary. 37 */ 38 static PetscErrorCode PCSetUp_SAVIENNACL(PC pc) 39 { 40 PC_SAVIENNACL *sa = (PC_SAVIENNACL*)pc->data; 41 PetscBool flg = PETSC_FALSE; 42 PetscErrorCode ierr; 43 Mat_SeqAIJViennaCL *gpustruct; 44 45 PetscFunctionBegin; 46 ierr = PetscObjectTypeCompare((PetscObject)pc->pmat,MATSEQAIJVIENNACL,&flg);CHKERRQ(ierr); 47 if (!flg) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"Currently only handles ViennaCL matrices"); 48 if (pc->setupcalled != 0) { 49 try { 50 delete sa->SAVIENNACL; 51 } catch(char *ex) { 52 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex); 53 } 54 } 55 try { 56 #if defined(PETSC_USE_COMPLEX) 57 gpustruct = NULL; 58 SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"No support for complex arithmetic in SAVIENNACL preconditioner"); 59 #else 60 ierr = MatViennaCLCopyToGPU(pc->pmat);CHKERRQ(ierr); 61 gpustruct = (Mat_SeqAIJViennaCL*)(pc->pmat->spptr); 62 63 viennacl::linalg::amg_tag amg_tag_sa_pmis; 64 amg_tag_sa_pmis.set_coarsening_method(viennacl::linalg::AMG_COARSENING_METHOD_MIS2_AGGREGATION); 65 amg_tag_sa_pmis.set_interpolation_method(viennacl::linalg::AMG_INTERPOLATION_METHOD_SMOOTHED_AGGREGATION); 66 ViennaCLAIJMatrix *mat = (ViennaCLAIJMatrix*)gpustruct->mat; 67 sa->SAVIENNACL = new viennacl::linalg::amg_precond<viennacl::compressed_matrix<PetscScalar> >(*mat, amg_tag_sa_pmis); 68 sa->SAVIENNACL->setup(); 69 #endif 70 } catch(char *ex) { 71 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex); 72 } 73 PetscFunctionReturn(0); 74 } 75 76 /* -------------------------------------------------------------------------- */ 77 /* 78 PCApply_SAVIENNACL - Applies the SAVIENNACL preconditioner to a vector. 79 80 Input Parameters: 81 . pc - the preconditioner context 82 . x - input vector 83 84 Output Parameter: 85 . y - output vector 86 87 Application Interface Routine: PCApply() 88 */ 89 static PetscErrorCode PCApply_SAVIENNACL(PC pc,Vec x,Vec y) 90 { 91 PC_SAVIENNACL *sac = (PC_SAVIENNACL*)pc->data; 92 PetscErrorCode ierr; 93 PetscBool flg1,flg2; 94 viennacl::vector<PetscScalar> const *xarray=NULL; 95 viennacl::vector<PetscScalar> *yarray=NULL; 96 97 PetscFunctionBegin; 98 /*how to apply a certain fixed number of iterations?*/ 99 ierr = PetscObjectTypeCompare((PetscObject)x,VECSEQVIENNACL,&flg1);CHKERRQ(ierr); 100 ierr = PetscObjectTypeCompare((PetscObject)y,VECSEQVIENNACL,&flg2);CHKERRQ(ierr); 101 if (!(flg1 && flg2)) SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP, "Currently only handles ViennaCL vectors"); 102 if (!sac->SAVIENNACL) { 103 ierr = PCSetUp_SAVIENNACL(pc);CHKERRQ(ierr); 104 } 105 ierr = VecSet(y,0.0);CHKERRQ(ierr); 106 ierr = VecViennaCLGetArrayRead(x,&xarray);CHKERRQ(ierr); 107 ierr = VecViennaCLGetArrayWrite(y,&yarray);CHKERRQ(ierr); 108 try { 109 #if defined(PETSC_USE_COMPLEX) 110 111 #else 112 *yarray = *xarray; 113 sac->SAVIENNACL->apply(*yarray); 114 #endif 115 } catch(char * ex) { 116 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex); 117 } 118 ierr = VecViennaCLRestoreArrayRead(x,&xarray);CHKERRQ(ierr); 119 ierr = VecViennaCLRestoreArrayWrite(y,&yarray);CHKERRQ(ierr); 120 ierr = PetscObjectStateIncrease((PetscObject)y);CHKERRQ(ierr); 121 PetscFunctionReturn(0); 122 } 123 /* -------------------------------------------------------------------------- */ 124 /* 125 PCDestroy_SAVIENNACL - Destroys the private context for the SAVIENNACL preconditioner 126 that was created with PCCreate_SAVIENNACL(). 127 128 Input Parameter: 129 . pc - the preconditioner context 130 131 Application Interface Routine: PCDestroy() 132 */ 133 static PetscErrorCode PCDestroy_SAVIENNACL(PC pc) 134 { 135 PC_SAVIENNACL *sac = (PC_SAVIENNACL*)pc->data; 136 PetscErrorCode ierr; 137 138 PetscFunctionBegin; 139 if (sac->SAVIENNACL) { 140 try { 141 delete sac->SAVIENNACL; 142 } catch(char * ex) { 143 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex); 144 } 145 } 146 147 /* 148 Free the private data structure that was hanging off the PC 149 */ 150 ierr = PetscFree(pc->data);CHKERRQ(ierr); 151 PetscFunctionReturn(0); 152 } 153 154 static PetscErrorCode PCSetFromOptions_SAVIENNACL(PetscOptionItems *PetscOptionsObject,PC pc) 155 { 156 PetscErrorCode ierr; 157 158 PetscFunctionBegin; 159 ierr = PetscOptionsHead(PetscOptionsObject,"SAVIENNACL options");CHKERRQ(ierr); 160 ierr = PetscOptionsTail();CHKERRQ(ierr); 161 PetscFunctionReturn(0); 162 } 163 164 /* -------------------------------------------------------------------------- */ 165 166 167 /*MC 168 PCSAViennaCL - A smoothed agglomeration algorithm that can be used via the CUDA, OpenCL, and OpenMP backends of ViennaCL 169 170 Level: advanced 171 172 .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC 173 174 M*/ 175 176 PETSC_EXTERN PetscErrorCode PCCreate_SAVIENNACL(PC pc) 177 { 178 PC_SAVIENNACL *sac; 179 PetscErrorCode ierr; 180 181 PetscFunctionBegin; 182 /* 183 Creates the private data structure for this preconditioner and 184 attach it to the PC object. 185 */ 186 ierr = PetscNewLog(pc,&sac);CHKERRQ(ierr); 187 pc->data = (void*)sac; 188 189 /* 190 Initialize the pointer to zero 191 Initialize number of v-cycles to default (1) 192 */ 193 sac->SAVIENNACL = 0; 194 195 /* 196 Set the pointers for the functions that are provided above. 197 Now when the user-level routines (such as PCApply(), PCDestroy(), etc.) 198 are called, they will automatically call these functions. Note we 199 choose not to provide a couple of these functions since they are 200 not needed. 201 */ 202 pc->ops->apply = PCApply_SAVIENNACL; 203 pc->ops->applytranspose = 0; 204 pc->ops->setup = PCSetUp_SAVIENNACL; 205 pc->ops->destroy = PCDestroy_SAVIENNACL; 206 pc->ops->setfromoptions = PCSetFromOptions_SAVIENNACL; 207 pc->ops->view = 0; 208 pc->ops->applyrichardson = 0; 209 pc->ops->applysymmetricleft = 0; 210 pc->ops->applysymmetricright = 0; 211 PetscFunctionReturn(0); 212 } 213 214