std::vector<Point<2> > &a) const
{
const unsigned int dim = 2;
- std::vector<Point<dim> > line_points (2);
+ const std::array<double, 2> interior_gl_points
+ {{0.5 - 0.5*std::sqrt(1.0/5.0), 0.5 + 0.5*std::sqrt(1.0/5.0)}};
+
// loop over each of the lines, and if it is at the boundary, then first get
// the boundary description and second compute the points on it. if not at
{
// first get the normal vectors at the two vertices of this line
// from the boundary description
- const Boundary<dim> &boundary
- = line->get_triangulation().get_boundary(line->boundary_id());
+ const Manifold<dim> &manifold = line->get_manifold();
- Boundary<dim>::FaceVertexNormals face_vertex_normals;
- boundary.get_normals_at_vertices (line, face_vertex_normals);
+ Manifold<dim>::FaceVertexNormals face_vertex_normals;
+ manifold.get_normals_at_vertices (line, face_vertex_normals);
// then transform them into interpolation points for a cubic
// polynomial
-face_vertex_normals[1][0] * std::sin(alpha)))
-2*c;
- QGaussLobatto<1> quad_points(4);
- const double t1 = quad_points.point(1)[0];
- const double t2 = quad_points.point(2)[0];
+ const double t1 = interior_gl_points[0];
+ const double t2 = interior_gl_points[1];
const double s_t1 = (((-b-c)*t1+b)*t1+c)*t1;
const double s_t2 = (((-b-c)*t2+b)*t2+c)*t2;
};
}
else
- // not at boundary
+ // not at boundary, so just use scaled Gauss-Lobatto points (i.e., use
+ // plain straight lines).
{
- static const StraightBoundary<dim> straight_boundary;
- straight_boundary.get_intermediate_points_on_line (line, line_points);
- a.insert (a.end(), line_points.begin(), line_points.end());
+ // Note that the zeroth Gauss-Lobatto point is a boundary point, so
+ // we push back mapped versions of the first and second.
+ a.push_back((1.0 - interior_gl_points[0])*line->vertex(0) +
+ (interior_gl_points[0]*line->vertex(1)));
+ a.push_back((1.0 - interior_gl_points[1])*line->vertex(0) +
+ (interior_gl_points[1]*line->vertex(1)));
}
}
}
#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/grid_reordering.h>
#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/manifold_lib.h>
#include <deal.II/grid/tria.h>
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/grid/tria_iterator.h>
-#include <deal.II/grid/tria_boundary_lib.h>
#include <deal.II/grid/intergrid_map.h>
#include <deal.II/distributed/tria.h>
// once and then re-arrange all
// interior nodes so that the mesh is
// the least distorted
- HyperShellBoundary<3> boundary (p);
+ SphericalManifold<3> boundary (p);
Triangulation<3> tmp;
hyper_shell (tmp, p, inner_radius, outer_radius, 12);
- tmp.set_boundary(0, boundary);
- tmp.set_boundary(1, boundary);
+ tmp.set_manifold(0, boundary);
+ tmp.set_manifold(1, boundary);
tmp.refine_global (1);
// let's determine the distance at
return true;
}
};
+
+
+
+ template <int dim, int spacedim>
+ const Manifold<dim, spacedim> &
+ get_default_flat_manifold()
+ {
+ static const FlatManifold<dim, spacedim> flat_manifold;
+ return flat_manifold;
+ }
}
}
//if we have found an entry, return it
return *(it->second);
}
- else
- {
- //if we have not found an entry connected with number, we return
- //straight_boundary
- return straight_boundary;
- }
+
+ // if we have not found an entry connected with number, we return
+ // the default (flat) manifold
+ return internal::Triangulation::get_default_flat_manifold<dim,spacedim>();
}