cubic solver
09-08-2023, 02:26 PM (This post was last modified: 03-03-2024 11:58 PM by Albert Chan.)
Post: #1
 Albert Chan Senior Member Posts: 2,562 Joined: Jul 2018
cubic solver
This program solve x^3 = a*x + b, for x

Code:
cubic_ab(a,b) :=  BEGIN LOCAL d; a, b := a/3, b/2;    /* x^3 = (3a)*x + (2b) */ d := sqrt(b^2-a^3) * (-1)^(sign(re(b))==-1); b := surd(b+d, 3); return [a/b, b];     /* alpha, beta */ END;
Code:
cubic(a,b) := BEGIN LOCAL w := exp(2*pi/3*i); [[1,1],[w,conj(w)],[conj(w),w]] * cubic_ab(a,b); END;

Michael Penn's example

Let x = sin(θ) + cos(θ), it reduced to cubic: x^3 = 3*x - 11/8

Cas> r := cubic(3, -11/8)
Cas> simplify(r)                 → [1/2,    (3*√5-1)/4,       (-3*√5-1)/4]
Cas> float(Ans)                  → [0.5, 1.42705098312, −1.92705098312]
Cas> float(r)                      → [1.42705098312, −1.92705098312, 0.5]

We have Cas simplify bug here! symbolic and approximate numbers don't match.
Luckily, the algorithm doesn't care which cube roots is used (no need for principle root).

But, this bug should be fixed.

Cas> b := cubic_ab(3, -11/8) [2]     → (1/256*(-√15*48*i-176))^(1/3)
Cas> polar(float(b))                         → 1., −0.77627903074
Cas> polar(simplify(b))                    → 1 , atan(√15)

Simplify bug had phase angle off by 2*pi/3:

Cas> atan(√15) - 2*pi/3.0               → −0.77627903074

Update: cubic_ab(a,b) code explained here
09-08-2023, 02:44 PM
Post: #2
 Albert Chan Senior Member Posts: 2,562 Joined: Jul 2018
RE: cubic solver
Root simplify returned wrong (not principle) root bug existed for XCas 1.5.0

But, it is worse on XCas 1.9.0-31 (mingw32 win32)
Result does not match phase *and* magnitude.

XCas> b := (1/256*(-√15*48*i-176))^(1/3)
XCas> polar(float(b))

1.0, -0.77627903074

XCas> simplify(polar(b))                     // good

$$\displaystyle 1\;,\;\frac{-\pi +\arctan \left(\frac{3}{11} \sqrt{15}\right)}{3}$$

XCas> simplify(polar(simplify(b)))       // bad

$$256\;,\;-\pi +\arctan \left(\frac{3}{11} \sqrt{15}\right)$$
09-08-2023, 06:56 PM
Post: #3
 parisse Senior Member Posts: 1,311 Joined: Dec 2013
RE: cubic solver
Yes, there is an additional check in recent versions of Xcas, and here it failed because an argument was not reduced mod 2*pi.
https://github.com/geogebra/giac/commit/...ea8fd94c9d
 « Next Oldest | Next Newest »

User(s) browsing this thread: 1 Guest(s)