Hello again everyone,
I have been promising myself that I would do some work on a library for large numbers since college. I read about them in a textbook by Binstock and Rex, but I didn't like the static nature of their implementation and I had a few other issues so I wrote my own.
Here is an excerpt from the main driver, test is currently set to 9999:
The program kept up with the system's graphical calculator (xcalc) until about 9^35 when the calculator produced an error message. The loop above produced 9^10000, it is in the file that I have attached. It has about the right number of digits 9494, so I'm rather confident with it.Code:a = bigintread("9");
b = bigintread("9");
for(i = 0; i < test; i++, bigintcpy(c, b))
{
bigintmult(a, b, c);
}
printf( "------------------------\n" );
printf( "for( i = 0; i < test; i++ ) c = c*a\n");
bigintprint(stdout, a, 0);
printf( "\n" );
bigintprint(stdout, b, 0);
printf( "\n" );
bigintprint(stdout, c, 0);
printf( "\n" );
Does anyone have that figure handy to verify it?
The program was running in the debugger for about seven hours before producing that figure, one factor of nine at a time lol. I think maybe it would be prudent to build a database for previously performed calculations. Then use a separate thread to find them if any while waiting for a calculation thread to finish, and using the first that came back successful of course.
The code for multiplication is here:
Other alternative implementations also came to mind. Like using a queue rather than a table format for the number datastructure, but it made more sense to me that it grow a 'page' at a time rather than double in size as a queue for a stream might according to the suggestion in the algorithm analysis textbook. If a carry had more than an allocation unit of figures, the instrinsic data type would have already overflowed is my justification. What sort of advantages might I have missed from using a queue for representation?Code:/*** bigintmult.c **********************************************
*
* Multiply two big integers a and b. Place the result into c.
*
*/
#include "bigint.h"
bigint *bigintmult( const bigint *a, const bigint *b, bigint *c )
{
bigint *r; // pointer to the current register (product)
bigint *swp; //
bigint *tmp;
short i;
short j;
short k;
short h;
int carry;
// verify/arrange for the register
if( NULL == c )
{
r = bigintcrt();
}
else
{
r = bigintzero(c);
}
// test that the operands are sane
if(not_a_number == a->s || not_a_number == b->s)
{
// NaN * b = NaN ^ a * NaN = NaN ^ NaN * NaN = NaN
r->s = not_a_number;
return r;
}
if(infinite_positive == a->s || infinite_negative == a->s)
{
// infinite * b = infinite
r->s = a->s;
return r;
}
if(infinite_positive == b->s || infinite_negative == b->s)
{
// a * infinite = infinite
r->s = b->s;
return r;
}
if(zero == a->s || zero == b->s)
{
return r;
}
// a != 0 ^ b != 0, a,b non-ideal integers
// clear and resize the register in case of overflow
tmp = bigintcrt();
bigintrealloc(tmp, 2*MAX(a->u, b->u) + 1);
swp = bigintcrt();
// perform multiplication
for(i = 0; i < b->mlen; ++i)
{
if(ZERO == b->m[i])
{
continue;
}
for(j = 0; j < a->mlen; j++)
{
tmp->m[i+j] = a->m[j]*b->m[i];
if(zero == tmp->s && ZERO != tmp->m[j+i])
{
tmp->s = positive;
}
}
// -- normalize tmp for addition to the result register --
for(k = 0; k < tmp->u*ARBITRARYMATH_ALLOCATION_UNIT; ++k)
{
if(9 < tmp->m[k])
{
carry = tmp->m[k]/10;
tmp->m[k] = tmp->m[k];
for(h = 1;
0 < carry
&& (k+h) < tmp->u*ARBITRARYMATH_ALLOCATION_UNIT;
++h, carry /= 10) // high_order_truncation will occur here
{
tmp->m[k+h] += (carry);
if(h+k+1 == tmp->u*ARBITRARYMATH_ALLOCATION_UNIT
&& (tmp->u+1)*ARBITRARYMATH_ALLOCATION_UNIT <= SHRT_MAX)
{
tmp->mlen = tmp->u*ARBITRARYMATH_ALLOCATION_UNIT;
bigintrealloc(tmp, tmp->u+1);
}
}
}
}
for(tmp->mlen = tmp->u*ARBITRARYMATH_ALLOCATION_UNIT;
0 < tmp->mlen; --tmp->mlen)
{
if(ZERO != tmp->m[tmp->mlen-1])
{
break;
}
}
bigintadd(r, tmp, swp);
bigintcpy(swp, r);
memset(tmp->m, ZERO, tmp->u*ARBITRARYMATH_ALLOCATION_UNIT*sizeof(bigdigit));
tmp->mlen = 1;
tmp->s = zero;
}
// select the correct sign
if( a->s != b->s )
{
r->s = negative;
}
else
{
r->s = positive;
}
// clean up
bigintdes(tmp);
bigintdes(swp);
bigintlength(r);
if(ceil(((float)r->mlen)/ARBITRARYMATH_ALLOCATION_UNIT) < r->u)
{
bigintrealloc(r, ceil(((float)r->mlen)/ARBITRARYMATH_ALLOCATION_UNIT));
}
return r;
}
Also, in order to wring out the last bit of speed from this Pentium 4, I think that the point of parallelism is at the intermediate products before summing. That doesn't help with the situation when multiplication occurs one unit at a time. Perhaps partitioning the figure by place value is more appropriate? Tips here would be great too. I had a great parallel algorithm analysis textbook while I was in college, but a friend indeed needed one.
Thank you for your time and effort and it is great to post with you again.
Best Regards,
New Ink -- Henry
Enclosure: Attachment 11274