Somethin' ain't right with cardinal splines

[mweatherley]

The issue

There is a mismatch between what the documentation says and how CubicCardinalSpline actually works. Namely, this:

Interpolation

The curve passes through every control point.

As implemented, this is not true. Instead, the curve passes through only the interior control points — that is, it doesn't actually go through the points at the ends.

Example

Here is a complete example (coutesy of @afonsolage):

use bevy::math::prelude::*;
use bevy::math::vec2;

fn main() {
    let control_points = [
        vec2(-1.0, 50.0),
        vec2(0.3, 100.0),
        vec2(0.4, 150.0),
        vec2(1.0, 150.0),
    ];

    let curve = CubicCardinalSpline::new(0.1, control_points).to_curve();

    for position in curve.iter_positions(10) {
        println!("Position: {position:?}");
    }
}

This produces the following output:

Position: Vec2(0.3, 100.0)
Position: Vec2(0.31351003, 102.165)
Position: Vec2(0.32608, 106.32)
Position: Vec2(0.33777002, 111.955)
Position: Vec2(0.34864, 118.56)
Position: Vec2(0.35875, 125.625)
Position: Vec2(0.36816, 132.64)
Position: Vec2(0.37693, 139.095)
Position: Vec2(0.38511997, 144.48)
Position: Vec2(0.39279, 148.28499)
Position: Vec2(0.39999998, 150.0)

As one can see, the curve only extends between 0.3 and 0.4, where the interior control points are defined. If you want it to actually pass through all those points, you need to extend your control sequence to a length of six, so that your desired endpoints become interior.

Background

Cardinal splines connect a series of points P_i with a sequence of cubic segments between P_i and P_{i+1} so that:

1. The curve segment B_i from P_i to P_{i+1} starts at P_i and ends at P_{i+1}
2. The derivatives of B_i and B_{i+1} at P_{i+1} coincide

Property (1) happens by using P_i and P_{i+1} as the outermost control points of the Bézier curve B_i. Property (2) is ensured by making the third control point of B_i, the second control point of B_{i+1}, and P_{i+1} all collinear; the point of doing this is that it makes the union of cubic segments differentiable (i.e. it doesn't have "kinks").

The general construction is such that the derivative at P_i is given by the difference between P_{i-1} and P_{i+1} — i.e., the tangent to the curve at P_i is parallel to the line between the points preceding and following it. I believe this part is implemented at present; however, this aspect of the construction leaves one thing open: what do we do at the endpoints (where P_{i+1} or P_{i-1} doesn't exist)?

A standard choice is to just use the line from P_0 to P_1 to determine the tangent at P_0 (and analogously for the other endpoint). Of course, there are other choices that one could make as well (e.g. just requiring the tangents for those to be specified manually).

Presently, either by accident or on purpose, we just exclude those points. Either way, we are currently misleading about it.

Which way, spline implementor?

There are basically two paths forward here:

1. Update the documentation so that it is truthful about how the endpoints are only used to determined the tangents for the interior points.
2. Update the cardinal spline implementation so that it actually passes through all the control points.

I think it's pretty clear that I would have not written such a long issue if I did not prefer (2) (which I believe does a good job of meeting expectations and providing ergonomics for curve interpolation between a number of points), but there is a downside: in principle, the current implementation allows the user to indirectly specify the tangent at the endpoints by choosing their control points carefully. (2) is also a breaking change.

I am unsure if this is a bug or if the implementation is this way intentionally, so I have written it as a docs issue.

[afonsolage]

Additional context which originated the issue: on discord help channel.

[NthTensor]

Looks like a bug to me. Docs state the classic properties of cardinals. I will need to review the history to see if this issue existed before I reworded the spline docs, it's possible I introduced this.

Let's fix the implementation.