zoukankan      html  css  js  c++  java
  • 「NOI2016」循环之美(小性质+min_25筛)

    传送门.

    题解


    感觉这题最难的是第一个结论。

    数值上不等,所以x/y首先要互质,然后如果在10进制是纯循环小数,不难想到y不是2、5的倍数就好了。

    因为十进制下除以2和5是除得尽的。

    必然会多出来的什么东西。

    如果是k进制,可以类比得(gcd(y,k)=1)

    证明:

    假设纯循环的位数是(l)

    (x*k^lequiv x(mod~y))

    (k^lequiv 1(mod~y))

    要存在(l)的话,就必须有(gcd(k,y)=1),反过来一样。

    反演:

    (Ans=sum_{i=1}^nsum_{j=1}^m[(i,j)=1]*[(j,k)=1]))

    (sum_{d=1}^{min(n,m)} mu(d)*[gcd(d,k)=1]*(n/d)*sum_{j=1}^{m/d}[gcd(j,k)=1])

    如果能解决前面的前缀和问题,后面的分块后容斥或者说预处理O(1)都不是问题。

    然后你发现前面是个积性函数,而且是可以min_25筛的积性函数,不过要筛个两遍。

    就没了……

    Code:


    #include<bits/stdc++.h>
    #define fo(i, x, y) for(int i = x, B = y; i <= B; i ++)
    #define ff(i, x, y) for(int i = x, B = y; i <  B; i ++)
    #define fd(i, x, y) for(int i = x, B = y; i >= B; i --)
    #define ll long long
    #define pp printf
    #define hh pp("
    ")
    using namespace std;
    
    int k;
    
    int d[20], d0;
    
    void fen(int x) {
    	for(int i = 2; i * i <= x; i ++) if(x % i == 0) {
    		d[++ d0] = i;
    		while(x % i == 0) x /= i;
    	}
    	if(x > 1) d[++ d0] = x;
    }
    
    namespace m25 {
    const int N = 2e5 + 5;
    int n, sq;
    int w[N], w0, i1[N], i2[N];
    int p[N], p0, bz[N];
    
    void sieve(int n) {
    	fo(i, 2, n) {
    		if(!bz[i]) p[++ p0] = i;
    		for(int j = 1; i * p[j] <= n; j ++) {
    			int k = i * p[j];
    			bz[k] = 1; if(i % p[j] == 0) break;
    		}
    	}
    }
    
    int g[N];
    #define num(x) ((x) <= sq ? i1[x] : i2[n / (x)])
    
    int sp[N];
    
    void build() {
    	sq = sqrt(n); sieve(sq);
    	for(int i = 1, j; i <= n; i = j + 1) {
    		j = n / (n / i);
    		w[++ w0] = n / i;
    		if(w[w0] <= sq) i1[w[w0]] = w0; else i2[i] = w0;
    		g[w0] = w[w0] - 1;
    	}
    	fo(j, 1, p0) for(int i = 1; p[j] * p[j] <= w[i] && i <= w0; i ++) {
    		int k = num(w[i] / p[j]);
    		g[i] -= (g[k] - (j - 1));
    	}
    	int l = 0;
    	fd(i, w0, 1) {
    		while(l < d0 && d[l + 1] <= w[i]) l ++;
    		g[i] -= l;
    	}
    	fo(i, 1, w0) g[i] = -g[i];
    	fo(i, 1, p0) sp[i] = sp[i - 1] - (k % p[i] != 0);
    	for(int j = p0; j >= 1; j --) if(k % p[j]) for(int i = 1; p[j] * p[j] <= w[i] && i <= w0; i ++) {
    		int k = num(w[i] / p[j]);
    		g[i] += (g[k] - sp[j]) * -1;
    	}
    	fo(i, 1, w0) g[i] ++;
    }
    
    ll q1(int x) {
    	return g[num(x)];
    }
    }
    using m25 :: q1;
    
    namespace m26 {
    const int N = 2e5 + 5;
    int n, sq;
    int w[N], w0, i1[N], i2[N];
    int p[N], p0, bz[N];
    
    void sieve(int n) {
    	fo(i, 2, n) {
    		if(!bz[i]) p[++ p0] = i;
    		for(int j = 1; i * p[j] <= n; j ++) {
    			int k = i * p[j];
    			bz[k] = 1; if(i % p[j] == 0) break;
    		}
    	}
    }
    
    int g[N];
    #define num(x) ((x) <= sq ? i1[x] : i2[n / (x)])
    
    int sp[N];
    
    void build() {
    	sq = sqrt(n); sieve(sq);
    	for(int i = 1, j; i <= n; i = j + 1) {
    		j = n / (n / i);
    		w[++ w0] = n / i;
    		if(w[w0] <= sq) i1[w[w0]] = w0; else i2[i] = w0;
    		g[w0] = w[w0] - 1;
    	}
    	fo(j, 1, p0) for(int i = 1; p[j] * p[j] <= w[i] && i <= w0; i ++) {
    		int k = num(w[i] / p[j]);
    		g[i] -= (g[k] - (j - 1));
    	}
    	int l = 0;
    	fd(i, w0, 1) {
    		while(l < d0 && d[l + 1] <= w[i]) l ++;
    		g[i] -= l;
    	}
    	fo(i, 1, w0) g[i] = -g[i];
    	fo(i, 1, p0) sp[i] = sp[i - 1] - (k % p[i] != 0);
    	for(int j = p0; j >= 1; j --) if(k % p[j]) for(int i = 1; p[j] * p[j] <= w[i] && i <= w0; i ++) {
    		int k = num(w[i] / p[j]);
    		g[i] += (g[k] - sp[j]) * -1;
    	}
    	fo(i, 1, w0) g[i] ++;
    }
    ll q2(int x) {
    	return g[num(x)];
    }
    }
    using m26 :: q2;
    
    int n, m;
    
    struct P {
    	int x, y;
    } t[1024]; int t0;
    
    void dg(int x, int y, int z) {
    	if(x > d0) {
    		t[++ t0] = (P) {y, z};
    		return;
    	}
    	dg(x + 1, y, z);
    	dg(x + 1, y * d[x], -z);
    }
    
    int cmp(P a, P b) { return a.x < b.x;}
    
    void query() {
    	ll ans = 0;
    	int lc = 0, c = 0;
    	for(int i = 1, j; i <= min(n, m); i = j + 1) {
    		int m1 = n / (n / i), m2 = m / (m / i);
    		if(m1 < m2) {
    			j = m1, c = 0;
    		} else {
    			j = m2, c = 1;
    		}
    		ll s = !c ? q1(j) : q2(j);
    		if(i > 1) s -= !lc ? q1(i - 1) : q2(i - 1);
    		lc = c;
    		ll s2 = 0;
    		fo(u, 1, t0) if(t[u].x <= m / i) s2 += t[u].y * (m / i / t[u].x);
    		ans += s * (n / i) * s2;
    	}
    	pp("%lld", ans);
    }
    
    int main() {
    	scanf("%d %d %d", &n, &m, &k);
    	fen(k);
    	m25 :: n = n;
    	m25 :: build();
    	m26 :: n = m;
    	m26 :: build();
    	dg(1, 1, 1);
    	sort(t + 1, t + t0 + 1, cmp);
    	query();
    }
    
  • 相关阅读:
    hello
    1234566i 还是规范
    萨嘎给
    DRF 商城项目
    DRF 商城项目
    DRF 商城项目
    Windows 上连接本地 Linux虚拟机上的 mysql 数据库
    xadmin 数据添加报错: IndexError: list index out of range
    python 实现聊天室
    xadmin 组件拓展自定义使用
  • 原文地址:https://www.cnblogs.com/coldchair/p/11390716.html
Copyright © 2011-2022 走看看