[petsc-users] Problem with interpolating projected field of quadratic order

[Matthew Knepley]

On Sat, Jun 18, 2016 at 9:26 AM, Matthew Knepley <knepley at gmail.com> wrote:

> On Sat, Jun 18, 2016 at 8:56 AM, Sander Arens <Sander.Arens at ugent.be>
> wrote:
>
>> Hello all,
>>
>> Recently I've been trying to run a fem problem using DMPlex, PetscFE,
>> etc. using a mesh created with gmsh. I used DMPlexReverseCell to get rid of
>> the negative jacobian determinants. However, I noticed another problem:
>> when I projected the coordinates on a quadratic fem-field and interpolated
>> them they didn't match the interpolated values of the original coordinate
>> field anymore.
>>
>> I attached a simple test to reproduce this. With linear interpolation I
>> see no problem, it's only when I start using quadratic interpolation.
>>
>> I have no idea why this is giving problems. Am I missing something simple
>> here or is this a bug?
>>
>
> Tests for this already exist
>
>
> https://bitbucket.org/petsc/petsc/src/996a00f1a1693e2d85044ed633fe1848cf357f2a/src/dm/impls/plex/examples/tests/ex3.c?at=master&fileviewer=file-view-default
>
> with run parameters here:
>
>
> https://bitbucket.org/petsc/petsc/src/996a00f1a1693e2d85044ed633fe1848cf357f2a/config/builder.py?at=master&fileviewer=file-view-default#builder.py-58
>
> and specifically for P2 triangles
>
>
> https://bitbucket.org/petsc/petsc/src/996a00f1a1693e2d85044ed633fe1848cf357f2a/config/builder.py?at=master&fileviewer=file-view-default#builder.py-72
>
> 1) Does the test run for you?
>
>   ./config/builder2.py check src/dm/impls/plex/examples/tests/ex3.c
>
> 2) If so, then could it be that you assume something about the order of
> functions?
>

I am having a hard time understanding what you are doing here:

  /* View the projected coordinates field */
  ierr = VecGetSubVector(u, fields[0], &subvec);CHKERRQ(ierr);
  ierr = PetscObjectSetName((PetscObject) subvec, fieldNames[0]);CHKERRQ(ierr);
  ierr = VecViewFromOptions(subvec, NULL,
"-projectedcoordinates_vec_view");CHKERRQ(ierr);
  ierr = VecRestoreSubVector(u, fields[0], &subvec);CHKERRQ(ierr);

  /* View the interpolated coordinates field */
  ierr = VecGetSubVector(u, fields[1], &subvec);CHKERRQ(ierr);
  ierr = PetscObjectSetName((PetscObject) subvec, fieldNames[1]);CHKERRQ(ierr);
  ierr = VecViewFromOptions(subvec, NULL,
"-interpolatedcoordinates_vec_view");CHKERRQ(ierr);
  ierr = VecRestoreSubVector(u, fields[1], &subvec);CHKERRQ(ierr);


So the values in 'u' are the FEM coefficients for the coordinate function
{x, y}. Why are you calling 'x' projected and
'y' interpolated?

  Thanks,

     Matt


>   Thanks,
>
>      Matt
>
>
>> Thanks,
>> Sander
>>
>
>
>
> --
> 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
>



-- 
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/20160618/1205193b/attachment.html>