1、 加特技
水题，树上背包dp。
#include
#include
#include
using namespace std;
int w[105],v[105],d[105],fa[105],f[105][505];
int node[105],next[105],a[105],vis[105];
int n,m,tot,i,j,k,ans;
void add(int x,int y)
{ node[++tot]=y; next[tot]=a[x]; a[x]=tot; }
void dp(int x){
if (w[x]>m) return;
f[x][w[x]]=v[x];
for (int i=a[x];i;i=next[i]){
int y = node[i];
dp(y);
for (int j=m;j>0;j--)
for (int k=j;k>=0;k--)
f[x][j]=max(f[x][j],f[x][j-k]+f[y][k]);
}
}
int main(){
freopen("duang.in","r",stdin);
freopen("duang.out","w",stdout);
scanf("%d%d",&n,&m);
for (i=1;i<=n;i++) scanf("%d",w+i);
for (i=1;i<=n;i++) scanf("%d",v+i);
for (i=1;i<=n;i++) scanf("%d",d+i);
for (i=0;i<=n;i++) fa[i]=i;
for (i=1;i<=n;i++)
if (!vis[i]){
for (j=i;j>0&&!vis[j];vis[j]=i,j=d[j]);
if (j==0||vis[j]
2、花
可以推出第一步答案为(n!)^2的约数个数
第二步其实要求一个long long范围内的质因数分解，Pollard_rho解决。
#include
#include
#include
#include
using namespace std;
typedef long long LL;
const int Mod=1000000007;
const int Maxn=1000000;
const int PRIME[9]={2,3,5,7,11,13,17,19,23};
int prime[Maxn],check[Maxn];
int n,i,j,x,tot,t,cnt,s[Maxn],dig[20];
LL ans,maxx;
void init(){
for (i=2;i
=Maxn) break; check[i*prime[j]]=j; if (i%prime[j]==0) break; } } } LL mult(LL a,LL b,LL Mod){ LL ret = 0; for (;b>0;b>>=1){ if (b&1) ret = (ret+a)%Mod; a = (a+a)%Mod; } return ret; } LL qck(LL a,LL b,LL Mod){ LL ret = 1; for (;b>0;b>>=1){ if (b&1) ret = mult(ret,a,Mod); a = mult(a,a,Mod); } return ret; } bool WITNESS(LL a,LL n){ LL t = 0, x = n-1; while (!(x&1)) x>>=1, t++; LL d = qck(a,x,n), dd; while (t--){ dd = mult(d,d,n); if (dd==1 && d!=1 && d!=n-1) return 1; d = dd; } return (d!=1); } bool Miller_Rabin(LL n){ if (n<=1) return 0; if (n==2) return 1; if (!(n&1)) return 0; for (i=0;i<9&&PRIME[i]
=n) p = pollard_rho(n,rand()%(n-1)+1); calc(p); calc(n/p); } #define sqr(a) ((a)*(a)) int main(){ freopen("flower.in","r",stdin); freopen("flower.out","w",stdout); scanf("%d%d",&n,&x); init(); for (i=1;i<=n;i++){ for (j=i;j>1;j/=prime[check[j]]) s[ check[j] ]+=2; } for (i=1,ans=1;i<=tot;i++) ans=ans*(LL)(s[i]+1)%Mod; if (x==0) dig[ cnt=1 ] = 0; else for (i=x;i>0;i/=10) dig[ ++cnt ] = i%10; for (i=cnt;i>0;i--) ans=ans*10+(LL)dig[i]; maxx=0; calc(ans); printf("%lld\n",maxx); return 0; }
3、连电线
任意重连一条树边，求修改后树直径的最小值和修改方案。
比较繁琐的树dp。大概就是想办法搞出子树里的直径，子树外的直径，那么割去子树连向外的边，添边的最优方案就是加在两条直径的中间。
开了20来个数组。。。。
#include
#include
#include
using namespace std;
const int Maxn=300005;
int dw1[Maxn],dw2[Maxn],dw3[Maxn],e1[Maxn],e2[Maxn],e3[Maxn];
int fa[Maxn],q[Maxn],node[Maxn<<1],next[Maxn<<1],a[Maxn],up[Maxn];
int lp[Maxn],lp1[Maxn],lp2[Maxn],le1[Maxn],le2[Maxn],lpu[Maxn],dep[Maxn];
int n,i,x,y,tot,l,r,lef,rig,ans,ans1,ans2,ans3,ans4,maxx,maxp;
void add(int x,int y){
node[++tot]=y; next[tot]=a[x]; a[x]=tot;
node[++tot]=x; next[tot]=a[y]; a[y]=tot;
}
void dp(){
for (q[l=r=1]=1;l<=r;l++)
for (i=a[q[l]];i;i=next[i])
if (fa[q[l]]!=node[i])
fa[ q[++r]=node[i] ] = q[l];
for (;r>0;r--){
x = q[r]; y = fa[x];
lp[x]=max(dw1[x]+dw2[x],lp1[x]);
if (y==0) continue;
if (dw1[x]+1>dw1[y]){
dw3[y] = dw2[y];
e3[y] = e2[y];
dw2[y] = dw1[y];
e2[y] = e1[y];
dw1[y] = dw1[x]+1;
e1[y] = x;
} else
if (dw1[x]+1>dw2[y]){
dw3[y] = dw2[y];
e3[y] = e2[y];
dw2[y] = dw1[x]+1;
e2[y] = x;
} else
if (dw1[x]+1>dw3[y]){
dw3[y] = dw1[x]+1;
e3[y] = x;
}
if (lp[x]>lp1[y]){
lp2[y] = lp1[y];
le2[y] = le1[y];
lp1[y] = lp[x];
le1[y] = x;
} else
if (lp[x]>lp2[y]){
lp2[y] = lp[x];
le2[y] = x;
}
}
ans=1000000000;
for (r=1;r<=n;r++){
x = q[r];
for (i=a[x];i;i=next[i])
if (node[i]!=fa[x]){
y = node[i];
lef = lp[y];
if (le1[x]==y) rig = max( lp2[x], lpu[x] );
else rig = max( lp1[x], lpu[x] );
if (e1[x]==y)
rig = max( rig, dw2[x]+max(dw3[x],up[x]) );
else if (e2[x]==y)
rig = max( rig, dw1[x]+max(dw3[x],up[x]) );
else
rig = max( rig, dw1[x]+max(dw2[x],up[x]) );
int tmp = max( max(lef,rig), (lef+1)/2+(rig+1)/2+1 );
if (ans>tmp){
ans = tmp;
ans1 = x; ans2 = y;
}
lpu[y] = rig;
if (e1[x]==y) up[y] = max(dw2[x],up[x])+1;
else up[y] = max(dw1[x],up[x])+1;
}
}
}
int calc(int X,int Y){
maxx=-1;
memset(fa,0,sizeof(fa));
for (q[l=r=1]=X,dep[X]=0;l<=r;l++){
if (maxx
4、拆机房
比较经典的网络流，二分时间time，把所有合法电脑连向门，原点连向电脑，门连向汇点容量为time（其实这里利用了贪心的思想：如果距离为k的电脑从这扇门出去，那么路径上至少有k个电脑会选择这扇门），跑最大流判断是不是所有电脑都有所归属。
trick：考试的时候没有想到贪心策略，无脑地把每个时刻的每个门暴力拆点居然过了。。。
#include
#include
#include
using namespace std;
const int dx[4]={1,-1,0,0};
const int dy[4]={0,0,1,-1};
const int INF=(1e9);
char s[30][30];
int p[30][30],d[405][80],q[85*400];
int a[85*400],c[85*400],js[85*400];
int node[85*8000],flow[85*8000],next[85*8000];
int n,m,i,j,k,cnt,num,ans,fl,minx,l,r,L,R,Mid,tot,N;
bool flag;
void bfs(){
memset(d,127/2,sizeof(d));
for (i=1;i<=n;i++)
for (j=1;j<=m;j++)
if (s[i][j]=='.')
p[i][j]=++num;
for (i=1;i<=n;i++)
for (j=1;j<=m;j++)
if (s[i][j]=='D'){
d[0][++cnt]=0;
for (q[l=r=1]=i*100+j;l<=r;l++){
int x = q[l]/100, y = q[l]%100;
for (k=0;k<4;k++)
if (s[x+dx[k]][y+dy[k]]=='.'){
if (d[ p[x+dx[k]][y+dy[k]] ][cnt] <= d[ p[x][y] ][cnt]+1) continue;
q[++r] = (x+dx[k])*100+(y+dy[k]);
d[ p[x+dx[k]][y+dy[k]] ][cnt] = d[ p[x][y] ][cnt]+1;
}
}
}
}
void add(int x,int y,int z){
node[++tot]=y; next[tot]=a[x]; a[x]=tot; flow[tot]=z;
node[++tot]=x; next[tot]=a[y]; a[y]=tot; flow[tot]=0;
}
void init(){
memset(js,0,sizeof(js));
memset(c,10,sizeof(c));
for (q[l=r=1]=N,js[c[N]=0]=1;l<=r;l++)
for (i=a[q[l]];i;i=next[i])
if ((i&1) && c[node[i]]>c[q[l]]+1)
js[ c[ q[++r] = node[i] ] = c[q[l]]+1 ]++;
}
void sap(int x){
if (x==N){ fl+=minx; flag=1; return; }
int minn=minx, e=N-1, i;
for (i=a[x];i;i=next[i])
if (flow[i]>0){
if (c[node[i]]+1==c[x]){
minx = min(minx,flow[i]);
sap(node[i]);
if (c[1]>=N) return;
if (flag) break; else minx=minn;
}
e = min(e, c[node[i]]);
}
if (flag) {flow[i]-=minx;flow[i^1]+=minx;return;}
if ( !(--js[c[x]]) ) c[1]=N;
js[ c[x]=e+1 ]++;
}
bool Judge(){
N = Mid*cnt+num+1+1;
memset(a,0,sizeof(a));
for (i=1,tot=1;i<=num;i++){
add(1,i+1,1);
for (j=1;j<=cnt;j++)
if (d[i][j]<=Mid)
add(i+1,num+1+(j-1)*Mid+d[i][j],1);
}
for (i=1;i<=cnt;i++){
for (j=1;j
>1;
if (Judge()) R=Mid-1, ans=Mid;
else L = Mid+1;
}
if (ans==INF) puts("impossible");
else printf("%d\n",ans);
return 0;
}
5、字符串
暴力好题。。。。
linux下最好不要用gets读入，后果自负
#include
#include
#include
using namespace std;
const int MaxP=(1e9)+7;
const int Maxn=2000005;
char T[Maxn],S[30];
int n,m,N,i,j,tmp,ans;
int len[30],g[30],H[Maxn],pow[Maxn];
bool f[Maxn];
#define gethash(i,j) (H[j]-pow[(j)-(i)+1]*H[(i)-1])
int main(){
freopen("string.in","r",stdin);
freopen("string.out","w",stdout);
scanf("%d%d\n",&n,&m);
for (i=1,pow[0]=1;i<=2000000;i++)
pow[i] = pow[i-1]*MaxP;
for (i=1;i<=n;i++){
//scanf("%s",S+1);
gets(S+1);
len[i]=strlen(S+1);
for (j=1;j<=len[i];j++)
g[i]=g[i]*MaxP+S[j];
}
while (m--){
//scanf("%s",T+1);
gets(T+1);
N=strlen(T+1);
for (i=1;i<=N;i++)
H[i]=H[i-1]*MaxP+T[i];
//memset(f,0,sizeof(f));
f[0]=1;
for (i=1,ans=0,tmp=0;i<=N;i++){
f[i]=0;
for (j=1;j<=n;j++)
if (len[j]<=i && f[i-len[j]] && gethash(i-len[j]+1,i)==g[j])
{ f[i] = 1; break; }
if (f[i]) ans = i, tmp=0;
else tmp++;
if (tmp>10) break;
}
printf("%d\n",ans);
}
return 0;
}
6、string
预处理出在询问某些位后仍然不能用的字符串（比较trick的方法），由于这张图是可拓扑的，所以直接推就可以了。
#include
#include
#include
using namespace std;
#define two(x) (1LL<<(LL)(x))
typedef unsigned UN;
typedef long long LL;
typedef unsigned long long UL;
int n,m,i,j,k,tmp,cnt,p1,p2,g[55][55];
LL f[1<<20];
char S[55][55];
double ans,opt[1<<20];
int calc(LL d){
int ret = 0;
for (;d;d&=(d-1)) ret++;
return ret;
}
int main(){
freopen("string.in","r",stdin);
freopen("string.out","w",stdout);
while (~scanf("%d\n",&n)){
for (i=0;i
0;i--){
for (j=0;j
7、block
假设这个山头向右可以积水，那么在它右边一定有不比他低的山头。预处理出向右的第一个不低于他的山头，这之间的积水很容易用前缀和求出。这样这就成了一个森林的样子。用倍增优化查询。向左的可以同样处理。
#include
#include
#include
using namespace std;
const int Mod=(1e9)+7;
const int Maxn=300005, INF=2*(1e9);
typedef long long LL;
int Case,n,m,i,j,l,r,h[Maxn],s[Maxn],q[Maxn],ans;
int lf[Maxn][20],rf[Maxn][20],lw[Maxn][20],rw[Maxn][20];
void init(){
h[0]=INF;
q[l=r=1]=0;
for (i=1;i<=n;i++){
while (l<=r && h[q[r]]
0;i--){
while (l
'9') ch=getchar(); while (ch>='0' && ch<='9') {ret=ret*10+ch-'0'; ch=getchar();} return ret; } void work(){ ans=0; while (m--){ //scanf("%d%d",&l,&r); l=read(); r=read(); for (i=19;i>=0;i--){ if (rf[l][i]==0 || rf[l][i]>r) continue; ans = (ans+rw[l][i])%Mod; l = rf[l][i]; } for (i=19;i>=0;i--){ if (lf[r][i]==0 || lf[r][i]
8、tree
发现我们只需要知道这个点的颜色和子树里染色的个数奇偶性，用f[i][0/1][0/1]表示染色i号点为0/1色，子树下染色0偶/1奇数个节点。
不考虑重复，很容易求出正序和逆序的染色方案数。
现在考虑如何去重。观察可能被重复计数的序列，不难发现这样的序列一定满足子树里的染色节点全为偶数或者全奇数，对应下来颜色就权威0/1或者01相间。
这样很好dp了。
#include
#include
#include
#include
using namespace std;
typedef long long LL;
const int Maxn=100005;
const int Mod=(1e9)+7;
int f[Maxn][2][2],x,tot,i,n;
vector
e[Maxn]; void dp(int x){ f[x][0][1] = f[x][1][1] = 1; f[x][0][0] = f[x][1][0] = 0; if (e[x].empty()) return; sort(e[x].begin(), e[x].end()); int len = e[x].size(); //正序 for (int i=0;i
=0;i--){ int y = e[x][i]; int tmp00 = k00; int tmp01 = k01; int tmp10 = k10; int tmp11 = k11; k00 = ( (tmp00 + (LL)tmp00*f[y][0][0]%Mod)%Mod + (LL)tmp01*f[y][1][1]%Mod )%Mod; k01 = ( (tmp01 + (LL)tmp00*f[y][0][1]%Mod)%Mod + (LL)tmp01*f[y][1][0]%Mod )%Mod; k10 = ( (tmp10 + (LL)tmp10*f[y][1][0]%Mod)%Mod + (LL)tmp11*f[y][0][1]%Mod )%Mod; k11 = ( (tmp11 + (LL)tmp10*f[y][1][1]%Mod)%Mod + (LL)tmp11*f[y][0][0]%Mod )%Mod; } f[x][0][0] = (f[x][0][0] + k00)%Mod; f[x][0][1] = (f[x][0][1] + k01)%Mod; f[x][1][0] = (f[x][1][0] + k10)%Mod; f[x][1][1] = (f[x][1][1] + k11)%Mod; //去重 int g0 = 1, g1 = 1; int h00 = 0, h01 = 1, h10 = 1, h11 = 0; for (int i=0;i
9、bkb
题目里保证n=p^a*q^b，想了很久不知道这是干神马用的。。。
题解里说会相互影响的只会是恰好p*q个组合里的扇叶，也就是说如果这里p*q的扇叶中有缺失，p平衡和q平衡里是能保留一种，由于每组之间互不影响，贪心保留哪一种即可。
#include
#include
#include
#include
using namespace std;
const int Maxn=2000005;
bool d[Maxn], flag;
int ans1,ans2,s[Maxn];
int check[Maxn],prime[Maxn];
int n,m,nn,x,tot,i,j,k,m1,m2,ans,p,q,N;
void init(){
for (i=2;i<=n;i++){
if (!check[i]) check[i]=i, prime[++tot]=i;
for (j=1;j<=tot;j++){
if (i*prime[j]>n) break;
check[ i*prime[j] ] = prime[j];
if (i%prime[j]==0) break;
}
}
}
void work0(){
for (i=0;i
m2) ans-=m1; else ans-=m2; } printf("%d\n",ans); } int main(){ freopen("bkb.in","r",stdin); freopen("bkb.out","w",stdout); scanf("%d%d",&n,&m); if (n==1) {printf("0\n"); return 0;} init(); nn = n; p = check[nn]; while (nn%p==0) nn/=p; q = check[nn]; if (q==1) work0(); else work1(); return 0; }
10、notamotua
给一个只有一个初始状态和一个终止状态的自动机，每次修改串中的一个位置上的字符，询问有多少次能被自动机接受。
匹配这个东东有可合并的性质，利用线段树保留信息即可。
#include
#include
#include
using namespace std;
const int Maxn=100005;
int change,stringLength,modifications;
int n,m,i,j,ans,go[Maxn*4][20],t[20][20];
int x[Maxn],c[Maxn],x0[Maxn],c0[Maxn],xa[Maxn],ca[Maxn];
void init(){
int initElements = change;
int letters = m;
for (int i = 0 ; i < initElements; i++){
x[i] = x0[i] % stringLength ;
c[i] = c0[i] % letters ;
}
for (int i = initElements; i < modifications; i++){
x[i] = 0; c[i] = 0;
for (int j = i - initElements; j
>1;
build(p*2+1,l,mid);
build(p*2+2,mid+1,r);
for (int i=0;i
x || x>r) return; if (l==r && l==x){ for (int i=0;i
>1; modify(p*2+1,l,mid,x,c); modify(p*2+2,mid+1,r,x,c); for (int i=0;i
11、双向桥
求由两条边组成的“桥”。具体方法参见我的另一篇博客：dzy loves Chinese。再用map计数即可
#include
12、seq
算法比较长，又要贴公式，所以博主就懒懒地送上题解
最后的转化和巧妙，值得学习。
另外关于分段那部分其实就是d对于k*d~(k+1)*d-1的贡献是相通的，这样枚举倍数可以做到O(nlogn)算出所有的F。
过去利用分段的性质都是对于单独的数n和若干d，这样的复杂度为O(n^1.5)，这里是一个d若干n，不可一概而论，两种方法都需掌握。
#include
#include
#include
using namespace std;
#define lowbit(x) ((x)&(-(x)))
const int Mod = (1e9)+7;
const int Maxn=100005;
typedef long long LL;
int check[Maxn],mu[Maxn],prime[Maxn],a[Maxn];
int g[Maxn],ny[Maxn],s1[Maxn<<1],s2[Maxn<<1];
int n,m,Case,ans,N,i,j,k,tot,t,tn;
int qck(int a,int b){
int ret = 1;
for (;b>0;b>>=1){
if (b&1) ret = (LL)ret*a%Mod;
a = (LL)a*a%Mod;
}
return ret;
}
void init(){
mu[1]=1;
for (i=2;i<=100000;i++){
if (!check[i]) mu[i]=-1, check[i]=i, prime[++tot]=i;
for (j=1;j<=tot;j++){
if (i*prime[j]>100000) break;
check[i*prime[j]]=prime[j];
if (i%prime[j]==0) {mu[i*prime[j]]=0;break;}
else mu[i*prime[j]]=-mu[i];
}
}
for (i=2,ny[1]=1;i<=100000;i++)
if (check[i]==i) ny[i] = qck(i,Mod-2);
else ny[i] = (LL)ny[i/check[i]]*ny[check[i]]%Mod;
for (i=1,ny[1]=1;i<=100000;i++){
j = i; g[i] = mu[i];
while (j>1){
t = 1;
tn = check[j];
while (j%tn==0){
j /= tn; t *= tn;
g[i] += mu[ i/t ];
}
}
if (g[i]<0) g[i] += Mod;
}
}
int main(){
freopen("seq.in","r",stdin);
freopen("seq.out","w",stdout);
init();
scanf("%d",&Case);
while (Case--){
scanf("%d",&n);
N = 1; m = 0;
memset(s1,0,sizeof(s1));
memset(s2,0,sizeof(s2));
for (i=1;i<=n;i++){
scanf("%d",&a[i]);
m = max(m, a[i]);
N = (LL)N*a[i]%Mod;
s1[a[i]] = (s1[a[i]]+ny[a[i]])%Mod;
s2[a[i]] = (s2[a[i]]+(LL)ny[a[i]]*ny[a[i]]%Mod)%Mod;
}
for (i=1;i<=2*m;i++){
s1[i] = (s1[i]+s1[i-1])%Mod;
s2[i] = (s2[i]+s2[i-1])%Mod;
}
ans = 0;
for (i=1;i<=m;i++){
int x1 = 0, x2 = 0;
for (j=1;i*j<=m;j++){
x1 = ( x1 + (LL)j * ( (s1[(j+1)*i-1]-s1[i*j-1]+Mod)%Mod ) %Mod ) %Mod;
x2 = ( x2 + (LL)j*j%Mod * ( (s2[(j+1)*i-1]-s2[i*j-1]+Mod)%Mod ) %Mod ) %Mod;
}
int tmp = ( (LL)x1*x1%Mod - x2 + Mod )%Mod;
tmp = (LL)tmp * g[i] %Mod;
ans = (ans + tmp) %Mod;
}
ans = (LL)ans * N %Mod * ny[2] %Mod;
printf("%d\n",ans);
}
return 0;
}
13、迷宫
裸高斯消元求期望，加些特技即可
。表示被高斯消元的误差吓傻了。。。。Eps调到1e-20才比较好。。。
#include
#include
#include
using namespace std;
const double Eps = 1e-20;
const int Maxn = 2555;
const int dx[4]={1,0,-1,0};
const int dy[4]={0,1,0,-1};
long double ans[Maxn],a[Maxn][Maxn],t;
double K[4][55][55];
int n,m,i,j,k,ii,jj,N;
double fabs(double a)
{ if (a<0) return -a; return a; }
void Gauss(){
for (i=1;i<=N;i++){
if (fabs(a[i][i])
0;i--){
for (j=i+1;j<=N;j++)
a[i][0] -= a[i][j]*ans[j];
ans[i] = a[i][0]/a[i][i];
}
}
int main(){
freopen("maze.in","r",stdin);
freopen("maze.out","w",stdout);
scanf("%d%d",&n,&m);
for (i=1;i<=n;i++)
for (j=1;j<=m;j++)
scanf("%lf",&K[0][i][j]);
for (i=1;i<=n;i++)
for (j=1;j<=m;j++)
scanf("%lf",&K[1][i][j]);
for (i=1;i<=n;i++)
for (j=1;j<=m;j++)
scanf("%lf",&K[2][i][j]);
for (i=1;i<=n;i++)
for (j=1;j<=m;j++)
scanf("%lf",&K[3][i][j]);
N = n*m;
for (i=1;i<=n;i++)
for (j=1;j<=m;j++){
if (i==n && j==m)
{
a[N][N] = 1; a[N][0] = 0; continue;
}
a[(i-1)*m+j][(i-1)*m+j] = -1;
a[(i-1)*m+j][0] = -1;
for (k=0;k<4;k++){
ii = i+dx[k]; jj = j+dy[k];
if (ii<1 || ii>n || jj<1 || jj>m)
continue;
a[(i-1)*m+j][(ii-1)*m+jj] = K[k][i][j];
}
}
Gauss();
printf("%.3lf\n",(double)ans[1]);
return 0;
}
14、树袋熊
抄袭noi兔农给差评，而且蒟蒻还不会做。。。
利用mod k下为1的作为小循环节，因为此时会让这一位标为0，假设前一位为d，下面的序列就会变为d,d,2*d,3*d......可以发现前面的系数为Fibonacci数列，我们要求第一个Fib满足Fib*d=1(mod k)，用逆元很容易搞出。最终我们可以找到大循环节，用矩乘加速即可。
#include
#include
#include
using namespace std;
const int Maxn=1000005;
typedef long long LL;
LL n,len[Maxn],LEN,ny[Maxn],l,r;
int Fib[Maxn*10],p[Maxn],v[Maxn];
int P,K,p0,p1,p2,k0,k1,k2,cnt,tmp0,tmp1,i;
struct Matrix
{
int a[3][3];
Matrix operator *(const Matrix &x)const
{
Matrix ret;
for (int i=0;i<3;i++)
for (int j=0;j<3;j++){
ret.a[i][j]=0;
for (int k=0;k<3;k++)
ret.a[i][j] = ( ret.a[i][j] + (LL)a[i][k]*x.a[k][j]%P ) %P;
}
return ret;
}
} Wx, W, A[Maxn], B;
int exgcd(int a,int b,LL &x,LL &y){
if (b==0){
x = 1; y = 0;
return a;
}
int ret = exgcd(b,a%b,x,y);
LL tmp = x;
x = y; y = tmp - (LL)(a/b)*y;
return ret;
}
void init(){
Fib[1] = Fib[2] = 1;
ny[1] = ny[2] = 1;
p[1] = 1;
for (i=1;i
K) Fib[i] -= K;
if (ny[Fib[i]]!=-1 && p[ny[Fib[i]]]==0)
p[ny[Fib[i]]]=i;
}
}
Matrix qck(Matrix A,LL p){
Matrix ret;
memset(ret.a,0,sizeof(ret.a));
ret.a[0][0] = ret.a[1][1] = ret.a[2][2] = 1;
for (;p>0;p>>=1){
if (p&1) ret = ret*A;
A = A*A;
}
return ret;
}
void work(){
W.a[0][1]=W.a[2][2]=1;
W.a[1][0]=W.a[1][1]=1;
Wx.a[0][1]=Wx.a[2][2]=1;
Wx.a[1][0]=Wx.a[1][1]=1;
Wx.a[1][2]=P-1;
k0 = k1 = 1;
p0 = p1 = 1;
for (i=3;i<=10*K&&i<=n;i++){
k2 = k0+k1;
p2 = p0+p1;
if (k2>K) k2 -= K;
if (p2>P) p2 -= P;
k0 = k1; k1 = k2;
p0 = p1; p1 = p2;
if (k1==1) {k1--;p1--;break;}
}
if (i>n){
printf("%d\n",p1);
return;
}
n = n-i;
if (p1<0) p1 += P;
memset(v,-1,sizeof(v));
bool flag = 0;
while (true){
if (p[k0]==0) {flag=1;break;}
int t = p[k0];
A[++cnt] = Wx * qck( W, t-1 );
k0 = (LL)Fib[t-1]*k0%K;
len[cnt] = t;
if (v[k0]>=0) break;
v[k0] = cnt+1;
}
if (flag){
for (i=1;i<=cnt&&n>=len[i];i++){
tmp0 = ( ( (LL)p0*A[i].a[0][0]%P+(LL)p1*A[i].a[0][1]%P )%P + A[i].a[0][2] ) %P;
tmp1 = ( ( (LL)p0*A[i].a[1][0]%P+(LL)p1*A[i].a[1][1]%P )%P + A[i].a[1][2] ) %P;
p0 = tmp0; p1 = tmp1;
n = n-len[i];
}
B = qck( W, n );
tmp0 = ( ( (LL)p0*B.a[0][0]%P+(LL)p1*B.a[0][1]%P )%P + B.a[0][2] ) %P;
tmp1 = ( ( (LL)p0*B.a[1][0]%P+(LL)p1*B.a[1][1]%P )%P + B.a[1][2] ) %P;
p0 = tmp0; p1 = tmp1;
printf("%d\n",p1);
return;
}
l = v[k0]; r = cnt+1;
for (i=1;i
len[i];i++){ tmp0 = ( ( (LL)p0*A[i].a[0][0]%P+(LL)p1*A[i].a[0][1]%P )%P + A[i].a[0][2] ) %P; tmp1 = ( ( (LL)p0*A[i].a[1][0]%P+(LL)p1*A[i].a[1][1]%P )%P + A[i].a[1][2] ) %P; p0 = tmp0; p1 = tmp1; n = n-len[i]; } if (i
=len[i];n-=len[i++]) B = A[i]*B; while (n--) B = W*B; tmp0 = ( ( (LL)p0*B.a[0][0]%P+(LL)p1*B.a[0][1]%P )%P + B.a[0][2] ) %P; tmp1 = ( ( (LL)p0*B.a[1][0]%P+(LL)p1*B.a[1][1]%P )%P + B.a[1][2] ) %P; p0 = tmp0; p1 = tmp1; printf("%d\n",p1); } int main(){ freopen("koala.in","r",stdin); freopen("koala.out","w",stdout); scanf("%lld%d%d",&n,&K,&P); init(); work(); return 0; }
15、mc
求不是回文子串的回文字序列。
考虑用暴力的话可以枚举中点k求有多少i满足S[k-i]=S[k+i]，假设为t个，那么子序列数为(2^t-1)，最后减去子串数。
优化暴力。不难发现求t值其实就是一个卷积的形式，用fft加速。
#include
#include
#include
#include
using namespace std;
const int Maxn=100005;
const int Mod=(1e9)+7;
int col[Maxn],two[Maxn<<3],s[Maxn<<3];
int num[Maxn],len[Maxn],son[Maxn][2],fail[Maxn];
int n,m,k,N,i,ans,st,c,p,t,x,last;
struct CP
{
double x, y;
CP operator +(const CP &a)const
{ return (CP){x+a.x, y+a.y}; }
CP operator -(const CP &a)const
{ return (CP){x-a.x, y-a.y}; }
CP operator *(const CP &a)const
{ return (CP){x*a.x-y*a.y, y*a.x+x*a.y}; }
} A[Maxn<<3];
void FFT(CP A[],int N,int flag){
for (int i=1, j=0;j
>=1)); if (i
>1)<=k;N<<=1); for (i=1,two[0]=1;i
扫一扫
674
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
