Actual source code: stsles.c
2: /*
3: The ST (spectral transformation) interface routines related to the
4: KSP object associated to it.
5: */
7: #include src/st/stimpl.h
11: /*
12: STAssociatedKSPSolve - Solves the linear system of equations associated
13: to the spectral transformation.
15: Input Parameters:
16: . st - the spectral transformation context
17: . b - right hand side vector
19: Output Parameter:
20: . x - computed solution
21: */
22: PetscErrorCode STAssociatedKSPSolve(ST st,Vec b,Vec x)
23: {
25: PetscInt its;
26: KSPConvergedReason reason;
32: if (!st->ksp) { SETERRQ(PETSC_ERR_SUP,"ST has no associated KSP"); }
33: KSPSolve(st->ksp,b,x);
34: KSPGetConvergedReason(st->ksp,&reason);
35: if (reason<0) { SETERRQ1(0,"Warning: KSP did not converge (%d)",reason); }
36: KSPGetIterationNumber(st->ksp,&its);
37: st->lineariterations += its;
38: PetscInfo1(st,"Linear solve iterations=%d\n",its);
39: return(0);
40: }
44: /*
45: STAssociatedKSPSolveTranspose - Solves the transpose of the linear
46: system of equations associated to the spectral transformation.
48: Input Parameters:
49: . st - the spectral transformation context
50: . b - right hand side vector
52: Output Parameter:
53: . x - computed solution
54: */
55: PetscErrorCode STAssociatedKSPSolveTranspose(ST st,Vec b,Vec x)
56: {
58: PetscInt its;
59: KSPConvergedReason reason;
65: if (!st->ksp) { SETERRQ(PETSC_ERR_SUP,"ST has no associated KSP"); }
66: KSPSolveTranspose(st->ksp,b,x);
67: KSPGetConvergedReason(st->ksp,&reason);
68: if (reason<0) { SETERRQ1(0,"Warning: KSP did not converge (%d)",reason); }
69: KSPGetIterationNumber(st->ksp,&its);
70: st->lineariterations += its;
71: PetscInfo1(st,"Linear solve iterations=%d\n",its);
72: return(0);
73: }
77: /*@
78: STSetKSP - Sets the KSP object associated with the spectral
79: transformation.
81: Not collective
83: Input Parameters:
84: + st - the spectral transformation context
85: - ksp - the linear system context
87: Level: advanced
89: @*/
90: PetscErrorCode STSetKSP(ST st,KSP ksp)
91: {
98: if (st->ksp) {
99: KSPDestroy(st->ksp);
100: }
101: st->ksp = ksp;
102: PetscObjectReference((PetscObject)st->ksp);
103: return(0);
104: }
108: /*@
109: STGetKSP - Gets the KSP object associated with the spectral
110: transformation.
112: Not collective
114: Input Parameter:
115: . st - the spectral transformation context
117: Output Parameter:
118: . ksp - the linear system context
120: Notes:
121: On output, the value of ksp can be PETSC_NULL if the combination of
122: eigenproblem type and selected transformation does not require to
123: solve a linear system of equations.
124:
125: Level: intermediate
127: @*/
128: PetscErrorCode STGetKSP(ST st,KSP* ksp)
129: {
132: if (!st->type_name) { SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call STSetType first"); }
133: if (ksp) *ksp = st->ksp;
134: return(0);
135: }
139: /*@
140: STGetOperationCounters - Gets the total number of operator applications,
141: inner product operations and linear iterations used by the ST object.
143: Not Collective
145: Input Parameter:
146: . st - the spectral transformation context
148: Output Parameter:
149: + ops - number of operator applications
150: . dots - number of inner product operations
151: - lits - number of linear iterations
153: Notes:
154: Any output parameter may be PETSC_NULL on input if not needed.
155:
156: Level: intermediate
158: .seealso: STResetOperationCounters()
159: @*/
160: PetscErrorCode STGetOperationCounters(ST st,int* ops,int* dots,int* lits)
161: {
164: if (ops) *ops = st->applys;
165: if (dots) *dots = st->innerproducts;
166: if (lits) *lits = st->lineariterations;
167: return(0);
168: }
172: /*@
173: STResetOperationCounters - Resets the counters for operator applications,
174: inner product operations and total number of linear iterations used by
175: the ST object.
177: Collective on ST
179: Input Parameter:
180: . st - the spectral transformation context
182: Level: intermediate
184: .seealso: STGetOperationCounters()
185: @*/
186: PetscErrorCode STResetOperationCounters(ST st)
187: {
190: st->lineariterations = 0;
191: st->applys = 0;
192: st->innerproducts = 0;
193: return(0);
194: }
198: PetscErrorCode STCheckNullSpace_Default(ST st,int n,const Vec V[])
199: {
201: int i,c;
202: PetscReal norm;
203: Vec *T,w;
204: Mat A;
205: PC pc;
206: MatNullSpace nullsp;
207:
209: PetscMalloc(n*sizeof(Vec),&T);
210: KSPGetPC(st->ksp,&pc);
211: PCGetOperators(pc,&A,PETSC_NULL,PETSC_NULL);
212: MatGetVecs(A,PETSC_NULL,&w);
213: c = 0;
214: for (i=0;i<n;i++) {
215: MatMult(A,V[i],w);
216: VecNorm(w,NORM_2,&norm);
217: if (norm < 1e-8) {
218: PetscInfo2(st,"Vector %i norm=%g\n",i,norm);
219: T[c] = V[i];
220: c++;
221: }
222: }
223: VecDestroy(w);
224: if (c>0) {
225: MatNullSpaceCreate(st->comm,PETSC_FALSE,c,T,&nullsp);
226: KSPSetNullSpace(st->ksp,nullsp);
227: MatNullSpaceDestroy(nullsp);
228: }
229: PetscFree(T);
230: return(0);
231: }
235: /*@
236: STCheckNullSpace - Given a set of vectors, this function tests each of
237: them to be a nullspace vector of the coefficient matrix of the associated
238: KSP object. All these nullspace vectors are passed to the KSP object.
240: Collective on ST
242: Input Parameters:
243: + st - the spectral transformation context
244: . n - number of vectors
245: - V - vectors to be checked
247: Note:
248: This function allows to handle singular pencils and to solve some problems
249: in which the nullspace is important (see the users guide for details).
250:
251: Level: developer
253: .seealso: EPSAttachDeflationSpace()
254: @*/
255: PetscErrorCode STCheckNullSpace(ST st,int n,const Vec V[])
256: {
260: if (n>0 && st->checknullspace) {
261: (*st->checknullspace)(st,n,V);
262: }
263: return(0);
264: }