二階命題論理の存在量化子を全称量化子で書き換える

TaPL24章のネタです。まずStatementがcoqで正しく書けているか怪しいですが…

Theorem exists_forall_second: forall T: Prop -> Prop, 
((exists X, T X) <-> (forall Y: Prop, (forall X, T X -> Y) -> Y)).
intros. split.
- intros. destruct H. now apply H0 in H.
- intros. apply H. intros. now exists X. Qed.  

ARC152D Halftree

D - Halftree

意外とこの解き方している人が少ないように見えたので。

例えば入力が7 3の場合は
(0,1)(1,2)(2,3)
と辺を張れば0~6の線グラフが得られることが分かると思います。
9 3の場合は頂点7, 8が残ってしまいますが、(7,5)と辺を張ることで、線グラフ+αの形で木が構成できます。
このように、できるだけ線グラフで構築して、残った頂点を適当に線グラフにつなげれば良いです。

提出

2023年の振り返りと2024年目標

修士を修了

特に問題なく、スムーズに卒業できたし、得られた結果も個人的には満足できるものでした。国際会議(といっても査読がすごく厳しいところではない)に出したり、arXivに投稿したりして、当時の自分はそれでいっぱいいっぱいになっていたけど、今見返してみると、アラがあったりして、他人に引用してもらえるようなものではなかったのかなと反省しています。

入社

外資系ソフトウェアエンジニアになりました。一応業務内容はこなせている気がしますが、コミュニケーション能力だったり、バグを報告/修正するというエンジニアとしてのマインドはまだ足りていないと痛感しています。あと同僚の思考スピードについていけていない気がするので、知識を体系化するか、または物事の考え方について考えを改めないといけないかもなと思っています。

趣味

4-6月に圏論、6月からは型システム入門をCoqで進めていました。圏論はAlg-dさんのサイト 圏論 | 壱大整域で学びました。あんまり学部/修士時代、プログラミングの基礎論に興味がなかったのですが、いざやってみると結構面白かったです。引き続き続けていきたいと思います。

目標

仕事: "エンジニアマインド"を身に着けるにします。そうすれば自然とアウトプットが増えるんじゃないかなと願っています。後シンプルに社会人として無理のないまともな生活を送りたいです。
趣味: HaskellとRustに入門するにします。

yukicoder: No.1062 素敵なスコア

https://yukicoder.me/problems/no/1062

 A\lt Bとします。
 K以下の要素を黒丸、 Kよりも大きい要素を白丸で表すことにすると以下のような状況になります。


このうち操作を行った後も同じ要素になるindexは以下の黒枠で囲った部分になります。



このように P,Q両方でswapされずに同じ要素になるindexと、両方でswapされて同じ要素になるindexがあることがわかります。前者の総数は解説にあるように (X^2+Y^2+Z^2) (N-1)!となります。
後者の総数を Sとして、 Sを求めることを考えます。
各index (1 \leq i \leq X) S_i=( P_i,Q_i両方で白丸で、操作で黒丸とswapされて P_i=Q_iとなる数列の総数)が求まれば、 S=2\sum_{i=1}^X S_iとなります。係数が2の理由ですが、○○がswapされる先は●●なので、 P_i,Q_i両方で黒丸のindexも同数だけあるからです。

 S_iの状況は以下のようになります。



これを、


こうじゃ。

なので、 S_i=(後半に●●i個、前半に○●Y個と●●X-i個を好きなように並べた列に○○Z-1個を挿入した列の総数) =X!Y!Z! \times {}_{Y+X-i}C_{Y}\times {}_{N-1}C_{Z-1}となります。よって、

 S=2X!Y!Z! \times {}_{N-1}C_{Z-1} \sum_{i=1}^X {}_{Y+X-i}C_{Y}=2X!Y!Z! \times {}_{N-1}C_{Z-1} \sum_{i=0}^{X-1} {}_{Y+i}C_{Y}

ですが、経路数を少し考えると

 \sum_{i=0}^{X-1} {}_{Y+i}C_{Y}={}_{X+Y}C_{X-1}

であることがわかるので、 S=2X!Y!Z! \times {}_{N-1}C_{Z-1} \times {}_{X+Y}C_{X-1}=\frac{2XZ}{Y+1}(N-1)!となります。
前半パートと合わせれば答えは (X^2+Y^2+Z^2+\frac{2XZ}{Y+1})(N-1)!です。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod = 998244353;

int N, A, B;

int main() {
    cin >> N >> A >> B;
    if(A > B) swap(A, B);
    ll x = A, y = B - A, z = N - B;
    ll fac = 1, nfac = 1;
    for(int i = 1; i <= N - 1; i++) {
        fac = fac * i % mod;
        if(i != y + 1) nfac = nfac * i % mod;
    }
    cout << ((x * x % mod + y * y % mod + z * z % mod) * fac 
            + 2 * x * z % mod * nfac) % mod << endl;
}

