Talk:ITP method

How to tune kappa parameters?
One of the weak spots in this article is how to tune the $$ \kappa $$ parameters, particularly $$ \kappa_1 $$. The original paper discusses the question a bit, but doesn't itself give any clear guidance. After a bit of experimentation (mostly on inverse arc length problems), I'm pretty confident that $$ \kappa_2 = 2 $$ is a reasonable default choice, suitable in most cases unless experiment or analysis strongly suggests a different value. But $$ \kappa_1 $$ is a bit of a different question. The current example in this article uses a value of $$ 0.1 $$, which is the same numerical value as the examples in the paper, but I think that fact is a bit misleading. The examples in the paper have $$b - a = 2$$, while the example here has $$b - a = 1$$, and there is scaling involved. If $$ \kappa_2 = 2 $$, then I think the heuristic that captures the intent of the paper is $$ \kappa_1 = 0.2 / (b - a) $$, and this is what I have documented in the kurbo implementation. Also note, if the example were to adopt this recommendation, then $$ n_0 = 0 $$ would work with fast convergence; I think you can say that you need more "slack" as provided by $$ n_0 $$ when $$ \kappa_1 $$ is lowballed.

On that last point, I feel like I might be skating a little close to the edge of WP:NOR here. I personally consider an open source project such as kurbo a reasonably reliable source, but can see how others may feel differently. Ideally there would be some survey paper or textbook chapter that we could cite, but, partly due to the technique being so new, there's precious little out there. I'm going to shamelessly self-promote in the interest of providing useful resources to readers, but open to hearing better ways to handle this. Raph Levien (talk) 16:19, 5 April 2021 (UTC)
 * , I have added a small note to link to your library (I *feel* that's compliant to NOR). Cheers 🦀🍺. --Artoria2e5 🌉 12:56, 4 July 2021 (UTC)
 * Oh, while we are at it, I feel the description has a gap in the description of $$n_0$$. The paper (and this Wikipedia article) only mentions it's &ge; 0, while your documentation seems to suggest a smaller useful range of [0,1] further constrained by the optimization to use usize? --Artoria2e5 🌉 13:13, 4 July 2021 (UTC)

Observations from tuning of epsilon, n0 and k1
As the WP:NOR policy does not apply to Talk pages, I think it's OK to share a some yet-unpublished recommendations here: I have conducted a quasi Monte Carlo parametric study of epsilon and $$n_0$$, using two alternative monotonic curves of the form $$y=x^{((1/c)-1)}$$ and $$y=1-(1-x)^{(c/(1-c))}$$, where $$c$$ is a ranomized constant in the range $$c \in (0,1)$$, and $$x$$ is in the range $$x \in (0,1)$$. Simulations were conducted in Visual Basic (machine floating-point epsilon = 1E-15). $$\kappa_2 = 2$$ was kept constant for reasons given above. My findings are as follows:

Tuning of epsilon

 * I tried values 1E-16 $$\le \epsilon \le$$ 0.1, in steps of factor 10.
 * As $$\epsilon$$ is reduced, the number of iterations increases and precision improves. However this effect levels out for small values of $$\epsilon \le$$ 1E-10.
 * For large values of $$\epsilon$$: 1E-9 to 0.1, there is little benefit of this algorithm compared to basic Bisection, both in terms of number of iterations, and precision of $$x_{ITP}$$.
 * For small values of $$\epsilon$$: 1E-16 to 1E-10, there is negligible difference in performance, both in terms of number of iterations and precision of $$x_{ITP}$$. It levels out, with median 8 iterations.
 * Therefore I recommend to set $$\epsilon$$ equal to the machine epsilon (1E-15 in my case), as there are no drawbacks compared to using lower precision.

Tuning of n0

 * I tried values $$0 \le n_0 \le 20$$ in steps of +1, and from 20 up to 50 in steps of +10. This was repeated for a range of $$\kappa_1$$ values; see below. Results here apply for a wide-range of $$\kappa_1$$ values.
 * $$n_0 = 0$$ is identical to Bisection.
 * $$n_0 = 1$$ halves the number of iterations compared to Bisection, and improved precision of $$x_{ITP}$$ (median 100000 times more accurate, and 3rd quartile 30 times more accurate), assuming small $$\epsilon \le $$1E-10.
 * Increasing $$n_0$$ further (i.e. $$n_0=2$$ and higher) leads to further improvements in precision (median and Q3 statistics) with the added benefit of reducing the Mean number of iterations (median and Q3 remain static). The effect levels out exponentially, with no further improvement over $$17 \le n_0$$.
 * The maximum possible number of iterations increases with $$n_0$$: Max iterations = $$n_{bisection} + n_0$$. However, the likelihood of this hapening sinks for higher $$n_0$$, so the Mean number of iterations actually sinks for higher $$n_0$$, but the effect levels out exponentially, with a minimum in the range $$n_0$$ = 14 to 17, as mentioned above.
 * Thus I recommend $$n_0 = 14$$ as an all-round value. I see little benefit in developing adaptive algorithms that automatically tune $$n_0$$, as I found little difference in the range $$11 \le n_0 \le 20$$.

Tuning of k1
Comment written by Peter.schild (talk) 15:28, 2 January 2023 (UTC)
 * I tried scaling factors in the range $$0.01 <= f <= 1.00$$ in increments of 0.01, in calculating $$\kappa_1 = f/(x_b-x_a)$$.
 * The arithmentic mean number iterations reduces exponentially as $$f$$ increases from zero towards approx 0.33, above which the number of iterations flattens out in the region $$0.33 < f <= 0.58$$, so the optimal value of $$f$$ may be anywhere in this region. Above $$0.58 \le f$$, the number of iterations increaces again gradually.
 * The median number of iterations is minimum in the range $$0.08 <= f <= 0.33$$, and the third-quartile is minimum in the range $$0.33 <= f <= 0.58$$.
 * Thus I recommend $$\kappa_1 = 0.33/(x_b-x_a)$$ as an all-round value, together with $$n_0=14$$ and $$\epsilon=$$1E-15, which worked well for the large range of monotonic function shapes that I tested. This gives a mean number of iterations (9.8) that is a third that of Bisection (30). This finding is obviously influenced by the distribution of monotonic function shapes that I tried.

Projection Step
In the computation of $$\rho_k \equiv \min\left\{\epsilon 2^{n_{1/2}+n_0-j} - \frac{b-a}{2},|x_t-x_{1/2}|\right\} $$, the term $$\epsilon 2^{n_{1/2}+n_0-j} - \frac{b-a}{2}$$ must be non-negative. Round-off errors could result in a negative value for this term, and this can cause the algorithm to be unstable. This is explained in the ITP paper, but it is not shown in the Algorithm in the ITP paper. So it should be $$\rho_k \equiv \min\left\{\max\left\{0,\epsilon 2^{n_{1/2}+n_0-j} - \frac{b-a}{2}\right\},|x_t-x_{1/2}|\right\} $$. One should also add r=max(0,r) in the pseudo-code. One should also subtract the unit of least precision before computing max(0,r), so the full expression should be r=max(0,r-ulp(r)). This is also in the ITP paper. One last interesting fact is that if r is ever equal to 0, then the ITP algorithm is the same as the bisection algorithm, so one can check for this case and switch to the bisection algorithm. KennethBaclawski (talk) 16:34, 22 March 2023 (UTC)


 * P.S. Yes, this is explained in Appendix B in the original paper. For programming languages that don't have the ULP function, one can use $$\rho =\max \{0,\rho_k -\epsilon|\rho_k| \}$$ Peter.schild (talk) 06:42, 16 May 2023 (UTC)