HDU 6395 – Sequence (分块+快速幂)

题意

计算$F_n \mod 10^9+7$,定义如下:$$\begin{cases} F_1 &=& A \\ F_2 &=& B \\ F_n &=& C\cdot{}F_{n-2}+D\cdot{}F_{n-1}+\left\lfloor\frac{P}{n}\right\rfloor \end{cases}$$

数据范围

  • $1 \leq T \leq 20$
  • $0 \leq A, B, C, D \leq 10^9$
  • $1 \leq P, n \leq 10^9$

思路

这题太硬核了我脑子里只有$O(n)$的算法肯定跑不过嘛

如果$\left\lfloor\frac{P}{n}\right\rfloor$的值是固定的那么肯定可以快速幂解决了,但是难就难在这个值会随着$n$变化。聪明的徐臣发现了当$n$很大的时候,这个值几乎是不变的,有很长一段区间可以用快速幂解决。

那么$n$多大的时候这样计算才不会浪费呢?其实根本不会浪费因为快速幂肯定跑的快啊!

对于开始值$i$,这一段快速幂计算区间为$\left[ i, P / (P / i) \right]$(注意除法会消掉余数,所以右端点正好是保持$\lfloor\frac{P}{N}\rfloor$值不变的最大值)。用矩阵快速幂疯狂向前推就可以了。

代码

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll modulo = (ll) 1e9 + 7;

class matrix{
 public:
  ll s[3][3];

  matrix(){
    memset(s, 0, sizeof(s));
  }
  matrix(ll Fn, ll Fn1, ll E) {
    s[0][0] = Fn;
    s[0][1] = Fn1;
    s[0][2] = E;
    s[1][0] = s[1][1] = s[1][2] = 0;
    s[2][0] = s[2][1] = s[2][2] = 0;
  }
  matrix(ll C, ll D) {
    s[0][0] = D;
    s[1][0] = C;
    s[0][1] = s[2][0] = s[2][2] = 1;
    s[0][2] = s[1][1] = s[1][2] = s[2][1] = 0;
  }
};

matrix mul(matrix &a, matrix &b) {
  matrix ret;
  for (int i = 0; i < 3; ++i) {
    for (int j = 0; j < 3; ++j) {
      for (int k = 0; k < 3; ++k) {
        ret.s[i][j] = (ret.s[i][j] + a.s[i][k] * b.s[k][j]) % modulo;
      }
    }
  }
  return ret;
}

ll A = 0, B = 0, C = 0, D = 0, P = 0, n = 0;

void xmul(ll &Fn, ll &Fn1, ll E, ll t) {
  matrix ret(Fn, Fn1, E);
  matrix b(C, D);
  while (t) {
    if (t & 1) {
      ret = mul(ret, b);
    }
    b = mul(b, b);
    t >>= 1;
  }
  Fn = ret.s[0][0];
  Fn1 = ret.s[0][1];
}

ll solve() {
  if (n == 1) return A;
  if (n == 2) return B;
  ll Fn1 = A, Fn = B, nex = 0, PN = 0;
  matrix res;
  for (ll i = 3; i <= n; ) {
    PN = P / i;
    if (PN == 0) {
      nex = n;
    } else {
      nex = min(n, P / PN);
    }
    xmul(Fn, Fn1, PN, nex - i + 1);
    i = nex + 1;
  }
  return Fn;
}

int main() {
  int T = 0;
  scanf("%d", &T);
  while (T--) {
    scanf("%lld%lld%lld%lld%lld%lld", &A, &B, &C, &D, &P, &n);
    printf("%lld\n", solve());
  }
}