簡單的數學題

題目連接

https://www.luogu.org/problemnew/show/P3768

題目描述

輸入一個正整數$n,n\le 10^{10}$和$p,p\le 1.1\times 10^9$.且$p$為質數.

計算${\sum }_{i=1}^n{\sum }_{j=1}^{n}ijgcd(i,j)$對$p$取模.

題解

方法一.暴力推公式

考慮把枚舉$gcd(i,j)$
原式=${\sum }_{d=1}^n{d}^3{\sum }_{i=1}^{n/d}{\sum }_{j=1}^{n/d}ij\times [gcd(i,j)=1]$

莫比烏斯反演處理后邊的式子
設$f(x)={\sum }_{i=1}^{n/d}{\sum }_{j=1}^{n/d}ij\times [gcd(i,j)=x]$
$g(x)={\sum }_{x|d}f(d)=x^2{\sum }_{i=1}^{n/dx}{\sum }_{j=1}^{n/dx}ij$
故$f(x)={\sum }_{x|p}\mu (p/x)g(p)={\sum }_{x|p}\mu (p/x)p^2{\sum }_{i=1}^{n/dp}{\sum }_{j=1}^{n/dp}ij$
因此$f(1)={\sum }_{p=1}^n\mu (p)p^2{\sum }_{i=1}^{n/dp}{\sum }_{j=1}^{n/dp}ij$

合并兩部分
原式=${\sum }_{d=1}^n{d}^3{\sum }_{p=1}^n\mu (p)p^2{\sum }_{i=1}^{n/dp}{\sum }_{j=1}^{n/dp}ij$
轉而枚舉$a=dp$.我們得到:
原式$={\sum }_{a=1}^n{\sum }_{d|a}{a}^{2}d\mu (a/d){\sum }_{i=1}^{n/a}i{\sum }_{j=1}^{n/a}j$
$={\sum }_{a=1}^n(\frac{?\frac{n}{a}?(?\frac{n}{a}?+1)}{2}{)}^2{a}^2{\sum }_{d|a}d\mu (a/d)$
我們發現$Id?\mu =?$是常見的狄利克雷卷積.
因此原式$={\sum }_{a=1}^n(\frac{?\frac{n}{a}?(?\frac{n}{a}?+1)}{2}{)}^2{a}^2?(a)$
化簡到這一步的時候,前面部分可以分塊計算,后面的部分$a^2?(a)$要快速計算前綴和.
于是我們想到了杜教篩:
記$f(x)=x^2?(x)$,找到積性函數$g(x)=x^2$與它做卷積使得$g(x)$的前綴和可以快速求出,而$(f?g)(x)={\sum }_{d|x}{\frac{x}{d}}^2\times d^2?(d)=x^2{\sum }_{d|x}?(x)=x^3$的前綴和也可以快速求出.
記$S(x)={\sum }_{i=1}^{x}f(x)$
根據杜教篩公式$g(1)S(n)={\sum }_{x=1}^n(f?g)(x)?{\sum }_{x=2}^{n}g(x)S(\frac{n}{x})$.
$S(n)={\sum }_{x=1}^n{x}^3?{\sum }_{x=2}^n{x}^{2}S(\frac{n}{x})$
寫到這就完了.

方法二.$?$卷積

根據公式:

$gcd(i,j)={\sum }_{d|gcd(i,j)}?(d)={\sum }_{d|i,d|j}?(d)$

因此

${\sum }_{i=1}^n{\sum }_{j=1}^{n}ijgcd(i,j)={\sum }_{i=1}^n{\sum }_{j=1}^{n}ij{\sum }_{d|i,d|j}?(d)={\sum }_{d=1}^n?(d)d^2{\sum }_{i=1}^{n/d}{\sum }_{j=1}^{n/d}ij={\sum }_{d=1}^n?(d)d^2(\frac{?\frac{n}{d}?(?\frac{n}{d}?+1)}{2}{)}^2$

同樣也得到了我們的式子,是不是推導方法簡單多了?

遇見$gcd(i,j)$的時候,多想想是不是可以用$?$卷積來做?

可能會減少很多不必要的推導過程.

代碼

#include <iostream> #include <algorithm> #include <cstring> #include <unordered_map> #define pr(x) std::cout << #x << ':' << x << std::endl #define rep(i,a,b) for(LL i = a;i <= b;++i) const int N = 1e7; typedef long long LL; LL n,p; LL phi[N+10]; int prime[N+10],zhi[N+10],low[N+10],pcnt; LL mod_pow(LL x,LL n) {LL res = 1;while(n) {if(n&1) res = res * x % p;x = x * x % p;n >>= 1;}return res; } LL inv6,inv2,inv4; void sieve() {pcnt = 0;low[1] = phi[1] = zhi[1] = 1;rep(i,2,N) {if(!zhi[i]) {phi[i] = i-1;prime[pcnt++] = i;low[i] = i;}for(LL j = 0;j < pcnt && prime[j]*i <= N;++j) {zhi[i*prime[j]] = 1;if(i % prime[j] == 0) {low[i*prime[j]] = low[i] * prime[j];if(i == low[i]) {phi[i*prime[j]] = phi[i]*prime[j];}else {phi[i*prime[j]] = phi[i/low[i]] * phi[low[i]*prime[j]];}break;}else{low[i*prime[j]] = prime[j];phi[i*prime[j]] = phi[i] * phi[prime[j]];}}} } LL sum3(LL n) {n %= p;LL res = n*(n+1)%p*(2*n+1)%p*inv6%p;return res; } LL sum4(LL n) {n %= p;LL x = n*(n+1)%p;return x*x%p*inv4%p; } std::unordered_map<LL,LL> vis,rec; LL F(LL n) {if(n <= N) return phi[n];if(vis[n]) return rec[n];LL res = sum4(n);for(LL x = 2,last;x <= n;x = last) {last = n/(n/x)+1;res = (res - ((sum3(last-1)-sum3(x-1)+p)%p*F(n/x)%p) + p) % p;}vis[n] = 1;return rec[n] = res; }signed main() {sieve();std::ios::sync_with_stdio(false);std::cin >> p >> n;inv6 = mod_pow(6,p-2);inv4 = mod_pow(4,p-2);inv2 = mod_pow(2,p-2);rep(i,1,N) phi[i] = ((phi[i]*i%p*i%p) + phi[i-1]) % p;LL last,ans = 0;for(LL x = 1;x <= n;x = last+1) {last = n/(n/x);LL y = n/x%p;LL z = (1+y)*(y)/2%p;ans = (ans + (z*z%p*((F(last)-F(x-1)+p)%p))) % p;}std::cout << ans << std::endl;return 0; }

總結

以上是生活随笔為你收集整理的P3768 简单的数学题 [狄利克雷卷积,杜教筛,莫比乌斯反演]的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。