Codeforces Round #616 (Div. 1)

https://codeforces.com/contest/1290

ダメダメだったため、撤退…あとでちゃんと復習します。

A
勘違いして変な条件で考えてしまいました…詰めるのもめちゃくちゃ遅かった。例によって条件を整理しきらずに紙に頼ってしまった。

B
adhocだと思えば解けたかもしれない…

C
解説みたらライトを辺と置き換える言い換えすらもできてなかった。(ライトを頂点とした2部グラフで見てた。)
でも解けそうな気はしたが…

codeforcesの問題、ちゃんと考察にadhocな点が1つ2つあって、実装もまあまあ重い感じなんだよな。
AtCoderはadhocに振り切れればいいけど。失敗しているときはだいたい考察が足りてないときではある気がする。頭を動かしましょう。

AOJ2445: MinimumCostPath

http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2445

N<400の時は愚直に求めることにします。
N>=400の時は、
左下、右上の100x100領域に注目します。
(0,0)->(左下の領域)->(真ん中の領域)->(右上の領域)->(N-1,N-1)と移動していくことになりますが、
真ん中の領域では必ず右または上のみの移動になります。

これを証明するためには、任意の経路を経路長を伸ばすことなく、真ん中の領域で右または上のみ移動する経路に変換できることを示せばよいです。
左下の領域から初めて出た場所を(x,100)とします。対称性より(100,y)の時も同様です。
0<=y<100において、障害物がない行が必ず存在します。その行を通って左下の領域から出ると、右または上のみの移動で、右上の領域に入ることが出来ます。
このとき経路長は伸びません。

この事実により、左下の領域から出た座標と右上の領域に入った座標について、最短経路の本数を求めれば良くなり、これは
https://atcoder.jp/contests/dp/tasks/dp_y
です。最後に通った障害物で包除する感じです。

int dx[4] = {1, -1, 0, 0};
int dy[4] = {0, 0, 1, -1};

int N, M;
int X[110], Y[110];
int B[510][510];
ll D[510][510];

void pre(int l, int r, int x, int y) {
	int n = r - l;
	rep(i, 0, n) {
		rep(j, 0, n) {
			B[i][j] = inf;
		}
	}
	rep(i, 0, M) {
		if(l <= X[i] && X[i] < r && l <= Y[i] && Y[i] < r) {
			B[X[i] - l][Y[i] - l] = -1;
		}
	}
	memset(D, 0, sizeof(D));
	x -= l;
	y -= l;
	queue<pi> que;
	que.push(pi(x, y)); 
	B[x][y] = 0;
	D[x][y] = 1;
	while(!que.empty()) {
		pi p = que.front(); que.pop();
		rep(i, 0, 4) {
			int nx = p.fst + dx[i], ny = p.sec + dy[i];
			if(0 <= nx && nx < n && 0 <= ny && ny < n && B[nx][ny] != -1) {
				if(B[nx][ny] >= B[p.fst][p.sec] + 1) {
					ADD(D[nx][ny], D[p.fst][p.sec]);
					if(B[nx][ny] > B[p.fst][p.sec] + 1) {
						B[nx][ny] = B[p.fst][p.sec] + 1;
						que.push(pi(nx, ny));
					}
				}
			}
		}
	}
}

ll qre(int x1, int y1, int x2, int y2) {
	vector<pi> vec;
	// debug(x1, y1, x2, y2);
	rep(i, 0, M) {
		if(x1 == X[i] && y1 == Y[i]) return 0;
		if(x2 == X[i] && y2 == Y[i]) return 0;
		if(x1 <= X[i] && X[i] <= x2 && y1 <= Y[i] && Y[i] <= y2) vec.pb(pi(X[i], Y[i]));
	}
	vec.pb(pi(x1, y1));
	vector<ll> dp(sz(vec), 0);
	sort(all(vec));
	rer(i, sz(vec), 0) {
		int a = vec[i].fst, b = vec[i].sec;
		dp[i] = C(x2 - a + y2 - b, x2 - a);
		rep(j, i + 1, sz(vec)) {
			int c = vec[j].fst, d = vec[j].sec;
			ADD(dp[i], mod - dp[j] * C(c - a + d - b, d - b) % mod);
		}
	}
	return dp[0];
}

