文章目录
- 普通floor_sum
- 洛谷P5170 【模板】类欧几里得算法
- 万能欧几里得算法
- 求 ∑ i = 1 n A i B ⌊ a i + b c ⌋ \sum_{i=1}^{n}A^iB^{\lfloor \frac{ai+b}{c} \rfloor} ∑i=1nAiB⌊cai+b⌋
- 求 ∑ i = 0 n ⌊ a i + b c ⌋ \sum_{i=0}^n \lfloor\frac{ai+b}{c}\rfloor ∑i=0n⌊cai+b⌋(简单版)
- 求1~n除M模R的数的二进制1的数量
- 求 ∑ x = 0 n x k 1 ⌊ a x + b c ⌋ k 2 \sum_{x=0}^n x^{k_1}\lfloor \frac{ax+b}{c} \rfloor^{k_2} ∑x=0nxk1⌊cax+b⌋k2
- 求 a x + b y ≤ n ax+by \leq n ax+by≤n的对数 ( x , y ) (x,y) (x,y)
- 其他例题
- CF1098E
- CF1182F
- abc313G
- abc402G
- 2020ICPC沈阳I
- LOJ6344
类欧几里得求解问题有如下形式: f ( a , b , c , n ) = ∑ i = 0 n − 1 ⌊ a i + b c ⌋ f(a,b,c,n)=\sum_{i=0}^{n-1}\lfloor \frac{ai+b}{c} \rfloor f(a,b,c,n)=∑i=0n−1⌊cai+b⌋
普通floor_sum
若
a
≥
c
a \geq c
a≥c或者
b
≥
c
b \geq c
b≥c,令
a
′
=
a
m
o
d
c
,
b
′
=
b
m
o
d
c
a'=a \mod c,b'=b \mod c
a′=amodc,b′=bmodc
f
(
a
,
b
,
c
,
n
)
=
∑
i
=
0
n
−
1
⌊
a
c
⌋
i
+
⌊
b
c
⌋
+
⌊
a
′
i
+
b
′
c
⌋
=
n
(
n
−
1
)
2
⌊
a
c
⌋
+
n
⌊
b
c
⌋
+
f
(
a
′
,
b
′
,
c
,
n
)
f(a,b,c,n)=\sum_{i=0}^{n-1} \lfloor \frac{a}{c} \rfloor i+\lfloor \frac{b}{c} \rfloor+\lfloor \frac{a'i+b'}{c} \rfloor=\frac{n(n-1)}{2}\lfloor \frac{a}{c} \rfloor+n\lfloor \frac{b}{c} \rfloor+f(a',b',c,n)
f(a,b,c,n)=i=0∑n−1⌊ca⌋i+⌊cb⌋+⌊ca′i+b′⌋=2n(n−1)⌊ca⌋+n⌊cb⌋+f(a′,b′,c,n)
就转化成了
a
<
c
,
b
<
c
a \lt c,b \lt c
a<c,b<c的情况
考虑
a
<
c
,
b
<
c
a \lt c,b \lt c
a<c,b<c的情况:
f
(
a
,
b
,
c
,
n
)
=
∑
i
=
0
n
−
1
∑
j
=
0
⌊
a
i
+
b
c
⌋
−
1
[
1
]
f(a,b,c,n)=\sum_{i=0}^{n-1}\sum_{j=0}^{\lfloor\frac{ai+b}{c}\rfloor-1}[1]
f(a,b,c,n)=i=0∑n−1j=0∑⌊cai+b⌋−1[1]
令
m
=
⌊
a
(
n
−
1
)
+
b
c
⌋
−
1
m=\lfloor\frac{a(n-1)+b}{c} \rfloor-1
m=⌊ca(n−1)+b⌋−1
f
(
a
,
b
,
c
,
n
)
=
∑
j
=
0
m
∑
i
=
0
n
−
1
[
j
<
⌊
a
i
+
b
c
⌋
]
=
∑
j
=
0
m
−
1
∑
i
=
0
n
−
1
[
c
(
j
+
1
)
≤
a
i
+
b
]
=
∑
j
=
0
m
−
1
∑
i
=
0
n
−
1
[
i
>
⌊
c
j
+
c
−
b
−
1
a
⌋
]
=
∑
j
=
0
m
−
1
[
n
−
1
−
⌊
c
j
+
c
−
b
−
1
a
⌋
]
=
(
n
−
1
)
m
−
f
(
c
,
c
−
b
−
1
,
a
,
m
)
f(a,b,c,n)=\sum_{j=0}^m\sum_{i=0}^{n-1}[j \lt \lfloor \frac{ai+b}{c} \rfloor]\\ =\sum_{j=0}^{m-1}\sum_{i=0}^{n-1}[c(j+1) \leq ai+b]\\ =\sum_{j=0}^{m-1}\sum_{i=0}^{n-1}[i \gt \lfloor \frac{cj+c-b-1}{a}\rfloor]\\ =\sum_{j=0}^{m-1}[n-1-\lfloor \frac{cj+c-b-1}{a}\rfloor]\\ =(n-1)m-f(c,c-b-1,a,m)
f(a,b,c,n)=j=0∑mi=0∑n−1[j<⌊cai+b⌋]=j=0∑m−1i=0∑n−1[c(j+1)≤ai+b]=j=0∑m−1i=0∑n−1[i>⌊acj+c−b−1⌋]=j=0∑m−1[n−1−⌊acj+c−b−1⌋]=(n−1)m−f(c,c−b−1,a,m)
复杂度为
O
(
log
max
(
a
,
c
)
)
O(\log \max(a,c))
O(logmax(a,c))
从几何意义讲,表示的是 y = a x + b c y=\frac{ax+b}{c} y=cax+b下的整数点
类似地, g ( a , b , c , n ) = ∑ i = 0 n − 1 i ⌊ a i + b c ⌋ , h ( a , b , c , n ) = ∑ i = 0 n − 1 ⌊ a i + b c ⌋ 2 g(a,b,c,n)=\sum_{i=0}^{n-1}i\lfloor \frac{ai+b}{c} \rfloor,h(a,b,c,n)=\sum_{i=0}^{n-1}\lfloor \frac{ai+b}{c} \rfloor^2 g(a,b,c,n)=∑i=0n−1i⌊cai+b⌋,h(a,b,c,n)=∑i=0n−1⌊cai+b⌋2
a ≥ c , b ≥ c a \geq c,b \geq c a≥c,b≥c的情况:
g ( a , b , c , n ) = n ( n − 1 ) ( 2 n − 1 ) 6 ⌊ a c ⌋ + n ( n − 1 ) 2 ⌊ b c ⌋ + g ( a ′ , b ′ , c , n ) g(a,b,c,n)=\frac{n(n-1)(2n-1)}{6}\lfloor \frac{a}{c} \rfloor+\frac{n(n-1)}{2}\lfloor\frac{b}{c}\rfloor+g(a',b',c,n) g(a,b,c,n)=6n(n−1)(2n−1)⌊ca⌋+2n(n−1)⌊cb⌋+g(a′,b′,c,n)
h ( a , b , c , n ) = h ( a ′ , b ′ , c , n ) + 2 ⌊ a c ⌋ g ( a ′ , b ′ , c , n ) + 2 ⌊ b c ⌋ f ( a ′ , b ′ , c , n ) + n ( n − 1 ) ( 2 n − 1 ) 6 ⌊ a c ⌋ 2 + n ( n − 1 ) ⌊ a c ⌋ ⌊ b c ⌋ + n ⌊ b c ⌋ 2 h(a,b,c,n)=h(a',b',c,n)+2\lfloor \frac{a}{c} \rfloor g(a',b',c,n)+2\lfloor\frac{b}{c}\rfloor f(a',b',c,n)+\frac{n(n-1)(2n-1)}{6}\lfloor \frac{a}{c} \rfloor^2+n(n-1)\lfloor \frac{a}{c} \rfloor \lfloor \frac{b}{c} \rfloor+n\lfloor\frac{b}{c}\rfloor^2 h(a,b,c,n)=h(a′,b′,c,n)+2⌊ca⌋g(a′,b′,c,n)+2⌊cb⌋f(a′,b′,c,n)+6n(n−1)(2n−1)⌊ca⌋2+n(n−1)⌊ca⌋⌊cb⌋+n⌊cb⌋2
a < c , b < c a \lt c,b \lt c a<c,b<c的情况:
g ( a , b , c , n ) = 1 2 ( m n ( n − 1 ) − h ( c , c − b − 1 , a , m ) − f ( c , c − b − 1 , a , m ) ) g(a,b,c,n)=\frac{1}{2}(mn(n-1)-h(c,c-b-1,a,m)-f(c,c-b-1,a,m)) g(a,b,c,n)=21(mn(n−1)−h(c,c−b−1,a,m)−f(c,c−b−1,a,m))
h ( a , b , c , n ) = ( n − 1 ) m 2 − 2 g ( c , c − b − 1 , a , m ) − f ( c , c − b − 1 , a , m ) h(a,b,c,n)=(n-1)m^2-2g(c,c-b-1,a,m)-f(c,c-b-1,a,m) h(a,b,c,n)=(n−1)m2−2g(c,c−b−1,a,m)−f(c,c−b−1,a,m)
洛谷P5170 【模板】类欧几里得算法
struct FloorSum{
static constexpr int MOD=998244353,inv2=499122177,inv6=166374059;
struct Node{
int f=0,g=0,h=0;
};
int mo(int x){return x>=MOD?x-MOD:x;}
int sqr(int x){return 1ll*x*x%MOD;}
Node f(int a,int b,int c,int n){
Node res;
if(!a){
res.f=1ll*b/c*n%MOD;
res.g=1ll*n*(n-1)/2%MOD*(b/c)%MOD;
res.h=1ll*sqr(b/c)*n%MOD;
}
else if(a>=c||b>=c){
Node o=f(a%c,b%c,c,n);
res.f=(o.f+1ll*n*(n-1)/2%MOD*(a/c)+1ll*b/c*n)%MOD;
res.g=(o.g+1ll*a/c*n%MOD*(n-1)%MOD*(n*2-1)%MOD*inv6+1ll*n*(n-1)/2%MOD*(b/c))%MOD;
res.h=(o.h+2ll*(a/c)*o.g+2ll*(b/c)*o.f+1ll*n*(n-1)%MOD*(n*2-1)%MOD*inv6%MOD*sqr(a/c)+1ll*n*(n-1)%MOD*(a/c)%MOD*(b/c)+1ll*n*sqr(b/c))%MOD;
}
else{
int m=(1ll*a*(n-1)+b)/c;
Node o=f(c,c-b-1,a,m);
res.f=(1ll*(n-1)*m-o.f+MOD)%MOD;
res.g=((1ll*m*n%MOD*(n-1)-o.h-o.f)%MOD+MOD)*inv2%MOD;
res.h=(1ll*(n-1)*m%MOD*m-mo(o.g*2)-o.f+MOD+MOD)%MOD;
}
return res;
}
};
void solve(){
FloorSum f;
int a,b,c,n;
cin>>n>>a>>b>>c;
auto res=f.f(a,b,c,n+1);
cout<<res.f<<" "<<res.h<<" "<<res.g<<"\n";
}
万能欧几里得算法
具体内容参考这篇博客:https://www.cnblogs.com/apjifengc/p/17492207.html
下面只讲做法推导
记如果经过网格的一条横线,
y
=
c
(
c
是大于等于
0
的任意正整数
)
y=c(c是大于等于0的任意正整数)
y=c(c是大于等于0的任意正整数),就乘上矩阵
U
U
U,如果经过
x
=
c
(
c
是大于等于
0
的任意正整数
)
x=c(c是大于等于0的任意正整数)
x=c(c是大于等于0的任意正整数),乘上矩阵
R
R
R,如果同时经过则先
U
U
U再
R
R
R
例子: y = ⌊ a i + b c ⌋ y=\lfloor \frac{ai+b}{c} \rfloor y=⌊cai+b⌋,我们维护向量 [ y , ∑ y , 1 ] [y,\sum y,1] [y,∑y,1]
则经过
U
U
U时,
y
+
1
y+1
y+1,经过
R
R
R,
∑
y
+
y
\sum y+y
∑y+y
U
=
(
1
0
0
0
1
0
1
0
1
)
,
R
=
(
1
1
0
0
1
0
0
0
1
)
U=\begin{pmatrix} 1 & 0 & 0\\ 0 & 1 & 0\\ 1 & 0 & 1 \end{pmatrix}, R=\begin{pmatrix} 1 & 1 & 0\\ 0 & 1 & 0\\ 0 & 0 & 1 \end{pmatrix}
U=
101010001
,R=
100110001
模板如下,复杂度为
O
(
T
log
max
(
a
,
c
)
)
O(T \log \max(a,c))
O(Tlogmax(a,c))
template<typename T>
struct FloorSum{
T qpow(T a,int b){
T res;
while(b){
if(b&1)res=res*a;
a=a*a;
b>>=1;
}
return res;
}
T f(int a,int b,int c,int n,T U,T R){
if(!n)return T();
if(b>=c)return qpow(U,b/c)*f(a,b%c,c,n,U,R);
if(a>=c)return f(a%c,b,c,n,U,qpow(U,a/c)*R);
int m=((__int128)a*n+b)/c;
if(!m)return qpow(R,n);
return qpow(R,(c-b-1)/a)*U*f(c,(c-b-1)%a,a,m-1,R,U)*qpow(R,n-((__int128)c*m-b-1)/a);
}
};
考虑模板题,要求的是 ∑ y , ∑ y 2 , ∑ x y \sum y,\sum y^2 ,\sum xy ∑y,∑y2,∑xy,如何合并两个区间
定义
X
l
/
r
,
Y
l
/
r
X_{l/r},Y_{l/r}
Xl/r,Yl/r表示左/右区间的
c
n
t
r
,
c
n
t
u
cntr,cntu
cntr,cntu的值(经过多少个
R
,
U
R,U
R,U)
∑
y
=
∑
y
l
+
∑
(
y
r
+
Y
l
)
=
∑
y
l
+
∑
y
r
+
X
r
Y
l
∑
y
2
=
∑
y
l
2
+
∑
(
y
r
+
Y
l
)
2
=
∑
y
l
2
+
∑
y
r
2
+
2
Y
l
∑
y
r
+
X
r
Y
l
2
∑
x
y
=
∑
x
l
y
l
+
∑
(
x
r
+
X
l
)
(
y
r
+
Y
l
)
=
∑
x
l
y
l
+
∑
x
r
y
r
+
Y
l
∑
x
r
+
X
l
∑
y
r
+
X
r
X
l
Y
l
∑
x
=
∑
x
l
+
∑
(
x
r
+
X
l
)
=
∑
x
l
+
∑
x
r
+
X
l
X
r
\sum y=\sum y_l+\sum (y_r+Y_l)=\sum y_l+\sum y_r+X_rY_l\\ \sum y^2=\sum y_l^2+\sum(y_r+Y_l)^2=\sum y_l^2+\sum y_r^2+2Y_l\sum y_r+X_rY_l^2\\ \sum xy=\sum x_ly_l+\sum (x_r+X_l)(y_r+Y_l)=\sum x_ly_l+\sum x_ry_r+Y_l\sum x_r+X_l\sum y_r+X_rX_lY_l\\ \sum x=\sum x_l+\sum (x_r+X_l)=\sum x_l+\sum x_r+X_lX_r
∑y=∑yl+∑(yr+Yl)=∑yl+∑yr+XrYl∑y2=∑yl2+∑(yr+Yl)2=∑yl2+∑yr2+2Yl∑yr+XrYl2∑xy=∑xlyl+∑(xr+Xl)(yr+Yl)=∑xlyl+∑xryr+Yl∑xr+Xl∑yr+XrXlYl∑x=∑xl+∑(xr+Xl)=∑xl+∑xr+XlXr
总结:
- 左边求和的数量为 X l X_l Xl,右边求和的数量为 X r X_r Xr
- 合并时,右边的 x , y x,y x,y加上 X l , Y l X_l,Y_l Xl,Yl
- U , R U,R U,R初始化时, U U U初始化 y y y, R R R初始化 x x x以及其他和 x x x相关的答案
写出如下结构体
struct Info{
int x,y,sumx,sumy,sumxy,sumy2;
Info():x(0),y(0),sumx(0),sumy(0),sumxy(0),sumy2(0){}
Info operator * (const Info &r){
Info res;
res.x=(x+r.x)%mod;
res.y=(y+r.y)%mod;
res.sumx=(sumx+r.sumx+x*r.x%mod)%mod;
res.sumy=(sumy+r.sumy+y*r.x%mod)%mod;
res.sumy2=(sumy2+r.sumy2+2*y%mod*r.sumy%mod+r.x*y%mod*y%mod)%mod;
res.sumxy=(sumxy+r.sumxy+y*r.sumx%mod+x*r.sumy%mod+x*y%mod*r.x%mod)%mod;
return res;
}
};
初始化 U , R U,R U,R,求解的是 [ 1 , n ] [1,n] [1,n]的范围,需要加上 i = 0 i=0 i=0的答案
void solve(){
int n,a,b,c;
cin>>n>>a>>b>>c;
Info U,R;
U.y=R.x=R.sumx=1;
FloorSum<Info> f;
auto ans=f.f(a,b,c,n,U,R);
ans.sumy=(ans.sumy+b/c)%mod;
ans.sumy2=(ans.sumy2+(b/c)*(b/c))%mod;
cout<<ans.sumy<<" "<<ans.sumy2<<" "<<ans.sumxy<<"\n";
}
求 ∑ i = 1 n A i B ⌊ a i + b c ⌋ \sum_{i=1}^{n}A^iB^{\lfloor \frac{ai+b}{c} \rfloor} ∑i=1nAiB⌊cai+b⌋
维护
∑
A
x
B
y
,
A
X
,
B
Y
\sum A^xB^y,A^X,B^Y
∑AxBy,AX,BY
∑
A
x
B
y
=
∑
A
x
l
B
y
l
+
∑
A
X
l
+
x
r
B
Y
l
+
y
r
=
∑
A
x
l
B
y
l
+
A
X
l
∑
A
x
r
B
y
r
B
Y
l
\sum A^xB^y=\sum A^{x_l}B^{y_l}+\sum A^{X_l+x_r}B^{Y_l+y_r}=\sum A^{x_l}B^{y_l}+A^{X_l}\sum A^{x_r}B^{y_r}B^{Y_l}
∑AxBy=∑AxlByl+∑AXl+xrBYl+yr=∑AxlByl+AXl∑AxrByrBYl
struct Matrix{
static constexpr int maxn=22;
int m[maxn][maxn];
Matrix(){
memset(m,0,sizeof m);
}
Matrix(int n){
init(n);
}
void init(int n){
memset(m,0,sizeof m);
for(int i=0;i<n;i++)m[i][i]=1;
}
Matrix operator*(const Matrix &a){
Matrix res;
for(int k=0;k<maxn;k++){
for(int i=0;i<maxn;i++){
for(int j=0;j<maxn;j++){
res.m[i][j]=(res.m[i][j]+m[i][k]*a.m[k][j]%mod)%mod;
}
}
}
return res;
}
Matrix operator+(const Matrix &a){
Matrix res;
for(int i=0;i<maxn;i++){
for(int j=0;j<maxn;j++){
res.m[i][j]=(m[i][j]+a.m[i][j])%mod;
}
}
return res;
}
};
struct Info{
Matrix sum,x,y;
Info():sum(),x(22),y(22){}
Info operator * (const Info &r){
Info res;
res.x=x*r.x;
res.y=y*r.y;
res.sum=sum+x*r.sum*y;
return res;
}
};
template<typename T>
struct FloorSum{
T qpow(T a,int b){
T res;
while(b){
if(b&1)res=res*a;
a=a*a;
b>>=1;
}
return res;
}
T f(int a,int b,int c,int n,T U,T R){
if(!n)return T();
if(b>=c)return qpow(U,b/c)*f(a,b%c,c,n,U,R);
if(a>=c)return f(a%c,b,c,n,U,qpow(U,a/c)*R);
int m=((__int128)a*n+b)/c;
if(!m)return qpow(R,n);
return qpow(R,(c-b-1)/a)*U*f(c,(c-b-1)%a,a,m-1,R,U)*qpow(R,n-((__int128)c*m-b-1)/a);
}
};
void solve(){
int p,q,r,l,n;
cin>>p>>q>>r>>l>>n;
Info R,U;
for(int i=0;i<n;i++){
for(int j=0;j<n;j++)cin>>R.x.m[i][j],R.sum.m[i][j]=R.x.m[i][j];
}
for(int i=0;i<n;i++){
for(int j=0;j<n;j++)cin>>U.y.m[i][j];
}
FloorSum<Info> f;
auto ans=f.f(p,r,q,l,U,R);
for(int i=0;i<n;i++){
for(int j=0;j<n;j++)cout<<ans.sum.m[i][j]<<" \n"[j==n-1];
}
}
求 ∑ i = 0 n ⌊ a i + b c ⌋ \sum_{i=0}^n \lfloor\frac{ai+b}{c}\rfloor ∑i=0n⌊cai+b⌋(简单版)
int FloorSum(int a,int b,int c,int n){// \sum_{i=0}^{n} (ai+b)/c
if(!a)return b/c*(n+1);
if(a<0){
b+=n*a;
a=-a;
}
if(b<0){
int d=(a-b-1)/a;
n-=d;
b+=a*d;
if(n<0)return 0;
}
//======== 负数计入答案 ======================
// if(b<0){
// int k=(c-b-1)/c;
// return FloorSum(a,b+k*c,c,n)-k*(n+1);
// }
//===========================================
if(a>=c||b>=c)return n*(n+1)/2*(a/c)+(b/c)*(n+1)+FloorSum(a%c,b%c,c,n);
int m=(a*n+b)/c;
return m*n-FloorSum(c,c-b-1,a,m-1);
}
求1~n除M模R的数的二进制1的数量
结论:对于一个数 A A A二进制第 k k k位为 1 1 1,有 ⌊ A + 2 k 2 k + 1 ⌋ − ⌊ A 2 k + 1 ⌋ = 1 \lfloor \frac{A+2^k}{2^{k+1}}\rfloor-\lfloor \frac{A}{2^{k+1}} \rfloor=1 ⌊2k+1A+2k⌋−⌊2k+1A⌋=1,另一种写法是 ⌊ A 2 k ⌋ − 2 ⌊ A 2 k + 1 ⌋ = 1 \lfloor \frac{A}{2^k} \rfloor-2\lfloor \frac{A}{2^{k+1}} \rfloor=1 ⌊2kA⌋−2⌊2k+1A⌋=1
可以用类欧求解
求 ∑ x = 0 n x k 1 ⌊ a x + b c ⌋ k 2 \sum_{x=0}^n x^{k_1}\lfloor \frac{ax+b}{c} \rfloor^{k_2} ∑x=0nxk1⌊cax+b⌋k2
∑ x k 1 y k 2 = ∑ x l k 1 y l k 2 + ∑ ( x r + X l ) k 1 ( y r + Y l ) k 2 ∑ ( x r + X l ) k 1 ( y r + Y l ) k 2 = ∑ ∑ i = 0 k 1 ∑ j = 0 k 2 ( k 1 i ) X l i x r k 1 − i ( k 2 j ) Y l j y r k 2 − j = ∑ i = 0 k 1 ∑ j = 0 k 2 ( k 1 i ) ( k 2 j ) X l i Y l j ∑ x r k 1 − i y r k 2 − j \sum x^{k_1}y^{k_2}=\sum x_l^{k_1}y_l^{k_2}+\sum(x_r+X_l)^{k_1}(y_r+Y_l)^{k_2}\\ \sum(x_r+X_l)^{k_1}(y_r+Y_l)^{k_2}=\sum \sum_{i=0}^{k_1}\sum_{j=0}^{k_2}\binom{k_1}{i}X_l^ix_r^{k_1-i}\binom{k_2}{j}Y_l^jy_r^{k_2-j}=\sum_{i=0}^{k_1}\sum_{j=0}^{k_2}\binom{k_1}{i}\binom{k_2}{j}X_l^iY_l^j \sum x_r^{k_1-i}y_r^{k_2-j} ∑xk1yk2=∑xlk1ylk2+∑(xr+Xl)k1(yr+Yl)k2∑(xr+Xl)k1(yr+Yl)k2=∑i=0∑k1j=0∑k2(ik1)Xlixrk1−i(jk2)Yljyrk2−j=i=0∑k1j=0∑k2(ik1)(jk2)XliYlj∑xrk1−iyrk2−j
因此启发我们维护 x , y , ∑ x i y j x,y,\sum x^iy^j x,y,∑xiyj
struct Info{
int x,y;
vector<vector<int>> sum;
Info():sum(vector (k1+1,vector<int>(k2+1))),x(0),y(0){}
Info operator * (const Info &r){
Info res;
res.x=(x+r.x)%mod;
res.y=(y+r.y)%mod;
int cx=x%mod,cy=y%mod;
for(int u=0;u<=k1;u++){
for(int v=0;v<=k2;v++){
res.sum[u][v]=sum[u][v];
for(int i=0,X=1;i<=u;i++,X=X*cx%mod){
for(int j=0,Y=1;j<=v;j++,Y=Y*cy%mod){
res.sum[u][v]=(res.sum[u][v]+comb.binom(u,i)*comb.binom(v,j)%mod*X%mod*Y%mod*r.sum[u-i][v-j]%mod)%mod;
}
}
}
}
return res;
}
};
void solve(){
int n,a,b,c;
cin>>n>>a>>b>>c>>k1>>k2;
FloorSum<Info> f;
Info U,R;
U.y=R.x=1;
for(int i=0;i<=k1;i++)R.sum[i][0]=1;
comb.init(15);
auto ans=f.f(a,b,c,n,U,R);
int res=ans.sum[k1][k2];
if(k1==0)res=(res+comb.qpow(b/c,k2))%mod;//k1=0,0^0=1,答案记得加上
cout<<res<<"\n";
}
求 a x + b y ≤ n ax+by \leq n ax+by≤n的对数 ( x , y ) (x,y) (x,y)
∑ i = 0 n a ⌊ n − i a + b b ⌋ \sum_{i=0}^{\frac{n}{a}}\lfloor \frac{n-ia+b}{b} \rfloor ∑i=0an⌊bn−ia+b⌋
void solve(){
int n,a,b;
cin>>n>>a>>b;
cout<<FloorSum(-a,n+b,b,n/a)<<'\n';
}
其他例题
CF1098E
给定长度为 n n n的序列 a i a_i ai,计算出每个 gcd ( a l , . . . . a r ) \gcd(a_l,....a_r) gcd(al,....ar)放进 b b b数组,从大到小排序;再计算 ∑ b [ l , . . . r ] \sum b[l,...r] ∑b[l,...r]放进数组 c c c,然后求出数组 c c c的中位数, 1 ≤ n ≤ 5 ∗ 10 4 , 1 ≤ a i ≤ 10 5 1 \leq n \leq 5*10^4, 1\leq a_i \leq 10^5 1≤n≤5∗104,1≤ai≤105
思路:
首先计算
b
b
b数组,以
i
i
i为右端点,向左gcd最多迭代
log
\log
log次,st表预处理区间gcd,然后二分跳
log
\log
log次,处理出
b
b
b每个数字的出现次数,这部分为
O
(
n
log
3
n
)
O(n \log^3 n)
O(nlog3n)
接下来二分答案mid,然后是推式子
- 如果区间都是相同的数,计算出最多容纳的数字
mx=min(cnt[i],mid/nums[i])
,然后就是等差数列求和ans+=mx*(2*cnt[i]-mx+1)/2;
- 然后是有不同的数字的情况,假设当前区间为 [ l , r ] [l,r] [l,r],那么记 s u m = b [ l + 1 , . . . . . r − 1 ] sum=b[l+1,.....r-1] sum=b[l+1,.....r−1],记 A = n u m s l , S A = c n t l , B = n u m s r , S B = c n t r A=nums_l,SA=cnt_l,B=nums_r,SB=cnt_r A=numsl,SA=cntl,B=numsr,SB=cntr
- 如果 A ∗ S A + s u m + B ∗ S B ≤ m i d A*SA+sum+B*SB \leq mid A∗SA+sum+B∗SB≤mid那么 S A SA SA个左端点和 S B SB SB个右端点任选,答案为 S A ∗ S B SA*SB SA∗SB
- 否则答案为 ∑ i = 1 S A ∑ j = 1 S B s u m + A i + B j ≤ m i d \sum_{i=1}^{SA}\sum_{j=1}^{SB}sum+Ai+Bj \leq mid ∑i=1SA∑j=1SBsum+Ai+Bj≤mid
j
≤
m
i
d
−
s
u
m
−
A
i
B
∑
i
=
1
S
A
∑
j
=
1
S
B
s
u
m
+
A
i
+
B
j
≤
m
i
d
=
∑
i
=
1
S
A
min
(
max
(
0
,
m
i
d
−
s
u
m
−
A
i
B
)
,
S
B
)
因为
min
(
a
,
b
)
=
a
−
m
a
x
(
0
,
a
−
b
)
=
∑
i
=
1
S
A
max
(
0
,
m
i
d
−
s
u
m
−
A
i
B
)
−
∑
i
=
1
S
A
max
(
0
,
max
(
0
,
m
i
d
−
s
u
m
−
A
i
B
)
−
S
B
)
=
∑
i
=
1
S
A
max
(
0
,
m
i
d
−
s
u
m
−
A
i
B
)
−
∑
i
=
1
S
A
max
(
0
,
max
(
−
S
B
,
m
i
d
−
s
u
m
−
A
i
−
B
∗
S
B
B
)
)
=
∑
i
=
1
S
A
max
(
0
,
m
i
d
−
s
u
m
−
A
i
B
)
−
∑
i
=
1
S
A
max
(
0
,
m
i
d
−
s
u
m
−
A
i
−
B
∗
S
B
B
)
j \leq \frac{mid-sum-Ai}{B}\\ \sum_{i=1}^{SA}\sum_{j=1}^{SB}sum+Ai+Bj \leq mid=\sum_{i=1}^{SA}\min(\max(0, \frac{mid-sum-Ai}{B}),SB)\\ \text{因为}\min(a,b)=a-max(0,a-b)\\ =\sum_{i=1}^{SA}\max(0, \frac{mid-sum-Ai}{B})-\sum_{i=1}^{SA}\max(0,\max(0, \frac{mid-sum-Ai}{B})-SB)=\\ \sum_{i=1}^{SA}\max(0, \frac{mid-sum-Ai}{B})-\sum_{i=1}^{SA}\max(0,\max(-SB, \frac{mid-sum-Ai-B*SB}{B}))\\ =\sum_{i=1}^{SA}\max(0, \frac{mid-sum-Ai}{B})-\sum_{i=1}^{SA}\max(0, \frac{mid-sum-Ai-B*SB}{B})
j≤Bmid−sum−Aii=1∑SAj=1∑SBsum+Ai+Bj≤mid=i=1∑SAmin(max(0,Bmid−sum−Ai),SB)因为min(a,b)=a−max(0,a−b)=i=1∑SAmax(0,Bmid−sum−Ai)−i=1∑SAmax(0,max(0,Bmid−sum−Ai)−SB)=i=1∑SAmax(0,Bmid−sum−Ai)−i=1∑SAmax(0,max(−SB,Bmid−sum−Ai−B∗SB))=i=1∑SAmax(0,Bmid−sum−Ai)−i=1∑SAmax(0,Bmid−sum−Ai−B∗SB)
然后就可以套板子了
可以双指针维护左右端点
int FloorSum(int a,int b,int c,int n){// \sum_{i=0}^{n} (ai+b)/c
if(!a)return b/c*(n+1);
if(a<0){
b+=n*a;
a=-a;
}
if(b<0){
int d=(a-b-1)/a;
n-=d;
b+=a*d;
}
if(n<0)return 0;
if(a>=c||b>=c)return n*(n+1)/2*(a/c)+(b/c)*(n+1)+FloorSum(a%c,b%c,c,n);
int m=(a*n+b)/c;
return m*n-FloorSum(c,c-b-1,a,m-1);
}
template<typename T,typename F>
struct SpraseTable{
int n;
vector<vector<T>> st;
F op;
void init(int n,const vector<T>& a){
this->n=n;
st.assign(n+1,vector<T>(__lg(n)+1));
for(int i=1;i<=n;i++)st[i][0]=a[i];
for(int j=1;1<<j<=n;j++){
for(int i=1;i+(1<<j)-1<=n;i++){
st[i][j]=op(st[i][j-1],st[i+(1<<j-1)][j-1]);
}
}
}
T query(int l,int r){
int k=__lg(r-l+1);
return op(st[l][k],st[r-(1<<k)+1][k]);
}
};
SpraseTable<int,decltype([](int a,int b){return __gcd(a,b);})> st;
void solve(){
int n;
cin>>n;
vector<int> a(n+1);
for(int i=1;i<=n;i++)cin>>a[i];
st.init(n,a);
map<int,int> mp;
for(int i=1;i<=n;i++){
int x=i;
while(x>=1){
int l=1,r=x;
int g=st.query(x,i);
while(l<r){
int mid=l+r>>1;
if(st.query(mid,i)==g)r=mid;
else l=mid+1;
}
mp[g]+=x-l+1;
x=l-1;
}
}
vector<int> nums,cnt;
for(auto &[x,y]:mp){
nums.push_back(x);
cnt.push_back(y);
}
int C=n*(n+1)/2;
C=C*(C+1)/2;
if(C&1)C=(C+1)/2;
else C=C/2;
auto check=[&](int mid){
int ans=0;
for(int i=0;i<nums.size();i++){
int mx=min(cnt[i],mid/nums[i]);
ans+=mx*(2*cnt[i]-mx+1)/2;
}
if(nums.size()<=1)return ans>=C;
auto get=[&](int l,int r,int sum){
if(r>=nums.size())return 0ll;
int a=nums[l],b=nums[r],sa=cnt[l],sb=cnt[r];
if(a*sa+sum+b*sb<=mid)return sa*sb;
return FloorSum(-a,mid-sum,b,sa)-max(0ll,(mid-sum)/b)-FloorSum(-a,mid-sum-b*sb,b,sa-1)+max(0ll,(mid-sum-b*sb)/b);
};
for(int i=0,j=1,k=0,sum=0,tot=0;i<nums.size();i++){
if(nums[i]>mid)break;
if(i>=j)j=i+1,sum=tot=0;
ans+=tot*cnt[i];
while(j<nums.size()&&sum+cnt[j]*nums[j]<=mid){
ans+=get(i,j,sum);
sum+=cnt[j]*nums[j];
tot+=cnt[j];
j++;
}
if(j<nums.size())ans+=get(i,j,sum);
if(i+1<j){
sum-=cnt[i+1]*nums[i+1];
tot-=cnt[i+1];
}
}
return ans>=C;
};
int l=1,r=1e15;
while(l<r){
int mid=l+r>>1;
if(check(mid))r=mid;
else l=mid+1;
}
cout<<l<<"\n";
}
CF1182F
题意:
求
[
a
,
b
]
[a,b]
[a,b]内
f
(
x
)
=
∣
sin
p
q
π
x
∣
f(x)=|\sin \frac{p}{q}\pi x|
f(x)=∣sinqpπx∣的最大值对应的最小整数
x
x
x,
0
≤
a
,
b
≤
10
9
,
1
≤
p
,
q
≤
10
9
0 \leq a,b \leq 10^9,1 \leq p,q \leq 10^9
0≤a,b≤109,1≤p,q≤109
思路:
- p x m o d q q \frac{px \mod q}{q} qpxmodq接近 1 2 \frac{1}{2} 21更优,二分这个距离,即 2 p x m o d 2 q ∈ [ q − m i d , q + m i d ] 2px \mod 2q \in [q-mid,q+mid] 2pxmod2q∈[q−mid,q+mid]
- 结论: [ p x m o d q ∈ [ l , r ] ] = [ ⌊ p x − l q ⌋ − ⌊ p x − r − 1 q ⌋ > 0 ] [px \mod q \in [l,r]]=[\lfloor \frac{px-l}{q} \rfloor-\lfloor \frac{px-r-1}{q} \rfloor \gt 0] [pxmodq∈[l,r]]=[⌊qpx−l⌋−⌊qpx−r−1⌋>0]
- 因此只需求和判断是否存在即可
- 求出最小距离后,用exgcd还原求解
int FloorSum(int a,int b,int c,int n){// \sum_{i=0}^{n} (ai+b)/c
if(!a)return b/c*(n+1);
if(a<0){
b+=n*a;
a=-a;
}
// if(b<0){
// int d=(a-b-1)/a;
// n-=d;
// b+=a*d;
// if(n<0)return 0;
// }
//======== 负数计入答案 ======================
if(b<0){
int k=(c-b-1)/c;
return FloorSum(a,b+k*c,c,n)-k*(n+1);
}
//===========================================
if(a>=c||b>=c)return n*(n+1)/2*(a/c)+(b/c)*(n+1)+FloorSum(a%c,b%c,c,n);
int m=(a*n+b)/c;
return m*n-FloorSum(c,c-b-1,a,m-1);
}
int exgcd(int a,int b,int &x,int &y){
if(!b){
x=1,y=0;
return a;
}
int g=exgcd(b,a%b,y,x);
y-=a/b*x;
return g;
}
void solve(){
int a,b,p,q;
cin>>a>>b>>p>>q;
int P=p<<1,Q=q<<1;
auto check=[&](int mid){
int l=q-mid,r=q+mid;
return FloorSum(P,-l,Q,b)-FloorSum(P,-l,Q,a-1)-FloorSum(P,-r-1,Q,b)+FloorSum(P,-r-1,Q,a-1)>0;
};
int l=0,r=1e9;
while(l<r){
int mid=l+r>>1;
if(check(mid))r=mid;
else l=mid+1;
}
int x,y;
int g=exgcd(P,Q,x,y);
int ans=1e9;
int dx=Q/g,dy=P/g;
if((q-l)%g==0){
int val=(q-l)/g*x;
val+=(a-val)/dx*dx;
while(val>=a)val-=dx;
while(val<a)val+=dx;
ans=min(ans,val);
}
if((q+l)%g==0){
int val=(q+l)/g*x;
val+=(a-val)/dx*dx;
while(val>=a)val-=dx;
while(val<a)val+=dx;
ans=min(ans,val);
}
cout<<ans<<"\n";
}
abc313G
abc402G
2020ICPC沈阳I
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
∣
2
π
(
i
M
+
j
)
H
M
−
2
π
j
M
∣
≤
2
π
A
H
M
]
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
∣
i
M
+
j
−
H
j
∣
≤
A
]
H
M
−
(
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
M
+
j
−
H
j
>
A
]
+
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
M
+
j
−
H
j
<
−
A
]
)
\sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[|\frac{2 \pi(iM+j)}{HM}-\frac{2 \pi j}{M}| \leq \frac{2 \pi A}{HM}]\\ \sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[|iM+j-Hj| \leq A]\\ HM-(\sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[iM+j-Hj\gt A]+\sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[iM+j-Hj\lt -A])
i=0∑H−1j=0∑M−1[∣HM2π(iM+j)−M2πj∣≤HM2πA]i=0∑H−1j=0∑M−1[∣iM+j−Hj∣≤A]HM−(i=0∑H−1j=0∑M−1[iM+j−Hj>A]+i=0∑H−1j=0∑M−1[iM+j−Hj<−A])
对于前半部分:
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
M
+
j
−
H
j
>
A
]
=
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
>
⌊
A
+
(
H
−
1
)
j
M
⌋
]
=
∑
j
=
0
M
−
1
H
−
1
−
⌊
A
+
(
H
−
1
)
j
M
⌋
=
M
(
H
−
1
)
−
∑
j
=
0
M
−
1
⌊
A
+
(
H
−
1
)
j
M
⌋
\sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[iM+j-Hj\gt A]=\\ \sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[i \gt \lfloor\frac{A+(H-1)j}{M}\rfloor]=\\ \sum_{j=0}^{M-1}H-1-\lfloor\frac{A+(H-1)j}{M}\rfloor=\\ M(H-1)-\sum_{j=0}^{M-1}\lfloor\frac{A+(H-1)j}{M}\rfloor
i=0∑H−1j=0∑M−1[iM+j−Hj>A]=i=0∑H−1j=0∑M−1[i>⌊MA+(H−1)j⌋]=j=0∑M−1H−1−⌊MA+(H−1)j⌋=M(H−1)−j=0∑M−1⌊MA+(H−1)j⌋
对于后半部分:
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
M
+
j
−
H
j
<
−
A
]
=
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
<
⌊
(
H
−
1
)
j
−
A
M
⌋
]
=
∑
i
=
0
H
−
1
∑
j
=
0
M
−
1
[
i
≤
⌊
(
H
−
1
)
j
−
A
−
1
M
⌋
]
=
∑
j
=
0
M
−
1
⌊
(
H
−
1
)
j
−
A
−
1
M
⌋
+
1
=
∑
j
=
0
M
−
1
⌊
(
H
−
1
)
j
−
A
−
1
+
M
M
⌋
\sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[iM+j-Hj\lt -A]=\\ \sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[i \lt \lfloor\frac{(H-1)j-A}{M}\rfloor]=\\ \sum_{i=0}^{H-1}\sum_{j=0}^{M-1}[i \leq \lfloor\frac{(H-1)j-A-1}{M}\rfloor]=\\ \sum_{j=0}^{M-1} \lfloor\frac{(H-1)j-A-1}{M}\rfloor+1=\\ \sum_{j=0}^{M-1} \lfloor\frac{(H-1)j-A-1+M}{M}\rfloor
i=0∑H−1j=0∑M−1[iM+j−Hj<−A]=i=0∑H−1j=0∑M−1[i<⌊M(H−1)j−A⌋]=i=0∑H−1j=0∑M−1[i≤⌊M(H−1)j−A−1⌋]=j=0∑M−1⌊M(H−1)j−A−1⌋+1=j=0∑M−1⌊M(H−1)j−A−1+M⌋
似乎要计算负数部分和特判HM
int FloorSum(int a,int b,int c,int n){// \sum_{i=0}^{n} (ai+b)/c
if(!a)return b/c*(n+1);
if(a<0){
b+=n*a;
a=-a;
}
// if(b<0){
// int d=(a-b-1)/a;
// n-=d;
// b+=a*d;
// if(n<0)return 0;
// }
//======== 负数计入答案 ======================
if(b<0){
int k=(c-b-1)/c;
return FloorSum(a,b+k*c,c,n)-k*(n+1);
}
//===========================================
if(a>=c||b>=c)return n*(n+1)/2*(a/c)+(b/c)*(n+1)+FloorSum(a%c,b%c,c,n);
int m=(a*n+b)/c;
return m*n-FloorSum(c,c-b-1,a,m-1);
}
void solve(){
int H,M,A;
cin>>H>>M>>A;
int ans=H*M;
ans-=M*(H-1)-FloorSum(H-1,A,M,M-1);
ans-=FloorSum(H-1,M-A-1,M,M-1);
cout<<min(H*M,ans)<<"\n";
}
补一发易理解的正解,设时间为 T T T,时针角速度 ω 1 = 2 π H M \omega_1=\frac{2\pi}{HM} ω1=HM2π,分针角速度 ω 2 = 2 π M \omega_2=\frac{2\pi}{M} ω2=M2π
因此 θ = T 2 π ( H − 1 ) H M m o d 2 π \theta=T\frac{2 \pi (H-1)}{HM} \mod 2\pi θ=THM2π(H−1)mod2π, min ( θ , 2 π − θ ) ≤ 2 π A H M \min(\theta,2\pi-\theta) \leq \frac{2 \pi A}{HM} min(θ,2π−θ)≤HM2πA
答案就是 [ T ( H − 1 ) m o d H M ≤ A ] + [ T ( H − 1 ) m o d H M ≥ H M − A ] [T(H-1) \mod HM \leq A] + [T(H-1) \mod HM \geq HM-A] [T(H−1)modHM≤A]+[T(H−1)modHM≥HM−A]
当 A = H M 2 A=\frac{HM}{2} A=2HM时会有重合,此时相当于角度小于等于 π \pi π,则任意时刻都满足,答案为 H M HM HM
否则来看一下式子 T ( H − 1 ) m o d H M ≤ A T(H-1) \mod HM \leq A T(H−1)modHM≤A,令 G = gcd ( H − 1 , H M ) G=\gcd(H-1,HM) G=gcd(H−1,HM)
如果 G = 1 G=1 G=1,因为 T ∈ [ 0 , H M − 1 ] T \in [0,HM-1] T∈[0,HM−1]是 H M HM HM的完全剩余系,所以 T ( H − 1 ) T(H-1) T(H−1)也是 H M HM HM的完全剩余系,答案就是 A + 1 A+1 A+1
否则,考虑 T H − 1 G m o d H M G ≤ ⌊ A G ⌋ T\frac{H-1}{G} \mod \frac{HM}{G} \leq \lfloor \frac{A}{G} \rfloor TGH−1modGHM≤⌊GA⌋,变成上面的情况,答案就是 ⌊ A G ⌋ + 1 \lfloor \frac{A}{G} \rfloor+1 ⌊GA⌋+1,相当于将 [ 0 , H M − 1 ] [0,HM-1] [0,HM−1]切分成 G G G段,因此总答案为 G ( ⌊ A G ⌋ + 1 ) G(\lfloor \frac{A}{G} \rfloor+1) G(⌊GA⌋+1)
考虑 T ( H − 1 ) m o d H M ≥ H M − A T(H-1) \mod HM \geq HM-A T(H−1)modHM≥HM−A
T H − 1 G m o d H M G ≥ ⌈ H M − A G ⌉ T \frac{H-1}{G} \mod \frac{HM}{G} \geq \lceil \frac{HM-A}{G} \rceil TGH−1modGHM≥⌈GHM−A⌉,答案就是 G ( H M − 1 G − ⌊ H M − A + G − 1 G ⌋ + 1 ) G(\frac{HM-1}{G}-\lfloor \frac{HM-A+G-1}{G} \rfloor+1) G(GHM−1−⌊GHM−A+G−1⌋+1)
void solve(){
int H,M,A;
cin>>H>>M>>A;
if(2*A==H*M){
cout<<H*M<<"\n";
return;
}
int G=__gcd(H-1,H*M);
cout<<G*(A/G+1+(H*M-1)/G-(H*M-A+G-1)/G+1)<<"\n";
}
LOJ6344
按位考虑一个数第
k
k
k位为
1
1
1,等价于
⌊
x
2
k
⌋
−
2
⌊
x
2
k
+
1
⌋
\lfloor \frac{x}{2^k} \rfloor -2\lfloor \frac{x}{2^{k+1}} \rfloor
⌊2kx⌋−2⌊2k+1x⌋
定义 a n s i ans_i ansi表示 ∑ i = 1 n ⌊ n − i ⌊ n i ⌋ 2 k ⌋ \sum_{i=1}^n \lfloor \frac{n-i\lfloor \frac{n}{i} \rfloor}{2^k}\rfloor ∑i=1n⌊2kn−i⌊in⌋⌋
对于小于等于 n \sqrt n n的数暴力做,这些块比较密集;大于的用整除分块+类欧
int FloorSum(int a, int b, int c, int n) { // \sum_{i=0}^{n} (ai+b)/c
if (!a)
return b / c * (n + 1);
if (a < 0) {
b += n * a;
a = -a;
}
if (b < 0) {
int d = (a - b - 1) / a;
n -= d;
b += a * d;
if (n < 0)
return 0;
}
//======== 负数计入答案 ======================
// if(b<0){
// int k=(c-b-1)/c;
// return FloorSum(a,b+k*c,c,n)-k*(n+1);
// }
//===========================================
if (a >= c || b >= c)
return n * (n + 1) / 2 * (a / c) + (b / c) * (n + 1) + FloorSum(a % c, b % c, c, n);
int m = (a * n + b) / c;
return m * n - FloorSum(c, c - b - 1, a, m - 1);
}
void solve() {
int n;
cin >> n;
int c = 0;
for (int i = 1; i <= min(n, (int)1e6); i++)
c ^= n % i;
vector<int> ans(50);
for (int l = 1e6 + 1, r; l <= n; l = r + 1) {
int v = n / l;
r = n / (n / l);
for (int i = 0; i < 40; i++) {
ans[i] += FloorSum(-v, n - l * v, 1ll << i, r - l);
}
}
int res = 0;
for (int i = 0; i < 39; i++) {
res += ((ans[i] - 2 * ans[i + 1]) & 1) * (1ll << i);
}
cout << (c ^ res) << "\n";
}