#include<stdio.h>
#include<stdlib.h>
#include<time.h>
#include<sys/time.h>
#include<cuda.h>
#include<cuda_runtime.h>
#define MAXlid 5
#define MAXst 40
#define MAXntr 180
#define MAXtd 7
#define MAXtid 6000
#define MAXtm4test 1200
#define MAXtm 1200
#define MAXtvtm 7
#define MAXstPline 20
int Nst;
int Ntid;
int maxDist;
short int trainTBnst[MAXtid];
short int trainTB[MAXtid][MAXstPline][3];
int outroute[MAXtm][3][2];
void input() {
int i,j,k;
scanf("%d", &Nst);
int nlid;
scanf("%d", &nlid);
int yobi;
for(i=0; i<nlid; i++) {
int m;
scanf("%d", &m);
for(j=0; j<m; j++) scanf("%d", &yobi);
}
for(i=0; i<Nst; i++) for(j=0; j<Nst; j++) scanf("%d", &yobi);
scanf("%d", &Ntid);
int tid;
for(tid=0; tid<Ntid; tid++) {
int m;
scanf("%d", &m);
trainTBnst[tid] = m;
int j;
for(j=0; j<m; j++) {
int st , tm, dist;
scanf("%d %d %d", &st, &tm, &dist);
trainTB[tid][j][0] = st;
trainTB[tid][j][1] = tm;
trainTB[tid][j][2] = dist;
}
}
for(i=0; i<MAXtm; i++) for(j=0; j<3; j++) for(k=0; k<2; k++) outroute[i][j][k] = -1;
{
void *p;
cudaError_t err;
err = cudaMalloc(&p, 4096);
if (err != cudaSuccess) {
fprintf(stderr, "cudaMalloc in input() failed!\n");
exit(1);
}
err = cudaFree(p);
if (err != cudaSuccess) {
fprintf(stderr, "cudaFree in input() failed!\n");
exit(1);
}
}
}
void output() {
int i,j;
printf("##########\n");
printf("maxDist: %d\n",maxDist);
for(i=0; i<MAXtm; i++) {
for(j=0; j<3; j++) printf("%2d %2d ", outroute[i][j][0], outroute[i][j][1]);
printf("\n");
}
}
int M,V;
struct edge{
int from,to;
short int cost,cap;
int rev;
bool type;
short int id;
};
struct edge2{
int from,to;
short int cost,cap;
int rev;
};
const int INF=30000;
int source,sink;
const int BASE=100;
const int MAX_V=100000,MAX_E=MAXstPline*MAXtid*2+MAXst*MAXtm*4+MAXst*2+BASE;
struct edge es[MAX_E];
struct edge2 es2[MAX_E];
struct edge *g[MAX_V][22];
int deg[MAX_V];
void prep(){
int i,j;
M=0;
V=Nst*MAXtm*2+5;
source=V-3;
sink=V-2;
for(i=0;i<Ntid;++i){
int back=-1,backtime;
for(j=0;j<trainTBnst[i];++j){
int sid=trainTB[i][j][0],tm=trainTB[i][j][1],dist=trainTB[i][j][2];
if(back!=-1){
es[M]=(edge){(backtime*Nst+back)*2+1,(sid+tm*Nst)*2,dist,1,M+1,1,i};++M;
es[M]=(edge){(sid+tm*Nst)*2,(backtime*Nst+back)*2+1,-dist,0,M-1,0,i};++M;
}
backtime=tm;back=sid;
}
}
for(i=0;i<Nst*MAXtm;++i){
es[M]=(edge){i*2,i*2+1,0,1,M+1,1,-1};++M;
es[M]=(edge){i*2+1,i*2,0,0,M-1,0,-1};++M;
}
for(i=0;i<Nst;++i){
es[M]=(edge){source,i*2,0,1,M+1,1,-1};++M;
es[M]=(edge){i*2,source,0,0,M+1,0,-1};++M;
es[M]=(edge){(i+(MAXtm-1)*Nst)*2+1,sink,0,1,M+1,1,-1};++M;
es[M]=(edge){sink,(i+(MAXtm-1)*Nst)*2+1,0,0,M-1,0,-1};++M;
}
for(i=0;i<Nst*(MAXtm-1);++i){
es[M]=(edge){i*2+1,(i+Nst)*2,0,1,M+1,1,-1};++M;
es[M]=(edge){(i+Nst)*2,i*2+1,0,0,M-1,0,-1};++M;
}
for(i=0;i<M;++i) g[es[i].from][deg[es[i].from]++]=&es[i];
}
short int cost[MAX_V];
int prevv[MAX_V],preve[MAX_V];
void max_cost_flow_1(int f=3){
int i,j;
int first_loop=1;
while(f>0){
for(i=0;i<V;++i){
prevv[i]=preve[i]=cost[i]=-INF;
}
cost[source]=0;
int changed=1;
int loop=0;
while(changed){
++loop;
changed=0;
int first=0;
for(i=source;i<V;(i==source?i=0:++i)){
if(i==source){
if(!first) first=1;
else i=sink;
}
for(j=0;j<deg[i];++j){
struct edge* e=g[i][j];
if(e->cap && cost[e->to]<cost[e->from]+e->cost){
cost[e->to]=cost[e->from]+e->cost;
prevv[e->to]=e->from;
preve[e->to]=j;
}
}
}
if(first_loop==1){
first_loop=0;
break;
}
}
if(cost[sink]<0){
return;
}
int cur=sink;
maxDist+=cost[sink];
while(cur!=source){
g[prevv[cur]][preve[cur]]->cap=0;
es[g[prevv[cur]][preve[cur]]->rev].cap=1;
cur=prevv[cur];
}
--f;
}
}
__device__ short int costDev[MAX_V];
__device__ unsigned int preveDev[MAX_V];
__global__ void initDev(){
costDev[blockIdx.x*blockDim.x+threadIdx.x]=SHRT_MIN;
}
__global__ void prep(int source){
costDev[source]=0;
}
__global__ void doMain(struct edge2* es,bool* changedFlag){
struct edge2* e=&es[blockIdx.x*blockDim.x+threadIdx.x];
if(e->cap && costDev[e->to]<costDev[e->from]+e->cost){
costDev[e->to]=costDev[e->from]+e->cost;
preveDev[e->to]=blockIdx.x*blockDim.x+threadIdx.x;
changedFlag[0]=true;
}
}
__global__ void doMain2(struct edge2* es){
struct edge2* e=&es[blockIdx.x*blockDim.x+threadIdx.x];
if(e->cap && costDev[e->to]<costDev[e->from]+e->cost){
costDev[e->to]=costDev[e->from]+e->cost;
preveDev[e->to]=blockIdx.x*blockDim.x+threadIdx.x;
}
}
__global__ void set(bool* changedFlag){
changedFlag[0]=false;
}
__global__ void makepath2(int sink,int* maxDist,int source,struct edge2* es){
maxDist[0]+=costDev[sink];
int cur=sink;
while(cur!=source){
es[preveDev[cur]].cap=0;
es[es[preveDev[cur]].rev].cap=1;
cur=es[preveDev[cur]].from;
}
}
void max_cost_flow_2(int f=3){
while(M%BASE){
es[M++]=(edge){0,0,0,0,0,0,-1};
}
int i;
for(i=0;i<M;++i){
es2[i].from=es[i].from;
es2[i].to= es[i].to;
es2[i].cost=es[i].cost;
es2[i].rev= es[i].rev;
es2[i].cap= es[i].cap;
}
struct edge2* pos;
bool* bpos;
bool btmp[1];btmp[0]=false;
cudaMalloc((void **)&bpos,sizeof(bool));
cudaMalloc((void **)&pos,sizeof(struct edge2)*M);
cudaMemcpy(pos,es2,sizeof(struct edge2)*M,cudaMemcpyHostToDevice);
int* maxDistDev;
cudaMalloc((void **)&maxDistDev,sizeof(int));
int tmp[1]={maxDist};
cudaMemcpy(maxDistDev,tmp,sizeof(int),cudaMemcpyHostToDevice);
int lastloop=0;
while(f>0){
initDev<<<1000,100>>>();
prep<<<1,1>>>(source);
cudaDeviceSynchronize();
int changed=V+1;
int loop=0;
while(changed--){
++loop;
if(loop>=lastloop/2 && loop%30==15){
set<<<1,1>>>(bpos);
doMain<<<M/BASE,BASE>>>(pos,bpos);
cudaMemcpy(btmp,bpos,sizeof(bool),cudaMemcpyDeviceToHost);
if(btmp[0]==false) break;
}else{
doMain2<<<M/BASE,BASE>>>(pos);
}
}
lastloop=loop;
makepath2<<<1,1>>>(sink,maxDistDev,source,pos);
--f;
}
cudaMemcpy(tmp,maxDistDev,sizeof(int),cudaMemcpyDeviceToHost);
maxDist=tmp[0];
cudaMemcpy(es2,pos,sizeof(struct edge2)*M,cudaMemcpyDeviceToHost);
cudaFree(maxDistDev);
cudaFree(pos);
cudaFree(bpos);
for(i=0;i<M;++i){
es[i].cap=es2[i].cap;
}
}
void makepath(){
int i,j,k;
for(i=0;i<3;++i){
int cur=source;
for(j=0;j<deg[cur];++j) if(g[cur][j]->cap==0){
g[cur][j]->type=0;g[cur][j]->cap=1;
int to=g[cur][j]->to;
cur=to;
}
while(cur!=sink){
int time=cur/2/Nst,st=cur/2%Nst;
for(j=0;j<deg[cur];++j) if(g[cur][j]->cap==0 && g[cur][j]->type==1){
g[cur][j]->type=0;g[cur][j]->cap=1;
cur=g[cur][j]->to;
break;
}
for(j=0;j<deg[cur];++j) if(g[cur][j]->cap==0 && g[cur][j]->type==1){
g[cur][j]->cap=1;g[cur][j]->type=0;
int to=g[cur][j]->to;
if(to==sink){
outroute[time][i][0]=st;
outroute[time][i][1]=-1;
cur=sink;
break;
}
if(to/2%Nst==st && time+1==to/2/Nst){
outroute[time][i][0]=st;
outroute[time][i][1]=-1;
cur=to;
break;
}
int nxtTime=to/2/Nst;
outroute[time][i][0]=st;
outroute[time][i][1]=g[cur][j]->id;
for(k=time+1;k<nxtTime;++k){
outroute[k][i][0]=-1;
outroute[k][i][1]=g[cur][j]->id;
}
cur=to;
break;
}
}
}
}
int main() {
struct timeval tstart, tlast;
input();
gettimeofday(&tstart, NULL);
prep();
max_cost_flow_1(1);
max_cost_flow_2(2);
makepath();
gettimeofday(&tlast, NULL);
printf("...end of execution... time %f\n",
(double)(tlast.tv_sec - tstart.tv_sec)+
(double)(tlast.tv_usec - tstart.tv_usec)/ 1000000);
output();
return 0;
}
JOI夏セミ 最小カット乱択
#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cstring>
#include<cmath>
#include<vector>
#define REP(i,m) for(int i=0;i<m;++i)
#define REPN(i,m,in) for(int i=in;i<m;++i)
#define ALL(t) (t).begin(),(t).end()
#define fr first
#define sc second
#define mp make_pair
#define pb push_back
#define dump(x) cerr<<#x<<" = "<<x<<endl
#define prl cerr<<"LINE:"<<__LINE__<<" is called"<<endl
using namespace std;
template<class T> void debug(T a,T b){ for(;a!=b;++a) cerr<<*a<<' ' ;cerr<<endl;}
const int INF=5e8;
typedef pair<int,int> pi;
static const int MAX_V=1005;
int n;
int mat[MAX_V][MAX_V],temp[MAX_V][MAX_V];
int memo[30][MAX_V][MAX_V];
int pos[MAX_V],sum[MAX_V];
int rec(int N,int dep){
if(N==2){
return mat[0][1];
}
REP(i,N) REP(j,N) memo[dep][i][j]=mat[i][j];
int res=INF;
REP(hoge,2){
REP(i,N) REP(j,N) mat[i][j]=memo[dep][i][j];
int h=N/sqrt(2.0),esnum=0;
REP(i,N){
REP(j,N) esnum+=mat[i][j];
sum[i]=esnum;
}
REP(i,N) pos[i]=i;
for(int i=N;i>h;--i){
if(esnum==0) return 0;
int e=rand()%esnum;
int s=INF,t=INF;
REP(j,i) if(sum[pos[j]]>e){
s=j;
break;
}
if(s>0) e-=sum[pos[s-1]];
int addit=0,sub;
REP(j,i){
if(mat[pos[s]][pos[j]]>e){
sub=mat[pos[s]][pos[j]];
t=j;break;
}else e-=mat[pos[s]][pos[j]];
}
if(s>t) swap(s,t);
REP(j,i) if(j!=s && j!=t){
mat[pos[j]][pos[s]]+=mat[pos[j]][pos[t]];
mat[pos[s]][pos[j]]+=mat[pos[t]][pos[j]];
addit+=mat[pos[t]][pos[j]];
}
addit-=sub;
sub*=2;
REPN(j,i-1,t) pos[j]=pos[j+1];
REPN(j,t,s) sum[pos[j]]+=addit;
REPN(j,i-1,t) sum[pos[j]]-=sub;
esnum-=sub;
}
REP(i,N) REP(j,N) temp[i][j]=mat[i][j];
REP(i,h) REP(j,h) mat[i][j]=temp[pos[i]][pos[j]];
res=min(res,rec(h,dep+1));
}
return res;
}
#include<ctime>
int firstState[MAX_V][MAX_V];
int main(){
clock_t start=clock(),end;
scanf("%d",&n);
REP(i,n) REP(j,n) scanf("%d",&mat[i][j]);
int res=INF,logN=0;
while((1<<logN)<n) ++logN;
REP(i,n) REP(j,n) firstState[i][j]=mat[i][j];
REP(hoge,logN){
REP(i,n) REP(j,n) mat[i][j]=firstState[i][j];
int temp=rec(n,0);
printf("%d\n",temp);
res=min(res,temp);
}
printf("%d\n",res);
end=clock();
printf("time is %.5f\n",(end-start)/(double)CLOCKS_PER_SEC);
return 0;
}
