Arbitrary-precision arithmetic Explanation

Question

I'm trying to learn C and have come across the inability to work with REALLY big numbers (i.e., 100 digits, 1000 digits, etc.). I am aware that there exist libraries to do this, but I want to attempt to implement it myself.

I just want to know if anyone has or can provide a very detailed, dumbed down explanation of arbitrary-precision arithmetic.

Answer 1

While re-inventing the wheel is extremely good for your personal edification and learning, its also an extremely large task. I don't want to dissuade you as its an important exercise and one that I've done myself, but you should be aware that there are subtle and complex issues at work that larger packages address.

For example, multiplication. Naively, you might think of the 'schoolboy' method, i.e. write one number above the other, then do long multiplication as you learned in school. example:

      123
    x  34
    -----
      492
+    3690
---------
     4182

but this method is extremely slow (O(n^2), n being the number of digits). Instead, modern bignum packages use either a discrete Fourier transform or a Numeric transform to turn this into an essentially O(n ln(n)) operation.

And this is just for integers. When you get into more complicated functions on some type of real representation of number (log, sqrt, exp, etc.) things get even more complicated.

If you'd like some theoretical background, I highly recommend reading the first chapter of Yap's book, "Fundamental Problems of Algorithmic Algebra". As already mentioned, the gmp bignum library is an excellent library. For real numbers, I've used mpfr and liked it.

Answer 2

You can do it with high school level of mathematics. Though more advanced algorithms are used in reality. So for example to add two 1024-byte numbers :

unsigned char first[1024], second[1024], result[1025];
unsigned char carry = 0;
unsigned int  sum   = 0;

for(size_t i = 0; i < 1024; i++)
{
    sum = first[i] + second[i] + carry;
    carry = sum - 255;
}

result will have to be bigger by one place in case of addition to take care of maximum values. Look at this :

9
   +
9
----
18

TTMath is a great library if you want to learn. It is built using C++. The above example was silly one, but this is how addition and subtraction is done in general!

A good reference about the subject is Computational complexity of mathematical operations. It tells you how much space is required for each operation you want to implement. For example, If you have two N-digit numbers, then you need 2N digits to store the result of multiplication.

As Mitch said, it is by far not an easy task to implement! I recommend you take a look at TTMath if you know C++.

Answer 3

One of the ultimate references (IMHO) is Knuth's TAOCP Volume II. It explains lots of algorithms for representing numbers and arithmetic operations on these representations.

@Book{Knuth:taocp:2,
   author    = {Knuth, Donald E.},
   title     = {The Art of Computer Programming},
   volume    = {2: Seminumerical Algorithms, second edition},
   year      = {1981},
   publisher = {\Range{Addison}{Wesley}},
   isbn      = {0-201-03822-6},
}

Answer 4

You do it in basically the same way you do with pencil and paper...

• The number is to be represented in a buffer (array) able to take on an arbitrary size (which means using malloc and realloc) as needed
• you implement basic arithmetic as much as possible using language supported structures, and deal with carries and moving the radix-point manually
• you scour numeric analysis texts to find efficient arguments for dealing by more complex function
• you only implement as much as you need.

Typically you will use as you basic unit of computation

• bytes containing with 0-99 or 0-255
• 16 bit words contaning wither 0-9999 or 0--65536
• 32 bit words containing...
• ...

as dictated by your architecture.

The choice of binary or decimal base depends on you desires for maximum space efficiency, human readability, and the presence of absence of Binary Coded Decimal (BCD) math support on your chip.

Answer 5

Assuming that you wish to write a big integer code yourself, this can be surprisingly simple to do, spoken as someone who did it recently (though in MATLAB.) Here are a few of the tricks I used:

• I stored each individual decimal digit as a double number. This makes many operations simple, especially output. While it does take up more storage than you might wish, memory is cheap here, and it makes multiplication very efficient if you can convolve a pair of vectors efficiently. Alternatively, you can store several decimal digits in a double, but beware then that convolution to do the multiplication can cause numerical problems on very large numbers.

• Store a sign bit separately.

• Addition of two numbers is mainly a matter of adding the digits, then check for a carry at each step.

• Multiplication of a pair of numbers is best done as convolution followed by a carry step, at least if you have a fast convolution code on tap.

• Even when you store the numbers as a string of individual decimal digits, division (also mod/rem ops) can be done to gain roughly 13 decimal digits at a time in the result. This is much more efficient than a divide that works on only 1 decimal digit at a time.

• To compute an integer power of an integer, compute the binary representation of the exponent. Then use repeated squaring operations to compute the powers as needed.

• Many operations (factoring, primality tests, etc.) will benefit from a powermod operation. That is, when you compute mod(a^p,N), reduce the result mod N at each step of the exponentiation where p has been expressed in a binary form. Do not compute a^p first, and then try to reduce it mod N.

Answer 6

Here's a simple ( naive ) example I did in PHP.

I implemented "Add" and "Multiply" and used that for an exponent example.

http://adevsoft.com/simple-php-arbitrary-precision-integer-big-num-example/

Code snip

// Add two big integers
function ba($a, $b)
{
    if( $a === "0" ) return $b;
    else if( $b === "0") return $a;

    $aa = str_split(strrev(strlen($a)>1?ltrim($a,"0"):$a), 9);
    $bb = str_split(strrev(strlen($b)>1?ltrim($b,"0"):$b), 9);
    $rr = Array();

    $maxC = max(Array(count($aa), count($bb)));
    $aa = array_pad(array_map("strrev", $aa),$maxC+1,"0");
    $bb = array_pad(array_map("strrev", $bb),$maxC+1,"0");

    for( $i=0; $i<=$maxC; $i++ )
    {
        $t = str_pad((string) ($aa[$i] + $bb[$i]), 9, "0", STR_PAD_LEFT);

        if( strlen($t) > 9 )
        {
            $aa[$i+1] = ba($aa[$i+1], substr($t,0,1));
            $t = substr($t, 1);
        }

        array_unshift($rr, $t);
     }

     return implode($rr);
}