Mathematician Finds Easier Way to Solve Quadratic Equations

[Namir]

Ok folks .... any one can talk about fractional derivatives ... would like to see us move a bit forward in math topics!!!!

Math is fun .... but trivial math can be a bit boring!!!

Namir

[Albert Chan]

"Easier" way is really for educationally purpose, instead of memorizing quadratic formula.

Let roots = m ± g (all parabola look like this)

(x - (m+g)) * (x - (m-g)) = x² + (-2m) x + (m²-g²)

Match this to x² + (b/a) x + (c/a) = 0

b/a = -2m
c/a = m²-g²

LHS is just constants. Solve for m, then for g

We can use same pattern trick for cubic roots. Solve for y, then z = p/(-3y)

(07-20-2021 10:12 PM)Albert Chan Wrote:  x³ + y³ + z³ − 3xyz = (x+y+z) * (x+yω+z/ω) * (x+y/ω+zω), where ω = e^(i*2*pi/3)

LHS is a depressed cubic: x³ + p*x + q , where p = -3yz, q = y³ + z³

We then setup a quadratic of t, with roots (y³, z³)

(t - y³)*(t - z³) = t² - q*t - (p/3)³       → t = q/2 ± √Δ, where Δ = (q/2)² + (p/3)³

https://math.stackexchange.com/questions...2ab-omega2

[Gerson W. Barbosa]

(04-24-2024 12:20 AM)Thomas Klemm Wrote:  

The lyrics are available in the video description. I’ll place the text in a code box for space reasons. Anyway, human languages are technically codes. Google translate is not perfect, but those of us non-German speakers can have a good idea of what is being sung. I know this hasn’t been addressed to all, but thanks for posting!

Code:


Liedtext zu pq-Formel (Mathe-Song)

Hast du schon mal versucht, eine Gleichung mit Quadraten drin zu lösen?
Wenn ja, dann weißt du sicher, dabei darf man nicht dösen,
denn, ob es eine Lösung, keine Lösung, zwei Lösungen gibt,
weiß man oft erst, wenn sich das dann durch eine Rechnung ergibt,
doch das ist gar nicht mal so schwer, spitz deine Ohren, hör' jetzt gut zu:
Du machst da einfach 0=x^2+px+q
draus und schon hast du damit fast die ganze Arbeit getan,
denn für den Rest eignest du dir nur noch die Lösungsformel an:

x = -p/2 ± √((p/2)^2 - q)

Du fragst dich vielleicht: Wie kommt man da drauf und warum geht das immer auf?
Nun ja, ich zeig's dir und schreib dir mal die Ausgangsgleichung auf.
Wir rechnen jetzt erst +(p/2)^2 und dann noch -q,
was das bringt, erklär ich dir jetzt im Nu:
Nimm die erste Binomische Formel herbei,
ersetze a mit x und b mit p/2
und den nächsten Schritt, den machst du jetzt,
indem du einfach nur das hier mit dem da ersetzt.

Jetzt ist das x schon mal vereinzelt, was noch stört ist das Quadrat,
doch nur das Wurzelziehen hat die Lösung auch noch nicht parat,
denn auch aus „Minus mal Minus" wird dann ja letztlich noch ein Plus
und deshalb ist an dieser Stelle hier ein Plusminus ein Muss.
Den allerletzten Schritt sollte jeder gleich verstehen
und normalerweise auch auf den ersten Blick sehen,
denn jetzt wird bei der Gleichung noch p/2 subtrahiert
und wir sind fertig, denn jetzt ist das x isoliert.

x = -p/2 ± √((p/2)^2 - q)

Since we’re at it I decided to take a look at the work you did fourteen years ago (last messages in this old thread) in order to try to fully understand it.

Given the second degree equation

\(x^2+px+q=0\)

we can find the first root, as we have seen, with help of a this very compact formula:

\(x_{1}=-\frac{p}{2}+\sqrt{\left({-\frac{p}{2}}\right)^2-q}\)

This is suitable for minimizing the number of programming steps, but like the traditional formula it is subject to cancellation errors in some extreme cases. Hence your motivation for another method which avoids this problem.
I’ll try to reproduce your steps to arrive at the formula you presented then.

