LOJ3403 「2020-2021 集训队作业」function

其实这个东西可以 OEIS 到,不过还是正经写一下吧。

考虑如何求出 P(x)P(x) ,我们对它质因数分解成 pici\prod p_i^{c_i} ,然后会发现这个东西可以对每个质因数分开考虑,最后使用 CRT 合并即可。考虑对于一个 picip_i^{c_i} ,它的阶是 φ(pici)=(pi1)pici1\varphi(p_i^{c_i}) = (p_i-1)p_i^{c_i-1} 。首先讨论 pi3p_i \ge 3 ,设原根为 ω\omega ,那么有 ω(pi1)pici1=1\omega^{(p_i-1)p_i^{c_i-1}} = 1 。讨论一下,如果 3(pi1)pici13 \not | (p_i - 1)p_i^{c_i-1} ,那么唯一的解就是 11 ,这个非常显然。否则,可以得到三个解,分别是原根的 13(pi1)pici1\frac 1 3 (p_i-1)p_i^{c_i-1}23(pi1)pici1\frac23 (p_i-1)p_i^{c_i-1}11 。也就是说,这种情况存在三个解的前提是 pi1(mod3)p_i \equiv 1 \pmod 3 或者 pi=3,ci>1p_i = 3 , c_i > 1

然后考虑 pi=2p_i = 2 ,显然解是个奇数。因此考虑设非 11 解可以表示成 y×2t+1y\times 2^t + 1 的形式,其中 yy 是奇数。于是它做三次方后是 (y2t)3+3(y2t)2+3y22t+1(y2^t)^3+3(y2^t)^2+3y^22^t+1 ,由于 3y23y^2 是个奇数,它对 2ci2^{c_i} 取模后后面两项一定会被保留,因此不可能等于 11 。于是我们证明了对于 2c2^c 的形式它一定只有一个解。

可以发现答案一定是 3k3^k 的形式。由于 ω(n)\omega(n)n2×1010n \le 2\times 10^{10} 最多只有 88 ,我们不妨筛出小于等于 nn 的答案为 30,31,,383^0,3^1,\dots ,3^8 的所有答案。分别考虑 Min_25 筛的两个部分。第一部分由于求得是质数,直接丢掉 33 ,就是求出小于等于 nnmod3=1\bmod 3 = 1 的质数的个数。这个东西看起来并不是积性的,根据 ZJK 教育的方法,设 f(x)f'(x) ,在 xmod3=1x \bmod 3 = 1 的时候为 11 ,在 xmod3=2x \bmod 3 = 2 的时候为 1-1 ,其他时候为 00 ,这个东西就具有积性,可以直接对这个东西做 Min_25 筛第一部分,再做一下 1n1 \sim n 的质数个数,这两个一起就可以算出 mod3=1\bmod 3 = 1 的质数个数。

第二部分直接按照 Min_25 的做法做即可。

复杂度 n34lognmax(ω(n))\frac{n^{\frac 3 4}}{\log n} \max(\omega(n))

#include "iostream" 
#include "algorithm"
#include "cstring"
#include "cstdio"
#include "cmath"
#include "vector"
#include "map"
#include "set"
#include "queue"
using namespace std;
#define MAXN 500006
//#define int long long
#define rep(i, a, b) for (int i = (a), i##end = (b); i <= i##end; ++i)
#define per(i, a, b) for (int i = (a), i##end = (b); i >= i##end; --i)
#define pii pair<int,int>
#define fi first
#define se second
#define mp make_pair
#define pb push_back
#define eb emplace_back
#define vi vector<int>
#define all(x) (x).begin() , (x).end()
#define mem( a ) memset( a , 0 , sizeof a )
typedef long long ll;
ll n;
int m;
ll A[MAXN];

struct poi {
	ll A[9];
	poi( ) { mem( A ); }
	poi( ll a , ll b ) {
		mem( A );
		A[0] = b , A[1] = a;
	}
	ll& operator [] ( const int r ) { return A[r]; }
	poi operator + ( poi x ) {
		poi S;
		rep( i , 0 , 8 ) S[i] = A[i] + x[i];
		return S;
	}
	poi operator - ( poi x ) {
		poi S;
		rep( i , 0 , 8 ) S[i] = A[i] - x[i];
		return S;
	}
};

int pri[MAXN] , sp[MAXN] , sp2[MAXN] , en;
int B;
void sieve( ) {
	B = sqrt( n + 0.5 );
	rep( i , 2 , B ) {
		if( !pri[i] ) pri[++ en] = i , sp[en] = sp[en - 1] + ( i % 3 == 1 ? 1 : i != 3 ? -1 : 0 ) , sp2[en] = sp2[en - 1] + ( i % 3 == 1 ? 1 : 0 );
		for( int j = 1 ; j <= en && pri[j] * i <= B ; ++ j ) {
			pri[i * pri[j]] = 1;
			if( i % pri[j] == 0 ) break;
		}
	}
}

int ff( ll x ) {
	return x % 3 == 1 ? 1 : x % 3 == 2 ? -1 : 0;
}

int cc( ll x ) {
	return ( x % 3 == 1 ? 0 : x % 3 == 2 ? -1 : -1 );
}
int A1[MAXN] , A2[MAXN];
int po( ll x ) {
	return x <= B ? A1[x] : A2[n / x];
}
ll g[MAXN] , gp[MAXN];
void getG( ) {
	int tot = 0;
	for( ll l = 1 , r ; l <= n ; l = r + 1 ) {
		r = n / ( n / l );
		++ tot;
		A[tot] = n / l;
		if( A[tot] <= B ) A1[A[tot]] = tot; else A2[n / A[tot]] = tot;
		g[tot] = cc( A[tot] ) , gp[tot] = A[tot] - 1;
	}
	rep( i , 1 , en ) 
		for( int j = 1 ; pri[i] * 1ll * pri[i] <= A[j] ; ++ j ) {
			g[j] = g[j] - ff( pri[i] ) * ( g[po( A[j] / pri[i] )] - sp[i - 1] );
			gp[j] = gp[j] - ( gp[po( A[j] / pri[i] )] - ( i - 1 ) );
		}
}

poi gg( int x ) {
	return poi( g[x] + gp[x] - 1 >> 1 , gp[x] - g[x] + 1 >> 1 );
}

bool f( int p , int a ) {
	return ( ( p == 3 && a > 1 ) || ( p % 3 == 1 ) );
}

poi S( ll n , int k ) {
	if( n <= pri[k] ) return poi();
	poi re = ( gg( po( n ) ) - poi( sp2[k] , k - sp2[k] ) );
	for( int i = k + 1 ; pri[i] * 1ll * pri[i] <= n && i <= en ; ++ i ) {
		ll cur = pri[i];
		for( int e = 1 ; cur <= n ; ++ e ) {
			poi pr = S( n / cur , i ) + poi( 0 , e != 1 );
			if( f( pri[i] , e ) ) per( k , 8 , 1 ) re[k] = ( re[k] + pr[k - 1] );
			else rep( k , 0 , 8 ) re[k] = ( re[k] + pr[k] );
			cur *= pri[i];
		}
	}
	return re;
}

void solve() {
	cin >> n >> m;
	++ m;
	int t = 0;
	while( m % 3 == 0 ) m /= 3 , ++ t;
	if( m != 1 ) return void( puts("0") );
	sieve( ) , getG( );
	cout << S( n , 0 )[t] + ( t == 0 ) << endl;
}

signed main() {
//    int T;cin >> T;while( T-- ) solve();
    solve();
}

\