© 2008 W.Ehrhardt Last update Oct. 03, 2008

Basic information about the MPArith package (V1.7.23)

Contents



Introduction
The MPArith (formerly mpint) package contains Pascal source for multi precision integer, rational, and floating point arithmetic; the units can be compiled with the usual Pascal versions that allow const parameters. MPArith does not have the most efficient code, but it is small, understandable, open source, and can be compiled with BP7. My Pascal routines use many public resources (source code libraries, books, articles), links are given in the references section.

One of the main deviations to the C versions is the usage of procedures and mp_error and/or exceptions instead of functions that return the error codes. This is analogous to IO procedures/IOResult, error handling is described below.


Unit/source file information
The MPArith package includes the following essential files, additionally there are other test programs, batch files, data files, etc...
mp_base.pas   MP integer arithmetic basic routines
mp_calc.pas   Parse and evaluate mp_int expressions
mp_conf.inc   Include file: config definitions/options MP library
mp_numth.pas  MP integer modular arithmetic and number theoretic functions
mp_prng.pas   MP interface to (C)PRNG, functions for PRNG generation
mp_ratio.pas  MP rational arithmetic routines
mp_rcalc.pas  Parse and evaluate mp_float expressions
mp_rconp.pas  124000 bit tables for pi, ln(2), ln(10)
mp_real.pas   MP floating point arithmetic routines
mp_rsa.pas    MP functions for basic RSA based public-key cryptography
mp_supp.pas   Supplemental routines for MPArith
mp_types.pas  MPArith type definitions and constants
std.inc       Standard include file with definitions and options
t_calc.pas    Simple interactive calculator for mp_int expressions
t_check.pas   Main test program for the integer  routines
t_chkfp.pas   Main test program for the floating point routines
t_chkrat.pas  Main test program for the rational routines
t_rcalc.pas   Simple interactive calculator for mp_float expressions
Hint for Delphi IDE users: If a t_*.pas test program cannot be compiled in the IDE, rename or copy the .pas file to a .dpr file and try again.


MPArith definitions
The basic mp_int type is a record with a pointer to a dynamically allocated array of mp_digits (plus some extras), a rational number is represented by a record of two mp_ints, and a real number is defined by an mp_int mantissa, an exponent, and a bit precision:

type
  TDigitArray = packed array[0..MaxDigits+63] of mp_digit; {the digits of an mp_int}
  PDigitArray = ^TDigitArray;                              {pointer to digit array}

type
  mp_int    = record                   {MP integer number           }
                pdigits : PDigitArray; {pointer to digit array      }
                alloc   : word;        {allocated digits in pdigits }
                used    : word;        {used digits in pdigits      }
                sign    : word;        {sign: MP_ZPOS or MP_NEG     }
                magic   : word;        {set to MP_MAGIC by mp_init  }
              end;

  pmp_int   = ^mp_int;                 {pointer to an MP integer    }
  pmp_digit = ^mp_digit;               {pointer to an MP digit      }

  mp_rat    = record                   {MP rational number          }
                num: mp_int;           {numerator                   }
                den: mp_int;           {denominator > 0             }
              end;
  pmp_rat   = ^mp_rat;                 {pointer to an MP rational   }

  mp_float  = record                   {MP floating point number    }
                mantissa: mp_int;      {mantissa of an mp_float     }
                exponent: longint;     {exponent of an mp_float     }
                bitprec : longint;     {bit precision=max bitsize of mantissa}
              end;
  pmp_float = ^mp_float;               {pointer to an MP float}
The mp_digit type can be configured with the mp_conf.inc file. Default: MP_32BIT if compiler supports int64, MP_16BIT otherwise.

{$ifdef MP_32BIT}
type
  mp_digit  = cardinal;           {type that holds an MP digit   }
  mp_word   = int64;              {type that holds two MP digits }
const
  DIGIT_BIT = 30;                 {number of bits of an MP digit }
  MAXDigits = 16000;              {max number of mp_int digits   }
{$else}
type
  mp_digit  = word;               {type that holds an MP digit   }
  mp_word   = longint;            {type that holds two MP digits }
const
  DIGIT_BIT = 15;                 {number of bits of an MP digit }
  MAXDigits = 32000;              {max number of mp_int digits   }
{$endif}
mp_digits are defined in mp_types by DIGIT_BIT, ie the number of bits used in the base type. An mp_digit must be able to hold DIGIT_BIT + 1 bits, an mp_word must be able to hold 2*DIGIT_BIT + 1 bits. DIGIT_BIT must be at least 8. MAXDigits is 16000 for MP_32BIT and 32000 otherwise, resulting in about 144444 decimal digits (with the 15/16 bit config) as shown by the simple example (Note: 479999 = 32000*15 - 1).


Other configuration parameters
In mp_types.pas there a few other parameters, that can be used to control certain features of the library (the default values are assigned in the initialization part of the unit):

var
  mp_mul_cutoff : word;    {default    32: Karatsuba multiplication cutoff    }
  mp_sqr_cutoff : word;    {default    32: Karatsuba square cutoff            }
  mp_clearzero  : boolean; {default false: clear memory to 0 before freemem.  }
                           {Note: 32 bit realloc may NOT clear the old memory!}

var
  mp_show_plus  : boolean; {default false: if true include "+" in ascii string}
  mp_uppercase  : boolean; {default false: uppercase chars for radix > 10     }
  mp_roundfloat : boolean; {default false: round float representation of mp_rat}

Two mp_ints are multiplied with the Karatsuba algorithm, if both factors have at least mp_mul_cutoff  mp_digits, for squaring the threshold is mp_sqr_cutoff. The Karatsuba algorithms are asymptotically faster than the standard ones, but because of the expensive structure the better performance is observed only above certain limits. The cutoff values depend on the processor, frequency, generated code etc. There are two test programs for tuning the Karatsuba cutoffs.

Error handling, conditional defines
In the file mp_conf.inc there are a few conditional defines related to error handling, the defaults are are listed below. During development it is recommended to use the debug settings or at least define MPC_HaltOnError.

If MPC_HaltOnError (and MPC_HaltOnArgCheck) are not defined, errors are returned via mp_error. This normally is an integer (thread) variable, but can be defined as a function (MPC_ErrFunc) for debugging. MPC_NOHALT overrides the MPC_HaltOnError and HaltOnArgCheck definitions.

MPC_UseExceptions is defined for Delphi 2+, VirtualPascal, and FPC 2.0.2+ but must be activated by setting other options (e.g. MPC_HaltOnError). If active, errors are signaled via exceptions otherwise via run time errors. (Bugs in older FPC versions do not allow to configure exceptions dynamically via include file or command line options.)

With a few exceptions the library routines do nothing if mp_error<>MP_OKAY (BUT mp_clear DOES deallocation, ...). If no RTEs or exceptions are enabled, mp_error should be checked and reset via set_mp_error(0). The function mp_result returns and resets mp_error (like IOResult for IO errors).

