SCUT - 354 - CC的简单多项式 - 杜教筛

Posted 2022-01-25 yinku

https://scut.online/p/354
跟多项式一点关系都没有。
注意到其实两个多项式在1处求值，那么就是他们的系数加起来。

列一列发现系数就是n以内两两求gcd的值，还自动把0去掉了。

那么就是
\(\sum\limits_i=1^n\sum\limits_i=1^ngcd(i^2,j^2)\)

这种情况就要枚举g但是为了方便我们也是枚举g而不是g平方
\(\sum\limits_g=1^ng^2\sum\limits_i=1^n\sum\limits_i=1^n[gcd(i^2,j^2)==g^2]\)

列一列gcd的分解式发现其实可以把平方约分掉
\(\sum\limits_g=1^ng^2\sum\limits_i=1^n\sum\limits_i=1^n[gcd(i,j)==g]\)

二话不说除以g
\(\sum\limits_g=1^ng^2\sum\limits_i=1^\lfloor\fracng\rfloor\sum\limits_i=1^\lfloor\fracng\rfloor[gcd(i,j)==1]\)

套上反演
\(\sum\limits_g=1^ng^2\sum\limits_i=1^\lfloor\fracng\rfloor\sum\limits_i=1^\lfloor\fracng\rfloor\sum\limits_k|gcd(i,j)\mu(k)\)
\(\sum\limits_g=1^ng^2\sum\limits_k=1^n\mu(k)\lfloor\fracngk\rfloor^2\)

枚举T
\(\sum\limits_T=1^n\lfloor\fracnT\rfloor^2\sum\limits_g|Tg^2\mu(\fracTg)\)

假如搞得出后面那个狄利克雷卷积的前缀和，那么可以分块回答，看看复杂度刚刚够的样子。

后面那个是
\(\sum\limits_g|Tg^2\mu(\fracTg)\)

也就是
\(id^2*\mu\)

嗷神提示卷一个东西，恒等函数\(I(n)=1\)
\((id^2*\mu)*I\)

结合律
\(id^2*(\mu*I)\)

后面那个就是
\(id^2*(\epsilon)\)

这个东西展开就是
\(\sum\limits_g|Tg^2[\fracTg==1]\)
\(T^2\)

也就是说卷积之后的前缀和是\(s_2\)，而恒等函数的前缀和就是\(s_0\)

上一波杜教筛就可以了

测试发现unorder_map和cc_hash_table差异非常小。

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;

inline ll read() 
    ll x = 0;
    //int f = 0;
    char c;
    do 
        c = getchar();
        /*if(c == '-')
            f = 1;*/
     while(c < '0' || c > '9');
    do 
        x = (x << 3) + (x << 1) + c - '0';
        c = getchar();
     while(c >= '0' && c <= '9');
    //return f ? -x : x;
    return x;


inline void _write(int x) 
    if(x > 9)
        _write(x / 10);
    putchar(x % 10 + '0');


inline void write(int x) 
    if(x < 0) 
        putchar('-');
        x = -x;
    
    _write(x);
    putchar('\n');


/*void TestCase(int ti);

int main() 
#ifdef Yinku
    freopen("Yinku.in", "r", stdin);
    //freopen("Yinku.out","w",stdout);
#endif // Yinku
    int T = 1;
    for(int ti = 1; ti <= T; ti++)
        TestCase(ti);
*/

/*---  ---*/

const int mod = 998244353;
const int inv2 = (mod + 1) >> 1;
const int MAXN = 1.5e7;

int pri[MAXN + 1];
int &pritop = pri[0];
int f[MAXN + 1];
int pk[MAXN + 1];

void sieve(int n = MAXN) 
    pk[1] = 1;
    f[1] = 1;
    for(int i = 2; i <= n; i++) 
        if(!pri[i]) 
            pri[++pritop] = i;
            pk[i] = i;
            f[i] = (1ll * i * i - 1ll) % mod;
        
        for(int j = 1; j <= pritop; j++) 
            int &p = pri[j];
            int t = i * p;
            if(t > n)
                break;
            pri[t] = 1;
            if(i % p) 
                pk[t] = p;
                f[t] = 1ll * f[i] * f[p] % mod;
             else 
                pk[t] = pk[i] * p;
                if(pk[t] == t) 
                    f[t] = 1ll * f[i] * p % mod * p % mod;
                 else 
                    f[t] = 1ll * f[t / pk[t]] * f[pk[t]] % mod;
                
                break;
            
        
    
    for(int i = 1; i <= n; i++) 
        f[i] = f[i - 1] + f[i];
        if(f[i] >= mod)
            f[i] -= mod;
    


/*inline int qpow(ll x, int n) 
    ll res = 1;
    while(n) 
        if(n & 1) 
            res *= x;
            if(res >= mod)
                res %= mod;
        
        x *= x;
        if(x >= mod)
            x %= mod;
        n >>= 1;
    
    return res;
*/

const int inv6 = 166374059; //qpow(6, mod - 2);

inline int s2(ll n) 
    if(n >= mod)
        n %= mod;
    ll tmp = n * (n + 1);
    if(tmp >= mod)
        tmp %= mod;
    tmp *= n * 2 + 1 ;
    if(tmp >= mod)
        tmp %= mod;
    tmp *= inv6;
    if(tmp >= mod)
        tmp %= mod;
    return tmp;


unordered_map<ll, int> Sf;

inline int F(ll n) 
    if(n <= MAXN)
        return f[n];
    if(Sf.count(n))
        return Sf[n];
    ll ret = s2(n);
    for(ll l = 2, r; l <= n; l = r + 1) 
        ll t = n / l;
        r = n / t;
        ret -= (r - l + 1) % mod * F(t) % mod;
        if(ret < 0)
            ret += mod;
    
    return Sf[n] = ret;


