Actual source code: gsnes.c

  1: #ifdef PETSC_RCS_HEADER
  2: static char vcid[] = "$Id: gsnes.c,v 1.11 2000/07/16 05:40:27 knepley Exp $";
  3: #endif

 5:  #include src/snes/snesimpl.h
 6:  #include src/gvec/gvecimpl.h
 7:  #include gsolver.h

  9: /*
 10:   This file provides routines for grid SNES objects. They support solution of
 11:   systems of nonlinear equations generated from a mesh and discretization
 12:   specified by a Grid object.
 13: */

 15: #undef  __FUNCT__
 17: int GSNESCreateStructures_Private(GSNES snes, void *funcCtx, void *jacCtx)
 18: {
 19:   Grid       grid;
 20:   GMat       J;
 21:   GVec       f;
 22:   PetscTruth isMatrixFree;
 23:   int        ierr;

 26:   GSNESGetGrid(snes, &grid);
 27:   GVecCreateConstrained(grid, &f);
 28:   GridIsMatrixFree(grid, &isMatrixFree);
 29:   if (isMatrixFree == PETSC_TRUE) {
 30:     if (grid->explicitConstraints == PETSC_TRUE) {
 31:       GMatCreateMF(grid, grid->constraintOrder, grid->constraintOrder, &J);
 32:     } else {
 33:       GMatCreateMF(grid, grid->order, grid->order, &J);
 34:     }
 35:   } else {
 36:     if (grid->explicitConstraints == PETSC_TRUE) {
 37:       GMatCreateRectangular(grid, grid->constraintOrder, grid->constraintOrder, &J);
 38:     } else {
 39:       GMatCreate(grid, &J);
 40:     }
 41:   }
 42:   SNESSetFunction(snes, f, (int (*)(SNES, Vec, Vec, void *)) GSNESEvaluateRhs, funcCtx);
 43:   SNESSetJacobian(snes, J, J, (int (*)(SNES, Vec, Mat *, Mat *, MatStructure *, void *)) GSNESEvaluateJacobian, jacCtx);
 44:   PetscLogObjectParent(grid, f);
 45:   PetscLogObjectParent(grid, J);
 46:   return(0);
 47: }

 49: #undef  __FUNCT__
 51: int GSNESDestroyStructures_Private(GSNES snes)
 52: {
 53:   GMat J;
 54:   GVec f;
 55:   int  ierr;

 58:   SNESGetFunction(snes, &f, PETSC_NULL, PETSC_NULL);
 59:   SNESGetJacobian(snes, &J, PETSC_NULL, PETSC_NULL, PETSC_NULL);
 60:   GMatDestroy(J);
 61:   GVecDestroy(f);
 62:   return(0);
 63: }

 65: #undef  __FUNCT__
 67: /*@C
 68:   GSNESDestroy - Destroys a grid SNES.

 70:   Collective on GSNES

 72:   Input Parameter:
 73: . snes - The nonlinear solver context

 75:   Level: beginner

 77: .keywords: grid SNES, destroy
 78: .seealso: SNESDestroy(), GSNESDuplicate()
 79: @*/
 80: int GSNESDestroy(GSNES snes)
 81: {

 86:   if (--snes->refct > 0)
 87:     return(0);
 88:   GSNESDestroyStructures_Private(snes);
 89:   SNESDestroy(snes);
 90:   return(0);
 91: }

 93: #undef  __FUNCT__
 95: /*@ 
 96:   GSNESView - Views a grid SNES.

 98:   Collective on GSNES

100:   Input Parameters:
101: + snes   - The grid SNES
102: - viewer - [Optional] A visualization context

104:   Options Database Key:
105: $ -snes_view : calls SNESView() at end of SNESSolve()

107:   Notes:
108:   The available visualization contexts include
109: $   VIEWER_STDOUT_SELF - standard output (default)
110: $   VIEWER_STDOUT_WORLD - synchronized standard
111: $     output where only the first processor opens
112: $     the file.  All other processors send their 
113: $     data to the first processor to print. 

115:   The user can open alternative vistualization contexts with
116: $   PetscViewerFileOpenASCII() - output vector to a specified file

118:   Level: beginner

120: .keywords: grid SNES, view, visualize
121: .seealso: SNESView()
122: @*/
123: int GSNESView(GSNES snes, PetscViewer viewer)
124: {
125:   Grid       grid;
126:   FILE      *fd;
127:   int        numActiveFields;
128:   int        field, func, op;
129:   PetscTruth isascii;
130:   int        ierr;

134:   if (!viewer) {
135:     viewer = PETSC_VIEWER_STDOUT_SELF;
136:   } else {
138:   }

140:   GSNESGetGrid(snes, &grid);
141:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &isascii);
142:   if (isascii == PETSC_TRUE) {
143:     GridView(grid, viewer);
144:     PetscViewerFlush(viewer);
145:     PetscViewerASCIIGetPointer(viewer, &fd);
146:     PetscFPrintf(grid->comm, fd, "GSNES Object:n");
147:     GridGetNumActiveFields(grid, &numActiveFields);
148:     if (numActiveFields == 1) {
149:       PetscFPrintf(snes->comm, fd, "  %d active field:n", numActiveFields);
150:     } else {
151:       PetscFPrintf(snes->comm, fd, "  %d active fields:n", numActiveFields);
152:     }
153:     for(field = 0; field < grid->numFields; field++) {
154:       if (grid->fields[field].isActive == PETSC_FALSE) continue;
155:       if (grid->fields[field].name) {
156:         PetscFPrintf(grid->comm, fd, "  %s field with %d components:n", grid->fields[field].name, grid->fields[field].numComp);
157:       } else {
158:         PetscFPrintf(grid->comm, fd, "  field %d with %d components:n", field, grid->fields[field].numComp);
159:       }
160:       for(func = 0; func < grid->numRhsFuncs; func++) {
161:         if (grid->rhsFuncs[func].field == field) {
162:           PetscFPrintf(grid->comm, fd, "    Rhs function with scalar multiplier %gn", PetscRealPart(grid->rhsFuncs[func].alpha));
163:         }
164:       }
165:       for(op = 0; op < grid->numRhsOps; op++) {
166:         if (grid->rhsOps[op].field == field) {
167:           if (grid->rhsOps[op].nonlinearOp != PETSC_NULL) {
168:             PetscFPrintf(grid->comm, fd, "    Rhs nonlinear operator with scalar multiplier %gn",
169:                          PetscRealPart(grid->rhsOps[field].alpha));
170:           } else {
171:             PetscFPrintf(grid->comm, fd, "    Rhs operator %d with scalar multiplier %gn",
172:                          grid->rhsOps[op].op, PetscRealPart(grid->rhsOps[op].alpha));
173:           }
174:         }
175:       }
176:       for(op = 0; op < grid->numMatOps; op++) {
177:         if (grid->matOps[op].field == field) {
178:           PetscFPrintf(grid->comm, fd, "    Jacobian operator %d with scalar multiplier %gn",
179:                        grid->matOps[op].op, PetscRealPart(grid->matOps[op].alpha));
180:         }
181:       }
182:     }
183:     SNESView(snes, viewer);
184:   }