If MPC_UseISAAC is defined, mp_prng uses the ISAAC cryptographic random number generator, otherwise the Taus88 generator is used. If MPC_Randomize is defined, the random number generator is initialized via randomize. If MPC_UseTSC is defined, the bits from the TSC are used together with randomize to seed the PRNG.

{.$define MPC_ArgCheck}         {Check arguments: initialized, range etc}
                                {RTE/exception if check failed          }
{$define MPC_Assert}            {Use assert for critical parts          }
{.$define MPC_Diagnostic}       {include diagnostic code/functions      }
{$define MPC_ECM_Primetable}    {Prime table in mp_ecm_factor (BIT32)   }
{.$define MPC_ErrFunc}          {mp_error is function not integer var   }
{$define MPC_HaltOnError}       {RTE/exceptions on errors vs. mp_error  }
{$define MPC_HaltOnArgCheck}    {RTE/exceptions on arg check failures   }
{.$define MPC_NOHALT}           {override MPC_HaltOnError/HaltOnArgCheck}
{.$define MPC_NOISAAC}          {override MPC_UseISAAC, no ISAAC        }
{.$define MPC_Randomize}        {randomize the PRN generator at startup }
                                {undef if predictable 'random numbers'  }
                                {are needed for debugging or testing    }
{.$define MPC_UseTSC}           {use TSC and randomize for PRNG seeds   }
{$define MPC_Ln2Ln10Tab}        {use tables for ln(2), ln(10) constants }
{$define MPC_Reduce_2k}         {use mp_reduce_2k in mp_exptmod         }
{$define MPC_UseKONG}           {Use Kong if sqrtmethod=0 and p=9 mod 16}
{$define MPC_UseGCD32}          {use 32 bit GCD in MP GCD routines      }
{$define MPC_UseISAAC}          {use ISAAC random number generator      }


{$ifdef debug}
  {define config options for debug mode}
  {$ifndef MPC_ArgCheck}        {$define MPC_ArgCheck}       {$endif}
  {$ifndef MPC_Assert}          {$define MPC_Assert}         {$endif}
  {$ifndef MPC_ErrFunc}         {$define MPC_ErrFunc}        {$endif}
  {$ifndef MPC_HaltOnError}     {$define MPC_HaltOnError}    {$endif}
  {$ifndef MPC_HaltOnArgCheck}  {$define MPC_HaltOnArgCheck} {$endif}
  {$ifndef MPC_Diagnostic}      {$define MPC_Diagnostic}     {$endif}
{$endif}


Simple complete example, general information
The MPArith units can be compiled with the usual Pascal versions that allow const parameters (BP 7.0, VP 2.1, FPC 1.0/2.0/2.2, and Delphi versions 1..7/9/10). As an example, here is the small t_intro.pas demo program, with Delphi use the -cc command line switch or load the t_ddemo.dpr Delphi GUI demo (it calculates Mersenne primes):

uses
  mp_types, mp_base, mp_numth;
var
  a: mp_int;
begin
  mp_init(a);
  mp_fib(271,a);
  writeln('Fibonacci(271) = ',mp_decimal(a));
  mp_2expt(a,479999);
  writeln('Number of decimal digits of 2^479999 = ', mp_radix_size(a,10));
  mp_clear(a);
end.

Fibonacci(271) = 193270471243015279782059101964580241188515112465021394429
Number of decimal digits of 2^479999 = 144496
Start with including mp_types in the uses statement (and the other units as needed) and declare the mp_int variables. Before working with mp_int variables, they must be created with one of the mp_init_... procedures. The allocated memory grows with the size of the numbers and is normally never decreased. If the variables are no longer needed they should be destroyed (the allocated memory is released) with one of the mp_clear_... routines; this should is essentially a must for local data.

Conversion from/to strings can be done with mp_read_decimal or mp_read_radix and mp_decimal or mp_radix_str and some other functions.

Unit mp_ratio.pas implements multi precision rational arithmetic. Most functions expect normalized rational input: the denominator must be greater than zero and relative prime to the numerator (zero=0/1). If the two fields are manipulated separately the s_mpr_normalize procedure can be used for normalization.

Unit mp_real.pas implements multi precision floating point arithmetic. Each mp_float number has its own precision (but normally most variables will have the same precision). Almost all functions expect normalized input: the bitsize of the mantissa is equal to the bit precision of the number. If the fields are manipulated separately the s_mpf_normalize procedures can be used for normalization. Although there is the function mpf_is_eq, equality of mp_floats is a fuzzy concept; mpf_is_eq_rel checks if the relative deviation is not greater than 1/2^(bitprec-1). mpf_reldev returns the value of the relative deviation and is used to detect failures and bad precisions in the t_chkfp test program.

Unit mp_calc.pas contains functions to handle mp_int expressions; the t_pfde.pas program tries to factor entered expressions. The t_calc.pas program is a simple interactive calculator (scriptable via redirection), it reads strings from stdin, parses, and evaluates them (and optionally assigns the values to one of the x,y,z variables). Here is the output of the ? help command:

Operators:  +  -  *  /  div  mod  ^  !  # (primorial)  % (same as mod)

Functions:  abs(a)       and(a,b)        binomial(a,b)  fermat(a)
            fib(a)       gcd(a,b)        invmod(a,b)    ispprime(a)
            jacobi(a,b)  kronecker(a,b)  lcm(a,b)       luc(a,b)
            max(a,b)     mersenne(a)     min(a,b)       nextprime(a)
            or(a,b)      prevprime(a)    random(a)      root(a,b)
            sqr(a)       sqrt(a)         sqrtmod(a,b)   xor(a,b)

Variables:  x  y  z,  Syntax: Var=expression[;] or Var=_[;]

Other    :  line terminator ";" shows only bitsize/chksum/time of result
            bin, dec, hex: set radix for result display
            "_<enter>" re-displays last result
            ".<enter>" displays time for last calculation

Unit mp_rcalc.pas and the demo calculator t_rcalc.pas handle mp_float expressions:

Operators:  +  -  *  /  ^           Constant:  pi

Functions:  abs(a)      acosh(a)    acosh1p(a)     agm(a,b)    arccos(a)
            arcsin(a)   arctan(a)   arctan2(y,x)   asinh(a)    atanh(a)
            CE(a)       CK(a)       cos(a)         cosh(a)     exp(a)
            expm1(a)    frac(a)     int(a)         ln(a)       ln1p(a)
            log10(a)    log2(a)     log(a,b)       max(a,b)    min(a,b)
            random(a)   sin(a)      sinh(a)        sqr(a)      sqrt(a)
            tan(a)      tanh(a)

Variables:  x  y  z,  Syntax: Var=expression[;] or Var=_[;]

Other    :  line terminator ";" shows only ldx/chksum/time of result
            bin, dec, hex: set radix for result display
            prec [n]: display/set bit precision
            "_<enter>" re-displays last result
            ".<enter>" displays time for last calculation


