Logo Search packages:      
Sourcecode: eli version File versions

strmath.c

static char rcsid[] = "$Id: strmath.c,v 3.8 1998/10/14 14:42:51 mjung Exp $";
/*    strmath.c: string mathematical routines.
 *
 *    Revision:   1.0  
 *    Author:   Bruce K. Haddon  
 *    Log:   C:/src/strmath/vcs/strmath.c_v  
 *
 *    Rev 1.0   20 Oct 1990 14:16:30   Bruce K. Haddon
 * Initial revision.
 *
 *    Release: 2.0 Beta 
 */
/* Copyright 1990, Bruce K. Haddon */

/* This file is part of the Eli Module Library.

The Eli Module Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Library General Public License as
published by the Free Software Foundation; either version 2 of the
License, or (at your option) any later version.

The Eli Module Library is distributed in the hope that it will be
useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
Library General Public License for more details.

You should have received a copy of the GNU Library General Public
License along with the Eli Module Library; see the file COPYING.LIB.
If not, write to the Free Software Foundation, Inc., 59 Temple Place -
Suite 330, Boston, MA 02111-1307, USA.  */

/* As a special exception, when this file is copied by Eli into the
   directory resulting from a :source derivation, you may use that
   created file as a part of that directory without restriction. */

#define _MSC_VER 600
#define UNSIGNED (0)


#if defined(__STDC__) || STDDEF > 0
#include <stddef.h>
#endif
#if defined(__STDC__) || STDARG > 0 || defined(__cplusplus)
#include <stdarg.h>
#else
#include <varargs.h>
#endif
#include <string.h>
#include <ctype.h>
#include <errno.h>
#include "strmath.h"

#define SIGNIFICAND_SIZE (2 * ARITH_SIZE + 2)

#define SIGN(digit, base) ( 2 * (int)(digit) >= (base) ? -1 :\
            ( 2 * (int)(digit) == (base) - 1 ? 0 : 1 ))
#define ABS(x) ( (x) >= 0 ? (x) : -(x) )
#define MIN(x,y) ((x) <= (y)? (x) : (y) )
#define MAX(x,y) ((x) >= (y)? (x) : (y) )
#define DIM(x) ((int)(sizeof(x)/sizeof(x[0])))
#define SETERROR(x) ( error ? 0 : (error = (int)(x), 0 ) )

#define ADX(X,N) {(X)->x += (N);}
#define INTXP(X) ((X)->x + 1)
#define SETXP(X,N) {(X)->x = (N)-1;}

typedef char digit_t;

typedef struct number_tag
      {
            long x;
            digit_t *r;
            int significand_size;
            char inexact;
      } number;

/***  Format of string numbers
 *~
 *    The strings accepted by this package are essentially those
 *    of most higher level languages, i.e.,
 *
 *>         [+/-][d*][.][d*][e[+/-]d*]
 *~
 *    where [] indicate optional parts, +/- indicates that a sign may 
 *    be present, |d indicates digits in the chosen base, * indicates 
 *    repetition, `.' is a period representing the separation between 
 *    the integer and fractional parts of the number, and `e' stands 
 *    for any one of the exponent symbols given in the control string 
 *    EXP_SYMBOL.  Strings generated by the package have the same 
 *    formats.  The actual characters used to represent digits, signs,
 *    fractional separators and exponent symbols are determined by
 *    default or |strmath settings (see function |strmath). 
 ***/

/**   The control structure contains the values that control 
 *    interpretation and generation of strings (see |strmath).
 **/

00099 static struct strmath_t {
      char *digits;
      char *exp_symbols;
      char *signs;
      char *separators;
      int exp_base;
      int integer_size;
      int round_size;
      int rounding;
      int denormalize;
      int inexact;
      int check_digits;
      int ignore_case;
      } control =
      {
            STRM_DIGITS_DEFAULT,
            STRM_EXP_SYMBOLS_DEFAULT,
            STRM_SIGNS_DEFAULT,
            STRM_SEPARATORS_DEFAULT,
            STRM_EXP_BASE_DEFAULT,
            STRM_INTEGER_SIZE_DEFAULT,
            STRM_ROUND_SIZE_DEFAULT,
            STRM_ROUNDING_DEFAULT,
            STRM_DENORMALIZE_DEFAULT,
            STRM_INEXACT_DEFAULT,
            STRM_CHECK_DIGITS_DEFAULT,
            STRM_IGNORE_CASE_DEFAULT
      };

static digit_t sm[SIGNIFICAND_SIZE];
static number m = {0L, sm, DIM(sm), '\0'};
static digit_t sn[SIGNIFICAND_SIZE];
static number n = {0L, sn, DIM(sn), '\0'};
static digit_t sr[6 * SIGNIFICAND_SIZE];
static number r = {0L, sr, DIM(sr), '\0'};
static digit_t ss[SIGNIFICAND_SIZE];
static number s = {0L, ss, DIM(ss), '\0'};
static digit_t sw[SIGNIFICAND_SIZE];
static number w = {0L, sw, DIM(sw), '\0'};
static digit_t sz[SIGNIFICAND_SIZE];
static number z = {0L, sz, DIM(sz), '\0'};

static int error;
static char *in;

#if defined(__STDC__) || defined(__cplusplus) || PROTOTYPES > 0
static int convert(number *, int, number *, int);
static int normalize(number *, int);
static int round(number *, int, int, int);
static int subtract(number *, number *, int);
static int multiply(int, number *, number *, number *, int);
static int divide(int, number *, number *, number *, int, int);
static int multbyd(number *, number *, int, int, int);
static int divby2(number *, int, int);
static int setvalue(number *, int, long, int);
static int loadval(number *, char *, int, int, int);
static char *unloadval(number *, char *, int);
static char *ltoab(long, int, int);
static int subber(number *, number *, int);
static int clear(number *);
static int negate(number *, int);
static int copy(number *, number *);
static int is_notless(number *, number *);
static int is_integer(number *);
static int shift(number *, long, int);
static int is_zero(number *, int);
static int ctod(int);
static char dtoc(int);
#else
static int convert();
static int normalize();
static int round();
static int subtract();
static int multiply();
static int divide();
static int multbyd();
static int divby2();
static int setvalue();
static int loadval();
static char *unloadval();
static char *ltoab();
static int subber();
static int clear();
static int negate();
static int copy();
static int is_notless();
static int is_integer();
static int shift();
static int is_zero();
static int ctod();
static char dtoc();
#endif


static char rstr[6 * ARITH_SIZE + 5]; /* Result string */

#if defined(__STDC__) || defined(__cplusplus) || STDARG > 0
/***  strmath --- set control values for strmath.
 *                                        \+1*/
