Programming Challenge: a classic trigonometry problem

[Albert Chan]

(03-09-2023 04:19 AM)Albert Chan Wrote:  Note: X ≈ C if C is big; X ≈ 1 if C is tiny.      → X = [C, C+1]

X guesses have numerical issues if C = ε.

F(X) = (X-C) - 1/X^3

X1 = ε      → Y1 = 0 - 1/ε^3, which is huge
X2 = 1+ε  → Y2 = 1 - 1/(1+ε)^3 ≈ 0, numerically has order of machine epsilon.

Secant line slope = (Y2-Y1)/(X2-X1) ≈ (0-(-1/ε^3)) / ((1+ε)-ε) = 1/ε^3

newX = X1 - Y1/slope ≈ ε - (-1/ε^3) / (1/ε^3) = 1+ε = X2

It is unable to improve. Worse, it assumed convergence, wrongly report X2 as the answer.

lua> C = 1e-6
lua> S.secant(fn'X: (X-C)-X^-3', C, C+1, 1e-9, true) -- BAD!!!
1.000001      1.000001
1.000001

We wanted (X1, X2) closer, so that (Y1, Y2) about same magnitude.
Here is my solution, X1, X2 = C+1, C+exp(-C)

If C is tiny, C + exp(-C) ≈ 1 + C²/2

This make (X1,X2) gap small, exactly what we wanted.

lua> C = 1e-6
lua> S.secant(fn'X: (X-C)-X^-3', C+1, C+exp(-C), 1e-9, true) -- GOOD
1.0000000000005        1.000000250000375
1.000000250000375      1.0000002500000937
1.0000002500000937

If C is huge, C + exp(-C) ≈ C (slightly bigger than C, likely a better guess)

Redo previous post example, for X

lua> a, b, c = 40, 30, 15
lua> D = sqrt(a*a-b*b)
lua> C = 4*c / D
lua> S.secant(fn'X: (X-C)-X^-3', C+1, C+exp(-C), 1e-9, true)
2.37132791792656       2.3441836383192522
2.3441836383192522     2.345306973289791
2.345306973289791      2.3453047639451214
2.3453047639451214     2.345304763754418
2.345304763754418