186:   return(0);
187: }

189: #undef  __FUNCT__
191: /*@C
192:   GSNESDuplicate - Duplicates a grid SNES.

194:   Collective on GSNES

196:   Input Parameter:
197: . snes    - The GSNES

199:   Output Parameter:
200: . newsnes - The new GSNES

202:   Level: beginner

204: .keywords: grid SNES, duplicate
205: .seealso: SNESDuplicate(), GSNESDestroy()
206: @*/
207: int GSNESDuplicate(GSNES snes, GSNES *newsnes)
208: {
212:   SETERRQ(PETSC_ERR_SUP, " ");
213: #if 0
214:   return(0);
215: #endif
216: }

218: #undef  __FUNCT__
220: /*@C
221:   GSNESDuplicateMonitors - Duplicates the monitors of a grid SNES in another.

223:   Collective on GSNES

225:   Input Parameter:
226: . snes    - The GSNES to copy from

228:   Output Parameter:
229: . newsnes - The GSNES to copy to

231:   Level: intermediate

233: .keywords: grid snes, duplicate, monitor
234: .seealso: GSNESDuplicate(), GSNESDestroy()
235: @*/
236: int GSNESDuplicateMonitors(GSNES snes, GSNES newsnes)
237: {
238:   int mon, mon2;


245:   for(mon = 0; mon < snes->numbermonitors; mon++)
246:   {
247:     for(mon2 = 0; mon2 < newsnes->numbermonitors; mon2++)
248:       if (newsnes->monitor[mon] == snes->monitor[mon])
249:         break;
250:     if (mon2 == newsnes->numbermonitors) {
251:       SNESSetMonitor(newsnes, snes->monitor[mon], snes->monitorcontext[mon], PETSC_NULL);
252:     }
253:   }
254:   return(0);
255: }

