This is actually not a bug. A spherical manifold is undefined on the center of the sphere/ball. You cannot attach set_all_manifold_ids of the spherical manifold and attach a spherical manifold to it, since the central cell (the one where the cell->center() coincides with the ball center) would not know what to do.
This is documented in the SphericalManfold documentation. See the phrase "In order for this manifold to work correctly, it cannot be attached to cells containing the center of the coordinate system. This point is a singular point of the coordinate transformation, and there taking averages does not make any sense.” Best, L. > On 25 Apr 2017, at 13:42, Alexander Gary Zimmerman > <[email protected]> wrote: > > When I refine a Triangulation that has been initialized with > GridGenerator::hyper_ball, and has a SphericalManifold attached (and having > set_all_manifold_ids on the Triangulation), some nan (not a number) values > are produced. > > On the other hand, I get the expected behavior when instead using a > GridGenerator::hyper_shell. > > I've documented the issue on my repo here: > https://github.com/alexanderzimmerman/nsb-pcm/issues/21 > > I'm using deal.II-v.8.5.0. I haven't tried reproducing this on any other > versions. > > These are in the repo, but here's the test that passes: > > #include <iostream> > #include <fstream> > #include <string> > > #include <deal.II/grid/tria.h> > #include <deal.II/grid/tria_accessor.h> > #include <deal.II/grid/tria_iterator.h> > #include <deal.II/grid/grid_generator.h> > #include <deal.II/grid/grid_refinement.h> > #include <deal.II/grid/grid_out.h> > #include <deal.II/grid/manifold_lib.h> > #include <deal.II/base/numbers.h> > #include <deal.II/base/signaling_nan.h> > > using namespace dealii; > > const unsigned int dim = 2; > > > void validate_vertices(Triangulation<dim> &tria) > { > Triangulation<dim>::cell_iterator cell; > > unsigned int nv = GeometryInfo<dim>::vertices_per_cell; > > for (cell = tria.begin(); cell != tria.end(); ++cell) > { > for (unsigned int v = 0; v < nv; ++v) > { > for (unsigned int i = 0; i < dim; ++i) > { > assert(!numbers::is_nan(cell->vertex(v)[i])); > } > > } > > } > > } > > > int main(int /*argc*/, char** /*argv*/) > { > Triangulation<dim> triangulation; > > GridGenerator::hyper_shell(triangulation, Point<dim>(), 0.1, 0.5); > > const unsigned int manifold_id = 0; > > triangulation.set_all_manifold_ids(manifold_id); > > SphericalManifold<dim> spherical_manifold; > > triangulation.set_manifold(manifold_id, spherical_manifold); > > GridOut grid_writer; > > > for (unsigned int r = 0; r < 2; ++r) > { > std::string file_name = "grid_"+std::to_string(r)+".vtk"; > > std::ofstream file(file_name); > > grid_writer.write_vtk(triangulation, file); > > validate_vertices(triangulation); > > triangulation.refine_global(1); > } > > triangulation.set_manifold(0); > > return 0; > } > > > and here's the test that fails: > > #include <iostream> > #include <fstream> > #include <string> > > #include <deal.II/grid/tria.h> > #include <deal.II/grid/tria_accessor.h> > #include <deal.II/grid/tria_iterator.h> > #include <deal.II/grid/grid_generator.h> > #include <deal.II/grid/grid_refinement.h> > #include <deal.II/grid/grid_out.h> > #include <deal.II/grid/manifold_lib.h> > #include <deal.II/base/numbers.h> > #include <deal.II/base/signaling_nan.h> > > > using namespace dealii; > > const unsigned int dim = 2; > > > void validate_vertices(Triangulation<dim> &tria) > { > Triangulation<dim>::cell_iterator cell; > > unsigned int nv = GeometryInfo<dim>::vertices_per_cell; > > for (cell = tria.begin(); cell != tria.end(); ++cell) > { > for (unsigned int v = 0; v < nv; ++v) > { > for (unsigned int i = 0; i < dim; ++i) > { > assert(!numbers::is_nan(cell->vertex(v)[i])); > } > > } > > } > > } > > > int main(int /*argc*/, char** /*argv*/) > { > Triangulation<dim> triangulation; > > GridGenerator::hyper_ball(triangulation); > > const unsigned int manifold_id = 0; > > triangulation.set_all_manifold_ids(manifold_id); > > SphericalManifold<dim> spherical_manifold; > > triangulation.set_manifold(manifold_id, spherical_manifold); > > GridOut grid_writer; > > > for (unsigned int r = 0; r < 2; ++r) > { > std::string file_name = "grid_"+std::to_string(r)+".vtk"; > > std::ofstream file(file_name); > > grid_writer.write_vtk(triangulation, file); > > validate_vertices(triangulation); > > triangulation.refine_global(1); > } > > triangulation.set_manifold(0); > > return 0; > } > > > Any ideas? > > Thanks, > > Alex > > P.S. for future reference, if I have an issue formulated like this, should I > raise it on the deal.II repo or here in the mailing list? > > -- > 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 [email protected]. > For more options, visit https://groups.google.com/d/optout. -- 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 [email protected]. For more options, visit https://groups.google.com/d/optout.
