[BACK]Return to nd.c CVS log [TXT][DIR] Up to [local] / OpenXM_contrib2 / asir2000 / engine

Diff for /OpenXM_contrib2/asir2000/engine/nd.c between version 1.64 and 1.65

version 1.64, 2003/09/12 01:12:40 version 1.65, 2003/09/12 08:01:51
Line 1 
Line 1 
 /* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.63 2003/09/11 09:03:53 noro Exp $ */  /* $OpenXM: OpenXM_contrib2/asir2000/engine/nd.c,v 1.64 2003/09/12 01:12:40 noro Exp $ */
   
 #include "ca.h"  #include "ca.h"
 #include "inline.h"  #include "inline.h"
Line 15 
Line 15 
 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 348  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,UINT *r,int *ind);  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 413  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 479  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 488  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 555  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 622  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 633  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 838  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 904  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 916  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 926  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 992  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 1003  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 1011  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 1077  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 1088  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 3746  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,UINT *r,int *ind)  void nm_ind_pair_to_vect_compress(int mod,UINT *s0,int n,NM_ind_pair pair,int *ind)
 {  {
         NM m;          NM m;
         NMV mr;          NMV mr;
Line 3762  void nm_ind_pair_to_vect_compress(int mod,UINT *s0,int
Line 3803  void nm_ind_pair_to_vect_compress(int mod,UINT *s0,int
         for ( i = j = 0, s = s0, mr = BDY(p); 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);                  ndl_add(d,DL(mr),t);
                 for ( ; !ndl_equal(t,s); s += nd_wpd, i++ );                  for ( ; !ndl_equal(t,s); s += nd_wpd, i++ );
                 r[j] = CM(mr);  
                 ind[j] = i;                  ind[j] = i;
         }          }
 }  }
   
   
   void ndv_reduce_vect(int m,UINT *svect,int **imat,NODE rp0)
   {
           int i,j,k,len,pos;
           UINT c,c1,c2,c3;
           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];
                   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]; DMAR(c1,c,c2,m,c3);
                                   svect[pos] = c3;
                           }
                   }
           }
   }
   
   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,*ivect;          UINT *s0vect,*svect,*p,*v;
         int *colstat;          int *colstat;
         int sugar,pos;          int **imat;
           int *rhead;
           int spcol,sprow;
           int sugar;
         PGeoBucket bucket;          PGeoBucket bucket;
         int t_0,t_1,t_2,t_3,t_4,t_5,t_c,t_e,t_20,t_21;  
   
         if ( !m )          if ( !m )
                 error("nd_f4 : not implemented");                  error("nd_f4 : not implemented");
Line 3800  NODE nd_f4(int m)
Line 3949  NODE nd_f4(int m)
         while ( d ) {          while ( d ) {
                 l = nd_minsugarp(d,&d);                  l = nd_minsugarp(d,&d);
                 sugar = SG(l);                  sugar = SG(l);
                 sp0 = 0;  
                 bucket = create_pbucket();                  bucket = create_pbucket();
                 t_0 = clock();                  if ( !(sp0 = nd_sp_f4(m,l,bucket))
                 for ( t = l, nsp = 0; t; t = NEXT(t) ) {                          || !(col = nd_symbolic_preproc(bucket,&s0vect,&rp0)) ) {
                         stat = nd_sp(m,0,t,&spol);                          for ( t = l; NEXT(t); t = NEXT(t) );
                         if ( !stat ) goto again;                          NEXT(t) = d; d = l;
                         if ( spol ) {                          d = nd_reconstruct(m,0,d);
                                 NEXTNODE(sp0,sp); BDY(sp) = (pointer)nd_dup(spol);                          continue;
                                 add_pbucket_symbolic(bucket,spol);  
                                 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++;  
                         }  
                 }                  }
                 t_1 = clock();  
                 if ( s0 ) NEXT(s) = 0;                  /* elimination (1st step) */
                 for ( i = 0, s = s0; s; s = NEXT(s), i++ );  
                 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,imat,rp0);
                 t_2 = clock();                          for ( i = 0; i < col; i++ ) if ( svect[i] ) break;
                 t_c = 0;                          if ( i < col ) {
                 t_e = 0;                                  spmat[sprow] = v = (UINT *)MALLOC_ATOMIC(spcol*sizeof(UINT));
                 rvect = (UINT *)ALLOCA(col*sizeof(UINT));                                  for ( j = k = 0; j < col; j++ )
                 ivect = (int *)ALLOCA(col*sizeof(int));                                          if ( !rhead[j] ) v[k++] = svect[j];
                 for ( rp = rp0; rp; rp = NEXT(rp) ) {                                  sprow++;
                         t_20 = clock();  
                         nm_ind_pair_to_vect_compress(m,s0vect,col,(NM_ind_pair)BDY(rp),rvect,ivect);  
                         t_21 = clock();  
                         t_c += t_21-t_20;  
                         k = ivect[0];  
                         len = LEN(nd_ps[((NM_ind_pair)BDY(rp))->index]);  
                         for ( i = 0; i < nsp; i++ ) {  
                                 svect = spmat[i];  
                                 if ( c = svect[k] )  
                                         for ( j = 0; j < len; j++ ) {  
                                                 pos = ivect[j];  
                                                 c1 = m-rvect[j];  
                                                 c2 = svect[pos]; DMAR(c1,c,c2,m,c3);  
                                                 svect[pos] = c3;  
                                         }  
                         }                          }
                         t_e +=  clock()-t_21;  
                 }                  }
                 t_4 = clock();                  /* free index arrays */
                 colstat = (int *)ALLOCA(col*sizeof(int));                  for ( i = 0; i < nred; i++ ) GC_free(imat[i]);
                 rank = generic_gauss_elim_mod(spmat,nsp,col,m,colstat);  
                 t_5 = clock();  
                 fprintf(asir_out,"sugar=%d,rank=%d ",sugar,rank);  
                 fprintf(asir_out,"symb=%f,conv1=%f,conv2=%f,elim1=%f,elim2=%f\n",  
                 (t_1-t_0)/(double)CLOCKS_PER_SEC,  
                 (t_2-t_1)/(double)CLOCKS_PER_SEC,  
                 (t_c)/(double)CLOCKS_PER_SEC,  
                 (t_e)/(double)CLOCKS_PER_SEC,  
                 (t_5-t_4)/(double)CLOCKS_PER_SEC);  
                 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);                  fprintf(asir_out,"sugar=%d,rank=%d\n",sugar,rank); fflush(asir_out);
   
                   /* adding new bases */
                   for ( i = sprow-1; i >= rank; i-- ) GC_free(spmat[i]);
                   for ( ; i >= 0; 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 ( 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;
Legend:
Removed from v.1.64  
changed lines
  Added in v.1.65
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>