version 1.3, 2000/08/22 05:04:05 |
version 1.7, 2020/10/04 03:14:08 |
|
|
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, |
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, |
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. |
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. |
* |
* |
* $OpenXM: OpenXM_contrib2/asir2000/engine/lmi.c,v 1.2 2000/08/21 08:31:27 noro Exp $ |
* $OpenXM: OpenXM_contrib2/asir2000/engine/lmi.c,v 1.6 2018/03/29 01:32:52 noro Exp $ |
*/ |
*/ |
#include "ca.h" |
#include "ca.h" |
#include "base.h" |
#include "base.h" |
#include "inline.h" |
#include "inline.h" |
|
|
typedef struct oGEN_LM { |
typedef struct oGEN_LM { |
int id; |
int id; |
N dense; |
N dense; |
N sparse; |
N sparse; |
int bit; |
int bit; |
unsigned int lower; |
unsigned int lower; |
} *GEN_LM; |
} *GEN_LM; |
|
|
void pwrlm0(N,N,N *); |
void pwrlm0(N,N,N *); |
void gen_simpn(N,N*); |
void gen_simpn(N,N*); |
void gen_simpn_force(N,N*); |
void gen_simpn_force(N,N*); |
|
void setmod_lf(N p); |
|
|
int lm_lazy; |
int lm_lazy; |
|
|
static GEN_LM current_mod_lm; |
static GEN_LM current_mod_lm; |
|
|
void random_lm(r) |
void random_lm(LM *r) |
LM *r; |
|
{ |
{ |
N n,t; |
N n,t; |
|
|
if ( !current_mod_lm ) |
if ( !current_mod_lm ) |
error("random_lm : current_mod_lm is not set"); |
error("random_lm : current_mod_lm is not set"); |
randomn(current_mod_lm->bit+1,&n); |
randomn(current_mod_lm->bit+1,&n); |
gen_simpn_force(n,&t); |
gen_simpn_force(n,&t); |
MKLM(t,*r); |
MKLM(t,*r); |
} |
} |
|
|
void ntosparsen(p,bits) |
void ntosparsen(N p,N *bits) |
N p; |
|
N *bits; |
|
{ |
{ |
int l,i,j,nz; |
int l,i,j,nz; |
N r; |
N r; |
unsigned int *pb,*rb; |
unsigned int *pb,*rb; |
if ( !p ) |
if ( !p ) |
*bits = 0; |
*bits = 0; |
else { |
else { |
l = n_bits(p); |
l = n_bits(p); |
pb = BD(p); |
pb = BD(p); |
for ( i = l-1, nz = 0; i >= 0; i-- ) |
for ( i = l-1, nz = 0; i >= 0; i-- ) |
if ( pb[i>>5]&(1<<i&31) ) nz++; |
if ( pb[i>>5]&(1<<i&31) ) nz++; |
*bits = r = NALLOC(nz); PL(r) = nz; |
*bits = r = NALLOC(nz); PL(r) = nz; |
rb = BD(r); |
rb = BD(r); |
for ( i = l-1, j = 0; i >= 0; i-- ) |
for ( i = l-1, j = 0; i >= 0; i-- ) |
if ( pb[i>>5]&(1<<i&31) ) |
if ( pb[i>>5]&(1<<i&31) ) |
rb[j++] = i; |
rb[j++] = i; |
} |
} |
} |
} |
|
|
void setmod_lm(p) |
void setmod_lm(N p) |
N p; |
|
{ |
{ |
int i; |
int i; |
|
Q q; |
|
|
if ( !current_mod_lm ) |
if ( !current_mod_lm ) |
current_mod_lm = (GEN_LM)MALLOC(sizeof(struct oGEN_LM)); |
current_mod_lm = (GEN_LM)MALLOC(sizeof(struct oGEN_LM)); |
bzero((char *)current_mod_lm,sizeof(struct oGEN_LM)); |
bzero((char *)current_mod_lm,sizeof(struct oGEN_LM)); |
current_mod_lm->dense = p; |
current_mod_lm->dense = p; |
ntosparsen(p,¤t_mod_lm->sparse); |
ntosparsen(p,¤t_mod_lm->sparse); |
current_mod_lm->bit = n_bits(p); |
current_mod_lm->bit = n_bits(p); |
current_mod_lm->id = UP2_DENSE; /* use ordinary rep. by default */ |
current_mod_lm->id = UP2_DENSE; /* use ordinary rep. by default */ |
if ( !(current_mod_lm->bit%32) ) { |
if ( !(current_mod_lm->bit%32) ) { |
for ( i = PL(p)-1; i >= 1; i-- ) |
for ( i = PL(p)-1; i >= 1; i-- ) |
if ( BD(p)[i] != 0xffffffff ) |
if ( BD(p)[i] != 0xffffffff ) |
break; |
break; |
if ( !i ) { |
if ( !i ) { |
current_mod_lm->lower = (unsigned int)((((UL)1)<<32)-((UL)BD(p)[0])); |
current_mod_lm->lower = (unsigned int)((((UL)1)<<32)-((UL)BD(p)[0])); |
current_mod_lm->id = UP2_SPARSE; /* XXX */ |
current_mod_lm->id = UP2_SPARSE; /* XXX */ |
} |
} |
} |
} |
|
setmod_lf(p); |
} |
} |
|
|
void getmod_lm(p) |
void getmod_lm(N *p) |
N *p; |
|
{ |
{ |
if ( !current_mod_lm ) |
if ( !current_mod_lm ) |
*p = 0; |
*p = 0; |
else |
else |
*p = current_mod_lm->dense; |
*p = current_mod_lm->dense; |
} |
} |
|
|
void simplm(n,r) |
void simplm(LM n,LM *r) |
LM n; |
|
LM *r; |
|
{ |
{ |
N rem; |
N rem; |
|
|
if ( !n ) |
if ( !n ) |
*r = 0; |
*r = 0; |
else if ( NID(n) != N_LM ) |
else if ( NID(n) != N_LM ) |
*r = n; |
*r = n; |
else { |
else { |
gen_simpn(n->body,&rem); |
gen_simpn(n->body,&rem); |
MKLM(rem,*r); |
MKLM(rem,*r); |
} |
} |
} |
} |
|
|
void qtolm(q,l) |
void simplm_force(LM n,LM *r) |
Q q; |
|
LM *l; |
|
{ |
{ |
N rn; |
N rem; |
LM a,b,c; |
|
|
|
if ( !q || (OID(q)==O_N && ((NID(q) == N_LM) || (NID(q) == N_GFPN))) ) { /* XXX */ |
if ( !n ) |
*l = (LM)q; |
*r = 0; |
} else if ( OID(q) == O_N && NID(q) == N_Q ) { |
else if ( NID(n) != N_LM ) |
gen_simpn(NM(q),&rn); MKLM(rn,a); |
*r = n; |
if ( SGN(q) < 0 ) { |
else { |
chsgnlm(a,&c); a = c; |
gen_simpn_force(n->body,&rem); |
} |
MKLM(rem,*r); |
if ( DN(q) ) { |
} |
gen_simpn_force(DN(q),&rn); MKLM(rn,b); |
|
if ( !b ) |
|
error("qtolm : invalid denominator"); |
|
divlm(a,b,l); |
|
} else |
|
*l = a; |
|
} else |
|
error("qtolm : invalid argument"); |
|
} |
} |
|
|
|
|
|
void qtolm(Q q,LM *l) |
|
{ |
|
N rn; |
|
LM a,b,c; |
|
|
|
if ( !q || (OID(q)==O_N && ((NID(q) == N_LM) || (NID(q) == N_GFPN))) ) { /* XXX */ |
|
*l = (LM)q; |
|
} else if ( OID(q) == O_N && NID(q) == N_Q ) { |
|
gen_simpn(NM(q),&rn); MKLM(rn,a); |
|
if ( SGN(q) < 0 ) { |
|
chsgnlm(a,&c); a = c; |
|
} |
|
if ( DN(q) ) { |
|
gen_simpn_force(DN(q),&rn); MKLM(rn,b); |
|
if ( !b ) |
|
error("qtolm : invalid denominator"); |
|
divlm(a,b,l); |
|
} else |
|
*l = a; |
|
} else |
|
error("qtolm : invalid argument"); |
|
} |
|
|
#define NZLM(a) ((a)&&(NID(a)==N_LM)) |
#define NZLM(a) ((a)&&(NID(a)==N_LM)) |
|
|
void addlm(a,b,c) |
void addlm(LM a,LM b,LM *c) |
LM a,b; |
|
LM *c; |
|
{ |
{ |
int s; |
N t,t1; |
N t,t1; |
LM z; |
LM z; |
|
|
|
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
if ( !a ) |
if ( !a ) |
*c = b; |
*c = b; |
else if ( !b ) |
else if ( !b ) |
*c = a; |
*c = a; |
else { |
else { |
addn(a->body,b->body,&t); |
addn(a->body,b->body,&t); |
gen_simpn(t,&t1); |
gen_simpn(t,&t1); |
MKLM(t1,*c); |
MKLM(t1,*c); |
} |
} |
} |
} |
|
|
void sublm(a,b,c) |
void sublm(LM a,LM b,LM *c) |
LM a,b; |
|
LM *c; |
|
{ |
{ |
int s; |
int s; |
N t,t1,mod; |
N t; |
LM z; |
LM z; |
|
|
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
if ( !b ) |
if ( !b ) |
*c = a; |
*c = a; |
else if ( !a ) |
else if ( !a ) |
chsgnlm(b,c); |
chsgnlm(b,c); |
else { |
else { |
s = subn(a->body,b->body,&t); |
s = subn(a->body,b->body,&t); |
if ( !s ) |
if ( !s ) |
*c = 0; |
*c = 0; |
else { |
else { |
MKLM(t,z); |
MKLM(t,z); |
if ( s < 0 ) |
if ( s < 0 ) |
chsgnlm(z,c); |
chsgnlm(z,c); |
else |
else |
*c = z; |
*c = z; |
} |
} |
} |
} |
} |
} |
|
|
void mullm(a,b,c) |
void mullm(LM a,LM b,LM *c) |
LM a,b; |
|
LM *c; |
|
{ |
{ |
N t,q,r; |
N t,r; |
LM z; |
LM z; |
|
|
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
if ( !a || !b ) |
if ( !a || !b ) |
*c = 0; |
*c = 0; |
else { |
else { |
kmuln(a->body,b->body,&t); |
kmuln(a->body,b->body,&t); |
gen_simpn(t,&r); |
gen_simpn(t,&r); |
MKLM(r,*c); |
MKLM(r,*c); |
} |
} |
} |
} |
|
|
void divlm(a,b,c) |
void divlm(LM a,LM b,LM *c) |
LM a,b; |
|
LM *c; |
|
{ |
{ |
LM r; |
LM r; |
Q t,m,i; |
Q t,m,i; |
LM z; |
LM z; |
|
|
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
if ( !b ) |
if ( !b ) |
error("divlm: division by 0"); |
error("divlm: division by 0"); |
else if ( !a ) |
else if ( !a ) |
*c = 0; |
*c = 0; |
else { |
else { |
NTOQ(b->body,1,t); NTOQ(current_mod_lm->dense,1,m); |
NTOQ(b->body,1,t); NTOQ(current_mod_lm->dense,1,m); |
invl(t,m,&i); |
invl(t,m,&i); |
MKLM(NM(i),r); |
MKLM(NM(i),r); |
mullm(a,r,c); |
mullm(a,r,c); |
} |
} |
} |
} |
|
|
void chsgnlm(a,c) |
void chsgnlm(LM a,LM *c) |
LM a,*c; |
|
{ |
{ |
LM t; |
LM t; |
N s,u; |
N s,u; |
|
|
if ( !a ) |
if ( !a ) |
*c = a; |
*c = a; |
else { |
else { |
qtolm((Q)a,&t); a = t; |
qtolm((Q)a,&t); a = t; |
gen_simpn_force(a->body,&s); |
gen_simpn_force(a->body,&s); |
subn(current_mod_lm->dense,s,&u); |
subn(current_mod_lm->dense,s,&u); |
MKLM(u,*c); |
MKLM(u,*c); |
} |
} |
} |
} |
|
|
void pwrlm(a,b,c) |
void pwrlm(LM a,Q b,LM *c) |
LM a; |
|
Q b; |
|
LM *c; |
|
{ |
{ |
LM t; |
LM t; |
N s; |
N s; |
|
|
if ( !b ) { |
if ( !b ) { |
MKLM(ONEN,*c); |
MKLM(ONEN,*c); |
} else if ( !a ) |
} else if ( !a ) |
*c = 0; |
*c = 0; |
else { |
else { |
qtolm((Q)a,&t); a = t; |
qtolm((Q)a,&t); a = t; |
pwrlm0(a->body,NM(b),&s); |
pwrlm0(a->body,NM(b),&s); |
MKLM(s,*c); |
MKLM(s,*c); |
} |
} |
} |
} |
|
|
void pwrlm0(a,n,c) |
void pwrlm0(N a,N n,N *c) |
N a,n; |
|
N *c; |
|
{ |
{ |
N n1,t,t1,t2,r1; |
N n1,t,t1,t2,r1; |
int r; |
int r; |
|
|
if ( !n ) |
if ( !n ) |
*c = ONEN; |
*c = ONEN; |
else if ( UNIN(n) ) |
else if ( UNIN(n) ) |
*c = a; |
*c = a; |
else { |
else { |
r = divin(n,2,&n1); |
r = divin(n,2,&n1); |
pwrlm0(a,n1,&t); |
pwrlm0(a,n1,&t); |
kmuln(t,t,&t1); gen_simpn(t1,&r1); |
kmuln(t,t,&t1); gen_simpn(t1,&r1); |
if ( r ) { |
if ( r ) { |
kmuln(r1,a,&t2); gen_simpn(t2,&r1); |
kmuln(r1,a,&t2); gen_simpn(t2,&r1); |
} |
} |
*c = r1; |
*c = r1; |
} |
} |
} |
} |
|
|
int cmplm(a,b) |
int cmplm(LM a,LM b) |
LM a,b; |
|
{ |
{ |
LM z; |
LM z; |
|
|
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
qtolm((Q)a,&z); a = z; qtolm((Q)b,&z); b = z; |
if ( !a ) |
if ( !a ) |
if ( !b ) |
if ( !b ) |
return 0; |
return 0; |
else |
else |
return -1; |
return -1; |
else if ( !b ) |
else if ( !b ) |
return 1; |
return 1; |
else |
else |
return cmpn(a->body,b->body); |
return cmpn(a->body,b->body); |
} |
} |
|
|
void remn_special(N,N,int,unsigned int ,N *); |
void remn_special(N,N,int,unsigned int ,N *); |
|
|
void gen_simpn(a,b) |
void gen_simpn(N a,N *b) |
N a; |
|
N *b; |
|
{ |
{ |
if ( !current_mod_lm ) |
if ( !current_mod_lm ) |
error("gen_simpn: current_mod_lm is not set"); |
error("gen_simpn: current_mod_lm is not set"); |
if ( lm_lazy ) |
if ( lm_lazy ) |
*b = a; |
*b = a; |
else if ( current_mod_lm->id == UP2_SPARSE ) |
else if ( current_mod_lm->id == UP2_SPARSE ) |
remn_special(a,current_mod_lm->dense,current_mod_lm->bit,current_mod_lm->lower,b); |
remn_special(a,current_mod_lm->dense,current_mod_lm->bit,current_mod_lm->lower,b); |
else |
else |
remn(a,current_mod_lm->dense,b); |
remn(a,current_mod_lm->dense,b); |
} |
} |
|
|
void gen_simpn_force(a,b) |
void gen_simpn_force(N a,N *b) |
N a; |
|
N *b; |
|
{ |
{ |
if ( !current_mod_lm ) |
if ( !current_mod_lm ) |
error("gen_simpn_force: current_mod_lm is not set"); |
error("gen_simpn_force: current_mod_lm is not set"); |
else if ( current_mod_lm->id == UP2_SPARSE ) |
else if ( current_mod_lm->id == UP2_SPARSE ) |
remn_special(a,current_mod_lm->dense,current_mod_lm->bit,current_mod_lm->lower,b); |
remn_special(a,current_mod_lm->dense,current_mod_lm->bit,current_mod_lm->lower,b); |
else |
else |
remn(a,current_mod_lm->dense,b); |
remn(a,current_mod_lm->dense,b); |
} |
} |
|
|