--- /dev/null
+//-------------------------------------------------------------------
+// Copyright (C) 2017 by the deal.II authors.
+//
+// This file is subject to LGPL and may not be distributed
+// without copyright and license information. Please refer
+// to the file deal.II/doc/license.html for the text and
+// further information on this license.
+//
+//-------------------------------------------------------------------
+
+
+// Test that transfinite interpolation manifold returns valid results on the
+// same geometry as transfinite_manifold_05 (ball inside square, curved
+// spherical surface). In an initial version, the line search in Newton for
+// the transfinite interpolation would eagerly search too far outside the
+// valid chart domain, leading to failures in the spherical manifold.
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/manifold_lib.h>
+
+
+int main ()
+{
+ initlog();
+
+ const int dim = 3;
+ Triangulation<dim> tria1, tria2, tria;
+ GridGenerator::hyper_shell(tria1, Point<dim>(), 0.4, std::sqrt(dim), 6);
+ GridGenerator::hyper_ball(tria2, Point<dim>(), 0.4);
+ GridGenerator::merge_triangulations(tria1, tria2, tria);
+ tria.set_all_manifold_ids(0);
+ for (typename Triangulation<dim>::cell_iterator cell = tria.begin();
+ cell != tria.end(); ++cell)
+ {
+ for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+ {
+ bool face_at_sphere_boundary = true;
+ for (unsigned int v=0; v<GeometryInfo<dim-1>::vertices_per_cell; ++v)
+ if (std::abs(cell->face(f)->vertex(v).norm()-0.4) > 1e-12)
+ face_at_sphere_boundary = false;
+ if (face_at_sphere_boundary)
+ cell->face(f)->set_all_manifold_ids(1);
+ }
+ }
+ static const SphericalManifold<dim> spherical_manifold;
+ tria.set_manifold(1, spherical_manifold);
+ static TransfiniteInterpolationManifold<dim> transfinite;
+ transfinite.initialize(tria);
+ tria.set_manifold(0, transfinite);
+
+ tria.refine_global(1);
+
+ deallog.precision(10);
+ deallog << "Cell centers" << std::endl;
+ for (auto cell : tria.cell_iterators())
+ deallog << cell->id() << " has center " << cell->center(/*respect_manifold*/true) << std::endl;
+
+ return 0;
+}
--- /dev/null
+
+DEAL::Cell centers
+DEAL::0_0: has center -1.249000903e-15 -3.108624469e-15 -0.7000000000
+DEAL::1_0: has center 0.7000000000 -2.164934898e-15 -5.107025913e-15
+DEAL::2_0: has center -1.335084820e-12 -7.147060721e-15 0.7000000000
+DEAL::3_0: has center -0.7000000000 9.436895709e-16 -1.340427769e-12
+DEAL::4_0: has center -3.719247132e-15 -0.7000000000 -1.335043187e-12
+DEAL::5_0: has center -6.689093723e-15 0.7000000000 -1.346367462e-12
+DEAL::6_0: has center 0.000000000 2.081668171e-17 2.081668171e-17
+DEAL::7_0: has center 5.551115123e-17 5.377642776e-17 -0.2422649731
+DEAL::8_0: has center 0.2422649731 5.724587471e-17 -7.028735233e-13
+DEAL::9_0: has center 5.724587471e-17 -7.028735233e-13 0.2422649731
+DEAL::10_0: has center -0.2422649731 5.724587471e-17 5.898059818e-17
+DEAL::11_0: has center 5.724587471e-17 -0.2422649731 -7.028735233e-13
+DEAL::12_0: has center 5.898059818e-17 0.2422649731 5.724587471e-17
+DEAL::0_1:0 has center -0.4117230837 -0.4117230837 -0.8354565986
+DEAL::0_1:1 has center 0.4117230837 -0.4117230837 -0.8354565986
+DEAL::0_1:2 has center -0.4117230837 0.4117230837 -0.8354565986
+DEAL::0_1:3 has center 0.4117230837 0.4117230837 -0.8354565986
+DEAL::0_1:4 has center -0.2348269975 -0.2348269975 -0.5064930967
+DEAL::0_1:5 has center 0.2348269975 -0.2348269975 -0.5064930967
+DEAL::0_1:6 has center -0.2348269975 0.2348269975 -0.5064930967
+DEAL::0_1:7 has center 0.2348269975 0.2348269975 -0.5064930967
+DEAL::1_1:0 has center 0.8354565986 -0.4117230837 -0.4117230837
+DEAL::1_1:1 has center 0.8354565986 0.4117230837 -0.4117230837
+DEAL::1_1:2 has center 0.5064930967 -0.2348269975 -0.2348269975
+DEAL::1_1:3 has center 0.5064930967 0.2348269975 -0.2348269975
+DEAL::1_1:4 has center 0.8354565986 -0.4117230837 0.4117230837
+DEAL::1_1:5 has center 0.8354565986 0.4117230837 0.4117230837
+DEAL::1_1:6 has center 0.5064930967 -0.2348269975 0.2348269975
+DEAL::1_1:7 has center 0.5064930967 0.2348269975 0.2348269975
+DEAL::2_1:0 has center -0.4117230837 -0.4117230837 0.8354565986
+DEAL::2_1:1 has center 0.4117230837 -0.4117230837 0.8354565986
+DEAL::2_1:2 has center -0.2348269975 -0.2348269975 0.5064930967
+DEAL::2_1:3 has center 0.2348269975 -0.2348269975 0.5064930967
+DEAL::2_1:4 has center -0.4117230837 0.4117230837 0.8354565986
+DEAL::2_1:5 has center 0.4117230837 0.4117230837 0.8354565986
+DEAL::2_1:6 has center -0.2348269975 0.2348269975 0.5064930967
+DEAL::2_1:7 has center 0.2348269975 0.2348269975 0.5064930967
+DEAL::3_1:0 has center -0.8354565986 -0.4117230837 -0.4117230837
+DEAL::3_1:1 has center -0.5064930967 -0.2348269975 -0.2348269975
+DEAL::3_1:2 has center -0.8354565986 0.4117230837 -0.4117230837
+DEAL::3_1:3 has center -0.5064930967 0.2348269975 -0.2348269975
+DEAL::3_1:4 has center -0.8354565986 -0.4117230837 0.4117230837
+DEAL::3_1:5 has center -0.5064930967 -0.2348269975 0.2348269975
+DEAL::3_1:6 has center -0.8354565986 0.4117230837 0.4117230837
+DEAL::3_1:7 has center -0.5064930967 0.2348269975 0.2348269975
+DEAL::4_1:0 has center -0.4117230837 -0.8354565986 -0.4117230837
+DEAL::4_1:1 has center 0.4117230837 -0.8354565986 -0.4117230837
+DEAL::4_1:2 has center -0.2348269975 -0.5064930967 -0.2348269975
+DEAL::4_1:3 has center 0.2348269975 -0.5064930967 -0.2348269975
+DEAL::4_1:4 has center -0.4117230837 -0.8354565986 0.4117230837
+DEAL::4_1:5 has center 0.4117230837 -0.8354565986 0.4117230837
+DEAL::4_1:6 has center -0.2348269975 -0.5064930967 0.2348269975
+DEAL::4_1:7 has center 0.2348269975 -0.5064930967 0.2348269975
+DEAL::5_1:0 has center -0.4117230837 0.8354565986 -0.4117230837
+DEAL::5_1:1 has center -0.2348269975 0.5064930967 -0.2348269975
+DEAL::5_1:2 has center 0.4117230837 0.8354565986 -0.4117230837
+DEAL::5_1:3 has center 0.2348269975 0.5064930967 -0.2348269975
+DEAL::5_1:4 has center -0.4117230837 0.8354565986 0.4117230837
+DEAL::5_1:5 has center -0.2348269975 0.5064930967 0.2348269975
+DEAL::5_1:6 has center 0.4117230837 0.8354565986 0.4117230837
+DEAL::5_1:7 has center 0.2348269975 0.5064930967 0.2348269975
+DEAL::6_1:0 has center -0.04226497308 -0.04226497308 -0.04226497308
+DEAL::6_1:1 has center 0.04226497308 -0.04226497308 -0.04226497308
+DEAL::6_1:2 has center -0.04226497308 0.04226497308 -0.04226497308
+DEAL::6_1:3 has center 0.04226497308 0.04226497308 -0.04226497308
+DEAL::6_1:4 has center -0.04226497308 -0.04226497308 0.04226497308
+DEAL::6_1:5 has center 0.04226497308 -0.04226497308 0.04226497308
+DEAL::6_1:6 has center -0.04226497308 0.04226497308 0.04226497308
+DEAL::6_1:7 has center 0.04226497308 0.04226497308 0.04226497308
+DEAL::7_1:0 has center -0.1205751096 -0.1205751096 -0.2776255832
+DEAL::7_1:1 has center 0.1205751096 -0.1205751096 -0.2776255832
+DEAL::7_1:2 has center -0.1205751096 0.1205751096 -0.2776255832
+DEAL::7_1:3 has center 0.1205751096 0.1205751096 -0.2776255832
+DEAL::7_1:4 has center -0.06842181353 -0.06842181353 -0.1488540582
+DEAL::7_1:5 has center 0.06842181353 -0.06842181353 -0.1488540582
+DEAL::7_1:6 has center -0.06842181353 0.06842181353 -0.1488540582
+DEAL::7_1:7 has center 0.06842181353 0.06842181353 -0.1488540582
+DEAL::8_1:0 has center 0.2776255832 -0.1205751096 -0.1205751096
+DEAL::8_1:1 has center 0.2776255832 0.1205751096 -0.1205751096
+DEAL::8_1:2 has center 0.1488540582 -0.06842181353 -0.06842181353
+DEAL::8_1:3 has center 0.1488540582 0.06842181353 -0.06842181353
+DEAL::8_1:4 has center 0.2776255832 -0.1205751096 0.1205751096
+DEAL::8_1:5 has center 0.2776255832 0.1205751096 0.1205751096
+DEAL::8_1:6 has center 0.1488540582 -0.06842181353 0.06842181353
+DEAL::8_1:7 has center 0.1488540582 0.06842181353 0.06842181353
+DEAL::9_1:0 has center -0.1205751096 -0.1205751096 0.2776255832
+DEAL::9_1:1 has center 0.1205751096 -0.1205751096 0.2776255832
+DEAL::9_1:2 has center -0.06842181353 -0.06842181353 0.1488540582
+DEAL::9_1:3 has center 0.06842181353 -0.06842181353 0.1488540582
+DEAL::9_1:4 has center -0.1205751096 0.1205751096 0.2776255832
+DEAL::9_1:5 has center 0.1205751096 0.1205751096 0.2776255832
+DEAL::9_1:6 has center -0.06842181353 0.06842181353 0.1488540582
+DEAL::9_1:7 has center 0.06842181353 0.06842181353 0.1488540582
+DEAL::10_1:0 has center -0.2776255832 -0.1205751096 -0.1205751096
+DEAL::10_1:1 has center -0.1488540582 -0.06842181353 -0.06842181353
+DEAL::10_1:2 has center -0.2776255832 0.1205751096 -0.1205751096
+DEAL::10_1:3 has center -0.1488540582 0.06842181353 -0.06842181353
+DEAL::10_1:4 has center -0.2776255832 -0.1205751096 0.1205751096
+DEAL::10_1:5 has center -0.1488540582 -0.06842181353 0.06842181353
+DEAL::10_1:6 has center -0.2776255832 0.1205751096 0.1205751096
+DEAL::10_1:7 has center -0.1488540582 0.06842181353 0.06842181353
+DEAL::11_1:0 has center -0.1205751096 -0.2776255832 -0.1205751096
+DEAL::11_1:1 has center 0.1205751096 -0.2776255832 -0.1205751096
+DEAL::11_1:2 has center -0.06842181353 -0.1488540582 -0.06842181353
+DEAL::11_1:3 has center 0.06842181353 -0.1488540582 -0.06842181353
+DEAL::11_1:4 has center -0.1205751096 -0.2776255832 0.1205751096
+DEAL::11_1:5 has center 0.1205751096 -0.2776255832 0.1205751096
+DEAL::11_1:6 has center -0.06842181353 -0.1488540582 0.06842181353
+DEAL::11_1:7 has center 0.06842181353 -0.1488540582 0.06842181353
+DEAL::12_1:0 has center -0.1205751096 0.2776255832 -0.1205751096
+DEAL::12_1:1 has center -0.06842181353 0.1488540582 -0.06842181353
+DEAL::12_1:2 has center 0.1205751096 0.2776255832 -0.1205751096
+DEAL::12_1:3 has center 0.06842181353 0.1488540582 -0.06842181353
+DEAL::12_1:4 has center -0.1205751096 0.2776255832 0.1205751096
+DEAL::12_1:5 has center -0.06842181353 0.1488540582 0.06842181353
+DEAL::12_1:6 has center 0.1205751096 0.2776255832 0.1205751096
+DEAL::12_1:7 has center 0.06842181353 0.1488540582 0.06842181353