[petsc-users] Memory violation in PetscFVLeastSquaresPseudoInverseSVD_Static

[Matthew Knepley]

On Thu, Oct 1, 2020 at 4:18 AM Pierre Seize <Pierre.Seize at onera.fr> wrote:

> Sure I'll try,
>
> I was thinking of using ls->B as the A matrix for dgelss, ls->work as
> work and ls->Binv as B. The result is then stored in ls->Binv, but in
> column-major format.
>
> Right now, the column-major result is transposed in
> PetscFVLeastSquaresPseudoInverseSVD_Static, and the row-major result is
> copied in the output in PetscFVComputeGradient_LeastSquares. I think
> it's because PetscFVLeastSquaresPseudoInverse_Static gives the result in
> row-major format.
>
> Would it be alright if I changed
> PetscFVLeastSquaresPseudoInverseSVD_Static so that the result would
> still be in column-major format ? I could include the result recopy in
> the if statement for example.
>
> Moreover, this would be to keep the compatibility with
> PetscFVLeastSquaresPseudoInverse_Static, but right now it is manually
> disabled (with useSVD = PETSC_TRUE), so I am worrying for nothing ?
>

As long as the layout is documented in the function manpage, I think that
change is fine. Nothing else uses
the code except the reconstruction right now.

  Thanks,

    Matt


> Pierre
>
>
> On 30/09/20 20:38, Jed Brown wrote:
> > Pierre Seize <Pierre.Seize at onera.fr> writes:
> >
> >> Hi,
> >>
> >> In PetscFVLeastSquaresPseudoInverseSVD_Static, there is
> >>     Brhs = work;
> >>     maxmn = PetscMax(m,n);
> >>     for (j=0; j<maxmn; j++) {
> >>       for (i=0; i<maxmn; i++) Brhs[i + j*maxmn] = 1.0*(i == j);
> >>     }
> >> where on the calling function, PetscFVComputeGradient_LeastSquares, we
> >> set the arguments m <= numFaces, n <= dim and work <= ls->work. The size
> >> of the work array is computed in PetscFVLeastSquaresSetMaxFaces_LS as:
> >>     ls->maxFaces = maxFaces;
> >>     m       = ls->maxFaces;
> >>     n       = dim;
> >>     nrhs    = ls->maxFaces;
> >>     minwork = 3*PetscMin(m,n) + PetscMax(2*PetscMin(m,n),
> >> PetscMax(PetscMax(m,n), nrhs)); /* required by LAPACK */
> > It's totally buggy because this formula is for the argument to dgelss,
> but the array is being used for a different purpose (to place Brhs).
> >
> >             WORK
> >
> >                       WORK is DOUBLE PRECISION array, dimension
> (MAX(1,LWORK))
> >                       On exit, if INFO = 0, WORK(1) returns the optimal
> LWORK.
> >
> >             LWORK
> >
> >                       LWORK is INTEGER
> >                       The dimension of the array WORK. LWORK >= 1, and
> also:
> >                       LWORK >= 3*min(M,N) + max( 2*min(M,N), max(M,N),
> NRHS )
> >                       For good performance, LWORK should generally be
> larger.
> >
> >                       If LWORK = -1, then a workspace query is assumed;
> the routine
> >                       only calculates the optimal size of the WORK
> array, returns
> >                       this value as the first entry of the WORK array,
> and no error
> >                       message related to LWORK is issued by XERBLA.
> >
> > There should be a separate allocation for Brhs and the work argument
> should be passed through to dgelss.
> >
> > The current code passes
> >
> >    tmpwork = Ainv;
> >
> > along to dgelss, but we don't know that it's the right size either.
> >
> >
> > Would you be willing to submit a merge request with your best attempt at
> fixing this.  I can help review and we'll get it into the 3.14.1 release?
> >
> >>     ls->workSize = 5*minwork; /* We can afford to be extra generous */
> >>
> >> In my example, the used size (maxmn * maxmn) is 81, and the actual size
> >> (ls->workSize) is 75, and therefore valgrind complains.
> >> Is it because I am missing something, or is it a bug ?
> >>
> >> Thanks
> >>
> >> Pierre Seize
>
>

-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20201001/4ca5d326/attachment.html>