HDU 4947 GCD Array 反演+树状数组维护

【题目大意】

有一个长度为l的数组，编号从1到l。两种操作。1、给定n，d，v，所有满足gcd（x，n）== d 的 ax 都加上v。2、询问前x个元素的和。

【思路】

如果n%d！=0，操作1无意义；如果n%d==0，利用反演，操作1可以视为对所有ai，加上 v* [ gcd（i/d，n/d）==1  ]，（[ ]符号表示，里面判断为真，值为1，否则值为0）可以变化为


     即：     


我们可以考虑维护一个pd数组，具体来说，就是对于操作1，枚举（n/d）的所有因子，使得dp[ 因子*d ] += v*u（p）。

有了pd数组，则 ai = ∑pd[ i的所有因子 ]。那么 ∑ai = （i/1）*pd[1] + （i/2）*pd[2] + （i/3）*pd[3] + ... + （i/i）*pd[i]

这个可以分段处理，每段的和可以用树状数组处理。这样的复杂度是 sqrt(x)*log(x)，这题很卡时间...这样不能过...

对于i <= sqrt（x*log（x）），我们暴力求，效率sqrt（x*log（x））；后面一部分，我们分段求，因为组多只有x / sqrt（x*log（x））段，每段用树状数组有一个log（x），综合效率还是sqrt（x*log（x））

//#pragma comment(linker, "/STACK:102400000,102400000")
#include<cstdio>
#include<cstring>
#include<vector>
#include<queue>
#include<cmath>
#include<cctype>
#include<string>
#include<algorithm>
#include<iostream>
#include<ctime>
#include<map>
#include<set>
using namespace std;
#define MP(x,y) make_pair((x),(y))
#define PB(x) push_back(x)
typedef long long LL;
//typedef unsigned __int64 ULL;
/* ****************** */
const LL INF = 1LL<<55;
const double INFF = 1e100;
const double eps = 1e-8;
const LL mod = 10000000007LL;
const int NN = 50010;
const int MM = 5000010;
/* ****************** */

bool vis[NN];
int prime[NN], mu[NN];
vector<int>yinzi[NN*4];
LL c[NN], pd[NN];

int getint()
{
char c;
while((c = getchar()), (c<'0' || c>'9'));
int d = c - '0';
while((c = getchar()), (c>='0' && c<='9'))
d = d*10 + c -'0';
return d;
}

void init()
{
int i, j, cnt = 0;
memset(vis, false, sizeof(vis));

mu[1] = 1;
for(i = 2; i <= 50000; i ++)
{
if(!vis[i])
{
prime[cnt++] = i;
mu[i] = -1;
}
for(j = 0; j < cnt && i*prime[j] <= 50000; j ++)
{
vis[i*prime[j]] = true;
if(i%prime[j])
mu[i*prime[j]] = -mu[i];
else
{
mu[i*prime[j]] = 0;
break;
}
}
}

for(i = 1; i <= 50000; i ++)
if(mu[i])
for(j = i; j <= 200000; j += i)
yinzi[j].PB(i);
}

int lowbit(int x)
{
return x&(-x);
}
void modify(int x,int n,LL ad)
{
for ( ; x <= n; x += lowbit(x))
c[x] += ad;
}
LL get_sum(int x)
{
LL sum = 0;
for( ; x > 0; x -= lowbit(x))
sum += c[x];
return sum;
}

int main()
{
init( );
int ee = 0;
int n, m;
LL t, pre, now;
int op, nn, d, x, limit, i, j, si, jj, temp;
while (scanf("%d%d", &n, &m) != EOF)
{
if (n==0 && m==0) break;
printf("Case #%d:\n", ++ee);

memset(c, 0, sizeof(c));
memset(pd, 0, sizeof(pd));

// limit = sqrt(n*log(n+0.0));
// cout<<"limit=="<<limit<<endl;

limit = 400;

for (i = 1; i <= m; i ++)
{
//scanf("%d", &op);
op = getint();
if (op == 1)
{
// scanf("%d%d%I64d", &nn, &d, &t);
nn = getint();
d = getint();
t = getint();
if (nn%d == 0)
{
si = yinzi[nn/d].size();
for (j = 0; j < si; j ++)
{
temp = yinzi[nn/d][j]*d;
9eee

if(temp <= n)
pd[temp] += t*mu[temp/d];
modify (temp, n, t*mu[temp/d]);
}
}
}
else
{
//scanf("%d", &x);
x = getint();
t = 0;
//limit = sqrt(x * log(x+0.0));
if (x <= limit)
{
for ( j = 1; j <= x; j ++)
{
t += (x/j)*pd[j];
}
}
else
{
pre = 0;
for ( j = 1; j <= limit; j ++)
{
t += x/j*pd[j];
pre += pd[j];
}
for( ; j <= x; j = jj + 1)
{
jj = x/(x/j);
now = get_sum(jj);
t += (x/j)*(now-pre);
pre = now;
}
}
printf("%I64d\n",t);
}
}
}
return 0;
}