Actual source code: vi.c

petsc-3.5.4 2015-05-23
Report Typos and Errors
  1: #include <petsc-private/snesimpl.h>  /*I "petscsnes.h" I*/
  2: #include <petscdm.h>

  6: /*@C
  7:    SNESVISetComputeVariableBounds - Sets a function  that is called to compute the variable bounds

  9:    Input parameter
 10: +  snes - the SNES context
 11: -  compute - computes the bounds

 13:    Level: advanced

 15: .seealso:   SNESVISetVariableBounds()

 17: @*/
 18: PetscErrorCode SNESVISetComputeVariableBounds(SNES snes, PetscErrorCode (*compute)(SNES,Vec,Vec))
 19: {
 20:   PetscErrorCode ierr,(*f)(SNES,PetscErrorCode (*)(SNES,Vec,Vec));

 24:   PetscObjectQueryFunction((PetscObject)snes,"SNESVISetComputeVariableBounds_C",&f);
 25:   if (!f) {SNESSetType(snes,SNESVINEWTONRSLS);}
 26:   PetscUseMethod(snes,"SNESVISetComputeVariableBounds_C",(SNES,PetscErrorCode (*)(SNES,Vec,Vec)),(snes,compute));
 27:   return(0);
 28: }

 32: PetscErrorCode SNESVISetComputeVariableBounds_VI(SNES snes,SNESVIComputeVariableBoundsFunction compute)
 33: {
 35:   snes->ops->computevariablebounds = compute;
 36:   return(0);
 37: }

 41: /*
 42:    SNESVIComputeInactiveSetIS - Gets the global indices for the bogus inactive set variables

 44:    Input parameter
 45: .  snes - the SNES context
 46: .  X    - the snes solution vector

 48:    Output parameter
 49: .  ISact - active set index set

 51:  */
 52: PetscErrorCode SNESVIComputeInactiveSetIS(Vec upper,Vec lower,Vec X,Vec F,IS *inact)
 53: {
 54:   PetscErrorCode    ierr;
 55:   const PetscScalar *x,*xl,*xu,*f;
 56:   PetscInt          *idx_act,i,nlocal,nloc_isact=0,ilow,ihigh,i1=0;

 59:   VecGetLocalSize(X,&nlocal);
 60:   VecGetOwnershipRange(X,&ilow,&ihigh);
 61:   VecGetArrayRead(X,&x);
 62:   VecGetArrayRead(lower,&xl);
 63:   VecGetArrayRead(upper,&xu);
 64:   VecGetArrayRead(F,&f);
 65:   /* Compute inactive set size */
 66:   for (i=0; i < nlocal; i++) {
 67:     if (((PetscRealPart(x[i]) > PetscRealPart(xl[i]) + 1.e-8 || (PetscRealPart(f[i]) < 0.0)) && ((PetscRealPart(x[i]) < PetscRealPart(xu[i]) - 1.e-8) || PetscRealPart(f[i]) > 0.0))) nloc_isact++;
 68:   }

 70:   PetscMalloc1(nloc_isact,&idx_act);

 72:   /* Set inactive set indices */
 73:   for (i=0; i < nlocal; i++) {
 74:     if (((PetscRealPart(x[i]) > PetscRealPart(xl[i]) + 1.e-8 || (PetscRealPart(f[i]) < 0.0)) && ((PetscRealPart(x[i]) < PetscRealPart(xu[i]) - 1.e-8) || PetscRealPart(f[i]) > 0.0))) idx_act[i1++] = ilow+i;
 75:   }

 77:   /* Create inactive set IS */
 78:   ISCreateGeneral(PetscObjectComm((PetscObject)upper),nloc_isact,idx_act,PETSC_OWN_POINTER,inact);

 80:   VecRestoreArrayRead(X,&x);
 81:   VecRestoreArrayRead(lower,&xl);
 82:   VecRestoreArrayRead(upper,&xu);
 83:   VecRestoreArrayRead(F,&f);
 84:   return(0);
 85: }

 87: /* --------------------------------------------------------------------------------------------------------*/

 91: PetscErrorCode  SNESMonitorVI(SNES snes,PetscInt its,PetscReal fgnorm,void *dummy)
 92: {
 93:   PetscErrorCode    ierr;
 94:   PetscViewer       viewer = dummy ? (PetscViewer) dummy : PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));
 95:   const PetscScalar *x,*xl,*xu,*f;
 96:   PetscInt          i,n,act[2] = {0,0},fact[2],N;
 97:   /* Number of components that actually hit the bounds (c.f. active variables) */
 98:   PetscInt  act_bound[2] = {0,0},fact_bound[2];
 99:   PetscReal rnorm,fnorm;
