Rational B-spline manifold
Non-uniform rational basis spline (NURBS) is also supported in BasicBSpline.jl package.
Rational B-spline manifold
Rational B-spline manifold is a parametric representation of a shape.
For given $d$-dimensional B-spline basis functions $B_{(i^1,p^1,k^1)} \otimes \cdots \otimes B_{(i^d,p^d,k^d)}$, given points $\bm{a}_{i^1 \dots i^d} \in V$ and real numbers $w_{i^1 \dots i^d} > 0$, rational B-spline manifold is defined by the following equality:
\[\bm{p}(t^1,\dots,t^d; \bm{a}_{i^1 \dots i^d}, w_{i^1 \dots i^d}) =\sum_{i^1,\dots,i^d} \frac{(B_{(i^1,p^1,k^1)} \otimes \cdots \otimes B_{(i^d,p^d,k^d)})(t^1,\dots,t^d) w_{i^1 \dots i^d}} {\sum\limits_{j^1,\dots,j^d}(B_{(j^1,p^1,k^1)} \otimes \cdots \otimes B_{(j^d,p^d,k^d)})(t^1,\dots,t^d) w_{j^1 \dots j^d}} \bm{a}_{i^1 \dots i^d}\]
Where $\bm{a}_{i^1,\dots,i^d}$ are called control points, and $w_{i^1 \dots i^d}$ are called weights.
BasicBSpline.RationalBSplineManifold
— TypeConstruct Rational B-spline manifold from given control points, weights and B-spline spaces.
Examples
julia> using StaticArrays, LinearAlgebra
julia> P = BSplineSpace{2}(KnotVector([0,0,0,1,1,1]))
BSplineSpace{2, Int64, KnotVector{Int64}}(KnotVector([0, 0, 0, 1, 1, 1]))
julia> w = [1, 1/√2, 1]
3-element Vector{Float64}:
1.0
0.7071067811865475
1.0
julia> a = [SVector(1,0), SVector(1,1), SVector(0,1)]
3-element Vector{SVector{2, Int64}}:
[1, 0]
[1, 1]
[0, 1]
julia> M = RationalBSplineManifold(a,w,P); # 1/4 arc
julia> M(0.3)
2-element SVector{2, Float64} with indices SOneTo(2):
0.8973756499953727
0.4412674277525845
julia> norm(M(0.3))
1.0
Properties
Similar to BSplineManifold
, RationalBSplineManifold
supports the following methods and properties.
- currying
refinement
- Affine commutativity