标签：lfloor tmp frac 求和 EZEC sum rfloor int P6825

P6825 「EZEC-4」求和

求：

\[\sum_{i=1}^{n}\sum_{j=1}^{n}\gcd(i,j)^{i+j} \]

先化简原式，有：

\[\begin{aligned} &\sum_{i=1}^{n}\sum_{j=1}^{n}\gcd(i,j)^{i+j}\\ &=\sum_{d=1}^{n}\sum_{i=1}^{n}\sum_{j=1}^{n}d^{i+j}[\gcd(i,j)=d]\\ &=\sum_{d=1}^{n}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}d^{d\cdot(i+j)}[\gcd(i,j)=1]\\ &=\sum_{d=1}^{n}\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}d^{d\cdot(i+j)}\sum_{p|\gcd(i,j)}\mu(p)\\ &=\sum_{d=1}^{n}\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}\mu(p)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor}d^{d\cdot(i+j)}[p|i][p|j]\\ &=\sum_{d=1}^{n}\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}\mu(p)\sum_{i=1}^{\lfloor\frac{n}{pd}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{pd}\rfloor}d^{dp\cdot(i+j)}\\ &=\sum_{d=1}^{n}\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}\mu(p)\sum_{i=1}^{\lfloor\frac{n}{T}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{T}\rfloor}d^{T\cdot(i+j)} & T=dp\\ &=\sum_{d=1}^{n}\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}\mu(p)\sum_{i=1}^{\lfloor\frac{n}{T}\rfloor}(d^T)^i\sum_{j=1}^{\lfloor\frac{n}{T}\rfloor}(d^T)^j\\ &=\sum_{d=1}^{n}\sum_{p=1}^{\lfloor\frac{n}{d}\rfloor}\mu(p)\left(\sum_{i=1}^{\lfloor\frac{n}{T}\rfloor}(d^T)^i\right)^2 \end{aligned} \]

再看后面那部分，令 \(x=\frac{T}{p}^T\)。

设 \(S_n\) 为等比数列前 \(n\) 项和，则有：

\[S_n=S_{\lfloor\frac{n}{2}\rfloor}(1+a^{\lfloor\frac{n}{2}\rfloor})+a^n[2 \nmid n] \]

即可 \(\mathcal O(\log n)\) 求解。

时间复杂度 \(\mathcal O(n \log n)\)。

#include <bits/stdc++.h>

using namespace std;

#define ll long long

#define pll pair<ll, ll>

#define mpr make_pair

#define fi first

#define se second

const int _ = 1.5e6 + 10;

int mod;

int read()
{
    int x = 0, f = 1;
    char c = getchar();
    while (c < '0' || c > '9')
    {
        if (c == '-')
            f = -1;
        c = getchar();
    }
    while (c >= '0' && c <= '9')
    {
        x = x * 10 + c - '0';
        c = getchar();
    }
    return x * f;
}

int qpow(int x, int y)
{
    int res = 1;
    while (y)
    {
        if (y & 1)
            res = 1ll * res * x % mod;
        x = 1ll * x * x % mod;
        y >>= 1;
    }
    return res;
}

pll qsum(int a, int n)
{
    if (n == 0)
        return mpr(1, 1);
    if (n == 1)
        return mpr(a, a);
    pll tmp = qsum(a, n / 2);
    if (n & 1)
        return mpr((tmp.fi * (tmp.se + 1) % mod + tmp.se * tmp.se % mod * a % mod) % mod, tmp.se * tmp.se % mod * a % mod);
    else
        return mpr(tmp.fi * (tmp.se + 1) % mod, tmp.se * tmp.se % mod);
}

int T, n, m, k;

bitset<_> vis;

int cnt, pri[_], mu[_];

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++)
    {
        if (!vis[i])
            pri[++cnt] = i, mu[i] = -1;
        for (int j = 1; j <= cnt && i * pri[j] <= n; j++)
        {
            vis[i * pri[j]] = 1;
            if (i % pri[j] == 0)
                break;
            mu[i * pri[j]] = -mu[i];
        }
    }
}

int solve(int n)
{
    int ans = 0;
    for (int d = 1; d <= n; d++)
    {
        int a = qpow(d, d);
        int b = a;
        for (int p = 1; p <= n / d; p++)
        {
            pll tmp = qsum(b, n / (d * p));
            ans = (ans + tmp.fi * tmp.fi % mod * (mu[p] + mod) % mod) % mod;
            b = 1ll * b * a % mod;
        }
    }
    return ans;
}

signed main()
{
    T = read();
    init(_ - 10);
    while (T--)
    {
        n = read(), mod = read();
        printf("%d\n", solve(n));
    }
    return 0;
}

标签：lfloor,tmp,frac,求和,EZEC,sum,rfloor,int,P6825
来源： https://www.cnblogs.com/orzz/p/15808138.html

本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享；
2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关；
3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关；
4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除；
5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。