===================================================================
RCS file: /home/cvs/OpenXM_contrib2/asir2000/builtin/array.c,v
retrieving revision 1.69
retrieving revision 1.75
diff -u -p -r1.69 -r1.75
--- OpenXM_contrib2/asir2000/builtin/array.c	2015/09/03 23:05:35	1.69
+++ OpenXM_contrib2/asir2000/builtin/array.c	2017/09/17 02:34:02	1.75
@@ -45,7 +45,7 @@
  * DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE,
  * PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE.
  *
- * $OpenXM: OpenXM_contrib2/asir2000/builtin/array.c,v 1.68 2015/08/14 13:51:54 fujimoto Exp $
+ * $OpenXM: OpenXM_contrib2/asir2000/builtin/array.c,v 1.74 2017/09/15 01:52:51 noro Exp $
 */
 #include "ca.h"
 #include "base.h"
@@ -68,7 +68,7 @@
 extern int DP_Print; /* XXX */
 
 
-void Pnewvect(), Pnewmat(), Psepvect(), Psize(), Pdet(), Pleqm(), Pleqm1(), Pgeninvm();
+void Pnewvect(), Pnewmat(), Psepvect(), Psize(), Pdet(), Pleqm(), Pleqm1(), Pgeninvm(), Ptriangleq();
 void Pinvmat();
 void Pnewbytearray(),Pmemoryplot_to_coord();
 
@@ -151,6 +151,7 @@ struct ftab array_tab[] = {
 	{"mat_col",Pmat_col,2},
 	{"lusolve_prep",Plusolve_prep,1},
 	{"lusolve_main",Plusolve_main,1},
+	{"triangleq",Ptriangleq,1},
 	{0,0,0},
 };
 
@@ -158,8 +159,12 @@ typedef struct _ent { int j; unsigned int e; } ent;
 
 ent *get_row(FILE *,int *l);
 void put_row(FILE *out,int l,ent *a);
-int lu_elim(int *l,ent **a,int k,int i,int mul,int mod);
+void lu_elim(int *l,ent **a,int k,int i,int mul,int mod);
+void lu_append(int *,ent **,int *,int,int,int);
+void solve_l(int *,ent **,int,int *,int);
+void solve_u(int *,ent **,int,int *,int);
 
+
 static int *ul,*ll;
 static ent **u,**l;
 static int modulus;
@@ -224,7 +229,7 @@ ent *get_row(FILE *in,int *l)
 	return a;
 }
 
