Codeforces 553E Kyoya and Train

解説を読んで実装した。FFTの練習。

  • 問題概要

愚直DPを高速化してください

  • 解法

DPの遷移の係数を求める式を見ると、DP列とある列との畳み込みになっている。FFTして一件落着!と行きたいところだが、DP[i]を求めないとDP[i+1]が求まらないので、オンラインでFFTを求める必要がある。

解説のように、うまいこと領域を分割するとlog T回FFTをすればいいことがわかる。
配列サイズ足りると思ったのに足りなかった。なんでだろう。doubleは遅いのでfloatにしてやっとAC。segtreeと書いてあるものがsegtreeじゃないのはご愛敬。(実装方針が途中で変わった)

#include<stdio.h>
#include<vector>
#include<algorithm>
#include<math.h>
#include<stdlib.h>
using namespace std;
class cd
{
public:
	float x, y;
	cd() { x = 0, y = 0; }
	cd(double x, double y) :x(x), y(y) {}//&#65533;&#65533;&#65533;&#1794;&#65533;&#65533;&#514;&#65533;
	cd operator + (cd p)
	{
		return cd(x + p.x, y + p.y);
	}
	cd operator - (cd p)
	{
		return cd(x - p.x, y - p.y);
	}
	cd operator * (cd p)
	{
		return cd(x*p.x - y*p.y, x*p.y + y*p.x);
	}
};
typedef pair<vector<cd>, vector<cd> >pvcd;
double pi = 3.1415926535897932384626433;
#define LOG 16
class FFT
{
public:
	cd dat[2][1 << LOG];
	cd omega[1 << LOG], romega[1 << LOG];
	int seg[1 << LOG];
	void init(int lg)
	{
		for (int i = 0; i < (1 << lg); i++)omega[i] = cd(cos(pi*i * 2 / double(1 << lg)), sin(pi*i * 2 / double(1 << lg)));
		for (int i = 0; i < (1 << lg); i++)romega[i] = cd(cos(pi*i * 2 / double(1 << lg)), -sin(pi*i * 2 / double(1 << lg)));
		for (int i = 0; i < (1 << lg); i++)
		{
			int t = 0;
			for (int j = 0; j < lg; j++)if (i&(1 << j))t += 1 << (lg - j - 1);
			seg[i] = t;
		}
	}
	vector<cd>transform(vector<cd>v, int lg, bool pal)//&#65533;T&#65533;C&#65533;Y2^lg, pal&#65533;&#847;&#65533;&#65533;&#970;&#65533;&#65533;&#770;&#386;&#65533;true, &#65533;t&#65533;&#970;&#65533;&#65533;&#770;&#386;&#65533;false
	{
		int cur = 0;
		v.resize(1 << lg);
		for (int i = 0; i < (1 << lg); i++)dat[0][seg[i]] = v[i];
		for (int i = lg - 1; i >= 0; i--)
		{
			for (int j = 0; j < (1 << lg); j += (1 << (lg - i)))
			{
				int s = 1 << (lg - i - 1);
				if (pal)
				{
					for (int k = 0; k < s; k++)dat[1 - cur][j + k] = dat[cur][j + k] + dat[cur][j + s + k] * omega[k << i];
					for (int k = 0; k < s; k++)dat[1 - cur][j + s + k] = dat[cur][j + k] - dat[cur][j + s + k] * omega[k << i];
				}
				else
				{
					for (int k = 0; k < s; k++)dat[1 - cur][j + k] = dat[cur][j + k] + dat[cur][j + s + k] * romega[k << i];
					for (int k = 0; k < s; k++)dat[1 - cur][j + s + k] = dat[cur][j + k] - dat[cur][j + s + k] * romega[k << i];
				}
			}
			cur = 1 - cur;
		}
		vector<cd>ret;
		for (int i = 0; i < (1 << lg); i++)ret.push_back(dat[cur][i]);
		return ret;
	}
	pvcd transform_2(vector<double>va, vector<double>vb, int lg, bool pal)//2&#65533;&#136;&#65533;C&#65533;&#65533;transform &#65533;&#65533;&#65533;&#65533;&#65533;&#65533;&#770;&#65533;
	{
		va.resize(1 << lg);
		vb.resize(1 << lg);
		vector<cd>r;
		r.resize(1 << lg);
		for (int i = 0; i < (1 << lg); i++)r[i] = cd(va[i], vb[i]);
		r = transform(r, lg, pal);
		vector<cd>ra, rb;
		ra.resize(1 << lg);
		rb.resize(1 << lg);
		for (int i = 0; i < (1 << lg); i++)
		{
			int j = ((1 << lg) - i)&((1 << lg) - 1);
			ra[i] = cd((r[i] + r[j]).x, (r[i] - r[j]).y) * cd(0.5, 0);
			rb[i] = cd((r[i] + r[j]).y, -(r[i] - r[j]).x) * cd(0.5, 0);
			//printf("%d %lf %lf    %lf %lf  %lf %lf\n", i, r[i].x, r[i].y, ra[i].x, ra[i].y, rb[i].x, rb[i].y);
		}
		return make_pair(ra, rb);
	}
};
#define SIZE 32768
FFT fft[17];
class segtree
{
public:
	double seg[1 << LOG];
	vector<cd>tim[20];
	double now[1 << LOG];
	void init(vector<double>t,int lg)
	{
		int pt = 0;
		for (int d = 0; d <= lg; d++)
		{
			vector<cd>z;
			for (int i = pt; i < (1 << d); i++)z.push_back(cd(t[i], 0));
			pt = 1 << d;
			tim[d] = fft[d].transform(z, d, true);
			//for (int i = 0; i < tim[d].size(); i++)if(fabs(tim[d][i].x)>1e15||fabs(tim[d][i].y)>1e15)printf("%d %d  %d %d  %lf %lf\n", d,i,z.size(),1<<d,tim[d][i].x,tim[d][i].y);
		}
	}
	double add(int a, double t)
	{
		now[a] += t*tim[0][0].x;
		seg[a] = t;
		int s = 1;
		for (int d = 0;; d++)
		{
			vector<cd>z;
			for (int i = a - s + 1; i <= a; i++)z.push_back(cd(seg[i],0));
			//printf("\n\n\n\n");
			//for (int i = 0; i < z.size(); i++)printf("%lf ", z[i]); printf("\n");
			z = fft[d + 1].transform(z, d + 1, true);
			//for (int i = 0; i < z.size(); i++)printf("%lf ", z[i]); printf("\n");
			for (int i = 0; i < z.size(); i++)z[i] = z[i] * tim[d + 1][i];
			//for (int i = 0; i < z.size(); i++)printf("%lf ", z[i]); printf("\n");
			z = fft[d + 1].transform(z, d + 1, false);
			//for (int i = 0; i < z.size(); i++)printf("%lf ", z[i]); printf("\n");
			for (int i = 0; i < z.size(); i++)now[i + a + 1] += z[i].x / double(1 << (d + 1));
			if ((a&(1 << d)) == 0)break;
			s <<= 1;
		}
		return now[a];
	}
};
segtree tree[100];
int dist[50][50];
double cst[100][20001];
double rui[100][20001];
int frm[100], to[100], len[100];
double dp[20001][50];
int main()
{
	/*int num;
	scanf("%d", &num);
	vector<double>va, vb;
	va.push_back(0);
	vb.push_back(0);
	for (int i = 0; i < num; i++)
	{
		int za, zb;
		scanf("%d%d", &za, &zb);
		va.push_back(za);
		vb.push_back(zb);
	}
	fft[0].init(18);
	pvcd r = fft[0].transform_2(va, vb, 18, true);
	for (int i = 0; i < r.first.size(); i++)r.first[i] = r.first[i] * r.second[i];
	r.first = fft[0].transform(r.first, 18, false);
	for (int i = 1; i <= num + num; i++)printf("%d\n", int(floor(r.first[i].x / (1 << 18) + 0.5)));*/
	int num, way, gen, pen;
	scanf("%d%d%d%d", &num, &way, &gen, &pen);
	for (int i = 0; i < num; i++)for (int j = 0; j < num; j++)dist[i][j] = 1000000000;
	for (int i = 0; i < num; i++)dist[i][i] = 0;
	for (int i = 0; i < 17; i++)fft[i].init(i);
	for (int i = 0; i < way; i++)
	{
		int za, zb, zc;
		scanf("%d%d%d", &za, &zb, &zc);
		za--, zb--;
		dist[za][zb] = min(dist[za][zb], zc);
		frm[i] = za, to[i] = zb, len[i] = zc;
		for (int j = 1; j <= gen; j++)scanf("%lf", &cst[i][j]), cst[i][j] /= 100000, rui[i][j] = rui[i][j - 1] + cst[i][j];
		vector<double>z;
		z.resize(1 << 16);
		for (int j = 0; j < gen; j++)z[j] = cst[i][j + 1];
		tree[i].init(z, 16);
	}
	for (int i = 0; i < num; i++)for (int j = 0; j < num; j++)for (int k = 0; k < num; k++)dist[j][k] = min(dist[j][k], dist[j][i] + dist[i][k]);
	for (int i = 0; i < num; i++)dp[0][i] = (double)dist[i][num - 1] + double(pen);
	dp[0][num - 1] = 0;
	for (int x = 0; x < gen; x++)
	{
		for (int i = 0; i < num; i++)dp[x + 1][i] = (double)dist[i][num - 1] + double(pen);
		dp[x + 1][num - 1] = 0;
		for (int i = 0; i < way; i++)
		{
			double z = tree[i].add(x, dp[x][to[i]]);
			dp[x + 1][frm[i]] = min(dp[x + 1][frm[i]], (double)len[i] + z + (double)(rui[i][gen] - rui[i][x + 1])*(dist[to[i]][num - 1] + pen));
		}
		//for (int i = 0; i < num; i++)printf("%lf ", dp[x + 1][i]); printf("\n");///
	}
	printf("%.10lf\n", dp[gen][0]);
}

