From 3645938c6efe68585b3b32cfe4be6b73484824de Mon Sep 17 00:00:00 2001 From: Rene Gassmoeller Date: Fri, 19 Aug 2016 17:16:36 -0600 Subject: [PATCH] Reduce use of Quadrature objects in Manifolds. --- include/deal.II/grid/manifold.h | 174 +++++++++++++++---- include/deal.II/grid/manifold_lib.h | 10 +- source/fe/mapping_manifold.cc | 14 +- source/fe/mapping_q_generic.cc | 9 +- source/grid/manifold.cc | 103 ++++++++--- source/grid/manifold_lib.cc | 15 +- source/grid/tria.cc | 3 +- source/grid/tria_accessor.cc | 9 +- tests/manifold/chart_manifold_03.cc | 3 +- tests/manifold/chart_manifold_03_embedded.cc | 3 +- tests/manifold/chart_manifold_09.cc | 2 +- tests/manifold/composition_manifold_01.cc | 2 +- tests/manifold/composition_manifold_02.cc | 2 +- tests/manifold/composition_manifold_03.cc | 2 +- tests/manifold/composition_manifold_04.cc | 4 +- tests/manifold/flat_manifold_03.cc | 3 +- tests/manifold/flat_manifold_07.cc | 2 +- tests/manifold/function_manifold_03.cc | 2 +- tests/manifold/polar_manifold_08.cc | 2 +- tests/manifold/spherical_manifold_01.cc | 10 +- tests/manifold/tensor_product_manifold_01.cc | 2 +- 21 files changed, 278 insertions(+), 98 deletions(-) diff --git a/include/deal.II/grid/manifold.h b/include/deal.II/grid/manifold.h index 5f3b7a8ee8..52646e054d 100644 --- a/include/deal.II/grid/manifold.h +++ b/include/deal.II/grid/manifold.h @@ -79,7 +79,54 @@ namespace Manifolds template Quadrature get_default_quadrature(const MeshIteratorType &iterator, - const bool with_laplace = false); + const bool with_laplace = false) DEAL_II_DEPRECATED; + + /** + * Given a general mesh iterator, construct vectors of quadrature points and + * weights that contain the following points: + * - If the iterator points to a line, then the quadrature points + * are the two vertices of the line. This results in a point vector + * with two points. + * - If the iterator points to a quad, then the quadrature points + * are the vertices and line mid-points. This results in a point vector + * with eight (4+4) points. + * - If the iterator points to a hex, then the quadrature points + * are the vertices, the line mid-points, and the face mid-points. + * This results in a points vector with 26 (8+12+6) points. + * + * The quadrature weights for these points are either chosen identically + * and equal to one over the number of quadrature points (if @p with_laplace + * is @p false), or in a way that gives points closer to the cell center + * (measured on the reference cell) a higher weight. These weights correspond + * to solving a Laplace equation and evaluating the solution at the quadrature + * points (if @p with_laplace is @p true). + * + * The function is primarily used to construct the input argument + * for the Manifold::get_new_point() function, which computes a new + * point on a manifold based on a weighted average of "surrounding" + * points represented by the quadrature points and weights stored in the + * returned pair of vectors. This function creates such an object based on + * the points that "surround" a cell, face, or edge, and weights + * are chosen in a way appropriate for computing the new "mid-point" + * of the object pointed to. An example of where this is necessary + * is for mesh refinement, where (using the 2d situation as an example) + * we need to first create new edge mid-points, and then a new cell-point. + * + * @param[in] iterator A mesh iterator that points to either a line, quad, + * or hex. + * @param[in] with_laplace Whether or not to compute the quadrature weights + * by solving a Laplace equation, as discussed above. + * @tparam MeshIteratorType An iterator type that corresponds to either + * Triangulation::cell_iterator (or variants such as + * Triangulation::active_cell_iterator or DoFHandler::cell_iterator) or + * that is the result of statements such as + * cell-@>face(f) or cell-@>line(l). + */ + template + std::pair >, + std::vector > + get_default_points_and_weights(const MeshIteratorType &iterator, + const bool with_laplace = false); } @@ -318,7 +365,28 @@ public: */ virtual Point - get_new_point (const Quadrature &quad) const; + get_new_point (const Quadrature &quad) const DEAL_II_DEPRECATED; + + /** + * Return the point which shall become the new vertex surrounded by the + * given points @p surrounding_points. @p weights contains appropriate + * weights for the surrounding points according to which the manifold + * determines the new point's position. + * + * In its default implementation it uses a pair-wise reduction of + * the points by calling the function get_intermediate_point() on the first + * two points, then on the resulting point and the next, until all points in + * the vector have been taken into account. User classes can get away by + * simply implementing the get_intermediate_point() function. Notice that + * by default the get_intermediate_point() function calls the + * project_to_manifold() function with the convex combination of its + * arguments. For simple situations you may get away by implementing + * only the project_to_manifold() function. + */ + virtual + Point + get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const; /** * Given a point which lies close to the given manifold, it modifies it and @@ -589,7 +657,33 @@ public: */ virtual Point - get_new_point(const Quadrature &quad) const; + get_new_point(const Quadrature &quad) const DEAL_II_DEPRECATED; + + /** + * Let the new point be the average sum of surrounding vertices. + * + * This particular implementation constructs the weighted average of the + * surrounding points, and then calls internally the function + * project_to_manifold(). The reason why we do it this way, is to allow lazy + * programmers to implement only the project_to_manifold() function for their + * own Manifold classes which are small (or trivial) perturbations of a flat + * manifold. This is the case whenever the coarse mesh is a decent + * approximation of the manifold geometry. In this case, the middle point of + * a cell is close to true middle point of the manifold, and a projection + * may suffice. + * + * For most simple geometries, it is possible to get reasonable results by + * deriving your own Manifold class from FlatManifold, and write a new + * interface only for the project_to_manifold function. You will have good + * approximations also with large deformations, as long as in the coarsest + * mesh size you are trying to refine, the middle point is not too far from + * the manifold mid point, i.e., as long as the coarse mesh size is small + * enough. + */ + virtual + Point + get_new_point(const std::vector > &surrounding_points, + const std::vector &weights) const; /** @@ -784,7 +878,16 @@ public: */ virtual Point - get_new_point(const Quadrature &quad) const; + get_new_point(const Quadrature &quad) const DEAL_II_DEPRECATED; + + /** + * Refer to the general documentation of this class and the documentation of + * the base class for more information. + */ + virtual + Point + get_new_point(const std::vector > &surrounding_points, + const std::vector &weights) const; /** * Pull back the given point in spacedim to the Euclidean chartdim @@ -953,12 +1056,25 @@ namespace Manifolds Quadrature get_default_quadrature(const MeshIteratorType &iterator, const bool with_laplace) + { + const std::pair >, + std::vector > points_and_weights = get_default_points_and_weights(iterator, + with_laplace); + return Quadrature(points_and_weights.first, + points_and_weights.second); + } + + template + std::pair >, + std::vector > + get_default_points_and_weights(const MeshIteratorType &iterator, + const bool with_laplace) { const int spacedim = MeshIteratorType::AccessorType::space_dimension; const int dim = MeshIteratorType::AccessorType::structure_dimension; - std::vector > sp; - std::vector wp; + std::pair >, + std::vector > points_weights; // note that the exact weights are chosen such as to minimize the @@ -968,32 +1084,32 @@ namespace Manifolds switch (dim) { case 1: - sp.resize(2); - wp.resize(2); - sp[0] = iterator->vertex(0); - wp[0] = .5; - sp[1] = iterator->vertex(1); - wp[1] = .5; + points_weights.first.resize(2); + points_weights.second.resize(2); + points_weights.first[0] = iterator->vertex(0); + points_weights.second[0] = .5; + points_weights.first[1] = iterator->vertex(1); + points_weights.second[1] = .5; break; case 2: - sp.resize(8); - wp.resize(8); + points_weights.first.resize(8); + points_weights.second.resize(8); for (unsigned int i=0; i<4; ++i) { - sp[i] = iterator->vertex(i); - sp[4+i] = ( iterator->line(i)->has_children() ? - iterator->line(i)->child(0)->vertex(1) : - iterator->line(i)->get_manifold().get_new_point_on_line(iterator->line(i)) ); + points_weights.first[i] = iterator->vertex(i); + points_weights.first[4+i] = ( iterator->line(i)->has_children() ? + iterator->line(i)->child(0)->vertex(1) : + iterator->line(i)->get_manifold().get_new_point_on_line(iterator->line(i)) ); } if (with_laplace) { - std::fill(wp.begin(), wp.begin()+4, 1.0/16.0); - std::fill(wp.begin()+4, wp.end(), 3.0/16.0); + std::fill(points_weights.second.begin(), points_weights.second.begin()+4, 1.0/16.0); + std::fill(points_weights.second.begin()+4, points_weights.second.end(), 3.0/16.0); } else - std::fill(wp.begin(), wp.end(), 1.0/8.0); + std::fill(points_weights.second.begin(), points_weights.second.end(), 1.0/8.0); break; case 3: { @@ -1003,9 +1119,9 @@ namespace Manifolds GeometryInfo::vertices_per_cell+ GeometryInfo::lines_per_cell+ GeometryInfo::faces_per_cell; - sp.resize(np); - wp.resize(np); - std::vector > *sp3 = reinterpret_cast > *>(&sp); + points_weights.first.resize(np); + points_weights.second.resize(np); + std::vector > *sp3 = reinterpret_cast > *>(&points_weights.first); unsigned int j=0; @@ -1016,33 +1132,33 @@ namespace Manifolds for (unsigned int i=0; i::vertices_per_cell; ++i, ++j) { (*sp3)[j] = hex->vertex(i); - wp[j] = 1.0/128.0; + points_weights.second[j] = 1.0/128.0; } for (unsigned int i=0; i::lines_per_cell; ++i, ++j) { (*sp3)[j] = (hex->line(i)->has_children() ? hex->line(i)->child(0)->vertex(1) : hex->line(i)->get_manifold().get_new_point_on_line(hex->line(i))); - wp[j] = 7.0/192.0; + points_weights.second[j] = 7.0/192.0; } for (unsigned int i=0; i::faces_per_cell; ++i, ++j) { (*sp3)[j] = (hex->quad(i)->has_children() ? hex->quad(i)->isotropic_child(0)->vertex(3) : hex->quad(i)->get_manifold().get_new_point_on_quad(hex->quad(i))); - wp[j] = 1.0/12.0; + points_weights.second[j] = 1.0/12.0; } // Overwrite the weights with 1/np if we don't want to use // laplace vectors. if (with_laplace == false) - std::fill(wp.begin(), wp.end(), 1.0/np); + std::fill(points_weights.second.begin(), points_weights.second.end(), 1.0/np); } break; default: Assert(false, ExcInternalError()); break; } - return Quadrature(sp,wp); + return points_weights; } } diff --git a/include/deal.II/grid/manifold_lib.h b/include/deal.II/grid/manifold_lib.h index 40a383aea4..f9583582bd 100644 --- a/include/deal.II/grid/manifold_lib.h +++ b/include/deal.II/grid/manifold_lib.h @@ -289,7 +289,15 @@ public: * the base class for a detailed description of what this function does. */ virtual Point - get_new_point(const Quadrature &quad) const; + get_new_point(const Quadrature &quad) const DEAL_II_DEPRECATED; + + /** + * Compute new points on the CylindricalManifold. See the documentation of + * the base class for a detailed description of what this function does. + */ + virtual Point + get_new_point(const std::vector > &surrounding_points, + const std::vector &weights) const; protected: /** diff --git a/source/fe/mapping_manifold.cc b/source/fe/mapping_manifold.cc index 481155db61..e0b7e4e1ec 100644 --- a/source/fe/mapping_manifold.cc +++ b/source/fe/mapping_manifold.cc @@ -219,7 +219,7 @@ transform_unit_to_real_cell (const typename Triangulation::cell_it vertices.push_back(cell->vertex(v)); weights.push_back(GeometryInfo::d_linear_shape_function(p,v)); } - return cell->get_manifold().get_new_point(Quadrature(vertices, weights)); + return cell->get_manifold().get_new_point(vertices, weights); } @@ -403,8 +403,8 @@ namespace internal for (unsigned int point=0; point - get_new_point(Quadrature(data.vertices, - data.cell_manifold_quadrature_weights[point+data_set])); + get_new_point(data.vertices, + data.cell_manifold_quadrature_weights[point+data_set]); } } } @@ -440,8 +440,8 @@ namespace internal // And get its image on the manifold: const Point P = data.manifold-> - get_new_point(Quadrature(data.vertices, - data.cell_manifold_quadrature_weights[point+data_set])); + get_new_point(data.vertices, + data.cell_manifold_quadrature_weights[point+data_set]); // To compute the Jacobian, we choose dim points aligned // with the dim reference axes, which are still in the @@ -469,8 +469,8 @@ namespace internal data.vertex_weights[j] = GeometryInfo::d_linear_shape_function(np, j); const Point NP= - data.manifold->get_new_point(Quadrature(data.vertices, - data.vertex_weights)); + data.manifold->get_new_point(data.vertices, + data.vertex_weights); const Tensor<1,spacedim> T = data.manifold->get_tangent_vector(P, NP); diff --git a/source/fe/mapping_q_generic.cc b/source/fe/mapping_q_generic.cc index 317dc0b5c1..5de96793cc 100644 --- a/source/fe/mapping_q_generic.cc +++ b/source/fe/mapping_q_generic.cc @@ -3484,8 +3484,7 @@ namespace const double x = line_support_points.point(i+1)[0]; w[1] = x; w[0] = (1-x); - Quadrature quadrature(surrounding_points, w); - points[i] = manifold.get_new_point(quadrature); + points[i] = manifold.get_new_point(surrounding_points, w); } break; } @@ -3511,8 +3510,7 @@ namespace for (unsigned int l=0; l<4; ++l) w[l] = GeometryInfo<2>::d_linear_shape_function(p, l); - Quadrature quadrature(surrounding_points, w); - points[c]=manifold.get_new_point(quadrature); + points[c]=manifold.get_new_point(surrounding_points, w); } } break; @@ -3542,8 +3540,7 @@ namespace for (unsigned int l=0; l<8; ++l) w[l] = GeometryInfo<3>::d_linear_shape_function(p, l); - Quadrature quadrature(surrounding_points, w); - points[c]=manifold.get_new_point(quadrature); + points[c]=manifold.get_new_point(surrounding_points, w); } } } diff --git a/source/grid/manifold.cc b/source/grid/manifold.cc index b318665216..d4f9a9fe9f 100644 --- a/source/grid/manifold.cc +++ b/source/grid/manifold.cc @@ -59,29 +59,67 @@ template Point Manifold:: get_new_point (const Quadrature &quad) const +{ + return get_new_point(quad.get_points(),quad.get_weights()); +} + +template +Point +Manifold:: +get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const { const double tol = 1e-10; + const unsigned int n_points = surrounding_points.size(); - Assert(quad.size() > 0, - ExcMessage("Quadrature should have at least one point.")); + Assert(n_points > 0, + ExcMessage("There should be at least one point.")); - Assert(std::abs(std::accumulate(quad.get_weights().begin(), quad.get_weights().end(), 0.0)-1.0) < tol, + Assert(n_points == weights.size(), + ExcMessage("There should be as many surrounding points as weights given.")); + + Assert(std::abs(std::accumulate(weights.begin(), weights.end(), 0.0)-1.0) < tol, ExcMessage("The weights for the individual points should sum to 1!")); - QSorted sorted_quad(quad); - Point p = sorted_quad.point(0); - double w = sorted_quad.weight(0); + // The number of points is small (at most 26 in 3D), therefore try to minize + // overhead and use a very simple search of O(N^2) to find the order in which + // to loop over points. + std::vector permutation(n_points); + std::vector found(n_points,false); + + unsigned int min_index = numbers::invalid_unsigned_int; + for (unsigned int i=0; i < n_points; ++i) + { + double min_weight = std::numeric_limits::max(); + + for (unsigned int j=0; j < n_points; ++j) + { + if ((!found[j]) && (weights[j] < min_weight)) + { + min_weight = weights[j]; + min_index = j; + } + } + permutation[i] = min_index; + found[min_index] = true; + } + + // Now loop over points in the order of their weights. This is done to + // produce unique points even if get_intermediate_points is not + // associative (as for the SphericalManifold). + Point p = surrounding_points[permutation[0]]; + double w = weights[permutation[0]]; - for (unsigned int i=1; i Manifold:: get_new_point_on_line (const typename Triangulation::line_iterator &line) const { - return get_new_point (get_default_quadrature(line)); + const std::pair >, std::vector > points_weights(get_default_points_and_weights(line)); + return get_new_point (points_weights.first,points_weights.second); } @@ -284,7 +323,8 @@ Point Manifold:: get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const { - return get_new_point (get_default_quadrature(quad)); + const std::pair >, std::vector > points_weights(get_default_points_and_weights(quad)); + return get_new_point (points_weights.first,points_weights.second); } @@ -407,7 +447,8 @@ Point<3> Manifold<3,3>:: get_new_point_on_hex (const Triangulation<3, 3>::hex_iterator &hex) const { - return get_new_point(get_default_quadrature(hex, true)); + const std::pair >, std::vector > points_weights(get_default_points_and_weights(hex,true)); + return get_new_point (points_weights.first,points_weights.second); } @@ -427,7 +468,7 @@ Manifold::get_tangent_vector(const Point &x1, w.push_back(epsilon); w.push_back(1.0-epsilon); - const Tensor<1,spacedim> neighbor_point = get_new_point (Quadrature(q, w)); + const Tensor<1,spacedim> neighbor_point = get_new_point (q, w); return (neighbor_point-x1)/epsilon; } @@ -447,11 +488,19 @@ Point FlatManifold:: get_new_point (const Quadrature &quad) const { - Assert(std::abs(std::accumulate(quad.get_weights().begin(), quad.get_weights().end(), 0.0)-1.0) < 1e-10, - ExcMessage("The weights for the individual points should sum to 1!")); + return get_new_point(quad.get_points(),quad.get_weights()); +} - const std::vector > &surrounding_points = quad.get_points(); - const std::vector &weights = quad.get_weights(); + + +template +Point +FlatManifold:: +get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const +{ + Assert(std::abs(std::accumulate(weights.begin(), weights.end(), 0.0)-1.0) < 1e-10, + ExcMessage("The weights for the individual points should sum to 1!")); Tensor<1,spacedim> minP = periodicity; @@ -554,15 +603,23 @@ Point ChartManifold:: get_new_point (const Quadrature &quad) const { - const std::vector > &surrounding_points = quad.get_points(); - const std::vector &weights = quad.get_weights(); + return get_new_point(quad.get_points(),quad.get_weights()); +} + + + +template +Point +ChartManifold:: +get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const +{ std::vector > chart_points(surrounding_points.size()); for (unsigned int i=0; i chart_quad(chart_points, weights); - const Point p_chart = sub_manifold.get_new_point(chart_quad); + const Point p_chart = sub_manifold.get_new_point(chart_points,weights); return push_forward(p_chart); } diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc index 0b0fd94ae0..e278616f97 100644 --- a/source/grid/manifold_lib.cc +++ b/source/grid/manifold_lib.cc @@ -308,11 +308,20 @@ Point CylindricalManifold:: get_new_point (const Quadrature &quad) const { - const std::vector > &surrounding_points = quad.get_points(); - const std::vector &weights = quad.get_weights(); + return get_new_point(quad.get_points(),quad.get_weights()); +} + + +template +Point +CylindricalManifold:: +get_new_point (const std::vector > &surrounding_points, + const std::vector &weights) const +{ // compute a proposed new point - Point middle = flat_manifold.get_new_point(quad); + const Point middle = flat_manifold.get_new_point(surrounding_points, + weights); double radius = 0; Tensor<1,spacedim> on_plane; diff --git a/source/grid/tria.cc b/source/grid/tria.cc index a70616c3c6..0dc044a918 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -4098,9 +4098,8 @@ namespace internal ps[1] = cell->face(GeometryInfo ::opposite_face[boundary_face]) ->child(0)->vertex(1); - Quadrature qs(ps,ws); triangulation.vertices[next_unused_vertex] - = cell->get_manifold().get_new_point(qs); + = cell->get_manifold().get_new_point(ps,ws); } } } diff --git a/source/grid/tria_accessor.cc b/source/grid/tria_accessor.cc index 594fede818..4f3b53a30c 100644 --- a/source/grid/tria_accessor.cc +++ b/source/grid/tria_accessor.cc @@ -1050,8 +1050,10 @@ namespace else { TriaRawIterator > it(obj); - Quadrature quadrature = Manifolds::get_default_quadrature(it, use_laplace); - return obj.get_manifold().get_new_point(quadrature); + const std::pair >, + std::vector > points_and_weights = Manifolds::get_default_points_and_weights(it, use_laplace); + return obj.get_manifold().get_new_point(points_and_weights.first, + points_and_weights.second); } } } @@ -1216,8 +1218,7 @@ TriaAccessor::intermediate_point (const Point quadrature(p, w); - return this->get_manifold().get_new_point(quadrature); + return this->get_manifold().get_new_point(p, w); } diff --git a/tests/manifold/chart_manifold_03.cc b/tests/manifold/chart_manifold_03.cc index 3e992b01c6..31c35e91b6 100644 --- a/tests/manifold/chart_manifold_03.cc +++ b/tests/manifold/chart_manifold_03.cc @@ -119,8 +119,7 @@ void test(unsigned int ref=1) for (unsigned int i=0; i(ps[i],ws); - middle = manifold.get_new_point(quad); + middle = manifold.get_new_point(ps[i],ws); deallog << "P0: " << ps[i][0] << " , P1: " << ps[i][1] << " , Middle: " << middle << std::endl; } diff --git a/tests/manifold/chart_manifold_03_embedded.cc b/tests/manifold/chart_manifold_03_embedded.cc index bc59316482..a1c4bd2b9d 100644 --- a/tests/manifold/chart_manifold_03_embedded.cc +++ b/tests/manifold/chart_manifold_03_embedded.cc @@ -127,8 +127,7 @@ void test(unsigned int ref=1) for (unsigned int i=0; i(ps[i],ws); - middle = manifold.get_new_point(quad); + middle = manifold.get_new_point(ps[i],ws); deallog << "P0: " << ps[i][0] << " , P1: " << ps[i][1] << " , Middle: " << middle << std::endl; } diff --git a/tests/manifold/chart_manifold_09.cc b/tests/manifold/chart_manifold_09.cc index 66fdebd283..c994caf9cf 100644 --- a/tests/manifold/chart_manifold_09.cc +++ b/tests/manifold/chart_manifold_09.cc @@ -60,7 +60,7 @@ main() weights[0] = (double)i/((double)n_intermediates-1); weights[1] = 1.0-weights[0]; - deallog << manifold.get_new_point(Quadrature<2>(points, weights)) << std::endl; + deallog << manifold.get_new_point(points, weights) << std::endl; } } diff --git a/tests/manifold/composition_manifold_01.cc b/tests/manifold/composition_manifold_01.cc index b6bfc236e6..4351b5cb75 100644 --- a/tests/manifold/composition_manifold_01.cc +++ b/tests/manifold/composition_manifold_01.cc @@ -59,7 +59,7 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(Quadrature(sp, w)); + Point ip = manifold.get_new_point(sp, w); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); diff --git a/tests/manifold/composition_manifold_02.cc b/tests/manifold/composition_manifold_02.cc index 6b91128bab..98c40d82fa 100644 --- a/tests/manifold/composition_manifold_02.cc +++ b/tests/manifold/composition_manifold_02.cc @@ -70,7 +70,7 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(Quadrature(sp, w)); + Point ip = manifold.get_new_point(sp, w); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); diff --git a/tests/manifold/composition_manifold_03.cc b/tests/manifold/composition_manifold_03.cc index bfbd11361d..2766d5e863 100644 --- a/tests/manifold/composition_manifold_03.cc +++ b/tests/manifold/composition_manifold_03.cc @@ -74,7 +74,7 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(Quadrature(sp, w)); + Point ip = manifold.get_new_point(sp, w); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); diff --git a/tests/manifold/composition_manifold_04.cc b/tests/manifold/composition_manifold_04.cc index a442601691..6ee8081701 100644 --- a/tests/manifold/composition_manifold_04.cc +++ b/tests/manifold/composition_manifold_04.cc @@ -75,7 +75,7 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(Quadrature(sp, w)); + Point ip = manifold.get_new_point(sp, w); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); @@ -95,7 +95,7 @@ int main () w[1] = 1.0 - w[0]; Point ip = manifold. - pull_back(manifold.get_new_point(Quadrature(sp, w))); + pull_back(manifold.get_new_point(sp, w)); ip[0] = w[1]; diff --git a/tests/manifold/flat_manifold_03.cc b/tests/manifold/flat_manifold_03.cc index 946142ac56..52a4f0612c 100644 --- a/tests/manifold/flat_manifold_03.cc +++ b/tests/manifold/flat_manifold_03.cc @@ -78,8 +78,7 @@ void test(unsigned int ref=1) for (unsigned int i=0; i(ps[i],ws); - middle = manifold.get_new_point(quad); + middle = manifold.get_new_point(ps[i],ws); deallog << "P0: " << ps[i][0] << " , P1: " << ps[i][1] << " , Middle: " << middle << std::endl; } diff --git a/tests/manifold/flat_manifold_07.cc b/tests/manifold/flat_manifold_07.cc index 2b7e0ad662..623eae5acf 100644 --- a/tests/manifold/flat_manifold_07.cc +++ b/tests/manifold/flat_manifold_07.cc @@ -60,7 +60,7 @@ main() weights[0] = (double)i/((double)n_intermediates-1); weights[1] = 1.0-weights[0]; - deallog << manifold.get_new_point(Quadrature<2>(points, weights)) << std::endl; + deallog << manifold.get_new_point(points, weights) << std::endl; } } diff --git a/tests/manifold/function_manifold_03.cc b/tests/manifold/function_manifold_03.cc index a698f7ca1b..136c8d767f 100644 --- a/tests/manifold/function_manifold_03.cc +++ b/tests/manifold/function_manifold_03.cc @@ -69,7 +69,7 @@ void test(unsigned int ref=1) w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(Quadrature(p, w)); + Point ip = manifold.get_new_point(p, w); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, p[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, p[1]); diff --git a/tests/manifold/polar_manifold_08.cc b/tests/manifold/polar_manifold_08.cc index 4638a037df..84c64c718e 100644 --- a/tests/manifold/polar_manifold_08.cc +++ b/tests/manifold/polar_manifold_08.cc @@ -66,7 +66,7 @@ main() weights[0] = (double)i/((double)n_intermediates-1); weights[1] = 1.0-weights[0]; - deallog << manifold.get_new_point(Quadrature<2>(points, weights)) << std::endl; + deallog << manifold.get_new_point(points, weights) << std::endl; } } diff --git a/tests/manifold/spherical_manifold_01.cc b/tests/manifold/spherical_manifold_01.cc index 1cef855d29..ecba0318cc 100644 --- a/tests/manifold/spherical_manifold_01.cc +++ b/tests/manifold/spherical_manifold_01.cc @@ -158,13 +158,9 @@ main() weights[1] = 1.0/3.0; weights[2] = 1.0/3.0; - Quadrature<3> quad1(points1, weights); - Quadrature<3> quad2(points2, weights); - Quadrature<3> quad3(points3, weights); - - Point<3> Q = manifold.get_new_point(quad1); - Point<3> S = manifold.get_new_point(quad2); - Point<3> T = manifold.get_new_point(quad3); + Point<3> Q = manifold.get_new_point(points1, weights); + Point<3> S = manifold.get_new_point(points2, weights); + Point<3> T = manifold.get_new_point(points3, weights); Point<3> P5(0.707107, 0.707107, 0.0); Point<3> P4(0.0, 0.0, 1.0); diff --git a/tests/manifold/tensor_product_manifold_01.cc b/tests/manifold/tensor_product_manifold_01.cc index a1ad09f494..bde990c3c7 100644 --- a/tests/manifold/tensor_product_manifold_01.cc +++ b/tests/manifold/tensor_product_manifold_01.cc @@ -63,7 +63,7 @@ void test1() w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(Quadrature(sp, w)); + Point ip = manifold.get_new_point(sp, w); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); -- 2.39.5