version 1.3, 2018/09/24 22:26:43 |
version 1.5, 2018/09/28 08:20:28 |
Line 353 void pwrz(Z n1,Z n,Z *nr) |
|
Line 353 void pwrz(Z n1,Z n,Z *nr) |
|
} else if ( !smallz(n) ) { |
} else if ( !smallz(n) ) { |
error("exponent too big."); *nr = 0; |
error("exponent too big."); *nr = 0; |
} else if ( n1->z && mpz_sgn(BDY((Z)n))>0 ) { |
} else if ( n1->z && mpz_sgn(BDY((Z)n))>0 ) { |
mpz_init(z); mpz_pow_ui(z,BDY(n1),QTOS(n)); MPZTOZ(z,*nr); |
mpz_init(z); mpz_pow_ui(z,BDY(n1),ZTOS(n)); MPZTOZ(z,*nr); |
} else { |
} else { |
MPZTOMPQ(BDY(n1),q); MPQTOQ(q,r); |
MPZTOMPQ(BDY(n1),q); MPQTOQ(q,r); |
pwrq(r,(Q)n,&p); *nr = (Z)p; |
pwrq(r,(Q)n,&p); *nr = (Z)p; |
Line 453 void gcdvz_estimate(VECT v,Z *q) |
|
Line 453 void gcdvz_estimate(VECT v,Z *q) |
|
else addz(s,b[i],&u); |
else addz(s,b[i],&u); |
s = u; |
s = u; |
} |
} |
for ( i = 0, t = 0; i < n; i++ ) { |
for ( t = 0; i < n; i++ ) { |
if ( b[i] && mpz_sgn(BDY(b[i]))<0 ) subz(t,b[i],&u); |
if ( b[i] && mpz_sgn(BDY(b[i]))<0 ) subz(t,b[i],&u); |
else addz(t,b[i],&u); |
else addz(t,b[i],&u); |
t = u; |
t = u; |
Line 461 void gcdvz_estimate(VECT v,Z *q) |
|
Line 461 void gcdvz_estimate(VECT v,Z *q) |
|
gcdz(s,t,q); |
gcdz(s,t,q); |
} |
} |
|
|
|
void gcdv_mpz_estimate(mpz_t g,mpz_t *b,int n) |
|
{ |
|
int m,m2,i,j; |
|
mpz_t s,t; |
|
|
|
mpz_init(g); |
|
for ( i = 0, m = 0; i < n; i++ ) |
|
if ( mpz_sgn(b[i]) ) m++; |
|
if ( !m ) { |
|
mpz_set_ui(g,0); |
|
return; |
|
} |
|
if ( m == 1 ) { |
|
for ( i = 0, m = 0; i < n; i++ ) |
|
if ( mpz_sgn(b[i]) ) break; |
|
if ( mpz_sgn(b[i])<0 ) mpz_neg(g,b[i]); |
|
else mpz_set(g,b[i]); |
|
return ; |
|
} |
|
m2 = m/2; |
|
mpz_init_set_ui(s,0); |
|
for ( i = j = 0; j < m2; i++ ) { |
|
if ( mpz_sgn(b[i]) ) { |
|
if ( mpz_sgn(b[i])<0 ) |
|
mpz_sub(s,s,b[i]); |
|
else |
|
mpz_add(s,s,b[i]); |
|
j++; |
|
} |
|
} |
|
mpz_init_set_ui(t,0); |
|
for ( ; i < n; i++ ) { |
|
if ( mpz_sgn(b[i]) ) { |
|
if ( mpz_sgn(b[i])<0 ) |
|
mpz_sub(t,t,b[i]); |
|
else |
|
mpz_add(t,t,b[i]); |
|
} |
|
} |
|
mpz_gcd(g,s,t); |
|
} |
|
|
|
|
void factorialz(unsigned int n,Z *nr) |
void factorialz(unsigned int n,Z *nr) |
{ |
{ |
mpz_t a; |
mpz_t a; |
Line 605 void pwrq(Q n1,Q n,Q *nr) |
|
Line 648 void pwrq(Q n1,Q n,Q *nr) |
|
} else if ( !smallz((Z)n) ) { |
} else if ( !smallz((Z)n) ) { |
error("exponent too big."); *nr = 0; |
error("exponent too big."); *nr = 0; |
} else { |
} else { |
e = QTOS(n); |
e = ZTOS(n); |
if ( e < 0 ) { |
if ( e < 0 ) { |
e = -e; |
e = -e; |
if ( n1->z ) { |
if ( n1->z ) { |
Line 662 void mkbc(int n,Z *t) |
|
Line 705 void mkbc(int n,Z *t) |
|
Z c,d,iq; |
Z c,d,iq; |
|
|
for ( t[0] = ONE, i = 1; i <= n/2; i++ ) { |
for ( t[0] = ONE, i = 1; i <= n/2; i++ ) { |
STOQ(n-i+1,c); mulz(t[i-1],c,&d); |
STOZ(n-i+1,c); mulz(t[i-1],c,&d); |
STOQ(i,iq); divsz(d,iq,&t[i]); |
STOZ(i,iq); divsz(d,iq,&t[i]); |
} |
} |
for ( ; i <= n; i++ ) |
for ( ; i <= n; i++ ) |
t[i] = t[n-i]; |
t[i] = t[n-i]; |
Line 1326 void isqrtz(Z a,Z *r) |
|
Line 1369 void isqrtz(Z a,Z *r) |
|
else { |
else { |
k = z_bits((Q)a); /* a <= 2^k-1 */ |
k = z_bits((Q)a); /* a <= 2^k-1 */ |
bshiftz(ONE,-((k>>1)+(k&1)),&x); /* a <= x^2 */ |
bshiftz(ONE,-((k>>1)+(k&1)),&x); /* a <= x^2 */ |
STOQ(2,two); |
STOZ(2,two); |
while ( 1 ) { |
while ( 1 ) { |
pwrz(x,two,&t); |
pwrz(x,two,&t); |
if ( cmpz(t,a) <= 0 ) { |
if ( cmpz(t,a) <= 0 ) { |
Line 1450 init_eg(&eg_mul1); init_eg(&eg_mul2); |
|
Line 1493 init_eg(&eg_mul1); init_eg(&eg_mul2); |
|
w = (int **)almat(row,col); |
w = (int **)almat(row,col); |
for ( ind = 0; ; ind++ ) { |
for ( ind = 0; ; ind++ ) { |
md = get_lprime(ind); |
md = get_lprime(ind); |
STOQ(md,mdq); |
STOZ(md,mdq); |
for ( i = 0; i < row; i++ ) |
for ( i = 0; i < row; i++ ) |
for ( j = 0, ai = a0[i], wi = w[i]; j < col; j++ ) |
for ( j = 0, ai = a0[i], wi = w[i]; j < col; j++ ) |
wi[j] = remqi((Q)ai[j],md); |
wi[j] = remqi((Q)ai[j],md); |
Line 1598 int generic_gauss_elim_hensel_dalg(MAT mat,DP *mb,MAT |
|
Line 1641 int generic_gauss_elim_hensel_dalg(MAT mat,DP *mb,MAT |
|
w = (int **)almat(row,col); |
w = (int **)almat(row,col); |
for ( ind = 0; ; ind++ ) { |
for ( ind = 0; ; ind++ ) { |
md = get_lprime(ind); |
md = get_lprime(ind); |
STOQ(md,mdq); |
STOZ(md,mdq); |
for ( i = 0; i < row; i++ ) |
for ( i = 0; i < row; i++ ) |
for ( j = 0, ai = a0[i], wi = w[i]; j < col; j++ ) |
for ( j = 0, ai = a0[i], wi = w[i]; j < col; j++ ) |
wi[j] = remqi((Q)ai[j],md); |
wi[j] = remqi((Q)ai[j],md); |
Line 1624 int generic_gauss_elim_hensel_dalg(MAT mat,DP *mb,MAT |
|
Line 1667 int generic_gauss_elim_hensel_dalg(MAT mat,DP *mb,MAT |
|
a = (Z **)almat_pointer(rank,rank); /* lhs mat */ |
a = (Z **)almat_pointer(rank,rank); /* lhs mat */ |
MKMAT(bmat,rank,col-rank); b = (Z **)bmat->body; /* lhs mat */ |
MKMAT(bmat,rank,col-rank); b = (Z **)bmat->body; /* lhs mat */ |
for ( j = li = ri = 0; j < col; j++ ) |
for ( j = li = ri = 0; j < col; j++ ) |
if ( cinfo[j] ) { |
if ( cinfo[j] > 0 ) { |
/* the column is in lhs */ |
/* the column is in lhs */ |
for ( i = 0; i < rank; i++ ) { |
for ( i = 0; i < rank; i++ ) { |
w[i][li] = w[i][j]; |
w[i][li] = w[i][j]; |
a[i][li] = a0[rinfo[i]][j]; |
a[i][li] = a0[rinfo[i]][j]; |
} |
} |
li++; |
li++; |
} else { |
} else if ( !cinfo[j] ) { |
/* the column is in rhs */ |
/* the column is in rhs */ |
for ( i = 0; i < rank; i++ ) |
for ( i = 0; i < rank; i++ ) |
b[i][ri] = a0[rinfo[i]][j]; |
b[i][ri] = a0[rinfo[i]][j]; |