SRM710 Div1 hard Hyperboxes

  • 問題概要

{1,2,...,N}^K の点を頂点とした K 次元超直方体を重ならないように M 個置く方法は何通りあるか。

  • 解法

超直方体同士が重ならない <-> どこかの次元でその超直方体が占める区間がdisjoint なので、1次元の場合に区間がどういう重なり方(M <= 6 個の区間のそれぞれのペアに対し重なるか重ならないかの高々 2^15 通り)をするものが何通りあるかを前計算した後、例の A and B に足す convolution で二分累乗する。
後半パートは 2^15 * 15 * log 10^9 となる。

後半パートが本質だと思ってなめてかかったら前半パートが計算量のよくわからないいくつかの方針の中で正しいものを選ばないといけないゲームで、そっちに手間取った。
迷走した痕跡を敢えて残してコードを貼っておきます。

#include<stdio.h>
#include<vector>
#include<algorithm>
using namespace std;
typedef long long ll;
ll mod=998244353;
#define LOG 15
class Andadd
{
public:
	ll a[LOG+1][1<<LOG],b[LOG+1][1<<LOG],ret[LOG+1][1<<LOG];
	vector<ll>conv(vector<ll>va,vector<ll>vb)
	{
		for(int i=0;i<va.size();i++)a[0][i]=va[i];
		for(int i=0;i<vb.size();i++)b[0][i]=vb[i];
		for(int i=0;i<LOG;i++)
		{
			for(int j=0;j<(1<<LOG);j+=(1<<(LOG-i)))
			{
				int s=1<<(LOG-i-1);
				for(int k=0;k<s;k++)a[i+1][j+k]=(a[i][j+k]+a[i][j+s+k])%mod,a[i+1][j+s+k]=a[i][j+s+k];
				for(int k=0;k<s;k++)b[i+1][j+k]=(b[i][j+k]+b[i][j+s+k])%mod,b[i+1][j+s+k]=b[i][j+s+k];
			}
		}
		for(int i=0;i<(1<<LOG);i++)ret[LOG][i]=a[LOG][i]*b[LOG][i]%mod;
		for(int i=LOG-1;i>=0;i--)
		{
			for(int j=0;j<(1<<LOG);j+=(1<<(LOG-i)))
			{
				int s=1<<(LOG-i-1);
				for(int k=0;k<s;k++)ret[i][j+k]=(ret[i+1][j+k]+mod-ret[i+1][j+s+k])%mod;
				for(int k=0;k<s;k++)ret[i][j+s+k]=ret[i+1][j+s+k];
			}
		}
		vector<ll>r;
		for(int i=0;i<(1<<LOG);i++)r.push_back(ret[0][i]);
		return r;
	}
};
ll po(ll a, ll b)
{
	if (b == 0)return 1;
	ll z = po(a, b / 2);
	z = z*z%mod;
	if (b & 1)z = z*a%mod;
	return z;
}
ll iiiinv(ll a)
{
	return po(a, mod - 2);
}
ll inv[100];
ll com(ll a,ll b)
{
	if(b<0||b>a)return 0;
	ll r=1;
	for(int i=1;i<=b;i++)r=r*(a-i+1)%mod*inv[i]%mod;
	return r;
}
ll dat[1<<15];
Andadd ad;
vector<ll>ppooww(vector<ll>a,ll b)
{
	if(b==1)return a;
	vector<ll>z=ppooww(a,b/2);
	z=ad.conv(z,z);
	if(b%2==1)z=ad.conv(z,a);
	return z;
}
#include<map>
map<ll,ll>dp[13];
typedef pair<ll,ll>pii;
int dig[6][6];
class Hyperboxes
{
public:
	int findCount(int len,int num,int dim)
	{
		for(int i=1;i<100;i++)inv[i]=iiiinv(i);
		/*vector<int>v;
		for(int i=0;i<num;i++)v.push_back(i+1),v.push_back(i+1);
		for(;;)
		{
			int now=0;
			bool fff=true;
			for(int i=0;i<v.size();i++)
			{
				if(now+1<v[i])fff=false;
				now=max(now,v[i]);
			}
		//	if(fff)
			{
				int d[10];
				fill(d,d+10,0);
				int st[10],go[10];
				vector<int>x;
				for(int i=0;i<v.size();i++)
				{
					if(d[v[i]]==0)d[v[i]]=1,x.push_back(v[i]),st[v[i]]=i;
					else x.push_back(-v[i]),go[v[i]]=i;
				}
				ll dp[13][13][13];
				for(int i=0;i<13;i++)for(int j=0;j<13;j++)for(int k=0;k<13;k++)dp[i][j][k]=0;
				dp[0][0][1]=1;
				for(int i=1;i<x.size();i++)
				{
					int t=0;
					if(x[i-1]<0&&x[i]>0)t=1;
					else if(x[i-1]>0&&x[i]>0&&v[i-1]<v[i])t=1;
					else if(x[i-1]<0&&x[i]<0&&v[i-1]<v[i])t=1;
					if(t==0)
					{
						for(int k=0;k<=num+num;k++)
						{
							for(int j=i-1;j>=0;j--)
							{
								if(v[i]==v[j])break;
								dp[i][j][k]+=dp[i-1][j][k];
							}
						}
					}
					for(int k=0;k<num+num;k++)
					{
						for(int j=0;j<i;j++)
						{
							dp[i][i][k+1]+=dp[i-1][j][k];
						}
					}
				}
				ll sum=0;
				for(int k=0;k<=num+num;k++)
				{
					ll ss=0;
					for(int j=0;j<=num+num;j++)ss+=dp[num+num-1][j][k];
					//printf(" %lld\n",ss);
					sum+=ss*com(len,k);
					sum%=mod;
				}
				int mask=0;
				int pt=0;
				for(int i=1;i<=num;i++)
				{
					for(int j=i+1;j<=num;j++)
					{
						if(!(go[i]<st[j]||go[j]<st[i]))mask+=(1<<pt);
						pt++;
					}
				}
				dat[mask]=(dat[mask]+sum)%mod;
				//printf("%d %lld\n",mask,sum);
			}
			if(!next_permutation(v.begin(),v.end()))break;
		}
	//	for(int i=0;i<2;i++)printf("%lld\n",dat[i]);*/
		int pppt=0;
		for(int i=0;i<num;i++)
		{
			for(int j=i+1;j<num;j++)
			{
				dig[i][j]=dig[j][i]=pppt;
				pppt++;
			}
		}
		dp[0][0]=1;
		for(int i=0;i<=num+num;i++)
		{
			map<ll,ll>::iterator it=dp[i].begin();
			for(;;)
			{
				if(it==dp[i].end())break;
				pii zz=*it;
				it++;
			//	printf("%d %lld %lld %lld\n",i,zz.first>>15,zz.first%(1<<15),zz.second);
				int x[6],mask=zz.first%(1<<15);
				zz.first/=(1<<15);
				for(int j=0;j<num;j++)x[j]=zz.first%3,zz.first/=3;
				vector<int>v0,v1;
				for(int j=0;j<num;j++)
				{
					if(x[j]==0)v0.push_back(j);
					else if(x[j]==1)v1.push_back(j);
				}
				if(v0.size()+v1.size()==0)
				{
					dat[mask]+=zz.second*com(len,i);
					dat[mask]%=mod;
					continue;
				}
				for(int p=0;p<(1<<v0.size());p++)
				{
					int y[6];
					int m=mask;
					for(int j=0;j<num;j++)y[j]=x[j];
					for(int j=0;j<v0.size();j++)
					{
						if(p&(1<<j))
						{
							for(int k=0;k<num;k++)if(y[k]==1)m|=1<<dig[v0[j]][k];
							y[v0[j]]=1;
						}
					}
					for(int q=0;q<(1<<v1.size());q++)
					{
						if(p==0&&q==0)continue;
						int z[6];
						for(int j=0;j<num;j++)z[j]=y[j];
						for(int j=0;j<v1.size();j++)if(q&(1<<j))z[v1[j]]=2;
						ll t=0;
						for(int j=num-1;j>=0;j--)t*=3,t+=z[j];
						t=(t<<15)+m;
						dp[i+1][t]=(dp[i+1][t]+zz.second)%mod;
					}
				}
			}
		}
		vector<ll>ans;
	//	for(int i=0;i<4;i++)printf("%lld\n",dat[i]);
		for(int i=0;i<(1<<15);i++)ans.push_back(dat[i]);
		vector<ll>r=ppooww(ans,dim);
		ll aaa=r[0];
	//	for(int i=1;i<=num;i++)aaa=aaa*i%mod;
		return int(aaa);
	}
};

