version 1.24, 2001/10/09 01:36:05 |
version 1.29, 2003/06/09 16:18:09 |
|
|
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, |
* DEVELOPER SHALL HAVE NO LIABILITY IN CONNECTION WITH THE USE, |
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. |
* PERFORMANCE OR NON-PERFORMANCE OF THE SOFTWARE. |
* |
* |
* $OpenXM: OpenXM_contrib2/asir2000/builtin/array.c,v 1.23 2001/10/01 01:58:01 noro Exp $ |
* $OpenXM: OpenXM_contrib2/asir2000/builtin/array.c,v 1.28 2003/05/29 16:44:59 saito Exp $ |
*/ |
*/ |
#include "ca.h" |
#include "ca.h" |
#include "base.h" |
#include "base.h" |
Line 64 void Pnewvect(), Pnewmat(), Psepvect(), Psize(), Pdet( |
|
Line 64 void Pnewvect(), Pnewmat(), Psepvect(), Psize(), Pdet( |
|
void Pinvmat(); |
void Pinvmat(); |
void Pnewbytearray(); |
void Pnewbytearray(); |
|
|
|
void Pgeneric_gauss_elim(); |
void Pgeneric_gauss_elim_mod(); |
void Pgeneric_gauss_elim_mod(); |
|
|
void Pmat_to_gfmmat(),Plu_gfmmat(),Psolve_by_lu_gfmmat(); |
void Pmat_to_gfmmat(),Plu_gfmmat(),Psolve_by_lu_gfmmat(); |
void Pgeninvm_swap(), Premainder(), Psremainder(), Pvtol(); |
void Pgeninvm_swap(), Premainder(), Psremainder(), Pvtol(); |
|
void Pgeninv_sf_swap(); |
void sepvect(); |
void sepvect(); |
void Pmulmat_gf2n(); |
void Pmulmat_gf2n(); |
void Pbconvmat_gf2n(); |
void Pbconvmat_gf2n(); |
Line 80 void Pirredpoly_up2(); |
|
Line 82 void Pirredpoly_up2(); |
|
void Pnbpoly_up2(); |
void Pnbpoly_up2(); |
void Pqsort(); |
void Pqsort(); |
void Pexponent_vector(); |
void Pexponent_vector(); |
|
void Pmat_swap_row_destructive(); |
|
void Pmat_swap_col_destructive(); |
|
void Pvect(); |
|
void Pmat(); |
|
void Pmatc(); |
|
|
struct ftab array_tab[] = { |
struct ftab array_tab[] = { |
{"solve_by_lu_gfmmat",Psolve_by_lu_gfmmat,4}, |
{"solve_by_lu_gfmmat",Psolve_by_lu_gfmmat,4}, |
{"lu_gfmmat",Plu_gfmmat,2}, |
{"lu_gfmmat",Plu_gfmmat,2}, |
{"mat_to_gfmmat",Pmat_to_gfmmat,2}, |
{"mat_to_gfmmat",Pmat_to_gfmmat,2}, |
|
{"generic_gauss_elim",Pgeneric_gauss_elim,1}, |
{"generic_gauss_elim_mod",Pgeneric_gauss_elim_mod,2}, |
{"generic_gauss_elim_mod",Pgeneric_gauss_elim_mod,2}, |
{"newvect",Pnewvect,-2}, |
{"newvect",Pnewvect,-2}, |
|
{"vect",Pvect,-99999999}, |
{"vector",Pnewvect,-2}, |
{"vector",Pnewvect,-2}, |
{"exponent_vector",Pexponent_vector,-99999999}, |
{"exponent_vector",Pexponent_vector,-99999999}, |
{"newmat",Pnewmat,-3}, |
{"newmat",Pnewmat,-3}, |
{"matrix",Pnewmat,-3}, |
{"matrix",Pnewmat,-3}, |
|
{"mat",Pmat,-99999999}, |
|
{"matr",Pmat,-99999999}, |
|
{"matc",Pmatc,-99999999}, |
{"newbytearray",Pnewbytearray,-2}, |
{"newbytearray",Pnewbytearray,-2}, |
{"sepmat_destructive",Psepmat_destructive,2}, |
{"sepmat_destructive",Psepmat_destructive,2}, |
{"sepvect",Psepvect,2}, |
{"sepvect",Psepvect,2}, |
Line 103 struct ftab array_tab[] = { |
|
Line 115 struct ftab array_tab[] = { |
|
{"leqm1",Pleqm1,2}, |
{"leqm1",Pleqm1,2}, |
{"geninvm",Pgeninvm,2}, |
{"geninvm",Pgeninvm,2}, |
{"geninvm_swap",Pgeninvm_swap,2}, |
{"geninvm_swap",Pgeninvm_swap,2}, |
|
{"geninv_sf_swap",Pgeninv_sf_swap,1}, |
{"remainder",Premainder,2}, |
{"remainder",Premainder,2}, |
{"sremainder",Psremainder,2}, |
{"sremainder",Psremainder,2}, |
{"mulmat_gf2n",Pmulmat_gf2n,1}, |
{"mulmat_gf2n",Pmulmat_gf2n,1}, |
Line 113 struct ftab array_tab[] = { |
|
Line 126 struct ftab array_tab[] = { |
|
{"x962_irredpoly_up2",Px962_irredpoly_up2,2}, |
{"x962_irredpoly_up2",Px962_irredpoly_up2,2}, |
{"irredpoly_up2",Pirredpoly_up2,2}, |
{"irredpoly_up2",Pirredpoly_up2,2}, |
{"nbpoly_up2",Pnbpoly_up2,2}, |
{"nbpoly_up2",Pnbpoly_up2,2}, |
|
{"mat_swap_row_destructive",Pmat_swap_row_destructive,3}, |
|
{"mat_swap_col_destructive",Pmat_swap_col_destructive,3}, |
{0,0,0}, |
{0,0,0}, |
}; |
}; |
|
|
Line 347 void Pnewvect(NODE arg,VECT *rp) |
|
Line 362 void Pnewvect(NODE arg,VECT *rp) |
|
*rp = vect; |
*rp = vect; |
} |
} |
|
|
|
void Pvect(NODE arg,VECT *rp) { |
|
int len,i,r; |
|
VECT vect; |
|
pointer *vb; |
|
NODE tn; |
|
|
|
if ( !arg ) { |
|
*rp =0; |
|
return; |
|
} |
|
|
|
for (len = 0, tn = arg; tn; tn = NEXT(tn), len++); |
|
if ( len == 1 ) { |
|
if ( ARG0(arg) != 0 ) { |
|
switch ( OID(ARG0(arg)) ) { |
|
case O_VECT: |
|
*rp = ARG0(arg); |
|
return; |
|
case O_LIST: |
|
for ( len = 0, tn = ARG0(arg); tn; tn = NEXT(tn), len++ ); |
|
MKVECT(vect,len-1); |
|
for ( i = 0, tn = BDY((LIST)ARG0(arg)), vb =BDY(vect); |
|
tn; i++, tn = NEXT(tn) ) |
|
vb[i] = (pointer)BDY(tn); |
|
*rp=vect; |
|
return; |
|
} |
|
} |
|
} |
|
MKVECT(vect,len); |
|
for ( i = 0, tn = arg, vb = BDY(vect); tn; i++, tn = NEXT(tn) ) |
|
vb[i] = (pointer)BDY(tn); |
|
*rp = vect; |
|
} |
|
|
void Pexponent_vector(NODE arg,DP *rp) |
void Pexponent_vector(NODE arg,DP *rp) |
{ |
{ |
nodetod(arg,rp); |
nodetod(arg,rp); |
Line 429 void Pnewmat(NODE arg,MAT *rp) |
|
Line 479 void Pnewmat(NODE arg,MAT *rp) |
|
*rp = m; |
*rp = m; |
} |
} |
|
|
|
void Pmat(NODE arg, MAT *rp) |
|
{ |
|
int row,col; |
|
int i; |
|
MAT m; |
|
pointer **mb; |
|
pointer *ent; |
|
NODE tn, sn; |
|
VECT v; |
|
|
|
if ( !arg ) { |
|
*rp =0; |
|
return; |
|
} |
|
|
|
for (row = 0, tn = arg; tn; tn = NEXT(tn), row++); |
|
if ( OID(ARG0(arg)) == O_VECT ) { |
|
v = ARG0(arg); |
|
col = v->len; |
|
} else if ( OID(ARG0(arg)) == O_LIST ) { |
|
for (col = 0, tn = BDY((LIST)ARG0(arg)); tn ; tn = NEXT(tn), col++); |
|
} |
|
|
|
MKMAT(m,row,col); |
|
for (row = 0, tn = arg, mb = BDY(m); tn; tn = NEXT(tn), row++) { |
|
if ( BDY(tn) == 0 ) { |
|
error("mat : invalid argument"); |
|
} else if ( OID(BDY(tn)) == O_VECT ) { |
|
v = tn->body; |
|
ent = BDY(v); |
|
for (i = 0; i < v->len; i++ ) mb[row][i] = (Obj)ent[i]; |
|
} else if ( OID(BDY(tn)) == O_LIST ) { |
|
for (col = 0, sn = BDY((LIST)BDY(tn)); sn; col++, sn = NEXT(sn) ) |
|
mb[row][col] = (pointer)BDY(sn); |
|
} else { |
|
error("mat : invalid argument"); |
|
} |
|
} |
|
*rp = m; |
|
} |
|
|
|
void Pmatc(NODE arg, MAT *rp) |
|
{ |
|
int row,col; |
|
int i; |
|
MAT m; |
|
pointer **mb; |
|
pointer *ent; |
|
NODE tn, sn; |
|
VECT v; |
|
|
|
if ( !arg ) { |
|
*rp =0; |
|
return; |
|
} |
|
|
|
for (col = 0, tn = arg; tn; tn = NEXT(tn), col++); |
|
if ( OID(ARG0(arg)) == O_VECT ) { |
|
v = ARG0(arg); |
|
row = v->len; |
|
} else if ( OID(ARG0(arg)) == O_LIST ) { |
|
for (row = 0, tn = BDY((LIST)ARG0(arg)); tn ; tn = NEXT(tn), row++); |
|
} |
|
|
|
MKMAT(m,row,col); |
|
for (col = 0, tn = arg, mb = BDY(m); tn; tn = NEXT(tn), col++) { |
|
if ( BDY(tn) == 0 ) { |
|
error("matc : invalid argument"); |
|
} else if ( OID(BDY(tn)) == O_VECT ) { |
|
v = tn->body; |
|
ent = BDY(v); |
|
for (i = 0; i < v->len; i++ ) mb[i][col] = (Obj)ent[i]; |
|
} else if ( OID(BDY(tn)) == O_LIST ) { |
|
for (row = 0, sn = BDY((LIST)BDY(tn)); sn; row++, sn = NEXT(sn) ) |
|
mb[row][col] = (pointer)BDY(sn); |
|
} else { |
|
error("matc : invalid argument"); |
|
} |
|
} |
|
*rp = m; |
|
} |
|
|
void Pvtol(NODE arg,LIST *rp) |
void Pvtol(NODE arg,LIST *rp) |
{ |
{ |
NODE n,n1; |
NODE n,n1; |
Line 635 void Pinvmat(NODE arg,LIST *rp) |
|
Line 767 void Pinvmat(NODE arg,LIST *rp) |
|
B[I] <-> x_{R[I]}+B[I][0]x_{C[0]}+B[I][1]x_{C[1]}+... |
B[I] <-> x_{R[I]}+B[I][0]x_{C[0]}+B[I][1]x_{C[1]}+... |
*/ |
*/ |
|
|
|
void Pgeneric_gauss_elim(NODE arg,LIST *rp) |
|
{ |
|
NODE n0; |
|
MAT m,nm; |
|
int *ri,*ci; |
|
VECT rind,cind; |
|
Q dn,q; |
|
int i,j,k,l,row,col,t,rank; |
|
|
|
asir_assert(ARG0(arg),O_MAT,"generic_gauss_elim"); |
|
m = (MAT)ARG0(arg); |
|
row = m->row; col = m->col; |
|
rank = generic_gauss_elim(m,&nm,&dn,&ri,&ci); |
|
t = col-rank; |
|
MKVECT(rind,rank); |
|
MKVECT(cind,t); |
|
for ( i = 0; i < rank; i++ ) { |
|
STOQ(ri[i],q); |
|
BDY(rind)[i] = (pointer)q; |
|
} |
|
for ( i = 0; i < t; i++ ) { |
|
STOQ(ci[i],q); |
|
BDY(cind)[i] = (pointer)q; |
|
} |
|
n0 = mknode(4,nm,dn,rind,cind); |
|
MKLIST(*rp,n0); |
|
} |
|
|
|
/* |
|
input : a row x col matrix A |
|
A[I] <-> A[I][0]*x_0+A[I][1]*x_1+... |
|
|
|
output : [B,R,C] |
|
B : a rank(A) x col-rank(A) matrix |
|
R : a vector of length rank(A) |
|
C : a vector of length col-rank(A) |
|
B[I] <-> x_{R[I]}+B[I][0]x_{C[0]}+B[I][1]x_{C[1]}+... |
|
*/ |
|
|
void Pgeneric_gauss_elim_mod(NODE arg,LIST *rp) |
void Pgeneric_gauss_elim_mod(NODE arg,LIST *rp) |
{ |
{ |
NODE n0; |
NODE n0; |
Line 949 int generic_gauss_elim_hensel(MAT mat,MAT *nmmat,Q *dn |
|
Line 1120 int generic_gauss_elim_hensel(MAT mat,MAT *nmmat,Q *dn |
|
} else |
} else |
wi[j] = 0; |
wi[j] = 0; |
|
|
rank = find_lhs_and_lu_mod(w,row,col,md,&rinfo,&cinfo); |
rank = find_lhs_and_lu_mod((unsigned int **)w,row,col,md,&rinfo,&cinfo); |
a = (Q **)almat_pointer(rank,rank); /* lhs mat */ |
a = (Q **)almat_pointer(rank,rank); /* lhs mat */ |
MKMAT(bmat,rank,col-rank); b = (Q **)bmat->body; /* lhs mat */ |
MKMAT(bmat,rank,col-rank); b = (Q **)bmat->body; /* lhs mat */ |
for ( j = li = ri = 0; j < col; j++ ) |
for ( j = li = ri = 0; j < col; j++ ) |
Line 1894 void mat_to_gfmmat(MAT m,unsigned int md,GFMMAT *rp) |
|
Line 2065 void mat_to_gfmmat(MAT m,unsigned int md,GFMMAT *rp) |
|
TOGFMMAT(row,col,wmat,*rp); |
TOGFMMAT(row,col,wmat,*rp); |
} |
} |
|
|
void Pgeninvm_swap(NODE arg,LIST *rp) |
void Pgeninvm_swap(arg,rp) |
|
NODE arg; |
|
LIST *rp; |
{ |
{ |
MAT m; |
MAT m; |
pointer **mat; |
pointer **mat; |
Line 1940 void Pgeninvm_swap(NODE arg,LIST *rp) |
|
Line 2113 void Pgeninvm_swap(NODE arg,LIST *rp) |
|
} |
} |
} |
} |
|
|
int gauss_elim_geninv_mod_swap(unsigned int **mat,int row,int col, |
gauss_elim_geninv_mod_swap(mat,row,col,md,invmatp,indexp) |
unsigned int md,unsigned int ***invmatp,int **indexp) |
unsigned int **mat; |
|
int row,col; |
|
unsigned int md; |
|
unsigned int ***invmatp; |
|
int **indexp; |
{ |
{ |
int i,j,k,inv,a,n,m; |
int i,j,k,inv,a,n,m; |
unsigned int *t,*pivot,*s; |
unsigned int *t,*pivot,*s; |
Line 1991 int gauss_elim_geninv_mod_swap(unsigned int **mat,int |
|
Line 2168 int gauss_elim_geninv_mod_swap(unsigned int **mat,int |
|
return 0; |
return 0; |
} |
} |
|
|
|
void Pgeninv_sf_swap(NODE arg,LIST *rp) |
|
{ |
|
MAT m; |
|
GFS **mat,**tmat; |
|
Q *tvect; |
|
GFS q; |
|
int **wmat,**invmat; |
|
int *index; |
|
unsigned int t; |
|
int i,j,row,col,status; |
|
MAT mat1; |
|
VECT vect1; |
|
NODE node1,node2; |
|
|
|
asir_assert(ARG0(arg),O_MAT,"geninv_sf_swap"); |
|
m = (MAT)ARG0(arg); |
|
row = m->row; col = m->col; mat = (GFS **)m->body; |
|
wmat = (int **)almat(row,col+row); |
|
for ( i = 0; i < row; i++ ) { |
|
bzero((char *)wmat[i],(col+row)*sizeof(int)); |
|
for ( j = 0; j < col; j++ ) |
|
if ( q = (GFS)mat[i][j] ) |
|
wmat[i][j] = FTOIF(CONT(q)); |
|
wmat[i][col+i] = _onesf(); |
|
} |
|
status = gauss_elim_geninv_sf_swap(wmat,row,col,&invmat,&index); |
|
if ( status > 0 ) |
|
*rp = 0; |
|
else { |
|
MKMAT(mat1,col,col); |
|
for ( i = 0, tmat = (GFS **)mat1->body; i < col; i++ ) |
|
for ( j = 0; j < col; j++ ) |
|
if ( t = invmat[i][j] ) { |
|
MKGFS(IFTOF(t),tmat[i][j]); |
|
} |
|
MKVECT(vect1,row); |
|
for ( i = 0, tvect = (Q *)vect1->body; i < row; i++ ) |
|
STOQ(index[i],tvect[i]); |
|
MKNODE(node2,vect1,0); MKNODE(node1,mat1,node2); MKLIST(*rp,node1); |
|
} |
|
} |
|
|
|
int gauss_elim_geninv_sf_swap(int **mat,int row,int col, |
|
int ***invmatp,int **indexp) |
|
{ |
|
int i,j,k,inv,a,n,m,u; |
|
int *t,*pivot,*s; |
|
int *index; |
|
int **invmat; |
|
|
|
n = col; m = row+col; |
|
*indexp = index = (int *)MALLOC_ATOMIC(row*sizeof(int)); |
|
for ( i = 0; i < row; i++ ) |
|
index[i] = i; |
|
for ( j = 0; j < n; j++ ) { |
|
for ( i = j; i < row && !mat[i][j]; i++ ); |
|
if ( i == row ) { |
|
*indexp = 0; *invmatp = 0; return 1; |
|
} |
|
if ( i != j ) { |
|
t = mat[i]; mat[i] = mat[j]; mat[j] = t; |
|
k = index[i]; index[i] = index[j]; index[j] = k; |
|
} |
|
pivot = mat[j]; |
|
inv = _invsf(pivot[j]); |
|
for ( k = j; k < m; k++ ) |
|
if ( pivot[k] ) |
|
pivot[k] = _mulsf(pivot[k],inv); |
|
for ( i = j+1; i < row; i++ ) { |
|
t = mat[i]; |
|
if ( a = t[j] ) |
|
for ( k = j, a = _chsgnsf(a); k < m; k++ ) |
|
if ( pivot[k] ) { |
|
u = _mulsf(pivot[k],a); |
|
t[k] = _addsf(u,t[k]); |
|
} |
|
} |
|
} |
|
for ( j = n-1; j >= 0; j-- ) { |
|
pivot = mat[j]; |
|
for ( i = j-1; i >= 0; i-- ) { |
|
t = mat[i]; |
|
if ( a = t[j] ) |
|
for ( k = j, a = _chsgnsf(a); k < m; k++ ) |
|
if ( pivot[k] ) { |
|
u = _mulsf(pivot[k],a); |
|
t[k] = _addsf(u,t[k]); |
|
} |
|
} |
|
} |
|
*invmatp = invmat = (int **)almat(col,col); |
|
for ( i = 0; i < col; i++ ) |
|
for ( j = 0, s = invmat[i], t = mat[i]; j < col; j++ ) |
|
s[j] = t[col+index[j]]; |
|
return 0; |
|
} |
|
|
void _addn(N,N,N); |
void _addn(N,N,N); |
int _subn(N,N,N); |
int _subn(N,N,N); |
void _muln(N,N,N); |
void _muln(N,N,N); |
Line 2185 void Pirredpoly_up2(NODE arg,GF2N *rp) |
|
Line 2459 void Pirredpoly_up2(NODE arg,GF2N *rp) |
|
*rp = 0; |
*rp = 0; |
} |
} |
|
|
|
void Pmat_swap_row_destructive(NODE arg, MAT *m) |
|
{ |
|
int i1,i2; |
|
pointer *t; |
|
MAT mat; |
|
|
|
asir_assert(ARG0(arg),O_MAT,"mat_swap_row_destructive"); |
|
asir_assert(ARG1(arg),O_N,"mat_swap_row_destructive"); |
|
asir_assert(ARG2(arg),O_N,"mat_swap_row_destructive"); |
|
mat = (MAT)ARG0(arg); |
|
i1 = QTOS((Q)ARG1(arg)); |
|
i2 = QTOS((Q)ARG2(arg)); |
|
if ( i1 < 0 || i2 < 0 || i1 >= mat->row || i2 >= mat->row ) |
|
error("mat_swap_row_destructive : Out of range"); |
|
t = mat->body[i1]; |
|
mat->body[i1] = mat->body[i2]; |
|
mat->body[i2] = t; |
|
*m = mat; |
|
} |
|
|
|
void Pmat_swap_col_destructive(NODE arg, MAT *m) |
|
{ |
|
int j1,j2,i,n; |
|
pointer *mi; |
|
pointer t; |
|
MAT mat; |
|
|
|
asir_assert(ARG0(arg),O_MAT,"mat_swap_col_destructive"); |
|
asir_assert(ARG1(arg),O_N,"mat_swap_col_destructive"); |
|
asir_assert(ARG2(arg),O_N,"mat_swap_col_destructive"); |
|
mat = (MAT)ARG0(arg); |
|
j1 = QTOS((Q)ARG1(arg)); |
|
j2 = QTOS((Q)ARG2(arg)); |
|
if ( j1 < 0 || j2 < 0 || j1 >= mat->col || j2 >= mat->col ) |
|
error("mat_swap_col_destructive : Out of range"); |
|
n = mat->row; |
|
for ( i = 0; i < n; i++ ) { |
|
mi = mat->body[i]; |
|
t = mi[j1]; mi[j1] = mi[j2]; mi[j2] = t; |
|
} |
|
*m = mat; |
|
} |
/* |
/* |
* f = type 'type' normal polynomial of degree m if exists |
* f = type 'type' normal polynomial of degree m if exists |
* IEEE P1363 A.7.2 |
* IEEE P1363 A.7.2 |