Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]

Groups > comp.programming > #313 > unrolled thread

Fourth root of integer

Back to article view | Back to comp.programming

Contents

  Fourth root of integer Old Wolf <oldwolf@inspire.net.nz> - 2011-05-15 20:30 -0700
    Re: Fourth root of integer kym@kymhorsell.com - 2011-05-16 04:10 +0000
      Re: Fourth root of integer Jacob <asdf@hsu.edu> - 2011-05-16 02:49 -0400
    Re: Fourth root of integer pete <pfiland@mindspring.com> - 2011-05-16 09:46 -0400
      Re: Fourth root of integer Charles Richmond <frizzle@tx.rr.com> - 2011-05-16 22:39 -0500
        Re: Fourth root of integer pete <pfiland@mindspring.com> - 2011-05-17 09:17 -0400
    Re: Fourth root of integer rossum <rossum48@coldmail.com> - 2011-05-16 16:41 +0100
    Re: Fourth root of integer Ralph Boland <rpboland@gmail.com> - 2011-05-16 11:06 -0700
    Re: Fourth root of integer Dann Corbit <dcorbit@connx.com> - 2011-05-16 13:28 -0700
      Re: Fourth root of integer Christopher Henrich <chenrich@monmouth.com> - 2011-05-18 23:41 -0400
        Re: Fourth root of integer rossum <rossum48@coldmail.com> - 2011-05-19 09:54 +0100
          Re: Fourth root of integer Christopher Henrich <chenrich@monmouth.com> - 2011-05-20 15:16 -0400
        Re: Fourth root of integer - fourth-root.mws [1/1] Dann Corbit <dcorbit@connx.com> - 2011-05-20 12:42 -0700
          Re: Fourth root of integer - fourth-root.mws [1/1] Christopher Henrich <chenrich@monmouth.com> - 2011-05-20 16:10 -0400
      Re: Fourth root of integer James Waldby <not@valid.invalid> - 2011-05-19 07:01 +0000
        Re: Fourth root of integer James Waldby <not@valid.invalid> - 2011-05-21 00:37 +0000

#313 — Fourth root of integer

Looking for an algorithm to find the fourth root of
a large integer (I have a multi-precision integer
library already). If there is no exact fourth root then
the floor or ceiling is OK, I don't need a 'remainder'
or anything.

Any ideas?

[toc] | [next] | [standalone]

#314

Old Wolf <oldwolf@inspire.net.nz> wrote:
> Looking for an algorithm to find the fourth root of
> a large integer (I have a multi-precision integer
> library already). If there is no exact fourth root then
> the floor or ceiling is OK, I don't need a 'remainder'
> or anything.
> Any ideas?

NR? 
Binary search?

-- 
[Help, my pattern matcher has developed socialistic tendencies!]
Mr. Robot is just another left wing alarmist.    
  -- Ed <ed.carpenter@ymail.com>, 19 Feb 2011 07:48 -0800 (PST)

[toc] | [prev] | [next] | [standalone]

#315

On 16 May 2011 04:10:58 GMT, kym@kymhorsell.com wrote:

>Old Wolf <oldwolf@inspire.net.nz> wrote:
>> Looking for an algorithm to find the fourth root of
>> a large integer (I have a multi-precision integer
>> library already). If there is no exact fourth root then
>> the floor or ceiling is OK, I don't need a 'remainder'
>> or anything.
>> Any ideas?
>
>NR? 
>Binary search?

Does the library have a square root function?

If the number in question doesn't have an exact 4th root, it doesn't
have an exact square root, either.

Will this work?
Take successive floors of square roots. Then, if needed, raise the
second square root and increments thereof to the 4th power to bracket
the original number. Could make incrementing incrementing less
expensive by expanding (n+1)^4, as n^4 would need to be calculated
directly only once.

[toc] | [prev] | [next] | [standalone]

#316

Old Wolf wrote:
> 
> Looking for an algorithm to find the fourth root of
> a large integer (I have a multi-precision integer
> library already). If there is no exact fourth root then
> the floor or ceiling is OK, I don't need a 'remainder'
> or anything.
> 
> Any ideas?

Why not just the square root of the square root?

unsigned fourth_root(unsigned n)
{
    return isqrt(isqrt(n));
}

