数学专题学习笔记
开坑于-2022.1.24
update1-2022.1.26
update2-2022.2.4
update3-2022.2.14
update4-2022.2.16
这是一个从约数开始的数学专题,至于能写多长我也不知道,有时间就学一点。
约数
约数应该没有什么好有疑问的地方,无非就是求一个数的约数,或者将题目转换成与约数有关的内容。
感觉大多数时是被其他算法利用的一个工具。
一些例题
P2926 [USACO08DEC]Patting Heads
一道比较简单的题目,考虑将问题的答案全部预处理出。
由于 \(\text{n}\) 只有 \(1e6\) ,所以直接统计每个数的出现次数,然后直接枚举,预处理所有的答案即可。
复杂度为 \(O(n \sqrt n )\)。
点击查看代码
#include <bits/stdc++.h>
using namespace std;
int n , a[100010] , cnt[1000010] , ans[1000010];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
int main()
{
n = read();
for(int i = 1;i <= n;i++)
a[i] = read() , cnt[a[i]]++;
for(int i = 1;i <= 1e6 + 10;i++)
{
for(int j = i;j <= 1e6 + 10;j += i)
{
ans[j] += cnt[i];
}
}
for(int i = 1;i <= n;i++)
cout << ans[a[i]] - 1 << endl;
return 0;
}
根据题意,我们发现直接求答案十分麻烦。
考虑将式子转换一下。
由于 \(\text{y}\) 为正整数。
所以
又因为 \(\text{x}\) , \(\text{y}\) 均为正整数。
所以 \(x > n!\)。
所以问题就变成了求 \(n!^2\) 的约数个数。
故直接求就可以了。
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int mod = 1e9 + 7;
const int maxn = 1e6 + 10;
int n , ans , cnt;
int prime[maxn] , isprime[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void init()
{
for(int i = 2;i <= n;i++)
{
if(!isprime[i]) prime[cnt++] = i;
for(int j = 0;prime[j] * i <= n;j++)
{
isprime[prime[j] * i] = 1;
if(i % prime[j] == 0) break;
}
}
}
signed main()
{
n = read() , ans = 1;
init();
for(int i = 0;i < cnt;i++)
{
int x = prime[i] , y = 0;
for(int j = n;j;j /= x) y += j / x;
ans = ans * (2 * y + 1) % mod;
}
cout << ans << endl;
return 0;
}
这里用了一个小技巧。
即求大数阶乘的约数。
我们可以将其质因数分解,这样如何求约数就解决了。
至于如何将其质因数分解,我们很容易就能知道,在 \(\text{n}\) 个数里,有 \(\lfloor\frac{n}{k}\rfloor\) 个是 \(\text{k}\) 的倍数。
那么就可以质因数分解了。
由于 \(n \leq 2 \times 10^9\) 。
所以需要考虑的质数只有九个,而最大的指数只有三十。
我们发现就可以直接暴搜求出答案。
点击查看代码
#include <bits/stdc++.h>
using namespace std;
int n , ans , num;
int prime[10] = {0 , 2 , 3 , 5 , 7 , 11 , 13 , 17 , 19 , 23};
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void dfs(int u , int last , int sum , int s)
{
if(s > num || (sum < ans && num == s))
{
num = s;
ans = sum;
}
if(u == 10) return;
for(int i = 1;i <= last;i++)
{
if((long long)sum * prime[u] > n) break;
sum *= prime[u];
dfs(u + 1 , i , sum , s * (i + 1));
}
}
int main()
{
n = read();
dfs(1 , 30 , 1 , 1);
cout << ans;
return 0;
}
欧拉函数
欧拉函数是小于等于 \(\text{x}\) 的数中与 \(\text{x}\) 互质的数的数目。
其中 \(a_{i}\) 为 \(\text{x}\) 的质因子,\(\text{n}\) 为 \(\text{x}\) 的因子个数。
几个性质
-
若 \(\text{x}\) 为质数,显然 \(φ(x)=x-1\)
-
若 \(m,n\) 互质,则 \(φ(mn)=φ(m)\timesφ(n)\)
-
若 \(x=2n~,~n \mod 2 = 1\),则 \(φ(x)=φ(n)\)
-
若 \(\text{x}\) 为质数,则 \(φ(x^k)= x^k-x^{k-1}\)
-
若 \(x>2\),则 \(φ(x) \mod 2 = 0\)
-
小于 \(n\) 且与 \(n\) 互质的数的总和为 \(φ(n) \times n \div 2(n>1)\)
-
\(n=\sum_{d|n}φ(d)\),即 \(\text{n}\) 包括 \(\text{n}\) 在内的所有因数的欧拉函数之和等于 \(\text{n}\),这一条也被称为欧拉反演。
求法
有两种求法。
- \(o(\sqrt n)\) 的时间复杂度求出单个欧拉函数。
点击查看代码
inline int eular(int x)
{
int res = x;
for(int i = 2;i * i <= x;i++)
if(x % i == 0)
{
res = res / i * (i - 1);
while(x % i == 0) x /= i;
}
if(x > 1) res = res / x * (x - 1);
return res;
}
- \(O(n\log(n))\) 的时间复杂度求出\(1-n\)欧拉函数。
点击查看代码
inline int eular(int x)
{
for(int i = 1;i <= x;i++) e[i] = i;
for(int i = 2;i <= x;i++)
if(e[i] == i)
for(int j = i;j < x;j++)
e[j] = e[j] / i * (i - 1);
}
一些例题
有空了再写题。
高斯消元
高斯消元主要用于求解线性方程组,方程组主要形式为:
如何实现
这里不介绍朴素的高斯消元了。
介绍一下高斯——约旦消元。
我们以样例为例,模拟一下消元步骤。
首先,列出增广矩阵。
将此时遍历到的一列中最大的数设为主元。
然后对其他所有的行,进行消元。
不断循环上步骤,最终将矩阵变成下矩阵。
格式化的结果
最终,输出答案,答案为 \(\frac{b_i}{a_i}\)
点击查看代码
#include <bits/stdc++.h>
using namespace std;
int n;
double a[110][110];
inline int read()
{
int s = 0 , f = 1; char ch;
while(!isdigit(ch = getchar())) if(ch == '-') f = -1;
while(isdigit(ch)) s = s * 10 + ch - '0' , ch = getchar();
return s * f;
}
int main()
{
n = read();
for(int i = 1;i <= n;i++)
{
for(int j = 1;j <= n + 1;j++)
a[i][j] = read();
}
for(int i = 1;i <= n;i++)
{
int m = i;
for(int j = i + 1;j <= n;j++)
{
if(fabs(a[j][i]) > fabs(a[m][i]))
m = j;
}
for(int j = 1;j <= n + 1;j++) swap(a[i][j] , a[m][j]);
if(int(a[i][i]) == 0)
{
puts("No Solution");
return 0;
}
for(int j = 1;j <= n;j++)
{
if(j != i)
{
double d = a[j][i] / a[i][i];
for(int k = i + 1;k <= n + 1;k++)
{
a[j][k] -= a[i][k] * d;
}
}
}
}
for(int i = 1;i <= n;i++)
printf("%.2lf\n" , a[i][i]);
return 0;
}
一些例题
一道十分套路的题。
特判一下无解与无穷解。
特判方法可以根据方程的特性,若是无穷解,则出现零是,最右边必然全是零,否则无解。
其余的就是板子了。
点击查看代码
#include <bits/stdc++.h>
using namespace std;
int n , cnt;
double a[100][100];
inline int read()
{
int s = 0 , f = 1; char ch;
while(!isdigit(ch = getchar())) if(ch == '-') f = -1;
while(isdigit(ch)) s = s * 10 + ch - '0' , ch = getchar();
return s * f;
}
int main()
{
n = read() , cnt = 1;
for(int i = 1;i <= n;i++)
for(int j = 1;j <= n + 1;j++)
a[i][j] = read();
for(int i = 1;i <= n;i++)
{
int m = cnt;
for(int j = cnt + 1;j <= n;j++)
if(fabs(a[j][i]) > fabs(a[m][i])) m = j;
if(a[m][i] == 0) continue;
for(int j = 1;j <= n + 1;j++)
swap(a[cnt][j] , a[m][j]);
for(int j = 1;j <= n;j++)
if(cnt != j)
{
double w = a[j][i] / a[cnt][i];
for(int k = i + 1;k <= n + 1;k++)
a[j][k] -= a[cnt][k] * w;
}
++cnt;
}
if(cnt <= n)
{
while(cnt <= n)
{
// cout << a[cnt][n] << endl;
if(a[cnt++][n + 1] != 0)
{
puts("-1");
return 0;
}
}
puts("0");
return 0;
}
for(int i = 1;i <= n;i++)
printf("x%d=%.2lf\n" , i , (int(a[i][n + 1] / a[i][i] * 100) == 0 ? 0 : a[i][n + 1] / a[i][i]));
return 0;
}
这道题主要就是列出线性方程,其余的就可以直接套板子。
由于一个球体的球心到球面的任意一个点距离相等,因此,我们需要求的 \((x_1 , x_2 , x_3 , x_4 ... x_n)\) 需要满足:
我们可以相邻两个作差。
化成矩阵即:
最后,就可以直接套板子了。
点击查看代码
#include <bits/stdc++.h>
using namespace std;
int n;
double a[20][20] , b[20][20];
inline int read()
{
int s = 0 , f = 1; char ch;
while(!isdigit(ch = getchar())) if(ch == '-') f = -1;
while(isdigit(ch)) s = s * 10 + ch - '0' , ch = getchar();
return s * f;
}
int main()
{
n = read();
for(int i = 1;i <= n + 1;i++)
for(int j = 1;j <= n;j++)
cin >> a[i][j];
for(int i = 1;i <= n;i++)
for(int j = 1;j <= n;j++)
b[i][j] = 2 * (a[i][j] - a[i + 1][j]),
b[i][n + 1] += (a[i][j] * a[i][j] - a[i + 1][j] * a[i + 1][j]);
for(int i = 1;i <= n;i++)
{
int m = i;
for(int j = i + 1;j <= n;j++)
if(fabs(b[j][i]) > fabs(b[m][i])) m = j;
for(int j = 1;j <= n + 1;j++)
swap(b[i][j] , b[m][j]);
for(int j = 1;j <= n;j++)
if(i != j)
{
double w = b[j][i] / b[i][i];
for(int k = i + 1;k <= n + 1;k++)
b[j][k] -= w * b[i][k];
}
}
for(int i = 1;i <= n;i++)
printf("%.3lf " , b[i][n + 1] / b[i][i]);
return 0;
}
欧几里得算法
非常简单的一个算法。
拓展欧几里得算法
扩展欧几里得算法是用来在已知 \((a,b)\) 时,求解一组 \((p,q)\) ,使得 \(p \times a+q \times b=\gcd(a,b)\)。
可以简单的阐述一下这个算法:
因为:
所以:
程序只需要不断递归操作,最后当 \(\text{b}\) 减少到 \(\text{0}\) 时,就可以得出 \(p=1,q=0\) ,然后递归回去就可以求出 \(\text{p,q}\) 了。
点击查看代码
inline int exgcd(int a , int b , int &x , int &y)
{
if(!b)
{
x = 1, y = 0;
return a;
}
int res = exgcd(b , a % b , x , y) , tmp = x;
x = y , y = tmp - a / b * y;
return res;
}
一些例题
模板题。
首先,必然可知,当且仅当 \(\gcd(a,b)|c\) 时,方程有整数解
然后,我们可以让方程变为拓展欧几里德算法可以求解的方程。
我们先求出:
那么
即为原方程的一组特解。
而整数通解为:
其中 \(x,y\) 的最小增减幅度分别为 \(\displaystyle\frac{b}{\gcd(a,b)}, \displaystyle\frac{a}{\gcd(a,b)}\) 即 \(b',a'\)。
接着就是,实现细节了。
1.无解
直接判断就可以了。
2.无正整数解
根据二元一次方程的通解可得,\(x,y\) 的增减性相反,所以若存在正整数解,则当 \(x\) 取得最小正整数解时, \(y\) 取得最大正整数解,若此时 \(y\le 0\),则说明该方程无正整数解。
3.有正整数解
同样,根据二元一次方程通解,求出 \(x,y\) 的最大值和最小值。
而正整数解的数量就为 \(\displaystyle\frac{x_{max}-x_{min}}{b'}+1\) 或 \(\displaystyle\frac{y_{max}-y_{min}}{a'}+1\)
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
int t;
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline int exgcd(int a , int b , int &x , int &y)
{
if(!b)
{
x = 1, y = 0;
return a;
}
int res = exgcd(b , a % b , x , y) , tmp = x;
x = y , y = tmp - a / b * y;
return res;
}
signed main()
{
t = read();
while(t--)
{
int a = read() , b = read() , c = read();
int x = 0 , y = 0 , z = __gcd(a , b);
if(c % z != 0)
{
puts("-1");
continue;
}
exgcd(a , b , x , y);
x = x * c / z , y = y * c / z;
a /= z , b /= z;
int t , xmax , xmin , ymax , ymin;
t = ceil((double)(1 - x) / b);
xmin = x + t * b , ymax = y - t * a;
t = ceil((double)(1 - y) / a);
xmax = x - t * b , ymin = y + t * a;
if(ymax <= 0) cout << xmin << " " << ymin << endl;
else cout << (xmax - xmin) / b + 1 << " " << xmin << " " << ymin << " " << xmax << " " << ymax << endl;
}
return 0;
}
几乎与上一题一模一样,只需要将式子改一下就可以了。
求出最小的 \(t\) 即可。
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
int t;
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline int exgcd(int a , int b , int &x , int &y)
{
if(!b)
{
x = 1, y = 0;
return a;
}
int res = exgcd(b , a % b , x , y) , tmp = x;
x = y , y = tmp - a / b * y;
return res;
}
signed main()
{
int a = read() , b = read() , m = read() , n = read() , l = read();
int x = 0 , y = 0 , z;
if((a - b) % (z = exgcd(n - m , l , x , y)) != 0)
puts("Impossible");
else cout << ((a - b) / z % (l / z) * x % (l / z) + l) % (l / z);
return 0;
}
BSGS算法
BSGS算法,全称Baby-Step-Giant-Step算法,又称大步小步算法。
它用于计算离散对数。
即一个最小的非负整数 \(l\),满足 \(a^l \equiv n \pmod p\)
我们可以运用整除分块的思想,将 \(l\) 设为 \(v\sqrt{p}-w\)。
显然有:
那么:
那么我们可以把右式的值存到一个哈希表里,每次遍历 \(v\) 时直接查表就可以了。
复杂度 \(O(\sqrt{p})\)
模板题
点击查看代码
#include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>
#define int long long
using namespace std;
using namespace __gnu_pbds;
int n , b , p , s;
gp_hash_table<int , int> q;
inline int read()
{
int s = 0 , f = 1; char ch;
while(!isdigit(ch = getchar())) if(ch == '-') f = -1;
while(isdigit(ch)) s = s * 10 + ch - '0' , ch = getchar();
return s * f;
}
inline int power(int x , int y , int mod)
{
int res = 1;
while(y)
{
if(y & 1) res = res * x % mod;
x = x * x % mod , y /= 2;
}
return res % mod;
}
signed main()
{
p = read() , b = read() , n = read() , s = ceil(sqrt(p));
for(int i = 1;i <= s;i++)
q[(n * power(b , i , p) % p)] = i;
for(int i = 1;i <= s;i++)
{
int res = power(b , i * s , p);
if(q.find(res) != q.end())
{
// cout << i << " " << s << " " << res << " " << q[res] << endl;
cout << i * s - q[res];
return 0;
}
}
cout << "no solution";
return 0;
}
exBSGS
拓展 \(BSGS\) 算法,即 \(ExBSGS\)。
在 \(BSGS\) 中,存在限制 \(\gcd(p, a) = 1\),如果两者不互质的话需要使用拓展 \(BSGS\) 算法。
继续考虑 \(\gcd(p,a)\),设其等于 \(d\)。
我们在原方程 \(a^x ≡ b \pmod p\) 中除掉一个 \(d\),变成\(\frac{a}{d}a^{x-1} \equiv \frac{b}{d} \pmod {\frac{p}{d}}\)。
重复这个过程,直到 \(gcd\) 为 1。或者再这个过程中 \(b\) 不能整除 \(d\),则无解。
于是问题变成了:
此时 \(gcd\) 为 1,可以使用 \(BSGS\) 求解。
这道题卡常属实有点狠。
直接套板子就可以了。
但要注意,因为这道题太过毒瘤,快速幂要写成递推式。
点击查看代码
#include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>
#define int long long
using namespace std;
using namespace __gnu_pbds;
int n , b , p , s , k , l;
gp_hash_table<int , int> q;
inline int read()
{
int s = 0 , f = 1; char ch;
while(!isdigit(ch = getchar())) if(ch == '-') f = -1;
while(isdigit(ch)) s = s * 10 + ch - '0' , ch = getchar();
return s * f;
}
inline int power(int x , int y , int mod)
{
int res = 1;
while(y)
{
if(y & 1) res = res * x % mod;
x = x * x % mod , y /= 2;
}
return res % mod;
}
inline void exbsgs()
{
q.clear() , b %= p , n %= p;
if(n == 1 || p == 1)
{
puts("0");
return;
}
l = 1 , k = 0;
while(1)
{
int res = __gcd(b , p);
if(res == 1) break;
if(n % res != 0)
{
puts("No Solution");
return;
}
p /= res , n /= res , l = (l * (b / res)) % p , k++;
if(l == n)
{
printf("%lld\n" , k);
return;
}
}
s = ceil(sqrt(p));
int ret = 1;
for(int i = 1;i <= s;i++)
{
ret = ret * b % p;
q[(n * ret % p)] = i;
}
ret = power(b , s , p);
int sumres = 1;
for(int i = 1;i <= s;i++)
{
sumres = sumres * ret % p;
int res = (sumres * l) % p;
if(q.find(res) != q.end())
{
printf("%lld\n" , i * s - q[res] + k);
return;
}
}
puts("No Solution");
}
signed main()
{
b = read() , p = read() , n = read();
while(b || p || n) exbsgs() , b = read() , p = read() , n = read();
return 0;
}
至于例题大多都直接套板子就可以做出来,这里就不多放题了。
中国剩余定理 CRT
求解方程组
但原版的 \(CRT\) 限制巨大,需要 \(m\) 全都互质,所以这里就不再介绍了,直接就介绍 \(exCRT\)
拓展中国剩余定理 exCRT
此时不保证所有的 \(m_i\) 互质。那么显然最终的答案就是在模 \(lcm(m_1, m_2, ..., m_n)\) 意义下进行的。
考虑如何合并两个方程的答案。
假设就是合并:
首先肯定要合并模数为 \(lcm(m_1, m_2)\) 。
显然有:
\(x = k_1m_1 + a_1 = k_2m_2 + a_2\)。
化简后就是 \(k_1m_1 ? k_2m_2 = a_2 ? a_1\)
\(m_1, m_2, a_1, a_2\) 都是常数,这玩意就是一个 \(exgcd\) 求解了。
如果 \(exgcd\) 无解显然整个方程组都无解,否则求解后得到 \(x\),合并模
数之后就得到了一个新的方程,然后接着把所有方程都合并就好了。
一些例题
点击查看代码
#include <bits/stdc++.h>
#define int __int128
using namespace std;
const int maxn = 100010;
int n , ans1 , ans2 , a[maxn] , b[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline int lcm(int a , int b)
{
return a / __gcd(a , b) * b;
}
inline void exgcd(int a , int b , int &x , int &y)
{
if(!b)
{
x = 1 , y = 0;
return;
}
exgcd(b , a % b , x , y);
int tmp = x;
x = y , y = tmp - a / b * y;
}
signed main()
{
n = read();
for(int i = 1;i <= n;i++) a[i] = read() , b[i] = read() % a[i];
ans1 = a[1] , ans2 = b[1];
for(int i = 2;i <= n;i++)
{
int x = 0 , y = 0 , gcd = __gcd(ans1 , a[i]);
exgcd(ans1 , a[i] , x , y);
int c = b[i] - ans2;
x = x * c / gcd % (a[i] / gcd);
if(x < 0) x += a[i] / gcd;
int mod = lcm(ans1 , a[i]);
ans2 = (ans1 * x + ans2) % mod;
if(ans2 < 0) ans2 += mod;
ans1 = mod;
}
cout << (long long)(ans2 % ans1);
return 0;
}
莫比乌斯反演
这个东西包含的内容很多,会时不时更新,并不能一次性学完,可能开始的内容及其基础,后面可能才会更加深入。
莫比乌斯函数:
设 \(n = \prod^m_{i=1}a_i^{p_i}\) 。
若 \(\forall i, p_i > 1\),那么 \(μ(n) = 0\)
否则 \(μ(n) = (-1)^m\)。
莫比乌斯反演:
\(g(n)=\sum\limits_{d|n}f(d)\)
\(f(n)=\sum\limits_{d|n}g(d)\mu(\frac{n}{d})\)
这个东西最主要的是多写题,多运用,熟悉的就更快。
一些例题
可以推一下式子:
根据函数意义:
这个式子就可以用整除分块来做了。
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int maxn = 50010;
int n , ans , tot , mu[maxn] , sum[maxn] , prime[maxn];
bool vis[maxn] , isprime[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void get_mu()
{
int s = 50000;
mu[1] = 1;
for(int i = 2;i <= s;i++)
{
if(!vis[i]) mu[i] = -1 , vis[i] = 1 , prime[++tot] = i , isprime[i] = 1;
for(int j = 1;j <= tot && i * prime[j] <= s;j++)
{
vis[prime[j] * i] = 1;
if(i % prime[j] == 0) break;
else mu[prime[j] * i] = -mu[i];
}
}
for(int i = 1;i <= s;i++) sum[i] = sum[i - 1] + mu[i];
}
signed main()
{
get_mu() , n = read();
for(int i = 1;i <= n;i++)
{
int a = read() , b = read() , d = read() , x = min(a , b);
ans = 0;
for(int l = 1 , r = 0;l <= x;l = r + 1)
r = min(a / (a / l) , b / (b / l)),
ans += (a / (l * d)) * (b / (l * d)) * (sum[r] - sum[l - 1]);
cout << ans << endl;
}
return 0;
}
虽然又是一道莫比乌斯反演模板题,但这个蒟蒻还是没能独立做出来。
但至少还是做出来了。
同样的,推一下式子:
我们设:\(n<m\)
我们知道莫比乌斯函数的性质:
所以:
枚举 \(k\) :
设 \(t = k\times d\):
前面一堆直接整除分块,后面的预处理前缀和即可。
复杂度\(O(t\sqrt{n})\)
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int maxn = 10000010;
int t , n , m , cnt , ans , h[maxn] , mu[maxn] , sum[maxn] , prime[maxn];
bool vis[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void get_mu()
{
int s = 1e7; mu[1] = 1;
for(int i = 2;i <= s;i++)
{
if(!vis[i]) mu[i] = -1 , prime[++cnt] = i , vis[i] = 1;
for(int j = 1;j <= cnt && i * prime[j] <= s;j++)
{
vis[i * prime[j]] = 1;
if(i % prime[j] == 0) break;
mu[i * prime[j]] = -mu[i];
}
}
for(int i = 1;i <= cnt;i++)
for(int j = 1;j * prime[i] <= s;j++) h[j * prime[i]] += mu[j];
for(int i = 1;i <= s;i++) sum[i] = sum[i - 1] + h[i];
}
signed main()
{
t = read() , get_mu();
while(t--)
{
n = read() , m = read() , ans = 0;
if(n > m) swap(n , m);
for(int l = 1 , r = 0;l <= n;l = r + 1)
{
r = min(n / (n / l) , m / (m / l));
ans += (n / l) * (m / l) * (sum[r] - sum[l - 1]);
}
printf("%lld\n" , ans);
}
return 0;
}
这一道题可以看做第一题的拓展:
那么根据容斥:
直接整除分块求解即可。
点击查看代码
#include <bits/stdc++.h>
using namespace std;
const int maxn = 50010;
int n , cnt , mu[maxn] , sum[maxn] , prime[maxn];
bool vis[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void get_mu()
{
int s = 50000; mu[1] = 1;
for(int i = 2;i <= s;i++)
{
if(!vis[i]) mu[i] = -1 , prime[++cnt] = i , vis[i] = 1;
for(int j = 1;j <= cnt && i * prime[j] <= s;j++)
{
vis[i * prime[j]] = 1;
if(i % prime[j] == 0) break;
mu[i * prime[j]] = -mu[i];
}
}
for(int i = 1;i <= s;i++) sum[i] = sum[i - 1] + mu[i];
}
inline int solve(int k , int x , int y)
{
int ans = 0;
if(x > y) swap(x , y);
for(int l = 1 , r = 0;l <= x;l = r + 1)
{
r = min(x / (x / l) , y / (y / l));
ans += (sum[r] - sum[l - 1]) * (x / (l * k)) * (y / (l * k));
}
return ans;
}
int main()
{
n = read() , get_mu();
while(n--)
{
int a = read() , b = read() , c = read() , d = read() , k = read();
cout << solve(k , b , d) - solve(k , a - 1 , d) - solve(k , b , c - 1) + solve(k , a - 1 , c - 1) << endl;
}
return 0;
}
推式子,推式子(〃'▽'〃):
令:
设:
则:
整除分块 \(O(\sqrt n)\) 求解即可。
带回原式:
太棒了,你会发现又是一个整除分块。
时间复杂度: \(O(n)\)
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int mod = 20101009;
const int maxn = 10000010;
int n , m , cnt , sum , inv4 , mu[maxn] , sum1[maxn] , sum2[maxn] , prime[maxn];
bool vis[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void get_mu()
{
int s = 1e7; mu[1] = 1;
for(int i = 2;i <= s;i++)
{
if(!vis[i]) mu[i] = -1 , vis[i] = 1 , prime[++cnt] = i;
for(int j = 1;j <= cnt && prime[j] * i <= s;j++)
{
vis[i * prime[j]] = 1;
if(i % prime[j] == 0) break;
mu[i * prime[j]] = -mu[i];
}
}
for(int i = 1;i <= s;i++) sum2[i] = (sum2[i - 1] + i) % mod;
for(int i = 1;i <= s;i++) sum1[i] = ((sum1[i - 1] + mu[i] * i * i % mod) % mod + mod) % mod;
}
inline int power(int x , int y)
{
int res = 1;
while(y)
{
if(y & 1) res = res * x % mod;
x = x * x % mod , y /= 2;
}
return res;
}
inline int inv(int x)
{
return power(x , mod - 2);
}
inline int calc(int x , int y)
{
return x * (x + 1) % mod * y % mod * (y + 1) % mod * inv4 % mod;
}
inline int solve(int x , int y)
{
int ans = 0;
if(x > y) swap(x , y);
for(int l = 1 , r = 0;l <= x;l = r + 1)
{
r = min(x / (x / l) , y / (y / l));
ans = (ans + ((sum1[r] - sum1[l - 1]) % mod + mod) % mod * calc(x / l , y / l) % mod) % mod;
}
return ans;
}
signed main()
{
n = read() , m = read();
inv4 = inv(4) , get_mu();
if(n > m) swap(n , m);
for(int l = 1 , r = 0;l <= n;l = r + 1)
{
r = min(n / (n / l) , m / (m / l));
sum = (sum + (sum2[r] - sum2[l - 1]) * ((solve(n / l , m / l) % mod + mod) % mod) % mod) % mod;
}
cout << (sum % mod + mod);
return 0;
}
上来就直接推式子:
首先,必然会有:
所以:
预处理:
所以原式:
好像现在这样就可以直接前缀和+整除分块了。
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int maxn = 50010;
int t , n , m , cnt , ans , g[maxn] , mu[maxn] , sum[maxn] , prime[maxn];
bool vis[maxn];
inline int read()
{
int asd = 0 , qwe = 1; char zxc;
while(!isdigit(zxc = getchar())) if(zxc == '-') qwe = -1;
while(isdigit(zxc)) asd = asd * 10 + zxc - '0' , zxc = getchar();
return asd * qwe;
}
inline void get_mu()
{
int s = 50000; mu[1] = 1;
for(int i = 2;i <= s;i++)
{
if(!vis[i]) mu[i] = -1 , prime[++cnt] = i , vis[i] = 1;
for(int j = 1;j <= cnt && prime[j] * i <= s;j++)
{
vis[i * prime[j]] = 1;
if(i % prime[j] == 0) break;
mu[i * prime[j]] = -mu[i];
}
}
for(int i = 1;i <= s;i++) sum[i] = sum[i - 1] + mu[i];
for(int i = 1;i <= s;i++) for(int l = 1 , r = 0;l <= i;l = r + 1)
r = i / (i / l) , g[i] += (i / l) * (r - l + 1);
}
signed main()
{
t = read() , get_mu();
while(t--)
{
n = read() , m = read() , ans = 0;
if(n > m) swap(n , m);
for(int l = 1 , r = 0;l <= n;l = r + 1)
{
r = min(n / (n / l) , m / (m / l));
ans += (sum[r] - sum[l - 1]) * g[n / l] * g[m / l];
}
cout << ans << endl;
}
return 0;
}