version 1.63, 2003/09/11 09:03:53 |
version 1.66, 2003/09/12 08:26:19 |
|
|
/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.62 2003/09/11 01:52:25 noro Exp $ */ |
/* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.65 2003/09/12 08:01:51 noro Exp $ */ |
|
|
#include "ca.h" |
#include "ca.h" |
#include "inline.h" |
#include "inline.h" |
|
#include <time.h> |
|
|
#if defined(__GNUC__) |
#if defined(__GNUC__) |
#define INLINE inline |
#define INLINE inline |
|
|
typedef unsigned int UINT; |
typedef unsigned int UINT; |
|
|
#define USE_GEOBUCKET 1 |
#define USE_GEOBUCKET 1 |
|
#define USE_UNROLL 1 |
|
|
#define REDTAB_LEN 32003 |
#define REDTAB_LEN 32003 |
|
|
Line 347 P ndvtop(int mod,VL vl,VL dvl,NDV p); |
|
Line 349 P ndvtop(int mod,VL vl,VL dvl,NDV p); |
|
NDV ndtondv(int mod,ND p); |
NDV ndtondv(int mod,ND p); |
ND ndvtond(int mod,NDV p); |
ND ndvtond(int mod,NDV p); |
int nm_ind_pair_to_vect(int m,UINT *s0,int n,NM_ind_pair pair,UINT *r); |
int nm_ind_pair_to_vect(int m,UINT *s0,int n,NM_ind_pair pair,UINT *r); |
|
void nm_ind_pair_to_vect_compress(int m,UINT *s0,int n,NM_ind_pair pair,int *ind); |
int nd_to_vect(int mod,UINT *s0,int n,ND d,UINT *r); |
int nd_to_vect(int mod,UINT *s0,int n,ND d,UINT *r); |
|
|
void nd_free_private_storage() |
void nd_free_private_storage() |
Line 411 INLINE int ndl_reducible(UINT *d1,UINT *d2) |
|
Line 414 INLINE int ndl_reducible(UINT *d1,UINT *d2) |
|
int i,j; |
int i,j; |
|
|
if ( TD(d1) < TD(d2) ) return 0; |
if ( TD(d1) < TD(d2) ) return 0; |
|
#if USE_UNROLL |
switch ( nd_bpe ) { |
switch ( nd_bpe ) { |
#if 0 |
|
case 3: |
case 3: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u1 = d1[i]; u2 = d2[i]; |
u1 = d1[i]; u2 = d2[i]; |
Line 477 INLINE int ndl_reducible(UINT *d1,UINT *d2) |
|
Line 480 INLINE int ndl_reducible(UINT *d1,UINT *d2) |
|
if ( d1[i] < d2[i] ) return 0; |
if ( d1[i] < d2[i] ) return 0; |
return 1; |
return 1; |
break; |
break; |
#endif |
|
default: |
default: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u1 = d1[i]; u2 = d2[i]; |
u1 = d1[i]; u2 = d2[i]; |
Line 486 INLINE int ndl_reducible(UINT *d1,UINT *d2) |
|
Line 488 INLINE int ndl_reducible(UINT *d1,UINT *d2) |
|
} |
} |
return 1; |
return 1; |
} |
} |
|
#else |
|
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
|
u1 = d1[i]; u2 = d2[i]; |
|
for ( j = 0; j < nd_epw; j++ ) |
|
if ( (u1&nd_mask[j]) < (u2&nd_mask[j]) ) return 0; |
|
} |
|
return 1; |
|
#endif |
} |
} |
|
|
/* |
/* |
Line 553 void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) |
|
Line 563 void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) |
|
UINT t1,t2,u,u1,u2; |
UINT t1,t2,u,u1,u2; |
int i,j,l; |
int i,j,l; |
|
|
|
#if USE_UNROLL |
switch ( nd_bpe ) { |
switch ( nd_bpe ) { |
#if 0 |
|
case 3: |
case 3: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u1 = d1[i]; u2 = d2[i]; |
u1 = d1[i]; u2 = d2[i]; |
Line 620 void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) |
|
Line 630 void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) |
|
d[i] = u1>u2?u1:u2; |
d[i] = u1>u2?u1:u2; |
} |
} |
break; |
break; |
#endif |
|
default: |
default: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u1 = d1[i]; u2 = d2[i]; |
u1 = d1[i]; u2 = d2[i]; |
Line 631 void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) |
|
Line 640 void ndl_lcm(UINT *d1,unsigned *d2,UINT *d) |
|
} |
} |
break; |
break; |
} |
} |
|
#else |
|
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
|
u1 = d1[i]; u2 = d2[i]; |
|
for ( j = 0, u = 0; j < nd_epw; j++ ) { |
|
t1 = (u1&nd_mask[j]); t2 = (u2&nd_mask[j]); u |= t1>t2?t1:t2; |
|
} |
|
d[i] = u; |
|
} |
|
#endif |
TD(d) = ndl_weight(d); |
TD(d) = ndl_weight(d); |
if ( nd_blockmask ) ndl_weight_mask(d); |
if ( nd_blockmask ) ndl_weight_mask(d); |
} |
} |
Line 836 int ndl_disjoint(UINT *d1,UINT *d2) |
|
Line 854 int ndl_disjoint(UINT *d1,UINT *d2) |
|
UINT t1,t2,u,u1,u2; |
UINT t1,t2,u,u1,u2; |
int i,j; |
int i,j; |
|
|
|
#if USE_UNROLL |
switch ( nd_bpe ) { |
switch ( nd_bpe ) { |
#if 0 |
|
case 3: |
case 3: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u1 = d1[i]; u2 = d2[i]; |
u1 = d1[i]; u2 = d2[i]; |
Line 902 int ndl_disjoint(UINT *d1,UINT *d2) |
|
Line 920 int ndl_disjoint(UINT *d1,UINT *d2) |
|
if ( d1[i] && d2[i] ) return 0; |
if ( d1[i] && d2[i] ) return 0; |
return 1; |
return 1; |
break; |
break; |
#endif |
|
default: |
default: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u1 = d1[i]; u2 = d2[i]; |
u1 = d1[i]; u2 = d2[i]; |
Line 914 int ndl_disjoint(UINT *d1,UINT *d2) |
|
Line 931 int ndl_disjoint(UINT *d1,UINT *d2) |
|
return 1; |
return 1; |
break; |
break; |
} |
} |
|
#else |
|
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
|
u1 = d1[i]; u2 = d2[i]; |
|
for ( j = 0; j < nd_epw; j++ ) { |
|
if ( (u1&nd_mask0) && (u2&nd_mask0) ) return 0; |
|
u1 >>= nd_bpe; u2 >>= nd_bpe; |
|
} |
|
} |
|
return 1; |
|
#endif |
} |
} |
|
|
int ndl_check_bound2(int index,UINT *d2) |
int ndl_check_bound2(int index,UINT *d2) |
Line 924 int ndl_check_bound2(int index,UINT *d2) |
|
Line 951 int ndl_check_bound2(int index,UINT *d2) |
|
|
|
d1 = nd_bound[index]; |
d1 = nd_bound[index]; |
ind = 0; |
ind = 0; |
|
#if USE_UNROLL |
switch ( nd_bpe ) { |
switch ( nd_bpe ) { |
#if 0 |
|
case 3: |
case 3: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u2 = d2[i]; |
u2 = d2[i]; |
Line 990 int ndl_check_bound2(int index,UINT *d2) |
|
Line 1017 int ndl_check_bound2(int index,UINT *d2) |
|
if ( d1[i]+d2[i]<d1[i] ) return 1; |
if ( d1[i]+d2[i]<d1[i] ) return 1; |
return 0; |
return 0; |
break; |
break; |
#endif |
|
default: |
default: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u2 = d2[i]; |
u2 = d2[i]; |
Line 1001 int ndl_check_bound2(int index,UINT *d2) |
|
Line 1027 int ndl_check_bound2(int index,UINT *d2) |
|
return 0; |
return 0; |
break; |
break; |
} |
} |
|
#else |
|
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
|
u2 = d2[i]; |
|
k = (nd_epw-1)*nd_bpe; |
|
for ( j = 0; j < nd_epw; j++, k -= nd_bpe ) |
|
if ( d1[ind++]+((u2>>k)&nd_mask0) > nd_mask0 ) return 1; |
|
} |
|
return 0; |
|
#endif |
} |
} |
|
|
int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
Line 1009 int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
|
Line 1044 int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
|
int i,j,ind,k; |
int i,j,ind,k; |
|
|
ind = 0; |
ind = 0; |
|
#if USE_UNROLL |
switch ( nd_bpe ) { |
switch ( nd_bpe ) { |
#if 0 |
|
case 3: |
case 3: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u2 = d2[i]; |
u2 = d2[i]; |
Line 1075 int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
|
Line 1110 int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
|
if ( d1[i]+d2[i]<d1[i] ) return 1; |
if ( d1[i]+d2[i]<d1[i] ) return 1; |
return 0; |
return 0; |
break; |
break; |
#endif |
|
default: |
default: |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
u2 = d2[i]; |
u2 = d2[i]; |
Line 1086 int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
|
Line 1120 int ndl_check_bound2_direct(UINT *d1,UINT *d2) |
|
return 0; |
return 0; |
break; |
break; |
} |
} |
|
#else |
|
for ( i = nd_exporigin; i < nd_wpd; i++ ) { |
|
u2 = d2[i]; |
|
k = (nd_epw-1)*nd_bpe; |
|
for ( j = 0; j < nd_epw; j++, k -= nd_bpe ) |
|
if ( d1[ind++]+((u2>>k)&nd_mask0) > nd_mask0 ) return 1; |
|
} |
|
return 0; |
|
#endif |
} |
} |
|
|
INLINE int ndl_hash_value(UINT *d) |
INLINE int ndl_hash_value(UINT *d) |
Line 3744 int nm_ind_pair_to_vect(int mod,UINT *s0,int n,NM_ind_ |
|
Line 3787 int nm_ind_pair_to_vect(int mod,UINT *s0,int n,NM_ind_ |
|
return i; |
return i; |
} |
} |
|
|
|
void nm_ind_pair_to_vect_compress(int mod,UINT *s0,int n, |
|
NM_ind_pair pair,int *ind) |
|
{ |
|
NM m; |
|
NMV mr; |
|
UINT *d,*t,*s; |
|
NDV p; |
|
int i,j,len; |
|
|
|
m = pair->mul; |
|
d = DL(m); |
|
p = nd_ps[pair->index]; |
|
len = LEN(p); |
|
t = (UINT *)ALLOCA(nd_wpd*sizeof(UINT)); |
|
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++ ); |
|
ind[j] = i; |
|
} |
|
} |
|
|
|
|
|
void ndv_reduce_vect(int m,UINT *svect,int col,int **imat,NODE rp0) |
|
{ |
|
int i,j,k,len,pos; |
|
UINT c,c1,c2,c3,up,lo,dmy; |
|
UINT *ivect; |
|
NDV redv; |
|
NMV mr0,mr; |
|
NODE rp; |
|
|
|
for ( rp = rp0, i = 0; rp; i++, rp = NEXT(rp) ) { |
|
ivect = imat[i]; |
|
k = ivect[0]; |
|
svect[k] %= m; |
|
if ( c = svect[k] ) { |
|
c = m-c; |
|
redv = nd_ps[((NM_ind_pair)BDY(rp))->index]; |
|
len = LEN(redv); |
|
mr0 = BDY(redv); |
|
for ( j = 0, mr = mr0; j < len; j++, NMV_ADV(mr) ) { |
|
pos = ivect[j]; |
|
c1 = CM(mr); |
|
c2 = svect[pos]; |
|
DMA(c1,c,c2,up,lo); |
|
if ( up ) { |
|
DSAB(m,up,lo,dmy,c3); svect[pos] = c3; |
|
} else |
|
svect[pos] = lo; |
|
} |
|
} |
|
} |
|
for ( i = 0; i < col; i++ ) |
|
if ( svect[i] >= (UINT)m ) svect[i] %= m; |
|
} |
|
|
|
NDV vect_to_ndv(UINT *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(nmv_adv*len); |
|
mr = mr0; |
|
p = s0vect; |
|
for ( j = k = 0; j < col; j++, p += nd_wpd ) |
|
if ( !rhead[j] ) { |
|
if ( c = vect[k++] ) { |
|
ndl_copy(p,DL(mr)); CM(mr) = c; NMV_ADV(mr); |
|
} |
|
} |
|
MKNDV(nd_nvar,mr0,len,r); |
|
return r; |
|
} |
|
} |
|
|
|
NODE nd_sp_f4(int m,ND_pairs l,PGeoBucket bucket) |
|
{ |
|
ND_pairs t; |
|
NODE sp0,sp; |
|
int stat; |
|
ND spol; |
|
|
|
sp0 = 0; |
|
for ( t = l; t; t = NEXT(t) ) { |
|
stat = nd_sp(m,0,t,&spol); |
|
if ( !stat ) return 0; |
|
if ( spol ) { |
|
NEXTNODE(sp0,sp); BDY(sp) = (pointer)nd_dup(spol); |
|
add_pbucket_symbolic(bucket,spol); |
|
} |
|
} |
|
return sp0; |
|
} |
|
|
|
int nd_symbolic_preproc(PGeoBucket bucket,UINT **s0vect,NODE *r) |
|
{ |
|
NODE rp0,rp; |
|
NM mul,head,s0,s; |
|
int index,col,i; |
|
RHist h; |
|
UINT *s0v,*p; |
|
NM_ind_pair pair; |
|
ND red; |
|
|
|
s0 = 0; rp0 = 0; col = 0; |
|
while ( 1 ) { |
|
head = remove_head_pbucket_symbolic(bucket); |
|
if ( !head ) break; |
|
if ( !s0 ) s0 = head; |
|
else NEXT(s) = head; |
|
s = head; |
|
index = ndl_find_reducer(DL(head)); |
|
if ( index >= 0 ) { |
|
h = nd_psh[index]; |
|
NEWNM(mul); |
|
ndl_sub(DL(head),DL(h),DL(mul)); |
|
if ( ndl_check_bound2(index,DL(mul)) ) return 0; |
|
MKNM_ind_pair(pair,mul,index); |
|
red = ndv_mul_nm_symbolic(mul,nd_ps[index]); |
|
add_pbucket_symbolic(bucket,nd_remove_head(red)); |
|
NEXTNODE(rp0,rp); BDY(rp) = (pointer)pair; |
|
} |
|
col++; |
|
} |
|
NEXT(rp) = 0; NEXT(s) = 0; |
|
s0v = (UINT *)MALLOC_ATOMIC(col*nd_wpd*sizeof(UINT)); |
|
for ( i = 0, p = s0v, s = s0; i < col; |
|
i++, p += nd_wpd, s = NEXT(s) ) ndl_copy(DL(s),p); |
|
*s0vect = s0v; |
|
*r = rp0; |
|
return col; |
|
} |
|
|
NODE nd_f4(int m) |
NODE nd_f4(int m) |
{ |
{ |
int i,nh,stat,index; |
int i,nh,stat,index; |
NODE r,g; |
NODE r,g; |
ND_pairs d,l,t; |
ND_pairs d,l,t; |
ND spol,red; |
ND spol,red; |
NDV nf; |
NDV nf,redv; |
NM_ind_pair pair,pair1,pair2; |
|
NM s0,s; |
NM s0,s; |
NODE sp0,sp,rp0,rp; |
NODE sp0,sp,rp0,rp; |
RHist h; |
int nsp,nred,col,rank,len,k,j,a; |
int nsp,nred,col,rank,len,k,j; |
UINT c; |
NMV mr0,mr; |
|
UINT c,c1,c2,c3; |
|
NM mul,head; |
|
UINT **spmat; |
UINT **spmat; |
UINT *s0vect,*rvect,*svect,*p; |
UINT *s0vect,*svect,*p,*v; |
int *colstat; |
int *colstat; |
|
int **imat; |
|
int *rhead; |
|
int spcol,sprow; |
int sugar; |
int sugar; |
PGeoBucket bucket; |
PGeoBucket bucket; |
|
struct oEGT eg0,eg1,eg_f4; |
|
|
if ( !m ) |
if ( !m ) |
error("nd_f4 : not implemented"); |
error("nd_f4 : not implemented"); |
Line 3774 NODE nd_f4(int m) |
|
Line 3956 NODE nd_f4(int m) |
|
g = update_base(g,i); |
g = update_base(g,i); |
} |
} |
while ( d ) { |
while ( d ) { |
|
get_eg(&eg0); |
l = nd_minsugarp(d,&d); |
l = nd_minsugarp(d,&d); |
sugar = SG(l); |
sugar = SG(l); |
fprintf(asir_out,"%d",sugar); |
|
sp0 = 0; |
|
bucket = create_pbucket(); |
bucket = create_pbucket(); |
for ( t = l, nsp = 0; t; t = NEXT(t) ) { |
if ( !(sp0 = nd_sp_f4(m,l,bucket)) |
stat = nd_sp(m,0,t,&spol); |
|| !(col = nd_symbolic_preproc(bucket,&s0vect,&rp0)) ) { |
if ( !stat ) goto again; |
for ( t = l; NEXT(t); t = NEXT(t) ); |
if ( spol ) { |
NEXT(t) = d; d = l; |
NEXTNODE(sp0,sp); BDY(sp) = (pointer)nd_dup(spol); |
d = nd_reconstruct(m,0,d); |
add_pbucket_symbolic(bucket,spol); |
continue; |
nsp++; |
|
} |
|
} |
} |
if ( sp0 ) NEXT(sp) = 0; |
|
s0 = 0; |
nsp = length(sp0); nred = length(rp0); spcol = col-nred; |
rp0 = 0; |
imat = (int **)MALLOC(nred*sizeof(int *)); |
while ( 1 ) { |
rhead = (int *)MALLOC_ATOMIC(col*sizeof(int)); |
head = remove_head_pbucket_symbolic(bucket); |
for ( i = 0; i < col; i++ ) rhead[i] = 0; |
if ( !head ) break; |
|
if ( !s0 ) s0 = head; |
/* construction of index arrays */ |
else NEXT(s) = head; |
for ( rp = rp0, i = 0; rp; i++, rp = NEXT(rp) ) { |
s = head; |
redv = nd_ps[((NM_ind_pair)BDY(rp))->index]; |
index = ndl_find_reducer(DL(head)); |
imat[i] = (int *)MALLOC_ATOMIC(LEN(redv)*sizeof(int)); |
if ( index >= 0 ) { |
nm_ind_pair_to_vect_compress(m,s0vect,col, |
h = nd_psh[index]; |
(NM_ind_pair)BDY(rp),imat[i]); |
NEWNM(mul); |
rhead[imat[i][0]] = 1; |
ndl_sub(DL(head),DL(h),DL(mul)); |
|
if ( ndl_check_bound2(index,DL(mul)) ) goto again; |
|
MKNM_ind_pair(pair,mul,index); |
|
red = ndv_mul_nm_symbolic(mul,nd_ps[index]); |
|
add_pbucket_symbolic(bucket,nd_remove_head(red)); |
|
NEXTNODE(rp0,rp); BDY(rp) = (pointer)pair; |
|
nred++; |
|
} |
|
} |
} |
if ( s0 ) NEXT(s) = 0; |
|
for ( i = 0, s = s0; s; s = NEXT(s), i++ ); |
/* elimination (1st step) */ |
col = i; |
|
s0vect = (UINT *)MALLOC(col*nd_wpd*sizeof(UINT)); |
|
for ( i = 0, p = s0vect, s = s0; i < col; |
|
i++, p += nd_wpd, s = NEXT(s) ) ndl_copy(DL(s),p); |
|
spmat = (UINT **)MALLOC(nsp*sizeof(UINT *)); |
spmat = (UINT **)MALLOC(nsp*sizeof(UINT *)); |
for ( i = 0, sp = sp0; i < nsp; i++, sp = NEXT(sp) ) { |
svect = (UINT *)MALLOC_ATOMIC(col*sizeof(UINT)); |
spmat[i] = (UINT *)MALLOC(col*sizeof(UINT)); |
for ( a = sprow = 0, sp = sp0; a < nsp; a++, sp = NEXT(sp) ) { |
nd_to_vect(m,s0vect,col,BDY(sp),spmat[i]); |
nd_to_vect(m,s0vect,col,BDY(sp),svect); |
} |
ndv_reduce_vect(m,svect,col,imat,rp0); |
rvect = (UINT *)ALLOCA(col*sizeof(UINT)); |
for ( i = 0; i < col; i++ ) if ( svect[i] ) break; |
for ( rp = rp0; rp; rp = NEXT(rp) ) { |
if ( i < col ) { |
k = nm_ind_pair_to_vect(m,s0vect,col,(NM_ind_pair)BDY(rp),rvect); |
spmat[sprow] = v = (UINT *)MALLOC_ATOMIC(spcol*sizeof(UINT)); |
for ( i = 0; i < nsp; i++ ) { |
for ( j = k = 0; j < col; j++ ) |
svect = spmat[i]; |
if ( !rhead[j] ) v[k++] = svect[j]; |
if ( c = svect[k] ) |
sprow++; |
for ( j = k; j < col; j++ ) { |
|
if ( c1 = rvect[j] ) { |
|
c1 = m-c1; c2 = svect[j]; DMAR(c1,c,c2,m,c3); |
|
svect[j] = c3; |
|
} |
|
} |
|
} |
} |
} |
} |
colstat = (int *)ALLOCA(col*sizeof(int)); |
/* free index arrays */ |
rank = generic_gauss_elim_mod(spmat,nsp,col,m,colstat); |
for ( i = 0; i < nred; i++ ) GC_free(imat[i]); |
printf("rank = %d\n",rank); |
|
if ( rank ) |
|
for ( j = 0, i = 0; j < col; j++ ) { |
|
if ( colstat[j] ) { |
|
for ( k = j, len = 0; k < col; k++ ) |
|
if ( spmat[i][k] ) len++; |
|
mr0 = (NMV)MALLOC_ATOMIC(nmv_adv*len); |
|
mr = mr0; |
|
p = s0vect+nd_wpd*j; |
|
for ( k = j; k < col; k++, p += nd_wpd ) |
|
if ( spmat[i][k] ) { |
|
ndl_copy(p,DL(mr)); |
|
CM(mr) = spmat[i][k]; |
|
NMV_ADV(mr); |
|
} |
|
MKNDV(nd_nvar,mr0,len,nf); |
|
SG(nf) = sugar; |
|
ndv_removecont(m,nf); |
|
nh = ndv_newps(nf,0); |
|
d = update_pairs(d,g,nh); |
|
g = update_base(g,nh); |
|
i++; |
|
} |
|
} |
|
continue; |
|
|
|
again: |
/* elimination (2nd step) */ |
for ( t = l; NEXT(t); t = NEXT(t) ); |
colstat = (int *)ALLOCA(spcol*sizeof(int)); |
NEXT(t) = d; d = l; |
rank = generic_gauss_elim_mod(spmat,sprow,spcol,m,colstat); |
d = nd_reconstruct(m,0,d); |
|
NEWNM(mul); |
get_eg(&eg1); init_eg(&eg_f4); add_eg(&eg_f4,&eg0,&eg1); |
|
fprintf(asir_out,"sugar=%d,nsp=%d,nred=%d,spmat=(%d,%d),rank=%d ", |
|
sugar,nsp,nred,sprow,spcol,rank); |
|
fprintf(asir_out,"%fsec\n",eg_f4.exectime+eg_f4.gctime); |
|
|
|
/* adding new bases */ |
|
for ( i = 0; i < rank; i++ ) { |
|
nf = vect_to_ndv(spmat[i],spcol,col,rhead,s0vect); |
|
SG(nf) = sugar; |
|
ndv_removecont(m,nf); |
|
nh = ndv_newps(nf,0); |
|
d = update_pairs(d,g,nh); |
|
g = update_base(g,nh); |
|
GC_free(spmat[i]); |
|
} |
|
for ( ; i < sprow; i++ ) GC_free(spmat[i]); |
} |
} |
for ( r = g; r; r = NEXT(r) ) BDY(r) = (pointer)nd_ps[(int)BDY(r)]; |
for ( r = g; r; r = NEXT(r) ) BDY(r) = (pointer)nd_ps[(int)BDY(r)]; |
return g; |
return g; |