<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Nov 13, 2013 at 1:33 PM, Eric Chamberland <span dir="ltr"><<a href="mailto:Eric.Chamberland@giref.ulaval.ca" target="_blank">Eric.Chamberland@giref.ulaval.ca</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hi all,<br>
<br>
we are in the process of adding new functionalities to our parallel finite-element code and are facing new challenges in computing the number of non-zeros of a sparse matrix.  Since we use Petsc, we want to be able to give this information before the first assembly (to pass to Mat{MPI,SEQ}<u></u>AIJSetPreallocation).<br>

<br>
Willing to rewrite the nonzeros counting algorithm from scratch if necessary, we asked ourselves a number of questions:<br>
<br>
1) Do algorithms in general count the *exact* number of non-zeros?  or do they just count a maximum and let Petsc release the memory after the first assembly?<br>
<br>
2) How bad is it to reserve too-much non-zeros before the first assembly?<br>
<br>
3) Where could we find some hints/articles examples or information about that kind of algorithm (working in parallel of course) ?<br>
<br>
Here is the context that we have up to now: we can manage a mixture of interpolations (P1, P2, discontinuous interpolations, etc) and of dimensions (scalar, 3D, etc) for all unknowns in a global matrix (example: a P2-P1 [Taylor-Hood] mixed formulation for Stokes problem).<br>

<br>
Since these parameters are determined at run-time (in fact they are user parameters), we have to think of a general algorithm to compute the "mostly" exact number of non-zeros for each matrix lines.<br>
<br>
Right now, we have a working, but memory consuming algorithm which computes the *exact* number on non-zeros per line (we exploded the memory for a problem with a little more than 2^31 unknowns with it -- using 64 bits indices of course).  In fact, we manage to store the number of neighbors for each type of mesh components (ex: for a vertice, we have the following number of neighbors: vertices, edges, faces, and each type of elements).  It is costly and moreover will not handle the new functionalities we want in the code (XFEM for example).<br>
</blockquote><div><br></div><div>DMPlex will do it already.</div><div><br></div><div>If you are opposed to using that, you can use Jed's simple method which is to "fake" assemble without values,</div><div>putting all the (i,j) pairs into a hash table.</div>
<div><br></div><div>  Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Thanks a lot!<span class="HOEnZb"><font color="#888888"><br>
<br>
Eric<br>
</font></span></blockquote></div><br><br clear="all"><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener
</div></div>