Multiply and divide the expression inside the surd by q and take the multiplying factor in the numerator out of it:

\(x_{1}=-\frac{p}{2}+\sqrt{q}\sqrt{\frac{\left({-\frac{p}{2}}\right)^2-q}{q}}\)

Multiply and divide the first term by √q:

\(x_{1}=\sqrt{q}\left({-\frac{p}{2\sqrt{q}}}\right)+\sqrt{q}\sqrt{\left({-\frac{p}{2\sqrt{q}}}\right)^2-1}\)

\(x_{1}=\sqrt{q}\left({-\frac{p}{2\sqrt{q}}+\sqrt{\left({-\frac{p}{2\sqrt{q}}}\right)^2-1}}\right)\)

Make

\(z=-\frac{p}{2\sqrt{q}}\)

Use the identity

\(\cosh^{-1}{(z)}=\ln\left({z+\sqrt{z^2-1}}\right)\)

or

\(e^{\cosh^{-1}{(z)}}=z+\sqrt{z^2-1}\)

Then

\(x_{1}=\sqrt{q}\times e^{\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

The second root can be computed as

\(x_{1} x_{2}= q\)

\(x_{2}=\frac{q}{x_{1}}\)

Done totally by hand. Hopefully no mistakes or typos.

HP-42 program

Code:


00 { 26-Byte Prgm }
01▸LBL "Q"
02 X<> ST Z
03 STO÷ ST Z
04 STO+ ST X
05 ÷
06 +/-
07 RCL ST Y
08 SQRT
09 X<>Y
10 RCL÷ ST Y
11 ACOSH
12 E↑X
13 ×
14 STO÷ ST Y
15 END

Example:

\(10^{-13}x^2-2x+1=0\)

E +/- 13 ENTER 2 +/- ENTER XEQ Q ->

T: previous X
Z: previous X
Y: 0.50000000001
X: 1.9999999996E13

Best regards,

Gerson.

[Albert Chan]

(05-01-2024 04:57 PM)Gerson W. Barbosa Wrote:  Given the second degree equation

\(x^2+px+q=0\)

... Then

\(x_{1}=\sqrt{q}\times e^{\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

The second root can be computed as

\(x_{1} x_{2}= q\)

\(x_{2}=\frac{q}{x_{1}}\)

Nice! And, we can make formula more compact, 2 roots with ±

(08-01-2021 03:08 PM)Albert Chan Wrote:  We can also solve for z² - 2*m*z + n² = 0

XCas> m,n := (-18+5i)/-2., sqrt(45.-15i)
XCas> proot([1, -2m, n*n])                                    → [3.0, 15.0-5.0*i]
XCas> n*tan(asin(n/m)/2), n/tan(asin(n/m)/2)        → [3.0, 15.0-5.0*i]
XCas> n/exp(acosh(m/n)), n*exp(acosh(m/n))        → [3.0, 15.0-5.0*i]

Again, combine both roots, z = n*exp(±acosh(m/n))

tan/asin quadratic formula based from tan half-angle formula:

sin(x) = 2t/(1+t^2), where t = tan(x/2)
t^2 - 2*(1/sin(x))*t + 1 = 0

exp/acosh quadratic formula based from cosh definition.

cosh(x) = (y + 1/y)/2, where y = e^x
y^2 - 2*cosh(x)*y + 1 = 0

Above patterns matched with scaled z² - 2*m*z + n² = 0

(z/n)^2 - 2*(m/n)*(z/n) + 1 = 0

[Thomas Klemm]

(05-01-2024 04:57 PM)Gerson W. Barbosa Wrote:  Since we’re at it I decided to take a look at the work you did fourteen years ago (last messages in this old thread) in order to try to fully understand it.

My first post in this forum was this reply: Re: Short quadratic solver (HP-15C)

(05-01-2024 04:57 PM)Gerson W. Barbosa Wrote:  (…) but like the traditional formula it is subject to cancellation errors in some extreme cases. Hence your motivation for another method which avoids this problem.

Palmer's program uses a different approach: A "Cadillac" Quadratic Solver for the hp 33s and HP35s
My article attempts to explain how it works: Solve the real quadratic equation \(c-2bz+az^2=0\) for real or complex roots.

[Gerson W. Barbosa]

(05-01-2024 06:51 PM)Albert Chan Wrote:  

(05-01-2024 04:57 PM)Gerson W. Barbosa Wrote:  Given the second degree equation

\(x^2+px+q=0\)

... Then

\(x_{1}=\sqrt{q}\times e^{\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

The second root can be computed as

\(x_{1} x_{2}= q\)

\(x_{2}=\frac{q}{x_{1}}\)

Nice! And, we can make formula more compact, 2 roots with ±

Yes,

\(x_{1,2}=\sqrt{q}\times e^{±\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

but I am mainly interested in program optimization. Anyway we can save one step with it, but at the cost of an extra couple of bytes:

Code:


00 { 28-Byte Prgm }
01▸LBL "Q"
02 RCL÷ ST Z
03 SQRT
04 X<> ST Z
05 +/-
06 STO+ ST X
07 ÷
08 RCL÷ ST Y
09 ACOSH
10 E↑X
11 RCL× ST Y
12 X<>Y
13 RCL÷ ST L
14 END

MEM says I have 1647 bytes left on the HP-42S and 433569792 on Free42. I don’t think memory should be a concern.

[Albert Chan]

(05-01-2024 11:08 PM)Gerson W. Barbosa Wrote:  \(x_{1,2}=\sqrt{q}\times e^{±\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

exp/cosh quadratic formula avoided subtraction canellation, but introduced another.
Exponential function slope is also an exponential function!

e^(x*(1-ε)) ≈ e^x - (x*e^x)*ε

RelErr(x) = ε      ⇒ RelErr(e^x) = (exact - approx) / exact ≈ (x*e^x)*ε / e^x = x*ε

Fortunately, acosh is similar to log, keeping size of x in check.
Still, Trig (tan/asin) quadratic version is perhaps more accurate.

Code:

00 { 24-Byte Prgm }
01▸LBL "T"          ; c     b   a
02 RCL÷ ST Z        ; q     b   a
03 SQRT             
04 X<> ST Z         ; a     b   √q
05 +/-              
06 ENTER            
07 +                ; -2a   b   √q   √q
08 ÷                
09 ÷                ; -2√q/p    √q   √q
10 ASIN
11 2
12 ÷
13 TAN
14 STO× ST Z
15 ÷
16 END

1E-13 Enter -2 Enter 1 XEQ "Q"
Y : 19,999,999,999,999.49999999999998749992      // error = 8 ULP
X : 5.000000000000125000000000006250021e-1     // error = -21 ULP

1E-13 Enter -2 Enter 1 XEQ "T"
Y: 5.000000000000125000000000006250001e-1       // error = -1 ULP
X: 19,999,999,999,999.4999999999999875              // error = 0 ULP

There is another issue. acosh(x) has real domain x≥1
Trig quadratic version does not have this problem. (tan/asin both odd functions)

1E-13 Enter 2 Enter 1 XEQ "Q"
Y: -19,999,999,999,999.49999999999998749992-2.31605661201243998216520079425921e-21i
X: -0.500000000000012500000000000625002+5.790141530031389462489503558740034e-35i

1E-13 Enter 2 Enter 1 XEQ "T"
Y: -5.000000000000125000000000006250001e-1
X: -19,999,999,999,999.4999999999999875

Update: here is e^asinh version. Determinant Δ = (p/2)² - q
It has the same RelErr(e^x) issue ... not recommended.

\(
x_{1,2}=
\sqrt{q}× e^{±\sinh^{-1} \sqrt{\frac{Δ}{q}} }
\)

[Gerson W. Barbosa]

(05-02-2024 02:00 PM)Albert Chan Wrote:  

(05-01-2024 11:08 PM)Gerson W. Barbosa Wrote:  \(x_{1,2}=\sqrt{q}\times e^{±\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

exp/cosh quadratic formula avoided subtraction canellation, but introduced another.
Exponential function slope is also an exponential function!

e^(x*(1-ε)) ≈ e^x - (x*e^x)*ε

RelErr(x) = ε      ⇒ RelErr(e^x) = (exact - approx) / exact ≈ (x*e^x)*ε / e^x = x*ε

Fortunately, acosh is similar to log, keeping size of x in check.
Still, Trig (tan/asin) quadratic version is perhaps more accurate.

Code:

00 { 24-Byte Prgm }
01▸LBL "T"          ; c     b   a
02 RCL÷ ST Z        ; q     b   a
03 SQRT             
04 X<> ST Z         ; a     b   √q
05 +/-              
06 ENTER            
07 +                ; -2a   b   √q   √q me
08 ÷                
09 ÷                ; -2√q/p    √q   √q
10 ASIN
11 2
12 ÷
13 TAN
14 STO× ST Z
15 ÷
16 END

Very nice! Thanks for this alternative method, useful for both Free42 and the HP-42S.

The following is longer, but preserves the previous stack register X for later use. Some optimization is still possible, I think.

Code:


00 { 31-Byte Prgm }
01▸LBL "T"
02 X<> ST Z
03 STO÷ ST Z
04 +/-
05 STO+ ST X
06 ÷
07 X<>Y
08 SQRT
09 STO ST Z
10 X<>Y
11 ÷
12 ASIN
13 2
14 ÷
15 TAN
16 RCL× ST Y
17 X<>Y
18 RCL÷ ST L
19 END

On the real HP-42S, in DEG mode:

E 13 +/- ENTER 2 +/- ENTER 1 XEQ T ->

Y: 4.99999999999E-1
Y: 2.00000000000E13

Same answers in RAD mode. RAD mode appears to be slightly more accurate. No reason to change the angular mode just for that reason, I think. (Mine is most of the time in DEG mode).

This is my very first calculator program, written for the SHARP EL-5813 back in 1981:

2ndF LRN ( ( ( K₂ +/- + ( K₂ x² - 4 K₁ K₃ ) 2ndF √ ) ÷ 2 K₁ ) 2ndF LOOK + K₂ ÷ K₁ ) +/- 2ndF LRN

I don’t have it anymore, but I remember the usage was

a STO K₁
b STO K₂
c STO K₃
COMP

Quadratic equations were often needed to be solved during my first college years, Physics then EE, so I always kept a quadratic solver in the programmable calculator I was using at the moment (TI-51-III, TI-59, HP-15C, HP-28S and HP-48GX). No Cadillac solvers needed for those trivial equations, though. My Beetle solvers did just fine! :-)

Edited to fix a typo and an orthographic error.

[Namir]

(05-02-2024 02:00 PM)Albert Chan Wrote:  

(05-01-2024 11:08 PM)Gerson W. Barbosa Wrote:  \(x_{1,2}=\sqrt{q}\times e^{±\cosh^{-1}{(-\frac{p}{2\sqrt{q}}})}\)

exp/cosh quadratic formula avoided subtraction canellation, but introduced another.
Exponential function slope is also an exponential function!

e^(x*(1-ε)) ≈ e^x - (x*e^x)*ε

RelErr(x) = ε      ⇒ RelErr(e^x) = (exact - approx) / exact ≈ (x*e^x)*ε / e^x = x*ε

Interesting formula but it uses more CPU power to calculate invers cosh and exp(x) in addition to the square root!

Namir

[Albert Chan]

Trivia: with signed zero, e^acosh(z) is an odd function!

acosh(z) = acos(z) * ±i       // sign to ensure acosh non-negative real part.
acos(z) + acos(-z) = pi       // acos has non-negative real part too: 0 .. pi

acosh(z) * ∓i + acosh(-z) * ±i = pi      // im(acos(z)) = - im(acos(-z))
acosh(z) - acosh(-z) = ±pi*i                // multiply both side by ±i
e^acosh(z) / e^acosh(-z) = -1            // Euler Identity, e^(pi*i) = -1

e^acosh(z) = - e^acosh(-z)

Identity might still hold without signed zero, but I have not checked edge cases.

---

Another proof for acosh(z) - acosh(-z) = ±pi*i

versin(x) = 1 - cos(x) = 2*sin(x/2)^2      → acos(1-x) = 2*asin(sqrt(x/2))

acos(z) = 2*asin(sqrt((1-z)/2)) = 2*atan(sqrt((1-z)/(1+z)))
acosh(z) = 2*asinh(sqrt((z-1)/2)) = 2*atanh(sqrt((z-1)/(z+1)))

acosh(z) - acosh(-z)
= 2*atanh(u) - 2*atanh(1/u))      , where u = sqrt((z-1)/(z+1))
= ln(v) - ln(-v)                            , where v = (1+u)/(1-u)
= ±pi*i

[Albert Chan]

(05-01-2024 04:57 PM)Gerson W. Barbosa Wrote:  Use the identity

\(\cosh^{-1}{(z)}=\ln\left({z+\sqrt{z^2-1}}\right)\)

or

\(e^{\cosh^{-1}{(z)}}=z+\sqrt{z^2-1}\)

Above identity does not hold, if sign(Re(z)) = -1

acosh(z) = ln(z + √(z+1) * √(z-1))      // 2 square roots to ensure acosh non-negative real part

Example: z = -2

ln(-2 + √(-1) * √(-3)) = ln(-2 − √3) ≈ +1.31696 + pi*i
ln(-2 + √((-2)^2 - 1)) = ln(-2 + √3) ≈ −1.31696 + pi*i

For acosh √(z²-1) form, we need to match sign, to ensure non-negative real part.
Again, assume sign(x) = ±1, never 0

\(e^{\cosh^{-1}{(z)}}=z + sign(\Re(z)) \;\sqrt{z^2-1}\)

It is now trivial to show e^acosh(z) is an odd function:

e^acosh(z) = - ((-z) + sign(Re(-z)) * √((-z)²-1)) = - e^acosh(-z)

Update: Perhaps better proof

e^z = (e^z+1/e^z)/2 + (e^z-1/e^z)/2 = cosh(z) + sinh(z)

e^acosh(z) = cosh(acosh(z)) + sinh(acosh(z)) = z + sign(Re(z)) * √(z²-1)

[Gerson W. Barbosa]

(05-03-2024 02:02 PM)Gerson W. Barbosa Wrote:  The following is longer, but preserves the previous stack register X for later use. Some optimization is still possible, I think.

Code:


00 { 31-Byte Prgm }
01▸LBL "T"
02 X<> ST Z
03 STO÷ ST Z
04 +/-
05 STO+ ST X
06 ÷
07 X<>Y
08 SQRT
09 STO ST Z
10 X<>Y
11 ÷
12 ASIN
13 2
14 ÷
15 TAN
16 RCL× ST Y
17 X<>Y
18 RCL÷ ST L
19 END

On the real HP-42S, in DEG mode:

E 13 +/- ENTER 2 +/- ENTER 1 XEQ T ->

Y: 4.99999999999E-1
Y: 2.00000000000E13

This saves four bytes and four steps:

Code:


00 { 27-Byte Prgm }
01▸LBL "T"
02 X<> ST Z
03 STO÷ ST Z
04 STO+ ST X
05 +/-
06 ÷
07 RCL ST Y
08 SQRT
09 RCL÷ ST Y
10 ASIN
11 COS
12 RCL× ST Y
13 +
14 STO÷ ST Y
15 END

On the real HP-42S, in DEG mode:

E 13 +/- ENTER 2 +/- ENTER 1 XEQ T ->

Y: 0.5
Y: 2.E13

On Free42, in DEG mode:

E 13 +/- ENTER 2 +/- ENTER 1 XEQ T ->

Y: 0.500000000000012500000000000625
X: 19,999,999,999,999.4999999999999875 

\(x_{1}=\left({\cos\left({\arcsin{\frac{\sqrt{q}}{-\frac{p}{2}}}}\right)+1}\right)\times-\frac{p}{2}\)

That is equivalent to Albert Chan’s formula.

\(x_{2}=\frac{q}{x_{1}}\)

[Albert Chan]

(05-04-2024 07:37 PM)Gerson W. Barbosa Wrote:  \(x_{1}=\left({\cos\left({\arcsin{\frac{\sqrt{q}}{-\frac{p}{2}}}}\right)+1}\right)\times-\frac{p}{2}\)

Going for x²+p*x+q=0 big root ... I like it!

cos(asin(z)) = √(1-z²)      // both sides have non-negative real part

x1 = (√(1 - q/(-p/2)²) + 1) * -p/2
x2 = q/x1

Another way to see this:

x^2 + p*x + q = 0            // x = (-p/2)*y
y^2 - 2y + sin(θ)^2 = 0    // sin(θ)^2 = q/(-p/2)^2
(y-1)^2 = cos(θ)^2
y = 1 ± cos(θ)

[Thomas Klemm]

The following geometric approach uses the right triangle altitude theorem.



Given \(c = p + q\) and \(h^2 = p \cdot q\) we get:

\(
\begin{align}
0
&=(x - p)(x - q) \\
&= x^2 - (p+q)x + p \cdot q \\
&= x^2 - c x + h^2 \\
\end{align}
\)

But then also:

\(
\begin{align}
c &= 2r \\
r^2 &= h^2 + d^2
\end{align}
\)

And therefore:

\(
\begin{align}
r &= \frac{c}{2} \\
d &= \sqrt{r^2 - h^2} \\
\\
p &= r + d \\
q &= r - d \\
\end{align}
\)

Example

\(
\begin{align}
x^2 - 5x + 6 = 0 \\
\end{align}
\)

\(
\begin{align}
c &= 5 \\
h^2 &= 6 \\
\\
r &= 2.5 \\
d
&= \sqrt{6.25 - 6} \\
&= \sqrt{0.25} \\
&= 0.5 \\
\\
p &= 2.5 + 0.5 = 3 \\
q &= 2.5 - 0.5 = 2 \\
\end{align}
\)

[Namir]

Let's keep in mind the KISS principle (keep it simple stupid), or what Einstein said, "Make it simple, but no simpler". Using fancy functions like cosh, exp, sin, and so on is fancy, but introduces a bit of complicatrion in my very humble opinion.

Wat's next? Solving systems of linear equations using Cramer's rule???

Namir

[carey]

(05-05-2024 02:24 AM)Namir Wrote:  Wat's next? Solving systems of linear equations using Cramer's rule???

Namir

And what’s wrong with using Cramer’s rule for linear systems? :)

[Thomas Klemm]

(05-05-2024 02:24 AM)Namir Wrote:  Wat's next? Solving systems of linear equations using Cramer's rule???

How about using the HP-67 to evaluate the determinant of a square matrix by repeated applications of the process of pivotal condensation?

[Maximilian Hohmann]

Hello!

(05-05-2024 02:24 AM)Namir Wrote:  Let's keep in mind the KISS principle (keep it simple stupid), or what Einstein said, "Make it simple, but no simpler".

Yes, it's amazing. A simple statement "Mathematician Finds Easier Way to Solve Quadratic Equations" has turned into "how many HP calculator nerds does it take to change a lightbulb?" ;-)

But apart from that I keep wondering why the quadratic equation and it's memorized solution is attributed such a high significance by schools. I have never since finishing school come across a quadratic equation in real life. Things are mostly simple and can be handled with linear calculations ("if my plane needs 1400 pounds of fuel in one hour, how much must I take along for a 2:25h flight?") or they are really complicated and require higher order polynomials, transcendental funcions or numeric solutions. But a quadratic equation? Never ever.

Regards
Max

[Namir]

(05-05-2024 06:31 AM)carey Wrote:  

(05-05-2024 02:24 AM)Namir Wrote:  Wat's next? Solving systems of linear equations using Cramer's rule???

Namir

And what’s wrong with using Cramer’s rule for linear systems?

In the late 80s I attended a programmers' conference in Boston. A speaker mentioned that he encountered a math application that was running WAY TOO SLOW! When he checked the source code, he found that the code author was using Cranmer's method to solve linear equations. You should have heard the attendees reaction .. a resounding ouuuuuuuh!!!

Just sharing!

[Namir]

(05-05-2024 11:07 AM)Thomas Klemm Wrote:  

(05-05-2024 02:24 AM)Namir Wrote:  Wat's next? Solving systems of linear equations using Cramer's rule???

How about using the HP-67 to evaluate the determinant of a square matrix by repeated applications of the process of pivotal condensation?

I don't see Cranmer's rule mentioned much in Numerical Analysis books when they discuss the topic of solving linear equations.`I guess these authors maintain a reasonable bar for algorithms used to solve linear equations.

Namir