Codeforces 329E Evil

見た目が解けるものに見えないが、実は解ける。見た目とのギャップがすごい。

  • 問題概要

平面上にN点与えられるので、N点全てを一回ずつ使うサイクルを作って、隣り合う2点間のマンハッタン距離の総和を最大化せよ。

  • 解法

x座標の中央値とy座標の中央値で切って領域を4つに分ける。全部の辺がこの対角線の領域同士を結んでいれば、無条件で最大。

それができない場合、まずNの偶奇で場合分けする。

  • Nが偶数のとき

左側同士を結ぶ辺と右側同士を結ぶ辺を1本ずつ、あるいは上側同士を結ぶ辺と下側同士を結ぶ辺を1本ずつはる必要がある。両方の場合について試してminをとる。「この2点を結ぶとどれくらい損するか」みたいなのを考えると、どこに辺をはればいいかは割と簡単にわかる。

  • Nが奇数で、(x座標の中央値,y座標の中央値)という点がある場合

上と大体同様だが、左側 or 右側 or 上側 or 下側同士を結ぶ辺のうちひとつをはればいい。

  • Nが奇数で、(x座標の中央値,y座標の中央値)という点がない場合

上記の上界を達成するサイクルが構成できる。

以上の場合分けをすると通る。面白い。

#include<stdio.h>
#include<vector>
#include<algorithm>
#include<math.h>
#include<stdlib.h>
using namespace std;
typedef long long ll;
typedef pair<ll, ll>pii;
#define INF 1000000000000000000LL
int main()
{
	int num;
	scanf("%d", &num);
	vector<pii>vx, vy;
	for (int i = 0; i < num; i++)
	{
		int za, zb;
		scanf("%d%d", &za, &zb);
		vx.push_back(make_pair(za, zb));
		vy.push_back(make_pair(zb, za));
	}
	sort(vx.begin(), vx.end());
	sort(vy.begin(), vy.end());
	ll mx = vx[vx.size() / 2].first, my = vy[vy.size() / 2].first;
	ll sx = 0, sy = 0;
	for (int i = 0; i < vx.size(); i++)sx += abs(vx[i].first - mx);
	for (int i = 0; i < vy.size(); i++)sy += abs(vy[i].first - my);
	if (num % 2 == 0)
	{
		int cnt = 0;
		for (int i = 0; i < num / 2; i++)if (make_pair(vx[i].second, vx[i].first) < vy[num / 2])cnt++;
		if (cnt == 0 || cnt == num / 2)
		{
			printf("%I64d\n", (sx + sy) * 2);
		}
		else
		{
			ll m1 = INF, m2 = INF, m3 = INF, m4 = INF;
			for (int i = 0; i < num / 2; i++)m1 = min(m1, mx - vx[i].first);
			for (int i = num / 2; i < num; i++)m2 = min(m2, vx[i].first - mx);
			for (int i = 0; i < num / 2; i++)m3 = min(m3, my - vy[i].first);
			for (int i = num / 2; i < num; i++)m4 = min(m4, vy[i].first - my);
			printf("%I64d\n", (sx + sy) * 2 - min(m1 + m2, m3 + m4) * 2);
		}
	}
	else
	{
		if (make_pair(vx[num / 2].second, vx[num / 2].first) != vy[num / 2])
		{
			printf("%I64d\n", (sx + sy) * 2);
		}
		else
		{
			int cnt = 0;
			for (int i = 0; i < num / 2; i++)if (make_pair(vx[i].second, vx[i].first) < vy[num / 2])cnt++;
			if (cnt == 0 || cnt == num / 2)
			{
				printf("%I64d\n", (sx + sy) * 2);
			}
			else
			{
				ll m1 = INF, m2 = INF, m3 = INF, m4 = INF;
				for (int i = 0; i < num / 2; i++)m1 = min(m1, mx - vx[i].first);
				for (int i = num / 2 + 1; i < num; i++)m2 = min(m2, vx[i].first - mx);
				for (int i = 0; i < num / 2; i++)m3 = min(m3, my - vy[i].first);
				for (int i = num / 2 + 1; i < num; i++)m4 = min(m4, vy[i].first - my);
				printf("%I64d\n", (sx + sy) * 2 - min(min(m1, m2), min(m3, m4)) * 2);
			}
		}
	}
}