int
strmath(int select, ...)
/*+~
 *    The argument |select of this function is one of the following 
 *    selection values.  A second value is required for each, to set a 
 *    control value for the operation of the strmath functions.
 *
 *    STRM_DIGITS, <string> sets the elements of the string to 
 *    represent the digits of the numbers represented by the strings. 
 *    e.g., "0123456789ABCDEF" for bases 2-16.  The matching of 
 *    alphabetic characters in input strings is further controlled by 
 *    STRM_IGNORE_CASE.  Output strings are build using strict 
 *    correspondence to the positions in this string.  The release 
 *    value is the string consisting of the digits '0'-'9', the upper 
 *    case alphabetic characters, the lower case alphabetic 
 *    characters, and the characters '%' and '$' (64 characters for 
 *    base 64).  Unpredictable results will be obtained if any character 
 *    is repeated. 
 *
 *    STRM_EXP_SYMBOLS, <string> sets of characters to be accepted as 
 *    representations of the exponent symbol.  Any character in the 
 *    string is acceptable in input strings; in output strings the 
 *    first member of the set is used.  The release value of this 
 *    string is "^".  Unpredictable results will be obtained if any 
 *    character is repeated or is used as a digit. 
 *
 *    STRM_SIGNS, <string> sets the characters to be accepted as 
 *    representations of negative and positive signs.  The first 
 *    member of the string is the character representing a negation, 
 *    and the remainder are all taken as representations of a positive 
 *    value. In output strings the first member is used to show 
 *    negation, and the second value to represent a positive sign (if 
 *    one is requested).  The release value of this string is "-+". 
 *    Unpredictable results will be obtained if any character is 
 *    repeated or is used as an exponent symbol or digit. 
 *
 *    STRM_SEPARATORS,<string> sets the characters to be accepted as 
 *    representations of the fractional separator.  Any are acceptable 
 *    in input strings; in output strings the first member of the set 
 *    is used.  The release value of this string is ".". Unpredictable 
 *    results will be obtained if any character is repeated or is used 
 *    as a sign, an exponent symbol, or digit.
 *
 *    STRM_EXP_BASE, <integer> defines the numeric base of the 
 *    exponents of numbers.  If this value is zero, the base of an 
 *    exponent is taken to be the same as the number itself (and is 
 *    supplied on each call to the strmath functions).  If non-zero, 
 *    then this value is used as the base of the exponent, independent 
 *    of the base of the number itself.
 *
 *    STRM_INTEGER_SIZE, <integer> is the maximum number of digits 
 *    that a representation of an integer may have.  This value, plus 
 *    STRM_ROUND_SIZE, may be no greater than ARITH_SIZE.  Values with 
 *    more significant digits than STRM_INTEGER_SIZE will be returned 
 *    represented as values with exponents, even if integral in value.
 *    The release value is ARITH_SIZE - 2.
 * 
 *    STRM_ROUND_SIZE, <integer> is the number of additional digits a 
 *    value may have if not represented as an integer.  This value,
 *    plus STRM_INTEGER_SIZE, may be no greater than ARITH_SIZE.
 *    Values with more digits than STRM_INTEGER_SIZE + STRM_ROUND_SIZE
 *    will be rounded on output to values with no more than that many
 *    digits.  The release value is 2.
 *
 *    STRM_ROUNDING, <mode> selects the rounding mode to be applied to 
 *    output values.  The modes are represented by one of the 
 *    following rounding selections:
 *          
 *    STRM_EVEN_ROUND: the output value is rounded to the nearest value 
 *          representable in the given number of digits.  If two 
 *          values are equally near, the one with an even last 
 *          digit will be returned.  This the release setting.
 *    STRM_ZERO_ROUND: the output value is rounded to the value
 *          representable in the given number of digits whose value 
 *          is closest to and no greater in magnitude.      
 *    STRM_UP_ROUND: the output value is rounded to the value
 *          representable in the given number of digits whose value 
 *          is closest to but no less than in magnitude.    
 *    STRM_DOWN_ROUND: the output value is rounded to the value
 *          representable in the given number of digits whose value 
 *          is closest to but no greater than in magnitude. 
 *    STRM_HAND_ROUND: the output value is rounded to the nearest value 
 *          representable in the given number of digits.  If two 
 *          values are equally near, the one with the greater 
 *          magnitude will be returned (as is normally done when
 *          computing by hand).
 * 
 *    STRM_DENORMALIZE, <0/1> selects, if 1, to allow values with 
 *    exponents that would normally be -1 to be instead represented on 
 *    output with a 0-valued digit before the fractional separator.  
 *    If set to 0, then the initial digit of an output value will 
 *    always be nonzero.  The release setting is 1.
 *    
 *    STRM_INEXACT, <0/1> selects, if 1, that input numbers with a 
 *    fractional separator are to be taken to represent inexact 
 *    values, and that, inexact values (values that have lost digits 
 *    in the course of computation, or result from inexact values) are 
 *    to be output with a fractional point included, followed by a 
 *    zero digit if there are no other digits.  If the 0 value is 
 *    selected, no fractional point will be included if all digits 
 *    after the fractional point are zeros. The release setting is 1.
 *      
 *    STRM_CHECK_DIGITS, <0/1> selects, if 1, that all characters of a 
 *    string be checked to see if they are valid digits (as defined by 
 *    the STRM_DIGITS string, and are     legitimate digits in the current 
 *    base.  If this is not so, an error indication is returned.  With 
 *    the 0 value selected, the number scan ends with the first 
 *    character that is not legitimately part of the number, and no 
 *    error is signalled.  The release setting is 1.
 *    
 *    STRM_IGNORE_CASE, <0/1> selects, if 1, that the case of 
 *    characters in input values, and the case of characters in the 
 *    STRM_DIGITS string, is to be ignored when matching and defining 
 *    the value of elements of an input number.  The 0 value indicates 
 *    that case is significant.  The release setting is 1.
 *    
 *    |strmath returns the value 1 if the requested setting succeeded,
 *    otherwise 0 is returned, and no selected value has been changed.
 *    The caller is responsible for ensuring that the contents of the
 *    string values are disjoint sets as described above.
 *
 ***/ 
#else
int
strmath(va_alist)
va_dcl
#endif
{ 
      va_list argptr;
      register char *temp1;
      register int temp2; 
      register int result = 0;

#if defined(__STDC__) || defined(__cplusplus) || STDARG > 0
      va_start(argptr, select);
#else
      int select;

      va_start(argptr);

      select = va_arg(argptr, int);
#endif

      switch( select )
      {
      case STRM_DIGITS:
            if( (temp1 = va_arg(argptr, char *)) )
            {
                  control.digits = temp1;
                  result = 1;
            }
            break;

      case STRM_EXP_SYMBOLS:
            if( (temp1 = va_arg(argptr, char *)) )
            {
                  control.exp_symbols = temp1;
                  result = 1;
            }
            break;

      case STRM_SIGNS:
            if( (temp1 = va_arg(argptr, char *)) )
            {
                  control.signs = temp1;
                  result = 1;
            }
            break;

      case STRM_SEPARATORS:
            if( (temp1 = va_arg(argptr, char *)) )
            {
                  control.separators = temp1;
                  result = 1;
            }
            break;

      case STRM_EXP_BASE:
            if( (temp2 = va_arg(argptr, int)) == 0 ||
                        (temp2 >= 2 && temp2 <= 64) )
            {
                  control.exp_base = temp2;
                  result = 1;
            }
            break;

      case STRM_INTEGER_SIZE:
            if( (temp2 = va_arg(argptr, int)) > 0 &&
                        temp2 + control.round_size <= ARITH_SIZE)
            {
                  control.integer_size = temp2;
                  result = 1;
            }
            break;

      case STRM_ROUND_SIZE:
            if( (temp2 = va_arg(argptr, int)) >= 0 &&
                        temp2 + control.integer_size <= ARITH_SIZE)
            {
                  control.round_size = temp2;
                  result = 1;
            }
            break;

      case STRM_ROUNDING:
            if( (temp2 = va_arg(argptr, int)) >= 0 &&
                        temp2 <= STRM_HAND_ROUND )
            {
                  control.rounding = temp2;
                  result = 1;
            }
            break;

      case STRM_DENORMALIZE:
            if( (temp2 = va_arg(argptr, int)) >= 0 && temp2 <= 1 )
            {
                  control.denormalize = temp2;
                  result = 1;
            }
            break;

      case STRM_INEXACT:
            if( (temp2 = va_arg(argptr, int)) >= 0 && temp2 <= 1 )
            {
                  control.inexact = temp2;
                  result = 1;
            }
            break;

      case STRM_CHECK_DIGITS:
            if( (temp2 = va_arg(argptr, int)) >= 0 && temp2 <= 1 )
            {
                  control.check_digits = temp2;
                  result = 1;
            }
            break;

      case STRM_IGNORE_CASE:
            if( (temp2 = va_arg(argptr, int)) >= 0 && temp2 <= 1 )
            {
                  control.ignore_case = temp2;
                  result = 1;
            }
            break;

      default:
            break;
      }
#if !defined(__STDC__) && STDARG == 0
      va_end(ap);
#endif

      return result;
      
} /* end math */

/***  Error values
 *
 *    In general (with the exception of |strnumb) the functions of
 *    |strmath return a NULL pointer when an error is found in the
 *    input values, or occurs during computation (e.g., overflow).
 *    When a null value is returned, the global |errno is set to 
 *    indicate the type of error.  If multiple errors occur, only the
 *    first is reported.
 *
 *    The error codes, defined in <errno.h>  used are:
 *    EINVAL: the was in error in the format of an input value,
 *          making it invalid.  (The function |strnumb can be used
 *          to assist in diagnozing these problems.
 *    EDOM: the given value falls outside the domain of the operation,
 *          e.g., it is too large, or negative when a positive
 *          value is required.
 *    ERANGE: the result of the computation cannot be adequately 
 *          represented (e.g., divide by zero, overflow).
 * 
 *    The usual rules for the use of |errno apply.  It is not 
 *    explicitly cleared.  If desired, this must be done before 
 *    calling the desired function.  It is set by these functions only
 *    if an error occurs.
 ***/

