<html><head><meta http-equiv="content-type" content="text/html; charset=utf-8"></head><body style="overflow-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;"><div><br></div>  I am concerned this is not good advice being provided.<div><br></div><div>  Let's back up and look more closely at your use case.</div><div><br></div><div>  * What is the ratio of new nonzero locations added compared to the initial number of nonzeros for your code, for each of your iterations?</div><div><br></div><div>  * Is it possible for many iterations, no or very few new nonzeros are being added?</div><div><br></div><div>  * Are many previous nonzero values becoming zero (or unneeded) later? Again as a ratio compared to the initial number of nonzeros?</div><div><br></div><div>  * Can you quantify the difference in time between initially filling the matrix and then refilling it using the reset preallocation and not using the reset preallocation?</div><div><br></div><div>The effect you report that resetting the preallocation results in slower code is possible if relatively few additional nonzero locations are being created.</div><div><br><div><br></div><div>  After a matrix is assembled with a given nonzero structure (regardless of how it was filled it, using preallocation or not), setting nonzero values into new locations will be slow due to needing to do possibly many mallocs() (as much as one for each new nonzero introduced). Resetting the initially provided preallocation removes the need for all the new mallocs(), but at the expense of needing additional bookkeeping while setting the values. If you did not preallocate originally, then there is no way to prevent the additional mallocs() the second time through, so if you never preallocate but need to add many new nonzero locations adding the new nonzeros will be time consuming; hence in that situation providing the initial preallocation (taking into account all future nonzeros appearing) will pay off in possibly a big way.</div><div><br></div><div>I will look at the bug you report for when <span style="font-family: monospace;">MatResetPreallocation()</span>is called before the first matrix assembly and see if it can be fixed.</div><div><br></div><div>Barry<br><div><br></div><div><br><div><br><blockquote type="cite"><div>On Jun 18, 2023, at 9:01 AM, Junchao Zhang <junchao.zhang@gmail.com> wrote:</div><br class="Apple-interchange-newline"><div><div dir="ltr">yeah,  see a C example at <a href="https://gitlab.com/petsc/petsc/-/blob/main/src/mat/tests/ex259.c">https://gitlab.com/petsc/petsc/-/blob/main/src/mat/tests/ex259.c</a><div><br></div><div>I guess you can code in this outline with petsc-3.19</div><div><br></div><div><code>MatCreate</code><br></div><div><code><code>MatSetSizes</code><br></code></div><div><code><code><code>MatSetFromOptions</code><br></code></code></div><div><p><font face="monospace">iteration-loop start</font></p><p><font face="monospace">        MatResetPreallocation(A,...)</font></p><p><font face="monospace">        Fill Matrix A with MatSetValues</font></p><p><font face="monospace">        MatAssemblyXX(A_WAVE,MAT_FINAL_ASSEMBLY,IERR)</font></p><p><font face="monospace"> iteration-loop end</font><br></p><div><br clear="all"><div><div dir="ltr" class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr">--Junchao Zhang</div></div></div><br></div></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Sun, Jun 18, 2023 at 7:37 AM Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div dir="ltr">On Sun, Jun 18, 2023 at 4:26 AM Karsten Lettmann <<a href="mailto:karsten.lettmann@uni-oldenburg.de" target="_blank">karsten.lettmann@uni-oldenburg.de</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
  
    
  
  <div><p>Dear Matthew,</p><p><br>
    </p><p>thanks for you help.</p><p><br>
    </p><p>1) I tested your suggestion to pass NULL as the arguments for the
      MatXAIJSetPreallocation.</p><p>So the old version was:</p><p>CALL
MatCreateMPIAIJ(MPI_GROUP,N_local,N_local,N_global,N_global,0,DNNZ,0,ONNZ,A,IERR)<br>
    </p><p>And after you suggestion it is now:</p><p>  CALL MatCreate(MPI_GROUP,A,IERR)<br>
        CALL MatSetType(A,MATAIJ,IERR)<br>
        CALL MatSetSizes(A,N_local,N_local,N_global,N_global,IERR)<br>
        CALL
