Cubic Bézier curve arc length parametrization reversal: find t given a length

Question

I am following this paper Approximate Arc Length Parametrization, M. Walter & A. Fournier, 1996 and have succesfully implemented the direct solution, as in finding the length $s(t)$ given $t$. This algorithm calculates the length of a curve given by $B(t)=f(P_0, P_1, P_2, P_3, t)$ at three different points ($t_1=\frac{1}{3}$, $t_2=\frac{2}{3}$, $t_3=1$), namely $s_1=s(\frac{1}{3})$, $s_2=s(\frac{2}{3})$ and $s_3=s(1)$. The length is approximated using Gauss-Legendre integration along the curve. Using these three points (and $t_0=0\Rightarrow s_0=s(0)=0$) we construct a polynomial of degree $n=3$ of the form $\hat{s}(t)=a_tt^3+b_tt^2+c_tt$ that approximates $s(t)$. This is used to lookup curve lengths at any point along the curve.

However, ideally I'd like to reverse this lookup, which the paper says is fairly easy. Just swap $t$ and $s(t)$ in the algorithm and that's it. We construct a polynomial ($n=3$) of the form $\hat{t}(s)=a_ss^3+b_ss^2+c_ss$. I found two ways this can be implemented:

1. Pick points at certain $s$

Calculate $s_3=s(1)$ (ie. the total curve length) and set $s_1=\frac{1}{3}s_3$ and $s_2=\frac{2}{3}s_3$. Now calculate the corresponding $t_1$, $t_2$, and $t_3=1$. The solution will be:

$\begin{bmatrix}t_1\\t_2\\1\end{bmatrix}=\begin{bmatrix}s_1^3&s_1^2&s_1\\s_2^3&s_2^2&s_2\\s_3^3&s_3^2&s_3\end{bmatrix}\begin{bmatrix}a_s\\b_s\\c_s\end{bmatrix}=\begin{bmatrix}\frac{1}{27}s_3^3&\frac{1}{9}s_3^2&\frac{1}{3}s_3\\\frac{8}{27}s_3^3&\frac{4}{9}s_3^2&\frac{2}{3}s_3\\s_3^3&s_3^2&s_3\end{bmatrix}\begin{bmatrix}a_s\\b_s\\c_s\end{bmatrix}$

$\Rightarrow\begin{bmatrix}a_s\\b_s\\c_s\end{bmatrix}= \left(\begin{bmatrix}\frac{1}{27}s_3^3&\frac{1}{9}s_3^2&\frac{1}{3}s_3\\\frac{8}{27}s_3^3&\frac{4}{9}s_3^2&\frac{2}{3}s_3\\s_3^3&s_3^2&s_3\end{bmatrix}\right)^{-1}\begin{bmatrix}t_1\\t_2\\1\end{bmatrix}$

This gives a fairly nice determinant for the matrix inverse, and we can obtain $t_1$ and $t_2$ using the bisection method on the Gauss-Legendre quadrature.

2. Pick points at certain $t$

Calculate $s_1=s(\frac{1}{3})$, $s_2=s(\frac{2}{3})$ and $s_3=s(1)$. Then solve for $a_s$, $b_s$ and $c_s$:

$\begin{bmatrix}t_1\\t_2\\t_3\end{bmatrix}=\begin{bmatrix}\frac{1}{3}\\\frac{2}{3}\\1\end{bmatrix}=\begin{bmatrix}s_1^3&s_1^2&s_1\\s_2^3&s_2^2&s_2\\s_3^3&s_3^2&s_3\end{bmatrix}\begin{bmatrix}a_s\\b_s\\c_s\end{bmatrix}$

$\Rightarrow\begin{bmatrix}a_s\\b_s\\c_s\end{bmatrix}=\begin{bmatrix}s_1^3&s_1^2&s_1\\s_2^3&s_2^2&s_2\\s_3^3&s_3^2&s_3\end{bmatrix}^{-1}\begin{bmatrix}\frac{1}{3}\\\frac{2}{3}\\1\end{bmatrix}$

This should in my opinion give valid results, but unless my implementation is wrong (I double checked most of it!), I'm wondering if there is a misstep in thought?

EDIT

Looks like I was implementing the inverse matrix wrong...seems to work now. I've chosen approach 1 as this will give an equal spacing along the known variable $s$ as the interpolation points. Approach 2 could, as I presume, in interpolation point that are too close to each other to realistically approach the real $t(s)$ function.