Actual source code: taolinearsolver_petsc.c
1: #include "tao_general.h"
3: #ifdef TAO_USE_PETSC
5: #include "taolinearsolver_petsc.h"
7: #include "../matrix/taomat_petsc.h"
8: #include "../vector/taovec_petsc.h"
10: TaoLinearSolverPetsc::TaoLinearSolverPetsc(KSP S):TaoLinearSolver(){
11: linear_its=0;
12: ksp=0;
13: this->pkspviewer=PETSC_VIEWER_STDOUT_WORLD;
14: this->SetKSP(S);
15: return;
16: }
18: TaoLinearSolverPetsc::~TaoLinearSolverPetsc(){
19: if (ksp){
20: KSPDestroy(ksp);
21: }
22: return;
23: }
28: /*@C
29: TaoWrapKSP - Create a new TaoLinearSolver object using PETSc KSP.
31: Input Parameter:
32: . S - a KSP
34: Output Parameter:
35: . SS - new TaoMat
37: Note:
38: A TaoLinearSolverPetsc is an object with the methods of an abstract
39: TaoLinearSolver object. A TaoLinearSolverPetsc contains an implementation
40: of the TaoLinearSolver methods. Routines using these vectors should
41: declare a pointer to a TaoLinearSolver, assign this pointer to the address
42: of a TaoLinearSolver object, use the pointer to invoke methods on the
43: object, and use this pointer as an argument when calling other routines.
44: This usage is different from the usage of a PETSc KSP. In PETSc,
45: applications will typically declare a KSP, and pass it as an argument into
46: routines. That is, applications will typically declare a pointer to a
47: TaoLinearSolver and use the pointer, or declare a KSP and use it directly.
49: Level: developer
51: @*/
52: int TaoWrapKSP( KSP S, TaoLinearSolverPetsc ** SS){
55: *SS = new TaoLinearSolverPetsc(S);
56: return(0);
57: }
61: /*@C
62: TaoLinearSolverGetKSP - If the TaoLinearSolver is of the TaoLinearSolverPetsc type, it gets the underlying PETSc KSP
64: Input Parameter:
65: + TS - the TaoLinearSolver
66: - S - the address of KSP
68: Output Parameter:
69: . S - address of the underlying PETSc KSP object
72: Note:
73: This routine does not create a KSP. It sets a pointer
74: to the location of an existing KSP.
76: Level: advanced
78: .seealso TaoVecGetPetscVec(), TaoMatGetPetscMat(), TaoAppGetKSP()
79: @*/
80: int TaoLinearSolverGetKSP( TaoLinearSolver *TS, KSP * S){
82: if (TS){
83: *S=((TaoLinearSolverPetsc *)TS)->GetKSP();
85: } else {
86: *S=0;
87: }
88: return(0);
89: }
94: int TaoLinearSolverPetsc::SetKSP(KSP ksp2){
95: int info;
97: if (ksp2){
99: PetscObjectReference((PetscObject)ksp2);
100: }
101: if (ksp){
102: info=KSPDestroy(ksp); CHKERRQ(info);
103: }
104: ksp=ksp2;
105: return(0);
106: }
110: int TaoLinearSolverPetsc::GetKSP(KSP *ksp2){
112: if (ksp2){
113: *ksp2=this->ksp;
114: }
115: return(0);
116: }
120: int TaoLinearSolverPetsc::PreSolve(TaoMat* m1)
121: {
122: TaoMatPetsc *MM = dynamic_cast <TaoMatPetsc *> (m1);
123: Mat mm, mm_pre;
124: MatStructure preflag;
125: int info;
128: info = MM->GetMatrix(&mm, &mm_pre, &preflag); CHKERRQ(info);
129: info = KSPSetOperators(ksp, mm, mm_pre, preflag); CHKERRQ(info);
130: return(0);
131: }
135: int TaoLinearSolverPetsc::SetTrustRadius(double rad)
136: {
137: KSPType ktype;
138: int info;
139: PetscTruth flg;
143: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
145: info = PetscStrcmp((char *)ktype, KSPSTCG, &flg); CHKERRQ(info);
146: if (flg == PETSC_TRUE) {
147: info = KSPSTCGSetRadius(ksp, rad); CHKERRQ(info);
148: }
150: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
151: if (flg == PETSC_TRUE) {
152: info = KSPGLTRSetRadius(ksp, rad); CHKERRQ(info);
153: }
155: return(0);
156: }
160: int TaoLinearSolverPetsc::GetNormDirection(double *norm_d)
161: {
162: KSPType ktype;
163: int info;
164: PetscTruth flg;
168: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
170: info = PetscStrcmp((char *)ktype, KSPSTCG, &flg); CHKERRQ(info);
171: if (flg == PETSC_TRUE) {
172: info = KSPSTCGGetNormD(ksp, norm_d); CHKERRQ(info);
173: }
175: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
176: if (flg == PETSC_TRUE) {
177: info = KSPGLTRGetNormD(ksp, norm_d); CHKERRQ(info);
178: }
180: return(0);
181: }
185: int TaoLinearSolverPetsc::GetObjFcn(double *o_fcn)
186: {
187: KSPType ktype;
188: int info;
189: PetscTruth flg;
193: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
195: info = PetscStrcmp((char *)ktype, KSPSTCG, &flg); CHKERRQ(info);
196: if (flg == PETSC_TRUE) {
197: info = KSPSTCGGetObjFcn(ksp, o_fcn); CHKERRQ(info);
198: }
200: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
201: if (flg == PETSC_TRUE) {
202: info = KSPGLTRGetObjFcn(ksp, o_fcn); CHKERRQ(info);
203: }
205: return(0);
206: }
210: int TaoLinearSolverPetsc::GetMinEig(double *e_min)
211: {
212: KSPType ktype;
213: int info;
214: PetscTruth flg;
218: *e_min = 0.0;
220: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
221: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
222: if (flg == PETSC_TRUE) {
223: info = KSPGLTRGetMinEig(ksp, e_min); CHKERRQ(info);
224: }
226: return(0);
227: }
231: int TaoLinearSolverPetsc::GetLambda(double *lambda)
232: {
233: KSPType ktype;
234: int info;
235: PetscTruth flg;
239: *lambda = 0.0;
241: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
242: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
243: if (flg == PETSC_TRUE) {
244: info = KSPGLTRGetLambda(ksp, lambda); CHKERRQ(info);
245: }
247: return(0);
248: }
252: int TaoLinearSolverPetsc::Solve(TaoVec* tv, TaoVec* tw, TaoTruth *flag)
253: {
254: TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
255: TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);
256: int info,its;
260: info = KSPSolve(ksp,pv->GetVec(),pw->GetVec()); CHKERRQ(info);
261: info = KSPGetIterationNumber(ksp, &its); CHKERRQ(info);
263: this->linear_its=PetscMax(its,-its);
264: if (its>0) *flag=TAO_TRUE;
265: else *flag=TAO_FALSE;
267: return(0);
268: }
272: int TaoLinearSolverPetsc::SolveTrustRegion(TaoVec *b, TaoVec *x,
273: double r, TaoTruth *flg)
274: {
275: int info;
278: info = this->SetTrustRadius(r); CHKERRQ(info);
279: info = this->Solve(b, x, flg); CHKERRQ(info);
280: return(0);
281: }
285: int TaoLinearSolverPetsc::GetNumberIterations(int * iters){
287: *iters=linear_its;
288: return(0);
289: }
294: int TaoLinearSolverPetsc::SetTolerances(double rtol, double atol,
295: double dtol, int maxits) {
296: int info;
298: info = KSPSetTolerances(ksp, rtol, atol, dtol, maxits); CHKERRQ(info);
299: return(0);
300: }
304: int TaoLinearSolverPetsc::View(){
305: int info;
307: info=KSPView(this->ksp,this->pkspviewer);CHKERRQ(info);
308: return(0);
309: }
313: int TaoLinearSolverPetsc::SetOptions(){
314: int info;
316: info=KSPSetFromOptions(ksp); CHKERRQ(info);
317: return(0);
318: }
321: #endif