Actual source code: chowiluviennacl.cxx
2: /* -------------------------------------------------------------------- */
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
11: #include <petsc/private/pcimpl.h>
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>
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;
25: /* -------------------------------------------------------------------------- */
26: /*
27: PCSetUp_CHOWILUVIENNACL - Prepares for the use of the CHOWILUVIENNACL preconditioner
28: by setting data structures and options.
30: Input Parameter:
31: . pc - the preconditioner context
33: Application Interface Routine: PCSetUp()
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: Mat_SeqAIJViennaCL *gpustruct;
45: PetscObjectTypeCompare((PetscObject)pc->pmat,MATSEQAIJVIENNACL,&flg);
47: if (pc->setupcalled != 0) {
48: try {
49: delete ilu->CHOWILUVIENNACL;
50: } catch(char *ex) {
51: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex);
52: }
53: }
54: try {
55: #if defined(PETSC_USE_COMPLEX)
56: gpustruct = NULL;
57: SETERRQ(PetscObjectComm((PetscObject)pc),PETSC_ERR_SUP,"No support for complex arithmetic in CHOWILUVIENNACL preconditioner");
58: #else
59: MatViennaCLCopyToGPU(pc->pmat);
60: gpustruct = (Mat_SeqAIJViennaCL*)(pc->pmat->spptr);
62: viennacl::linalg::chow_patel_tag ilu_tag;
63: ViennaCLAIJMatrix *mat = (ViennaCLAIJMatrix*)gpustruct->mat;
64: ilu->CHOWILUVIENNACL = new viennacl::linalg::chow_patel_ilu_precond<viennacl::compressed_matrix<PetscScalar> >(*mat, ilu_tag);
65: #endif
66: } catch(char *ex) {
67: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex);
68: }
69: return 0;
70: }
72: /* -------------------------------------------------------------------------- */
73: /*
74: PCApply_CHOWILUVIENNACL - Applies the CHOWILUVIENNACL preconditioner to a vector.
76: Input Parameters:
77: . pc - the preconditioner context
78: . x - input vector
80: Output Parameter:
81: . y - output vector
83: Application Interface Routine: PCApply()
84: */
85: static PetscErrorCode PCApply_CHOWILUVIENNACL(PC pc,Vec x,Vec y)
86: {
87: PC_CHOWILUVIENNACL *ilu = (PC_CHOWILUVIENNACL*)pc->data;
88: PetscBool flg1,flg2;
89: viennacl::vector<PetscScalar> const *xarray=NULL;
90: viennacl::vector<PetscScalar> *yarray=NULL;
92: /*how to apply a certain fixed number of iterations?*/
93: PetscObjectTypeCompare((PetscObject)x,VECSEQVIENNACL,&flg1);
94: PetscObjectTypeCompare((PetscObject)y,VECSEQVIENNACL,&flg2);
96: if (!ilu->CHOWILUVIENNACL) {
97: PCSetUp_CHOWILUVIENNACL(pc);
98: }
99: VecSet(y,0.0);
100: VecViennaCLGetArrayRead(x,&xarray);
101: VecViennaCLGetArrayWrite(y,&yarray);
102: try {
103: #if defined(PETSC_USE_COMPLEX)
105: #else
106: *yarray = *xarray;
107: ilu->CHOWILUVIENNACL->apply(*yarray);
108: #endif
109: } catch(char * ex) {
110: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex);
111: }
112: VecViennaCLRestoreArrayRead(x,&xarray);
113: VecViennaCLRestoreArrayWrite(y,&yarray);
114: PetscObjectStateIncrease((PetscObject)y);
115: return 0;
116: }
117: /* -------------------------------------------------------------------------- */
118: /*
119: PCDestroy_CHOWILUVIENNACL - Destroys the private context for the CHOWILUVIENNACL preconditioner
120: that was created with PCCreate_CHOWILUVIENNACL().
122: Input Parameter:
123: . pc - the preconditioner context
125: Application Interface Routine: PCDestroy()
126: */
127: static PetscErrorCode PCDestroy_CHOWILUVIENNACL(PC pc)
128: {
129: PC_CHOWILUVIENNACL *ilu = (PC_CHOWILUVIENNACL*)pc->data;
131: if (ilu->CHOWILUVIENNACL) {
132: try {
133: delete ilu->CHOWILUVIENNACL;
134: } catch(char *ex) {
135: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"ViennaCL error: %s", ex);
136: }
137: }
139: /*
140: Free the private data structure that was hanging off the PC
141: */
142: PetscFree(pc->data);
143: return 0;
144: }
146: static PetscErrorCode PCSetFromOptions_CHOWILUVIENNACL(PetscOptionItems *PetscOptionsObject,PC pc)
147: {
148: PetscOptionsHead(PetscOptionsObject,"CHOWILUVIENNACL options");
149: PetscOptionsTail();
150: return 0;
151: }
153: /* -------------------------------------------------------------------------- */
155: /*MC
156: PCCHOWILUViennaCL - A smoothed agglomeration algorithm that can be used via the CUDA, OpenCL, and OpenMP backends of ViennaCL
158: Level: advanced
160: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC
162: M*/
164: PETSC_EXTERN PetscErrorCode PCCreate_CHOWILUVIENNACL(PC pc)
165: {
166: PC_CHOWILUVIENNACL *ilu;
168: /*
169: Creates the private data structure for this preconditioner and
170: attach it to the PC object.
171: */
172: PetscNewLog(pc,&ilu);
173: pc->data = (void*)ilu;
175: /*
176: Initialize the pointer to zero
177: Initialize number of v-cycles to default (1)
178: */
179: ilu->CHOWILUVIENNACL = 0;
181: /*
182: Set the pointers for the functions that are provided above.
183: Now when the user-level routines (such as PCApply(), PCDestroy(), etc.)
184: are called, they will automatically call these functions. Note we
185: choose not to provide a couple of these functions since they are
186: not needed.
187: */
188: pc->ops->apply = PCApply_CHOWILUVIENNACL;
189: pc->ops->applytranspose = 0;
190: pc->ops->setup = PCSetUp_CHOWILUVIENNACL;
191: pc->ops->destroy = PCDestroy_CHOWILUVIENNACL;
192: pc->ops->setfromoptions = PCSetFromOptions_CHOWILUVIENNACL;
193: pc->ops->view = 0;
194: pc->ops->applyrichardson = 0;
195: pc->ops->applysymmetricleft = 0;
196: pc->ops->applysymmetricright = 0;
197: return 0;
198: }