Как распараллелить последовательность Фибоначчи до 10 ^ 5 терма

Question

Мне нужно распараллелить (с помощью openmp) последовательность Фибоначчи из этого последовательного кода, чтобы вычислить 10 ⁵ -й член последовательности, но я застрял на 3 недели без какой-либо хорошей идеи, кто-тоесть идеи или советы о том, как это сделать?

Вот последовательный код в C:

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

#define MAX 100010
#define LEN 25001

char seq[MAX][LEN];

void add(int a, int b) {
    int i, aux, s;

    for (i = 0, aux = 0; seq[a][i] != '\0' && seq[b][i] != '\0'; i++) {
        s = seq[a][i] + seq[b][i] + aux - '0' - '0';
        aux = s / 10;
        seq[a + 1][i] = s % 10 + '0';
    }

    while (seq[a][i] != '\0') {
        s = seq[a][i] + aux - '0';
        aux = s / 10;
        seq[a + 1][i] = s % 10 + '0';
        i++;
    }

    while (seq[b][i] != '\0') {
        s = seq[b][i] + aux - '0';
        aux = s / 10;
        seq[a + 1][i] = s % 10 + '0';
        i++;
    }

    if (aux != 0)
        seq[a + 1][i++] = aux + '0';

    seq[a + 1][i] = '\0';
}

int main() {
    int n, i, len;

    seq[0][0] = '0';
    seq[0][1] = '\0';
    seq[1][0] = '1';
    seq[1][1] = '\0';

    for (i = 2; i < MAX; i++)
        add(i - 1, i - 2);

    scanf("%d", &n);

    len = strlen(seq[n]);
    for (i = 0; i <= len - 1; i++)
        printf("%c", seq[n][len - 1 - i]);
    printf("\n");
    fflush(stdout);

    return 0;
}

Answer 1

Вместо того, чтобы пытаться распараллелить сложение bignum, что сложно, вы можете попробовать и вычислить несколько терминов параллельно:

F(n+1) = F(n) + F(n-1)
F(n+2) = F(n+1) + F(n) = 2*F(n) + F(n-1)
F(n+3) = F(n+2) + F(n+1) = 2*F(n+1) + F(n) = 2*F(n) + 2*F(n-1) + F(n) = 3*F(n) + 2*F(n-1)
...

Обратите также внимание, что вы должны вычислять блоки цифр за раз: 8или 9 базовых-10 цифр могут быть вычислены с использованием 32-битных элементов массива.

Вот модифицированная версия с несколькими улучшениями:

• она вычисляет блоки по 8 цифр за раз
• он может принимать аргументы командной строки
• он использует намного меньше памяти
• он может обрабатывать гораздо большие значения
• он намного эффективнее (20x)

Вы сможете легко распараллелить его.

/* Parallelisable bignum Fibonacci computation by chqrlie */
#include <stdio.h>
#include <stdint.h>
#include <inttypes.h>
#include <stdlib.h>
#include <string.h>

#if 0

/* 2500ms for fib(100000) */
#define DIGIT  10
#define NDIGIT 1
#define FMT "d"
typedef unsigned char digit_t;

#elif 1

/* 279ms for fib(100000) */
#define DIGIT  100000000
#define NDIGIT 8
#define FMT PRIu32
typedef uint32_t digit_t;

#else

/* 720ms for fib(100000) */
#define DIGIT  1000000000000000000
#define NDIGIT 18
#define FMT PRIu64
typedef uint64_t digit_t;

#endif

int add1(digit_t *c, const digit_t *a, int alen, const digit_t *b, int blen) {
    digit_t aux, s;
    int i;

    /* assuming alen >= blen */
    for (i = 0, aux = 0; i < alen; i++) {
        s = a[i] + b[i] + aux;
        aux = s / DIGIT;
        c[i] = s % DIGIT;
    }
    if (aux != 0) {
        c[i++] = (digit_t)aux;
    }
    return i;
}

int add2(digit_t *c, const digit_t *a, int alen, const digit_t *b, int blen) {
    digit_t aux, s;
    int i;

    /* assuming alen >= blen */
    for (i = 0, aux = 0; i < alen; i++) {
        s = 2 * a[i] + b[i] + aux;
        aux = s / DIGIT;
        c[i] = s % DIGIT;
    }
    if (aux != 0) {
        c[i++] = (digit_t)aux;
    }
    return i;
}

int add3(digit_t *c, const digit_t *a, int alen, const digit_t *b, int blen) {
    digit_t aux, s;
    int i;

    /* assuming alen >= blen */
    for (i = 0, aux = 0; i < alen; i++) {
        s = 3 * a[i] + 2 * b[i] + aux;
        aux = s / DIGIT;
        c[i] = s % DIGIT;
    }
    if (aux != 0) {
        c[i++] = (digit_t)aux;
    }
    return i;
}

int add4(digit_t *c, const digit_t *a, int alen, const digit_t *b, int blen) {
    digit_t aux, s;
    int i;

    /* assuming alen >= blen */
    for (i = 0, aux = 0; i < alen; i++) {
        s = 5 * a[i] + 3 * b[i] + aux;
        aux = s / DIGIT;
        c[i] = s % DIGIT;
    }
    if (aux != 0) {
        c[i++] = (digit_t)aux;
    }
    return i;
}