/***  strnumb --- scan an string representing a value.
 *                                        \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strnumb(char *a, int *flag, int base)
#else
strnumb(a, flag, base)
char *a;
int *flag;
int base;
#endif
/*+~
 *    The string |a is scanned and interpreted as a value, in the 
 *    given base.  If an error is found, the appropriate error value 
 *    (as described above) is returned in the variable |flag, 
 *    otherwise |flag will be set to zero.  A pointer is returned, 
 *    pointing into the string |a, indicating the last character 
 *    examined.  It will be the character causing the error if an 
 *    error value is returned, otherwise it will be the character that 
 *    terminated the scan of the number.  The termination of the scan 
 *    will be controlled by the setting of STRM_CHECK_DIGITS. 
 ***/
{
      (void) loadval(&r, a, 1, 0, base);
      *flag = error;
      return in;
      
} /* end numb */


/***    strsqrt --- find sqrt of value represented by a string.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strsqrt(char *a, int base)
#else
strsqrt(a, base)
char *a;
int base;
#endif
/*+~
 *      The square root of the value represented by the string |a
 *      is found, returning the answer as a (static) string,
 *    or NULL if the value represented is negative, or there are
 *    errors in the string.
 ***/
{     /* This algorithm is taken from "Software Manual for the Elementary
       * Functions", by W. J. Cody and W. M. Waite (Prentice-Hall, 1980
       * ISBN 0-13-822064-6).  The implementation for non-binary floating
       * point machines was used.  Comments refer to numbers in the flow
       * chart on page 18.
       */

      int N, count;

      /* 2, 3 */
      if( loadval(&r, a, 1, 0, base) < 0 || error || ABS(r.x) > MAX_EXP  )
      {
            SETERROR(EDOM);
            return (char *) NULL;
      }

      /* 1 */
      (void) copy(&m, &r);
      (void) normalize(&m, UNSIGNED);
      if (is_zero(&m, 0)) return unloadval(&m, (char *)NULL, base);

      /* 4 */
      N = INTXP(&m);

      /* 5 */
      if (N & 01) {
            int approx = 0, temp = (int)m.r[0];
            SETXP(&m, 1); N--;
            while( (approx + 1) * approx < temp ) ++approx;
            (void)setvalue(&n, approx, 1L, base);
      } else  {
            SETXP(&m, 0);

      /* 6 */
            copy(&r, &m);
            loadval(&n, ".175241", 1, 0, 10);
              if (base != 10) {
                      (void)convert(&w, base, &n, 10);
                      copy(&n, &w);
              }
            (void)negate(&n, base);
            (void)subtract(&r, &n, base);
            loadval(&n, ".086462", 1, 0, 10);
              if (base != 10) {
                      (void)convert(&w, base, &n, 10);
                      copy(&n, &w);
              }
            (void)divide(1, &s, &n, &r,
                  control.round_size + control.integer_size + 2, base);
            loadval(&n, ".580661", 1, 0, 10);
              if (base != 10) {
                      (void)convert(&w, base, &n, 10);
                      copy(&n, &w);
              }
            copy(&r, &m);
            (void)divby2(&r,
                  control.round_size + control.integer_size + 2, base);
            (void)negate(&r, base);
            (void)subtract(&n, &r, base);
            (void)subtract(&n, &s, base);
      }

      /* 7 */
      for (count = 1;
           count < control.round_size + control.integer_size + 2;
           count *= 4) {
            copy(&r, &m); copy(&s, &n);
            (void)divide(-1, &z, &r, &s,
                  control.round_size + control.integer_size + 2, base);
            (void)subtract(&n, &z, base);
      
            copy(&r, &m); copy(&s, &n);
            (void)divide(-1, &z, &r, &s,
                  control.round_size + control.integer_size + 2, base);
            loadval(&s, ".25", 1, 0, 10);
            if (base != 10) {
                  (void)convert(&w, base, &s, 10);
                  copy(&s, &w);
            }
            (void)multiply(1, &r, &s, &n, base);
            (void)subtract(&r, &z, base);
            copy(&n, &r);
      }

      /* 10 */
      ADX(&n, N/2);

      return unloadval(&n, (char *) NULL, base);

} /* end sqrt */


/***    strpow --- find power of value represented by a string.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strpow(char *a, char *b, int base)
#else
strpow(a, b, base)
char *a;
char *b;
int base;
#endif
/*+~
 *      An integer power of the value represented by the string |a
 *      is found, the power being given by string |b, returning the
 *    answer as a (static) string, or NULL if |b is not an integer,
 *    or there are other errors in the strings, or the result 
 *    overflows.
 ***/
{

      int signflag1 = loadval(&r, a, 0, 0, base);
      int signflag2 = loadval(&m, b, 0, error, base);
      int zero = is_zero(&r, 0);

      /* Error if the value to be raised is zero and the power to
         which to raise it zero or negative. */
      if( !is_integer(&m) || (zero &&
                  (is_zero(&m, 0) || signflag2 < 0)) )
            SETERROR(EDOM);

      if( error || zero )
            return unloadval(&r, (char * )NULL, base);

      /* If the power value is even, then the sign of the result
         must be positive, otherwise the sign of the result will be
         the sign of the value to be raised.  */
      (void) copy(&s, &m);
      if( !divby2(&s, (int) s.x + 1, base) )
            signflag1 = 1;

      (void) setvalue(&n, 1, 1L, base);

      while( !is_zero(&m, 0) )
      {
            if( divby2(&m, (int) m.x + 1, base) )
            {
                  (void) copy(&s, &n);
                  if( !multiply(0, &n, &s, &r, base) )
                        SETERROR(ERANGE);
            }

            (void) copy(&s, &r);
            (void) copy(&z, &r);
            if( !multiply(0, &r, &s, &z, base) )
                  SETERROR(ERANGE);
      }

      /* If sign of power was negative, then invert the result,
         i.e., divide into 1 */
      if( signflag2 < 0 )
      {
            (void) setvalue(&s, 1, 1L, base);
            (void) copy(&z, &n);
            (void) divide(0, &n, &s, &z, control.round_size +
                              control.integer_size + 2, base);
      }

      if( signflag1 < 0 )
            (void) negate(&n, base);

      return unloadval(&n, (char *) NULL, base);

} /* end pow */


/***    stradd --- add values represented by strings.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
stradd(char *a, char *b, int base)
#else
stradd(a, b, base)
char *a;
char *b;
int base;
#endif
/*+~
 *      The values represented by the strings |a and |b are added,
 *      returning the answer as a (static) string, or NULL if there
 *    are errors in the strings, or the result overflows.
 ***/
{
      (void) loadval(&m, a, 1, 0, base);
      (void) loadval(&n, b, -1, error, base);
      (void) subtract(&m, &n, base);
      return unloadval(&m, (char *) NULL, base);

} /* end add */


/***    strsub --- subtract values represented by strings.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strsub(char *a, char *b, int base)
#else
strsub(a, b, base)
char *a;
char *b;
int base;
#endif
/*+~
 *      The values represented by the strings |a and |b are
 *      subtracted, returning the answer as a (static) string,
 *    or NULL if there are errors in the strings, or the result
 *    overflows.
 ***/
{
      (void) loadval(&m, a, 1, 0, base);
      (void) loadval(&n, b, 1, error, base);
      (void) subtract(&m, &n, base);
      return unloadval(&m, (char *) NULL, base);

} /* end sub */


/***    strmult --- multiply values represented as strings.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strmult(char *a, char *b, int base)
#else
strmult(a, b, base)
char *a, *b;
int base;
#endif
/*+~
 *      The strings |a and |b are multiplied, and the answer is
 *      returned as a (static) string, or NULL if there are errors
 *    in the strings, or the result overflows.
 ***/
{
      int signflag;

      signflag = loadval(&m, a, 0, 0, base);
      signflag *= loadval(&n, b, 0, error, base);
      (void) multiply(signflag, &r, &n, &m, base);
      return unloadval(&r, (char *) NULL, base);

} /* end mult */


