[petsc-users] Memory corruption with two-dimensional array and PetscMemzero

[Matthew Knepley]

On Tue, May 20, 2014 at 7:12 AM, Christophe Ortiz <
christophe.ortiz at ciemat.es> wrote:

> I found another problem when using two-dimensional arrays defined using
> pointers of pointers.
>
> When I use a "classical" two-dimensional array defined by
>
> PetscScalar  array[dof][dof];
>

This declaration will use contiguous memory since its on the stack.


> and then build the Jacobian using
>
> ierr = MatSetValuesBlocked(*Jpre,1,&row,1,&col,&array[0][0],INSERT_VALUES);
>
> It works fine.
>
> The problem comes when I define the two-dimensional array as follows:
>
> PetscScalar  **array;
>

This one uses non-contiguous memory since its on the heap.


>  array = (PetscScalar**)malloc(sizeof(PetscScalar*) * dof);
>  for (k = 0; k < dof; k++){
>    array[k] = (PetscScalar*)malloc(sizeof(PetscScalar) * dof);
>   }
>
> When I pass it to MatSetValuesBlocked() there is a problem. Either Petsc
> complains because I am not passing it the right way or when it accepts it,
> wrong data is passed because the solution is not correct. Maybe Petsc
> expect dof*dof values and only dof are passed ?
>

You can only pass contiguous memory to MatSetValues().

   Matt


> How a two-dimensional array declared with pointers of pointers should be
> passed to MatSetValuesBlocked() ?
>
> Many thanks in advance.
> Christophe
>
>
>
>
> On Thu, May 15, 2014 at 2:08 AM, Jed Brown <jed at jedbrown.org> wrote:
>
>> Christophe Ortiz <christophe.ortiz at ciemat.es> writes:
>>
>> > Hi all,
>> >
>> > I am experiencing some problems of memory corruption with
>> PetscMemzero().
>> >
>> > I set the values of the Jacobian by blocks using MatSetValuesBlocked().
>> To
>> > do so, I use some temporary two-dimensional arrays[dof][dof] that I must
>> > reset at each loop.
>> >
>> > Inside FormIJacobian, for instance, I declare the following
>> two-dimensional
>> > array:
>> >
>> >    PetscScalar  diag[dof][dof];
>> >
>> > and then, to zero the array diag[][] I do
>> >
>> >    ierr = PetscMemzero(diag,dof*dof*sizeof(PetscScalar));
>>
>> Note that this can also be spelled
>>
>>   PetscMemzero(diag,sizeof diag);
>>
>> > Then, inside main(), once dof is determined, I allocate memory for diag
>> as
>> > follows:
>> >
>> >  diag = (PetscScalar**)malloc(sizeof(PetscScalar*) * dof);
>> >
>> >  for (k = 0; k < dof; k++){
>> >   diag[k] = (PetscScalar*)malloc(sizeof(PetscScalar) * dof);
>> >  }
>> > That is, the classical way to allocate memory using the pointer
>> notation.
>>
>> Note that you can do a contiguous allocation by creating a Vec, then use
>> VecGetArray2D to get 2D indexing of it.
>>
>
>


-- 
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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140520/b85bf640/attachment-0001.html>