99 ^ 99 in C

[asit]

I know it is really tough to calculate 99^99 in C. Its not only
difficult but also time consuming to get the o/p.

But how perl interpreter calculates it( $x = 99 ** 99; print 4x
within seconds???/

which algorithm it uses ???

 

[vippstar]

I know it is really tough to calculate 99^99 in C. Its not only
difficult but also time consuming to get the o/p.

You're mistaken.
99^99 is 0 in C. (^ is XOR, not some expontation operator)

But how perl interpreter calculates it( $x = 99 ** 99; print 4x
within seconds???/
which algorithm it uses ???

Which perl implementation? (Don't answer)
The implementation is allowed to calculate it however it wants.
Even s/99 **
99/369729637649726772657187905628805440595668764281741102430259972423552570455277523421410650010128232727940978889548326540119429996769494359451621570193644014418071060667659301384999779999159200499899/
;-)

For a serious answer, you might want to look into bignum libraries in
C.
GMP is one of the best and free.
<http://gmplib.org/>

To calculate 99 ** 99, you'd have:
mpz_t x;
mpz_init(x);
mpz_ui_powm_ui(x, 99, 99); /* x = 99 ** 99 */
/* ... */
mpz_clear(x);


However, it's important to note you shouldn't ask further questions
about gmp lib, since it's off-topic here.

 

[Flash Gordon]

asit wrote, On 27/09/08 13:55:

I know it is really tough to calculate 99^99 in C. Its not only
difficult but also time consuming to get the o/p.

How time consuming it is depends on how you do it. As evidence you have
what you wrote below...

But how perl interpreter calculates it( $x = 99 ** 99; print 4x
within seconds???/

which algorithm it uses ???

Look at the source code for Perl. Possibly try asking on an appropriate
perl group or mailing list.

 

[Bartc]

You're mistaken.
99^99 is 0 in C. (^ is XOR, not some expontation operator)

This always comes up. In English, ^ if often used to mean 'raised to the
power of'. The rest of the context makes this clear. In fact only in C code
is ^ used for exclusive-or.

Which perl implementation? (Don't answer)
The implementation is allowed to calculate it however it wants.
Even s/99 **
99/369729637649726772657187905628805440595668764281741102430259972423552570455277523421410650010128232727940978889548326540119429996769494359451621570193644014418071060667659301384999779999159200499899/
;-)

How accurate an answer is needed? In C I can do pow(99.0,99) with no trouble
and I can do it ten million times in half a second. About 3.7e197 (or
3.7x10^197).

 

[vippstar]

This always comes up. In English, ^ if often used to mean 'raised to the
power of'. The rest of the context makes this clear. In fact only in C code
is ^ used for exclusive-or.

You obviously missed the part of my post after that saying
"For a serious answer, [...]"

 

[Bartc]

I know it is really tough to calculate 99^99 in C. Its not only
difficult but also time consuming to get the o/p.

But how perl interpreter calculates it( $x = 99 ** 99; print 4x
within seconds???/

Within /seconds/? You mean microseconds?

Assuming you want an exact answer:

I calculated 99**99 on a big int library of mine that works with individual
decimal digits and runs under an interpreted language. You can't get much
slower. But it only took 1 millisecond.

which algorithm it uses ???

The algorithm, if you can call it that, consists of multiplying by 99, 98
times. There are better ways of doing it.

If you don't want to use an existing library, a function in C, working with
fixed precision numbers (some sort of array of ints), to multiply like this
wouldn't be that difficult to program (but neither trivial enough to present
here..).

The problem is displaying the answer, which will likely involve division,
rather more difficult.

 

[Keith Thompson]

I know it is really tough to calculate 99^99 in C. Its not only
difficult but also time consuming to get the o/p.

But how perl interpreter calculates it( $x = 99 ** 99; print 4x
within seconds???/

which algorithm it uses ???

<OT>
According to Perl's documentation ("perldoc perlop"), its "**"
operator simply uses the C pow() function. And no, it doesn't take
"seconds" to execute.
</OT>

 

[Ben Pfaff]

This always comes up. In English, ^ if often used to mean 'raised to
the power of'. The rest of the context makes this clear. In fact only
in C code is ^ used for exclusive-or.

I noticed yesterday that even otherwise respected book authors
get this wrong. Steve McConnell's _Code Complete_ contains the
following code on page 194 labeled as a "C Example":

/* Compute roots of a quadratic equation.
This assumes that (b^2-4*a*c) is positive. */

Temp = sqrt( b^2 - 4*a*c );
root[0] = ( -b + Temp ) / ( 2 * a );
root[1] = ( -b - Temp ) / ( 2 * a );

It will, of course, compute nothing of the sort and in fact b^2
is a constraint violation if 'b' is declared as a floating-point
type.

 

[Antoninus Twink]

To calculate 99 ** 99, you'd have:
mpz_t x;
mpz_init(x);
mpz_ui_powm_ui(x, 99, 99); /* x = 99 ** 99 */

Wrong. You need mpz_ui_pow_ui().

(The *powm* functions do exponentiation modulo some integer.)

 

[Antoninus Twink]

asit said:

Nope. Assuming you mean "to the power of", it's
13C8DB009372548E1D1EE21004B41F9D0A2F1C7A4\
CA716024FF714032153A5AB531817C13EF1F2F7F4\
01B48535A8C22BBDC1EE93B772F0DD094A8FD9727\
DAE16CEB8E4A25451BBA34A32099615746B57BC8BB

and it took my C program about a hundredth of a second to find this
out.

That seems incredibly slow unless you're running it on a PDP-11.

Is that using your hacked together "do everything with integers
represented as strings" library that you seem to be very proud of? Why
not use an off-the-shelf professional-grade multiprecision library, of
the sort you're obviously incapable of writing yourself?

 

[Bartc]

Malcolm McLean said:


Yes, it is. Given that high performance was not a design goal, I don't
think it does too shabbily.


<shrug> Numbers are numbers. Conversion to decimal is available, but takes
a little longer.

Well, I thought it would. It seems I was mistaken. Adding the conversion
appears to have stripped a couple of ms off the runtime!

369729637649726772657187905628805440595668764281\
741102430259972423552570455277523421410650010128\
232727940978889548326540119429996769494359451621\
570193644014418071060667659301384999779999159200\
499899

real 0m0.008s
user 0m0.000s
sys 0m0.000s

I suspect this may be faster than you think. How long does it take to
calculate 99**99 one thousand times, not including printing the result?

 

[Peter Nilsson]

I know it is really tough to calculate 99^99 in C. Its not
only difficult but also time consuming to get the o/p.
...

Not really. Even my old 1980's 1MHz Hitachi Peach managed to
do this sort of stuff in a few seconds with interpreted Basic.

% type 99.c
#include <stdio.h>
#include <stdlib.h>
#include <string.h>

const char digit[] = "0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZ";

int main(int argc, char **argv)
{
static char ans[1024] = "1";
int i, acc, carry;
long b = 0;
char *p;

if (argc > 1) b = strtol(argv[1], 0, 10);
if (b < 2 || b > 36) b = 10;

for (i = 0; i < 99; i++)
{
carry = 0;
for (p = ans; *p; p++)
{
acc = carry + 99 * (strchr(digit, *p) - digit);
*p = digit[acc % b];
carry = acc / b;
}

while (carry)
{
*p++ = digit[carry % b];
carry /= b;
}

*p = 0;
}

for (i = 0, p = ans + strlen(ans); ans < p; )
{
putchar(*--p);
if (++i % 40 == 0) putchar('\n');
}
if (i % 40) putchar('\n');

return 0;
}

% acc 99.c

% a
3697296376497267726571879056288054405956
6876428174110243025997242355257045527752
3421410650010128232727940978889548326540
1194299967694943594516215701936440144180
71060667659301384999779999159200499899

% a 16
13C8DB009372548E1D1EE21004B41F9D0A2F1C7A
4CA716024FF714032153A5AB531817C13EF1F2F7
F401B48535A8C22BBDC1EE93B772F0DD094A8FD9
727DAE16CEB8E4A25451BBA34A32099615746B57
BC8BB

% a 11
56114575632052945390991A9741129022605573
0373168542389887382353926641A47423956A45
6377959452500000000000000000000000000000
0000000000000000000000000000000000000000
000000000000000000000000000000

 

[Willem]

Richard Heathfield wrote:
) Simpler to program, without increasing the number of multiplications:
)
) u = 99 * 99 * 99
) v = u * u
) w = v * v
) x = w * w
) y = 99 * x * x
) z = 99 * y * y

