File: [local] / OpenXM_contrib2 / asir2000 / engine / gfspn.c (download)
Revision 1.2, Mon Sep 3 01:06:01 2001 UTC (22 years, 9 months ago) by noro
Branch: MAIN
Changes since 1.1: +2 -0
lines
Updated version number.
The previous commit includes getopt() (no arguments).
|
/* OpenXM */
#include "ca.h"
#include "base.h"
UM current_mod_gfspn;
void setmod_gfspn(p)
UM p;
{
current_mod_gfspn = p;
}
void getmod_gfspn(up)
UM *up;
{
*up = current_mod_gfspn;
}
void simpgfspn(n,r)
GFSPN n;
GFSPN *r;
{
UM t,q;
if ( !n )
*r = 0;
else if ( NID(n) != N_GFSPN )
*r = n;
else {
t = UMALLOC(DEG(BDY(n)));
q = W_UMALLOC(DEG(BDY(n)));
cpyum(BDY(n),t);
divsfum(t,current_mod_gfspn,q);
MKGFSPN(t,*r);
}
}
#define NZGFSPN(a) ((a)&&(OID(a)==O_N)&&(NID(a)==N_GFSPN))
void ntogfspn(a,b)
Obj a;
GFSPN *b;
{
UM t;
GFS c;
if ( !a || (OID(a)==O_N && NID(a) == N_GFSPN) )
*b = (GFSPN)a;
else if ( OID(a) == O_N && (NID(a) == N_GFS || NID(a) == N_Q) ) {
ntogfs((Num)a,&c);
if ( !b )
*b = 0;
else {
t = UMALLOC(0); ptosfum((P)c,t);
MKGFSPN(t,*b);
}
} else
error("ntogfspn : invalid argument");
}
void addgfspn(a,b,c)
GFSPN a,b;
GFSPN *c;
{
UM t,q;
GFSPN z;
int d;
ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
if ( !a )
*c = b;
else if ( !b )
*c = a;
else {
d = MAX(DEG(BDY(a)),DEG(BDY(b)));
t = UMALLOC(d);
q = W_UMALLOC(d);
addsfum(BDY(a),BDY(b),t);
divsfum(t,current_mod_gfspn,q);
MKGFSPN(t,*c);
}
}
void subgfspn(a,b,c)
GFSPN a,b;
GFSPN *c;
{
UM t,q;
GFSPN z;
int d;
ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
if ( !a )
chsgngfspn(b,c);
else if ( !b )
*c = a;
else {
d = MAX(DEG(BDY(a)),DEG(BDY(b)));
t = UMALLOC(d);
q = W_UMALLOC(d);
subsfum(BDY(a),BDY(b),t);
divsfum(t,current_mod_gfspn,q);
MKGFSPN(t,*c);
}
}
extern int up_lazy;
void mulgfspn(a,b,c)
GFSPN a,b;
GFSPN *c;
{
UM t,q;
GFSPN z;
int d;
ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
if ( !a || !b )
*c = 0;
else {
d = DEG(BDY(a))+DEG(BDY(b));
t = UMALLOC(d);
q = W_UMALLOC(d);
mulsfum(BDY(a),BDY(b),t);
divsfum(t,current_mod_gfspn,q);
MKGFSPN(t,*c);
}
}
void divgfspn(a,b,c)
GFSPN a,b;
GFSPN *c;
{
GFSPN z;
int d;
UM wb,wc,wd,we,t,q;
ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
if ( !b )
error("divgfspn: division by 0");
else if ( !a )
*c = 0;
else {
wb = W_UMALLOC(DEG(BDY(b))); cpyum(BDY(b),wb);
d = DEG(current_mod_gfspn);
wc = W_UMALLOC(d); cpyum(current_mod_gfspn,wc);
wd = W_UMALLOC(2*d); we = W_UMALLOC(2*d);
/* wd*wb+we*wc=1 */
eucsfum(wb,wc,wd,we);
d = DEG(BDY(a))+DEG(wd);
t = UMALLOC(d);
q = W_UMALLOC(d);
mulsfum(a,wd,t);
divsfum(t,current_mod_gfspn,q);
MKGFSPN(t,*c);
}
}
void invgfspn(b,c)
GFSPN b;
GFSPN *c;
{
GFSPN z;
int d;
UM wb,wc,wd,we,t;
ntogfspn((Obj)b,&z); b = z;
if ( !b )
error("divgfspn: division by 0");
else {
wb = W_UMALLOC(DEG(BDY(b))); cpyum(BDY(b),wb);
d = DEG(current_mod_gfspn);
wc = W_UMALLOC(d); cpyum(current_mod_gfspn,wc);
wd = W_UMALLOC(2*d); we = W_UMALLOC(2*d);
/* wd*wb+we*wc=1 */
eucsfum(wb,wc,wd,we);
d = DEG(wd);
t = UMALLOC(d);
cpyum(wd,t);
MKGFSPN(t,*c);
}
}
void chsgngfspn(a,c)
GFSPN a,*c;
{
GFSPN z;
int d;
struct oUM zero;
UM t;
ntogfspn((Obj)a,&z); a = z;
if ( !a )
*c = 0;
else {
d = DEG(BDY(a));
t = UMALLOC(d);
DEG(&zero) = -1;
subsfum(&zero,BDY(a),t);
MKGFSPN(t,*c);
}
}
void pwrgfspn(a,b,c)
GFSPN a;
Q b;
GFSPN *c;
{
GFSPN z;
UM t,x,y,q;
int d,k;
N e;
ntogfspn((Obj)a,&z); a = z;
if ( !b ) {
t = UMALLOC(0); DEG(t) = 0; COEF(t)[0] = _onesf();
MKGFSPN(t,*c);
} else if ( !a )
*c = 0;
else {
d = DEG(current_mod_gfspn);
/* y = 1 */
y = UMALLOC(d); DEG(y) = 0; COEF(y)[0] = _onesf();
t = W_UMALLOC(2*d); q = W_UMALLOC(2*d);
/* x = simplify(a) */
x = W_UMALLOC(DEG(BDY(a))); cpyum(BDY(a),x);
divsfum(x,current_mod_gfspn,q);
if ( DEG(x) < 0 ) {
*c = 0;
} else {
e = NM(b);
for ( k = n_bits(e)-1; k >= 0; k-- ) {
mulsfum(y,y,t);
divsfum(t,current_mod_gfspn,q);
cpyum(t,y);
if ( e->b[k/32] & (1<<(k%32)) ) {
mulsfum(y,x,t);
divsfum(t,current_mod_gfspn,q);
cpyum(t,y);
}
}
MKGFSPN(y,*c);
}
}
}
int cmpgfspn(a,b)
GFSPN a,b;
{
GFSPN z;
ntogfspn((Obj)a,&z); a = z; ntogfspn((Obj)b,&z); b = z;
if ( !a )
if ( !b )
return 0;
else
return -1;
else if ( !b )
return 1;
else
return compsfum(BDY(a),BDY(b));
}
void randomgfspn(r)
GFSPN *r;
{
int i,d;
UM t;
if ( !current_mod_gfspn )
error("randomgfspn : current_mod_gfspn is not set");
d = DEG(current_mod_gfspn);
t = UMALLOC(d-1);
randsfum(d,t);
degum(t,d-1);
if ( DEG(t) < 0 )
*r = 0;
else {
MKGFSPN(t,*r);
}
}