unsigned isqrt(unsigned n)
{
    unsigned power_of_4, root;

    power_of_4 = ((unsigned)-1 >> 2 - (unsigned)-1 % 3) + 1;
    root = 0;
    do {
        if (n >= root + power_of_4) {
            n -= root + power_of_4;
            root += power_of_4 << 1;
        }
        root >>= 1;
        power_of_4 >>= 2;
    } while (power_of_4 != 0);
    return root;
}

-- 
pete

[toc] | [prev] | [next] | [standalone]

#320

On 5/16/11 8:46 AM, pete wrote:
> Old Wolf wrote:
>>
>> Looking for an algorithm to find the fourth root of
>> a large integer (I have a multi-precision integer
>> library already). If there is no exact fourth root then
>> the floor or ceiling is OK, I don't need a 'remainder'
>> or anything.
>>
>> Any ideas?
>
> Why not just the square root of the square root?
>

The fourth root of (2143/22) gives PI to about 8 decimal places.
(Using your hand calculator, take the square root twice.)


-- 
+----------------------------------------+
|     Charles and Francis Richmond       |
|                                        |
|  plano dot net at aquaporin4 dot com   |
+----------------------------------------+

[toc] | [prev] | [next] | [standalone]

#321

Charles Richmond wrote:
> 
> On 5/16/11 8:46 AM, pete wrote:
> > Old Wolf wrote:
> >>
> >> Looking for an algorithm to find the fourth root of
> >> a large integer (I have a multi-precision integer
> >> library already). If there is no exact fourth root then
> >> the floor or ceiling is OK, I don't need a 'remainder'
> >> or anything.
> >>
> >> Any ideas?
> >
> > Why not just the square root of the square root?
> >
> 
> The fourth root of (2143/22) gives PI to about 8 decimal places.
> (Using your hand calculator, take the square root twice.)

The line in the original post which mentions
"no exact fourth root then the floor or ceiling is OK", 
suugests that OP is talking about integer math.

However, according to my calculator
                       (2143/22) == 97.409090909090909090909090909091
            square root(2143/22) == 9.8696043947612666756665932864345
square root(square root(2143/22))== 3.141592652582646125206037179644

So, I guess that I don't see your point.

-- 
pete

[toc] | [prev] | [next] | [standalone]

#317

On Sun, 15 May 2011 20:30:50 -0700 (PDT), Old Wolf
<oldwolf@inspire.net.nz> wrote:

>Looking for an algorithm to find the fourth root of
>a large integer (I have a multi-precision integer
>library already). If there is no exact fourth root then
>the floor or ceiling is OK, I don't need a 'remainder'
>or anything.
>
>Any ideas?
Try Newton-Raphson.

If x_n is an approximation to the fourth root of N then

 x_(n+1) = x_n - ((x_n)^4 - N) / (4 * (x_n)^3)

is a better approximation.

Iterate until the difference between x_n and x_(n+1) is small enough
for you.

You may be able to speed things up if you start with a reasonably
accurate approximation for x_0.  I have not tested it but 1 followed
by a quarter the number of digits in N worth of zeros might be a good
place to start.

rossum

[toc] | [prev] | [next] | [standalone]

#318

On May 15, 9:30 pm, Old Wolf <oldw...@inspire.net.nz> wrote:
> Looking for an algorithm to find the fourth root of
> a large integer (I have a multi-precision integer
> library already). If there is no exact fourth root then
> the floor or ceiling is OK, I don't need a 'remainder'
> or anything.
>
> Any ideas?

If you google "square root algorithm" you will find good numerical
approximation
algorithms and also an algorithm that computes the result exactly and
efficiently.
Applying these algorithms twice solves your problem.

The exact algorithm works in any integer base and is also
generalizable to
kth roots for positive integer k;  sorry I don't have references here
but I learned
it in a course so it must be right right?

If I were to implement the exact algorithm I would implement it in
binary as that
is easier.

Ralph Boland

[toc] | [prev] | [next] | [standalone]

#319

In article <40fdbc1f-8ca8-47bf-b701-5ca18bb9c3c9
@z13g2000prk.googlegroups.com>, oldwolf@inspire.net.nz says...
> 
> Looking for an algorithm to find the fourth root of
> a large integer (I have a multi-precision integer
> library already). If there is no exact fourth root then
> the floor or ceiling is OK, I don't need a 'remainder'
> or anything.
> 
> Any ideas?

