From: Daniel Arndt Date: Sun, 5 Nov 2017 19:23:36 +0000 (+0100) Subject: Use PolarManifold in SphericalManifold for 2D X-Git-Tag: v9.0.0-rc1~818^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=55471e19722e5c73288a146484bb2193239a5e33;p=dealii.git Use PolarManifold in SphericalManifold for 2D --- diff --git a/include/deal.II/grid/manifold_lib.h b/include/deal.II/grid/manifold_lib.h index bb808be3de..6a0294afe5 100644 --- a/include/deal.II/grid/manifold_lib.h +++ b/include/deal.II/grid/manifold_lib.h @@ -248,6 +248,13 @@ public: * The center of the spherical coordinate system. */ const Point center; + +private: + + /** + * A manifold description to be used for get_new_point in 2D. + **/ + const PolarManifold polar_manifold; }; diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc index 5256a5ead9..94d25fc02e 100644 --- a/source/grid/manifold_lib.cc +++ b/source/grid/manifold_lib.cc @@ -245,7 +245,8 @@ PolarManifold::push_forward_gradient(const Point &spheri template SphericalManifold::SphericalManifold(const Point center): - center(center) + center(center), + polar_manifold(center) {} @@ -403,8 +404,10 @@ get_new_point (const ArrayView> &vertices, rho /= total_weights; } - if (spacedim<2) - return center + rho*candidate; + // If not in 3D, just use the implementation from PolarManifold + // after we verified that the candidate is not the center. + if (spacedim<3) + return polar_manifold.get_new_point(vertices, weights); // Step 2: // Do Newton-style iterations to improve the estimate.