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