Actual source code: iguess.c

  1: #include <petsc/private/kspimpl.h>

  3: PetscFunctionList KSPGuessList = NULL;
  4: static PetscBool  KSPGuessRegisterAllCalled;

  6: /*
  7:   KSPGuessRegister -  Adds a method for initial guess computation in Krylov subspace solver package.

  9:    Not Collective

 11:    Input Parameters:
 12: +  name_solver - name of a new user-defined solver
 13: -  routine_create - routine to create method context

 15:    Notes:
 16:    KSPGuessRegister() may be called multiple times to add several user-defined solvers.

 18:    Sample usage:
 19: .vb
 20:    KSPGuessRegister("my_initial_guess",MyInitialGuessCreate);
 21: .ve

 23:    Then, it can be chosen with the procedural interface via
 24: $     KSPSetGuessType(ksp,"my_initial_guess")
 25:    or at runtime via the option
 26: $     -ksp_guess_type my_initial_guess

 28:    Level: advanced

 30: .seealso: `KSPGuess`, `KSPGuessRegisterAll()`

 32: @*/
 33: PetscErrorCode KSPGuessRegister(const char sname[], PetscErrorCode (*function)(KSPGuess))
 34: {
 35:   KSPInitializePackage();
 36:   PetscFunctionListAdd(&KSPGuessList, sname, function);
 37:   return 0;
 38: }

 40: /*
 41:   KSPGuessRegisterAll - Registers all KSPGuess implementations in the KSP package.

 43:   Not Collective

 45:   Level: advanced

 47: .seealso: `KSPRegisterAll()`, `KSPInitializePackage()`
 48: */
 49: PetscErrorCode KSPGuessRegisterAll(void)
 50: {
 51:   if (KSPGuessRegisterAllCalled) return 0;
 52:   KSPGuessRegisterAllCalled = PETSC_TRUE;
 53:   KSPGuessRegister(KSPGUESSFISCHER, KSPGuessCreate_Fischer);
 54:   KSPGuessRegister(KSPGUESSPOD, KSPGuessCreate_POD);
 55:   return 0;
 56: }

 58: /*@
 59:     KSPGuessSetFromOptions - Sets the options for a KSPGuess from the options database

 61:     Collective on guess

 63:     Input Parameter:
 64: .    guess - KSPGuess object

 66:    Level: intermediate

 68: .seealso: `KSPGuess`, `KSPGetGuess()`, `KSPSetGuessType()`, `KSPGuessType`
 69: @*/
 70: PetscErrorCode KSPGuessSetFromOptions(KSPGuess guess)
 71: {
 73:   PetscTryTypeMethod(guess, setfromoptions);
 74:   return 0;
 75: }

 77: /*@
 78:     KSPGuessSetTolerance - Sets the relative tolerance used in either eigenvalue (POD) or singular value (Fischer type 3) calculations. Ignored by the first and second Fischer types.

 80:     Collective on guess

 82:     Input Parameter:
 83: .    guess - KSPGuess object

 85:    Level: intermediate

 87: .seealso: `KSPGuess`, `KSPGuessType`, `KSPGuessSetFromOptions()`
 88: @*/
 89: PetscErrorCode KSPGuessSetTolerance(KSPGuess guess, PetscReal tol)
 90: {
 92:   PetscTryTypeMethod(guess, settolerance, tol);
 93:   return 0;
 94: }

 96: /*@
 97:    KSPGuessDestroy - Destroys KSPGuess context.

 99:    Collective on kspGuess

101:    Input Parameter:
102: .  guess - initial guess object

104:    Level: beginner

106: .seealso: `KSPGuessCreate()`, `KSPGuess`, `KSPGuessType`
107: @*/
108: PetscErrorCode KSPGuessDestroy(KSPGuess *guess)
109: {
110:   if (!*guess) return 0;
112:   if (--((PetscObject)(*guess))->refct > 0) {
113:     *guess = NULL;
114:     return 0;
115:   }
116:   PetscTryTypeMethod((*guess), destroy);
117:   MatDestroy(&(*guess)->A);
118:   PetscHeaderDestroy(guess);
119:   return 0;
120: }

122: /*@C
123:    KSPGuessView - View the KSPGuess object

125:    Logically Collective on guess

127:    Input Parameters:
128: +  guess  - the initial guess object for the Krylov method
129: -  viewer - the viewer object

131:    Notes:

133:   Level: intermediate

135: .seealso: `KSP`, `KSPGuess`, `KSPGuessType`, `KSPGuessRegister()`, `KSPGuessCreate()`, `PetscViewer`
136: @*/
137: PetscErrorCode KSPGuessView(KSPGuess guess, PetscViewer view)
138: {
139:   PetscBool ascii;

142:   if (!view) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)guess), &view);
145:   PetscObjectTypeCompare((PetscObject)view, PETSCVIEWERASCII, &ascii);
146:   if (ascii) {
147:     PetscObjectPrintClassNamePrefixType((PetscObject)guess, view);
148:     PetscViewerASCIIPushTab(view);
149:     PetscTryTypeMethod(guess, view, view);
150:     PetscViewerASCIIPopTab(view);
151:   }
152:   return 0;
153: }