/***    strdivf --- (floating) divide values represented as strings.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strdivf(char *a, char *b, int base)
#else
strdivf(a, b, base)
char *a, *b;
int base;
#endif
/*+~
 *      The string |a is divided by |b, and the answer is
 *    returned as a (static) string, representing the result of
 *    a full divide (with fractional part), or NULL if there
 *    are errors in the strings, or the result overflows.
 ***/
{
      int signflag;

      signflag = loadval(&m, a, 0, 0, base);
      signflag *= loadval(&n, b, 0, error, base);
      (void )divide(signflag, &r, &m, &n, control.integer_size
            + control.round_size + 2, base);
      return unloadval(&r, (char *)
      NULL, base);

} /* end divf */


/***    strdivi --- (integer) divide values represented by strings.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strdivi(char *a, char *b, int base)
#else
strdivi(a, b, base)
char *a, *b;
int base;
#endif
/*+~
 *      The value |b is divided into |a, and the answer is returned
 *    as a (static) string, representing the result of an integer
 *    divide, i.e., how many integral times the divisor |b will go
 *    into the dividend |a, or NULL if there are errors in the
 *    strings, or the result overflows.
 ***/
{
      int signflag;
      long positions;

      signflag = loadval(&m, a, 0, 0, base);
      signflag *= loadval(&n, b, 0, error, base);

      if( (positions = m.x - n.x + 1) <= MAX_EXP )
            (void) divide(signflag, &r, &m, &n, (int) positions, base);
      else
            SETERROR(EDOM);

      r.inexact = 0;
      return unloadval(&r, (char *) NULL, base);

} /* end divi */


/***    strrem --- remainder of values represented by strings.
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strrem(char *a, char *b, int base)
#else
strrem(a, b, base)
char *a, *b;
int base;
#endif
/*+~
 *      The value |b is divided into |a, and the answer is returned
 *    as a (static) string being the remainder of an integer division
 *    (the fractional part of the result of the division), or NULL if
 *    there are errors in the strings, or the result overflows.
 ***/
{
      int signflag;
      long positions;

      signflag = loadval(&m, a, 0, 0, base);
      (void) loadval(&n, b, 0, error, base);

      if( is_zero(&n, 0 ) )
      {
            if( signflag < 0 )
                  (void) negate(&m, base);
      } else if( (positions = m.x - n.x + 1) <= MAX_EXP )
            (void) divide(signflag, &r, &m, &n, (int) positions, base);
      else
            SETERROR(EDOM);

      return unloadval(&m, (char *) NULL, base);

} /* end rem */


/***  strnorm --- normalize a value represented by a string
 *                                              \+1*/
char *
#if defined(__STDC__) || defined(__cplusplus)
strnorm(char *a, int oldbase, int newbase, char *symbols)
#else
strnorm(a, oldbase, newbase, symbols)
char *a;
int oldbase;
int newbase;
char *symbols;
#endif
/*+~
 *    The value represented by string |a is converted from radix
 *    |oldbase to radix |newbase and normalized.  The result is
 *    returned as a (static) string.
 *
 *    If |symbols is NULL then the result will be formatted according
 *    to the normal output rules (i.e., with the defined fractional
 *    separator if the result cannot be represented as an integral
 *    value, and with an exponent marker and exponent if the value
 *    cannot be represented as an integer_size integer).
 *
 *    If |symbols is the empty string (i.e., the first character is
 *    '\0') and the result cannot be represented as an integer within
 *    integer_size characters), then a NULL string will be returned.
 *
 *    Otherwise (if |symbols is a non-empty string) the result is
 *    guaranteed to contain a fractional separator (with a single zero
 *    character if the fractional part of the value is zero).  In this
 *    case, any exponent marker will be the first character of
 *    |symbols, and then, if further characters are present in the
 *    |symbols string, they will be used, respectively, as the
 *    fractional separator, the minus sign, and a plus sign.  If these
 *    additional characters are not provided, the normally defined
 *    characters will be used.
 ***/
{
      int signflag;

      if( oldbase == newbase )
            /* This is just a normalization without a change of
               base, so load directly to the final number. */
            signflag = loadval(&n, a, 0, 0, oldbase);
      else
      {
            /* The base is to be converted, so load the original
               number. */
            signflag = loadval(&r, a, 0, 0, oldbase);
            if( !convert(&n, newbase, &r, oldbase) )
                  SETERROR(ERANGE);
      }


      if( signflag < 0 )
            (void) negate(&n, newbase);

      if( symbols && !*symbols && !is_integer(&n) )
            SETERROR(ERANGE);

      return unloadval(&n, symbols, newbase);

}

/**   convert --- from one base to another
 *                                        \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
convert(number *reg1, int newbase, number *reg2, int oldbase)
#else
convert(reg1, newbase, reg2, oldbase)
number *reg1;
int newbase;
number *reg2;
int oldbase;
#endif
/*+~
 *    The value in |reg2, represented in |oldbase, is converted to
 *    a new representation in |reg1, in |newbase. The value 1 is
 *    returned if this is done successfully, 0 otherwise.
 *
 *    If |inexact is present in |reg2, or any accuracy loss
 *    occurs during the conversion, then |inexact in |reg1 will
 *    be set.
 **/
{
      int carry;
      int limit;
      int i;
      long k;
      long exponent = 0L;

      (void) clear(reg1);
      reg1->inexact = reg2->inexact;

      /* If there is an error in the number, no more work is needed,
         so return immediately with the output value. */

      if( error || reg2->x > MAX_EXP )
            return 0;

      /* If the input is zero, then return now, with reg1 cleared. */
      if( is_zero(reg2, 0) )
            return 1;


      /* If the value is purely fractional, then multiply it by the
         new base until there is a least one non-fraction digit,
         keeping count so to be able to adjust the exponent of the
         result. */

      while( reg2->x < 0L )
      {
            carry = multbyd(reg2, reg2, newbase, 0, oldbase);
            while( carry )
            {
                  (void) shift(reg2, 1L, UNSIGNED);
                  reg2->r[0] = (digit_t) (carry % oldbase);
                  carry /= oldbase;
            }
            --exponent;
      }

      /* The new value will now start with an integer part.  The
         fractional point is at the right-hand end of the number. */
      reg1->x = reg1->significand_size - 1;

      /* Set |limit to be the index one past the last actual
         non-fractional digit to be processed (i.e. the first
         fractional digit if there is one. */

      limit = ( reg2->x > (long) reg2->significand_size ?
                  reg2->significand_size : (int) reg2->x + 1 );

      /* Trailing implied zeros will be treated in the next step:
         trailing real zeros can be treated the same way. */
      while( reg2->r[limit - 1] == 0 )
            --limit;

      /* For each actual non-fraction digit in the number, use
         Horner's method to do the conversion to the new number.
         If the new number overflows, treat remaining digits as
         zeros.  */

      for( i = 0; i < limit; ++i )
      {
            carry = (int) reg2->r[i];

            /* Zero out digits as they are used. */
            reg2->r[i] = 0;

            if( (carry =
                  multbyd(reg1, reg1, oldbase, carry, newbase)) )
            {
                  int j;
                  /* Conversion space has run out, so put back
                     the overflow, and treat any remaining digits
                     as though they were zeros.  */
                  while( carry )
                  {
                        (void) shift(reg1, 1L, UNSIGNED);
                        reg1->r[0] = (digit_t) (carry % newbase);
                        carry /= newbase;
                  }
                  ++i;
                  for( j = i; j < limit; ++j )
                        reg2->r[j] = 0;
                  reg1->inexact = 1;
                  break;
            }
      }

      /* For each non-fraction digit in the number (actual or implied),
         do the conversion to the new number.  */

      for( k = (long) i; k <= reg2->x; ++k )
      {
            carry = multbyd(reg1, reg1, oldbase, 0, newbase);
            while( carry )
            {
                  (void) shift(reg1, 1L, UNSIGNED);
                  reg1->r[0] = (digit_t) (carry % newbase);
                  carry /= newbase;
            }
      }

      if( !normalize(reg1, UNSIGNED) )
            return 0;

      /* Now the conversion of the fractional part may start:
         Firstly the fractional part of the original number is
         normalized: the integer part has been already removed, and
         then any implicit zeros at the beginning of the fractional
         part of made explicit by shifting by the negative of the
         exponent. */

      if( !is_zero(reg2, 0) )
      {
            (void) normalize(reg2, UNSIGNED);
            if( reg2->x < -1L )
                  (void) shift(reg2, -(reg2->x + 1), UNSIGNED);


            for(k = reg1->x + 1L;
                        k < (long) reg1->significand_size; ++k )
                  reg1->r[(int) k] =
                        (digit_t) multbyd(reg2, reg2, 
                                    newbase, 0, oldbase);
      }

      /* Put reg1 back into signed form */
      if( SIGN(reg1->r[0], newbase) <= 0 )
            (void) shift(reg1, 1L, UNSIGNED);

      /* Adjust the exponent for the initial "normalization". */
      reg1->x += exponent;

      return 1;

} /* end convert */


