Annotation of OpenXM_contrib2/asir2000/builtin/pf.c, Revision 1.1
1.1 ! noro 1: /* $OpenXM: OpenXM/src/asir99/builtin/pf.c,v 1.1.1.1 1999/11/10 08:12:26 noro Exp $ */
! 2: #include "ca.h"
! 3: #include "math.h"
! 4: #include "parse.h"
! 5: #if 0
! 6: #include <alloca.h>
! 7: #endif
! 8:
! 9: double const_pi(),const_e();
! 10:
! 11: void make_ihyp(void);
! 12: void make_hyp(void);
! 13: void make_itri(void);
! 14: void make_tri(void);
! 15: void make_exp(void);
! 16: void simplify_pow(PFINS,Obj *);
! 17:
! 18: void Pfunctor(),Pargs(),Pfunargs(),Pvtype(),Pcall(),Pdeval();
! 19: void Pregister_handler();
! 20:
! 21: struct ftab puref_tab[] = {
! 22: {"functor",Pfunctor,1},
! 23: {"args",Pargs,1},
! 24: {"funargs",Pfunargs,1},
! 25: {"register_handler",Pregister_handler,1},
! 26: {"call",Pcall,2},
! 27: {"vtype",Pvtype,1},
! 28: {"deval",Pdeval,1},
! 29: {0,0,0},
! 30: };
! 31:
! 32: #if PARI
! 33: int p_pi(),p_e();
! 34: int p_log(),p_exp(),p_pow();
! 35: int p_sin(),p_cos(),p_tan(),p_asin(),p_acos(),p_atan();
! 36: int p_sinh(),p_cosh(),p_tanh(),p_asinh(),p_acosh(),p_atanh();
! 37: #else
! 38: int p_pi,p_e;
! 39: int p_log,p_exp,p_pow;
! 40: int p_sin,p_cos,p_tan,p_asin,p_acos,p_atan;
! 41: int p_sinh,p_cosh,p_tanh,p_asinh,p_acosh,p_atanh;
! 42: #endif
! 43:
! 44: static V *uarg,*darg;
! 45: static P x,y;
! 46: static PF pidef,edef;
! 47: static PF logdef,expdef,powdef;
! 48: static PF sindef,cosdef,tandef;
! 49: static PF asindef,acosdef,atandef;
! 50: static PF sinhdef,coshdef,tanhdef;
! 51: static PF asinhdef,acoshdef,atanhdef;
! 52:
! 53: #define OALLOC(p,n) ((p)=(Obj *)CALLOC((n),sizeof(Obj)))
! 54:
! 55: double const_pi() { return 3.14159265358979323846264338327950288; }
! 56: double const_e() { return 2.718281828459045235360287471352662497; }
! 57:
! 58: void pf_init() {
! 59: uarg = (V *)CALLOC(1,sizeof(V));
! 60: uarg[0] = &oVAR[26]; MKV(uarg[0],x);
! 61:
! 62: darg = (V *)CALLOC(2,sizeof(V));
! 63: darg[0] = &oVAR[26];
! 64: darg[1] = &oVAR[27]; MKV(darg[1],y);
! 65:
! 66: mkpf("@pi",0,0,0,(int (*)())p_pi,const_pi,0,&pidef);
! 67: mkpf("@e",0,0,0,(int (*)())p_e,const_e,0,&edef);
! 68:
! 69: mkpf("log",0,1,uarg,(int (*)())p_log,log,0,&logdef);
! 70: mkpf("exp",0,1,uarg,(int (*)())p_exp,exp,0,&expdef);
! 71: mkpf("pow",0,2,darg,(int (*)())p_pow,pow,(int (*)())simplify_pow,&powdef);
! 72:
! 73: mkpf("sin",0,1,uarg,(int (*)())p_sin,sin,0,&sindef);
! 74: mkpf("cos",0,1,uarg,(int (*)())p_cos,cos,0,&cosdef);
! 75: mkpf("tan",0,1,uarg,(int (*)())p_tan,tan,0,&tandef);
! 76: mkpf("asin",0,1,uarg,(int (*)())p_asin,asin,0,&asindef);
! 77: mkpf("acos",0,1,uarg,(int (*)())p_acos,acos,0,&acosdef);
! 78: mkpf("atan",0,1,uarg,(int (*)())p_atan,atan,0,&atandef);
! 79:
! 80: mkpf("sinh",0,1,uarg,(int (*)())p_sinh,sinh,0,&sinhdef);
! 81: mkpf("cosh",0,1,uarg,(int (*)())p_cosh,cosh,0,&coshdef);
! 82: mkpf("tanh",0,1,uarg,(int (*)())p_tanh,tanh,0,&tanhdef);
! 83: #if !defined(VISUAL)
! 84: mkpf("asinh",0,1,uarg,(int (*)())p_asinh,asinh,0,&asinhdef);
! 85: mkpf("acosh",0,1,uarg,(int (*)())p_acosh,acosh,0,&acoshdef);
! 86: mkpf("atanh",0,1,uarg,(int (*)())p_atanh,atanh,0,&atanhdef);
! 87: #endif
! 88: make_exp();
! 89: make_tri();
! 90: make_itri();
! 91: make_hyp();
! 92: #if !defined(VISUAL)
! 93: make_ihyp();
! 94: #endif
! 95: }
! 96:
! 97: void make_exp() {
! 98: V v;
! 99: P u,vexp,vlog,vpow;
! 100: Obj *args;
! 101:
! 102: mkpfins(expdef,uarg,&v); MKV(v,vexp);
! 103: mkpfins(powdef,darg,&v); MKV(v,vpow);
! 104: mkpfins(logdef,uarg,&v); MKV(v,vlog);
! 105:
! 106: /* d/dx(log(x)) = 1/x */
! 107: OALLOC(logdef->deriv,1); divr(CO,(Obj)ONE,(Obj)x,&logdef->deriv[0]);
! 108:
! 109: /* d/dx(exp(x)) = exp(x) */
! 110: OALLOC(expdef->deriv,1); expdef->deriv[0] = (Obj)vexp;
! 111:
! 112: /* d/dy(x^y) = log(x)*x^y */
! 113: OALLOC(powdef->deriv,2); mulp(CO,vpow,vlog,(P *)&powdef->deriv[1]);
! 114:
! 115: /* d/dx(x^y) = y*x^(y-1) */
! 116: args = (Obj *)ALLOCA(2*sizeof(Obj));
! 117: args[0] = (Obj)x; subp(CO,y,(P)ONE,(P *)&args[1]);
! 118: _mkpfins(powdef,args,&v); MKV(v,u);
! 119: mulr(CO,(Obj)u,(Obj)y,&powdef->deriv[0]);
! 120: }
! 121:
! 122: void make_tri() {
! 123: V v;
! 124: P vcos,vsin,vtan,t;
! 125:
! 126: mkpfins(cosdef,uarg,&v); MKV(v,vcos);
! 127: mkpfins(sindef,uarg,&v); MKV(v,vsin);
! 128: mkpfins(tandef,uarg,&v); MKV(v,vtan);
! 129:
! 130: /* d/dx(sin(x)) = cos(x) */
! 131: OALLOC(sindef->deriv,1); sindef->deriv[0] = (Obj)vcos;
! 132:
! 133: /* d/dx(cos(x)) = -sin(x) */
! 134: OALLOC(cosdef->deriv,1); chsgnp(vsin,(P *)&cosdef->deriv[0]);
! 135:
! 136: /* d/dx(tan(x)) = 1+tan(x)^2 */
! 137: OALLOC(tandef->deriv,1);
! 138: mulr(CO,(Obj)vtan,(Obj)vtan,(Obj *)&t); addp(CO,(P)ONE,t,(P *)&tandef->deriv[0]);
! 139: }
! 140:
! 141: void make_itri() {
! 142: P t,xx;
! 143: Q mtwo;
! 144: V v;
! 145: Obj *args;
! 146:
! 147: /* d/dx(asin(x)) = (1-x^2)^(-1/2) */
! 148: OALLOC(asindef->deriv,1);
! 149: args = (Obj *)ALLOCA(2*sizeof(Obj));
! 150: mulp(CO,x,x,&xx); subp(CO,(P)ONE,xx,(P *)&args[0]);
! 151: STOQ(-2,mtwo); divq(ONE,mtwo,(Q *)&args[1]);
! 152: _mkpfins(powdef,args,&v); MKV(v,t);
! 153: asindef->deriv[0] = (Obj)t;
! 154:
! 155: /* d/dx(acos(x)) = -(1-x^2)^(-1/2) */
! 156: OALLOC(acosdef->deriv,1); chsgnp((P)asindef->deriv[0],(P *)&acosdef->deriv[0]);
! 157:
! 158: /* d/dx(atan(x)) = 1/(x^2+1) */
! 159: OALLOC(atandef->deriv,1);
! 160: addp(CO,(P)ONE,xx,&t); divr(CO,(Obj)ONE,(Obj)t,&atandef->deriv[0]);
! 161: }
! 162:
! 163: void make_hyp() {
! 164: V v;
! 165: P vcosh,vsinh,vtanh,t;
! 166:
! 167: mkpfins(coshdef,uarg,&v); MKV(v,vcosh);
! 168: mkpfins(sinhdef,uarg,&v); MKV(v,vsinh);
! 169: mkpfins(tanhdef,uarg,&v); MKV(v,vtanh);
! 170:
! 171: /* d/dx(sinh(x)) = cosh(x) */
! 172: OALLOC(sinhdef->deriv,1); sinhdef->deriv[0] = (Obj)vcosh;
! 173:
! 174: /* d/dx(cosh(x)) = sinh(x) */
! 175: OALLOC(coshdef->deriv,1); coshdef->deriv[0] = (Obj)vsinh;
! 176:
! 177: /* d/dx(tanh(x)) = 1-tanh(x)^2 */
! 178: OALLOC(tanhdef->deriv,1);
! 179: mulr(CO,(Obj)vtanh,(Obj)vtanh,(Obj *)&t); subp(CO,(P)ONE,t,(P *)&tanhdef->deriv[0]);
! 180: }
! 181:
! 182: void make_ihyp() {
! 183: P t,xx;
! 184: Q mtwo;
! 185: V v;
! 186: Obj *args;
! 187:
! 188: /* d/dx(asinh(x)) = (1+x^2)^(-1/2) */
! 189: OALLOC(asinhdef->deriv,1);
! 190: args = (Obj *)ALLOCA(2*sizeof(Obj));
! 191: mulp(CO,x,x,&xx); addp(CO,(P)ONE,xx,(P *)&args[0]);
! 192: STOQ(-2,mtwo); divq(ONE,mtwo,(Q *)&args[1]);
! 193: _mkpfins(powdef,args,&v); MKV(v,t);
! 194: asinhdef->deriv[0] = (Obj)t;
! 195:
! 196: /* d/dx(acosh(x)) = (x^2-1)^(-1/2) */
! 197: OALLOC(acoshdef->deriv,1);
! 198: subp(CO,xx,(P)ONE,(P *)&args[0]);
! 199: _mkpfins(powdef,args,&v); MKV(v,t);
! 200: acoshdef->deriv[0] = (Obj)t;
! 201:
! 202: /* d/dx(atanh(x)) = 1/(1-x^2) */
! 203: OALLOC(atanhdef->deriv,1);
! 204: subp(CO,(P)ONE,xx,&t); divr(CO,(Obj)ONE,(Obj)t,&atanhdef->deriv[0]);
! 205: }
! 206:
! 207: void mkpow(vl,a,e,r)
! 208: VL vl;
! 209: Obj a;
! 210: Obj e;
! 211: Obj *r;
! 212: {
! 213: PFINS ins;
! 214: PFAD ad;
! 215:
! 216: ins = (PFINS)CALLOC(1,sizeof(PF)+2*sizeof(struct oPFAD));
! 217: ins->pf = powdef; ad = ins->ad;
! 218: ad[0].d = 0; ad[0].arg = a; ad[1].d = 0; ad[1].arg = e;
! 219: simplify_ins(ins,r);
! 220: }
! 221:
! 222: void simplify_pow(ins,rp)
! 223: PFINS ins;
! 224: Obj *rp;
! 225: {
! 226: PF pf;
! 227: PFAD ad;
! 228: Obj a0,a1;
! 229: V v;
! 230: P t;
! 231:
! 232: pf = ins->pf; ad = ins->ad; a0 = ad[0].arg; a1 = ad[1].arg;
! 233: if ( !a1 )
! 234: *rp = (Obj)ONE;
! 235: else if ( !a0 )
! 236: *rp = 0;
! 237: else if ( NUM(a1) && INT(a1) )
! 238: arf_pwr(CO,a0,a1,rp);
! 239: else {
! 240: instov(ins,&v); MKV(v,t); *rp = (Obj)t;
! 241: }
! 242: }
! 243:
! 244: #define ISPFINS(p)\
! 245: (p)&&(ID(p) == O_P)&&((int)VR((P)p)->attr!=V_PF)&&\
! 246: UNIQ(DEG(DC((P)p)))&&UNIQ(COEF(DC((P)p)))
! 247:
! 248: void Pfunctor(arg,rp)
! 249: NODE arg;
! 250: P *rp;
! 251: {
! 252: P p;
! 253: FUNC t;
! 254: PF pf;
! 255: PFINS ins;
! 256:
! 257: p = (P)ARG0(arg);
! 258: if ( !ISPFINS(p) )
! 259: *rp = 0;
! 260: else {
! 261: ins = (PFINS)VR(p)->priv; pf = ins->pf;
! 262: t = (FUNC)MALLOC(sizeof(struct oFUNC));
! 263: t->name = pf->name; t->id = A_PURE; t->argc = pf->argc;
! 264: t->f.puref = pf;
! 265: makesrvar(t,rp);
! 266: }
! 267: }
! 268:
! 269: void Pargs(arg,rp)
! 270: NODE arg;
! 271: LIST *rp;
! 272: {
! 273: P p;
! 274: PF pf;
! 275: PFAD ad;
! 276: PFINS ins;
! 277: NODE n,n0;
! 278: int i;
! 279:
! 280: p = (P)ARG0(arg);
! 281: if ( !ISPFINS(p) )
! 282: *rp = 0;
! 283: else {
! 284: ins = (PFINS)VR(p)->priv; ad = ins->ad; pf = ins->pf;
! 285: for ( i = 0, n0 = 0; i < pf->argc; i++ ) {
! 286: NEXTNODE(n0,n); BDY(n) = (pointer)ad[i].arg;
! 287: }
! 288: if ( n0 )
! 289: NEXT(n) = 0;
! 290: MKLIST(*rp,n0);
! 291: }
! 292: }
! 293:
! 294: void Pfunargs(arg,rp)
! 295: NODE arg;
! 296: LIST *rp;
! 297: {
! 298: P p;
! 299: P f;
! 300: FUNC t;
! 301: PF pf;
! 302: PFINS ins;
! 303: PFAD ad;
! 304: NODE n,n0;
! 305: int i;
! 306:
! 307: p = (P)ARG0(arg);
! 308: if ( !ISPFINS(p) )
! 309: *rp = 0;
! 310: else {
! 311: ins = (PFINS)VR(p)->priv; ad = ins->ad; pf = ins->pf;
! 312: t = (FUNC)MALLOC(sizeof(struct oFUNC));
! 313: t->name = pf->name; t->id = A_PURE; t->argc = pf->argc;
! 314: t->f.puref = pf;
! 315: makesrvar(t,&f);
! 316: n0 = 0; NEXTNODE(n0,n); BDY(n) = (pointer)f;
! 317: for ( i = 0; i < pf->argc; i++ ) {
! 318: NEXTNODE(n0,n); BDY(n) = (pointer)ad[i].arg;
! 319: }
! 320: NEXT(n) = 0;
! 321: MKLIST(*rp,n0);
! 322: }
! 323: }
! 324:
! 325: void Pvtype(arg,rp)
! 326: NODE arg;
! 327: Q *rp;
! 328: {
! 329: P p;
! 330:
! 331: p = (P)ARG0(arg);
! 332: if ( !p || ID(p) != O_P )
! 333: *rp = 0;
! 334: else
! 335: STOQ((int)VR(p)->attr,*rp);
! 336: }
! 337:
! 338: extern FUNC registered_handler;
! 339:
! 340: void Pregister_handler(arg,rp)
! 341: NODE arg;
! 342: Q *rp;
! 343: {
! 344: P p;
! 345: V v;
! 346: FUNC func;
! 347:
! 348: p = (P)ARG0(arg);
! 349: if ( !p )
! 350: registered_handler = 0;
! 351: else if ( OID(p) != 2 )
! 352: error("register_hanlder : invalid argument");
! 353: v = VR(p);
! 354: if ( (int)v->attr != V_SR )
! 355: error("register_hanlder : no such function");
! 356: else {
! 357: func = (FUNC)v->priv;
! 358: if ( func->argc )
! 359: error("register_hanlder : the function must be with no argument");
! 360: else {
! 361: registered_handler = func;
! 362: *rp = ONE;
! 363: }
! 364: }
! 365: }
! 366:
! 367: void Pcall(arg,rp)
! 368: NODE arg;
! 369: Obj *rp;
! 370: {
! 371: P p;
! 372: V v;
! 373:
! 374: p = (P)ARG0(arg);
! 375: if ( !p || OID(p) != 2 )
! 376: error("call : invalid argument");
! 377: v = VR(p);
! 378: if ( (int)v->attr != V_SR )
! 379: error("call : no such function");
! 380:
! 381: else
! 382: *rp = (Obj)bevalf((FUNC)v->priv,BDY((LIST)ARG1(arg)));
! 383: }
! 384:
! 385: void Pdeval(arg,rp)
! 386: NODE arg;
! 387: Obj *rp;
! 388: {
! 389: asir_assert(ARG0(arg),O_R,"deval");
! 390: devalr(CO,(Obj)ARG0(arg),rp);
! 391: }
! 392:
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>