Return to b1num.c CVS log

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

Annotation of 43BSD/contrib/B/src/bint/b1num.c, revision 1.1

1.1     ! root        1: /* Copyright (c) Stichting Mathematisch Centrum, Amsterdam, 1985. */
        !             2: 
        !             3: /*
        !             4:   $Header: b1num.c,v 1.4 85/08/22 16:51:59 timo Exp $
        !             5: */
        !             6: 
        !             7: /* B numbers, basic external interface */
        !             8: 
        !             9: #include "b.h"
        !            10: #include "b0con.h"
        !            11: #include "b1obj.h"
        !            12: #include "b1num.h"
        !            13: 
        !            14: /*
        !            15:  * This file contains operations on numbers that are not predefined
        !            16:  * B functions (the latter are defined in `bfun.c').
        !            17:  * This includes conversion of numeric B values to C `int's and
        !            18:  * `double's, (numval() and intval()),
        !            19:  * but also utilities for comparing numeric values and hashing a numeric
        !            20:  * value to something usable as initialization for the random generator
        !            21:  * without chance of overflow (so numval(v) is unusable).
        !            22:  * It is also possible to build numbers of all types using mk_integer,
        !            23:  * mk_exact (or mk_rat) and mk_approx.  Note the rather irregular
        !            24:  * (historical!) argument structure for these: mk_approx has a
        !            25:  * `double' argument, where mk_exact and mk_rat have two values
        !            26:  * which must be `integer' (not `int's!) as arguments.
        !            27:  * The random number generator used by the DRAW and CHOOSE statements
        !            28:  * is also in this file.
        !            29:  */
        !            30: 
        !            31: /*
        !            32:  * ival is used internally by intval() and large().
        !            33:  * It converts an integer to double precision, yielding
        !            34:  * a huge value when overflow occurs (but giving no error).
        !            35:  */
        !            36: 
        !            37: Hidden double ival(u) integer u; {
        !            38:        double x = 0;
        !            39:        register int i;
        !            40: 
        !            41:        if (IsSmallInt(u)) return SmallIntVal(u);
        !            42:        for (i = Length(u)-1; i >= 0; --i) {
        !            43:                if (x >= Maxreal/BASE)
        !            44:                        return Msd(u) < 0 ? -Maxreal : Maxreal;
        !            45:                x = x*BASE + Digit(u, i);
        !            46:        }
        !            47: 
        !            48:        return x;
        !            49: }
        !            50: 
        !            51: 
        !            52: /* Determine if a value may be converted to an int */
        !            53: 
        !            54: Visible bool large(v) value v; {
        !            55:        double r;
        !            56:        if (!Is_number(v) || !integral(v)) {
        !            57:                error(MESS(1300, "number not an integer"));
        !            58:                return No;
        !            59:        }
        !            60:        if (Rational(v)) v= (value) Numerator((rational)v);
        !            61:        r= ival((integer)v);
        !            62:        if (r > Maxint) return Yes;
        !            63:        return No;
        !            64: }
        !            65: 
        !            66: /* return the C `int' value of a B numeric value, if it exists. */
        !            67: 
        !            68: Visible int intval(v) value v; {
        !            69:        /* v must be an Integral number or a Rational with Denominator==1
        !            70:            which may result from n round x [via mk_exact]!. */
        !            71:        double i;
        !            72:        if (IsSmallInt(v)) i= SmallIntVal(v);
        !            73:        else {
        !            74:                if (!Is_number(v)) syserr(MESS(1301, "intval on non-number"));
        !            75:                if (!integral(v)) {
        !            76:                        error(MESS(1302, "number not an integer"));
        !            77:                        return 0;
        !            78:                }
        !            79:                if (Rational(v)) v= (value) Numerator((rational)v);
        !            80:                i= ival((integer)v);
        !            81:        }
        !            82:        if (i > Maxint || i < -Maxint) {
        !            83:                error(MESS(1303, "exceedingly large integer"));
        !            84:                return 0;
        !            85:        }
        !            86:        return (int) i;
        !            87: }
        !            88: 
        !            89: 
        !            90: /* convert an int to an intlet */
        !            91: 
        !            92: Visible intlet propintlet(i) int i; {
        !            93:        if (i > Maxintlet || i < -Maxintlet) {
        !            94:                error(MESS(1304, "exceedingly large integer"));
        !            95:                return 0;
        !            96:        }
        !            97:        return i;
        !            98: }
        !            99: 
        !           100: 
        !           101: /*
        !           102:  * determine if a number is integer
        !           103:  */
        !           104: 
        !           105: Visible bool integral(v) value v; {
        !           106:        if (Integral(v) || (Rational(v) && Denominator((rational)v) == int_1))
        !           107:                return Yes;
        !           108:        else return No;
        !           109: }
        !           110: 
        !           111: /*
        !           112:  * mk_exact makes an exact number out of two integers.
        !           113:  * The third argument is the desired number of digits after the decimal
        !           114:  * point when the number is printed (see comments in `bfun.c' for
        !           115:  * `round' and code in `bnuC.c').
        !           116:  * This printing size (if positive) is hidden in the otherwise nearly
        !           117:  * unused * 'Size' field of the value struct
        !           118:  * (cf. definition of Rational(v) etc.).
        !           119:  */
        !           120: 
        !           121: Visible value mk_exact(p, q, len) integer p, q; int len; {
        !           122:        rational r = mk_rat(p, q, len);
        !           123: 
        !           124:        if (Denominator(r) == int_1 && len <= 0) {
        !           125:                integer i = (integer) Copy(Numerator(r));
        !           126:                release((value) r);
        !           127:                return (value) i;
        !           128:        }
        !           129: 
        !           130:        return (value) r;
        !           131: }
        !           132: 
        !           133: #define uSMALL 1
        !           134: #define uINT 2
        !           135: #define uRAT 4
        !           136: #define uAPP 8
        !           137: #define vSMALL 16
        !           138: #define vINT 32
        !           139: #define vRAT 64
        !           140: #define vAPP 128
        !           141: 
        !           142: Visible relation numcomp(u, v) value u, v; {
        !           143:        int tu, tv; relation s;
        !           144: 
        !           145:        if (IsSmallInt(u)) tu = uSMALL;
        !           146:        else if (Length(u) >= 0) tu = uINT;
        !           147:        else if (Length(u) <= -2) tu = uRAT;
        !           148:        else tu = uAPP;
        !           149:        if (IsSmallInt(v)) tv = vSMALL;
        !           150:        else if (Length(v) >= 0) tv = vINT;
        !           151:        else if (Length(v) <= -2) tv = vRAT;
        !           152:        else tv = vAPP;
        !           153: 
        !           154:        switch (tu|tv) {
        !           155: 
        !           156:        case uSMALL|vSMALL: return SmallIntVal(u) - SmallIntVal(v);
        !           157: 
        !           158:        case uSMALL|vINT:
        !           159:        case uINT|vSMALL:
        !           160:        case uINT|vINT: return int_comp((integer)u, (integer)v);
        !           161: 
        !           162:        case uSMALL|vRAT:
        !           163:        case uINT|vRAT:
        !           164:                u = (value) mk_rat((integer)u, int_1, 0);
        !           165:                s = rat_comp((rational)u, (rational)v);
        !           166:                release(u);
        !           167:                return s;
        !           168: 
        !           169:        case uRAT|vRAT:
        !           170:                return rat_comp((rational)u, (rational)v);
        !           171: 
        !           172:        case uRAT|vSMALL:
        !           173:        case uRAT|vINT:
        !           174:                v = (value) mk_rat((integer)v, int_1, 0);
        !           175:                s = rat_comp((rational)u, (rational)v);
        !           176:                release(v);
        !           177:                return s;
        !           178: 
        !           179:        case uSMALL|vAPP:
        !           180:        case uINT|vAPP:
        !           181:        case uRAT|vAPP:
        !           182:                u = approximate(u);
        !           183:                s = app_comp((real)u, (real)v);
        !           184:                release(u);
        !           185:                return s == 0 ? -1 : s; /* u < ~u = v */
        !           186: 
        !           187:        case uAPP|vAPP:
        !           188:                return app_comp((real)u, (real)v);
        !           189: 
        !           190:        case uAPP|vSMALL:
        !           191:        case uAPP|vINT:
        !           192:        case uAPP|vRAT:
        !           193:                v = approximate(v);
        !           194:                s = app_comp((real)u, (real)v);
        !           195:                release(v);
        !           196:                return s == 0 ? 1 : s; /* u = ~v > v */
        !           197: 
        !           198:        default: syserr(MESS(1305, "num_comp")); /* NOTREACHED */
        !           199: 
        !           200:        }
        !           201: }
        !           202: 
        !           203: 
        !           204: /*
        !           205:  * Deliver 10**n, where n is a (maybe negative!) C integer.
        !           206:  * The result is a value (integer or rational, actually).
        !           207:  */
        !           208: 
        !           209: Visible value tento(n) int n; {
        !           210:        if (n < 0) {
        !           211:                integer i= int_tento(-n);
        !           212:                value v= (value) mk_exact(int_1, i, 0);
        !           213:                release((value) i);
        !           214:                return v;
        !           215:        }
        !           216:        return (value) int_tento(n);
        !           217: }
        !           218: 
        !           219: 
        !           220: /*
        !           221:  * numval returns the numerical value of any numeric B value
        !           222:  * as a C `double'.
        !           223:  */
        !           224: 
        !           225: Visible double numval(u) value u; {
        !           226:        double expo, frac;
        !           227: 
        !           228:        if (!Is_number(u)) {
        !           229:                error(MESS(1306, "value not a number"));
        !           230:                return 0.0;
        !           231:        }
        !           232:        u = approximate(u);
        !           233:        expo = Expo((real)u), frac = Frac((real)u);
        !           234:        release(u);
        !           235:        if (expo > Maxexpo) {
        !           236:                error(MESS(1307, "approximate number too large to be handled"));
        !           237:                return 0.0;
        !           238:        }
        !           239:        if(expo < Minexpo) return 0.0;
        !           240:        return ldexp(frac, (int)expo);
        !           241: }
        !           242: 
        !           243: 
        !           244: /*
        !           245:  * Random numbers
        !           246:  */
        !           247: 
        !           248: 
        !           249: /*
        !           250:  * numhash produces a `double' number that depends on the value of
        !           251:  * v, useful for initializing the random generator.
        !           252:  * Needs rewriting, so that for instance numhash(n) never equals n.
        !           253:  * The magic numbers here are chosen at random.
        !           254:  */
        !           255: 
        !           256: Visible double numhash(v) value v; {
        !           257:        if (Integral(v)) {
        !           258:                double d = 0;
        !           259:                register int i;
        !           260: 
        !           261:                if (IsSmallInt(v)) return SmallIntVal(v);
        !           262:                for (i = Length(v) - 1; i >= 0; --i) {
        !           263:                        d *= 2;
        !           264:                        d += Digit((integer)v, i);
        !           265:                }
        !           266: 
        !           267:                return d;
        !           268:        }
        !           269: 
        !           270:        if (Rational(v))
        !           271:                return .777 * numhash((value) Numerator((rational)v)) +
        !           272:                       .123 * numhash((value) Denominator((rational)v));
        !           273: 
        !           274:        return numval(v); /* Which fails for HUGE reals.  Never mind. */
        !           275: }
        !           276: 
        !           277: 
        !           278: /* Initialize the random generator */
        !           279: 
        !           280: double lastran;
        !           281: 
        !           282: Hidden Procedure setran (seed) double seed; {
        !           283:        double x;
        !           284: 
        !           285:        x = seed >= 0 ? seed : -seed;
        !           286:        while (x >= 1) x /= 10;
        !           287:        lastran = floor(67108864.0*x);
        !           288: }
        !           289: 
        !           290: Visible Procedure set_random(v) value v; {
        !           291:        setran((double) hash(v));
        !           292: }
        !           293: 
        !           294: 
        !           295: /* Return a random number in [0, 1). */
        !           296: 
        !           297: Visible value random() {
        !           298:        double p;
        !           299: 
        !           300:        p = 26353589.0 * lastran + 1;
        !           301:        lastran = p - 67108864.0*floor(p/67108864.0);
        !           302: 
        !           303:        return (value) mk_approx(lastran / 67108864.0, 0.0);
        !           304: }
        !           305: 
        !           306: Visible Procedure initnum() {
        !           307:        rat_init();
        !           308:        setran((double) SEED);
        !           309: }
        !           310: 
        !           311: Visible Procedure endnum() {
        !           312:        endrat();
        !           313: }

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