100:   double    tmp;

103:   VecGetLocalSize(snes->vec_sol,&n);
104:   VecGetSize(snes->vec_sol,&N);
105:   VecGetArrayRead(snes->xl,&xl);
106:   VecGetArrayRead(snes->xu,&xu);
107:   VecGetArrayRead(snes->vec_sol,&x);
108:   VecGetArrayRead(snes->vec_func,&f);

110:   rnorm = 0.0;
111:   for (i=0; i<n; i++) {
112:     if (((PetscRealPart(x[i]) > PetscRealPart(xl[i]) + 1.e-8 || (PetscRealPart(f[i]) < 0.0)) && ((PetscRealPart(x[i]) < PetscRealPart(xu[i]) - 1.e-8) || PetscRealPart(f[i]) > 0.0))) rnorm += PetscRealPart(PetscConj(f[i])*f[i]);
113:     else if (PetscRealPart(x[i]) <= PetscRealPart(xl[i]) + 1.e-8 && PetscRealPart(f[i]) >= 0.0) act[0]++;
114:     else if (PetscRealPart(x[i]) >= PetscRealPart(xu[i]) - 1.e-8 && PetscRealPart(f[i]) <= 0.0) act[1]++;
115:     else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"Can never get here");
116:   }

118:   for (i=0; i<n; i++) {
119:     if (PetscRealPart(x[i]) <= PetscRealPart(xl[i]) + 1.e-8) act_bound[0]++;
120:     else if (PetscRealPart(x[i]) >= PetscRealPart(xu[i]) - 1.e-8) act_bound[1]++;
121:   }
122:   VecRestoreArrayRead(snes->vec_func,&f);
123:   VecRestoreArrayRead(snes->xl,&xl);
124:   VecRestoreArrayRead(snes->xu,&xu);
125:   VecRestoreArrayRead(snes->vec_sol,&x);
126:   MPI_Allreduce(&rnorm,&fnorm,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)snes));
127:   MPI_Allreduce(act,fact,2,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)snes));
128:   MPI_Allreduce(act_bound,fact_bound,2,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)snes));
129:   fnorm = PetscSqrtReal(fnorm);

131:   PetscViewerASCIIAddTab(viewer,((PetscObject)snes)->tablevel);
132:   if (snes->ntruebounds) tmp = ((double)(fact[0]+fact[1]))/((double)snes->ntruebounds);
133:   else tmp = 0.0;
134:   PetscViewerASCIIPrintf(viewer,"%3D SNES VI Function norm %14.12e Active lower constraints %D/%D upper constraints %D/%D Percent of total %g Percent of bounded %g\n",its,(double)fnorm,fact[0],fact_bound[0],fact[1],fact_bound[1],((double)(fact[0]+fact[1]))/((double)N),tmp);

136:   PetscViewerASCIISubtractTab(viewer,((PetscObject)snes)->tablevel);
137:   return(0);
138: }