/**     normalize --- normalize the content of a number.
 *                                        \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
normalize(number *reg, int base)
#else
normalize(reg, base)
number *reg;
int base;
#endif
/*+~
 *      The given number, |reg, is normalized in the given base, or,
 *    if the base is zero, |reg is normalized as an unsigned value.
 **/
{
      register int i = 0;
      register digit_t *r0 = reg->r;
      digit_t sgn = (digit_t)( !base ? 0 : 
                  ( SIGN(r0[0], base) < 0 ? base - 1 : 0 ));

      /* If the exponent is too large, then this is an error. */
      if( ABS(reg->x) > MAX_EXP )
            return 0;

      /* Look along the number, looking for the first digit
         that is not the same as the sign in the current base. */

      while( i < reg->significand_size && (sgn == *r0++) )
                         ++i;

      /* If the entirety of the number is zero, then the value
         is zero, and this is normally represented with the exponent
         also zero, otherwise move the significand to the left,
         correspondingly decreasing the exponent. */

      if( !sgn && i == reg->significand_size )
      {
            reg->x = 0;
            return 1;
      }
      

      /* If the digit that is different to the sign represents
         the same sign, then it may be moved to the zeroth digit
         position, otherwise move the digit before it to the
         digit position.  (Note, if i == 0, then the sign represented
         by the digit will always be the same as the sign of the value) */

      if( base && (SIGN(sgn, base) != SIGN(*--r0, base)) )
            --i;

      if( i != 0 )
      {
            int j = i;  
            register digit_t *rd = reg->r;
      
            r0 = &reg->r[i];
            reg->x -= i;
            
            /* Move the register up the correct number of places. */
            for(i = reg->significand_size - 1 - j; i >= 0 ; --i )
                  *rd++ = *r0++;
                  
            /* Fill the remainder with 0's. */
            for(i = j - 1; i >= 0 ; --i )
                  *rd++ = (digit_t) 0;
      }

      return 1;

} /* end normalize */


/**     round --- round a number.
 *                                        \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
round(number *reg, int round_option, int round_size, int base)
#else
round(reg, round_option, round_size, base)
number *reg;
int round_option;
int round_size;
int base;
#endif
/*+~
 *    The given number, |reg is rounded at the digit specified by
 *    |round_size.  The actual style of rounding is specified by
 *    |round_option, and may be one of the five allowed types of
 *    rounding.  If the number overflows during rounding, the 
 *    number is right-shifted one place, the sign restored, and the
 *    exponend incrmented.  If digits are lost in the rounding, then
 *    |inexact is set.
 **/
{
      register int j;
      register digit_t *r0 = reg->r;
      int operate = 0;

      /* If the first digit serves only as a sign, then move the
         rounding position over by one to accomodate this fact. */

      if( r0[0] == 0 || r0[0] == (digit_t)(base - 1) )
            ++round_size;

      /* If the rounding changes anything (i.e., there are digits
         at the rounding position or beyond, then set the inexact
         flag on the number.  */

      if( !is_zero(reg, round_size) )
            reg->inexact = 1;

      switch( round_option )
      {
            /* Round to nearest; 0.5 goes away from zero. */
      case STRM_HAND_ROUND:
            operate = ( 2 * ((int) r0[round_size] -
                  (SIGN(r0[0], base) < 0 ? 1 : 0 )) >= base );
            break;
            
            /* Round to nearest; 0.5 goes to nearest even value. */
            
            /* Candidate to round if round digit >= 0. 5  0.5 rounds
               up only if the next digit is odd (i.e. rounding up will
               make that digit even). This means give up if next digit
               is even */
      
      default:
      case STRM_EVEN_ROUND:
            operate = ( 2 * ((int) r0[round_size]) >= base ) &&   
                    !(round_size > 0 && 
                    2 * ((int) r0[round_size]) ==     base &&
                    !((int) r0[round_size - 1] & 01));
            break;

            /* round towards zero, i.e., truncate absolute value. */          
      case STRM_ZERO_ROUND:
            if( SIGN(reg->r[0], base ) >= 0 )
                  break;
            /* Fall into UP_ROUND */
            
            /* Round towards plus infinity, i.e., if non-zero
               fraction, truncate and add 1. */                                
      case STRM_UP_ROUND:
            operate = !is_zero(reg, round_size);
            break;

            /* Round towards minus infinity, i.e. truncate. */
      case STRM_DOWN_ROUND:
            break;      
      }
      
      if( operate )
      {
            register digit_t *z0 = &z.r[(j = round_size - 1)];
            
            /* Create value -1 to required number of digits 
               in register z, then subtract it to cause rounding. */
            (void) clear(&z);       
            for( ; j >= 0; --j )
                  *z0-- = (digit_t)(base - 1);
            (void) subber(reg, &z, base);
      }
            
      r0 = &r0[(j = reg->significand_size - 1)];      
      for( ; j >= round_size; --j )
            *r0-- = 0;
            
      return 1;

} /* end round */


/**   subtract --- |reg1 = |reg1 - |reg2.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
subtract(number *reg1, number *reg2, int base)
#else
subtract(reg1, reg2, base)
number *reg1;
number *reg2;
int base;
#endif
/*+~
 *    The value in |reg2 is subtracted from the value in |reg1
 *    and the value is returned in |reg1.  The value of |reg2 is
 *    (in principle) destroyed.
 **/
{
      if ( !is_zero(reg2, 0) ) {
      /* Shift the significand of the number with smaller exponent to
         line up with the significand of the number with larger exponent,
         then run the subtractor. */
            if( is_zero(reg1, 0) )
                  reg1->x = reg2->x;
            else if( reg1->x < reg2->x ) {
                  (void) shift(reg1, (long) reg2->x - reg1->x, base);
                  reg1->x = reg2->x;
            } else
                  (void) shift(reg2, (long) reg1->x - reg2->x, base);
      
            (void) subber(reg1, reg2, base);
            reg1->inexact |= reg2->inexact;
      }

      /* Normalize the result. */
      return normalize(reg1, base);

} /* end subtract */


/**   multiply --- reg1 = reg2 * reg3
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
multiply(int signflag, number *reg1, number *reg2, number *reg3, int base)
#else
multiply(signflag, reg1, reg2, reg3, base)
int signflag;
number *reg1, *reg2, *reg3;
int base;
#endif
/*+~
 *    The positive value in |reg2 is multiplied by the positive value
 *    in |reg3, and the result is placed in |reg1.  The result is given
 *    the sign indicated by |signflag (negative if 1).
 *
 *    At exit, the contents of |reg2 have been destroyed, |reg3 is
 *    left unchanged,   and the result is in |reg1.  The value 1 is
 *    returned if all has gone well.
 *
 *    The number |wreg is used as a temporary working number.
 **/
{
      int i;
      digit_t *r3 = reg3->r;
      long exponent;

      (void) negate(reg2, base);
      exponent = reg2->x + reg3->x + 1;
      (void) clear(reg1);
      if( ABS(exponent) > MAX_EXP ||
            ABS(reg2->x) > MAX_EXP ||
                  ABS(reg3->x) > MAX_EXP )
      {
            reg1->x = MAX_EXP + 1;
            return 0;
      }

      /* For each digit in the first multiplier, shift the (negated)
         second multiplier right one position, multiply it by that
         digit (into wreg), and subtract it (wreg) from the
         (accumulating) result. */
      for( i = 0; i < control.round_size + control.integer_size + 1; ++i )
      {
            (void) shift(reg2, 1L, base); 

            (void) multbyd(&w, reg2, (int) r3[i], 0, base);

            (void) subber(reg1, &w, base);
      }
      reg1->x = exponent;
      reg1->inexact |= reg2->inexact | reg3->inexact;

      if( signflag < 0 )
            (void) negate(reg1, base);

      return normalize(reg1, base);
      
} /* end multiply */


/**   divide --- reg1 = reg2 / reg3
 *                                        \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
divide(int signflag, number *reg1, number *reg2, number *reg3, int positions, int base)
#else
divide(signflag, reg1, reg2, reg3, positions, base)
int signflag;
number *reg1;
number *reg2;
number *reg3;
int positions;
int base;
#endif
/*+~
 *    The positive value in |reg2 is divided by the value in |reg3, to
 *    give a quotient in |reg1, and a remainder in |reg2.  |Positions
 *    determines the number of digits for which the division proceeds.
 *    The quotient and the remainder are both given a sign determined
 *    by |signflag (1 = negative).  The value in |reg3 is destroyed.
 *
 *    The value 1 is returned if no errors occured in the division.
 **/
{
      int i;
      long exponent;

      (void) clear(reg1);
      exponent = reg2->x - reg3->x + 1;

      if( ABS(reg2->x) > MAX_EXP || ABS(exponent) > MAX_EXP ||
                  ABS(reg3->x) > MAX_EXP || is_zero(reg3, 0) )
      {
            reg1->x = MAX_EXP + 1;
            return 0;
      }

      (void) shift(reg2, 1L, UNSIGNED);
      for( i = 0; i <= positions; ++i )
      {
            while( is_notless(reg2, reg3) )
            {
                  ++reg1->r[i];
                  (void) subber(reg2, reg3, base);
            }
            (void) shift(reg3, 1L, UNSIGNED);
      }
      reg1->x = exponent;
      if( signflag < 0 )
      {
            (void) negate(reg1, base);
            (void) negate(reg2, base);
      }
      reg1->inexact |= reg2->inexact | reg3->inexact;

      (void) normalize(reg2, base);
      return normalize(reg1, base);

} /* end divide */


/**   multbyd --- |reg1 = |reg2 * d
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
multbyd(number *reg1, number *reg2, int d, int carry_in, int base)
#else
multbyd(reg1, reg2, d, carry_in, base)
number *reg1;
number *reg2;
int d;
int carry_in;
int base;
#endif
/*+~
 *    The contents of |reg2 is taken to be an unsigned number.
 *    The value returned in |reg1 is |d times the value in |reg2.
 *    |reg1 and |reg2 may be the same number.  If not, |reg2 is
 *    returned unchanged.  |Carry_in is the initial carry into the
 *    multiplication.
 *    
 *    The value returned by |multbyd is the value of the carry 
 *    from the top of |reg1.
 **/
{
      register int i = MIN(reg1->significand_size - 1,
                        reg2->significand_size - 1);
      register digit_t *r1 = &reg1->r[i];
      register digit_t *r2 = &reg2->r[i];
      register int carry = carry_in;

      /* The additional nested loops are an optimization validated
         by profiling the execution of this code.  */
      while( i-- >= 0 )
      {
            register int temp = (int) *r2-- * d + carry;
      
            /* If temp is small enough, there is no need to
               do the division and the remainder.  */
            if( temp < base )
            {
                  carry = 0;
                  *r1-- = (digit_t) temp;
                  
            } else
            {
                  carry = temp / base;
                  *r1-- = (digit_t)(temp % base);
            }
            
            /* If the next digit is zero, then only the 
               carry needs to be propogated.  */
            while(i >= 0 && !*r2 )
            {
                  --i, --r2;
                  temp = carry;
                  carry = temp / base;
                  *r1-- = (digit_t)(temp % base);
                  /* If the carry is zero, and the digits in
                     the register remain zero, then only zeros
                     need to be inserted in the result.  */             
                  if( !carry )
                  {
                        while( i >= 0 && !*r2 )
                        {
                              --i, --r2;
                              *r1-- = 0;
                        }
                  }
            }
      }

      return carry;

} /* end multbyd */


/**   divby2 --- |reg = |reg / 2
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
divby2(number *reg, int positions, int base)
#else
divby2(reg, positions, base)
number *reg;
int positions;
int base;
#endif
/*+~
 *    The given (signed) number is divided by 2.  The value return
 *    is the remainder of that division (i.e., is 1 if the number 
 *    being divided was odd).
 **/
{
      register int i = 0;
      register digit_t *r0 = reg->r;
      register int carry = ( SIGN(r0[0], base) < 0 ? 1 : 0 );
      
      for( ; i < positions; ++i )
      {
            register unsigned int temp = (unsigned int)
                        ((int) *r0 + carry * base);
            *r0++ = (digit_t)(temp >> 1);
            carry = temp & 1;
      }
      (void) normalize(reg, base);

      return carry;

} /* end divby2 */


/**   setvalue --- |reg = unit^(exponent - 1)
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
setvalue(number *reg, int unit, long exponent, int base)
#else
setvalue(reg, unit, exponent, base)
number *reg;
int unit;
long exponent;
int base;
#endif
/*+~
 *    The number |reg is cleared, its exponent set to the given
 *    |exponent. and the value |unit set into the digit 1 position
 *    of the number.  It is then normalized.
 **/
{
      (void) clear(reg);
      reg->x = exponent;
      reg->r[1] =  (digit_t) unit;
      return normalize(reg, base);
}


/**     loadval --- load a value to a number and exponent
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
loadval(number *reg, char *instr, int flag, int errorval, int base)
#else
loadval(reg, instr, flag, errorval, base)
number *reg;
char *instr;
int flag;
int errorval;
int base;
#endif
/*+~
 *      The string |in is scanned in the given |base for a number, and
 *      the significand is returned in the number |reg, represented
 *      in the same |base, and the exponent in the given variable |px
 *      as an integral value.  |Flag determines the handling of the
 *      sign: if positive, the number is returned correctly signed,
 *      and if negative, the number is negated.  If the |flag is zero,
 *      then the number returned as its absolute value.  In all cases,
 *      the value returned by |loadval is 0 if the value should be
 *      zero or positive, and 1 if the value should be negative.
 **/
{
      int j = 0;
      int point = 0;
      int signflag = 1;
      short separator_seen = 0;
      short digit_seen = 0;

      error = errorval;

      (void) clear(reg);
      if( !(in = instr) || base < 2 || base > 64 )
      {
            SETERROR(EINVAL);
            return 0;
      }

      /* If a set of signs is defined, and the character is one of
         the signs, then step over it.  The value is negative if the
         character being stepped over is the first defined sign.  */

      if( control.signs && *in && strchr(control.signs, *in) )
            signflag = ( *in++ == control.signs[0] ? -1 : 1);

      /* Ignore leading zeros after the sign.  Zero is represented
         by the zero entry in the digits string. */
      while( *in == control.digits[0] )
      {
            digit_seen = 1;
            ++in;
      }

      /* If the first significant character is the fractional part
         separator, then note that the separator has been seen,
         then count off the zeros.  */
      if( *in && strchr(control.separators, *in) )
      {
            separator_seen = 1;
            point = j;

            /* Move separator point left for each zero after the
               fractional separator. */
            while( *++in == control.digits[0] )
            {
                  digit_seen = 1;
                  --point;
            }
      }

      /* Now assemble the significant digits of the number.
         If the fractional separator has not been seen, then
         note where it is if it is seen.  Excess digits are
         ignored, and, if non-zero, represent accuracy loss. */
      for( ; *in && !strchr(control.exp_symbols, *in); ++in )
      {
            int v;

            if( strchr(control.separators, *in) )
                  if( !separator_seen )
                  {
                        point = j;
                        separator_seen = 1;
                        reg->inexact = (char)
                              control.inexact;
                        continue;
                  } else
                        break;

            if( (v = ctod(*in)) >= base || v < 0 )
                  break;

            digit_seen = 1;
            if( j < reg->significand_size )
                  reg->r[j] = (digit_t) v;
            else if( v != 0 )
                  reg->inexact = 1;
            ++j;
      }

      /* Insert an implicit point at end if one not seen */
      if( !separator_seen ) point = j;

      /* If there is an exponent, pick it up, and then adjust it
         by where the point was seen  */
      if( *in && strchr(control.exp_symbols, *in) )
      {
            int v;
            int esignflag = 1;
            int expbase = ( control.exp_base ?
                        control.exp_base : base );

            ++in;
            if( control.signs && *in && strchr(control.signs, *in) )
                  esignflag = ( *in++ == control.signs[0] ? -1 : 1 );
            for( ; *in; ++in )
            {
                  if( (v = ctod(*in)) >= expbase || v < 0 )
                        break;

                  /* If there is an exponent, and no digits seen,
                     then create the implied `1'. */
                  if( !digit_seen )
                  {
                        reg->r[j++] = 1;
                        ++point;
                  }
                  digit_seen = 1;
                  if( (reg->x = reg->x * (long) expbase +
                              (long) v) > MAX_EXP )
                  {
                        SETERROR(ERANGE);
                        break;
                  }
            }
            if( esignflag < 0 )
                  reg->x = -reg->x;
      }

      if( !digit_seen || (control.check_digits && *in &&
                  (strchr(control.digits, *in) ||
                   strchr(control.exp_symbols, *in) ||
                   strchr(control.separators, *in))) )
            SETERROR(EINVAL);

      reg->x += (long)(point - 1);

      if( SIGN(reg->r[0], base) <= 0 )
            (void) shift(reg, 1L, UNSIGNED);

      if( (flag > 0 && signflag < 0) || (flag < 0 && signflag > 0 ) )
            (void) negate(reg, base);

      /* Normalize the result */
      if( !normalize(reg, base) )
            SETERROR(EDOM);

      return signflag;

} /* end loadval */


