chiark - git - mdw - catacomb/blob - math/mp.h

   1 /* -*-c-*-
   2  *
   3  * Simple multiprecision arithmetic
   4  *
   5  * (c) 1999 Straylight/Edgeware
   6  */
   7
   8 /*----- Licensing notice --------------------------------------------------*
   9  *
  10  * This file is part of Catacomb.
  11  *
  12  * Catacomb is free software; you can redistribute it and/or modify
  13  * it under the terms of the GNU Library General Public License as
  14  * published by the Free Software Foundation; either version 2 of the
  15  * License, or (at your option) any later version.
  16  *
  17  * Catacomb is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  20  * GNU Library General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU Library General Public
  23  * License along with Catacomb; if not, write to the Free
  24  * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
  25  * MA 02111-1307, USA.
  26  */
  27
  28 #ifndef CATACOMB_MP_H
  29 #define CATACOMB_MP_H
  30
  31 #ifdef __cplusplus
  32   extern "C" {
  33 #endif
  34
  35 /*----- Header files ------------------------------------------------------*/
  36
  37 #include <assert.h>
  38 #include <string.h>
  39
  40 #include <mLib/sub.h>
  41
  42 #ifndef CATACOMB_MPW_H
  43 #  include "mpw.h"
  44 #endif
  45
  46 #ifndef CATACOMB_ARENA_H
  47 #  include "arena.h"
  48 #endif
  49
  50 #ifndef CATACOMB_MPARENA_H
  51 #  include "mparena.h"
  52 #endif
  53
  54 #ifndef CATACOMB_MPX_H
  55 #  include "mpx.h"
  56 #endif
  57
  58 /*----- Data structures ---------------------------------------------------*/
  59
  60 /* --- A multiprecision integer --- */
  61
  62 typedef struct mp {
  63   mpw *v, *vl;                          /* Vector of digits, current limit */
  64   size_t sz;                            /* Size of digit buffer in words */
  65   mparena *a;                           /* Arena for buffer allocation */
  66   unsigned f;                           /* Flags (see below) */
  67   unsigned ref;                         /* Reference counter */
  68 } mp;
  69
  70 #define MP_NEG 1u                       /* Negative (signed magnitude) */
  71 #define MP_BURN 2u                      /* Secret (viral flag) */
  72 #define MP_CONST 4u                     /* Uses strange memory allocation */
  73 #define MP_UNDEF 8u                     /* Contains nothing interesting */
  74 #define MP_DESTROYED 16u                /* Has been destroyed */
  75
  76 /* --- A factor for simultaneous exponentation --- *
  77  *
  78  * Used by the Montgomery and Barrett exponentiators.
  79  */
  80
  81 typedef struct mp_expfactor {
  82   mp *base;
  83   mp *exp;
  84 } mp_expfactor;
  85
  86 /*----- Useful constants --------------------------------------------------*/
  87
  88 extern mp mp_const[];
  89
  90 #define MP_ZERO  (&mp_const[0])
  91 #define MP_ONE   (&mp_const[1])
  92 #define MP_TWO   (&mp_const[2])
  93 #define MP_THREE (&mp_const[3])
  94 #define MP_FOUR  (&mp_const[4])
  95 #define MP_FIVE  (&mp_const[5])
  96 #define MP_TEN   (&mp_const[6])
  97 #define MP_256   (&mp_const[7])
  98 #define MP_MONE  (&mp_const[8])
  99
 100 #define MP_NEW ((mp *)0)
 101 #define MP_NEWSEC (&mp_const[9])
 102
 103 /*----- Trivial macros ----------------------------------------------------*/
 104
 105 /* --- @MP_LEN@ --- *
 106  *
 107  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 108  *
 109  * Returns:     Length of the integer, in words.
 110  */
 111
 112 #define MP_LEN(m) ((m)->vl - ((m)->v))
 113
 114 /*----- Memory management and reference counting --------------------------*/
 115
 116 /* --- @mp_new@ --- *
 117  *
 118  * Arguments:   @size_t sz@ = size of vector required
 119  *              @unsigned f@ = flags to set
 120  *
 121  * Returns:     Pointer to a new MP structure.
 122  *
 123  * Use:         Allocates a new multiprecision integer.  The data space is
 124  *              allocated from either the standard global or secret arena,
 125  *              depending on the initial flags requested.
 126  */
 127
 128 extern mp *mp_new(size_t /*sz*/, unsigned /*f*/);
 129
 130 /* --- @mp_create@ --- *
 131  *
 132  * Arguments:   @size_t sz@ = size of vector required
 133  *
 134  * Returns:     Pointer to pristine new MP structure with enough memory
 135  *              bolted onto it.
 136  *
 137  * Use:         Creates a new multiprecision integer with indeterminate
 138  *              contents.  The integer has a single reference.
 139  */
 140
 141 extern mp *mp_create(size_t /*sz*/);
 142
 143 /* --- @mp_createsecure@ --- *
 144  *
 145  * Arguments:   @size_t sz@ = size of vector required
 146  *
 147  * Returns:     Pointer to pristine new MP structure with enough memory
 148  *              bolted onto it.
 149  *
 150  * Use:         Creates a new multiprecision integer with indeterminate
 151  *              contents.  The integer has a single reference.  The integer's
 152  *              data space is allocated from the secure arena.  Its burn flag
 153  *              is set.
 154  */
 155
 156 extern mp *mp_createsecure(size_t /*sz*/);
 157
 158 /* --- @mp_build@ --- *
 159  *
 160  * Arguments:   @mp *m@ = pointer to an MP block to fill in
 161  *              @mpw *v@ = pointer to a word array
 162  *              @mpw *vl@ = pointer just past end of array
 163  *
 164  * Returns:     ---
 165  *
 166  * Use:         Creates a multiprecision integer representing some smallish
 167  *              number.  You must provide storage for the number and dispose
 168  *              of it when you've finished with it.  The number is marked as
 169  *              constant while it exists.
 170  */
 171
 172 extern void mp_build(mp */*m*/, mpw */*v*/, mpw */*vl*/);
 173
 174 /* --- @mp_destroy@ --- *
 175  *
 176  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 177  *
 178  * Returns:     ---
 179  *
 180  * Use:         Destroys a multiprecision integer. The reference count isn't
 181  *              checked.  Don't use this function if you don't know what
 182  *              you're doing: use @mp_drop@ instead.
 183  */
 184
 185 extern void mp_destroy(mp */*m*/);
 186
 187 /* --- @mp_copy@ --- *
 188  *
 189  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 190  *
 191  * Returns:     A copy of the given multiprecision integer.
 192  *
 193  * Use:         Copies the given integer.  In fact you just get another
 194  *              reference to the same old one again.
 195  */
 196
 197 extern mp *mp_copy(mp */*m*/);
 198
 199 #define MP_COPY(m) ((m)->ref++, (m))
 200
 201 /* --- @mp_drop@ --- *
 202  *
 203  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 204  *
 205  * Returns:     ---
 206  *
 207  * Use:         Drops a reference to an integer which isn't wanted any more.
 208  *              If there are no more references, the integer is destroyed.
 209  */
 210
 211 extern void mp_drop(mp */*m*/);
 212
 213 #define MP_DROP(m) do {                                                 \
 214   mp *_mm = (m);                                                        \
 215   _mm->ref--;                                                           \
 216   if (_mm->ref == 0 && !(_mm->f & MP_CONST))                            \
 217     mp_destroy(_mm);                                                    \
 218 } while (0)
 219
 220 /* --- @mp_split@ --- *
 221  *
 222  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 223  *
 224  * Returns:     A reference to the same integer, possibly with a different
 225  *              address.
 226  *
 227  * Use:         Splits off a modifiable version of the integer referred to.
 228  */
 229
 230 extern mp *mp_split(mp */*m*/);
 231
 232 #define MP_SPLIT(m) do {                                                \
 233   mp *_m = (m);                                                         \
 234   if ((_m->f & MP_CONST) || _m->ref > 1) {                              \
 235     size_t _len = MP_LEN(_m);                                           \
 236     mp *_mm = mp_new(_len, _m->f);                                      \
 237     if (!(_m->f & MP_UNDEF))                                            \
 238       memcpy(_mm->v, _m->v, MPWS(_len));                                \
 239     _m->ref--;                                                          \
 240     _m = _mm;                                                           \
 241   }                                                                     \
 242   (m) = _m;                                                             \
 243 } while (0)
 244
 245 /* --- @mp_resize@ --- *
 246  *
 247  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 248  *              @size_t sz@ = new size
 249  *
 250  * Returns:     ---
 251  *
 252  * Use:         Resizes the vector containing the integer's digits.  The new
 253  *              size must be at least as large as the current integer's
 254  *              length.  This isn't really intended for client use.
 255  */
 256
 257 extern void mp_resize(mp */*m*/, size_t /*sz*/);
 258
 259 #define MP_RESIZE(m, ssz) do {                                          \
 260   mp *_m = (m);                                                         \
 261   size_t _sz = (ssz);                                                   \
 262   mparena *_a = (_m->f & MP_BURN) ? MPARENA_SECURE : MPARENA_GLOBAL;    \
 263   mpw *_v;                                                              \
 264   size_t _len = MP_LEN(_m);                                             \
 265   assert(((void)"can't make size less than length", _sz >= _len));      \
 266   _v = mpalloc(_a, _sz);                                                \
 267   if (!(_m->f & MP_UNDEF))                                              \
 268     memcpy(_v, _m->v, MPWS(_len));                                      \
 269   if (_m->f & MP_BURN)                                                  \
 270     memset(_m->v, 0, MPWS(_m->sz));                                     \
 271   mpfree(_m->a, _m->v);                                                 \
 272   _m->a = _a;                                                           \
 273   _m->v = _v;                                                           \
 274   _m->vl = _v + _len;                                                   \
 275 } while (0)
 276
 277 /* --- @mp_ensure@ --- *
 278  *
 279  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 280  *              @size_t sz@ = required size
 281  *
 282  * Returns:     ---
 283  *
 284  * Use:         Ensures that the integer has enough space for @sz@ digits.
 285  *              The value is not changed.
 286  */
 287
 288 extern void mp_ensure(mp */*m*/, size_t /*sz*/);
 289
 290 #define MP_ENSURE(m, ssz) do {                                          \
 291   mp *_m = (m);                                                         \
 292   size_t _ssz = (ssz);                                                  \
 293   size_t _len = MP_LEN(_m);                                             \
 294   if (_ssz >= _len) {                                                   \
 295     if (_ssz > _m->sz)                                                  \
 296       mp_resize(_m, _ssz);                                              \
 297     if (!(_m->f & MP_UNDEF) && _ssz > _len)                             \
 298       memset(_m->vl, 0, MPWS(_ssz - _len));                             \
 299     _m->vl = _m->v + _ssz;                                              \
 300   }                                                                     \
 301 } while (0)
 302
 303 /* --- @mp_dest@ --- *
 304  *
 305  * Arguments:   @mp *m@ = a suggested destination integer
 306  *              @size_t sz@ = size required for result, in digits
 307  *              @unsigned f@ = various flags
 308  *
 309  * Returns:     A pointer to an appropriate destination.
 310  *
 311  * Use:         Converts a suggested destination into a real destination with
 312  *              the required properties.  If the real destination is @d@,
 313  *              then the following properties will hold:
 314  *
 315  *                * @d@ will have exactly one reference.
 316  *
 317  *                * If @m@ is not @MP_NEW@, then the contents of @m@ will not
 318  *                  change, unless @f@ has the @MP_UNDEF@ flag set.
 319  *
 320  *                * If @m@ is not @MP_NEW@, then he reference count of @m@ on
 321  *                  entry is equal to the sum of the counts of @d@ and @m@ on
 322  *                  exit.
 323  *
 324  *                * The size of @d@ will be at least @sz@.
 325  *
 326  *                * If @f@ has the @MP_BURN@ flag set, then @d@ will be
 327  *                  allocated from @MPARENA_SECURE@.
 328  *
 329  *              Understanding this function is crucial to using Catacomb's
 330  *              multiprecision integer library effectively.
 331  */
 332
 333 extern mp *mp_dest(mp */*m*/, size_t /*sz*/, unsigned /*f*/);
 334
 335 #define MP_DEST(m, ssz, f) do {                                         \
 336   mp *_m = (m);                                                         \
 337   size_t _ssz = (ssz);                                                  \
 338   unsigned _f = (f);                                                    \
 339   _m = mp_dest(_m, _ssz, _f);                                           \
 340   (m) = _m;                                                             \
 341 } while (0)
 342
 343 /*----- Size manipulation -------------------------------------------------*/
 344
 345 /* --- @mp_shrink@ --- *
 346  *
 347  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 348  *
 349  * Returns:     ---
 350  *
 351  * Use:         Reduces the recorded length of an integer.  This doesn't
 352  *              reduce the amount of memory used, although it can improve
 353  *              performance a bit.  To reduce memory, use @mp_minimize@
 354  *              instead.  This can't change the value of an integer, and is
 355  *              therefore safe to use even when there are multiple
 356  *              references.
 357  */
 358
 359 extern void mp_shrink(mp */*m*/);
 360
 361 #define MP_SHRINK(m) do {                                               \
 362   mp *_mm = (m);                                                        \
 363   MPX_SHRINK(_mm->v, _mm->vl);                                          \
 364   if (MP_ZEROP(_mm))                                                    \
 365     _mm->f &= ~MP_NEG;                                                  \
 366 } while (0)
 367
 368 /* --- @mp_minimize@ --- *
 369  *
 370  * Arguments:   @mp *m@ = pointer to a multiprecision integer
 371  *
 372  * Returns:     ---
 373  *
 374  * Use:         Reduces the amount of memory an integer uses.  It's best to
 375  *              do this to numbers which aren't going to change in the
 376  *              future.
 377  */
 378
 379 extern void mp_minimize(mp */*m*/);
 380
 381 /*----- Bit scanning ------------------------------------------------------*/
 382
 383 #ifndef CATACOMB_MPSCAN_H
 384 #  include "mpscan.h"
 385 #endif
 386
 387 /* --- @mp_scan@ --- *
 388  *
 389  * Arguments:   @mpscan *sc@ = pointer to bitscanner block
 390  *              @const mp *m@ = pointer to a multiprecision integer
 391  *
 392  * Returns:     ---
 393  *
 394  * Use:         Initializes a bitscanner on a multiprecision integer.
 395  */
 396
 397 extern void mp_scan(mpscan */*sc*/, const mp */*m*/);
 398
 399 #define MP_SCAN(sc, m) do {                                             \
 400   const mp *_mm = (m);                                                  \
 401   mpscan *_sc = (sc);                                                   \
 402   MPSCAN_INITX(_sc, _mm->v, _mm->vl);                                   \
 403 } while (0)
 404
 405 /* --- @mp_rscan@ --- *
 406  *
 407  * Arguments:   @mpscan *sc@ = pointer to bitscanner block
 408  *              @const mp *m@ = pointer to a multiprecision integer
 409  *
 410  * Returns:     ---
 411  *
 412  * Use:         Initializes a reverse bitscanner on a multiprecision
 413  *              integer.
 414  */
 415
 416 extern void mp_rscan(mpscan */*sc*/, const mp */*m*/);
 417
 418 #define MP_RSCAN(sc, m) do {                                            \
 419   const mp *_mm = (m);                                                  \
 420   mpscan *_sc = (sc);                                                   \
 421   MPSCAN_RINITX(_sc, _mm->v, _mm->vl);                                  \
 422 } while (0)
 423
 424 /* --- Other bitscanning aliases --- */
 425
 426 #define mp_step mpscan_step
 427 #define mp_bit mpscan_bit
 428 #define mp_rstep mpscan_rstep
 429 #define mp_rbit mpscan_rbit
 430
 431 #define MP_STEP MPSCAN_STEP
 432 #define MP_BIT MPSCAN_BIT
 433 #define MP_RSTEP MPSCAN_RSTEP
 434 #define MP_RBIT MPSCAN_RBIT
 435
 436 /*----- Loading and storing -----------------------------------------------*/
 437
 438 /* --- @mp_octets@ --- *
 439  *
 440  * Arguments:   @const mp *m@ = a multiprecision integer
 441  *
 442  * Returns:     The number of octets required to represent @m@.
 443  *
 444  * Use:         Calculates the external storage required for a multiprecision
 445  *              integer.
 446  */
 447
 448 extern size_t mp_octets(const mp */*m*/);
 449
 450 /* --- @mp_octets2c@ --- *
 451  *
 452  * Arguments:   @const mp *m@ = a multiprecision integer
 453  *
 454  * Returns:     The number of octets required to represent @m@.
 455  *
 456  * Use:         Calculates the external storage required for a multiprecision
 457  *              integer represented as two's complement.
 458  */
 459
 460 extern size_t mp_octets2c(const mp */*m*/);
 461
 462 /* --- @mp_bits@ --- *
 463  *
 464  * Arguments:   @const mp *m@ = a multiprecision integer
 465  *
 466  * Returns:     The number of bits required to represent @m@.
 467  *
 468  * Use:         Calculates the external storage required for a multiprecision
 469  *              integer.
 470  */
 471
 472 extern unsigned long mp_bits(const mp */*m*/);
 473
 474 /* --- @mp_loadl@ --- *
 475  *
 476  * Arguments:   @mp *d@ = destination
 477  *              @const void *pv@ = pointer to source data
 478  *              @size_t sz@ = size of the source data
 479  *
 480  * Returns:     Resulting multiprecision number.
 481  *
 482  * Use:         Loads a multiprecision number from an array of octets.  The
 483  *              first byte in the array is the least significant.  More
 484  *              formally, if the bytes are %$b_0, b_1, \ldots, b_{n-1}$%
 485  *              then the result is %$N = \sum_{0 \le i < n} b_i 2^{8i}$%.
 486  */
 487
 488 extern mp *mp_loadl(mp */*d*/, const void */*pv*/, size_t /*sz*/);
 489
 490 /* --- @mp_storel@ --- *
 491  *
 492  * Arguments:   @const mp *m@ = source
 493  *              @void *pv@ = pointer to output array
 494  *              @size_t sz@ = size of the output array
 495  *
 496  * Returns:     ---
 497  *
 498  * Use:         Stores a multiprecision number in an array of octets.  The
 499  *              first byte in the array is the least significant.  If the
 500  *              array is too small to represent the number, high-order bits
 501  *              are truncated; if the array is too large, high order bytes
 502  *              are filled with zeros.  More formally, if the number is
 503  *              %$N = \sum{0 \le i} b_i 2^{8i}$% where %$0 \le b_i < 256$%,
 504  *              then the array is %$b_0, b_1, \ldots, b_{n-1}$%.
 505  */
 506
 507 extern void mp_storel(const mp */*m*/, void */*pv*/, size_t /*sz*/);
 508
 509 /* --- @mp_loadb@ --- *
 510  *
 511  * Arguments:   @mp *d@ = destination
 512  *              @const void *pv@ = pointer to source data
 513  *              @size_t sz@ = size of the source data
 514  *
 515  * Returns:     Resulting multiprecision number.
 516  *
 517  * Use:         Loads a multiprecision number from an array of octets.  The
 518  *              last byte in the array is the least significant.  More
 519  *              formally, if the bytes are %$b_{n-1}, b_{n-2}, \ldots, b_0$%
 520  *              then the result is %$N = \sum_{0 \le i < n} b_i 2^{8i}$%.
 521  */
 522
 523 extern mp *mp_loadb(mp */*d*/, const void */*pv*/, size_t /*sz*/);
 524
 525 /* --- @mp_storeb@ --- *
 526  *
 527  * Arguments:   @const mp *m@ = source
 528  *              @void *pv@ = pointer to output array
 529  *              @size_t sz@ = size of the output array
 530  *
 531  * Returns:     ---
 532  *
 533  * Use:         Stores a multiprecision number in an array of octets.  The
 534  *              last byte in the array is the least significant.  If the
 535  *              array is too small to represent the number, high-order bits
 536  *              are truncated; if the array is too large, high order bytes
 537  *              are filled with zeros.  More formally, if the number is
 538  *              %$N = \sum{0 \le i} b_i 2^{8i}$% where %$0 \le b_i < 256$%,
 539  *              then the array is %$b_{n-1}, b_{n-2}, \ldots, b_0$%.
 540  */
 541
 542 extern void mp_storeb(const mp */*m*/, void */*pv*/, size_t /*sz*/);
 543
 544 /* --- @mp_loadl2c@ --- *
 545  *
 546  * Arguments:   @mp *d@ = destination
 547  *              @const void *pv@ = pointer to source data
 548  *              @size_t sz@ = size of the source data
 549  *
 550  * Returns:     Resulting multiprecision number.
 551  *
 552  * Use:         Loads a multiprecision number from an array of octets as
 553  *              two's complement.  The first byte in the array is the least
 554  *              significant.
 555  */
 556
 557 extern mp *mp_loadl2c(mp */*d*/, const void */*pv*/, size_t /*sz*/);
 558
 559 /* --- @mp_storel2c@ --- *
 560  *
 561  * Arguments:   @const mp *m@ = source
 562  *              @void *pv@ = pointer to output array
 563  *              @size_t sz@ = size of the output array
 564  *
 565  * Returns:     ---
 566  *
 567  * Use:         Stores a multiprecision number in an array of octets as two's
 568  *              complement.  The first byte in the array is the least
 569  *              significant.  If the array is too small to represent the
 570  *              number, high-order bits are truncated; if the array is too
 571  *              large, high order bytes are sign-extended.
 572  */
 573
 574 extern void mp_storel2c(const mp */*m*/, void */*pv*/, size_t /*sz*/);
 575
 576 /* --- @mp_loadb2c@ --- *
 577  *
 578  * Arguments:   @mp *d@ = destination
 579  *              @const void *pv@ = pointer to source data
 580  *              @size_t sz@ = size of the source data
 581  *
 582  * Returns:     Resulting multiprecision number.
 583  *
 584  * Use:         Loads a multiprecision number from an array of octets as
 585  *              two's complement.  The last byte in the array is the least
 586  *              significant.
 587  */
 588
 589 extern mp *mp_loadb2c(mp */*d*/, const void */*pv*/, size_t /*sz*/);
 590
 591 /* --- @mp_storeb2c@ --- *
 592  *
 593  * Arguments:   @const mp *m@ = source
 594  *              @void *pv@ = pointer to output array
 595  *              @size_t sz@ = size of the output array
 596  *
 597  * Returns:     ---
 598  *
 599  * Use:         Stores a multiprecision number in an array of octets, as
 600  *              two's complement.  The last byte in the array is the least
 601  *              significant.  If the array is too small to represent the
 602  *              number, high-order bits are truncated; if the array is too
 603  *              large, high order bytes are sign-extended.
 604  */
 605
 606 extern void mp_storeb2c(const mp */*m*/, void */*pv*/, size_t /*sz*/);
 607
 608 /*----- Bit operations ----------------------------------------------------*/
 609
 610 /* --- @mp_not@ --- *
 611  *
 612  * Arguments:   @mp *d@ = destination
 613  *              @mp *a@ = source
 614  *
 615  * Returns:     The bitwise complement of the source.
 616  */
 617
 618 extern mp *mp_not(mp */*d*/, mp */*a*/);
 619
 620 /* --- @mp_bitop@ --- *
 621  *
 622  * Arguments:   @mp *d@ = destination
 623  *              @mp *a, *b@ = sources
 624  *
 625  * Returns:     The result of the given bitwise operation.  These functions
 626  *              don't handle negative numbers at all sensibly.  For that, use
 627  *              the @...2c@ variants.  The functions are named after the
 628  *              truth tables they generate:
 629  *
 630  *                      a:      0011
 631  *                      b:      0101
 632  *                      @mpx_bitXXXX@
 633  */
 634
 635 #define MP_BITDECL(string)                                              \
 636   extern mp *mp_bit##string(mp */*d*/, mp */*a*/, mp */*b*/);
 637 MPX_DOBIN(MP_BITDECL)
 638
 639 /* --- @mp_[n]and@, @mp_[n]or@, @mp_[n]xor@, @mp_not@ --- *
 640  *
 641  * Synonyms for the commonly-used functions.
 642  */
 643
 644 #define mp_and  mp_bit0001
 645 #define mp_or   mp_bit0111
 646 #define mp_nand mp_bit1110
 647 #define mp_nor  mp_bit1000
 648 #define mp_xor  mp_bit0110
 649
 650 /* --- @mp_testbit@ --- *
 651  *
 652  * Arguments:   @mp *x@ = a large integer
 653  *              @unsigned long n@ = which bit to test
 654  *
 655  * Returns:     Nonzero if the bit is set, zero if not.
 656  */
 657
 658 extern int mp_testbit(mp */*x*/, unsigned long /*n*/);
 659
 660 /* --- @mp_setbit@, @mp_clearbit@ --- *
 661  *
 662  * Arguments:   @mp *d@ = a destination
 663  *              @mp *x@ = a large integer
 664  *              @unsigned long n@ = which bit to modify
 665  *
 666  * Returns:     The argument @x@, with the appropriate bit set or cleared.
 667  */
 668
 669 extern mp *mp_setbit(mp */*d*/, mp */*x*/, unsigned long /*n*/);
 670 extern mp *mp_clearbit(mp */*d*/, mp */*x*/, unsigned long /*n*/);
 671
 672 /* --- @mp_lsl@, @mp_lslc@, @mp_lsr@ --- *
 673  *
 674  * Arguments:   @mp *d@ = destination
 675  *              @mp *a@ = source
 676  *              @size_t n@ = number of bits to move
 677  *
 678  * Returns:     Result, @a@ shifted left or right by @n@.
 679  *
 680  * Use:         Bitwise shift operators.  @mp_lslc@ fills the bits introduced
 681  *              on the right with ones instead of zeroes: it's used
 682  *              internally by @mp_lsl2c@, though it may be useful on its
 683  *              own.
 684  */
 685
 686 extern mp *mp_lsl(mp */*d*/, mp */*a*/, size_t /*n*/);
 687 extern mp *mp_lslc(mp */*d*/, mp */*a*/, size_t /*n*/);
 688 extern mp *mp_lsr(mp */*d*/, mp */*a*/, size_t /*n*/);
 689
 690 /* --- @mp_not2c@ --- *
 691  *
 692  * Arguments:   @mp *d@ = destination
 693  *              @mp *a@ = source
 694  *
 695  * Returns:     The sign-extended complement of the argument.
 696  */
 697
 698 extern mp *mp_not2c(mp */*d*/, mp */*a*/);
 699
 700 /* --- @mp_bitop2c@ --- *
 701  *
 702  * Arguments:   @mp *d@ = destination
 703  *              @mp *a, *b@ = sources
 704  *
 705  * Returns:     The result of the given bitwise operation.  Negative numbers
 706  *              are treated as two's complement, sign-extended infinitely to
 707  *              the left.  The functions are named after the truth tables
 708  *              they generate:
 709  *
 710  *                      a:      0011
 711  *                      b:      0101
 712  *                      @mpx_bitXXXX@
 713  */
 714
 715 #define MP_BIT2CDECL(string)                                            \
 716   extern mp *mp_bit##string##2c(mp */*d*/, mp */*a*/, mp */*b*/);
 717 MPX_DOBIN(MP_BIT2CDECL)
 718
 719 /* --- @mp_[n]and@, @mp_[n]or@, @mp_[n]xor@, @mp_not@ --- *
 720  *
 721  * Synonyms for the commonly-used functions.
 722  */
 723
 724 #define mp_and2c  mp_bit00012c
 725 #define mp_or2c   mp_bit01112c
 726 #define mp_nand2c mp_bit11102c
 727 #define mp_nor2c  mp_bit10002c
 728 #define mp_xor2c  mp_bit01102c
 729
 730 /* --- @mp_lsl2c@, @mp_lsr2c@ --- *
 731  *
 732  * Arguments:   @mp *d@ = destination
 733  *              @mp *a@ = source
 734  *              @size_t n@ = number of bits to move
 735  *
 736  * Returns:     Result, @a@ shifted left or right by @n@.  Handles the
 737  *              pretence of sign-extension for negative numbers.
 738  */
 739
 740 extern mp *mp_lsl2c(mp */*d*/, mp */*a*/, size_t /*n*/);
 741 extern mp *mp_lsr2c(mp */*d*/, mp */*a*/, size_t /*n*/);
 742
 743 /* --- @mp_testbit2c@ --- *
 744  *
 745  * Arguments:   @mp *x@ = a large integer
 746  *              @unsigned long n@ = which bit to test
 747  *
 748  * Returns:     Nonzero if the bit is set, zero if not.  Fakes up two's
 749  *              complement representation.
 750  */
 751
 752 extern int mp_testbit2c(mp */*x*/, unsigned long /*n*/);
 753
 754 /* --- @mp_setbit2c@, @mp_clearbit2c@ --- *
 755  *
 756  * Arguments:   @mp *d@ = a destination
 757  *              @mp *x@ = a large integer
 758  *              @unsigned long n@ = which bit to modify
 759  *
 760  * Returns:     The argument @x@, with the appropriate bit set or cleared.
 761  *              Fakes up two's complement representation.
 762  */
 763
 764 extern mp *mp_setbit2c(mp */*d*/, mp */*x*/, unsigned long /*n*/);
 765 extern mp *mp_clearbit2c(mp */*d*/, mp */*x*/, unsigned long /*n*/);
 766
 767 /*----- Comparisons -------------------------------------------------------*/
 768
 769 /* --- @mp_eq@ --- *
 770  *
 771  * Arguments:   @const mp *a, *b@ = two numbers
 772  *
 773  * Returns:     Nonzero if the numbers are equal.
 774  */
 775
 776 extern int mp_eq(const mp */*a*/, const mp */*b*/);
 777
 778 #define MP_EQ(a, b)                                                     \
 779   ((((a)->f ^ (b)->f) & MP_NEG) == 0 &&                                 \
 780    mpx_ueq((a)->v, (a)->vl, (b)->v, (b)->vl))
 781
 782 /* --- @mp_cmp@ --- *
 783  *
 784  * Arguments:   @const mp *a, *b@ = two numbers
 785  *
 786  * Returns:     Less than, equal to or greater than zero, according to
 787  *              whether @a@ is less than, equal to or greater than @b@.
 788  */
 789
 790 extern int mp_cmp(const mp */*a*/, const mp */*b*/);
 791
 792 #define MP_CMP(a, op, b) (mp_cmp((a), (b)) op 0)
 793
 794 /* --- Other handy macros --- */
 795
 796 #define MP_NEGP(x) ((x)->f & MP_NEG)
 797 #define MP_ZEROP(x) (!MP_LEN(x))
 798 #define MP_POSP(x) (!MP_NEGP(x) && !MP_ZEROP(x))
 799 #define MP_ODDP(x) (!MP_ZEROP(x) && ((x)->v[0] & 1u))
 800 #define MP_EVENP(x) (!MP_ODDP(x))
 801
 802 /*----- Arithmetic operations ---------------------------------------------*/
 803
 804 /* --- @mp_neg@ --- *
 805  *
 806  * Arguments:   @mp *d@ = destination
 807  *              @mp *a@ = argument
 808  *
 809  * Returns:     The negation of the argument.
 810  *
 811  * Use:         Negates its argument.
 812  */
 813
 814 extern mp *mp_neg(mp */*d*/, mp */*a*/);
 815
 816 /* --- @mp_add@ --- *
 817  *
 818  * Arguments:   @mp *d@ = destination
 819  *              @mp *a, *b@ = sources
 820  *
 821  * Returns:     Result, @a@ added to @b@.
 822  */
 823
 824 extern mp *mp_add(mp */*d*/, mp */*a*/, mp */*b*/);
 825
 826 /* --- @mp_sub@ --- *
 827  *
 828  * Arguments:   @mp *d@ = destination
 829  *              @mp *a, *b@ = sources
 830  *
 831  * Returns:     Result, @b@ subtracted from @a@.
 832  */
 833
 834 extern mp *mp_sub(mp */*d*/, mp */*a*/, mp */*b*/);
 835
 836 /* --- @mp_mul@ --- *
 837  *
 838  * Arguments:   @mp *d@ = destination
 839  *              @mp *a, *b@ = sources
 840  *
 841  * Returns:     Result, @a@ multiplied by @b@.
 842  */
 843
 844 extern mp *mp_mul(mp */*d*/, mp */*a*/, mp */*b*/);
 845
 846 /* --- @mp_sqr@ --- *
 847  *
 848  * Arguments:   @mp *d@ = destination
 849  *              @mp *a@ = source
 850  *
 851  * Returns:     Result, @a@ squared.
 852  */
 853
 854 extern mp *mp_sqr(mp */*d*/, mp */*a*/);
 855
 856 /* --- @mp_div@ --- *
 857  *
 858  * Arguments:   @mp **qq, **rr@ = destination, quotient and remainder
 859  *              @mp *a, *b@ = sources
 860  *
 861  * Use:         Calculates the quotient and remainder when @a@ is divided by
 862  *              @b@.
 863  */
 864
 865 extern void mp_div(mp **/*qq*/, mp **/*rr*/, mp */*a*/, mp */*b*/);
 866
 867 /* --- @mp_exp@ --- *
 868  *
 869  * Arguments:   @mp *d@ = fake destination
 870  *              @mp *a@ = base
 871  *              @mp *e@ = exponent
 872  *
 873  * Returns:     Result, %$a^e$%.
 874  */
 875
 876 extern mp *mp_exp(mp */*d*/, mp */*a*/, mp */*e*/);
 877
 878 /* --- @mp_odd@ --- *
 879  *
 880  * Arguments:   @mp *d@ = pointer to destination integer
 881  *              @mp *m@ = pointer to source integer
 882  *              @size_t *s@ = where to store the power of 2
 883  *
 884  * Returns:     An odd integer integer %$t$% such that %$m = 2^s t$%.
 885  *
 886  * Use:         Computes a power of two and an odd integer which, when
 887  *              multiplied, give a specified result.  This sort of thing is
 888  *              useful in number theory quite often.
 889  */
 890
 891 extern mp *mp_odd(mp */*d*/, mp */*m*/, size_t */*s*/);
 892
 893 /*----- More advanced algorithms ------------------------------------------*/
 894
 895 /* --- @mp_sqrt@ --- *
 896  *
 897  * Arguments:   @mp *d@ = pointer to destination integer
 898  *              @mp *a@ = (nonnegative) integer to take square root of
 899  *
 900  * Returns:     The largest integer %$x$% such that %$x^2 \le a$%.
 901  *
 902  * Use:         Computes integer square roots.
 903  *
 904  *              The current implementation isn't very good: it uses the
 905  *              Newton-Raphson method to find an approximation to %$a$%.  If
 906  *              there's any demand for a better version, I'll write one.
 907  */
 908
 909 extern mp *mp_sqrt(mp */*d*/, mp */*a*/);
 910
 911 /* --- @mp_gcd@ --- *
 912  *
 913  * Arguments:   @mp **gcd, **xx, **yy@ = where to write the results
 914  *              @mp *a, *b@ = sources (must be nonzero)
 915  *
 916  * Returns:     ---
 917  *
 918  * Use:         Calculates @gcd(a, b)@, and two numbers @x@ and @y@ such that
 919  *              @ax + by = gcd(a, b)@.  This is useful for computing modular
 920  *              inverses.  Neither @a@ nor @b@ may be zero.
 921  */
 922
 923 extern void mp_gcd(mp **/*gcd*/, mp **/*xx*/, mp **/*yy*/,
 924                    mp */*a*/, mp */*b*/);
 925
 926 /* -- @mp_modinv@ --- *
 927  *
 928  * Arguments:   @mp *d@ = destination
 929  *              @mp *x@ = argument
 930  *              @mp *p@ = modulus
 931  *
 932  * Returns:     The inverse %$x^{-1} \bmod p$%.
 933  *
 934  * Use:         Computes a modular inverse.    An assertion fails if %$p$%
 935  *              has no inverse.
 936  */
 937
 938 extern mp *mp_modinv(mp */*d*/, mp */*x*/, mp */*p*/);
 939
 940 /* --- @mp_jacobi@ --- *
 941  *
 942  * Arguments:   @mp *a@ = an integer
 943  *              @mp *n@ = another integer
 944  *
 945  * Returns:     @-1@, @0@ or @1@ -- the Jacobi symbol %$J(a, n)$%.
 946  *
 947  * Use:         Computes the Kronecker symbol %$\jacobi{a}{n}$%.  If @n@ is
 948  *              prime, this is the Legendre symbol and is equal to 1 if and
 949  *              only if @a@ is a quadratic residue mod @n@.  The result is
 950  *              zero if and only if @a@ and @n@ have a common factor greater
 951  *              than one.
 952  *
 953  *              If @n@ is composite, then this computes the Kronecker symbol
 954  *
 955  *                %$\jacobi{a}{n}=\jacobi{a}{u}\prod_i\jacobi{a}{p_i}^{e_i}$%
 956  *
 957  *              where %$n = u p_0^{e_0} \ldots p_{n-1}^{e_{n-1}}$% is the
 958  *              prime factorization of %$n$%.  The missing bits are:
 959  *
 960  *                * %$\jacobi{a}{1} = 1$%;
 961  *                * %$\jacobi{a}{-1} = 1$% if @a@ is negative, or 1 if
 962  *                  positive;
 963  *                * %$\jacobi{a}{0} = 0$%;
 964  *                * %$\jacobi{a}{2}$ is 0 if @a@ is even, 1 if @a@ is
 965  *                  congruent to 1 or 7 (mod 8), or %$-1$% otherwise.
 966  *
 967  *              If %$n$% is positive and odd, then this is the Jacobi
 968  *              symbol.  (The Kronecker symbol is a consistant domain
 969  *              extension; the Jacobi symbol was implemented first, and the
 970  *              name stuck.)
 971  */
 972
 973 extern int mp_jacobi(mp */*a*/, mp */*n*/);
 974
 975 /* --- @mp_modsqrt@ --- *
 976  *
 977  * Arguments:   @mp *d@ = destination integer
 978  *              @mp *a@ = source integer
 979  *              @mp *p@ = modulus (must be prime)
 980  *
 981  * Returns:     If %$a$% is a quadratic residue, a square root of %$a$%; else
 982  *              a null pointer.
 983  *
 984  * Use:         Returns an integer %$x$% such that %$x^2 \equiv a \pmod{p}$%,
 985  *              if one exists; else a null pointer.  This function will not
 986  *              work if %$p$% is composite: you must factor the modulus, take
 987  *              a square root mod each factor, and recombine the results
 988  *              using the Chinese Remainder Theorem.
 989  *
 990  *              We guarantee that the square root returned is the smallest
 991  *              one (i.e., the `positive' square root).
 992  */
 993
 994 extern mp *mp_modsqrt(mp */*d*/, mp */*a*/, mp */*p*/);
 995
 996 /* --- @mp_modexp@ --- *
 997  *
 998  * Arguments:   @mp *d@ = fake destination
 999  *              @mp *x@ = base of exponentiation
1000  *              @mp *e@ = exponent
1001  *              @mp *n@ = modulus (must be positive)
1002  *
1003  * Returns:     The value %$x^e \bmod n$%.
1004  */
1005
1006 extern mp *mp_modexp(mp */*d*/, mp */*x*/, mp */*e*/, mp */*n*/);
1007
1008 /*----- Test harness support ----------------------------------------------*/
1009
1010 #include <mLib/testrig.h>
1011
1012 #ifndef CATACOMB_MPTEXT_H
1013 #  include "mptext.h"
1014 #endif
1015
1016 extern const test_type type_mp;
1017
1018 /*----- That's all, folks -------------------------------------------------*/
1019
1020 #ifdef __cplusplus
1021   }
1022 #endif
1023
1024 #endif