155: /*@
156:    KSPGuessCreate - Creates the default KSPGuess context.

158:    Collective

160:    Input Parameter:
161: .  comm - MPI communicator

163:    Output Parameter:
164: .  guess - location to put the KSPGuess context

166:    Notes:
167:    The default KSPGuess type is XXX

169:    Level: beginner

171: .seealso: `KSPSolve()`, `KSPGuessDestroy()`, `KSPGuess`, `KSPGuessType`, `KSP`
172: @*/
173: PetscErrorCode KSPGuessCreate(MPI_Comm comm, KSPGuess *guess)
174: {
175:   KSPGuess tguess;

178:   *guess = NULL;
179:   KSPInitializePackage();
180:   PetscHeaderCreate(tguess, KSPGUESS_CLASSID, "KSPGuess", "Initial guess for Krylov Method", "KSPGuess", comm, KSPGuessDestroy, KSPGuessView);
181:   tguess->omatstate = -1;
182:   *guess            = tguess;
183:   return 0;
184: }

186: /*@C
187:    KSPGuessSetType - Sets the type of a KSPGuess

189:    Logically Collective on guess

191:    Input Parameters:
192: +  guess - the initial guess object for the Krylov method
193: -  type  - a known KSPGuess method

195:    Options Database Key:
196: .  -ksp_guess_type  <method> - Sets the method; use -help for a list
197:     of available methods

199:    Notes:

201:   Level: intermediate

203: .seealso: `KSP`, `KSPGuess`, `KSPGuessType`, `KSPGuessRegister()`, `KSPGuessCreate()`

205: @*/
206: PetscErrorCode KSPGuessSetType(KSPGuess guess, KSPGuessType type)
207: {
208:   PetscBool match;
209:   PetscErrorCode (*r)(KSPGuess);


214:   PetscObjectTypeCompare((PetscObject)guess, type, &match);
215:   if (match) return 0;

217:   PetscFunctionListFind(KSPGuessList, type, &r);
219:   PetscTryTypeMethod(guess, destroy);
220:   guess->ops->destroy = NULL;

222:   PetscMemzero(guess->ops, sizeof(struct _KSPGuessOps));
223:   PetscObjectChangeTypeName((PetscObject)guess, type);
224:   (*r)(guess);
225:   return 0;
226: }

228: /*@C
229:    KSPGuessGetType - Gets the KSPGuess type as a string from the KSPGuess object.

231:    Not Collective

233:    Input Parameter:
234: .  guess - the initial guess context

236:    Output Parameter:
237: .  name - name of KSPGuess method

239:    Level: intermediate

241: .seealso: `KSPGuessSetType()`
242: @*/
243: PetscErrorCode KSPGuessGetType(KSPGuess guess, KSPGuessType *type)
244: {
247:   *type = ((PetscObject)guess)->type_name;
248:   return 0;
249: }

251: /*@
252:     KSPGuessUpdate - Updates the guess object with the current solution and rhs vector

254:    Collective on guess

256:    Input Parameters:
257: +  guess - the initial guess context
258: .  rhs   - the corresponding rhs
259: -  sol   - the computed solution

261:    Level: intermediate

263: .seealso: `KSPGuessCreate()`, `KSPGuess`
264: @*/
265: PetscErrorCode KSPGuessUpdate(KSPGuess guess, Vec rhs, Vec sol)
266: {
270:   PetscTryTypeMethod(guess, update, rhs, sol);
271:   return 0;
272: }

274: /*@
275:     KSPGuessFormGuess - Form the initial guess

277:    Collective on guess

279:    Input Parameters:
280: +  guess - the initial guess context
281: .  rhs   - the current rhs vector
282: -  sol   - the initial guess vector

284:    Level: intermediate

286: .seealso: `KSPGuessCreate()`, `KSPGuess`
287: @*/
288: PetscErrorCode KSPGuessFormGuess(KSPGuess guess, Vec rhs, Vec sol)
289: {
293:   PetscTryTypeMethod(guess, formguess, rhs, sol);
294:   return 0;
295: }

297: /*@
298:     KSPGuessSetUp - Setup the initial guess object

300:    Collective on guess

302:    Input Parameter:
303: -  guess - the initial guess context

305:    Level: intermediate

307: .seealso: `KSPGuessCreate()`, `KSPGuess`
308: @*/
309: PetscErrorCode KSPGuessSetUp(KSPGuess guess)
310: {
311:   PetscObjectState matstate;
312:   PetscInt         oM = 0, oN = 0, M, N;
313:   Mat              omat = NULL;
314:   PC               pc;
315:   PetscBool        reuse;

318:   if (guess->A) {
319:     omat = guess->A;
320:     MatGetSize(guess->A, &oM, &oN);
321:   }
322:   KSPGetOperators(guess->ksp, &guess->A, NULL);
323:   KSPGetPC(guess->ksp, &pc);
324:   PCGetReusePreconditioner(pc, &reuse);
325:   PetscObjectReference((PetscObject)guess->A);
326:   MatGetSize(guess->A, &M, &N);
327:   PetscObjectStateGet((PetscObject)guess->A, &matstate);
328:   if (M != oM || N != oN) {
329:     PetscInfo(guess, "Resetting KSPGuess since matrix sizes have changed (%" PetscInt_FMT " != %" PetscInt_FMT ", %" PetscInt_FMT " != %" PetscInt_FMT ")\n", oM, M, oN, N);
330:   } else if (!reuse && (omat != guess->A || guess->omatstate != matstate)) {
331:     PetscInfo(guess, "Resetting KSPGuess since %s has changed\n", omat != guess->A ? "matrix" : "matrix state");
332:     PetscTryTypeMethod(guess, reset);
333:   } else if (reuse) {
334:     PetscInfo(guess, "Not resettting KSPGuess since reuse preconditioner has been specified\n");
335:   } else {
336:     PetscInfo(guess, "KSPGuess status unchanged\n");
337:   }
338:   PetscTryTypeMethod(guess, setup);
339:   guess->omatstate = matstate;
340:   MatDestroy(&omat);
341:   return 0;
342: }