Assuming you can get any sort of decent estimate of the root, you can 
implement this root polisher in multiple precision:

/*
xprimeB :=

            16        12          8  2        3  4       4    4
        x (x   + 120 x   a + 546 x  a  + 336 a  x  + 21 a ) (x  - a)
    x - -------------------------------------------------------------
         20        16           12  2         8  3        4  4      5
        x   + 203 x   a + 1554 x   a  + 1918 x  a  + 413 x  a  + 7 a

*/
double xprimeB(double x, double a)
{
    double x04 = x * x *x * x;
    double x08 = x04 * x04;
    double x12 = x08 * x04;
    double x16 = x12 * x04;
    double x20 = x16 * x04;

    double a2 = a * a;
    double a3 = a2 * a;
    double a4 = a3 * a;
    double a5 = a4 * a;

    return x - x * (x16 + 120.0 * x12 * a + 546.0 * x08 * a2 + 336.0 * 
a3 * x04 + 21.0 * a4) * (x04-a)/(x20 + 203.0  *  x16 * a  +  1554.0  *  
x12  * a2  +  1918.0  *  x08 * a3 + 413.0 * x04  * a4  +  7.0 * a5);
}

#ifdef UNIT_TEST
#include <math.h>
#include <stdio.h>
void polish(double x, double a)
{
	  int i;
		puts("root polisher iterations:");
    for (i = 0; i < 10; i++)
    {
    	x = xprimeB(x, a);
    	printf("estimate = %22.17g, for 4th root of %22.17g\n", x, a);
    }
}
int main(void)
{
	double test0 = 7166.0 * 7166.0 * 7166.0 * 7166.0 ;
	double test1 = 7166.5 * 7166.5 * 7166.5 * 7166.5 ;
	double test2 = 0.00001;
	double test3 = 1;
	polish(6000, test0);
	polish(6000, test1);
	polish(0.1, test2);
	polish(0.9, test3);
	return 0;
}
/* Possible output:
root polisher iterations:
estimate =     7166.0061387691758, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
estimate =                   7166, for 4th root of       
2636982303621136
root polisher iterations:
estimate =     7166.5061577672914, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
estimate =                 7166.5, for 4th root of       
2637718351152645
root polisher iterations:
estimate =   0.056395451781820163, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
estimate =   0.056234132519034911, for 4th root of 1.0000000000000001e-
005
root polisher iterations:
estimate =     1.0000000149948693, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
estimate =                      1, for 4th root of                      
1
*/
#endif

[toc] | [prev] | [next] | [standalone]

#336

In article <MPG.283b28c6fb000ad5989777@news.eternal-september.org>,
 Dann Corbit <dcorbit@connx.com> wrote:

> In article <40fdbc1f-8ca8-47bf-b701-5ca18bb9c3c9
> @z13g2000prk.googlegroups.com>, oldwolf@inspire.net.nz says...
> > 
> > Looking for an algorithm to find the fourth root of
> > a large integer (I have a multi-precision integer
> > library already). If there is no exact fourth root then
> > the floor or ceiling is OK, I don't need a 'remainder'
> > or anything.
> > 
> > Any ideas?
> 
> Assuming you can get any sort of decent estimate of the root, you can 
> implement this root polisher in multiple precision:
> 
> /*
> xprimeB :=
> 
>             16        12          8  2        3  4       4    4
>         x (x   + 120 x   a + 546 x  a  + 336 a  x  + 21 a ) (x  - a)
>     x - -------------------------------------------------------------
>          20        16           12  2         8  3        4  4      5
>         x   + 203 x   a + 1554 x   a  + 1918 x  a  + 413 x  a  + 7 a
> 
> */

A brief experiment with Mathematica reveals that this is a 7-th order 
approximation, in the sense that if x is close to the fourth root of a, 
then the relative error of xprimeB is the seventh power of the relative 
error of x. Roughly speaking, the number of valid significant figures in 
xprimeB is seven times that for x.

Gosh.

There has to be a general theory, for different orders and for different 
n-th roots. I would be grateful for references.

-- 
Christopher J. Henrich
chenrich@monmouth.com
http://www.mathinteract.com
"A bad analogy is like a leaky screwdriver." -- Boon

