Fibonacci Ateşi#12483

Fibonacci Ateşi (Editoryal)

değeri 10^{18}'e kadar çıkabildiğinden Fibonacci sayılarını tek tek dolaşmak mümkün değildir. Matris yöntemi de uygun değildir: -ıncı kuvvetlerin toplamı için yaklaşık boyutunda bir matris gerekirdi, ise 10^5'e kadar çıkabilir. Bu nedenle Fibonacci sayıları için açık bir formül — Binet formülü kullanırız:

f_i = \frac{\varphi^i - \psi^i}{\sqrt{5}},

burada

\varphi = \frac{1 + \sqrt{5}}{2}, \quad \psi = \frac{1 - \sqrt{5}}{2}

ve sayıları x^2 = x + 1 denkleminin kökleridir ve Vieta formüllerine göre \varphi\psi = -1'dir.

Kuvvetin açılımı. Binet formülünün -ıncı kuvvetini alalım ve Newton binom açılımını uygulayalım:

f_i^k = \frac{(\varphi^i - \psi^i)^k}{(\sqrt{5})^k} = \frac{1}{(\sqrt{5})^k}\sum_{j=0}^{k}(-1)^j C_k^j\, \varphi^{i(k-j)}\psi^{ij} = \frac{1}{(\sqrt{5})^k}\sum_{j=0}^{k}(-1)^j C_k^j\, \left(\varphi^{k-j}\psi^j\right)^i

Şu gösterimi kullanalım:

t_j = \varphi^{k-j}\psi^j, \quad j = 0, 1, \ldots, k

üzerinden 'den 'ye kadar toplayalım ve toplama sırasını değiştirelim:

S = \sum_{i=1}^{n} f_i^k = \frac{1}{(\sqrt{5})^k}\sum_{j=0}^{k}(-1)^j C_k^j \sum_{i=1}^{n} t_j^{\,i}

Böylece terimli toplam, k + 1 terimli bir toplama indirgenmiş olur.

Geometrik dizi. İçteki toplam, ortak çarpanı t = t_j olan bir geometrik dizinin toplamıdır:

G(t) = \sum_{i=1}^{n} t^i = \begin{cases} \dfrac{t\,(t^n - 1)}{t - 1}, & t \neq 1, \\[2mm] n, & t = 1 \end{cases}

t = 1 durumu gerçekten ortaya çıkar ve ayrıca ele alınmalıdır (aksi halde sıfıra bölme oluşur). Örneğin sayısı 'e bölünüyorsa, j = k / 2 için

t_{k/2} = (\varphi\psi)^{k/2} = (-1)^{k/2} = 1

elde ederiz.

t^n kuvveti ikili (hızlı) üs alma yöntemiyle O(\log n) sürede hesaplanır; bu yüzden her terim n = 10^{18} için bile hızlıca bulunur.

Ortak çarpanların ardışık hesaplanması. Her için ve 'nin kuvvetini yeniden hesaplamaya gerek yoktur. Ardışık ortak çarpanlar şu bağıntıyla ilişkilidir:

t_{j+1} = t_j \cdot \frac{\psi}{\varphi},

dolayısıyla t_0 = \varphi^k değerini bir kez hesaplamak ve ardından her iterasyonda mevcut değeri ile çarpmak yeterlidir.

Modüler hesaplamalar. Kayan noktalı aritmetik burada uygun değildir: p = 10^9 + 7 modunda kesin cevap gereklidir. a + b\sqrt{5} biçimindeki sayılarla çalışacağız; burada ve , moduna göre kalanlardır. Bu tür sayıların toplanması ve çarpılması şu kurallara göre yapılır:

(a + b\sqrt{5}) + (c + d\sqrt{5}) = (a + c) + (b + d)\sqrt{5}, \\ (a + b\sqrt{5})(c + d\sqrt{5}) = (ac + 5bd) + (ad + bc)\sqrt{5}

a + b\sqrt{5} biçimindeki sayılar \mathbb{Z}_p[\sqrt{5}] cismini oluşturur. Bu cisimde sıfırdan farklı her elemanın tersi vardır:

(a + b\sqrt{5})^{-1} = \frac{a - b\sqrt{5}}{a^2 - 5b^2},