Codeforces 429E Points and Segments

前4問は何だったのかという感じの良問。

問題概要: 区間がN個与えられる。この区間たちに 1 or -1 を割り当てて、どの座標でもそれを含む区間の値の和の絶対値が 1 以下になるようなものを構成せよ。

解法: まず全座標がdistinctなときにできるなら適当に順番を決めればdistinctじゃなくてもできる。

全座標がdistinctなとき、条件は区間が偶数個重なっているところ全てで和が 0 であることと同値。座標を前から見て、区間が消えたり現れたりするのを 2 つごとに見ると、

であることが必要十分なことがわかる。

これを満たす割り当てがあるかどうかだが、関係があるところにそれっぽく辺をはると(区間の始点/終点のどちらを見てはられた辺かを考えると)これは偶サイクルっぽい感じ(性格には、「違うパリティ」の辺が偶数本)なことがわかるので、必ず割り当てがある。解けた。

#include<stdio.h>
#include<vector>
#include<algorithm>
using namespace std;
#define SIZE 200000
class unionfind
{
public:
	int par[SIZE];
	int ran[SIZE];
	int ren[SIZE];
	void init()
	{
		for (int i = 0; i<SIZE; i++)
		{
			par[i] = i;
			ran[i] = 0;
			ren[i] = 1;
		}
	}
	int find(int a)
	{
		if (a == par[a])return a;
		else return par[a] = find(par[a]);
	}
	void unite(int a, int b)
	{
		a = find(a);
		b = find(b);
		if (a == b)return;
		if (ran[a]>ran[b])
		{
			par[b] = a;
			ren[a] += ren[b];
		}
		else
		{
			par[a] = b;
			ren[b] += ren[a];
		}
		if (ran[a] == ran[b])ran[b]++;
	}
};
unionfind uf;
typedef pair<int, int>pii;
typedef pair<pii, int>pi3;
bool flag[200000];
int main()
{
	int num;
	scanf("%d", &num);
	vector<pi3>vec;
	for (int i = 0; i < num; i++)
	{
		int za, zb;
		scanf("%d%d", &za, &zb);
		vec.push_back(make_pair(make_pair(za, 0), i));
		vec.push_back(make_pair(make_pair(zb, 1), i));
	}
	sort(vec.begin(), vec.end());
	uf.init();
	for (int i = 0; i < num; i++)
	{
		if (vec[i * 2].first.second == vec[i * 2 + 1].first.second)
		{
			uf.unite(vec[i * 2].second * 2, vec[i * 2 + 1].second * 2 + 1);
			uf.unite(vec[i * 2].second * 2 + 1, vec[i * 2 + 1].second * 2);
		}
		else
		{
			uf.unite(vec[i * 2].second * 2, vec[i * 2 + 1].second * 2);
			uf.unite(vec[i * 2].second * 2 + 1, vec[i * 2 + 1].second * 2 + 1);
		}
	}
	for (int i = 0; i < num; i++)
	{
		if (i * 2 == uf.find(i * 2))
		{
			if (!flag[uf.find(i * 2 + 1)])
			{
				flag[i * 2] = true;
			}
		}
		if (i * 2 + 1 == uf.find(i * 2 + 1))
		{
			if (!flag[uf.find(i * 2)])
			{
				flag[i * 2 + 1] = true;
			}
		}
		//printf("%d %d\n", uf.find(i * 2), uf.find(i * 2 + 1));
	}
	for (int i = 0; i < num; i++)
	{
		if (flag[uf.find(i * 2)])printf("0 ");
		else printf("1 ");
	}
	printf("\n");
}

