<div dir="ltr">Hi Dmitry,<div><br></div><div>OK, good to hear we're seeing the same behavior for the example.</div><div><br></div><div>Regarding this comment:</div><div><br></div><div><br></div><div class="gmail_extra"><div class="gmail_quote"><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><div>libMesh needs to change the way configure extracts PETSc information -- configuration data were moved:</div><div>conf --> lib/petsc-conf</div><div>${PETSC_ARCH}/conf --> ${PETSC_ARCH}/lib/petsc-conf</div><div><br></div><div>At one point I started looking at m4/petsc.m4, but that got put on the back burner.  For now making the relevant symlinks by hand lets you configure and build libMesh with petsc/master.</div></div></div></blockquote><div><br></div><div><br></div><div>So do you suggest that the next step here is to build libmesh against petsc/master so that we can try your PETSc pull request that implements MatReset() to see if that gets this example working?</div><div><br></div><div>David</div><div><br></div><div><br></div><div><br></div><div> </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"><div dir="ltr"><div><div><div><div class="gmail_quote">On Mon Feb 23 2015 at 9:15:44 AM David Knezevic <<a href="mailto:david.knezevic@akselos.com" target="_blank">david.knezevic@akselos.com</a>> 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">Hi Dmitry,<div><br></div><div>Thanks very much for testing out the example.</div><div><br></div><div><span style="font-size:12.8000001907349px">examples/systems_of_equations/</span><span style="font-size:12.8000001907349px"><u></u>ex8 works fine for me in serial, but it fails for me if I run with more than 1 MPI process. Can you try it with, say, 2 or 4 MPI processes?</span><br></div><div><span style="font-size:12.8000001907349px"><br></span></div><div><span style="font-size:12.8000001907349px">If we need access to MatReset in libMesh to get this to work, I'll be happy to work on a libMesh pull request for that.</span></div><div><span style="font-size:12.8000001907349px"><br></span></div><div><span style="font-size:12.8000001907349px">David</span></div><div><span style="font-size:12.8000001907349px"><br></span></div></div><div class="gmail_extra"><br clear="all"><div><div><div dir="ltr"><span style="color:rgb(0,0,0);font-family:Calibri,sans-serif;font-size:12px"><pre cols="72">--</pre><pre cols="72">David J. Knezevic | CTO
Akselos | 17 Bay State Road | Boston, MA | 02215
Phone (office): <a href="tel:%2B1-857-265-2238" value="+18572652238" target="_blank">+1-857-265-2238</a><br>Phone (mobile): <a href="tel:%2B1-617-599-4755" value="+16175994755" target="_blank">+1-617-599-4755</a>
Web: <a href="http://www.akselos.com" target="_blank">http://www.akselos.com</a></pre></span></div></div></div></div><div class="gmail_extra">
<br><div class="gmail_quote">On Mon, Feb 23, 2015 at 10:08 AM, Dmitry Karpeyev <span dir="ltr"><<a href="mailto:karpeev@mcs.anl.gov" target="_blank">karpeev@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 dir="ltr">David,<br><br><div>What code are you running when you encounter this error?  I'm trying to reproduce it and</div><div>I tried examples/systems_of_equations/<u></u>ex8, but it ran for me, ostensibly to completion.</div><div><br></div><div>I have a small PETSc pull request that implements MatReset(), which passes a small PETSc test,</div><div>but libMesh needs some work to be able to build against petsc/master because of some recent</div><div>changes to PETSc. </div><span><font color="#888888"><div><br></div><div>Dmitry.</div></font></span><div><div><br><div class="gmail_quote">On Mon Feb 23 2015 at 7:17:06 AM David Knezevic <<a href="mailto:david.knezevic@akselos.com" target="_blank">david.knezevic@akselos.com</a>> 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">Hi Barry, hi Dmitry,<div><br></div><div>I set the matrix to BAIJ and back to AIJ, and the code got a bit further. But I now run into the error pasted below (Note that I'm now using "--with-debugging=1"):</div><div><br></div><div><div>PETSC ERROR: --------------------- Error Message ------------------------------<u></u>------------------------------<u></u>--</div></div></div><div dir="ltr"><div><div>PETSC ERROR: Petsc has generated inconsistent data</div></div></div><div dir="ltr"><div><div>PETSC ERROR: MPIAIJ Matrix was assembled but is missing garray</div></div></div><div dir="ltr"><div><div>PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html" target="_blank">http://www.mcs.anl.gov/petsc/<u></u>documentation/faq.html</a> for trouble shooting.</div><div>PETSC ERROR: Petsc Release Version 3.5.2, Sep, 08, 2014 </div><div>PETSC ERROR: ./example-dbg on a arch-linux2-c-debug named david-Lenovo by dknez Mon Feb 23 08:05:44 2015</div><div>PETSC ERROR: Configure options --with-shared-libraries=1 --with-debugging=1 --download-suitesparse=1 --download-parmetis=1 --download-blacs=1 --download-scalapack=1 --download-mumps=1 --download-metis --download-superlu_dist --prefix=/home/dknez/software/<u></u>libmesh_install/dbg_real/petsc --download-hypre</div><div>PETSC ERROR: #1 MatCreateColmap_MPIAIJ_<u></u>Private() line 361 in /home/dknez/software/petsc-3.<u></u>5.2/src/mat/impls/aij/mpi/<u></u>mpiaij.c</div><div>PETSC ERROR: #2 MatSetValues_MPIAIJ() line 538 in /home/dknez/software/petsc-3.<u></u>5.2/src/mat/impls/aij/mpi/<u></u>mpiaij.c</div><div>PETSC ERROR: #3 MatSetValues() line 1136 in /home/dknez/software/petsc-3.<u></u>5.2/src/mat/interface/matrix.c</div><div>PETSC ERROR: #4 add_matrix() line 765 in /home/dknez/software/libmesh-<u></u>src/src/numerics/petsc_matrix.<u></u>C</div><div>------------------------------<u></u>------------------------------<u></u>--------------</div></div><div><br></div><div>This occurs when I try to set some entries of the matrix. Do you have any suggestions on how I can resolve this?</div><div><br></div><div>Thanks!</div></div><div dir="ltr"><div>David</div></div><div dir="ltr"><div><br></div><div><br></div><div><br></div><div><br></div><div>On Sun, Feb 22, 2015 at 10:22 PM, Dmitry Karpeyev <span dir="ltr"><<a href="mailto:dkarpeev@gmail.com" target="_blank">dkarpeev@gmail.com</a>></span> wrote:<br></div><div class="gmail_extra"><div class="gmail_quote"><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><br><div class="gmail_quote"><span>On Sun Feb 22 2015 at 9:15:22 PM Barry Smith <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> 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"><br>
> On Feb 22, 2015, at 9:09 PM, David Knezevic <<a href="mailto:david.knezevic@akselos.com" target="_blank">david.knezevic@akselos.com</a>> wrote:<br>
><br>
> Hi Dmitry,<br>
><br>
> Thanks for the suggestion. I tried MatSetType(mat,MATMPIAIJ) followed by MatXAIJSetPreallocation(...), but unfortunately this still gives me the same error as before: "nnz cannot be greater than row length: local row 168 value 24 rowlength 0".<br>
><br>
> I gather that the idea here is that MatSetType builds a new matrix object, and then I should be able to pre-allocate for that new matrix however I like, right? Was I supposed to clear the matrix object somehow before calling MatSetType? (I didn't do any sort of clear operation.)<br>
<br>
  If the type doesn't change then MatSetType() won't do anything. You can try setting the type to BAIJ and then setting the type back to AIJ. This may/should clear out the matrix.<br></blockquote></span><div>Ah, yes.  If the type is the same as before it does quit early, but changing the type and then back will clear out and rebuild the matrix.   We need </div><div><span style="font-size:13.1999998092651px">something like MatReset() to do the equivalent thing. </span></div><span><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">
<br>
><br>
> As I said earlier, I'll make a dbg PETSc build, so hopefully that will help shed some light on what's going wrong for me.<br></blockquote></span><div>I think it's always a good idea to have a dbg build of PETSc when you doing things like these. </div><span><font color="#888888"><div><br></div><div>Dmitry. </div></font></span><div><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">
<br>
   Don't bother, what I suggested won't work.<br>
<br>
  Barry<br>
<br>
<br>
><br>
> Thanks,<br>
> David<br>
><br>
><br>
><br>
><br>
> On Sun, Feb 22, 2015 at 6:02 PM, Dmitry Karpeyev <<a href="mailto:dkarpeev@gmail.com" target="_blank">dkarpeev@gmail.com</a>> wrote:<br>
> David,<br>
> It might be easier to just rebuild the whole matrix from scratch: you would in effect be doing all that with disassembling and resetting the preallocation.<br>
> MatSetType(mat,MATMPIAIJ)<br>
> or<br>
> PetscObjectGetType((<u></u>PetscObjec<u></u>t)mat,&type);<br>
> MatSetType(mat,type);<br>
> followed by<br>
> MatXAIJSetPreallocation(...);<br>
> should do.<br>
> Dmitry.<br>
><br>
><br>
> On Sun Feb 22 2015 at 4:45:46 PM Barry Smith <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> wrote:<br>
><br>
>  Do not call for SeqAIJ matrix. Do not call before the first time you have preallocated and put entries in the matrix and done the MatAssemblyBegin/End()<br>
><br>
>   If it still crashes you'll need to try the debugger<br>
><br>
>   Barry<br>
><br>
> > On Feb 22, 2015, at 4:09 PM, David Knezevic <<a href="mailto:david.knezevic@akselos.com" target="_blank">david.knezevic@akselos.com</a>> wrote:<br>
> ><br>
> > Hi Barry,<br>
> ><br>
> > Thanks for your help, much appreciated.<br>
> ><br>
> > I added a prototype for MatDisAssemble_MPIAIJ:<br>
> > PETSC_INTERN PetscErrorCode MatDisAssemble_MPIAIJ(Mat);<br>
> ><br>
> > and I added a call to MatDisAssemble_MPIAIJ before MatMPIAIJSetPreallocation. However, I get a segfault on the call to MatDisAssemble_MPIAIJ. The segfault occurs in both serial and parallel.<br>
> ><br>
> > FYI, I'm using Petsc 3.5.2, and I'm not using a non-debug build (though I could rebuild PETSc in debug mode if you think that would help figure out what's happening here).<br>
> ><br>
> > Thanks,<br>
> > David<br>
> ><br>
> ><br>
> ><br>
> > On Sun, Feb 22, 2015 at 1:13 PM, Barry Smith <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> wrote:<br>
> ><br>
> >   David,<br>
> ><br>
> >    This is an obscure little feature of MatMPIAIJ,   each time you change the sparsity pattern before you call the MatMPIAIJSetPreallocation you need to call  MatDisAssemble_MPIAIJ(Mat mat).    This is a private PETSc function so you need to provide your own prototype for it above the function you use it in.<br>
> ><br>
> >   Let us know if this resolves the problem.<br>
> ><br>
> >    Barry<br>
> ><br>
> > We never really intended that people would call MatMPIAIJSetPreallocation() AFTER they had already used the matrix.<br>
> ><br>
> ><br>
> > > On Feb 22, 2015, at 6:50 AM, David Knezevic <<a href="mailto:david.knezevic@akselos.com" target="_blank">david.knezevic@akselos.com</a>> wrote:<br>
> > ><br>
> > > Hi all,<br>
> > ><br>
> > > I've implemented a solver for a contact problem using SNES. The sparsity pattern of the jacobian matrix needs to change at each nonlinear iteration (because the elements which are in contact can change), so I tried to deal with this by calling MatSeqAIJSetPreallocation and MatMPIAIJSetPreallocation during each iteration in order to update the preallocation.<br>
> > ><br>
> > > This seems to work fine in serial, but with two or more MPI processes I run into the error "nnz cannot be greater than row length", e.g.:<br>
> > > nnz cannot be greater than row length: local row 528 value 12 rowlength 0<br>
> > ><br>
> > > This error is from the call to<br>
> > > MatSeqAIJSetPreallocation(b-><u></u>B<u></u>,o_nz,o_nnz); in MatMPIAIJSetPreallocation_<u></u>MPIA<u></u>IJ.<br>
> > ><br>
> > > Any guidance on what the problem might be would be most appreciated. For example, I was wondering if there is a problem with calling SetPreallocation on a matrix that has already been preallocated?<br>
> > ><br>
> > > Some notes:<br>
> > > - I'm using PETSc via libMesh<br>
> > > - The code that triggers this issue is available as a PR on the libMesh github repo, in case anyone is interested: <a href="https://github.com/libMesh/libmesh/pull/460/" target="_blank">https://github.com/libMesh/<u></u>lib<u></u>mesh/pull/460/</a><br>
> > > - I can try to make a minimal pure-PETSc example that reproduces this error, if that would be helpful.<br>
> > ><br>
> > > Many thanks,<br>
> > > David<br>
> > ><br>
> ><br>
> ><br>
><br>
><br>
<br>
</blockquote></div></div></div></div>
</blockquote></div><br></div></div></blockquote></div></div></div></div>
</blockquote></div><br></div></blockquote></div></div></div></div></div>
</blockquote></div></div></div>