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