真っ白な、幻の世界

最近、おかしな夢をよく見る。居眠りをしているときに、よく見る夢だ。

 

その夢は、それを見たのかどうかすらも分からないほどに抽象的で、一切の物理的な動きを、五感の働きを、焦りを、緊迫を含まない。一切の具体的な印象は、残らない。ふと目が覚めると、異世界の冒険から記憶を消して戻ってくることを選択したファンタジーの主人公が、日常の再開の瞬間に一瞬覚えるような、何か強く大きな世界がそのまま消滅したかのような違和感と、そこに大切な何かを置き忘れたかのような虚無感のみが、そこに残されている。

 

人はそれを思い出そうと、記憶にとどめようと懸命になる。しかしそう決意するほどの理性が復活した頃には、それはすでに記憶の彼方に飛んで行っており、些細な手がかりすらも、最早残されてはいない。確実に先ほどまでは存在したであろう「答え」を探したくて、それを言葉にしたくて、でもそう思ったころにはすでにどこを探せばいいのかすら分からなくて。そして、あまりのどうしようもなさから、「答え」など始めからそんなところにはなかったのだと、すべては気のせいで、そんなことは時間の無駄なのだと、決めつけて考えるのをやめてしまうことが、一切の矛盾をはらまなくて。

 

そして矛盾するようだが、今日はその「夢」の中身を、少しだけこちら側へ持ち帰ることができたので、それを紹介しようと思う。

 