ayrıca a^2 - 5b^2 paydası moduna göre sıfırdan farklıdır. Bu, çözümdeki tüm bölmelerin geçerli olmasını garanti eder. ve formüllerindeki 'ye bölme, 2^{-1} \bmod p = (p + 1) / 2 ile çarpma ile değiştirilir.

2'nin tersi. moduna göre 'nin tersi, 2x \equiv 1 \pmod{p} koşulunu sağlayan sayısıdır.

x = (p + 1) / 2 sayısının bu koşulu sağladığını kontrol edelim:

Demek ki (p + 1) / 2 tam olarak 2^{-1} \bmod p'dir.

Aranan toplam bir tam sayıdır; bu nedenle son a + b\sqrt{5} değerinde katsayısı mutlaka sıfır olur ve cevap 'dır.

Örnek

İlk örneği ele alalım: , \boldsymbol{k = 1}. k = 1 için formül şu hâli alır:

S = \frac{1}{\sqrt{5}}\left(\sum_{i=1}^{4}\varphi^i - \sum_{i=1}^{4}\psi^i\right)

\varphi^2 = \varphi + 1 bağıntısını kullanarak 'nin kuvvetlerini hesaplayalım:

\varphi = \frac{1 + \sqrt{5}}{2}, \quad \varphi^2 = \frac{3 + \sqrt{5}}{2}, \quad \varphi^3 = \frac{4 + 2\sqrt{5}}{2}, \quad \varphi^4 = \frac{7 + 3\sqrt{5}}{2}

Bunların toplamı \dfrac{15 + 7\sqrt{5}}{2}'dir. 'nin kuvvetleri \sqrt{5} yerine -\sqrt{5} yazılarak elde edilir, bu yüzden onların toplamı \dfrac{15 - 7\sqrt{5}}{2}'dir. Dolayısıyla

S = \frac{1}{\sqrt{5}}\left(\frac{15 + 7\sqrt{5}}{2} - \frac{15 - 7\sqrt{5}}{2}\right) = \frac{7\sqrt{5}}{\sqrt{5}} = 7,

bu da cevapla örtüşür: 1 + 1 + 2 + 3 = 7.

Bir örnek ele alalım: \boldsymbol{k = 3}. Farkın küpünü açalım:

f_i^3 = \frac{(\varphi^i - \psi^i)^3}{(\sqrt{5})^3} = \frac{\varphi^{3i} - 3\varphi^{2i}\psi^i + 3\varphi^i\psi^{2i} - \psi^{3i}}{5\sqrt{5}}

Dizilerin ortak çarpanları şunlardır (\varphi\psi = -1 eşitliğini kullanıyoruz):

t_0 = \varphi^3, \quad t_1 = \varphi^2\psi = -\varphi, \quad t_2 = \varphi\psi^2 = -\psi, \quad t_3 = \psi^3

O halde

\sum_{i=1}^{n} f_i^3 = \frac{1}{5\sqrt{5}}\left(\sum_{i=1}^{n}(\varphi^3)^i - 3\sum_{i=1}^{n}(-\varphi)^i + 3\sum_{i=1}^{n}(-\psi)^i - \sum_{i=1}^{n}(\psi^3)^i\right)

Dört toplamın her biri geometrik dizi toplam formülüyle hesaplanır. Bunu n = 3 için kontrol edelim. 'nin kuvvetleri şunlardır:

\varphi^3 = 2 + \sqrt{5}, \quad \varphi^6 = 9 + 4\sqrt{5}, \quad \varphi^9 = 38 + 17\sqrt{5}

O halde

\sum_{i=1}^{3}(\varphi^3)^i = \varphi^3 + \varphi^6 + \varphi^9 = (2 + \sqrt{5}) + (9 + 4\sqrt{5}) + (38 + 17\sqrt{5}) = 49 + 22\sqrt{5}, \\ \sum_{i=1}^{3}(-\varphi)^i = -\varphi + \varphi^2 - \varphi^3 = \frac{-(1 + \sqrt{5}) + (3 + \sqrt{5}) - (4 + 2\sqrt{5})}{2} = -1 - \sqrt{5}

için toplamlar \sqrt{5} yerine -\sqrt{5} yazılarak elde edilir:

\sum_{i=1}^{3}(-\psi)^i = -1 + \sqrt{5}, \quad \sum_{i=1}^{3}(\psi^3)^i = 49 - 22\sqrt{5}

Bulunan değerleri yerine koyalım:

\sum_{i=1}^{3} f_i^3 = \frac{(49 + 22\sqrt{5}) - 3(-1 - \sqrt{5}) + 3(-1 + \sqrt{5}) - (49 - 22\sqrt{5})}{5\sqrt{5}} = \frac{50\sqrt{5}}{5\sqrt{5}} = 10

Sonuç doğrudan hesaplamayla örtüşür:

f_1^3 + f_2^3 + f_3^3 = 1^3 + 1^3 + 2^3 = 1 + 1 + 8 = 10

\boldsymbol{k = 4} için bir örnek ele alalım (özel durum). Burada

t_0 = \varphi^4, \quad t_1 = \varphi^3\psi = -\varphi^2, \quad t_2 = \varphi^2\psi^2 = 1, \quad t_3 = \varphi\psi^3 = -\psi^2, \quad t_4 = \psi^4

ve

\sum_{i=1}^{n} f_i^4 = \frac{1}{25}\Biggl(\sum_{i=1}^{n}(\varphi^4)^i - 4\sum_{i=1}^{n}(-\varphi^2)^i + 6n - 4\sum_{i=1}^{n}(-\psi^2)^i + \sum_{i=1}^{n}(\psi^4)^i\Biggr)

Ortadaki terim t_2 = 1 ortak çarpanına karşılık gelir: geometrik dizi formülü bu terime uygulanamaz ve 1 + 1 + \ldots + 1 toplamı 'ye eşittir. Programda ayrıca ele alınan durum tam olarak budur. Örneğin n = 3 için şunu elde ederiz:

f_1^4 + f_2^4 + f_3^4 = 1^4 + 1^4 + 2^4 = 1 + 1 + 16 = 18

Algoritmanın gerçekleştirilmesi

Modül için sabitini tanımlayalım.

#define MOD 1000000007

Binom katsayılarının değerlerini dizisinde saklayacağız.

vector<long long> Cnk;

Z_{MOD} [\sqrt{5}] cisminde aritmetiği gerçekleştirmek için yapısını tanımlayalım. Bu yapı x + y \sqrt{5} biçimindeki sayıları temsil eder.

struct Num
{
  long long x, y;
  Num(long long x = 0, long long y = 0) : x(x% MOD), y(y% MOD) {}
};

add fonksiyonu Z_{MOD} [\sqrt{5}] cisminde toplama işlemini gerçekleştirir.

Num add(const Num& a, const Num& b)
{
  return Num((a.x + b.x) % MOD, (a.y + b.y) % MOD);
}

mult fonksiyonu Z_{MOD} [\sqrt{5}] cisminde çarpma işlemini gerçekleştirir.

Num mult(const Num& a, const Num& b)
{
  return Num(
    (a.x * b.x + 5 * a.y % MOD * b.y) % MOD,
    (a.x * b.y + a.y * b.x) % MOD
  );
}

pw fonksiyonu Z_{MOD} [\sqrt{5}] cisminde ikili üs almayı gerçekleştirir; burada üssü bir tam sayıdır.

Num pw(Num base, long long exp)
{
  Num res(1, 0);
  while (exp)
  {
    if (exp & 1) res = mult(res, base);
    base = mult(base, base);
    exp >>= 1;
  }
  return res;
}

pw fonksiyonu modunda a^b değerini hızlı ikili üs alma ile hesaplar.

long long pw(long long a, long long b)
{
  long long res = 1;
  a %= MOD;
  while (b)
  {
    if (b & 1) res = (res * a) % MOD;
    a = a * a % MOD;
    b >>= 1;
  }
  return res;
}

inv fonksiyonu Z_{MOD} [\sqrt{5}] cisminde elemanının çarpımsal tersini hesaplar.

Eğer

a = x + y \sqrt{5},

ise, tersi şu şekilde verilir:

a^{-1} = \frac{1}{x+y\sqrt{5}} = \frac{x-y\sqrt{5}}{x^2-5y^2}
Num inv(const Num& a)
{
  long long den = (a.x * a.x - 5 * a.y * a.y) % MOD;
  if (den < 0) den += MOD;

  long long invDen = pw(den, MOD - 2);

  return Num(
    a.x * invDen % MOD,
    (MOD - a.y) * invDen % MOD
  );
}