void printseq(const digit_t *s, int len) {
    printf("%"FMT, s[len - 1]);
    for (int i = 1; i < len; i++)
        printf("%.*"FMT, NDIGIT, s[len - 1 - i]);
    printf("\n");
}

int main(int argc, char *argv[]) {
    int MIN, i, LEN, MAX;

    if (argc > 1) {
        MAX = MIN = strtol(argv[1], NULL, 0);
        if (argc > 2)
            MAX = strtol(argv[2], NULL, 0);
    } else {
        scanf("%d", &MIN);
        MAX = MIN;
    }

    /* length if fib(n) is less than n*log10(phi)+2 */
    LEN = (MAX * 20910ULL) / 100000 / NDIGIT + 2;
    /* allocate 8 bignums */
    int *slen = calloc(sizeof(*slen), 8);
    digit_t (*seq)[LEN] = calloc(sizeof(*seq), 8);

    if (slen == NULL || seq == NULL) {
        fprintf(stderr, "memory allocation error\n");
        return 1;
    }

    seq[0][0] = 0;
    slen[0] = 1;
    if (0 >= MIN) printseq(seq[0], slen[0]);
    seq[1][0] = 1;
    slen[1] = 1;
    if (1 >= MIN) printseq(seq[1], slen[1]);

    for (i = 2; i <= MAX && (MAX + 1 - i) % 4 != 0; i++) {
        slen[i] = add1(seq[i], seq[i - 1], slen[i - 1], seq[i - 2], slen[i - 2]);
        if (i >= MIN) printseq(seq[i], slen[i]);
    }
    for (; i <= MAX; i += 4) {
        int im2 = (i - 2) & 7;
        int im1 = (i - 1) & 7;
        int i0 = (i + 0) & 7;
        int i1 = (i + 1) & 7;
        int i2 = (i + 2) & 7;
        int i3 = (i + 3) & 7;
        /* the next 4 calls can be parallelised */
        slen[i0] = add1(seq[i0], seq[im1], slen[im1], seq[im2], slen[im2]);
        slen[i1] = add2(seq[i1], seq[im1], slen[im1], seq[im2], slen[im2]);
        slen[i2] = add3(seq[i2], seq[im1], slen[im1], seq[im2], slen[im2]);
        slen[i3] = add4(seq[i3], seq[im1], slen[im1], seq[im2], slen[im2]);
        /* the print calls must be called sequentially */
        if (i + 0 >= MIN) printseq(seq[i0], slen[i0]);
        if (i + 1 >= MIN) printseq(seq[i1], slen[i1]);
        if (i + 2 >= MIN) printseq(seq[i2], slen[i2]);
        if (i + 3 >= MIN) printseq(seq[i3], slen[i3]);
    }
    free(slen);
    free(seq);
    return 0;
}

Answer 2

Как вы уже знаете, F (100000) - астрономически огромное число.И чтобы вычислить это значение, вы должны сложить два других действительно огромных числа F (99999) и F (99998).

Вот мой совет:

У вас есть два числа, которые много(тысячи) цифр длинных и N процессоров.Вы можете разделить сложение на несколько потоков.Например:

F(167) == 35600075545958458963222876581316753
F(168) == 57602132235424755886206198685365216

Чтобы вычислить F (169), вы должны добавить оба этих числа выше.Но давайте рассматривать это как 4 отдельных сложения по 9 цифр в каждом.

          (A)         (B)         (C)        (D)
F(167) == 035600075   545958458   963222876  581316753
F(168) == 057602132   235424755   886206198  685365216
======================================================= + 
           93202207   781383213  1849429074 1266681969
                                 ^          ^

Итак, теперь у нас есть 4 суммирования.Два из них, C и D, имеют операцию переноса.Так что нам просто нужно скорректировать результат слева от каждого на + 1

           93202207   781383213   849429074  266681969
                             +1          +1
======================================================= + 
F(169) ==  93202207   781383214   849429075 266681969

Итак, ваш алгоритм для вычисления Fib (100000) выглядит примерно так:Где BigNumber - это структура, представляющая ваши цифры.Вы используете массив символов, что также приемлемо.

void fib(int stop)
{
    BigNumber  f0 = 0;
    BigNumber  f1 = 1;
    BigNumber  f2 = 1;

    BigNumber* pF0 = &f0;
    BigNumber* pF1 = &f1;
    BigNumber* pF2 = &f2;

    for (int i = 2; i <= stop; i++)
    {
        ParallelAdd(pF0, pF1, pF2); // *pF2 = *pF1 + *pF0

        // shift via pointers
        pF0 = pF1;
        pF1 = pF2;
        pF2 = pF0;
    }
    Print(f2);
}

main()
{
    fib(100000);
}

Ваш ParallelAdd разбит числа, переданные через указатель, f1 и f0, на N групп по K цифр в каждой, где N - этоколичество процессоров, которые вы имеете в наличии или хотите использовать.Затем каждый процессор вычисляет сложение каждого из них, используя уже имеющийся код.После того, как эти N операций завершены, отсканируйте результирующий набор, чтобы увидеть, какие дополнения содержали результат длиной K + 1 цифра, а затем скорректируйте с помощью логики +1, как описано выше.Объедините в одну строку и вставьте обратно в адрес указателя, на который ссылается pF2.