Actual source code: dmsnes.c

  1: #include <petsc/private/snesimpl.h>
  2: #include <petsc/private/dmimpl.h>

  4: static PetscErrorCode DMSNESDestroy(DMSNES *kdm)
  5: {
  6:   if (!*kdm) return 0;
  8:   if (--((PetscObject)(*kdm))->refct > 0) {*kdm = NULL; return 0;}
  9:   if ((*kdm)->ops->destroy) ((*kdm)->ops->destroy)(*kdm);
 10:   PetscHeaderDestroy(kdm);
 11:   return 0;
 12: }

 14: PetscErrorCode DMSNESLoad(DMSNES kdm,PetscViewer viewer)
 15: {
 16:   PetscViewerBinaryRead(viewer,&kdm->ops->computefunction,1,NULL,PETSC_FUNCTION);
 17:   PetscViewerBinaryRead(viewer,&kdm->ops->computejacobian,1,NULL,PETSC_FUNCTION);
 18:   return 0;
 19: }

 21: PetscErrorCode DMSNESView(DMSNES kdm,PetscViewer viewer)
 22: {
 23:   PetscBool      isascii,isbinary;

 25:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 26:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
 27:   if (isascii) {
 28: #if defined(PETSC_SERIALIZE_FUNCTIONS) && defined(PETSC_SERIALIZE_FUNCTIONS_VIEW)
 29:     const char *fname;

 31:     PetscFPTFind(kdm->ops->computefunction,&fname);
 32:     if (fname) {
 33:       PetscViewerASCIIPrintf(viewer,"Function used by SNES: %s\n",fname);
 34:     }
 35:     PetscFPTFind(kdm->ops->computejacobian,&fname);
 36:     if (fname) {
 37:       PetscViewerASCIIPrintf(viewer,"Jacobian function used by SNES: %s\n",fname);
 38:     }
 39: #endif
 40:   } else if (isbinary) {
 41:     struct {
 42:       PetscErrorCode (*func)(SNES,Vec,Vec,void*);
 43:     } funcstruct;
 44:     struct {
 45:       PetscErrorCode (*jac)(SNES,Vec,Mat,Mat,void*);
 46:     } jacstruct;
 47:     funcstruct.func = kdm->ops->computefunction;
 48:     jacstruct.jac   = kdm->ops->computejacobian;
 49:     PetscViewerBinaryWrite(viewer,&funcstruct,1,PETSC_FUNCTION);
 50:     PetscViewerBinaryWrite(viewer,&jacstruct,1,PETSC_FUNCTION);
 51:   }
 52:   return 0;
 53: }

 55: static PetscErrorCode DMSNESCreate(MPI_Comm comm,DMSNES *kdm)
 56: {
 57:   SNESInitializePackage();
 58:   PetscHeaderCreate(*kdm, DMSNES_CLASSID,  "DMSNES", "DMSNES", "DMSNES", comm, DMSNESDestroy, DMSNESView);
 59:   return 0;
 60: }

 62: /* Attaches the DMSNES to the coarse level.
 63:  * Under what conditions should we copy versus duplicate?
 64:  */
 65: static PetscErrorCode DMCoarsenHook_DMSNES(DM dm,DM dmc,void *ctx)
 66: {
 67:   DMCopyDMSNES(dm,dmc);
 68:   return 0;
 69: }

 71: /* This could restrict auxiliary information to the coarse level.
 72:  */
 73: static PetscErrorCode DMRestrictHook_DMSNES(DM dm,Mat Restrict,Vec rscale,Mat Inject,DM dmc,void *ctx)
 74: {
 75:   return 0;
 76: }

 78: /* Attaches the DMSNES to the subdomain. */
 79: static PetscErrorCode DMSubDomainHook_DMSNES(DM dm,DM subdm,void *ctx)
 80: {
 81:   DMCopyDMSNES(dm,subdm);
 82:   return 0;
 83: }

 85: /* This could restrict auxiliary information to the coarse level.
 86:  */
 87: static PetscErrorCode DMSubDomainRestrictHook_DMSNES(DM dm,VecScatter gscat,VecScatter lscat,DM subdm,void *ctx)
 88: {
 89:   return 0;
 90: }

 92: static PetscErrorCode DMRefineHook_DMSNES(DM dm,DM dmf,void *ctx)
 93: {
 94:   DMCopyDMSNES(dm,dmf);
 95:   return 0;
 96: }

 98: /* This could restrict auxiliary information to the coarse level.
 99:  */
100: static PetscErrorCode DMInterpolateHook_DMSNES(DM dm,Mat Interp,DM dmf,void *ctx)
101: {
102:   return 0;
103: }

105: /*@C
106:    DMSNESCopy - copies the information in a DMSNES to another DMSNES

108:    Not Collective

110:    Input Parameters:
111: +  kdm - Original DMSNES
112: -  nkdm - DMSNES to receive the data, should have been created with DMSNESCreate()

114:    Level: developer

116: .seealso: DMSNESCreate(), DMSNESDestroy()
117: @*/
118: PetscErrorCode DMSNESCopy(DMSNES kdm,DMSNES nkdm)
119: {
122:   nkdm->ops->computefunction  = kdm->ops->computefunction;
123:   nkdm->ops->computejacobian  = kdm->ops->computejacobian;
124:   nkdm->ops->computegs        = kdm->ops->computegs;
125:   nkdm->ops->computeobjective = kdm->ops->computeobjective;
126:   nkdm->ops->computepjacobian = kdm->ops->computepjacobian;
127:   nkdm->ops->computepfunction = kdm->ops->computepfunction;
128:   nkdm->ops->destroy          = kdm->ops->destroy;
129:   nkdm->ops->duplicate        = kdm->ops->duplicate;

131:   nkdm->functionctx  = kdm->functionctx;
132:   nkdm->gsctx        = kdm->gsctx;
133:   nkdm->pctx         = kdm->pctx;
134:   nkdm->jacobianctx  = kdm->jacobianctx;
135:   nkdm->objectivectx = kdm->objectivectx;
136:   nkdm->originaldm   = kdm->originaldm;

138:   /*
139:   nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
140:   nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
141:   nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
142:   */

144:   /* implementation specific copy hooks */
145:   if (kdm->ops->duplicate) (*kdm->ops->duplicate)(kdm,nkdm);
146:   return 0;
147: }

149: /*@C
150:    DMGetDMSNES - get read-only private DMSNES context from a DM

152:    Not Collective

154:    Input Parameter:
155: .  dm - DM to be used with SNES

157:    Output Parameter:
158: .  snesdm - private DMSNES context

160:    Level: developer

162:    Notes:
163:    Use DMGetDMSNESWrite() if write access is needed. The DMSNESSetXXX API should be used wherever possible.

165: .seealso: DMGetDMSNESWrite()
166: @*/
167: PetscErrorCode DMGetDMSNES(DM dm,DMSNES *snesdm)
168: {
170:   *snesdm = (DMSNES) dm->dmsnes;
171:   if (!*snesdm) {
172:     PetscInfo(dm,"Creating new DMSNES\n");
173:     DMSNESCreate(PetscObjectComm((PetscObject)dm),snesdm);

175:     dm->dmsnes            = (PetscObject) *snesdm;
176:     (*snesdm)->originaldm = dm;
177:     DMCoarsenHookAdd(dm,DMCoarsenHook_DMSNES,DMRestrictHook_DMSNES,NULL);
178:     DMRefineHookAdd(dm,DMRefineHook_DMSNES,DMInterpolateHook_DMSNES,NULL);
179:     DMSubDomainHookAdd(dm,DMSubDomainHook_DMSNES,DMSubDomainRestrictHook_DMSNES,NULL);
180:   }
181:   return 0;
182: }

184: /*@C
185:    DMGetDMSNESWrite - get write access to private DMSNES context from a DM

187:    Not Collective

189:    Input Parameter:
190: .  dm - DM to be used with SNES

192:    Output Parameter:
193: .  snesdm - private DMSNES context

195:    Level: developer

197: .seealso: DMGetDMSNES()
198: @*/
199: PetscErrorCode DMGetDMSNESWrite(DM dm,DMSNES *snesdm)
200: {
201:   DMSNES         sdm;

204:   DMGetDMSNES(dm,&sdm);
206:   if (sdm->originaldm != dm) {  /* Copy on write */
207:     DMSNES oldsdm = sdm;
208:     PetscInfo(dm,"Copying DMSNES due to write\n");
209:     DMSNESCreate(PetscObjectComm((PetscObject)dm),&sdm);
210:     DMSNESCopy(oldsdm,sdm);
211:     DMSNESDestroy((DMSNES*)&dm->dmsnes);
212:     dm->dmsnes = (PetscObject)sdm;
213:     sdm->originaldm = dm;
214:   }
215:   *snesdm = sdm;
216:   return 0;
217: }

219: /*@C
220:    DMCopyDMSNES - copies a DM context to a new DM

222:    Logically Collective

224:    Input Parameters:
225: +  dmsrc - DM to obtain context from
226: -  dmdest - DM to add context to

228:    Level: developer

230:    Note:
231:    The context is copied by reference. This function does not ensure that a context exists.

233: .seealso: DMGetDMSNES(), SNESSetDM()
234: @*/
235: PetscErrorCode DMCopyDMSNES(DM dmsrc,DM dmdest)
236: {
239:   if (!dmdest->dmsnes) DMSNESCreate(PetscObjectComm((PetscObject) dmdest), (DMSNES *) &dmdest->dmsnes);
240:   DMSNESCopy((DMSNES) dmsrc->dmsnes, (DMSNES) dmdest->dmsnes);
241:   DMCoarsenHookAdd(dmdest,DMCoarsenHook_DMSNES,NULL,NULL);
242:   DMRefineHookAdd(dmdest,DMRefineHook_DMSNES,NULL,NULL);
243:   DMSubDomainHookAdd(dmdest,DMSubDomainHook_DMSNES,DMSubDomainRestrictHook_DMSNES,NULL);
244:   return 0;
245: }