Note the bit pattern of 99: 1100011

a = 99 (^1) 1
b = 99 * a * a (^3) 1
c = b * b (^6) 0
d = c * c (^12) 0
e = d * d (^24) 0
f = 99 * e * e (^49) 1
g = 99 * f * f (^99) 1

Which leads to the general O(log N) algorithm:

double pow(double fac, unsigned int exp)
{
double res = 1.0;
unsigned int mask = 1;
while (mask <= (exp >> 1)) {
mask <<= 1;
}
while (mask > 0) {
if (exp & mask) {
res = res * fac;
}
res = res * res;
mask >>= 1;
}
return res;
}

NB: This does one multiply too many,
which could be optimized by making (exp == 0) a special case.


SaSW, Willem
--
Disclaimer: I am in no way responsible for any of the statements
made in the above text. For all I know I might be
drugged or something..
No I'm not paranoid. You all think I'm paranoid, don't you !
#EOT

 

[user923005]

This always comes up. In English, ^ if often used to mean 'raised to
the power of'. The rest of the context makes this clear. In fact only
in C code is ^ used for exclusive-or.

I noticed yesterday that even otherwise respected book authors
get this wrong.  Steve McConnell's _Code Complete_ contains the
following code on page 194 labeled as a "C Example":

    /* Compute roots of a quadratic equation.
       This assumes that (b^2-4*a*c) is positive. */

    Temp    = sqrt( b^2 - 4*a*c );
    root[0] = ( -b + Temp ) / ( 2 * a );
    root[1] = ( -b - Temp ) / ( 2 * a );