257: #undef  __FUNCT__
259: /*@C
260:   GSNESReallocate - Reallocates the storage for a grid SNES.

262:   Collective on GSNES

264:   Input Parameter:
265: . snes - The GSNES

267:   Level: advanced

269: .keywords: grid snes, reallocate
270: .seealso: GridReform(), GSNESDuplicate(), GSNESDestroy()
271: @*/
272: int GSNESReallocate(GSNES snes)
273: {
274:   void *funcCtx;
275:   void *jacCtx;
276:   int   ierr;

280:   /* Preserve contexts */
281:   SNESGetFunction(snes, PETSC_NULL, &funcCtx, PETSC_NULL);
282:   SNESGetJacobian(snes, PETSC_NULL, PETSC_NULL, &jacCtx, PETSC_NULL);
283:   GSNESDestroyStructures_Private(snes);
284:   GSNESCreateStructures_Private(snes, funcCtx, jacCtx);
285:   return(0);
286: }

288: #undef  __FUNCT__
290: /*@
291:   GSNESGetGrid - Gets the grid from a grid SNES.

293:   Not collective

295:   Input Parameter:
296: . snes - The grid SNES

298:   Output Parameter:
299: . grid - The grid context

301:   Level: intermediate

303: .keywords: grid SNES, grid, get
304: .seealso: GVecGetGrid(), GMatGetGrid()
305: @*/
306: int GSNESGetGrid(GSNES snes, Grid *grid)
307: {


314:   PetscObjectQuery((PetscObject) snes, "Grid", (PetscObject *) grid);
316:   return(0);
317: }

319: #undef  __FUNCT__
321: /*@
322:   GSNESEvaluateRhs - Constructs the Rhs vector for Newton's equation.

324:   Collective on GSNES

326:   Input Parameters:
327: + snes - The grid SNES
328: . x    - The current iterate
329: - ctx  - The optional user context

331:   Output Parameter:
332: . f    - The function value

334:   Note:
335:   This function initializes the vector f to zero before calculation.

337:   Level: advanced

339: .keywords: grid SNES, rhs
340: .seealso: GridEvaluateRhs()
341: @*/
342: int GSNESEvaluateRhs(GSNES snes, GVec x, GVec f, PetscObject ctx)
343: {
344:   Grid        grid;
345:   PetscScalar zero = 0.0;
346:   int         ierr;

350:   GSNESGetGrid(snes, &grid);
351:   /* Initialize vector */
352:   VecSet(&zero, f);
353:   /* Form function */
354:   GridEvaluateRhs(grid, x, f, ctx);
355:   return(0);
356: }

358: #undef  __FUNCT__
360: /*@
361:   GSNESEvaluateRhsFunction - Constructs only the constant portion of the nonlinear equation,
362:   that is the portion arising from functions independent of the solution variable x.

364:   Collective on GSNES

366:   Input Parameters:
367: + snes - The grid SNES
368: . x    - The current iterate
369: - ctx  - The optional user context

371:   Output Parameter:
372: . f    - The function value

374:   Notes:
375:   This function actually constructs -f since f is usually written on the right, and the user
376:   cannot just multiply by -1 after boundary conditions are imposed.

378:   Level: advanced

380: .keywords: grid SNES, rhs, function
381: .seealso: GSNESEvaluateRhsOperator(), GSNESEvaluateRhs(), GSNESEvaluateJacobian()
382: @*/
383: int GSNESEvaluateRhsFunction(GSNES snes, GVec x, GVec f, void *ctx)
384: {
385:   Grid        grid;
386:   PetscScalar zero = 0.0;
387:   int         ierr;

391:   GSNESGetGrid(snes, &grid);
392:   /* Initialize vector */
393:   VecSet(&zero, f);
394:   /* Form function */
395:   GridScaleRhs(grid, -1.0);
396:   GridEvaluateRhsFunction(grid, x, f, ctx);
397:   GridScaleRhs(grid, -1.0);
398:   return(0);
399: }