divide fonksiyonu Z_{MOD} [\sqrt{5}] cisminde bölme işlemini gerçekleştirir:

\frac{a}{b} = a \cdot b^{-1}
Num divide(const Num& a, const Num& b)
{
  return mult(a, inv(b));
}

Programın ana kısmı. Giriş verilerini okuyalım.

scanf("%lld %lld", &n, &k);

C_k^0, C_k^1, C_k^2, ..., C_k^k binom katsayılarını hesaplayıp dizisine kaydedelim. Şunu belirtelim:

Cnk[i] = C_k^i = C_k^{i-1} \cdot \frac{k - i + 1}{i}
Cnk.resize(k + 1);
Cnk[0] = 1;
for (i = 1; i <= k; i++)
{
  Cnk[i] = Cnk[i - 1] * (k - i + 1) % MOD;
  Cnk[i] = Cnk[i] * pw(i, MOD - 2) % MOD;
}

inv2 = 2^{-1} \mod MOD = (MOD + 1) / 2 değerini hesaplayalım.

long long inv2 = (MOD + 1) / 2;

Sabitleri başlatalım: phi = \phi = \frac{1 + \sqrt{5}}{2} ve psi = \psi = \frac{1 - \sqrt{5}}{2}.

Num phi(inv2, inv2);
Num psi(inv2, MOD - inv2);

ratio = \psi / \phi değerini hesaplayalım.

Num ratio = divide(psi, phi);

Z_{MOD} [\sqrt{5}] cisminde ve sabitlerini başlatalım.

Num one(1, 0);
Num minus_one(MOD - 1, 0);

Cevap değişkeninde biriktirilecektir.

Num res(0, 0);

t = \phi^k değerini başlatalım. Ardından döngüde bu değeri ile çarparak şu değerleri elde edeceğiz:

\phi^k, \phi^{k-1}\psi, \phi^{k-2}\psi^2, ..., \phi\psi^{k-1}, \psi^k
Num t = pw(phi, k);

değişkeninde şu toplamı hesaplıyoruz:

\sum_{i=1}^{n} f_i^k = \frac{1}{(\sqrt{5})^k} \sum_{j=0}^{k} (-1)^j C_k^j \sum_{i=1}^{n} \left(\phi^{k-j}\psi^j\right)^i

Aranan toplam k + 1 terim içerir. Bu nedenle k + 1 iterasyon yaparız.

for (j = 0; j <= k; j++)
{
  Num term;

Bu noktada t = \phi^{k-j} \psi^j'dir. Geometrik dizinin toplamını hesaplayalım:

term = \sum_{i=1}^{n} t^i = \frac{t(t^n-1)}{t-1}

Ancak t = 1 durumu ayrıca ele alınmalıdır, çünkü geometrik dizi toplam formülünde 'a bölme oluşur. t = 1 olduğunda toplam 1 + 1 + 1 + ... + 1 ( kez) biçimini alır ve 'ye eşittir.

  if (t.x == 1 && t.y == 0)
    term = Num(n % MOD, 0);
  else
  {
    Num tn = pw(t, n);
    Num num = add(tn, minus_one);     // t^n - 1
    term = mult(t, num);              // t * (t^n - 1)
    Num den = add(t, minus_one);      // t - 1
    term = divide(term, den);         // t * (t^n - 1) / (t - 1)
  }

term = (-1)^j C_k^j \cdot term değerini hesaplayalım.

  long long coef = Cnk[j];
  if (j & 1) coef = MOD - coef;

  term.x = term.x * coef % MOD;
  term.y = term.y * coef % MOD;

Sıradaki terimini sonucuna ekleyelim.

  res = add(res, term);

'yi ile çarparak t = \phi^{k-j-1} \psi^{j + 1} elde edelim.

  t = mult(t, ratio);
}

Geriye 'i (\sqrt{5})^k'ya bölmek kalır.

Num sqrt5(0, 1);
Num sqrt5k = pw(sqrt5, k);
res = divide(res, sqrt5k);

Cevabı yazdıralım.

printf("%lld\n", res.x);