<div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">On Wed, Oct 16, 2013 at 4:22 PM, Christophe Ortiz <span dir="ltr"><<a href="mailto:christophe.ortiz@ciemat.es" target="_blank">christophe.ortiz@ciemat.es</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><div class="gmail_extra"><br><br><div class="gmail_quote">
<div><div>On Wed, Oct 16, 2013 at 4:10 PM, Christophe Ortiz <span dir="ltr"><<a href="mailto:christophe.ortiz@ciemat.es" target="_blank">christophe.ortiz@ciemat.es</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div dir="ltr"><br><div class="gmail_extra"><br><div class="gmail_quote">
<div>On Wed, Oct 16, 2013 at 3:36 PM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div>Christophe Ortiz <<a href="mailto:christophe.ortiz@ciemat.es" target="_blank">christophe.ortiz@ciemat.es</a>> writes:<br>
> I found where is the problem. I had<br>
> set TSARKIMEXSetFullyImplicit(ts,PETSC_TRUE);<br>
> I found it produces wrong solutions, either with ARKIMEX3, 1BEE or A2.<br>
> Solution becomes negative. Without this option, all ARKIMEX types work very<br>
> well and very fast.<br>
<br>
</div>I'm glad they work as IMEX methods (please check that the solutions are<br>
correct), but I would like to find out why they are not working in fully<br>
implicit mode. Is your problem similar to ex25? </blockquote><div><br></div></div><div>Yes. The only difference is that I put everything (diffusion and reaction terms) under IFunction and IJacobian.</div><div><br></div>
<div>My system is</div>
<div><br></div><div>u_t - alpha.u_xx + k.u.v = 0</div><div>v_t - beta.v_xx + k.u.v = 0</div><div> </div><div>with Dirichlet boundary conditions</div><div>u(x=0,t) = u(x=L,t) = cste</div><div>v(x=0,t) = v(x=L,t) = cste</div>
<div>
<div><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"> (If I remember<br>
correctly, that example uses a DAE formulation for boundary conditions<br>
and thus has trouble with some methods with explicit stages. There may<br>
be more going on.)<br></blockquote><div><br></div></div><div>Now you say it, I'm having some troubles. After checking that it worked very well with only diffusion terms with ARKIMEX, I introduced the reaction term. With ARKIMEX 1BEE/A2, it goes very fast up to t=0.27 s, then, steps start being rejected and then it goes very slowly, with small timesteps. See output below.</div>
<div>Do you think it is due to the fact that I use IFunction to define boundary conditions, which produces troubles with ARKIMEX methods ?</div><div>What would be the solution ?</div><div><br></div><div>Christophe</div><div>
<br></div></div></div></div></blockquote><div><br></div></div></div><div>I just tried with TSARKIMEXSetFullyImplicit(ts,PETSC_TRUE); It goes much better and reaches final time....but solution is clearly not correct. The initial gaussian distributions do not change. Like if diffusion was null.</div>
</div></div></div></blockquote><div><br></div><div>Hi,</div><div><br></div><div>I've tried something else. Since ARKIMEX is able to converge up to t~0.20s and gives the correct solution for this time, and then comes into trouble, I thought that it could come from the SNES.</div>
<div>Then, instead of using SNESNEWTONLS, I tried SNESNEWTONTR. I summarize my findings:</div><div><br></div><div>- With NEWTONTR, ARKIMEX (3, 1BEE or A2) converges.</div><div><br></div><div>-With NEWTONTR, ARKIMEX 3 converges much faster than 1BEE or A2 (much less timesteps).</div>
<div><br></div><div>- No conflicts between NEWTONTR and ARKIMEXFullyImplicit as with NEWTONLS. Give the same results.</div><div><br></div><div><br></div><div>Then I came back to SNESNEWTONLS</div><div><br></div><div><div>
- With NEWTONLS, ARKIMEX does not converges.</div><div><br></div><div>- When ARKIMEXFullyImplicit is used with NEWTONLS, it converges, but solution clearly incorrect. It seems it does nothing, initial distributions do not change.</div>
</div><div><br></div><div><br></div><div>And I played a bit with the tolerances:</div><div><br></div><div>- No influence of rtol (increasing or decreasing it). It stops at same point.</div><div><br></div><div>- Changing stol from 1e-8 to 1e-7 seems to ease as it goes up to longer time (6.98 s instead of 0.20) but then stops with an error:</div>
<div><br></div><div><div>[0]PETSC ERROR: --------------------- Error Message ------------------------------------</div><div>[0]PETSC ERROR: Floating point exception!</div><div>[0]PETSC ERROR: Vec entry at local location 450 is not-a-number or infinite at end of function: Parameter number 3!</div>
<div>...</div></div><div>...</div><div><div>[0]PETSC ERROR: VecValidValues() line 30 in src/vec/vec/interface/rvector.c</div><div>[0]PETSC ERROR: SNESComputeFunction() line 1998 in src/snes/interface/snes.c</div><div>[0]PETSC ERROR: SNESSolve_NEWTONLS() line 162 in src/snes/impls/ls/ls.c</div>
<div>[0]PETSC ERROR: SNESSolve() line 3636 in src/snes/interface/snes.c</div><div>[0]PETSC ERROR: TSStep_ARKIMEX() line 765 in src/ts/impls/arkimex/arkimex.c</div><div>[0]PETSC ERROR: TSStep() line 2458 in src/ts/interface/ts.c</div>
<div>[0]PETSC ERROR: TSSolve() line 2583 in src/ts/interface/ts.c</div><div>[0]PETSC ERROR: main() line 392 in src/ts/examples/tutorials/diffusion.c</div></div><div><br></div><div>Maybe all this will give you a hint of what's going on...</div>
<div>Saludos,</div><div>Christophe</div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div></div></div></div>