ref: 32b537dbdb56584c3dfdb9f95f1d1cec0359d7d7
dir: /3rd/mp/test.c/
#include "platform.h"
#include "mp.h"
#include "ieee754.h"
double D_PNAN, D_NNAN, D_PINF, D_NINF;
float F_PNAN, F_NNAN, F_PINF, F_NINF;
static int loops = 1;
static char str[16][8192];
static int istr = -1;
static char *
MPB(int base, mpint *m)
{
char *s, *b;
bool minus;
istr = (istr+1) % nelem(str);
b = str[istr];
s = mptoa(m, base, b+2, sizeof(str[istr])-2);
if(base == 10)
return s;
minus = s[0] == '-';
if(minus)
s++;
if(base == 2){
*--s = 'b';
*--s = '0';
}else if(base == 8){
*--s = '0';
}else if(base == 16){
*--s = 'x';
*--s = '0';
}
if(minus)
*--s = '-';
return s;
}
static char *
MP(mpint *m)
{
char *b;
istr = (istr+1) % nelem(str);
b = str[istr];
return mptoa(m, 10, b, sizeof(str[istr]));
}
static int64_t
nsec(void)
{
return 0;
}
static void
prng(uint8_t *p, int n)
{
while(n-- > 0)
*p++ = rand();
}
static void
testconv(char *str)
{
int i, base[] = {2,8,10,16/*,32,64*/};
mpint *b;
char *p;
printf("testconv \"%s\":\n", str);
b = strtomp(str, nil, 16, nil);
for(i=0; i < nelem(base); i++){
p = mptoa(b, base[i], nil, 0);
if(p == nil){
fprintf(stderr, "mptoa(base=%d) -> nil\n", base[i]);
exit(2);
}
printf("base%d: %s = ", base[i], p);
if(strtomp(p, nil, base[i], b) == nil){
fprintf(stderr, "A strtomp(%s) -> nil\n", p);
exit(2);
}
free(p);
printf("%s\n", MPB(base[i], b));
switch(base[i]){
case 2:
case 8:
case 10:
case 16:
asprintf(&p, "%s", MPB(base[i], b));
printf("# %s = ", p);
if(strtomp(p, nil, 0, b) == nil){
fprintf(stderr, "B strtomp(%s) -> nil\n", p);
exit(2);
}
free(p);
printf("%s\n", MPB(base[i], b));
break;
}
}
mpfree(b);
}
static void
testshift(char *str)
{
mpint *b1, *b2;
int i;
b1 = strtomp(str, nil, 16, nil);
b2 = mpnew(0);
for(i = 0; i < 64; i++){
mpleft(b1, i, b2);
printf("%2.2d %s\n", i, MP(b2));
}
for(i = 0; i < 64; i++){
mpright(b2, i, b1);
printf("%2.2d %s\n", i, MP(b1));
}
mpfree(b1);
mpfree(b2);
}
static void
testaddsub(char *str)
{
mpint *b1, *b2;
int i;
b1 = strtomp(str, nil, 16, nil);
b2 = mpnew(0);
for(i = 0; i < 16; i++){
mpadd(b1, b2, b2);
printf("%2.2d %s\n", i, MP(b2));
}
for(i = 0; i < 16; i++){
mpsub(b2, b1, b2);
printf("%2.2d %s\n", i, MP(b2));
}
mpfree(b1);
mpfree(b2);
}
static void
testvecdigmuladd(char *str, mpdigit d)
{
mpint *b, *b2;
int i;
int64_t now;
b = strtomp(str, nil, 16, nil);
b2 = mpnew(0);
mpbits(b2, (b->top+1)*Dbits);
now = nsec();
for(i = 0; i < loops; i++){
memset(b2->p, 0, b2->top*Dbytes);
mpvecdigmuladd(b->p, b->top, d, b2->p);
}
if(loops > 1)
printf("%"PRId64" ns for a %d*%d vecdigmul\n", (nsec()-now)/loops, b->top*Dbits, Dbits);
mpnorm(b2);
printf("0 + %s * %"PRIx64" = %s\n", MP(b), (uint64_t)d, MP(b2));
mpfree(b);
mpfree(b2);
}
static void
testvecdigmulsub(char *str, mpdigit d)
{
mpint *b, *b2;
int i;
int64_t now;
b = strtomp(str, nil, 16, nil);
b2 = mpnew(0);
mpbits(b2, (b->top+1)*Dbits);
now = nsec();
for(i = 0; i < loops; i++){
memset(b2->p, 0, b2->top*Dbytes);
mpvecdigmulsub(b->p, b->top, d, b2->p);
}
if(loops > 1)
printf("%"PRId64" ns for a %d*%d vecdigmul\n", (nsec()-now)/loops, b->top*Dbits, Dbits);
mpnorm(b2);
printf("0 - %s * %"PRIx64" = %s\n", MP(b), (uint64_t)d, MP(b2));
mpfree(b);
mpfree(b2);
}
static void
testmul(char *str)
{
mpint *b, *b1, *b2;
int64_t now;
int i;
b = strtomp(str, nil, 16, nil);
b1 = mpcopy(b);
b2 = mpnew(0);
now = nsec();
for(i = 0; i < loops; i++)
mpmul(b, b1, b2);
if(loops > 1)
printf("%"PRId64" µs for a %d*%d mult\n", (nsec()-now)/(loops*1000),
b->top*Dbits, b1->top*Dbits);
printf("%s * %s = %s\n", MP(b), MP(b1), MP(b2));
mpfree(b);
mpfree(b1);
mpfree(b2);
}
static void
testmul2(mpint *b, mpint *b1)
{
mpint *b2;
int64_t now;
int i;
b2 = mpnew(0);
now = nsec();
for(i = 0; i < loops; i++)
mpmul(b, b1, b2);
if(loops > 1)
printf("%"PRId64" µs for a %d*%d mult\n", (nsec()-now)/(loops*1000), b->top*Dbits, b1->top*Dbits);
printf("%s * ", MP(b));
printf("%s = ", MP(b1));
printf("%s\n", MP(b2));
mpfree(b2);
}
static void
testdigdiv(char *str, mpdigit d)
{
mpint *b;
mpdigit q;
int i;
int64_t now;
b = strtomp(str, nil, 16, nil);
now = nsec();
for(i = 0; i < loops; i++)
mpdigdiv(b->p, d, &q);
if(loops > 1)
printf("%"PRId64" ns for a %d / %d div\n", (nsec()-now)/loops, 2*Dbits, Dbits);
printf("%s / %"PRIx64" = %"PRIx64"\n", MP(b), (uint64_t)d, (uint64_t)q);
mpfree(b);
}
static void
testdiv(mpint *x, mpint *y)
{
mpint *b2, *b3;
int64_t now;
int i;
b2 = mpnew(0);
b3 = mpnew(0);
now = nsec();
for(i = 0; i < loops; i++)
mpdiv(x, y, b2, b3);
if(loops > 1)
printf("%"PRId64" µs for a %d/%d div\n", (nsec()-now)/(1000*loops),
x->top*Dbits, y->top*Dbits);
printf("%s / %s = %s %s\n", MP(x), MP(y), MP(b2), MP(b3));
mpfree(b2);
mpfree(b3);
}
static void
testmod(mpint *x, mpint *y)
{
mpint *r;
int64_t now;
int i;
r = mpnew(0);
now = nsec();
for(i = 0; i < loops; i++)
mpmod(x, y, r);
if(loops > 1)
printf("%"PRId64" µs for a %d/%d mod\n", (nsec()-now)/(1000*loops),
x->top*Dbits, y->top*Dbits);
printf("%s mod %s = %s\n", MP(x), MP(y), MP(r));
mpfree(r);
}
static void
testinvert(mpint *x, mpint *y)
{
mpint *r, *d1, *d2;
int64_t now;
int i;
r = mpnew(0);
d1 = mpnew(0);
d2 = mpnew(0);
now = nsec();
mpextendedgcd(x, y, r, d1, d2);
mpdiv(x, r, x, d1);
mpdiv(y, r, y, d1);
for(i = 0; i < loops; i++)
mpinvert(x, y, r);
if(loops > 1)
printf("%"PRId64" µs for a %d in %d invert\n", (nsec()-now)/(1000*loops),
x->top*Dbits, y->top*Dbits);
printf("%s**-1 mod %s = %s\n", MP(x), MP(y), MP(r));
mpmul(r, x, d1);
mpmod(d1, y, d2);
printf("%s*%s mod %s = %s\n", MP(x), MP(r), MP(y), MP(d2));
mpfree(r);
mpfree(d1);
mpfree(d2);
}
static void
testsub1(char *a, char *b)
{
mpint *b1, *b2, *b3;
b1 = strtomp(a, nil, 16, nil);
b2 = strtomp(b, nil, 16, nil);
b3 = mpnew(0);
mpsub(b1, b2, b3);
printf("%s - %s = %s\n", MP(b1), MP(b2), MP(b3));
}
static void
testmul1(char *a, char *b)
{
mpint *b1, *b2, *b3;
b1 = strtomp(a, nil, 16, nil);
b2 = strtomp(b, nil, 16, nil);
b3 = mpnew(0);
mpmul(b1, b2, b3);
printf("%s * %s = %s\n", MP(b1), MP(b2), MP(b3));
}
static void
testexp(char *base, char *exp, char *mod)
{
mpint *b, *e, *m, *res;
int i;
uint64_t now;
b = strtomp(base, nil, 16, nil);
e = strtomp(exp, nil, 16, nil);
res = mpnew(0);
if(mod != nil)
m = strtomp(mod, nil, 16, nil);
else
m = nil;
now = nsec();
for(i = 0; i < loops; i++)
mpexp(b, e, m, res);
if(loops > 1)
printf("%"PRIu64"µs for a %d to the %d bit exp\n", (nsec()-now)/(loops*1000),
b->top*Dbits, e->top*Dbits);
if(m != nil)
printf("%s ^ %s mod %s == %s\n", MP(b), MP(e), MP(m), MP(res));
else
printf("%s ^ %s == %s\n", MP(b), MP(e), MP(res));
mpfree(b);
mpfree(e);
mpfree(res);
if(m != nil)
mpfree(m);
}
static void
testgcd(void)
{
mpint *a, *b, *d, *x, *y, *t1, *t2;
int i;
uint64_t now, then;
uint64_t etime;
d = mpnew(0);
x = mpnew(0);
y = mpnew(0);
t1 = mpnew(0);
t2 = mpnew(0);
etime = 0;
a = strtomp("4EECAB3E04C4E6BC1F49D438731450396BF272B4D7B08F91C38E88ADCD281699889AFF872E2204C80CCAA8E460797103DE539D5DF8335A9B20C0B44886384F134C517287202FCA914D8A5096446B40CD861C641EF9C2730CB057D7B133F4C2B16BBD3D75FDDBD9151AAF0F9144AAA473AC93CF945DBFE0859FB685D5CBD0A8B3", nil, 16, nil);
b = strtomp("C41CFBE4D4846F67A3DF7DE9921A49D3B42DC33728427AB159CEC8CBBDB12B5F0C244F1A734AEB9840804EA3C25036AD1B61AFF3ABBC247CD4B384224567A863A6F020E7EE9795554BCD08ABAD7321AF27E1E92E3DB1C6E7E94FAAE590AE9C48F96D93D178E809401ABE8A534A1EC44359733475A36A70C7B425125062B1142D", nil, 16, nil);
mpextendedgcd(a, b, d, x, y);
printf("gcd %s*%s+%s*%s = %s?\n", MP(a), MP(x), MP(b), MP(y), MP(d));
mpfree(a);
mpfree(b);
for(i = 0; i < loops; i++){
a = mprand(2048, prng, nil);
b = mprand(2048, prng, nil);
then = nsec();
mpextendedgcd(a, b, d, x, y);
now = nsec();
etime += now-then;
mpmul(a, x, t1);
mpmul(b, y, t2);
mpadd(t1, t2, t2);
if(mpcmp(d, t2) != 0)
printf("%d gcd %s*%s+%s*%s != %s\n", i, MP(a), MP(x), MP(b), MP(y), MP(d));
// else
// printf("%d euclid %s*%s+%s*%s == %s\n", i, MP(a), MP(x), MP(b), MP(y), MP(d));
mpfree(a);
mpfree(b);
}
mpfree(x);
mpfree(y);
mpfree(d);
mpfree(t1);
mpfree(t2);
if(loops > 1)
printf("binary %"PRIu64"\n", etime);
}
int
main(int argc, char **argv)
{
mpint *x, *y;
if(argc == 3 && strcmp(argv[1], "-n") == 0)
loops = atoi(argv[2]);
D_PNAN = D_NNAN = strtod("+NaN", nil);
D_PINF = D_NINF = strtod("+Inf", nil);
union ieee754_double *d;
d = (union ieee754_double *)&D_NNAN;
d->ieee.negative = 1;
d = (union ieee754_double *)&D_NINF;
d->ieee.negative = 1;
srand(0);
mpsetminbits(2*Dbits);
testshift("1111111111111111");
testaddsub("fffffffffffffffff");
testdigdiv("1234567812345678", 0x76543218);
testdigdiv("1ffff", 0xffff);
testdigdiv("ffff", 0xffff);
testdigdiv("fff", 0xffff);
testdigdiv("effffffff", 0xffff);
testdigdiv("ffffffff", 0x1);
testdigdiv("ffffffff", 0);
testdigdiv("200000000", 2);
testdigdiv("ffffff00fffffff1", 0xfffffff1);
testvecdigmuladd("fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", 2);
testconv("0");
testconv("-abc0123456789abcedf");
testconv("abc0123456789abcedf");
testconv("ffffffff");
testconv("aaaaaaaaaaaaaaaaa");
testconv("1111111111111111");
testconv("33333333333333333333333333333333");
testvecdigmulsub("fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", 2);
testsub1("1FFFFFFFE00000000", "FFFFFFFE00000001");
testmul1("ffffffff", "f");
testmul("ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff");
testmul1("100000000000000000000000000000000000000000000000000000002000000000000000000000000000000000000000000000000000000030000000000000000000000000000000000000000000000000000000400000000000000000000000000000000000000000000000000000004FFFFFFFFFFFFFFFE0000000200000000000000000000000000000003FFFFFFFFFFFFFFFE0000000200000000000000000000000000000002FFFFFFFFFFFFFFFE0000000200000000000000000000000000000001FFFFFFFFFFFFFFFE0000000200000000000000000000000000000000FFFFFFFFFFFFFFFE0000000200000000FFFFFFFE00000001", "FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF");
testmul1("1000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001", "1000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001");
testmul1("1000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001", "1000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001000000000000000000000001");
testmul1("FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000", "FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000");
x = mprand(256, prng, nil);
y = mprand(128, prng, nil);
testdiv(x, y);
x = mprand(2048, prng, nil);
y = mprand(1024, prng, nil);
testdiv(x, y);
x = mprand(4*1024, prng, nil);
y = mprand(4*1024, prng, nil);
testmul2(x, y);
testsub1("677132C9", "-A26559B6");
testgcd();
x = mprand(512, prng, nil);
x->sign = -1;
y = mprand(256, prng, nil);
testdiv(x, y);
testmod(x, y);
x->sign = 1;
testinvert(y, x);
testexp("111111111", "222", "1000000000000000000000");
testexp("ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", "ffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffff", "100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000");
return 0;
}