Actual source code: epreconditioner.c
1: /*
2: Wrappers for PETSc PC ESI implementation
3: */
5: #include esi/petsc/preconditioner.h
7: esi::petsc::Preconditioner<double,int>::Preconditioner(PC ipc)
8: {
10: this->pc = ipc;
11: this->pobject = (PetscObject)this->pc;
12: PetscObjectGetComm((PetscObject)this->pc,&this->comm);if (ierr) return;
13: PetscObjectReference((PetscObject)ipc);if (ierr) return;
14: }
17: esi::petsc::Preconditioner<double,int>::~Preconditioner()
18: {
20: PetscObjectDereference((PetscObject)this->pc);if (ierr) return;
21: }
23: esi::ErrorCode esi::petsc::Preconditioner<double,int>::getInterface(const char* name, void *& iface)
24: {
25: PetscTruth flg;
27: if (!PetscStrcmp(name,"esi::Object",&flg),flg){
28: iface = (void *) (esi::Object *) this;
29: } else if (!PetscStrcmp(name,"esi::Operator",&flg),flg){
30: iface = (void *) (esi::Operator<double,int> *) this;
31: } else if (!PetscStrcmp(name,"esi::Preconditioner",&flg),flg){
32: iface = (void *) (esi::Preconditioner<double,int> *) this;
33: } else if (!PetscStrcmp(name,"PC",&flg),flg){
34: iface = (void *) this->pc;
35: } else if (!PetscStrcmp(name,"esi::petsc::Preconditioner",&flg),flg){
36: iface = (void *) (esi::petsc::Preconditioner<double,int> *) this;
37: } else {
38: iface = 0;
39: }
40: return 0;
41: }
43: esi::ErrorCode esi::petsc::Preconditioner<double,int>::getInterfacesSupported(esi::Argv * list)
44: {
45: list->appendArg("esi::Object");
46: list->appendArg("esi::Operator");
47: list->appendArg("esi::Preconditioner");
48: list->appendArg("esi::petsc::Preconditioner");
49: list->appendArg("PC");
50: return 0;
51: }
54: esi::ErrorCode esi::petsc::Preconditioner<double,int>::apply( esi::Vector<double,int> &xx,esi::Vector<double,int> &yy)
55: {
57: Vec py,px;
59: yy.getInterface("Vec",reinterpret_cast<void*&>(py));
60: xx.getInterface("Vec",reinterpret_cast<void*&>(px));
62: PCSetVector(this->pc,px);
63: return PCApply(this->pc,px,py);
64: }
66: esi::ErrorCode esi::petsc::Preconditioner<double,int>::solve( esi::Vector<double,int> &xx,esi::Vector<double,int> &yy)
67: {
68: return this->apply(xx,yy);
69: }
71: esi::ErrorCode esi::petsc::Preconditioner<double,int>::solveLeft( esi::Vector<double,int> &xx,esi::Vector<double,int> &yy)
72: {
74: Vec py,px;
76: yy.getInterface("Vec",reinterpret_cast<void*&>(py));
77: xx.getInterface("Vec",reinterpret_cast<void*&>(px));
79: return PCApplySymmetricLeft(this->pc,px,py);
80: }
82: esi::ErrorCode esi::petsc::Preconditioner<double,int>::solveRight( esi::Vector<double,int> &xx,esi::Vector<double,int> &yy)
83: {
85: Vec py,px;
87: yy.getInterface("Vec",reinterpret_cast<void*&>(py));
88: xx.getInterface("Vec",reinterpret_cast<void*&>(px));
90: return PCApplySymmetricRight(this->pc,px,py);
91: }
93: esi::ErrorCode esi::petsc::Preconditioner<double,int>::applyB( esi::Vector<double,int> &xx,esi::Vector<double,int> &yy)
94: {
95: int ierr;
96: Vec py,px,work;
97: PCSide iside;
99: yy.getInterface("Vec",reinterpret_cast<void*&>(py));
100: xx.getInterface("Vec",reinterpret_cast<void*&>(px));
101: VecDuplicate(py,&work);
102: if (this->side == esi::PRECONDITIONER_LEFT) iside = PC_LEFT;
103: if (this->side == esi::PRECONDITIONER_RIGHT) iside = PC_RIGHT;
104: if (this->side == esi::PRECONDITIONER_TWO_SIDED) iside = PC_SYMMETRIC;
105: PCApplyBAorAB(this->pc,iside,px,py,work);
106: VecDestroy(work);
107: return 0;
108: }
110: esi::ErrorCode esi::petsc::Preconditioner<double,int>::setup()
111: {
112: return 0;
113: }
115: esi::ErrorCode esi::petsc::Preconditioner<double,int>::setPreconditionerSide(esi::PreconditionerSide iside)
116: {
117: this->side = iside;
118: return 0;
119: }
121: esi::ErrorCode esi::petsc::Preconditioner<double,int>::getPreconditionerSide(esi::PreconditionerSide & iside)
122: {
123: iside = this->side;
124: return 0;
125: }
127: esi::ErrorCode esi::petsc::Preconditioner<double,int>::setOperator( esi::Operator<double,int> &op)
128: {
129: /*
130: For now require Operator to be a PETSc Mat
131: */
132: Mat A;
133: int op.getInterface("Mat",reinterpret_cast<void*&>(A));
134: PCSetOperators(this->pc,A,A,DIFFERENT_NONZERO_PATTERN);
135: return 0;
136: }