Skip to content

Commit

Permalink
Merge pull request #112 from oberbichler/feature/span-tools
Browse files Browse the repository at this point in the history 
Extend span functions
  • Loading branch information
@oberbichler
oberbichler authored Apr 24, 2020
2 parents 3876fe4 + 6d49c61 commit 47dff0e
Show file tree
Hide file tree
Showing 3 changed files with 129 additions and 21 deletions.
87 changes: 73 additions & 14 deletions include/anurbs/Algorithm/BrepFaceIntegrationPoints.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,8 +16,7 @@ class BrepFaceIntegrationPoints
using Vector = linear_algebra::Vector<3>;

public: // static methods
static IntegrationPointList<2> get(const BrepFace& face,
const double tolerance)
static IntegrationPointList<2> get(const BrepFace& face, const double tolerance)
{
TrimmedSurfaceClipping clipper(tolerance * 10, tolerance / 10);

Expand Down Expand Up @@ -48,12 +47,10 @@ class BrepFaceIntegrationPoints
}

if (clipper.span_trim_type(i, j) == TrimTypes::Full) {
const auto points = IntegrationPoints::get(degree_u + 1,
degree_v + 1, clipper.span_u(i), clipper.span_v(j));
const auto points = IntegrationPoints::get(degree_u + 1, degree_v + 1, clipper.span_u(i), clipper.span_v(j));

for (const auto [u, v, weight] : points) {
const auto ders =
surface_geometry.derivatives_at(u, v, 1);
const auto ders = surface_geometry.derivatives_at(u, v, 1);

const auto n = cross(ders[1], ders[2]);

Expand All @@ -63,17 +60,14 @@ class BrepFaceIntegrationPoints
const auto polygons = clipper.span_polygons(i, j);

for (const auto& polygon : polygons) {
const auto points =
PolygonIntegrationPoints::get(degree, polygon);
const auto points = PolygonIntegrationPoints::get(degree, polygon);

for (const auto [u, v, weight] : points) {
const auto ders =
surface_geometry.derivatives_at(u, v, 1);
const auto ders = surface_geometry.derivatives_at(u, v, 1);

const auto n = cross(ders[1], ders[2]);

integration_points.emplace_back(u, v,
weight * norm(n));
integration_points.emplace_back(u, v, weight * norm(n));
}
}
}
Expand All @@ -83,6 +77,71 @@ class BrepFaceIntegrationPoints
return integration_points;
}

static std::vector<std::tuple<Index, Index, IntegrationPointList<2>>> get_with_spans(const BrepFace& face, const double tolerance)
{
TrimmedSurfaceClipping clipper(tolerance * 10, tolerance / 10);

std::vector<std::tuple<Index, Index, IntegrationPointList<2>>> integration_points;

for (const auto loop : face.loops()) {
clipper.begin_loop();

for (const auto trim : loop->trims()) {
clipper.add_curve(*trim->curve_2d());
}

clipper.end_loop();
}

const auto& surface_geometry = *face.surface_geometry().data();
clipper.compute(surface_geometry.spans_u(), surface_geometry.spans_v());

const auto degree_u = surface_geometry.degree_u();
const auto degree_v = surface_geometry.degree_v();

const auto degree = std::max(degree_u, degree_v) + 1;

for (Index i = 0; i < clipper.nb_spans_u(); i++) {
for (Index j = 0; j < clipper.nb_spans_v(); j++) {
if (clipper.span_trim_type(i, j) == TrimTypes::Empty) {
continue;
}

IntegrationPointList<2> span_integration_points;

if (clipper.span_trim_type(i, j) == TrimTypes::Full) {
const auto points = IntegrationPoints::get(degree_u + 1, degree_v + 1, clipper.span_u(i), clipper.span_v(j));

for (const auto [u, v, weight] : points) {
const auto ders = surface_geometry.derivatives_at(u, v, 1);

const auto n = cross(ders[1], ders[2]);

span_integration_points.emplace_back(u, v, weight * norm(n));
}
} else {
const auto polygons = clipper.span_polygons(i, j);

for (const auto& polygon : polygons) {
const auto points = PolygonIntegrationPoints::get(degree, polygon);

for (const auto [u, v, weight] : points) {
const auto ders = surface_geometry.derivatives_at(u, v, 1);

const auto n = cross(ders[1], ders[2]);

span_integration_points.emplace_back(u, v, weight * norm(n));
}
}
}

integration_points.push_back({i, j, span_integration_points});
}
}

return integration_points;
}

public: // python
static void register_python(pybind11::module& m)
{
Expand All @@ -91,8 +150,8 @@ class BrepFaceIntegrationPoints

using Type = BrepFaceIntegrationPoints;

m.def("integration_points", [](const BrepFace& face, double tolerance) {
return Type::get(face, tolerance); }, "face"_a, "tolerance"_a);
m.def("integration_points", [](const BrepFace& face, double tolerance) { return Type::get(face, tolerance); }, "face"_a, "tolerance"_a);
m.def("integration_points_with_spans", [](const BrepFace& face, double tolerance) { return Type::get_with_spans(face, tolerance); }, "face"_a, "tolerance"_a);
}
};

Expand Down
31 changes: 25 additions & 6 deletions include/anurbs/Geometry/NurbsCurveGeometry.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -291,8 +291,7 @@ struct NurbsCurveGeometry : public CurveBase<TDimension>
return point;
}

