Solve Quartic with Casio FX-115MS

[Albert Chan]

y^4 + A y^2 + B y + C = 0
y^4 + 2D y^2 + D^2 = (2D - A) y^2 - B y + (D^2 - C)

For some value of D, RHS can be made a perfect square:

(y^2 + D)^2 = (2D - A) * (y - (B/2)/(2D-A))^2

Let E = √(2D-A), F = -B/(2E), Quartic => (y^2 + E y + D+F)(y^2 - E y + D-F) = 0

Example: f(x) = x^4 + 2 x^3 + 3 x^2 + 4 x + 5 = 0

Depress the Quartic (remove cubic term) by substitute x = y + k, where k = -b/(4a)
We get depressed quartic: y^4 + A y^2 + B y + C

k = -b/(4a) = -2/4 = -0.5
C = f(k) = 3.5625
B = f(k−1) = 4.0625
A = f(k+1) = 8.0625
B = (A - B) / 2 = 2
A = A - 1 - B - C = 1.5

Find D such that RHS is perfect square, i.e. discriminant = B^2 - 4(2D - A)(D^2 - C) = 0
We do not want complex values for E, so 2D ≥ A, D^2 ≥ C:

D ≥ max(A/2, √C) ≈ max(0.75, 1.904) = 1.904

A trick to force Casio FX-115MS to solve for D: 0ABC + B^2 - 4(2D - A)(D^2 - C)

D = 1.99130646 (Casio SOLVE, with guessed D = 2)
E = √(2D-A) = 1.575630959
F = -B/(2E) = -0.6346663818

Quadratic Solver, enter a=1, b=+E, c=D+F, roots = −0.78781548 ± 0.857896758 i
Quadratic Solver, enter a=1, b=−E, c=D-F , roots = +0.78781548 ± 1.416093080 i

x roots = y roots + k = -1.28781548 ± 0.857896758 i, 0.28781548 ± 1.416093080 i

For calculator without Quad Solver, we can use this:

y^2 + E y + (D+F) = 0 ==> y = −E/2 ± √(−F - D/2 - A/4)
y^2 − E y + (D−F) = 0 ==> y = +E/2 ± √(+F - D/2 - A/4)

Since C=(D+F)(D-F), we might derive 1 from the other, hopefully better accuracy too.

[Albert Chan]

(11-29-2018 05:28 PM)Albert Chan Wrote:  A trick to force Casio FX-115MS to solve for D: 0ABC + B^2 - 4(2D - A)(D^2 - C)

Discover a simpler way to solve for D, without the 0ABC + ... hack

B^2 - 4(2D - A)(D^2 - C) [SOLVE]

Enter the values like before, but scroll back to D= line before pressing [SOLVE]

If A, B, C were values already stored, scroll directly to D=, enter guess, scroll up, then [SOLVE]

If calculated D too small, below max(A/2, √C), press [SOLVE] again with higher guess.

This was vaguely (if followed the keystrokes) stated in the manual ... Oh well

Edit: although D can be solved with Cubic Solver, I worried about its bugs.
Calculated D must be very accurate, otherwise all 4 roots will be bad.

Without the bugs, Cubic Solver is easier. Just enter the coefficients:
a = -8, b = 4A, c = 8C, d = B^2-4AC, we get x1 = 1.99130646

[Albert Chan]

(11-29-2018 05:28 PM)Albert Chan Wrote:  Example: f(x) = x^4 + 2 x^3 + 3 x^2 + 4 x + 5 = 0

Depress the Quartic (remove cubic term) by substitute x = y + k, where k = -b/(4a)
We get depressed quartic: y^4 + A y^2 + B y + C

Using Casio CALC, evaluate f(x) at y = 0, ±1 ...

Another way to calculate shifted polynomial coefficients is repeated synthetic divisions.
Shift x-axis to k = -2/(4*1) = -0.5

k> 1 2 3 4 5
k> 1 1.5 2.25 2.875 3.5625 <-- C
k> 1 1.0 1.75 2 <-- B
k> 1 0.5 1.5 <-- A

Repeated differentiations also work:

C = f(k) = k^4 + 2 k^3 + 3 k^2 + 4 k + 5 = 3.5625
B = f'(k) = 4k^3 + 6k^2 + 6k + 4 = 2
A = f''(k)/2! = (12k^2 + 12k + 6) / 2! = 1.5

As expected: y^3 coef. = (24k + 12) / 3! = 0, y^4 coef. = (24) / 4! = 1

[Albert Chan]

(11-29-2018 05:28 PM)Albert Chan Wrote:  Find D such that RHS is perfect square, i.e. discriminant = B^2 - 4(2D - A)(D^2 - C) = 0
We do not want complex values for E, so 2D ≥ A, D^2 ≥ C:

