On Thu, Nov 24, 2011 at 2:36 PM, Jed Brown <jedbrown at mcs.anl.gov> wrote: > On Thu, Nov 24, 2011 at 07:33, Dominik Szczerba <dominik at itis.ethz.ch> > wrote: >> >> MatCreateMPIAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, NTOT, >> NTOT, 500, PETSC_NULL, 500, PETSC_NULL, &A); >> >> (with 500 a big overshoot) >> >> then I do in a loop: >> >> 1. ierr = MatMPIAIJSetPreallocation(A, 500, PETSC_NULL, 500, >> PETSC_NULL); CHKERRQ(ierr); >> 2. ierr = MatZeroEntries(A); CHKERRQ(ierr); >> >> Line 1. does not matter the first time in the loop, info.mallocs is >> always 0, but if it is commented out, info.mallocs is non-zero the >> second time loop is executed. >> So I am confused, is the preallocation in MatCreateMPIAIJ not >> sufficient? Does it need to be refreshed? > > You must be changing the number of nonzeros in each iteration of the loop.
Yes, that is true, but 1) within the overshot preallocation margin 2) I do it at each execution of the loop, so I do not see why only the second time there are mallocs. I display mallocs directly before KSPSolve(). > I suggest assembling all possible entries (even the ones that happen to be > zero) in the first assembly. Then there will be space and you can reuse the > data structure. This will require a bigger change in my code. Is always calling MatMPIAIJSetPreallocation a bad (inefficient) alternative? Thanks Dominik
