Return to B1fun.c CVS log

Up to [CSRG BSD Unix] / 43BSD / contrib / B / src / bsmall

Annotation of 43BSD/contrib/B/src/bsmall/B1fun.c, revision 1.1

1.1     ! root        1: /* Copyright (c) Stichting Mathematisch Centrum, Amsterdam, 1984. */
        !             2: /* $Header: B1fun.c,v 1.1 84/06/28 00:48:54 timo Exp $ */
        !             3: 
        !             4: /* B functions */
        !             5: #include "b.h"
        !             6: #include "b1obj.h"
        !             7: #include "b2sem.h"
        !             8: #include "B1num.h"
        !             9: 
        !            10: #define Maxlen(x, y) (Length(x) > Length(y) ? Length(x) : Length(y))
        !            11: #define Sumlen(x, y) (Length(x) + Length(y))
        !            12: 
        !            13: value sum(x, y) register value x, y; {
        !            14:        value r, z;
        !            15:        Checknum(x); Checknum(y);
        !            16:        if (!Exact(x) || !Exact(y)) return mk_approx(Numval(x) + Numval(y));
        !            17:        z= mk_exact(Denominator(x), Denominator(y), 0);
        !            18:        r= mk_exact(Numerator(x)*Denominator(z)+Numerator(y)*Numerator(z),
        !            19:                        Denominator(x)*Denominator(z), Maxlen(x, y));
        !            20:        release(z);
        !            21:        return r;
        !            22: }
        !            23: 
        !            24: value negated(x) register value x; {
        !            25:        Checknum(x);
        !            26:        if (!Exact(x)) return mk_approx(-Numval(x));
        !            27:        return mk_exact(-Numerator(x), Denominator(x), Length(x));
        !            28: }
        !            29: 
        !            30: value diff(x, y) register value x, y; {
        !            31:        value r, my;
        !            32:        r= sum(x, my= negated(y));
        !            33:        release(my);
        !            34:        return r;
        !            35: }
        !            36: 
        !            37: value inverse(x) value x;{
        !            38:        Checknum(x);
        !            39:        if (Numval(x) == 0) error("in x/y, y is zero");
        !            40:        if (!Exact(x)) return mk_approx(1.0/Numval(x));
        !            41:        return mk_exact(Denominator(x), Numerator(x), Length(x));
        !            42: }
        !            43: 
        !            44: value prod(x, y) register value x, y; {
        !            45:        value a, b, r;
        !            46:        Checknum(x); Checknum(y);
        !            47:        if (!Exact(x) || !Exact(y)) return mk_approx(Numval(x) * Numval(y));
        !            48:        a= mk_exact(Numerator(x), Denominator(y), 0);
        !            49:        b= mk_exact(Numerator(y), Denominator(x), 0);
        !            50:        r= mk_exact(Numerator(a)*Numerator(b), Denominator(a)*Denominator(b),
        !            51:                                                Sumlen(x, y));
        !            52:        release(a); release(b);
        !            53:        return r;
        !            54: }
        !            55: 
        !            56: value quot(x, y) register value x, y; {
        !            57:        value r, iy;
        !            58:        r= prod(x, iy= inverse(y));
        !            59:        release(iy);
        !            60:        return r;
        !            61: }
        !            62: 
        !            63: #define Even(x) ((x) == Two*floor((x)/Two))
        !            64: value power(x, y) register value x, y; {
        !            65:        Checknum(x); Checknum(y);
        !            66:        if (Exact(y)) {
        !            67:                integer py= Numerator(y), qy= Denominator(y);
        !            68:                if (Integral(y) && Exact(x)) {
        !            69:                        integer px, qx, ppx, pqx, Ppx, Pqx;
        !            70:                        if (py == Zero) return mk_int(One);
        !            71:                        if (py > Zero) {
        !            72:                                px= Numerator(x);
        !            73:                                qx= Denominator(x);
        !            74:                        } else {
        !            75:                                py= -py;
        !            76:                                px= Denominator(x);
        !            77:                                qx= Numerator(x);
        !            78:                        }
        !            79:                        ppx= pqx= One;
        !            80:                        Ppx= px; Pqx= qx;
        !            81:                        while (py >= Two) {
        !            82:                                if (!Even(py)) {
        !            83:                                        ppx*= Ppx; pqx*= Pqx;
        !            84:                                }
        !            85:                                Ppx*= Ppx; Pqx*= Pqx;
        !            86:                                py= floor(py/Two);
        !            87:                        }
        !            88:                        ppx*= Ppx; pqx*= Pqx;
        !            89:                        return mk_exact(ppx, pqx, 0);
        !            90:                } /* END Integral(y) && Exact(x) */
        !            91:                else {
        !            92:                        double vx= Numval(x);
        !            93:                        short sx= vx < 0 ? -1 : vx == 0 ? 0 : 1;
        !            94:                        if (sx < 0 && Even(qy))
        !            95:                            error("in x**(p/q), x is negative and q is even");
        !            96:                        if (sx == 0 && py < Zero)
        !            97:                            error("0**y with negative y");
        !            98:                        if (sx < 0 && Even(py)) sx= 1;
        !            99:                        return mk_approx(sx * pow(fabs(vx), py/qy));
        !           100:                }
        !           101:        } /* END Exact(y) */
        !           102:        else {
        !           103:                double vx= Numval(x), vy= Approxval(y);
        !           104:                if (vy == 0) return mk_approx(1.0);
        !           105:                if (vx < 0)
        !           106:                        error("in x**y, x is negative and y is not exact");
        !           107:                if (vx == 0 && vy < 0)
        !           108:                        error("0E0**y with negative y");
        !           109:                return mk_approx(pow(vx, vy));
        !           110:        }
        !           111: }
        !           112: 
        !           113: value root2(n, x) register value n, x; {
        !           114:        value r, in;
        !           115:        Checknum(n);
        !           116:        if (Numval(n) == 0) error("in x root y, x is zero");
        !           117:        r= power(x, in= inverse(n));
        !           118:        release(in);
        !           119:        return r;
        !           120: }
        !           121: 
        !           122: value absval(x) register value x; {
        !           123:        Checknum(x);
        !           124:        if (!Exact(x)) return mk_approx(fabs(Numval(x)));
        !           125:        return mk_exact((integer) fabs((double) Numerator(x)), Denominator(x), Length(x));
        !           126: }
        !           127: 
        !           128: value signum(x) register value x; {
        !           129:        double v= numval(x);
        !           130:        return mk_int(v < 0 ? -One : v == 0 ? Zero : One);
        !           131: }
        !           132: 
        !           133: value floorf(x) register value x; {
        !           134:        return mk_int(floor(numval(x)));
        !           135: }
        !           136: 
        !           137: value ceilf(x) register value x; {
        !           138:        return mk_int(ceil(numval(x)));
        !           139: }
        !           140: 
        !           141: value round1(x) register value x; {
        !           142:        return mk_int(floor(numval(x) + .5));
        !           143: }
        !           144: 
        !           145: value round2(n, x) register value n, x; {
        !           146:        value ten, tenp, xtenp, r0, r;
        !           147:        Checknum(n);
        !           148:        if (!Integral(n)) error("in n round x, n is not an integer");
        !           149:        ten= mk_integer(10);
        !           150:        tenp= power(ten, n);
        !           151:        xtenp= prod(x, tenp);
        !           152:        r0= round1(xtenp);
        !           153:        r= mk_exact(Numerator(r0), Numerator(tenp), propintlet((int) Numerator(n)));
        !           154:        release(ten); release(tenp); release(xtenp); release(r0);
        !           155:        return r;
        !           156: }
        !           157: 
        !           158: value mod(a, n) register value a, n; {
        !           159:        value f, p, d;
        !           160:        Checknum(a); Checknum(n);
        !           161:        f= mk_int(floor(Numval(a) / Numval(n)));
        !           162:        p= prod(n, f);
        !           163:        d= diff(a, p);
        !           164:        release(f); release(p);
        !           165:        return d;
        !           166: }
        !           167: 
        !           168: double lastran;
        !           169: 
        !           170: setran (seed) double seed;
        !           171: {double x;
        !           172:  x= seed >= 0 ? seed : -seed;
        !           173:  while (x >= 1) x/= 10;
        !           174:  lastran= floor(67108864.0 * x);
        !           175: }
        !           176: 
        !           177: set_random(v) value v; {
        !           178:        setran((double) hash(v));
        !           179: }
        !           180: 
        !           181: value random() /* 0 <= r < 1 */
        !           182: {double p;
        !           183:  p= 26353589.0 * lastran + 1;
        !           184:  lastran= p - 67108864.0 * floor (p / 67108864.0);
        !           185:  return mk_approx(lastran / 67108864.0);
        !           186: }
        !           187: 
        !           188: value root1(v) value v; {
        !           189:        value two= mk_integer(2);
        !           190:        v= root2(two, v);
        !           191:        release(two);
        !           192:        return(v);
        !           193: }
        !           194: 
        !           195: value pi() { return mk_approx(3.141592653589793238462); }
        !           196: value e() { return mk_approx(exp(1.0)); }
        !           197: 
        !           198: value sin1(v) value v; { return mk_approx(sin(numval(v))); }
        !           199: value cos1(v) value v; { return mk_approx(cos(numval(v))); }
        !           200: value tan1(v) value v; { return mk_approx(tan(numval(v))); }
        !           201: value atn1(v) value v; { return mk_approx(atan(numval(v))); }
        !           202: value exp1(v) value v; { return mk_approx(exp(numval(v))); }
        !           203: value log1(v) value v; { return mk_approx(log(numval(v))); }
        !           204: 
        !           205: value log2(u, v) value u, v;{
        !           206:        return mk_approx(log(numval(v)) / log(numval(u)));
        !           207: }
        !           208: 
        !           209: value atn2(u, v) value u, v; {
        !           210:        return mk_approx(atan2(numval(v), numval(u)));
        !           211: }

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.