virtual
Point<spacedim>
get_new_point(const std::vector<Point<spacedim> > &surrounding_points,
- const std::vector<double> &weights) const = 0;
+ const std::vector<double> &weights) const = 0;
/**
* Return equally spaced intermediate points between the given
*/
virtual
void
- get_intermediate_points(std::vector<Point<spacedim> > &points,
- const std::vector<Point<spacedim> > &surrounding_points) const;
+ get_intermediate_points(const std::vector<Point<spacedim> > &surrounding_points,
+ std::vector<Point<spacedim> > &points) const;
/**
* Given a point which lies close to the given manifold, it modifies
* normal_vector. Derived classes can overload that function.
*/
virtual
- void get_normals_at_points(std::vector<Point<spacedim> > &normals,
- const std::vector<Point<spacedim> > &points) const;
+ void get_normals_at_points(const std::vector<Point<spacedim> > &points,
+ std::vector<Point<spacedim> > &normals) const;
/**
* Compute the normal to the manifold at the given point. The
* those points. Input arguments must be of the same size.
*/
virtual
- void get_normals_at_points(std::vector<Point<spacedim> > &normals,
- const std::vector<Point<spacedim> > &points) const;
+ void get_normals_at_points(const std::vector<Point<spacedim> > &points,
+ std::vector<Point<spacedim> > &normals) const;
};
ps[0] = line->vertex(0);
ps[1] = line->vertex(1);
- line->get_manifold().get_intermediate_points (line_points, ps);
+ line->get_manifold().get_intermediate_points (ps, line_points);
if (dim==3)
{
vertices.resize(4);
for(unsigned int i=0;i<4; ++i)
vertices[i] = face->vertex(i);
- face->get_manifold().get_intermediate_points(quad_points, vertices);
+ face->get_manifold().get_intermediate_points(vertices, quad_points);
// in 3D, the orientation, flip and rotation of the face might not
// match what we expect here, namely the standard orientation. thus
// reorder points accordingly. since a Mapping uses the same shape
vertices.resize(4);
for(unsigned int i=0;i<4; ++i)
vertices[i] = face->vertex(i);
- face->get_manifold().get_intermediate_points (quad_points, vertices);
+ face->get_manifold().get_intermediate_points (vertices, quad_points);
// in 3D, the orientation, flip and rotation of the face might
// not match what we expect here, namely the standard
// orientation. thus reorder points accordingly. since a Mapping
for(unsigned int i=0; i<4; ++i)
vertices[i] = cell->vertex(i);
- cell->get_manifold().get_intermediate_points (quad_points, vertices);
+ cell->get_manifold().get_intermediate_points (vertices, quad_points);
for (unsigned int i=0; i<quad_points.size(); ++i)
a.push_back(quad_points[i]);
}
template <int spacedim>
void
-Manifold<spacedim>::get_normals_at_points(std::vector<Point<spacedim> > &normals,
- const std::vector<Point<spacedim> > &points) const
+Manifold<spacedim>::get_normals_at_points(const std::vector<Point<spacedim> > &points,
+ std::vector<Point<spacedim> > &normals) const
{
AssertDimension(normals.size(), points.size());
for(unsigned int i=0; i<normals.size(); ++i)
template <int spacedim>
void
-Manifold<spacedim>::get_intermediate_points(std::vector<Point<spacedim> > &points,
- const std::vector<Point<spacedim> > &surrounding_points) const
+Manifold<spacedim>::get_intermediate_points(const std::vector<Point<spacedim> > &surrounding_points,
+ std::vector<Point<spacedim> > &points) const
{
Assert(surrounding_points.size() >= 2, ExcMessage("At least 2 surrounding points are required"));
const unsigned int n=points.size();
}
-// template<int dim>
-// CellData<dim>::CellData()
-// {
-// // Set the face to be interior
-// boundary_id = numbers::internal_face_boundary_id;
-// // And the manifold to be invalid
-// manifold_id = numbers::invalid_manifold_id;
-// }
-
-
bool
SubCellData::check_consistency (const unsigned int dim) const
{
template <int dim, int spacedim>
void
-StraightBoundary<dim, spacedim>::get_normals_at_points(std::vector<Point<spacedim> > &normals,
- const std::vector<Point<spacedim> > &points) const
+StraightBoundary<dim, spacedim>::get_normals_at_points(const std::vector<Point<spacedim> > &points,
+ std::vector<Point<spacedim> > &normals) const
{
AssertDimension(normals.size(), points.size());
switch(dim)