Square Constraints
求 ([0,2N-1]) 的排列 (P) 数量,满足 (Nleq P_i^2+i^2 leq (2N)^2, forall 0leq i leq 2N-1) 。
(1leq N leq 250) 。
对于每一个数 (i) 可以填的位置是一个区间 ([l_i,r_i]) ,问题等价于给每个 (i) 找一个合法位置,使得数和位置一一对应,其中
注意到当 (l_i,r_i) 任取的时候问题等价于求二分图完美匹配的方案数,这是一个 NPC 问题,所以要充分利用到题目的性质。
性质1.
性质2.
注意到当 (l_i = 0) 时,方案数等价于将 (r_i) 从小到大排序后的 (prod_{i=0}^{2N-1} r_i+1-(i-1))
性质3.
也就是说对于 (i < N) 的部分,如果某个数强制违反了下界,那么一定会在性质2的排序过程中处于没有违反下界的 (i < N) 的部分前面。
有了上面三个性质,我们考虑容斥,设强制 (k) 个位置违反下界后,所有位置满足上界的方案数为 (F_k) ,根据容斥原理,恰好 (0) 个位置不违反下界,所有位置满足上界的方案数为 (sum_{i=0}^N (-1)^i F_i) 。
假设 (k) 确定了,我们求 (F_k) ,对于 (i<N) 的部分以 (l_i-1) 为第一关键字,(r_i) 为第二关键字,对于 (igeq N) 的部分,以 (r_i) 为第一关键字,(l_i-1) 为第二关键字从小到大排序。
记 (dp[i][j]) 表示前 (i) 个位置强制违反了 (j) 个下界对于性质2点式子的贡献,前 (i) 个数中 (< N) 的数的数量为 (t) , 那么有
对于每一个 (K) 都计算一遍总复杂度 (mathcal O(n^3)) 。
理解这个式子的关键在于后面贡献的计算方法,对于 (i geq N) 的部分,前面所有 (i geq N) 的部分都会提供贡献,前面所有 (i < N) 并违反下界的部分也会提供贡献。
对于 (i < N) 的部分,如果不违反下界,根据性质3,所有 (i geq N) 的部分以及违反下界的 (i < N) 的部分会提供贡献;如果违反下界,(i geq N) 在其前面的部分会提供贡献,(i < N) 的部分在其前面并违反下界的部分会提供贡献。
这题可以这样做的原因在于存在一种考虑顺序满足每个数的贡献仅与前面违反下界的数的个数有关。
code
/*program by mangoyang*/
#pragma GCC optimize("Ofast", "inline")
#include<bits/stdc++.h>
#define inf (0x3f3f3f3f)
#define Max(a, b) ((a) > (b) ? (a) : (b))
#define Min(a, b) ((a) < (b) ? (a) : (b))
typedef long long ll;
using namespace std;
template <class T>
inline void read(T &x){
int ch = 0, f = 0; x = 0;
for(; !isdigit(ch); ch = getchar()) if(ch == '-') f = 1;
for(; isdigit(ch); ch = getchar()) x = x * 10 + ch - 48;
if(f) x = -x;
}
const int N = 505;
int dp[N][N], n, mod, ans;
inline void up(int &x, int y){
x = x + y >= mod ? x + y - mod : x + y;
}
struct Node{ int id, l, r; } a[N];
inline bool cmp(Node A, Node B){
int a1 = A.id < n ? A.l - 1 : A.r;
int a2 = A.id < n ? A.r : A.l - 1;
int b1 = B.id < n ? B.l - 1: B.r;
int b2 = B.id < n ? B.r : B.l - 1;
return a1 != b1 ? a1 < b1 : a2 != b2 ? a2 < b2 : A.id > B.id;
}
inline int solve(int k){
memset(dp, 0, sizeof(dp));
int sum = 0;
dp[0][0] = 1;
for(int i = 1; i <= 2 * n; i++){
for(int j = 0; j <= k; j++){
if(a[i].id >= n){
int c = a[i].r + 1 - (i - 1 - (sum - j));
if(c > 0) up(dp[i][j], 1ll * dp[i-1][j] * c % mod);
}
else{
int c = a[i].r + 1 - (2 * n - a[i].id - 1 + k - j);
if(c > 0) up(dp[i][j], 1ll * dp[i-1][j] * c % mod);
if(j){
c = a[i].l - (i - 1 - (sum - j + 1));
if(c > 0) up(dp[i][j], 1ll * dp[i-1][j-1] * c % mod);
}
}
}
sum += a[i].id < n;
}
return dp[2*n][k];
}
int main(){
read(n), read(mod);
for(int i = 0; i < 2 * n; i++){
a[i+1].id = i;
if(i < n) a[i+1].l = ceil(sqrt(n * n - i * i));
a[i+1].r = floor(sqrt(4 * n * n - i * i));
a[i+1].r = min(a[i+1].r, 2 * n - 1);
}
sort(a + 1, a + 2 * n + 1, cmp);
for(int i = 0; i <= n; i++){
if(i & 1) up(ans, mod - solve(i));
else up(ans, solve(i));
}
cout << ans << endl;
return 0;
}