答：欧几里得+......+
int gcd(int a,int b)
{
return b==0?a:gcd(b,a%b);
}
这段gcd代码大家一定很熟悉^_^。gcd嘛！又叫欧几里得算法，这和我们今天要说的莫比乌斯有(非常^n)密切关系！
这两天做了李总挂的一套莫比乌斯的题目深有感触。那我们赶紧看看这两个人是如何快乐的玩耍的！
1.莫比乌斯函数定义
,其中每个p[i]都是素数。
也就是说一个数x，它的素因子分解中只要有一个素因子的幂次超过1那么它的莫比乌斯函数值直接为0，否则就按素因子的个数r，如果r是奇数那就是-1，偶数就是1，此外，一般 。
这是个比较简单的函数，取值只有三种-1,1,0.先来算单个莫比乌斯函数值(这个基本不用)
int Mobius(int x)
{
if(x==1) return 1;
int cnt=0;
for(int i=2;i*i<=x;i++)
if(x%i==0)
{
cnt++;
x/=i;
if(x%i==0) return 0;//i这个素因子已经超过一个
}
if(x>1) cnt++;
return cnt%2?-1:1;
}
常用的是利用线性筛法去打表，O(n)算出一定范围内的莫比乌斯函数值。
务必理解一下线性筛，它价值不至于筛素数，而且基于线性筛可以快速算各种各样的函数值。记住线性筛素数时，每个数总是被它的最小素因子筛掉(自己手动模拟加深理解)
const int N=1e7+5;
int pri[N],tot=0,mb[N];
bool vis[N];
void init()
{
memset(vis,false,sizeof vis);
mb[1]=1;//care for it
for(int i=2;i<N;i++)
{
if(!vis[i])
{
pri[tot++]=i;
mb[i]=-1;//i 是素数
}
for(int j=0;j<tot && i*pri[j]<N;j++)
{
vis[i*pri[j]]=true;
//这里需要一种递推思想
if(i%pri[j]) mb[i*pri[j]]=-mb[i];
//素数pri[j]是i没有的，所以增加一个素数就在mb[i]基础上乘mb[pri[j]]=-1;
else {
mb[i*pri[j]]=0;//i*pri[j]至少含有两个pri[j]
break;
}
}
}
}
这个模板很常用。单纯的莫比乌斯函数就这些。
2.莫比乌斯函数的由来及简单应用
最初，莫比乌斯函数就是容斥原理在数论方面的一种应用，从代码上说就是容斥的递归形式改写成了循环的形式。
来看一个简单的问题：给n个素数p[i]，求有多少个a,满足1<=a<=m && a不是任意一个p[i]的倍数(这里不考虑数据范围，假设都很少(我知道小数据可以暴力))
显然用容斥的思想很容易：设f[i]表示[1,m]内是i的倍数的个数，f[i]=m/i,然后容斥，不妨假设n=3,p[1]=2,p[2]=3,p[3]=5;
那么,当数量较多时显然用递归写比较方便嘛！
int ans=0;
void dfs(int pos,int choose,int tag)
{
if(pos>n)
{
ans+=tag*f(choose);
return;
}
//每个数选或者不选
dfs(pos+1,choose*p[pos],-tag);
dfs(pos+1,choose,tag)
}
我们把注意力放在tag变量上，也就是系数，只有1或-1，有些要加，有些要减，仔细观察发现当你选了奇数个素因子时是-1，偶数的话就是1，这正好就对应莫比乌斯函数，而且for循环的话，莫比乌斯函数值的0正好把无关的没计算在内，它作为系数而存在。而这个问题p[i]是离散的用莫比乌斯函数写比容斥更加麻烦，然而这都不是重点，至少你得先理解清楚莫比乌斯函数，做过容斥原理专题的，不妨看看有没有可以转换成莫比乌斯函数去做的，虽然有点牵强，因为后面的内容比较抽象、枯燥，所以务必理解莫比乌斯函数与容斥的关系。
既然莫比乌斯函数是容斥原理在数论上的一种应用，当它遇到了gcd时，就青出于蓝了！！！
3.莫比乌斯反演
在讲这个前，我先列举下，莫比乌斯反演通常能解决什么样的问题：
1.求有多少对a,b满足gcd(a,b)=1,其中1<=a<=n,1<=b<=m;都考虑有序对，(2,1)和(1,2)是不一样的,不对(1,1)是唯一的。
2.
求有多少对a,b满足gcd(a,b)是素数,其中1<=a<=n,1<=b<=m;
3.求有多少对a,b满足gcd(a,b)=k,其中1<=a<=n,1<=b<=m;
4.求有多少对a,b满足gcd(a,b)是完全平方数，其中1<=a<=n,1<=b<=m;
5.a,b的范围是一些输入的数，即离散的集合。
6.求
7........
总的来说就是求一个集合(不一定每次都是[1,n])里gcd满足一些乱七八糟的条件的对数，这是计数，显然计数都可以用来求和或者求乘积，因为无非就是算一个数出现多少次，然后算贡献值,这个很重要。
好了，开始讲啦！
先来看看反演：
设f(x)表示gcd(a,b)=x的对数(a,b都是有一定范围的)，这个函数一般就是我们直接要去求的，这个一般很好设出来。
然后设：，(新手注意里面有个数论符号 d|x 表示x是d的倍数，专业术语叫d整除x ,比如F(4)=f(1)+f(2)+f(4) )先不要管F(x)的含义，这是设F(x)的套路，也就说，f(x)是根据题目变化的，F(x)和f(x)关系确实不变的。当然设完之后就要考虑F(x)的含义了。一般通过F(x)的含义可以很容易算出F(x),后面具体题目会讲的.那么现在问题是:我们要求的是f(x),而F(x)很好求，这个时候就用反演，用好求的F(x)，来表示f(x);
,第一次看这个反演公式一般都是一脸懵逼的+......+,所以必须解释下！
不过这个反演公式似乎不常用，还有一个：
，此处对应的F(x)和上面不一样：,d是有上限的，注意两个公式区别，下面主要解释这个反演公式。
来个例题：
求有多少对a,b满足gcd(a,b)=1,其中1<=a<=n,1<=b<=m;最基础最简单的一个
设f(x)表示gcd(a,b)=1的对数，1<=a<=n,1<=b<=m
设F(x):，当然d肯定不超过min(n,m)，先看F(x)的含义,首先f(d)是gcd=d对数，而d是x的倍数,那么F(x)就是gcd(a,b)是x倍数的对数,只要一对(a,b)的gcd是x的倍数就统计在内。
求F(x):既然gcd(a,b)是x的倍数，那么a=x*i,b=x*j，只要a,b不越界，i,j可以任意取值，i,j可以看成是x的倍数，记住我们是要计数，有多少对(i,j),这个时候条件已经没有什么限制，[1,n]肯定有个i满足条件,同理[1,m]是个j
先组合就是任选i,j就是两个相乘，所以,一下就算出来了
用反演公式：，然后要求的是f(1)值
仔细观察这个式子，回想前面说的容斥与莫比乌斯函数的关系，这个和式就是容斥的循环版本：
比如令min(n,m)=3吧，d=1就是ans+=F(1)，结合F(x)的含义，就是加所有gcd是1的倍数的数，然后d=2是mb[2]=-1;
ans+=-F(2),减去gcd是2的倍数的,最后减去gcd是3的倍数的。。哈哈，是容斥吧！还不理解就手动写几个看看。
当然这里F(d)里面是一个整除，所以循环可以分块加速，复杂度O(sqrt(n))，具体方法查看分块加速
bzoj 2818
求gcd为素数的对数，a,b,都是[1,n]，n<=1e7.
这个上一题一样的反演公式，得到反演公式后就是枚举当x为素数时的答案在求和就是行了
p是素数,这样会tle，所以需要求和技巧：交换求和次序。这里我们枚举p，然后对p的每个倍数d去求和，反过来考虑一个d被运算了多少次，就是d有多少个素因子，d的每个素因子都会使得d运算一次，所以：似乎看起来没什么区别，重点在于我们枚举的对象不一样，现在我们枚举的是d，然后p是素数同时也是d的因子。这样枚举的好处是让d/p相对于一个固定的d而已是独立的，我们单独处理第二个和式，令 ,这就是个函数，d确定了函数值就定了，我们枚举一个d正好符合函数,注意F(d)不能放进来，它还与n有关.
最终，既然把T[d]来出来，显然我是要打表咯，这里线性筛就发挥作用了。
其实还是挺麻烦的。
就和筛莫比乌斯函数一样，代码中我们只要考虑三个地方就可以了，i是素数，i%pri[j]!=0 和i%pri[j]==0的情况。
最好自己先想想。
1.i是素数很简单。能得到T[i]的值
2.if(i%pri[j])我们算T[i*pri[j] ],按照求和公式需要考虑它（注意不是i而是i*pri[j]）每个素因子，不妨设为pi,显然pri[j]本身就是某一个pi,既然是求和，当pi==pri[j]时，否则把i/pi看成整体，利用莫比乌斯函数的积性：
u[i*pri[j]]=u[i]*u[pri[j]]=-u[i](注意前提是i%pri[j]!=0),所以，最后两种加起来
T[i*pri[j]]=u[i]-T[i];
3.if(i%pri[j]==0)同上讨论pi,当pi==pri[j]时与分母抵消一个pri[j]，就等于u[i],其他的就简单了pri[j]抵消不了，那i*pri[j]的pri[j]的幂次超过1那就结果直接是0,所以u[i*pri[j]]=u[i].
当然这个显然也要求T[i]的前缀和，因为要对F[d]分块
#include<cstdio>
#include<iostream>
#include<cstring>
#include<vector>
#include<algorithm>
using namespace std;
const int N=1e7+9;
typedef long long ll;
int pri[N],tot=0;
int mb[N];
bool vis[N];
int sum[N];
void init()
{
memset(vis,false,sizeof vis);
mb[1]=1;
sum[0]=0;//这里偷个懒，sum[i]筛素数时就是刚刚说的T[i]，然后在求一个前缀和
sum[1]=0;
for(int i=2; i<N; i++)
{
if(!vis[i])
{
pri[tot++]=i;
mb[i]=-1;
sum[i]=1;
}
for(int j=0; j<tot && i*pri[j]<N; j++)
{
vis[i*pri[j]]=true;
if(i%pri[j])
{
mb[i*pri[j]]=-mb[i];
sum[i*pri[j]]=mb[i]-sum[i];
}
else
{
mb[i*pri[j]]=0;
sum[i*pri[j]]=mb[i];
break;
}
}
}
for(int i=1; i<N; i++) //前缀和
sum[i]+=sum[i-1];
}
int main()
{
int n;
init();
scanf("%d",&n);
ll ans=0;
for(int i=2,la=0; i<=n; i=la+1)
{
la=n/(n/i);
ans+=(ll)(sum[la]-sum[i-1])*(n/i)*(n/i);
}
printf("%lld\n",ans);
return 0;
}
再来一道同类型的 HDU 5663，和上一题类似，枚举平方数然后方法都一样
为了继续讲下去先要学个定理
看起来就是个很简单的公式，只有n是1时那么莫比乌斯和才为1，其他都是0.这公式用来装换gcd
再看一个和式小标的写法：
,不难理解gcd(a,b)是d的倍数，那么a,b一定都是d倍数
bzoj 2154 这题求
我们知道lcm(a,b)=a*b/gcd(a,b);
设i=a*gcd(i,j),j=b*gcd(i,j)
那么我们枚举gcd的话，假设枚举gcd=d，此时gcd(i,j)=d的条件显然就是gcd(i/d,j/d)=1,我们在和式后面乘上一个布尔表达式使之正确：
,干掉布尔表达式
，
d已经只与第一个和式有关了，提到外面，然后枚举t
,把t提出来，里面其实就是枚举t的倍数i,j(这和最初的那个变量不一样),因为上一步的a,b都是t的倍数。里面几个变量已经独立了。
,最后两个和式就是等差求和嘛，然后对于每一个d,把n/d，m/d看成整体，那么可以对t进行分块(当然d也可以分块嘛)，只要有整除t的就行了，有分块就要有前缀和,这里t*t要和u[t]合在一起打表
#pragma comment(linker, "/STACK:10240000,10240000")
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<queue>
#include<set>
#include<vector>
#include<map>
#include<stack>
#include<cmath>
#include<algorithm>
using namespace std;
const double R=0.5772156649015328606065120900;
const int N=1e7+5;
const int mod=20101009;
const int INF=0x3f3f3f3f;
const double eps=1e-8;
const double pi=acos(-1.0);
typedef long long ll;
int pri[N],tot=0,mb[N];
bool vis[N];
ll sum[N];
void init(int n)
{
n++;
memset(vis,false,(n+2)*sizeof(vis[0]));
sum[0]=0;
mb[1]=1;
sum[1]=1;
for(int i=2;i<n;i++)
{
if(!vis[i]) {
mb[i]=-1;
pri[tot++]=i;
}
for(int j=0;j<tot && i*pri[j]<n;j++)
{
vis[i*pri[j]]=true;
if(i%pri[j]) mb[i*pri[j]]=-mb[i];
else {
mb[i*pri[j]]=0;
break;
}
}
sum[i]=sum[i-1]+(ll)mb[i]*i*i;//打表
sum[i]%=mod;
}
}
ll ff(ll n)//等差求和
{
return (ll)n*(n+1)/2%mod;
}
ll get(int n,int m)
{
ll ans=0;
for(int i=1,la=0;i<=n;i=la+1)//枚举t，分块
{
la=min(n/(n/i),m/(m/i));
ans=(ans+((sum[la]-sum[i-1])*ff(n/i)%mod)*ff(m/i))%mod;
}
return ans;
}
int main()
{
int n,m;
scanf("%d%d",&n,&m);
if(n>m) n^=m^=n^=m;
init(n);
ll ans=0;
for(ll i=1,la=0;i<=n;i=la+1)//枚举d
{
la=min(n/(n/i),m/(m/i));
ll t=get(n/i,m/i);//每一个d都要去算里面的和式
ans=(ans+t*((la+i)*(la-i+1)/2%mod))%mod;//这里d的一个前缀就是等差嘛
}
ans+=mod;
ans%=mod;
printf("%lld\n",ans);
return 0;
}
bzoj 3994 &codeforces 235E
这两题要用同一个定理
d(n)表示n的约数个数
给个题解连接吧 235E 3994
235E代码：(代码很丑，建议理解题解，然后自己写)
#include<cstdio>
#include<iostream>
#include<cstring>
#include<vector>
using namespace std;
const int N=2005;
const int mod=1073741824;
typedef long long ll;
vector<int>g[N];
int dp[N][N];
int gcd(int a,int b)
{
if(b==0) return a;
if(dp[a][b]!=-1) return dp[a][b];
return dp[a][b]=gcd(b,a%b);
}
int mb[N];
int pri[N],tot=0;
bool vis[N];
void init()
{
memset(vis,false,sizeof vis);
mb[1]=1;
for(int i=2;i<N;i++)
{
if(!vis[i]) {
pri[tot++]=i;
mb[i]=-1;
}
for(int j=0;j<tot && i*pri[j]<N;j++)
{
vis[i*pri[j]]=true;
if(i%pri[j]) mb[i*pri[j]]=-mb[i];
else {
mb[i*pri[j]]=0;
break;
}
}
}
memset(dp,-1,sizeof dp);
for(int i=1;i<N;i++)
{
for(int j=1;j<N;j++)
if(gcd(i,j)==1) g[i].push_back(j);//预处理因子
}
}
int main()
{
int a,b,c;
init();
scanf("%d%d%d",&a,&b,&c);
ll ans=0;
int up=min(b,c);
//其实就直接三重循环枚举
for(int i=1;i<=a;i++)//one
{
int n=g[i].size();
ll sum=0;
for(int ii=0;ii<n && g[i][ii]<=up;ii++)//two
{
int d=g[i][ii];
if(mb[d]==0) continue;
ll t1=0;
int up1=b/d;
ll tmp=mb[d];
for(int jj=0;jj<n && g[i][jj]<=up1;jj++)//three
t1=(t1+up1/g[i][jj])%mod;
tmp=tmp*t1%mod;
int up2=c/d;
ll t2=0;
for(int kk=0;kk<n && g[i][kk]<=up2;kk++)//three
t2=(t2+up2/g[i][kk])%mod;
tmp=tmp*t2%mod;
sum=(sum+tmp)%mod;
}
ans=(ans+a/i*sum)%mod;
}
printf("%I64d\n",ans);
return 0;
}
HDU 4947
这题主要是要用树状数组维护。
对于每一个更新操作可以看做a[x]+=v*[gcd(x,n)==d),令x=d*i,n=d*j;所以gcd(i,j)==1,依旧是把比尔表示换成和式
,变量比较多，n,d看成常量，我们要更新的是a[x]，那么看x满足的条件：首先就是对于每个(n/d)的因子t，x必须是t*d的倍数，也就是说满足这个条件的x都要更新，不过我们又要枚举t又要枚举x复杂度很高
我们还是枚举t,但是我想一次就把所有满足条件的x更新完，我对a[x]进行拆分：
,f[j]是自己构造的一个函数，如果对于每个t*d，我更新f[t*d]自然在求和的过程中，每个符合的a[x]都会把更新的f[t*d]加上去，至于复杂度我们还要继续推到，现在我们看每次query的答案等于什么：
,看似复杂，我交换求和变量：
,每次更新f[i]的值，用树状数组维护，后面有整除又可以分块。
#pragma comment(linker, "/STACK:10240000,10240000")
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<queue>
#include<set>
#include<vector>
#include<map>
#include<stack>
#include<cmath>
#include<algorithm>
using namespace std;
const double R=0.5772156649015328606065120900;
const int N=1e5+5;
const int mod=1e9+7;
const int INF=0x3f3f3f3f;
const double eps=1e-8;
const double pi=acos(-1.0);
typedef long long ll;
const int M=2e5+5;
int pri[M],tot=0,mb[M];
bool vis[M];
vector<int>g[M];
void init()
{
memset(vis,false,sizeof vis);
mb[1]=1;
for(int i=2;i<M;i++)
{
if(!vis[i]){
mb[i]=-1;
pri[tot++]=i;
}
for(int j=0;j<tot && i*pri[j]<M;j++)
{
vis[i*pri[j]]=true;
if(i%pri[j]) mb[i*pri[j]]=-mb[i];
else {
mb[i*pri[j]]=0;
break;
}
}
}
for(int i=1;i<M;i++)
{
for(int j=i;j<M;j+=i)
g[j].push_back(i);
}
}
ll f[N];
ll sum(int i)
{
ll s=0;
while(i>0) {
s+=f[i];
i-=i&-i;
}
return s;
}
int ri;
void add(int i,ll x)
{
while(i<=ri)
{
f[i]+=x;
i+=i&-i;
}
}
void update(int n,int d,ll v)
{
if(n%d) return;
n=n/d;
int re=g[n].size();
for(int i=0;i<re;i++)//枚举每个n/d的因子t去更新f[t*d]
{
int j=g[n][i];
if(j*d>ri) break;
add(j*d,v*mb[j]);
}
}
ll query(int n)
{
ll ans=0;
for(int i=1,la=0;i<=n;i=la+1)
{
la=n/(n/i);
ans+=(ll)(sum(la)-sum(i-1))*(n/i);
}
return ans;
}
int main()
{
init();
int q;
int kase=1;
while(~scanf("%d%d",&ri,&q) && ri+q)
{
printf("Case #%d:\n",kase++);
memset(f,0,(ri+10)*sizeof(f[0]));
while(q--)
{
int typ;
scanf("%d",&typ);
int n,d,v;
if(typ==1){
scanf("%d%d%d",&n,&d,&v);
update(n,d,v);
}else {
scanf("%d",&n);
printf("%I64d\n",query(n));
}
}
}
return 0;
}
/*
*/
现在我们来看一些离散的问题
面对这类问题一定要明确一个事实：n的因子个数不会很多，比如n<=1e5是，最多也就128个，自己打表试试吧！
这类问题总是要把集合里的数按照因子倍数分类，即f[i]表示i的倍数的个数，然后用容斥的思想，我一般写法是莫比乌斯函数
codeforces 548E
题意就是说给出n个数a[1]...a[n]，刚开始集合是空的，然后有q个查询，每个查询输入一个下标x，如果集合中已经有第x个元素a[x]了就把它删掉，否则就放进去，求每次操作后集合元素中右多少对互质的。
有种dp思想，如果每次对集合元素求对数的话很难，那就考虑每次操作的那个数对答案的影响，最初是0，插入一个数就看这个元素能增加多少对互质的数，删除也是要更新被删掉数的影响的，那么问题就变成了一个数与一个集合里的数有多少互质，我们保存集合里的元素的f[i]值：是i的倍数的有多少个，利用f[i]值容斥一下就可以了，枚举a[x]的因子，对 f[ a[x]的因子 ] 容斥，当然莫比乌斯函数就可以了。
#pragma comment(linker, "/STACK:10240000,10240000")
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<queue>
#include<set>
#include<vector>
#include<map>
#include<stack>
#include<cmath>
#include<algorithm>
using namespace std;
const double R=0.5772156649015328606065120900;
const int N=5e5+5;
const int mod=1e9+7;
const int INF=0x3f3f3f3f;
const double eps=1e-8;
const double pi=acos(-1.0);
typedef long long ll;
vector<int>g[N];
int f[N];
int a[N];
int pri[N],tot=0,mb[N];
bool vis[N],has[200009];
void init()
{
memset(vis,false,sizeof vis);
mb[1]=1;
for(int i=2;i<N;i++)
{
if(!vis[i])
{
mb[i]=-1;
pri[tot++]=i;
}
for(int j=0;j<tot && i*pri[j]<N;j++)
{
vis[i*pri[j]]=true;
if(i%pri[j]) mb[i*pri[j]]=-mb[i];
else {
mb[i*pri[j]]=0;
break;
}
}
}
}
void resolve(int x)
{
if(g[x].size()!=0) return;//已经分解了的数不用再分解了
int up=sqrt(1.0*x);
for(int i=1;i<=up;i++)
if(x%i==0)
{
g[x].push_back(i);
if(i*i!=x) g[x].push_back(x/i);
}
}
void Insert(int x,ll& ans)
{
resolve(x);
int re=g[x].size();
for(int i=0;i<re;i++)
{
int d=g[x][i];
ans+=f[d]*mb[d];
f[d]++;//更新f[d]
}
}
void Delete(int x,ll& ans)
{
resolve(x);
int re=g[x].size();
for(int i=0;i<re;i++)
{
int d=g[x][i];
f[d]--;
ans-=f[d]*mb[d];//和Insert类似
}
}
int main()
{
init();
int n,q;
scanf("%d%d",&n,&q);
for(int i=1;i<=n;i++)
scanf("%d",&a[i]);
memset(f,0,sizeof f);
memset(has,false,sizeof has);
ll ans=0;//持续更新
while(q--)
{
int x;
scanf("%d",&x);
if(has[x]) {
Delete(a[x],ans);
has[x]=false;//标记是否存在
} else {
Insert(a[x],ans);
has[x]=true;
}
printf("%I64d\n",ans);
}
return 0;
}
/*
*/
是
