[petsc-dev] broken DMGetInterpolation_DA !!!!

[Jed Brown]

On Thu, Jul 21, 2011 at 14:48, Barry Smith <bsmith at mcs.anl.gov> wrote:

>
>
>  Whoever wrote the DMGetInterpolation_DA_Q1() for when coordinates are
> provided in DA
>
>   1) it is completely broken for the periodic case!
>

There is a serious semantic problem with periodic non-regular coordinates:
we don't know how big the wrap element is. There is no way to determine this
automatically, and there is no way to store it.


>
>   2) it doesn't actually use the coordinates provided by the user in the DA
> ; it uses uniform coordinates it cooks up on the fly
>
> /* compute local coordinate arrays */
>    nxi  = ratioi + 1;
>    neta = ratioj + 1;
>    ierr = PetscMalloc(sizeof(PetscScalar)*nxi,&xi);CHKERRQ(ierr);
>    ierr = PetscMalloc(sizeof(PetscScalar)*neta,&eta);CHKERRQ(ierr);
>    for (li=0; li<nxi; li++) {
>      xi[li] = -1.0 + (PetscScalar)li*(2.0/(PetscScalar)(nxi-1));
>    }
>    for (lj=0; lj<neta; lj++) {
>      eta[lj] = -1.0 + (PetscScalar)lj*(2.0/(PetscScalar)(neta-1));
>    }
>

These are reference coordinates. They get mapped using a Q1 local coordinate
mapping. Dave wrote this code, but it was tested to give second order
accuracy for non-uniform mapped grids.

Are you using non-uniform periodic grids? If you have regular periodic
grids, but you still want to store coordinates explicitly, then we could
write some special case to detect that you set coordinates that happened to
be uniform, therefore we can guess that the wrap cell should also be
uniform, and the code for no provided coordinates would be correct. What do
you want to do?
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20110721/033ad421/attachment.html>