Алгоритм умножения целых чисел произвольной точности (bignum)

Я пишу небольшую бигнум-библиотеку для домашнего задания. Я должен реализовать умножение Карацубы, но перед этим я хотел бы написать наивную процедуру умножения.

Я следую руководству, написанному Полом Циммерманом под названием «Современная компьютерная арифметика», которое бесплатно доступен в Интернете.

На странице 4 есть описание алгоритма BasecaseMultiply, который выполняет умножение в начальной школе.

Я понимаю шаги 2, 3, где B^j - сдвиг цифры 1, j раз. Но я не понимаю шагов 1 и 3, где у нас A*b_j. Как должно выполняться это умножение, если умножение bignum еще не определено?

Будет ли операция «*» в этом алгоритме просто методом повторного сложения?

Вот части, которые я написал до сих пор. Я провел их модульное тестирование, поэтому по большей части они кажутся правильными:

Структура, которую я использую для своего бигнума, выглядит следующим образом:

#define BIGNUM_DIGITS 2048
typedef uint32_t u_hw; // halfword
typedef uint64_t u_w; // word

typedef struct {
    unsigned int sign; // 0 or 1
    unsigned int n_digits;
    u_hw digits[BIGNUM_DIGITS];
} bn;

Доступные на данный момент подпрограммы:

bn *bn_add(bn *a, bn *b); // returns a+b as a newly allocated bn
void bn_lshift(bn *b, int d); // shifts d digits to the left, retains sign
int bn_cmp(bn *a, bn *b); // returns 1 if a>b, 0 if a=b, -1 if a<b

c algorithm math multiplication bignum
person Community    schedule 02.05.2010    source источник
comment
Вы можете попробовать сделать процедуру умножения больших чисел рекурсивной, поскольку вам необходимо выполнить A * b_j. Я не уверен, как это будет определено, так что вряд ли это ответ, но это, по крайней мере, идея. Хотя я уверен, что есть лучшее решение. :)   -  person Dustin    schedule 03.05.2010

Ответы (2)


arrow_upward
2
arrow_downward

Я написал алгоритм умножения некоторое время назад, и у меня есть этот комментарий вверху. Если у вас есть два числа x и y одинакового размера (одинаковые n_digits), вы должны умножить их таким образом, чтобы получить n, в котором будет вдвое больше цифр. Часть сложности алгоритма связана с определением того, какие биты не следует умножать, если n_digits не одинаковы для обоих входов.

Начиная справа, n0 равно x0*y0, и вы избавляетесь от переполнения. Теперь n1 ​​представляет собой сумму x1*y0 и y1*x0, а предыдущее переполнение сдвинуто на размер вашей цифры. Если вы используете 32-битные цифры в 64-битной математике, это означает, что n0 = low32(x0*y0) и вы несете high32(x0*y0) в качестве переполнения. Вы можете видеть, что если бы вы фактически использовали 32-битные цифры, вы не могли бы добавить центральные столбцы, не превысив 64 бит, поэтому вы, вероятно, используете 30 или 31-битные цифры.

Если у вас есть 30 бит на цифру, это означает, что вы можете умножить два 8-значных числа вместе. Сначала напишите этот алгоритм, чтобы он принимал два небольших буфера с n_digits до 8 и использовал собственную математику для арифметики. Затем реализуйте его снова, взяв n_digits произвольного размера и используя первую версию вместе с вашим методом сдвига и добавления, чтобы умножить фрагменты цифр 8x8 за раз.

/*
    X*Y = N

                          x0     y3
                            \   /  
                             \ /   
                              X    
                      x1     /|\     y2
                        \   / | \   /  
                         \ /  |  \ /   
                          X   |   X    
                  x2     /|\  |  /|\     y1
                    \   / | \ | / | \   /  
                     \ /  |  \|/  |  \ /   
                      X   |   X   |   X    
              x3     /|\  |  /|\  |  /|\     y0
                \   / | \ | / | \ | / | \   /
                 \ /  |  \|/  |  \|/  |  \ /
                  V   |   X   |   X   |   V
                  |\  |  /|\  |  /|\  |  /|
                  | \ | / | \ | / | \ | / |
                  |  \|/  |  \|/  |  \|/  |
                  |   V   |   X   |   V   |
                  |   |\  |  /|\  |  /|   |
                  |   | \ | / | \ | / |   |
                  |   |  \|/  |  \|/  |   |
                  |   |   V   |   V   |   |
                  |   |   |\  |  /|   |   |
                  |   |   | \ | / |   |   |
                  |   |   |  \|/  |   |   |
                  |   |   |   V   |   |   |
                  |   |   |   |   |   |   |
              n7  n6  n5  n4  n3  n2  n1  n0
*/
person Community    schedule 02.05.2010
comment
Я не уверен, что вы подразумеваете под превышением 64 бит с 32-битными цифрами? Как видите, у меня есть 32-битные цифры, и я планировал использовать 64-битный размер слова для управления полным диапазоном. - person snap; 03.05.2010
comment
Кстати, крутая диаграмма, очень удобно с умножением. - person snap; 03.05.2010
comment
Базовый вариант умножения позволяет избежать добавления столбцов, проходя по одной диагонали за раз и создавая результат по мере его продвижения, поэтому вы можете игнорировать эту проблему. Если вы сложите произведения двух 32-битных умножений, вы получите 65-битное значение. Вы либо должны справиться с этим переполнением (легко в ассемблере, хлопотно в C), либо вам нужно использовать меньше битов в каждой цифре. - person drawnonward; 03.05.2010
comment
Между прочим, дешевле выполнить c+=ab, чем просто выполнить c=ab, потому что в последнем случае вы должны инициализировать c нулем. - person drawnonward; 03.05.2010
comment
Вы не имеете в виду только обратное? т.е. в первом случае вы должны инициализировать c нулем, так как в более позднем ('=') вы просто полностью переопределяете c. - person Matthieu M.; 03.05.2010
comment
В этом случае c — это буфер размером 8 КБ. Алгоритм умножения добавляет промежуточные результаты к произвольным цифрам по мере его выполнения, поэтому все цифры числа c должны быть инициализированы. На практике c=ab вычисляется как c=0; с+=аб; поэтому, если c уже имеет значение, добавление бесплатно. - person drawnonward; 03.05.2010

arrow_upward
1
arrow_downward

Чтобы выполнить A*b_j, вам нужно выполнить школьное умножение бигнума на однозначную цифру. В конечном итоге вам придется сложить вместе несколько двузначных продуктов:

bn *R = ZERO;
for(int i = 0; i < n; i++) {
  bn S = {0, 2};
  S.digits[0] = a[i] * b_j;
  S.digits[1] = (((u_w)a[i]) * b_j) >> 32;  // order depends on endianness
  bn_lshift(S, i);
  R = bn_add(R, S);
}

Конечно, это очень неэффективно.

person Community    schedule 02.05.2010