/**     unloadval --- unload the value from a number.
 *                                              \+1*/
static char *
#if defined(__STDC__) || defined(__cplusplus)
unloadval(number *reg, char *symbols, int base)
#else
unloadval(reg, symbols, base)
number *reg;
char *symbols;
int base;
#endif
/*+~
 *      A pointer to a (static) string is returned, which is the
 *      representation of the number in the number, |reg, whose |base
 *      is as given, and with exponent |*px.
 *
 **/
{
      int lastplace;

      int initialplace = 0;
      char *out = rstr;
      char exp_symbol = control.exp_symbols[0];
      char fract_symbol = control.separators[0];
      char minus_symbol = control.signs[0];
      char plus_symbol = '\0';

      /* Assume result will be in floating form i.e. with exponent. */
      int floating = 1;

      /* Assume fractional part is optional. */
      int fract_part = control.inexact & reg->inexact;

      /* Are there substitute characters to use as punctuation in the
         output string.  If yes, then make the substitions, and
         insist that a fractional part is produced in the output.  */
      if( symbols && symbols[0] )
      {
            fract_part = 1;
            exp_symbol = symbols[0];
            if( symbols[1] )
            {
                  fract_symbol = symbols[1];
                  if( symbols[2] )
                  {
                        minus_symbol = symbols[2];
                        if( symbols[3] )
                              plus_symbol = symbols[3];
                  }

            }

      }

      (void) round(reg, control.rounding, control.integer_size +
                  control.round_size, base);


      if( error || ABS(reg->x) > MAX_EXP )
      {
            SETERROR(ERANGE);
            errno = error;
            return (char *) 0;
      }


      if( SIGN(reg->r[0], base) < 0 )
      {
            *out++ = minus_symbol;
            (void) negate(reg, base);
      } else if( plus_symbol )
            *out++ = plus_symbol;

      /* If the point place is in the number (i.e., the
         exponent is in the range [0..integer_size]), and
         there are no digits after the point, then convert as
         an integer.  */
      if( !fract_part && is_integer(reg) )
            floating = 0;


      /* If there is no initial digit (i.e., it is zero) and
         the next digit is significant (i.e., it is not zero)
         then use that digit as the starting place, and
         decrement the exponent by one to compensate.  */
      if( !reg->r[0] && reg->r[1] )
      {
            ++initialplace;
            --reg->x;
      }

      if( floating )
      {
            /* Special case of exponent = -1: this may be
               converted in the denormalized form
               "0.xyz..", with no exponent value, if
               requested.  */
            if( reg->x == -1 && control.denormalize )
            {
                  *out++ = dtoc(0);
                  reg->x = 0;
            } else
                  /* ... else convert in the form "x.yz.."
                     with an exponent if needed. */
                  *out++ = dtoc((int) reg->r[initialplace++]);

            /* Scan backwards to find the last significant digit
               to include in the conversion. The first digit has
               already been generated.  */
            lastplace = control.round_size + control.integer_size - 1;
            while( !reg->r[lastplace] && --lastplace >= initialplace );


            /* If fractional part not going to be produced, but
               is required, make it happen.  */

            if( initialplace > lastplace && fract_part )
                  lastplace = initialplace;

            /* Only if there is something after the point should
               the point be included. */
            if( initialplace <= lastplace )
                  *out++ = fract_symbol;
      }
      else /* not floating */
            lastplace =  initialplace + (int) reg->x;

      /* Generate the rest of the digits.  */
      while( initialplace <= lastplace )
            *out++ = dtoc((int) reg->r[initialplace++]);

      *out = '\0';            /* Terminate the string.  */

      /* If floating format, and the exponent is non-zero,
         then add the exponent to the string.  */
      if( floating && reg->x )
      {
            *out++ = exp_symbol;
            (void) strcpy(out, ltoab((long) reg->x, 0,
                  ( control.exp_base ? control.exp_base : base )));
      }

      return rstr;

} /* end unloadval */


/**     ltoab --- long integer to ascii, with base.
 *                                              \+1*/
static char *
#if defined(__STDC__) || defined(__cplusplus)
ltoab(long i, int digits, int base)
#else
ltoab(i, digits, base)
long i;
int digits;
int base;
#endif
/*+~
 *      Generate a string of digits in the given base |base,
 *    representing the number |i, containing a minimum of |digits
 *    digits, or else the number of digits needed to represent
 *    |i.  The returned string is allocated statically, so must
 *    be saved if needed.
 **/
{
      static char result[sizeof(i) * 8 + 2];
      char *p = &result[sizeof(result)/sizeof(result[0]) - 1];
      /* |signflag is the sign of the result. */
      int signflag = ( i >= 0 ? 0 : (i = -i, 1) );

      *p = '\0';
      do
      {
             *--p = dtoc((int)(i % base));
            i /= base;
            --digits;
      } while( i || digits > 0);

      if( signflag ) *--p = '-';

      return p;

} /* end ltoab */


/**     subber --- subtract numbers.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
subber(number *reg1, number *reg2, int base)
#else
subber(reg1, reg2, base)
number *reg1;
number *reg2;
int base;
#endif
/*+~
 *      Subtract the |number |reg2 from the |number |reg1 in the
 *      given |base, assuming that both numbers have equal exponents.
 *    If the result overflows, |reg1 will be shifted right one place,
 *    and its exponent incremented by 1.
 **/
{
      register int i = MIN(reg1->significand_size - 1, 
                        reg2->significand_size - 1);
      register digit_t *r1 = &reg1->r[i];
      register digit_t *r2 = &reg2->r[i];
      register int borrow = 0;
      int sr1 = SIGN(reg1->r[0], base);     /* Initial sign of reg1 */

      /* To perform the subtraction, the sum of the reg2 digit
         and the previous borrow (initially zero) are subtracted
         from the reg1 digit.  If this digit goes negative, it is
         increased by the value of the base, and borrow set,
         otherwise the borrow is zero. */

      for( ; i >= 0; --i )
      {
            register int v = (int) *r1;
            if( (v -= (int) *r2-- + borrow) < 0 )
            {
                  *r1-- = (digit_t)(v + base);
                  borrow = 1;
            } else
            {
                  *r1-- = (digit_t) v;
                  borrow = 0;
            }
      }

      /* If reg1 and reg2 are initially the same sign, then the
         subtraction cannot overflow.  If they are different, then
         overflow has occurred if the result does not have the
         same sign as reg1.  [The reason for this latter is that
         if reg1 is initially negative, then reg2 must be positive
         when it does not have the same sign, and so the
         subtraction of reg2 should make reg1 more negative, i.e.,
         is should have the same sign as it started with.  If it
         doesn't, it has overflowed.  Similarly for the case of
         reg1 being initially positive. */

      if( sr1 != SIGN(reg2->r[0], base) &&
                  sr1 != SIGN(reg1->r[0], base) )
      {
            /* The shift sign-extends.  But the sign is known
               to be wrong, even though the remainder of the
               digits are correct, so the sign is corrected by
               making it the same as the initial value of reg1. */

            if( shift(reg1, 1L, UNSIGNED) )
                  reg1->r[0] = (digit_t) ( sr1 > 0 ? 0 : base - 1 );
            else
                  return 0;
      }
      return 1;

} /* end subber */