401: #undef  __FUNCT__
403: /*@
404:   GSNESEvaluateRhsOperator - Constructs only the portion of the nonlinear equation dependent of the solution variable x.

406:   Collective on GSNES

408:   Input Parameters:
409: + snes - The grid SNES
410: . x    - The current iterate
411: - ctx  - The optional user context

413:   Output Parameter:
414: . f    - The function value

416:   Level: advanced

418: .keywords: grid SNES, rhs, operator
419: .seealso: GSNESEvaluateRhsLinearOperator(), GSNESEvaluateRhsNonlinearOperator(), GSNESEvaluateRhsFunction(),
420:           GSNESEvaluateRhs(), GSNESEvaluateJacobian()
421: @*/
422: int GSNESEvaluateRhsOperator(GSNES snes, GVec x, GVec f, void *ctx)
423: {
424:   Grid        grid;
425:   PetscTruth  reduceElement;
426:   PetscScalar zero = 0.0;
427:   int         ierr;

431:   GSNESGetGrid(snes, &grid);
432:   /* Initialize vector */
433:   VecSet(&zero, f);
434:   /* Turn off boundary values */
435:   GridGetReduceElement(grid, &reduceElement);
436:   GridSetReduceElement(grid, PETSC_FALSE);
437:   /* Apply operator */
438:   GridEvaluateRhsOperator(grid, x, f, PETSC_TRUE, PETSC_TRUE, ctx);
439:   /* Restore boundary values */
440:   GridSetReduceElement(grid, reduceElement);
441:   return(0);
442: }

444: #undef  __FUNCT__
446: /*@
447:   GSNESEvaluateRhsLinearOperator - Constructs only the portion of the nonlinear equation
448:   linearly dependent on the solution variable x.

450:   Collective on GSNES

452:   Input Parameters:
453: + snes - The grid SNES
454: . x    - The current iterate
455: - ctx  - The optional user context

457:   Output Parameter:
458: . f    - The function value

460:   Level: advanced

462: .keywords: grid SNES, rhs, operator, linear
463: .seealso: GSNESEvaluateRhsOperator(), GSNESEvaluateRhsFunction(), GSNESEvaluateRhs(), GSNESEvaluateJacobian()
464: @*/
465: int GSNESEvaluateRhsLinearOperator(GSNES snes, GVec x, GVec f, void *ctx)
466: {
467:   Grid        grid;
468:   PetscTruth  reduceElement;
469:   PetscScalar zero = 0.0;
470:   int         ierr;

474:   GSNESGetGrid(snes, &grid);
475:   /* Initialize vector */
476:   VecSet(&zero, f);
477:   /* Turn off boundary values */
478:   GridGetReduceElement(grid, &reduceElement);
479:   GridSetReduceElement(grid, PETSC_FALSE);
480:   /* Apply operator */
481:   GridEvaluateRhsOperator(grid, x, f, PETSC_TRUE, PETSC_FALSE, ctx);
482:   /* Restore boundary values */
483:   GridSetReduceElement(grid, reduceElement);
484:   return(0);
485: }

