Задача
Последовательность Фибоначчи определяется рекурсивным правилом:
Fn = Fn−1 + Fn−2, где F1 = 1 и F2 = 1.
Таким образом, первые 12 членов последовательности равны:
F1 = 1
F2 = 1
F3 = 2
F4 = 3
F5 = 5
F6 = 8
F7 = 13
F8 = 21
F9 = 34
F10 = 55
F11 = 89
F12 = 144
Двенадцатый член F12 - первый член последовательности, который содержит три цифры.
Каков порядковый номер первого члена последовательности Фибоначчи, содержащего 1000 цифр?
Решение
Сразу скажу, что задачу можно решить в лоб, просто складывая числа, пока длина очередного числа не станет равна 1000. Для этого нужна будет только арифметика больших чисел, которая уже неоднократно применялась в предыдущих задачах. Это будет быстро.
Но сначала мне казалось, что лобовое решение здесь уже не пройдёт, поэтому я обратился к более продвинутым способам вычисления n-го числа Фибоначчи. Их два: один это возведение матриц в степень, а другой это вычисления золотого сечения, где тоже используется возведение в степень, а так как мы всё равно работаем с большими числами, то особой разницы между этими способами нет.
Я выбрал матрицы, потому что я их понимаю немножко больше. Три последних числа Фибоначчи, например, 8, 5, 3, можно представить в виде матрицы 2x2:
8 5
5 3
Теперь если эту матрицу умножить на матрицy
1 1
1 0
То правила умножения матриц приведут к тому, что 8 банально сложится с 5, 5 с 3, а 3 станет 5:
13 8
8 5
Т.е. все числа, которые были в матрице, получили "апгрейд" в виде следующего числа Фибоначчи.
Почему это работает, я не могу чётко понять, просто вот такие правила, если им следовать, то будет получаться именно так, и мне придётся пока с этим жить.
Матричные степени
Однако умножать матрицу на [1,1 | 1,0], чтобы получить следующее число, бесперспективно. С таким же успехом можно просто вычислять следующее число без всяких матриц.
Что делает матрицы полезными, это возведение их в степень. Для этого берём матрицу с допустим 6-м номером Фибоначчи 8: [8,5 | 5,3], умножаем на саму себя и получаем первое число новой матрицы (8*8+5*5) = 89, это уже 11-й номер Фибоначчи. Т.е. мы получили не номер n+1, а номер n*2-1. Если полученную матрицу умножить на саму себя ещё раз, мы получим n*4, далее n*8, n*16 и т.д.
Это сильно ускоряет процесс и мы можем достичь нужного номера за время O(log2 n). Как это работает, то есть почему степень матриц заменяет суммирование, я понимаю ещё меньше. Но что делать.
Стратегия решения
Начиная с первой матрицы [1,1 | 1,0], умножать матрицу на саму себя, пока не получим в её первом числе переполнение более 1000 символов.
После этого откатимся на предыдущую матрицу и будем просто домножать её на [1,1 | 1,0], получая по одному следующему числу, пока не дойдём до длины ровно 1000.
В этом месте также возможна оптимизация. Так, если матрицу нельзя умножить на саму себя из-за переполнения, её можно было бы умножить на какую-то более младшую степень самой себя. Но для этого нужно хранить промежуточные степени, что мне делать было уже лень.
Программа получилась большой.
Арифметика больших чисел
Базируется на предыдущих выпусках, но здесь она доработана, чтобы работать с операндами произвольной длины, и также добавлена проверка на переполнение. Кроме того, добавлены сервисные функции для обнуления, присвоения значения и копирования.
Так как промежуточные результаты вычисления держать негде, использование больших чисел становится очень похожим на ассемблер. Операцию типа
e = a * b + с * d
нужно разбивать на отдельные элементарные операции. Это выглядит так:
e = a
e *= b
tmp = c
tmp *= d
e += tmp
Память
Мне нужны две матрицы, в каждой по 4 числа длиной 1000 знаков, и ещё одно 1000-значное число для промежуточных операций. Итого нужно 9 * 1000 байт памяти.
#define LENGTH 1000
static char num_buffer[LENGTH * 9];
Создаю матрицы, состоящие из 4-х структур типа Number, у каждой структуры свой указатель на буфер:
Перемножение матриц:
Здесь я немного смухлевал. Так как в Фибоначчи-матрицах два числа по диагонали одинаковые, я вычисляю только одно из них, а второе просто копирую. А ещё лучше отказаться от матриц как таковых и хранить тупо три числа. Математические правила от этого не изменятся, просто нужно знать, что на что умножать и с чем складывать. Но эта мысль пришла ко мне слишком опосля :)
Также я смухлевал с продвижением матрицы на одно число дальше. Здесь не надо ничего умножать, а просто два раза сложить и два раза скопировать:
Ну и главный код:
В этом месте инициализируется первая матрица [1,1 | 1,0], счётчик степеней exp, и два указателя dst_ptr и src_ptr. Они будут использоваться, чтобы менять местами матрицы при вычислениях. То есть сначала одна матрица умножается на себя и результат попадает в другую матрицу, затем другая матрица умножается на себя и результат попадает в первую матрицу. Так они постоянно меняются друг с другом. А точнее, меняются указатели на них.
Как только умножение матриц дало сбой (возникло переполнение), в указателе dst_ptr будет недоумноженный результат, а в src_ptr – последняя вычисленная матрица. Эту матрицу добиваем до нужной длины уже простыми итерациями mat_next().
Ну и как я говорил, можно убрать перемножение матриц и оставить просто итерации, по скорости я не заметил никакой разницы :)
Ссылка на онлайн-компилятор языка C с текстом программы
Подбока всех задач: