15
\$\begingroup\$

For my mathematical research, I am to deal with very large numbers and I know that the standard C types are not going to stand a chance against them due to their limited capacity, so I decided to write a super calculator in C language in hope that I won't have to resort to a super computer.

I started with an add function that I want to be reviewed for bugs if any, and possible improvement in performance and general programming aspects. Note that, on my system, the program was able to add a number consisting of 200 million digits of 9 to the same number in 0.98 seconds, which is pretty impressive, but still, if there is anything in performance that can be enhanced, let me know. That is because the multiplication function is going to be dependent on this function, and things will get slower because multiplication is repeated addition.

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>

char *add(char *x, char *y);

int main(void)
{
    clock_t t1 = clock();
    /* ............................................................................................................................. */
    size_t nd = 200 * 1000 * 1000; // 2 hundred million digits
    char *x = (char *)malloc(nd + 1);
    char *y = (char *)malloc(nd + 1);
    if (!x || !y) {
        fputs("error: memory allocation failed.\n", stderr);
        goto dismantle;
    }
    /* assign x and y */
    for (size_t i = 0; i < nd; i++) {
        x[i] = '9';
        y[i] = '9';
    }
    /* NUL-terminate the arrays */
    x[nd] = '\0';
    y[nd] = '\0';
    /* add x and y */
    char *z = add(x, y);
    if (z) {
        //printf("x + y = %s\n", z);
        free(z);
    }
dismantle:
    if (x) free(x);
    if (y) free(y);
    /* ............................................................................................................................. */
    clock_t t2 = clock();
    fprintf(stderr, "time elapsed: %.4f\n", (double)(t2 - t1) / CLOCKS_PER_SEC);
    fprintf(stderr, "Press any key to continue . . . ");
    getchar();
    return 0;
}

void ReverseArray(char *buf)
{
    char tmp, *ptr = buf + strlen(buf) - 1;
    while (buf < ptr) {
        tmp = *buf;
        *buf++ = *ptr;
        *ptr-- = tmp;
    }
}

char *add(char *x, char *y)
{
    /* store number of digits in size_t variables*/
    size_t xlen = strlen(x);
    size_t ylen = strlen(y);
    /* essential variables */
    char *z = NULL;
    int r = -1;
    int val, var = (xlen >= ylen);
    size_t n = 0;
    /* pointers to x and y */
    char *p1, *p2;
    char *t1, *t2;
    /* assign pointers according to the value of var */
    if (var) {
        p1 = x + xlen - 1;  t1 = x;
        p2 = y + ylen - 1;  t2 = y;
        z = (char *)malloc(xlen + 2);
    }
    else {
        p1 = y + ylen - 1;  t1 = y;
        p2 = x + xlen - 1;  t2 = x;
        z = (char *)malloc(ylen + 2);
    }
    /* check for allocation failure */
    if (!z) {
        fputs("error: memory allocation failed.\n", stderr);
        return NULL;
    }
    /* stage 1 */
    while (p2 >= t2)
    {
        if (r == -1) {
            val = (*p1-- - '0') + (*p2-- - '0');
        }
        else {
            val = (*p1-- - '0') + (*p2-- - '0') + 1;
        }
        if (val > 9) {
            r = val - 10;
            z[n++] = r + '0';
        }
        else {
            z[n++] = val + '0';
            r = -1;
        }
    }
    /* stage 2 */
    while (p1 >= t1)
    {
        if (r == -1) {
            z[n++] = *p1--;
        }
        else {
            val = (*p1-- - '0') + 1;
            if (val > 9) {
                r = val - 10;
                z[n++] = r + '0';
            }
            else {
                z[n++] = val + '0';
                r = -1;
            }
        }
    }
    /* stage 3 */
    if (r != -1) {
        z[n++] = 1 + '0';
    }
    z[n] = '\0';
    ReverseArray(z);
    return z;
}
\$\endgroup\$
  • 10
    \$\begingroup\$ Did you consider to use one of the many (freely available) big number libraries, instead of rolling that stuff on your own? \$\endgroup\$ – πάντα ῥεῖ Jan 22 '17 at 19:03
  • \$\begingroup\$ @πάνταῥεῖ Actually, no. The idea did not cross my mind. \$\endgroup\$ – machine_1 Jan 22 '17 at 20:11
  • 3
    \$\begingroup\$ 200 * 1000 * 1000 --> Careful with using narrow constants to build wide values. (size_t)200 * 1000 * 1000 is more portable. \$\endgroup\$ – chux - Reinstate Monica Feb 1 '17 at 3:31
12
\$\begingroup\$

You are using textual strings to represent numbers. This is not ideal for the processor. It is better to store the numbers as arrays of 32 or 64-bit unsigned integers, depending on the type of CPU. The only issue of course is that you might have to convert from and to the textual representation for your input and output, but if you have to perform many operations on your numbers, this is a small cost to pay.

There are existing libraries that can handle very large numbers (usually called "bignums"). One of them is the GNU Multiple Precision Arithmetic library. You can use that one if you want, or look at its source code to see how they handle large numbers.

\$\endgroup\$
  • \$\begingroup\$ If I were to store numbers as arrays of 32 or 64 bit unsigned integers, then my program will not be able to handle large numbers because that will eat up my memory. \$\endgroup\$ – machine_1 Jan 27 '17 at 13:29
  • 7
    \$\begingroup\$ @machine_1: You will need less memory if you use binary encoding in unsigned integers than if you use decimal characters. \$\endgroup\$ – G. Sliepen Jan 27 '17 at 20:14
  • \$\begingroup\$ +1 for GMP! I have a simple calculator using GMP here \$\endgroup\$ – luser droog Jan 28 '17 at 6:07
5
+50
\$\begingroup\$

You can pack 18 digits into a uint64_t. If you do it right, you will obtain a substantial speed-up for your addition operation:

typedef struct big_integer {
    uint64_t* data;
    size_t data_length;
} big_integer;

static size_t DIGITS_PER_UINT64 = 18;
static uint64_t MODULO = 1000000000000000000L;

big_integer* big_integer_create(const char* number_text)
{
    size_t integer_text_length = strlen(number_text);
    size_t data_length = integer_text_length / DIGITS_PER_UINT64 +
                        (integer_text_length % DIGITS_PER_UINT64 != 0 ? 1 : 0);

    uint64_t* data = calloc(data_length, sizeof(uint64_t));

    size_t chars_loaded = 0;
    size_t data_index = 0;
    uint64_t sum = 0;
    uint64_t factor = 1;

    for (int i = integer_text_length - 1; i >= 0; --i)
    {
        const char character = number_text[i];

        sum += factor * (character - '0');
        factor *= 10;

        if (++chars_loaded % DIGITS_PER_UINT64 == 0)
        {
            data[data_index++] = sum;
            sum = 0;
            factor = 1;
        }
    }

    if (chars_loaded % DIGITS_PER_UINT64 != 0)
    {
        data[data_index] = sum;
    }

    big_integer* result = malloc(sizeof *result);
    result->data = data;
    result->data_length = data_length;
    return result;
}

void big_integer_free(big_integer* big_int)
{
    free(big_int->data);
    free(big_int);
}

void big_integer_print(big_integer* big_int, FILE* file)
{
    int start_index = (int) big_int->data_length - 1;

    if (big_int->data[start_index] == 0)
    {
        --start_index;

        if (start_index >= 0)
        {
            fprintf(file, "%llu", big_int->data[start_index]);
        }

        for (int i = start_index - 1; i >= 0; --i)
        {
            fprintf(file, "%018llu", big_int->data[i]);
        }
    }
    else
    {
        fprintf(file, "%llu", big_int->data[start_index--]);

        for (int i = start_index; i >= 0; --i)
        {
            fprintf(file, "%018llu", big_int->data[i]);
        }
    }
}

big_integer* big_integer_copy(big_integer* big_int)
{
    big_integer* ret = malloc(sizeof *ret);
    uint64_t* data = malloc(big_int->data_length * sizeof(uint64_t));
    ret->data_length = big_int->data_length;
    ret->data = data;
    memcpy(data, big_int->data, big_int->data_length * sizeof(uint64_t));
    return ret;
}

big_integer* big_integer_add(big_integer* a, big_integer* b)
{
    size_t a_length = a->data_length;
    size_t b_length = b->data_length;
    size_t data_length = MAX(a_length, b_length) + 1;
    uint64_t* data = calloc(data_length, sizeof(uint64_t));

    size_t end_index = MIN(a_length, b_length);
    int index = 0;
    uint64_t carry = 0;
    big_integer* result = malloc(sizeof *result);

    while (index != end_index)
    {
        uint64_t tmp = a->data[index] + b->data[index] + carry;

        if (tmp >= MODULO)
        {
            carry = 1;
            data[index] = tmp % MODULO;
        }
        else
        {
            carry = 0;
            data[index] = tmp;
        }

        ++index;
    }

    while (index < a_length)
    {
        uint64_t tmp = a->data[index] + carry;

        if (tmp >= MODULO)
        {
            carry = 1;
            data[index++] = tmp % MODULO;
        }
        else
        {
            carry = 0;
            data[index++] = tmp;
        }
    }

    while (index < b_length)
    {
        uint64_t tmp = b->data[index] + carry;

        if (tmp >= MODULO)
        {
            carry = 1;
            data[index++] = tmp % MODULO;
        }
        else
        {
            carry = 0;
            data[index++] = tmp;
        }
    }

    if (carry > 0) {
        data[data_length - 1] = 1;
    }

    result->data = data;
    result->data_length = data_length;
    return result;
}

(The entire demonstration program here.)

Now, when I run my demo program (see the Gist above) with -O3 optimization flag, I get something like

time elapsed: 1.1515, total time: 1.4847
coderodde time: 0.2157, total time: 1.2360

Hope that helps.

Note I have updated the code once again in the linked gist. Make sure you are up-to-date.

\$\endgroup\$
  • \$\begingroup\$ I'm a bit confused as to how you can store 18 digits in a uint_64. Each number between 0 and 9 requires at least 4 bits to represent, so surely the means 64 bits can represent 16 different numbers, not 18. What am I missing here? \$\endgroup\$ – DemonessJess Jan 30 '17 at 17:34
  • 2
    \$\begingroup\$ @psychedelic_alex I use each uint64_t to hold a huge "digit" from 0 to 999...99. Now, if I add, say, 2 to the maximum "digit", I get 100....1. Then, I take the modulo of 100..0, and arrive to 1 with a carry flag for the next uint64_t. I cannot explain this formally, sorry. :-( \$\endgroup\$ – coderodde Jan 30 '17 at 17:44
  • 1
    \$\begingroup\$ @psychedelic_alex numbers are never encoded as individual digits then appended together like you propose because that is wasteful: those 4 bytes can represent [0-15] but you'd only use [0-9] so you're wasting 40% of expressivity on that digit. With 16 such digits in 64bits you're wasting 99.98% of expressivity! Instead, the numbers are entirely represented in binary base (or rather in octal base [0-9ABCDE]). So a 64bit number can represent 0-18446744073709551616 which means up to 18 digits you're fine. \$\endgroup\$ – MrBrushy Jan 31 '17 at 12:51
  • 2
    \$\begingroup\$ Using uint64_t as the base type is problematic with multiplication and division (depending on how it is chained) as typically intermediate 2x types like int128_t are needed - which may not exist. For that reason, suggest using (u)int32_t as the work-horse type or even uint16_t is that is the processor's native size. \$\endgroup\$ – chux - Reinstate Monica Feb 1 '17 at 3:27
2
\$\begingroup\$

Calculations naturally start from the least significant digits, and if they are stored at the lowest address then the computation can start there without the need for strlen(). Removing strlen() would then save time.

\$\endgroup\$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.