SDOI2013 R1 Day2
目录
2018.3.25 Test
时间：7:30~11:30 (最后半小时不做了)
期望得分：50+100+20=170
实际得分：40+44+20=104
总结
T1：善用容斥。
T2：要求输出小数当然有它的道理。。
T3：行列分开考虑。看好怎么取模。如要求最大的答案，然后对答案取模，那么在计算过程中不能取模，最后求答案时才可以！
T1 BZOJ.3129.[SDOI2013]方程(扩展Lucas 容斥)
如果n个未知数是正整数且没有限制，那么方案数是C(m-1,n-1).(把m个数分割成n份)
不需要把x变为非负整数解，然后去求C(n-k+1,k).(当然无所谓吧。。按标程来)
\(\geq\) 的限制很好做，减去即可。那么对于 \(xi\leq Ai\) 的我们容斥，去掉x的限制，然后把 \(x>Ai\) 的减掉。
利用扩展Lucas计算组合数.复杂度 \(O(2^{n1}*Ex\_Lucas)\).
复杂度有点高了。。扩展Lucas预处理阶乘！(预处理后飞快啊，以前从来没预处理。。)
//1620kb 1364ms
//预处理前:840kb 10552ms
#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=1e4+23;
int n,P,n1,n2,m,A[10],cnt,pi[666],pk[666];
LL Ans,fac[100005]/*这个开多大还是不清楚*/;
inline int read()
{
int now=0;register char c=gc();
for(;!isdigit(c);c=gc());
for(;isdigit(c);now=now*10+c-'0',c=gc());
return now;
}
LL FP(LL x,int k,LL p)
{
LL t=1ll;
for(; k; k>>=1,x=x*x%p)
if(k&1) t=t*x%p;
return t;
}
void Exgcd(int a,int b,int &x,int &y)
{
if(!b) x=1,y=0;
else Exgcd(b,a%b,y,x),y-=a/b*x;
}
inline LL inv(int a,int p)
{
int x,y; Exgcd(a,p,x,y);
return (LL)((x%p+p)%p);
}
void Init_P(int p)
{
for(int i=2,lim=sqrt(p+1); i<=lim; ++i)
if(!(p%i))
{
int _pk=1;
while(!(p%i)) p/=i,_pk*=i;
pi[++cnt]=i, pk[cnt]=_pk;
if(p==1) break;
}
if(p!=1) pi[++cnt]=pk[cnt]=p;
}
LL Fact(int n,int pi,int pk)
{
if(!n) return 1ll;//边界！
LL res=1ll;
if(n>=pk) res=FP(fac[pk],n/pk,pk);
(res*=fac[n%pk])%=pk;
return res*Fact(n/pi,pi,pk)%pk;
}
LL C(int n,int m,int pi,int pk,int mod)
{
fac[0]=1;
for(int i=1; i<=pk; ++i)
if(i%pi) fac[i]=fac[i-1]*i%pk;
else fac[i]=fac[i-1];
LL a=Fact(n,pi,pk),b=Fact(m,pi,pk),c=Fact(n-m,pi,pk);int k=0;
for(int i=n; i; i/=pi) k+=i/pi;
for(int i=m; i; i/=pi) k-=i/pi;
for(int i=n-m; i; i/=pi) k-=i/pi;
LL res=a*inv(b,pk)%pk*inv(c,pk)%pk*FP(pi,k,pk)%pk;
return 1ll*res*(mod/pk)%mod*inv(mod/pk,pk)%mod;
}
LL Ex_Lucas(int n,int m)
{
if(n<m) return 0ll;
LL res=0;
for(int i=1; i<=cnt; ++i)
(res+=C(n,m,pi[i],pk[i],P))%=P;
return res;
}
void DFS(int p,int cnt,int sum)
{
if(p>n1)
if(cnt&1) (Ans-=Ex_Lucas(m-sum-1,n-1))%=P;//前边的数减去A[i]，即满足>A[i]。
else (Ans+=Ex_Lucas(m-sum-1,n-1))%=P;//C(m-sum-1,n-1) not C(n-sum-1,m-1)!
else DFS(p+1,cnt,sum),DFS(p+1,cnt+1,sum+A[p]);
}
int main()
{
int T=read(); P=read();
Init_P(P);
while(T--)
{
n=read(),n1=read(),n2=read(),m=read();
LL tot=0;
for(int i=1; i<=n1; ++i) A[i]=read();
for(int i=1; i<=n2; ++i) m-=read()-1;//去下界
if(m<=1) {puts("0"); continue;}
if(m==n) {puts("1"); continue;}
Ans=0ll, DFS(1,0,0), printf("%lld\n",(Ans+P)%P);//取模！
}
return 0;
}
T2 洛谷.3305.[SDOI2013]费用流(最大流ISAP 二分)
题目链接(BZOJ无SPJ)
题意即是求满足最大流条件下 使得流量最大的边尽量小。花费显然是最大的某条边上的流量*P。
但P是实数范围内！这说明流量也是实数，二分答案跑最大流检验。限制的话，直接重设流量。
贪心的话最大流过程中流量一直是整数，所以就GG了。。
#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
const double eps=1e-6;
const int N=218,M=2018;
int n,m,P,Max,src,des,Enum,H[N],cur[N],fr[M],to[M],nxt[M],cost[M],num[N],lev[N],pre[N],que[N];
double cap[M],Cap[M];
inline int read()
{
int now=0;register char c=gc();
for(;!isdigit(c);c=gc());
for(;isdigit(c);now=now*10+c-'0',c=gc());
return now;
}
inline void AddEdge(int u,int v,int w)
{
Max=std::max(Max,w);
to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, Cap[Enum]=cap[Enum]=w;
to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, Cap[Enum]=cap[Enum]=0;
}
bool BFS()
{
for(int i=src; i<des; ++i) lev[i]=des+1;
lev[des]=0, que[0]=des; int h=0,t=1;
while(h<t)
{
int x=que[h++];
for(int i=H[x]; i; i=nxt[i])
if(lev[to[i]]==des+1 && cap[i^1])
lev[to[i]]=lev[x]+1, que[t++]=to[i];
}
return lev[src]<=des;
}
double Augment()
{
double mn=1e8;
for(int i=des; i!=src; i=fr[pre[i]])
mn=std::min(mn,cap[pre[i]]);
for(int i=des; i!=src; i=fr[pre[i]])
cap[pre[i]]-=mn, cap[pre[i]^1]+=mn;
return /*fabs(mn)<eps?0:*/mn;
}
double ISAP()
{
if(!BFS()) return 0;
for(int i=src; i<=des; ++i) ++num[lev[i]],cur[i]=H[i];
int x=src; double res=0;
while(lev[src]<=des)
{
if(x==des) x=src,res+=Augment();
bool can=0;
for(int i=cur[x]; i; i=nxt[i])
if(lev[to[i]]==lev[x]-1 && cap[i])
{
can=1, cur[x]=i, pre[x=to[i]]=i;
break;
}
if(!can)
{
int mn=des;
for(int i=H[x]; i; i=nxt[i])
if(cap[i]) mn=std::min(mn,lev[to[i]]);
if(!--num[lev[x]]) break;
++num[lev[x]=mn+1];
cur[x]=H[x];
if(x!=src) x=fr[pre[x]];
}
}
return res;
}
void Init(double lim){
for(int i=2; i<=Enum; ++i) cap[i]=std::min(Cap[i],lim);
}
double Solve(double mxflow)
{
double l=1.0,r=Max,m;
while(r-l>1e-4)
if(Init(m=(l+r)*0.5),fabs(mxflow-ISAP())<eps) r=m;
else l=m;
return l;
}
int main()
{
n=read(),m=read(),P=read(),Enum=src=1,des=n;
for(int u,v,w,i=1; i<=m; ++i)
u=read(),v=read(),w=read(),AddEdge(u,v,w);//f&=(w==1);
double res1=ISAP();
printf("%.0lf\n%.4lf",res1,1.0*Solve(res1)*P);
return 0;
}
T3 BZOJ.3131.[SDOI2013]淘金(数位DP 堆)
40pts:
行列是独立的！即同在第i行的元素一定都会加到f(i)行，它们对所有列位置的贡献同其它行一样是固定的，只是加的行不同而已。列同理。而且行列的贡献是相同的。
所以求出(一维下)有多少个数能转移到某个位置i。找出两两相乘(二维)最大的K个。-> 这40pts我也不会写啊QAQ
其实是用堆 然后每个节点维护两个指针，拿一个就放两个？并不，这样在x==y时会有重。(其实也有办法 大概要判重之类？)
如果x是初始的s_x，则不断扩展 push(x,y+1)；如果x!=y，则push(x+1,y)(维护一个有序数对，感觉好6啊。。)
而且在x!=y时加两遍，因为行列是可以在反过来的！
(也可以先push进n个去再扩展吧。。)
100pts:
\(n\)很大无法直接算出某个位置的\(val[]\)。
在\(n\)的范围内\(f(x)\)的种类只有不是很多种(大概不到\(20000\)?)，且因子只会有\(2.3.5.7\)，可以通过数位DP算出\(y=2^i*3^j*5^k*7^l\)时有多少个x的数位乘积=\(y(x\leq n)\).
\(dp[len][0/1][i][j][k][l]\)表示当前是\(n\)的第\(len\)位，前面是否达到了上界，以及\(y\)的因子个数分别为\(i,j,k,l\)时的\(x\)个数。
最后统计所有\(y\)的\(val[]\)。
注意dp的细节: 可以有前导\(0\)，除了最高位都要计算前导\(0\)的贡献，即\(dp[now][0][0][0][0][0]\)++(不卡上界)；最高位初始化时是\(dp[now][1][0][0][0][0]=1\)(上界)。
（并没有自己写代码。。借鉴(抄)下xxy博客了）
另外这个取模是在取最多的情况下算出后取模！\(val[]\)运算过程中不能取模。(是这样吧)
//12124kb 740ms
#include <queue>
#include <cstdio>
#include <cctype>
#include <cstring>
#include <algorithm>
#define gc() getchar()
#define mod (1000000007)
typedef long long LL;
const int N=1e5+7;
int K,cnt,f[11][5]/*多一位数影响的数值*/,A[20];
//LL n,val[N],dp[2][2][40][26][18][15];
LL n,val[N],dp[2][2][41][27][19][16];
struct Node
{
int x,y; LL v;
Node(int _x,int _y): x(_x),y(_y),v(val[x]*val[y]) {}
bool operator <(const Node &a)const{
return v<a.v;//val升序
}
};
std::priority_queue<Node> q;
inline bool cmp(const int &a,const int &b){
return a>b;
}
inline void Pre()
{//0:2 1:3 2:5 3:7
f[2][0]=1, f[4][0]=2, f[6][0]=1, f[8][0]=3,
f[3][1]=1, f[6][1]=1, f[9][1]=2,
f[5][2]=1, f[7][3]=1;
}
inline LL FP(LL x,int k)
{
LL t=1;
for(; k; k>>=1,x=x*x)
if(k&1) t=t*x;
return t;
}
void DP()
{//2:40 3:26 5:18 7:15
Pre();
LL t=n; int len=0;
while(t) A[len++]=t%10,t/=10;//我一定是sb了。。下标从0开始还用A[0]做长度
std::reverse(A,A+len);
int now=0,nxt=1;
for(int p=0; p<len; ++p,std::swap(now,nxt))//神tm还拿l做变量。。
{
if(!p) dp[now][1][0][0][0][0]=1;
else ++dp[now][0][0][0][0][0];
memset(dp[nxt],0,sizeof dp[nxt]);
for(int s=0; s<=1; ++s)
for(int i=0; i<=40; ++i)
for(int j=0; j<=26; ++j)
for(int k=0; k<=18; ++k)
for(int l=0; l<=15; ++l)
if(dp[now][s][i][j][k][l])//怎么这么。。
for(int w=1,lim=s?A[p]:9; w<=lim; ++w)//不能取模 //另外不计算0！
dp[nxt][s&&w==lim][i+f[w][0]][j+f[w][1]][k+f[w][2]][l+f[w][3]]+=dp[now][s][i][j][k][l];
}
LL tmp2=1,sum;
for(int i=0; i<=40; ++i)
{
LL tmp3=tmp2;
for(int j=0; j<=26; ++j)
{
LL tmp5=tmp3;
for(int k=0; k<=18; ++k)
{
LL tmp7=tmp5;
for(int l=0; l<=15; ++l)
{
if(sum=(dp[now][0][i][j][k][l]+dp[now][1][i][j][k][l])%mod) val[++cnt]=sum;
if((tmp7*=7)>n) break;
}
if((tmp5*=5)>n) break;
}
if((tmp3*=3)>n) break;
}
if((tmp2<<=1)>n) break;
}
}
int main()
{
scanf("%lld%d",&n,&K);
DP();
std::sort(val+1,val+1+cnt,cmp);
LL res=0;
q.push(Node(1,1));
while(!q.empty())
{
Node rt=q.top(); q.pop();
(res+=rt.v)%=mod;
if(!--K) break;
if(rt.x!=rt.y){
(res+=rt.v)%=mod;
if(!--K) break;
q.push(Node(rt.x+1,rt.y));
}
if(rt.x==1) q.push(Node(rt.x,rt.y+1));
}
printf("%lld",res);
return 0;
}
考试代码
T1
恶心的Lucas+Ex_Lucas。
#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=2e4+23;
int n,n1,n2,m,A[10],B[10],fac[N],f[10][N],sum[N],cnt,pi[666],pk[666];
inline int read()
{
int now=0;register char c=gc();
for(;!isdigit(c);c=gc());
for(;isdigit(c);now=now*10+c-'0',c=gc());
return now;
}
LL FP(LL x,int k,LL p)
{
LL t=1ll;
for(; k; k>>=1,x=x*x%p)
if(k&1) t=t*x%p;
return t;
}
void Exgcd(int a,int b,int &x,int &y)
{
if(!b) x=1,y=0;
else Exgcd(b,a%b,y,x),y-=a/b*x;
}
inline LL inv(int a,int p)
{
int x,y; Exgcd(a,p,x,y);
return (LL)((x%p+p)%p);
}
LL Calc(int n,int m,int p)
{
if(n<m) return 0ll;
return 1ll*fac[n]*inv(fac[m],p)%p*inv(fac[n-m],p)%p;
}
LL Lucas(int n,int m,int p)//p=10007
{
// printf("Lucas:C(%d,%d) mod %d = ",n,m,p);
if(n<m) return 0ll;
LL res=1;
for(; m&&res; n/=p,m/=p)
(res*=Calc(n%p,m%p,p))%=p;
// printf("%d\n",res);
return res;
}
LL Fact(int n,int pi,int pk)
{
if(!n) return 1ll;//边界！
LL res=1ll;
if(n>=pk)
{
for(int i=2; i<pk; ++i)
if(i%pi) (res*=i)%=pk;
res=FP(res,res/pk,pk);//这怎么写成res/pk了。。
}
for(int i=2; i<=n%pk; ++i)
if(i%pi) (res*=i)%=pk;
return res*Fact(n/pi,pi,pk)%pk;
}
LL C(int n,int m,int pi,int pk,int mod)
{
LL a=Fact(n,pi,pk),b=Fact(m,pi,pk),c=Fact(n-m,pi,pk);int k=0;
for(int i=n; i; i/=pi) k+=i/pi;
for(int i=m; i; i/=pi) k-=i/pi;
for(int i=n-m; i; i/=pi) k-=i/pi;
LL res=a*inv(b,pk)%pk*inv(c,pk)%pk*FP(pi,k,pk)%pk;
return res*(mod/pk)%mod*inv(mod/pk,pk)%mod;
}
LL Ex_Lucas(int n,int m,int p)
{
// printf("Ex_Lucas:C(%d,%d) mod %d = ",n,m,p);
if(n<m) return 0ll;
if(n==m) return 1ll;
int res=0;
for(int i=1; i<=cnt; ++i)
(res+=C(n,m,pi[i],pk[i],p))%=p;
// printf("%d\n",res);
return (LL)res;
}
void Init_P(int p)
{
for(int i=2,m=sqrt(p+1); i<=m; ++i)
if(!(p%i))
{
int pk_=1;
while(!(p%i)) p/=i,pk_*=i;
pi[++cnt]=i, pk[cnt]=pk_;
if(p==1) break;
}
}
int main()
{
freopen("equation.in","r",stdin);
freopen("equation.out","w",stdout);
int T=read(),p=read();
fac[0]=1;
if(p==10007)
for(int x,y,i=1; i<p; ++i) fac[i]=fac[i-1]*i%p;
else Init_P(p);
while(T--)
{
n=read(),n1=read(),n2=read(),m=read()-n+n2;
LL tot=0;
for(int i=1; i<=n1; ++i) A[i]=read()-1,tot+=A[i];
for(int i=1; i<=n2; ++i) B[i]=read(),m-=B[i];
// printf("1: n:%d m:%d\n",n,m);
if(m<0) {puts("0"); continue;}
if(!m) {puts("1"); continue;}
if(m==1)
{
int res=n;
for(int i=1; i<=n1; ++i) if(!A[i]) --res;
printf("%d\n",res);
continue;
}
A[n1+1]=0;
for(int j=0; j<=A[1]; ++j) f[1][j]=1;
sum[0]=1;
for(int j=1; j<=A[1]+A[2]; ++j) sum[j]=sum[j-1]+f[1][j];
for(int tmp=A[1],tmp2=A[1]+A[2],i=2; i<=n1; ++i)//O(8m)
{
tmp+=A[i], f[i][0]=1;
for(int j=1; j<=tmp; ++j)
{
if(j<std::min(j,A[i])) f[i][j]=sum[j];
else f[i][j]=sum[j]-sum[j-std::min(j,A[i])-1], f[i][j]=(f[i][j]%p+p)%p;
// printf("f[%d][%d]:%d %d bef:%d %d\n",i,j,f[i][j],p,sum[j],sum[j-std::min(j-1,A[i])-1]);
}
tmp2+=A[i+1], sum[0]=1;
for(int j=1; j<=tmp2; ++j)
{
sum[j]=sum[j-1]+f[i][j];
if(sum[j]>=p) sum[j]-=p;
}
}
LL res=0; f[n1][0]=1;
// printf("2: n:%d m:%d\n",n,m);
tot=std::min(tot,(LL)m);
// for(int i=0; i<=tot; ++i) printf("f:%d:%d\n",i,f[n1][i]);
if(p==10007) for(int i=0; i<=tot; ++i) (res+=1ll*f[n1][i]*Lucas(n+m-i-n1-1,m-i,p))%=p;
else for(int i=0; i<=tot; ++i) (res+=1ll*f[n1][i]*Ex_Lucas(n+m-i-n1-1,m-i,p))%=p;
printf("%I64d\n",res);
}
return 0;
}
T2
答案都比标程正好大了一倍，气啊。
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
const int N=218,M=2018,INF=0x3f3f3f3f;
int n,m,P,Max,src,des,Enum,H[N],cur[N],fr[M],to[M],cap[M],nxt[M],cost[M],num[N],lev[N],pre[N],que[N];
inline int read()
{
int now=0;register char c=gc();
for(;!isdigit(c);c=gc());
for(;isdigit(c);now=now*10+c-'0',c=gc());
return now;
}
inline void AddEdge(int u,int v,int w)
{
to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, cap[Enum]=w;
to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, cap[Enum]=0;
}
bool BFS()
{
for(int i=src; i<des; ++i) lev[i]=des+1;
lev[des]=0, que[0]=des; int h=0,t=1;
while(h<t)
{
int x=que[h++];
for(int i=H[x]; i; i=nxt[i])
if(lev[to[i]]==des+1 && cap[i^1])
lev[to[i]]=lev[x]+1, que[t++]=to[i];
}
return lev[src]<=des;
}
int Augment()
{
int mn=INF;
for(int i=des; i!=src; i=fr[pre[i]])
mn=std::min(mn,cap[pre[i]]);
for(int i=des; i!=src; i=fr[pre[i]])
cap[pre[i]]-=mn, cap[pre[i]^1]+=mn;
// printf("Augment:%d\n",mn);
// Max=std::max(Max,mn);
return mn;
}
int ISAP()
{
if(!BFS()) return 0;
for(int i=src; i<=des; ++i) ++num[lev[i]],cur[i]=H[i];
int x=src,res=0;
while(lev[src]<=des)
{
if(x==des) x=src,res+=Augment();
int way=0,mn=0;
for(int i=cur[x]; i; i=nxt[i])
if(lev[to[i]]==lev[x]-1 && cap[i])
if(!mn) way=mn=i;
else if(cap[way]>cap[i]) way=i;
if(way) cur[x]=mn, pre[x=to[way]]=way;
else
{
int mn=des;
for(int i=H[x]; i; i=nxt[i])
if(cap[i]) mn=std::min(mn,lev[to[i]]);
if(!--num[lev[x]]) break;
++num[lev[x]=mn+1];
cur[x]=H[x];
if(x!=src) x=fr[pre[x]];
}
}
return res;
}
int main()
{
freopen("costflow.in","r",stdin);
freopen("costflow.out","w",stdout);
n=read(),m=read(),P=read(),Enum=src=1,des=n;
for(int u,v,w,i=1; i<=m; ++i)
u=read(),v=read(),w=read(),AddEdge(u,v,w);//f&=(w==1);
int res1=ISAP();
// if(f) {printf("%d\n%d.0000",res1,P); return 0;}
for(int i=3; i<=Enum; i+=2) Max=std::max(Max,cap[i]);//应是从最大流后全局的反向弧找Max
printf("%d\n%d.0000",res1,Max*P);
return 0;
}
T3
#include <queue>
#include <cstdio>
#include <cctype>
#include <cstring>
#include <algorithm>
#define gc() getchar()
#define mod (1000000007)
typedef long long LL;
const int N=1e6+7;
int K,f[N+3],sum[1005][1005];//size!
LL n;
std::priority_queue<int> q;
void Violence()
{
for(int i=1; i<=n; ++i)
for(int j=1; j<=n; ++j)
++sum[f[i]][f[j]];
for(int i=1; i<=n; ++i)
for(int j=1; j<=n; ++j)
if(sum[i][j]) q.push(sum[i][j]);
// for(int i=1; i<=n; ++i,putchar('\n'))
// for(int j=1; j<=n; ++j)
// printf("%d ",sum[i][j]);
LL res=0;
while(K-- && !q.empty()) res+=q.top(),q.pop();
printf("%I64d",res%mod);
}
int main()
{
freopen("gold.in","r",stdin);
freopen("gold.out","w",stdout);
scanf("%I64d%d",&n,&K);
f[0]=1;
for(int i=1; i<N; ++i)
if(i%10) f[i]=f[i/10]*(i%10);
else f[i]=0;
f[0]=0;
// for(n=100; n<=300; ++n,putchar('\n'))
// memset(sum,0,sizeof sum),printf("%I64d:\t",n),Violence();
if(n<=1000) {Violence(); return 0;}
printf("%d",rand()%mod);
return 0;
}
------------------------------------------------------------------------------------------------------------------------
很久以前的奇怪但现在依旧成立的签名
attack is our red sun $$\color{red}{\boxed{\color{red}{attack\ is\ our\ red\ sun}}}$$ ------------------------------------------------------------------------------------------------------------------------
很久以前的奇怪但现在依旧成立的签名
attack is our red sun $$\color{red}{\boxed{\color{red}{attack\ is\ our\ red\ sun}}}$$ ------------------------------------------------------------------------------------------------------------------------