It will, of course, compute nothing of the sort and in fact b^2
is a constraint violation if 'b' is declared as a floating-point
type.

Even if this were done:
Temp = sqrt( b*b - 4*a*c );
instead of this:
Temp = sqrt( b^2 - 4*a*c );
it's still a very, very poor way to calculate the quadratic formula
(if b*b is about equal to 4*a*c then up to half of the digits of
precision are lost).
See:
http://www.physics.ohio-state.edu/~dws/grouplinks/floating_point_math.pdf
pages 180-182

Still plenty of room for improvement with this method, but it is
better than the naive:

#include <stdio.h>
#include <math.h>

void rootcalc(double a, double b, double c, double *br1,
double *br2, double *gr1, double *gr2)
{

double four_a_c = 4.0 * a * c;
double temp = sqrt(b * b - four_a_c);

*br1 = (-b + temp) / (2.0 * a);
*br2 = (-b - temp) / (2.0 * a);

*gr1 = -2.0 * c / (b + temp);
*gr2 = -2.0 * c / (b - temp);

printf("y = %23.18gx^2 + %23.18gx + %23.18g\n", a, b, c);
printf("naive roots=%23.18g %23.18g\ncareful roots=%23.18g %23.18g
\n",
*br1, *br2, *gr1, *gr2);
printf("relative errors = %23.18g %23.18g, 4*a*c = %23.18g\n\n",
(*gr1 - *br1) / *gr1, (*br2 - *gr2) / *br2, four_a_c);
}

int main(void)
{
double a,
b,
c;
double br1=0,
br2=0;
double gr1=0,
gr2=0;

double temp;
double four_a_c;
int n;

a = 1.0;
b = 1.0;
c = 1.0;

for (n = 1; n <= 18; n++) {
c /= 10.0;
rootcalc(a, b, c, &br1, &br2, &gr1, &gr2);
}
rootcalc(1.0, 5000000.0, 0.00001, &br1, &br2, &gr1, &gr2);
rootcalc(6.0, 5.0, -4.0, &br1, &br2, &gr1, &gr2);
rootcalc(6.0e30, 5.0e30, -4.0e30, &br1, &br2, &gr1, &gr2);
rootcalc(1.0, -1.0e5, 1.0, &br1, &br2, &gr1, &gr2);
rootcalc(1.0, -4.0, 3.999999, &br1, &br2, &gr1, &gr2);
rootcalc(1.0e-30, -1.0e30, 1.0e30, &br1, &br2, &gr1, &gr2);
/* Linear degenerate case left as an exercise for the reader: */
rootcalc(0.0, 1.0, 1.0, &br1, &br2, &gr1, &gr2);
/* Constant degenerate case left as an exercise for the reader: */
rootcalc(0.0, 0.0, 1.0, &br1, &br2, &gr1, &gr2);
return 0;
}

 

[Pilcrow]

asit said:


Nope. Assuming you mean "to the power of", it's
13C8DB009372548E1D1EE21004B41F9D0A2F1C7A4\
CA716024FF714032153A5AB531817C13EF1F2F7F4\
01B48535A8C22BBDC1EE93B772F0DD094A8FD9727\
DAE16CEB8E4A25451BBA34A32099615746B57BC8BB

and it took my C program about a hundredth of a second to find this out.



I suppose it depends on how often you want to do it. A hundredth of a
second here, a fiftieth there, pretty soon you're talking serious time.


SECONDS? Sheesh, no wonder nobody uses Perl any more.


A slow one, it seems.

I just calculated 99**99 one million times with perl. It took less than
0.18 second, even with the program overhead. Here is the code, and the
result:
-------------------------------------------------------------
C:\perl
use Time::HiRes 'gettimeofday';
($sec, $usec) = Time::HiRes::gettimeofday;
for(0..1000000){$a = 99.**99}
($sec2, $usec2) = Time::HiRes::gettimeofday;
$tsec = $sec2 -$sec;
$tusec = $usec2 - $usec;
print "$tsec $tusec\n";
^D
0 179863
--------------------------------------------------------------

 