/**     clear --- clear number to zeroes.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
clear(number *reg)
#else
clear(reg)
number *reg;
#endif
/*+~
 *      The number |reg is cleared to all zeroes, and its exponent
 *    set to zero.
 **/

{
      register int i;
      register digit_t *r0 = reg->r;

      for( i = reg->significand_size - 1; i >= 0; --i )
            *r0++ = 0;

      reg->x = 0L;
      reg->inexact = '\0';

      return 1;

} /* end clear */


/**     negate --- negate a number.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
negate(number *reg, int base)
#else
negate(reg, base)
number *reg;
int base;
#endif
/*+~
 *      The given number, |reg, is subtracted from zero, in the given
 *      |base, to negate it.  If the result overflows, |reg will
 *      be shifted right one place. 
 **/
{
      int rbase = base - 1;
      int sr0 = SIGN(reg->r[0], base);     /* Initial sign of reg */
      register int i= reg->significand_size - 1;
      register digit_t *r0 = &reg->r[i];

      /* To perform the negation, trailing zeros are ignored. The
         first significant digit is subtracted from base to find
         its negation, or, as is done here, that digit decremented
         will be subtracted from |rbase (= |base - 1). After this
         digit is seen, the (virtual) borrow will always be set,
         so the subtraction is effectively from |base - 1 (as a
         borrow will always be carried forward). */

      for(; i >= 0 && !*r0; --r0, --i );

      if( i >= 0 )
            --*r0;

      for( ; i >= 0; --r0, --i )
            *r0 = (digit_t)(rbase - (int) *r0);

      /* There is overflow if the sign of the original value was
         negative and the sign of the result is also negative. */

      if( sr0 < 0 && SIGN(reg->r[0], base) < 0 )
      {
            /* The shift is done without sign extension, as the
               sign is known to be wrong, even though the remainder
               of the digits are correct, so the sign is corrected
               by making it +ve. */

            if( !shift(reg, 1L, UNSIGNED) )
                  return 0;
      }
      return 1;

} /* end negate */


/**     copy --- copy one number to another.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
copy(number *reg1, number *reg2)
#else
copy(reg1, reg2)
number *reg1, *reg2;
#endif
/*+~
 *      Number |reg2 is copied to |reg1.  If reg2 is longer, it
 *    is filled with zeroes, if shorter, the number is truncated to
 *    fit.
 **/
{
      register int i = MIN(reg1->significand_size - 1,
                              reg2->significand_size - 1);
      register digit_t *r1 = reg1->r;
      register digit_t *r2 = reg2->r;

      (void) clear(reg1);
      for(; i >= 0; --i )
            *r1++ = *r2++;
      reg1->x = reg2->x;
      reg1->inexact = reg2->inexact;
      
      return 1;
}


/**     is_notless --- compare magnitude of positive significands.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
is_notless(number *reg1, number *reg2)
#else
is_notless(reg1, reg2)
number *reg1, *reg2;
#endif
/*+~
 *      If the significand of |reg1 is not less (i.e., greater than or
 *    equal) to that of the significand of |reg2, return 1, otherwise
 *    return 0.
 *
 *      Both |reg1 and |reg2 are required to contain positive
 *    significands for correct results.
 **/
{
      register int i = 0;
      int limit = MIN(reg1->significand_size - 1,
                              reg2->significand_size - 1);
      register digit_t *r1 = reg1->r;
      register digit_t *r2 = reg2->r;

      while( i < limit && *r1 == *r2 ) ++r1, ++r2, ++i;

      return ( i == limit + 1 ? 1 : *r1 >= *r2 );

} /* end is_notless */


/**     is_integer --- is the number an integer.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
is_integer(number *reg)
#else
is_integer(reg)
number *reg;
#endif
/*+~
 *      A number is an integer if it is not inexact, has its fractional
 *    point in the number between the zeroth digit and the 
 *    |integer_size digit inclusive, and is zero to the right of the
 *    fractional point.
 **/
{
      return !reg->inexact &&
                  reg->x >= 0  && reg->x <= control.integer_size
                   && is_zero(reg, (int) reg->x + 1);
} /* end is_integer */


/**     shift --- right shift number.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
shift(number *reg, long dist, int base)
#else
shift(reg, dist, base)
number *reg;
long dist;
int base;
#endif
/*+~
 *      Shift the number |reg to the right, by |dist places.  |Dist
 *    must not be negative.  Movement beyond the end of the number
 *    clears the number.  If digits fall off the right hand end,
 *    |inexact is set, and the value is truncated towards
 *    negative infinity.  If UNSIGNED shift is requested (i.e., a
 *    base value = 0, then the value is zero filled on the left,
 *    otherwise the sign is extended at the left hand end of the
 *    number.  The exponent of the number is increased by |dist,
 *    or is set to zero if the number clears of all digits.
 **/
{
      register int i;
      register digit_t *r0;
      register digit_t *rd;
      int places;
      digit_t sgn = (digit_t) ( !base ? 0 :
                  ( SIGN(reg->r[0], base) < 0 ? base - 1 : 0 ));

      /* Is the shift going to move out all the digits.  If
         so, just set that many places to clear.  The exponent
         is biased so as to become zero.  Otherwise, the number
         of places to move is that given.  */
      if( dist >= reg->significand_size )
      {
            places = reg->significand_size;
            sgn = 0;
            reg->x = -dist;
      } else
            places = (int) dist;

      /* If real digits are moved out, set inexact.  */
      if( !is_zero(reg, (int) reg->significand_size - places) )
            reg->inexact = 1;

      /* Move the digits, and replace the ones moved with copies
         of the sign (set to 0 is a zero base was given. */
         
      r0 = &reg->r[(i = reg->significand_size - places - 1)];
      rd = &reg->r[i + places];
      for( ; i >= 0; --i )
            *rd-- = *r0--;

      rd = &reg->r[(i = places - 1)];
      for( ; i >= 0; --i )
            *rd-- = sgn;

      /* Adjust the exponent. */
      reg->x += dist;

      return 1;

} /* end shift */


/**     is_zero --- test if a part of a number has a value of zero.
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
is_zero(number *reg, int position)
#else
is_zero(reg, position)
number *reg;
int position;
#endif
/*+~
 *      The given number |reg is checked to see if all digits
 *    from the given position in its significand, to the end of
 *    the significand, are all zero: if so, the value 1 is
 *    returned, otherwise 0.
 **/
{
      int i = position;
      int limit = reg->significand_size;
      register digit_t *r0 = &reg->r[i];

      for( ; i < limit; ++i )
            if( *r0++ ) return 0;
      return 1;

} /* end is_zero */


/**     ctod --- character to digit value
 *                                              \+1*/
static int
#if defined(__STDC__) || defined(__cplusplus)
ctod(int c)
#else
ctod(c)
int c;
#endif
/*+~
 *      The given character is converted to the integer value
 *      it represents: the character must be a character in
 *    the |control.|digits string: otherwise a negative error
 *    indication is returned.
 **/
{
      int d, k, limit = (int) strlen(control.digits);

      if( control.digits )
      {
            if( control.ignore_case && islower(c) )
                  c = toupper(c);

            for( k = 0; k < limit; ++k )
            {
                  d = (int) control.digits[k];
                  if( control.ignore_case && islower(d) )
                        d = toupper(d);
                  if( c == d )
                        return k;
            }
      }

      return -1;

} /* end ctod */


/**     dtoc --- digit to character value
 *                                              \+1*/
static char
#if defined(__STDC__) || defined(__cplusplus)
dtoc(int i)
#else
dtoc(i)
int i;
#endif
/*+~
 *      The given integer value |i is converted to the character
 *      that represents it: the value must be in the range
 *      0..63: otherwise the routine may misfunction with no
 *      indication.  |dtoc is the inverse function to |ctod.
 **/
{
      return control.digits[i];

} /* end dtoc */

Generated by  Doxygen 1.6.0   Back to index