OpenXM_contrib2/asir2018/engine/nd.c - diff

Return to nd.c CVS log

Up to [local] / OpenXM_contrib2 / asir2018 / engine

Diff for /OpenXM_contrib2/asir2018/engine/nd.c between version 1.1 and 1.20

version 1.1, 2018/09/19 05:45:07

version 1.20, 2019/09/15 08:46:12

Line 1

/* $OpenXM$ */

/* $OpenXM: OpenXM_contrib2/asir2018/engine/nd.c,v 1.19 2019/09/04 05:32:10 noro Exp $ */

#include "nd.h"

struct oEGT eg_search;

int Nnd_add,Nf4_red;

struct oEGT eg_search,f4_symb,f4_conv,f4_elim1,f4_elim2;

int diag_period = 6;

int weight_check = 1;

int (*ndl_compare_function)(UINT *a1,UINT *a2);

/* for schreyer order */

int (*ndl_base_compare_function)(UINT *a1,UINT *a2);

int nd_dcomp;

int nd_rref2;

NM _nm_free_list;

Line 77 NDV pltondv(VL vl,VL dvl,LIST p);

Line 80 NDV pltondv(VL vl,VL dvl,LIST p);

void pltozpl(LIST l,Q *cont,LIST *pp);

void ndl_max(UINT *d1,unsigned *d2,UINT *d);

void nmtodp(int mod,NM m,DP *r);

void ndltodp(UINT *d,DP *r);

NODE reverse_node(NODE n);

P ndc_div(int mod,union oNDC a,union oNDC b);

P ndctop(int mod,union oNDC c);

Line 86 void parse_nd_option(NODE opt);

Line 90 void parse_nd_option(NODE opt);

void dltondl(int n,DL dl,UINT *r);

DP ndvtodp(int mod,NDV p);

DP ndtodp(int mod,ND p);

DPM ndvtodpm(int mod,NDV p);

NDV dpmtondv(int mod,DPM p);

int dpm_getdeg(DPM p,int *rank);

void dpm_ptozp(DPM p,Z *cont,DPM *r);

int compdmm(int nv,DMM a,DMM b);

void Pdp_set_weight(NODE,VECT *);

void Pox_cmo_rpc(NODE,Obj *);

Line 469 int ndl_weight(UINT *d)

Line 478 int ndl_weight(UINT *d)

for ( j = 0; j < nd_epw; j++, u>>=nd_bpe )

t += (u&nd_mask0);

}

if ( nd_module && current_module_weight_vector && MPOS(d) )

if ( nd_module && nd_module_rank && MPOS(d) )

t += current_module_weight_vector[MPOS(d)];

t += nd_module_weight[MPOS(d)-1];

for ( i = nd_exporigin; i < nd_wpd; i++ )

if ( d[i] && !t )

printf("afo\n");

return t;

}

Line 485 int ndl_weight2(UINT *d)

Line 497 int ndl_weight2(UINT *d)

u = GET_EXP(d,i);

t += nd_sugarweight[i]*u;

}

if ( nd_module && current_module_weight_vector && MPOS(d) )

if ( nd_module && nd_module_rank && MPOS(d) )

t += current_module_weight_vector[MPOS(d)];

t += nd_module_weight[MPOS(d)-1];

return t;

}

Line 566 int ndl_matrix_compare(UINT *d1,UINT *d2)

Line 578 int ndl_matrix_compare(UINT *d1,UINT *d2)

Z t1;

Z t,t2;

for ( j = 0; j < nd_nvar; j++ )

nd_work_vector[j] = GET_EXP(d1,j)-GET_EXP(d2,j);