MatXAIJSetPreallocation(A,1,DNNZ,ONNZ,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,IERR)<br>
    </p><p><br>
    </p><p>  Setting block-size = 1.<br>
    </p><p><br>
    </p><p>2) Concerning the error with MatResetPreallocation:</p><p>We have an iterative loop, in which the matrix A is filled very
      often with different non-zero structure.<br>
      Further, I read in the manual pages that due to performance
      issues, one should preallocate enough space, as operations as
      matsetvalues might be time consuming due to additional <br>
      on-demand allocations.</p><p><br>
    </p><p>So I did the following coding in principle:</p><p><br>
    </p><p>    Set matrix A the first time with preallocation</p><p><br>
    </p><p>    iteration-loop start</p><p>        MatResetPreallocation(A,...)</p><p>        MatZeroEntries (A)<br>
    </p><p>        Fill Matrix A</p><p>        MatAssemblyXX(A_WAVE,MAT_FINAL_ASSEMBLY,IERR)</p><p><br>
    </p><p>    iteration-loop end<br>
    </p><p><br>
    </p><p>With these settings, the code run with 2 CPU.<br>
      But with 1 CPU I got an error, which was in MatResetPreallocation.<br>
      I could not understand, why the above code works with 2 CPU but
      not with 1 CPU.</p><p>At the moment, I believe the reason for this error seems to be a
      pre-check, that is done in SeqAIJ but not in MPIAIJ fo a valid and
      present matrix A.</p><p>(Now, an image is included showing the codings of   :<br>
<a href="https://petsc.org/release/src/mat/impls/aij/seq/aij.c.html#MatResetPreallocation_SeqAIJ" target="_blank">https://petsc.org/release/src/mat/impls/aij/seq/aij.c.html#MatResetPreallocation_SeqAIJ</a><br>
<a href="https://petsc.org/release/src/mat/impls/aij/mpi/mpiaij.c.html#MatResetPreallocation_MPIAIJ" target="_blank">https://petsc.org/release/src/mat/impls/aij/mpi/mpiaij.c.html#MatResetPreallocation_MPIAIJ</a>  
      )    <br>
    </p><p><br>
    </p><p><span id="cid:188ce81459a2eb3091d1"><fgcfkbigjgelpmgo.png></span></p><p><br>
    </p><p><br>
    </p><p>So, it seems for me at the moment, that the first
      MatResetPreallocation (when the iteration loop is entered the
      first time) is done on an not-assembled matrix A.<br>
      So for one CPU I got an error, while 2 CPUs seem to have been more
      tolerant.<br>
      (I'm not sure, if this interpretation is correct.)<br>
    </p><p><br>
    </p><p>So, I changed the coding in that way, that I check the assembly
      status before the preallocation.</p><p><br>
    </p><p>Using the coding:</p><p>    CALL MatAssembled(A,A_assembled,ierr)<br>
          IF (A_assembled .eqv. PETSC_TRUE) then<br>
              CALL MatResetPreallocation(A,ierr)<br>
          ENDIF    <br>
    </p><p>then worked for 1 and 2 CPU.<br>
    </p><p><br>
    </p><p><br>
    </p><p>3) There was another finding, which hopefully is correct.<br>
    </p><p><br>
    </p><p>Actually, I did this MatResetPreallocation to have a better
      performance when filling the matrix A later, as was suggested on
      the manual pages.<br>
      <br>
    </p><p>However, I found (if I did nothing wrong) that this
      MatResetPreallocation was much more time consuming than the
      additional (and unwanted) allocations done during filling the
      matrix.<br>
      <br>
    </p><p>So, in the end, my code seems to be faster, when <b>not</b>
      doing this in the iteration loop:<br>
    </p><p>    CALL MatAssembled(A,A_assembled,ierr)<br>
          IF (A_assembled .eqv. PETSC_TRUE) then<br>
              CALL MatResetPreallocation(A,ierr)<br>
          ENDIF    </p><p><br>
    </p><p>As I told you, I'm a beginner to PETSC and I do not know, if I
      have done it correctly???</p></div></blockquote><div>I think this may be correct now. We have rewritten Mat so that inserting values is much more efficient, and</div><div>can be done online, so preallocation is not really needed anymore. It is possible that this default mechanism</div><div>is faster than the old preallocation.</div><div><br></div><div>I would try the code without preallocation, using the latest release, and see how it performs.</div><div><br></div><div>  Thanks,</div><div><br></div><div>    Matt </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div><p>Best, Karsten</p>
    <br>
    <blockquote type="cite">
      <div>
        <div dir="ltr">
          <div class="gmail_quote">
            <div>The arguments are a combination of the AIJ and SBAIJ
              arguments. You can just pass NULL for the SBAIJ args.<br>
            </div>
            <blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
              <div><p>Then I  ran into issues with Resetpreallocation, that
                  I do not understand.</p>
              </div>
            </blockquote>
            <div>Can you send the error?</div>
            <div><br>
            </div>
            <div>  Thanks,</div>
            <div><br>
            </div>
            <div>    Matt <br>
            </div>
            <blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
              <div><p>I want this, because we have an iterative procedure,
                  where the matrix A_wave and its non-zero structure are
                  changing very often.<br>
                </p><p>I try to find the reason for my problem.</p><p><br>
                </p><p><br>
                </p><p>I really thank you for your answer, that helped me to
                  understand things a bit.<br>
                </p><p><br>
                </p><p>I wish you all the best, Karsten<br>
                </p><p><br>
                </p><p><br>
                </p>
                <div>Am 15.06.23 um 16:51 schrieb Matthew Knepley:<br>
                </div>
                <blockquote type="cite">
                  <div style="border:1pt solid rgb(0,51,51);padding:2pt"><p style="line-height:10pt;background:rgb(213,234,255)"><span>ACHTUNG!</span><span> Diese E-Mail kommt von Extern!
                        <span style="color:rgb(255,0,0)">WARNING!</span>
                        This email originated off-campus.<br>
                      </span></p>
                  </div>
                  <div>
                    <div dir="ltr">
                      <div dir="ltr">On Thu, Jun 15, 2023 at 8:32 AM
                        Karsten Lettmann <<a href="mailto:karsten.lettmann@uni-oldenburg.de" target="_blank">karsten.lettmann@uni-oldenburg.de</a>>
                        wrote:<br>
                      </div>
                      <div class="gmail_quote">
                        <blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
                          Dear all,<br>
                          <br>
                          <br>
                          I'm quite new to PETSC. So I hope the
                          following questions are not too <br>
                          stupid.<br>
                          <br>
                          <br>
                          1) We have a (Fortran) code, that we want to
                          update from an older PETSC <br>
                          version (petsc.2.3.3-p16) to a newer version.<br>
                          <br>
                          Inside the old code, for creating matrices A,
                          there are function calls <br>
                          of the from:<br>
                          MatCreateMPIAIJ<br>
                          <br>
                          In the reference page for this old version it
                          says:<br>
                          When calling this routine with a single
                          process communicator, a matrix <br>
                          of type SEQAIJ is returned.<br>
                          <br>
                          So I assume the following behavior of this old
                          routine:<br>
                          - for N_proc == 1:<br>
                              a matrix of type SEQAIJ is returned.<br>
                          <br>
                          - for N_proc > 1:<br>
                              a matrix of type MPIAIJ is returned<br>
                          <br>
                          <br>
                          <br>
                          2a) So, in the new code, we want to have a
                          similar behavior.<br>
                          <br>
                          I found that this function is not present any
                          more in the newer PETSC <br>
                          versions.<br>
                          <br>
                          Instead, one might use: MatCreateAIJ(….)<br>
                          ( <a href="https://petsc.org/release/manualpages/Mat/MatCreateAIJ/" rel="noreferrer" target="_blank">
https://petsc.org/release/manualpages/Mat/MatCreateAIJ/</a> )<br>
                          <br>
                          If I understand the reference page of the
                          function correctly, then, <br>
                          actually, a similar behavior should be
                          expected:<br>
                          <br>
                          - for N_proc == 1:<br>
                              a matrix of type SEQAIJ is returned.<br>
                          <br>
                          - for N_proc > 1:<br>
                              a matrix of type MPIAIJ is returned<br>
                          <br>
                          <br>
                          2b) However, on the reference page, there is
                          the note:<br>
                          <br>
                          It is recommended that one use the
                          MatCreate(), MatSetType() and/or <br>
                          MatSetFromOptions(), MatXXXXSetPreallocation()
                          paradigm instead of this <br>
                          routine directly.<br>
                          <br>
                          So, if I want the behavior above, it is
                          recommended to code it like <br>
                          this, isn't it:<br>
                          <br>
                          If (N_Proc == 1)<br>
                          <br>
                               MatCreate(.. ,A ,...)<br>
                               MatSetType(…,A, MATSEQAIJ,..)<br>
                               MatSetSizes(…,A, ..)<br>
                               MatSeqAIJSetPreallocation(,...A,...)<br>
                          <br>
                          else<br>
                          <br>
                               MatCreate(.. ,A ,...)<br>
                               MatSetType(…,A, MATMPIAIJ,..)<br>
                               MatSetSizes(…,A, ..)<br>
                               MatMPIAIJSetPreallocation(,...A,...)<br>
                        </blockquote>
                        <div><br>
                        </div>
                        <div>You can use</div>
                        <div><br>
                        </div>
                        <div>  MatCreate(comm, &A);</div>
                        <div>  MatSetType(A, MATAIJ);</div>
                        <div>  MatSetSizes(A, ...);</div>
                        <div>  MatXAIJSetPreallocation(A, ...);</div>
                        <div><br>
                        </div>
                        <div>We recommend this because we would like to
                          get rid of the convenience functions that</div>
                        <div>wrap up exactly this code.</div>
                        <div><br>
                        </div>
                        <div>  Thanks,</div>
                        <div><br>
                        </div>
                        <div>     Matt</div>
                        <div> </div>
                        <blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
                          end<br>
                          <br>
                          <br>
                          <br>
                          3) So my questions are:<br>
                          <br>
                          - Is my present understanding correct?<br>
                          <br>
                          If  yes:<br>
                          <br>
                          - Why might using MatCreateAIJ(….) for my case
                          not be helpful?<br>
                          <br>
                          - So, why is it recommended to use the way 2b)
                          instead of this <br>
                          MatCreateAIJ(….) ?<br>
                          <br>
                          <br>
                          Best, Karsten<br>
                          <br>
                          <br>
                          <br>
                          <br>
                          -- <br>
                          ICBM<br>
                          Section: Physical Oceanography<br>
                          Universitaet Oldenburg<br>
                          Postfach 5634<br>
                          D-26046 Oldenburg<br>
                          Germany<br>
                          <br>
                          Tel:    +49 (0)441 798 4061<br>
                          email:  <a href="mailto:karsten.lettmann@uni-oldenburg.de" target="_blank">karsten.lettmann@uni-oldenburg.de</a><br>
                          <br>
                        </blockquote>
                      </div>
                      <br clear="all">
                      <div><br>
                      </div>
                      <span class="gmail_signature_prefix">-- </span><br>
                      <div dir="ltr" class="gmail_signature">
                        <div dir="ltr">
                          <div>
                            <div dir="ltr">
                              <div>
                                <div dir="ltr">
                                  <div>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><br>
                                  </div>
                                  <div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
                                  </div>
                                </div>
                              </div>
                            </div>
                          </div>
                        </div>
                      </div>
                    </div>
                  </div>
                </blockquote>
                <pre cols="72">-- 
