From: David Wells Date: Fri, 21 Jul 2017 13:31:43 +0000 (-0400) Subject: Break the Manifold interface for performance. X-Git-Tag: v9.0.0-rc1~1099^2~4 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b91365f526b6a81e080db8476cb7edf9a89f1cbd;p=dealii.git Break the Manifold interface for performance. This commit changes the interfaces of Manifolds::get_default_points_and_weights Manifold::project_to_manifold Manifold::get_new_point Manifold::add_new_points to use ArrayView instead of std::vector. In addition, the interface of add_new_points has been changed to populate the ArrayView argument instead of appending to the end of the array. This breaks the public interface of Manifold for the sake of improving performance by about 30%: profiling indicates that, when creating a grid with a Manifold, we spend about 30% of our time purely calling new and delete since we must create and destroy so many std::vectors. --- diff --git a/include/deal.II/grid/manifold.h b/include/deal.II/grid/manifold.h index 5dc6d96f64..aa41005992 100644 --- a/include/deal.II/grid/manifold.h +++ b/include/deal.II/grid/manifold.h @@ -20,6 +20,7 @@ /*---------------------------- manifold.h ---------------------------*/ #include +#include #include #include #include @@ -38,6 +39,25 @@ template class Table; */ namespace Manifolds { + /** + * A constexpr helper function that returns the number of + * default points for the structure type pointed to by the given + * MeshIteratorType. See the documentation of + * Manifolds::get_default_points_and_weights() for more information. + */ + template + inline + constexpr std::size_t n_default_points_per_cell() + { + // Note that in C++11 a constexpr function can only have a return + // statement, so we cannot alias the structure dimension + return GeometryInfo::vertices_per_cell + + GeometryInfo::lines_per_cell + + GeometryInfo::quads_per_cell + + GeometryInfo::hexes_per_cell + - 1; // don't count the cell itself, just the bounding objects + } + /** * Given a general mesh iterator, construct a quadrature object that * contains the following points: @@ -85,7 +105,7 @@ namespace Manifolds const bool with_laplace = false) DEAL_II_DEPRECATED; /** - * Given a general mesh iterator, construct vectors of quadrature points and + * Given a general mesh iterator, construct arrays 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 @@ -126,13 +146,13 @@ namespace Manifolds * cell-@>face(f) or cell-@>line(l). */ template - std::pair >, - std::vector > + std::pair, + n_default_points_per_cell()>, + std::array()> > get_default_points_and_weights(const MeshIteratorType &iterator, const bool with_laplace = false); } - /** * Manifolds are used to describe the geometry of boundaries of domains as * well as the geometry of the interior. Manifold objects are therefore @@ -373,20 +393,18 @@ public: */ virtual Point - get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const; + get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const; + /** - * Compute a new set of points that interpolate between the given points - * @p surrounding_points. @p weights is a table with as many columns as - * @p surrounding_points.size(). The number of rows in @p weights determines - * how many new points will be computed and appended to the last input - * argument @p new_points. After exit of this function, the size of - * @p new_points equals the size at entry plus the number of rows in - * @p weights. + * Compute a new set of points that interpolate between the given points @p + * surrounding_points. @p weights is a table with as many columns as @p + * surrounding_points.size(). The number of rows in @p weights must match + * the length of @p new_points. * * In its default implementation, this function simply calls get_new_point() - * on each row of @p weights and appends those points to the output vector + * on each row of @p weights and writes those points into the output array * @p new_points. However, this function is more efficient if multiple new * points need to be generated like in MappingQGeneric and the manifold does * expensive transformations between a chart space and the physical space, @@ -396,14 +414,14 @@ public: * by implementing only the get_new_point() function. * * The implementation does not allow for @p surrounding_points and - * @p new_points to point to the same vector, so make sure to pass different + * @p new_points to point to the same array, so make sure to pass different * objects into the function. */ virtual void - add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const; + add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const; /** * Given a point which lies close to the given manifold, it modifies it and @@ -419,8 +437,8 @@ public: * the default behavior should work out of the box. */ virtual - Point project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const; + Point project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const; /** * Backward compatibility interface. Return the point which shall become @@ -704,17 +722,15 @@ public: */ virtual Point - get_new_point(const std::vector > &surrounding_points, - const std::vector &weights) const; + get_new_point(const ArrayView> &surrounding_points, + const ArrayView &weights) const override; + /** - * Compute a new set of points that interpolate between the given points - * @p surrounding_points. @p weights is a table with as many columns as - * @p surrounding_points.size(). The number of rows in @p weights determines - * how many new points will be computed and appended to the last input - * argument @p new_points. After exit of this function, the size of - * @p new_points equals the size at entry plus the number of rows in - * @p weights. + * Compute a new set of points that interpolate between the given points @p + * surrounding_points. @p weights is a table with as many columns as @p + * surrounding_points.size(). The number of rows in @p weights must match + * the length of @p new_points. * * For this particular implementation, the interpolation of the * @p surrounding_points according to the @p weights is simply performed in @@ -722,9 +738,9 @@ public: */ virtual void - add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const; + add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const; /** * Project to FlatManifold. This is the identity function for flat, @@ -735,8 +751,8 @@ public: */ virtual Point - project_to_manifold (const std::vector > &points, - const Point &candidate) const; + project_to_manifold (const ArrayView> &points, + const Point &candidate) const; /** * Return a vector that, at $\mathbf x_1$, is tangential to @@ -922,17 +938,14 @@ public: */ virtual Point - get_new_point(const std::vector > &surrounding_points, - const std::vector &weights) const; + get_new_point(const ArrayView> &surrounding_points, + const ArrayView &weights) const override; /** - * Compute a new set of points that interpolate between the given points - * @p surrounding_points. @p weights is a table with as many columns as - * @p surrounding_points.size(). The number of rows in @p weights determines - * how many new points will be computed and appended to the last input - * argument @p new_points. After exit of this function, the size of - * @p new_points equals the size at entry plus the number of rows in - * @p weights. + * Compute a new set of points that interpolate between the given points @p + * surrounding_points. @p weights is a table with as many columns as @p + * surrounding_points.size(). The number of rows in @p weights must match + * the length of @p new_points. * * The implementation of this function first transforms the * @p surrounding_points to the chart space by calling pull_back(). Then, new @@ -951,10 +964,9 @@ public: */ virtual void - add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const; - + add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const; /** * Pull back the given point in spacedim to the Euclidean chartdim * dimensional space. @@ -1076,8 +1088,6 @@ private: }; - - /* -------------- declaration of explicit specializations ------------- */ #ifndef DOXYGEN @@ -1130,24 +1140,30 @@ namespace Manifolds 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); + const auto points_and_weights = get_default_points_and_weights(iterator, with_laplace); + static const int spacedim = MeshIteratorType::AccessorType::space_dimension; + return Quadrature + (std::vector>(points_and_weights.first.begin(), + points_and_weights.first.end()), + std::vector(points_and_weights.second.begin(), + points_and_weights.second.end())); } + + template - std::pair >, - std::vector > + std::pair, + n_default_points_per_cell()>, + std::array()> > 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; + const int spacedim = MeshIteratorType::AccessorType::space_dimension; + constexpr std::size_t points_per_cell = n_default_points_per_cell(); - std::pair >, - std::vector > points_weights; + std::pair, points_per_cell>, std::array > + points_weights; // note that the exact weights are chosen such as to minimize the @@ -1157,16 +1173,16 @@ namespace Manifolds switch (dim) { case 1: - points_weights.first.resize(2); - points_weights.second.resize(2); + Assert(points_weights.first.size() == 2, ExcInternalError()); + Assert(points_weights.second.size() == 2, ExcInternalError()); 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: - points_weights.first.resize(8); - points_weights.second.resize(8); + Assert(points_weights.first.size() == 8, ExcInternalError()); + Assert(points_weights.second.size() == 8, ExcInternalError()); for (unsigned int i=0; i<4; ++i) { @@ -1192,9 +1208,10 @@ namespace Manifolds GeometryInfo::vertices_per_cell+ GeometryInfo::lines_per_cell+ GeometryInfo::faces_per_cell; - points_weights.first.resize(np); - points_weights.second.resize(np); - std::vector > *sp3 = reinterpret_cast > *>(&points_weights.first); + Assert(points_weights.first.size() == np, ExcInternalError()); + Assert(points_weights.second.size() == np, ExcInternalError()); + auto *sp3 = reinterpret_cast, n_default_points_per_cell()> *> + (&points_weights.first); unsigned int j=0; diff --git a/include/deal.II/grid/manifold_lib.h b/include/deal.II/grid/manifold_lib.h index 80a7b05bbf..62c2200c54 100644 --- a/include/deal.II/grid/manifold_lib.h +++ b/include/deal.II/grid/manifold_lib.h @@ -22,6 +22,8 @@ #include #include +#include + DEAL_II_NAMESPACE_OPEN /** @@ -233,8 +235,8 @@ public: */ virtual Point - get_new_point (const std::vector > &vertices, - const std::vector &weights) const; + get_new_point (const ArrayView> &vertices, + const ArrayView &weights) const override; /** * The center of the spherical coordinate system. @@ -303,12 +305,12 @@ public: * 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 override; + get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const override; protected: /** - * A vector orthogonal to direcetion. + * A vector orthogonal to the normal direction. */ const Tensor<1,spacedim> normal_direction; @@ -617,11 +619,10 @@ private: * the case for PolarManifold but not for Spherical manifold, so be careful * when using the latter. In case the quality of the manifold is not good * enough, upon mesh refinement it may happen that the transformation to a - * chart inside the get_new_point() or add_new_points() methods produces points - * that are outside the unit cell. Then this class throws an exception of - * type Manifold@::ExcTransformationFailed. In that case, - * the mesh should be refined before attaching this class, as done in the - * following example: + * chart inside the get_new_point() or add_new_points() methods produces + * points that are outside the unit cell. Then this class throws an exception + * of type Mapping::ExcTransformationFailed. In that case, the mesh should be + * refined before attaching this class, as done in the following example: * * @code * SphericalManifold spherical_manifold; @@ -696,17 +697,14 @@ public: */ virtual Point - get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const; + get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const override; /** - * Compute a new set of points that interpolate between the given points - * @p surrounding_points. @p weights is a table with as many columns as - * @p surrounding_points.size(). The number of rows in @p weights determines - * how many new points will be computed and appended to the last input - * argument @p new_points. After exit of this function, the size of - * @p new_points equals the size at entry plus the number of rows in - * @p weights. + * Compute a new set of points that interpolate between the given points @p + * surrounding_points. @p weights is a table with as many columns as @p + * surrounding_points.size(). The number of columns in @p weights must match + * the length of @p new_points. * * The implementation in this class overrides the method in the base class * and computes the new point by a transfinite interpolation. The first step @@ -723,9 +721,9 @@ public: */ virtual void - add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const; + add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const; private: /** @@ -739,16 +737,18 @@ private: * 10 entries of a the indices cell->index(). */ std::array - get_possible_cells_around_points(const std::vector > &surrounding_points) const; + get_possible_cells_around_points(const ArrayView> &surrounding_points) const; /** - * Finalizes the identification of the correct chart and returns the location - * of the surrounding points on the chart. This method internally calls - * @p get_possible_cells_around_points(). + * Finalizes the identification of the correct chart and populates @p + * chart_points with the pullbacks of the surrounding points. This method + * internally calls @p get_possible_cells_around_points(). + * + * Returns an iterator to the cell on which the chart is defined. */ - std::pair::cell_iterator, - std::vector > > - compute_chart_points(const std::vector > &surrounding_points) const; + typename Triangulation::cell_iterator + compute_chart_points(const ArrayView> &surrounding_points, + ArrayView> chart_points) const; /** * Pull back operation into the unit coordinates on the given coarse cell. diff --git a/include/deal.II/opencascade/boundary_lib.h b/include/deal.II/opencascade/boundary_lib.h index 7055db2880..6f576b4c43 100644 --- a/include/deal.II/opencascade/boundary_lib.h +++ b/include/deal.II/opencascade/boundary_lib.h @@ -87,8 +87,8 @@ namespace OpenCASCADE * algorithms. */ virtual Point - project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const; + project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const; private: @@ -149,8 +149,8 @@ namespace OpenCASCADE * projection algorithms. */ virtual Point - project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const; + project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const; private: /** @@ -236,8 +236,8 @@ namespace OpenCASCADE * exception is thrown. */ virtual Point - project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const; + project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const; private: /** diff --git a/source/fe/mapping_manifold.cc b/source/fe/mapping_manifold.cc index 2b98672b4f..f827927517 100644 --- a/source/fe/mapping_manifold.cc +++ b/source/fe/mapping_manifold.cc @@ -213,14 +213,18 @@ MappingManifold:: transform_unit_to_real_cell (const typename Triangulation::cell_iterator &cell, const Point &p) const { - std::vector > vertices; - std::vector weights; + std::array, GeometryInfo::vertices_per_cell> vertices; + std::array::vertices_per_cell> weights; + for (unsigned int v=0; v::vertices_per_cell; ++v) { - vertices.push_back(cell->vertex(v)); - weights.push_back(GeometryInfo::d_linear_shape_function(p,v)); + vertices[v] = cell->vertex(v); + weights[v] = GeometryInfo::d_linear_shape_function(p, v); } - return cell->get_manifold().get_new_point(vertices, weights); + return cell->get_manifold().get_new_point(make_array_view(vertices.begin(), + vertices.end()), + make_array_view(weights.begin(), + weights.end())); } @@ -375,9 +379,9 @@ namespace internal { for (unsigned int point=0; point - get_new_point(data.vertices, - data.cell_manifold_quadrature_weights[point+data_set]); + quadrature_points[point] = data.manifold->get_new_point + (make_array_view(data.vertices), + make_array_view(data.cell_manifold_quadrature_weights[point+data_set])); } } } @@ -413,9 +417,9 @@ namespace internal const Point &p = data.quad.point(point+data_set); // And get its image on the manifold: - const Point P = data.manifold-> - get_new_point(data.vertices, - data.cell_manifold_quadrature_weights[point+data_set]); + const Point P = data.manifold->get_new_point + (make_array_view(data.vertices), + make_array_view(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 @@ -443,8 +447,8 @@ namespace internal data.vertex_weights[j] = GeometryInfo::d_linear_shape_function(np, j); const Point NP= - data.manifold->get_new_point(data.vertices, - data.vertex_weights); + data.manifold->get_new_point(make_array_view(data.vertices), + make_array_view(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 72ab60eaa9..4a3987ffb4 100644 --- a/source/fe/mapping_q_generic.cc +++ b/source/fe/mapping_q_generic.cc @@ -3810,11 +3810,21 @@ add_line_support_points (const typename Triangulation::cell_iterat } else { - tmp_points.resize(2); - tmp_points[0] = cell->vertex(GeometryInfo::line_to_cell_vertices(line_no, 0)); - tmp_points[1] = cell->vertex(GeometryInfo::line_to_cell_vertices(line_no, 1)); - manifold.add_new_points(tmp_points, - support_point_weights_perimeter_to_interior[0], a); + const std::array, 2> vertices + { + { + cell->vertex(GeometryInfo::line_to_cell_vertices(line_no, 0)), + cell->vertex(GeometryInfo::line_to_cell_vertices(line_no, 1)) + } + }; + + const std::size_t n_rows = support_point_weights_perimeter_to_interior[0].size(0); + a.resize(a.size() + n_rows); + auto a_view = make_array_view(a.end() - n_rows, a.end()); + manifold.add_new_points(make_array_view(vertices.begin(), + vertices.end()), + support_point_weights_perimeter_to_interior[0], + a_view); } } } @@ -3890,7 +3900,9 @@ add_quad_support_points(const Triangulation<3,3>::cell_iterator &cell, // extract the points surrounding a quad from the points // already computed. First get the 4 vertices and then the points on // the four lines - tmp_points.resize(4 + 4*(polynomial_degree-1)); + boost::container::small_vector, 200> + tmp_points(GeometryInfo<2>::vertices_per_cell + + GeometryInfo<2>::lines_per_cell*(polynomial_degree-1)); for (unsigned int v=0; v::vertices_per_cell; ++v) tmp_points[v] = a[GeometryInfo<3>::face_to_cell_vertices(face_no,v)]; if (polynomial_degree > 1) @@ -3900,9 +3912,14 @@ add_quad_support_points(const Triangulation<3,3>::cell_iterator &cell, a[GeometryInfo<3>::vertices_per_cell + (polynomial_degree-1)* GeometryInfo<3>::face_to_cell_lines(face_no,line) + i]; - face->get_manifold().add_new_points (tmp_points, + + const std::size_t n_rows = support_point_weights_perimeter_to_interior[1].size(0); + a.resize(a.size() + n_rows); + auto a_view = make_array_view(a.end() - n_rows, a.end()); + face->get_manifold().add_new_points (make_array_view(tmp_points.begin(), + tmp_points.end()), support_point_weights_perimeter_to_interior[1], - a); + a_view); } } } @@ -3924,9 +3941,10 @@ add_quad_support_points(const Triangulation<2,3>::cell_iterator &cell, } else { - std::vector > vertices; + std::array, GeometryInfo<2>::vertices_per_cell> vertices; for (unsigned int i=0; i::vertices_per_cell; ++i) - vertices.push_back(cell->vertex(i)); + vertices[i] = cell->vertex(i); + Table<2,double> weights(Utilities::fixed_power<2>(polynomial_degree-1), GeometryInfo<2>::vertices_per_cell); for (unsigned int q=0, q2=0; q2::cell_iterator &cell, weights(q,i) = GeometryInfo<2>::d_linear_shape_function(point, i); } // TODO: use all surrounding points once Boundary path is removed - cell->get_manifold().add_new_points(vertices, weights, a); + const std::size_t n_rows = weights.size(0); + a.resize(a.size() + n_rows); + auto a_view = make_array_view(a.end() - n_rows, a.end()); + cell->get_manifold().add_new_points(make_array_view(vertices.begin(), + vertices.end()), + weights, + a_view); } } @@ -3989,8 +4013,13 @@ compute_mapping_support_points(const typename Triangulation::cell_ if (all_manifold_ids_are_equal) { - std::vector > vertices(a); - cell->get_manifold().add_new_points(vertices, support_point_weights_cell, a); + const std::size_t n_rows = support_point_weights_cell.size(0); + a.resize(a.size() + n_rows); + auto a_view = make_array_view(a.end() - n_rows, a.end()); + cell->get_manifold().add_new_points(make_array_view(a.begin(), + a.end() - n_rows), + support_point_weights_cell, + a_view); } else switch (dim) @@ -4008,10 +4037,13 @@ compute_mapping_support_points(const typename Triangulation::cell_ add_quad_support_points(cell, a); else { - std::vector > tmp_points(a); - cell->get_manifold().add_new_points(tmp_points, + const std::size_t n_rows = support_point_weights_perimeter_to_interior[1].size(0); + a.resize(a.size() + n_rows); + auto a_view = make_array_view(a.end() - n_rows, a.end()); + cell->get_manifold().add_new_points(make_array_view(a.begin(), + a.end() - n_rows), support_point_weights_perimeter_to_interior[1], - a); + a_view); } break; @@ -4022,10 +4054,13 @@ compute_mapping_support_points(const typename Triangulation::cell_ // then compute the interior points { - std::vector > tmp_points(a); - cell->get_manifold().add_new_points(tmp_points, + const std::size_t n_rows = support_point_weights_perimeter_to_interior[2].size(0); + a.resize(a.size() + n_rows); + auto a_view = make_array_view(a.end() - n_rows, a.end()); + cell->get_manifold().add_new_points(make_array_view(a.begin(), + a.end() - n_rows), support_point_weights_perimeter_to_interior[2], - a); + a_view); } break; diff --git a/source/grid/manifold.cc b/source/grid/manifold.cc index c2d2238f51..1968c74bbc 100644 --- a/source/grid/manifold.cc +++ b/source/grid/manifold.cc @@ -29,27 +29,7 @@ DEAL_II_NAMESPACE_OPEN using namespace Manifolds; -// This structure is used as comparison function for std::sort when sorting -// points according to their weight. -struct CompareWeights -{ -public: - CompareWeights(const std::vector &weights) - : - compare_weights(weights) - {} - - bool operator() (unsigned int a, unsigned int b) const - { - return compare_weights[a] < compare_weights[b]; - } - -private: - const std::vector &compare_weights; -}; - /* -------------------------- Manifold --------------------- */ - template Manifold::~Manifold () {} @@ -59,7 +39,7 @@ Manifold::~Manifold () template Point Manifold:: -project_to_manifold (const std::vector > &, +project_to_manifold (const ArrayView> &, const Point &) const { Assert (false, ExcPureFunctionCalled()); @@ -75,10 +55,9 @@ get_intermediate_point (const Point &p1, const Point &p2, const double w) const { - std::vector > vertices; - vertices.push_back(p1); - vertices.push_back(p2); - return project_to_manifold(vertices, w * p2 + (1-w)*p1 ); + const std::array, 2> vertices {{p1, p2}}; + return project_to_manifold(make_array_view(vertices.begin(), vertices.end()), + w * p2 + (1-w)*p1); } @@ -86,8 +65,8 @@ get_intermediate_point (const Point &p1, template Point Manifold:: -get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const +get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const { const double tol = 1e-10; const unsigned int n_points = surrounding_points.size(); @@ -106,7 +85,11 @@ get_new_point (const std::vector > &surrounding_points, // associative (as for the SphericalManifold). boost::container::small_vector permutation(n_points); std::iota(permutation.begin(), permutation.end(), 0u); - std::sort(permutation.begin(), permutation.end(), CompareWeights(weights)); + std::sort(permutation.begin(), permutation.end(), + [&weights](const std::size_t a, const std::size_t b) + { + return weights[a] < weights[b]; + }); // Now loop over points in the order of their associated weight Point p = surrounding_points[permutation[0]]; @@ -133,22 +116,24 @@ get_new_point (const std::vector > &surrounding_points, template void Manifold:: -add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const +add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const { AssertDimension(surrounding_points.size(), weights.size(1)); - Assert(&surrounding_points != &new_points, - ExcMessage("surrounding_points and new_points cannot be the same " - "array")); - const unsigned int n_points = surrounding_points.size(); - std::vector local_weights(n_points); + const std::size_t n_points = surrounding_points.size(); + // TODO find a better dimension-dependent estimate for the size of this + // vector + boost::container::small_vector local_weights(n_points); for (unsigned int row=0; row Manifold:: get_new_point_on_line (const typename Triangulation::line_iterator &line) const { - const std::pair >, std::vector > points_weights(get_default_points_and_weights(line)); - return get_new_point (points_weights.first,points_weights.second); + const auto points_weights = get_default_points_and_weights(line); + return get_new_point (make_array_view(points_weights.first.begin(), + points_weights.first.end()), + make_array_view(points_weights.second.begin(), + points_weights.second.end())); } @@ -364,8 +352,11 @@ Point Manifold:: get_new_point_on_quad (const typename Triangulation::quad_iterator &quad) const { - const std::pair >, std::vector > points_weights(get_default_points_and_weights(quad)); - return get_new_point (points_weights.first,points_weights.second); + const auto points_weights = get_default_points_and_weights(quad); + return get_new_point (make_array_view(points_weights.first.begin(), + points_weights.first.end()), + make_array_view(points_weights.second.begin(), + points_weights.second.end())); } @@ -492,8 +483,11 @@ Point<3> Manifold<3,3>:: get_new_point_on_hex (const Triangulation<3, 3>::hex_iterator &hex) const { - const std::pair >, std::vector > points_weights(get_default_points_and_weights(hex,true)); - return get_new_point (points_weights.first,points_weights.second); + const auto points_weights = get_default_points_and_weights(hex, true); + return get_new_point (make_array_view(points_weights.first.begin(), + points_weights.first.end()), + make_array_view(points_weights.second.begin(), + points_weights.second.end())); } @@ -505,15 +499,12 @@ Manifold::get_tangent_vector(const Point &x1, { const double epsilon = 1e-8; - std::vector > q; - q.push_back(x1); - q.push_back(x2); - - std::vector w; - w.push_back(epsilon); - w.push_back(1.0-epsilon); - - const Tensor<1,spacedim> neighbor_point = get_new_point (q, w); + const std::array, 2> points {{x1, x2}}; + const std::array weights {{epsilon, 1.0 - epsilon}}; + const Point neighbor_point = get_new_point (make_array_view(points.begin(), + points.end()), + make_array_view(weights.begin(), + weights.end())); return (neighbor_point-x1)/epsilon; } @@ -533,8 +524,8 @@ FlatManifold::FlatManifold (const Tensor<1,spacedim> &periodicity, template Point FlatManifold:: -get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const +get_new_point (const ArrayView> &surrounding_points, + const ArrayView &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!")); @@ -588,15 +579,15 @@ get_new_point (const std::vector > &surrounding_points, template void FlatManifold:: -add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const +add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const { AssertDimension(surrounding_points.size(), weights.size(1)); if (weights.size(0) == 0) return; - const unsigned int n_points = surrounding_points.size(); + const std::size_t n_points = surrounding_points.size(); Tensor<1,spacedim> minP = periodicity; for (unsigned int d=0; d > &surrounding_points, // check whether periodicity shifts some of the points. Only do this if // necessary to avoid memory allocation const Point *surrounding_points_start = &surrounding_points[0]; - std::vector > modified_points; + + boost::container::small_vector, 200> modified_points; bool adjust_periodicity = false; for (unsigned int d=0; d 0) @@ -624,7 +616,9 @@ add_new_points (const std::vector > &surrounding_points, } if (adjust_periodicity == true) { - modified_points = surrounding_points; + modified_points.resize(surrounding_points.size()); + std::copy(surrounding_points.begin(), surrounding_points.end(), + modified_points.begin()); for (unsigned int d=0; d 0) for (unsigned int i=0; i > &surrounding_points, if (new_point[d] < 0) new_point[d] += periodicity[d]; - new_points.push_back(project_to_manifold(surrounding_points, new_point)); + // TODO should this use surrounding_points_start or surrounding_points? + // The older version used surrounding_points + new_points[row] = project_to_manifold(make_array_view(surrounding_points.begin(), + surrounding_points.end()), + new_point); } } - template Point -FlatManifold::project_to_manifold (const std::vector > &/*vertices*/, - const Point &candidate) const +FlatManifold::project_to_manifold +(const ArrayView> &/*vertices*/, + const Point &candidate) const { return candidate; } @@ -717,15 +715,19 @@ ChartManifold::ChartManifold (const Tensor<1,chartdim> &p template Point ChartManifold:: -get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const +get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const { - std::vector > chart_points(surrounding_points.size()); + const std::size_t n_points = surrounding_points.size(); + + boost::container::small_vector, 200> chart_points(n_points); - for (unsigned int i=0; i p_chart = sub_manifold.get_new_point(chart_points,weights); + const Point p_chart = sub_manifold.get_new_point + (make_array_view(chart_points.begin(), chart_points.end()), + weights); return push_forward(p_chart); } @@ -735,25 +737,28 @@ get_new_point (const std::vector > &surrounding_points, template void ChartManifold:: -add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const +add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView> new_points) const { Assert(weights.size(0) > 0, ExcEmptyObject()); AssertDimension(surrounding_points.size(), weights.size(1)); - const unsigned int n_points = surrounding_points.size(); + const std::size_t n_points = surrounding_points.size(); - std::vector > chart_points(n_points); - for (unsigned int i=0; i, 200> chart_points(n_points); + for (std::size_t i=0; i > new_points_on_chart; - new_points_on_chart.reserve(weights.size(0)); - sub_manifold.add_new_points(chart_points, weights, new_points_on_chart); + boost::container::small_vector, 200> new_points_on_chart(weights.size(0)); + sub_manifold.add_new_points(make_array_view(chart_points.begin(), + chart_points.end()), + weights, + make_array_view(new_points_on_chart.begin(), + new_points_on_chart.end())); - for (unsigned int row=0; row &p1, template Point SphericalManifold:: -get_new_point (const std::vector > &vertices, - const std::vector &weights) const +get_new_point (const ArrayView> &vertices, + const ArrayView &weights) const { const unsigned int n_points = vertices.size(); @@ -368,8 +368,8 @@ CylindricalManifold::CylindricalManifold(const Point &d template Point CylindricalManifold:: -get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const +get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const { // First check if the average in space lies on the axis. Point middle; @@ -536,11 +536,12 @@ DerivativeForm<1,chartdim, spacedim> FunctionManifold::push_forward_gradient(const Point &chart_point) const { DerivativeForm<1, chartdim, spacedim> DF; - std::vector > gradients(spacedim); - push_forward_function->vector_gradient(chart_point, gradients); for (unsigned int i=0; igradient(chart_point, i); + for (unsigned int j=0; j weights(GeometryInfo<2>::vertices_per_face); - std::vector > points(GeometryInfo<2>::vertices_per_face); + std::array::vertices_per_face> weights; + std::array, GeometryInfo<2>::vertices_per_face> points; + // note that the views are immutable, but the arrays are not + const auto weights_view = make_array_view(weights.begin(), weights.end()); + const auto points_view = make_array_view(points.begin(), points.end()); + for (unsigned int line=0; line::lines_per_cell; ++line) { const double my_weight = line%2 ? chart_point[line/2] : 1-chart_point[line/2]; @@ -748,7 +753,8 @@ namespace weights[0] = 1. - line_point; weights[1] = line_point; new_point += my_weight * - cell.line(line)->get_manifold().get_new_point(points, weights); + cell.line(line)->get_manifold().get_new_point(points_view, + weights_view); } } @@ -790,8 +796,12 @@ namespace weights_lines[line] = 0; // start with the contributions of the faces - std::vector weights; - std::vector > points; + std::array::vertices_per_cell> weights; + std::array, GeometryInfo<2>::vertices_per_cell> points; + // note that the views are immutable, but the arrays are not + const auto weights_view = make_array_view(weights.begin(), weights.end()); + const auto points_view = make_array_view(points.begin(), points.end()); + for (unsigned int face=0; face::faces_per_cell; ++face) { Point<2> quad_point(chart_point[(face/2+1)%3], chart_point[(face/2+2)%3]); @@ -817,15 +827,14 @@ namespace } else { - points.resize(GeometryInfo<2>::vertices_per_cell); - weights.resize(GeometryInfo<2>::vertices_per_cell); for (unsigned int v=0; v::vertices_per_cell; ++v) { points[v] = cell.vertex(GeometryInfo<3>::face_to_cell_vertices(face,v)); weights[v] = GeometryInfo<2>::d_linear_shape_function(quad_point, v); } new_point += my_weight * - cell.face(face)->get_manifold().get_new_point(points, weights); + cell.face(face)->get_manifold().get_new_point(points_view, + weights_view); } } @@ -860,14 +869,13 @@ namespace } else { - points.resize(2); - weights.resize(2); points[0] = cell.vertex(GeometryInfo<3>::line_to_cell_vertices(line,0)); points[1] = cell.vertex(GeometryInfo<3>::line_to_cell_vertices(line,1)); weights[0] = 1. - line_point; weights[1] = line_point; new_point -= my_weight * - cell.line(line)->get_manifold().get_new_point(points, weights); + cell.line(line)->get_manifold().get_new_point(points_view, + weights_view); } } @@ -995,7 +1003,7 @@ TransfiniteInterpolationManifold template std::array TransfiniteInterpolationManifold -::get_possible_cells_around_points(const std::vector > &points) const +::get_possible_cells_around_points(const ArrayView> &points) const { // The methods to identify cells around points in GridTools are all written // for the active cells, but we are here looking at some cells at the coarse @@ -1012,7 +1020,7 @@ TransfiniteInterpolationManifold typename Triangulation::cell_iterator cell = triangulation->begin(level_coarse), endc = triangulation->end(level_coarse); - std::vector > distances_and_cells; + boost::container::small_vector, 200> distances_and_cells; for ( ; cell != endc; ++cell) { // only consider cells where the current manifold is attached @@ -1069,14 +1077,14 @@ TransfiniteInterpolationManifold template -std::pair::cell_iterator, - std::vector > > - TransfiniteInterpolationManifold - ::compute_chart_points (const std::vector > &surrounding_points) const +typename Triangulation::cell_iterator +TransfiniteInterpolationManifold +::compute_chart_points (const ArrayView> &surrounding_points, + ArrayView> chart_points) const { - std::pair::cell_iterator, - std::vector > > chart_points; - chart_points.second.resize(surrounding_points.size()); + Assert(surrounding_points.size() == chart_points.size(), + ExcMessage("The chart points array view must be as large as the " + "surrounding points array view.")); std::array nearby_cells = get_possible_cells_around_points(surrounding_points); @@ -1092,11 +1100,11 @@ std::pair::cell_iterator, bool inside_unit_cell = true; for (unsigned int i=0; i::is_inside_unit_cell(chart_points.second[i], + if (GeometryInfo::is_inside_unit_cell(chart_points[i], 1e-6) == false) { inside_unit_cell = false; @@ -1105,16 +1113,14 @@ std::pair::cell_iterator, } if (inside_unit_cell == true) { - chart_points.first = cell; - return chart_points; + return cell; } } // a valid inversion should have returned a point above. AssertThrow(false, (typename Mapping::ExcTransformationFailed())); - chart_points.second.clear(); - return chart_points; + return typename Triangulation::cell_iterator(); } @@ -1122,16 +1128,17 @@ std::pair::cell_iterator, template Point TransfiniteInterpolationManifold -::get_new_point (const std::vector > &surrounding_points, - const std::vector &weights) const +::get_new_point (const ArrayView> &surrounding_points, + const ArrayView &weights) const { - const std::pair::cell_iterator, - std::vector > > chart_points = - compute_chart_points(surrounding_points); + boost::container::small_vector, 100> chart_points(surrounding_points.size()); + ArrayView> chart_points_view = make_array_view(chart_points.begin(), + chart_points.end()); + const auto cell = compute_chart_points(surrounding_points, chart_points_view); - const Point p_chart = chart_manifold.get_new_point(chart_points.second,weights); + const Point p_chart = chart_manifold.get_new_point (chart_points_view, weights); - return push_forward(chart_points.first, p_chart); + return push_forward(cell, p_chart); } @@ -1139,23 +1146,26 @@ TransfiniteInterpolationManifold template void TransfiniteInterpolationManifold:: -add_new_points (const std::vector > &surrounding_points, - const Table<2,double> &weights, - std::vector > &new_points) const +add_new_points (const ArrayView> &surrounding_points, + const Table<2,double> &weights, + ArrayView > new_points) const { Assert(weights.size(0) > 0, ExcEmptyObject()); AssertDimension(surrounding_points.size(), weights.size(1)); - const std::pair::cell_iterator, - std::vector > > chart_points = - compute_chart_points(surrounding_points); + boost::container::small_vector, 100> chart_points(surrounding_points.size()); + ArrayView> chart_points_view = make_array_view(chart_points.begin(), + chart_points.end()); + const auto cell = compute_chart_points(surrounding_points, chart_points_view); - std::vector > new_points_on_chart; - new_points_on_chart.reserve(weights.size(0)); - chart_manifold.add_new_points(chart_points.second, weights, new_points_on_chart); + boost::container::small_vector, 100> new_points_on_chart(weights.size(0)); + chart_manifold.add_new_points(chart_points_view, + weights, + make_array_view(new_points_on_chart.begin(), + new_points_on_chart.end())); for (unsigned int row=0; row > ps(2); - std::vector ws(2, 0.5); - ps[0] = cell->face(boundary_face) - ->child(0)->vertex(1); - ps[1] = cell->face(GeometryInfo - ::opposite_face[boundary_face]) - ->child(0)->vertex(1); + const std::array, 2> ps + { + { + cell->face(boundary_face)->child(0)->vertex(1), + cell->face(GeometryInfo::opposite_face[boundary_face]) + ->child(0)->vertex(1) + } + }; + const std::array ws {{0.5, 0.5}}; triangulation.vertices[next_unused_vertex] - = cell->get_manifold().get_new_point(ps,ws); + = cell->get_manifold().get_new_point(make_array_view(ps.begin(), + ps.end()), + make_array_view(ws.begin(), + ws.end())); } } } diff --git a/source/grid/tria_accessor.cc b/source/grid/tria_accessor.cc index 59f3f1fa23..106b9615a9 100644 --- a/source/grid/tria_accessor.cc +++ b/source/grid/tria_accessor.cc @@ -1051,10 +1051,11 @@ namespace else { TriaRawIterator > it(obj); - 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); + const auto points_and_weights = Manifolds::get_default_points_and_weights(it, use_laplace); + return obj.get_manifold().get_new_point(make_array_view(points_and_weights.first.begin(), + points_and_weights.first.end()), + make_array_view(points_and_weights.second.begin(), + points_and_weights.second.end())); } } } @@ -1206,8 +1207,8 @@ Point TriaAccessor::intermediate_point (const Point &coordinates) const { // Surrounding points and weights. - std::vector > p(GeometryInfo::vertices_per_cell); - std::vector w(GeometryInfo::vertices_per_cell); + std::array, GeometryInfo::vertices_per_cell> p; + std::array::vertices_per_cell> w; for (unsigned int i=0; i::vertices_per_cell; ++i) { @@ -1215,7 +1216,8 @@ TriaAccessor::intermediate_point (const Point::d_linear_shape_function(coordinates, i); } - return this->get_manifold().get_new_point(p, w); + return this->get_manifold().get_new_point(make_array_view(p.begin(), p.end()), + make_array_view(w.begin(), w.end())); } diff --git a/source/opencascade/boundary_lib.cc b/source/opencascade/boundary_lib.cc index 98284f8326..8125c5bfd2 100644 --- a/source/opencascade/boundary_lib.cc +++ b/source/opencascade/boundary_lib.cc @@ -90,8 +90,8 @@ namespace OpenCASCADE template Point NormalProjectionBoundary:: - project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const + project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const { (void)surrounding_points; #ifdef DEBUG @@ -120,8 +120,8 @@ namespace OpenCASCADE template Point DirectionalProjectionBoundary:: - project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const + project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const { (void)surrounding_points; #ifdef DEBUG @@ -151,8 +151,8 @@ namespace OpenCASCADE template Point NormalToMeshProjectionBoundary:: - project_to_manifold (const std::vector > &surrounding_points, - const Point &candidate) const + project_to_manifold (const ArrayView> &surrounding_points, + const Point &candidate) const { TopoDS_Shape out_shape; Tensor<1,3> average_normal; diff --git a/tests/manifold/chart_manifold_03.cc b/tests/manifold/chart_manifold_03.cc index 5ce140d19d..b7ee85c83f 100644 --- a/tests/manifold/chart_manifold_03.cc +++ b/tests/manifold/chart_manifold_03.cc @@ -117,7 +117,8 @@ void test(unsigned int ref=1) for (unsigned int i=0; i ip = manifold.get_new_point(sp, w); + Point ip = manifold.get_new_point(make_array_view(sp), + make_array_view(w)); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); @@ -70,4 +71,3 @@ int main () return 0; } - diff --git a/tests/manifold/composition_manifold_02.cc b/tests/manifold/composition_manifold_02.cc index 79dacbf933..2b9660fc32 100644 --- a/tests/manifold/composition_manifold_02.cc +++ b/tests/manifold/composition_manifold_02.cc @@ -68,7 +68,8 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(sp, w); + Point ip = manifold.get_new_point(make_array_view(sp), + make_array_view(w)); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); @@ -80,4 +81,3 @@ int main () return 0; } - diff --git a/tests/manifold/composition_manifold_03.cc b/tests/manifold/composition_manifold_03.cc index 986fb87b71..4fc292eaa9 100644 --- a/tests/manifold/composition_manifold_03.cc +++ b/tests/manifold/composition_manifold_03.cc @@ -72,7 +72,8 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(sp, w); + Point ip = manifold.get_new_point(make_array_view(sp), + make_array_view(w)); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); @@ -84,4 +85,3 @@ int main () return 0; } - diff --git a/tests/manifold/composition_manifold_04.cc b/tests/manifold/composition_manifold_04.cc index c38d7193e9..b4b0cd2fc0 100644 --- a/tests/manifold/composition_manifold_04.cc +++ b/tests/manifold/composition_manifold_04.cc @@ -73,7 +73,8 @@ int main () w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(sp, w); + Point ip = manifold.get_new_point(make_array_view(sp), + make_array_view(w)); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); @@ -93,7 +94,8 @@ int main () w[1] = 1.0 - w[0]; Point ip = manifold. - pull_back(manifold.get_new_point(sp, w)); + pull_back(manifold.get_new_point(make_array_view(sp), + make_array_view(w))); ip[0] = w[1]; diff --git a/tests/manifold/flat_manifold_03.cc b/tests/manifold/flat_manifold_03.cc index bde2e81ed8..58d7058dae 100644 --- a/tests/manifold/flat_manifold_03.cc +++ b/tests/manifold/flat_manifold_03.cc @@ -76,7 +76,8 @@ void test(unsigned int ref=1) for (unsigned int i=0; i ip = manifold.get_new_point(p, w); + Point ip = manifold.get_new_point(make_array_view(p), + make_array_view(w)); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, p[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, p[1]); @@ -89,4 +90,3 @@ int main () return 0; } - diff --git a/tests/manifold/polar_manifold_08.cc b/tests/manifold/polar_manifold_08.cc index 84c64c718e..0ab7a76dc8 100644 --- a/tests/manifold/polar_manifold_08.cc +++ b/tests/manifold/polar_manifold_08.cc @@ -66,12 +66,11 @@ main() weights[0] = (double)i/((double)n_intermediates-1); weights[1] = 1.0-weights[0]; - deallog << manifold.get_new_point(points, weights) << std::endl; + deallog << manifold.get_new_point(make_array_view(points), + make_array_view(weights)) + << std::endl; } } return 0; } - - - diff --git a/tests/manifold/projection_manifold_01.cc b/tests/manifold/projection_manifold_01.cc index 64d65dcc2a..b7d53db2f4 100644 --- a/tests/manifold/projection_manifold_01.cc +++ b/tests/manifold/projection_manifold_01.cc @@ -31,8 +31,8 @@ class MyManifold : public Manifold { public: Point - project_to_manifold (const std::vector > &vertices, - const Point &candidate) const + project_to_manifold (const ArrayView > &vertices, + const Point &candidate) const override { // Shift the y coordinate to 4*x*(1-x) Point p = candidate; @@ -73,4 +73,3 @@ int main () return 0; } - diff --git a/tests/manifold/spherical_manifold_01.cc b/tests/manifold/spherical_manifold_01.cc index ecba0318cc..93bbbbaaeb 100644 --- a/tests/manifold/spherical_manifold_01.cc +++ b/tests/manifold/spherical_manifold_01.cc @@ -158,9 +158,12 @@ main() weights[1] = 1.0/3.0; weights[2] = 1.0/3.0; - 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> Q = manifold.get_new_point(make_array_view(points1), + make_array_view(weights)); + Point<3> S = manifold.get_new_point(make_array_view(points2), + make_array_view(weights)); + Point<3> T = manifold.get_new_point(make_array_view(points3), + make_array_view(weights)); Point<3> P5(0.707107, 0.707107, 0.0); Point<3> P4(0.0, 0.0, 1.0); @@ -177,6 +180,3 @@ main() // Quadrature (const std::vector< Point< dim > > &points, const std::vector< double > &weights); return 0; } - - - diff --git a/tests/manifold/tensor_product_manifold_01.cc b/tests/manifold/tensor_product_manifold_01.cc index 5301ac4250..598f003d24 100644 --- a/tests/manifold/tensor_product_manifold_01.cc +++ b/tests/manifold/tensor_product_manifold_01.cc @@ -61,7 +61,8 @@ void test1() w[0] = 1.0-(double)i/((double)n_intermediates); w[1] = 1.0 - w[0]; - Point ip = manifold.get_new_point(sp, w); + Point ip = manifold.get_new_point(make_array_view(sp), + make_array_view(w)); Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); @@ -83,4 +84,3 @@ int main () return 0; } - diff --git a/tests/manifold/transfinite_manifold_01.cc b/tests/manifold/transfinite_manifold_01.cc index 7c83b2819f..1abfbecb70 100644 --- a/tests/manifold/transfinite_manifold_01.cc +++ b/tests/manifold/transfinite_manifold_01.cc @@ -43,14 +43,18 @@ void do_test(const Triangulation &tria) std::vector weights(2); weights[0] = 0.1; weights[1] = 0.9; - Point p = cell->get_manifold().get_new_point(points, weights); - Point pref = cell->face(face)->get_manifold().get_new_point(points, weights); + Point p = cell->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); + Point pref = cell->face(face)->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); deallog << "Distance between cell manifold and face manifold: " << (pref-p) << std::endl; weights[0] = 0.55; weights[1] = 0.45; - p = cell->get_manifold().get_new_point(points, weights); - pref = cell->face(face)->get_manifold().get_new_point(points, weights); + p = cell->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); + pref = cell->face(face)->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); deallog << "Distance between cell manifold and face manifold: " << (pref-p) << std::endl; } @@ -147,4 +151,3 @@ int main () return 0; } - diff --git a/tests/manifold/transfinite_manifold_02.cc b/tests/manifold/transfinite_manifold_02.cc index b3d47d4c07..32f2c0dc23 100644 --- a/tests/manifold/transfinite_manifold_02.cc +++ b/tests/manifold/transfinite_manifold_02.cc @@ -43,14 +43,18 @@ void do_test(const Triangulation &tria) std::vector weights(2); weights[0] = 0.1; weights[1] = 0.9; - Point p = cell->get_manifold().get_new_point(points, weights); - Point pref = cell->face(face)->get_manifold().get_new_point(points, weights); + Point p = cell->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); + Point pref = cell->face(face)->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); deallog << "Distance between cell manifold and face manifold: " << (pref-p) << std::endl; weights[0] = 0.55; weights[1] = 0.45; - p = cell->get_manifold().get_new_point(points, weights); - pref = cell->face(face)->get_manifold().get_new_point(points, weights); + p = cell->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); + pref = cell->face(face)->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); deallog << "Distance between cell manifold and face manifold: " << (pref-p) << std::endl; } @@ -156,4 +160,3 @@ int main () return 0; } - diff --git a/tests/manifold/transfinite_manifold_04.cc b/tests/manifold/transfinite_manifold_04.cc index c427246c88..37e635d82f 100644 --- a/tests/manifold/transfinite_manifold_04.cc +++ b/tests/manifold/transfinite_manifold_04.cc @@ -43,14 +43,18 @@ void do_test(const Triangulation &tria) std::vector weights(2); weights[0] = 0.1; weights[1] = 0.9; - Point p = cell->get_manifold().get_new_point(points, weights); - Point pref = cell->face(face)->get_manifold().get_new_point(points, weights); + Point p = cell->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); + Point pref = cell->face(face)->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); deallog << "Distance between cell manifold and face manifold: " << (pref-p) << std::endl; weights[0] = 0.55; weights[1] = 0.45; - p = cell->get_manifold().get_new_point(points, weights); - pref = cell->face(face)->get_manifold().get_new_point(points, weights); + p = cell->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); + pref = cell->face(face)->get_manifold().get_new_point(make_array_view(points), + make_array_view(weights)); deallog << "Distance between cell manifold and face manifold: " << (pref-p) << std::endl; } @@ -118,5 +122,3 @@ int main () return 0; } - -