[toc] | [prev] | [next] | [standalone]

#338

On Wed, 18 May 2011 23:41:20 -0400, Christopher Henrich
<chenrich@monmouth.com> wrote:

>There has to be a general theory, for different orders and for different 
>n-th roots. I would be grateful for references.
Newton-Raphson?

rossum

[toc] | [prev] | [next] | [standalone]

#343

In article <ulm9t69b52grsmf761g790d5qm39c7smv8@4ax.com>,
 rossum <rossum48@coldmail.com> wrote:

> On Wed, 18 May 2011 23:41:20 -0400, Christopher Henrich
> <chenrich@monmouth.com> wrote:
> 
> >There has to be a general theory, for different orders and for different 
> >n-th roots. I would be grateful for references.
> Newton-Raphson?
> 
> rossum
Well, Newton-Raphson is a method of order 2. What I would like to see is 
a generalization to methods of order k >= 2, such as the one Dann Corbit 
offered, of order 7.

Here's a thought about how to do it. 

Suppose we are trying to solve the equation F(x) = A, and we have an 
approximate solution x_0; that is, if y_0 = F(x_0), then y_0 is close to 
A. The graph of the equation y = F(x) is a curve C; we draw the line 
tangent to the curve at (x_0, y_0), and we easily find the value of x_1 
such that (x_1, A) is on the tangent line. This is a graphical way of 
describing the Newton-Raphson method.

That tangent line has contact of order 2 with the curve C; being a 
straight line, it can be described as the graph of the equation x = P(y) 
where P is a fairly simple function, namely a first-degree polynomial.
Generalize: let P be a polynomial of degree k-1, such that the graphs of 
"y = F(x)" and "x = P(y)" have contact of order k at (x_0, y_0). Then 
let x_1 = P(A). 

Generalize a bit further: instead of a polynomial, use a rational 
function: x = P(y)/Q(y) where P and Q are polynomials whose degrees add 
up to k-1.

These ideas are a confluent form of "inverse interpolation." This is 
generally classified as a topic in the "calculus of finite differences." 
I am sure that confluent inverse interpolation has been invented more 
than once, but I don't know any specific references. If Dann Corbit is 
following this thread, I hope he will provide some.

-- 
Christopher J. Henrich
chenrich@monmouth.com
http://www.mathinteract.com
"A bad analogy is like a leaky screwdriver." -- Boon

[toc] | [prev] | [next] | [standalone]

#344 — Re: Fourth root of integer - fourth-root.mws [1/1]

[Multipart message — attachments visible in raw view] — view raw

In article <chenrich-9D8F8E.23411918052011@news.eternal-september.org>, 
chenrich@monmouth.com says...
> 
> In article <MPG.283b28c6fb000ad5989777@news.eternal-september.org>,
>  Dann Corbit <dcorbit@connx.com> wrote:
> 
> > In article <40fdbc1f-8ca8-47bf-b701-5ca18bb9c3c9
> > @z13g2000prk.googlegroups.com>, oldwolf@inspire.net.nz says...
> > > 
> > > Looking for an algorithm to find the fourth root of
> > > a large integer (I have a multi-precision integer
> > > library already). If there is no exact fourth root then
> > > the floor or ceiling is OK, I don't need a 'remainder'
> > > or anything.
> > > 
> > > Any ideas?
> > 
> > Assuming you can get any sort of decent estimate of the root, you can 
> > implement this root polisher in multiple precision:
> > 
> > /*
> > xprimeB :=
> > 
> >             16        12          8  2        3  4       4    4
> >         x (x   + 120 x   a + 546 x  a  + 336 a  x  + 21 a ) (x  - a)
> >     x - -------------------------------------------------------------
> >          20        16           12  2         8  3        4  4      5
> >         x   + 203 x   a + 1554 x   a  + 1918 x  a  + 413 x  a  + 7 a
> > 
> > */
> 
> A brief experiment with Mathematica reveals that this is a 7-th order 
> approximation, in the sense that if x is close to the fourth root of a, 
> then the relative error of xprimeB is the seventh power of the relative 
> error of x. Roughly speaking, the number of valid significant figures in 
> xprimeB is seven times that for x.
> 
> Gosh.
> 
> There has to be a general theory, for different orders and for different 
> n-th roots. I would be grateful for references.

