version 1.17, 2020/01/09 01:47:41 |
version 1.21, 2021/10/08 07:27:05 |
|
|
/* $OpenXM: OpenXM_contrib2/asir2018/engine/Q.c,v 1.16 2019/05/22 02:41:16 noro Exp $ */ |
/* $OpenXM: OpenXM_contrib2/asir2018/engine/Q.c,v 1.20 2020/12/04 08:09:33 noro Exp $ */ |
#include "ca.h" |
#include "ca.h" |
#include "gmp.h" |
#include "gmp.h" |
#include "base.h" |
#include "base.h" |
Line 17 extern int DP_Print; |
|
Line 17 extern int DP_Print; |
|
|
|
void isqrtz(Z a,Z *r); |
void isqrtz(Z a,Z *r); |
void bshiftz(Z a,int n,Z *r); |
void bshiftz(Z a,int n,Z *r); |
|
int mpz_inttorat(mpz_t c,mpz_t m,mpz_t b,mpz_t nm,mpz_t dn); |
|
int generic_gauss_elim_hensel64(MAT mat,MAT *nmmat,Z *dn,int **rindp,int **cindp,DP *mb); |
|
int find_lhs_and_lu_mod64(mp_limb_t **a,int row,int col,mp_limb_t md,int **rinfo,int **cinfo); |
|
void solve_by_lu_mod64(mp_limb_t **a,int n,mp_limb_t md,mp_limb_signed_t **b,int l,int normalize); |
|
|
void *gc_realloc(void *p,size_t osize,size_t nsize) |
void *gc_realloc(void *p,size_t osize,size_t nsize) |
{ |
{ |
Line 93 void dnq(Q q,Z *r) |
|
Line 97 void dnq(Q q,Z *r) |
|
int sgnq(Q q) |
int sgnq(Q q) |
{ |
{ |
if ( !q ) return 0; |
if ( !q ) return 0; |
else if ( q->z ) return mpz_sgn(BDY((Z)q)); |
else if ( q->z ) { |
else return mpz_sgn(mpq_numref(BDY(q))); |
int sgn; |
|
sgn = mpz_sgn(BDY((Z)q)); |
|
return sgn; |
|
} else return mpz_sgn(mpq_numref(BDY(q))); |
} |
} |
|
|
Q mpqtozq(mpq_t a) |
Q mpqtozq(mpq_t a) |
Line 388 int cmpz(Z q1,Z q2) |
|
Line 395 int cmpz(Z q1,Z q2) |
|
else if ( sgn < 0 ) return -1; |
else if ( sgn < 0 ) return -1; |
else return 0; |
else return 0; |
} |
} |
|
/* XXX */ |
|
return 0; |
} |
} |
|
|
void gcdz(Z n1,Z n2,Z *nq) |
void gcdz(Z n1,Z n2,Z *nq) |
Line 665 void pwrq(Q n1,Q n,Q *nr) |
|
Line 674 void pwrq(Q n1,Q n,Q *nr) |
|
nm[0] = mpq_denref(BDY(n1))[0]; |
nm[0] = mpq_denref(BDY(n1))[0]; |
dn[0] = mpq_numref(BDY(n1))[0]; |
dn[0] = mpq_numref(BDY(n1))[0]; |
} |
} |
|
if ( mpz_sgn(dn)<0 ) { |
|
mpz_neg(nm,nm); |
|
mpz_neg(dn,dn); |
|
} |
} else { |
} else { |
if ( n1->z ) { |
if ( n1->z ) { |
nm[0] = BDY((Z)n1)[0]; |
nm[0] = BDY((Z)n1)[0]; |
Line 783 void lgp(P p,Z *g,Z *l) |
|
Line 796 void lgp(P p,Z *g,Z *l) |
|
|
|
if ( NUM(p) ) { |
if ( NUM(p) ) { |
if ( ((Q)p)->z ) { |
if ( ((Q)p)->z ) { |
MPZTOZ(BDY((Z)p),*g); |
absz((Z)p,g); |
*l = ONE; |
*l = ONE; |
} else { |
} else { |
MPZTOZ(mpq_numref(BDY((Q)p)),*g); |
MPZTOZ(mpq_numref(BDY((Q)p)),g1); absz(g1,g); |
MPZTOZ(mpq_denref(BDY((Q)p)),*l); |
MPZTOZ(mpq_denref(BDY((Q)p)),*l); |
} |
} |
} else { |
} else { |
Line 1232 int intmtoratm(MAT mat,Z md,MAT nm,Z *dn) |
|
Line 1245 int intmtoratm(MAT mat,Z md,MAT nm,Z *dn) |
|
return 0; |
return 0; |
row = mat->row; col = mat->col; |
row = mat->row; col = mat->col; |
bshiftz(md,1,&t); |
bshiftz(md,1,&t); |
isqrt(t,&s); |
isqrtz(t,&s); |
bshiftz(s,64,&b); |
bshiftz(s,64,&b); |
if ( !b ) b = ONE; |
if ( !b ) b = ONE; |
dn0 = ONE; |
dn0 = ONE; |
Line 2263 int generic_gauss_elim_hensel64(MAT mat,MAT *nmmat,Z * |
|
Line 2276 int generic_gauss_elim_hensel64(MAT mat,MAT *nmmat,Z * |
|
for ( j = 0, bi = b[i], wi = wc[i]; j < ri; j++ ) |
for ( j = 0, bi = b[i], wi = wc[i]; j < ri; j++ ) |
wi[j] = mpz_fdiv_ui(bi[j],md); |
wi[j] = mpz_fdiv_ui(bi[j],md); |
/* wc = A^(-1)wc; wc is not normalized */ |
/* wc = A^(-1)wc; wc is not normalized */ |
solve_by_lu_mod64(w,rank,md,wc,ri,0); |
solve_by_lu_mod64(w,rank,md,(mp_limb_signed_t **)wc,ri,0); |
/* x += q*wc */ |
/* x += q*wc */ |
for ( i = 0; i < rank; i++ ) |
for ( i = 0; i < rank; i++ ) |
for ( j = 0, wi = wc[i]; j < ri; j++ ) |
for ( j = 0, wi = wc[i]; j < ri; j++ ) |