problem with a projection using SphericalManifold

[Nicola Giuliani]

Dear all,

I have the following snippet of code.

  const unsigned int dim = 2;

  double factor = 1./std::sqrt(dim);

  Point<dim> p1test{factor, -factor};

  Point<dim> p2test{factor, factor};

  Point<dim> p3test{-factor, -factor};

  std::vector<Point<dim> > surrounding_points;

  surrounding_points.push_back(p1test);

  surrounding_points.push_back(p2test);

  surrounding_points.push_back(p3test);

  std::vector<double> surrounding_weights(3);

  surrounding_weights[0] = 0.5;

  surrounding_weights[1] = 0.25;

  surrounding_weights[2] = 0.25;  

  auto new_vertex = spherical_manifold.get_new_point(

            ArrayView<Point<spacedim>>(&surrounding_points[0],

                                       surrounding_points.size()),

            ArrayView<double>(&surrounding_weights[0],

                              surrounding_weights.size()));

I would expect the new point to coincide with the original one, but instead it collapses on the second. Looking into the code of spherical manifold the associated angles are: 7/4 pi, 1/4 pi and 5/4 pi. From this representation I would expect the projected point to be at 5/4 pi. 

I really don't understand this behavior of SphericalManifold, do you have an explanation?

Bests,

Nicola

[luca.heltai]

Nicola,

I think you are hitting a periodicity issue here.

The points with angle 7/4 pi and 5/4 pi are at a distance (from pi/4) of 6/4 pi = -pi/2 (this is the shorter distance, and this is what is used in the code), and pi = -pi (this could be either -pi or pi, and it is in fact not well defined to use points at exactly T/2 distance where T is the periodicity, since either T/2 or -T/2 are the same point) from 1/4 pi.

With your weights, you may get -pi/4 -pi/4 + pi/4 = - pi/4 or -pi/4 + pi/4 + pi/4 = pi/4 (which is what you observe). The only solution to this is to make sure that the initial grid does not contain a single cell, with a periodic boundary associated to the it. This produces a manifold which is not well defined, and does not behave as one would expect.

L.

> --
> The deal.II project is located at http://www.dealii.org/
> For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en
> ---
> You received this message because you are subscribed to the Google Groups "deal.II User Group" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to dealii+un...@googlegroups.com.
> To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/ca22bb00-481d-49ae-aa38-1572415d8340%40googlegroups.com.