How to fix problems solving for values of BSpline function?

Question

I generated a spline function called f using the BSplineFunction with domain {t, 0, 1}.

f = BSplineFunction[{{0., 5.81152}, {-0.909122, 5.73997}, {-1.86805, 5.5031}, 
  {-1.79586, 5.52709}, {-2.40811, 5.28912}, {-3.4109, 4.70527}, {-4.68131,   3.501}, 
  {-4.93131, 2.99958}, {-5.09697, 2.6673}, {-5.34697, 2.16588}, {-5.59697, 1.66446}, 
  {-5.84697, 1.16304}}, SplineWeights -> {10, 10, 10, 1, 1, 1, 1, 1, 1, 10, 10, 10}]

I then wanted to find the points at which the function intersects the piecewise linear function joining the spline's control points. Graphically, I know the intersection occurs within the domain of f.

k[t_] = {-4.93131 - 0.915666 t, 2.99958 - 1.83654 t}

However, running the code:

FindRoot[First[k[0]] == First[f[t]], {t, .5, 0, 1}]

where k[0] is the end of the line, produces a negative number. Is there a way to solve for the intersection points of those two functions?

What's f? Did you mean g? Plus there are more problems (in addition to the extra square bracket in your FindRoot). For instance, First[g[t]] just yields t. So your FindRoot statements reduces to FindRoot[-4.93131==t,{t,0.5}] which will return t->-4.9313. — rhomboidRhipper, Jul 01 '14 at 16:34
@rhomboidRhipper The syntax problems should be fixed. f[t] does not return a value of t...it yields x and y coordinates for the domain 0<t<1. FindRoot[First[k[0]] == First[f[t]], {t, .5, 0, 1}] should yield the value of t at which the x coordinate is equal to -4.93131. Is that correct or am I completely missing something? — Kaisey, Jul 02 '14 at 13:34

score 4 · Answer 1 · answered Sep 06 '17 at 00:55

This is yet another problem that can be solved by explicitly representing the B-spline curve with BSplineBasis[], as was done here:

pts = {{1, 1}, {2, 3}, {3, -1}, {4, 1}, {5, 0}};

(* piecewise linear interpolant *)
lin = Interpolation[pts, InterpolationOrder -> 1];

n = Min[Length[pts] - 1, 3]; (* B-spline degree *)
m = Length[pts];
(* clamped uniform knots for B-spline *)
knots = {ConstantArray[0, n + 1], Range[m - n - 1]/(m - n),
         ConstantArray[1, n + 1]} // Flatten;

{xu, yu} = Transpose[pts];

(* B-spline component functions *)
bf[t_] = xu.Table[BSplineBasis[{n, knots}, i - 1, t], {i, Length[pts]}];
bg[t_] = yu.Table[BSplineBasis[{n, knots}, i - 1, t], {i, Length[pts]}];

From here we use a trick similar to what was done in this answer; use the MeshFunctions option of Plot[] to find initial estimates for the intersection, and then polish those estimates with FindRoot[]:

crs = Cases[Normal[Plot[bg[t] - lin[bf[t]], {t, 0, 1},
                        Mesh -> {{0}}, MeshFunctions -> {#2 &}]], 
            Point[{x_, _}] :> (Through[{bf, bg}[\[FormalT]]] /. 
                               FindRoot[bg[\[FormalT]] - lin[bf[\[FormalT]]],
                                        {\[FormalT], x}]), ∞]
   {{2.3914022076831603, 1.4343911692673563}, {3.707060176019142, 0.4141203520382841}}

Visualize the geometry:

ParametricPlot[{bf[t], bg[t]}, {t, 0, 1}, 
               Epilog -> {{Directive[AbsoluteThickness[1.6], ColorData[97, 2]], 
                           Line[pts]},
                          {Directive[AbsolutePointSize[6], ColorData[97, 3]], 
                           Point[pts]},
                          {Directive[AbsolutePointSize[8], ColorData[97, 4]], 
                           Point[crs]}},
               PlotRange -> {Automatic, {-1, 3}}, PlotRangePadding -> Scaled[.05]]

B-spline and piecewise linear function, with intersections

score 1 · Answer 2 · answered Jan 26 '21 at 05:44

Define the spline for numerical values only and solve for x and y coordinate simultaneously.

g[t_?NumericQ] = f[t]
fr = FindRoot[{First[k[0]], a} == g[t], {t, .7}, {a, 1}]
(*   {t -> 0.674177, a -> 2.99958}   *)
{f[t], k[0]} /. fr
(*   {{-4.93131, 2.99958}, {-4.93131, 2.99958}}   *)

Michael E2 · Answer 3 · 2021-04-02T16:38:33.440

1

Since V10, one can use Indexed[.., 1] instead of First[..] or Part:

FindRoot[Indexed[k[0], 1] == Indexed[f[t], 1], {t, .5, 0, 1}]
(*  {t -> 0.674177}  *)

edited Apr 02 '21 at 16:38

answered Jan 26 '21 at 15:07

Michael E2

235,386
17
334
747

score 1 · Answer 4 · edited Jan 26 '21 at 04:57

Clear["Global`*"]
pts = {{1, 1}, {2, 3}, {3, -1}, {4, 1}, {5, 0}};
f1 = Interpolation[pts, InterpolationOrder -> 1];
f2 = BSplineFunction[pts];
Show[Graphics[{Red, Point[pts], PointSize[0.02], 
   Point[f2 /@ {0.2, 0.4}], Green, Line[pts]}, Axes -> True], 
   ParametricPlot[f2[t], {t, 0, 1}]]

First,determine approximate range of t, in my example tmin = 0.2, tmax = 0.4.

Then,use method of dichotomy to find the point we need.

dichotomy[{t1_, t2_}] :=
 Module[{x1, x2, y1, y2, tmid, x3, y3},
  {x1, x2} = First@f2[#] & /@ {t1, t2};
  {y1, y2} = f1 /@ {x1, x2};
  tmid = (t1 + t2)/2;
  {x3, y3} = f2[tmid];
  If[(y1 - Last@f2@t1)*(f1[x3] - Last@f2@tmid) < 0, {t1, tmid}, {tmid, t2}]
  ]
tmin = 0.2; tmax = 0.4; error = 10^-5.;
ans = First@NestWhile[dichotomy, {0.2, 0.4}, 
     With[{p = f2@First@#}, f1@First@p - Last@p > error] &]

0.306755

Last, check the answer:

f2[ans]
f1@First@f2[ans]

{2.3914, 1.4344}

1.4344

Thank you! @Chenminqi Although, with the shorthand, I'm having trouble understanding the dichotomy function. Would you mind explaining? — Kaisey, Jul 01 '14 at 16:46

How to fix problems solving for values of BSpline function?

4 Answers4