lib/bitvector/BitVector.c

   1 #ifndef MODULE_BIT_VECTOR
   2 #define MODULE_BIT_VECTOR
   3 /*****************************************************************************/
   4 /*  MODULE NAME:  BitVector.c                           MODULE TYPE:  (adt)  */
   5 /*****************************************************************************/
   6 /*  MODULE IMPORTS:                                                          */
   7 /*****************************************************************************/
   8 #include <stdlib.h>                                 /*  MODULE TYPE:  (sys)  */
   9 #include <limits.h>                                 /*  MODULE TYPE:  (sys)  */
  10 #include <string.h>                                 /*  MODULE TYPE:  (sys)  */
  11 #include <ctype.h>                                  /*  MODULE TYPE:  (sys)  */
  12 #include "ToolBox.h"                                /*  MODULE TYPE:  (dat)  */
  13 /*****************************************************************************/
  14 /*  MODULE INTERFACE:                                                        */
  15 /*****************************************************************************/
  16
  17 typedef enum
  18     {
  19         ErrCode_Ok = 0,    /* everything went allright                       */
  20
  21         ErrCode_Type,      /* types word and size_t have incompatible sizes  */
  22         ErrCode_Bits,      /* bits of word and sizeof(word) are inconsistent */
  23         ErrCode_Word,      /* size of word is less than 16 bits              */
  24         ErrCode_Long,      /* size of word is greater than size of long      */
  25         ErrCode_Powr,      /* number of bits of word is not a power of two   */
  26         ErrCode_Loga,      /* error in calculation of logarithm              */
  27
  28         ErrCode_Null,      /* unable to allocate memory                      */
  29
  30         ErrCode_Indx,      /* index out of range                             */
  31         ErrCode_Ordr,      /* minimum > maximum index                        */
  32         ErrCode_Size,      /* bit vector size mismatch                       */
  33         ErrCode_Pars,      /* input string syntax error                      */
  34         ErrCode_Ovfl,      /* numeric overflow error                         */
  35         ErrCode_Same,      /* operands must be distinct                      */
  36         ErrCode_Expo,      /* exponent must be positive                      */
  37         ErrCode_Zero       /* division by zero error                         */
  38     } ErrCode;
  39
  40 typedef wordptr *listptr;
  41
  42 /* ===> MISCELLANEOUS BASIC FUNCTIONS: <=== */
  43
  44 charptr BitVector_Error      (ErrCode error);  /* return string for err code */
  45
  46 ErrCode BitVector_Boot       (void);                 /* 0 = ok, 1..7 = error */
  47
  48 N_word  BitVector_Size       (N_int bits);  /* bit vector size (# of words)  */
  49 N_word  BitVector_Mask       (N_int bits);  /* bit vector mask (unused bits) */
  50
  51 /* ===> CLASS METHODS: <=== */
  52
  53 charptr BitVector_Version    (void);                /* return version string */
  54
  55 N_int   BitVector_Word_Bits  (void);     /* return # of bits in machine word */
  56 N_int   BitVector_Long_Bits  (void);    /* return # of bits in unsigned long */
  57
  58 /* ===> CONSTRUCTOR METHODS: <=== */
  59
  60 wordptr BitVector_Create     (N_int bits, boolean clear);          /* malloc */
  61 listptr BitVector_Create_List(N_int bits, boolean clear, N_int count);
  62
  63 wordptr BitVector_Resize     (wordptr oldaddr, N_int bits);       /* realloc */
  64
  65 wordptr BitVector_Shadow     (wordptr addr); /* make new same size but empty */
  66 wordptr BitVector_Clone      (wordptr addr);         /* make exact duplicate */
  67
  68 wordptr BitVector_Concat     (wordptr X, wordptr Y); /* return concatenation */
  69
  70 /* ===> DESTRUCTOR METHODS: <=== */
  71
  72 void    BitVector_Dispose            (charptr string);             /* string */
  73 void    BitVector_Destroy            (wordptr addr);               /* bitvec */
  74 void    BitVector_Destroy_List       (listptr list, N_int count);  /* list   */
  75
  76 /* ===> OBJECT METHODS: <=== */
  77
  78 /* ===> bit vector copy function: */
  79
  80 void    BitVector_Copy       (wordptr X, wordptr Y);              /* X = Y   */
  81
  82 /* ===> bit vector initialization: */
  83
  84 void    BitVector_Empty      (wordptr addr);                      /* X = {}  */
  85 void    BitVector_Fill       (wordptr addr);                      /* X = ~{} */
  86 void    BitVector_Flip       (wordptr addr);                      /* X = ~X  */
  87
  88 void    BitVector_Primes     (wordptr addr);
  89
  90 /* ===> miscellaneous functions: */
  91
  92 void    BitVector_Reverse    (wordptr X, wordptr Y);
  93
  94 /* ===> bit vector interval operations and functions: */
  95
  96 void    BitVector_Interval_Empty     (wordptr addr, N_int lower, N_int upper);
  97 void    BitVector_Interval_Fill      (wordptr addr, N_int lower, N_int upper);
  98 void    BitVector_Interval_Flip      (wordptr addr, N_int lower, N_int upper);
  99 void    BitVector_Interval_Reverse   (wordptr addr, N_int lower, N_int upper);
 100
 101 boolean BitVector_interval_scan_inc  (wordptr addr, N_int start,
 102                                       N_intptr min, N_intptr max);
 103 boolean BitVector_interval_scan_dec  (wordptr addr, N_int start,
 104                                       N_intptr min, N_intptr max);
 105
 106 void    BitVector_Interval_Copy      (wordptr X, wordptr Y, N_int Xoffset,
 107                                       N_int Yoffset, N_int length);
 108
 109 wordptr BitVector_Interval_Substitute(wordptr X, wordptr Y,
 110                                       N_int Xoffset, N_int Xlength,
 111                                       N_int Yoffset, N_int Ylength);
 112
 113 /* ===> bit vector test functions: */
 114
 115 boolean BitVector_is_empty   (wordptr addr);                  /* X == {} ?   */
 116 boolean BitVector_is_full    (wordptr addr);                  /* X == ~{} ?  */
 117
 118 boolean BitVector_equal      (wordptr X, wordptr Y);          /* X == Y ?    */
 119 Z_int   BitVector_Lexicompare(wordptr X, wordptr Y);          /* X <,=,> Y ? */
 120 Z_int   BitVector_Compare    (wordptr X, wordptr Y);          /* X <,=,> Y ? */
 121
 122 /* ===> bit vector string conversion functions: */
 123
 124 charptr BitVector_to_Hex     (wordptr addr);
 125 ErrCode BitVector_from_Hex   (wordptr addr, charptr string);
 126
 127 charptr BitVector_to_Bin     (wordptr addr);
 128 ErrCode BitVector_from_Bin   (wordptr addr, charptr string);
 129
 130 charptr BitVector_to_Dec     (wordptr addr);
 131 ErrCode BitVector_from_Dec   (wordptr addr, charptr string);
 132
 133 charptr BitVector_to_Enum    (wordptr addr);
 134 ErrCode BitVector_from_Enum  (wordptr addr, charptr string);
 135
 136 /* ===> bit vector bit operations, functions & tests: */
 137
 138 void    BitVector_Bit_Off    (wordptr addr, N_int index); /*  X = X \ {x}    */
 139 void    BitVector_Bit_On     (wordptr addr, N_int index); /*  X = X + {x}    */
 140 boolean BitVector_bit_flip   (wordptr addr, N_int index); /* (X+{x})\(X*{x}) */
 141
 142 boolean BitVector_bit_test   (wordptr addr, N_int index); /*  {x} in X ?     */
 143
 144 void    BitVector_Bit_Copy   (wordptr addr, N_int index, boolean bit);
 145
 146 /* ===> bit vector bit shift & rotate functions: */
 147
 148 void    BitVector_LSB                (wordptr addr, boolean bit);
 149 void    BitVector_MSB                (wordptr addr, boolean bit);
 150 boolean BitVector_lsb_               (wordptr addr);
 151 boolean BitVector_msb_               (wordptr addr);
 152 boolean BitVector_rotate_left        (wordptr addr);
 153 boolean BitVector_rotate_right       (wordptr addr);
 154 boolean BitVector_shift_left         (wordptr addr, boolean carry_in);
 155 boolean BitVector_shift_right        (wordptr addr, boolean carry_in);
 156 void    BitVector_Move_Left          (wordptr addr, N_int bits);
 157 void    BitVector_Move_Right         (wordptr addr, N_int bits);
 158
 159 /* ===> bit vector insert/delete bits: */
 160
 161 void    BitVector_Insert     (wordptr addr, N_int offset, N_int count,
 162                               boolean clear);
 163 void    BitVector_Delete     (wordptr addr, N_int offset, N_int count,
 164                               boolean clear);
 165
 166 /* ===> bit vector arithmetic: */
 167
 168 boolean BitVector_increment  (wordptr addr);                        /*  X++  */
 169 boolean BitVector_decrement  (wordptr addr);                        /*  X--  */
 170
 171 boolean BitVector_compute    (wordptr X, wordptr Y, wordptr Z, boolean minus,
 172                                                                boolean *carry);
 173 boolean BitVector_add        (wordptr X, wordptr Y, wordptr Z, boolean *carry);
 174 boolean BitVector_sub        (wordptr X, wordptr Y, wordptr Z, boolean *carry);
 175 boolean BitVector_inc        (wordptr X, wordptr Y);
 176 boolean BitVector_dec        (wordptr X, wordptr Y);
 177
 178 void    BitVector_Negate     (wordptr X, wordptr Y);
 179 void    BitVector_Absolute   (wordptr X, wordptr Y);
 180 Z_int   BitVector_Sign       (wordptr addr);
 181 ErrCode BitVector_Mul_Pos    (wordptr X, wordptr Y, wordptr Z, boolean strict);
 182 ErrCode BitVector_Multiply   (wordptr X, wordptr Y, wordptr Z);
 183 ErrCode BitVector_Div_Pos    (wordptr Q, wordptr X, wordptr Y, wordptr R);
 184 ErrCode BitVector_Divide     (wordptr Q, wordptr X, wordptr Y, wordptr R);
 185 ErrCode BitVector_GCD        (wordptr X, wordptr Y, wordptr Z);
 186 ErrCode BitVector_GCD2       (wordptr U, wordptr V, wordptr W,      /*   O   */
 187                                          wordptr X, wordptr Y);     /*   I   */
 188 ErrCode BitVector_Power      (wordptr X, wordptr Y, wordptr Z);
 189
 190 /* ===> direct memory access functions: */
 191
 192 void    BitVector_Block_Store(wordptr addr, charptr buffer, N_int length);
 193 charptr BitVector_Block_Read (wordptr addr, N_intptr length);
 194
 195 /* ===> word array functions: */
 196
 197 void    BitVector_Word_Store (wordptr addr, N_int offset, N_int value);
 198 N_int   BitVector_Word_Read  (wordptr addr, N_int offset);
 199
 200 void    BitVector_Word_Insert(wordptr addr, N_int offset, N_int count,
 201                               boolean clear);
 202 void    BitVector_Word_Delete(wordptr addr, N_int offset, N_int count,
 203                               boolean clear);
 204
 205 /* ===> arbitrary size chunk functions: */
 206
 207 void    BitVector_Chunk_Store(wordptr addr, N_int chunksize,
 208                               N_int offset, N_long value);
 209 N_long  BitVector_Chunk_Read (wordptr addr, N_int chunksize,
 210                               N_int offset);
 211
 212 /* ===> set operations: */
 213
 214 void    Set_Union            (wordptr X, wordptr Y, wordptr Z); /* X = Y + Z */
 215 void    Set_Intersection     (wordptr X, wordptr Y, wordptr Z); /* X = Y * Z */
 216 void    Set_Difference       (wordptr X, wordptr Y, wordptr Z); /* X = Y \ Z */
 217 void    Set_ExclusiveOr      (wordptr X, wordptr Y, wordptr Z); /*(Y+Z)\(Y*Z)*/
 218 void    Set_Complement       (wordptr X, wordptr Y);            /* X = ~Y    */
 219
 220 /* ===> set functions: */
 221
 222 boolean Set_subset           (wordptr X, wordptr Y);            /* X in Y ?  */
 223
 224 N_int   Set_Norm             (wordptr addr);                    /* = | X |   */
 225 N_int   Set_Norm2            (wordptr addr);                    /* = | X |   */
 226 N_int   Set_Norm3            (wordptr addr);                    /* = | X |   */
 227 Z_long  Set_Min              (wordptr addr);                    /* = min(X)  */
 228 Z_long  Set_Max              (wordptr addr);                    /* = max(X)  */
 229
 230 /* ===> matrix-of-booleans operations: */
 231
 232 void    Matrix_Multiplication(wordptr X, N_int rowsX, N_int colsX,
 233                               wordptr Y, N_int rowsY, N_int colsY,
 234                               wordptr Z, N_int rowsZ, N_int colsZ);
 235
 236 void    Matrix_Product       (wordptr X, N_int rowsX, N_int colsX,
 237                               wordptr Y, N_int rowsY, N_int colsY,
 238                               wordptr Z, N_int rowsZ, N_int colsZ);
 239
 240 void    Matrix_Closure       (wordptr addr, N_int rows, N_int cols);
 241
 242 void    Matrix_Transpose     (wordptr X, N_int rowsX, N_int colsX,
 243                               wordptr Y, N_int rowsY, N_int colsY);
 244
 245 /*****************************************************************************/
 246 /*  MODULE RESOURCES:                                                        */
 247 /*****************************************************************************/
 248
 249 #define bits_(BitVector) *(BitVector-3)
 250 #define size_(BitVector) *(BitVector-2)
 251 #define mask_(BitVector) *(BitVector-1)
 252
 253 #define  ERRCODE_TYPE  "sizeof(word) > sizeof(size_t)"
 254 #define  ERRCODE_BITS  "bits(word) != sizeof(word)*8"
 255 #define  ERRCODE_WORD  "bits(word) < 16"
 256 #define  ERRCODE_LONG  "bits(word) > bits(long)"
 257 #define  ERRCODE_POWR  "bits(word) != 2^x"
 258 #define  ERRCODE_LOGA  "bits(word) != 2^ld(bits(word))"
 259 #define  ERRCODE_NULL  "unable to allocate memory"
 260 #define  ERRCODE_INDX  "index out of range"
 261 #define  ERRCODE_ORDR  "minimum > maximum index"
 262 #define  ERRCODE_SIZE  "bit vector size mismatch"
 263 #define  ERRCODE_PARS  "input string syntax error"
 264 #define  ERRCODE_OVFL  "numeric overflow error"
 265 #define  ERRCODE_SAME  "result vector(s) must be distinct"
 266 #define  ERRCODE_EXPO  "exponent must be positive"
 267 #define  ERRCODE_ZERO  "division by zero error"
 268 #define  ERRCODE_OOPS  "unexpected internal error - please contact author"
 269
 270 const N_int BitVector_BYTENORM[256] =
 271 {
 272     0x00, 0x01, 0x01, 0x02,  0x01, 0x02, 0x02, 0x03,
 273     0x01, 0x02, 0x02, 0x03,  0x02, 0x03, 0x03, 0x04, /* 0x00 */
 274     0x01, 0x02, 0x02, 0x03,  0x02, 0x03, 0x03, 0x04,
 275     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05, /* 0x10 */
 276     0x01, 0x02, 0x02, 0x03,  0x02, 0x03, 0x03, 0x04,
 277     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05, /* 0x20 */
 278     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05,
 279     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06, /* 0x30 */
 280     0x01, 0x02, 0x02, 0x03,  0x02, 0x03, 0x03, 0x04,
 281     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05, /* 0x40 */
 282     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05,
 283     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06, /* 0x50 */
 284     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05,
 285     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06, /* 0x60 */
 286     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06,
 287     0x04, 0x05, 0x05, 0x06,  0x05, 0x06, 0x06, 0x07, /* 0x70 */
 288     0x01, 0x02, 0x02, 0x03,  0x02, 0x03, 0x03, 0x04,
 289     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05, /* 0x80 */
 290     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05,
 291     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06, /* 0x90 */
 292     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05,
 293     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06, /* 0xA0 */
 294     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06,
 295     0x04, 0x05, 0x05, 0x06,  0x05, 0x06, 0x06, 0x07, /* 0xB0 */
 296     0x02, 0x03, 0x03, 0x04,  0x03, 0x04, 0x04, 0x05,
 297     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06, /* 0xC0 */
 298     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06,
 299     0x04, 0x05, 0x05, 0x06,  0x05, 0x06, 0x06, 0x07, /* 0xD0 */
 300     0x03, 0x04, 0x04, 0x05,  0x04, 0x05, 0x05, 0x06,
 301     0x04, 0x05, 0x05, 0x06,  0x05, 0x06, 0x06, 0x07, /* 0xE0 */
 302     0x04, 0x05, 0x05, 0x06,  0x05, 0x06, 0x06, 0x07,
 303     0x05, 0x06, 0x06, 0x07,  0x06, 0x07, 0x07, 0x08  /* 0xF0 */
 304 };
 305
 306 /*****************************************************************************/
 307 /*  MODULE IMPLEMENTATION:                                                   */
 308 /*****************************************************************************/
 309
 310     /**********************************************/
 311     /* global implementation-intrinsic constants: */
 312     /**********************************************/
 313
 314 #define BIT_VECTOR_HIDDEN_WORDS 3
 315
 316     /*****************************************************************/
 317     /* global machine-dependent constants (set by "BitVector_Boot"): */
 318     /*****************************************************************/
 319
 320 static N_word BITS;     /* = # of bits in machine word (must be power of 2)  */
 321 static N_word MODMASK;  /* = BITS - 1 (mask for calculating modulo BITS)     */
 322 static N_word LOGBITS;  /* = ld(BITS) (logarithmus dualis)                   */
 323 static N_word FACTOR;   /* = ld(BITS / 8) (ld of # of bytes)                 */
 324
 325 static N_word LSB = 1;  /* = mask for least significant bit                  */
 326 static N_word MSB;      /* = mask for most significant bit                   */
 327
 328 static N_word LONGBITS; /* = # of bits in unsigned long                      */
 329
 330 static N_word LOG10;    /* = logarithm to base 10 of BITS - 1                */
 331 static N_word EXP10;    /* = largest possible power of 10 in signed int      */
 332
 333     /********************************************************************/
 334     /* global bit mask table for fast access (set by "BitVector_Boot"): */
 335     /********************************************************************/
 336
 337 static wordptr BITMASKTAB;
 338
 339     /*****************************/
 340     /* global macro definitions: */
 341     /*****************************/
 342
 343 #define BIT_VECTOR_ZERO_WORDS(target,count) \
 344     while (count-- > 0) *target++ = 0;
 345
 346 #define BIT_VECTOR_FILL_WORDS(target,fill,count) \
 347     while (count-- > 0) *target++ = fill;
 348
 349 #define BIT_VECTOR_FLIP_WORDS(target,flip,count) \
 350     while (count-- > 0) *target++ ^= flip;
 351
 352 #define BIT_VECTOR_COPY_WORDS(target,source,count) \
 353     while (count-- > 0) *target++ = *source++;
 354
 355 #define BIT_VECTOR_BACK_WORDS(target,source,count) \
 356     { target += count; source += count; while (count-- > 0) *--target = *--source; }
 357
 358 #define BIT_VECTOR_CLR_BIT(address,index) \
 359     *(address+(index>>LOGBITS)) &= NOT BITMASKTAB[index AND MODMASK];
 360
 361 #define BIT_VECTOR_SET_BIT(address,index) \
 362     *(address+(index>>LOGBITS)) |= BITMASKTAB[index AND MODMASK];
 363
 364 #define BIT_VECTOR_TST_BIT(address,index) \
 365     ((*(address+(index>>LOGBITS)) AND BITMASKTAB[index AND MODMASK]) != 0)
 366
 367 #define BIT_VECTOR_FLP_BIT(address,index,mask) \
 368     (mask = BITMASKTAB[index AND MODMASK]), \
 369     (((*(addr+(index>>LOGBITS)) ^= mask) AND mask) != 0)
 370
 371 #define BIT_VECTOR_DIGITIZE(type,value,digit) \
 372     value = (type) ((digit = value) / 10); \
 373     digit -= value * 10; \
 374     digit += (type) '0';
 375
 376     /*********************************************************/
 377     /* private low-level functions (potentially dangerous!): */
 378     /*********************************************************/
 379
 380 static N_word power10(N_word x)
 381 {
 382     N_word y = 1;
 383
 384     while (x-- > 0) y *= 10;
 385     return(y);
 386 }
 387
 388 static void BIT_VECTOR_zro_words(wordptr addr, N_word count)
 389 {
 390     BIT_VECTOR_ZERO_WORDS(addr,count)
 391 }
 392
 393 static void BIT_VECTOR_cpy_words(wordptr target, wordptr source, N_word count)
 394 {
 395     BIT_VECTOR_COPY_WORDS(target,source,count)
 396 }
 397
 398 static void BIT_VECTOR_mov_words(wordptr target, wordptr source, N_word count)
 399 {
 400     if (target != source)
 401     {
 402         if (target < source) BIT_VECTOR_COPY_WORDS(target,source,count)
 403         else                 BIT_VECTOR_BACK_WORDS(target,source,count)
 404     }
 405 }
 406
 407 static void BIT_VECTOR_ins_words(wordptr addr, N_word total, N_word count,
 408                                  boolean clear)
 409 {
 410     N_word length;
 411
 412     if ((total > 0) and (count > 0))
 413     {
 414         if (count > total) count = total;
 415         length = total - count;
 416         if (length > 0) BIT_VECTOR_mov_words(addr+count,addr,length);
 417         if (clear)      BIT_VECTOR_zro_words(addr,count);
 418     }
 419 }
 420
 421 static void BIT_VECTOR_del_words(wordptr addr, N_word total, N_word count,
 422                                  boolean clear)
 423 {
 424     N_word length;
 425
 426     if ((total > 0) and (count > 0))
 427     {
 428         if (count > total) count = total;
 429         length = total - count;
 430         if (length > 0) BIT_VECTOR_mov_words(addr,addr+count,length);
 431         if (clear)      BIT_VECTOR_zro_words(addr+length,count);
 432     }
 433 }
 434
 435 static void BIT_VECTOR_reverse(charptr string, N_word length)
 436 {
 437     charptr last;
 438     N_char  temp;
 439
 440     if (length > 1)
 441     {
 442         last = string + length - 1;
 443         while (string < last)
 444         {
 445             temp = *string;
 446             *string = *last;
 447             *last = temp;
 448             string++;
 449             last--;
 450         }
 451     }
 452 }
 453
 454 static N_word BIT_VECTOR_int2str(charptr string, N_word value)
 455 {
 456     N_word  length;
 457     N_word  digit;
 458     charptr work;
 459
 460     work = string;
 461     if (value > 0)
 462     {
 463         length = 0;
 464         while (value > 0)
 465         {
 466             BIT_VECTOR_DIGITIZE(N_word,value,digit)
 467             *work++ = (N_char) digit;
 468             length++;
 469         }
 470         BIT_VECTOR_reverse(string,length);
 471     }
 472     else
 473     {
 474         length = 1;
 475         *work++ = (N_char) '0';
 476     }
 477     return(length);
 478 }
 479
 480 static N_word BIT_VECTOR_str2int(charptr string, N_word *value)
 481 {
 482     N_word  length;
 483     N_word  digit;
 484
 485     *value = 0;
 486     length = 0;
 487     digit = (N_word) *string++;
 488     /* separate because isdigit() is likely a macro! */
 489     while (isdigit((int)digit) != 0)
 490     {
 491         length++;
 492         digit -= (N_word) '0';
 493         if (*value) *value *= 10;
 494         *value += digit;
 495         digit = (N_word) *string++;
 496     }
 497     return(length);
 498 }
 499
 500     /********************************************/
 501     /* routine to convert error code to string: */
 502     /********************************************/
 503
 504 charptr BitVector_Error(ErrCode error)
 505 {
 506     switch (error)
 507     {
 508         case ErrCode_Ok:   return( (charptr)     NULL     ); break;
 509         case ErrCode_Type: return( (charptr) ERRCODE_TYPE ); break;
 510         case ErrCode_Bits: return( (charptr) ERRCODE_BITS ); break;
 511         case ErrCode_Word: return( (charptr) ERRCODE_WORD ); break;
 512         case ErrCode_Long: return( (charptr) ERRCODE_LONG ); break;
 513         case ErrCode_Powr: return( (charptr) ERRCODE_POWR ); break;
 514         case ErrCode_Loga: return( (charptr) ERRCODE_LOGA ); break;
 515         case ErrCode_Null: return( (charptr) ERRCODE_NULL ); break;
 516         case ErrCode_Indx: return( (charptr) ERRCODE_INDX ); break;
 517         case ErrCode_Ordr: return( (charptr) ERRCODE_ORDR ); break;
 518         case ErrCode_Size: return( (charptr) ERRCODE_SIZE ); break;
 519         case ErrCode_Pars: return( (charptr) ERRCODE_PARS ); break;
 520         case ErrCode_Ovfl: return( (charptr) ERRCODE_OVFL ); break;
 521         case ErrCode_Same: return( (charptr) ERRCODE_SAME ); break;
 522         case ErrCode_Expo: return( (charptr) ERRCODE_EXPO ); break;
 523         case ErrCode_Zero: return( (charptr) ERRCODE_ZERO ); break;
 524         default:           return( (charptr) ERRCODE_OOPS ); break;
 525     }
 526 }
 527
 528     /*****************************************/
 529     /* automatic self-configuration routine: */
 530     /*****************************************/
 531
 532     /*******************************************************/
 533     /*                                                     */
 534     /*   MUST be called once prior to any other function   */
 535     /*   to initialize the machine dependent constants     */
 536     /*   of this package! (But call only ONCE, or you      */
 537     /*   will suffer memory leaks!)                        */
 538     /*                                                     */
 539     /*******************************************************/
 540
 541 ErrCode BitVector_Boot(void)
 542 {
 543     N_long longsample = 1L;
 544     N_word sample = LSB;
 545     N_word lsb;
 546
 547     if (sizeof(N_word) > sizeof(size_t)) return(ErrCode_Type);
 548
 549     BITS = 1;
 550     while (sample <<= 1) BITS++;    /* determine # of bits in a machine word */
 551
 552     if (BITS != (sizeof(N_word) << 3)) return(ErrCode_Bits);
 553
 554     if (BITS < 16) return(ErrCode_Word);
 555
 556     LONGBITS = 1;
 557     while (longsample <<= 1) LONGBITS++;  /* = # of bits in an unsigned long */
 558
 559     if (BITS > LONGBITS) return(ErrCode_Long);
 560
 561     LOGBITS = 0;
 562     sample = BITS;
 563     lsb = (sample AND LSB);
 564     while ((sample >>= 1) and (not lsb))
 565     {
 566         LOGBITS++;
 567         lsb = (sample AND LSB);
 568     }
 569
 570     if (sample) return(ErrCode_Powr);      /* # of bits is not a power of 2! */
 571
 572     if (BITS != (LSB << LOGBITS)) return(ErrCode_Loga);
 573
 574     MODMASK = BITS - 1;
 575     FACTOR = LOGBITS - 3;  /* ld(BITS / 8) = ld(BITS) - ld(8) = ld(BITS) - 3 */
 576     MSB = (LSB << MODMASK);
 577
 578     BITMASKTAB = (wordptr) malloc((size_t) (BITS << FACTOR));
 579
 580     if (BITMASKTAB == NULL) return(ErrCode_Null);
 581
 582     for ( sample = 0; sample < BITS; sample++ )
 583     {
 584         BITMASKTAB[sample] = (LSB << sample);
 585     }
 586
 587     LOG10 = (N_word) (MODMASK * 0.30103); /* = (BITS - 1) * ( ln 2 / ln 10 ) */
 588     EXP10 = power10(LOG10);
 589
 590     return(ErrCode_Ok);
 591 }
 592
 593 N_word BitVector_Size(N_int bits)           /* bit vector size (# of words)  */
 594 {
 595     N_word size;
 596
 597     size = bits >> LOGBITS;
 598     if (bits AND MODMASK) size++;
 599     return(size);
 600 }
 601
 602 N_word BitVector_Mask(N_int bits)           /* bit vector mask (unused bits) */
 603 {
 604     N_word mask;
 605
 606     mask = bits AND MODMASK;
 607     if (mask) mask = (N_word) ~(~0L << mask); else mask = (N_word) ~0L;
 608     return(mask);
 609 }
 610
 611 charptr BitVector_Version(void)
 612 {
 613     return((charptr)"6.4");
 614 }
 615
 616 N_int BitVector_Word_Bits(void)
 617 {
 618     return(BITS);
 619 }
 620
 621 N_int BitVector_Long_Bits(void)
 622 {
 623     return(LONGBITS);
 624 }
 625
 626 /********************************************************************/
 627 /*                                                                  */
 628 /*  WARNING: Do not "free()" constant character strings, i.e.,      */
 629 /*           don't call "BitVector_Dispose()" for strings returned  */
 630 /*           by "BitVector_Error()" or "BitVector_Version()"!       */
 631 /*                                                                  */
 632 /*  ONLY call this function for strings allocated with "malloc()",  */
 633 /*  i.e., the strings returned by the functions "BitVector_to_*()"  */
 634 /*  and "BitVector_Block_Read()"!                                   */
 635 /*                                                                  */
 636 /********************************************************************/
 637
 638 void BitVector_Dispose(charptr string)                      /* free string   */
 639 {
 640     if (string != NULL) free((voidptr) string);
 641 }
 642
 643 void BitVector_Destroy(wordptr addr)                        /* free bitvec   */
 644 {
 645     if (addr != NULL)
 646     {
 647         addr -= BIT_VECTOR_HIDDEN_WORDS;
 648         free((voidptr) addr);
 649     }
 650 }
 651
 652 void BitVector_Destroy_List(listptr list, N_int count)      /* free list     */
 653 {
 654     listptr slot;
 655
 656     if (list != NULL)
 657     {
 658         slot = list;
 659         while (count-- > 0)
 660         {
 661             BitVector_Destroy(*slot++);
 662         }
 663         free((voidptr) list);
 664     }
 665 }
 666
 667 wordptr BitVector_Create(N_int bits, boolean clear)         /* malloc        */
 668 {
 669     N_word  size;
 670     N_word  mask;
 671     N_word  bytes;
 672     wordptr addr;
 673     wordptr zero;
 674
 675     size = BitVector_Size(bits);
 676     mask = BitVector_Mask(bits);
 677     bytes = (size + BIT_VECTOR_HIDDEN_WORDS) << FACTOR;
 678     addr = (wordptr) malloc((size_t) bytes);
 679     if (addr != NULL)
 680     {
 681         *addr++ = bits;
 682         *addr++ = size;
 683         *addr++ = mask;
 684         if (clear)
 685         {
 686             zero = addr;
 687             BIT_VECTOR_ZERO_WORDS(zero,size)
 688         }
 689     }
 690     return(addr);
 691 }
 692
 693 listptr BitVector_Create_List(N_int bits, boolean clear, N_int count)
 694 {
 695     listptr list = NULL;
 696     listptr slot;
 697     wordptr addr;
 698     N_int   i;
 699
 700     if (count > 0)
 701     {
 702         list = (listptr) malloc(sizeof(wordptr) * count);
 703         if (list != NULL)
 704         {
 705             slot = list;
 706             for ( i = 0; i < count; i++ )
 707             {
 708                 addr = BitVector_Create(bits,clear);
 709                 if (addr == NULL)
 710                 {
 711                     BitVector_Destroy_List(list,i);
 712                     return(NULL);
 713                 }
 714                 *slot++ = addr;
 715             }
 716         }
 717     }
 718     return(list);
 719 }
 720
 721 wordptr BitVector_Resize(wordptr oldaddr, N_int bits)       /* realloc       */
 722 {
 723     N_word  bytes;
 724     N_word  oldsize;
 725     N_word  oldmask;
 726     N_word  newsize;
 727     N_word  newmask;
 728     wordptr newaddr;
 729     wordptr source;
 730     wordptr target;
 731
 732     oldsize = size_(oldaddr);
 733     oldmask = mask_(oldaddr);
 734     newsize = BitVector_Size(bits);
 735     newmask = BitVector_Mask(bits);
 736     if (oldsize > 0) *(oldaddr+oldsize-1) &= oldmask;
 737     if (newsize <= oldsize)
 738     {
 739         newaddr = oldaddr;
 740         bits_(newaddr) = bits;
 741         size_(newaddr) = newsize;
 742         mask_(newaddr) = newmask;
 743         if (newsize > 0) *(newaddr+newsize-1) &= newmask;
 744     }
 745     else
 746     {
 747         bytes = (newsize + BIT_VECTOR_HIDDEN_WORDS) << FACTOR;
 748         newaddr = (wordptr) malloc((size_t) bytes);
 749         if (newaddr != NULL)
 750         {
 751             *newaddr++ = bits;
 752             *newaddr++ = newsize;
 753             *newaddr++ = newmask;
 754             target = newaddr;
 755             source = oldaddr;
 756             newsize -= oldsize;
 757             BIT_VECTOR_COPY_WORDS(target,source,oldsize)
 758             BIT_VECTOR_ZERO_WORDS(target,newsize)
 759         }
 760         BitVector_Destroy(oldaddr);
 761     }
 762     return(newaddr);
 763 }
 764
 765 wordptr BitVector_Shadow(wordptr addr)     /* makes new, same size but empty */
 766 {
 767     return( BitVector_Create(bits_(addr),true) );
 768 }
 769
 770 wordptr BitVector_Clone(wordptr addr)               /* makes exact duplicate */
 771 {
 772     N_word  bits;
 773     wordptr twin;
 774
 775     bits = bits_(addr);
 776     twin = BitVector_Create(bits,false);
 777     if ((twin != NULL) and (bits > 0))
 778         BIT_VECTOR_cpy_words(twin,addr,size_(addr));
 779     return(twin);
 780 }
 781
 782 wordptr BitVector_Concat(wordptr X, wordptr Y)      /* returns concatenation */
 783 {
 784     /* BEWARE that X = most significant part, Y = least significant part! */
 785
 786     N_word  bitsX;
 787     N_word  bitsY;
 788     N_word  bitsZ;
 789     wordptr Z;
 790
 791     bitsX = bits_(X);
 792     bitsY = bits_(Y);
 793     bitsZ = bitsX + bitsY;
 794     Z = BitVector_Create(bitsZ,false);
 795     if ((Z != NULL) and (bitsZ > 0))
 796     {
 797         BIT_VECTOR_cpy_words(Z,Y,size_(Y));
 798         BitVector_Interval_Copy(Z,X,bitsY,0,bitsX);
 799         *(Z+size_(Z)-1) &= mask_(Z);
 800     }
 801     return(Z);
 802 }
 803
 804 void BitVector_Copy(wordptr X, wordptr Y)                           /* X = Y */
 805 {
 806     N_word  sizeX = size_(X);
 807     N_word  sizeY = size_(Y);
 808     N_word  maskX = mask_(X);
 809     N_word  maskY = mask_(Y);
 810     N_word  fill  = 0;
 811     wordptr lastX;
 812     wordptr lastY;
 813
 814     if ((X != Y) and (sizeX > 0))
 815     {
 816         lastX = X + sizeX - 1;
 817         if (sizeY > 0)
 818         {
 819             lastY = Y + sizeY - 1;
 820             if ( (*lastY AND (maskY AND NOT (maskY >> 1))) == 0 ) *lastY &= maskY;
 821             else
 822             {
 823                 fill = (N_word) ~0L;
 824                 *lastY |= NOT maskY;
 825             }
 826             while ((sizeX > 0) and (sizeY > 0))
 827             {
 828                 *X++ = *Y++;
 829                 sizeX--;
 830                 sizeY--;
 831             }
 832             *lastY &= maskY;
 833         }
 834         while (sizeX-- > 0) *X++ = fill;
 835         *lastX &= maskX;
 836     }
 837 }
 838
 839 void BitVector_Empty(wordptr addr)                        /* X = {}  clr all */
 840 {
 841     N_word size = size_(addr);
 842
 843     BIT_VECTOR_ZERO_WORDS(addr,size)
 844 }
 845
 846 void BitVector_Fill(wordptr addr)                         /* X = ~{} set all */
 847 {
 848     N_word size = size_(addr);
 849     N_word mask = mask_(addr);
 850     N_word fill = (N_word) ~0L;
 851
 852     if (size > 0)
 853     {
 854         BIT_VECTOR_FILL_WORDS(addr,fill,size)
 855         *(--addr) &= mask;
 856     }
 857 }
 858
 859 void BitVector_Flip(wordptr addr)                         /* X = ~X flip all */
 860 {
 861     N_word size = size_(addr);
 862     N_word mask = mask_(addr);
 863     N_word flip = (N_word) ~0L;
 864
 865     if (size > 0)
 866     {
 867         BIT_VECTOR_FLIP_WORDS(addr,flip,size)
 868         *(--addr) &= mask;
 869     }
 870 }
 871
 872 void BitVector_Primes(wordptr addr)
 873 {
 874     N_word  bits = bits_(addr);
 875     N_word  size = size_(addr);
 876     wordptr work;
 877     N_word  temp;
 878     N_word  i,j;
 879
 880     if (size > 0)
 881     {
 882         temp = 0xAAAA;
 883         i = BITS >> 4;
 884         while (--i > 0)
 885         {
 886             temp <<= 16;
 887             temp |= 0xAAAA;
 888         }
 889         i = size;
 890         work = addr;
 891         *work++ = temp XOR 0x0006;
 892         while (--i > 0) *work++ = temp;
 893         for ( i = 3; (j = i * i) < bits; i += 2 )
 894         {
 895             for ( ; j < bits; j += i ) BIT_VECTOR_CLR_BIT(addr,j)
 896         }
 897         *(addr+size-1) &= mask_(addr);
 898     }
 899 }
 900
 901 void BitVector_Reverse(wordptr X, wordptr Y)
 902 {
 903     N_word bits = bits_(X);
 904     N_word mask;
 905     N_word bit;
 906     N_word value;
 907
 908     if (bits > 0)
 909     {
 910         if (X == Y) BitVector_Interval_Reverse(X,0,bits-1);
 911         else if (bits == bits_(Y))
 912         {
 913 /*          mask = mask_(Y);  */
 914 /*          mask &= NOT (mask >> 1);  */
 915             mask = BITMASKTAB[(bits-1) AND MODMASK];
 916             Y += size_(Y) - 1;
 917             value = 0;
 918             bit = LSB;
 919             while (bits-- > 0)
 920             {
 921                 if ((*Y AND mask) != 0)
 922                 {
 923                     value |= bit;
 924                 }
 925                 if (not (mask >>= 1))
 926                 {
 927                     Y--;
 928                     mask = MSB;
 929                 }
 930                 if (not (bit <<= 1))
 931                 {
 932                     *X++ = value;
 933                     value = 0;
 934                     bit = LSB;
 935                 }
 936             }
 937             if (bit > LSB) *X = value;
 938         }
 939     }
 940 }
 941
 942 void BitVector_Interval_Empty(wordptr addr, N_int lower, N_int upper)
 943 {                                                  /* X = X \ [lower..upper] */
 944     N_word  bits = bits_(addr);
 945     N_word  size = size_(addr);
 946     wordptr loaddr;
 947     wordptr hiaddr;
 948     N_word  lobase;
 949     N_word  hibase;
 950     N_word  lomask;
 951     N_word  himask;
 952     N_word  diff;
 953
 954     if ((size > 0) and (lower < bits) and (upper < bits) and (lower <= upper))
 955     {
 956         lobase = lower >> LOGBITS;
 957         hibase = upper >> LOGBITS;
 958         diff = hibase - lobase;
 959         loaddr = addr + lobase;
 960         hiaddr = addr + hibase;
 961
 962         lomask = (N_word)   (~0L << (lower AND MODMASK));
 963         himask = (N_word) ~((~0L << (upper AND MODMASK)) << 1);
 964
 965         if (diff == 0)
 966         {
 967             *loaddr &= NOT (lomask AND himask);
 968         }
 969         else
 970         {
 971             *loaddr++ &= NOT lomask;
 972             while (--diff > 0)
 973             {
 974                 *loaddr++ = 0;
 975             }
 976             *hiaddr &= NOT himask;
 977         }
 978     }
 979 }
 980
 981 void BitVector_Interval_Fill(wordptr addr, N_int lower, N_int upper)
 982 {                                                  /* X = X + [lower..upper] */
 983     N_word  bits = bits_(addr);
 984     N_word  size = size_(addr);
 985     N_word  fill = (N_word) ~0L;
 986     wordptr loaddr;
 987     wordptr hiaddr;
 988     N_word  lobase;
 989     N_word  hibase;
 990     N_word  lomask;
 991     N_word  himask;
 992     N_word  diff;
 993
 994     if ((size > 0) and (lower < bits) and (upper < bits) and (lower <= upper))
 995     {
 996         lobase = lower >> LOGBITS;
 997         hibase = upper >> LOGBITS;
 998         diff = hibase - lobase;
 999         loaddr = addr + lobase;
1000         hiaddr = addr + hibase;
1001
1002         lomask = (N_word)   (~0L << (lower AND MODMASK));
1003         himask = (N_word) ~((~0L << (upper AND MODMASK)) << 1);
1004
1005         if (diff == 0)
1006         {
1007             *loaddr |= (lomask AND himask);
1008         }
1009         else
1010         {
1011             *loaddr++ |= lomask;
1012             while (--diff > 0)
1013             {
1014                 *loaddr++ = fill;
1015             }
1016             *hiaddr |= himask;
1017         }
1018         *(addr+size-1) &= mask_(addr);
1019     }
1020 }
1021
1022 void BitVector_Interval_Flip(wordptr addr, N_int lower, N_int upper)
1023 {                                                  /* X = X ^ [lower..upper] */
1024     N_word  bits = bits_(addr);
1025     N_word  size = size_(addr);
1026     N_word  flip = (N_word) ~0L;
1027     wordptr loaddr;
1028     wordptr hiaddr;
1029     N_word  lobase;
1030     N_word  hibase;
1031     N_word  lomask;
1032     N_word  himask;
1033     N_word  diff;
1034
1035     if ((size > 0) and (lower < bits) and (upper < bits) and (lower <= upper))
1036     {
1037         lobase = lower >> LOGBITS;
1038         hibase = upper >> LOGBITS;
1039         diff = hibase - lobase;
1040         loaddr = addr + lobase;
1041         hiaddr = addr + hibase;
1042
1043         lomask = (N_word)   (~0L << (lower AND MODMASK));
1044         himask = (N_word) ~((~0L << (upper AND MODMASK)) << 1);
1045
1046         if (diff == 0)
1047         {
1048             *loaddr ^= (lomask AND himask);
1049         }
1050         else
1051         {
1052             *loaddr++ ^= lomask;
1053             while (--diff > 0)
1054             {
1055                 *loaddr++ ^= flip;
1056             }
1057             *hiaddr ^= himask;
1058         }
1059         *(addr+size-1) &= mask_(addr);
1060     }
1061 }
1062
1063 void BitVector_Interval_Reverse(wordptr addr, N_int lower, N_int upper)
1064 {
1065     N_word  bits = bits_(addr);
1066     wordptr loaddr;
1067     wordptr hiaddr;
1068     N_word  lomask;
1069     N_word  himask;
1070
1071     if ((bits > 0) and (lower < bits) and (upper < bits) and (lower < upper))
1072     {
1073         loaddr = addr + (lower >> LOGBITS);
1074         hiaddr = addr + (upper >> LOGBITS);
1075         lomask = BITMASKTAB[lower AND MODMASK];
1076         himask = BITMASKTAB[upper AND MODMASK];
1077         for ( bits = upper - lower + 1; bits > 1; bits -= 2 )
1078         {
1079             if (((*loaddr AND lomask) != 0) XOR ((*hiaddr AND himask) != 0))
1080             {
1081                 *loaddr ^= lomask;  /* swap bits only if they differ! */
1082                 *hiaddr ^= himask;
1083             }
1084             if (not (lomask <<= 1))
1085             {
1086                 lomask = LSB;
1087                 loaddr++;
1088             }
1089             if (not (himask >>= 1))
1090             {
1091                 himask = MSB;
1092                 hiaddr--;
1093             }
1094         }
1095     }
1096 }
1097
1098 boolean BitVector_interval_scan_inc(wordptr addr, N_int start,
1099                                     N_intptr min, N_intptr max)
1100 {
1101     N_word  size = size_(addr);
1102     N_word  mask = mask_(addr);
1103     N_word  offset;
1104     N_word  bitmask;
1105     N_word  value;
1106     boolean empty;
1107
1108     if ((size == 0) or (start >= bits_(addr))) return(false);
1109
1110     *min = start;
1111     *max = start;
1112
1113     offset = start >> LOGBITS;
1114
1115     *(addr+size-1) &= mask;
1116
1117     addr += offset;
1118     size -= offset;
1119
1120     bitmask = BITMASKTAB[start AND MODMASK];
1121     mask = NOT (bitmask OR (bitmask - 1));
1122
1123     value = *addr++;
1124     if ((value AND bitmask) == 0)
1125     {
1126         value &= mask;
1127         if (value == 0)
1128         {
1129             offset++;
1130             empty = true;
1131             while (empty and (--size > 0))
1132             {
1133                 if ((value = *addr++)) empty = false; else offset++;
1134             }
1135             if (empty) return(false);
1136         }
1137         start = offset << LOGBITS;
1138         bitmask = LSB;
1139         mask = value;
1140         while (not (mask AND LSB))
1141         {
1142             bitmask <<= 1;
1143             mask >>= 1;
1144             start++;
1145         }
1146         mask = NOT (bitmask OR (bitmask - 1));
1147         *min = start;
1148         *max = start;
1149     }
1150     value = NOT value;
1151     value &= mask;
1152     if (value == 0)
1153     {
1154         offset++;
1155         empty = true;
1156         while (empty and (--size > 0))
1157         {
1158             if ((value = NOT *addr++)) empty = false; else offset++;
1159         }
1160         if (empty) value = LSB;
1161     }
1162     start = offset << LOGBITS;
1163     while (not (value AND LSB))
1164     {
1165         value >>= 1;
1166         start++;
1167     }
1168     *max = --start;
1169     return(true);
1170 }
1171
1172 boolean BitVector_interval_scan_dec(wordptr addr, N_int start,
1173                                     N_intptr min, N_intptr max)
1174 {
1175     N_word  size = size_(addr);
1176     N_word  mask = mask_(addr);
1177     N_word offset;
1178     N_word bitmask;
1179     N_word value;
1180     boolean empty;
1181
1182     if ((size == 0) or (start >= bits_(addr))) return(false);
1183
1184     *min = start;
1185     *max = start;
1186
1187     offset = start >> LOGBITS;
1188
1189     if (offset >= size) return(false);
1190
1191     *(addr+size-1) &= mask;
1192
1193     addr += offset;
1194     size = ++offset;
1195
1196     bitmask = BITMASKTAB[start AND MODMASK];
1197     mask = (bitmask - 1);
1198
1199     value = *addr--;
1200     if ((value AND bitmask) == 0)
1201     {
1202         value &= mask;
1203         if (value == 0)
1204         {
1205             offset--;
1206             empty = true;
1207             while (empty and (--size > 0))
1208             {
1209                 if ((value = *addr--)) empty = false; else offset--;
1210             }
1211             if (empty) return(false);
1212         }
1213         start = offset << LOGBITS;
1214         bitmask = MSB;
1215         mask = value;
1216         while (not (mask AND MSB))
1217         {
1218             bitmask >>= 1;
1219             mask <<= 1;
1220             start--;
1221         }
1222         mask = (bitmask - 1);
1223         *max = --start;
1224         *min = start;
1225     }
1226     value = NOT value;
1227     value &= mask;
1228     if (value == 0)
1229     {
1230         offset--;
1231         empty = true;
1232         while (empty and (--size > 0))
1233         {
1234             if ((value = NOT *addr--)) empty = false; else offset--;
1235         }
1236         if (empty) value = MSB;
1237     }
1238     start = offset << LOGBITS;
1239     while (not (value AND MSB))
1240     {
1241         value <<= 1;
1242         start--;
1243     }
1244     *min = start;
1245     return(true);
1246 }
1247
1248 void BitVector_Interval_Copy(wordptr X, wordptr Y, N_int Xoffset,
1249                              N_int Yoffset, N_int length)
1250 {
1251     N_word  bitsX = bits_(X);
1252     N_word  bitsY = bits_(Y);
1253     N_word  source = 0;        /* silence compiler warning */
1254     N_word  target = 0;        /* silence compiler warning */
1255     N_word  s_lo_base;
1256     N_word  s_hi_base;
1257     N_word  s_lo_bit;
1258     N_word  s_hi_bit;
1259     N_word  s_base;
1260     N_word  s_lower = 0;       /* silence compiler warning */
1261     N_word  s_upper = 0;       /* silence compiler warning */
1262     N_word  s_bits;
1263     N_word  s_min;
1264     N_word  s_max;
1265     N_word  t_lo_base;
1266     N_word  t_hi_base;
1267     N_word  t_lo_bit;
1268     N_word  t_hi_bit;
1269     N_word  t_base;
1270     N_word  t_lower = 0;       /* silence compiler warning */
1271     N_word  t_upper = 0;       /* silence compiler warning */
1272     N_word  t_bits;
1273     N_word  t_min;
1274     N_word  mask;
1275     N_word  bits;
1276     N_word  select;
1277     boolean ascending;
1278     boolean notfirst;
1279     wordptr Z = X;
1280
1281     if ((length > 0) and (Xoffset < bitsX) and (Yoffset < bitsY))
1282     {
1283         if ((Xoffset + length) > bitsX) length = bitsX - Xoffset;
1284         if ((Yoffset + length) > bitsY) length = bitsY - Yoffset;
1285
1286         ascending = (Xoffset <= Yoffset);
1287
1288         s_lo_base = Yoffset >> LOGBITS;
1289         s_lo_bit = Yoffset AND MODMASK;
1290         Yoffset += --length;
1291         s_hi_base = Yoffset >> LOGBITS;
1292         s_hi_bit = Yoffset AND MODMASK;
1293
1294         t_lo_base = Xoffset >> LOGBITS;
1295         t_lo_bit = Xoffset AND MODMASK;
1296         Xoffset += length;
1297         t_hi_base = Xoffset >> LOGBITS;
1298         t_hi_bit = Xoffset AND MODMASK;
1299
1300         if (ascending)
1301         {
1302             s_base = s_lo_base;
1303             t_base = t_lo_base;
1304         }
1305         else
1306         {
1307             s_base = s_hi_base;
1308             t_base = t_hi_base;
1309         }
1310         s_bits = 0;
1311         t_bits = 0;
1312         Y += s_base;
1313         X += t_base;
1314         notfirst = false;
1315         while (true)
1316         {
1317             if (t_bits == 0)
1318             {
1319                 if (notfirst)
1320                 {
1321                     *X = target;
1322                     if (ascending)
1323                     {
1324                         if (t_base == t_hi_base) break;
1325                         t_base++;
1326                         X++;
1327                     }
1328                     else
1329                     {
1330                         if (t_base == t_lo_base) break;
1331                         t_base--;
1332                         X--;
1333                     }
1334                 }
1335                 select = ((t_base == t_hi_base) << 1) OR (t_base == t_lo_base);
1336                 switch (select)
1337                 {
1338                     case 0:
1339                         t_lower = 0;
1340                         t_upper = BITS - 1;
1341                         t_bits = BITS;
1342                         target = 0;
1343                         break;
1344                     case 1:
1345                         t_lower = t_lo_bit;
1346                         t_upper = BITS - 1;
1347                         t_bits = BITS - t_lo_bit;
1348                         mask = (N_word) (~0L << t_lower);
1349                         target = *X AND NOT mask;
1350                         break;
1351                     case 2:
1352                         t_lower = 0;
1353                         t_upper = t_hi_bit;
1354                         t_bits = t_hi_bit + 1;
1355                         mask = (N_word) ((~0L << t_upper) << 1);
1356                         target = *X AND mask;
1357                         break;
1358                     case 3:
1359                         t_lower = t_lo_bit;
1360                         t_upper = t_hi_bit;
1361                         t_bits = t_hi_bit - t_lo_bit + 1;
1362                         mask = (N_word) (~0L << t_lower);
1363                         mask &= (N_word) ~((~0L << t_upper) << 1);
1364                         target = *X AND NOT mask;
1365                         break;
1366                 }
1367             }
1368             if (s_bits == 0)
1369             {
1370                 if (notfirst)
1371                 {
1372                     if (ascending)
1373                     {
1374                         if (s_base == s_hi_base) break;
1375                         s_base++;
1376                         Y++;
1377                     }
1378                     else
1379                     {
1380                         if (s_base == s_lo_base) break;
1381                         s_base--;
1382                         Y--;
1383                     }
1384                 }
1385                 source = *Y;
1386                 select = ((s_base == s_hi_base) << 1) OR (s_base == s_lo_base);
1387                 switch (select)
1388                 {
1389                     case 0:
1390                         s_lower = 0;
1391                         s_upper = BITS - 1;
1392                         s_bits = BITS;
1393                         break;
1394                     case 1:
1395                         s_lower = s_lo_bit;
1396                         s_upper = BITS - 1;
1397                         s_bits = BITS - s_lo_bit;
1398                         break;
1399                     case 2:
1400                         s_lower = 0;
1401                         s_upper = s_hi_bit;
1402                         s_bits = s_hi_bit + 1;
1403                         break;
1404                     case 3:
1405                         s_lower = s_lo_bit;
1406                         s_upper = s_hi_bit;
1407                         s_bits = s_hi_bit - s_lo_bit + 1;
1408                         break;
1409                 }
1410             }
1411             notfirst = true;
1412             if (s_bits > t_bits)
1413             {
1414                 bits = t_bits - 1;
1415                 if (ascending)
1416                 {
1417                     s_min = s_lower;
1418                     s_max = s_lower + bits;
1419                 }
1420                 else
1421                 {
1422                     s_max = s_upper;
1423                     s_min = s_upper - bits;
1424                 }
1425                 t_min = t_lower;
1426             }
1427             else
1428             {
1429                 bits = s_bits - 1;
1430                 if (ascending) t_min = t_lower;
1431                 else           t_min = t_upper - bits;
1432                 s_min = s_lower;
1433                 s_max = s_upper;
1434             }
1435             bits++;
1436             mask = (N_word) (~0L << s_min);
1437             mask &= (N_word) ~((~0L << s_max) << 1);
1438             if (s_min == t_min) target |= (source AND mask);
1439             else
1440             {
1441                 if (s_min < t_min) target |= (source AND mask) << (t_min-s_min);
1442                 else               target |= (source AND mask) >> (s_min-t_min);
1443             }
1444             if (ascending)
1445             {
1446                 s_lower += bits;
1447                 t_lower += bits;
1448             }
1449             else
1450             {
1451                 s_upper -= bits;
1452                 t_upper -= bits;
1453             }
1454             s_bits -= bits;
1455             t_bits -= bits;
1456         }
1457         *(Z+size_(Z)-1) &= mask_(Z);
1458     }
1459 }
1460
1461
1462 wordptr BitVector_Interval_Substitute(wordptr X, wordptr Y,
1463                                       N_int Xoffset, N_int Xlength,
1464                                       N_int Yoffset, N_int Ylength)
1465 {
1466     N_word Xbits = bits_(X);
1467     N_word Ybits = bits_(Y);
1468     N_word limit;
1469     N_word diff;
1470
1471     if ((Xoffset <= Xbits) and (Yoffset <= Ybits))
1472     {
1473         limit = Xoffset + Xlength;
1474         if (limit > Xbits)
1475         {
1476             limit = Xbits;
1477             Xlength = Xbits - Xoffset;
1478         }
1479         if ((Yoffset + Ylength) > Ybits)
1480         {
1481             Ylength = Ybits - Yoffset;
1482         }
1483         if (Xlength == Ylength)
1484         {
1485             if ((Ylength > 0) and ((X != Y) or (Xoffset != Yoffset)))
1486             {
1487                 BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1488             }
1489         }
1490         else /* Xlength != Ylength */
1491         {
1492             if (Xlength > Ylength)
1493             {
1494                 diff = Xlength - Ylength;
1495                 if (Ylength > 0) BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1496                 if (limit < Xbits) BitVector_Delete(X,Xoffset+Ylength,diff,false);
1497                 if ((X = BitVector_Resize(X,Xbits-diff)) == NULL) return(NULL);
1498             }
1499             else /* Ylength > Xlength  ==>  Ylength > 0 */
1500             {
1501                 diff = Ylength - Xlength;
1502                 if (X != Y)
1503                 {
1504                     if ((X = BitVector_Resize(X,Xbits+diff)) == NULL) return(NULL);
1505                     if (limit < Xbits) BitVector_Insert(X,limit,diff,false);
1506                     BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1507                 }
1508                 else /* in-place */
1509                 {
1510                     if ((Y = X = BitVector_Resize(X,Xbits+diff)) == NULL) return(NULL);
1511                     if (limit >= Xbits)
1512                     {
1513                         BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1514                     }
1515                     else /* limit < Xbits */
1516                     {
1517                         BitVector_Insert(X,limit,diff,false);
1518                         if ((Yoffset+Ylength) <= limit)
1519                         {
1520                             BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1521                         }
1522                         else /* overlaps or lies above critical area */
1523                         {
1524                             if (limit <= Yoffset)
1525                             {
1526                                 Yoffset += diff;
1527                                 BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1528                             }
1529                             else /* Yoffset < limit */
1530                             {
1531                                 Xlength = limit - Yoffset;
1532                                 BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Xlength);
1533                                 Yoffset = Xoffset + Ylength; /* = limit + diff */
1534                                 Xoffset += Xlength;
1535                                 Ylength -= Xlength;
1536                                 BitVector_Interval_Copy(X,Y,Xoffset,Yoffset,Ylength);
1537                             }
1538                         }
1539                     }
1540                 }
1541             }
1542         }
1543     }
1544     return(X);
1545 }
1546
1547 boolean BitVector_is_empty(wordptr addr)                    /* X == {} ?     */
1548 {
1549     N_word  size = size_(addr);
1550     boolean r = true;
1551
1552     if (size > 0)
1553     {
1554         *(addr+size-1) &= mask_(addr);
1555         while (r and (size-- > 0)) r = ( *addr++ == 0 );
1556     }
1557     return(r);
1558 }
1559
1560 boolean BitVector_is_full(wordptr addr)                     /* X == ~{} ?    */
1561 {
1562     N_word  size = size_(addr);
1563     N_word  mask = mask_(addr);
1564     boolean r = false;
1565     wordptr last;
1566
1567     if (size > 0)
1568     {
1569         r = true;
1570         last = addr + size - 1;
1571         *last |= NOT mask;
1572         while (r and (size-- > 0)) r = ( NOT *addr++ == 0 );
1573         *last &= mask;
1574     }
1575     return(r);
1576 }
1577
1578 boolean BitVector_equal(wordptr X, wordptr Y)               /* X == Y ?      */
1579 {
1580     N_word  size = size_(X);
1581     N_word  mask = mask_(X);
1582     boolean r = false;
1583
1584     if (bits_(X) == bits_(Y))
1585     {
1586         r = true;
1587         if (size > 0)
1588         {
1589             *(X+size-1) &= mask;
1590             *(Y+size-1) &= mask;
1591             while (r and (size-- > 0)) r = (*X++ == *Y++);
1592         }
1593     }
1594     return(r);
1595 }
1596
1597 Z_int BitVector_Lexicompare(wordptr X, wordptr Y)           /* X <,=,> Y ?   */
1598 {                                                           /*  unsigned     */
1599     N_word  bitsX = bits_(X);
1600     N_word  bitsY = bits_(Y);
1601     N_word  size  = size_(X);
1602     boolean r = true;
1603
1604     if (bitsX == bitsY)
1605     {
1606         if (size > 0)
1607         {
1608             X += size;
1609             Y += size;
1610             while (r and (size-- > 0)) r = (*(--X) == *(--Y));
1611         }
1612         if (r) return((Z_int) 0);
1613         else
1614         {
1615             if (*X < *Y) return((Z_int) -1); else return((Z_int) 1);
1616         }
1617     }
1618     else
1619     {
1620         if (bitsX < bitsY) return((Z_int) -1); else return((Z_int) 1);
1621     }
1622 }
1623
1624 Z_int BitVector_Compare(wordptr X, wordptr Y)               /* X <,=,> Y ?   */
1625 {                                                           /*   signed      */
1626     N_word  bitsX = bits_(X);
1627     N_word  bitsY = bits_(Y);
1628     N_word  size  = size_(X);
1629     N_word  mask  = mask_(X);
1630     N_word  sign;
1631     boolean r = true;
1632
1633     if (bitsX == bitsY)
1634     {
1635         if (size > 0)
1636         {
1637             X += size;
1638             Y += size;
1639             mask &= NOT (mask >> 1);
1640             if ((sign = (*(X-1) AND mask)) != (*(Y-1) AND mask))
1641             {
1642                 if (sign) return((Z_int) -1); else return((Z_int) 1);
1643             }
1644             while (r and (size-- > 0)) r = (*(--X) == *(--Y));
1645         }
1646         if (r) return((Z_int) 0);
1647         else
1648         {
1649             if (*X < *Y) return((Z_int) -1); else return((Z_int) 1);
1650         }
1651     }
1652     else
1653     {
1654         if (bitsX < bitsY) return((Z_int) -1); else return((Z_int) 1);
1655     }
1656 }
1657
1658 charptr BitVector_to_Hex(wordptr addr)
1659 {
1660     N_word  bits = bits_(addr);
1661     N_word  size = size_(addr);
1662     N_word  value;
1663     N_word  count;
1664     N_word  digit;
1665     N_word  length;
1666     charptr string;
1667
1668     length = bits >> 2;
1669     if (bits AND 0x0003) length++;
1670     string = (charptr) malloc((size_t) (length+1));
1671     if (string == NULL) return(NULL);
1672     string += length;
1673     *string = (N_char) '\0';
1674     if (size > 0)
1675     {
1676         *(addr+size-1) &= mask_(addr);
1677         while ((size-- > 0) and (length > 0))
1678         {
1679             value = *addr++;
1680             count = BITS >> 2;
1681             while ((count-- > 0) and (length > 0))
1682             {
1683                 digit = value AND 0x000F;
1684                 if (digit > 9) digit += (N_word) 'A' - 10;
1685                 else           digit += (N_word) '0';
1686                 *(--string) = (N_char) digit; length--;
1687                 if ((count > 0) and (length > 0)) value >>= 4;
1688             }
1689         }
1690     }
1691     return(string);
1692 }
1693
1694 ErrCode BitVector_from_Hex(wordptr addr, charptr string)
1695 {
1696     N_word  size = size_(addr);
1697     N_word  mask = mask_(addr);
1698     boolean ok = true;
1699     N_word  length;
1700     N_word  value;
1701     N_word  count;
1702     int     digit;
1703
1704     if (size > 0)
1705     {
1706         length = strlen((char *) string);
1707         string += length;
1708         while (size-- > 0)
1709         {
1710             value = 0;
1711             for ( count = 0; (ok and (length > 0) and (count < BITS)); count += 4 )
1712             {
1713                 digit = (int) *(--string); length--;
1714                 /* separate because toupper() is likely a macro! */
1715                 digit = toupper(digit);
1716                 if ((ok = (isxdigit(digit) != 0)))
1717                 {
1718                     if (digit >= (int) 'A') digit -= (int) 'A' - 10;
1719                     else                    digit -= (int) '0';
1720                     value |= (((N_word) digit) << count);
1721                 }
1722             }
1723             *addr++ = value;
1724         }
1725         *(--addr) &= mask;
1726     }
1727     if (ok) return(ErrCode_Ok);
1728     else    return(ErrCode_Pars);
1729 }
1730
1731 charptr BitVector_to_Bin(wordptr addr)
1732 {
1733     N_word  size = size_(addr);
1734     N_word  value;
1735     N_word  count;
1736     N_word  digit;
1737     N_word  length;
1738     charptr string;
1739
1740     length = bits_(addr);
1741     string = (charptr) malloc((size_t) (length+1));
1742     if (string == NULL) return(NULL);
1743     string += length;
1744     *string = (N_char) '\0';
1745     if (size > 0)
1746     {
1747         *(addr+size-1) &= mask_(addr);
1748         while (size-- > 0)
1749         {
1750             value = *addr++;
1751             count = BITS;
1752             if (count > length) count = length;
1753             while (count-- > 0)
1754             {
1755                 digit = value AND 0x0001;
1756                 digit += (N_word) '0';
1757                 *(--string) = (N_char) digit; length--;
1758                 if (count > 0) value >>= 1;
1759             }
1760         }
1761     }
1762     return(string);
1763 }
1764
1765 ErrCode BitVector_from_Bin(wordptr addr, charptr string)
1766 {
1767     N_word  size = size_(addr);
1768     N_word  mask = mask_(addr);
1769     boolean ok = true;
1770     N_word  length;
1771     N_word  value;
1772     N_word  count;
1773     int     digit;
1774
1775     if (size > 0)
1776     {
1777         length = strlen((char *) string);
1778         string += length;
1779         while (size-- > 0)
1780         {
1781             value = 0;
1782             for ( count = 0; (ok and (length > 0) and (count < BITS)); count++ )
1783             {
1784                 digit = (int) *(--string); length--;
1785                 switch (digit)
1786                 {
1787                     case (int) '0':
1788                         break;
1789                     case (int) '1':
1790                         value |= BITMASKTAB[count];
1791                         break;
1792                     default:
1793                         ok = false;
1794                         break;
1795                 }
1796             }
1797             *addr++ = value;
1798         }
1799         *(--addr) &= mask;
1800     }
1801     if (ok) return(ErrCode_Ok);
1802     else    return(ErrCode_Pars);
1803 }
1804
1805 charptr BitVector_to_Dec(wordptr addr)
1806 {
1807     N_word  bits = bits_(addr);
1808     N_word  length;
1809     N_word  digits;
1810     N_word  count;
1811     N_word  q;
1812     N_word  r;
1813     boolean loop;
1814     charptr result;
1815     charptr string;
1816     wordptr quot;
1817     wordptr rest;
1818     wordptr temp;
1819     wordptr base;
1820     Z_int   sign;
1821
1822     length = (N_word) (bits / 3.3);        /* digits = bits * ln(2) / ln(10) */
1823     length += 2; /* compensate for truncating & provide space for minus sign */
1824     result = (charptr) malloc((size_t) (length+1));    /* remember the '\0'! */
1825     if (result == NULL) return(NULL);
1826     string = result;
1827     sign = BitVector_Sign(addr);
1828     if ((bits < 4) or (sign == 0))
1829     {
1830         if (bits > 0) digits = *addr; else digits = (N_word) 0;
1831         if (sign < 0) digits = ((N_word)(-((Z_word)digits))) AND mask_(addr);
1832         *string++ = (N_char) digits + (N_char) '0';
1833         digits = 1;
1834     }
1835     else
1836     {
1837         quot = BitVector_Create(bits,false);
1838         if (quot == NULL)
1839         {
1840             BitVector_Dispose(result);
1841             return(NULL);
1842         }
1843         rest = BitVector_Create(bits,false);
1844         if (rest == NULL)
1845         {
1846             BitVector_Dispose(result);
1847             BitVector_Destroy(quot);
1848             return(NULL);
1849         }
1850         temp = BitVector_Create(bits,false);
1851         if (temp == NULL)
1852         {
1853             BitVector_Dispose(result);
1854             BitVector_Destroy(quot);
1855             BitVector_Destroy(rest);
1856             return(NULL);
1857         }
1858         base = BitVector_Create(bits,true);
1859         if (base == NULL)
1860         {
1861             BitVector_Dispose(result);
1862             BitVector_Destroy(quot);
1863             BitVector_Destroy(rest);
1864             BitVector_Destroy(temp);
1865             return(NULL);
1866         }
1867         if (sign < 0) BitVector_Negate(quot,addr);
1868         else           BitVector_Copy(quot,addr);
1869         digits = 0;
1870         *base = EXP10;
1871         loop = (bits >= BITS);
1872         do
1873         {
1874             if (loop)
1875             {
1876                 BitVector_Copy(temp,quot);
1877                 if (BitVector_Div_Pos(quot,temp,base,rest))
1878                 {
1879                     BitVector_Dispose(result); /* emergency exit */
1880                     BitVector_Destroy(quot);
1881                     BitVector_Destroy(rest);   /* should never occur */
1882                     BitVector_Destroy(temp);   /* under normal operation */
1883                     BitVector_Destroy(base);
1884                     return(NULL);
1885                 }
1886                 loop = not BitVector_is_empty(quot);
1887                 q = *rest;
1888             }
1889             else q = *quot;
1890             count = LOG10;
1891             while (((loop and (count-- > 0)) or ((not loop) and (q != 0))) and
1892                 (digits < length))
1893             {
1894                 if (q != 0)
1895                 {
1896                     BIT_VECTOR_DIGITIZE(N_word,q,r)
1897                 }
1898                 else r = (N_word) '0';
1899                 *string++ = (N_char) r;
1900                 digits++;
1901             }
1902         }
1903         while (loop and (digits < length));
1904         BitVector_Destroy(quot);
1905         BitVector_Destroy(rest);
1906         BitVector_Destroy(temp);
1907         BitVector_Destroy(base);
1908     }
1909     if ((sign < 0) and (digits < length))
1910     {
1911         *string++ = (N_char) '-';
1912         digits++;
1913     }
1914     *string = (N_char) '\0';
1915     BIT_VECTOR_reverse(result,digits);
1916     return(result);
1917 }
1918
1919 ErrCode BitVector_from_Dec(wordptr addr, charptr string)
1920 {
1921     ErrCode error = ErrCode_Ok;
1922     N_word  bits = bits_(addr);
1923     N_word  mask = mask_(addr);
1924     boolean init = (bits > BITS);
1925     boolean minus;
1926     boolean shift;
1927     boolean carry;
1928     wordptr term;
1929     wordptr base;
1930     wordptr prod;
1931     wordptr rank;
1932     wordptr temp;
1933     N_word  accu;
1934     N_word  powr;
1935     N_word  count;
1936     N_word  length;
1937     int     digit;
1938
1939     if (bits > 0)
1940     {
1941         length = strlen((char *) string);
1942         if (length == 0) return(ErrCode_Pars);
1943         digit = (int) *string;
1944         if ((minus = (digit == (int) '-')) or
1945                      (digit == (int) '+'))
1946         {
1947             string++;
1948             if (--length == 0) return(ErrCode_Pars);
1949         }
1950         string += length;
1951         term = BitVector_Create(BITS,false);
1952         if (term == NULL)
1953         {
1954             return(ErrCode_Null);
1955         }
1956         base = BitVector_Create(BITS,false);
1957         if (base == NULL)
1958         {
1959             BitVector_Destroy(term);
1960             return(ErrCode_Null);
1961         }
1962         prod = BitVector_Create(bits,init);
1963         if (prod == NULL)
1964         {
1965             BitVector_Destroy(term);
1966             BitVector_Destroy(base);
1967             return(ErrCode_Null);
1968         }
1969         rank = BitVector_Create(bits,init);
1970         if (rank == NULL)
1971         {
1972             BitVector_Destroy(term);
1973             BitVector_Destroy(base);
1974             BitVector_Destroy(prod);
1975             return(ErrCode_Null);
1976         }
1977         temp = BitVector_Create(bits,false);
1978         if (temp == NULL)
1979         {
1980             BitVector_Destroy(term);
1981             BitVector_Destroy(base);
1982             BitVector_Destroy(prod);
1983             BitVector_Destroy(rank);
1984             return(ErrCode_Null);
1985         }
1986         BitVector_Empty(addr);
1987         *base = EXP10;
1988         shift = false;
1989         while ((not error) and (length > 0))
1990         {
1991             accu = 0;
1992             powr = 1;
1993             count = LOG10;
1994             while ((not error) and (length > 0) and (count-- > 0))
1995             {
1996                 digit = (int) *(--string); length--;
1997                 /* separate because isdigit() is likely a macro! */
1998                 if (isdigit(digit) != 0)
1999                 {
2000                     accu += ((N_word) digit - (N_word) '0') * powr;
2001                     powr *= 10;
2002                 }
2003                 else error = ErrCode_Pars;
2004             }
2005             if (not error)
2006             {
2007                 if (shift)
2008                 {
2009                     *term = accu;
2010                     BitVector_Copy(temp,rank);
2011                     error = BitVector_Mul_Pos(prod,temp,term,false);
2012                 }
2013                 else
2014                 {
2015                     *prod = accu;
2016                     if ((not init) and ((accu AND NOT mask) != 0)) error = ErrCode_Ovfl;
2017                 }
2018                 if (not error)
2019                 {
2020                     carry = false;
2021                     BitVector_compute(addr,addr,prod,false,&carry);
2022                     /* ignores sign change (= overflow) but not */
2023                     /* numbers too large (= carry) for resulting bit vector */
2024                     if (carry) error = ErrCode_Ovfl;
2025                     else
2026                     {
2027                         if (length > 0)
2028                         {
2029                             if (shift)
2030                             {
2031                                 BitVector_Copy(temp,rank);
2032                                 error = BitVector_Mul_Pos(rank,temp,base,false);
2033                             }
2034                             else
2035                             {
2036                                 *rank = *base;
2037                                 shift = true;
2038                             }
2039                         }
2040                     }
2041                 }
2042             }
2043         }
2044         BitVector_Destroy(term);
2045         BitVector_Destroy(base);
2046         BitVector_Destroy(prod);
2047         BitVector_Destroy(rank);
2048         BitVector_Destroy(temp);
2049         if (not error and minus)
2050         {
2051             BitVector_Negate(addr,addr);
2052             if ((*(addr + size_(addr) - 1) AND mask AND NOT (mask >> 1)) == 0)
2053                 error = ErrCode_Ovfl;
2054         }
2055     }
2056     return(error);
2057 }
2058
2059 charptr BitVector_to_Enum(wordptr addr)
2060 {
2061     N_word  bits = bits_(addr);
2062     N_word  sample;
2063     N_word  length;
2064     N_word  digits;
2065     N_word  factor;
2066     N_word  power;
2067     N_word  start;
2068     N_word  min;
2069     N_word  max;
2070     charptr string;
2071     charptr target;
2072     boolean comma;
2073
2074     if (bits > 0)
2075     {
2076         sample = bits - 1;  /* greatest possible index */
2077         length = 2;         /* account for index 0 and terminating '\0' */
2078         digits = 1;         /* account for intervening dashes and commas */
2079         factor = 1;
2080         power = 10;
2081         while (sample >= (power-1))
2082         {
2083             length += ++digits * factor * 6;  /* 9,90,900,9000,... (9*2/3 = 6) */
2084             factor = power;
2085             power *= 10;
2086         }
2087         if (sample > --factor)
2088         {
2089             sample -= factor;
2090             factor = (N_word) ( sample / 3 );
2091             factor = (factor << 1) + (sample - (factor * 3));
2092             length += ++digits * factor;
2093         }
2094     }
2095     else length = 1;
2096     string = (charptr) malloc((size_t) length);
2097     if (string == NULL) return(NULL);
2098     start = 0;
2099     comma = false;
2100     target = string;
2101     while ((start < bits) and BitVector_interval_scan_inc(addr,start,&min,&max))
2102     {
2103         start = max + 2;
2104         if (comma) *target++ = (N_char) ',';
2105         if (min == max)
2106         {
2107             target += BIT_VECTOR_int2str(target,min);
2108         }
2109         else
2110         {
2111             if (min+1 == max)
2112             {
2113                 target += BIT_VECTOR_int2str(target,min);
2114                 *target++ = (N_char) ',';
2115                 target += BIT_VECTOR_int2str(target,max);
2116             }
2117             else
2118             {
2119                 target += BIT_VECTOR_int2str(target,min);
2120                 *target++ = (N_char) '-';
2121                 target += BIT_VECTOR_int2str(target,max);
2122             }
2123         }
2124         comma = true;
2125     }
2126     *target = (N_char) '\0';
2127     return(string);
2128 }
2129
2130 ErrCode BitVector_from_Enum(wordptr addr, charptr string)
2131 {
2132     ErrCode error = ErrCode_Ok;
2133     N_word  bits = bits_(addr);
2134     N_word  state = 1;
2135     N_word  token;
2136     N_word  index;
2137     N_word  start = 0;         /* silence compiler warning */
2138
2139     if (bits > 0)
2140     {
2141         BitVector_Empty(addr);
2142         while ((not error) and (state != 0))
2143         {
2144             token = (N_word) *string;
2145             /* separate because isdigit() is likely a macro! */
2146             if (isdigit((int)token) != 0)
2147             {
2148                 string += BIT_VECTOR_str2int(string,&index);
2149                 if (index < bits) token = (N_word) '0';
2150                 else error = ErrCode_Indx;
2151             }
2152             else string++;
2153             if (not error)
2154             switch (state)
2155             {
2156                 case 1:
2157                     switch (token)
2158                     {
2159                         case (N_word) '0':
2160                             state = 2;
2161                             break;
2162                         case (N_word) '\0':
2163                             state = 0;
2164                             break;
2165                         default:
2166                             error = ErrCode_Pars;
2167                             break;
2168                     }
2169                     break;
2170                 case 2:
2171                     switch (token)
2172                     {
2173                         case (N_word) '-':
2174                             start = index;
2175                             state = 3;
2176                             break;
2177                         case (N_word) ',':
2178                             BIT_VECTOR_SET_BIT(addr,index)
2179                             state = 5;
2180                             break;
2181                         case (N_word) '\0':
2182                             BIT_VECTOR_SET_BIT(addr,index)
2183                             state = 0;
2184                             break;
2185                         default:
2186                             error = ErrCode_Pars;
2187                             break;
2188                     }
2189                     break;
2190                 case 3:
2191                     switch (token)
2192                     {
2193                         case (N_word) '0':
2194                             if (start < index)
2195                                 BitVector_Interval_Fill(addr,start,index);
2196                             else if (start == index)
2197                                 BIT_VECTOR_SET_BIT(addr,index)
2198                             else error = ErrCode_Ordr;
2199                             state = 4;
2200                             break;
2201                         default:
2202                             error = ErrCode_Pars;
2203                             break;
2204                     }
2205                     break;
2206                 case 4:
2207                     switch (token)
2208                     {
2209                         case (N_word) ',':
2210                             state = 5;
2211                             break;
2212                         case (N_word) '\0':
2213                             state = 0;
2214                             break;
2215                         default:
2216                             error = ErrCode_Pars;
2217                             break;
2218                     }
2219                     break;
2220                 case 5:
2221                     switch (token)
2222                     {
2223                         case (N_word) '0':
2224                             state = 2;
2225                             break;
2226                         default:
2227                             error = ErrCode_Pars;
2228                             break;
2229                     }
2230                     break;
2231             }
2232         }
2233     }
2234     return(error);
2235 }
2236
2237 void BitVector_Bit_Off(wordptr addr, N_int index)           /* X = X \ {x}   */
2238 {
2239     if (index < bits_(addr)) BIT_VECTOR_CLR_BIT(addr,index)
2240 }
2241
2242 void BitVector_Bit_On(wordptr addr, N_int index)            /* X = X + {x}   */
2243 {
2244     if (index < bits_(addr)) BIT_VECTOR_SET_BIT(addr,index)
2245 }
2246
2247 boolean BitVector_bit_flip(wordptr addr, N_int index)   /* X=(X+{x})\(X*{x}) */
2248 {
2249     N_word mask;
2250
2251     if (index < bits_(addr)) return( BIT_VECTOR_FLP_BIT(addr,index,mask) );
2252     else                     return( false );
2253 }
2254
2255 boolean BitVector_bit_test(wordptr addr, N_int index)       /* {x} in X ?    */
2256 {
2257     if (index < bits_(addr)) return( BIT_VECTOR_TST_BIT(addr,index) );
2258     else                     return( false );
2259 }
2260
2261 void BitVector_Bit_Copy(wordptr addr, N_int index, boolean bit)
2262 {
2263     if (index < bits_(addr))
2264     {
2265         if (bit) BIT_VECTOR_SET_BIT(addr,index)
2266         else     BIT_VECTOR_CLR_BIT(addr,index)
2267     }
2268 }
2269
2270 void BitVector_LSB(wordptr addr, boolean bit)
2271 {
2272     if (bits_(addr) > 0)
2273     {
2274         if (bit) *addr |= LSB;
2275         else     *addr &= NOT LSB;
2276     }
2277 }
2278
2279 void BitVector_MSB(wordptr addr, boolean bit)
2280 {
2281     N_word size = size_(addr);
2282     N_word mask = mask_(addr);
2283
2284     if (size-- > 0)
2285     {
2286         if (bit) *(addr+size) |= mask AND NOT (mask >> 1);
2287         else     *(addr+size) &= NOT mask OR (mask >> 1);
2288     }
2289 }
2290
2291 boolean BitVector_lsb_(wordptr addr)
2292 {
2293     if (size_(addr) > 0) return( (*addr AND LSB) != 0 );
2294     else                 return( false );
2295 }
2296
2297 boolean BitVector_msb_(wordptr addr)
2298 {
2299     N_word size = size_(addr);
2300     N_word mask = mask_(addr);
2301
2302     if (size-- > 0)
2303         return( (*(addr+size) AND (mask AND NOT (mask >> 1))) != 0 );
2304     else
2305         return( false );
2306 }
2307
2308 boolean BitVector_rotate_left(wordptr addr)
2309 {
2310     N_word  size = size_(addr);
2311     N_word  mask = mask_(addr);
2312     N_word  msb;
2313     boolean carry_in;
2314     boolean carry_out = false;
2315
2316     if (size > 0)
2317     {
2318         msb = mask AND NOT (mask >> 1);
2319         carry_in = ((*(addr+size-1) AND msb) != 0);
2320         while (size-- > 1)
2321         {
2322             carry_out = ((*addr AND MSB) != 0);
2323             *addr <<= 1;
2324             if (carry_in) *addr |= LSB;
2325             carry_in = carry_out;
2326             addr++;
2327         }
2328         carry_out = ((*addr AND msb) != 0);
2329         *addr <<= 1;
2330         if (carry_in) *addr |= LSB;
2331         *addr &= mask;
2332     }
2333     return(carry_out);
2334 }
2335
2336 boolean BitVector_rotate_right(wordptr addr)
2337 {
2338     N_word  size = size_(addr);
2339     N_word  mask = mask_(addr);
2340     N_word  msb;
2341     boolean carry_in;
2342     boolean carry_out = false;
2343
2344     if (size > 0)
2345     {
2346         msb = mask AND NOT (mask >> 1);
2347         carry_in = ((*addr AND LSB) != 0);
2348         addr += size-1;
2349         *addr &= mask;
2350         carry_out = ((*addr AND LSB) != 0);
2351         *addr >>= 1;
2352         if (carry_in) *addr |= msb;
2353         carry_in = carry_out;
2354         addr--;
2355         size--;
2356         while (size-- > 0)
2357         {
2358             carry_out = ((*addr AND LSB) != 0);
2359             *addr >>= 1;
2360             if (carry_in) *addr |= MSB;
2361             carry_in = carry_out;
2362             addr--;
2363         }
2364     }
2365     return(carry_out);
2366 }
2367
2368 boolean BitVector_shift_left(wordptr addr, boolean carry_in)
2369 {
2370     N_word  size = size_(addr);
2371     N_word  mask = mask_(addr);
2372     N_word  msb;
2373     boolean carry_out = carry_in;
2374
2375     if (size > 0)
2376     {
2377         msb = mask AND NOT (mask >> 1);
2378         while (size-- > 1)
2379         {
2380             carry_out = ((*addr AND MSB) != 0);
2381             *addr <<= 1;
2382             if (carry_in) *addr |= LSB;
2383             carry_in = carry_out;
2384             addr++;
2385         }
2386         carry_out = ((*addr AND msb) != 0);
2387         *addr <<= 1;
2388         if (carry_in) *addr |= LSB;
2389         *addr &= mask;
2390     }
2391     return(carry_out);
2392 }
2393
2394 boolean BitVector_shift_right(wordptr addr, boolean carry_in)
2395 {
2396     N_word  size = size_(addr);
2397     N_word  mask = mask_(addr);
2398     N_word  msb;
2399     boolean carry_out = carry_in;
2400
2401     if (size > 0)
2402     {
2403         msb = mask AND NOT (mask >> 1);
2404         addr += size-1;
2405         *addr &= mask;
2406         carry_out = ((*addr AND LSB) != 0);
2407         *addr >>= 1;
2408         if (carry_in) *addr |= msb;
2409         carry_in = carry_out;
2410         addr--;
2411         size--;
2412         while (size-- > 0)
2413         {
2414             carry_out = ((*addr AND LSB) != 0);
2415             *addr >>= 1;
2416             if (carry_in) *addr |= MSB;
2417             carry_in = carry_out;
2418             addr--;
2419         }
2420     }
2421     return(carry_out);
2422 }
2423
2424 void BitVector_Move_Left(wordptr addr, N_int bits)
2425 {
2426     N_word count;
2427     N_word words;
2428
2429     if (bits > 0)
2430     {
2431         count = bits AND MODMASK;
2432         words = bits >> LOGBITS;
2433         if (bits >= bits_(addr)) BitVector_Empty(addr);
2434         else
2435         {
2436             while (count-- > 0) BitVector_shift_left(addr,0);
2437             BitVector_Word_Insert(addr,0,words,true);
2438         }
2439     }
2440 }
2441
2442 void BitVector_Move_Right(wordptr addr, N_int bits)
2443 {
2444     N_word count;
2445     N_word words;
2446
2447     if (bits > 0)
2448     {
2449         count = bits AND MODMASK;
2450         words = bits >> LOGBITS;
2451         if (bits >= bits_(addr)) BitVector_Empty(addr);
2452         else
2453         {
2454             while (count-- > 0) BitVector_shift_right(addr,0);
2455             BitVector_Word_Delete(addr,0,words,true);
2456         }
2457     }
2458 }
2459
2460 void BitVector_Insert(wordptr addr, N_int offset, N_int count, boolean clear)
2461 {
2462     N_word bits = bits_(addr);
2463     N_word last;
2464
2465     if ((count > 0) and (offset < bits))
2466     {
2467         last = offset + count;
2468         if (last < bits)
2469         {
2470             BitVector_Interval_Copy(addr,addr,last,offset,(bits-last));
2471         }
2472         else last = bits;
2473         if (clear) BitVector_Interval_Empty(addr,offset,(last-1));
2474     }
2475 }
2476
2477 void BitVector_Delete(wordptr addr, N_int offset, N_int count, boolean clear)
2478 {
2479     N_word bits = bits_(addr);
2480     N_word last;
2481
2482     if ((count > 0) and (offset < bits))
2483     {
2484         last = offset + count;
2485         if (last < bits)
2486         {
2487             BitVector_Interval_Copy(addr,addr,offset,last,(bits-last));
2488         }
2489         else count = bits - offset;
2490         if (clear) BitVector_Interval_Empty(addr,(bits-count),(bits-1));
2491     }
2492 }
2493
2494 boolean BitVector_increment(wordptr addr)                   /* X++           */
2495 {
2496     N_word  size  = size_(addr);
2497     N_word  mask  = mask_(addr);
2498     wordptr last  = addr + size - 1;
2499     boolean carry = true;
2500
2501     if (size > 0)
2502     {
2503         *last |= NOT mask;
2504         while (carry and (size-- > 0))
2505         {
2506             carry = (++(*addr++) == 0);
2507         }
2508         *last &= mask;
2509     }
2510     return(carry);
2511 }
2512
2513 boolean BitVector_decrement(wordptr addr)                   /* X--           */
2514 {
2515     N_word  size  = size_(addr);
2516     N_word  mask  = mask_(addr);
2517     wordptr last  = addr + size - 1;
2518     boolean carry = true;
2519
2520     if (size > 0)
2521     {
2522         *last &= mask;
2523         while (carry and (size-- > 0))
2524         {
2525             carry = (*addr == 0);
2526             --(*addr++);
2527         }
2528         *last &= mask;
2529     }
2530     return(carry);
2531 }
2532
2533 boolean BitVector_compute(wordptr X, wordptr Y, wordptr Z, boolean minus, boolean *carry)
2534 {
2535     N_word size = size_(X);
2536     N_word mask = mask_(X);
2537     N_word vv = 0;
2538     N_word cc;
2539     N_word mm;
2540     N_word yy;
2541     N_word zz;
2542     N_word lo;
2543     N_word hi;
2544
2545     if (size > 0)
2546     {
2547         if (minus) cc = (*carry == 0);
2548         else       cc = (*carry != 0);
2549         /* deal with (size-1) least significant full words first: */
2550         while (--size > 0)
2551         {
2552             yy = *Y++;
2553             if (minus) zz = (N_word) NOT ( Z ? *Z++ : 0 );
2554             else       zz = (N_word)     ( Z ? *Z++ : 0 );
2555             lo = (yy AND LSB) + (zz AND LSB) + cc;
2556             hi = (yy >> 1) + (zz >> 1) + (lo >> 1);
2557             cc = ((hi AND MSB) != 0);
2558             *X++ = (hi << 1) OR (lo AND LSB);
2559         }
2560         /* deal with most significant word (may be used only partially): */
2561         yy = *Y AND mask;
2562         if (minus) zz = (N_word) NOT ( Z ? *Z : 0 );
2563         else       zz = (N_word)     ( Z ? *Z : 0 );
2564         zz &= mask;
2565         if (mask == LSB) /* special case, only one bit used */
2566         {
2567             vv = cc;
2568             lo = yy + zz + cc;
2569             cc = (lo >> 1);
2570             vv ^= cc;
2571             *X = lo AND LSB;
2572         }
2573         else
2574         {
2575             if (NOT mask) /* not all bits are used, but more than one */
2576             {
2577                 mm = (mask >> 1);
2578                 vv = (yy AND mm) + (zz AND mm) + cc;
2579                 mm = mask AND NOT mm;
2580                 lo = yy + zz + cc;
2581                 cc = (lo >> 1);
2582                 vv ^= cc;
2583                 vv &= mm;
2584                 cc &= mm;
2585                 *X = lo AND mask;
2586             }
2587             else /* other special case, all bits are used */
2588             {
2589                 mm = NOT MSB;
2590                 lo = (yy AND mm) + (zz AND mm) + cc;
2591                 vv = lo AND MSB;
2592                 hi = ((yy AND MSB) >> 1) + ((zz AND MSB) >> 1) + (vv >> 1);
2593                 cc = hi AND MSB;
2594                 vv ^= cc;
2595                 *X = (hi << 1) OR (lo AND mm);
2596             }
2597         }
2598         if (minus) *carry = (cc == 0);
2599         else       *carry = (cc != 0);
2600     }
2601     return(vv != 0);
2602 }
2603
2604 boolean BitVector_add(wordptr X, wordptr Y, wordptr Z, boolean *carry)
2605 {
2606     return(BitVector_compute(X,Y,Z,false,carry));
2607 }
2608
2609 boolean BitVector_sub(wordptr X, wordptr Y, wordptr Z, boolean *carry)
2610 {
2611     return(BitVector_compute(X,Y,Z,true,carry));
2612 }
2613
2614 boolean BitVector_inc(wordptr X, wordptr Y)
2615 {
2616     boolean carry = true;
2617
2618     return(BitVector_compute(X,Y,NULL,false,&carry));
2619 }
2620
2621 boolean BitVector_dec(wordptr X, wordptr Y)
2622 {
2623     boolean carry = true;
2624
2625     return(BitVector_compute(X,Y,NULL,true,&carry));
2626 }
2627
2628 void BitVector_Negate(wordptr X, wordptr Y)
2629 {
2630     N_word  size  = size_(X);
2631     N_word  mask  = mask_(X);
2632     boolean carry = true;
2633
2634     if (size > 0)
2635     {
2636         while (size-- > 0)
2637         {
2638             *X = NOT *Y++;
2639             if (carry)
2640             {
2641                 carry = (++(*X) == 0);
2642             }
2643             X++;
2644         }
2645         *(--X) &= mask;
2646     }
2647 }
2648
2649 void BitVector_Absolute(wordptr X, wordptr Y)
2650 {
2651     N_word size = size_(Y);
2652     N_word mask = mask_(Y);
2653
2654     if (size > 0)
2655     {
2656         if (*(Y+size-1) AND (mask AND NOT (mask >> 1))) BitVector_Negate(X,Y);
2657         else                                             BitVector_Copy(X,Y);
2658     }
2659 }
2660
2661 Z_int BitVector_Sign(wordptr addr)
2662 {
2663     N_word  size = size_(addr);
2664     N_word  mask = mask_(addr);
2665     wordptr last = addr + size - 1;
2666     boolean r    = true;
2667
2668     if (size > 0)
2669     {
2670         *last &= mask;
2671         while (r and (size-- > 0)) r = ( *addr++ == 0 );
2672     }
2673     if (r) return((Z_int) 0);
2674     else
2675     {
2676         if (*last AND (mask AND NOT (mask >> 1))) return((Z_int) -1);
2677         else                                      return((Z_int)  1);
2678     }
2679 }
2680
2681 ErrCode BitVector_Mul_Pos(wordptr X, wordptr Y, wordptr Z, boolean strict)
2682 {
2683     N_word  mask;
2684     N_word  limit;
2685     N_word  count;
2686     Z_long  last;
2687     wordptr sign;
2688     boolean carry;
2689     boolean overflow;
2690     boolean ok = true;
2691
2692     /*
2693        Requirements:
2694          -  X, Y and Z must be distinct
2695          -  X and Y must have equal sizes (whereas Z may be any size!)
2696          -  Z should always contain the SMALLER of the two factors Y and Z
2697        Constraints:
2698          -  The contents of Y (and of X, of course) are destroyed
2699             (only Z is preserved!)
2700     */
2701
2702     if ((X == Y) or (X == Z) or (Y == Z)) return(ErrCode_Same);
2703     if (bits_(X) != bits_(Y)) return(ErrCode_Size);
2704     BitVector_Empty(X);
2705     if (BitVector_is_empty(Y)) return(ErrCode_Ok); /* exit also taken if bits_(Y)==0 */
2706     if ((last = Set_Max(Z)) < 0L) return(ErrCode_Ok);
2707     limit = (N_word) last;
2708     sign = Y + size_(Y) - 1;
2709     mask = mask_(Y);
2710     *sign &= mask;
2711     mask &= NOT (mask >> 1);
2712     for ( count = 0; (ok and (count <= limit)); count++ )
2713     {
2714         if ( BIT_VECTOR_TST_BIT(Z,count) )
2715         {
2716             carry = false;
2717             overflow = BitVector_compute(X,X,Y,false,&carry);
2718             if (strict) ok = not (carry or overflow);
2719             else        ok = not  carry;
2720         }
2721         if (ok and (count < limit))
2722         {
2723             carry = BitVector_shift_left(Y,0);
2724             if (strict)
2725             {
2726                 overflow = ((*sign AND mask) != 0);
2727                 ok = not (carry or overflow);
2728             }
2729             else ok = not carry;
2730         }
2731     }
2732     if (ok) return(ErrCode_Ok); else return(ErrCode_Ovfl);
2733 }
2734
2735 ErrCode BitVector_Multiply(wordptr X, wordptr Y, wordptr Z)
2736 {
2737     ErrCode error = ErrCode_Ok;
2738     N_word  bit_x = bits_(X);
2739     N_word  bit_y = bits_(Y);
2740     N_word  bit_z = bits_(Z);
2741     N_word  size;
2742     N_word  mask;
2743     N_word  msb;
2744     wordptr ptr_y;
2745     wordptr ptr_z;
2746     boolean sgn_x;
2747     boolean sgn_y;
2748     boolean sgn_z;
2749     boolean zero;
2750     wordptr A;
2751     wordptr B;
2752
2753     /*
2754        Requirements:
2755          -  Y and Z must have equal sizes
2756          -  X must have at least the same size as Y and Z but may be larger (!)
2757        Features:
2758          -  The contents of Y and Z are preserved
2759          -  X may be identical with Y or Z (or both!)
2760             (in-place multiplication is possible!)
2761     */
2762
2763     if ((bit_y != bit_z) or (bit_x < bit_y)) return(ErrCode_Size);
2764     if (BitVector_is_empty(Y) or BitVector_is_empty(Z))
2765     {
2766         BitVector_Empty(X);
2767     }
2768     else
2769     {
2770         A = BitVector_Create(bit_y,false);
2771         if (A == NULL) return(ErrCode_Null);
2772         B = BitVector_Create(bit_z,false);
2773         if (B == NULL) { BitVector_Destroy(A); return(ErrCode_Null); }
2774         size  = size_(Y);
2775         mask  = mask_(Y);
2776         msb   = (mask AND NOT (mask >> 1));
2777         sgn_y = (((*(Y+size-1) &= mask) AND msb) != 0);
2778         sgn_z = (((*(Z+size-1) &= mask) AND msb) != 0);
2779         sgn_x = sgn_y XOR sgn_z;
2780         if (sgn_y) BitVector_Negate(A,Y); else BitVector_Copy(A,Y);
2781         if (sgn_z) BitVector_Negate(B,Z); else BitVector_Copy(B,Z);
2782         ptr_y = A + size;
2783         ptr_z = B + size;
2784         zero = true;
2785         while (zero and (size-- > 0))
2786         {
2787             zero &= (*(--ptr_y) == 0);
2788             zero &= (*(--ptr_z) == 0);
2789         }
2790         if (*ptr_y > *ptr_z)
2791         {
2792             if (bit_x > bit_y)
2793             {
2794                 A = BitVector_Resize(A,bit_x);
2795                 if (A == NULL) { BitVector_Destroy(B); return(ErrCode_Null); }
2796             }
2797             error = BitVector_Mul_Pos(X,A,B,true);
2798         }
2799         else
2800         {
2801             if (bit_x > bit_z)
2802             {
2803                 B = BitVector_Resize(B,bit_x);
2804                 if (B == NULL) { BitVector_Destroy(A); return(ErrCode_Null); }
2805             }
2806             error = BitVector_Mul_Pos(X,B,A,true);
2807         }
2808         if ((not error) and sgn_x) BitVector_Negate(X,X);
2809         BitVector_Destroy(A);
2810         BitVector_Destroy(B);
2811     }
2812     return(error);
2813 }
2814
2815 ErrCode BitVector_Div_Pos(wordptr Q, wordptr X, wordptr Y, wordptr R)
2816 {
2817     N_word  bits = bits_(Q);
2818     N_word  mask;
2819     wordptr addr;
2820     Z_long  last;
2821     boolean flag;
2822     boolean copy = false; /* flags whether valid rest is in R (0) or X (1) */
2823
2824     /*
2825        Requirements:
2826          -  All bit vectors must have equal sizes
2827          -  Q, X, Y and R must all be distinct bit vectors
2828          -  Y must be non-zero (of course!)
2829        Constraints:
2830          -  The contents of X (and Q and R, of course) are destroyed
2831             (only Y is preserved!)
2832     */
2833
2834     if ((bits != bits_(X)) or (bits != bits_(Y)) or (bits != bits_(R)))
2835         return(ErrCode_Size);
2836     if ((Q == X) or (Q == Y) or (Q == R) or (X == Y) or (X == R) or (Y == R))
2837         return(ErrCode_Same);
2838     if (BitVector_is_empty(Y))
2839         return(ErrCode_Zero);
2840
2841     BitVector_Empty(R);
2842     BitVector_Copy(Q,X);
2843     if ((last = Set_Max(Q)) < 0L) return(ErrCode_Ok);
2844     bits = (N_word) ++last;
2845     while (bits-- > 0)
2846     {
2847         addr = Q + (bits >> LOGBITS);
2848         mask = BITMASKTAB[bits AND MODMASK];
2849         flag = ((*addr AND mask) != 0);
2850         if (copy)
2851         {
2852             BitVector_shift_left(X,flag);
2853             flag = false;
2854             BitVector_compute(R,X,Y,true,&flag);
2855         }
2856         else
2857         {
2858             BitVector_shift_left(R,flag);
2859             flag = false;
2860             BitVector_compute(X,R,Y,true,&flag);
2861         }
2862         if (flag) *addr &= NOT mask;
2863         else
2864         {
2865             *addr |= mask;
2866             copy = not copy;
2867         }
2868     }
2869     if (copy) BitVector_Copy(R,X);
2870     return(ErrCode_Ok);
2871 }
2872
2873 ErrCode BitVector_Divide(wordptr Q, wordptr X, wordptr Y, wordptr R)
2874 {
2875     ErrCode error = ErrCode_Ok;
2876     N_word  bits = bits_(Q);
2877     N_word  size = size_(Q);
2878     N_word  mask = mask_(Q);
2879     N_word  msb = (mask AND NOT (mask >> 1));
2880     boolean sgn_q;
2881     boolean sgn_x;
2882     boolean sgn_y;
2883     wordptr A;
2884     wordptr B;
2885
2886     /*
2887        Requirements:
2888          -  All bit vectors must have equal sizes
2889          -  Q and R must be two distinct bit vectors
2890          -  Y must be non-zero (of course!)
2891        Features:
2892          -  The contents of X and Y are preserved
2893          -  Q may be identical with X or Y (or both)
2894             (in-place division is possible!)
2895          -  R may be identical with X or Y (or both)
2896             (but not identical with Q!)
2897     */
2898
2899     if ((bits != bits_(X)) or (bits != bits_(Y)) or (bits != bits_(R)))
2900         return(ErrCode_Size);
2901     if (Q == R)
2902         return(ErrCode_Same);
2903     if (BitVector_is_empty(Y))
2904         return(ErrCode_Zero);
2905
2906     if (BitVector_is_empty(X))
2907     {
2908         BitVector_Empty(Q);
2909         BitVector_Empty(R);
2910     }
2911     else
2912     {
2913         A = BitVector_Create(bits,false);
2914         if (A == NULL) return(ErrCode_Null);
2915         B = BitVector_Create(bits,false);
2916         if (B == NULL) { BitVector_Destroy(A); return(ErrCode_Null); }
2917         size--;
2918         sgn_x = (((*(X+size) &= mask) AND msb) != 0);
2919         sgn_y = (((*(Y+size) &= mask) AND msb) != 0);
2920         sgn_q = sgn_x XOR sgn_y;
2921         if (sgn_x) BitVector_Negate(A,X); else BitVector_Copy(A,X);
2922         if (sgn_y) BitVector_Negate(B,Y); else BitVector_Copy(B,Y);
2923         if (not (error = BitVector_Div_Pos(Q,A,B,R)))
2924         {
2925             if (sgn_q) BitVector_Negate(Q,Q);
2926             if (sgn_x) BitVector_Negate(R,R);
2927         }
2928         BitVector_Destroy(A);
2929         BitVector_Destroy(B);
2930     }
2931     return(error);
2932 }
2933
2934 ErrCode BitVector_GCD(wordptr X, wordptr Y, wordptr Z)
2935 {
2936     ErrCode error = ErrCode_Ok;
2937     N_word  bits = bits_(X);
2938     N_word  size = size_(X);
2939     N_word  mask = mask_(X);
2940     N_word  msb = (mask AND NOT (mask >> 1));
2941     boolean sgn_a;
2942     boolean sgn_b;
2943     boolean sgn_r;
2944     wordptr Q;
2945     wordptr R;
2946     wordptr A;
2947     wordptr B;
2948     wordptr T;
2949
2950     /*
2951        Requirements:
2952          -  All bit vectors must have equal sizes
2953        Features:
2954          -  The contents of Y and Z are preserved
2955          -  X may be identical with Y or Z (or both)
2956             (in-place is possible!)
2957          -  GCD(0,z) == GCD(z,0) == z
2958          -  negative values are handled correctly
2959     */
2960
2961     if ((bits != bits_(Y)) or (bits != bits_(Z))) return(ErrCode_Size);
2962     if (BitVector_is_empty(Y))
2963     {
2964         if (X != Z) BitVector_Copy(X,Z);
2965         return(ErrCode_Ok);
2966     }
2967     if (BitVector_is_empty(Z))
2968     {
2969         if (X != Y) BitVector_Copy(X,Y);
2970         return(ErrCode_Ok);
2971     }
2972     Q = BitVector_Create(bits,false);
2973     if (Q == NULL)
2974     {
2975         return(ErrCode_Null);
2976     }
2977     R = BitVector_Create(bits,false);
2978     if (R == NULL)
2979     {
2980         BitVector_Destroy(Q);
2981         return(ErrCode_Null);
2982     }
2983     A = BitVector_Create(bits,false);
2984     if (A == NULL)
2985     {
2986         BitVector_Destroy(Q);
2987         BitVector_Destroy(R);
2988         return(ErrCode_Null);
2989     }
2990     B = BitVector_Create(bits,false);
2991     if (B == NULL)
2992     {
2993         BitVector_Destroy(Q);
2994         BitVector_Destroy(R);
2995         BitVector_Destroy(A);
2996         return(ErrCode_Null);
2997     }
2998     size--;
2999     sgn_a = (((*(Y+size) &= mask) AND msb) != 0);
3000     sgn_b = (((*(Z+size) &= mask) AND msb) != 0);
3001     if (sgn_a) BitVector_Negate(A,Y); else BitVector_Copy(A,Y);
3002     if (sgn_b) BitVector_Negate(B,Z); else BitVector_Copy(B,Z);
3003     while (not error)
3004     {
3005         if (not (error = BitVector_Div_Pos(Q,A,B,R)))
3006         {
3007             if (BitVector_is_empty(R)) break;
3008             T = A; sgn_r = sgn_a;
3009             A = B; sgn_a = sgn_b;
3010             B = R; sgn_b = sgn_r;
3011             R = T;
3012         }
3013     }
3014     if (not error)
3015     {
3016         if (sgn_b) BitVector_Negate(X,B); else BitVector_Copy(X,B);
3017     }
3018     BitVector_Destroy(Q);
3019     BitVector_Destroy(R);
3020     BitVector_Destroy(A);
3021     BitVector_Destroy(B);
3022     return(error);
3023 }
3024
3025 ErrCode BitVector_GCD2(wordptr U, wordptr V, wordptr W, wordptr X, wordptr Y)
3026 {
3027     ErrCode error = ErrCode_Ok;
3028     N_word  bits = bits_(U);
3029     N_word  size = size_(U);
3030     N_word  mask = mask_(U);
3031     N_word  msb = (mask AND NOT (mask >> 1));
3032     boolean minus;
3033     boolean carry;
3034     boolean sgn_q;
3035     boolean sgn_r;
3036     boolean sgn_a;
3037     boolean sgn_b;
3038     boolean sgn_x;
3039     boolean sgn_y;
3040     listptr L;
3041     wordptr Q;
3042     wordptr R;
3043     wordptr A;
3044     wordptr B;
3045     wordptr T;
3046     wordptr X1;
3047     wordptr X2;
3048     wordptr X3;
3049     wordptr Y1;
3050     wordptr Y2;
3051     wordptr Y3;
3052     wordptr Z;
3053
3054     /*
3055        Requirements:
3056          -  All bit vectors must have equal sizes
3057          -  U, V, and W must all be distinct bit vectors
3058        Features:
3059          -  The contents of X and Y are preserved
3060          -  U, V and W may be identical with X or Y (or both,
3061             provided that U, V and W are mutually distinct)
3062             (i.e., in-place is possible!)
3063          -  GCD(0,z) == GCD(z,0) == z
3064          -  negative values are handled correctly
3065     */
3066
3067     if ((bits != bits_(V)) or
3068         (bits != bits_(W)) or
3069         (bits != bits_(X)) or
3070         (bits != bits_(Y)))
3071     {
3072         return(ErrCode_Size);
3073     }
3074     if ((U == V) or (U == W) or (V == W))
3075     {
3076         return(ErrCode_Same);
3077     }
3078     if (BitVector_is_empty(X))
3079     {
3080         if (U != Y) BitVector_Copy(U,Y);
3081         BitVector_Empty(V);
3082         BitVector_Empty(W);
3083         *W = 1;
3084         return(ErrCode_Ok);
3085     }
3086     if (BitVector_is_empty(Y))
3087     {
3088         if (U != X) BitVector_Copy(U,X);
3089         BitVector_Empty(V);
3090         BitVector_Empty(W);
3091         *V = 1;
3092         return(ErrCode_Ok);
3093     }
3094     if ((L = BitVector_Create_List(bits,false,11)) == NULL)
3095     {
3096         return(ErrCode_Null);
3097     }
3098     Q  = L[0];
3099     R  = L[1];
3100     A  = L[2];
3101     B  = L[3];
3102     X1 = L[4];
3103     X2 = L[5];
3104     X3 = L[6];
3105     Y1 = L[7];
3106     Y2 = L[8];
3107     Y3 = L[9];
3108     Z  = L[10];
3109     size--;
3110     sgn_a = (((*(X+size) &= mask) AND msb) != 0);
3111     sgn_b = (((*(Y+size) &= mask) AND msb) != 0);
3112     if (sgn_a) BitVector_Negate(A,X); else BitVector_Copy(A,X);
3113     if (sgn_b) BitVector_Negate(B,Y); else BitVector_Copy(B,Y);
3114     BitVector_Empty(X1);
3115     BitVector_Empty(X2);
3116     *X1 = 1;
3117     BitVector_Empty(Y1);
3118     BitVector_Empty(Y2);
3119     *Y2 = 1;
3120     sgn_x = false;
3121     sgn_y = false;
3122     while (not error)
3123     {
3124         if ((error = BitVector_Div_Pos(Q,A,B,R)))
3125         {
3126             break;
3127         }
3128         if (BitVector_is_empty(R))
3129         {
3130             break;
3131         }
3132         sgn_q = sgn_a XOR sgn_b;
3133
3134         if (sgn_x) BitVector_Negate(Z,X2); else BitVector_Copy(Z,X2);
3135         if ((error = BitVector_Mul_Pos(X3,Z,Q,true)))
3136         {
3137             break;
3138         }
3139         minus = not (sgn_x XOR sgn_q);
3140         carry = 0;
3141         if (BitVector_compute(X3,X1,X3,minus,&carry))
3142         {
3143             error = ErrCode_Ovfl;
3144             break;
3145         }
3146         sgn_x = (((*(X3+size) &= mask) AND msb) != 0);
3147
3148         if (sgn_y) BitVector_Negate(Z,Y2); else BitVector_Copy(Z,Y2);
3149         if ((error = BitVector_Mul_Pos(Y3,Z,Q,true)))
3150         {
3151             break;
3152         }
3153         minus = not (sgn_y XOR sgn_q);
3154         carry = 0;
3155         if (BitVector_compute(Y3,Y1,Y3,minus,&carry))
3156         {
3157             error = ErrCode_Ovfl;
3158             break;
3159         }
3160         sgn_y = (((*(Y3+size) &= mask) AND msb) != 0);
3161
3162         T = A; sgn_r = sgn_a;
3163         A = B; sgn_a = sgn_b;
3164         B = R; sgn_b = sgn_r;
3165         R = T;
3166
3167         T = X1;
3168         X1 = X2;
3169         X2 = X3;
3170         X3 = T;
3171
3172         T = Y1;
3173         Y1 = Y2;
3174         Y2 = Y3;
3175         Y3 = T;
3176     }
3177     if (not error)
3178     {
3179         if (sgn_b) BitVector_Negate(U,B); else BitVector_Copy(U,B);
3180         BitVector_Copy(V,X2);
3181         BitVector_Copy(W,Y2);
3182     }
3183     BitVector_Destroy_List(L,11);
3184     return(error);
3185 }
3186
3187 ErrCode BitVector_Power(wordptr X, wordptr Y, wordptr Z)
3188 {
3189     ErrCode error = ErrCode_Ok;
3190     N_word  bits  = bits_(X);
3191     boolean first = true;
3192     Z_long  last;
3193     N_word  limit;
3194     N_word  count;
3195     wordptr T;
3196
3197     /*
3198        Requirements:
3199          -  X must have at least the same size as Y but may be larger (!)
3200          -  X may not be identical with Z
3201          -  Z must be positive
3202        Features:
3203          -  The contents of Y and Z are preserved
3204     */
3205
3206     if (X == Z) return(ErrCode_Same);
3207     if (bits < bits_(Y)) return(ErrCode_Size);
3208     if (BitVector_msb_(Z)) return(ErrCode_Expo);
3209     if ((last = Set_Max(Z)) < 0L)
3210     {
3211         if (bits < 2) return(ErrCode_Ovfl);
3212         BitVector_Empty(X);
3213         *X |= LSB;
3214         return(ErrCode_Ok);                             /* anything ^ 0 == 1 */
3215     }
3216     if (BitVector_is_empty(Y))
3217     {
3218         if (X != Y) BitVector_Empty(X);
3219         return(ErrCode_Ok);                    /* 0 ^ anything not zero == 0 */
3220     }
3221     T = BitVector_Create(bits,false);
3222     if (T == NULL) return(ErrCode_Null);
3223     limit = (N_word) last;
3224     for ( count = 0; ((!error) and (count <= limit)); count++ )
3225     {
3226         if ( BIT_VECTOR_TST_BIT(Z,count) )
3227         {
3228             if (first)
3229             {
3230                 first = false;
3231                 if (count) {             BitVector_Copy(X,T); }
3232                 else       { if (X != Y) BitVector_Copy(X,Y); }
3233             }
3234             else error = BitVector_Multiply(X,T,X); /* order important because T > X */
3235         }
3236         if ((!error) and (count < limit))
3237         {
3238             if (count) error = BitVector_Multiply(T,T,T);
3239             else       error = BitVector_Multiply(T,Y,Y);
3240         }
3241     }
3242     BitVector_Destroy(T);
3243     return(error);
3244 }
3245
3246 void BitVector_Block_Store(wordptr addr, charptr buffer, N_int length)
3247 {
3248     N_word  size = size_(addr);
3249     N_word  mask = mask_(addr);
3250     N_word  value;
3251     N_word  count;
3252
3253     /* provide translation for independence of endian-ness: */
3254     if (size > 0)
3255     {
3256         while (size-- > 0)
3257         {
3258             value = 0;
3259             for ( count = 0; (length > 0) and (count < BITS); count += 8 )
3260             {
3261                 value |= (((N_word) *buffer++) << count); length--;
3262             }
3263             *addr++ = value;
3264         }
3265         *(--addr) &= mask;
3266     }
3267 }
3268
3269 charptr BitVector_Block_Read(wordptr addr, N_intptr length)
3270 {
3271     N_word  size = size_(addr);
3272     N_word  value;
3273     N_word  count;
3274     charptr buffer;
3275     charptr target;
3276
3277     /* provide translation for independence of endian-ness: */
3278     *length = size << FACTOR;
3279     buffer = (charptr) malloc((size_t) ((*length)+1));
3280     if (buffer == NULL) return(NULL);
3281     target = buffer;
3282     if (size > 0)
3283     {
3284         *(addr+size-1) &= mask_(addr);
3285         while (size-- > 0)
3286         {
3287             value = *addr++;
3288             count = BITS >> 3;
3289             while (count-- > 0)
3290             {
3291                 *target++ = (N_char) (value AND 0x00FF);
3292                 if (count > 0) value >>= 8;
3293             }
3294         }
3295     }
3296     *target = (N_char) '\0';
3297     return(buffer);
3298 }
3299
3300 void BitVector_Word_Store(wordptr addr, N_int offset, N_int value)
3301 {
3302     N_word size = size_(addr);
3303
3304     if (size > 0)
3305     {
3306         if (offset < size) *(addr+offset) = value;
3307         *(addr+size-1) &= mask_(addr);
3308     }
3309 }
3310
3311 N_int BitVector_Word_Read(wordptr addr, N_int offset)
3312 {
3313     N_word size = size_(addr);
3314
3315     if (size > 0)
3316     {
3317         *(addr+size-1) &= mask_(addr);
3318         if (offset < size) return( *(addr+offset) );
3319     }
3320     return( (N_int) 0 );
3321 }
3322
3323 void BitVector_Word_Insert(wordptr addr, N_int offset, N_int count,
3324                            boolean clear)
3325 {
3326     N_word  size = size_(addr);
3327     N_word  mask = mask_(addr);
3328     wordptr last = addr+size-1;
3329
3330     if (size > 0)
3331     {
3332         *last &= mask;
3333         if (offset > size) offset = size;
3334         BIT_VECTOR_ins_words(addr+offset,size-offset,count,clear);
3335         *last &= mask;
3336     }
3337 }
3338
3339 void BitVector_Word_Delete(wordptr addr, N_int offset, N_int count,
3340                            boolean clear)
3341 {
3342     N_word  size = size_(addr);
3343     N_word  mask = mask_(addr);
3344     wordptr last = addr+size-1;
3345
3346     if (size > 0)
3347     {
3348         *last &= mask;
3349         if (offset > size) offset = size;
3350         BIT_VECTOR_del_words(addr+offset,size-offset,count,clear);
3351         *last &= mask;
3352     }
3353 }
3354
3355 void BitVector_Chunk_Store(wordptr addr, N_int chunksize, N_int offset,
3356                            N_long value)
3357 {
3358     N_word bits = bits_(addr);
3359     N_word mask;
3360     N_word temp;
3361
3362     if ((chunksize > 0) and (offset < bits))
3363     {
3364         if (chunksize > LONGBITS) chunksize = LONGBITS;
3365         if ((offset + chunksize) > bits) chunksize = bits - offset;
3366         addr += offset >> LOGBITS;
3367         offset &= MODMASK;
3368         while (chunksize > 0)
3369         {
3370             mask = (N_word) (~0L << offset);
3371             bits = offset + chunksize;
3372             if (bits < BITS)
3373             {
3374                 mask &= (N_word) ~(~0L << bits);
3375                 bits = chunksize;
3376             }
3377             else bits = BITS - offset;
3378             temp = (N_word) (value << offset);
3379             temp &= mask;
3380             *addr &= NOT mask;
3381             *addr++ |= temp;
3382             value >>= bits;
3383             chunksize -= bits;
3384             offset = 0;
3385         }
3386     }
3387 }
3388
3389 N_long BitVector_Chunk_Read(wordptr addr, N_int chunksize, N_int offset)
3390 {
3391     N_word bits = bits_(addr);
3392     N_word chunkbits = 0;
3393     N_long value = 0L;
3394     N_long temp;
3395     N_word mask;
3396
3397     if ((chunksize > 0) and (offset < bits))
3398     {
3399         if (chunksize > LONGBITS) chunksize = LONGBITS;
3400         if ((offset + chunksize) > bits) chunksize = bits - offset;
3401         addr += offset >> LOGBITS;
3402         offset &= MODMASK;
3403         while (chunksize > 0)
3404         {
3405             bits = offset + chunksize;
3406             if (bits < BITS)
3407             {
3408                 mask = (N_word) ~(~0L << bits);
3409                 bits = chunksize;
3410             }
3411             else
3412             {
3413                 mask = (N_word) ~0L;
3414                 bits = BITS - offset;
3415             }
3416             temp = (N_long) ((*addr++ AND mask) >> offset);
3417             value |= temp << chunkbits;
3418             chunkbits += bits;
3419             chunksize -= bits;
3420             offset = 0;
3421         }
3422     }
3423     return(value);
3424 }
3425
3426     /*******************/
3427     /* set operations: */
3428     /*******************/
3429
3430 void Set_Union(wordptr X, wordptr Y, wordptr Z)             /* X = Y + Z     */
3431 {
3432     N_word bits = bits_(X);
3433     N_word size = size_(X);
3434     N_word mask = mask_(X);
3435
3436     if ((size > 0) and (bits == bits_(Y)) and (bits == bits_(Z)))
3437     {
3438         while (size-- > 0) *X++ = *Y++ OR *Z++;
3439         *(--X) &= mask;
3440     }
3441 }
3442
3443 void Set_Intersection(wordptr X, wordptr Y, wordptr Z)      /* X = Y * Z     */
3444 {
3445     N_word bits = bits_(X);
3446     N_word size = size_(X);
3447     N_word mask = mask_(X);
3448
3449     if ((size > 0) and (bits == bits_(Y)) and (bits == bits_(Z)))
3450     {
3451         while (size-- > 0) *X++ = *Y++ AND *Z++;
3452         *(--X) &= mask;
3453     }
3454 }
3455
3456 void Set_Difference(wordptr X, wordptr Y, wordptr Z)        /* X = Y \ Z     */
3457 {
3458     N_word bits = bits_(X);
3459     N_word size = size_(X);
3460     N_word mask = mask_(X);
3461
3462     if ((size > 0) and (bits == bits_(Y)) and (bits == bits_(Z)))
3463     {
3464         while (size-- > 0) *X++ = *Y++ AND NOT *Z++;
3465         *(--X) &= mask;
3466     }
3467 }
3468
3469 void Set_ExclusiveOr(wordptr X, wordptr Y, wordptr Z)       /* X=(Y+Z)\(Y*Z) */
3470 {
3471     N_word bits = bits_(X);
3472     N_word size = size_(X);
3473     N_word mask = mask_(X);
3474
3475     if ((size > 0) and (bits == bits_(Y)) and (bits == bits_(Z)))
3476     {
3477         while (size-- > 0) *X++ = *Y++ XOR *Z++;
3478         *(--X) &= mask;
3479     }
3480 }
3481
3482 void Set_Complement(wordptr X, wordptr Y)                   /* X = ~Y        */
3483 {
3484     N_word size = size_(X);
3485     N_word mask = mask_(X);
3486
3487     if ((size > 0) and (bits_(X) == bits_(Y)))
3488     {
3489         while (size-- > 0) *X++ = NOT *Y++;
3490         *(--X) &= mask;
3491     }
3492 }
3493
3494     /******************/
3495     /* set functions: */
3496     /******************/
3497
3498 boolean Set_subset(wordptr X, wordptr Y)                    /* X subset Y ?  */
3499 {
3500     N_word size = size_(X);
3501     boolean r = false;
3502
3503     if ((size > 0) and (bits_(X) == bits_(Y)))
3504     {
3505         r = true;
3506         while (r and (size-- > 0)) r = ((*X++ AND NOT *Y++) == 0);
3507     }
3508     return(r);
3509 }
3510
3511 N_int Set_Norm(wordptr addr)                                /* = | X |       */
3512 {
3513     byteptr byte;
3514     N_word  bytes;
3515     N_int   n;
3516
3517     byte = (byteptr) addr;
3518     bytes = size_(addr) << FACTOR;
3519     n = 0;
3520     while (bytes-- > 0)
3521     {
3522         n += BitVector_BYTENORM[*byte++];
3523     }
3524     return(n);
3525 }
3526
3527 N_int Set_Norm2(wordptr addr)                               /* = | X |       */
3528 {
3529     N_word  size = size_(addr);
3530     N_word  w0,w1;
3531     N_int   n,k;
3532
3533     n = 0;
3534     while (size-- > 0)
3535     {
3536         k = 0;
3537         w1 = NOT (w0 = *addr++);
3538         while (w0 and w1)
3539         {
3540             w0 &= w0 - 1;
3541             w1 &= w1 - 1;
3542             k++;
3543         }
3544         if (w0 == 0) n += k;
3545         else         n += BITS - k;
3546     }
3547     return(n);
3548 }
3549
3550 N_int Set_Norm3(wordptr addr)                               /* = | X |       */
3551 {
3552     N_word  size  = size_(addr);
3553     N_int   count = 0;
3554     N_word  c;
3555
3556     while (size-- > 0)
3557     {
3558         c = *addr++;
3559         while (c)
3560         {
3561             c &= c - 1;
3562             count++;
3563         }
3564     }
3565     return(count);
3566 }
3567
3568 Z_long Set_Min(wordptr addr)                                /* = min(X)      */
3569 {
3570     boolean empty = true;
3571     N_word  size  = size_(addr);
3572     N_word  i     = 0;
3573     N_word  c     = 0;         /* silence compiler warning */
3574
3575     while (empty and (size-- > 0))
3576     {
3577         if ((c = *addr++)) empty = false; else i++;
3578     }
3579     if (empty) return((Z_long) LONG_MAX);                  /* plus infinity  */
3580     i <<= LOGBITS;
3581     while (not (c AND LSB))
3582     {
3583         c >>= 1;
3584         i++;
3585     }
3586     return((Z_long) i);
3587 }
3588
3589 Z_long Set_Max(wordptr addr)                                /* = max(X)      */
3590 {
3591     boolean empty = true;
3592     N_word  size  = size_(addr);
3593     N_word  i     = size;
3594     N_word  c     = 0;         /* silence compiler warning */
3595
3596     addr += size-1;
3597     while (empty and (size-- > 0))
3598     {
3599         if ((c = *addr--)) empty = false; else i--;
3600     }
3601     if (empty) return((Z_long) LONG_MIN);                  /* minus infinity */
3602     i <<= LOGBITS;
3603     while (not (c AND MSB))
3604     {
3605         c <<= 1;
3606         i--;
3607     }
3608     return((Z_long) --i);
3609 }
3610
3611     /**********************************/
3612     /* matrix-of-booleans operations: */
3613     /**********************************/
3614
3615 void Matrix_Multiplication(wordptr X, N_int rowsX, N_int colsX,
3616                            wordptr Y, N_int rowsY, N_int colsY,
3617                            wordptr Z, N_int rowsZ, N_int colsZ)
3618 {
3619     N_word i;
3620     N_word j;
3621     N_word k;
3622     N_word indxX;
3623     N_word indxY;
3624     N_word indxZ;
3625     N_word termX;
3626     N_word termY;
3627     N_word sum;
3628
3629   if ((colsY == rowsZ) and (rowsX == rowsY) and (colsX == colsZ) and
3630       (bits_(X) == rowsX*colsX) and
3631       (bits_(Y) == rowsY*colsY) and
3632       (bits_(Z) == rowsZ*colsZ))
3633   {
3634     for ( i = 0; i < rowsY; i++ )
3635     {
3636         termX = i * colsX;
3637         termY = i * colsY;
3638         for ( j = 0; j < colsZ; j++ )
3639         {
3640             indxX = termX + j;
3641             sum = 0;
3642             for ( k = 0; k < colsY; k++ )
3643             {
3644                 indxY = termY + k;
3645                 indxZ = k * colsZ + j;
3646                 if ( BIT_VECTOR_TST_BIT(Y,indxY) &&
3647                      BIT_VECTOR_TST_BIT(Z,indxZ) ) sum ^= 1;
3648             }
3649             if (sum) BIT_VECTOR_SET_BIT(X,indxX)
3650             else     BIT_VECTOR_CLR_BIT(X,indxX)
3651         }
3652     }
3653   }
3654 }
3655
3656 void Matrix_Product(wordptr X, N_int rowsX, N_int colsX,
3657                     wordptr Y, N_int rowsY, N_int colsY,
3658                     wordptr Z, N_int rowsZ, N_int colsZ)
3659 {
3660     N_word i;
3661     N_word j;
3662     N_word k;
3663     N_word indxX;
3664     N_word indxY;
3665     N_word indxZ;
3666     N_word termX;
3667     N_word termY;
3668     N_word sum;
3669
3670   if ((colsY == rowsZ) and (rowsX == rowsY) and (colsX == colsZ) and
3671       (bits_(X) == rowsX*colsX) and
3672       (bits_(Y) == rowsY*colsY) and
3673       (bits_(Z) == rowsZ*colsZ))
3674   {
3675     for ( i = 0; i < rowsY; i++ )
3676     {
3677         termX = i * colsX;
3678         termY = i * colsY;
3679         for ( j = 0; j < colsZ; j++ )
3680         {
3681             indxX = termX + j;
3682             sum = 0;
3683             for ( k = 0; k < colsY; k++ )
3684             {
3685                 indxY = termY + k;
3686                 indxZ = k * colsZ + j;
3687                 if ( BIT_VECTOR_TST_BIT(Y,indxY) &&
3688                      BIT_VECTOR_TST_BIT(Z,indxZ) ) sum |= 1;
3689             }
3690             if (sum) BIT_VECTOR_SET_BIT(X,indxX)
3691             else     BIT_VECTOR_CLR_BIT(X,indxX)
3692         }
3693     }
3694   }
3695 }
3696
3697 void Matrix_Closure(wordptr addr, N_int rows, N_int cols)
3698 {
3699     N_word i;
3700     N_word j;
3701     N_word k;
3702     N_word ii;
3703     N_word ij;
3704     N_word ik;
3705     N_word kj;
3706     N_word termi;
3707     N_word termk;
3708
3709   if ((rows == cols) and (bits_(addr) == rows*cols))
3710   {
3711     for ( i = 0; i < rows; i++ )
3712     {
3713         ii = i * cols + i;
3714         BIT_VECTOR_SET_BIT(addr,ii)
3715     }
3716     for ( k = 0; k < rows; k++ )
3717     {
3718         termk = k * cols;
3719         for ( i = 0; i < rows; i++ )
3720         {
3721             termi = i * cols;
3722             ik = termi + k;
3723             for ( j = 0; j < rows; j++ )
3724             {
3725                 ij = termi + j;
3726                 kj = termk + j;
3727                 if ( BIT_VECTOR_TST_BIT(addr,ik) &&
3728                      BIT_VECTOR_TST_BIT(addr,kj) )
3729                      BIT_VECTOR_SET_BIT(addr,ij)
3730             }
3731         }
3732     }
3733   }
3734 }
3735
3736 void Matrix_Transpose(wordptr X, N_int rowsX, N_int colsX,
3737                       wordptr Y, N_int rowsY, N_int colsY)
3738 {
3739     N_word  i;
3740     N_word  j;
3741     N_word  ii;
3742     N_word  ij;
3743     N_word  ji;
3744     N_word  addii;
3745     N_word  addij;
3746     N_word  addji;
3747     N_word  bitii;
3748     N_word  bitij;
3749     N_word  bitji;
3750     N_word  termi;
3751     N_word  termj;
3752     boolean swap;
3753
3754   /* BEWARE that "in-place" is ONLY possible if the matrix is quadratic!! */
3755
3756   if ((rowsX == colsY) and (colsX == rowsY) and
3757       (bits_(X) == rowsX*colsX) and
3758       (bits_(Y) == rowsY*colsY))
3759   {
3760     if (rowsY == colsY) /* in-place is possible! */
3761     {
3762         for ( i = 0; i < rowsY; i++ )
3763         {
3764             termi = i * colsY;
3765             for ( j = 0; j < i; j++ )
3766             {
3767                 termj = j * colsX;
3768                 ij = termi + j;
3769                 ji = termj + i;
3770                 addij = ij >> LOGBITS;
3771                 addji = ji >> LOGBITS;
3772                 bitij = BITMASKTAB[ij AND MODMASK];
3773                 bitji = BITMASKTAB[ji AND MODMASK];
3774                 swap = ((*(Y+addij) AND bitij) != 0);
3775                 if ((*(Y+addji) AND bitji) != 0)
3776                      *(X+addij) |=     bitij;
3777                 else
3778                      *(X+addij) &= NOT bitij;
3779                 if (swap)
3780                      *(X+addji) |=     bitji;
3781                 else
3782                      *(X+addji) &= NOT bitji;
3783             }
3784             ii = termi + i;
3785             addii = ii >> LOGBITS;
3786             bitii = BITMASKTAB[ii AND MODMASK];
3787             if ((*(Y+addii) AND bitii) != 0)
3788                  *(X+addii) |=     bitii;
3789             else
3790                  *(X+addii) &= NOT bitii;
3791         }
3792     }
3793     else /* rowsX != colsX, in-place is NOT possible! */
3794     {
3795         for ( i = 0; i < rowsY; i++ )
3796         {
3797             termi = i * colsY;
3798             for ( j = 0; j < colsY; j++ )
3799             {
3800                 termj = j * colsX;
3801                 ij = termi + j;
3802                 ji = termj + i;
3803                 addij = ij >> LOGBITS;
3804                 addji = ji >> LOGBITS;
3805                 bitij = BITMASKTAB[ij AND MODMASK];
3806                 bitji = BITMASKTAB[ji AND MODMASK];
3807                 if ((*(Y+addij) AND bitij) != 0)
3808                      *(X+addji) |=     bitji;
3809                 else
3810                      *(X+addji) &= NOT bitji;
3811             }
3812         }
3813     }
3814   }
3815 }
3816
3817 /*****************************************************************************/
3818 /*  VERSION:  6.4                                                            */
3819 /*****************************************************************************/
3820 /*  VERSION HISTORY:                                                         */
3821 /*****************************************************************************/
3822 /*                                                                           */
3823 /*    Version 6.4  03.10.04  Added C++ comp. directives. Improved "Norm()".  */
3824 /*    Version 6.3  28.09.02  Added "Create_List()" and "GCD2()".             */
3825 /*    Version 6.2  15.09.02  Overhauled error handling. Fixed "GCD()".       */
3826 /*    Version 6.1  08.10.01  Make VMS linker happy: _lsb,_msb => _lsb_,_msb_ */
3827 /*    Version 6.0  08.10.00  Corrected overflow handling.                    */
3828 /*    Version 5.8  14.07.00  Added "Power()". Changed "Copy()".              */
3829 /*    Version 5.7  19.05.99  Quickened "Div_Pos()". Added "Product()".       */
3830 /*    Version 5.6  02.11.98  Leading zeros eliminated in "to_Hex()".         */
3831 /*    Version 5.5  21.09.98  Fixed bug of uninitialized "error" in Multiply. */
3832 /*    Version 5.4  07.09.98  Fixed bug of uninitialized "error" in Divide.   */
3833 /*    Version 5.3  12.05.98  Improved Norm. Completed history.               */
3834 /*    Version 5.2  31.03.98  Improved Norm.                                  */
3835 /*    Version 5.1  09.03.98  No changes.                                     */
3836 /*    Version 5.0  01.03.98  Major additions and rewrite.                    */
3837 /*    Version 4.2  16.07.97  Added is_empty, is_full.                        */
3838 /*    Version 4.1  30.06.97  Added word-ins/del, move-left/right, inc/dec.   */
3839 /*    Version 4.0  23.04.97  Rewrite. Added bit shift and bool. matrix ops.  */
3840 /*    Version 3.2  04.02.97  Added interval methods.                         */
3841 /*    Version 3.1  21.01.97  Fixed bug on 64 bit machines.                   */
3842 /*    Version 3.0  12.01.97  Added flip.                                     */
3843 /*    Version 2.0  14.12.96  Efficiency and consistency improvements.        */
3844 /*    Version 1.1  08.01.96  Added Resize and ExclusiveOr.                   */
3845 /*    Version 1.0  14.12.95  First version under UNIX (with Perl module).    */
3846 /*    Version 0.9  01.11.93  First version of C library under MS-DOS.        */
3847 /*    Version 0.1  ??.??.89  First version in Turbo Pascal under CP/M.       */
3848 /*                                                                           */
3849 /*****************************************************************************/
3850 /*  AUTHOR:                                                                  */
3851 /*****************************************************************************/
3852 /*                                                                           */
3853 /*    Steffen Beyer                                                          */
3854 /*    mailto:sb@engelschall.com                                              */
3855 /*    http://www.engelschall.com/u/sb/download/                              */
3856 /*                                                                           */
3857 /*****************************************************************************/
3858 /*  COPYRIGHT:                                                               */
3859 /*****************************************************************************/
3860 /*                                                                           */
3861 /*    Copyright (c) 1995 - 2004 by Steffen Beyer.                            */
3862 /*    All rights reserved.                                                   */
3863 /*                                                                           */
3864 /*****************************************************************************/
3865 /*  LICENSE:                                                                 */
3866 /*****************************************************************************/
3867 /*                                                                           */
3868 /*    This library is free software; you can redistribute it and/or          */
3869 /*    modify it under the terms of the GNU Library General Public            */
3870 /*    License as published by the Free Software Foundation; either           */
3871 /*    version 2 of the License, or (at your option) any later version.       */
3872 /*                                                                           */
3873 /*    This library is distributed in the hope that it will be useful,        */
3874 /*    but WITHOUT ANY WARRANTY; without even the implied warranty of         */
3875 /*    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU       */
3876 /*    Library General Public License for more details.                       */
3877 /*                                                                           */
3878 /*    You should have received a copy of the GNU Library General Public      */
3879 /*    License along with this library; if not, write to the                  */
3880 /*    Free Software Foundation, Inc.,                                        */
3881 /*    59 Temple Place, Suite 330, Boston, MA 02111-1307 USA                  */
3882 /*                                                                           */
3883 /*    or download a copy from ftp://ftp.gnu.org/pub/gnu/COPYING.LIB-2.0      */
3884 /*                                                                           */
3885 /*****************************************************************************/
3886 #endif