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: const KSPType ktype;
138: int info;
139: PetscTruth flg;
143: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
145: info = PetscStrcmp((char *)ktype, KSPNASH, &flg); CHKERRQ(info);
146: if (flg == PETSC_TRUE) {
147: info = KSPNASHSetRadius(ksp, rad); CHKERRQ(info);
148: }
150: info = PetscStrcmp((char *)ktype, KSPSTCG, &flg); CHKERRQ(info);
151: if (flg == PETSC_TRUE) {
152: info = KSPSTCGSetRadius(ksp, rad); CHKERRQ(info);
153: }
155: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
156: if (flg == PETSC_TRUE) {
157: info = KSPGLTRSetRadius(ksp, rad); CHKERRQ(info);
158: }
160: return(0);
161: }
165: int TaoLinearSolverPetsc::GetNormDirection(double *norm_d)
166: {
167: const KSPType ktype;
168: int info;
169: PetscTruth flg;
173: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
175: info = PetscStrcmp((char *)ktype, KSPNASH, &flg); CHKERRQ(info);
176: if (flg == PETSC_TRUE) {
177: info = KSPNASHGetNormD(ksp, norm_d); CHKERRQ(info);
178: }
180: info = PetscStrcmp((char *)ktype, KSPSTCG, &flg); CHKERRQ(info);
181: if (flg == PETSC_TRUE) {
182: info = KSPSTCGGetNormD(ksp, norm_d); CHKERRQ(info);
183: }
185: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
186: if (flg == PETSC_TRUE) {
187: info = KSPGLTRGetNormD(ksp, norm_d); CHKERRQ(info);
188: }
190: return(0);
191: }
195: int TaoLinearSolverPetsc::GetObjFcn(double *o_fcn)
196: {
197: const KSPType ktype;
198: int info;
199: PetscTruth flg;
203: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
205: info = PetscStrcmp((char *)ktype, KSPNASH, &flg); CHKERRQ(info);
206: if (flg == PETSC_TRUE) {
207: info = KSPNASHGetObjFcn(ksp, o_fcn); CHKERRQ(info);
208: }
210: info = PetscStrcmp((char *)ktype, KSPSTCG, &flg); CHKERRQ(info);
211: if (flg == PETSC_TRUE) {
212: info = KSPSTCGGetObjFcn(ksp, o_fcn); CHKERRQ(info);
213: }
215: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
216: if (flg == PETSC_TRUE) {
217: info = KSPGLTRGetObjFcn(ksp, o_fcn); CHKERRQ(info);
218: }
220: return(0);
221: }
225: int TaoLinearSolverPetsc::GetMinEig(double *e_min)
226: {
227: const KSPType ktype;
228: int info;
229: PetscTruth flg;
233: *e_min = 0.0;
235: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
236: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
237: if (flg == PETSC_TRUE) {
238: info = KSPGLTRGetMinEig(ksp, e_min); CHKERRQ(info);
239: }
241: return(0);
242: }
246: int TaoLinearSolverPetsc::GetLambda(double *lambda)
247: {
248: const KSPType ktype;
249: int info;
250: PetscTruth flg;
254: *lambda = 0.0;
256: info = KSPGetType(ksp, &ktype); CHKERRQ(info);
257: info = PetscStrcmp((char *)ktype, KSPGLTR, &flg); CHKERRQ(info);
258: if (flg == PETSC_TRUE) {
259: info = KSPGLTRGetLambda(ksp, lambda); CHKERRQ(info);
260: }
262: return(0);
263: }
267: int TaoLinearSolverPetsc::Solve(TaoVec* tv, TaoVec* tw, TaoTruth *flag)
268: {
269: TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
270: TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);
271: int info;
272: PetscInt its;
276: info = KSPSolve(ksp,pv->GetVec(),pw->GetVec()); CHKERRQ(info);
277: info = KSPGetIterationNumber(ksp, &its); CHKERRQ(info);
279: this->linear_its=PetscMax(its,-its);
280: if (its>0) *flag=TAO_TRUE;
281: else *flag=TAO_FALSE;
283: return(0);
284: }
288: int TaoLinearSolverPetsc::SolveTrustRegion(TaoVec *b, TaoVec *x,
289: double r, TaoTruth *flg)
290: {
291: int info;
294: info = this->SetTrustRadius(r); CHKERRQ(info);
295: info = this->Solve(b, x, flg); CHKERRQ(info);
296: return(0);
297: }
301: int TaoLinearSolverPetsc::GetNumberIterations(int * iters){
303: *iters=linear_its;
304: return(0);
305: }
310: int TaoLinearSolverPetsc::SetTolerances(double rtol, double atol,
311: double dtol, int maxits) {
312: int info;
314: info = KSPSetTolerances(ksp, rtol, atol, dtol, maxits); CHKERRQ(info);
315: return(0);
316: }
320: int TaoLinearSolverPetsc::View(){
321: int info;
323: info=KSPView(this->ksp,this->pkspviewer);CHKERRQ(info);
324: return(0);
325: }
329: int TaoLinearSolverPetsc::SetOptions(){
330: int info;
332: info=KSPSetFromOptions(ksp); CHKERRQ(info);
333: return(0);
334: }
337: #endif