[petsc-users] Matrix assembly in SNES problem

[Barry Smith]

> On Jul 22, 2020, at 9:35 AM, Adolfo Rodriguez <adantra at gmail.com> wrote:
> 
> Hong,
> 
> Thanks for your email. My question now is how to handle situations where the matrix structure changes during the iteration process?

   The routine you are calling MatSeqAIJSetPreallocationCSR() is fine. It allows you to put in new CSR data structures with different nonzero patters. 

    See my previous email.

   Barry

> Adolfo 
> 
> On Wed, Jul 22, 2020 at 8:32 AM Zhang, Hong <hongzhang at anl.gov <mailto:hongzhang at anl.gov>> wrote:
> 
> 
> > On Jul 22, 2020, at 7:03 AM, Adolfo Rodriguez <adantra at gmail.com <mailto:adantra at gmail.com>> wrote:
> > 
> > I am trying to replace the non-linear solver in a flow simulation problem where the matrix sparsity can change during the iterations. I tried successfully to create the matrix within the FormJacobian function but I have a memory leak which I don't know how to fix.
> > 
> > My FormJacobian function is something like this:
> > 
> > .. FormJacobian( ..){
> >      matrixconstruction(); // returns csr matrix in c++ regular vectors
> >      MatSeqAIJSetPreallocationCSR ...;
> >      //loop over the entries and 
> >      MatSetValues;
> >      MatAssemblyBegin();
> >      MatAssemblyEnd();
> >      return 0;
> > }
> 
> Yes, the leak is caused by the new matrix created every time FormJacobian is called. The working matrix for SNES needs to be allocated only once and passed to SNES by using
> 
> SNESSetJacobian(SNES snes,Mat Amat,Mat Pmat,PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
> 
> In your FormJacobian(SNES snes,Vec x,Mat A,Mat B,void *ctx), you just need to copy the csr matrix to the matrix A, which has already been allocated.
> 
> The working matrix should be destroyed after SNESSolve().
> 
> Hong (Mr.)
> 
> 
> > This approach works for small problems but, when dealing with large ones I observe a memory leak. I am pretty sure that It is related to the fact that I am creating new matrices every time FormJacobian is called but I don't know how to destroy it between iterations.
> > 
> > Any suggestions?
> > 
> > Regards,
> > 
> > Adolfo
> 

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20200722/1c1201d1/attachment.html>