if ( nd_top_weight ) {

if ( OID(nd_top_weight) == O_VECT ) {

mat = (Z **)&BDY((VECT)nd_top_weight);

row = 1;

} else {

mat = (Z **)BDY((MAT)nd_top_weight);

row = ((MAT)nd_top_weight)->row;

}

for ( i = 0; i < row; i++ ) {

w = mat[i];

for ( j = 0, t = 0; j < nd_nvar; j++ ) {

STOQ(nd_work_vector[j],t1);

STOZ(nd_work_vector[j],t1);

mulz(w[j],t1,&t2);

addz(t,t2,&t1);

t = t1;

}

if ( t ) {

s = sgnz(t);

if ( s > 0 ) return 1;

else if ( s < 0 ) return -1;

}

for ( i = 0; i < nd_matrix_len; i++ ) {

v = nd_matrix[i];

for ( j = 0, s = 0; j < nd_nvar; j++ )

s += v[j]*nd_work_vector[j];

if ( s > 0 ) return 1;

else if ( s < 0 ) return -1;

}

for ( i = 0; i < nd_matrix_len; i++ ) {

v = nd_matrix[i];

for ( j = 0, s = 0; j < nd_nvar; j++ )

s += v[j]*nd_work_vector[j];

if ( s > 0 ) return 1;

else if ( s < 0 ) return -1;

}

if ( !ndl_equal(d1,d2) )

error("afo");

error("ndl_matrix_compare : invalid matrix");

return 0;

}

int ndl_composite_compare(UINT *d1,UINT *d2)

Line 704 int ndl_module_grlex_compare(UINT *d1,UINT *d2)

Line 716 int ndl_module_grlex_compare(UINT *d1,UINT *d2)

{

int i,c;

if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

// if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

if ( nd_ispot ) {

if ( nd_pot_nelim && MPOS(d1)>=nd_pot_nelim+1 && MPOS(d2) >= nd_pot_nelim+1 ) {

if ( TD(d1) > TD(d2) ) return 1;

else if ( TD(d1) < TD(d2) ) return -1;

if ( (c = ndl_lex_compare(d1,d2)) != 0 ) return c;

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

return 0;

}

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

}

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

}

if ( TD(d1) > TD(d2) ) return 1;

else if ( TD(d1) < TD(d2) ) return -1;

if ( (c = ndl_lex_compare(d1,d2)) != 0 ) return c;

Line 731 int ndl_module_glex_compare(UINT *d1,UINT *d2)

Line 743 int ndl_module_glex_compare(UINT *d1,UINT *d2)

{

int i,c;

if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

// if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

if ( nd_ispot ) {

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

Line 750 int ndl_module_lex_compare(UINT *d1,UINT *d2)

Line 762 int ndl_module_lex_compare(UINT *d1,UINT *d2)

{

int i,c;

if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

// if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

if ( nd_ispot ) {

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

Line 767 int ndl_module_block_compare(UINT *d1,UINT *d2)

Line 779 int ndl_module_block_compare(UINT *d1,UINT *d2)

{

int i,c;

if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

// if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

if ( nd_ispot ) {

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

Line 784 int ndl_module_matrix_compare(UINT *d1,UINT *d2)

Line 796 int ndl_module_matrix_compare(UINT *d1,UINT *d2)

{

int i,c;

if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

// if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

if ( nd_ispot ) {

if ( MPOS(d1) < MPOS(d2) ) return 1;

else if ( MPOS(d1) > MPOS(d2) ) return -1;

Line 801 int ndl_module_composite_compare(UINT *d1,UINT *d2)

Line 813 int ndl_module_composite_compare(UINT *d1,UINT *d2)

{

int i,c;

if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

// if ( nd_module_rank && (c = ndl_module_weight_compare(d1,d2)) ) return c;

if ( nd_ispot ) {

if ( MPOS(d1) > MPOS(d2) ) return 1;

else if ( MPOS(d1) < MPOS(d2) ) return -1;

Line 1125 int ndl_check_bound2(int index,UINT *d2)

Line 1137 int ndl_check_bound2(int index,UINT *d2)

INLINE int ndl_hash_value(UINT *d)

{

int i;

int r;

UINT r;

r = 0;

for ( i = 0; i < nd_wpd; i++ )

r = ((r<<16)+d[i])%REDTAB_LEN;

r = (r*1511+d[i]);

r %= REDTAB_LEN;

return r;

}

Line 1216 ND nd_add(int mod,ND p1,ND p2)

Line 1229 ND nd_add(int mod,ND p1,ND p2)

ND r;

NM m1,m2,mr0,mr,s;

Nnd_add++;

if ( !p1 ) return p2;

else if ( !p2 ) return p1;

else if ( mod == -1 ) return nd_add_sf(p1,p2);

Line 1412 ND nd_reduce2(int mod,ND d,ND g,NDV p,NM mul,NDC dn,Ob

Line 1426 ND nd_reduce2(int mod,ND d,ND g,NDV p,NM mul,NDC dn,Ob

}

*divp = (Obj)credp;

} else {

igcd_cofactor(HCQ(g),HCQ(p),&gcd,&cg,&cred);

igcd_cofactor(HCZ(g),HCZ(p),&gcd,&cg,&cred);

chsgnz(cg,&CQ(mul));

chsgnz(cg,&CZ(mul));

nd_mul_c_q(d,(P)cred); nd_mul_c_q(g,(P)cred);

if ( dn ) {

mulz(dn->z,cred,&tq); dn->z = tq;

Line 1424 ND nd_reduce2(int mod,ND d,ND g,NDV p,NM mul,NDC dn,Ob

Line 1438 ND nd_reduce2(int mod,ND d,ND g,NDV p,NM mul,NDC dn,Ob

}

/* ret=1 : success, ret=0 : overflow */

int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND *rp)

int nd_nf(int mod,ND d,ND g,NDV *ps,int full,ND *rp)

{

NM m,mrd,tail;

NM mul;

Line 1443 int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND

Line 1457 int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND

union oNDC hg;

P cont;

if ( dn ) {

if ( mod )

dn->m = 1;

else if ( nd_vc )

dn->r = (R)ONE;

else

dn->z = ONE;

}

if ( !g ) {

*rp = d;

return 1;

Line 1473 int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND

Line 1479 int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND

}

p = nd_demand ? ndv_load(index) : ps[index];

/* d+g -> div*(d+g)+mul*p */

g = nd_reduce2(mod,d,g,p,mul,dn,&div);

g = nd_reduce2(mod,d,g,p,mul,0,&div);

if ( nd_gentrace ) {

/* Trace=[div,index,mul,ONE] */

STOQ(index,iq);

STOZ(index,iq);

nmtodp(mod,mul,&dmul);

node = mknode(4,div,iq,dmul,ONE);

}

Line 1484 int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND

Line 1490 int nd_nf(int mod,ND d,ND g,NDV *ps,int full,NDC dn,ND

if ( !mod && g && !nd_vc && ((double)(p_mag(HCP(g))) > hmag) ) {

hg = HCU(g);

nd_removecont2(d,g);

if ( dn || nd_gentrace ) {

if ( nd_gentrace ) {

/* overwrite cont : Trace=[div,index,mul,cont] */

/* exact division */

cont = ndc_div(mod,hg,HCU(g));

if ( dn ) {

if ( nd_vc ) {

divr(nd_vc,(Obj)dn->r,(Obj)cont,&tr);

reductr(nd_vc,(Obj)tr,&tr1); dn->r = (R)tr1;

} else divz(dn->z,(Z)cont,&dn->z);

}

if ( nd_gentrace && !UNIQ(cont) ) ARG3(node) = (pointer)cont;

}

hmag = ((double)p_mag(HCP(g)))*nd_scale;

Line 1543 int nd_nf_pbucket(int mod,ND g,NDV *ps,int full,ND *rp

Line 1544 int nd_nf_pbucket(int mod,ND g,NDV *ps,int full,ND *rp

}

sugar = SG(g);

n = NV(g);

if ( !mod ) hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;

if ( !mod ) hmag = ((double)p_mag((P)HCZ(g)))*nd_scale;

bucket = create_pbucket();

add_pbucket(mod,bucket,g);

d = 0;

Line 1586 int nd_nf_pbucket(int mod,ND g,NDV *ps,int full,ND *rp

Line 1587 int nd_nf_pbucket(int mod,ND g,NDV *ps,int full,ND *rp

c1 = invm(HCM(p),mod); c2 = mod-HCM(g);

DMAR(c1,c2,0,mod,c); CM(mul) = c;

} else {

igcd_cofactor(HCQ(g),HCQ(p),&gcd,&cg,&cred);

igcd_cofactor(HCZ(g),HCZ(p),&gcd,&cg,&cred);

chsgnz(cg,&CQ(mul));

chsgnz(cg,&CZ(mul));

nd_mul_c_q(d,(P)cred);

mulq_pbucket(bucket,cred);

g = bucket->body[hindex];

gmag = (double)p_mag((P)HCQ(g));

gmag = (double)p_mag((P)HCZ(g));

}

red = ndv_mul_nm(mod,mul,p);

bucket->body[hindex] = nd_remove_head(g);

Line 1607 int nd_nf_pbucket(int mod,ND g,NDV *ps,int full,ND *rp

Line 1608 int nd_nf_pbucket(int mod,ND g,NDV *ps,int full,ND *rp

return 1;

}

nd_removecont2(d,g);

hmag = ((double)p_mag((P)HCQ(g)))*nd_scale;

hmag = ((double)p_mag((P)HCZ(g)))*nd_scale;

add_pbucket(mod,bucket,g);

}

} else if ( !full ) {

Line 1662 again:

Line 1663 again:

if ( m == -2 ) ndv_mod(m,r);

#endif

d = ndvtond(m,r);

stat = nd_nf(m,0,d,nd_ps,0,0,&nf);

stat = nd_nf(m,0,d,nd_ps,0,&nf);

if ( !stat ) {

nd_reconstruct(0,0);

goto again;

Line 1670 again:

Line 1671 again:

if ( nd_gentrace ) {

nd_tracelist = reverse_node(nd_tracelist);

MKLIST(list,nd_tracelist);

STOQ(i,q); s = mknode(2,q,list); MKLIST(list,s);

STOZ(i,q); s = mknode(2,q,list); MKLIST(list,s);

MKNODE(s,list,nd_alltracelist);

nd_alltracelist = s; nd_tracelist = 0;

}

Line 1868 int head_pbucket_q(PGeoBucket g)

Line 1869 int head_pbucket_q(PGeoBucket g)

if ( j < 0 ) {

j = i;

gj = g->body[j];

sum = HCQ(gj);

sum = HCZ(gj);

} else {

nv = NV(gi);

c = DL_COMPARE(HDL(gi),HDL(gj));

if ( c > 0 ) {

if ( sum ) HCQ(gj) = sum;

if ( sum ) HCZ(gj) = sum;

else g->body[j] = nd_remove_head(gj);

j = i;

gj = g->body[j];

sum = HCQ(gj);

sum = HCZ(gj);

} else if ( c == 0 ) {

addz(sum,HCQ(gi),&t);

addz(sum,HCZ(gi),&t);

sum = t;

g->body[i] = nd_remove_head(gi);

}

Line 1887 int head_pbucket_q(PGeoBucket g)

Line 1888 int head_pbucket_q(PGeoBucket g)

}

if ( j < 0 ) return -1;

else if ( sum ) {

HCQ(gj) = sum;

HCZ(gj) = sum;

return j;

} else

g->body[j] = nd_remove_head(gj);

Line 2012 void register_hcf(NDV p)

Line 2013 void register_hcf(NDV p)

int do_diagonalize(int sugar,int m)

{

int i,nh,stat;

NODE r,g,t;

ND h,nf,s,head;

NDV nfv;

Q q;

P nm,nmp,dn,mnp,dnp,cont,cont1;

union oNDC hc;

NODE node;

LIST l;

Z iq;

for ( i = nd_psn-1; i >= 0 && SG(nd_psh[i]) == sugar; i-- ) {

if ( nd_gentrace ) {

/* Trace = [1,index,1,1] */

STOQ(i,iq); node = mknode(4,ONE,iq,ONE,ONE);

STOZ(i,iq); node = mknode(4,ONE,iq,ONE,ONE);

MKLIST(l,node); MKNODE(nd_tracelist,l,0);

}

if ( nd_demand )

nfv = ndv_load(i);

else

nfv = nd_ps[i];

s = ndvtond(m,nfv);

s = nd_separate_head(s,&head);

stat = nd_nf(m,head,s,nd_ps,1,0,&nf);

if ( !stat ) return 0;

ndv_free(nfv);

hc = HCU(nf); nd_removecont(m,nf);

cont = ndc_div(m,hc,HCU(nf));

if ( nd_gentrace ) finalize_tracelist(i,cont);

nfv = ndtondv(m,nf);

nd_free(nf);

nd_bound[i] = ndv_compute_bound(nfv);

if ( !m ) register_hcf(nfv);

if ( nd_demand ) {

ndv_save(nfv,i);

ndv_free(nfv);

} else

nd_ps[i] = nfv;

}

return 1;

if ( nd_demand )

nfv = ndv_load(i);

else

nfv = nd_ps[i];

s = ndvtond(m,nfv);

s = nd_separate_head(s,&head);

stat = nd_nf(m,head,s,nd_ps,1,&nf);

if ( !stat ) return 0;

ndv_free(nfv);

hc = HCU(nf); nd_removecont(m,nf);

/* exact division */

cont = ndc_div(m,hc,HCU(nf));

if ( nd_gentrace ) finalize_tracelist(i,cont);

nfv = ndtondv(m,nf);

nd_free(nf);

nd_bound[i] = ndv_compute_bound(nfv);

if ( !m ) register_hcf(nfv);

if ( nd_demand ) {

ndv_save(nfv,i);

ndv_free(nfv);

} else

nd_ps[i] = nfv;

}

return 1;

}

LIST compute_splist()

Line 2069 LIST compute_splist()

Line 2071 LIST compute_splist()

}

for ( t = d, tn0 = 0; t; t = NEXT(t) ) {

NEXTNODE(tn0,tn);

STOQ(t->i1,i1); STOQ(t->i2,i2);

STOZ(t->i1,i1); STOZ(t->i2,i2);

node = mknode(2,i1,i2); MKLIST(l0,node);

BDY(tn) = l0;

}

Line 2081 LIST compute_splist()

Line 2083 LIST compute_splist()

NODE nd_gb(int m,int ishomo,int checkonly,int gensyz,int **indp)

{

int i,nh,sugar,stat;

NODE r,g,t;

ND_pairs d;

ND_pairs l;

ND h,nf,s,head,nf1;

NDV nfv;

Z q;

union oNDC dn,hc;

int diag_count = 0;

P cont;

LIST list;

g = 0; d = 0;

Nnd_add = 0;

for ( i = 0; i < nd_psn; i++ ) {

g = 0; d = 0;

d = update_pairs(d,g,i,gensyz);

for ( i = 0; i < nd_psn; i++ ) {

g = update_base(g,i);

d = update_pairs(d,g,i,gensyz);

}

g = update_base(g,i);

sugar = 0;

}

while ( d ) {

sugar = 0;

while ( d ) {

again:

l = nd_minp(d,&d);

if ( MaxDeg > 0 && SG(l) > MaxDeg ) break;

if ( SG(l) != sugar ) {

if ( ishomo ) {

diag_count = 0;

stat = do_diagonalize(sugar,m);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

}

sugar = SG(l);

if ( DP_Print ) fprintf(asir_out,"%d",sugar);

}

stat = nd_sp(m,0,l,&h);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

}

sugar = SG(l);

if ( DP_Print ) fprintf(asir_out,"%d",sugar);

}

stat = nd_sp(m,0,l,&h);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

}

#if USE_GEOBUCKET

stat = (m&&!nd_gentrace)?nd_nf_pbucket(m,h,nd_ps,!Top,&nf)

:nd_nf(m,0,h,nd_ps,!Top,0,&nf);

:nd_nf(m,0,h,nd_ps,!Top,&nf);

#else

stat = nd_nf(m,0,h,nd_ps,!Top,0,&nf);

stat = nd_nf(m,0,h,nd_ps,!Top,&nf);

#endif

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

} else if ( nf ) {

if ( checkonly || gensyz ) return 0;

if ( nd_newelim ) {

if ( nd_module ) {

if ( MPOS(HDL(nf)) > 1 ) return 0;

} else if ( !(HDL(nf)[nd_exporigin] & nd_mask[0]) ) return 0;

}

if ( DP_Print ) { printf("+"); fflush(stdout); }

hc = HCU(nf);

nd_removecont(m,nf);

if ( !m && nd_nalg ) {

nd_monic(0,&nf);

nd_removecont(m,nf);

}

if ( nd_gentrace ) {

/* exact division */

cont = ndc_div(m,hc,HCU(nf));

if ( m || !UNIQ(cont) ) {

t = mknode(4,NULLP,NULLP,NULLP,cont);

MKLIST(list,t); MKNODE(t,list,nd_tracelist);

nd_tracelist = t;

}

nfv = ndtondv(m,nf); nd_free(nf);

nh = ndv_newps(m,nfv,0,0);

if ( !m && (ishomo && ++diag_count == diag_period) ) {

diag_count = 0;

stat = do_diagonalize(sugar,m);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

d = update_pairs(d,g,nh,0);

g = update_base(g,nh);

FREENDP(l);

} else {

if ( nd_gentrace && gensyz ) {

nd_tracelist = reverse_node(nd_tracelist);

MKLIST(list,nd_tracelist);

STOQ(-1,q); t = mknode(2,q,list); MKLIST(list,t);

MKNODE(t,list,nd_alltracelist);

nd_alltracelist = t; nd_tracelist = 0;

}

if ( DP_Print ) { printf("."); fflush(stdout); }

nfv = ndtondv(m,nf); nd_free(nf);

FREENDP(l);

nh = ndv_newps(m,nfv,0,0);

if ( !m && (ishomo && ++diag_count == diag_period) ) {

diag_count = 0;

stat = do_diagonalize(sugar,m);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

conv_ilist(nd_demand,0,g,indp);

d = update_pairs(d,g,nh,0);

if ( !checkonly && DP_Print ) { printf("nd_gb done.\n"); fflush(stdout); }

g = update_base(g,nh);

FREENDP(l);

} else {

if ( nd_gentrace && gensyz ) {

nd_tracelist = reverse_node(nd_tracelist);

MKLIST(list,nd_tracelist);

STOZ(-1,q); t = mknode(2,q,list); MKLIST(list,t);

MKNODE(t,list,nd_alltracelist);

nd_alltracelist = t; nd_tracelist = 0;

}

if ( DP_Print ) { printf("."); fflush(stdout); }

FREENDP(l);

}

conv_ilist(nd_demand,0,g,indp);

if ( !checkonly && DP_Print ) { printf("nd_gb done. Number of nd_add=%d\n",Nnd_add); fflush(stdout); }

return g;

}

Line 2196 int check_splist(int m,NODE splist)

Line 2200 int check_splist(int m,NODE splist)

for ( d = 0, t = splist; t; t = NEXT(t) ) {

p = BDY((LIST)BDY(t));

NEXTND_pairs(d,r);

r->i1 = QTOS((Q)ARG0(p)); r->i2 = QTOS((Q)ARG1(p));

r->i1 = ZTOS((Q)ARG0(p)); r->i2 = ZTOS((Q)ARG1(p));

ndl_lcm(DL(nd_psh[r->i1]),DL(nd_psh[r->i2]),r->lcm);

SG(r) = TD(LCM(r)); /* XXX */

}

if ( d ) NEXT(r) = 0;

while ( d ) {

again:

l = nd_minp(d,&d);

stat = nd_sp(m,0,l,&h);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

}

stat = nd_nf(m,0,h,nd_ps,!Top,0,&nf);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

} else if ( nf ) return 0;

if ( DP_Print) { printf("."); fflush(stdout); }

}

stat = nd_nf(m,0,h,nd_ps,!Top,&nf);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(0,d);

goto again;

} else if ( nf ) return 0;

if ( DP_Print) { printf("."); fflush(stdout); }

}

if ( DP_Print) { printf("done.\n"); fflush(stdout); }

return 1;

}

Line 2227 again:

Line 2231 again:

int check_splist_f4(int m,NODE splist)

{

UINT *s0vect;

PGeoBucket bucket;

NODE p,rp0,t;

ND_pairs d,r,l,ll;

int col,stat;

for ( d = 0, t = splist; t; t = NEXT(t) ) {

p = BDY((LIST)BDY(t));

NEXTND_pairs(d,r);

r->i1 = QTOS((Q)ARG0(p)); r->i2 = QTOS((Q)ARG1(p));

r->i1 = ZTOS((Q)ARG0(p)); r->i2 = ZTOS((Q)ARG1(p));

ndl_lcm(DL(nd_psh[r->i1]),DL(nd_psh[r->i2]),r->lcm);

SG(r) = TD(LCM(r)); /* XXX */

}

if ( d ) NEXT(r) = 0;

while ( d ) {

l = nd_minsugarp(d,&d);

bucket = create_pbucket();

stat = nd_sp_f4(m,0,l,bucket);

if ( !stat ) {

for ( ll = l; NEXT(ll); ll = NEXT(ll) );

NEXT(ll) = d; d = l;

d = nd_reconstruct(0,d);

continue;

}

if ( bucket->m < 0 ) continue;

col = nd_symbolic_preproc(bucket,0,&s0vect,&rp0);

if ( !col ) {

for ( ll = l; NEXT(ll); ll = NEXT(ll) );

NEXT(ll) = d; d = l;

d = nd_reconstruct(0,d);

continue;

}

if ( nd_f4_red(m,l,0,s0vect,col,rp0,0) ) return 0;

}

return 1;

if ( bucket->m < 0 ) continue;

col = nd_symbolic_preproc(bucket,0,&s0vect,&rp0);

if ( !col ) {

for ( ll = l; NEXT(ll); ll = NEXT(ll) );

NEXT(ll) = d; d = l;

d = nd_reconstruct(0,d);

continue;

}

if ( nd_f4_red(m,l,0,s0vect,col,rp0,0) ) return 0;

}

return 1;

}

int do_diagonalize_trace(int sugar,int m)

{

int i,nh,stat;

NODE r,g,t;

ND h,nf,nfq,s,head;

NDV nfv,nfqv;

Q q,den,num;

union oNDC hc;

NODE node;

LIST l;

Z iq;

P cont,cont1;

for ( i = nd_psn-1; i >= 0 && SG(nd_psh[i]) == sugar; i-- ) {

if ( nd_gentrace ) {

/* Trace = [1,index,1,1] */

STOQ(i,iq); node = mknode(4,ONE,iq,ONE,ONE);

STOZ(i,iq); node = mknode(4,ONE,iq,ONE,ONE);

MKLIST(l,node); MKNODE(nd_tracelist,l,0);

}

/* for nd_ps */

s = ndvtond(m,nd_ps[i]);

s = nd_separate_head(s,&head);

stat = nd_nf_pbucket(m,s,nd_ps,1,&nf);

if ( !stat ) return 0;

nf = nd_add(m,head,nf);

ndv_free(nd_ps[i]);

nd_ps[i] = ndtondv(m,nf);

nd_free(nf);

/* for nd_ps_trace */

if ( nd_demand )

nfv = ndv_load(i);

else

nfv = nd_ps_trace[i];

s = ndvtond(0,nfv);

s = nd_separate_head(s,&head);

stat = nd_nf(0,head,s,nd_ps_trace,1,0,&nf);

stat = nd_nf(0,head,s,nd_ps_trace,1,&nf);

if ( !stat ) return 0;

ndv_free(nfv);

hc = HCU(nf); nd_removecont(0,nf);

/* exact division */

cont = ndc_div(0,hc,HCU(nf));

if ( nd_gentrace ) finalize_tracelist(i,cont);

nfv = ndtondv(0,nf);

nd_free(nf);

nd_bound[i] = ndv_compute_bound(nfv);

register_hcf(nfv);

if ( nd_demand ) {

ndv_save(nfv,i);

ndv_free(nfv);

} else

nd_ps_trace[i] = nfv;

}

return 1;

}

static struct oEGT eg_invdalg;

Line 2335 void nd_subst_vector(VL vl,P p,NODE subst,P *r)

Line 2340 void nd_subst_vector(VL vl,P p,NODE subst,P *r)

NODE nd_gb_trace(int m,int ishomo,int **indp)

{

int i,nh,sugar,stat;

NODE r,g,t;

ND_pairs d;

ND_pairs l;

ND h,nf,nfq,s,head;

NDV nfv,nfqv;

Z q,den,num;

P hc;

union oNDC dn,hnfq;

struct oEGT eg_monic,egm0,egm1;

int diag_count = 0;

P cont;

LIST list;

init_eg(&eg_monic);

init_eg(&eg_invdalg);

init_eg(&eg_le);

g = 0; d = 0;

for ( i = 0; i < nd_psn; i++ ) {

d = update_pairs(d,g,i,0);

g = update_base(g,i);

}

sugar = 0;

while ( d ) {

again:

l = nd_minp(d,&d);

if ( MaxDeg > 0 && SG(l) > MaxDeg ) break;

if ( SG(l) != sugar ) {

#if 1

if ( ishomo ) {

if ( DP_Print > 2 ) fprintf(asir_out,"|");

stat = do_diagonalize_trace(sugar,m);

if ( DP_Print > 2 ) fprintf(asir_out,"|");

diag_count = 0;

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

#endif

sugar = SG(l);

if ( DP_Print ) fprintf(asir_out,"%d",sugar);

}

stat = nd_sp(m,0,l,&h);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

#endif

sugar = SG(l);

if ( DP_Print ) fprintf(asir_out,"%d",sugar);

}

stat = nd_sp(m,0,l,&h);

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

#if USE_GEOBUCKET

stat = nd_nf_pbucket(m,h,nd_ps,!Top,&nf);

#else

stat = nd_nf(m,0,h,nd_ps,!Top,0,&nf);

stat = nd_nf(m,0,h,nd_ps,!Top,&nf);

#endif

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

} else if ( nf ) {

if ( nd_demand ) {

nfqv = ndv_load(nd_psn);

nfq = ndvtond(0,nfqv);

} else

nfq = 0;

if ( !nfq ) {

if ( !nd_sp(0,1,l,&h) || !nd_nf(0,0,h,nd_ps_trace,!Top,&nfq) ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

if ( nfq ) {

/* m|HC(nfq) => failure */

if ( nd_vc ) {

nd_subst_vector(nd_vc,HCP(nfq),nd_subst,&hc); q = (Z)hc;

} else

q = HCZ(nfq);

if ( !remqi((Q)q,m) ) return 0;

if ( DP_Print ) { printf("+"); fflush(stdout); }

hnfq = HCU(nfq);

if ( nd_nalg ) {

/* m|DN(HC(nf)^(-1)) => failure */

get_eg(&egm0);

if ( !nd_monic(m,&nfq) ) return 0;

get_eg(&egm1); add_eg(&eg_monic,&egm0,&egm1);

nd_removecont(0,nfq); nfqv = ndtondv(0,nfq); nd_free(nfq);

nfv = ndv_dup(0,nfqv); ndv_mod(m,nfv); nd_free(nf);

} else {

nd_removecont(0,nfq); nfqv = ndtondv(0,nfq); nd_free(nfq);

nd_removecont(m,nf); nfv = ndtondv(m,nf); nd_free(nf);

}

if ( nd_gentrace ) {

/* exact division */

cont = ndc_div(0,hnfq,HCU(nfqv));

if ( !UNIQ(cont) ) {

t = mknode(4,NULLP,NULLP,NULLP,cont);

MKLIST(list,t); MKNODE(t,list,nd_tracelist);

nd_tracelist = t;

}

nh = ndv_newps(0,nfv,nfqv,0);

if ( ishomo && ++diag_count == diag_period ) {

diag_count = 0;

if ( DP_Print > 2 ) fprintf(asir_out,"|");

stat = do_diagonalize_trace(sugar,m);

if ( DP_Print > 2 ) fprintf(asir_out,"|");

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

} else if ( nf ) {

}

if ( nd_demand ) {

nfqv = ndv_load(nd_psn);

nfq = ndvtond(0,nfqv);

} else

nfq = 0;

if ( !nfq ) {

if ( !nd_sp(0,1,l,&h) || !nd_nf(0,0,h,nd_ps_trace,!Top,0,&nfq) ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

if ( nfq ) {

/* m|HC(nfq) => failure */

if ( nd_vc ) {

nd_subst_vector(nd_vc,HCP(nfq),nd_subst,&hc); q = (Z)hc;

} else

q = HCQ(nfq);

if ( !remqi((Q)q,m) ) return 0;

if ( DP_Print ) { printf("+"); fflush(stdout); }

hnfq = HCU(nfq);

if ( nd_nalg ) {

/* m|DN(HC(nf)^(-1)) => failure */

get_eg(&egm0);

if ( !nd_monic(m,&nfq) ) return 0;

get_eg(&egm1); add_eg(&eg_monic,&egm0,&egm1);

nd_removecont(0,nfq); nfqv = ndtondv(0,nfq); nd_free(nfq);

nfv = ndv_dup(0,nfqv); ndv_mod(m,nfv); nd_free(nf);

} else {

nd_removecont(0,nfq); nfqv = ndtondv(0,nfq); nd_free(nfq);

nd_removecont(m,nf); nfv = ndtondv(m,nf); nd_free(nf);

}

if ( nd_gentrace ) {

cont = ndc_div(0,hnfq,HCU(nfqv));

if ( !UNIQ(cont) ) {

t = mknode(4,NULLP,NULLP,NULLP,cont);

MKLIST(list,t); MKNODE(t,list,nd_tracelist);

nd_tracelist = t;

}

nh = ndv_newps(0,nfv,nfqv,0);

if ( ishomo && ++diag_count == diag_period ) {

diag_count = 0;

if ( DP_Print > 2 ) fprintf(asir_out,"|");

stat = do_diagonalize_trace(sugar,m);

if ( DP_Print > 2 ) fprintf(asir_out,"|");

if ( !stat ) {

NEXT(l) = d; d = l;

d = nd_reconstruct(1,d);

goto again;

}

d = update_pairs(d,g,nh,0);

g = update_base(g,nh);

} else {

if ( DP_Print ) { printf("*"); fflush(stdout); }

}

} else {

if ( DP_Print ) { printf("."); fflush(stdout); }

}

FREENDP(l);

d = update_pairs(d,g,nh,0);

g = update_base(g,nh);

} else {

if ( DP_Print ) { printf("*"); fflush(stdout); }

}

} else {

if ( DP_Print ) { printf("."); fflush(stdout); }

}

if ( nd_nalg ) {

FREENDP(l);

if ( DP_Print ) {

}

print_eg("monic",&eg_monic);

if ( nd_nalg ) {

print_eg("invdalg",&eg_invdalg);

if ( DP_Print ) {

print_eg("le",&eg_le);

print_eg("monic",&eg_monic);

}

print_eg("invdalg",&eg_invdalg);

print_eg("le",&eg_le);

}

conv_ilist(nd_demand,1,g,indp);

if ( DP_Print ) { printf("nd_gb_trace done.\n"); fflush(stdout); }

return g;

}

int ndv_compare(NDV *p1,NDV *p2)

Line 2510 NODE ndv_reduceall(int m,NODE f)

Line 2516 NODE ndv_reduceall(int m,NODE f)

perm = (int *)MALLOC(n*sizeof(int));

if ( nd_gentrace ) {

for ( t = nd_tracelist, i = 0; i < n; i++, t = NEXT(t) )

perm[i] = QTOS((Q)ARG1(BDY((LIST)BDY(t))));

perm[i] = ZTOS((Q)ARG1(BDY((LIST)BDY(t))));

}

for ( i = 0; i < n; ) {

if ( nd_gentrace ) {

/* Trace = [1,index,1,1] */

STOQ(i,iq); node = mknode(4,ONE,iq,ONE,ONE);

STOZ(i,iq); node = mknode(4,ONE,iq,ONE,ONE);

MKLIST(l,node); MKNODE(nd_tracelist,l,0);

}

g = ndvtond(m,nd_ps[i]);

g = nd_separate_head(g,&head);

stat = nd_nf(m,head,g,nd_ps,1,0,&nf);

stat = nd_nf(m,head,g,nd_ps,1,&nf);

if ( !stat )

nd_reconstruct(0,0);

else {

Line 2529 NODE ndv_reduceall(int m,NODE f)

Line 2535 NODE ndv_reduceall(int m,NODE f)

hc = HCU(nf); nd_removecont(m,nf);

if ( nd_gentrace ) {

for ( t = nd_tracelist; t; t = NEXT(t) ) {

jq = ARG1(BDY((LIST)BDY(t))); j = QTOS(jq);

jq = ARG1(BDY((LIST)BDY(t))); j = ZTOS(jq);

STOQ(perm[j],jq); ARG1(BDY((LIST)BDY(t))) = jq;

STOZ(perm[j],jq); ARG1(BDY((LIST)BDY(t))) = jq;

}

/* exact division */

cont = ndc_div(m,hc,HCU(nf));

finalize_tracelist(perm[i],cont);

}

Line 2609 ND_pairs nd_newpairs( NODE g, int t )

Line 2616 ND_pairs nd_newpairs( NODE g, int t )

dl = DL(nd_psh[t]);

ts = SG(nd_psh[t]) - TD(dl);

if ( nd_module && nd_intersect && (MPOS(dl) > 1) ) return 0;

if ( nd_module && nd_intersect && (MPOS(dl) > nd_intersect) ) return 0;

for ( r0 = 0, h = g; h; h = NEXT(h) ) {

if ( nd_module && (MPOS(DL(nd_psh[(long)BDY(h)])) != MPOS(dl)) )

continue;

Line 2647 ND_pairs nd_ipairtospair(NODE ipair)

Line 2654 ND_pairs nd_ipairtospair(NODE ipair)

for ( r0 = 0, t = ipair; t; t = NEXT(t) ) {

NEXTND_pairs(r0,r);

tn = BDY((LIST)BDY(t));

r->i1 = QTOS((Q)ARG0(tn));

r->i1 = ZTOS((Q)ARG0(tn));

r->i2 = QTOS((Q)ARG1(tn));

r->i2 = ZTOS((Q)ARG1(tn));

ndl_lcm(DL(nd_psh[r->i1]),DL(nd_psh[r->i2]),r->lcm);

s1 = SG(nd_psh[r->i1])-TD(DL(nd_psh[r->i1]));

s2 = SG(nd_psh[r->i2])-TD(DL(nd_psh[r->i2]));

Line 2987 int ndv_newps(int m,NDV a,NDV aq,int f4)

Line 2994 int ndv_newps(int m,NDV a,NDV aq,int f4)

if ( nd_gentrace ) {

/* reverse the tracelist and append it to alltracelist */

nd_tracelist = reverse_node(nd_tracelist); MKLIST(l,nd_tracelist);

STOQ(nd_psn,iq); tn = mknode(2,iq,l); MKLIST(l,tn);

STOZ(nd_psn,iq); tn = mknode(2,iq,l); MKLIST(l,tn);

MKNODE(tn,l,nd_alltracelist); nd_alltracelist = tn; nd_tracelist = 0;

}

return nd_psn++;

Line 2998 int ndv_newps(int m,NDV a,NDV aq,int f4)

Line 3005 int ndv_newps(int m,NDV a,NDV aq,int f4)

int ndv_setup(int mod,int trace,NODE f,int dont_sort,int dont_removecont)

{

int i,j,td,len,max;

NODE s,s0,f0,tn;

UINT *d;

RHist r;

NDVI w;

NDV a,am;

union oNDC hc;

NODE node;

P hcp;

Z iq,jq;

LIST l;

nd_found = 0; nd_notfirst = 0; nd_create = 0;

/* initialize the tracelist */

nd_tracelist = 0;

for ( nd_psn = 0, s = f; s; s = NEXT(s) ) if ( BDY(s) ) nd_psn++;

w = (NDVI)MALLOC(nd_psn*sizeof(struct oNDVI));

for ( i = j = 0, s = f; s; s = NEXT(s), j++ )

if ( BDY(s) ) { w[i].p = BDY(s); w[i].i = j; i++; }

if ( !dont_sort ) {

/* XXX heuristic */

if ( !nd_ord->id && (nd_ord->ord.simple<2) )

qsort(w,nd_psn,sizeof(struct oNDVI),

(int (*)(const void *,const void *))ndvi_compare_rev);

else

qsort(w,nd_psn,sizeof(struct oNDVI),

(int (*)(const void *,const void *))ndvi_compare);

}

nd_pslen = 2*nd_psn;

nd_ps = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_ps_trace = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_ps_sym = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_ps_trace_sym = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_psh = (RHist *)MALLOC(nd_pslen*sizeof(RHist));

nd_bound = (UINT **)MALLOC(nd_pslen*sizeof(UINT *));

nd_hcf = 0;

if ( trace && nd_vc )

makesubst(nd_vc,&nd_subst);

else

nd_subst = 0;

qsort(w,nd_psn,sizeof(struct oNDVI),

(int (*)(const void *,const void *))ndvi_compare);

}

nd_pslen = 2*nd_psn;

nd_ps = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_ps_trace = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_ps_sym = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_ps_trace_sym = (NDV *)MALLOC(nd_pslen*sizeof(NDV));

nd_psh = (RHist *)MALLOC(nd_pslen*sizeof(RHist));

nd_bound = (UINT **)MALLOC(nd_pslen*sizeof(UINT *));

nd_hcf = 0;

if ( !nd_red )

if ( trace && nd_vc )

nd_red = (RHist *)MALLOC(REDTAB_LEN*sizeof(RHist));

makesubst(nd_vc,&nd_subst);

for ( i = 0; i < REDTAB_LEN; i++ ) nd_red[i] = 0;

else

for ( i = 0; i < nd_psn; i++ ) {

nd_subst = 0;

hc = HCU(w[i].p);

if ( trace ) {

if ( !nd_red )

if ( mod == -2 ) {

nd_red = (RHist *)MALLOC(REDTAB_LEN*sizeof(RHist));

/* over a large finite field */

for ( i = 0; i < REDTAB_LEN; i++ ) nd_red[i] = 0;

/* trace = small modulus */

for ( i = 0; i < nd_psn; i++ ) {

a = nd_ps_trace[i] = ndv_dup(-2,w[i].p);

hc = HCU(w[i].p);

ndv_mod(-2,a);

if ( trace ) {

if ( !dont_removecont) ndv_removecont(-2,a);

if ( mod == -2 ) {

am = nd_ps[i] = ndv_dup(trace,w[i].p);

/* over a large finite field */

ndv_mod(trace,am);

/* trace = small modulus */

if ( DL_COMPARE(HDL(am),HDL(a)) )

a = nd_ps_trace[i] = ndv_dup(-2,w[i].p);

return 0;

ndv_mod(-2,a);

ndv_removecont(trace,am);

if ( !dont_removecont) ndv_removecont(-2,a);

} else {

am = nd_ps[i] = ndv_dup(trace,w[i].p);

a = nd_ps_trace[i] = ndv_dup(0,w[i].p);

ndv_mod(trace,am);

if ( !dont_removecont) ndv_removecont(0,a);

if ( DL_COMPARE(HDL(am),HDL(a)) )

register_hcf(a);

return 0;

am = nd_ps[i] = ndv_dup(mod,a);

ndv_removecont(trace,am);

ndv_mod(mod,am);

} else {

if ( DL_COMPARE(HDL(am),HDL(a)) )

a = nd_ps_trace[i] = ndv_dup(0,w[i].p);

return 0;

if ( !dont_removecont) ndv_removecont(0,a);

ndv_removecont(mod,am);

register_hcf(a);

}

am = nd_ps[i] = ndv_dup(mod,a);

} else {

ndv_mod(mod,am);

a = nd_ps[i] = ndv_dup(mod,w[i].p);

if ( DL_COMPARE(HDL(am),HDL(a)) )

if ( mod || !dont_removecont ) ndv_removecont(mod,a);

return 0;

if ( !mod ) register_hcf(a);

ndv_removecont(mod,am);

}

if ( nd_gentrace ) {

} else {

STOQ(i,iq); STOQ(w[i].i,jq); node = mknode(3,iq,jq,ONE);

a = nd_ps[i] = ndv_dup(mod,w[i].p);

if ( mod || !dont_removecont ) ndv_removecont(mod,a);

if ( !mod ) register_hcf(a);

}

if ( nd_gentrace ) {

STOZ(i,iq); STOZ(w[i].i,jq); node = mknode(3,iq,jq,ONE);

/* exact division */

if ( !dont_removecont )

ARG2(node) = (pointer)ndc_div(trace?0:mod,hc,HCU(a));

MKLIST(l,node); NEXTNODE(nd_tracelist,tn); BDY(tn) = l;

}

NEWRHist(r); SG(r) = HTD(a); ndl_copy(HDL(a),DL(r));

nd_bound[i] = ndv_compute_bound(a);

nd_psh[i] = r;

if ( nd_demand ) {

if ( trace ) {

ndv_save(nd_ps_trace[i],i);

nd_ps_sym[i] = ndv_symbolic(mod,nd_ps_trace[i]);

nd_ps_trace_sym[i] = ndv_symbolic(mod,nd_ps_trace[i]);

nd_ps_trace[i] = 0;

} else {

ndv_save(nd_ps[i],i);

nd_ps_sym[i] = ndv_symbolic(mod,nd_ps[i]);

nd_ps[i] = 0;

}

if ( nd_gentrace && nd_tracelist ) NEXT(tn) = 0;

NEWRHist(r); SG(r) = HTD(a); ndl_copy(HDL(a),DL(r));

return 1;

nd_bound[i] = ndv_compute_bound(a);

nd_psh[i] = r;

if ( nd_demand ) {

if ( trace ) {

ndv_save(nd_ps_trace[i],i);

nd_ps_sym[i] = ndv_symbolic(mod,nd_ps_trace[i]);

nd_ps_trace_sym[i] = ndv_symbolic(mod,nd_ps_trace[i]);

nd_ps_trace[i] = 0;

} else {

ndv_save(nd_ps[i],i);

nd_ps_sym[i] = ndv_symbolic(mod,nd_ps[i]);

nd_ps[i] = 0;

}

if ( nd_gentrace && nd_tracelist ) NEXT(tn) = 0;

return 1;

}

struct order_spec *append_block(struct order_spec *spec,

Line 3227 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

Line 3235 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

int *perm;

EPOS oepos;

int obpe,oadv,ompos,cbpe;

VECT hvect;

nd_module = 0;

if ( !m && Demand ) nd_demand = 1;

Line 3270 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

Line 3279 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

for ( t = BDY(f), max = 1; t; t = NEXT(t) )

for ( tv = vv; tv; tv = NEXT(tv) ) {

if ( nd_module ) {

s = BDY((LIST)BDY(t));

if ( OID(BDY(t)) == O_DPM ) {

trank = length(s);

e = dpm_getdeg((DPM)BDY(t),&trank);

mrank = MAX(mrank,trank);

max = MAX(e,max);

for ( ; s; s = NEXT(s) ) {

mrank = MAX(mrank,trank);

e = getdeg(tv->v,(P)BDY(s));

} else {

max = MAX(e,max);

s = BDY((LIST)BDY(t));

trank = length(s);

mrank = MAX(mrank,trank);

for ( ; s; s = NEXT(s) ) {

e = getdeg(tv->v,(P)BDY(s));

max = MAX(e,max);

}

} else {

e = getdeg(tv->v,(P)BDY(t));

Line 3287 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

Line 3302 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

ishomo = 1;

for ( fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {

if ( nd_module ) {

if ( !m && !nd_gentrace ) pltozpl((LIST)BDY(t),&dmy,&zpl);

if ( OID(BDY(t)) == O_DPM ) {

else zpl = (LIST)BDY(t);

Z cont;

DPM zdpm;

if ( !m && !nd_gentrace ) dpm_ptozp((DPM)BDY(t),&cont,&zdpm);

else zdpm = (DPM)BDY(t);

b = (pointer)dpmtondv(m,zdpm);

} else {

if ( !m && !nd_gentrace ) pltozpl((LIST)BDY(t),&dmy,&zpl);

else zpl = (LIST)BDY(t);

b = (pointer)pltondv(CO,vv,zpl);

}

} else {

if ( !m && !nd_gentrace ) ptozp((P)BDY(t),1,&dmy,&zp);

else zp = (P)BDY(t);

Line 3337 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

Line 3361 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

if ( !x ) {

*rp = 0; return;

}

if ( nd_gentrace ) {

MKVECT(hvect,nd_psn);

for ( i = 0; i < nd_psn; i++ )

ndltodp(nd_psh[i]->dl,(DP *)&BDY(hvect)[i]);

}

if ( !ishomo && homo ) {

/* dehomogenization */

for ( t = x; t; t = NEXT(t) ) ndv_dehomogenize((NDV)BDY(t),ord);

Line 3346 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

Line 3375 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

nd_demand = 0;

if ( nd_module && nd_intersect ) {

for ( j = nd_psn-1, x = 0; j >= 0; j-- )

if ( MPOS(DL(nd_psh[j])) > 1 ) {

if ( MPOS(DL(nd_psh[j])) > nd_intersect ) {

MKNODE(xx,(pointer)((unsigned long)j),x); x = xx;

}

conv_ilist(nd_demand,0,x,0);

Line 3370 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

Line 3399 void nd_gr(LIST f,LIST v,int m,int homo,int retdp,int

nd_setup_parameters(nd_nvar,0);

FINAL:

for ( r0 = 0, t = x; t; t = NEXT(t) ) {

NEXTNODE(r0,r);

if ( nd_module ) BDY(r) = ndvtopl(m,CO,vv,BDY(t),mrank);

if ( nd_module ) {

else if ( retdp ) BDY(r) = ndvtodp(m,BDY(t));

if ( retdp ) BDY(r) = ndvtodpm(m,BDY(t));

else BDY(r) = ndvtop(m,CO,vv,BDY(t));

else BDY(r) = ndvtopl(m,CO,vv,BDY(t),mrank);

} else if ( retdp ) BDY(r) = ndvtodp(m,BDY(t));

else BDY(r) = ndvtop(m,CO,vv,BDY(t));

}

if ( r0 ) NEXT(r) = 0;

if ( !m && nd_nalg )

Line 3381 FINAL:

Line 3412 FINAL:

MKLIST(*rp,r0);

if ( nd_gentrace ) {

if ( f4 ) {

STOQ(16,bpe);

STOZ(16,bpe);

STOQ(nd_last_nonzero,last_nonzero);

STOZ(nd_last_nonzero,last_nonzero);

tr = mknode(5,*rp,(!ishomo&&homo)?ONE:0,BDY(nzlist),bpe,last_nonzero); MKLIST(*rp,tr);

tr = mknode(6,*rp,(!ishomo&&homo)?ONE:0,BDY(nzlist),bpe,last_nonzero,hvect); MKLIST(*rp,tr);

} else {

tl1 = reverse_node(tl1); tl2 = reverse_node(tl2);

tl3 = reverse_node(tl3);

Line 3392 FINAL:

Line 3422 FINAL:

for ( t = tl2; t; t = NEXT(t) ) {

/* s = [i,[*,j,*,*],...] */

s = BDY((LIST)BDY(t));

j = perm[QTOS((Q)ARG0(s))]; STOQ(j,jq); ARG0(s) = (pointer)jq;

j = perm[ZTOS((Q)ARG0(s))]; STOZ(j,jq); ARG0(s) = (pointer)jq;

for ( s = BDY((LIST)ARG1(s)); s; s = NEXT(s) ) {

j = perm[QTOS((Q)ARG1(BDY((LIST)BDY(s))))]; STOQ(j,jq);

j = perm[ZTOS((Q)ARG1(BDY((LIST)BDY(s))))]; STOZ(j,jq);

ARG1(BDY((LIST)BDY(s))) = (pointer)jq;

}

for ( j = length(x)-1, t = 0; j >= 0; j-- ) {

STOQ(perm[j],jq); MKNODE(s,jq,t); t = s;

STOZ(perm[j],jq); MKNODE(s,jq,t); t = s;

}

MKLIST(l1,tl1); MKLIST(l2,tl2); MKLIST(l3,t); MKLIST(l4,tl3);

MKLIST(l5,tl4);

STOQ(nd_bpe,bpe);

STOZ(nd_bpe,bpe);

tr = mknode(8,*rp,(!ishomo&&homo)?ONE:0,l1,l2,l3,l4,l5,bpe); MKLIST(*rp,tr);

tr = mknode(9,*rp,(!ishomo&&homo)?ONE:0,l1,l2,l3,l4,l5,bpe,hvect); MKLIST(*rp,tr);

}

#if 0

Line 3512 NDV recompute_trace(NODE ti,NDV *p,int mod)

Line 3542 NDV recompute_trace(NODE ti,NDV *p,int mod)

NODE sj;

NDV red;

Obj mj;

static int afo=0;

afo++;

mul = (NM)MALLOC(sizeof(struct oNM)+(nd_wpd-1)*sizeof(UINT));

CM(mul) = 1;

tail = 0;

for ( i = 0, d = r = 0; ti; ti = NEXT(ti), i++ ) {

sj = BDY((LIST)BDY(ti));

if ( ARG0(sj) ) {

red = p[QTOS((Q)ARG1(sj))];

red = p[ZTOS((Q)ARG1(sj))];

mj = (Obj)ARG2(sj);

if ( OID(mj) != O_DP ) ndl_zero(DL(mul));

else dltondl(nd_nvar,BDY((DP)mj)->dl,DL(mul));

Line 3590 void nd_gr_recompute_trace(LIST f,LIST v,int m,struct

Line 3618 void nd_gr_recompute_trace(LIST f,LIST v,int m,struct

break;

}

nd_init_ord(ord);

nd_bpe = QTOS((Q)ARG7(BDY(tlist)));

nd_bpe = ZTOS((Q)ARG7(BDY(tlist)));

nd_setup_parameters(nvar,0);

len = length(BDY(f));

Line 3610 void nd_gr_recompute_trace(LIST f,LIST v,int m,struct

Line 3638 void nd_gr_recompute_trace(LIST f,LIST v,int m,struct

trace = NEXT(permtrace);

for ( i = length(perm)-1, t = trace; t; t = NEXT(t) ) {

j = QTOS((Q)ARG0(BDY((LIST)BDY(t))));

j = ZTOS((Q)ARG0(BDY((LIST)BDY(t))));

if ( j > i ) i = j;

}

n = i+1;

pb = (NDV *)MALLOC(n*sizeof(NDV *));

for ( t = perm, i = 0; t; t = NEXT(t), i++ ) {

ti = BDY((LIST)BDY(t));

pb[QTOS((Q)ARG0(ti))] = db[QTOS((Q)ARG1(ti))];

pb[ZTOS((Q)ARG0(ti))] = db[ZTOS((Q)ARG1(ti))];

}

for ( t = trace; t; t = NEXT(t) ) {

ti = BDY((LIST)BDY(t));

pb[QTOS((Q)ARG0(ti))] = recompute_trace(BDY((LIST)ARG1(ti)),pb,m);

pb[ZTOS((Q)ARG0(ti))] = recompute_trace(BDY((LIST)ARG1(ti)),pb,m);

if ( !pb[QTOS((Q)ARG0(ti))] ) { *rp = 0; return; }

if ( !pb[ZTOS((Q)ARG0(ti))] ) { *rp = 0; return; }

if ( DP_Print ) {

fprintf(asir_out,"."); fflush(asir_out);

}

for ( t = intred; t; t = NEXT(t) ) {

ti = BDY((LIST)BDY(t));

pb[QTOS((Q)ARG0(ti))] = recompute_trace(BDY((LIST)ARG1(ti)),pb,m);

pb[ZTOS((Q)ARG0(ti))] = recompute_trace(BDY((LIST)ARG1(ti)),pb,m);

if ( !pb[QTOS((Q)ARG0(ti))] ) { *rp = 0; return; }

if ( !pb[ZTOS((Q)ARG0(ti))] ) { *rp = 0; return; }

if ( DP_Print ) {

fprintf(asir_out,"*"); fflush(asir_out);

}

for ( r0 = 0, t = ind; t; t = NEXT(t) ) {

NEXTNODE(r0,r);

b = pb[QTOS((Q)BDY(t))];

b = pb[ZTOS((Q)BDY(t))];

ndv_mul_c(m,b,invm(HCM(b),m));

#if 0

BDY(r) = ndvtop(m,CO,vv,pb[QTOS((Q)BDY(t))]);

BDY(r) = ndvtop(m,CO,vv,pb[ZTOS((Q)BDY(t))]);

#else

BDY(r) = ndvtodp(m,pb[QTOS((Q)BDY(t))]);

BDY(r) = ndvtodp(m,pb[ZTOS((Q)BDY(t))]);

#endif

}

if ( r0 ) NEXT(r) = 0;

Line 3650 void nd_gr_recompute_trace(LIST f,LIST v,int m,struct

Line 3678 void nd_gr_recompute_trace(LIST f,LIST v,int m,struct

if ( DP_Print ) fprintf(asir_out,"\n");

}

void nd_gr_trace(LIST f,LIST v,int trace,int homo,int f4,struct order_spec *ord,LIST *rp)

void nd_gr_trace(LIST f,LIST v,int trace,int homo,int retdp,int f4,struct order_spec *ord,LIST *rp)

{

VL tv,fv,vv,vc,av;

NODE fd,fd0,in0,in,r,r0,t,s,cand,alist;

Line 3673 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3701 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

int *perm;

int j,ret;

Z jq,bpe;

VECT hvect;

nd_module = 0;

nd_lf = 0;

Line 3727 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3756 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

for ( t = BDY(f), max = 1; t; t = NEXT(t) )

for ( tv = vv; tv; tv = NEXT(tv) ) {

if ( nd_module ) {

if ( OID(BDY(t)) == O_DPM ) {

e = dpm_getdeg((DPM)BDY(t),&trank);

max = MAX(e,max);

mrank = MAX(mrank,trank);

} else {

s = BDY((LIST)BDY(t));

trank = length(s);

mrank = MAX(mrank,trank);

Line 3734 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3768 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

e = getdeg(tv->v,(P)BDY(s));

max = MAX(e,max);

}

} else {

e = getdeg(tv->v,(P)BDY(t));

max = MAX(e,max);

Line 3744 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3779 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

ishomo = 1;

for ( in0 = 0, fd0 = 0, t = BDY(f); t; t = NEXT(t) ) {

if ( nd_module ) {

if ( !nd_gentrace ) pltozpl((LIST)BDY(t),&dmy,&zpl);

if ( OID(BDY(t)) == O_DPM ) {

else zpl = (LIST)BDY(t);

Z cont;

DPM zdpm;

if ( !nd_gentrace ) dpm_ptozp((DPM)BDY(t),&cont,&zdpm);

else zdpm = (DPM)BDY(t);

c = (pointer)dpmtondv(m,zdpm);

} else {

if ( !nd_gentrace ) pltozpl((LIST)BDY(t),&dmy,&zpl);

else zpl = (LIST)BDY(t);

c = (pointer)pltondv(CO,vv,zpl);

}

} else {

if ( !nd_gentrace ) ptozp((P)BDY(t),1,&dmy,&zp);

else zp = (P)BDY(t);

c = (pointer)ptondv(CO,vv,zp);

}

if ( ishomo )

ishomo = ishomo && ndv_ishomo(c);

Line 3790 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3834 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

else m = get_lprime(++mindex);

continue;

}

if ( nd_gentrace ) {

MKVECT(hvect,nd_psn);

for ( i = 0; i < nd_psn; i++ )

ndltodp(nd_psh[i]->dl,(DP *)&BDY(hvect)[i]);

}

if ( !ishomo && homo ) {

/* dehomogenization */

for ( t = cand; t; t = NEXT(t) ) ndv_dehomogenize((NDV)BDY(t),ord);

Line 3837 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3886 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

}

get_eg(&eg1); init_eg(&eg_check); add_eg(&eg_check,&eg0,&eg1);

if ( DP_Print )

fprintf(asir_out,"check=%.3fsec,",eg_check.exectime+eg_check.gctime);

fprintf(asir_out,"check=%.3fsec\n",eg_check.exectime);

/* dp->p */

nd_bpe = cbpe;

nd_setup_parameters(nd_nvar,0);

for ( r = cand; r; r = NEXT(r) ) {

if ( nd_module ) BDY(r) = ndvtopl(0,CO,vv,BDY(r),mrank);

if ( nd_module ) {

else BDY(r) = (pointer)ndvtop(0,CO,vv,BDY(r));

if ( retdp ) BDY(r) = ndvtodpm(0,BDY(r));

else BDY(r) = ndvtopl(0,CO,vv,BDY(r),mrank);

} else if ( retdp ) BDY(r) = ndvtodp(0,BDY(r));

else BDY(r) = (pointer)ndvtop(0,CO,vv,BDY(r));

}

if ( nd_nalg )

cand = postprocess_algcoef(av,alist,cand);

Line 3855 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

Line 3907 void nd_gr_trace(LIST f,LIST v,int trace,int homo,int

for ( t = tl2; t; t = NEXT(t) ) {

/* s = [i,[*,j,*,*],...] */

s = BDY((LIST)BDY(t));

j = perm[QTOS((Q)ARG0(s))]; STOQ(j,jq); ARG0(s) = (pointer)jq;

j = perm[ZTOS((Q)ARG0(s))]; STOZ(j,jq); ARG0(s) = (pointer)jq;

for ( s = BDY((LIST)ARG1(s)); s; s = NEXT(s) ) {

j = perm[QTOS((Q)ARG1(BDY((LIST)BDY(s))))]; STOQ(j,jq);

j = perm[ZTOS((Q)ARG1(BDY((LIST)BDY(s))))]; STOZ(j,jq);

ARG1(BDY((LIST)BDY(s))) = (pointer)jq;

}

for ( j = length(cand)-1, t = 0; j >= 0; j-- ) {

STOQ(perm[j],jq); MKNODE(s,jq,t); t = s;

STOZ(perm[j],jq); MKNODE(s,jq,t); t = s;

}

MKLIST(l1,tl1); MKLIST(l2,tl2); MKLIST(l3,t); MKLIST(l4,tl3);

MKLIST(l5,tl4);

STOQ(nd_bpe,bpe);

STOZ(nd_bpe,bpe);

tr = mknode(8,*rp,(!ishomo&&homo)?ONE:0,l1,l2,l3,l4,l5,bpe); MKLIST(*rp,tr);

tr = mknode(9,*rp,(!ishomo&&homo)?ONE:0,l1,l2,l3,l4,l5,bpe,hvect); MKLIST(*rp,tr);

}

Line 3931 void nmtodp(int mod,NM m,DP *r)

Line 3983 void nmtodp(int mod,NM m,DP *r)

*r = dp;

}

void ndltodp(UINT *d,DP *r)

{

DP dp;

MP mr;

NEWMP(mr);

mr->dl = ndltodl(nd_nvar,d);

mr->c = (Obj)ONE;

NEXT(mr) = 0; MKDP(nd_nvar,mr,dp); dp->sugar = mr->dl->td;

*r = dp;

}

void ndl_print(UINT *dl)

{

int n;

Line 3980 void nd_print_q(ND p)

Line 4044 void nd_print_q(ND p)

else {

for ( m = BDY(p); m; m = NEXT(m) ) {

printf("+");

printexpr(CO,(Obj)CQ(m));

printexpr(CO,(Obj)CZ(m));

printf("*");

ndl_print(DL(m));

}

Line 4014 void nd_removecont(int mod,ND p)

Line 4078 void nd_removecont(int mod,ND p)

w = (Z *)MALLOC(n*sizeof(Q));

v.len = n;

v.body = (pointer *)w;

for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ ) w[i] = CQ(m);

for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ ) w[i] = CZ(m);

removecont_array((P *)w,n,1);

for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ ) CQ(m) = w[i];

for ( m = BDY(p), i = 0; i < n; m = NEXT(m), i++ ) CZ(m) = w[i];

}

Line 4035 void nd_removecont2(ND p1,ND p2)

Line 4099 void nd_removecont2(ND p1,ND p2)

v.body = (pointer *)w;

i = 0;

if ( p1 )

for ( m = BDY(p1); i < n1; m = NEXT(m), i++ ) w[i] = CQ(m);

for ( m = BDY(p1); i < n1; m = NEXT(m), i++ ) w[i] = CZ(m);

if ( p2 )

for ( m = BDY(p2); i < n; m = NEXT(m), i++ ) w[i] = CQ(m);

for ( m = BDY(p2); i < n; m = NEXT(m), i++ ) w[i] = CZ(m);

removecont_array((P *)w,n,1);

i = 0;

if ( p1 )

for ( m = BDY(p1); i < n1; m = NEXT(m), i++ ) CQ(m) = w[i];

for ( m = BDY(p1); i < n1; m = NEXT(m), i++ ) CZ(m) = w[i];

if ( p2 )

for ( m = BDY(p2); i < n; m = NEXT(m), i++ ) CQ(m) = w[i];

for ( m = BDY(p2); i < n; m = NEXT(m), i++ ) CZ(m) = w[i];

}

void ndv_removecont(int mod,NDV p)

Line 4095 void ndv_homogenize(NDV p,int obpe,int oadv,EPOS oepos

Line 4159 void ndv_homogenize(NDV p,int obpe,int oadv,EPOS oepos

NMV m,mr0,mr,t;

len = p->len;

for ( m = BDY(p), i = 0, max = 1; i < len; NMV_OADV(m), i++ )

for ( m = BDY(p), i = 0, max = 0; i < len; NMV_OADV(m), i++ )

max = MAX(max,TD(DL(m)));

mr0 = nmv_adv>oadv?(NMV)REALLOC(BDY(p),len*nmv_adv):BDY(p);

m = (NMV)((char *)mr0+(len-1)*oadv);

Line 4103 void ndv_homogenize(NDV p,int obpe,int oadv,EPOS oepos

Line 4167 void ndv_homogenize(NDV p,int obpe,int oadv,EPOS oepos

t = (NMV)MALLOC(nmv_adv);

for ( i = 0; i < len; i++, NMV_OPREV(m), NMV_PREV(mr) ) {

ndl_homogenize(DL(m),DL(t),obpe,oepos,ompos,max);

CQ(mr) = CQ(m);

CZ(mr) = CZ(m);

ndl_copy(DL(t),DL(mr));

}

NV(p)++;

Line 4128 void ndv_dehomogenize(NDV p,struct order_spec *ord)

Line 4192 void ndv_dehomogenize(NDV p,struct order_spec *ord)

if ( newwpd != nd_wpd ) {

newadv = ROUND_FOR_ALIGN(sizeof(struct oNMV)+(newwpd-1)*sizeof(UINT));

for ( m = r = BDY(p), i = 0; i < len; NMV_ADV(m), NDV_NADV(r), i++ ) {

CQ(r) = CQ(m);

CZ(r) = CZ(m);

if ( nd_module ) pos = MPOS(DL(m));

for ( j = 0; j < newexporigin; j++ ) DL(r)[j] = DL(m)[j];

adj = nd_exporigin-newexporigin;

Line 4219 void removecont_array_q(Z *c,int n)

Line 4283 void removecont_array_q(Z *c,int n)

v.id = O_VECT; v.len = n; v.body = (pointer *)r;

gcdvz(&v,&d1);

gcdz(d0,d1,&gcd);

divz(d0,gcd,&a);

/* exact division */

divsz(d0,gcd,&a);

for ( i = 0; i < n; i++ ) {

mulz(a,q[i],&u);

if ( r[i] ) {

divz(r[i],gcd,&u1);

/* exact division */

divsz(r[i],gcd,&u1);

addz(u,u1,&q[i]);

} else

q[i] = u;

Line 4232 void removecont_array_q(Z *c,int n)

Line 4298 void removecont_array_q(Z *c,int n)

for ( i = 0; i < n; i++ ) c[i] = q[i];

}

void gcdv_mpz_estimate(mpz_t d0,mpz_t *c,int n);

void mpz_removecont_array(mpz_t *c,int n)

{

mpz_t d0,a,u,u1,gcd;

int i,j;

static mpz_t *q,*r;

static int c_len = 0;

for ( i = 0; i < n; i++ )

if ( mpz_sgn(c[i]) ) break;

if ( i == n ) return;

gcdv_mpz_estimate(d0,c,n);

if ( n > c_len ) {

q = (mpz_t *)MALLOC(n*sizeof(mpz_t));

r = (mpz_t *)MALLOC(n*sizeof(mpz_t));

c_len = n;

}

for ( i = 0; i < n; i++ ) {

mpz_init(q[i]); mpz_init(r[i]);

mpz_fdiv_qr(q[i],r[i],c[i],d0);

}

for ( i = 0; i < n; i++ )

if ( mpz_sgn(r[i]) ) break;

mpz_init(gcd); mpz_init(a); mpz_init(u); mpz_init(u1);

if ( i < n ) {

mpz_gcd(gcd,d0,r[i]);

for ( j = i+1; j < n; j++ ) mpz_gcd(gcd,gcd,r[j]);

mpz_div(a,d0,gcd);

for ( i = 0; i < n; i++ ) {

mpz_mul(u,a,q[i]);

if ( mpz_sgn(r[i]) ) {

mpz_div(u1,r[i],gcd);

mpz_add(q[i],u,u1);

} else

mpz_set(q[i],u);

}

for ( i = 0; i < n; i++ )

mpz_set(c[i],q[i]);

}

void nd_mul_c(int mod,ND p,int mul)

{

NM m;

Line 4646 int nd_sp(int mod,int trace,ND_pairs p,ND *rp)

Line 4754 int nd_sp(int mod,int trace,ND_pairs p,ND *rp)

divsp(nd_vc,HCP(p2),gp,&CP(m1));

divsp(nd_vc,HCP(p1),gp,&tp); chsgnp(tp,&CP(m2));

} else {

igcd_cofactor(HCQ(p1),HCQ(p2),&g,&t,&CQ(m1)); chsgnz(t,&CQ(m2));

igcd_cofactor(HCZ(p1),HCZ(p2),&g,&t,&CZ(m1)); chsgnz(t,&CZ(m2));

}

t1 = ndv_mul_nm(mod,m1,p1); t2 = ndv_mul_nm(mod,m2,p2);

*rp = nd_add(mod,t1,t2);

if ( nd_gentrace ) {

/* nd_tracelist is initialized */

STOQ(p->i1,iq); nmtodp(mod,m1,&d); node = mknode(4,ONE,iq,d,ONE);

STOZ(p->i1,iq); nmtodp(mod,m1,&d); node = mknode(4,ONE,iq,d,ONE);

MKLIST(hist,node); MKNODE(nd_tracelist,hist,0);

STOQ(p->i2,iq); nmtodp(mod,m2,&d); node = mknode(4,ONE,iq,d,ONE);

STOZ(p->i2,iq); nmtodp(mod,m2,&d); node = mknode(4,ONE,iq,d,ONE);

MKLIST(hist,node); MKNODE(node,hist,nd_tracelist);

nd_tracelist = node;

}

Line 4701 void ndv_mul_c_q(NDV p,Z mul)

Line 4809 void ndv_mul_c_q(NDV p,Z mul)

if ( !p ) return;

len = LEN(p);

for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {

mulz(CQ(m),mul,&c); CQ(m) = c;

mulz(CZ(m),mul,&c); CZ(m) = c;

}

Line 4767 void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta

Line 4875 void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta

} else if ( nd_vc )

mulp(nd_vc,CP(m0),CP(m1),&CP(m));

else

mulz(CQ(m0),CQ(m1),&CQ(m));

mulz(CZ(m0),CZ(m1),&CZ(m));

for ( i = 0; i < nd_wpd; i++ ) d[i] = 0;

homo = n&1 ? 1 : 0;

if ( homo ) {

Line 4824 void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta

Line 4932 void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta

} else if ( nd_vc )

mulp(nd_vc,CP(tab[u]),(P)q,&CP(tab[u]));

else {

mulz(CQ(tab[u]),q,&q1); CQ(tab[u]) = q1;

mulz(CZ(tab[u]),q,&q1); CZ(tab[u]) = q1;

}

Line 4838 void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta

Line 4946 void weyl_mul_nm_nmv(int n,int mod,NM m0,NMV m1,NM *ta

} else if ( nd_vc )

mulp(nd_vc,CP(tab[u]),(P)q,&CP(t));

else

mulz(CQ(tab[u]),q,&CQ(t));

mulz(CZ(tab[u]),q,&CZ(t));

*p = t;

}

Line 4978 ND nd_quo(int mod,PGeoBucket bucket,NDV d)

Line 5086 ND nd_quo(int mod,PGeoBucket bucket,NDV d)

DMAR(c1,c2,0,mod,c); CM(mq) = c;

CM(tm) = mod-c;

} else {

divsz(HCQ(p),HCQ(d),&CQ(mq));

divsz(HCZ(p),HCZ(d),&CZ(mq));

chsgnz(CQ(mq),&CQ(tm));

chsgnz(CZ(mq),&CZ(tm));

}

t = ndv_mul_nmv_trunc(mod,tm,d,HDL(d));

bucket->body[hindex] = nd_remove_head(p);

Line 5011 void ndv_realloc(NDV p,int obpe,int oadv,EPOS oepos)

Line 5119 void ndv_realloc(NDV p,int obpe,int oadv,EPOS oepos)

mr = (NMV)((char *)mr0+(len-1)*nmv_adv);

t = (NMV)MALLOC(nmv_adv);

for ( i = 0; i < len; i++, NMV_OPREV(m), NMV_PREV(mr) ) {

CQ(t) = CQ(m);

CZ(t) = CZ(m);

for ( k = 0; k < nd_wpd; k++ ) DL(t)[k] = 0;

ndl_reconstruct(DL(m),DL(t),obpe,oepos);

CQ(mr) = CQ(t);

CZ(mr) = CZ(t);

ndl_copy(DL(t),DL(mr));

}

BDY(p) = mr0;

Line 5032 NDV ndv_dup_realloc(NDV p,int obpe,int oadv,EPOS oepos

Line 5140 NDV ndv_dup_realloc(NDV p,int obpe,int oadv,EPOS oepos

for ( i = 0; i < len; i++, NMV_OADV(m), NMV_ADV(mr) ) {

ndl_zero(DL(mr));

ndl_reconstruct(DL(m),DL(mr),obpe,oepos);

CQ(mr) = CQ(m);

CZ(mr) = CZ(m);

}

MKNDV(NV(p),mr0,len,r);

SG(r) = SG(p);

Line 5052 NDV ndv_dup(int mod,NDV p)

Line 5160 NDV ndv_dup(int mod,NDV p)

m0 = m = (NMV)((mod>0 || mod==-1)?MALLOC_ATOMIC(len*nmv_adv):MALLOC(len*nmv_adv));

for ( t = BDY(p), i = 0; i < len; i++, NMV_ADV(t), NMV_ADV(m) ) {

ndl_copy(DL(t),DL(m));

CQ(m) = CQ(t);

CZ(m) = CZ(t);

}

MKNDV(NV(p),m0,len,d);

SG(d) = SG(p);

Line 5070 NDV ndv_symbolic(int mod,NDV p)

Line 5178 NDV ndv_symbolic(int mod,NDV p)

m0 = m = (NMV)((mod>0||mod==-1)?MALLOC_ATOMIC(len*nmv_adv):MALLOC(len*nmv_adv));

for ( t = BDY(p), i = 0; i < len; i++, NMV_ADV(t), NMV_ADV(m) ) {

ndl_copy(DL(t),DL(m));

CQ(m) = ONE;

CZ(m) = ONE;

}

MKNDV(NV(p),m0,len,d);

SG(d) = SG(p);

Line 5086 ND nd_dup(ND p)

Line 5194 ND nd_dup(ND p)

for ( m0 = 0, t = BDY(p); t; t = NEXT(t) ) {

NEXTNM(m0,m);

ndl_copy(DL(t),DL(m));

CQ(m) = CQ(t);

CZ(m) = CZ(t);

}

if ( m0 ) NEXT(m) = 0;

MKND(NV(p),m0,LEN(p),d);

Line 5135 void ndv_mod(int mod,NDV p)

Line 5243 void ndv_mod(int mod,NDV p)

nd_subst_vector(nd_vc,CP(t),nd_subst,&cp);

c = (Z)cp;

} else

c = CQ(t);

c = CZ(t);

r = remqi((Q)c,mod);

if ( r ) {

CM(d) = r;

Line 5157 NDV ptondv(VL vl,VL dvl,P p)

Line 5265 NDV ptondv(VL vl,VL dvl,P p)

void pltozpl(LIST l,Q *cont,LIST *pp)

{

NODE nd,nd1;

int n;

P *pl;

Q *cl;

int i;

P dmy;

Z dvr;

Z dvr,inv;

LIST r;

nd = BDY(l); n = length(nd);

pl = (P *)MALLOC(n*sizeof(P));

cl = (Q *)MALLOC(n*sizeof(P));

cl = (Q *)MALLOC(n*sizeof(Q));

for ( i = 0; i < n; i++, nd = NEXT(nd) )

for ( i = 0; i < n; i++, nd = NEXT(nd) ) {

ptozp((P)BDY(nd),1,&cl[i],&dmy);

qltozl(cl,n,&dvr);

}

nd = BDY(l);

qltozl(cl,n,&dvr);

for ( i = 0; i < n; i++, nd = NEXT(nd) ) {

divz(ONE,dvr,&inv);

divsp(CO,(P)BDY(nd),(P)dvr,&pl[i]);

nd = BDY(l);

}

for ( i = 0; i < n; i++, nd = NEXT(nd) )

nd = 0;

divsp(CO,(P)BDY(nd),(P)dvr,&pl[i]);

for ( i = n-1; i >= 0; i-- ) {

nd = 0;

MKNODE(nd1,pl[i],nd); nd = nd1;

for ( i = n-1; i >= 0; i-- ) {

}

MKNODE(nd1,pl[i],nd); nd = nd1;

MKLIST(r,nd);

}

*pp = r;

MKLIST(r,nd);

*pp = r;

}

/* (a1,a2,...,an) -> a1*e(1)+...+an*e(n) */

Line 5228 ND ptond(VL vl,VL dvl,P p)

Line 5337 ND ptond(VL vl,VL dvl,P p)

ndl_zero(DL(m));

if ( !INT((Q)p) )

error("ptond : input must be integer-coefficient");

CQ(m) = (Z)p;

CZ(m) = (Z)p;

NEXT(m) = 0;

MKND(nd_nvar,m,1,r);

SG(r) = 0;

Line 5249 ND ptond(VL vl,VL dvl,P p)

Line 5358 ND ptond(VL vl,VL dvl,P p)

} else {

NEWNM(m0); d = DL(m0);

for ( j = k-1, s = 0; j >= 0; j-- ) {

ndl_zero(d); e = QTOS(DEG(w[j])); PUT_EXP(d,i,e);

ndl_zero(d); e = ZTOS(DEG(w[j])); PUT_EXP(d,i,e);

TD(d) = MUL_WEIGHT(e,i);

if ( nd_blockmask) ndl_weight_mask(d);

if ( nd_module ) MPOS(d) = 0;

Line 5287 P ndvtop(int mod,VL vl,VL dvl,NDV p)

Line 5396 P ndvtop(int mod,VL vl,VL dvl,NDV p)

} else if ( mod == -2 ) {

c = (P)CZ(m);

} else if ( mod > 0 ) {

STOQ(CM(m),q); c = (P)q;

STOZ(CM(m),q); c = (P)q;

} else

c = CP(m);

d = DL(m);

for ( i = 0, t = c, tvl = dvl; i < n; tvl = NEXT(tvl), i++ ) {

MKV(tvl->v,r); e = GET_EXP(d,i); STOQ(e,q);

MKV(tvl->v,r); e = GET_EXP(d,i); STOZ(e,q);

pwrp(vl,r,q,&u); mulp(vl,t,u,&w); t = w;

}

addp(vl,s,t,&u); s = u;

Line 5326 LIST ndvtopl(int mod,VL vl,VL dvl,NDV p,int rank)

Line 5435 LIST ndvtopl(int mod,VL vl,VL dvl,NDV p,int rank)

if ( mod == -1 ) {

e = IFTOF(CM(m)); MKGFS(e,gfs); c = (P)gfs;

} else if ( mod ) {

STOQ(CM(m),q); c = (P)q;

STOZ(CM(m),q); c = (P)q;

} else

c = CP(m);

d = DL(m);

for ( i = 0, t = c, tvl = dvl; i < n; tvl = NEXT(tvl), i++ ) {

MKV(tvl->v,r); e = GET_EXP(d,i); STOQ(e,q);

MKV(tvl->v,r); e = GET_EXP(d,i); STOZ(e,q);

pwrp(vl,r,q,&u); mulp(vl,t,u,&w); t = w;

}

addp(vl,a[MPOS(d)],t,&u); a[MPOS(d)] = u;

Line 5363 NDV ndtondv(int mod,ND p)

Line 5472 NDV ndtondv(int mod,ND p)

#endif

for ( t = BDY(p), i = 0; t; t = NEXT(t), i++, NMV_ADV(m) ) {

ndl_copy(DL(t),DL(m));

CQ(m) = CQ(t);

CZ(m) = CZ(t);

}

MKNDV(NV(p),m0,len,d);

SG(d) = SG(p);

return d;

}

static int dmm_comp_nv;

int dmm_comp(DMM *a,DMM *b)

{

return -compdmm(dmm_comp_nv,*a,*b);

}

void dmm_sort_by_ord(DMM *a,int len,int nv)

{

dmm_comp_nv = nv;

qsort(a,len,sizeof(DMM),(int (*)(const void *,const void *))dmm_comp);

}

void dpm_sort(DPM p,DPM *rp)

{

DMM t,t1;

int len,i,n;

DMM *a;

DPM d;

if ( !p ) *rp = 0;

for ( t = BDY(p), len = 0; t; t = NEXT(t), len++ );

a = (DMM *)MALLOC(len*sizeof(DMM));

for ( i = 0, t = BDY(p); i < len; i++, t = NEXT(t) ) a[i] = t;

n = p->nv;

dmm_sort_by_ord(a,len,n);

t = 0;

for ( i = len-1; i >= 0; i-- ) {

NEWDMM(t1);

t1->c = a[i]->c;

t1->dl = a[i]->dl;

t1->pos = a[i]->pos;

t1->next = t;

t = t1;

}

MKDPM(n,t,d);

SG(d) = SG(p);

*rp = d;

}

int dpm_comp(DPM *a,DPM *b)

{

return compdpm(CO,*a,*b);

}

NODE dpm_sort_list(NODE l)

{

int i,len;

NODE t,t1;

DPM *a;

len = length(l);

a = (DPM *)MALLOC(len*sizeof(DPM));

for ( t = l, i = 0; i < len; i++, t = NEXT(t) ) a[i] = (DPM)BDY(t);

qsort(a,len,sizeof(DPM),(int (*)(const void *,const void *))dpm_comp);

t = 0;

for ( i = len-1; i >= 0; i-- ) {

MKNODE(t1,(pointer)a[i],t); t = t1;

}

return t;

}

int nmv_comp(NMV a,NMV b)

{

return -DL_COMPARE(a->dl,b->dl);

}

NDV dpmtondv(int mod,DPM p)

{

NDV d;

NMV m,m0;

DMM t;

DMM *a;

int i,len,n;

if ( !p ) return 0;

for ( t = BDY(p), len = 0; t; t = NEXT(t), len++ );

a = (DMM *)MALLOC(len*sizeof(DMM));

for ( i = 0, t = BDY(p); i < len; i++, t = NEXT(t) ) a[i] = t;

n = p->nv;

dmm_sort_by_ord(a,len,n);

if ( mod > 0 || mod == -1 )

m0 = m = (NMV)MALLOC_ATOMIC_IGNORE_OFF_PAGE(len*nmv_adv);

else

m0 = m = MALLOC(len*nmv_adv);

#if 0

ndv_alloc += nmv_adv*len;

#endif

for ( i = 0; i < len; i++, NMV_ADV(m) ) {

dltondl(n,a[i]->dl,DL(m));

MPOS(DL(m)) = a[i]->pos;

TD(DL(m)) = ndl_weight(DL(m));

CZ(m) = (Z)a[i]->c;

}

qsort(m0,len,nmv_adv,(int (*)(const void *,const void *))nmv_comp);

MKNDV(NV(p),m0,len,d);

SG(d) = SG(p);

return d;

}

ND ndvtond(int mod,NDV p)

{

ND d;

Line 5383 ND ndvtond(int mod,NDV p)

Line 5592 ND ndvtond(int mod,NDV p)

for ( t = BDY(p), i = 0; i < len; NMV_ADV(t), i++ ) {

NEXTNM(m0,m);

ndl_copy(DL(t),DL(m));

CQ(m) = CQ(t);

CZ(m) = CZ(t);

}

NEXT(m) = 0;

MKND(NV(p),m0,len,d);

Line 5412 DP ndvtodp(int mod,NDV p)

Line 5621 DP ndvtodp(int mod,NDV p)

return d;

}

DPM ndvtodpm(int mod,NDV p)

{

DMM m,m0;

DPM d;

NMV t;

int i,len;

if ( !p ) return 0;

m0 = 0;

len = p->len;

for ( t = BDY(p), i = 0; i < len; NMV_ADV(t), i++ ) {

NEXTDMM(m0,m);

m->dl = ndltodl(nd_nvar,DL(t));

m->c = (Obj)ndctop(mod,t->c);

m->pos = MPOS(DL(t));

}

NEXT(m) = 0;

MKDPM(nd_nvar,m0,d);

SG(d) = SG(p);

return d;

}

DP ndtodp(int mod,ND p)

{

MP m,m0;

Line 5460 void ndv_print_q(NDV p)

Line 5692 void ndv_print_q(NDV p)

len = LEN(p);

for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {

printf("+");

printexpr(CO,(Obj)CQ(m));

printexpr(CO,(Obj)CZ(m));

printf("*");

ndl_print(DL(m));

}

Line 5499 NODE ndv_reducebase(NODE x,int *perm)

Line 5731 NODE ndv_reducebase(NODE x,int *perm)

/* XXX incomplete */

extern int dpm_ordtype;

void nd_init_ord(struct order_spec *ord)

{

nd_module = (ord->id >= 256);

Line 5510 void nd_init_ord(struct order_spec *ord)

Line 5744 void nd_init_ord(struct order_spec *ord)

nd_poly_weight = ord->top_weight;

nd_module_rank = ord->module_rank;

nd_module_weight = ord->module_top_weight;

dpm_ordtype = ord->ispot;

}

nd_matrix = 0;

nd_matrix_len = 0;

Line 5742 void nd_nf_p(Obj f,LIST g,LIST v,int m,struct order_sp

Line 5977 void nd_nf_p(Obj f,LIST g,LIST v,int m,struct order_sp

/* dont sort, dont removecont */

ndv_setup(m,0,in0,1,1);

nd_scale=2;

stat = nd_nf(m,0,ndf,nd_ps,1,0,&nf);

stat = nd_nf(m,0,ndf,nd_ps,1,&nf);

if ( !stat )

error("nd_nf_p : exponent too large");

if ( nd_module ) *rp = (Obj)ndvtopl(m,CO,vv,ndtondv(m,nf),mrank);

Line 5765 int nd_to_vect(int mod,UINT *s0,int n,ND d,UINT *r)

Line 6000 int nd_to_vect(int mod,UINT *s0,int n,ND d,UINT *r)

return i;

}

#if defined(__GNUC__) && SIZEOF_LONG==8

#define MOD128(a,c,m) ((a)=(((c)!=0||((a)>=(m)))?(((((U128)(c))<<64)+(a))%(m)):(a)))

int nd_to_vect64(int mod,UINT *s0,int n,ND d,U64 *r)

{

NM m;

UINT *t,*s;

int i;

for ( i = 0; i < n; i++ ) r[i] = 0;

for ( i = 0, s = s0, m = BDY(d); m; m = NEXT(m) ) {

t = DL(m);

for ( ; !ndl_equal(t,s); s += nd_wpd, i++ );

r[i] = (U64)CM(m);

}

for ( i = 0; !r[i]; i++ );

return i;

}

#endif

int nd_to_vect_q(UINT *s0,int n,ND d,Z *r)

{

NM m;

Line 5796 int nd_to_vect_q(UINT *s0,int n,ND d,Z *r)

Line 6010 int nd_to_vect_q(UINT *s0,int n,ND d,Z *r)

for ( i = 0, s = s0, m = BDY(d); m; m = NEXT(m) ) {

t = DL(m);

for ( ; !ndl_equal(t,s); s += nd_wpd, i++ );

r[i] = CQ(m);

r[i] = CZ(m);

}

for ( i = 0; !r[i]; i++ );

return i;

Line 5874 Z *nm_ind_pair_to_vect(int mod,UINT *s0,int n,NM_ind_p

Line 6088 Z *nm_ind_pair_to_vect(int mod,UINT *s0,int n,NM_ind_p

for ( i = j = 0, s = s0, mr = BDY(p); j < len; j++, NMV_ADV(mr) ) {

ndl_add(d,DL(mr),t);

for ( ; !ndl_equal(t,s); s += nd_wpd, i++ );

r[i] = CQ(mr);

r[i] = CZ(mr);

}

return r;

}

IndArray nm_ind_pair_to_vect_compress(int trace,UINT *s0,int n,int *s0hash,NM_ind_pair pair)

IndArray nm_ind_pair_to_vect_compress(int trace,UINT *s0,int n,NM_ind_pair pair,int start)

{

NM m;

NMV mr;

UINT *d,*t,*s;

UINT *d,*t,*s,*u;

NDV p;

unsigned char *ivc;

unsigned short *ivs;

UINT *v,*ivi,*s0v;

int i,j,len,prev,diff,cdiff,h;

int i,j,len,prev,diff,cdiff,h,st,ed,md,c;

IndArray r;

struct oEGT eg0,eg1;

m = pair->mul;

d = DL(m);

Line 5902 struct oEGT eg0,eg1;

Line 6115 struct oEGT eg0,eg1;

len = LEN(p);

t = (UINT *)MALLOC(nd_wpd*sizeof(UINT));

v = (unsigned int *)MALLOC(len*sizeof(unsigned int));

get_eg(&eg0);

for ( prev = start, mr = BDY(p), j = 0; j < len; j++, NMV_ADV(mr) ) {

for ( i = j = 0, s = s0, mr = BDY(p); j < len; j++, NMV_ADV(mr) ) {

ndl_add(d,DL(mr),t);

st = prev;

h = ndl_hash_value(t);

ed = n;

for ( ; h != s0hash[i] || !ndl_equal(t,s); s += nd_wpd, i++ );

while ( ed > st ) {

v[j] = i;

md = (st+ed)/2;

u = s0+md*nd_wpd;

c = DL_COMPARE(u,t);

if ( c == 0 ) break;

else if ( c > 0 ) st = md;

else ed = md;

}

prev = v[j] = md;

}

get_eg(&eg1); add_eg(&eg_search,&eg0,&eg1);

r = (IndArray)MALLOC(sizeof(struct oIndArray));

r->head = v[0];

diff = 0;

Line 5952 void expand_array(Z *svect,Z *cvect,int n)

Line 6171 void expand_array(Z *svect,Z *cvect,int n)

if ( svect[i] ) svect[i] = cvect[j++];

}

#if 0

int ndv_reduce_vect_q(Z *svect,int trace,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

{

int i,j,k,len,pos,prev,nz;

Line 5966 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

Line 6186 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

int maxrs;

double hmag;

Z *cvect;

int l;

maxrs = 0;

for ( i = 0; i < col && !svect[i]; i++ );

Line 5980 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

Line 6201 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

redv = nd_demand?ndv_load(rp0[i]->index)

:(trace?nd_ps_trace[rp0[i]->index]:nd_ps[rp0[i]->index]);

len = LEN(redv); mr = BDY(redv);

igcd_cofactor(svect[k],CQ(mr),&gcd,&cs,&cr);

igcd_cofactor(svect[k],CZ(mr),&gcd,&cs,&cr);

chsgnz(cs,&mcs);

if ( !UNIQ(cr) ) {

for ( j = 0; j < col; j++ ) {

Line 5993 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

Line 6214 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

ivc = ivect->index.c;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivc[j]; prev = pos;

mulz(CQ(mr),mcs,&c2); addz(svect[pos],c2,&t); svect[pos] = t;

muladdtoz(CZ(mr),mcs,&svect[pos]);

}

break;

case 2:

ivs = ivect->index.s;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivs[j]; prev = pos;

mulz(CQ(mr),mcs,&c2); addz(svect[pos],c2,&t); svect[pos] = t;

muladdtoz(CZ(mr),mcs,&svect[pos]);

}

break;

case 4:

ivi = ivect->index.i;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivi[j]; prev = pos;

mulz(CQ(mr),mcs,&c2); addz(svect[pos],c2,&t); svect[pos] = t;

muladdtoz(CZ(mr),mcs,&svect[pos]);

}

break;

}

Line 6029 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

Line 6250 int ndv_reduce_vect_q(Z *svect,int trace,int col,IndAr

}

return maxrs;

}

#else

int ndv_reduce_vect(int m,UINT *svect,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

/* direct mpz version */

int ndv_reduce_vect_q(Z *svect0,int trace,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

{

int i,j,k,len,pos,prev;

UINT c,c1,c2,c3,up,lo,dmy;

mpz_t cs,cr,gcd;

IndArray ivect;

unsigned char *ivc;

unsigned short *ivs;

Line 6042 int ndv_reduce_vect(int m,UINT *svect,int col,IndArray

Line 6265 int ndv_reduce_vect(int m,UINT *svect,int col,IndArray

NMV mr;

NODE rp;

int maxrs;

double hmag;

int l;

static mpz_t *svect;

static int svect_len=0;

maxrs = 0;

for ( i = 0; i < col && !svect0[i]; i++ );

if ( i == col ) return maxrs;

hmag = p_mag((P)svect0[i])*nd_scale;

if ( col > svect_len ) {

svect = (mpz_t *)MALLOC(col*sizeof(mpz_t));

svect_len = col;

}

for ( i = 0; i < col; i++ ) {

mpz_init(svect[i]);

if ( svect0[i] )

mpz_set(svect[i],BDY(svect0[i]));

else

mpz_set_ui(svect[i],0);

}

mpz_init(gcd); mpz_init(cs); mpz_init(cr);

for ( i = 0; i < nred; i++ ) {

ivect = imat[i];

k = ivect->head; svect[k] %= m;

k = ivect->head;

if ( (c = svect[k]) != 0 ) {

if ( mpz_sgn(svect[k]) ) {

maxrs = MAX(maxrs,rp0[i]->sugar);

c = m-c; redv = nd_ps[rp0[i]->index];

redv = nd_demand?ndv_load(rp0[i]->index)

:(trace?nd_ps_trace[rp0[i]->index]:nd_ps[rp0[i]->index]);

len = LEN(redv); mr = BDY(redv);

svect[k] = 0; prev = k;

mpz_gcd(gcd,svect[k],BDY(CZ(mr)));

mpz_div(cs,svect[k],gcd);

mpz_div(cr,BDY(CZ(mr)),gcd);

mpz_neg(cs,cs);

if ( MUNIMPZ(cr) )

for ( j = 0; j < col; j++ ) mpz_neg(svect[j],svect[j]);

else if ( !UNIMPZ(cr) )

for ( j = 0; j < col; j++ ) {

if ( mpz_sgn(svect[j]) ) mpz_mul(svect[j],svect[j],cr);

}

mpz_set_ui(svect[k],0);

prev = k;

switch ( ivect->width ) {

case 1:

ivc = ivect->index.c;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivc[j]; c1 = CM(mr); prev = pos;

pos = prev+ivc[j]; prev = pos;

if ( c1 ) {

mpz_addmul(svect[pos],BDY(CZ(mr)),cs);

c2 = svect[pos];

DMA(c1,c,c2,up,lo);

if ( up ) { DSAB(m,up,lo,dmy,c3); svect[pos] = c3;

} else svect[pos] = lo;

}

break;

case 2:

ivs = ivect->index.s;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivs[j]; c1 = CM(mr);

pos = prev+ivs[j]; prev = pos;

prev = pos;

mpz_addmul(svect[pos],BDY(CZ(mr)),cs);

if ( c1 ) {

c2 = svect[pos];

DMA(c1,c,c2,up,lo);

if ( up ) { DSAB(m,up,lo,dmy,c3); svect[pos] = c3;

} else svect[pos] = lo;

}

break;

case 4:

ivi = ivect->index.i;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivi[j]; c1 = CM(mr);

pos = prev+ivi[j]; prev = pos;

prev = pos;

mpz_addmul(svect[pos],BDY(CZ(mr)),cs);

if ( c1 ) {

c2 = svect[pos];

DMA(c1,c,c2,up,lo);

if ( up ) { DSAB(m,up,lo,dmy,c3); svect[pos] = c3;

} else svect[pos] = lo;

}

break;

}

for ( j = k+1; j < col && !svect[j]; j++ );

if ( j == col ) break;

if ( hmag && ((double)mpz_sizeinbase(svect[j],2) > hmag) ) {

mpz_removecont_array(svect,col);

hmag = ((double)mpz_sizeinbase(svect[j],2))*nd_scale;

}

mpz_removecont_array(svect,col);

if ( DP_Print ) {

fprintf(asir_out,"-"); fflush(asir_out);

}

for ( i = 0; i < col; i++ )

if ( svect[i] >= (UINT)m ) svect[i] %= m;

if ( mpz_sgn(svect[i]) ) MPZTOZ(svect[i],svect0[i]);

else svect0[i] = 0;

return maxrs;

}

#endif

#if defined(__GNUC__) && SIZEOF_LONG==8

int ndv_reduce_vect(int m,UINT *svect,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

int ndv_reduce_vect64(int m,U64 *svect,U64 *cvect,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

{

int i,j,k,len,pos,prev;

U64 a,c,c1,c2;

UINT c,c1,c2,c3,up,lo,dmy;

IndArray ivect;

unsigned char *ivc;

unsigned short *ivs;

Line 6114 int ndv_reduce_vect64(int m,U64 *svect,U64 *cvect,int

Line 6361 int ndv_reduce_vect64(int m,U64 *svect,U64 *cvect,int

NODE rp;

int maxrs;

for ( i = 0; i < col; i++ ) cvect[i] = 0;

maxrs = 0;

for ( i = 0; i < nred; i++ ) {

ivect = imat[i];

k = ivect->head;

k = ivect->head; svect[k] %= m;

a = svect[k]; c = cvect[k];

MOD128(a,c,m);

svect[k] = a; cvect[k] = 0;

if ( (c = svect[k]) != 0 ) {

maxrs = MAX(maxrs,rp0[i]->sugar);

c = m-c; redv = nd_ps[rp0[i]->index];

Line 6132 int ndv_reduce_vect64(int m,U64 *svect,U64 *cvect,int

Line 6375 int ndv_reduce_vect64(int m,U64 *svect,U64 *cvect,int

ivc = ivect->index.c;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivc[j]; c1 = CM(mr); prev = pos;

if ( c1 ) {

c2 = svect[pos]+c1*c;

c2 = svect[pos];

if ( c2 < svect[pos] ) cvect[pos]++;

DMA(c1,c,c2,up,lo);

svect[pos] = c2;

if ( up ) { DSAB(m,up,lo,dmy,c3); svect[pos] = c3;

}

} else svect[pos] = lo;

}

break;

case 2:

ivs = ivect->index.s;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivs[j]; c1 = CM(mr); prev = pos;

pos = prev+ivs[j]; c1 = CM(mr);

if ( c1 ) {

prev = pos;

c2 = svect[pos]+c1*c;

if ( c1 ) {

if ( c2 < svect[pos] ) cvect[pos]++;

c2 = svect[pos];

svect[pos] = c2;

DMA(c1,c,c2,up,lo);

}

if ( up ) { DSAB(m,up,lo,dmy,c3); svect[pos] = c3;

} else svect[pos] = lo;

}

break;

case 4:

ivi = ivect->index.i;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivi[j]; c1 = CM(mr); prev = pos;

pos = prev+ivi[j]; c1 = CM(mr);

if ( c1 ) {

prev = pos;

c2 = svect[pos]+c1*c;

if ( c1 ) {

if ( c2 < svect[pos] ) cvect[pos]++;

c2 = svect[pos];

svect[pos] = c2;

DMA(c1,c,c2,up,lo);

}

if ( up ) { DSAB(m,up,lo,dmy,c3); svect[pos] = c3;

} else svect[pos] = lo;

}

break;

}

for ( i = 0; i < col; i++ ) {

for ( i = 0; i < col; i++ )

a = svect[i]; c = cvect[i]; MOD128(a,c,m); svect[i] = a;

if ( svect[i] >= (UINT)m ) svect[i] %= m;

}

return maxrs;

}

#endif

int ndv_reduce_vect_sf(int m,UINT *svect,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

{

Line 6428 NDV vect_to_ndv(UINT *vect,int spcol,int col,int *rhea

Line 6674 NDV vect_to_ndv(UINT *vect,int spcol,int col,int *rhea

}

#if defined(__GNUC__) && SIZEOF_LONG==8

NDV vect64_to_ndv(U64 *vect,int spcol,int col,int *rhead,UINT *s0vect)

{

int j,k,len;

UINT *p;

UINT c;

NDV r;

NMV mr0,mr;

for ( j = 0, len = 0; j < spcol; j++ ) if ( vect[j] ) len++;

if ( !len ) return 0;

else {

mr0 = (NMV)MALLOC_ATOMIC_IGNORE_OFF_PAGE(nmv_adv*len);

#if 0

ndv_alloc += nmv_adv*len;

#endif

mr = mr0;

p = s0vect;

for ( j = k = 0; j < col; j++, p += nd_wpd )

if ( !rhead[j] ) {

if ( (c = (UINT)vect[k++]) != 0 ) {

ndl_copy(p,DL(mr)); CM(mr) = c; NMV_ADV(mr);

}

MKNDV(nd_nvar,mr0,len,r);

return r;

}

#endif

NDV vect_to_ndv_2(unsigned long *vect,int col,UINT *s0vect)

{

int j,k,len;

Line 6484 NDV vect_to_ndv_2(unsigned long *vect,int col,UINT *s0

Line 6700 NDV vect_to_ndv_2(unsigned long *vect,int col,UINT *s0

NDV vect_to_ndv_q(Z *vect,int spcol,int col,int *rhead,UINT *s0vect)

{

int j,k,len;

UINT *p;

Z c;

NDV r;

NMV mr0,mr;

for ( j = 0, len = 0; j < spcol; j++ ) if ( vect[j] ) len++;

if ( !len ) return 0;

else {

mr0 = (NMV)MALLOC(nmv_adv*len);

#if 0

ndv_alloc += nmv_adv*len;

#endif

mr = mr0;

p = s0vect;

for ( j = k = 0; j < col; j++, p += nd_wpd )

for ( j = k = 0; j < col; j++, p += nd_wpd ) {

if ( !rhead[j] ) {

if ( (c = vect[k++]) != 0 ) {

if ( !INT(c) )

error("afo");

error("vect_to_ndv_q : components must be integers");

ndl_copy(p,DL(mr)); CQ(mr) = c; NMV_ADV(mr);

ndl_copy(p,DL(mr)); CZ(mr) = c; NMV_ADV(mr);

}

MKNDV(nd_nvar,mr0,len,r);

return r;

}

MKNDV(nd_nvar,mr0,len,r);

return r;

}

NDV vect_to_ndv_lf(mpz_t *vect,int spcol,int col,int *rhead,UINT *s0vect)

Line 6563 NDV plain_vect_to_ndv_q(Z *vect,int col,UINT *s0vect)

Line 6780 NDV plain_vect_to_ndv_q(Z *vect,int col,UINT *s0vect)

for ( j = k = 0; j < col; j++, p += nd_wpd, k++ )

if ( (c = vect[k]) != 0 ) {

if ( !INT(c) )

error("afo");

error("plain_vect_to_ndv_q : components must be integers");

ndl_copy(p,DL(mr)); CQ(mr) = c; NMV_ADV(mr);

ndl_copy(p,DL(mr)); CZ(mr) = c; NMV_ADV(mr);

}

MKNDV(nd_nvar,mr0,len,r);

return r;

Line 6598 int nd_symbolic_preproc(PGeoBucket bucket,int trace,UI

Line 6815 int nd_symbolic_preproc(PGeoBucket bucket,int trace,UI

NM_ind_pair pair;

ND red;

NDV *ps;

static int afo;

s0 = 0; rp0 = 0; col = 0;

if ( nd_demand )

Line 6667 NODE nd_f4(int m,int checkonly,int **indp)

Line 6883 NODE nd_f4(int m,int checkonly,int **indp)

PGeoBucket bucket;

struct oEGT eg0,eg1,eg_f4;

Z i1,i2,sugarq;

init_eg(&f4_symb); init_eg(&f4_conv); init_eg(&f4_conv); init_eg(&f4_elim1); init_eg(&f4_elim2);

#if 0

ndv_alloc = 0;

#endif

Nf4_red=0;

g = 0; d = 0;

for ( i = 0; i < nd_psn; i++ ) {

d = update_pairs(d,g,i,0);

Line 6686 NODE nd_f4(int m,int checkonly,int **indp)

Line 6905 NODE nd_f4(int m,int checkonly,int **indp)

if ( MaxDeg > 0 && sugar > MaxDeg ) break;

if ( nzlist_t ) {

node = BDY((LIST)BDY(nzlist_t));

sugarh = QTOS((Q)ARG0(node));

sugarh = ZTOS((Q)ARG0(node));

tn = BDY((LIST)ARG1(node));

if ( !tn ) {

nzlist_t = NEXT(nzlist_t);

Line 6711 NODE nd_f4(int m,int checkonly,int **indp)

Line 6930 NODE nd_f4(int m,int checkonly,int **indp)

d = nd_reconstruct(0,d);

continue;

}

get_eg(&eg1); init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg1);

get_eg(&eg1); init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg1); add_eg(&f4_symb,&eg0,&eg1);

if ( DP_Print )

fprintf(asir_out,"sugar=%d,symb=%.3fsec,",

sugar,eg_f4.exectime+eg_f4.gctime);

sugar,eg_f4.exectime);

nflist = nd_f4_red(m,nd_nzlist?lh:l,0,s0vect,col,rp0,nd_gentrace?&ll:0);

if ( checkonly && nflist ) return 0;

/* adding new bases */

Line 6740 NODE nd_f4(int m,int checkonly,int **indp)

Line 6959 NODE nd_f4(int m,int checkonly,int **indp)

if ( nd_gentrace ) {

for ( t = ll, tn0 = 0; t; t = NEXT(t) ) {

NEXTNODE(tn0,tn);

STOQ(t->i1,i1); STOQ(t->i2,i2);

STOZ(t->i1,i1); STOZ(t->i2,i2);

node = mknode(2,i1,i2); MKLIST(l0,node);

BDY(tn) = l0;

}

if ( tn0 ) NEXT(tn) = 0; MKLIST(l0,tn0);

STOQ(sugar,sugarq); node = mknode(2,sugarq,l0); MKLIST(l1,node);

STOZ(sugar,sugarq); node = mknode(2,sugarq,l0); MKLIST(l1,node);

MKNODE(node,l1,nzlist); nzlist = node;

}

if ( nd_nzlist ) nzlist_t = NEXT(nzlist_t);

Line 6760 NODE nd_f4(int m,int checkonly,int **indp)

Line 6979 NODE nd_f4(int m,int checkonly,int **indp)

#if 0

fprintf(asir_out,"ndv_alloc=%d\n",ndv_alloc);

#endif

if ( DP_Print ) {

fprintf(asir_out,"number of red=%d,",Nf4_red);

fprintf(asir_out,"symb=%.3fsec,conv=%.3fsec,elim1=%.3fsec,elim2=%.3fsec\n",

f4_symb.exectime,f4_conv.exectime,f4_elim1.exectime,f4_elim2.exectime);

}

conv_ilist(nd_demand,0,g,indp);

return g;

}

Line 6813 NODE nd_f4_trace(int m,int **indp)

Line 7037 NODE nd_f4_trace(int m,int **indp)

}

get_eg(&eg1); init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg1);

if ( DP_Print )

fprintf(asir_out,"sugar=%d,symb=%.3fsec,",

fprintf(asir_out,"\nsugar=%d,symb=%.3fsec,",

sugar,eg_f4.exectime+eg_f4.gctime);

sugar,eg_f4.exectime);

nflist = nd_f4_red(m,l,0,s0vect,col,rp0,&l0);

if ( !l0 ) continue;

l = l0;

Line 6841 NODE nd_f4_trace(int m,int **indp)

Line 7065 NODE nd_f4_trace(int m,int **indp)

for ( r = nflist; r; r = NEXT(r) ) {

nfqv = (NDV)BDY(r);

ndv_removecont(0,nfqv);

if ( !remqi((Q)HCQ(nfqv),m) ) return 0;

if ( !remqi((Q)HCZ(nfqv),m) ) return 0;

if ( nd_nalg ) {

ND nf1;

Line 6983 NODE nd_f4_red_2(ND_pairs sp0,UINT *s0vect,int col,NOD

Line 7207 NODE nd_f4_red_2(ND_pairs sp0,UINT *s0vect,int col,NOD

unsigned long *v;

get_eg(&eg0);

init_eg(&eg_search);

for ( sp = sp0, nsp = 0; sp; sp = NEXT(sp), nsp++ );

nred = length(rp0);

mat = alloc_matrix(nsp,col);

Line 7028 init_eg(&eg_search);

Line 7251 init_eg(&eg_search);

init_eg(&eg_elim2); add_eg(&eg_elim2,&eg1,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,elim2=%.3fsec,",

eg_elim1.exectime+eg_elim1.gctime,eg_elim2.exectime+eg_elim2.gctime);

eg_elim1.exectime,eg_elim2.exectime);

fflush(asir_out);

}

return r0;

Line 7038 init_eg(&eg_search);

Line 7261 init_eg(&eg_search);

NODE nd_f4_red(int m,ND_pairs sp0,int trace,UINT *s0vect,int col,NODE rp0,ND_pairs *nz)

{

IndArray *imat;

int nsp,nred,i;

int nsp,nred,i,start;

int *rhead;

NODE r0,rp;

ND_pairs sp;

NM_ind_pair *rvect;

UINT *s;

int *s0hash;

struct oEGT eg0,eg1,eg_conv;

if ( m == 2 && nd_rref2 )

return nd_f4_red_2(sp0,s0vect,col,rp0,nz);

init_eg(&eg_search);

for ( sp = sp0, nsp = 0; sp; sp = NEXT(sp), nsp++ );

nred = length(rp0);

imat = (IndArray *)MALLOC(nred*sizeof(IndArray));

Line 7057 init_eg(&eg_search);

Line 7280 init_eg(&eg_search);

for ( i = 0; i < col; i++ ) rhead[i] = 0;

/* construction of index arrays */

get_eg(&eg0);

if ( DP_Print ) {

fprintf(asir_out,"%dx%d,",nsp+nred,col);

fflush(asir_out);

}

rvect = (NM_ind_pair *)MALLOC(nred*sizeof(NM_ind_pair));

s0hash = (int *)MALLOC(col*sizeof(int));

for ( start = 0, rp = rp0, i = 0; rp; i++, rp = NEXT(rp) ) {

for ( i = 0, s = s0vect; i < col; i++, s += nd_wpd )

s0hash[i] = ndl_hash_value(s);

for ( rp = rp0, i = 0; rp; i++, rp = NEXT(rp) ) {

rvect[i] = (NM_ind_pair)BDY(rp);

imat[i] = nm_ind_pair_to_vect_compress(trace,s0vect,col,s0hash,rvect[i]);

imat[i] = nm_ind_pair_to_vect_compress(trace,s0vect,col,rvect[i],start);

rhead[imat[i]->head] = 1;

start = imat[i]->head;

}

get_eg(&eg1); init_eg(&eg_conv); add_eg(&eg_conv,&eg0,&eg1); add_eg(&f4_conv,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"conv=%.3fsec,",eg_conv.exectime);

fflush(asir_out);

}

if ( m > 0 )

#if defined(__GNUC__) && SIZEOF_LONG==8

#if SIZEOF_LONG==8

r0 = nd_f4_red_mod64_main(m,sp0,nsp,s0vect,col,rvect,rhead,imat,nred,nz);

#else

r0 = nd_f4_red_main(m,sp0,nsp,s0vect,col,rvect,rhead,imat,nred,nz);

Line 7081 init_eg(&eg_search);

Line 7309 init_eg(&eg_search);

r0 = nd_f4_red_lf_main(m,sp0,nsp,trace,s0vect,col,rvect,rhead,imat,nred);

else

r0 = nd_f4_red_q_main(sp0,nsp,trace,s0vect,col,rvect,rhead,imat,nred);

#if 0

if ( DP_Print ) print_eg("search",&eg_search);

#endif

return r0;

}

Line 7134 NODE nd_f4_red_main(int m,ND_pairs sp0,int nsp,UINT *s

Line 7359 NODE nd_f4_red_main(int m,ND_pairs sp0,int nsp,UINT *s

}

get_eg(&eg1); init_eg(&eg_f4_1); add_eg(&eg_f4_1,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime+eg_f4_1.gctime);

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime);

fflush(asir_out);

}

/* free index arrays */

Line 7159 NODE nd_f4_red_main(int m,ND_pairs sp0,int nsp,UINT *s

Line 7384 NODE nd_f4_red_main(int m,ND_pairs sp0,int nsp,UINT *s

get_eg(&eg2); init_eg(&eg_f4_2); add_eg(&eg_f4_2,&eg1,&eg2);

init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime+eg_f4_2.gctime);

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime);

fprintf(asir_out,"nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ",

nsp,nred,sprow,spcol,rank);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime+eg_f4.gctime);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime);

}

if ( nz ) {

for ( i = 0; i < rank-1; i++ ) NEXT(spactive[i]) = spactive[i+1];

Line 7175 NODE nd_f4_red_main(int m,ND_pairs sp0,int nsp,UINT *s

Line 7400 NODE nd_f4_red_main(int m,ND_pairs sp0,int nsp,UINT *s

return r0;

}

#if defined(__GNUC__) && SIZEOF_LONG==8

/* for Fp, 2^15=<p<2^29 */

NODE nd_f4_red_mod64_main(int m,ND_pairs sp0,int nsp,UINT *s0vect,int col,

NM_ind_pair *rvect,int *rhead,IndArray *imat,int nred,ND_pairs *nz)

{

int spcol,sprow,a;

int i,j,k,l,rank;

NODE r0,r;

ND_pairs sp;

ND spol;

U64 **spmat;

U64 *svect,*cvect;

U64 *v;

int *colstat;

struct oEGT eg0,eg1,eg2,eg_f4,eg_f4_1,eg_f4_2;

int maxrs;

int *spsugar;

ND_pairs *spactive;

spcol = col-nred;

get_eg(&eg0);

/* elimination (1st step) */

spmat = (U64 **)MALLOC(nsp*sizeof(U64 *));

svect = (U64 *)MALLOC(col*sizeof(U64));

cvect = (U64 *)MALLOC(col*sizeof(U64));

spsugar = (int *)MALLOC(nsp*sizeof(int));

spactive = !nz?0:(ND_pairs *)MALLOC(nsp*sizeof(ND_pairs));

for ( a = sprow = 0, sp = sp0; a < nsp; a++, sp = NEXT(sp) ) {

nd_sp(m,0,sp,&spol);

if ( !spol ) continue;

nd_to_vect64(m,s0vect,col,spol,svect);

maxrs = ndv_reduce_vect64(m,svect,cvect,col,imat,rvect,nred);

for ( i = 0; i < col; i++ ) if ( svect[i] ) break;

if ( i < col ) {

spmat[sprow] = v = (U64 *)MALLOC_ATOMIC(spcol*sizeof(U64));

for ( j = k = 0; j < col; j++ )

if ( !rhead[j] ) v[k++] = (UINT)svect[j];

spsugar[sprow] = MAX(maxrs,SG(spol));

if ( nz )

spactive[sprow] = sp;

sprow++;

}

nd_free(spol);

}

get_eg(&eg1); init_eg(&eg_f4_1); add_eg(&eg_f4_1,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime+eg_f4_1.gctime);

fflush(asir_out);

}

/* free index arrays */

for ( i = 0; i < nred; i++ ) GCFREE(imat[i]->index.c);

/* elimination (2nd step) */

colstat = (int *)MALLOC(spcol*sizeof(int));

rank = nd_gauss_elim_mod64(spmat,spsugar,spactive,sprow,spcol,m,colstat);

r0 = 0;

for ( i = 0; i < rank; i++ ) {

NEXTNODE(r0,r); BDY(r) =

(pointer)vect64_to_ndv(spmat[i],spcol,col,rhead,s0vect);

SG((NDV)BDY(r)) = spsugar[i];

GCFREE(spmat[i]);

}

if ( r0 ) NEXT(r) = 0;

for ( ; i < sprow; i++ ) GCFREE(spmat[i]);

get_eg(&eg2); init_eg(&eg_f4_2); add_eg(&eg_f4_2,&eg1,&eg2);

init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime+eg_f4_2.gctime);

fprintf(asir_out,"nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ",

nsp,nred,sprow,spcol,rank);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime+eg_f4.gctime);

}

if ( nz ) {

for ( i = 0; i < rank-1; i++ ) NEXT(spactive[i]) = spactive[i+1];

if ( rank > 0 ) {

NEXT(spactive[rank-1]) = 0;

*nz = spactive[0];

} else

*nz = 0;

}

return r0;

}

#endif

/* for small finite fields */

NODE nd_f4_red_sf_main(int m,ND_pairs sp0,int nsp,UINT *s0vect,int col,

Line 7305 NODE nd_f4_red_sf_main(int m,ND_pairs sp0,int nsp,UINT

Line 7445 NODE nd_f4_red_sf_main(int m,ND_pairs sp0,int nsp,UINT

}

get_eg(&eg1); init_eg(&eg_f4_1); add_eg(&eg_f4_1,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime+eg_f4_1.gctime);

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime);

fflush(asir_out);

}

/* free index arrays */

Line 7327 NODE nd_f4_red_sf_main(int m,ND_pairs sp0,int nsp,UINT

Line 7467 NODE nd_f4_red_sf_main(int m,ND_pairs sp0,int nsp,UINT

get_eg(&eg2); init_eg(&eg_f4_2); add_eg(&eg_f4_2,&eg1,&eg2);

init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime+eg_f4_2.gctime);

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime);

fprintf(asir_out,"nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ",

nsp,nred,sprow,spcol,rank);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime+eg_f4.gctime);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime);

}

if ( nz ) {

for ( i = 0; i < rank-1; i++ ) NEXT(spactive[i]) = spactive[i+1];

Line 7382 NODE nd_f4_red_lf_main(int m,ND_pairs sp0,int nsp,int

Line 7522 NODE nd_f4_red_lf_main(int m,ND_pairs sp0,int nsp,int

}

get_eg(&eg1); init_eg(&eg_f4_1); add_eg(&eg_f4_1,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime+eg_f4_1.gctime);

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime);

fflush(asir_out);

}

/* free index arrays */

Line 7417 NODE nd_f4_red_lf_main(int m,ND_pairs sp0,int nsp,int

Line 7557 NODE nd_f4_red_lf_main(int m,ND_pairs sp0,int nsp,int

get_eg(&eg2); init_eg(&eg_f4_2); add_eg(&eg_f4_2,&eg1,&eg2);

init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime+eg_f4_2.gctime);

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime);

fprintf(asir_out,"nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ",

nsp,nred,sprow,spcol,rank);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime+eg_f4.gctime);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime);

}

return r0;

}

Line 7464 NODE nd_f4_red_q_main(ND_pairs sp0,int nsp,int trace,U

Line 7604 NODE nd_f4_red_q_main(ND_pairs sp0,int nsp,int trace,U

}

get_eg(&eg1); init_eg(&eg_f4_1); add_eg(&eg_f4_1,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime+eg_f4_1.gctime);

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime);

fflush(asir_out);

}

/* free index arrays */

Line 7498 NODE nd_f4_red_q_main(ND_pairs sp0,int nsp,int trace,U

Line 7638 NODE nd_f4_red_q_main(ND_pairs sp0,int nsp,int trace,U

get_eg(&eg2); init_eg(&eg_f4_2); add_eg(&eg_f4_2,&eg1,&eg2);

init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime+eg_f4_2.gctime);

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime);

fprintf(asir_out,"nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ",

nsp,nred,sprow,spcol,rank);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime+eg_f4.gctime);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime);

}

return r0;

}

Line 7672 int nd_gauss_elim_mod(UINT **mat0,int *sugar,ND_pairs

Line 7812 int nd_gauss_elim_mod(UINT **mat0,int *sugar,ND_pairs

return rank;

}

#if defined(__GNUC__) && SIZEOF_LONG==8

int nd_gauss_elim_mod64(U64 **mat,int *sugar,ND_pairs *spactive,int row,int col,int md,int *colstat)

{

int i,j,k,l,rank,s;

U64 inv;

U64 a;

UINT c;

U64 *t,*pivot,*pk;

UINT *ck;

UINT **cmat;

UINT *ct;

ND_pairs pair;

cmat = (UINT **)MALLOC(row*sizeof(UINT *));

for ( i = 0; i < row; i++ ) {

cmat[i] = MALLOC_ATOMIC(col*sizeof(UINT));

bzero(cmat[i],col*sizeof(UINT));

}

for ( rank = 0, j = 0; j < col; j++ ) {

for ( i = rank; i < row; i++ ) {

a = mat[i][j]; c = cmat[i][j];

MOD128(a,c,md);

mat[i][j] = a; cmat[i][j] = 0;

}

for ( i = rank; i < row; i++ )

if ( mat[i][j] )

break;

if ( i == row ) {

colstat[j] = 0;

continue;

} else

colstat[j] = 1;

if ( i != rank ) {

t = mat[i]; mat[i] = mat[rank]; mat[rank] = t;

ct = cmat[i]; cmat[i] = cmat[rank]; cmat[rank] = ct;

s = sugar[i]; sugar[i] = sugar[rank]; sugar[rank] = s;

if ( spactive ) {

pair = spactive[i]; spactive[i] = spactive[rank];

spactive[rank] = pair;

}

/* column j is normalized */

s = sugar[rank];

inv = invm((UINT)mat[rank][j],md);

/* normalize pivot row */

for ( k = j, pk = mat[rank]+j, ck = cmat[rank]+j; k < col; k++, pk++, ck++ ) {

a = *pk; c = *ck; MOD128(a,c,md); *pk = (a*inv)%md; *ck = 0;

}

for ( i = rank+1; i < row; i++ ) {

if ( (a = mat[i][j]) != 0 ) {

sugar[i] = MAX(sugar[i],s);

red_by_vect64(md,mat[i]+j,cmat[i]+j,mat[rank]+j,(int)(md-a),col-j);

}

rank++;

}

for ( j = col-1, l = rank-1; j >= 0; j-- )

if ( colstat[j] ) {

for ( k = j, pk = mat[l]+j, ck = cmat[l]+j; k < col; k++, pk++, ck++ ) {

a = *pk; c = *ck; MOD128(a,c,md); *pk = a; *ck = 0;

}

s = sugar[l];

for ( i = 0; i < l; i++ ) {

a = mat[i][j]; c = cmat[i][j]; MOD128(a,c,md); mat[i][j] = a; cmat[i][j] = 0;

if ( a ) {

sugar[i] = MAX(sugar[i],s);

red_by_vect64(md,mat[i]+j,cmat[i]+j,mat[l]+j,(int)(md-a),col-j);

}

l--;

}

for ( i = 0; i < row; i++ ) GCFREE(cmat[i]);

GCFREE(cmat);

return rank;

}

#endif

int nd_gauss_elim_sf(UINT **mat0,int *sugar,int row,int col,int md,int *colstat)

{

int i,j,k,l,inv,a,rank,s;

Line 7813 int ndv_ishomo(NDV p)

Line 7875 int ndv_ishomo(NDV p)

h = TD(DL(m));

NMV_ADV(m);

for ( len--; len; len--, NMV_ADV(m) )

if ( TD(DL(m)) != h ) return 0;

if ( TD(DL(m)) != h ) {

return 0;

}

return 1;

}

Line 7842 void ndv_save(NDV p,int index)

Line 7906 void ndv_save(NDV p,int index)

write_int(s,(unsigned int *)&len);

for ( m = BDY(p), i = 0; i < len; i++, NMV_ADV(m) ) {

saveobj(s,(Obj)CQ(m));

saveobj(s,(Obj)CZ(m));

dl = DL(m);

td = TD(dl);

write_int(s,(unsigned int *)&td);

Line 7908 NDV ndv_load(int index)

Line 7972 NDV ndv_load(int index)

m0 = m = MALLOC(len*nmv_adv);

for ( i = 0; i < len; i++, NMV_ADV(m) ) {

loadobj(s,&obj); CQ(m) = (Z)obj;

loadobj(s,&obj); CZ(m) = (Z)obj;

dl = DL(m);

ndl_zero(dl);

read_int(s,(unsigned int *)&td); TD(dl) = td;

Line 8007 void nd_det(int mod,MAT f,P *rp)

Line 8071 void nd_det(int mod,MAT f,P *rp)

for ( j = 0; j < n; j++ ) {

if ( !m[i][j] ) continue;

lgp(m[i][j],&nm,&dn);

gcdz(dn0,dn,&gn); divz(dn0,gn,&qn); mulz(qn,dn,&dn0);

gcdz(dn0,dn,&gn); divsz(dn0,gn,&qn); mulz(qn,dn,&dn0);

}

if ( !UNIZ(dn0) ) {

ds = dn0;

Line 8144 ND ndv_mul_nmv_trunc(int mod,NMV m0,NDV p,UINT *d)

Line 8208 ND ndv_mul_nmv_trunc(int mod,NMV m0,NDV p,UINT *d)

}

} else {

q = CQ(m0);

q = CZ(m0);

for ( i = 0; i < len; i++, NMV_ADV(m) ) {

ndl_add(DL(m),d0,DL(tnm));

if ( ndl_reducible(DL(tnm),d) ) {

NEXTNM(mr0,mr);

mulz(CQ(m),q,&CQ(mr));

mulz(CZ(m),q,&CZ(mr));

ndl_copy(DL(tnm),DL(mr));

}

Line 8278 int nd_monic(int mod,ND *p)

Line 8342 int nd_monic(int mod,ND *p)

is_lc = 1;

while ( 1 ) {

NEWMP(mp0); mp = mp0;

mp->c = (Obj)CQ(m);

mp->c = (Obj)CZ(m);

mp->dl = nd_separate_d(DL(m),DL(ma));

NEWNM(mb);

for ( m = NEXT(m); m; m = NEXT(m) ) {

alg = nd_separate_d(DL(m),DL(mb));

if ( !ndl_equal(DL(ma),DL(mb)) )

break;

NEXTMP(mp0,mp); mp->c = (Obj)CQ(m); mp->dl = alg;

NEXTMP(mp0,mp); mp->c = (Obj)CZ(m); mp->dl = alg;

}

NEXT(mp) = 0;

MKDP(nd_nalg,mp0,nm);

Line 8318 int nd_monic(int mod,ND *p)

Line 8382 int nd_monic(int mod,ND *p)

/* l = lcm(denoms) */

l = ln;

for ( mr0 = 0, m = ma0; m; m = NEXT(m) ) {

divz(l,CA(m)->dn,&mul);

divsz(l,CA(m)->dn,&mul);

for ( mp = BDY(CA(m)->nm); mp; mp = NEXT(mp) ) {

NEXTNM(mr0,mr);

mulz((Z)mp->c,mul,&CQ(mr));

mulz((Z)mp->c,mul,&CZ(mr));

dl = mp->dl;

td = TD(DL(m));

ndl_copy(DL(m),DL(mr));

Line 8361 P ndc_div(int mod,union oNDC a,union oNDC b)

Line 8425 P ndc_div(int mod,union oNDC a,union oNDC b)

int inv,t;

if ( mod == -1 ) c.m = _mulsf(a.m,_invsf(b.m));

else if ( mod == -2 ) divlf(a.gz,b.gz,&c.gz);

else if ( mod == -2 ) divlf(a.z,b.z,&c.z);

else if ( mod ) {

inv = invm(b.m,mod);

DMAR(a.m,inv,0,mod,t); c.m = t;

Line 8381 P ndctop(int mod,union oNDC c)

Line 8445 P ndctop(int mod,union oNDC c)

if ( mod == -1 ) {

e = IFTOF(c.m); MKGFS(e,gfs); return (P)gfs;

} else if ( mod == -2 ) {

q = c.gz; return (P)q;

q = c.z; return (P)q;

} else if ( mod > 0 ) {

STOQ(c.m,q); return (P)q;

STOZ(c.m,q); return (P)q;

} else

return (P)c.p;

}

Line 8401 void finalize_tracelist(int i,P cont)

Line 8465 void finalize_tracelist(int i,P cont)

MKLIST(l,node); MKNODE(node,l,nd_tracelist);

nd_tracelist = node;

}

STOQ(i,iq);

STOZ(i,iq);

nd_tracelist = reverse_node(nd_tracelist);

MKLIST(l,nd_tracelist);

node = mknode(2,iq,l); MKLIST(l,node);

Line 8453 void parse_nd_option(NODE opt)

Line 8517 void parse_nd_option(NODE opt)

nd_gbblock = MALLOC((2*length(u)+1)*sizeof(int));

for ( i = 0; u; u = NEXT(u) ) {

p = BDY((LIST)BDY(u));

s = nd_gbblock[i++] = QTOS((Q)BDY(p));

s = nd_gbblock[i++] = ZTOS((Q)BDY(p));

nd_gbblock[i++] = s+QTOS((Q)BDY(NEXT(p)))-1;

nd_gbblock[i++] = s+ZTOS((Q)BDY(NEXT(p)))-1;

}

nd_gbblock[i] = -1;

} else

Line 8463 void parse_nd_option(NODE opt)

Line 8527 void parse_nd_option(NODE opt)

nd_newelim = value?1:0;

else if ( !strcmp(key,"intersect") )

nd_intersect = value?1:0;

else if ( !strcmp(key,"syzgen") )

nd_intersect = ZTOS((Q)value);

else if ( !strcmp(key,"lf") )

nd_lf = value?1:0;

else if ( !strcmp(key,"trace") ) {

if ( value ) {

u = BDY((LIST)value);

nd_nzlist = BDY((LIST)ARG2(u));

nd_bpe = QTOS((Q)ARG3(u));

nd_bpe = ZTOS((Q)ARG3(u));

}

} else if ( !strcmp(key,"f4red") ) {

nd_f4red = QTOS((Q)value);

nd_f4red = ZTOS((Q)value);

} else if ( !strcmp(key,"rank0") ) {

nd_rank0 = value?1:0;

} else if ( !strcmp(key,"splist") ) {

Line 8484 void parse_nd_option(NODE opt)

Line 8550 void parse_nd_option(NODE opt)

n = length(u);

nd_sugarweight = MALLOC(n*sizeof(int));

for ( i = 0; i < n; i++, u = NEXT(u) )

nd_sugarweight[i] = QTOS((Q)BDY(u));

nd_sugarweight[i] = ZTOS((Q)BDY(u));

}

Line 8508 ND mdptond(DP d)

Line 8574 ND mdptond(DP d)

else {

NEWNM(m);

dltondl(NV(d),BDY(d)->dl,DL(m));

CQ(m) = (Z)BDY(d)->c;

CZ(m) = (Z)BDY(d)->c;

NEXT(m) = 0;

MKND(NV(d),m,1,r);

}

Line 8573 ND *btog(NODE ti,ND **p,int nb,int mod)

Line 8639 ND *btog(NODE ti,ND **p,int nb,int mod)

s = BDY((LIST)BDY(t));

if ( ARG0(s) ) {

m = mdptond((DP)ARG2(s));

ptomp(mod,(P)HCQ(m),&c);

ptomp(mod,(P)HCZ(m),&c);

if ( (ci = ((MQ)c)->cont) != 0 ) {

HCM(m) = ci;

pi = p[QTOS((Q)ARG1(s))];

pi = p[ZTOS((Q)ARG1(s))];

for ( i = 0; i < nb; i++ ) {

tp = nd_mul_nm(mod,BDY(m),pi[i]);

add_pbucket(mod,r[i],tp);

Line 8614 ND *btog_lf(NODE ti,ND **p,int nb)

Line 8680 ND *btog_lf(NODE ti,ND **p,int nb)

s = BDY((LIST)BDY(t));

if ( ARG0(s) ) {

m = mdptond((DP)ARG2(s));

simp_ff((Obj)HCQ(m),(Obj *)&lm);

simp_ff((Obj)HCZ(m),(Obj *)&lm);

if ( lm ) {

lmtolf(lm,&lf); HCZ(m) = lf;

pi = p[QTOS((Q)ARG1(s))];

pi = p[ZTOS((Q)ARG1(s))];

for ( i = 0; i < nb; i++ ) {

tp = nd_mul_nm_lf(BDY(m),pi[i]);

add_pbucket(-2,r[i],tp);

Line 8649 ND btog_one(NODE ti,ND *p,int nb,int mod)

Line 8715 ND btog_one(NODE ti,ND *p,int nb,int mod)

s = BDY((LIST)BDY(t));

if ( ARG0(s) ) {

m = mdptond((DP)ARG2(s));

ptomp(mod,(P)HCQ(m),&c);

ptomp(mod,(P)HCZ(m),&c);

if ( (ci = ((MQ)c)->cont) != 0 ) {

HCM(m) = ci;

pi = p[j=QTOS((Q)ARG1(s))];

pi = p[j=ZTOS((Q)ARG1(s))];

if ( !pi ) {

pi = nd_load_mod(j);

tp = nd_mul_nm(mod,BDY(m),pi);

Line 8704 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

Line 8770 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

}

nd_init_ord(ord);

#if 0

nd_bpe = QTOS((Q)ARG7(BDY(tlist)));

nd_bpe = ZTOS((Q)ARG7(BDY(tlist)));

#else

nd_bpe = 32;

#endif

Line 8714 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

Line 8780 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

ind = BDY((LIST)ARG4(BDY(tlist)));

perm = BDY((LIST)BDY(permtrace)); trace =NEXT(permtrace);

for ( i = length(perm)-1, t = trace; t; t = NEXT(t) ) {

j = QTOS((Q)BDY(BDY((LIST)BDY(t))));

j = ZTOS((Q)BDY(BDY((LIST)BDY(t))));

if ( j > i ) i = j;

}

n = i+1;

Line 8722 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

Line 8788 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

p = (ND **)MALLOC(n*sizeof(ND *));

for ( t = perm, i = 0; t; t = NEXT(t), i++ ) {

pi = BDY((LIST)BDY(t));

pi0 = QTOS((Q)ARG0(pi)); pi1 = QTOS((Q)ARG1(pi));

pi0 = ZTOS((Q)ARG0(pi)); pi1 = ZTOS((Q)ARG1(pi));

p[pi0] = c = (ND *)MALLOC(nb*sizeof(ND));

ptomp(mod,(P)ARG2(pi),&inv);

((MQ)inv)->cont = invm(((MQ)inv)->cont,mod);

Line 8733 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

Line 8799 MAT nd_btog(LIST f,LIST v,int mod,struct order_spec *o

for ( t = trace,i=0; t; t = NEXT(t), i++ ) {

printf("%d ",i); fflush(stdout);

ti = BDY((LIST)BDY(t));

p[j=QTOS((Q)ARG0(ti))] = btog(BDY((LIST)ARG1(ti)),p,nb,mod);

p[j=ZTOS((Q)ARG0(ti))] = btog(BDY((LIST)ARG1(ti)),p,nb,mod);

}

for ( t = intred, i=0; t; t = NEXT(t), i++ ) {

printf("%d ",i); fflush(stdout);

ti = BDY((LIST)BDY(t));

p[j=QTOS((Q)ARG0(ti))] = btog(BDY((LIST)ARG1(ti)),p,nb,mod);

p[j=ZTOS((Q)ARG0(ti))] = btog(BDY((LIST)ARG1(ti)),p,nb,mod);

}

m = length(ind);

MKMAT(mat,nb,m);

for ( j = 0, t = ind; j < m; j++, t = NEXT(t) )

for ( i = 0, c = p[QTOS((Q)BDY(t))]; i < nb; i++ )

for ( i = 0, c = p[ZTOS((Q)BDY(t))]; i < nb; i++ )

BDY(mat)[i][j] = ndtodp(mod,c[i]);

return mat;

}

Line 8773 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

Line 8839 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

}

nd_init_ord(ord);

#if 0

nd_bpe = QTOS((Q)ARG7(BDY(tlist)));

nd_bpe = ZTOS((Q)ARG7(BDY(tlist)));

#else

nd_bpe = 32;

#endif

Line 8783 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

Line 8849 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

ind = BDY((LIST)ARG4(BDY(tlist)));

perm = BDY((LIST)BDY(permtrace)); trace =NEXT(permtrace);

for ( i = length(perm)-1, t = trace; t; t = NEXT(t) ) {

j = QTOS((Q)BDY(BDY((LIST)BDY(t))));

j = ZTOS((Q)BDY(BDY((LIST)BDY(t))));

if ( j > i ) i = j;

}

n = i+1;

Line 8791 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

Line 8857 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

p = (ND **)MALLOC(n*sizeof(ND *));

for ( t = perm, i = 0; t; t = NEXT(t), i++ ) {

pi = BDY((LIST)BDY(t));

pi0 = QTOS((Q)ARG0(pi)); pi1 = QTOS((Q)ARG1(pi));

pi0 = ZTOS((Q)ARG0(pi)); pi1 = ZTOS((Q)ARG1(pi));

p[pi0] = c = (ND *)MALLOC(nb*sizeof(ND));

simp_ff((Obj)ARG2(pi),(Obj *)&lm); lmtolf(lm,&lf); invz(lf,current_mod_lf,&inv);

u = ptond(CO,vv,(P)ONE);

Line 8801 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

Line 8867 MAT nd_btog_lf(LIST f,LIST v,struct order_spec *ord,LI

for ( t = trace,i=0; t; t = NEXT(t), i++ ) {

printf("%d ",i); fflush(stdout);

ti = BDY((LIST)BDY(t));

p[j=QTOS((Q)ARG0(ti))] = btog_lf(BDY((LIST)ARG1(ti)),p,nb);

p[j=ZTOS((Q)ARG0(ti))] = btog_lf(BDY((LIST)ARG1(ti)),p,nb);

}

for ( t = intred, i=0; t; t = NEXT(t), i++ ) {

printf("%d ",i); fflush(stdout);

ti = BDY((LIST)BDY(t));

p[j=QTOS((Q)ARG0(ti))] = btog_lf(BDY((LIST)ARG1(ti)),p,nb);

p[j=ZTOS((Q)ARG0(ti))] = btog_lf(BDY((LIST)ARG1(ti)),p,nb);

}

m = length(ind);

MKMAT(mat,nb,m);

for ( j = 0, t = ind; j < m; j++, t = NEXT(t) )

for ( i = 0, c = p[QTOS((Q)BDY(t))]; i < nb; i++ )

for ( i = 0, c = p[ZTOS((Q)BDY(t))]; i < nb; i++ )

BDY(mat)[i][j] = ndtodp(-2,c[i]);

return mat;

}

Line 8844 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

Line 8910 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

}

nd_init_ord(ord);

#if 0

nd_bpe = QTOS((Q)ARG7(BDY(tlist)));

nd_bpe = ZTOS((Q)ARG7(BDY(tlist)));

#else

nd_bpe = 32;

#endif

Line 8854 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

Line 8920 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

ind = BDY((LIST)ARG4(BDY(tlist)));

perm = BDY((LIST)BDY(permtrace)); trace =NEXT(permtrace);

for ( i = length(perm)-1, t = trace; t; t = NEXT(t) ) {

j = QTOS((Q)BDY(BDY((LIST)BDY(t))));

j = ZTOS((Q)BDY(BDY((LIST)BDY(t))));

if ( j > i ) i = j;

}

n = i+1;

Line 8862 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

Line 8928 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

p = (ND *)MALLOC(n*sizeof(ND *));

for ( t = perm, i = 0; t; t = NEXT(t), i++ ) {

pi = BDY((LIST)BDY(t));

pi0 = QTOS((Q)ARG0(pi)); pi1 = QTOS((Q)ARG1(pi));

pi0 = ZTOS((Q)ARG0(pi)); pi1 = ZTOS((Q)ARG1(pi));

if ( pi1 == pos ) {

ptomp(mod,(P)ARG2(pi),&inv);

((MQ)inv)->cont = invm(((MQ)inv)->cont,mod);

Line 8874 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

Line 8940 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

for ( t = trace,i=0; t; t = NEXT(t), i++ ) {

printf("%d ",i); fflush(stdout);

ti = BDY((LIST)BDY(t));

p[j=QTOS((Q)ARG0(ti))] = btog_one(BDY((LIST)ARG1(ti)),p,nb,mod);

p[j=ZTOS((Q)ARG0(ti))] = btog_one(BDY((LIST)ARG1(ti)),p,nb,mod);

if ( Demand ) {

nd_save_mod(p[j],j); nd_free(p[j]); p[j] = 0;

}

Line 8882 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

Line 8948 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

for ( t = intred, i=0; t; t = NEXT(t), i++ ) {

printf("%d ",i); fflush(stdout);

ti = BDY((LIST)BDY(t));

p[j=QTOS((Q)ARG0(ti))] = btog_one(BDY((LIST)ARG1(ti)),p,nb,mod);

p[j=ZTOS((Q)ARG0(ti))] = btog_one(BDY((LIST)ARG1(ti)),p,nb,mod);

if ( Demand ) {

nd_save_mod(p[j],j); nd_free(p[j]); p[j] = 0;

}

Line 8890 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

Line 8956 VECT nd_btog_one(LIST f,LIST v,int mod,struct order_sp

m = length(ind);

MKVECT(vect,m);

for ( j = 0, t = ind; j < m; j++, t = NEXT(t) ) {

u = p[QTOS((Q)BDY(t))];

u = p[ZTOS((Q)BDY(t))];

if ( !u ) {

u = nd_load_mod(QTOS((Q)BDY(t)));

u = nd_load_mod(ZTOS((Q)BDY(t)));

BDY(vect)[j] = ndtodp(mod,u);

nd_free(u);

} else

Line 9044 void nd_f4_lf_trace(LIST f,LIST v,int trace,int homo,s

Line 9110 void nd_f4_lf_trace(LIST f,LIST v,int trace,int homo,s

}

get_eg(&eg1); init_eg(&eg_check); add_eg(&eg_check,&eg0,&eg1);

if ( DP_Print )

fprintf(asir_out,"check=%.3fsec\n",eg_check.exectime+eg_check.gctime);

fprintf(asir_out,"check=%.3fsec\n",eg_check.exectime);

/* dp->p */

nd_bpe = cbpe;

nd_setup_parameters(nd_nvar,0);

Line 9104 NODE nd_f4_lf_trace_main(int m,int **indp)

Line 9170 NODE nd_f4_lf_trace_main(int m,int **indp)

}

get_eg(&eg1); init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg1);

if ( DP_Print )

fprintf(asir_out,"sugar=%d,symb=%.3fsec,",

fprintf(asir_out,"\nsugar=%d,symb=%.3fsec,",sugar,eg_f4.exectime);

sugar,eg_f4.exectime+eg_f4.gctime);

nflist = nd_f4_red(m,l,0,s0vect,col,rp0,&l0);

if ( !l0 ) continue;

l = l0;

Line 9144 NODE nd_f4_lf_trace_main(int m,int **indp)

Line 9209 NODE nd_f4_lf_trace_main(int m,int **indp)

conv_ilist(nd_demand,1,g,indp);

return g;

}

#if SIZEOF_LONG==8

NDV vect64_to_ndv(mp_limb_t *vect,int spcol,int col,int *rhead,UINT *s0vect)

{

int j,k,len;

UINT *p;

UINT c;

NDV r;

NMV mr0,mr;

for ( j = 0, len = 0; j < spcol; j++ ) if ( vect[j] ) len++;

if ( !len ) return 0;

else {

mr0 = (NMV)MALLOC_ATOMIC_IGNORE_OFF_PAGE(nmv_adv*len);

#if 0

ndv_alloc += nmv_adv*len;

#endif

mr = mr0;

p = s0vect;

for ( j = k = 0; j < col; j++, p += nd_wpd )

if ( !rhead[j] ) {

if ( (c = (UINT)vect[k++]) != 0 ) {

ndl_copy(p,DL(mr)); CM(mr) = c; NMV_ADV(mr);

}

MKNDV(nd_nvar,mr0,len,r);

return r;

}

int nd_to_vect64(int mod,UINT *s0,int n,ND d,mp_limb_t *r)

{

NM m;

UINT *t,*s,*u;

int i,st,ed,md,prev,c;

for ( i = 0; i < n; i++ ) r[i] = 0;

prev = 0;

for ( i = 0, m = BDY(d); m; m = NEXT(m) ) {

t = DL(m);

st = prev;

ed = n;

while ( ed > st ) {

md = (st+ed)/2;

u = s0+md*nd_wpd;

c = DL_COMPARE(u,t);

if ( c == 0 ) break;

else if ( c > 0 ) st = md;

else ed = md;

}

r[md] = (mp_limb_t)CM(m);

prev = md;

}

for ( i = 0; !r[i]; i++ );

return i;

}

#define MOD128(a,c,m) ((a)=(((c)!=0||((a)>=(m)))?(((((U128)(c))<<64)+(a))%(m)):(a)))

int ndv_reduce_vect64(int m,mp_limb_t *svect,mp_limb_t *cvect,int col,IndArray *imat,NM_ind_pair *rp0,int nred)

{

int i,j,k,len,pos,prev;

mp_limb_t a,c,c1,c2;

IndArray ivect;

unsigned char *ivc;

unsigned short *ivs;

unsigned int *ivi;

NDV redv;

NMV mr;

NODE rp;

int maxrs;

for ( i = 0; i < col; i++ ) cvect[i] = 0;

maxrs = 0;

for ( i = 0; i < nred; i++ ) {

ivect = imat[i];

k = ivect->head;

a = svect[k]; c = cvect[k];

MOD128(a,c,m);

svect[k] = a; cvect[k] = 0;

if ( (c = svect[k]) != 0 ) {

Nf4_red++;

maxrs = MAX(maxrs,rp0[i]->sugar);

c = m-c; redv = nd_ps[rp0[i]->index];

len = LEN(redv); mr = BDY(redv);

svect[k] = 0; prev = k;

switch ( ivect->width ) {

case 1:

ivc = ivect->index.c;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivc[j]; c1 = CM(mr); prev = pos;

c2 = svect[pos]+c1*c;

if ( c2 < svect[pos] ) cvect[pos]++;

svect[pos] = c2;

}

break;

case 2:

ivs = ivect->index.s;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivs[j]; c1 = CM(mr); prev = pos;

c2 = svect[pos]+c1*c;

if ( c2 < svect[pos] ) cvect[pos]++;

svect[pos] = c2;

}

break;

case 4:

ivi = ivect->index.i;

for ( j = 1, NMV_ADV(mr); j < len; j++, NMV_ADV(mr) ) {

pos = prev+ivi[j]; c1 = CM(mr); prev = pos;

c2 = svect[pos]+c1*c;

if ( c2 < svect[pos] ) cvect[pos]++;

svect[pos] = c2;

}

break;

}

for ( i = 0; i < col; i++ ) {

a = svect[i]; c = cvect[i]; MOD128(a,c,m); svect[i] = a;

}

return maxrs;

}

/* for Fp, 2^15=<p<2^29 */

NODE nd_f4_red_mod64_main(int m,ND_pairs sp0,int nsp,UINT *s0vect,int col,

NM_ind_pair *rvect,int *rhead,IndArray *imat,int nred,ND_pairs *nz)

{

int spcol,sprow,a;

int i,j,k,l,rank;

NODE r0,r;

ND_pairs sp;

ND spol;

mp_limb_t **spmat;

mp_limb_t *svect,*cvect;

mp_limb_t *v;

int *colstat;

struct oEGT eg0,eg1,eg2,eg_f4,eg_f4_1,eg_f4_2;

int maxrs;

int *spsugar;

ND_pairs *spactive;

spcol = col-nred;

get_eg(&eg0);

/* elimination (1st step) */

spmat = (mp_limb_t **)MALLOC(nsp*sizeof(mp_limb_t *));

svect = (mp_limb_t *)MALLOC(col*sizeof(mp_limb_t));

cvect = (mp_limb_t *)MALLOC(col*sizeof(mp_limb_t));

spsugar = (int *)MALLOC(nsp*sizeof(int));

spactive = !nz?0:(ND_pairs *)MALLOC(nsp*sizeof(ND_pairs));

for ( a = sprow = 0, sp = sp0; a < nsp; a++, sp = NEXT(sp) ) {

nd_sp(m,0,sp,&spol);

if ( !spol ) continue;

nd_to_vect64(m,s0vect,col,spol,svect);

maxrs = ndv_reduce_vect64(m,svect,cvect,col,imat,rvect,nred);

for ( i = 0; i < col; i++ ) if ( svect[i] ) break;

if ( i < col ) {

spmat[sprow] = v = (mp_limb_t *)MALLOC_ATOMIC(spcol*sizeof(mp_limb_t));

for ( j = k = 0; j < col; j++ )

if ( !rhead[j] ) v[k++] = (UINT)svect[j];

spsugar[sprow] = MAX(maxrs,SG(spol));

if ( nz )

spactive[sprow] = sp;

sprow++;

}

nd_free(spol);

}

get_eg(&eg1); init_eg(&eg_f4_1); add_eg(&eg_f4_1,&eg0,&eg1); add_eg(&f4_elim1,&eg0,&eg1);

if ( DP_Print ) {

fprintf(asir_out,"elim1=%.3fsec,",eg_f4_1.exectime);

fflush(asir_out);

}

/* free index arrays */

for ( i = 0; i < nred; i++ ) GCFREE(imat[i]->index.c);

/* elimination (2nd step) */

colstat = (int *)MALLOC(spcol*sizeof(int));

rank = nd_gauss_elim_mod64(spmat,spsugar,spactive,sprow,spcol,m,colstat);

r0 = 0;

for ( i = 0; i < rank; i++ ) {

NEXTNODE(r0,r); BDY(r) =

(pointer)vect64_to_ndv(spmat[i],spcol,col,rhead,s0vect);

SG((NDV)BDY(r)) = spsugar[i];

GCFREE(spmat[i]);

}

if ( r0 ) NEXT(r) = 0;

for ( ; i < sprow; i++ ) GCFREE(spmat[i]);

get_eg(&eg2); init_eg(&eg_f4_2); add_eg(&eg_f4_2,&eg1,&eg2); add_eg(&f4_elim2,&eg1,&eg2);

init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg2);

if ( DP_Print ) {

fprintf(asir_out,"elim2=%.3fsec,",eg_f4_2.exectime);

fprintf(asir_out,"nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ",

nsp,nred,sprow,spcol,rank);

fprintf(asir_out,"%.3fsec,",eg_f4.exectime);

}

if ( nz ) {

for ( i = 0; i < rank-1; i++ ) NEXT(spactive[i]) = spactive[i+1];

if ( rank > 0 ) {

NEXT(spactive[rank-1]) = 0;

*nz = spactive[0];

} else

*nz = 0;

}

return r0;

}

int nd_gauss_elim_mod64(mp_limb_t **mat,int *sugar,ND_pairs *spactive,int row,int col,int md,int *colstat)

{

int i,j,k,l,rank,s;

mp_limb_t inv;

mp_limb_t a;

UINT c;

mp_limb_t *t,*pivot,*pk;

UINT *ck;

UINT **cmat;

UINT *ct;

ND_pairs pair;

cmat = (UINT **)MALLOC(row*sizeof(UINT *));

for ( i = 0; i < row; i++ ) {

cmat[i] = MALLOC_ATOMIC(col*sizeof(UINT));

bzero(cmat[i],col*sizeof(UINT));

}

for ( rank = 0, j = 0; j < col; j++ ) {

for ( i = rank; i < row; i++ ) {

a = mat[i][j]; c = cmat[i][j];

MOD128(a,c,md);

mat[i][j] = a; cmat[i][j] = 0;

}

for ( i = rank; i < row; i++ )

if ( mat[i][j] )

break;

if ( i == row ) {

colstat[j] = 0;

continue;

} else

colstat[j] = 1;

if ( i != rank ) {

t = mat[i]; mat[i] = mat[rank]; mat[rank] = t;

ct = cmat[i]; cmat[i] = cmat[rank]; cmat[rank] = ct;

s = sugar[i]; sugar[i] = sugar[rank]; sugar[rank] = s;

if ( spactive ) {

pair = spactive[i]; spactive[i] = spactive[rank];

spactive[rank] = pair;

}

/* column j is normalized */

s = sugar[rank];

inv = invm((UINT)mat[rank][j],md);

/* normalize pivot row */

for ( k = j, pk = mat[rank]+j, ck = cmat[rank]+j; k < col; k++, pk++, ck++ ) {

a = *pk; c = *ck; MOD128(a,c,md); *pk = (a*inv)%md; *ck = 0;

}

for ( i = rank+1; i < row; i++ ) {

if ( (a = mat[i][j]) != 0 ) {

sugar[i] = MAX(sugar[i],s);

red_by_vect64(md,mat[i]+j,cmat[i]+j,mat[rank]+j,(int)(md-a),col-j);

Nf4_red++;

}

rank++;

}

for ( j = col-1, l = rank-1; j >= 0; j-- )

if ( colstat[j] ) {

for ( k = j, pk = mat[l]+j, ck = cmat[l]+j; k < col; k++, pk++, ck++ ) {

a = *pk; c = *ck; MOD128(a,c,md); *pk = a; *ck = 0;

}

s = sugar[l];

for ( i = 0; i < l; i++ ) {

a = mat[i][j]; c = cmat[i][j]; MOD128(a,c,md); mat[i][j] = a; cmat[i][j] = 0;

if ( a ) {

sugar[i] = MAX(sugar[i],s);

red_by_vect64(md,mat[i]+j,cmat[i]+j,mat[l]+j,(int)(md-a),col-j);

Nf4_red++;

}

l--;

}

for ( i = 0; i < row; i++ ) GCFREE(cmat[i]);

GCFREE(cmat);

return rank;

}

#endif

FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>