140: /*
141:      Checks if J^T F = 0 which implies we've found a local minimum of the norm of the function,
142:     || F(u) ||_2 but not a zero, F(u) = 0. In the case when one cannot compute J^T F we use the fact that
143:     0 = (J^T F)^T W = F^T J W iff W not in the null space of J. Thanks for Jorge More
144:     for this trick. One assumes that the probability that W is in the null space of J is very, very small.
145: */
148: PetscErrorCode SNESVICheckLocalMin_Private(SNES snes,Mat A,Vec F,Vec W,PetscReal fnorm,PetscBool *ismin)
149: {
150:   PetscReal      a1;
152:   PetscBool      hastranspose;

155:   *ismin = PETSC_FALSE;
156:   MatHasOperation(A,MATOP_MULT_TRANSPOSE,&hastranspose);
157:   if (hastranspose) {
158:     /* Compute || J^T F|| */
159:     MatMultTranspose(A,F,W);
160:     VecNorm(W,NORM_2,&a1);
161:     PetscInfo1(snes,"|| J^T F|| %g near zero implies found a local minimum\n",(double)(a1/fnorm));
162:     if (a1/fnorm < 1.e-4) *ismin = PETSC_TRUE;
163:   } else {
164:     Vec         work;
165:     PetscScalar result;
166:     PetscReal   wnorm;

168:     VecSetRandom(W,NULL);
169:     VecNorm(W,NORM_2,&wnorm);
170:     VecDuplicate(W,&work);
171:     MatMult(A,W,work);
172:     VecDot(F,work,&result);
173:     VecDestroy(&work);
174:     a1   = PetscAbsScalar(result)/(fnorm*wnorm);
175:     PetscInfo1(snes,"(F^T J random)/(|| F ||*||random|| %g near zero implies found a local minimum\n",(double)a1);
176:     if (a1 < 1.e-4) *ismin = PETSC_TRUE;
177:   }
178:   return(0);
179: }

181: /*
182:      Checks if J^T(F - J*X) = 0
183: */
186: PetscErrorCode SNESVICheckResidual_Private(SNES snes,Mat A,Vec F,Vec X,Vec W1,Vec W2)
187: {
188:   PetscReal      a1,a2;
190:   PetscBool      hastranspose;

193:   MatHasOperation(A,MATOP_MULT_TRANSPOSE,&hastranspose);
194:   if (hastranspose) {
195:     MatMult(A,X,W1);
196:     VecAXPY(W1,-1.0,F);

198:     /* Compute || J^T W|| */
199:     MatMultTranspose(A,W1,W2);
200:     VecNorm(W1,NORM_2,&a1);
201:     VecNorm(W2,NORM_2,&a2);
202:     if (a1 != 0.0) {
203:       PetscInfo1(snes,"||J^T(F-Ax)||/||F-AX|| %g near zero implies inconsistent rhs\n",(double)(a2/a1));
204:     }
205:   }
206:   return(0);
207: }

209: /*
210:   SNESConvergedDefault_VI - Checks the convergence of the semismooth newton algorithm.

212:   Notes:
213:   The convergence criterion currently implemented is
214:   merit < abstol
215:   merit < rtol*merit_initial
216: */
219: PetscErrorCode SNESConvergedDefault_VI(SNES snes,PetscInt it,PetscReal xnorm,PetscReal gradnorm,PetscReal fnorm,SNESConvergedReason *reason,void *dummy)
220: {


227:   *reason = SNES_CONVERGED_ITERATING;

229:   if (!it) {
230:     /* set parameter for default relative tolerance convergence test */
231:     snes->ttol = fnorm*snes->rtol;
232:   }
233:   if (fnorm != fnorm) {
234:     PetscInfo(snes,"Failed to converged, function norm is NaN\n");
235:     *reason = SNES_DIVERGED_FNORM_NAN;
236:   } else if (fnorm < snes->abstol) {
237:     PetscInfo2(snes,"Converged due to function norm %g < %g\n",(double)fnorm,(double)snes->abstol);
238:     *reason = SNES_CONVERGED_FNORM_ABS;
239:   } else if (snes->nfuncs >= snes->max_funcs) {
240:     PetscInfo2(snes,"Exceeded maximum number of function evaluations: %D > %D\n",snes->nfuncs,snes->max_funcs);
241:     *reason = SNES_DIVERGED_FUNCTION_COUNT;
242:   }

244:   if (it && !*reason) {
245:     if (fnorm < snes->ttol) {
246:       PetscInfo2(snes,"Converged due to function norm %g < %g (relative tolerance)\n",(double)fnorm,(double)snes->ttol);
247:       *reason = SNES_CONVERGED_FNORM_RELATIVE;
248:     }
249:   }
250:   return(0);
251: }


254: /* -------------------------------------------------------------------------- */
255: /*
256:    SNESVIProjectOntoBounds - Projects X onto the feasible region so that Xl[i] <= X[i] <= Xu[i] for i = 1...n.

258:    Input Parameters:
259: .  SNES - nonlinear solver context

261:    Output Parameters:
262: .  X - Bound projected X

264: */