The original reference would be:
Reference: "On Infinitely Many Algorithms For Solving Equations"
by Ernst Schroder, Translated by G. W. Stewart.

I made a generalization of his methods using hints he had in his paper 
for further experiments.  The 'B' family of solvers is extremely robust 
(IOW, it won't get stymied by mulitple roots and other things that in 
general can be problematic).

I have a generic Maple worksheet that can create equation solvers for 
anything that is differentiable.

You can find all the related stuff here (including lots of C and C++ 
samples that were generated):
http://cap.connx.com/schroder/

Also, attached is the specific worksheet I used to generate the polisher 
in question.

[toc] | [prev] | [next] | [standalone]

#345 — Re: Fourth root of integer - fourth-root.mws [1/1]

In article <MPG.284064176cd8327c989778@news.eternal-september.org>,
 Dann Corbit <dcorbit@connx.com> wrote:

> > There has to be a general theory, for different orders and for different 
> > n-th roots. I would be grateful for references.
> 
> The original reference would be:
> Reference: "On Infinitely Many Algorithms For Solving Equations"
> by Ernst Schroder, Translated by G. W. Stewart.
> 
> I made a generalization of his methods using hints he had in his paper 
> for further experiments.  The 'B' family of solvers is extremely robust 
> (IOW, it won't get stymied by mulitple roots and other things that in 
> general can be problematic).
> 
> I have a generic Maple worksheet that can create equation solvers for 
> anything that is differentiable.
> 
> You can find all the related stuff here (including lots of C and C++ 
> samples that were generated):
> http://cap.connx.com/schroder/
> 
> Also, attached is the specific worksheet I used to generate the polisher 
> in question.

Thank you!

-- 
Christopher J. Henrich
chenrich@monmouth.com
http://www.mathinteract.com
"A bad analogy is like a leaky screwdriver." -- Boon

[toc] | [prev] | [next] | [standalone]

#337

On Mon, 16 May 2011 13:28:29 -0700, Dann Corbit wrote:
>  oldwolf@inspire.net.nz says...
>> 
>> Looking for an algorithm to find the fourth root of a large integer (I
>> have a multi-precision integer library already). If there is no exact
>> fourth root then the floor or ceiling is OK, I don't need a 'remainder'
>> or anything.
...
> Assuming you can get any sort of decent estimate of the root, you can
> implement this root polisher in multiple precision:
> 
> /*
> xprimeB :=
> 
>             16        12          8  2        3  4       4    4
>         x (x   + 120 x   a + 546 x  a  + 336 a  x  + 21 a ) (x  - a)
>     x - -------------------------------------------------------------
>          20        16           12  2         8  3        4  4      5
>         x   + 203 x   a + 1554 x   a  + 1918 x  a  + 413 x  a  + 7 a
> 
> */
> double xprimeB(double x, double a)

[snip code and sample output]