inline int S(ll n) 
    ll res = 0;
    for(ll l = 1, r; l <= n; l = r + 1) 
        ll t = n / l;
        r = n / t;
        ll tmp = F(r) - F(l - 1);
        if(tmp < 0)
            tmp += mod;
        if(t >= mod)
            t %= mod;
        t *= t;
        if(t >= mod)
            t %= mod;
        tmp *= t;
        if(tmp >= mod)
            tmp %= mod;
        res += tmp;
        if(res >= mod)
            res -= mod;
    
    return res;


int main() 
#ifdef Yinku
    freopen("Yinku.in", "r", stdin);
    //freopen("Yinku.out","w",stdout);
#endif // Yinku
    sieve();
    int T = read();
    while(T--) 
        write(S(read()));
    

但是好像发现一个问题，别人都是欧拉函数的？
\(\sum\limits_g=1^ng^2\sum\limits_i=1^\lfloor\fracng\rfloor\sum\limits_i=1^\lfloor\fracng\rfloor[gcd(i,j)==1]\)

这里内部记为
\(S(n)=\sum\limits_i=1^n\sum\limits_i=1^n[gcd(i,j)==1]\)

n以内互质对的个数？那么枚举较大的那个，小的那个要和他互质，就是欧拉函数，大小互换多了一倍，其中(1,1)重复去掉一个
\(S(n)=-1+2\sum\limits_i=1^n\varphi(i)\)

所以原式就是
\(\sum\limits_g=1^ng^2S(\lfloor\fracng\rfloor)\)

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;

inline ll read() 
    ll x = 0;
    //int f = 0;
    char c;
    do 
        c = getchar();
        /*if(c == '-')
            f = 1;*/
     while(c < '0' || c > '9');
    do 
        x = (x << 3) + (x << 1) + c - '0';
        c = getchar();
     while(c >= '0' && c <= '9');
    //return f ? -x : x;
    return x;


inline void _write(int x) 
    if(x > 9)
        _write(x / 10);
    putchar(x % 10 + '0');


inline void write(int x) 
    if(x < 0) 
        putchar('-');
        x = -x;
    
    _write(x);
    putchar('\n');


/*void TestCase(int ti);

int main() 
#ifdef Yinku
    freopen("Yinku.in", "r", stdin);
    //freopen("Yinku.out","w",stdout);
#endif // Yinku
    int T = 1;
    for(int ti = 1; ti <= T; ti++)
        TestCase(ti);
*/

/*---  ---*/

const int mod = 998244353;
const int inv2 = (mod + 1) >> 1;
const int MAXN = 1.6e7;

int pri[MAXN + 1];
int &pritop = pri[0];
int phi[MAXN + 1];

void sieve(int n = MAXN) 
    phi[1] = 1;
    for(int i = 2; i <= n; i++) 
        if(!pri[i]) 
            pri[++pritop] = i;
            phi[i] = i - 1;
        
        for(int j = 1; j <= pritop; j++) 
            int &p = pri[j];
            int t = i * p;
            if(t > n)
                break;
            pri[t] = 1;
            if(i % p) 
                phi[t] = phi[i] * phi[p];
             else 
                phi[t] = phi[i] * p;
                break;
            
        
    
    for(int i = 1; i <= n; i++) 
        phi[i] += phi[i - 1];
        if(phi[i] >= mod)
            phi[i] -= mod;
    


/*inline int qpow(ll x, int n) 
    ll res = 1;
    while(n) 
        if(n & 1) 
            res *= x;
            if(res >= mod)
                res %= mod;
        
        x *= x;
        if(x >= mod)
            x %= mod;
        n >>= 1;
    
    return res;
*/


inline int s1(ll n) 
    if(n >= mod)
        n %= mod;
    return n * (n + 1) % mod * inv2 % mod;


const int inv6 = 166374059; //qpow(6, mod - 2);
inline int s2(ll n) 
    if(n >= mod)
        n %= mod;
    ll tmp = n * (n + 1);
    if(tmp >= mod)
        tmp %= mod;
    tmp *= n * 2 + 1 ;
    if(tmp >= mod)
        tmp %= mod;
    tmp *= inv6;
    if(tmp >= mod)
        tmp %= mod;
    return tmp;


unordered_map<ll, int> Sphi;

inline int Phi(ll n) 
    if(n <= MAXN)
        return phi[n];
    if(Sphi.count(n))
        return Sphi[n];
    int ret = s1(n);
    for(ll l = 2, r, t; l <= n; l = r + 1) 
        t = n / l;
        r = n / t;
        ret -= (r - l + 1) % mod * Phi(t) % mod;
        if(ret < 0)
            ret += mod;
    
    return Sphi[n] = ret;


inline int S(ll n) 
    return (2 * Phi(n) - 1) % mod;


inline int Ans(ll n) 
    int res = 0;
    for(ll l = 1, r, t; l <= n; l = r + 1) 
        t = n / l;
        r = n / t;
        ll tmp = s2(r) - s2(l - 1);
        if(tmp < 0)
            tmp += mod;
        tmp *= S(t);
        if(tmp >= mod)
            tmp %= mod;
        res += tmp;
        if(res >= mod)
            res -= mod;
    
    return res;


int main() 
#ifdef Yinku
    freopen("Yinku.in", "r", stdin);
    //freopen("Yinku.out","w",stdout);
#endif // Yinku
    sieve();
    int T = read();
    while(T--) 
        write(Ans(read()));