std::pair<std::vector<Index>, linear_algebra::MatrixXd> shape_functions_at(const double t,
const Index order) const
std::pair<std::vector<Index>, Eigen::MatrixXd> shape_functions_at(const double t, const Index order) const
{
NurbsCurveShapeFunction shape_function(degree(), order);

Expand Down Expand Up @@ -334,6 +333,24 @@ struct NurbsCurveGeometry : public CurveBase<TDimension>
return result;
}

Index span_at(const double t) const
{
return Nurbs::upper_span(degree(), knots(), t);
}

std::vector<Index> nonzero_pole_indices_at_span(const Index span) const
{
std::vector<Index> result(degree() + 1);

const Index first_nonzero_pole_index = span - degree() + 1;

for (Index i = 0; i < length(result); i++) {
result[i] = first_nonzero_pole_index + i;
}

return result;
}

double greville_point(Index index) const
{
double u = 0.0;
Expand Down Expand Up @@ -408,9 +425,9 @@ struct NurbsCurveGeometry : public CurveBase<TDimension>
py::class_<Type, Base, Holder>(m, name.c_str())
// constructors
.def(py::init<const Index, const Index, const bool>(), "degree"_a, "nb_poles"_a, "is_rational"_a)
.def(py::init<const Index, const Eigen::VectorXd, const Poles>(), "degree"_a, "knots"_a, "poles"_a)
.def(py::init<const Index, const Eigen::VectorXd, const Poles, const Eigen::VectorXd>(), "degree"_a, "knots"_a, "poles"_a, "weights"_a)
.def(py::init<const Index, const Eigen::VectorXd, const std::vector<ControlPoint>>(), "degree"_a, "knots"_a, "control_points"_a)
.def(py::init<const Index, const Knots, const Poles>(), "degree"_a, "knots"_a, "poles"_a)
.def(py::init<const Index, const Knots, const Poles, const Weights>(), "degree"_a, "knots"_a, "poles"_a, "weights"_a)
.def(py::init<const Index, const Knots, const std::vector<ControlPoint>>(), "degree"_a, "knots"_a, "control_points"_a)
// read-only properties
.def_property_readonly("is_rational", &Type::is_rational)
.def_property_readonly("nb_knots", &Type::nb_knots)
Expand All @@ -420,8 +437,10 @@ struct NurbsCurveGeometry : public CurveBase<TDimension>
.def_property("poles", py::overload_cast<>(&Type::poles), &Type::set_poles)
.def_property("weights", py::overload_cast<>(&Type::weights), &Type::set_weights)
// methods
.def("shape_functions_at", &Type::shape_functions_at, "t"_a, "order"_a)
.def("greville_point", &Type::greville_point, "index"_a)
.def("nonzero_pole_indices_at_span", &Type::nonzero_pole_indices_at_span, "span"_a)
.def("shape_functions_at", &Type::shape_functions_at, "t"_a, "order"_a)
.def("span_at", &Type::span_at, "t"_a)
// .def("clone", &Type::clone)
;

Expand Down
32 changes: 31 additions & 1 deletion include/anurbs/Geometry/NurbsSurfaceGeometry.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -414,6 +414,34 @@ class NurbsSurfaceGeometry : public SurfaceBase<TDimension> {
m_weights = value;
}

//

std::pair<Index, Index> span_at(const double u, const double v) const
{
const auto span_u = Nurbs::upper_span(degree_u(), knots_u(), u);
const auto span_v = Nurbs::upper_span(degree_v(), knots_v(), v);

return {span_u, span_v};
}

std::vector<Index> nonzero_pole_indices_at_span(const Index span_u, const Index span_v) const
{
std::vector<Index> indices((degree_u() + 1) * (degree_v() + 1));

const Index first_nonzero_pole_index_u = span_u - degree_u() + 1;
const Index first_nonzero_pole_index_v = span_v - degree_v() + 1;

auto it = indices.begin();

for (Index i = 0; i < degree_u() + 1; i++) {
for (Index j = 0; j < degree_v() + 1; j++) {
*(it++) = to_single_index(first_nonzero_pole_index_u + i, first_nonzero_pole_index_v + j);
}
}

return indices;
}

template <typename TValue, typename TValues>
TValue evaluate_at(TValues values, const double u, const double v) const
{
Expand Down Expand Up @@ -691,7 +719,9 @@ class NurbsSurfaceGeometry : public SurfaceBase<TDimension> {
.def("shape_functions_at", &Type::shape_functions_at, "u"_a, "v"_a, "order"_a)
.def("weight", (double (Type::*)(const Index) const) & Type::weight, "index"_a)
.def("weight", (double (Type::*)(const Index, const Index) const) & Type::weight, "index_u"_a, "index_v"_a)
.def("greville_point", &Type::greville_point, "index_u"_a, "index_v"_a);
.def("greville_point", &Type::greville_point, "index_u"_a, "index_v"_a)
.def("nonzero_pole_indices_at_span", &Type::nonzero_pole_indices_at_span, "span_u"_a, "span_v"_a)
.def("span_at", &Type::span_at, "u"_a, "v"_a);

m.def("add", [](Model& model, Holder data) { model.add<Type>(data); });

Expand Down

0 comments on commit 47dff0e

Please sign in to comment.