-int lu_gauss(int *ul,ent **u,int *ll,ent **l,int n,int mod)
+void lu_gauss(int *ul,ent **u,int *ll,ent **l,int n,int mod)
 {
 	int i,j,k,s,mul;
 	unsigned int inv;
@@ -247,7 +252,7 @@ int lu_gauss(int *ul,ent **u,int *ll,ent **l,int n,int
 
 #define INITLEN 10
 
-lu_append(int *l,ent **a,int *l2,int k,int i,int mul)
+void lu_append(int *l,ent **a,int *l2,int k,int i,int mul)
 {
 	int len;
 	ent *p;
@@ -270,7 +275,7 @@ lu_append(int *l,ent **a,int *l2,int k,int i,int mul)
 
 /* a[k] = a[k]-mul*a[i] */ 
 
-int lu_elim(int *l,ent **a,int k,int i,int mul,int mod)
+void lu_elim(int *l,ent **a,int k,int i,int mul,int mod)
 {
 	ent *ak,*ai,*w;
 	int lk,li,j,m,p,q,r,s,t,j0;
@@ -318,7 +323,7 @@ int lu_elim(int *l,ent **a,int k,int i,int mul,int mod
 	l[k] = j;
 }
 
-int solve_l(int *ll,ent **l,int n,int *rhs,int mod)
+void solve_l(int *ll,ent **l,int n,int *rhs,int mod)
 {
 	int j,k,s,len;
 	ent *p;
@@ -332,7 +337,7 @@ int solve_l(int *ll,ent **l,int n,int *rhs,int mod)
 	}
 }
 
-int solve_u(int *ul,ent **u,int n,int *rhs,int mod)
+void solve_u(int *ul,ent **u,int n,int *rhs,int mod)
 {
 	int j,k,s,len,inv;
 	ent *p;
@@ -2351,6 +2356,24 @@ void red_by_vect(int m,unsigned int *p,unsigned int *r
 		}
 }
 
+#if defined(__GNUC__) && SIZEOF_LONG==8
+/* 64bit vector += UNIT vector(normalized) */
+
+void red_by_vect64(int m, U64 *p,unsigned int *c,U64 *r,unsigned int hc,int len)
+{
+  U64 t;
+
+  /* (p[0],c[0]) is normalized */
+  *p++ = 0; *c++ = 0; r++; len--;
+  for ( ; len; len--, r++, p++, c++ )
+    if ( *r ) {
+      t = (*p)+(*r)*hc;
+      if ( t < *p ) (*c)++;
+      *p = t;
+    }
+}
+#endif
+
 void red_by_vect_sf(int m,unsigned int *p,unsigned int *r,unsigned int hc,int len)
 {
 	*p++ = 0; r++; len--;
@@ -2359,6 +2382,22 @@ void red_by_vect_sf(int m,unsigned int *p,unsigned int
 			*p = _addsf(_mulsf(*r,hc),*p);
 }
 
+extern GZ current_mod_lf;
+extern int current_mod_lf_size;
+
+void red_by_vect_lf(mpz_t *p,mpz_t *r,mpz_t hc,int len)
+{
+	mpz_set_ui(*p++,0); r++; len--;
+	for ( ; len; len--, r++, p++ ) {
+       mpz_addmul(*p,*r,hc);
+#if 0
+       if ( mpz_size(*p) > current_mod_lf_size )
+         mpz_mod(*p,*p,BDY(current_mod_lf));
+#endif
+    }
+}
+
+
 extern unsigned int **psca;
 
 void reduce_sp_by_red_mod_compress (int *sp,CDP *redmat,int *ind,
@@ -3070,9 +3109,7 @@ void mat_to_gfmmat(MAT m,unsigned int md,GFMMAT *rp)
 	TOGFMMAT(row,col,wmat,*rp);
 }
 
-void Pgeninvm_swap(arg,rp)
-NODE arg;
-LIST *rp;
+void Pgeninvm_swap(NODE arg,LIST *rp)
 {
 	MAT m;
 	pointer **mat;
@@ -3118,12 +3155,8 @@ LIST *rp;
 	}
 }
 
-gauss_elim_geninv_mod_swap(mat,row,col,md,invmatp,indexp)
-unsigned int **mat;
-int row,col;
-unsigned int md;
-unsigned int ***invmatp;
-int **indexp;
+int gauss_elim_geninv_mod_swap(unsigned int **mat,int row,int col,unsigned int md,
+    unsigned int ***invmatp,int **indexp)
 {
 	int i,j,k,inv,a,n,m;
 	unsigned int *t,*pivot,*s;
@@ -3830,4 +3863,45 @@ void Pmat_col(NODE arg,VECT *rp)
 		BDY(vect)[i] = BDY(mat)[i][j];
 	}
 	*rp = vect;
+}
+
+NODE triangleq(NODE e)
+{
+  int n,i,k;
+  V v;
+  VL vl;
+  P *p;
+  NODE r,r1;
+
+  n = length(e);
+  p = (P *)MALLOC(n*sizeof(P));
+  for ( i = 0; i < n; i++, e = NEXT(e) ) p[i] = (P)BDY(e);
+  i = 0;
+  while ( 1 ) {
+    for ( ; i < n && !p[i]; i++ );
+    if ( i == n ) break;
+    if ( OID(p[i]) == O_N ) return 0;
+    v = p[i]->v;
+    for ( k = i+1; k < n; k++ )
+      if ( p[k] ) {
+        if ( OID(p[k]) == O_N ) return 0;
+        if ( p[k]->v == v ) p[k] = 0;
+      }
+    i++;
+  }
+  for ( r = 0, i = 0; i < n; i++ ) {
+    if ( p[i] ) {
+      MKNODE(r1,p[i],r); r = r1;
+    }
+  }
+  return r;
+}
+
+void Ptriangleq(NODE arg,LIST *rp)
+{
+  NODE ret;
+
+  asir_assert(ARG0(arg),O_LIST,"sparseleq");
+  ret = triangleq(BDY((LIST)ARG0(arg)));
+  MKLIST(*rp,ret);
 }