<div dir="ltr">Hi Barry,<div><br></div><div>I compiled with mpich configured using --enable-g=meminit to get rid of MPI errors in Valgrind. Doing this reduced the number of errors to 2. I have attached the Valgrind output.</div><div><br></div><div>I'm using GAMG+GMRES for in each linear iteration of SNES. The linear solver converges with CONVERGED_RTOL for the first 6 iterations and with CONVERGED_STEP_LENGTH after that. I'm still very confused about why this is happening. Any thoughts/ideas?</div><div><br></div><div>Thanks,</div><div>Harshad</div></div><div class="gmail_extra"><br><div class="gmail_quote">On Thu, Sep 8, 2016 at 11:26 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"><br>
  Install your MPI with --download-mpich as a PETSc ./configure option, this will eliminate all the MPICH valgrind errors. Then send as an attachment the resulting valgrind file.<br>
<br>
  I do not 100 % trust any code that produces such valgrind errors.<br>
<span class="HOEnZb"><font color="#888888"><br>
   Barry<br>
</font></span><div class="HOEnZb"><div class="h5"><br>
<br>
<br>
> On Sep 8, 2016, at 10:12 PM, Harshad Sahasrabudhe <<a href="mailto:hsahasra@purdue.edu">hsahasra@purdue.edu</a>> wrote:<br>
><br>
> Hi Barry,<br>
><br>
> Thanks for the reply. My code is in C. I ran with Valgrind and found many "Conditional jump or move depends on uninitialized value(s)", "Invalid read" and "Use of uninitialized value" errors. I think all of them are from the libraries I'm using (LibMesh, Boost, MPI, etc.). I'm not really sure what I'm looking for in the Valgrind output. At the end of the file, I get:<br>
><br>
> ==40223== More than 10000000 total errors detected.  I'm not reporting any more.<br>
> ==40223== Final error counts will be inaccurate.  Go fix your program!<br>
> ==40223== Rerun with --error-limit=no to disable this cutoff.  Note<br>
> ==40223== that errors may occur in your program without prior warning from<br>
> ==40223== Valgrind, because errors are no longer being displayed.<br>
><br>
> Can you give some suggestions on how I should proceed?<br>
><br>
> Thanks,<br>
> Harshad<br>
><br>
> On Thu, Sep 8, 2016 at 1:59 PM, Barry Smith <<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>> wrote:<br>
><br>
>    This is very odd. CONVERGED_STEP_LENGTH for KSP is very specialized and should never occur with GMRES.<br>
><br>
>    Can you run with valgrind to make sure there is no memory corruption? <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/<wbr>documentation/faq.html#<wbr>valgrind</a><br>
><br>
>    Is your code fortran or C?<br>
><br>
>    Barry<br>
><br>
> > On Sep 8, 2016, at 10:38 AM, Harshad Sahasrabudhe <<a href="mailto:hsahasra@purdue.edu">hsahasra@purdue.edu</a>> wrote:<br>
> ><br>
> > Hi,<br>
> ><br>
> > I'm using GAMG + GMRES for my Poisson problem. The solver converges with KSP_CONVERGED_STEP_LENGTH at a residual of 9.773346857844e-02, which is much higher than what I need (I need a tolerance of at least 1E-8). I am not able to figure out which tolerance I need to set to avoid convergence due to CONVERGED_STEP_LENGTH.<br>
> ><br>
> > Any help is appreciated! Output of -ksp_view and -ksp_monitor:<br>
> ><br>
> >     0 KSP Residual norm 3.121347818142e+00<br>
> >     1 KSP Residual norm 9.773346857844e-02<br>
> >   Linear solve converged due to CONVERGED_STEP_LENGTH iterations 1<br>
> > KSP Object: 1 MPI processes<br>
> >   type: gmres<br>
> >     GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement<br>
> >     GMRES: happy breakdown tolerance 1e-30<br>
> >   maximum iterations=10000, initial guess is zero<br>
> >   tolerances:  relative=1e-08, absolute=1e-50, divergence=10000<br>
> >   left preconditioning<br>
> >   using PRECONDITIONED norm type for convergence test<br>
> > PC Object: 1 MPI processes<br>
> >   type: gamg<br>
> >     MG: type is MULTIPLICATIVE, levels=2 cycles=v<br>
> >       Cycles per PCApply=1<br>
> >       Using Galerkin computed coarse grid matrices<br>
> >   Coarse grid solver -- level ------------------------------<wbr>-<br>
> >     KSP Object:    (mg_coarse_)     1 MPI processes<br>
> >       type: preonly<br>
> >       maximum iterations=1, initial guess is zero<br>
> >       tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>
> >       left preconditioning<br>
> >       using NONE norm type for convergence test<br>
> >     PC Object:    (mg_coarse_)     1 MPI processes<br>
> >       type: bjacobi<br>
> >         block Jacobi: number of blocks = 1<br>
> >         Local solve is same for all blocks, in the following KSP and PC objects:<br>
> >         KSP Object:        (mg_coarse_sub_)         1 MPI processes<br>
> >           type: preonly<br>
> >           maximum iterations=1, initial guess is zero<br>
> >           tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>
> >           left preconditioning<br>
> >           using NONE norm type for convergence test<br>
> >         PC Object:        (mg_coarse_sub_)         1 MPI processes<br>
> >           type: lu<br>
> >             LU: out-of-place factorization<br>
> >             tolerance for zero pivot 2.22045e-14<br>
> >             using diagonal shift on blocks to prevent zero pivot [INBLOCKS]<br>
> >             matrix ordering: nd<br>
> >             factor fill ratio given 5, needed 1.91048<br>
> >               Factored matrix follows:<br>
> >                 Mat Object:                 1 MPI processes<br>
> >                   type: seqaij<br>
> >                   rows=284, cols=284<br>
> >                   package used to perform factorization: petsc<br>
> >                   total: nonzeros=7726, allocated nonzeros=7726<br>
> >                   total number of mallocs used during MatSetValues calls =0<br>
> >                     using I-node routines: found 133 nodes, limit used is 5<br>
> >           linear system matrix = precond matrix:<br>
> >           Mat Object:           1 MPI processes<br>
> >             type: seqaij<br>
> >             rows=284, cols=284<br>
> >             total: nonzeros=4044, allocated nonzeros=4044<br>
> >             total number of mallocs used during MatSetValues calls =0<br>
> >               not using I-node routines<br>
> >       linear system matrix = precond matrix:<br>
> >       Mat Object:       1 MPI processes<br>
> >         type: seqaij<br>
> >         rows=284, cols=284<br>
> >         total: nonzeros=4044, allocated nonzeros=4044<br>
> >         total number of mallocs used during MatSetValues calls =0<br>
> >           not using I-node routines<br>
> >   Down solver (pre-smoother) on level 1 ------------------------------<wbr>-<br>
> >     KSP Object:    (mg_levels_1_)     1 MPI processes<br>
> >       type: chebyshev<br>
> >         Chebyshev: eigenvalue estimates:  min = 0.195339, max = 4.10212<br>
> >       maximum iterations=2<br>
> >       tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>
> >       left preconditioning<br>
> >       using nonzero initial guess<br>
> >       using NONE norm type for convergence test<br>
> >     PC Object:    (mg_levels_1_)     1 MPI processes<br>
> >       type: sor<br>
> >         SOR: type = local_symmetric, iterations = 1, local iterations = 1, omega = 1<br>
> >       linear system matrix = precond matrix:<br>
> >       Mat Object:      ()       1 MPI processes<br>
> >         type: seqaij<br>
> >         rows=9036, cols=9036<br>
> >         total: nonzeros=192256, allocated nonzeros=192256<br>
> >         total number of mallocs used during MatSetValues calls =0<br>
> >           not using I-node routines<br>
> >   Up solver (post-smoother) same as down solver (pre-smoother)<br>
> >   linear system matrix = precond matrix:<br>
> >   Mat Object:  ()   1 MPI processes<br>
> >     type: seqaij<br>
> >     rows=9036, cols=9036<br>
> >     total: nonzeros=192256, allocated nonzeros=192256<br>
> >     total number of mallocs used during MatSetValues calls =0<br>
> >       not using I-node routines<br>
> ><br>
> > Thanks,<br>
> > Harshad<br>
><br>
><br>
<br>
</div></div></blockquote></div><br></div>