version 1.13, 2019/03/18 07:00:33 |
version 1.18, 2020/10/04 03:14:09 |
|
|
/* $OpenXM: OpenXM_contrib2/asir2018/engine/Q.c,v 1.12 2019/03/17 02:48:13 noro Exp $ */ |
/* $OpenXM: OpenXM_contrib2/asir2018/engine/Q.c,v 1.17 2020/01/09 01:47:41 noro Exp $ */ |
#include "ca.h" |
#include "ca.h" |
#include "gmp.h" |
#include "gmp.h" |
#include "base.h" |
#include "base.h" |
#include "inline.h" |
#include "inline.h" |
|
|
mpz_t ONEMPZ; |
mpz_t ONEMPZ; |
Z ONE; |
extern Z ONE; |
int lf_lazy; |
int lf_lazy; |
Z current_mod_lf; |
Z current_mod_lf; |
int current_mod_lf_size; |
int current_mod_lf_size; |
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 558 void subq(Q n1,Q n2,Q *nr) |
|
Line 562 void subq(Q n1,Q n2,Q *nr) |
|
|
|
if ( !n1 ) { |
if ( !n1 ) { |
if ( !n2 ) *nr = 0; |
if ( !n2 ) *nr = 0; |
else if ( n1->z ) chsgnz((Z)n1,(Z *)nr); |
else if ( n2->z ) chsgnz((Z)n2,(Z *)nr); |
else { |
else { |
mpq_init(t); mpq_neg(t,BDY(n2)); MPQTOQ(t,*nr); |
mpq_init(t); mpq_neg(t,BDY(n2)); MPQTOQ(t,*nr); |
} |
} |
Line 762 void lgp(P p,Z *g,Z *l); |
|
Line 766 void lgp(P p,Z *g,Z *l); |
|
|
|
void ptozp(P p,int sgn,Q *c,P *pr) |
void ptozp(P p,int sgn,Q *c,P *pr) |
{ |
{ |
Z nm,dn; |
Z nm,dn,nm1; |
|
|
if ( !p ) { |
if ( !p ) { |
*c = 0; *pr = 0; |
*c = 0; *pr = 0; |
} else { |
} else { |
lgp(p,&nm,&dn); |
lgp(p,&nm,&dn); |
|
if ( sgn < 0 ) { |
|
chsgnz(nm,&nm1); nm = nm1; |
|
} |
divz(nm,dn,(Z *)c); |
divz(nm,dn,(Z *)c); |
divsp(CO,p,(P)*c,pr); |
divsp(CO,p,(P)*c,pr); |
} |
} |
Line 1229 int intmtoratm(MAT mat,Z md,MAT nm,Z *dn) |
|
Line 1236 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 1337 int mpz_inttorat(mpz_t c,mpz_t m,mpz_t b,mpz_t nm,mpz_ |
|
Line 1344 int mpz_inttorat(mpz_t c,mpz_t m,mpz_t b,mpz_t nm,mpz_ |
|
} |
} |
if ( mpz_cmp(v1,b) >= 0 ) return 0; |
if ( mpz_cmp(v1,b) >= 0 ) return 0; |
else { |
else { |
|
mpz_gcd(t,v1,v2); |
|
if ( UNIMPZ(t) ) |
|
mpz_set_ui(r1,0); |
|
else { |
|
/* v1 /= t, v2 /= t, t=c*v1-v2, r1=t%m */ |
|
mpz_divexact(v1,v1,t); mpz_divexact(v2,v2,t); |
|
mpz_mul(t,c,v1); mpz_sub(t,t,v2); mpz_mod(r1,t,m); |
|
} |
|
if ( mpz_sgn(r1) ) return 0; |
if ( mpz_sgn(v1)<0 ) { |
if ( mpz_sgn(v1)<0 ) { |
mpz_neg(dn,v1); mpz_neg(nm,v2); |
mpz_neg(dn,v1); mpz_neg(nm,v2); |
} else { |
} else { |
Line 1348 int mpz_inttorat(mpz_t c,mpz_t m,mpz_t b,mpz_t nm,mpz_ |
|
Line 1364 int mpz_inttorat(mpz_t c,mpz_t m,mpz_t b,mpz_t nm,mpz_ |
|
|
|
int inttorat(Z c,Z m,Z b,Z *nmp,Z *dnp) |
int inttorat(Z c,Z m,Z b,Z *nmp,Z *dnp) |
{ |
{ |
Z qq,t,u1,v1,r1; |
Z qq,t,s,r,u1,v1,r1; |
Z q,u2,v2,r2; |
Z q,u2,v2,r2; |
|
|
u1 = 0; v1 = ONE; u2 = m; v2 = c; |
u1 = 0; v1 = ONE; u2 = m; v2 = c; |
Line 1358 int inttorat(Z c,Z m,Z b,Z *nmp,Z *dnp) |
|
Line 1374 int inttorat(Z c,Z m,Z b,Z *nmp,Z *dnp) |
|
} |
} |
if ( cmpz(v1,b) >= 0 ) return 0; |
if ( cmpz(v1,b) >= 0 ) return 0; |
else { |
else { |
if ( mpz_sgn(BDY(v1))<0 ) { |
/* reduction and check */ |
chsgnz(v1,dnp); chsgnz(v2,nmp); |
/* v2/v1 = u2/u1, c*u1-u2 = 0 mod m? */ |
|
gcdz(v1,v2,&t); |
|
if ( UNIZ(t) ) { |
|
u1 = v1; u2 = v2; r = 0; |
} else { |
} else { |
*dnp = v1; *nmp = v2; |
divsz(v1,t,&u1); divsz(v2,t,&u2); |
|
mulz(c,u1,&t); subz(t,u2,&s); remz(s,m,&r); |
|
} |
|
if ( r ) return 0; |
|
if ( mpz_sgn(BDY(u1))<0 ) { |
|
chsgnz(u1,dnp); chsgnz(u2,nmp); |
|
} else { |
|
*dnp = u1; *nmp = u2; |
} |
} |
return 1; |
return 1; |
} |
} |