<div dir="ltr"><div><div>Thanks Barry and Matthew.<br><br></div>@Barry: I'm following the same procedure as you've mentioned - PetscOptionsSetValue() precede SNESSetFromOptions. Here's the snippet for my code:<br><br>-----------------------------------------------------------------------------------------------------------<br><br> error = SNESCreate(PETSC_COMM_WORLD,&snes);<br> checkPETScError(error,<br> "SNESCreate failed.");<br> <br> error = SNESSetType(snes, SNESQN);<br> checkPETScError(error,<br> "SNESSetType failed.");<br> <br> error = SNESQNSetType(snes, SNES_QN_LBFGS);<br> checkPETScError(error,<br> "SNESQNSetType failed.");<br> <br> error = SNESQNSetScaleType(snes, SNES_QN_SCALE_SHANNO);<br> checkPETScError(error,<br> "SNESQNSetScaleType failed.");<br> <br> error = SNESQNSetRestartType(snes, SNES_QN_RESTART_PERIODIC);<br> checkPETScError(error,<br> "SNESQNSetRestartType failed.");<br> <br> error = PetscOptionsSetValue("-snes_qn_m","500");<br> checkPETScError(error,<br> "PETScOptionsSetValue failed.");<br><br> SNESLineSearch linesearch;<br> error = SNESGetLineSearch(snes,&linesearch);<br> checkPETScError(error,<br> "SNESGetLineSearch failed.");<br> <br> error = SNESLineSearchSetType(linesearch,SNESLINESEARCHCP);<br> checkPETScError(error,<br> "SNESLineSearchSetType failed.");<br> <br> error = PetscOptionsSetValue("-snes_linesearch_max_it", "1");<br> checkPETScError(error,<br> "PetscOptionsSetValue failed.");<br><br> error = SNESLineSearchView(linesearch, PETSC_VIEWER_STDOUT_WORLD);<br> checkPETScError(error,<br> "SNESLineSearchView failed.");<br><br> error =SNESLineSearchSetMonitor(linesearch, <br> PETSC_TRUE);<br> checkPETScError(error,<br> "SNESLineSearchSet Monitor failed.");<br><br> error = SNESLineSearchSetFromOptions(linesearch);<br> checkPETScError(error,<br> "SNESLineSearchSetFromOptions failed.");<br><br> SNESLineSearchReason lineSearchReason;<br> error = SNESLineSearchGetReason(linesearch, &lineSearchReason);<br> checkPETScError(error,<br> "SNESLineSearchGetReason failed.");<br><br> error = SNESSetFunction(snes,r,FormFunction,&petscData);<br> checkPETScError(error,<br> "SNESSetFunction failed.");<br> <br> //<br> // Customize KSP <br> //<br> error = SNESGetKSP(snes,&ksp);<br> checkPETScError(error,<br> "SNESGetKSP failed.");<br> <br> error = KSPSetType(ksp,KSPGMRES);<br> checkPETScError(error,<br> "KSPSetType failed.");<br> <br> error = KSPGMRESSetRestart(ksp,300);<br> checkPETScError(error,<br> "KSPGMRESSetRestart failed.");<br> <br> error = KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);<br> checkPETScError(error,<br> "KSPSetInitialGuessNonzero failed.");<br> <br> error = KSPGetPC(ksp,&pc);<br> checkPETScError(error,<br> "KSPGetPC failed.");<br> <br> error = PCSetType(pc,PCJACOBI); <br> checkPETScError(error,<br> "PCSetType failed.");<br> <br> error = PCSetReusePreconditioner(pc,PETSC_TRUE);<br> checkPETScError(error,<br> "PCSetReusePreconditioner failed.");<br> <br> error = KSPSetTolerances(ksp,<br> PETSC_DEFAULT,<br> 1e-15,<br> 1e7,<br> 10000); <br> checkPETScError(error,<br> "KSPSetTolerances failed.");<br><br> error = KSPSetFromOptions(ksp);<br> checkPETScError(error,<br> "Call to KSPSetFromOptions() failed.");<br> <br> //<br> //get reason for non-convergence<br> //<br> KSPConvergedReason kspReason;<br> error = KSPGetConvergedReason(ksp, &kspReason);<br> checkPETScError(error,<br> "Call to KSPGetConvergedReason() failed.");<br><br> if(kspReason < 0)<br> {<br> if(debugLevel != 0)<br> std::cout<<"Other kind of divergence in SNES-KSP : "<< kspReason <<std::endl;<br> <br> } <br><br> PetscInt lag = 1;<br> error = SNESSetLagPreconditioner(snes,<br> lag);<br> checkPETScError(error,<br> "Call to SNESSetLagPreconditioner() failed.");<br> <br> PetscInt maxFails = 2;<br> error = SNESSetMaxLinearSolveFailures(snes,maxFails);<br> checkPETScError(error,<br> "Call to SNESSetMaxLinearSolveFailures() failed.");<br> <br> PetscReal abstol = 1e-13; // absolute convergence tolerance<br> PetscInt maxit = 100000;<br> error = SNESSetTolerances(snes,<br> abstol,<br> PETSC_DEFAULT,<br> PETSC_DEFAULT,<br> maxit,<br> maxit);<br> checkPETScError(error,<br> "SNESSetTolerances failed.");<br> <br> error = SNESView(snes,<br> PETSC_VIEWER_STDOUT_WORLD);<br> checkPETScError(error,<br> "Call to SNESView() failed.");<br><br> error = SNESMonitorSet(snes,SNESMonitorDefault,PETSC_NULL,PETSC_NULL);<br> checkPETScError(error,<br> "Call to SNESMonitorSet() failed.");<br><br> error = SNESSetFromOptions(snes);<br> checkPETScError(error,<br> "Call to SNESSetFromOptions() failed.");<br><br><br> //<br> // Solve the system<br> //<br> error = SNESSolve(snes,PETSC_NULL,x);<br> checkPETScError(error,<br> "Call to SNESSolve() failed.");<br><br> SNESConvergedReason reason;<br> error = SNESGetConvergedReason(snes,&reason);<br> checkPETScError(error,<br> "Call to SNESGetConvergedReason() failed.");<br><br>------------------------------------------------------------------------------------------------------------------------------------<br><br></div><div>Also, I didn't find any SNESQN examples in my snes/examples folder (using petsc-3.6.3). <br>Moreover, the Powell descent condition seems to be only declared and then assigned a value through the PetscOptionsReal call. Beyond that I didn't find any other mention of it. I was grepping for powell_downhill variable. (Note: powell_downhill features in 3.6.3 and not in 3.7 version). <br><br></div><div>Thanks,<br></div><div>Bikash<br></div><div><br></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Wed, Feb 14, 2018 at 7:02 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><span class="">On Wed, Feb 14, 2018 at 6:43 PM, Smith, Barry F. <span dir="ltr"><<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><br>
Hmm,<br>
<br>
1) make sure you call PetscOptionsSetValue() before you call to SNESSetFromOptions()<br>
<br>
2) make sure you call SNESSetFromOptions()<br>
<br>
3) did you add a prefix to the SNES object? If so make sure you include it in the PetscOptionsSetValue() call.<br>
<br>
I can't see a reason why it won't work. Does it work with the PETSc examples for you or not?<br>
<br>
Regarding the Powell descent option, I'm afraid you'll need to examine the code for exact details. src/snes/impls/qn/qn.c</blockquote><div><br></div></span><div>Here is the description</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/939b553f045c5ba32242d0d49e80e4934ed3bf76/src/snes/impls/qn/qn.c?at=master&fileviewer=file-view-default#qn.c-451" target="_blank">https://bitbucket.org/petsc/<wbr>petsc/src/<wbr>939b553f045c5ba32242d0d49e80e4<wbr>934ed3bf76/src/snes/impls/qn/<wbr>qn.c?at=master&fileviewer=<wbr>file-view-default#qn.c-451</a></div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><span class=""><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><span class="m_3822503704025141643gmail-HOEnZb"><font color="#888888"><br>
Barry<br>
</font></span><div class="m_3822503704025141643gmail-HOEnZb"><div class="m_3822503704025141643gmail-h5"><br>
<br>
> On Feb 14, 2018, at 5:25 PM, Bikash Kanungo <<a href="mailto:bikash@umich.edu" target="_blank">bikash@umich.edu</a>> wrote:<br>
><br>
> Hi,<br>
><br>
> I'm using the L-BFGS QN solver. In order to set the number of past states (also the restart size if I use Periodic restart), to say 50, I'm using PetscOptionsSetValue("-snes_qn<wbr>_m", "50"). However while running, it still shows "Stored subspace size: 10", i.e., the default value of 10 is not overwritten.<br>
><br>
> Additionally, I would like to know more about the the -snes_qn_powell_descent option. For Powell restart, one uses a gamma parameter which I believe is defined by the -snes_qn_powell_gamma option. What exactly does the descent condition do? It would be useful if there are good references to it.<br>
><br>
> Thanks,<br>
> Biksah<br>
><br>
> --<br>
> Bikash S. Kanungo<br>
> PhD Student<br>
> Computational Materials Physics Group<br>
> Mechanical Engineering<br>
> University of Michigan<br>
><br>
<br>
</div></div></blockquote></span></div><span class="HOEnZb"><font color="#888888"><br><br clear="all"><div><br></div>-- <br><div class="m_3822503704025141643gmail_signature"><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.caam.rice.edu/~mk51/" target="_blank">https://www.cse.buffalo.edu/~<wbr>knepley/</a><br></div></div></div></div></div>
</font></span></div></div>
</blockquote></div><br><br clear="all"><br>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div><div><div><div><font color="#666666">Bikash S. Kanungo<br></font></div><font color="#666666">PhD Student<br></font></div><font color="#666666">Computational Materials Physics Group<br></font></div><font color="#666666">Mechanical Engineering <br></font></div><font color="#666666">University of Michigan<br><br></font></div></div>
</div>