D ≥ max(1.5/2, √3.5625) = max(0.75, 1.904) = 1.904

A trick to force Casio FX-115MS to solve for D: 0ABC + B^2 - 4(2D - A)(D^2 - C)

D = 1.99130646 (Casio SOLVE, with guessed D = 2)

For many calculators, [SOLVE] might not be available.
We can setup Newtons iterations to get D

2 =
Ans - (B² - 4(2 Ans - A)(Ans² - C)) ÷ (8(C + Ans(A - 3 Ans
= → 1.99137931
= → 1.991306465
= → 1.99130646

Note: for Casio FX-115MS, it is not necessary to put bracket before 8
However, without the bracket, formula is ambiguous: https://plus.maths.org/content/pemdas-paradox

[Albert Chan]

An old algorithm for solving Quartics: Brown's method
(source: W.S. Brown, "Solution of Biquadratic Equations", Aircraft Engineering, 16,14 (1944)

Quartics: x^4 + A x^3 + B x^2 + C x + D = 0

Step 1: get Y, the largest real root of resolvent cubics

y^3 - B y^2 + (AC-4D) y + D(4B-A²)-C² = 0

it will work for any root of cubics, but later quadratics might involve complex numbers.

Step 2: solve 2 quadratics

c^2 - A c + (B-Y) = 0
d^2 - Y d + D = 0

Step 3: make sure roots, c1 d2 + c2 d1 = C. If not, swap d1, d2.
This may required calculation of c1 d1 + c2 d2 as well, due to numerical imprecision.

I have a slight improvement in this step, see example below

Step 4: solve 2 more quadratics, the 4 roots is the roots of Quartics:

x^2 + c1 x + d1 = 0
x^2 + c2 x + d2 = 0

Example: solve x^4 + 2 x^3 + 3 x^2 + 4 x + 5 = 0

A = 2
B = 3
C = 4
D = 5
E = AC-4D = -12
F = D(4B-A²)-C² = 24

y^3 - B y^2 + E y + F = 0 → Y = 4.482612919

c^2 - A c + (B-Y) = 0        → [c1,c2] = [-0.575630959, 2.575630959]
d^2 - Y d + D = 0             → [d1,d2] = [2.088157364, 2.394455555]

if B and Y nearly equal, we may use another form for B-Y

Y^2 (Y-B) + E Y + F = 0 → B-Y = (E Y + F) / Y²

Calculate AY/2
If AY/2 > C, make sure roots of c,d sorted the same way
If AY/2 < C, make sure roots of c,d sorted opposite way

For this case, AY/2 = Y > C, thus c,d roots sorted the same way.

x^2 + c1 x + d1 = 0 → x = 0.2878154795 ± 1.416093081 i
x^2 + c2 x + d2 = 0 → x = -1.287815480 ± 0.8578967571 i

[Albert Chan]

(08-30-2019 02:02 PM)Albert Chan Wrote:  An old algorithm for solving Quartics: Brown's method ...

c^2 - A c + (B-Y) = 0
d^2 - Y d + D = 0

Step 3: make sure roots, c1 d2 + c2 d1 = C. If not, swap d1, d2 ...

Calculate AY/2
If AY/2 > C, make sure roots of c,d sorted the same way
If AY/2 < C, make sure roots of c,d sorted opposite way

To show AY/2 trick work, let roots, [c1,c2] = [p+q, p-q], [d1,d2] = [r+s,r-s]

C = c1 d2 + c2 d1
= (p+q)(r-s) + (p-q)(r+s)
= (pr - ps + qr - qs) + (pr + ps - qr - qs)
= 2pr - 2qs
= 2(A/2)(Y/2) - 2qs
= AY/2 - 2qs

If AY/2 = C, then 2qs = 0, sort order does not matter.
If AY/2 > C, then 2qs > 0, thus c,d roots sorted the same way.
If AY/2 < C, then 2qs < 0, thus c,d roots sorted opposite way.

[Albert Chan]

Found a Quartic Solver from https://thomasokken.com/free42/

It uses the same resolvent cubics as Brown's method, without sorting the roots.

(08-30-2019 02:02 PM)Albert Chan Wrote:  Example: solve x^4 + 2 x^3 + 3 x^2 + 4 x + 5 = 0

A = 2
B = 3
C = 4
D = 5
E = AC-4D = -12
F = D(4B-A²)-C² = 24

y^3 - B y^2 + E y + F = 0 → Y = 4.482612919

c^2 - A c + (B-Y) = 0        → [c1,c2] = [-0.575630959, 2.575630959]
d^2 - Y d + D = 0             → [d1,d2] = [2.088157364, 2.394455555]

if B and Y nearly equal, we may use another form for B-Y

Y^2 (Y-B) + E Y + F = 0 → B-Y = (E Y + F) / Y²

Calculate AY/2
If AY/2 > C, make sure roots of c,d sorted the same way
If AY/2 < C, make sure roots of c,d sorted opposite way

For this case, AY/2 = Y > C, thus c,d roots sorted the same way.

x^2 + c1 x + d1 = 0 → x = 0.2878154795 ± 1.416093081 i
x^2 + c2 x + d2 = 0 → x = -1.287815480 ± 0.8578967571 i

Upto calculations of Y is identical

E = √(A²/4 - B + Y) = 1.575630959
F = √(Y²/4 - D) = 0.153149095

→ [c1,c2] = A/2 ± E = [2.575630959, -0.575630959]
→ [d1,d2] = Y/2 ∓ F = [2.088157364, 2.394455555]

Unfortunately, the roots are in wrong order !

Solving (x^2 + c1 x + d1)(x^2 + c2 x + d2) = 0, gives all wrong roots.

Can anyone try the actual code, or am I mistaken somewhere ?
Source in https://thomasokken.com/free42/42progs/, named G4G

[Albert Chan]

(02-01-2020 12:27 PM)Albert Chan Wrote:  Unfortunately, the roots are in wrong order !

Solving (x^2 + c1 x + d1)(x^2 + c2 x + d2) = 0, gives all wrong roots.

Can anyone try the actual code, or am I mistaken somewhere ?
Source in https://thomasokken.com/free42/42progs/, named G4G

Tried above program just now, for Quartic x^4 + 2 x^3 + 3 x^2 + 4 x + 5 = 0

2 Enter 3 Enter 4 Enter 5 XEQ "G4G"

Output "X=I,Z=I" ⇒ roots = Y ± X*I, T ± Z*I ≈ -1.2878 ± -0.6555*I, 0.2878 ± -1.5204*I

Correct root should be -1.2878 ± 0.8579*I, 0.2878 ± 1.416*I

G4G code indeed gives the wrong roots

[SlideRule]

(08-24-2019 12:47 PM)Albert Chan Wrote:  Note: … without the bracket, formula is ambiguous: https://plus.maths.org/content/pemdas-paradox

Most excellent article! I deployed 'GEMA' in my technical math classes, placing an additional emphasis on 'equality' { = ≠ ≡ ≈ < ≤ > ≥ := etc. } notation.
The article parallels (for me anyway) the many aspects relating fallacy and ambiguity in language as documented in Fallacies Arising from Ambiguity (APPLIED LOGIC SERIES), DOUGLAS WALTON, Department of Philosophy, University of Winnipeg, Canada, © 1996, ISBN 978-94-015-8632-0 (eBook).

SlideRule

[Albert Chan]

(08-30-2019 02:02 PM)Albert Chan Wrote:  Example: solve x^4 + 2 x^3 + 3 x^2 + 4 x + 5 = 0

A = 2
B = 3
C = 4
D = 5
E = AC-4D = -12
F = D(4B-A²)-C² = 24

y^3 - B y^2 + E y + F = 0 → Y = 4.482612919

c^2 - A c + (B-Y) = 0        → [c1,c2] = [-0.575630959, 2.575630959]
d^2 - Y d + D = 0             → [d1,d2] = [2.088157364, 2.394455555]

Calculate AY/2
If AY/2 > C, make sure roots of c,d sorted the same way
If AY/2 < C, make sure roots of c,d sorted opposite way

For this case, AY/2 = Y > C, thus c,d roots sorted the same way.

x^2 + c1 x + d1 = 0 → x = 0.2878154795 ± 1.416093081 i
x^2 + c2 x + d2 = 0 → x = -1.287815480 ± 0.8578967571 i

We solve quadratic 4 times, 1 for c's, 1 for d's, 2 more for 4 roots.
We can combine them to 2, although this does not save anything.

Note: sqrt always give non-negative real part, sorting only mean 1 swap, c's or d's.
Below trick swap c if needed, see (32S) Polynomial Root Finder from 32SII manual

c1, c2 = A/2 ± sign(A*Y/2-C) * √((A/2)²-B+Y)
d1, d2 = Y/2 ± √((Y/2)²-D)

x^4 + A*x^3 + B*x^2 + C*x + D = (x^2 + c1*x + d1) * (x^2 + c2*x + d2)

My preference is solve for c's and d's first, then swap c if needed.

This removed the need to have A,B,C,D real
This removed the need to have Y be largest real. (no need to solve all cubic roots!)
We check for C == c1*d2 + c2*d1, swap c's if not true, then solve the 2 quadratics.