<div dir="ltr">I've attached two modified versions of ex1:<div><br></div><div>ex1_powell.c uses the Powell restart</div><div>ex1_none.c uses no restart</div><div><br></div><div>For the default initial guess (x0 = [0.5,0.5]), both converge just fine. However, for the initial guess x0 = [3.,3.], the Powell solution fails to converge, while None and Periodic both still converge. This is with the "easy" equation set (run without -hard).<br></div><div><br></div><div>Interestingly enough, the Powell restart still "finishes" in a reasonable number of iterations (7 iterations), but the residual is very large (on the order of 1e254).</div></div><div class="gmail_extra"><br><div class="gmail_quote">On Thu, Jul 14, 2016 at 5:50 PM, Barry Smith <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:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div class="HOEnZb"><div class="h5"><br>
> On Jul 14, 2016, at 6:18 PM, Andrew Ho <<a href="mailto:andrewh0@uw.edu">andrewh0@uw.edu</a>> wrote:<br>
><br>
> I am trying to solve a simple ionization/recombination ODE using PETSc's quasi-newton SNES.<br>
><br>
> This is a basic non-linear coupled ODE system:<br>
><br>
> delta = -a u^2 + b u v<br>
> d_t u = delta<br>
> d_t v = -delta<br>
><br>
> a and b are constants.<br>
><br>
> I wrote a backwards Euler root finding function (yes, I know the TS module has BE implemented, but this is more of a learning exercise).<br>
><br>
> Here is the function evaluation:<br>
><br>
> struct ion_rec_ctx<br>
> {<br>
> PetscScalar rate_a, rate_b;<br>
> PetscScalar dt;<br>
> };<br>
> PetscErrorCode bdf1(SNES snes, Vec x, Vec f, void *ctx)<br>
> {<br>
> const PetscScalar *xx;<br>
> PetscScalar *ff;<br>
> ion_rec_ctx& params = *reinterpret_cast<ion_rec_ctx*>(ctx);<br>
> CHKERRQ(VecGetArrayRead(x, &xx));<br>
> CHKERRQ(VecGetArray(f,&ff));<br>
> auto delta = (-params.rate_a*xx[0]*xx[0]+params.rate_b*xx[1]*xx[0]);<br>
> ff[0] = xx[0]-params.dt*delta;<br>
> ff[1] = xx[1]-params.dt*-delta;<br>
> CHKERRQ(VecRestoreArrayRead(x,&xx));<br>
> CHKERRQ(VecRestoreArray(f,&ff));<br>
> return 0;<br>
> }<br>
><br>
> To setup the solver and solve one time step:<br>
><br>
> // q0, q1, and res are Vec's previously initialized<br>
> // initial conditions: q0 = [1e19,1e19]<br>
> SNES solver;<br>
> CHKERRQ(SNESCreate(comm, &solver));<br>
> CHKERRQ(SNESSetType(solver, SNESQN));<br>
> CHKERRQ(SNESQNSetType(solver, SNES_QN_LBFGS));<br>
> ion_rec_ctx params = {9.59e-16, 1.15e-19, 1.};<br>
> CHKERRQ(SNESSetFunction(solver, res, &bdf1, ¶ms));<br>
> CHKERRQ(SNESSolve(solver, q0, q1));<br>
><br>
> When I run this, the solver fails to converge to a solution for this rather large time step.<br>
> The solution produced when the SNES module finally gives up is:<br>
><br>
> q1 = [-2.72647e142, 2.72647e142]<br>
><br>
> For reference, when I disable the scale and restart types, I get these values:<br>
><br>
> q1 = [1.0279e17, 1.98972e19]<br>
><br>
> This is only a problem when I use the SNES_QN_RESTART_POWELL restart type (seems to be regardless of the scale type type). I get reasonable answers for other combinations of restart/scale type. I've tried every combination of restart type/scale type except for SNES_QN_SCALE_JACOBIAN (my ultimate application doesn't have an available Jacobian), and only cases using SNES_QN_RESTART_POWELL are failing.<br>
><br>
> I'm unfamiliar with Powell's restart criterion, but is it suppose to work reasonably well with Quasi-Newton methods? I tried it on the simple problem given in this example: <a href="http://www.mcs.anl.gov/petsc/petsc-current/src/snes/examples/tutorials/ex1.c.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/petsc-current/src/snes/examples/tutorials/ex1.c.html</a><br>
><br>
> And Powell restarts also fails to converge to a meaningful solution (solving for f(x) = [1,1], for x0 = [1,1]), but the other restart methods do converge properly.<br>
<br>
</div></div> Could you please send the exact options you are using for the ex1.c that both fail and work and we'll see if there is some problem with the Powell restart.<br>
<br>
Thanks<br>
<span class="HOEnZb"><font color="#888888"><br>
Barry<br>
</font></span><div class="HOEnZb"><div class="h5"><br>
><br>
> Software information:<br>
><br>
> PETSc version 3.7.2 (built from git maint branch)<br>
> PETSc arch: arch-linux2-c-opt<br>
> OS: Ubuntu 15.04 x64<br>
> Compiler: gcc 4.9.2<br>
<br>
</div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr">Andrew Ho</div></div>
</div>