xref: /petsc/src/ksp/pc/impls/saviennacl/saviennacl.cxx (revision ccb4e88a40f0b86eaeca07ff64c64e4de2fae686)
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 #define PETSC_SKIP_IMMINTRIN_H_CUDAWORKAROUND 1
10 #include <petsc/private/pcimpl.h>   /*I "petscpc.h" I*/
11 #include <../src/mat/impls/aij/seq/aij.h>
12 #include <../src/vec/vec/impls/dvecimpl.h>
13 #include <../src/mat/impls/aij/seq/seqviennacl/viennaclmatimpl.h>
14 #include <../src/vec/vec/impls/seq/seqviennacl/viennaclvecimpl.h>
15 #include <viennacl/linalg/amg.hpp>
16 
17 /*
18    Private context (data structure) for the SAVIENNACL preconditioner.
19 */
20 typedef struct {
21   viennacl::linalg::amg_precond<viennacl::compressed_matrix<PetscScalar> > *SAVIENNACL;
22 } PC_SAVIENNACL;
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 = VecViennaCLGetArrayRead(x,&xarray);CHKERRQ(ierr);
106   ierr = VecViennaCLGetArrayWrite(y,&yarray);CHKERRQ(ierr);
107   try {
108 #if !defined(PETSC_USE_COMPLEX)
109     *yarray = *xarray;
110     sac->SAVIENNACL->apply(*yarray);
111 #endif
112   } catch(char * ex) {
113     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex);
114   }
115   ierr = VecViennaCLRestoreArrayRead(x,&xarray);CHKERRQ(ierr);
116   ierr = VecViennaCLRestoreArrayWrite(y,&yarray);CHKERRQ(ierr);
117   ierr = PetscObjectStateIncrease((PetscObject)y);CHKERRQ(ierr);
118   PetscFunctionReturn(0);
119 }
120 /* -------------------------------------------------------------------------- */
121 /*
122    PCDestroy_SAVIENNACL - Destroys the private context for the SAVIENNACL preconditioner
123    that was created with PCCreate_SAVIENNACL().
124 
125    Input Parameter:
126 .  pc - the preconditioner context
127 
128    Application Interface Routine: PCDestroy()
129 */
130 static PetscErrorCode PCDestroy_SAVIENNACL(PC pc)
131 {
132   PC_SAVIENNACL  *sac = (PC_SAVIENNACL*)pc->data;
133   PetscErrorCode ierr;
134 
135   PetscFunctionBegin;
136   if (sac->SAVIENNACL) {
137     try {
138       delete sac->SAVIENNACL;
139     } catch(char * ex) {
140       SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex);
141     }
142   }
143 
144   /*
145       Free the private data structure that was hanging off the PC
146   */
147   ierr = PetscFree(pc->data);CHKERRQ(ierr);
148   PetscFunctionReturn(0);
149 }
150 
151 static PetscErrorCode PCSetFromOptions_SAVIENNACL(PetscOptionItems *PetscOptionsObject,PC pc)
152 {
153   PetscErrorCode ierr;
154 
155   PetscFunctionBegin;
156   ierr = PetscOptionsHead(PetscOptionsObject,"SAVIENNACL options");CHKERRQ(ierr);
157   ierr = PetscOptionsTail();CHKERRQ(ierr);
158   PetscFunctionReturn(0);
159 }
160 
161 /* -------------------------------------------------------------------------- */
162 
163 /*MC
164      PCSAViennaCL  - A smoothed agglomeration algorithm that can be used via the CUDA, OpenCL, and OpenMP backends of ViennaCL
165 
166    Level: advanced
167 
168 .seealso:  PCCreate(), PCSetType(), PCType (for list of available types), PC
169 
170 M*/
171 
172 PETSC_EXTERN PetscErrorCode PCCreate_SAVIENNACL(PC pc)
173 {
174   PC_SAVIENNACL  *sac;
175   PetscErrorCode ierr;
176 
177   PetscFunctionBegin;
178   /*
179      Creates the private data structure for this preconditioner and
180      attach it to the PC object.
181   */
182   ierr     = PetscNewLog(pc,&sac);CHKERRQ(ierr);
183   pc->data = (void*)sac;
184 
185   /*
186      Initialize the pointer to zero
187      Initialize number of v-cycles to default (1)
188   */
189   sac->SAVIENNACL = 0;
190 
191   /*
192       Set the pointers for the functions that are provided above.
193       Now when the user-level routines (such as PCApply(), PCDestroy(), etc.)
194       are called, they will automatically call these functions.  Note we
195       choose not to provide a couple of these functions since they are
196       not needed.
197   */
198   pc->ops->apply               = PCApply_SAVIENNACL;
199   pc->ops->applytranspose      = 0;
200   pc->ops->setup               = PCSetUp_SAVIENNACL;
201   pc->ops->destroy             = PCDestroy_SAVIENNACL;
202   pc->ops->setfromoptions      = PCSetFromOptions_SAVIENNACL;
203   pc->ops->view                = 0;
204   pc->ops->applyrichardson     = 0;
205   pc->ops->applysymmetricleft  = 0;
206   pc->ops->applysymmetricright = 0;
207   PetscFunctionReturn(0);
208 }
209 
210