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: {
 24:   PetscErrorCode ierr;
 25:   int            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:   PetscLogInfo(st,"ST: linear solve iterations=%d\n",its);
 39:   return(0);
 40: }

 44: /*@
 45:    STSetKSP - Sets the KSP object associated with the spectral 
 46:    transformation.

 48:    Not collective

 50:    Input Parameters:
 51: +  st   - the spectral transformation context
 52: -  ksp  - the linear system context

 54:    Level: advanced

 56: @*/
 57: PetscErrorCode STSetKSP(ST st,KSP ksp)
 58: {
 59:   PetscErrorCode ierr;

 65:   if (st->ksp) {
 66:     KSPDestroy(st->ksp);
 67:   }
 68:   st->ksp = ksp;
 69:   PetscObjectReference((PetscObject)st->ksp);
 70:   return(0);
 71: }

 75: /*@
 76:    STGetKSP - Gets the KSP object associated with the spectral
 77:    transformation.

 79:    Not collective

 81:    Input Parameter:
 82: .  st - the spectral transformation context

 84:    Output Parameter:
 85: .  ksp  - the linear system context

 87:    Notes:
 88:    On output, the value of ksp can be PETSC_NULL if the combination of 
 89:    eigenproblem type and selected transformation does not require to 
 90:    solve a linear system of equations.
 91:    
 92:    Level: intermediate

 94: @*/
 95: PetscErrorCode STGetKSP(ST st,KSP* ksp)
 96: {
 99:   if (!st->type_name) { SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Must call STSetType first"); }
100:   if (ksp)  *ksp = st->ksp;
101:   return(0);
102: }

106: /*@
107:    STGetNumberLinearIterations - Gets the total number of linear iterations
108:    used by the ST object.

110:    Not Collective

112:    Input Parameter:
113: .  st - the spectral transformation context

115:    Output Parameter:
116: .  lits - number of linear iterations

118:    Level: intermediate

120: .seealso: STResetNumberLinearIterations()
121: @*/
122: PetscErrorCode STGetNumberLinearIterations(ST st,int* lits)
123: {
127:   *lits = st->lineariterations;
128:   return(0);
129: }

133: /*@
134:    STResetNumberLinearIterations - Resets the counter for total number of 
135:    linear iterations used by the ST object.

137:    Collective on ST

139:    Input Parameter:
140: .  st - the spectral transformation context

142:    Level: intermediate

144: .seealso: STGetNumberLinearIterations()
145: @*/
146: PetscErrorCode STResetNumberLinearIterations(ST st)
147: {
150:   st->lineariterations = 0;
151:   return(0);
152: }

156: PetscErrorCode STCheckNullSpace_Default(ST st,int n,Vec* V)
157: {
158:   PetscErrorCode ierr;
159:   int            i,c;
160:   PetscReal      norm;
161:   Vec            *T,w;
162:   Mat            A;
163:   PC             pc;
164:   MatNullSpace   nullsp;
165: 
167:   PetscMalloc(n*sizeof(Vec),&T);
168:   KSPGetPC(st->ksp,&pc);
169:   PCGetOperators(pc,&A,PETSC_NULL,PETSC_NULL);
170:   MatGetVecs(A,PETSC_NULL,&w);
171:   c = 0;
172:   for (i=0;i<n;i++) {
173:     MatMult(A,V[i],w);
174:     VecNorm(w,NORM_2,&norm);
175:     if (norm < 1e-8) {
176:       PetscLogInfo(st,"STCheckNullSpace: vector %i norm=%g\n",i,norm);
177:       T[c] = V[i];
178:       c++;
179:     }
180:   }
181:   VecDestroy(w);
182:   if (c>0) {
183:     MatNullSpaceCreate(st->comm,PETSC_FALSE,c,T,&nullsp);
184:     KSPSetNullSpace(st->ksp,nullsp);
185:     MatNullSpaceDestroy(nullsp);
186:   }
187:   PetscFree(T);
188:   return(0);
189: }

193: /*@C
194:    STCheckNullSpace - Given a set of vectors, this function test each of
195:    them to be a nullspace vector of the coefficient matrix of the associated
196:    KSP object. All these nullspace vectors are passed to the KSP object.

198:    Collective on ST

200:    Input Parameters:
201: +  st - the spectral transformation context
202: .  n  - number of vectors
203: -  V  - vectors to be checked

205:    Note:
206:    This function allows to handle singular pencils and to solve some problems
207:    in which the nullspace is important (see the users guide for details).
208:    
209:    Level: developer

211: .seealso: EPSAttachDeflationSpace()
212: @*/
213: PetscErrorCode STCheckNullSpace(ST st,int n,Vec* V)
214: {
215:   PetscErrorCode ierr;

218:   if (n>0 && st->checknullspace) {
219:     (*st->checknullspace)(st,n,V);
220:   }
221:   return(0);
222: }