ICBM
Section: Physical Oceanography
Universitaet Oldenburg
Postfach 5634
D-26046 Oldenburg
Germany

Tel:    +49 (0)441 798 4061
email:  <a href="mailto:karsten.lettmann@uni-oldenburg.de" target="_blank">karsten.lettmann@uni-oldenburg.de</a></pre>
              </div>
            </blockquote>
          </div>
          <br clear="all">
          <div><br>
          </div>
          <span class="gmail_signature_prefix">-- </span><br>
          <div dir="ltr" class="gmail_signature">
            <div dir="ltr">
              <div>
                <div dir="ltr">
                  <div>
                    <div dir="ltr">
                      <div>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><br>
                      </div>
                      <div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
                      </div>
                    </div>
                  </div>
                </div>
              </div>
            </div>
          </div>
        </div>
      </div>
    </blockquote>
    <pre cols="72">-- 
ICBM
Section: Physical Oceanography
Universitaet Oldenburg
Postfach 5634
D-26046 Oldenburg
Germany

Tel:    +49 (0)441 798 4061
email:  <a href="mailto:karsten.lettmann@uni-oldenburg.de" target="_blank">karsten.lettmann@uni-oldenburg.de</a></pre>
  </div>

</blockquote></div><br clear="all"><div><br></div><span class="gmail_signature_prefix">-- </span><br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>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><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>
</blockquote></div>
</div></blockquote></div><br></div></div></div></body></html>