268: PetscErrorCode SNESVIProjectOntoBounds(SNES snes,Vec X)
269: {
270:   PetscErrorCode    ierr;
271:   const PetscScalar *xl,*xu;
272:   PetscScalar       *x;
273:   PetscInt          i,n;

276:   VecGetLocalSize(X,&n);
277:   VecGetArray(X,&x);
278:   VecGetArrayRead(snes->xl,&xl);
279:   VecGetArrayRead(snes->xu,&xu);

281:   for (i = 0; i<n; i++) {
282:     if (PetscRealPart(x[i]) < PetscRealPart(xl[i])) x[i] = xl[i];
283:     else if (PetscRealPart(x[i]) > PetscRealPart(xu[i])) x[i] = xu[i];
284:   }
285:   VecRestoreArray(X,&x);
286:   VecRestoreArrayRead(snes->xl,&xl);
287:   VecRestoreArrayRead(snes->xu,&xu);
288:   return(0);
289: }


294: /*
295:    SNESVIGetActiveSetIndices - Gets the global indices for the active set variables

297:    Input parameter
298: .  snes - the SNES context
299: .  X    - the snes solution vector
300: .  F    - the nonlinear function vector

302:    Output parameter
303: .  ISact - active set index set
304:  */
305: PetscErrorCode SNESVIGetActiveSetIS(SNES snes,Vec X,Vec F,IS *ISact)
306: {
307:   PetscErrorCode    ierr;
308:   Vec               Xl=snes->xl,Xu=snes->xu;
309:   const PetscScalar *x,*f,*xl,*xu;
310:   PetscInt          *idx_act,i,nlocal,nloc_isact=0,ilow,ihigh,i1=0;

313:   VecGetLocalSize(X,&nlocal);
314:   VecGetOwnershipRange(X,&ilow,&ihigh);
315:   VecGetArrayRead(X,&x);
316:   VecGetArrayRead(Xl,&xl);
317:   VecGetArrayRead(Xu,&xu);
318:   VecGetArrayRead(F,&f);
319:   /* Compute active set size */
320:   for (i=0; i < nlocal;i++) {
321:     if (!((PetscRealPart(x[i]) > PetscRealPart(xl[i]) + 1.e-8 || (PetscRealPart(f[i]) < 0.0)) && ((PetscRealPart(x[i]) < PetscRealPart(xu[i]) - 1.e-8) || PetscRealPart(f[i]) > 0.0))) nloc_isact++;
322:   }

324:   PetscMalloc1(nloc_isact,&idx_act);

326:   /* Set active set indices */
327:   for (i=0; i < nlocal; i++) {
328:     if (!((PetscRealPart(x[i]) > PetscRealPart(xl[i]) + 1.e-8 || (PetscRealPart(f[i]) < 0.0)) && ((PetscRealPart(x[i]) < PetscRealPart(xu[i]) - 1.e-8) || PetscRealPart(f[i]) > 0.0))) idx_act[i1++] = ilow+i;
329:   }

331:   /* Create active set IS */
332:   ISCreateGeneral(PetscObjectComm((PetscObject)snes),nloc_isact,idx_act,PETSC_OWN_POINTER,ISact);

334:   VecRestoreArrayRead(X,&x);
335:   VecRestoreArrayRead(Xl,&xl);
336:   VecRestoreArrayRead(Xu,&xu);
337:   VecRestoreArrayRead(F,&f);
338:   return(0);
339: }

343: PetscErrorCode SNESVICreateIndexSets_RS(SNES snes,Vec X,Vec F,IS *ISact,IS *ISinact)
344: {
346:   PetscInt       rstart,rend;

349:   SNESVIGetActiveSetIS(snes,X,F,ISact);
350:   VecGetOwnershipRange(X,&rstart,&rend);
351:   ISComplement(*ISact,rstart,rend,ISinact);
352:   return(0);
353: }