References
The MPArith package started as a Pascal port of the C libraries from [1] and [2] with background information from [3] and [5]. Some test and verification cases are taken from [14] or calculated with [12] and [13].

The number theoretic functions are influenced by [4], [8], and [10]. The random number generators used in the package come from [16] and [17]. Expressions parsing and evaluation is described in [19], [20], [21], and basic rational functions can be found in [14] and [22].

For me the most valuable general references in the list are [3], [4], [5] (with the focus on cryptography), and [10].
 [1]  LibTomMath, Tom St Denis, http://math.libtomcrypt.com/download.html
 [2]  MPI, M.J. Fromberger, http://spinning-yarns.org/michael/sw/
 [3]  D.E. Knuth: The Art of computer programming, Volume 2,
      Seminumerical Algorithms, 3rd ed., 1998
      http://www-cs-faculty.stanford.edu/~knuth/taocp.html
 [4]  O. Forster: Algorithmische Zahlentheorie, 1996
      http://www.mathematik.uni-muenchen.de/~forster/books/azth/algzth.html
 [5]  (HAC) A. Menezes, P. van Oorschot, S. Vanstone: Handbook of
      Applied Cryptography, 1996, http://www.cacr.math.uwaterloo.ca/hac/
 [6]  R.P. Brent, Factor: an integer factorization program for
      the IBM PC, Report TR-CS-89-23, October 1989,
      http://wwwmaths.anu.edu.au/~brent/pub/pub117.html
      http://wwwmaths.anu.edu.au/~brent/ftp/rpb117/rpb117.exe
 [7]  P. Ribenboim: The New Book of Prime Number Records, 3rd ed., 1995.
 [8]  Marcel Martin: NX - Numerics. Open source library of multiprecision
      numbers for Turbo Delphi and Free Pascal
      http://www.ellipsa.net/public/nx/index.html
 [9]  Wei Dai: Lucas Sequences in Cryptography, http://www.weidai.com/lucas.html
      or http://www.eskimo.com/~weidai/lucas.html
[10]  R. Crandall, C. Pomerance: Prime Numbers, A Computational Perspective, 2nd ed., 2005
[11]  Fast Factorial Functions: Peter Luschny's Homepage of Factorial Algorithms,
      http://www.luschny.de/math/factorial/FastFactorialFunctions.htm
[12]  PARI/GP, http://pari.math.u-bordeaux.fr/
[13]  ARIBAS: O. Forster's open source interpreter for big integer arithmetic,
      http://www.mathematik.uni-muenchen.de/~forster/sw/aribas.html
[14]  IMath library, M.J. Fromberger, http://spinning-yarns.org/michael/sw/
[15]  The GNU Multiple Precision Arithmetic Library, http://gmplib.org/
[16]  P. L'Ecuyer, "Maximally Equidistributed Combined Tausworthe Generators",
      Mathematics of Computation 65, 213 (1996), 203-213. Online version with
      corrections: http://www.iro.umontreal.ca/~lecuyer/myftp/papers/tausme.ps
[17]  B. Jenkins, ISAAC: a fast cryptographic random number generator
      http://burtleburtle.net/bob/rand/isaacafa.html
[18]  RFC 2313 - PKCS #1: RSA Encryption Version 1.5
      http://www.faqs.org/rfcs/rfc2313.html and
      RFC 3447 - PKCS #1: RSA Encryption Version 2.1
      http://www.faqs.org/rfcs/rfc3447.html
[19]  T. Norvell: Parsing Expressions by Recursive Descent,
      http://www.engr.mun.ca/~theo/Misc/exp_parsing.htm
[20]  G. Toal's tutorial pages OperatorPrecedence.html,CompilersOneOhOne.html,
      GrahamToalsCompilerDemo.html at http://www.gtoal.com/software/
[21]  T.R. Nicely's parser.c in factor1.zip from http://www.trnicely.net
      derived from GMP demo pexpr.c (see [15])
[22]  LiDIA - A Library for Computational Number Theory,
      http://www.cdc.informatik.tu-darmstadt.de/TI/LiDIA/
[23]  J. Shallit, J. Sorenson, A Binary Algorithm for the Jacobi Symbol,
      SIGSAM Bulletin, 27(1), 4-11, 1993; available online at
      http://euclid.butler.edu/~sorenson/papers/binjac.ps or
      http://citeseer.ist.psu.edu/article/shallit93binary.html
[24]  H. Cohen, A Course in Computational Algebraic Number Theory, 4th printing, 2000
[25]  IEEE P1363/Draft. Standard Specifications for Public Key Cryptography.
      Annex A (informative). Number-Theoretic Background.
      http://grouper.ieee.org/groups/1363/P1363/draft.html
[26]  A. Adler, J.E. Coury: The Theory of Numbers, 1995
[27]  T. Papanikolaou: libF - Eine lange Gleitpunktarithmetik, Diplomarbeit 1995,
      http://www.cdc.informatik.tu-darmstadt.de/reports/reports/papa.diplom.ps.gz
[28]  J.P. Sorenson, An analysis of Lehmer's Euclidean GCD algorithm,
      ACM International Symposium on Symbolic and Algebraic Computation, 1995
      http://euclid.butler.edu/~sorenson/papers/lehmer.pdf
[29]  V. Shoup, A Computational Introduction to Number Theory and Algebra, Version 2, 2008
      http://shoup.net/ntb/

Implemented functions
The MPArith package implements the following functions:

General integer functions

