Skip to main content

How to deal with the condition $u_i=u_{i+1}$ in B-Spline basis function?


The defintion of B-Spline basis function as shown below:


Let $\vec{U}=\{u_0,u_1,\ldots,u_m\}$ a nondecreasing sequence of real numbers,i.e, $u_i\leq u_{i+1}\quad i=0,1,2\ldots m-1$


$$N_{i,0}(u)= \begin{cases} 1 & u_i\leq u

Although I know that Mathematica owns a built-in function BSplineBasis, however, I would like to write my auxiliay function $N_{i,p}(u)$ to learn the NURBS theory and Mathematica programming.


NBSpline


Alogrithm:




enter image description here



 (*=======================Caculate N[i,0](u)================================*)
NBSpline[i_Integer, 0, u_Symbol, U : {Sequence[_] ..}] /;
i <= Length[U] - 2 :=
Piecewise[
{{1, U[[i + 1]] <= u < U[[i + 2]]},
{0, u < U[[i + 1]] || u >= U[[i + 2]]}}]


(*=======================Caculate N[i,p](u)================================*)

NBSpline[i_Integer, p_Integer, u_Symbol, U : {Sequence[_] ..}?OrderedQ] /;
p > 0 && i + p <= Length[U] - 2 :=
Module[{ini},
ini = Table[NBSpline[j, 0, u, U], {j, i, i + p}];
First@Simplify@
Nest[
Dot @@@
(Thread@

{Partition[#, 2, 1],
With[{m = i + p - Length@# + 1},
Table[
{(u - U[[k + 1]])/(U[[k + m + 1]] - U[[k + 1]]),
(U[[k + m + 2]] - u)/(U[[k + m + 2]] - U[[k + 2]])}, {k, i, i + Length@# - 2}]]}) &,
ini, p]
]

Test


NBSpline[1, 3, u, {1, 2, 3, 4, 5, 7}] // TraditionalForm



enter image description here



In my function NBSpline I avoid the condition $u_i=u_{i+1}$, because it will occured the case $\frac{0}{0}$


In the book "The NURBS book", it defines this quotient $\frac{0}{0}$ to be zero.


Question


How to deal with the condition $\frac{0}{0}$ that I sometimes need to set it to 0 ? Namely, How to deal with the condition $u_i=u_{i+1}$ in B-Spline basis function?



Answer



Here is one way to deal with repeated entries in U. One can define a function to compute the coefficient, using one rule when $u_i = u_j$ and the general formula otherwise. One might put extra conditions on the patterns in coeff below, but if the function is called only within NBSpline, then one might assume the conditions are met.



ClearAll[coeff];
coeff[u_, i_, j_, U_] /; U[[i]] == U[[j]] := 0;
coeff[u_, i_, j_, U_] := (u - U[[i]])/(U[[j]] - U[[i]])

Then change the definition of NBSpline for p != 0 as follows.


NBSpline[i_Integer, p_Integer, u_Symbol, 
U : {Sequence[_] ..}?OrderedQ] /; p > 0 && i + p <= Length[U] - 2 :=
Module[{ini}, ini = Table[NBSpline[j, 0, u, U], {j, i, i + p}];
First@Simplify@
Nest[Dot @@@ (Thread@{Partition[#, 2, 1],

With[{m = i + p - Length@# + 1},
Table[{
coeff[u, k + 1, k + m + 1, U],
coeff[u, k + m + 2, k + 2, U]},
{k, i, i + Length@# - 2}]]}) &, ini, p]]

Example:


NBSpline[1, 3, u, {1, 2, 2, 4, 5, 7}]

Mathematica graphics



The output of NBSpline[1, 3, u, {1, 2, 3, 4, 5, 7}] agrees with the output in the question.


P.S. The pattern U : {Sequence[_] ..}?OrderedQ is equivalent to U_List?OrderedQ. You might want a check that restricts U to be a list of numbers, since an ordered list of symbols such as {a, b, c} passes the OrderedQ test. The pattern U_?(VectorQ[#, NumericQ] && OrderedQ[#] &) is one way.


Comments

Popular posts from this blog

plotting - Filling between two spheres in SphericalPlot3D

Manipulate[ SphericalPlot3D[{1, 2 - n}, {θ, 0, Pi}, {ϕ, 0, 1.5 Pi}, Mesh -> None, PlotPoints -> 15, PlotRange -> {-2.2, 2.2}], {n, 0, 1}] I cant' seem to be able to make a filling between two spheres. I've already tried the obvious Filling -> {1 -> {2}} but Mathematica doesn't seem to like that option. Is there any easy way around this or ... Answer There is no built-in filling in SphericalPlot3D . One option is to use ParametricPlot3D to draw the surfaces between the two shells: Manipulate[ Show[SphericalPlot3D[{1, 2 - n}, {θ, 0, Pi}, {ϕ, 0, 1.5 Pi}, PlotPoints -> 15, PlotRange -> {-2.2, 2.2}], ParametricPlot3D[{ r {Sin[t] Cos[1.5 Pi], Sin[t] Sin[1.5 Pi], Cos[t]}, r {Sin[t] Cos[0 Pi], Sin[t] Sin[0 Pi], Cos[t]}}, {r, 1, 2 - n}, {t, 0, Pi}, PlotStyle -> Yellow, Mesh -> {2, 15}]], {n, 0, 1}]

plotting - Plot 4D data with color as 4th dimension

I have a list of 4D data (x position, y position, amplitude, wavelength). I want to plot x, y, and amplitude on a 3D plot and have the color of the points correspond to the wavelength. I have seen many examples using functions to define color but my wavelength cannot be expressed by an analytic function. Is there a simple way to do this? Answer Here a another possible way to visualize 4D data: data = Flatten[Table[{x, y, x^2 + y^2, Sin[x - y]}, {x, -Pi, Pi,Pi/10}, {y,-Pi,Pi, Pi/10}], 1]; You can use the function Point along with VertexColors . Now the points are places using the first three elements and the color is determined by the fourth. In this case I used Hue, but you can use whatever you prefer. Graphics3D[ Point[data[[All, 1 ;; 3]], VertexColors -> Hue /@ data[[All, 4]]], Axes -> True, BoxRatios -> {1, 1, 1/GoldenRatio}]

plotting - Adding a thick curve to a regionplot

Suppose we have the following simple RegionPlot: f[x_] := 1 - x^2 g[x_] := 1 - 0.5 x^2 RegionPlot[{y < f[x], f[x] < y < g[x], y > g[x]}, {x, 0, 2}, {y, 0, 2}] Now I'm trying to change the curve defined by $y=g[x]$ into a thick black curve, while leaving all other boundaries in the plot unchanged. I've tried adding the region $y=g[x]$ and playing with the plotstyle, which didn't work, and I've tried BoundaryStyle, which changed all the boundaries in the plot. Now I'm kinda out of ideas... Any help would be appreciated! Answer With f[x_] := 1 - x^2 g[x_] := 1 - 0.5 x^2 You can use Epilog to add the thick line: RegionPlot[{y < f[x], f[x] < y < g[x], y > g[x]}, {x, 0, 2}, {y, 0, 2}, PlotPoints -> 50, Epilog -> (Plot[g[x], {x, 0, 2}, PlotStyle -> {Black, Thick}][[1]]), PlotStyle -> {Directive[Yellow, Opacity[0.4]], Directive[Pink, Opacity[0.4]],