version 1.3, 2000/04/20 02:20:15 |
version 1.7, 2000/08/22 05:03:56 |
|
|
/* $OpenXM: OpenXM_contrib2/asir2000/builtin/array.c,v 1.2 2000/03/14 05:25:43 noro Exp $ */ |
/* |
|
* Copyright (c) 1994-2000 FUJITSU LABORATORIES LIMITED |
|
* All rights reserved. |
|
* |
|
* FUJITSU LABORATORIES LIMITED ("FLL") hereby grants you a limited, |
|
* non-exclusive and royalty-free license to use, copy, modify and |
|
* redistribute, solely for non-commercial and non-profit purposes, the |
|
* computer program, "Risa/Asir" ("SOFTWARE"), subject to the terms and |
|
* conditions of this Agreement. For the avoidance of doubt, you acquire |
|
* only a limited right to use the SOFTWARE hereunder, and FLL or any |
|
* third party developer retains all rights, including but not limited to |
|
* copyrights, in and to the SOFTWARE. |
|
* |
|
* (1) FLL does not grant you a license in any way for commercial |
|
* purposes. You may use the SOFTWARE only for non-commercial and |
|
* non-profit purposes only, such as academic, research and internal |
|
* business use. |
|
* (2) The SOFTWARE is protected by the Copyright Law of Japan and |
|
* international copyright treaties. If you make copies of the SOFTWARE, |
|
* with or without modification, as permitted hereunder, you shall affix |
|
* to all such copies of the SOFTWARE the above copyright notice. |
|
* (3) An explicit reference to this SOFTWARE and its copyright owner |
|
* shall be made on your publication or presentation in any form of the |
|
* results obtained by use of the SOFTWARE. |
|
* (4) In the event that you modify the SOFTWARE, you shall notify FLL by |
|
* e-mail at risa-admin@sec.flab.fujitsu.co.jp of the detailed specification |
|
* for such modification or the source code of the modified part of the |
|
* SOFTWARE. |
|
* |
|
* THE SOFTWARE IS PROVIDED AS IS WITHOUT ANY WARRANTY OF ANY KIND. FLL |
|
* MAKES ABSOLUTELY NO WARRANTIES, EXPRESSED, IMPLIED OR STATUTORY, AND |
|
* EXPRESSLY DISCLAIMS ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS |
|
* FOR A PARTICULAR PURPOSE OR NONINFRINGEMENT OF THIRD PARTIES' |
|
* RIGHTS. NO FLL DEALER, AGENT, EMPLOYEES IS AUTHORIZED TO MAKE ANY |
|
* MODIFICATIONS, EXTENSIONS, OR ADDITIONS TO THIS WARRANTY. |
|
* UNDER NO CIRCUMSTANCES AND UNDER NO LEGAL THEORY, TORT, CONTRACT, |
|
* OR OTHERWISE, SHALL FLL BE LIABLE TO YOU OR ANY OTHER PERSON FOR ANY |
|
* DIRECT, INDIRECT, SPECIAL, INCIDENTAL, PUNITIVE OR CONSEQUENTIAL |
|
* DAMAGES OF ANY CHARACTER, INCLUDING, WITHOUT LIMITATION, DAMAGES |
|
* ARISING OUT OF OR RELATING TO THE SOFTWARE OR THIS AGREEMENT, DAMAGES |
|
* FOR LOSS OF GOODWILL, WORK STOPPAGE, OR LOSS OF DATA, OR FOR ANY |
|
* DAMAGES, EVEN IF FLL SHALL HAVE BEEN INFORMED OF THE POSSIBILITY OF |
|
* SUCH DAMAGES, OR FOR ANY CLAIM BY ANY OTHER PARTY. EVEN IF A PART |
|
* OF THE SOFTWARE HAS BEEN DEVELOPED BY A THIRD PARTY, THE THIRD PARTY |
|
* 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.6 2000/08/21 08:31:18 noro Exp $ |
|
*/ |
#include "ca.h" |
#include "ca.h" |
#include "base.h" |
#include "base.h" |
#include "parse.h" |
#include "parse.h" |
#include "inline.h" |
#include "inline.h" |
/* |
|
|
#if 0 |
#undef DMAR |
#undef DMAR |
#define DMAR(a1,a2,a3,d,r) (r)=dmar(a1,a2,a3,d); |
#define DMAR(a1,a2,a3,d,r) (r)=dmar(a1,a2,a3,d); |
*/ |
#endif |
|
|
extern int Print; /* XXX */ |
extern int Print; /* XXX */ |
|
|
|
|
|
|
asir_assert(ARG0(arg),O_N,"newvect"); |
asir_assert(ARG0(arg),O_N,"newvect"); |
len = QTOS((Q)ARG0(arg)); |
len = QTOS((Q)ARG0(arg)); |
if ( len <= 0 ) |
if ( len < 0 ) |
error("newvect : invalid size"); |
error("newvect : invalid size"); |
MKVECT(vect,len); |
MKVECT(vect,len); |
if ( argc(arg) == 2 ) { |
if ( argc(arg) == 2 ) { |
|
|
asir_assert(ARG0(arg),O_N,"newmat"); |
asir_assert(ARG0(arg),O_N,"newmat"); |
asir_assert(ARG1(arg),O_N,"newmat"); |
asir_assert(ARG1(arg),O_N,"newmat"); |
row = QTOS((Q)ARG0(arg)); col = QTOS((Q)ARG1(arg)); |
row = QTOS((Q)ARG0(arg)); col = QTOS((Q)ARG1(arg)); |
if ( row <= 0 || col <= 0 ) |
if ( row < 0 || col < 0 ) |
error("newmat : invalid size"); |
error("newmat : invalid size"); |
MKMAT(m,row,col); |
MKMAT(m,row,col); |
if ( argc(arg) == 3 ) { |
if ( argc(arg) == 3 ) { |
|
|
return -1; |
return -1; |
} |
} |
|
|
struct oEGT eg_mod,eg_elim,eg_chrem,eg_gschk,eg_intrat,eg_symb; |
struct oEGT eg_mod,eg_elim,eg_elim1,eg_elim2,eg_chrem,eg_gschk,eg_intrat,eg_symb; |
|
|
int generic_gauss_elim(mat,nm,dn,rindp,cindp) |
int generic_gauss_elim(mat,nm,dn,rindp,cindp) |
MAT mat; |
MAT mat; |
|
|
return 1; |
return 1; |
} |
} |
|
|
|
#define ONE_STEP1 if ( zzz = *s ) { DMAR(zzz,hc,*tj,md,*tj) } tj++; s++; |
|
|
|
void reduce_reducers_mod(mat,row,col,md) |
|
int **mat; |
|
int row,col; |
|
int md; |
|
{ |
|
int i,j,k,l,hc,zzz; |
|
int *t,*s,*tj,*ind; |
|
|
|
/* reduce the reducers */ |
|
ind = (int *)ALLOCA(row*sizeof(int)); |
|
for ( i = 0; i < row; i++ ) { |
|
t = mat[i]; |
|
for ( j = 0; j < col && !t[j]; j++ ); |
|
/* register the position of the head term */ |
|
ind[i] = j; |
|
for ( l = i-1; l >= 0; l-- ) { |
|
/* reduce mat[i] by mat[l] */ |
|
if ( hc = t[ind[l]] ) { |
|
/* mat[i] = mat[i]-hc*mat[l] */ |
|
j = ind[l]; |
|
s = mat[l]+j; |
|
tj = t+j; |
|
hc = md-hc; |
|
k = col-j; |
|
for ( ; k >= 64; k -= 64 ) { |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
ONE_STEP1 ONE_STEP1 ONE_STEP1 ONE_STEP1 |
|
} |
|
for ( ; k >= 0; k-- ) { |
|
if ( zzz = *s ) { DMAR(zzz,hc,*tj,md,*tj) } tj++; s++; |
|
} |
|
} |
|
} |
|
} |
|
} |
|
|
|
/* |
|
mat[i] : reducers (i=0,...,nred-1) |
|
spolys (i=nred,...,row-1) |
|
mat[0] < mat[1] < ... < mat[nred-1] w.r.t the term order |
|
1. reduce the reducers |
|
2. reduce spolys by the reduced reducers |
|
*/ |
|
|
|
void pre_reduce_mod(mat,row,col,nred,md) |
|
int **mat; |
|
int row,col,nred; |
|
int md; |
|
{ |
|
int i,j,k,l,hc,inv; |
|
int *t,*s,*tk,*ind; |
|
|
|
#if 1 |
|
/* reduce the reducers */ |
|
ind = (int *)ALLOCA(row*sizeof(int)); |
|
for ( i = 0; i < nred; i++ ) { |
|
/* make mat[i] monic and mat[i] by mat[0],...,mat[i-1] */ |
|
t = mat[i]; |
|
for ( j = 0; j < col && !t[j]; j++ ); |
|
/* register the position of the head term */ |
|
ind[i] = j; |
|
inv = invm(t[j],md); |
|
for ( k = j; k < col; k++ ) |
|
if ( t[k] ) |
|
DMAR(t[k],inv,0,md,t[k]) |
|
for ( l = i-1; l >= 0; l-- ) { |
|
/* reduce mat[i] by mat[l] */ |
|
if ( hc = t[ind[l]] ) { |
|
/* mat[i] = mat[i]-hc*mat[l] */ |
|
for ( k = ind[l], hc = md-hc, s = mat[l]+k, tk = t+k; |
|
k < col; k++, tk++, s++ ) |
|
if ( *s ) |
|
DMAR(*s,hc,*tk,md,*tk) |
|
} |
|
} |
|
} |
|
/* reduce the spolys */ |
|
for ( i = nred; i < row; i++ ) { |
|
t = mat[i]; |
|
for ( l = nred-1; l >= 0; l-- ) { |
|
/* reduce mat[i] by mat[l] */ |
|
if ( hc = t[ind[l]] ) { |
|
/* mat[i] = mat[i]-hc*mat[l] */ |
|
for ( k = ind[l], hc = md-hc, s = mat[l]+k, tk = t+k; |
|
k < col; k++, tk++, s++ ) |
|
if ( *s ) |
|
DMAR(*s,hc,*tk,md,*tk) |
|
} |
|
} |
|
} |
|
#endif |
|
} |
|
/* |
|
mat[i] : reducers (i=0,...,nred-1) |
|
mat[0] < mat[1] < ... < mat[nred-1] w.r.t the term order |
|
*/ |
|
|
|
void reduce_sp_by_red_mod(sp,redmat,ind,nred,col,md) |
|
int *sp,**redmat; |
|
int *ind; |
|
int nred,col; |
|
int md; |
|
{ |
|
int i,j,k,hc,zzz; |
|
int *t,*s,*tj; |
|
|
|
/* reduce the spolys by redmat */ |
|
for ( i = nred-1; i >= 0; i-- ) { |
|
/* reduce sp by redmat[i] */ |
|
if ( hc = sp[ind[i]] ) { |
|
/* sp = sp-hc*redmat[i] */ |
|
j = ind[i]; |
|
hc = md-hc; |
|
s = redmat[i]+j; |
|
tj = sp+j; |
|
for ( k = col-j; k >= 0; k-- ) { |
|
if ( zzz = *s ) { DMAR(zzz,hc,*tj,md,*tj) } tj++; s++; |
|
} |
|
} |
|
} |
|
} |
|
|
|
#define ONE_STEP2 if ( zzz = *pk ) { DMAR(zzz,a,*tk,md,*tk) } pk++; tk++; |
|
|
int generic_gauss_elim_mod(mat,row,col,md,colstat) |
int generic_gauss_elim_mod(mat,row,col,md,colstat) |
int **mat; |
int **mat; |
int row,col,md; |
int row,col,md; |
int *colstat; |
int *colstat; |
{ |
{ |
int i,j,k,l,inv,a,rank; |
int i,j,k,l,inv,a,rank,zzz; |
int *t,*pivot; |
int *t,*pivot,*pk,*tk; |
|
|
for ( rank = 0, j = 0; j < col; j++ ) { |
for ( rank = 0, j = 0; j < col; j++ ) { |
for ( i = rank; i < row && !mat[i][j]; i++ ); |
for ( i = rank; i < row && !mat[i][j]; i++ ); |
|
|
} |
} |
pivot = mat[rank]; |
pivot = mat[rank]; |
inv = invm(pivot[j],md); |
inv = invm(pivot[j],md); |
for ( k = j; k < col; k++ ) |
for ( k = j, pk = pivot+k; k < col; k++, pk++ ) |
if ( pivot[k] ) { |
if ( *pk ) { |
DMAR(pivot[k],inv,0,md,pivot[k]) |
DMAR(*pk,inv,0,md,*pk) |
} |
} |
for ( i = rank+1; i < row; i++ ) { |
for ( i = rank+1; i < row; i++ ) { |
t = mat[i]; |
t = mat[i]; |
if ( a = t[j] ) |
if ( a = t[j] ) { |
for ( k = j, a = md - a; k < col; k++ ) |
a = md - a; pk = pivot+j; tk = t+j; |
if ( pivot[k] ) { |
k = col-j; |
DMAR(pivot[k],a,t[k],md,t[k]) |
for ( ; k >= 64; k -= 64 ) { |
} |
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
} |
|
for ( ; k >= 0; k -- ) { |
|
if ( zzz = *pk ) { DMAR(zzz,a,*tk,md,*tk) } pk++; tk++; |
|
} |
|
} |
} |
} |
rank++; |
rank++; |
} |
} |
|
|
pivot = mat[l]; |
pivot = mat[l]; |
for ( i = 0; i < l; i++ ) { |
for ( i = 0; i < l; i++ ) { |
t = mat[i]; |
t = mat[i]; |
if ( a = t[j] ) |
if ( a = t[j] ) { |
for ( k = j, a = md-a; k < col; k++ ) |
a = md-a; pk = pivot+j; tk = t+j; |
if ( pivot[k] ) { |
k = col-j; |
DMAR(pivot[k],a,t[k],md,t[k]) |
for ( ; k >= 64; k -= 64 ) { |
} |
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
ONE_STEP2 ONE_STEP2 ONE_STEP2 ONE_STEP2 |
|
} |
|
for ( ; k >= 0; k -- ) { |
|
if ( zzz = *pk ) { DMAR(zzz,a,*tk,md,*tk) } pk++; tk++; |
|
} |
|
} |
} |
} |
l--; |
l--; |
} |
} |