背景が真っ白な世界に、どこまでも平行に、決して交わらずにまっすぐに伸びる、鉄道の路線が走っている。ちょうど人間の身体と同じサイズの電車が、望んだときに、望んだ場所に現れる。夢の中の僕は、そのとても速い電車に、好きなところで乗って、好きなだけ進んで、好きなところで降りる。道中も、降りた先も、ただどこまでも壁紙のように真っ白な、無音の世界で、幸せというどこまでも抽象的なものでさえ、具体的な形を帯びているようで。その様子をはるか上から俯瞰する僕が、感覚を同じくする存在として、一切の引いた目線を持たずに、存在していて。

 

夢の中で僕は、なんだってできて。電車はすぐに、目的地に到着して。でもその世界は真っ白だから、具体的な何かをすることは絶対になくて。その世界には具体的なものは何一つ存在しないのだから、矛盾をはらむことすらあり得なくて。ただまっすぐに電車を乗り継いで、その抽象的な、そしておそらく、感情そのものが実体をもって現れたようなその背景を、幸せが手でつかめるかのようなその世界を、ただ無条件に味わって。

 

夢とは所詮、自らの脳の見せる錯覚である。自らの頭の中にないものは、決して現れない。しかし我々は、その夢の中の我々の想像力にしばしば感心する。はっきりした頭では決して考えられないような、支離滅裂で、前提もめちゃくちゃで、しかし面白い世界を、夢は見せてくれる。まるでそこには自分をはるかに超越した存在がいて、昼間の我々が必死で探し求める「答え」を持っており、その片鱗を見せて、それに必死で触れようとする我々を、からかって遊んでいるかのように。

 

それゆえ人は、夢を思い出そうと躍起になる。根拠はないし、おそらくはただ寝ぼけているだけなのに、そこに存在したような気がする「答え」を、取り返そうとする。取り返すこと、それは往々にして言葉にすることだ。そして言葉にしようと準備をすることすら、その脆い世界を粉々にして、跡形もなく消滅させるのに、十分すぎる衝撃となる。

 

人は、正当化を行う生き物だ。そして、人は正当化が上手だ。人を殺すことすらも、納得のできる形に昇華させてしまうほどに。

 

それゆえ人は、答えなどそこには存在しなかったのだ、すべてはただの錯覚だと、納得することにする。そもそも答えが、自分の頭の中にすでに存在したなんて、そんなことはあるわけがないと。そして後には、この居眠りは、そして期待は、何の役にも立たなかったのだという虚無感のみが残る。それでも夢は同じように、人に答えの片割れを、変わらず見せ続けてくる。

 

今日の僕は、その答えの片鱗を、夢から奪ってくることができた。それは答えと呼べるほどはっきりとしたものではなくて。そこで見えたものは、現実の僕の行動を変えるような教えを、何一つ含まなくて。しかし少なくとも、普段の虚無感を、味わわずに済んだ。言葉にできない、記憶の中にすらも存在できないように思えたその世界を堪能した自分を、そこで覚えた、どこまでも強くて、そしてどこまでも脆い感情の欠片を、持って帰ってくることができた。

 

夢とは不確実なものだ。古今東西、あらゆる人が夢の強さに惹かれ、夢を操ろうとした。所詮自分の頭の中で起こっていることなのだ。頭の中で、完結する世界なのだ。意のままに操れても、何らおかしくはないと。それでもその試みには、まだ誰も成功していない。

 

だから僕が、この幾度となく見て、そして初めて持って帰ってくることのできたこの夢を、もう一度見ることがあるのかどうかはわからない。だが少なくとも、これだけは確実に言える。夢の片割れを、一度、ほんの少しでも心に残しておくことができた僕には、もうあの正当化を行う必要も、権利も、残されていないということだ。

 

僕は加藤和也にはなれない

何かに熱心になること、そのことしか考えられないほど熱中すること、そして、純粋に好きという気持ちだけでそれを続けられること。

 

それは素晴らしいことだと思うし、皆そう思っている。寝るのも忘れてバットを振る野球選手は美しいし、ただひたすらに作品を作り続ける伝統工芸の職人は尊敬される。

 

そのあまりにも明確で、強い気持ちは、人を魅了する。そして純粋な羨望の対象となる。それゆえ人はこう言う。「本当に好きな何かを見つけなさい」と。一切の悪意を含まず、純粋に、本当に楽しい人生を送るためのアドバイスとして。それを見つけた自分の覚えるであろう真の楽しさを、究極の喜びを夢想して。

 

もちろん僕もその一人である。僕は、周りより数学ができた。数学オリンピックの日本代表にもなった。今でも、大学の細分化された学問体系の、僕が専攻したいと思っているとある分野では、少なくとも国内の自分と同学年の学生の中では、誰にも負けていないと思っている。

 

そして、僕は数学が好きだ。人が10分で考えるのをやめるところを、3時間考えていられる程度には。そして悲しいことに、普通の人から見ると、それはいつまででも考え続けられることと、区別がつかないのだ。

 

熱心になれること、好きでいられること、それは究極の喜びを夢想する人間にとっての規範である。究極の喜びを得るために、人はそれを好きでなければならない。本の証明が理解できなくても、疲れていたとしても、それをやり続けなければならない。朝から晩まで、そのことばかり考えていなければならない。好きなのだから。好きだと決めたのだから。

 

好きでいるかどうか、熱心でいられるかどうか。それは気持ちの問題だ。それゆえ人は熱心になれない自分を責める。起きてから寝るまで16時間、どこかで必ず気が抜ける自分を責める。休み休みに6時間しか考えられない自分を、たとえ6時間考えられること自体が一般には特別なことであったとしても、心が弱いと責める。16時間、考えられる人は確実に存在するのだから。そういう人についての話を聞くのだから。そして、自分は、本当に好きであるはずのことについて、才能なんかではなくて、気持ちで負けているのだから。

 

加藤和也という数学者がいる。数学に熱中するあまり駅を半裸で歩いて補導されたとか、そういうエピソードのある数学者だ。我々はそのエピソードを聞いて、とりあえず変な人だと笑う。しかし確実に、それを格好いいと感じる。そうなりたいと望む。補導されたいとは思わないにせよ、自分が半裸で外を歩いているという明らかな異常性にすら気付かないほどに、何かに熱中してみたいと思う。そして、自分がそれができる人間であることを示すような、武勇伝を欲しがる。人はその武勇伝のために、自分の姿を見ようと開きかけている目を、強制的に塞ぐことさえする。

 

人は自分を見る目を塞いでまで、熱心になろうともがく。好きでないことは不誠実だ。その対象に、そして自分自身に誠実であるために、現実の自分の姿に目をつぶる。自分をその理想像に投影し、そしてその理想を省みることをしない。理想と現実の乖離を、頑なに認めない。

 

スクリーンには、屈託のない笑顔を浮かべた未来の自分が、何かすごそうな賞を受賞しているのが映っている。その映画館の座席では、作り笑いを浮かべた青年が、ただポップコーンをつまんでいる。本当はポップコーンなどつまんでいる暇はないはずなのに。本当にスクリーンに映っているべきは、すごそうな賞を受賞している自分などではなく、その好きなことを、普段通りに楽しんでいる、普段の自分であるべきなのに。

Codeforces 739D Recover a functional graph

変な問題だけど解法はそこそこきれい。

パス長とサイクル長のペアたちが満たすべき条件は、

  • 各サイクル長に対し、パス長0の点の数はサイクル長の倍数(サイクルをいくつつくるかの条件)
  • 各サイクル長に対し、パス長xの点が存在しないならパス長x+1の点も存在しない

である。この条件を満たすときの構成は容易。

まず、サイクル長が不定でパス長が定まっているような頂点であって、サイクル長もパス長も定まっているどの頂点のパス長よりもパス長が長いものたちは、全部同じサイクルにくっつけるのがよいことがわかる。このときのサイクル長を全探索する。

すると、条件は、「まだパス長やサイクル長が定まっていない頂点たちから、パス長xでサイクル長yの点が何個欲しい」という形で書ける。これを列挙し、パス長やサイクル長を定めていない頂点から適切に辺をはって最大流を流す。欲しい頂点と使える頂点が全部対応づいていればあとは適切に構成すればOK。

パス長もサイクル長も定まっている頂点は何をしてもOKなのでわざわざ辺をはらずに足りないときにとってきて使うことにすると辺の数が減り、計算量は全体でO(N^3)となる。

最初はもうちょっとgreedyでいけると思ってもう少し小さくフローを使おうとしていたが、嘘なことに気付き、そのあと全体をフローに突っ込んでしまえばいいことに気付いた。迷走である。ICPCでやったら破門。

パス長0,正,不定、サイクル長0,不定の6通りを全部違う変数で持とうとしたら実装が大変なことになったのでまとめて書き直した。これでも大変だけどだいぶましになってると思う。

本番でAC1人とかでもまあ時間かければ解けるっぽいけど、時間をかけないと解けないのでダメ。4時間くらいかかった。嘘に走ったり実装方針を立てるのに失敗せずにスパッと最短距離で実装できる能力が欲しいなぁ。

