=================================================================== RCS file: /home/cvs/OpenXM_contrib2/asir2018/engine/nd.c,v retrieving revision 1.28 retrieving revision 1.30 diff -u -p -r1.28 -r1.30 --- OpenXM_contrib2/asir2018/engine/nd.c 2020/06/30 01:52:17 1.28 +++ OpenXM_contrib2/asir2018/engine/nd.c 2020/07/03 03:37:59 1.30 @@ -1,4 +1,4 @@ -/* $OpenXM: OpenXM_contrib2/asir2018/engine/nd.c,v 1.27 2020/06/25 02:53:31 noro Exp $ */ +/* $OpenXM: OpenXM_contrib2/asir2018/engine/nd.c,v 1.29 2020/07/02 09:24:17 noro Exp $ */ #include "nd.h" @@ -50,7 +50,7 @@ static NDV *nd_ps_trace; static NDV *nd_ps_sym; static NDV *nd_ps_trace_sym; static RHist *nd_psh; -static int nd_psn,nd_pslen; +static int nd_psn,nd_pslen,nd_nbase; static RHist *nd_red; static int *nd_work_vector; static int **nd_matrix; @@ -72,6 +72,7 @@ static int nd_splist; static int *nd_sugarweight; static int nd_f4red,nd_rank0,nd_last_nonzero; static DL *nd_sba_hm; +static NODE *nd_sba_pos; NumberField get_numberfield(); UINT *nd_det_compute_bound(NDV **dm,int n,int j); @@ -1250,7 +1251,7 @@ INLINE int ndl_find_reducer_s(UINT *dg,SIG sig) quo->pos = nd_psh[i]->sig->pos; ret = comp_sig(sig,quo); if ( ret > 0 ) { singular = 0; break; } - if ( ret == 0 ) { singular = 1; } + if ( ret == 0 ) { fprintf(asir_out,"s"); fflush(asir_out); singular = 1; } } } if ( singular ) return -1; @@ -2483,8 +2484,8 @@ print_eg("update",&eg_update); return g; } -ND_pairs update_pairs_s(ND_pairs d,NODE g,int t,NODE syz); -ND_pairs nd_newpairs_s( NODE g, int t ,NODE syz); +ND_pairs update_pairs_s(ND_pairs d,int t,NODE *syz); +ND_pairs nd_newpairs_s(int t ,NODE *syz); int nd_nf_pbucket_s(int mod,ND g,NDV *ps,int full,ND *nf); int nd_nf_s(int mod,ND d,ND g,NDV *ps,int full,ND *nf); @@ -2543,10 +2544,13 @@ ND_pairs remove_spair_s(ND_pairs d,SIG sig) return (ND_pairs)root.next; } +int _dl_redble_ext(DL,DL,DL,int); + int small_lcm(ND_pairs l) { SIG sig; int i; + NODE t; static DL lcm,mul,quo; static int nvar; @@ -2555,11 +2559,10 @@ int small_lcm(ND_pairs l) } sig = l->sig; _ndltodl(l->lcm,lcm); +#if 0 for ( i = 0; i < nd_psn; i++ ) { if ( sig->pos == nd_psh[i]->sig->pos && - _dl_redble(DL(nd_psh[i]->sig),DL(sig),nd_nvar) ) { - _copydl(nd_nvar,DL(sig),quo); - _subfromdl(nd_nvar,DL(nd_psh[i]->sig),quo); + _dl_redble_ext(DL(nd_psh[i]->sig),DL(sig),quo,nd_nvar) ) { _ndltodl(DL(nd_psh[i]),mul); _addtodl(nd_nvar,quo,mul); if ( (*cmpdl)(nd_nvar,lcm,mul) > 0 ) @@ -2568,6 +2571,19 @@ int small_lcm(ND_pairs l) } if ( i < nd_psn ) return 1; else return 0; +#else + for ( t = nd_sba_pos[sig->pos]; t; t = t->next ) { + i = (long)BDY(t); + if ( _dl_redble_ext(DL(nd_psh[i]->sig),DL(sig),quo,nd_nvar) ) { + _ndltodl(DL(nd_psh[i]),mul); + _addtodl(nd_nvar,quo,mul); + if ( (*cmpdl)(nd_nvar,lcm,mul) > 0 ) + break; + } + } + if ( t ) return 1; + else return 0; +#endif } ND_pairs remove_large_lcm(ND_pairs d) @@ -2590,10 +2606,12 @@ ND_pairs remove_large_lcm(ND_pairs d) struct oEGT eg_create,eg_newpairs,eg_merge; +NODE conv_ilist_s(int demand,int trace,int **indp); + NODE nd_sba_buch(int m,int ishomo,int **indp) { - int i,nh,sugar,stat; - NODE r,g,t; + int i,j,nh,sugar,stat; + NODE r,t,g; ND_pairs d; ND_pairs l; ND h,nf,s,head,nf1; @@ -2603,25 +2621,31 @@ NODE nd_sba_buch(int m,int ishomo,int **indp) P cont; LIST list; SIG sig; - NODE syzlist; + NODE *syzlist; int Nredundant; DL lcm,quo,mul; - struct oEGT eg1,eg2,eg_update,eg_remove; + struct oEGT eg1,eg2,eg_update,eg_remove,eg_large; init_eg(&eg_remove); - syzlist = 0; + syzlist = (NODE *)MALLOC(nd_psn*sizeof(NODE)); Nsyz = 0; Nnd_add = 0; Nredundant = 0; - g = 0; d = 0; + d = 0; for ( i = 0; i < nd_psn; i++ ) { - d = update_pairs_s(d,g,i,0); - g = append_one(g,i); + d = update_pairs_s(d,i,syzlist); } + for ( i = 0; i < nd_psn; i++ ) + for ( j = i+1; j < nd_psn; j++ ) { + NEWSIG(sig); sig->pos = j; + _copydl(nd_nvar,nd_sba_hm[i],sig->dl); + syzlist[sig->pos] = insert_sig(syzlist[sig->pos],sig); + } sugar = 0; NEWDL(lcm,nd_nvar); NEWDL(quo,nd_nvar); NEWDL(mul,nd_nvar); init_eg(&eg_create); init_eg(&eg_merge); +init_eg(&eg_large); while ( d ) { again: if ( DP_Print ) { @@ -2667,20 +2691,20 @@ again: nfv = ndtondv(m,nf); nd_free(nf); nh = ndv_newps(m,nfv,0); - d = update_pairs_s(d,g,nh,syzlist); - g = append_one(g,nh); + d = update_pairs_s(d,nh,syzlist); + nd_sba_pos[sig->pos] = append_one(nd_sba_pos[sig->pos],nh); FREENDP(l); } else { // syzygy get_eg(&eg1); d = remove_spair_s(d,sig); get_eg(&eg2); add_eg(&eg_remove,&eg1,&eg2); - syzlist = insert_sig(syzlist,sig); + syzlist[sig->pos] = insert_sig(syzlist[sig->pos],sig); if ( DP_Print ) { printf("."); fflush(stdout); } FREENDP(l); } } - conv_ilist(nd_demand,0,g,indp); + g = conv_ilist_s(nd_demand,0,indp); if ( DP_Print ) { printf("\nnd_sba done. nd_add=%d,Nsyz=%d,Nredundant=%d\n",Nnd_add,Nsyz,Nredundant); fflush(stdout); @@ -3111,14 +3135,14 @@ ND_pairs update_pairs( ND_pairs d, NODE /* of index */ ND_pairs merge_pairs_s(ND_pairs d,ND_pairs d1); -ND_pairs update_pairs_s( ND_pairs d, NODE /* of index */ g, int t,NODE syz) +ND_pairs update_pairs_s( ND_pairs d, int t,NODE *syz) { ND_pairs d1; struct oEGT eg1,eg2,eg3; - if ( !g ) return d; + if ( !t ) return d; get_eg(&eg1); - d1 = nd_newpairs_s(g,t,syz); + d1 = nd_newpairs_s(t,syz); get_eg(&eg2); add_eg(&eg_create,&eg1,&eg2); d = merge_pairs_s(d,d1); get_eg(&eg3); add_eg(&eg_merge,&eg2,&eg3); @@ -3162,7 +3186,6 @@ ND_pairs nd_newpairs( NODE g, int t ) return r0; } - int comp_sig(SIG s1,SIG s2) { #if 0 @@ -3323,27 +3346,44 @@ ND_pairs insert_pair_s(ND_pairs l,ND_pairs s) } } -ND_pairs nd_newpairs_s( NODE g, int t, NODE syz) +INLINE int __dl_redble(DL d1,DL d2,int nvar) { + int i; + + if ( d1->td > d2->td ) + return 0; + for ( i = nvar-1; i >= 0; i-- ) + if ( d1->d[i] > d2->d[i] ) + break; + if ( i >= 0 ) + return 0; + else + return 1; +} + +ND_pairs nd_newpairs_s(int t, NODE *syz) +{ NODE h,s; UINT *dl; int ts,ret,i; ND_pairs r,r0,_sp,sp; - SIG _sig1,_sig2,spsig,tsig; + SIG spsig,tsig; + static int nvar; + static SIG _sig1,_sig2; struct oEGT eg1,eg2,eg3,eg4; NEWND_pairs(_sp); - NEWSIG(_sig1); NEWSIG(_sig2); + if ( !_sig1 || nvar != nd_nvar ) { + nvar = nd_nvar; NEWSIG(_sig1); NEWSIG(_sig2); + } r0 = 0; for ( i = 0; i < t; i++ ) { ret = _create_spair_s(i,t,_sp,_sig1,_sig2); -// for ( h = g; h; h = NEXT(h) ) { -// ret = _create_spair_s((long)BDY(h),t,_sp,_sig1,_sig2); if ( ret ) { spsig = _sp->sig; - for ( s = syz; s; s = s->next ) { + for ( s = syz[spsig->pos]; s; s = s->next ) { tsig = (SIG)s->body; - if ( tsig->pos == spsig->pos && _dl_redble(DL(tsig),DL(spsig),nd_nvar) ) + if ( _dl_redble(DL(tsig),DL(spsig),nd_nvar) ) break; } if ( s == 0 ) { @@ -3853,6 +3893,11 @@ int ndv_setup(int mod,int trace,NODE f,int dont_sort,i NEWDL(nd_sba_hm[i],nd_nvar); _ndltodl(DL(nd_psh[i]),nd_sba_hm[i]); } + nd_sba_pos = (NODE *)MALLOC(nd_psn*sizeof(NODE)); + for ( i = 0; i < nd_psn; i++ ) { + j = nd_psh[i]->sig->pos; + nd_sba_pos[j] = append_one(nd_sba_pos[j],i); + } } if ( nd_gentrace && nd_tracelist ) NEXT(tn) = 0; return 1; @@ -4210,6 +4255,7 @@ void nd_sba(LIST f,LIST v,int m,int homo,int retdp,int int *perm; EPOS oepos; int obpe,oadv,ompos,cbpe; + struct oEGT eg0,eg1,egconv; nd_module = 0; nd_demand = 0; @@ -4281,6 +4327,7 @@ void nd_sba(LIST f,LIST v,int m,int homo,int retdp,int x = ndv_reducebase(x,perm); x = ndv_reduceall(m,x); nd_setup_parameters(nd_nvar,0); + get_eg(&eg0); for ( r0 = 0, t = x; t; t = NEXT(t) ) { NEXTNODE(r0,r); if ( retdp ) BDY(r) = ndvtodp(m,BDY(t)); @@ -4288,6 +4335,8 @@ void nd_sba(LIST f,LIST v,int m,int homo,int retdp,int } if ( r0 ) NEXT(r) = 0; MKLIST(*rp,r0); + get_eg(&eg1); init_eg(&egconv); add_eg(&egconv,&eg0,&eg1); + print_eg("conv",&egconv); fprintf(asir_out,"\n"); } void nd_gr_postproc(LIST f,LIST v,int m,struct order_spec *ord,int do_check,LIST *rp) @@ -5530,6 +5579,91 @@ ND_pairs nd_reconstruct(int trace,ND_pairs d) return s0; } +void nd_reconstruct_s(int trace,ND_pairs *d) +{ + int i,obpe,oadv,h; + static NM prev_nm_free_list; + static ND_pairs prev_ndp_free_list; + RHist mr0,mr; + RHist r; + RHist *old_red; + ND_pairs s0,s,t; + EPOS oepos; + + obpe = nd_bpe; + oadv = nmv_adv; + oepos = nd_epos; + if ( obpe < 2 ) nd_bpe = 2; + else if ( obpe < 3 ) nd_bpe = 3; + else if ( obpe < 4 ) nd_bpe = 4; + else if ( obpe < 5 ) nd_bpe = 5; + else if ( obpe < 6 ) nd_bpe = 6; + else if ( obpe < 8 ) nd_bpe = 8; + else if ( obpe < 10 ) nd_bpe = 10; + else if ( obpe < 16 ) nd_bpe = 16; + else if ( obpe < 32 ) nd_bpe = 32; + else error("nd_reconstruct_s : exponent too large"); + + nd_setup_parameters(nd_nvar,0); + prev_nm_free_list = _nm_free_list; + prev_ndp_free_list = _ndp_free_list; + _nm_free_list = 0; + _ndp_free_list = 0; + for ( i = nd_psn-1; i >= 0; i-- ) { + ndv_realloc(nd_ps[i],obpe,oadv,oepos); + ndv_realloc(nd_ps_sym[i],obpe,oadv,oepos); + } + if ( trace ) + for ( i = nd_psn-1; i >= 0; i-- ) { + ndv_realloc(nd_ps_trace[i],obpe,oadv,oepos); + ndv_realloc(nd_ps_trace_sym[i],obpe,oadv,oepos); + } + + for ( i = 0; i < nd_nbase; i++ ) { + s0 = 0; + for ( t = d[i]; t; t = NEXT(t) ) { + NEXTND_pairs(s0,s); + s->i1 = t->i1; + s->i2 = t->i2; + s->sig = t->sig; + SG(s) = SG(t); + ndl_reconstruct(LCM(t),LCM(s),obpe,oepos); + } + d[i] = s0; + } + + old_red = (RHist *)MALLOC(REDTAB_LEN*sizeof(RHist)); + for ( i = 0; i < REDTAB_LEN; i++ ) { + old_red[i] = nd_red[i]; + nd_red[i] = 0; + } + for ( i = 0; i < REDTAB_LEN; i++ ) + for ( r = old_red[i]; r; r = NEXT(r) ) { + NEWRHist(mr); + mr->index = r->index; + SG(mr) = SG(r); + ndl_reconstruct(DL(r),DL(mr),obpe,oepos); + h = ndl_hash_value(DL(mr)); + NEXT(mr) = nd_red[h]; + nd_red[h] = mr; + mr->sig = r->sig; + } + for ( i = 0; i < REDTAB_LEN; i++ ) old_red[i] = 0; + old_red = 0; + for ( i = 0; i < nd_psn; i++ ) { + NEWRHist(r); SG(r) = SG(nd_psh[i]); + ndl_reconstruct(DL(nd_psh[i]),DL(r),obpe,oepos); + r->sig = nd_psh[i]->sig; + nd_psh[i] = r; + } + if ( s0 ) NEXT(s) = 0; + prev_nm_free_list = 0; + prev_ndp_free_list = 0; +#if 0 + GC_gcollect(); +#endif +} + void ndl_reconstruct(UINT *d,UINT *r,int obpe,EPOS oepos) { int n,i,ei,oepw,omask0,j,s,ord_l,l; @@ -9400,6 +9534,25 @@ void conv_ilist(int demand,int trace,NODE g,int **indp if ( indp ) *indp = ind; } +NODE conv_ilist_s(int demand,int trace,int **indp) +{ + int n,i,j; + int *ind; + NODE g0,g; + + n = nd_psn; + ind = (int *)MALLOC(n*sizeof(int)); + g0 = 0; + for ( i = 0; i < n; i++ ) { + ind[i] = i; + NEXTNODE(g0,g); + BDY(g) = (pointer)(demand?ndv_load(i):(trace?nd_ps_trace[i]:nd_ps[i])); + } + if ( g0 ) NEXT(g) = 0; + if ( indp ) *indp = ind; + return g0; +} + void parse_nd_option(NODE opt) { NODE t,p,u; @@ -10540,7 +10693,7 @@ NODE nd_f4_red_mod64_main_s(int m,ND_pairs sp0,int nsp for ( a = sprow = 0, sp = sp0; a < nsp; a++, sp = NEXT(sp) ) { nd_sp(m,0,sp,&spol); if ( !spol ) { - *syzlistp = insert_sig(*syzlistp,sp->sig); + syzlistp[sp->sig->pos] = insert_sig(syzlistp[sp->sig->pos],sp->sig); continue; } svect = (mp_limb_t *)MALLOC(col*sizeof(mp_limb_t)); @@ -10553,7 +10706,7 @@ NODE nd_f4_red_mod64_main_s(int m,ND_pairs sp0,int nsp spsig[sprow] = sp->sig; sprow++; } else { - *syzlistp = insert_sig(*syzlistp,sp->sig); + syzlistp[sp->sig->pos] = insert_sig(syzlistp[sp->sig->pos],sp->sig); } nd_free(spol); } @@ -10580,7 +10733,7 @@ NODE nd_f4_red_mod64_main_s(int m,ND_pairs sp0,int nsp for ( ; i < sprow; i++ ) { GCFREE(spmat[i]); - *syzlistp = insert_sig(*syzlistp,spsig[i]); + syzlistp[spsig[i]->pos] = insert_sig(syzlistp[spsig[i]->pos],spsig[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); @@ -10754,18 +10907,19 @@ NODE nd_sba_f4(int m,int **indp) int *rhead; int spcol,sprow; int sugar,sugarh; + int nbase; PGeoBucket bucket; struct oEGT eg0,eg1,eg_f4; Z i1,i2,sugarq; - NODE syzlist; + NODE *syzlist; Nf4_red=0; - g = 0; d = 0; - syzlist = 0; + d = 0; + syzlist = (NODE *)MALLOC(nd_psn*sizeof(NODE)); for ( i = 0; i < nd_psn; i++ ) { - d = update_pairs_s(d,g,i,syzlist); - g = append_one(g,i); + d = update_pairs_s(d,i,syzlist); } + nd_nbase = nd_psn; f4red = 1; while ( d ) { l = nd_minsugarp_s(d,&d); @@ -10790,7 +10944,7 @@ NODE nd_sba_f4(int m,int **indp) } if ( DP_Print ) fprintf(asir_out,"sugar=%d,",sugar); l = remove_large_lcm(l); - nflist = nd_f4_red_s(m,l,0,s0vect,col,rp0,&syzlist); + nflist = nd_f4_red_s(m,l,0,s0vect,col,rp0,syzlist); /* adding new bases */ for ( r = nflist; r; r = NEXT(r) ) { ND tmp,tmpred; @@ -10807,14 +10961,15 @@ NODE nd_sba_f4(int m,int **indp) nf = ndtondv(m,tmpred); ndv_removecont(m,nf); nh = ndv_newps(m,nf,0); - d = update_pairs_s(d,g,nh,syzlist); - g = append_one(g,nh); + d = update_pairs_s(d,nh,syzlist); + nd_sba_pos[sig->pos] = append_one(nd_sba_pos[sig->pos],nh); } else { - syzlist = insert_sig(syzlist,sig); + syzlist[sig->pos] = insert_sig(syzlist[sig->pos],sig); } } - for ( r = syzlist; r; r = NEXT(r) ) - d = remove_spair_s(d,(SIG)BDY(r)); + for ( i = 0; i < nd_nbase; i++ ) + for ( r = syzlist[i]; r; r = NEXT(r) ) + d = remove_spair_s(d,(SIG)BDY(r)); if ( DP_Print ) { fprintf(asir_out,"f4red=%d,gblen=%d\n",f4red,length(g)); fflush(asir_out); } @@ -10825,6 +10980,6 @@ NODE nd_sba_f4(int m,int **indp) if ( DP_Print ) { fprintf(asir_out,"number of red=%d,",Nf4_red); } - conv_ilist(nd_demand,0,g,indp); + g = conv_ilist_s(nd_demand,0,indp); return g; }