357: PetscErrorCode SNESVIComputeInactiveSetFnorm(SNES snes,Vec F,Vec X, PetscReal *fnorm)
358: {
359:   PetscErrorCode    ierr;
360:   const PetscScalar *x,*xl,*xu,*f;
361:   PetscInt          i,n;
362:   PetscReal         rnorm;

365:   VecGetLocalSize(X,&n);
366:   VecGetArrayRead(snes->xl,&xl);
367:   VecGetArrayRead(snes->xu,&xu);
368:   VecGetArrayRead(X,&x);
369:   VecGetArrayRead(F,&f);
370:   rnorm = 0.0;
371:   for (i=0; i<n; i++) {
372:     if (((PetscRealPart(x[i]) > PetscRealPart(xl[i]) + 1.e-8 || (PetscRealPart(f[i]) < 0.0)) && ((PetscRealPart(x[i]) < PetscRealPart(xu[i]) - 1.e-8) || PetscRealPart(f[i]) > 0.0))) rnorm += PetscRealPart(PetscConj(f[i])*f[i]);
373:   }
374:   VecRestoreArrayRead(F,&f);
375:   VecRestoreArrayRead(snes->xl,&xl);
376:   VecRestoreArrayRead(snes->xu,&xu);
377:   VecRestoreArrayRead(X,&x);
378:   MPI_Allreduce(&rnorm,fnorm,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)snes));
379:   *fnorm = PetscSqrtReal(*fnorm);
380:   return(0);
381: }

385: PetscErrorCode SNESVIDMComputeVariableBounds(SNES snes,Vec xl, Vec xu)
386: {

390:   DMComputeVariableBounds(snes->dm, xl, xu);
391:   return(0);
392: }


395: /* -------------------------------------------------------------------------- */
396: /*
397:    SNESSetUp_VI - Does setup common to all VI solvers -- basically makes sure bounds have been properly set up
398:    of the SNESVI nonlinear solver.

400:    Input Parameter:
401: .  snes - the SNES context
402: .  x - the solution vector

404:    Application Interface Routine: SNESSetUp()

406:    Notes:
407:    For basic use of the SNES solvers, the user need not explicitly call
408:    SNESSetUp(), since these actions will automatically occur during
409:    the call to SNESSolve().
410:  */
413: PetscErrorCode SNESSetUp_VI(SNES snes)
414: {
416:   PetscInt       i_start[3],i_end[3];

419:   SNESSetWorkVecs(snes,1);
420:   SNESSetUpMatrices(snes);

422:   if (!snes->ops->computevariablebounds && snes->dm) {
423:     PetscBool flag;
424:     DMHasVariableBounds(snes->dm, &flag);

426:     snes->ops->computevariablebounds = SNESVIDMComputeVariableBounds;
427:   }
428:   if (!snes->usersetbounds) {
429:     if (snes->ops->computevariablebounds) {
430:       if (!snes->xl) {VecDuplicate(snes->vec_sol,&snes->xl);}
431:       if (!snes->xu) {VecDuplicate(snes->vec_sol,&snes->xu);}
432:       (*snes->ops->computevariablebounds)(snes,snes->xl,snes->xu);
433:     } else if (!snes->xl && !snes->xu) {
434:       /* If the lower and upper bound on variables are not set, set it to -Inf and Inf */
435:       VecDuplicate(snes->vec_sol, &snes->xl);
436:       VecSet(snes->xl,PETSC_NINFINITY);
437:       VecDuplicate(snes->vec_sol, &snes->xu);
438:       VecSet(snes->xu,PETSC_INFINITY);
439:     } else {
440:       /* Check if lower bound, upper bound and solution vector distribution across the processors is identical */
441:       VecGetOwnershipRange(snes->vec_sol,i_start,i_end);
442:       VecGetOwnershipRange(snes->xl,i_start+1,i_end+1);
443:       VecGetOwnershipRange(snes->xu,i_start+2,i_end+2);
444:       if ((i_start[0] != i_start[1]) || (i_start[0] != i_start[2]) || (i_end[0] != i_end[1]) || (i_end[0] != i_end[2]))
445:         SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Distribution of lower bound, upper bound and the solution vector should be identical across all the processors.");
446:     }
447:   }
448:   return(0);
449: }
450: /* -------------------------------------------------------------------------- */
453: PetscErrorCode SNESReset_VI(SNES snes)
454: {

458:   VecDestroy(&snes->xl);
459:   VecDestroy(&snes->xu);
460:   snes->usersetbounds = PETSC_FALSE;
461:   return(0);
462: }