487: #undef  __FUNCT__
489: /*@
490:   GSNESEvaluateRhsNonlinearOperator - Constructs only the portion of the nonlinear equation
491:   nonlinearly dependent on the solution variable x.

493:   Collective on GSNES

495:   Input Parameters:
496: + snes - The grid SNES
497: . n    - The number of input vectors
498: . vecs - The input vectors
499: - ctx  - The optional user context

501:   Output Parameter:
502: . f    - The function value

504:   Level: advanced

506:   Note:
507:   If there is only one input argument, it will be used as every argument for
508:   multilinear functions.

510: .keywords: grid SNES, rhs, operator, nonlinear
511: .seealso: GSNESEvaluateRhsOperator(), GSNESEvaluateRhsFunction(), GSNESEvaluateRhs(), GSNESEvaluateJacobian()
512: @*/
513: int GSNESEvaluateRhsNonlinearOperator(GSNES snes, int n, GVec vecs[], GVec f, void *ctx)
514: {
515:   Grid              grid;
516:   NonlinearOperator op;
517:   int               field;
518:   PetscScalar       alpha;
519:   PetscTruth        isALE;
520:   PetscTruth        reduceElement;
521:   PetscScalar       zero = 0.0;
522:   int               ierr;

526:   GSNESGetGrid(snes, &grid);
527:   /* Initialize vector */
528:   VecSet(&zero, f);
529:   /* Turn off boundary values */
530:   GridGetReduceElement(grid, &reduceElement);
531:   GridSetReduceElement(grid, PETSC_FALSE);
532:   /* Apply operator */
533:   if (n == 1) {
534:     GridEvaluateRhsOperator(grid, vecs[0], f, PETSC_FALSE, PETSC_TRUE, ctx);
535:   } else {
536:     GridGetNonlinearOperatorStart(grid, &op, &field, &alpha, &isALE);
537:     while(op != PETSC_NULL) {
538:       GVecEvaluateNonlinearOperatorGalerkin(f, n, vecs, 1, &field, op, alpha, isALE, ctx);
539:       GridGetNonlinearOperatorNext(grid, &op, &field, &alpha, &isALE);
540:     }
541:   }
542:   /* Restore boundary values */
543:   GridSetReduceElement(grid, reduceElement);
544:   return(0);
545: }

547: #undef  __FUNCT__
549: /*@
550:   GSNESEvaluateJacobian - Constructs the system matrix for Newton's equation

552:   Collective on GSNES

554:   Input Parameters:
555: + snes - The grid SNES
556: . x    - The current iterate
557: . flag - The matrix structure flag
558: - ctx  - The optional user context

560:   Output Parameters:
561: + J    - The Jacobian matrix
562: - M    - The preconditioner for the Jacobian matrix, usually the same as J

564:   Level: advanced

566: .keywords: grid SNES, jacobian
567: .seealso: GSNESEvaluateJacobianMF(), GSNESEvaluateRhs(), GridEvaluateSystemMatrix()
568: @*/
569: int GSNESEvaluateJacobian(GSNES snes, GVec x, GMat *J, GMat *M, MatStructure *flag, PetscObject ctx)
570: {
571:   Grid grid;
572:   int  ierr;

576:   GSNESGetGrid(snes, &grid);
577:   /* Initialize matrix */
578:   MatZeroEntries(*J);
579:   /* Form function */
580:   GridEvaluateSystemMatrix(grid, x, J, M, flag, ctx);
581:   /* Apply boundary conditions */
582:   GMatSetBoundary(*J, 1.0, ctx);
583:   return(0);
584: }

586: #undef  __FUNCT__
588: /*@
589:   GSNESEvaluateJacobianMF - Constructs the application of the system matrix for Newton's equation

591:   Collective on GSNES

593:   Input Parameters:
594: + snes - The grid SNES
595: . x    - The current iterate
596: . y    - The vector to which to apply the Jacobian
597: - ctx  - The optional user context

599:   Output Parameter:
600: . f    - The vector J(x) y

602:   Level: advanced

604: .keywords: grid SNES< jacobian, mf
605: .seealso: GSNESEvaluateJacobian(), GSNESEvaluateRhs(), GridEvaluateSystemMatrix(), GVecEvaluateJacobian()
606: @*/
607: int GSNESEvaluateJacobianMF(GSNES snes, GVec x, GVec y, GVec f, void *ctx)
608: {
609:   Grid        grid;
610:   PetscScalar zero = 0.0;
611:   int         ierr;

615:   GSNESGetGrid(snes, &grid);
616:   /* Initialize vector  */
617:   VecSet(&zero, f);
618:   /* Form function */
619:   GVecEvaluateJacobian(f, x, y, ctx);
620:   return(0);
621: }

623: #undef __FUNCT__  
625: /*@
626:   GSNESRhsBC - This sets zero boundary conditions for the Rhs.

628:   Collective on GSNES

630:   Input Parameters:
631: + snes - The SNES
632: - ctx  - The user-supplied context

634:   Output Parameter:
635: . rhs  - The modified Rhs

637:   Level: advanced

639: .keywords: grid SNES, set, boundary conditions
640: .seealso GSNESSetSolutionBC()
641: @*/
642: int GSNESRhsBC(GSNES snes, GVec rhs, void *ctx)
643: {

649:   GVecSetBoundaryZero(rhs, ctx);
650:   return(0);
651: }

653: #undef __FUNCT__  
655: /*@
656:   GSNESSolutionBC - This sets the boundary conditions registered with the grid
657:   on the solution vector.

659:   Collective on GSNES

661:   Input Parameters:
662: + snes - The SNES
663: - ctx  - The user-supplied context

665:   Output Parameter:
666: . sol  - The modified solution

668:   Level: advanced

670: .keywords: GSNES, set, boundary conditions
671: .seealso: GSNESSetRhsBC()
672: @*/
673: int GSNESSolutionBC(GSNES snes, GVec sol, void *ctx)
674: {

680:   GVecSetBoundary(sol, ctx);
681:   return(0);
682: }

684: #undef __FUNCT__  
686: /*@
687:   GSNESUpdate - This is a general purpose function which is called at the beginning of every iterate.

689:   Collective on GSNES

691:   Input Parameters:
692: + snes - The nonlinear solver context
693: - step - The current step

695:   Level: advanced

697: .keywords: GSNES, update
698: .seealso: SNESSolve()
699: @*/
700: int GSNESUpdate(GSNES snes, int step)
701: {
702:   Grid       grid;
703:   Mesh       mesh;
704:   MeshMover  mover;
705:   PetscTruth isMoving;
706:   int        ierr;

709:   /* Calculate mesh velocity */
710:   GSNESGetGrid(snes, &grid);
711:   GridGetMesh(grid, &mesh);
712:   MeshGetMovement(mesh, &isMoving);
713:   if (isMoving == PETSC_TRUE) {
714:     MeshGetMover(mesh, &mover);
715:     MeshMoverCalcNodeVelocities(mover, PETSC_TRUE);
716:   }
717:   return(0);
718: }

720: #undef  __FUNCT__
722: /*@
723:   GSNESCreate - Creates a grid SNES associated with a particular discretization context.

725:   Collective on Grid

727:   Input Parameter:
728: + grid - The discretization context
729: - ctx  - The optional user context

731:   Output Parameter:
732: . snes - The nonlinear solver context

734:   Level: beginner

736: .keywords: GSNES, create
737: .seealso: GSNESDestroy()
738: @*/
739: int GSNESCreate(Grid grid, void *ctx, GSNES *snes)
740: {

746: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
747:   GSolverInitializePackage(PETSC_NULL);
748: #endif

750:   /* Create the SNES */
751:   SNESCreate(grid->comm,snes);
752:   PetscObjectCompose((PetscObject) *snes, "Grid", (PetscObject) grid);
753:   (*snes)->isGSNES = PETSC_TRUE;

755:   /* Set boundary conditions */
756:   SNESSetRhsBC(*snes, GSNESRhsBC);
757:   SNESSetSolutionBC(*snes, GSNESSolutionBC);
758:   if (grid->reduceSystem == PETSC_TRUE) {
759:     /* Save space and time by reducing the system right in construction rather than using an explicit matvec and subtraction */
760:     grid->reduceElement = PETSC_TRUE;
761:     /* Since we solve J x = -F, we need a BC multiplier of -1 */
762:     GridSetBCMultiplier(grid, -1.0);
763:   }

765:   /* Set update function */
766:   SNESSetUpdate(*snes, GSNESUpdate);

768:   /* Setup SNES data structures */
769:   GSNESCreateStructures_Private(*snes, ctx, ctx);
770:   return(0);
771: }

773: #undef  __FUNCT__
775: int GSNESDefaultComputeJacobianWithColoring(SNES snes, GVec x, GMat *J, GMat *M, MatStructure *flag, void *ctx)
776: {
777:   SETERRQ(PETSC_ERR_SUP, " ");
778: #if 0
779:   Grid grid;
780:   int  ierr;

787:   GMatGetGrid(*M, &grid);
789:   if (!grid->fdcoloring) {
790:     GMatFDColoringCreate(*M);
791:   }
792:   SNESDefaultComputeJacobianWithColoring(snes, x, J, M, flag, grid->fdcoloring);
793:   return(0);
794: #endif
795: }