mp_2expt                computes a = 2^b
mp_abs                  absolute value, b = |a|
mp_add                  high level addition (handles signs)
mp_addmod               calculate d = a + b (mod c)
mp_add_d                single digit addition
mp_add_int              calculate c = a + b
mp_and                  calculate c = a and b
mp_binomial             calculate the binomial coefficient a = (n over k)
mp_bitsize              return the number of bits in a (index of highest bit), 0 if no bit is set
mp_checksum             return a checksum for a, -1 if mp_error<>MP_OKAY, -2 if not initialized
mp_chs                  change sign, b = -a
mp_clamp                trim unused digits
mp_clear[x]             clear [x] mp_ints (x=1..9)
mp_clear_multi          clear a vector of mp_ints
mp_clear_multi_p        clear a list of mp_ints given as a ptr vector
mp_clrbit               clear bit n of a, no action if out of range, (1 = bit 0)
mp_cmp                  compare a with a longint, return sign(a-b)
mp_cmp_d                compare a with a longint, return sign(a-b)
mp_cmp_int              compare a with a longint, return sign(a-b)
mp_cmp_mag              compare magnitude of two mp_ints (unsigned), return sign(|a|-|b|)
mp_cmp_mag_d            compare |a| with a digit, return sign(|a|-b)
mp_cnt_lsb              count the number of least significant bits which are zero
mp_copy                 copy an mp_int, b: = a
mp_cornacchia           Solve x^2 + d*y^2 = p, p prime, 0 < d <= p with Cornacchia's algorithm
mp_cornacchia4          Solve x^2 + |d|*y^2 = 4p, p prime, -4p < d < 0 with the modified Cornacchia algorithm
mp_coshmult             Forster's coshmult function, d=mod_coshmult(a,b,c)
mp_crt_setup            Calculate CRT coefficients c[i] for pairwise co-prime moduli m[i], i=0..n-1.
mp_crt_setupf           Calculate CRT coefficients, return true if c[i] are successfully calculated.
mp_crt_single           Calculate x with x mod m[i] = v[i], pairwise co-prime moduli m[i], i=0..n-1.
mp_crt_solve            Calculate x with x mod m[i] = v[i], coefficients c[i] must be precalulated with mp_crt_setup
mp_dec                  decrement an mp_int by 1
mp_dec_int              calculate a = a - b
mp_dfact                calculate double factorial a = n!!, a=n*(n-2)...
mp_div                  integer signed division, c^ = a div b, d^ = a mod b
mp_div_2                divide by 2, b = a/2
mp_div_2k               divide by 2^b; quotient in c, optional remainder in pd^, sign(pd^)=sign(a)
mp_div_d                single digit division, pc^ = a div b, d = a mod b
mp_div_int              integer signed division, pc^ = a div b, d = a mod b
mp_ecm_brent            find a factor f of n by Brent's ECM method, f=0 if no factor found
mp_ecm_factor           find a factor f of n by inversionless ECM method, f=0 if no factor found
mp_ecm_simple           simplified version of mp_ecm_factor with CMax=0, C1=ECM_C1Min
mp_exch                 exchange two mp_ints
mp_expt                 calculate c = a^b, b>=0
mp_exptmod              computes d=a^b mod c, c>0
mp_exptmod_d            computes d=a^b mod c, c>0, b longint
mp_exptmod_win          computes y=g^|e| mod p, internal sliding windows
mp_expt_d               calculate c = a^b
mp_expt_int             calculate c = a^b, b>=0
mp_fact                 calculate a = factorial(n) using Recursive Split, error if n > MaxFact
mp_fermat               return nth Fermat number, fn = 2^(2^n)+1 (MP_RANGE error for n>MaxFermat)
mp_fib                  calculate Fibonacci number fn=fib(n), fib(-n)=(-1)^(n-1)*fib(n)
mp_fib2                 calculate two Fibonacci numbers fn=fib(n), f1=fib(n-1), n>=0
mp_gcd                  calculate c = gcd(a,b) using the binary method
mp_gcd1                 calculate c = gcd(a,b) using the binary method, return true if c=1 and no error
mp_gcd_euclid           calculate c = gcd(a,b), non optimized Euclid
mp_gcd_ml               calculate u = gcd(a,b) using the Sorensen's Modified Lehmer method
mp_get_int              get the lower signed 31 bits of an mp_int
mp_grow                 grow an mp_int to a given size (new part is zerofilled)
mp_gr_mod               reduces x to x mod N, N > 1, using generalized reciprocal iteration.
mp_gr_setup             calculate the generalized reciprocal for N>1, @RN<>@N
mp_inc                  increment an mp_int by 1
mp_inc_int              calculate a = a + b
mp_init[x]              initialize [x] mp_ints (x=1..9)
mp_init_copy            creates a then copies b into it
mp_init_multi           initialize a vector of mp_ints
mp_init_multi_p         initialize a list of mp_ints given as a ptr vector
mp_init_set             initialize and set a digit
mp_init_set_int         initialize and set a to a longint
mp_init_size            initialize to a size of digits, rounded up to multiple of mp_allocprec
mp_invmod               compute c = a^-1 (mod b), b>0, via mp_xgcd, MP_UNDEF error if there is no inverse
mp_invmodf              compute c = a^-1 (mod b), b>0, via mp_xgcd, return true if inverse exists
mp_is0                  initialized and a = 0
mp_is1                  initialized and a = 1
mp_is1a                 initialized and abs(a) = 1
mp_isbit                test if bit n of a is set, (1 = bit 0)
mp_iseven               initialized and even
mp_isMersennePrime      Lucas-Lehmer test for Mersenne number m=2^p-1
mp_isodd                initialized and odd
mp_iszero               initialized and zero
mp_is_eq                return a = b
mp_is_ge                return a >= b
mp_is_gt                return a > b
mp_is_le                return a <= b
mp_is_lt                return a < b
mp_is_ne                return a <> b
mp_is_pow2              check if |a| is a power of 2, if true, return n with |a|=2^n
mp_is_pow2_d            check if d is power of 2, if true, return n with d=2^n
mp_is_power             Calculate smallest prime p with a=b^p; p=1,b=a if a is no power
mp_is_pprime            Test if a is prime (BPSW pseudo prime if a>2^32)
mp_is_pprime_ex         Test if a is prime (BPSW pseudo prime if a>2^32); trial division up to smax
mp_is_primepower        return true if a=b^k, b prime, k>1, otherwise false and a=b^k, k=1 if no power
mp_is_slpsp             strong Lucas pseudo prime test for a
mp_is_spsp              strong pseudo prime test of n to base of a>1
mp_is_spsp_d            strong pseudo prime test of n to mp_digit base of a>1
mp_is_square            test if a is square
mp_is_square2           test if a is square, return sqrt(a) if a is a square and psqrt<>nil
mp_jacobi               computes the Jacobi/Legendre symbol (a|n), n: odd and > 2
mp_jacobi_lm            computes the Jacobi/Legendre symbol (a|n), n: odd and > 2
mp_jacobi_ml            computes the Jacobi/Legendre symbol (a|n), n: odd and > 2
mp_karatsuba_mul        calculate c = |a| * |b| using Karatsuba Multiplication
mp_karatsuba_sqr        Karatsuba squaring, computes b = a*a using three half size squarings
mp_kronecker            computes the Kronecker symbol (a|n)
mp_lcm                  computes least common multiple as |a*b|/(a, b)
mp_lshd                 shift left a certain amount of digits
mp_lucas                calculate Lucas number lk=luc(k), luc(-k)=(-1)^k*luc(k)
mp_lucas2               calculate two Lucas numbers lk=luc(k), l1=luc(k-1), k>=0
mp_lucasv               calculate v[k] of Lucas V sequence for p,q, p^2-4q <>0
mp_lucasv2              calculate v=v[k], w=v[k+1] of Lucas V sequence
mp_lucasv2p             calculate v[k], v[k+1] of Lucas V sequence, pointer version
mp_lucasvmod            calculate v[k] mod n of Lucas V sequence for p,q
mp_makeodd              return b,s with a = 2^s*b if a<>0, b=0,s=-1 otherwise
mp_mersenne             return nth Mersenne number, mn = 2^n-1, MP_RANGE err for n>MaxMersenne
mp_miller_rabin         Miller-Rabin test of n, security parameter t, from HAC p. 139 Alg.4.24
mp_mod                  calculate c = a mod b, 0 <= c < b
mp_mod_2k               calculate c = a mod 2^b, 0 <= c < 2^b
mp_mod_d                calculate c = a mod b, 0 <= c < b (digit version)
mp_mod_int              calculate c = a mod b
mp_montgomery_calcnorm  calculate R=B^n, n=number of digits in m, B=2^DIGIT_BIT (=mp_int radix)
mp_montgomery_reduce    computes xR^-1 == x (mod N) via Montgomery Reduction
mp_montgomery_setup     calculates rho for Montgomery reduction
mp_mul                  high level multiplication (handles sign)
mp_mulmod               calculate d = a * b (mod c)
mp_mul_2                multiply by 2, b = a*2
mp_mul_2k               shift left by a certain bit count [synonym for mp_shl
mp_mul_d                multiply by a digit
mp_mul_int              multiply by a 32 bit integer
mp_nextprime            next prime >= n, 2 if n<=2
mp_n_root               calculate the n'th root of a, a must >= 0 if n is even. b=0 if n<1
mp_OddProd              calculate p=prod(2*i+1),i=a+1...b;  p=1 if a<=b
mp_or                   calculate c = a or b
mp_pell                 calculate the smallest solution of x^2 - d*y^2 = 1
mp_pell4                calculate the smallest solution of x^2 - d*y^2 = r, r in [-4,+4,-1,+1]
mp_pollard_pm1          find a factor f of N with p-1 method, f=0 if no factor found
mp_pollard_rho          find a factor f of N, f=0 if no factor found
mp_popcount             get population count = number of 1-bits in a
mp_ppexpo               product of primes B0 < p <= B1 and integers B0 < n*n <= B1
mp_prevprime            previous prime <= n, 0 if n<2
mp_primorial            Primorial of n;  a = n# = product of primes <= n
mp_prod_int             calculate a = product of first n elements of longint array b
mp_radix_size           return size of ASCII representation
mp_rand                 makes a pseudo-random int of a given size
mp_rand_bits            make a pseudo-random mp_int of a given bit size
mp_rand_bits_ex         make pseudo-random a with bitsize <= bits, if sethi highest bit is set
mp_rand_ex              make a pseudo-random mp_int, if not sethi then a[digits-1] my be zero
mp_rand_prime           generate random prime of bitsize, pt: prime type of p
mp_rand_radix           makes a pseudo-random int of size radix^digits
mp_read_decimal         read an mp_int from a decimal ASCII pchar
mp_read_decimal_str     read an mp_int from a decimal ASCII string[255]
mp_read_radix           read an mp_int from a ASCII pchar in given radix
mp_read_radix_str       read an mp_int from a ASCII string[255] in given radix
mp_read_signed_bin      read signed bin, big endian, first byte is 0=positive or 1=negative
mp_read_unsigned_bin    reads a unsigned mp_int, assumes the msb is stored first [big endian]
mp_reduce               reduces x mod m, assumes x < m^2, mu is precomputed via mp_reduce_setup
mp_reduce_2k            reduces a mod n where n is of the form 2^p-d
mp_reduce_2k_setup      determine setup value d for unrestricted diminished radix reduction, a>=0
mp_reduce_is_2k         determine if mp_reduce_2k can be used
mp_reduce_setup         pre-calculate the value required for Barrett reduction
mp_reverse              reverse an array of char, used for radix code
mp_rshd                 shift right a certain amount of digits
mp_set                  set a to digit b
mp_setbit               set bit n of a, error if n<0 or n>MP_MAXBIT (1 = bit 0)
mp_set_int              set a to a longint
mp_set_pow              set a to b^c, a=0 for c<0
mp_shl                  Shift left a, c = a*2^b; c=a if b<=0
mp_shr                  Shift right a, c = a/2^b; c=a if b<=0
mp_shr1                 divide a by 2, a = a/2
mp_sign                 return sign(a): -1 if a<0, 0 if a=0, +1 if a>0
mp_signed_bin_size      get the size in bytes for an signed equivalent
mp_small_factor         calculate small digit factor or 0, start with f0, f will be < $FFFF
mp_sqr                  computes b = a*a
mp_sqrmod               calculate d = a * a (mod c)
mp_sqrt                 computes b = sqrt(a)
mp_sqrtmod              calculate square root b of a with b*b = a mod p, p prime, with Jacobi check
mp_sqrtmod2k            calculate square root b of an integer a with b*b = a mod 2^k
mp_sqrtmodp2            calculate square root b of a with b*b = a mod p^2, p odd prime
mp_sqrtmodpk            calculate square root b of a with b*b = a mod p^k, p odd prime
mp_sqrtmodpq            calculate square roots +x,-x,+y,-y of a mod (pq), p,q different primes
mp_sqrtmod_ex           calculate square root b of a with b*b = a mod p, optional Jacobi check
mp_sub                  high level subtraction (handles signs)
mp_submod               calculate d = a - b (mod c)
mp_sub_d                single digit subtraction
mp_sub_int              calculate c = a - b
mp_todouble             convert a to double, +-inf if too large
mp_todouble_ex          convert a*2^x to double, +-inf if too large
mp_toextended           convert a to extended, +-inf if too large
mp_toextended_ex        convert a*2^x to extended, +-inf if too large
mp_toradix              stores mp_int as a ASCII string in a given radix, better use mp_toradix_n
mp_toradix_n            stores a mp_int as a ASCII string in a given radix (2..MAXRadix)
mp_to_signed_bin_n      store in signed big-endian format, max n bytes, return no. of bytes stored
mp_to_unsigned_bin_n    store in unsigned big-endian format, max n bytes, return no. of bytes stored
mp_unsigned_bin_size    get the size in bytes for an unsigned equivalent
mp_williams_pp1         find a factor f of N with William's p+1 method, f=0 if no success
mp_xgcd                 extended gcd algorithm, calculate  a*p1^ + b*p2^ = p3^ = gcd(a,b)
mp_xgcd_bin             extended binary gcd, calculate a*p1^ + b*p2^ = p3^ = gcd(a,b)
mp_xor                  calculate c = a xor b
mp_zero                 set a to zero

s_mp_add                low level addition c=a+b, based on HAC pp.594, algorithm 14.7
s_mp_add_d              single digit addition, no init check, b<>0
s_mp_binom_l            internal binomial coefficient for small k, no init check
s_mp_checksum           update Adler32 checksum with Msg data; init with adler=1
s_mp_chs                change sign, assumes that a is initialized
s_mp_div                integer signed division, no init check
s_mp_expt_dl            calculate c = a^b, returns 0 for b<0
s_mp_expt_wl            calculate c = a^b, returns 0 for b<0
s_mp_mod_2k             calculate c = a mod 2^b, -(|a| mod 2^b) if a < 0
s_mp_mul_digs           multiplies |a| * |b| and only computes up to digs digits of result
s_mp_mul_high_digs      multiplies |a| * |b| and does not compute the lower digs digits
s_mp_mul_int            multiply by a 32 bit integer, c=a*b, tmp is an initialized temporary
s_mp_read_radix         read an ASCII pchar in a given radix into a, breaks on sep and #0, no init check
s_mp_sqr                low level squaring, b = a*a, HAC pp.596-597, algorithm 14.16
s_mp_sqrtmod2k          calculate unique square root b of an odd integer a with b*b = a mod 2^k
s_mp_sub                low level subtraction (assumes |a| > |b|), HAC pp.595 algorithm 14.9
s_mp_sub_d              single digit subtraction, no init check, b<>0
s_mp_toradix_n          convert an mp_int to ASCII for a given radix, plus=show '+', no init check
s_mp_write_radix        write radix representation to file tf


RSA related functions

mp_i2osp                Convert a nonnegative mp_int to an octet string of a specified length
mp_i2pchar              Convert a nonnegative mp_int to a pchar with a specified max length
mp_os2ip                Convert an octet string of length ilen to a nonnegative mp_int
mp_pkcs1v15_decode      EME-PKCS1-v1_5 decoding; true if decoding is successful, false otherwise
mp_pkcs1v15_decrypt     Decrypt a message using RSA and EME-PKCS1-v1_5 padding
mp_pkcs1v15_decrypt2    Decrypt a message using RSA/CRT and EME-PKCS1-v1_5 padding
mp_pkcs1v15_encode      EME-PKCS1-v1_5 encoding; true if decoding is successful, false otherwise
mp_pkcs1v15_encrypt     Encrypt a message using RSA and EME-PKCS1-v1_5 padding
mp_pkcs1v15_maxlen      Maximum message length for RSA modulus n using PKCS1-v1_5 encoding
mp_rsadp                Basic RSA decryption operation, m=c^d mod n
mp_rsadp2               Basic RSA decryption operation for private key CRT record
mp_rsaep                Basic RSA encryption operation, c=m^e mod n
mp_rsa_calc_d           Get RSA decryption exp. d from private key record and encryption exp. e
mp_rsa_calc_npq         Generate RSA primes p,q; q<p, n=p*q, osize: octet size of n
mp_rsa_calc_private     Calculate remaining fields of private RSA/CRT key from e,p,q
mp_rsa_clear_private    Clear fields of private RSA/CRT key
mp_rsa_init_private     Initialize fields of private RSA/CRT key
mp_rsa_keygen1          Basic RSA private key pair (n, d) generation; osize: octet size of n
mp_rsa_keygen2          Generate private RSA/CRT key prk and modulus n; osize: octet size of n



mp_calc/mp_rcalc related functions

mp_calculate            parse and evaluate string psz
mp_calc_errorstr        translates known error codes
mp_clear_eval           clear the mp_ints of evr
mp_clear_expr           release memory used by e
mp_eval                 evaluate expression tree e, result in evr
mp_eval_mod             evaluate expression tree e mod m, result in evr
mp_init_eval            initialize the mp_ints of evr
mp_parse                parse string psz into expression tree e
mpf_calculate           parse and evaluate string psz
mpf_calc_errorstr       translate known error codes
mpf_clear_eval          clear the mp_floats of evr
mpf_clear_expr          release memory used by e
mpf_eval                evaluate expression tree e, result in evr
mpf_init_eval           initialize the mp_floats of evr
mpf_parse               parse string psz into expression



Functions specific to I/O or Pascal

bigalloc                allocate heap > 64K, returns nil if error, no diagnostics
mp_adecimal             convert to decimal ansistring, max 65000 digits
mp_ahex                 convert to hex ansistring, max 65000 digits
mp_alloc                allocate and zero heap, returns nil if error
mp_arctanw              compute sum := mul*arctan(1/x), to prec Radix digits, word version
mp_atanhw               compute sum := mul*atanh(1/x), to prec Radix digits
mp_decimal              convert to decimal, max 255 digits
mp_div_w                divide a by a single word b
mp_dumpa                dump all fields of a
mp_dumpu                dump used fields of a
mp_dump_diagctr         dump diagnostic counters
mp_dump_meminfo         writes mp_memstat to output (if MPC_Diagnostic is defined)
mp_dump_memused         writes mp_memused to output (if MPC_Diagnostic is defined)
mp_freemem              deallocate heap if p<>nil, p will be set to nil
mp_getmem               allocate heap, returns nil if error
mp_get_allocprec        returns current value of mp_allocprec
mp_hex                  convert to hex string, max 255 digits
mp_init_prim            initialize to a size of digits, rounded up to multiple of mp_allocprec
mp_is_longint           test if a fits into longint, if true set b := a
mp_memused              return total allocated memory
mp_mod_w                calculate r = a mod b for a single word b
mp_mul_w                multiply by a word
mp_not_init             sanity check if a is initialized, does not catch all cases!
mp_not_init_multi       sanity check if all elements of a are initialized, does not catch all cases!
mp_output_decimal       write decimal representation to output
mp_output_radix         write radix representation to output
mp_radix_astr           convert to radix representation ansistring, max 65000 digits
mp_radix_str            convert to radix representation, max 255 digits
mp_random_byte          returns a random byte
mp_random_digit         returns a random mp_digit
mp_random_int           returns a random signed longint
mp_random_long          returns a random positive longint
mp_random_radix         returns a random word in the range 0..radix-1
mp_random_randomize     initialize PRNG via randomize
mp_random_read          read len bytes from the PRNG to dest
mp_random_seed          initialize PRNG with array of longint
mp_random_seed1         initialize PRNG with a longint
mp_random_word          returns a random word
mp_realloc              reallocate heap to new size, if newsize>oldsize the new allocated space is zerofilled
mp_reset_diagctr        resets diagnostic counters to 0
mp_result               return and reset mp_error
mp_set_allocprec        set new alloc prec 8..64, will be rounded up to power of 2
mp_set_progress         makes PP the new progress proc
mp_set_short            set a to a shortint
mp_set_w                set a to a word
mp_writeln              writeln a to stdout with leading msg
mp_write_decimal        write decimal representation to file tf
mp_write_radix          write radix representation to file tf
set_mp_error            Set error variable


FPU and 32 bit functions

add32_ovr               add z=x+y with overflow detection
bitsize32               return the number of bits in a (index of highest bit), 0 if no bit is set
DblNaN                  return double NaN (Not a Number)
DblNegInf               return negative double infinity
DblPosInf               return positive double infinity
FindFirstPrime32        Find first prime >= n and initialize ctx for FindNextPrime32
FindNextPrime32         Find next 32 bit prime (DWORD interpretation, see note), success if ctx.prime<>0
frexpd                  returns m,e with d=m*2^e and 0.5 <= abs(m) < 1  (double)
frexpx                  returns m,e with d=m*2^e and 0.5 <= abs(m) < 1  (extended)
GCD32                   calculate GCD of two longints
GCD32U                  calculate GCD of two longints (DWORD interpretation)
invmod32                return a^-1 mod b, b>1. Result is 0 if gcd(a,b><>1 or b<2
IsPow2_w                check if w is power of 2, if true, return n with w=2^n
IsPrime16               test if N is prime
IsPrime32               test if longint N is prim
isqrt32                 returns floor(sqrt(abs(a))
is_spsp32               Strong pseudo prime test for N with base A
is_spsp32A              Strong pseudo prime test for N with a number of bases
jacobi32                computes the Jacobi/Legendre symbol (a|b), b: odd and > 2
ldexpd                  returns d*2^e (double)
ldexpx                  returns d*2^e (extended)
nextprime32             next 32 bit prime >= n (DWORD interpretation, see note)
nextprime32_array       fill an array with the next 32 bit primes >= n (DWORD interpretation)
popcount16              get population count = number of 1-bits in a word
popcount32              get population count = number of 1-bits in a longint
prevprime32             previous 32 bit prime <= n, prevprime32(0)=0, (DWORD interpretation)



Rational arithmetic functions

mpr_abs                 absolute value, b = |a|
mpr_add                 add two mp_rats: c = a+b
mpr_add_mpi             add mp_int to mp_rat: c = a+b
mpr_adecimal            convert to radix representation ansistring, max 65000 digits
mpr_ceil                return b := ceil(a)
mpr_checksum            return a checksum for a, -1 if mp_error<>MP_OKAY, -2 if not initialized
mpr_chs                 change sign, b = -a
mpr_clear               free an mp_rat
mpr_clear[x]            clear [x] mp_rats
mpr_clear_multi         clear a vector of mp_rats
mpr_clear_multi_p       clear a list of mp_rats given as a ptr vector
mpr_cmp                 compare two mp_rats (signed), return sign(a-b)
mpr_cmp_mag             compare magnitude of two mp_rats (unsigned), return sign(|a|-|b|)
mpr_copy                copy an mp_rat, b: = a
mpr_decimal             convert a to decimal, max 255 chars
mpr_div                 divide two mp_rats: c = a/b, b<>0
mpr_div_mpi             divide mp_rat by mp_int: c = a/b, b<>0
mpr_exch                exchange two mp_rats
mpr_expt                calculate c = a^b, a<>0 for b<0, 0^0=1
mpr_floor               return b := floor(a)
mpr_frac                return b := frac(a) = a - trunc(a)
mpr_init                initialize an mp_rat
mpr_init[x]             initialize [x] mp_rats
mpr_init_copy           creates a then copies b into it
mpr_init_multi          initialize a vector of mp_rats.
mpr_init_multi_p        initialize a list of mp_rats given as a ptr vector.
mpr_init_size           initialize an mp_rat to given number of digits
mpr_inv                 invert an mp_rat, b = 1/a, a<>0
mpr_is0                 return a=0
mpr_is_eq               return a = b
mpr_is_ge               return a >= b
mpr_is_gt               return a > b
mpr_is_le               return a <= b
mpr_is_lt               return a < b
mpr_is_mp               return true if a is initialized and a.den=1
mpr_is_ne               return a <> b
mpr_mul                 multiply two mp_rats: c = a*b
mpr_mul_mpi             multiply mp_rat and mp_int: c = a*b
mpr_not_init            sanity check if a is initialized, does not catch all cases!
mpr_output_decimal      write decimal representation to output
mpr_output_radix        write radix representation to output
mpr_radix_astr          convert to radix representation ansistring, max 65000 digits
mpr_radix_size          return size of ASCII representation (incl. sign and #0)
mpr_radix_str           convert to radix representation, max 255 digits
mpr_read_decimal        read a ASCII decimal string into a. str may contain a single '/'
mpr_read_double         convert d to an mp_rat
mpr_read_float_decimal  read a ASCII float decimal string into a. str may contain a single '.'
mpr_read_float_radix    read a ASCII float radix string into a. str may contain a single '.'
mpr_read_radix          read a ASCII radix string into a. str may contain a single '/'
mpr_round               return b := round(a), round(-a)=-round(a), round(n+0.5)=n+1
mpr_set                 set a to n/d, d<>0
mpr_set1                set a to n/1
mpr_set_int             set a to n/d, d<>0
mpr_sub                 subtract two mp_rats: c = a-b
mpr_sub_mpi             subtract mp_int from mp_rat: c = a-b
mpr_todouble            convert a to double, +-inf if too large, 0 if to small
mpr_tofloat_astr        convert to float representation with prec, max 65000 digits
mpr_tofloat_n           convert to float format for a given radix, prec digits after '.'
mpr_tofloat_str         convert to float representation, prec digits after '.', max 255 chars
mpr_toradix_n           convert an mp_rat to an ASCII string for a given radix (2..MAXRadix)
mpr_trunc               return b := trunc(a)
mpr_write_decimal       write decimal representation to file tf
mpr_write_radix         write radix representation to file tf
mpr_zero                set a to zero

s_mpr_add_sub           add or subtract two mp_rats
s_mpr_normalize         normalizes a, assumes a is initialized
s_mpr_normsign          make denominator positive, assumes a is initialized
s_mpr_qr                calculate q := abs(num) div abs(den); r := abs(num) mod abs(den), no init check
s_mpr_tofloat_n         convert to float format for a given radix, prec digits after '.', assumes a is initialized



Floating point functions

mpf_10expt              calculate b = 10^a
mpf_2expt               calculate b = 2^a
mpf_abs                 absolute value, b = |a|
mpf_acosh               calculate b = acosh(a), acosh(-a)=acosh(a), |a| >= 1. Note: for a near 1 the function acosh1p(a-1) should be used
mpf_acosh1p             calculate b = acosh(1+a), a>=0
mpf_add                 calculate c = a+b
mpf_add_ext             calculate c = a+b
mpf_add_mpi             calculate c = a+b
mpf_adecimal            convert to decimal representation with ndd decimal digits, max 65000 digits
mpf_agm                 calculate c = AGM(|a|,|b|)
mpf_arccos              calculate b = arccos(a), |a| <= 1
mpf_arcsin              calculate b = arcsin(a), |a| <= 1
mpf_arctan              calculate b = arctan(a)
mpf_arctan2             calculate a = arctan(y/x) with special treatment for zero x or y. a is the principal value of arg(x + i*y)
mpf_asinh               calculate b = asinh(a)
mpf_atanh               calculate b = atanh(a), |a| < 1
mpf_ccell1              Complementary complete elliptic integrals of the first kind
mpf_ccell12             Complementary complete elliptic integrals of the 1st and 2nd kind
mpf_ccell2              Complementary complete elliptic integral of the 2nd kind CE = CE(k)
mpf_checksum            return a checksum for a, -1 if mp_error<>MP_OKAY, -2 if not initialized
mpf_chg_prec            change bitprec of a to newprec
mpf_chs                 change sign, b = -a
mpf_clear[x]            clear [x] mp_floats
mpf_cmp                 compare two mp_floats, return sign(a-b)
mpf_cmp_mag             compare magnitude of two mp_floats, return sign(|a|-|b|)
mpf_copy                copy a to b with b.bitprec=a.bitprec
mpf_copyp               copy a to b, preserve b.bitprec
mpf_cos                 calculate b = cos(a)
mpf_cosh                calculate b = cosh(a), a < 2^31 * ln(2)
mpf_decimal             convert to decimal representation with ndd>0 decimal digits, max 255 chars
mpf_div                 calculate c = a/b
mpf_div_d               calculate c = a/b
mpf_div_ext             calculate c = a/b
mpf_div_int             calculate c = a/b
mpf_div_mpi             calculate c = a/b
mpf_exch                exchange two mp_floats (including bitprec!!)
mpf_exp                 calculate b = exp(a)
mpf_expm1               calculate b = exp(a)-1, a < 2^31 * ln(2); special version for small a
mpf_expt                calculate c = a^b, a>0
mpf_expt_int            calculate c = a^b
mpf_frac                set b to the fractional part of a; frac(x)=x-int(x)
mpf_get_default_prec    return current default precision, initial=240
mpf_iexpt               calculate c = a^b, a>0
mpf_init                initialize an mp_float with default precision
mpf_initp               initialize an mp_float with bit precision prec
mpf_initp[x]            initialize [x] mp_floats with bit precision prec
mpf_initp_multi_p       initialize with bit precision prec a list of mp_floats given as a pointer
mpf_int                 set b to the integer part of a; i.e. is b rounded toward zero
mpf_inv                 calculate b = 1/a
mpf_is0                 return true if a=0
mpf_is1                 return true if a=1
mpf_is1a                return true if abs(a)=1
mpf_is_eq               return a = b
mpf_is_eq_rel           Check if |a-b| <= r*2^(1-b.bitprec);  r=1 if b=0, r=|b| otherwise
mpf_is_ge               return a >= b
mpf_is_gt               return a > b
mpf_is_le               return a <= b
mpf_is_lt               return a < b
mpf_is_ne               return a <> b
mpf_ln                  calculate b=ln(a), a>0
mpf_ln1p                calculate b = ln(1+a), a>-1; special version for small a
mpf_log10               calculate b = log10(a), a>0
mpf_log2                calculate b = log2(a), a>0
mpf_mul                 calculate c = a*b
mpf_mul_2k              calculate b = a*2^k
mpf_mul_d               multiply by a digit
mpf_mul_ext             calculate c = a*b
mpf_mul_int             multiply by a 32 bit integer
mpf_mul_mpi             calculate c = a*b
mpf_not_init            sanity check if a is initialized, does not catch all cases!
mpf_output_decimal      write an mp_float to output using decimal scientific notation
mpf_output_radix        write an mp_float to output using radix scientific notation
mpf_random              set to a random number uniformly distributed in [0,1)
mpf_read_decimal        read a from ASCII float decimal string
mpf_read_hex            read a from ASCII float hexadecimal string
mpf_read_radix          read a from ASCII float radix string
mpf_set0                set a=0, a.bitprec is preserved
mpf_set1                set a=1, a.bitprec is preserved
mpf_set_default_prec    set new default precision
mpf_set_ext             set a to an extended. Error if x = NAN or INF
mpf_set_int             set a to a longint
mpf_set_ln10            set a to ln(10), preserve a.bitprec
mpf_set_ln10p           set a to ln(10) with bit precision prec
mpf_set_ln10p2k         set a to ln(10)*2^k with bit precision prec
mpf_set_ln2             set a to ln(2), preserve a.bitprec
mpf_set_ln2p            set a to ln(2) with bit precision prec
mpf_set_ln2p2k          set a to ln(2)*2^k with bit precision prec
mpf_set_mpi             set a to an mp_int
mpf_set_mpi2k           set a to m*2^e (build a from mantissa and exponent)
mpf_set_pi              set a to pi, preserve a.bitprec
mpf_set_pi2k            set a to pi*2^k, preserve a.bitprec
mpf_set_pip             set a to pi with bit precision prec
mpf_set_pip2k           set a to pi*2^k with bit precision prec
mpf_sin                 calculate b = sin(a)
mpf_sinh                calculate b = sinh(a), a < 2^31 * ln(2)
mpf_sqr                 calculate b = a*a
mpf_sqrt                calculate b = sqrt(a)
mpf_sub                 calculate c = a-b
mpf_sub_mpi             calculate c = a-b
mpf_tan                 calculate b = tan(a)
mpf_tanh                calculate b = tanh(a)
mpf_todecimal_n         convert an mp_float to decimal scientific notation
mpf_todouble            convert a to double, +-inf if too large
mpf_toextended          convert a to extended, +-inf if too large
mpf_tohex_n             convert an mp_float to hexadecimal scientific notation
mpf_toradix_n           convert an mp_float to radix scientific notation
mpf_trig                calculate pc^=cos(a), ps^=sin(a), pt^=tan(a) if pointers are <> nil
mpf_trunc               truncate mp_float to mp_int
mpf_writeln             writeln a to output with leading msg
mpf_write_decimal       write an mp_float to file tf using decimal scientific notation
mpf_write_radix         write an mp_float to file tf using radix scientific notation

s_mpf_abs               absolute value of an mp_float, no init check
s_mpf_add1              calculate b = a+1; no init checks
s_mpf_agm               calculate c = AGM(|a|,|b|), if ps<>nil set ps^=sum(2^k*c_k^2); ps<>@c, no init check
s_mpf_ccell12           Complementary complete elliptic integrals of the 1st and 2nd kind using AGM algorithm; with init checks
s_mpf_chs               change sign of an mp_float, no init check
s_mpf_cmp_mag           compare magnitude of two mp_floats, return sign(|a|-|b|), no init check
s_mpf_dec1              calculate a = a-1; no init check
s_mpf_inc               calculate x = x+y
s_mpf_inc1              calculate a = a+1; no init check
s_mpf_incexp            increment a.exponent by x, do nothing if a.mantissa=0
s_mpf_incf              calculate x = x+y; return true if x is changed
s_mpf_is0               return true if a=0, no init checks
s_mpf_is_ge0            return true if a>=0, no init checks
s_mpf_is_le0            return true if a<=0, no init checks
s_mpf_is_neg            return true if a<0, no init checks
s_mpf_ldx               calculate ldx(a) = a.exponent + mp_bitsize(a.mantissa) = floor(log2(|a|))+1 for a<>0
s_mpf_mod_pi2k          calculate b = a mod pi*2^k, c in [0,pi*2^k)
s_mpf_normalize         normalize an mp_float
s_mpf_normalizep        normalize an mp_float to bit precision prec
s_mpf_toradix_n         convert an mp_float to radix scientific notation