块速幂/光速幂

最近队友拿来一个蛮有意思的东西

给定一个类似于斐波那契数列的东西,让你去求它的第\(n\)\((n<=10^{18})\),由于模数为1e9+7,可以先降一波幂,相当于\(n<=10^9+7\),要么去推下特征方程要么矩阵快速幂,但是这题把带log的做法卡了,引出了一个被叫做光速幂的东西,但我感觉叫块速幂更正常点...

就叫块速幂吧,实际上很简单,如果求\(p^q\)\(p\)已经给定,\(q\)的上限也已经确定是在\(1e12\)左右或者以下,设\(sq=\sqrt{q}\)那么我们线性预处理处理出\(p \cdots p^{sq}\),再处理出\(p^{sq},p^{2sq}\cdots\)

那看到这个式子应该很显然了\(p^q=p^{q/sq}p^{sq}p^{q\%sq}\)

对于所有\(n<=1e12\)的询问我们都能\(o(1)\)计算

类似板子的东西(这玩意真的有必要存板子么

struct KSM
{
    int fac[MAXN];
    int facp[MAXN];
    ll p;               //给定p
    ll maxq = 1e9 + 10; //已知的p的幂次上限(maxq<=1e12)
    int sqrtq;
    void init(ll x)
    {
        sqrtq = int(sqrt(maxq)) + 1;
        p = x;
        //预处理到p^sqrtq
        fac[0] = 1;
        for (int i = 1; i <= sqrtq; i++)
            fac[i] = 1LL * fac[i - 1] * p % MOD;
        //再处理(p^sqtr)^sqtr
        facp[0] = 1;
        for (int i = 1; i <= sqrtq; i++)
            facp[i] = 1LL * facp[i - 1] * fac[sqrtq] % MOD;
    }
    int cal(int x)
    {
        return 1LL * facp[x / sqrtq] * fac[x % sqrtq] % MOD;
    }
}

队友拿来的那道题

传送门

解题过程

\(a_n=233a_{n-1}+666a_{n-2}\)

去网上搜下斐波那契通项公式怎么求就知道了

这种数列基本可以令\(a_n=p^n\)次将递推公式化成\(p^2=233p+666\)

解得\(p_1=\frac{233+\sqrt{56953}}{2}\), \(p_2=\frac{233-\sqrt{56953}}{2}\)

\(A,B\)为任意实数,\(a_n=Ap_1^n+Bp_2^n\)

\(a_1=1,a_2=233\)分别代入

如果把数值一开始带进去除了会受到233和666的毒害以外毫无意义..先把\(p_1,p_2\)放在里面推一波

可以得到\(B=\frac{233-p_1}{p_2(p_2-p_1)},A=\frac{233-p_2}{p_1(p_1-p_2)}\)

所以\(B=\frac{\frac{233-\sqrt{56953}}{2}}{p_2(p_2-p_1)}=\frac{p_2}{p_2(p_2-p_1)}=\frac{1}{p_2-p_1}=\frac{1}{\sqrt{56953}}\)

同理\(A=\frac{-1}{\sqrt{56953}}\)

所以\(a_n=\frac{(233+\sqrt{56953})^n-(233-\sqrt{56953})^n}{2^n\sqrt{56953}}\)

掏出二次剩余板子找到\(\sqrt{56953}\)的二次剩余是\(188305837\)

所以在膜1e9+7条件下

\(a_n=\frac{(\frac{233+188305837}{2})^n-(\frac{233-188305837}{2})^n}{188305837}\)

\(=\frac{(94153035)^n-(905847205)^n}{188305837}\)

好像算好了 还被卡了一手常数

#include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/hash_policy.hpp>
#include <ext/pb_ds/tree_policy.hpp>
#include <ext/pb_ds/trie_policy.hpp>
using namespace __gnu_pbds;
using namespace std;
// freopen("k.in", "r", stdin);
// freopen("k.out", "w", stdout);
// clock_t c1 = clock();
// std::cerr << "Time:" << clock() - c1 <<"ms" << std::endl;
//#pragma comment(linker, "/STACK:1024000000,1024000000")
mt19937 rnd(time(NULL));
#define de(a) cout << #a << " = " << a << endl
#define rep(i, a, n) for (int i = a; i <= n; i++)
#define per(i, a, n) for (int i = n; i >= a; i--)
#define ls ((x) << 1)
#define rs ((x) << 1 | 1)
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> PII;
typedef pair<double, double> PDD;
typedef pair<char, char> PCC;
typedef pair<ll, ll> PLL;
typedef vector<int> VI;
#define inf 0x3f3f3f3f
const ll INF = 0x3f3f3f3f3f3f3f3f;
const ll MAXN = 1e5 + 7;
const ll MAXM = 4e5 + 7;
const ll MOD = 1e9 + 7;
const double eps = 1e-7;
ll qpw(ll a, ll b)
{
    ll ans = 1;
    while (b)
    {
        if (b & 1)
            ans = (ans * a) % MOD;
        (a *= a) %= MOD;
        b >>= 1;
    }
    return ans;
}
struct KSM
{
    int fac[MAXN];
    int facp[MAXN];
    ll p;               //给定p
    ll maxq = 1e9 + 10; //已知的p的幂次上限(maxq<=1e12)
    int sqrtq;
    void init(ll x)
    {
        sqrtq = int(sqrt(maxq)) + 1;
        p = x;
        //预处理到p^sqrtq
        fac[0] = 1;
        for (int i = 1; i <= sqrtq; i++)
            fac[i] = 1LL * fac[i - 1] * p % MOD;
        //再处理(p^sqtr)^sqtr
        facp[0] = 1;
        for (int i = 1; i <= sqrtq; i++)
            facp[i] = 1LL * facp[i - 1] * fac[sqrtq] % MOD;
    }
    int cal(int x)
    {
        return 1LL * facp[x / sqrtq] * fac[x % sqrtq] % MOD;
    }
} a, b;

namespace Mker
{
    unsigned long long SA, SB, SC;
    void init() { scanf("%llu%llu%llu", &SA, &SB, &SC); }
    unsigned long long rand()
    {
        SA ^= SA << 32, SA ^= SA >> 13, SA ^= SA << 1;
        unsigned long long t = SA;
        SA = SB, SB = SC, SC ^= t ^ SA;
        return SC;
    }
} // namespace Mker
int del(int a, int b)
{
    int tmp = a - b;
    if (tmp < 0)
        tmp += MOD;
    return tmp;
}
int main()
{
    int t;
    a.init(94153035), b.init(905847205);
    int inv = qpw(188305837, MOD - 2);
    int ans = 0;
    scanf("%d", &t);
    Mker::init();
    while (t--)
    {
        int n = Mker::rand() % (MOD - 1);
        ans ^= 1LL * del(a.cal(n), b.cal(n)) * inv % MOD;
    }
    printf("%d\n", ans);
    return 0;
}

原文链接: https://www.cnblogs.com/graytido/p/13917063.html

欢迎关注

微信关注下方公众号,第一时间获取干货硬货;公众号内回复【pdf】免费获取数百本计算机经典书籍

    块速幂/光速幂

原创文章受到原创版权保护。转载请注明出处:https://www.ccppcoding.com/archives/204705

非原创文章文中已经注明原地址,如有侵权,联系删除

关注公众号【高性能架构探索】,第一时间获取最新文章

转载文章受原作者版权保护。转载请注明原作者出处!

(0)
上一篇 2023年2月12日 下午9:57
下一篇 2023年2月12日 下午9:57

相关推荐