[petsc-users] petsc4py DA with dof>1 confusion

[Lisandro Dalcin]

On 4 July 2011 14:25, Matthew Emmett <memmett at unc.edu> wrote:
> Hi everyone,
>
> I am trying to use petsc4py in a fairly straight forward finite volume
> shallow-water code (the PDE has three unknowns). I have a structured
> grid and am using a DA to communicate ghost cells.  When I use dof=1,
> everything works as I would expect.  However, when I try to use dof=3,
> I am confused by the results.  I have probably mis-understood
> something -- any suggestions would be appreciated.
>
> For example (see attached script), first I create the DA ('dof' and
> 'width' are defined):
>
>  da = PETSc.DA().create(dim=2, sizes=[8, 8], dof=dof,
>                         boundary_type='periodic',
>                         stencil_width=width,
>                         stencil_type='box')
>
> Next, I create global and local vectors:
>
>  gvec = da.createGlobalVec()
>  lvec = da.createLocalVector()
>
> Then, put something into the local vectors, and scatter:
>
>  with da.getVecArray(lvec) as va:
>    if dof == 1:
>      va.array[width:-width,width:-width] = MPI.COMM_WORLD.rank+1
>    else:
>      va.array[width:-width,width:-width,0] = MPI.COMM_WORLD.rank+1
>  da.localToGlobal(lvec, gvec)
>
> Finally, get the local vector again and print:
>
>  da.globalToLocal(gvec, lvec)
>  with da.getVecArray(lvec) as va:
>    if dof == 1:
>      print va.array[:,:]
>    else:
>      print va.array[:,:,0]
>
> With 'dof=1' I get the following on the second processor:
>
> [[ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]
>  [ 1.  2.  2.  2.  2.  1.]]
>
> With 'dof=3' I get the following on the second processor:
>
> [[ 0.  0.  0.  0.  0.  0.]
>  [ 0.  0.  0.  0.  2.  0.]
>  [ 0.  0.  0.  0.  2.  0.]
>  [ 0.  0.  0.  2.  2.  0.]
>  [ 0.  0.  0.  2.  2.  0.]
>  [ 0.  0.  0.  2.  2.  0.]
>  [ 0.  0.  0.  2.  2.  0.]
>  [ 0.  0.  0.  2.  2.  2.]
>  [ 0.  0.  0.  2.  2.  2.]
>  [ 0.  0.  0.  0.  0.  2.]]
>
> I was expecting to see the same as the 'dof=1' case, since I am
> setting (and printing) va.array[:,:,0].
>
>
> Again, any suggestions would be appreciated.  Thanks,
> Matt
>
> PS, thanks to the developers of PETSc and petsc4py for all your hard
> work.  I had the pleasure of meeting Jed and Matt at KAUST recently.
>

Mmm, it seems that this functionality is broken for the case of
boundary_type='periodic'. I'll take a look.


-- 
Lisandro Dalcin
---------------
CIMEC (INTEC/CONICET-UNL)
Predio CONICET-Santa Fe
Colectora RN 168 Km 472, Paraje El Pozo
3000 Santa Fe, Argentina
Tel: +54-342-4511594 (ext 1011)
Tel/Fax: +54-342-4511169