void solve() {
	int bound = 100;
	cin >> N >> M;
	C_init(2 * N);
	rep(i, 0, M) {
		cin >> X[i] >> Y[i];
		X[i]--; Y[i]--;
	}
	if(N < 400) {
		pre(0, N, 0, 0);
		cout << D[N - 1][N - 1] << "\n";
	}
	else {
		vector<pair<pi, pl>> vl, vr;
		pre(0, bound, 0, 0);
		rep(i, 0, bound) {
			if(B[i][bound - 1] != -1)
			vl.pb(make_pair(pi(i, bound), pl(B[i][bound - 1] + 1, D[i][bound - 1])));
			if(B[bound - 1][i] != -1)
			vl.pb(make_pair(pi(bound, i), pl(B[bound - 1][i] + 1, D[bound - 1][i])));
		}
		pre(N - bound, N, N - 1, N - 1);
		rep(i, 0, bound) {
			if(B[i][0] != -1)
			vr.pb(make_pair(pi(N - bound + i, N - bound - 1), pl(B[i][0] + 1, D[i][0])));
			if(B[0][i] != -1)
			vr.pb(make_pair(pi(N - bound - 1, N - bound + i), pl(B[0][i] + 1, D[0][i])));
		}
		ll md = inf;
		rep(i, 0, sz(vl)) {
			rep(j, 0, sz(vr)) {
				ll v = vl[i].sec.fst + vr[j].sec.fst + (vr[j].fst.sec - vl[i].fst.sec) + (vr[j].fst.fst - vl[i].fst.fst);
				MIN(md, v);
				// debug(v, vl[i], vr[j]);
			}
		}
		// debug(md, vl, vr);
		ll res = 0;
		rep(i, 0, sz(vl)) {
			rep(j, 0, sz(vr)) {
				if(md == vl[i].sec.fst + vr[j].sec.fst + (vr[j].fst.sec - vl[i].fst.sec) + (vr[j].fst.fst - vl[i].fst.fst)) {
					ADD(res, vl[i].sec.sec * vr[j].sec.sec % mod 
							* qre(vl[i].fst.fst, vl[i].fst.sec, vr[j].fst.fst, vr[j].fst.sec) % mod);
				}
			}
		}
		cout << res << "\n";
	}
}

2パートあるだけで実装辛くなる…

AOJ2327: Sky Jump

http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=2327

すげえ〜〜〜〜〜

x座標がXに到達した時のYの最大値と最小値を求めてしまえば、そこの間の値はtを連続的に動かすことにより取ることができるのでYes or Noの判定が出来ます。

最小値はロケットのうちどれか1つだけ使えば達成できます。

最大値を考えます。
max f(t_1,....,t_n) = (Σw_i-1/2 g t_i^2)
where
h(t_1,...,t_n)=0, t_1>=0, t_2>=0, ... ,t_n>=0
です。

fにマイナスを掛けて最小値問題にした後、KKT条件は以下のようになります。

∇(-f)+\mu_1 ∇(-t1)+...+\mu_n ∇(-tn)+\lambda ∇(h)=0
\mu_i t_i = 0
\mu_i >= 0 t_i >= 0
h=0

∇をt_iごとに展開すると、
-w_i+g t_i - \mu_i + \lambda v_i = 0
となります。

t_i,\mu_iを削除することを目標にすると、以下のように変形できます。

Σ(v_i^2/g \mu_i - v_i^2/g \lambda) = X - Σv_i w_i / g
\mu_i = max(v_i \lambda - w_i, 0)
t_i = max((w_i - v_i \lambda)/g,0)

これで\lambdaで条件を書き表すことが出来ました。
u(\lambda)=Σ(v_i^2/g \mu_i - v_i^2/g \lambda)とすると、最初の式はu(\lambda) = X - Σv_i w_i / g となります。
さらにu(\lambda)は単調減少であり、\lambda-> -∞ のときu(\lambda) -> -∞, \lambda -> +∞のときu(\lambda) = -Σv_i w_i / g となるので、条件を満たすような\lambdaはただ1つだけ存在します。
単調性より\lambdaは二分探索で求めることができます。boundですが[-20000, 1000]とすれば十分です。
\lambdaが求まれば\t_i,さらにはfもわかるので最大値が求められます。

int N;
double X, Y;
double V[1010], W[1010];

double f(double lamb) {
	double res = 0;
	rep(i, 0, N) {
		res += max(V[i] * lamb - W[i], 0.0) * V[i] / 9.8 - V[i] * V[i] / 9.8 * lamb;
	}
	return res;
}

void solve() {
	while(true) {
		cin >> N;
		if(N == 0) break;
		rep(i, 0, N) cin >> V[i] >> W[i];
		cin >> X >> Y;
		double sum = X;
		rep(i, 0, N) sum -= V[i] * W[i] / 9.8;
		double ok = -20000, ng = 1000;
		rep(_, 0, 60) {
			double m = (ok + ng) / 2.0;
			if(f(m) >= sum) ok = m;
			else ng = m;
		}
		double mf = 0.0;
		rep(i, 0, N) {
			double t = max((W[i] - V[i] * ok) / 9.8, 0.0);
			mf += W[i] * t - 9.8 * 0.5 * t * t;
		}
		// debug(mf, ok);
		bool ans = true;
		if(mf + eps < Y) ans = false;
		double found = false;
		rep(i, 0, N) {
			double t = X / V[i];
			if(W[i] * t - 9.8 * 0.5 * t * t <= Y) found = true;
		}
		if(!found) ans = false;
		if(ans) cout << "Yes\n";
		else cout << "No\n";
	}
}