「EC Final 2019」狄利克雷 k 次根
题意:给出一个函数的前 n n n项 g ( 1 ) , g ( 2 ) . . . . g ( n ) g(1),g(2)....g(n) g(1),g(2)....g(n),求它在狄利克雷卷积意义下的 k k k次根 f f f的前 n n n项,这里的狄利克雷卷积是指 h ( n ) = ∑ d ∣ n f ( d ) g ( n d ) m o d p h(n) = \sum_{d|n} f(d)g(\frac nd)\bmod p h(n)=∑d∣nf(d)g(dn)modp。
n
≤
1
e
5
n \leq 1e5
n≤1e5做法:
容易发现
f
p
=
ϵ
f^p = \epsilon
fp=ϵ,也就是狄利克雷卷积的单位元。
但是证明并不显然
证明:
f
p
f^p
fp的组合含义是,我们选择了
p
p
p个数他们的乘积等于
n
n
n,这
p
p
p个数假设为
a
1
,
a
2
.
.
.
a
p
a_1,a_2...a_p
a1,a2...ap,则
f
p
(
n
)
=
∑
a
1
+
a
2
+
.
.
a
p
=
n
∏
i
=
1
p
f
(
a
i
)
f^p(n) = \sum_{a_1+a_2+..a_p=n} \prod_{i=1}^p f(a_i)
fp(n)=∑a1+a2+..ap=n∏i=1pf(ai)
考虑一个
a
1
,
a
2
.
.
.
a
p
a_1,a_2...a_p
a1,a2...ap的方案,我们可以把两个不同的
a
i
,
a
j
a_i,a_j
ai,aj值交换,仍然可以得到一个
∏
i
=
1
p
f
(
a
i
)
\prod_{i=1}^p f(a_i)
∏i=1pf(ai)相等的方案,假设说
a
1
,
a
2
.
.
.
a
p
a_1,a_2...a_p
a1,a2...ap有
k
k
k种不同的数,每种有
b
1
,
b
2
.
.
b
k
b_1,b_2..b_k
b1,b2..bk个,那么本质和
a
1
,
a
2
.
.
.
a
p
a_1,a_2...a_p
a1,a2...ap相同的方案数就有
p
!
∏
i
=
1
k
b
i
!
\frac {p!}{\prod_{i=1}^k b_i!}
∏i=1kbi!p!种,
根据库默尔定理,当
p
p
p是质数时
p
!
∏
i
=
1
k
b
i
!
m
o
d
p
≢
0
\frac {p!}{\prod_{i=1}^k b_i!} \bmod p \not \equiv 0
∏i=1kbi!p!modp≡0
当且仅当在
p
p
p进制下将
b
i
b_i
bi全部加起来的时候不发生进位。
但是因为
b
i
b_i
bi的和刚好为
p
p
p,所以不发生进位的情况只有一种也就是只有一个
b
b
b,
b
1
=
p
b_1 = p
b1=p。
也就是在选择的
p
p
p个数完全相同为
a
1
a_1
a1时对
a
1
p
a_1^p
a1p有
f
(
a
1
)
p
f(a_1)^p
f(a1)p的贡献,其他时候都为
0
0
0。
可以推出
f
p
(
1
)
=
1
f^p(1)=1
fp(1)=1,也就是选择
p
p
p个
1
1
1,至于第二小的情况选择
p
p
p个
2
2
2,在这题里面因为
p
=
998244353
p = 998244353
p=998244353,所以对于
1
<
x
<
2
998244353
1\lt x\lt 2^{998244353}
1<x<2998244353的所有
x
x
x都有
f
p
(
x
)
m
o
d
p
=
0
f^p(x) \bmod p = 0
fp(x)modp=0。
所以对于
n
≥
2
p
n \geq 2^p
n≥2p的情况可以说我们这个式子是不正确的。
所以求出
k
m
o
d
p
k \bmod p
kmodp的逆元
b
b
b后即可直接
g
b
=
f
k
b
=
f
a
p
+
1
=
f
g^b = f^{kb} = f^{ap+1} = f
gb=fkb=fap+1=f
狄利克雷卷积快速幂即可。
C
o
d
e
Code
Code
#include<bits/stdc++.h>
#define maxn 1000005
#define mod 998244353
using namespace std;
int n,K,a[maxn],b[maxn];
int Pow(int b,int k){ int r=1;for(;k;k>>=1,b=1ll*b*b%mod) if(k&1) r=1ll*r*b%mod; return r; }
void mul(int *a,int *b,int *c){
static int d[maxn];
for(int i=1;i<=n;i++) d[i] = 0;
for(int i=1;i<=n;i++) for(int j=1,k=i;k<=n;k+=i,j++)
d[k] = (d[k] + 1ll * a[i] * b[j]) % mod;
for(int i=1;i<=n;i++) c[i] = d[i];
}
int main(){
scanf("%d%d",&n,&K);
for(int i=1;i<=n;i++) scanf("%d",&a[i]);
b[1] = 1;K = Pow(K , mod-2);
for(;K;K>>=1,mul(a,a,a)) if(K&1) mul(a,b,b);
for(int i=1;i<=n;i++) printf("%d%c",b[i]," \n"[i==n]);
}
n
≤
1
e
6
做
法
:
n \leq 1e6做法:
n≤1e6做法:
考虑狄利克雷生成函数
F
(
x
)
=
∑
i
=
1
n
f
(
i
)
i
x
F(x) = \sum_{i=1}^n f(i)i^x
F(x)=∑i=1nf(i)ix
狄利克雷卷积就是他们的生成函数相乘。
对于
g
g
g的生成函数求
ln
\ln
ln再
÷
k
\div k
÷k再
exp
\exp
exp即可。
有
ln
F
(
x
)
=
C
+
∫
F
(
x
)
′
F
(
x
)
\ln F(x) = C + \int \frac {F(x)'}{F(x)}
lnF(x)=C+∫F(x)F(x)′
其中
C
C
C是一个与
x
x
x无关的常数,是因为对于常数
a
a
a显然
ln
a
≠
∫
a
′
a
\ln a \neq \int \frac {a'}a
lna=∫aa′,所以需要一个
C
=
ln
a
C=\ln a
C=lna来补全。
但是这个题它的常数是
f
(
1
)
1
x
f(1)1^x
f(1)1x,这就很为难了,你不能不说它是常数,但是他的确看上去和
x
x
x有关。
就这样吧,
C
=
ln
f
(
1
)
1
x
=
0
C = \ln f(1)1^x = 0
C=lnf(1)1x=0
F
(
x
)
′
=
∑
i
=
1
n
f
(
i
)
i
x
ln
i
F(x)' = \sum_{i=1}^n f(i)i^x\ln i
F(x)′=∑i=1nf(i)ixlni
同时可以得到
∫
F
(
x
)
=
∑
i
=
2
n
f
(
i
)
i
x
ln
i
+
C
\int F(x) = \frac {\sum_{i=2}^n f(i)i^x}{\ln i} + C
∫F(x)=lni∑i=2nf(i)ix+C
exp
F
(
x
)
′
=
exp
F
(
x
)
×
F
(
x
)
′
\exp F(x)' = \exp F(x) \times F(x)'
expF(x)′=expF(x)×F(x)′
exp
F
(
x
)
=
∫
exp
F
(
x
)
×
F
(
x
)
′
\exp F(x) = \int \exp F(x) \times F(x)'
expF(x)=∫expF(x)×F(x)′
因为
F
(
x
)
′
F(x)'
F(x)′不存在
1
x
1^x
1x项所以可以递推出
exp
F
(
x
)
\exp F(x)
expF(x)。
最大的问题是如何处理
ln
i
\ln i
lni。
需要重新定义导数和
ln
\ln
ln,而积分和
exp
\exp
exp则分别是他们的逆运算。
新的定义是一种线性变换满足链式法则(
ln
F
(
x
)
=
∫
F
(
x
)
′
F
(
x
)
\ln F(x) = \int \frac {F(x)'}{F(x)}
lnF(x)=∫F(x)F(x)′)
对于
ln
\ln
ln需要满足
exp
(
ln
k
F
(
x
)
)
=
F
(
x
)
k
\exp(\ln kF(x)) = F(x)^k
exp(lnkF(x))=F(x)k,所以需要满足
ln
a
b
=
ln
a
+
ln
b
\ln ab = \ln a + \ln b
lnab=lna+lnb,
ln
1
=
0
\ln 1 = 0
ln1=0
对于此题令
ln
a
\ln a
lna为
a
a
a的质因子个数,
F
(
x
)
′
=
∑
i
=
1
n
f
(
i
)
i
x
ln
i
F(x)' = \sum_{i=1}^n f(i)i^x\ln i
F(x)′=∑i=1nf(i)ixlni即可。
A C C o d e \mathcal AC \ Code AC Code
#include<bits/stdc++.h>
#define maxn 1000005
#define mod 998244353
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
using namespace std;
int n,K,a[maxn],b[maxn];
int pr[maxn],vis[maxn],sg[maxn],cnt_pr;
int inv[maxn];
void ln(int *a){
static int b[maxn];
rep(i,1,n) b[i] = 1ll * a[i] * sg[i] % mod;
rep(i,1,n){
for(int j=2,k=2*i;k<=n;k+=i,j++)
b[k] = (b[k] - 1ll * b[i] * a[j]) % mod;
b[i] = 1ll * b[i] * inv[sg[i]] % mod;
}
rep(i,1,n) a[i] = b[i];
}
void exp(int *a){
static int b[maxn];
rep(i,1,n) b[i] = 1ll * a[i] * sg[i] % mod , a[i] = 0;
a[1] = 1;
rep(i,1,n){
if(i>1) a[i] = 1ll * a[i] * inv[sg[i]] % mod;
for(int j=2,k=2*i;k<=n;k+=i,j++)
a[k] = (a[k] + 1ll * a[i] * b[j]) % mod;
}
}
int Pow(int b,int k){ int r=1;for(;k;k>>=1,b=1ll*b*b%mod) if(k&1) r=1ll*r*b%mod;return r; }
int main(){
scanf("%d%d",&n,&K);
for(int i=1;i<=n;i++) scanf("%d",&a[i]);
inv[0] = inv[1] = 1;
rep(i,2,n){
if(!vis[i]) pr[cnt_pr++] = i , sg[i] = 1;
for(int j=0;pr[j] * i <= n;j++){
vis[i * pr[j]] = 1;
sg[i * pr[j]] = sg[i] + 1;
if(i % pr[j] == 0) break;
}
}
rep(i,2,20) inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod;
ln(a);K = Pow(K , mod-2);
rep(i,1,n) a[i] = 1ll * a[i] * K % mod;
exp(a);
rep(i,1,n) printf("%d%c",a[i]," \n"[i==n]);
}
2018 ACM-ICPC World Finals H .Single Cut of Failure
给出一个 w × h w\times h w×h的矩形,并给出若干条端点在矩形边上且不和矩形的边重合的线段,给出最少用几条端点在矩形边上的线段可以切割所有线段的方案。
发现每条线段都会被主对角线和副对角线之一切割,所以答案最多是
2
2
2。
将第
i
i
i条线段的两个端点都标号为
i
i
i,对于所有端点极角排序后求是否有一个长度为
n
n
n的连续子序列中没有相同标号的点即可。
A C C o d e \mathcal AC \ Code AC Code
#include<bits/stdc++.h>
#define maxn 2000006
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
#define db double
using namespace std;
int n,w,h;
int c[maxn],cnt[maxn],A[maxn],B[maxn];
db ang[maxn];
bool cmp(const int &u,const int &v){ return ang[u] < ang[v]; }
int main(){
scanf("%d%d%d",&n,&w,&h);
rep(i,1,n){
int a,b,c,d;
scanf("%d%d%d%d",&a,&b,&c,&d);
::c[i*2-1] = i*2-1 , ::c[i*2] = i*2;
ang[i*2-1] = atan2(b-h/2.0,a-w/2.0);
A[i*2-1] = a , B[i*2-1] = b;
ang[i*2] = atan2(d-h/2.0,c-w/2.0);
A[i*2] = c , B[i*2] = d;
}
sort(c+1,c+1+2*n,cmp);
c[0] = c[2*n] , c[2*n+1] = c[1];
int j=1;
rep(i,1,2*n){
for(;cnt[(c[i]+1)/2] == 1;)
cnt[(c[j++]+1)/2] = 0;
cnt[(c[i]+1)/2] = 1;
if(i - j + 1 == n){
printf("%d\n",1);
if(A[c[i]] == w) printf("%lf %lf ",(db)A[c[i]],B[c[i]]+0.5);
else if(B[c[i]] == h) printf("%lf %lf ",A[c[i]]-0.5,(db)B[c[i]]);
else if(A[c[i]] == 0) printf("%lf %lf ",(db)A[c[i]],B[c[i]]-0.5);
else printf("%lf %lf ",A[c[i]]+0.5,(db)B[c[i]]);
i = j-1;
if(A[c[i]] == w) printf("%lf %lf\n",(db)A[c[i]],B[c[i]]+0.5);
else if(B[c[i]] == h) printf("%lf %lf\n",A[c[i]]-0.5,(db)B[c[i]]);
else if(A[c[i]] == 0) printf("%lf %lf\n",(db)A[c[i]],B[c[i]]-0.5);
else printf("%lf %lf\n",A[c[i]]+0.5,(db)B[c[i]]);
return 0;
}
}
printf("%d\n",2);
printf("0.5 0 %lf %lf\n",w-0.5,(db)h);
printf("%lf %lf %lf %lf\n",0.5,(db)h,w-0.5,0.0);
}
2018 ACM-ICPC World Finals I. Triangles
给出一个图,数等边三角形个数。
维护往左,往左上,往右上的最长能走多远,
然后对于每一排横着扫,用BIT在
(
x
,
y
)
(x,y)
(x,y)处
+
1
+1
+1,移动到
(
x
+
D
,
y
)
(x+D,y)
(x+D,y)处
−
1
-1
−1,其中
D
D
D是
(
x
,
y
)
(x,y)
(x,y)往右上能走的最长距离,以此求出以
(
a
,
b
)
(a,b)
(a,b)为右下角的三角形个数。
注意还需要处理以
(
a
,
b
)
(a,b)
(a,b)为左上角的三角形个数。
A C C o d e \mathcal AC \ Code AC Code
#include<bits/stdc++.h>
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
#define LL long long
#define pb push_back
using namespace std;
int r,c;
char s[6005][12005];
int U[6005][12005],D[6005][12005];
int tr[24005];
vector<int>G[24005];
void upd(int u,int v){ for(;u<=4*c;u+=u&-u) tr[u] += v; }
int qry(int u){ int r=0;for(;u;u-=u&-u) r+=tr[u];return r; }
int main(){
scanf("%d%d\n",&r,&c);
rep(i,1,2*r-1)
fgets(s[i]+1,2*c+2,stdin);
LL ans = 0;
queue<int>q;
for(int i=1;i<=2*r-1;i+=2){
for(int j=(i/2&1?3:1);j<=2*c-1;j+=4){
if(s[i][j-1] == ' '){
for(;!q.empty();q.pop()){
int x = q.front();
if(x + D[i][x] >= j)
upd(2 * c - x , -1) , G[x+D[i][x]].clear();
}
}
if(s[i-1][j-1] != ' ' && s[i-1][j-1])
U[i][j] = U[i-2][j-2] + 1;
else
U[i][j] = 0;
if(s[i-1][j+1] != ' ' && s[i-1][j+1])
D[i][j] = D[i-2][j+2] + 4;
else
D[i][j] = 0;
ans += qry(2*c-(j-4*U[i][j]));
q.push(j);
upd(2*c-j,1);
G[j+D[i][j]].pb(j);
for(;!G[j].empty();)
upd(2*c-G[j].back(),-1),G[j].pop_back();
}
for(;!q.empty();q.pop()){
int x = q.front();
if(x + D[i][x] > 2*c-1)
upd(2 * c - x , -1) , G[x+D[i][x]].clear();
}
}
for(int i=2*r-1;i>=1;i-=2){
for(int j=(i/2&1?3:1);j<=2*c-1;j+=4){
if(s[i][j-1] == ' '){
for(;!q.empty();q.pop()){
int x = q.front();
if(x + D[i][x] >= j)
upd(2 * c - x , -1) , G[x+D[i][x]].clear();
}
}
if(s[i+1][j-1] != ' ' && s[i+1][j-1])
U[i][j] = U[i+2][j-2] + 1;
else
U[i][j] = 0;
if(s[i+1][j+1] != ' ' && s[i+1][j+1])
D[i][j] = D[i+2][j+2] + 4;
else
D[i][j] = 0;
ans += qry(2*c-(j-4*U[i][j]));
q.push(j);
upd(2*c-j,1);
G[j+D[i][j]].pb(j);
for(;!G[j].empty();)
upd(2*c-G[j].back(),-1),G[j].pop_back();
}
for(;!q.empty();q.pop()){
int x = q.front();
if(x + D[i][x] > 2*c-1)
upd(2 * c - x , -1) , G[x+D[i][x]].clear();
}
}
printf("%lld\n",ans);
}
「ICPC World Finals 2018」征服世界
一棵树上给出一些军队的位置和每个位置需要的军队数,求移动军队的距离之和的最小值。
发现这是个裸的费用流问题,这个费用流的实质是最小权完备匹配,于是在树上从下至上枚举
l
c
a
lca
lca,然后两个点
x
,
y
x,y
x,y的匹配就是
d
e
p
x
+
d
e
p
y
−
2
∗
d
e
p
l
c
a
dep_x+dep_y-2*dep_{lca}
depx+depy−2∗deplca,因为需要的军队必须被满足,所以我们把被需要的点的权值减去
∞
\infty
∞,最后加上即可,这样我们直接做最小权匹配得到的答案也一定是完备匹配。
那么每次就是用堆维护寻找
v
a
l
x
+
v
a
l
y
−
2
∗
d
e
p
l
c
a
<
0
val_x+val_y-2*dep_{lca}\lt 0
valx+valy−2∗deplca<0的点对
x
,
y
x,y
x,y,然后将他们匹配,因为是求最小权匹配,只有费用为负我们才需要考虑加入。
考虑费用流会有将一个点换一个点匹配的操作,那么当
x
,
y
x,y
x,y匹配时我们将
x
x
x的权值改为放弃
y
y
y的收益即
2
∗
d
e
p
l
c
a
−
v
a
l
y
2*dep_{lca}-val_y
2∗deplca−valy再放入堆中即可。
如果之前已经在
z
z
z处匹配了的一对点
x
,
y
x,y
x,y在
z
z
z的祖先处同时换了点匹配,那么
z
z
z的父边会被经过两次一定不优,也就是说每次的增广路长度是
O
(
1
)
O(1)
O(1)的,增广
O
(
n
)
O(n)
O(n)次所以总的堆操作次数也是
O
(
n
)
O(n)
O(n)的,如果是在
z
z
z处同时换了点匹配,可以贪心的证明这种情况不会存在。
A
C
C
o
d
e
\mathcal AC\ Code
AC Code
#include<bits/stdc++.h>
#define maxn 250005
#define maxp 4000006
#define LL long long
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
#define inf ((LL)(1e12))
using namespace std;
char cb[1<<16],*cs=cb,*ct=cb;
#define getc() (cs==ct&&(ct=(cs=cb)+fread(cb,1,1<<16,stdin),cs==ct)?0:*cs++)
void read(int &res){
char ch;bool f=0;
for(;!isdigit(ch=getc());) if(ch=='-') f=1;
for(res=ch-'0';isdigit(ch=getc());res=res*10+ch-'0');
(f) && (res=-res);
}
int n;
int info[maxn],Prev[maxn<<1],to[maxn<<1],cst[maxn<<1],cnt_e;
void Node(int u,int v,int c){ Prev[++cnt_e]=info[u],info[u]=cnt_e,to[cnt_e]=v,cst[cnt_e]=c; }
LL dep[maxn],v[maxp];
int rt[maxn][2],x[maxn],y[maxn],h[maxp],lc[maxp],rc[maxp],tot;
void merge(int &u,int l,int r){
if(!l || !r) return (void)(u = l+r);
if(v[l] < v[r]) u = l , merge(rc[u] , rc[l] , r);
else u = r , merge(rc[u] , l , rc[r]);
if(h[lc[u]] < h[rc[u]]) swap(lc[u],rc[u]);
h[u] = h[rc[u]] + 1;
}
LL ans;
void dfs(int u,int ff){
for(;x[u] > y[u];x[u]--)
v[++tot] = dep[u] , merge(rt[u][0],rt[u][0],tot);
for(;x[u] < y[u];x[u]++)
v[++tot] = dep[u] - inf , merge(rt[u][1],rt[u][1],tot);
for(int i=info[u],v;i;i=Prev[i]) if((v=to[i])^ff){
dep[v] = dep[u] + cst[i];
dfs(v,u);
merge(rt[u][0],rt[u][0],rt[v][0]);
merge(rt[u][1],rt[u][1],rt[v][1]);
}
for(;rt[u][0] && rt[u][1] && v[rt[u][0]] + v[rt[u][1]] - dep[u] * 2 < 0;){
ans += v[rt[u][0]] + v[rt[u][1]] - dep[u] * 2;
int a = rt[u][0] , b = rt[u][1];
LL va = v[a] , vb = v[b];
merge(rt[u][0],lc[rt[u][0]],rc[rt[u][0]]);
merge(rt[u][1],lc[rt[u][1]],rc[rt[u][1]]);
v[a] = dep[u] * 2 - vb;
lc[a] = rc[a] = 0;
merge(rt[u][0],rt[u][0],a);
v[b] = dep[u] * 2 - va;
lc[b] = rc[b] = 0;
merge(rt[u][1],rt[u][1],b);
}
}
int main(){
read(n);
rep(i,1,n-1){
int u,v,w;read(u),read(v),read(w);
Node(u,v,w),Node(v,u,w);
}
int sm = 0;
rep(i,1,n){
read(x[i]),read(y[i]);
if(x[i] < y[i]) sm += y[i] - x[i];
}
dfs(1,0);
printf("%lld\n",ans+inf*sm);
}
「ICPC World Finals 2018」绿宝石之岛
一开始 n n n个人每人都有一颗宝石,在 d d d个夜晚中每个晚上都会等概率选出一个宝石将其变成两个,求宝石最多的 r r r个人的宝石总数的期望值。
假设最后第
i
i
i个人多拿了
a
i
a_i
ai个宝石。
那么在这多出来的
∑
i
a
i
=
d
\sum_{i}a_i = d
∑iai=d个宝石中,按顺序分配给这
i
i
i个人的方案数为
d
!
∏
a
i
!
\frac {d!}{\prod a_i!}
∏ai!d!。
但是对于第
i
i
i个人多拿到了他的第
k
k
k个宝石的时候,有
k
k
k个宝石可以用来分裂,也就是说方案要乘上一个
∏
a
i
!
\prod {a_i!}
∏ai!
所以对于任何一个合法的
{
a
i
}
\{a_i\}
{ai}方案数都为
d
!
d!
d!。
所以答案就是所有合法的
{
a
i
}
\{a_i\}
{ai}的前
r
r
r大之和
÷
\div
÷合法的
{
a
i
}
\{a_i\}
{ai}方案数。
设前者为
s
u
m
i
,
j
sum_{i,j}
sumi,j表示多出来的前
i
i
i个宝石分给
j
j
j个人的所有方案中前
r
r
r大之和,
后者为
c
n
t
i
,
j
cnt_{i,j}
cnti,j表示多出来的前
i
i
i个宝石分给
j
j
j个人的方案数。
而对于
s
u
m
sum
sum的计算有一个巧妙的
d
p
dp
dp,那就是我们每次选一些人把他们的宝石数
+
1
+1
+1,后一次选出的人一定要包含所有的前一次选出的人,可以发现一个宝石分配方案和一个操作方案一一对应。
这样转化的话选出的人数如果大于
r
r
r,则我们会给前
r
r
r大之和
+
r
+r
+r,否则选出的人都是前
r
r
r大的。
对于
c
n
t
cnt
cnt的
d
p
dp
dp可能有很多但是就照着
s
u
m
sum
sum的
d
p
dp
dp写比较省事。
A C C o d e \mathcal AC \ Code AC Code
#include<bits/stdc++.h>
#define maxn 505
#define db double
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
using namespace std;
db C[maxn][maxn],f[maxn][maxn],c[maxn][maxn];
int n,d,r;
int main(){
scanf("%d%d%d",&n,&d,&r);
rep(i,C[0][0]=1,n) rep(j,C[i][0]=1,i)
C[i][j] = C[i-1][j-1] + C[i-1][j];
rep(i,1,n) c[0][i] = 1;
rep(i,1,d) rep(j,1,n) rep(k,1,min(i,j))
c[i][j] += C[j][k] * c[i-k][k],
f[i][j] += C[j][k] * (f[i-k][k] + min(k,r) * c[i-k][k]);
printf("%.10lf\n",f[d][n] / c[d][n] + (db)r);
}
「ICPC World Finals 2018」跳过罪恶
二维网格图,图中每一个格子都是一幢楼,楼都有高度,每次可以从一个楼顶部的中心跳到另外一个楼顶部的中心,轨迹是个斜抛运动, g = 9.80665 m / s g = 9.80665 m/s g=9.80665m/s,初始速度需要固定都为 v v v,需要中间不会碰到别的楼才可以跳,给出起点楼,求每个楼最多蹦几次可以到。
连边即可跑
B
F
S
BFS
BFS,对于两幢楼
a
,
b
a,b
a,b。
设他们的高度之差为
h
h
h,水平距离为
d
d
d。
那么设水平速度为
x
x
x,竖直速度为
y
y
y(初始时的),滞空时长为
t
t
t。
有
x
2
+
y
2
=
v
2
x^2+y^2 = v^2
x2+y2=v2
y
t
−
g
t
2
2
=
h
yt-\frac {gt^2}2 = h
yt−2gt2=h
x
t
=
d
xt = d
xt=d
对于第一个式子
×
t
2
\times t^2
×t2然后把别的式子化为
x
t
,
y
t
xt,yt
xt,yt的形式带入到这个式子中可以得到:
g
2
4
t
4
+
(
h
g
−
v
2
)
t
2
+
d
2
+
h
2
=
0
\frac {g^2}4t^4 + (hg-v^2)t^2+d^2+h^2=0
4g2t4+(hg−v2)t2+d2+h2=0
然后就可以解出
t
t
t,注意我们需要飞的尽量高,这样才更能避开中间的楼,所以需要
t
t
t最大的解。
然后在楼房交界处计算一下到这里时的高度和楼房高度比较一下即可。
A C C o d e \mathcal AC \ Code AC Code
#include<bits/stdc++.h>
#define maxn 22
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
#define db double
#define mp make_pair
using namespace std;
int dx,dy,lx,ly;
const db G = 9.80665;
db w,v,h[maxn][maxn];
db sqr(db a){ return a * a; }
int c[maxn][maxn][maxn][maxn];
int dis[maxn][maxn];
int main(){
scanf("%d%d%lf%lf%d%d",&dx,&dy,&w,&v,&lx,&ly);
rep(i,1,dy) rep(j,1,dx) scanf("%lf",&h[i][j]);
rep(i,1,dy) rep(j,1,dx) rep(a,1,dy) rep(b,1,dx) if(a!=i || b!=j){
db A = sqr(G) / 4 , B = G * (h[a][b] - h[i][j]) - sqr(v) , C = sqr(h[a][b] - h[i][j]) + sqr((i-a) * w) + sqr((j-b) * w);
if(B*B-4*A*C >= 0){
db t = (-B + sqrt(B*B - 4*A*C)) / (2 * A);
if(t >= 0){
t = sqrt(t);
db d = sqrt(sqr((i-a) * w) + sqr((j-b) * w)) , xv = sqrt(sqr(d / t) / sqr(d) * sqr((j-b) * w)) , yv = sqrt(sqr(d / t) / sqr(d) * sqr((i-a) * w));
if(a < i) yv = -yv;
if(b < j) xv = -xv;
db H = h[a][b] - h[i][j] , zv = (H + G * sqr(t) / 2) / t;
bool flg = 0;
rep(k,min(j,b),max(j,b)-1){
t = ((k-j)*w+0.5*w) / xv;
db ny = (i-1)*w + 0.5*w + t * yv;
int p = ceil(ny / w);
if(max(h[p][k] , h[p][k+1]) > zv * t - G * sqr(t) / 2 + h[i][j]){
flg = 1;
break;
}
if(ny / w == p){
if(max(h[p+1][k] , h[p+1][k+1]) > zv * t - G * sqr(t) / 2 + h[i][j]){
flg = 1;
break;
}
}
}
rep(k,min(i,a),max(i,a)-1){
t = ((k-i)*w+0.5*w) / yv;
db nx = (j-1)*w + 0.5*w + t * xv;
int p = ceil(nx / w);
if(max(h[k][p] , h[k+1][p]) > zv * t - G * sqr(t) / 2 + h[i][j]){
flg = 1;
break;
}
if(nx / w == p){
if(max(h[k][p+1] , h[k+1][p+1]) > zv * t - G * sqr(t) / 2 + h[i][j]){
flg = 1;
break;
}
}
}
if(!flg)
c[i][j][a][b] = 1;
}
}
}
queue<pair<int,int> >q;
memset(dis,-1,sizeof dis);
q.push(mp(ly,lx)),dis[ly][lx] = 0;
for(int x,y;!q.empty();){
x = q.front().first , y = q.front().second , q.pop();
rep(i,1,dy) rep(j,1,dx) if(c[x][y][i][j] && dis[i][j] == -1)
dis[i][j] = dis[x][y] + 1 , q.push(mp(i,j));
}
rep(i,1,dy) rep(j,1,dx){
if(dis[i][j] >= 0)
printf("%d%c",dis[i][j]," \n"[j==dx]);
else
printf("X%c"," \n"[j==dx]);
}
}
「ICPC World Finals 2018」熊猫保护区
给出一个多边形,求多边形内部的所有点中和距离其最近的顶点之间的最远距离。
对于每个点 x x x,其他所有点和 x x x的中垂线的一边离 x x x最近,那么将这些半平面求交再求交中在多边形中的部分离 x x x最远的点的距离,对于所有的 x x x取 max \max max即可。(这就是 O ( n 2 log n ) O(n^2\log n) O(n2logn)求 voronoi \texttt{voronoi} voronoi图的方法。)
如何求半平面交在非凸多边形内离一点最远的距离?
此题可以暴力对于半平面交的每条边和多边形的每条边求交点,因为
voronoi
\texttt{voronoi}
voronoi图是
Delaunay
\texttt{Delaunay}
Delaunay三角剖分的对偶图,点数是
O
(
n
)
O(n)
O(n),边数也是
O
(
n
)
O(n)
O(n)的。
总复杂度是
O
(
n
2
log
n
+
n
2
)
O(n^2\log n + n^2)
O(n2logn+n2)的。
A
C
C
o
d
e
\mathcal AC \ Code
AC Code
#include<bits/stdc++.h>
#define maxn 2055
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
#define eps 1e-10
#define db double
#define Ct const
#define Pt Point
using namespace std;
int n;
db sqr(db a){ return a * a; }
struct Pt{
db x,y;
Pt(Ct db &x=0,Ct db &y=0):x(x),y(y){}
Pt operator +(Ct Pt &B)Ct{ return Pt(x+B.x,y+B.y); }
Pt operator -(Ct Pt &B)Ct{ return Pt(x-B.x,y-B.y); }
db operator *(Ct Pt &B)Ct{ return x * B.y - y * B.x; }
Pt operator *(Ct db &B)Ct{ return Pt(x * B , y * B); }
db dist(Ct Pt &B)Ct{ return sqrt(sqr(B.y-y) + sqr(B.x-x)); }
}P[maxn],qp[maxn];
int dcmp(db a){ return a < -eps ? -1 : a > eps ? 1 : 0; }
struct Ln{
Pt S,T;
db ang;
Ln(Ct Pt &S=0,Ct Pt &T=0):S(S),T(T){ang = atan2(T.y-S.y,T.x-S.x);}
bool operator <(Ct Ln &B)Ct{ return dcmp(ang-B.ang) ? ang < B.ang : (B.T-S) * (T-S) < 0; }
}L[maxn],q[maxn];
int cnt = 0;
db ans = 0;
Pt Ipt(Ct Pt &p1,Ct Pt &v1,Ct Pt &p2,Ct Pt &v2){
return p1 + (v1 * (((p2 - p1)*v2) / (v1*v2)));
}
bool chk1(Ct Pt &u){
static db ag = 0.1234 , c = cos(ag) , s = sin(ag);
Pt v = u + Pt(c * 1e7 , s * 1e7);
int cnt = 0;
rep(i,1,n){
if(u.dist(P[i]) + u.dist(P[i%n+1]) == P[i].dist(P[i%n+1])) return 1;
if(((P[i%n+1] - u) * (v-u)) * ((P[i] - u) * (v-u)) <= 0
&& ((u - P[i]) * (P[i%n+1]-P[i])) * ((v - P[i]) * (P[i%n+1] - P[i])) <= 0)
cnt ^= 1;
}
return cnt;
}
int main(){
scanf("%d",&n);
rep(i,1,n) scanf("%lf%lf",&P[i].x,&P[i].y);
rep(i,1,n){
cnt = 0;
L[++cnt] = Ln(Pt(1e4,-1e4),Pt(-1e4,-1e4));
L[++cnt] = Ln(Pt(1e4,1e4),Pt(1e4,-1e4));
L[++cnt] = Ln(Pt(-1e4,1e4),Pt(1e4,1e4));
L[++cnt] = Ln(Pt(-1e4,-1e4),Pt(-1e4,1e4));
rep(j,1,n) if(i^j){
Pt A = (P[i] + P[j]) * 0.5 , B = (P[j] - P[i]);
swap(B.x,B.y),B.y=-B.y;
L[++cnt] = Ln(A , A+B);
}
sort(L+1,L+1+cnt);
int ql=0,qr=-1;
rep(j,1,cnt) if(j == 1 || dcmp(L[j].ang - L[j-1].ang)){
for(;qr-ql>=1 && (qp[qr-1] - L[j].S) * (L[j].T - L[j].S) < 0;qr--);
for(;qr-ql>=1 && (qp[ql] - L[j].S) * (L[j].T - L[j].S) < 0;ql++);
q[++qr] = L[j];
if(qr-ql>=1) qp[qr-1] = Ipt(q[qr].S,q[qr].T-q[qr].S,q[qr-1].S,q[qr-1].T-q[qr-1].S);
}
for(;qr-ql>=2 && (qp[qr-1] - q[ql].S) * (q[ql].T - q[ql].S) < 0;qr--);
qp[qr] = Ipt(q[qr].S,q[qr].T-q[qr].S,q[ql].S,q[ql].T-q[ql].S);
qp[qr+1] = qp[ql];
rep(j,ql,qr){
if(chk1(qp[j]))
ans = max(ans , qp[j].dist(P[i]));
rep(k,1,n)
if(
((qp[j]-P[k])*(P[k%n+1]-P[k])) * ((qp[j+1]-P[k])*(P[k%n+1]-P[k])) < 0
&& ((P[k]-qp[j])*(qp[j+1]-qp[j])) * ((P[k%n+1]-qp[j])*(qp[j+1]-qp[j])) < 0
){
ans = max(ans , P[i].dist(Ipt(qp[j],qp[j+1]-qp[j],P[k],P[k%n+1]-P[k])));
}
}
}
printf("%.10lf\n",ans);
}
2019 ICPC Asia-East Continent Final B.black and white
给出一个 n × m n\times m n×m的网格,你从 ( 0 , 0 ) (0,0) (0,0)走到 ( n , m ) (n,m) (n,m),每次只能让横纵坐标之一加一,这个网格是黑白相间染色的, ( 0 , 0 ) (0,0) (0,0)为左下角的格子是白色,求你走出的路径的左边的白色格子数 − - −黑色格子数等于 k k k的方案数。
统计横着走一步后加上上面的白色格子数减去上面的黑色格子数,
那么每次加上的数只有
0
0
0和
±
1
\pm1
±1。
为了解决
±
1
\pm 1
±1我们一次走两步(这需要
2
∣
n
+
m
2|n+m
2∣n+m,如不满足只需要枚举最后一步怎么走即可),那么如果这两步方向相同,则答案不变(竖着走则没有贡献,横着走如果上面的格子数是偶数则没有贡献,否则是一列白比黑多一,另一列少一,刚好抵消。)
我们把方向相同的步数删去之后假设每次都是第一步横着第二步竖着,一共走了
p
p
p组,因为
(
0
,
0
)
(0,0)
(0,0)一定是白色,所以这样走的贡献是
⌊
p
+
1
2
⌋
\lfloor \frac {p+1}2 \rfloor
⌊2p+1⌋,可以通过画图讨论得到。
那么把
p
p
p组中一些(
q
q
q个)换成第一步竖着第二步横着,那么实际上就是在边界上把原来在轮廓左侧的格子放到了右侧,可以发现这些格子都是白色的,所以实际上的总和就是
⌊
p
+
q
+
1
2
⌋
−
q
\lfloor \frac {p+q+1}2 \rfloor-q
⌊2p+q+1⌋−q,所以我们枚举
p
+
q
p+q
p+q即可算出
q
q
q然后组合数一下即可。
A C C o d e \mathcal AC \ Code AC Code
#include<bits/stdc++.h>
#define maxn 100005
#define mod 998244353
#define rep(i,j,k) for(int i=(j),LIM=(k);i<=LIM;i++)
#define per(i,j,k) for(int i=(j),LIM=(k);i>=LIM;i--)
using namespace std;
int fac[maxn],inv[maxn],invf[maxn];
int f(int n,int m,int K){
int ret = 0;
rep(p,0,min(n,m)){
int q = (p+1)/2-K;
if((n-p) % 2 == 0 && (m-p) % 2 == 0 && n >= p && m >= p && p >= q && q >= 0)
ret = (ret + fac[(n+m)/2] * 1ll * invf[p-q] % mod * invf[q] % mod * invf[(n-p)/2] % mod * invf[(m-p)/2]) % mod;
}
return ret;
}
int main(){
fac[0] = fac[1] = inv[0] = inv[1] = invf[0] = invf[1] = 1;
rep(i,2,maxn-1)
fac[i] = 1ll * fac[i-1] * i % mod,
inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod,
invf[i] = 1ll * invf[i-1] * inv[i] % mod;
int T;
scanf("%d",&T);
for(int n,m,K;T--;){
scanf("%d%d%d",&n,&m,&K);
if(n+m&1) printf("%d\n",(f(n,m-1,K) + f(n-1,m,K+(m&1))) % mod);
else printf("%d\n",f(n,m,K));
}
}