If this is a 7th-order approximation (as suggested in Christopher 
Henrich's followup) it won't require many iterations, even for numbers
with many digits.  However, in multiple precision arithmetic, even a 
single iteration of the above may be slower than just finding the 
square root twice.  Typical MP packages use efficient square root 
computations that take slightly longer than a single divide.

For example, on my computer using mpfr (gmp) and with default 
precision 256 bits, one iteration of xprimeB takes ~ 3.8 times as 
long as 2 mpfr_sqrt calls.  Here are some times, in seconds per 
half-million passes.  DeltaI is t_I-t_1.

I   Item          t_I
1. Empty loop:  0.247415   Passes: 500000
2. 2-Sqrt:      1.015866   Delta2:  0.768451
3. xprimeB:     3.172006   Delta3:  2.924591
4. divide:      0.500479   Delta4:  0.253064
(t_3-t_1)/(t_2-t_1):    3.805826

Code for loop 2 that took about .77 seconds longer than the 
empty loop consists of following, where RM is MPFR_RNDU rounding 
mode, RSEED=13343, ir is a 32-bit unsigned long int, and ttime() 
is as shown below sig:

  srandom(RSEED);
  t0 = ttime(0);
  for (i=0; i<REPS; ++i) {
    ir = random ();
    mpfr_set_ui (q, ir, RM);
    mpfr_sqr  (q2, q, RM);
    mpfr_sqr  (q4, q2, RM);
    mpfr_sqrt (r2, q4, RM);
    mpfr_sqrt (r,  r2, RM);
  }
  t2 = ttime(t0);

The "empty loop" is like loop 2 but without mpfr_sqrt() 
calls.

Loop 3 is like loop 2 but with 
    xprimeB (q, q4, pkb);
in place of the two mpfr_sqrt() calls.   xprimeB() was
coded with mpfr_mul(), mpfr_add(), mpfr_div() etc, using
global mpfr_t variables initialized when main() began.

Loop 4 is like loop 2 but has
    mpfr_div (pkb, q4, q2, RM);
in place of the two mpfr_sqrt() calls.  For numbers as in 
this test, 2 sqrts take about as long as 3 divides.

-- 
jiw

double ttime(double base) {
  struct timeval tod;
  gettimeofday(&tod, NULL);
  return tod.tv_sec + tod.tv_usec/1e6 - base;
}

[toc] | [prev] | [next] | [standalone]

#346

On Thu, 19 May 2011 07:01:14 +0000, James Waldby wrote:
> On Mon, 16 May 2011 13:28:29 -0700, Dann Corbit wrote:
>>  oldwolf@inspire.net.nz says...
>>> 
>>> Looking for an algorithm to find the fourth root of a large integer (I
>>> have a multi-precision integer library already). If there is no exact
>>> fourth root then the floor or ceiling is OK, I don't need a
>>> 'remainder' or anything.
> ...
>> Assuming you can get any sort of decent estimate of the root, you can
>> implement this root polisher in multiple precision:
[snip equation]
> [snip code and sample output]
> [...] in multiple precision arithmetic, even a
> single iteration of the above may be slower than just finding the square
> root twice.  Typical MP packages use efficient square root computations
> that take slightly longer than a single divide.
> 
> For example, on my computer using mpfr (gmp) and with default precision
> 256 bits, one iteration of xprimeB takes ~ 3.8 times as long as 2
> mpfr_sqrt calls.  Here are some times, in seconds per half-million
> passes.  DeltaI is t_I-t_1.
> 
> I   Item          t_I
> 1. Empty loop:  0.247415   Passes: 500000
> 2. 2-Sqrt:      1.015866   Delta2:  0.768451
> 3. xprimeB:     3.172006   Delta3:  2.924591
> 4. divide:      0.500479   Delta4:  0.253064
> (t_3-t_1)/(t_2-t_1):   3.805826
> 
[snip explanations of loops]

Since OP wants 4th roots of integers, I recoded program from mpfr to
gmp (ie, using calls like mpz_sqrt, mpz_div, mpz_mul rather than 
mpfr_sqrt, mpfr_div, mpfr_mul) and got the following timings, for
the cases where the 4th powers are 128-bit numbers (like in the 
mpfr test results above) and 256-bit numbers:

I   Item          t_I      Bits:    128
1. Empty loop:  0.062849   Passes:  500000
2. 2-Sqrt:      0.348615   Delta2:  0.285766
3. xprimeB:     1.949296   Delta3:  1.886447
4. divide:      0.136300   Delta4:  0.073451
(t_3-t_1)/(t_2-t_1):    6.601371
(t_2-t_1)/(t_4-t_1):    3.890576

I   Item          t_I      Bits:    256
1. Empty loop:  0.126367   Passes:  500000
2. 2-Sqrt:      0.575044   Delta2:  0.448677
3. xprimeB:     4.153543   Delta3:  4.027176
4. divide:      0.286537   Delta4:  0.160170
(t_3-t_1)/(t_2-t_1):    8.975670
(t_2-t_1)/(t_4-t_1):    2.801253

In the latter case, in each of 4 loops, the same 500000 random 
64-bit numbers were each raised to the 4th power, and treated 
as follows:  Loop 1: nothing else; Loop 2: 2 square roots;  
Loop 3: 1 pass of xprimeB from Dann Corbit's post; Loop 4: 
4th power divided by 2nd power.  The two ratios show that 1 
pass of xprimeB took about 9 times as long as extracting 2 
square roots, and that  extracting 2 square roots took about
3 times as long as 1 divide.

-- 
jiw

[toc] | [prev] | [standalone]

Back to top | Article view | comp.programming

csiph-web