Saturday, November 7, 2020

HP Prime: Approximate of a Quartic Root

 HP Prime: Approximate of a Quartic Root


Introduction


In the article "Square Root & Cube Root Algorithms" (see source below), teacher and educator Dave  Elgin wrote an article on how his Advanced Higher Applied Mathematics class developed algorithms to estimate the square root and cubic roots of numbers.  The algorithms are first based off of Netwon's Method, then uses selected initial guesses and solving for linear systems.


Elgin's Derivation 


Square Root


Suppose x^2 = k.   To set up for Newton's Method, let f(x) = x^2 - k, with df/dx = 2x.  Then


x_n+1 = x_n - (x_n^2 - k) / (2 * x_n) = 1/2 * (x_n + k / x_n^2)


(The article uses N for the number to find the root of, but to eliminate confusion, I use k.  Completely a choice of labels.)


The class used a pattern to determine guesses for the square root of k:


k / (x_n),  (k^2) / (x_n^3),  (k^3) / (x_n^5), and so on.   


Let g(x) be an iterative function where x_n+1 = g(x_n) and for a second-order approximation:


g(x) = a1 * x + a2 * k / x + a3 * k^2 / x^3


Two derivatives of g(x) are taken:


g'(x) = a1 - a2 * k / x^2 - 3 * a3 * k^2 / x^4, with g'(x_n) = 0


g''(x) = 2 * a2 * k / x^3 + 12 * a3 * k^2 / x^5 with g''(x_n) = 0


This requires the following system of equations to be solved for a1, a2, and a3 (after substituting x = √k):


1 = a1 + a2 + a3

0 = a1 - a2 - 3 * a3

0 = 2 * a2 + 12 * a3


leading to the solutions a1 = 3/8, a2 = 3/4, and a3 = -1/8


Hence the second-order recursive function for the square root is (after simplifying):


g(x) = 3/8 * x + (3 * k) / (4 * x) -  k^2 / (8 * x^3)


which translates to 


x_n+1 = 3/8 * x_n + (3 * k) / (4 * x_n) -  k^2 / (8 * x_n^3)


Cube Root


They repeat the same process for the cubic root, which I will briefly outline here:


x^3 = k,   f(x) = x^3 - k,  f'(x) = 3x^2


x_n+1 = x_n - (x^3 - k) / (3x^2) = 1/3 * (2 * x_n - k / (x_n^2))


With guess of k / (x_n^2) and k^2 / (x_n^5) used, the iterative function is set up as:


g(x) = a1 * x + a2 * k / x^2 + a3 * k^2 / x^5


and


g'(x) = a1 - 2 * a2 * k / x^3 - 5 * a3 * k^2 / x^6


g''(x) = 6 * a2 * k / x^4 + 30 * a3 * k^2 / x^7


and with g(x_n) = x_n+1, g'(x_n) = 0, g''(x_n) = 0 and substituting x = k^1/3, the system becomes:


1 = a1 + a2 + a3

0 = a1 - 2 * a2 - 5 * a3

0 = 6 * a2 + 30 * a3


with the solutions a1 = 5/9, a2 = 5/9, and a3 = -1/9, giving the second-order recursive function:


g(x) = 5/9 * x + (5 * k) / (9 * x^2) - k^2 / (9 * x^5)


The article shows the derivation of a third-order recursive function for both square and cube root. 


Deriving a Second-Order Algorithm for Quartic Roots


Let's use a similar approach used in Elgin's article to develop an algorithm to calculate the fourth (quartic) root:  


k^1/4 = x


Let f(x) = x^4 - k,  then f'(x) = 4*x^3, and


x_n+1 = x_n - (x_n^4 - k) / (4 * x_n^3) = 3/4 * x_n - k / (4 * x_n^3)


Use guesses x_n, k / (x_n^3), k^2 / (x_n^7), we set up the equations:


g(x) = a1 * x + a2 * k / x^3 + a3 * k^2 / x^7


g'(x) = a1 - 3 * a2 * k / x^4 - 7 * a3 * k^2 / x^8


g''(x) = 12 * a2 * k / x^5 + 56 * a3 * k^2 / x^9



Setting g(x) = x^1/4, g'(x) = 0, g''(x) = 0, and setting g(k^1/4), we get the system:


1 = a1 + a2 + a3

0 = a1 - 3 * a2 - 7 * a3

0 = 12 * a2 + 56 * a3


The solutions to above system:  a1 = 21/32, a2 = 7/16, a3 = -3/32, which gives the second order recursive  equation:


g(x) = 21/32 * x + (7 * k) / (32 * x^3) - (3 * k) / (32 * x^7)

 

The program FTHROOT use the recursive equation to approximate the quartic root. 



HP Prime Program:  FTHROOT


EXPORT FTHROOT(k)

BEGIN

// EWS 2020-10-21

// Approx 4th Root

LOCAL r,r0,r1,ri;

r:=k^0.25;

r0:=0; 

r1:=√k;

ri:=0;

WHILE ABS(r0-r1)>1ᴇ−10 DO

ri:=ri+1;

r0:=r1;

r1:=(21*r0)/32+(7*k)/(16*r0^3)-(3*k^2)/(32*r0^7);

END;

PRINT();

PRINT("4√"+PRINT(k));

PRINT("Root = "+STRING(r));

PRINT("------");

PRINT("Approximation: "+

STRING(r1));

PRINT("Iterations: "+STRING(ri));


END;


The choice of a good first guess is necessary with any iterative root finding process. The program FTHROOT chooses the square root of k for an initial guess.   The goal is to seek a positive root.


Examples


Each example is followed by a set of screen shots, which include setting up Sequences and their graphs on the HP Prime.  


Example 1


k = 176.4

Result:  3.64438831256 (algorithm took 7 iterations with initial guess √176.4)



Example 2

k = 5525
Result:  8.62150472576 (algorithm took 9 iterations with initial guess √5525)




Source

Elgin, Dave.  "Square Root & Cube Root Algorithms"  The Mathematical Association.  Mathematics in School, Jan. 2006, Vol. 35, No. 1 pp. 30-31.  https://www.jstor.org/stable/30215863

Eddie

All original content copyright, © 2011-2020.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author.