464: /*
465:    SNESDestroy_VI - Destroys the private SNES_VI context that was created
466:    with SNESCreate_VI().

468:    Input Parameter:
469: .  snes - the SNES context

471:    Application Interface Routine: SNESDestroy()
472:  */
475: PetscErrorCode SNESDestroy_VI(SNES snes)
476: {

480:   PetscFree(snes->data);

482:   /* clear composed functions */
483:   PetscObjectComposeFunction((PetscObject)snes,"SNESLineSearchSet_C",NULL);
484:   PetscObjectComposeFunction((PetscObject)snes,"SNESLineSearchSetMonitor_C",NULL);
485:   return(0);
486: }

490: /*@
491:    SNESVISetVariableBounds - Sets the lower and upper bounds for the solution vector. xl <= x <= xu.

493:    Input Parameters:
494: .  snes - the SNES context.
495: .  xl   - lower bound.
496: .  xu   - upper bound.

498:    Notes:
499:    If this routine is not called then the lower and upper bounds are set to
500:    PETSC_INFINITY and PETSC_NINFINITY respectively during SNESSetUp().

502:    Level: advanced

504: @*/
505: PetscErrorCode SNESVISetVariableBounds(SNES snes, Vec xl, Vec xu)
506: {
507:   PetscErrorCode ierr,(*f)(SNES,Vec,Vec);

513:   PetscObjectQueryFunction((PetscObject)snes,"SNESVISetVariableBounds_C",&f);
514:   if (!f) {SNESSetType(snes,SNESVINEWTONRSLS);}
515:   PetscUseMethod(snes,"SNESVISetVariableBounds_C",(SNES,Vec,Vec),(snes,xl,xu));
516:   snes->usersetbounds = PETSC_TRUE;
517:   return(0);
518: }

522: PetscErrorCode SNESVISetVariableBounds_VI(SNES snes,Vec xl,Vec xu)
523: {
524:   PetscErrorCode    ierr;
525:   const PetscScalar *xxl,*xxu;
526:   PetscInt          i,n, cnt = 0;

529:   SNESGetFunction(snes,&snes->vec_func,NULL,NULL);
530:   if (!snes->vec_func) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() or SNESSetDM() first");
531:   {
532:     PetscInt xlN,xuN,N;
533:     VecGetSize(xl,&xlN);
534:     VecGetSize(xu,&xuN);
535:     VecGetSize(snes->vec_func,&N);
536:     if (xlN != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector lengths lower bound = %D solution vector = %D",xlN,N);
537:     if (xuN != N) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Incompatible vector lengths: upper bound = %D solution vector = %D",xuN,N);
538:   }
539:   SNESSetType(snes,SNESVINEWTONRSLS);
540:   PetscObjectReference((PetscObject)xl);
541:   PetscObjectReference((PetscObject)xu);
542:   VecDestroy(&snes->xl);
543:   VecDestroy(&snes->xu);
544:   snes->xl = xl;
545:   snes->xu = xu;
546:   VecGetLocalSize(xl,&n);
547:   VecGetArrayRead(xl,&xxl);
548:   VecGetArrayRead(xu,&xxu);
549:   for (i=0; i<n; i++) cnt += ((xxl[i] != PETSC_NINFINITY) || (xxu[i] != PETSC_INFINITY));

551:   MPI_Allreduce(&cnt,&snes->ntruebounds,1,MPIU_INT,MPI_SUM,PetscObjectComm((PetscObject)snes));
552:   VecRestoreArrayRead(xl,&xxl);
553:   VecRestoreArrayRead(xu,&xxu);
554:   return(0);
555: }

559: PetscErrorCode SNESSetFromOptions_VI(SNES snes)
560: {
562:   PetscBool      flg;
563:   SNESLineSearch linesearch;

566:   PetscOptionsHead("SNES VI options");
567:   PetscOptionsBool("-snes_vi_monitor","Monitor all non-active variables","None",PETSC_FALSE,&flg,0);
568:   if (flg) {
569:     SNESMonitorSet(snes,SNESMonitorVI,0,0);
570:   }
571:   if (!snes->linesearch) {
572:     SNESGetLineSearch(snes, &linesearch);
573:     SNESLineSearchSetType(linesearch, SNESLINESEARCHBT);
574:     SNESLineSearchBTSetAlpha(linesearch, 0.0);
575:   }
576:   PetscOptionsTail();
577:   return(0);
578: }