797: #undef  __FUNCT__
799: /*@
800:   GSNESSolutionMonitor - Monitors solution at each GSNES iteration.

802:   Collective on GSNES

804:   Input Parameters:
805: + snes  - The nonlinear solver context
806: . it    - The number of iterations so far
807: . rnorm - The current (approximate) residual norm
808: - ctx   - The viewer

810:   Level: advanced

812: .keywords: grid SNES, monitor, solution
813: .seealso: GSNESResidualMonitor(), GSNESErrorMonitor(), SNESDefaultMonitor(), SNESSetMonitor()
814: @*/
815: int GSNESSolutionMonitor(GSNES snes, int it, PetscReal rnorm, void *ctx)
816: {
817:   PetscViewer viewer = (PetscViewer) ctx;
818:   Vec         x;
819:   int         ierr;

823:   SNESGetSolution(snes, &x);
824:   VecView(x, viewer);
825:   return(0);
826: }

828: #undef  __FUNCT__
830: /*@
831:   GSNESResidualMonitor - Monitors residual at each SNES iteration.

833:   Collective on GSNES

835:   Input Parameters:
836: + snes  - The nonlinear solver context
837: . it    - The number of iterations so far
838: . rnorm - The current (approximate) residual norm
839: - ctx   - The viewer

841:   Level: advanced

843: .keywords: grid SNES, monitor, residual
844: .seealso: GSNESSolutionMonitor(), GSNESErrorMonitor(), SNESDefaultMonitor(), SNESSetMonitor()
845: @*/
846: int GSNESResidualMonitor(GSNES snes, int it, PetscReal rnorm, void *ctx)
847: {
848:   PetscViewer viewer = (PetscViewer) ctx;
849:   Vec         x;
850:   int         ierr;

854:   SNESGetFunction(snes, &x, PETSC_NULL, PETSC_NULL);
855:   VecView(x, viewer);
856:   return(0);
857: }

859: #undef  __FUNCT__
861: /*@
862:   GSNESErrorMonitor - Displays the error at each iteration.

864:   Collective on GSNES

866:   Input Parameters:
867: + snes   - The nonlinear context
868: . it     - The number of iterations so far
869: . rnorm  - The current (approximate) residual norm
870: - monCtx - A GSNESErrorMonitorCtx

872:   Notes: 
873:   The final argument to SNESSetMonitor() with this routine must be a pointer to a GSNESErrorMonitorCtx.

875:   Level: advanced

877: .keywords: grid SNES, monitor, residual
878: .seealso: GSNESSolutionMonitor(), GSNESResidualMonitor(), SNESDefaultMonitor(), SNESSetMonitor()
879: @*/
880: int GSNESErrorMonitor(GSNES snes, int it, PetscReal rnorm, void *monCtx)
881: {
882:   GSNESErrorMonitorCtx *errorCtx = (GSNESErrorMonitorCtx *) monCtx;
883:   void                 *ctx      = errorCtx->ctx;
884:   PetscScalar           minusOne = -1.0;
885:   GVec                  x, e;
886:   PetscReal             norm_2, norm_max, sol_norm_2;
887:   int                   ierr;

891:   if (it == 0) {
892:     /* Build solution at the beginning */
893:     (*errorCtx->solutionFunc)(errorCtx->solution, ctx);
894:   }
895:   SNESGetSolution(snes, &x);
896:   GVecGetWorkGVec(errorCtx->solution, &e);
897:   VecWAXPY(&minusOne, x, errorCtx->solution, e);
898:   VecView(e, errorCtx->error_viewer);

900:   /* Compute 2-norm and max-norm of error */
901:   if (errorCtx->norm_error_viewer) {
902:     VecNorm(e, NORM_2, &norm_2);
903:     VecNorm(e, NORM_MAX, &norm_max);
904:     VecNorm(errorCtx->solution, NORM_2, &sol_norm_2);
905:     PetscViewerASCIIPrintf(errorCtx->norm_error_viewer,
906:                            "Iteration %d residual norm %g error 2-norm %g error max-norm %g exact sol norm-2 %gn",
907:                            it, rnorm, norm_2, norm_max, sol_norm_2);
908:   }
909:   GVecRestoreWorkGVec(errorCtx->solution, &e);
910:   return(0);
911: }