[Pilcrow]

I just calculated 99**99 one million times with perl. It took less than
0.18 second, even with the program overhead. Here is the code, and the
result:
-------------------------------------------------------------
C:\perl
use Time::HiRes 'gettimeofday';
($sec, $usec) = Time::HiRes::gettimeofday;
for(0..1000000){$a = 99.**99}
($sec2, $usec2) = Time::HiRes::gettimeofday;
$tsec = $sec2 -$sec;
$tusec = $usec2 - $usec;
print "$tsec $tusec\n";
^D
0 179863
--------------------------------------------------------------

oops, it only calculated it once, then optimisation converted it to a
constant. but it went thru the loop a million times. hard to know how
long it took for the calc, but probably a usec or so.

 

[Ben Bacarisse]

Pilcrow said:

Well, it's easy to find out - just omit the loop, and re-time.

I did this, and got 0.004 seconds. I was suspicious, so I made the thing
print its result, getting 0.017 seconds, which suggests that the 0.004 for
just the calc was probably right. Trouble is, the print gave me this:

3.69729637649727e+197

which is too high by
2273428120943711945594043312357182588975697400275764474295447224765\
78589349989871767272059021110451673459880570003230505640548378429806\
355985581928939332340698615000220000840799500101

If the answer doesn't have to be right, I can get it in 0 seconds.

Even if the answer has to be exact you can get it in 0 seconds. There
is no need to spend any time at all calculating a know constant.

 

[Bartc]

oops, it only calculated it once, then optimisation converted it to a
constant. but it went thru the loop a million times. hard to know how
long it took for the calc, but probably a usec or so.

Optimisers are a pain when you want to do benchmarks.

Suppose you just repeat the $a=99.**99 statement one million times in the
source? (Better start with just 1000 times.) It may not be clever enough to
optimise that.

But that "." in "99." is worrying; I think this is supposed to be an integer
calculation.

 

[Bartc]

Bartc said:


I removed the . and ran the script again, but it made no difference to the
result, which is still wrong by a huge amount (albeit a vanishingly small
percentage).

I tried this (after a quick course in perl):

use Math::BigInt;
sub bint { Math::BigInt->new(shift); }

$x = bint(99)**bint(99);

print $x;

With the $x= line repeated 1000 times, it took about a second on my low-end
Windows machine.

And gave the right answer.

 

[Willem]

Bartc wrote:
) I tried this (after a quick course in perl):
)
) use Math::BigInt;
) sub bint { Math::BigInt->new(shift); }
)
) $x = bint(99)**bint(99);
)
) print $x;
)
) With the $x= line repeated 1000 times, it took about a second on my low-end
) Windows machine.
)
) And gave the right answer.

This is comp.lang.c after all:

(Don't mind the scary memory leaks... )

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

char * bmul(const char *a, const char *b)
{
size_t i, j, al, bl;
char *res, s, c;
al = strlen(a);
bl = strlen(b);
res = calloc(al + bl + 1, 1);

for (i = 0; i < al; i++) {
c = 0;
for (j = 0; j < bl; j++) {
s = c + (a - '0') * (b[j] - '0');
if (res[i + j]) s += res[i + j] - '0';
c = s/10;
s = s%10;
res[i + j] = s + '0';
}
while (c) {
s = c;
if (res[i + j]) s += res[i + j] - '0';
c = s/10;
s = s%10;
res[i + j] = s + '0';
j++;
}
}
return res;
}

char *bpow(const char *num, int exp)
{
char *res = "1";
int mask = 1;
while (mask <= (exp >> 1)) {
mask <<= 1;
}
while (mask) {
res = bmul(res, res);
if (exp & mask) {
res = bmul(num, res);
}
mask >>= 1;
}
return res;
}

char *rev(const char *str)
{
size_t l, i;
char *res;
l = strlen(str);
res = malloc(l + 1);
for (i = 0; i < l; i++) {
res[l-i-1] = str;
}
res[l] = 0;
return res;
}

int main(int argc, char *argv[])
{
int i;
if (argc < 3) return 1;
for (i = 0; i < 1000; i++) {
bpow(rev(argv[1]), atoi(argv[2]));
}
printf("%s ^ %d = %s\n", argv[1], atoi(argv[2]),
rev(bpow(rev(argv[1]), atoi(argv[2]))));
return 0;
}


Which took about 0.4 seconds for those iterations.

(Then I noticed that I forgot to turn on optimization,
which cut the time to less than 0.2 seconds.)


SaSW, Willem
--
Disclaimer: I am in no way responsible for any of the statements
made in the above text. For all I know I might be
drugged or something..
No I'm not paranoid. You all think I'm paranoid, don't you !
#EOT