<div dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">On Sun, Aug 24, 2014 at 9:40 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=""><br>
On Aug 24, 2014, at 9:06 PM, Derek Gaston <<a href="mailto:friedmud@gmail.com">friedmud@gmail.com</a>> wrote:<br>
<br>
> There must be something I'm missing about SNESSetFunctionDomainError().<br>
><br>
> I'm using a matrix free newton solve and I was hoping that calling SNESSetFunctionDomainError() would completely stop the nonlinear solver... but that doesn't appear to be the case.<br>
><br>
> Now - if I call that function during the initial residual computation it appears to work (the solve never starts).  However, if I call it during the _middle_ of a matrix free nonlinear solve during a residual calculation - it appears to do absolutely nothing.<br>


<br>
</div>   You are not missing anything. The way SNESSetFunctionDomainError() works is it merely sets a flag in the SNES object. The solver code then can check that flag after a function evaluation and generate an if that flag is set (or perhaps do some kind of recovery). So, for example, at the beginning of SNESSolve_NewtonLS() we have<br>


<br>
 } else {<br>
    if (!snes->vec_func_init_set) {<br>
      ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);<br>
      ierr = SNESGetFunctionDomainError(snes, &domainerror);CHKERRQ(ierr);<br>
      if (domainerror) {<br>
        snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN;<br>
        PetscFunctionReturn(0);</blockquote><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
Now the reason it is not erroring out for you is because<br>
<br>
1) our MatMult_MFFD() totally ignores the SNES domain error business because it doesn’t even know about SNES<br>
<br>
2) inside the various nonlinear solvers and many complicated line search codes we may not always be checking the domainerror flag.<br>
<br>
  Both of these things, of course, should be fixed. We need to add support to the MatMult_MFFD() to indicate domain errors and we need to check all the code to make sure we always handle the domain error flag.<br><br></blockquote>

<div>In particular, we check for domain errors inside the linesearch and after the linear update. Which means that, if after the linear update you are still hitting a domain error, the solver will quit.</div><div> </div>
<div>
This problem is somewhat related to the question of what to do inside SNESVI (bounded) MF solves: if a constraint is active, it's quite possible that differencing in a given direction will cause a domain error (this situation might arise, for example, in Relap-7 with its bounds on densities). For reduced space methods we might simply mask out the active components of a differencing direction, since the linear solve ignores those anyway.  It's less clear to me what to do with the bounds that are hit, but are inactive. Or even interior points that are within a differencing parameter of the boundary of the feasible set? The same goes for semismooth solvers (or, perhaps, other methods that we might want to eventually incorporate from Tao).</div>

<div><br></div><div>Also, while we are at it, we ought to unify MatMFFD with SNESComputeJacobianDefault(), which more or less duplicates the MFFD algorithms.</div><div><br></div><div>Dmitry.</div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

<span class="HOEnZb"><font color="#888888">
  Barry<br>
</font></span><div class="HOEnZb"><div class="h5"><br>
<br>
><br>
> Here's my stupid debugging code that's at the top of my residual callback method:<br>
><br>
>     static int count = 0;<br>
>     count++;<br>
><br>
>     if (count == 5)<br>
>     {<br>
>       std::cout<<"stopping solve!"<<std::endl;<br>
>       SNESSetFunctionDomainError(snes);<br>
>     }<br>
><br>
><br>
> During the solve "stopping solve!" gets printed out after computing the 3rd linear iteration of the first newton step... but the solve just continues without stopping...<br>
><br>
> Is there something more I need to do here?<br>
><br>
> Thanks!<br>
><br>
> Derek<br>
<br>
</div></div></blockquote></div><br></div></div>