913: EXTERN_C_BEGIN
914: #undef  __FUNCT__
916: int GSNESOptionsChecker_Private(GSNES snes)
917: {
918:   char       *prefix;
919:   int         loc[4], nmax;
920:   MPI_Comm    comm;
921:   PetscViewer viewer;
922:   PetscDraw   draw;
923:   PetscTruth  opt;
924:   int         ierr;

927:   SNESGetOptionsPrefix(snes, &prefix);
928:   PetscObjectGetComm((PetscObject) snes, &comm);

930:   nmax = 4;
931:   loc[0] = PETSC_DECIDE; loc[1] = PETSC_DECIDE; loc[2] = 300; loc[3] = 300;
932:   PetscOptionsGetIntArray(prefix, "-gvec_snes_solutionmonitor", loc, &nmax, &opt);
933:   if (opt) {
934:     PetscViewerDrawOpen(comm, 0, 0, loc[0], loc[1], loc[2], loc[3], &viewer);
935:     PetscViewerDrawGetDraw(viewer, 0, &draw);
936:     PetscDrawSetTitle(draw, "Approx. Solution");
937:     PetscLogObjectParent(snes, (PetscObject) viewer);
938:     SNESSetMonitor(snes, GSNESSolutionMonitor, (void *) viewer, PETSC_NULL);
939:   }
940:   nmax = 4;
941:   loc[0] = PETSC_DECIDE; loc[1] = PETSC_DECIDE; loc[2] = 300; loc[3] = 300;
942:   PetscOptionsGetIntArray(prefix, "-gvec_snes_residualmonitor", loc, &nmax, &opt);
943:   if (opt) {
944:     PetscViewerDrawOpen(comm, 0, 0, loc[0], loc[1], loc[2], loc[3], &viewer);
945:     PetscViewerDrawGetDraw(viewer, 0, &draw);
946:     PetscDrawSetTitle(draw, "Residual");
947:     PetscLogObjectParent(snes, (PetscObject) viewer);
948:     SNESSetMonitor(snes, GSNESResidualMonitor, (void *) viewer, PETSC_NULL);
949:   }
950: #if 0
951:   nmax = 4;
952:   loc[0] = PETSC_DECIDE; loc[1] = PETSC_DECIDE; loc[2] = 300; loc[3] = 300;
953:   PetscOptionsGetIntArray(prefix, "-gvec_snes_errormonitor", loc, &nmax, &opt);
954:   if (opt) {
955:     PetscViewerDrawOpen(comm, 0, 0, loc[0], loc[1], loc[2], loc[3], &viewer);
956:     PetscViewerDrawGetDraw(viewer, 0, &draw);
957:     PetscDrawSetTitle(draw, "Error");
958:     PetscLogObjectParent(snes, (PetscObject) viewer);
959:     SNESSetMonitor(snes, GSNESErrorMonitor, (void *) viewer, PETSC_NULL);
960:   }
961: #endif
962:   /*--------------------------------------------------------------------- */
963:   PetscOptionsHasName(prefix, "-gvec_snes_fd", &opt);
964:   if (opt) {
965:     GMat A,B;
966:     SNESGetJacobian(snes, &A, &B, PETSC_NULL, PETSC_NULL);
967:     SNESSetJacobian(snes, A, B, GSNESDefaultComputeJacobianWithColoring, 0);
968:     PetscLogInfo(snes, "GSNESSetFromOptions: Setting gvec finite difference Jacobian matrixn");
969:   }
970:   /*--------------------------------------------------------------------- */
971:   PetscOptionsHasName(PETSC_NULL, "-help", &opt);
972:   if (opt) {
973:     char *pprefix;
974:     int   len = 2;
975:     if (prefix != PETSC_NULL) {
976:       PetscStrlen(prefix, &len);
977:     }
978:     PetscMalloc((len+2) * sizeof(char), &pprefix);
979:     PetscStrcpy(pprefix, "-");
980:     if (prefix != PETSC_NULL)
981:       PetscStrcat(pprefix, prefix);
982:     PetscPrintf(comm," Additional SNES Monitor options for grid vectorsn");
983:     PetscPrintf(comm,"   %sgvec_snes_solutionmonitorn", pprefix);
984:     PetscPrintf(comm,"   %sgvec_snes_residualmonitorn", pprefix);
985:     PetscPrintf(comm,"   %sgvec_snes_errormonitor -- not working yetn", pprefix);
986:     PetscFree(pprefix);
987:   }
988:   return(0);
989: }
990: EXTERN_C_END