#include<stdio.h>
#include<vector>
#include<algorithm>
#include<string>
#include<iostream>
#include<stdlib.h>
using namespace std;
int get(string s)
{
	if (s == "?")return 400;
	int t = 0;
	for (int i = 0; i < s.size(); i++)t = t * 10 + s[i] - '0';
	return t;
}
typedef pair<int, int>pii;
vector<int>dat[500][500];
vector<int>fr;
vector<int>zdat[500][500];
vector<int>zfr;
int ans[300];
#define SIZE 1000
vector<int>pat[SIZE];
vector<int>cap[SIZE];
vector<int>rev[SIZE];
void adde(int a, int b, int c)
{
	pat[a].push_back(b);
	pat[b].push_back(a);
	cap[a].push_back(c);
	cap[b].push_back(0);
	rev[a].push_back(pat[b].size() - 1);
	rev[b].push_back(pat[a].size() - 1);
}
int frv[SIZE], fre[SIZE];
bool flag[SIZE];
void dfs(int node)
{
	if (flag[node])return;
	flag[node] = true;
	for (int i = 0; i<pat[node].size(); i++)
	{
		if (cap[node][i]>0 && (!flag[pat[node][i]]))
		{
			frv[pat[node][i]] = node;
			fre[pat[node][i]] = i;
			dfs(pat[node][i]);
		}
	}
}
int nod;
int maxflow(int st, int go)
{
	int ret = 0;
	for (;;)
	{
		fill(flag, flag + nod, false);
		dfs(st);
		if (!flag[go])return ret;
		int mini = 1000000000;
		int now = go;
		for (;;)
		{
			if (now == st)break;
			int t = fre[now];
			now = frv[now];
			mini = min(mini, cap[now][t]);
		}
		ret += mini;
		//flow++;
		now = go;
		for (;;)
		{
			if (now == st)break;
			int t = fre[now];
			int nx = now;
			now = frv[now];
			cap[now][t] -= mini;
			cap[nx][rev[now][t]] += mini;
		}
	}
}
int main()
{
	int num;
	scanf("%d", &num);
	int mmxx = 0;
	for (int i = 0; i < num; i++)
	{
		string za, zb;
		cin >> za >> zb;
		int a = get(za), b = get(zb);
		if (a < 400 || b < 400)zdat[b][a].push_back(i);
		else zfr.push_back(i);
		if (b < 400 && a < 400)mmxx = max(mmxx, a);
	}
	vector<pii>adv;
	for (int i = mmxx + 1; i <= num; i++)
	{
		for (int j = 0; j < zdat[400][i].size(); i++)
		{
			adv.push_back(make_pair(i, zdat[400][i][j]));
		}
		zdat[400][i].clear();
	}
	//for (int i = 0; i < adv.size(); i++)printf("....%d %d\n", adv[i].first,adv[i].second);
	for (int p = 1; p <= num; p++)
	{
		for (int i = 0; i < 500; i++)for (int j = 0; j < 500;j++)dat[i][j] = zdat[i][j];
		fr = zfr;
		int sum = 0;
		for (int i = 0; i < 500; i++)sum += dat[p][i].size();
		for (int i = 0; i < adv.size(); i++)dat[p][adv[i].first].push_back(adv[i].second);
		vector<pii>nee;
		for (int i = 1; i <= num; i++)
		{
			int s = 0;
			for (int j = 0; j < 500; j++)s += dat[i][j].size();
			if (s == 0)continue;
			int maxi = 0;
			for (int j = 1; j <= num; j++)if (dat[i][j].size() > 0)maxi = j;
			for (int j = 1; j <= maxi; j++)if (dat[i][j].size() == 0)nee.push_back(make_pair(i, j));
			int x = max(1, int(dat[i][0].size() + i - 1) / i)*i;
			for (int j = 0; j < x - dat[i][0].size(); j++)nee.push_back(make_pair(i, 0));
		}
		if (nee.size() > num)continue;
		for (int i = 0; i < SIZE; i++)
		{
			pat[i].clear();
			rev[i].clear();
			cap[i].clear();
		}
		nod = SIZE;
		for (int i = 0; i < nee.size(); i++)
		{
			adde(i, nee[i].first + 330, 1);
			adde(i, nee[i].second + 660, 1);
		}
		for (int i = 0; i < nee.size(); i++)adde(998, i, 1);
		for (int i = 0; i <= num; i++)adde(i + 330, 999, dat[i][400].size());
		for (int i = 0; i <= num; i++)adde(i + 660, 999, dat[400][i].size());
		int ret = maxflow(998, 999);
		if (ret + fr.size() < nee.size())continue;
		//printf("%d:\n", p);
		for (int i = 0; i < nee.size(); i++)
		{
			//printf("   %d %d\n", nee[i].first, nee[i].second);
			bool f = false;
			for (int j = 0; j < pat[i].size(); j++)
			{
				if (pat[i][j] != 998 && cap[i][j] == 0)
				{
					if (pat[i][j] < 660)
					{
						dat[nee[i].first][nee[i].second].push_back(dat[nee[i].first][400][dat[nee[i].first][400].size() - 1]);
						dat[nee[i].first][400].pop_back();
					}
					else
					{
						dat[nee[i].first][nee[i].second].push_back(dat[400][nee[i].second][dat[400][nee[i].second].size() - 1]);
						dat[400][nee[i].second].pop_back();
					}
					f = true;
					break;
				}
			}
			if (!f)
			{
				dat[nee[i].first][nee[i].second].push_back(fr[fr.size() - 1]);
				fr.pop_back();
			}
		}
		for (int i = 0; i <= num; i++)
		{
			for (int j = 0; j < dat[i][400].size(); j++)dat[i][1].push_back(dat[i][400][j]);
			if (i != 0)for (int j = 0; j < dat[400][i].size(); j++)dat[p][i].push_back(dat[400][i][j]);
			else for (int j = 0; j < dat[400][i].size(); j++)dat[1][i].push_back(dat[400][i][j]);
		}
		for (int i = 0; i < fr.size(); i++)dat[1][0].push_back(fr[i]);
		//for (int i = 0; i <= num; i++)for (int j = 0; j <= num; j++)for (int k = 0; k < dat[i][j].size();k++)printf("%d %d %d\n", i, j, dat[i][j][k]);
		for (int i = 1; i <= num; i++)
		{
			for (int j = 0; j < dat[i][0].size(); j++)
			{
				ans[dat[i][0][j]] = dat[i][0][(j + 1) % i + (j / i)*i];
			}
			for (int j = 1; j <= num; j++)
			{
				for (int k = 0; k < dat[i][j].size(); k++)
				{
					if (dat[i][j - 1].size() == 0)
					{
						goto l01;
					}
					//printf(" %d %d\n", dat[i][j][k], dat[i][j - 1][0]);
					ans[dat[i][j][k]] = dat[i][j - 1][0];
				}
			}
		}
		for (int i = 0; i < num; i++)printf("%d ", ans[i] + 1);
		printf("\n");
		return 0;
	l01:;
	}
	printf("-1\n");
}