247: /*@C
248:    DMSNESSetFunction - set SNES residual evaluation function

250:    Not Collective

252:    Input Parameters:
253: +  dm - DM to be used with SNES
254: .  f - residual evaluation function; see SNESFunction for details
255: -  ctx - context for residual evaluation

257:    Level: advanced

259:    Note:
260:    SNESSetFunction() is normally used, but it calls this function internally because the user context is actually
261:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
262:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

264: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESFunction
265: @*/
266: PetscErrorCode DMSNESSetFunction(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
267: {
268:   DMSNES         sdm;

271:   if (f || ctx) {
272:     DMGetDMSNESWrite(dm,&sdm);
273:   }
274:   if (f) sdm->ops->computefunction = f;
275:   if (ctx) sdm->functionctx = ctx;
276:   return 0;
277: }

279: /*@C
280:    DMSNESSetMFFunction - set SNES residual evaluation function used in applying the matrix-free Jacobian with -snes_mf_operator

282:    Logically Collective on dm

284:    Input Parameters:
285: +  dm - DM to be used with SNES
286: -  f - residual evaluation function; see SNESFunction for details

288:    Level: advanced

290: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESFunction, DMSNESSetFunction()
291: @*/
292: PetscErrorCode DMSNESSetMFFunction(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
293: {
294:   DMSNES         sdm;

297:   if (f || ctx) {
298:     DMGetDMSNESWrite(dm,&sdm);
299:   }
300:   if (f) sdm->ops->computemffunction = f;
301:   if (ctx) sdm->mffunctionctx = ctx;
302:   return 0;
303: }

305: /*@C
306:    DMSNESGetFunction - get SNES residual evaluation function

308:    Not Collective

310:    Input Parameter:
311: .  dm - DM to be used with SNES

313:    Output Parameters:
314: +  f - residual evaluation function; see SNESFunction for details
315: -  ctx - context for residual evaluation

317:    Level: advanced

319:    Note:
320:    SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
321:    associated with the DM.

323: .seealso: DMSNESSetContext(), DMSNESSetFunction(), SNESSetFunction(), SNESFunction
324: @*/
325: PetscErrorCode DMSNESGetFunction(DM dm,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx)
326: {
327:   DMSNES         sdm;

330:   DMGetDMSNES(dm,&sdm);
331:   if (f) *f = sdm->ops->computefunction;
332:   if (ctx) *ctx = sdm->functionctx;
333:   return 0;
334: }

336: /*@C
337:    DMSNESSetObjective - set SNES objective evaluation function

339:    Not Collective

341:    Input Parameters:
342: +  dm - DM to be used with SNES
343: .  obj - objective evaluation function; see SNESObjectiveFunction for details
344: -  ctx - context for residual evaluation

346:    Level: advanced

348: .seealso: DMSNESSetContext(), SNESGetObjective(), DMSNESSetFunction()
349: @*/
350: PetscErrorCode DMSNESSetObjective(DM dm,PetscErrorCode (*obj)(SNES,Vec,PetscReal*,void*),void *ctx)
351: {
352:   DMSNES         sdm;

355:   if (obj || ctx) {
356:     DMGetDMSNESWrite(dm,&sdm);
357:   }
358:   if (obj) sdm->ops->computeobjective = obj;
359:   if (ctx) sdm->objectivectx = ctx;
360:   return 0;
361: }

363: /*@C
364:    DMSNESGetObjective - get SNES objective evaluation function

366:    Not Collective

368:    Input Parameter:
369: .  dm - DM to be used with SNES

371:    Output Parameters:
372: +  obj- residual evaluation function; see SNESObjectiveFunction for details
373: -  ctx - context for residual evaluation

375:    Level: advanced

377:    Note:
378:    SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
379:    associated with the DM.

381: .seealso: DMSNESSetContext(), DMSNESSetObjective(), SNESSetFunction()
382: @*/
383: PetscErrorCode DMSNESGetObjective(DM dm,PetscErrorCode (**obj)(SNES,Vec,PetscReal*,void*),void **ctx)
384: {
385:   DMSNES         sdm;

388:   DMGetDMSNES(dm,&sdm);
389:   if (obj) *obj = sdm->ops->computeobjective;
390:   if (ctx) *ctx = sdm->objectivectx;
391:   return 0;
392: }

394: /*@C
395:    DMSNESSetNGS - set SNES Gauss-Seidel relaxation function

397:    Not Collective

399:    Input Parameters:
400: +  dm - DM to be used with SNES
401: .  f  - relaxation function, see SNESGSFunction
402: -  ctx - context for residual evaluation

404:    Level: advanced

406:    Note:
407:    SNESSetNGS() is normally used, but it calls this function internally because the user context is actually
408:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
409:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

411: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), DMSNESSetFunction(), SNESGSFunction
412: @*/
413: PetscErrorCode DMSNESSetNGS(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
414: {
415:   DMSNES         sdm;

418:   if (f || ctx) {
419:     DMGetDMSNESWrite(dm,&sdm);
420:   }
421:   if (f) sdm->ops->computegs = f;
422:   if (ctx) sdm->gsctx = ctx;
423:   return 0;
424: }

426: /*@C
427:    DMSNESGetNGS - get SNES Gauss-Seidel relaxation function

429:    Not Collective

431:    Input Parameter:
432: .  dm - DM to be used with SNES

434:    Output Parameters:
435: +  f - relaxation function which performs Gauss-Seidel sweeps, see SNESGSFunction
436: -  ctx - context for residual evaluation

438:    Level: advanced

440:    Note:
441:    SNESGetNGS() is normally used, but it calls this function internally because the user context is actually
442:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
443:    not. If DM took a more central role at some later date, this could become the primary method of setting the residual.

445: .seealso: DMSNESSetContext(), SNESGetNGS(), DMSNESGetJacobian(), DMSNESGetFunction(), SNESNGSFunction
446: @*/
447: PetscErrorCode DMSNESGetNGS(DM dm,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx)
448: {
449:   DMSNES         sdm;

452:   DMGetDMSNES(dm,&sdm);
453:   if (f) *f = sdm->ops->computegs;
454:   if (ctx) *ctx = sdm->gsctx;
455:   return 0;
456: }

458: /*@C
459:    DMSNESSetJacobian - set SNES Jacobian evaluation function

461:    Not Collective

463:    Input Parameters:
464: +  dm - DM to be used with SNES
465: .  J - Jacobian evaluation function
466: -  ctx - context for residual evaluation

468:    Level: advanced

470:    Note:
471:    SNESSetJacobian() is normally used, but it calls this function internally because the user context is actually
472:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
473:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.

475: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESGetJacobian(), SNESSetJacobian(), SNESJacobianFunction
476: @*/
477: PetscErrorCode DMSNESSetJacobian(DM dm,PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
478: {
479:   DMSNES         sdm;

482:   if (J || ctx) {
483:     DMGetDMSNESWrite(dm,&sdm);
484:   }
485:   if (J) sdm->ops->computejacobian = J;
486:   if (ctx) sdm->jacobianctx = ctx;
487:   return 0;
488: }

490: /*@C
491:    DMSNESGetJacobian - get SNES Jacobian evaluation function

493:    Not Collective

495:    Input Parameter:
496: .  dm - DM to be used with SNES

498:    Output Parameters:
499: +  J - Jacobian evaluation function; see SNESJacobianFunction for all calling sequence
500: -  ctx - context for residual evaluation

502:    Level: advanced

504:    Note:
505:    SNESGetJacobian() is normally used, but it calls this function internally because the user context is actually
506:    associated with the DM.  This makes the interface consistent regardless of whether the user interacts with a DM or
507:    not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.

509: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESJacobianFunction
510: @*/
511: PetscErrorCode DMSNESGetJacobian(DM dm,PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx)
512: {
513:   DMSNES         sdm;

516:   DMGetDMSNES(dm,&sdm);
517:   if (J) *J = sdm->ops->computejacobian;
518:   if (ctx) *ctx = sdm->jacobianctx;
519:   return 0;
520: }

522: /*@C
523:    DMSNESSetPicard - set SNES Picard iteration matrix and RHS evaluation functions.

525:    Not Collective

527:    Input Parameters:
528: +  dm - DM to be used with SNES
529: .  b - RHS evaluation function
530: .  J - Picard matrix evaluation function
531: -  ctx - context for residual evaluation

533:    Level: advanced

535: .seealso: SNESSetPicard(), DMSNESSetFunction(), DMSNESSetJacobian()
536: @*/
537: PetscErrorCode DMSNESSetPicard(DM dm,PetscErrorCode (*b)(SNES,Vec,Vec,void*),PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
538: {
539:   DMSNES         sdm;

542:   DMGetDMSNES(dm,&sdm);
543:   if (b) sdm->ops->computepfunction = b;
544:   if (J) sdm->ops->computepjacobian = J;
545:   if (ctx) sdm->pctx = ctx;
546:   return 0;
547: }

549: /*@C
550:    DMSNESGetPicard - get SNES Picard iteration evaluation functions

552:    Not Collective

554:    Input Parameter:
555: .  dm - DM to be used with SNES

557:    Output Parameters:
558: +  b - RHS evaluation function; see SNESFunction for details
559: .  J  - RHS evaluation function; see SNESJacobianFunction for detailsa
560: -  ctx - context for residual evaluation

562:    Level: advanced

564: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian()
565: @*/
566: PetscErrorCode DMSNESGetPicard(DM dm,PetscErrorCode (**b)(SNES,Vec,Vec,void*),PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx)
567: {
568:   DMSNES         sdm;

571:   DMGetDMSNES(dm,&sdm);
572:   if (b) *b = sdm->ops->computepfunction;
573:   if (J) *J = sdm->ops->computepjacobian;
574:   if (ctx) *ctx = sdm->pctx;
575:   return 0;
576: }