DELETED _lib/dataStructure/fenwickTree.cpp
Index: _lib/dataStructure/fenwickTree.cpp
==================================================================
--- _lib/dataStructure/fenwickTree.cpp
+++ _lib/dataStructure/fenwickTree.cpp
@@ -1,30 +0,0 @@
-//-------------------------------------------------------------
-// Fenwick Tree
-//   O(log N) per each operation
-//
-// Verified by
-//   - SRM424 Div1 LV3
-//-------------------------------------------------------------
-
-template<typename T = LL>
-struct FenwickTree
-{
-	vector<T> x;
-	FenwickTree(size_t n, const T& v = T()) : x(n, v) {}
-
-	void incr( int k, const T& a ) { // z[k] += a;
-		for(; k < x.size(); k|=k+1)
-			x[k] += a;
-	}
-
-	T sum(int i, int j) { // z[i]+...+z[j] : inclusive
-		if(i)
-			return sum(0, j) - sum(0, i-1);
-		else {
-			T v = T();
-			for(; j>=0; j=(j&(j+1))-1)
-				v += x[j];
-			return v;
-		}
-	}
-};

DELETED _lib/dataStructure/rmq.cpp
Index: _lib/dataStructure/rmq.cpp
==================================================================
--- _lib/dataStructure/rmq.cpp
+++ _lib/dataStructure/rmq.cpp
@@ -1,53 +0,0 @@
-//-------------------------------------------------------------
-// Range Minimum Query
-//   O(N log N) construction
-//   O(1) per each query
-//   returns the LEFTMOST/ALL minimum index (I hope)
-//
-// Verified by
-//   - SRM337 Div1 LV2
-//   - SRM431 Div1 LV3
-//-------------------------------------------------------------
-
-template<typename T>
-struct RMQ
-{
-	vector< vector<int> > rm;
-	vector<T> d;
-
-	RMQ( const vector<T>& d ) : d(d) {
-		int n = d.size();
-
-		// rm[k][x] = i s.t. d[i] is the minimum in [x, x+2^k)
-		rm.push_back( vector<int>(n) );
-		for(int x=0; x<n; ++x)
-			rm[0][x] = x;
-		for(int k=1; (1<<k)<=n; ++k) {
-			rm.push_back( rm[k-1] );
-			for(int x=0; x+(1<<k-1)<n; ++x)
-				if( d[rm[k][x]] > d[rm[k-1][x + (1<<k-1)]] )
-					rm[k][x] = rm[k-1][x + (1<<k-1)];
-		}
-	}
-
-	int operator()(int L, int R) const { // inclusive
-		int k=0;
-		for(; L+(1<<k) < R-(1<<k)+1; ++k) {}
-		LL i = rm[k][L];
-		LL j = rm[k][R-(1<<k)+1];
-		return (d[i]<=d[j] ? i : j);
-	}
-
-	vector<int> all(int L, int R) const {
-		vector<int> ans;
-		int minValue = d[(*this)(L, R)];
-		while( L <= R ) {
-			int C = (*this)(L, R);
-			if( minValue < d[C] )
-				break;
-			ans.push_back(C);
-			L = C+1;
-		}
-		return ans;
-	}
-};

DELETED _lib/doc/hitori-sharp.pdf
Index: _lib/doc/hitori-sharp.pdf
==================================================================
--- _lib/doc/hitori-sharp.pdf
+++ _lib/doc/hitori-sharp.pdf
cannot compute difference between binary files

DELETED _lib/doc/lb_ub.txt
Index: _lib/doc/lb_ub.txt
==================================================================
--- _lib/doc/lb_ub.txt
+++ _lib/doc/lb_ub.txt
@@ -1,21 +0,0 @@
-lb = lower_bound(k)
-  ... k-1] lb [k ...
-
-
-ub = upper_bound(k)
-  .., k] ub [k+1 ...
-
-
-なので任意の開閉区間をイテレータの [) 区間に直すには
-
-  [A, B]
-    = [ lower_bound(A), upper_bound(B) )
-
-  [A, B)
-    = [ lower_bound(A), lower_bound(B) )
-
-  (A, B]
-    = [ upper_bound(A), upper_bound(B) )
-
-  (A, B)
-    = [ upper_bound(A), lower_bound(B) )

DELETED _lib/doc/nim.txt
Index: _lib/doc/nim.txt
==================================================================
--- _lib/doc/nim.txt
+++ _lib/doc/nim.txt
@@ -1,38 +0,0 @@
-// SRM338 Div1 LV3
-
-定理のステートメント
-
-    対象となるゲーム
-      - 二人ゲーム
-      - 二人とも動かせる手の集合は同じ
-      - 動かせなくなった方が負け
-      - 無限ループしない
-
-    nim
-      - サイズ s1, s2, .., sk の山がある
-      - どれか一つ山を選んで 1 ～ 山サイズ の任意個取り除ける
-      - 最後の山を取った方が勝ち（山がなくなって打つ手が無くなった方が負け）
-
-    *n をサイズ n の山が1個だけの nim とする
-      - *0 は負け
-      - *n は勝ち n>=1
-
-    ゲーム
-      - 状態 G から打てる手で遷移する先を G1, ..., Gk とする
-        G = {G1, ..., Gk} と書く
-
-    等価
-      - 二つのゲームG, Fが等価なことを G ≡ F と書く
-        等価というのは勝ち負け一致よりもっと強い。二つのゲームの和を取ると必敗になること
-        *n と *n は同サイズの山2つのnimは必敗なので等価。同サイズでなければ必勝なので等価でない
-
-    定理
-      - G1≡*n1, ..., Gk≡*nk ならば、G={G1, ..., Gk} ≡ *mex(n1,...,nk)
-        where mex(S) = min(nat \setminus S)
-        帰納的に、全てのゲームはなんらかの *n と等価
-
-
-おまけ
-  @G を G のnim値とする。つまりG≡*@G
-  G = X + Y ならば @G = @X xor @Y
-  山が複数のnimの勝敗は全部xorとった1つ山のnimに等しいのと同じ原理。

DELETED _lib/gcj/!gcj-smpl.aml
Index: _lib/gcj/!gcj-smpl.aml
==================================================================
--- _lib/gcj/!gcj-smpl.aml
+++ _lib/gcj/!gcj-smpl.aml
@@ -1,40 +0,0 @@
-(****************************************************************************
- * Written in
- *    Alice 1.4
- *    http://www.ps.uni-saarland.de/alice/
- ****************************************************************************)
-
-fun solve (n : int, k : int) : bool =
-	let
-		open Word
-		infix 5 << andb
-		val mask = (0wx1 << fromInt n) - 0wx1
-	in
-		(fromInt k) andb mask = mask
-	end
-
-fun main () =
-	let
-		fun readLine () =
-			case TextIO.inputLine TextIO.stdIn of
-			| NONE    => []
-			| SOME(s) => map Int.fromString (String.tokens (fn x => x = #" ") s)
-
-		fun parseOne c =
-			case readLine() of
-			| [SOME n, SOME k] => (c, (n, k))
-			| _                => assert false
-
-		fun spawn solveOne (c, inp) =
-			(c, solve inp)
-
-		fun printOne (c, ans) =
-			print ("Case #" ^ Int.toString (c+1) ^ ": " ^ (if ans then "ON" else "OFF") ^ "\n")
-	in
-		case readLine () of
-		| [SOME t] => List.app printOne (List.map solveOne (List.tabulate (t, parseOne)))
-		| _        => assert false
-	end
-
-do main ()
-do OS.Process.exit 0

DELETED _lib/gcj/!gcj-smpl.as
Index: _lib/gcj/!gcj-smpl.as
==================================================================
--- _lib/gcj/!gcj-smpl.as
+++ _lib/gcj/!gcj-smpl.as
@@ -1,33 +0,0 @@
-/////////////////////////////////////////////////////////////////////////////
-// Written in
-//    ActionScript (using mtasc 1.14)
-//    http://www.mtasc.org/
-/////////////////////////////////////////////////////////////////////////////
-
-class A
-{
-	static function main(mc)
-	{
-		var input_string = "<paste the input data here>";
-		var input = input_string.split("\n");
-		var T = Number(input[0]);
-
-		var output_string = ""
-		for(var C=1; C<=T; ++C)
-		{
-			var theCase = input[C].split(" ");
-			var N = Number(theCase[0]);
-			var K = Number(theCase[1]);
-			output_string += "Case #" + C + ": " + (solve(N,K)?"ON":"OFF") + "\n"
-		}
-
-		_root.createTextField("tf",0,0,0,800,600);
-		_root.tf.text = output_string;
-	}
-
-	static function solve(N, K)
-	{
-		var mask = (1<<N) - 1;
-		return (K & mask) == mask;
-	}
-}

DELETED _lib/gcj/!gcj-smpl.bl
Index: _lib/gcj/!gcj-smpl.bl
==================================================================
--- _lib/gcj/!gcj-smpl.bl
+++ _lib/gcj/!gcj-smpl.bl
@@ -1,133 +0,0 @@
-#############################################################################
-# Written in
-#    Blue 1.7.5
-#    http://www.lechak.info/blue/
-#############################################################################
-
-global FROM_STRING = func {
-	arg s;
-
-	a = [];
-	i = 0;
-	loop {
-		a = a.append( s.substr(i,i+1).num() );
-		(i=i+1) >= s.length() ? return;
-	};
-	return a;
-};
-
-global TO_STRING = func {
-	arg x;
-	lexical s = "";
-	x.map( func{ s = s + this; } );
-	return s;
-};
-
-global LESS = func {
-	arg x;
-	arg y;
-
-	xl = x.length();
-	yl = y.length();
-	(xl < yl) ? return 1;
-	(xl > yl) ? return 0;
-	i = 0;
-	return loop {
-		(i == xl)     ? return 0;
-		(x[i] < y[i]) ? return 1;
-		(x[i] > y[i]) ? return 0;
-		i = i+1;
-	};
-};
-
-global SUB = func {
-	arg x;
-	arg y;
-
-	x.length() > y.length() ? (y = [].resize(x.length()-y.length(),0).merge(y));
-	z = [].resize(x.length(), 0);
-	i = x.length()-1;
-	carry = 0;
-	loop {
-		z[i]  = x[i] - y[i] - carry;
-		carry = z[i] < 0 ? (z[i]=z[i]+10; 1) : 0;
-		(i=i-1) < 0 ? return;
-	};
-	head = 0;
-	loop {
-		head == z.length()-1 ? return;
-		z[head] != 0 ? return;
-		head = head + 1;
-	};
-	return z.slice(head, z.length());
-};
-
-global DBL = func {
-	arg x;
-
-	z = [].resize(x.length(), 0);
-	i = x.length()-1;
-	carry = 0;
-	loop {
-		z[i]  = x[i] + x[i] + carry;
-		carry = z[i] > 9 ? (z[i]=z[i]-10; 1) : 0;
-		((i=i-1) < 0) ? return;
-	};
-	return (carry==1 ? [1].merge(z) : z);
-};
-
-global DIF = func {
-	arg x;
-	arg y;
-	return LESS(x, y) ? SUB(y, x) : SUB(x, y);
-};
-
-global MOD = func {
-	arg x;
-	arg y;
-
-	LESS(x,y) ? return x;
-	z = MOD(x, DBL(y));
-	LESS(z,y) ? return z;
-	return SUB(z,y);
-};
-
-global ISZERO = func {
-	arg x;
-	return x[0] == 0;
-};
-
-global GCD = func {
-	arg x;
-	arg y;
-	return ISZERO(x) ? y : GCD(MOD(y,x),x);
-};
-
-#############################################################################
-
-solve = func {
-	arg N;
-	arg t;
-	t = t.map( func{return FROM_STRING(this);} );
-
-	lexical t0 = t[0];
-	lexical g  = DIF(t0, t[1]);
-	t.slice(2, t.length()).map(func{
-		g = GCD(g, DIF(t0, this));
-	});
-	r = MOD(t[0], g);
-	return TO_STRING( ISZERO(r) ? r : DIF(g,r) );
-};
-
-input = sys.library("streams").stdio().read(1000000000).split("\n").map(func{return this.split(" ");});
-C = input[0][0].num();
-
-CaseID = (args.length() >= 2 ? args[1].num() : 1);
-loop {
-	"Case #".print();
-	CaseID.print();
-	": ".print();
-	solve( input[CaseID][0].num(), input[CaseID].slice(1,input[CaseID].length()) ).print();
-	"\n".print();
-	((CaseID = CaseID+1) > C) ? return;
-};

DELETED _lib/gcj/!gcj-smpl.cl
Index: _lib/gcj/!gcj-smpl.cl
==================================================================
--- _lib/gcj/!gcj-smpl.cl
+++ _lib/gcj/!gcj-smpl.cl
@@ -1,150 +0,0 @@
-/////////////////////////////////////////////////////////////////////////////
-// Writte in
-//    Claire v3.3.37
-//    http://www.claire-language.com/
-/////////////////////////////////////////////////////////////////////////////
-
-D :: 1000000000
-
-long <: ephemeral_object(hi:integer, lo:integer)
-
-[long!(x: integer) : long ->
-	long(hi = x / D, lo = x mod D)
-]
-
-[+(x:long, y:long) : long ->
-	let lolo := x.lo - D + y.lo in
-	(
-		if ( lolo < 0 )
-			long(hi = x.hi + y.hi,     lo = lolo + D)
-		else
-			long(hi = x.hi + y.hi + 1, lo = lolo)
-	)
-]
-
-[-(x:long, y:long) : long ->
-	let lolo := x.lo - y.lo in
-	(
-		if ( lolo < 0 )
-			long(hi = x.hi - y.hi - 1, lo = lolo + D)
-		else
-			long(hi = x.hi - y.hi,     lo = lolo)
-	)
-]
-
-[*(x:long, y:integer) : long ->
-	let r := long!(0) in
-	let c := x in
-		(while (y > 0) (
-			(if (y mod 2 = 1) r :+ c),
-			c := c + c,
-			y :/ 2
-		),
-		r)
-]
-
-[+(x:long, y:integer) : long ->
-	x + long!(y)
-]
-
-[<(x:long, y:long) : boolean ->
-	if (x.hi < y.hi)
-		true
-	else if (x.hi > y.hi)
-		false
-	else
-		x.lo < y.lo
-]
-[>=(x:long, y:long) : boolean -> not(x < y)]
-[<=(x:long, y:long) : boolean -> not(y < x)]
-[>(x:long, y:long) : boolean -> (y < x)]
-
-[>=(x:long, y:integer) : boolean -> not(x < long!(y))]
-[<=(x:long, y:integer) : boolean -> not(long!(y) < x)]
-[<(x:long, y:integer) : boolean -> (x < long!(y))]
-[>(x:long, y:integer) : boolean -> (long!(y) < x)]
-
-[>=(x:integer, y:long) : boolean -> not(long!(x) < y)]
-[<=(x:integer, y:long) : boolean -> not(y < long!(x))]
-[<(x:integer, y:long) : boolean -> (long!(x) < y)]
-[>(x:integer, y:long) : boolean -> (y < long!(x))]
-
-
-[string!(x: long) : string ->
-	if (x.hi > 0)
-		let los := string!(x.lo) in
-			string!(x.hi) /+ make_string(9 - length(los), '0') /+ los
-	else
-		string!(x.lo)
-]
-
-/////////////////////////////////////////////////////////////////////////////
-
-[readLine() : list<integer> ->
-	list<integer>{integer!(s) | s in explode(freadline(stdin), " ")}
-]
-
-[solve(R:integer, k:integer, N:integer, g:list<integer>) ->
-	let dest := make_list(N, 0) in
-	let earn := make_list(N, long!(0)) in
-	(
-		for s in (0 .. N - 1) (
-			let ride := long!(0) in
-			let i    := 0        in
-			(
-				while (i < N)
-				(
-					let ride2 := ride + g[((s + i) mod N) + 1] in
-					(if (ride2 > k)
-						break()
-					else
-						(ride := ride2, i :+ 1))
-				),
-				earn[s + 1] := ride,
-				dest[s + 1] := ((s + i) mod N) + 1
-			)
-		),
-		let firstVisitTime := make_list(N,       -1) in
-		let firstVisitEarn := make_list(N, long!(0)) in
-		let q         := 1        in
-		let totalEarn := long!(0) in
-		let i         := 0        in
-		(
-			(while (i < R) (
-				if (firstVisitTime[q] = -1)
-				(
-					firstVisitTime[q] := i,
-					firstVisitEarn[q] := totalEarn,
-					totalEarn         :+ earn[q],
-					q                 := dest[q],
-					i                 :+ 1
-				)
-				else
-				(
-					let loopSize := i - firstVisitTime[q] in
-					let loopEarn := totalEarn - firstVisitEarn[q] in
-					let rest     := R - i in
-					(
-						totalEarn :+ loopEarn * (rest / loopSize),
-						// clear
-						firstVisitTime := make_list(N, -1),
-						i              := R - (rest mod loopSize)
-					)
-				)
-			)),
-			princ(string!(totalEarn))
-		)
-	)
-]
-
-[main() ->
-	let T := readLine()[1] in
-		for C in (1 .. T)
-		(
-			printf("Case #~A: ", C),
-			let RkN := readLine() in solve(RkN[1], RkN[2], RkN[3], readLine()),
-			princ("\n")
-		)
-]
-
-(main())

DELETED _lib/gcj/!gcj-smpl.cs
Index: _lib/gcj/!gcj-smpl.cs
==================================================================
--- _lib/gcj/!gcj-smpl.cs
+++ _lib/gcj/!gcj-smpl.cs
@@ -1,60 +0,0 @@
-using System;
-using System.Linq;
-
-class C
-{
-	static void Main()
-	{
-		long[] Ts = readLongArray();
-		long   T  = Ts[0];
-
-		for(long C=1; C<=T; ++C)
-		{
-			long[] RkN = readLongArray();
-			long   R = RkN[0];
-			long   k = RkN[1];
-			long   N = RkN[2];
-			long[] g = readLongArray();
-
-			Console.WriteLine("Case #{0}: {1}", C, solveSmall(R,k,N,g));
-		}
-	}
-
-	static long[] readLongArray()
-	{
-		return (from s in Console.ReadLine().Split(' ') select long.Parse(s)).ToArray();
-	}
-
-	static long solveSmall(long R, long k, long N, long[] g)
-	{
-		long totalEarn = 0;
-
-		for(int i=0; i<R; ++i)
-		{
-			long ride = 0;
-			for(int q=0; q<g.Length; ++q)
-				if( ride + g[q] > k )
-				{
-					g = rotate(g, q);
-					break;
-				}
-				else
-				{
-					ride += g[q];
-				}
-			totalEarn += ride;
-		}
-
-		return totalEarn;
-	}
-
-	static long[] rotate(long[] g, int s)
-	{
-		long[] g2 = new long[g.Length];
-		for(int i=s; i<g.Length; ++i)
-			g2[i-s] = g[i];
-		for(int i=0; i<s; ++i)
-			g2[i+g.Length-s] = g[i];
-		return g2;
-	}
-}

DELETED _lib/gcj/!gcj-smpl.vbs
Index: _lib/gcj/!gcj-smpl.vbs
==================================================================
--- _lib/gcj/!gcj-smpl.vbs
+++ _lib/gcj/!gcj-smpl.vbs
@@ -1,39 +0,0 @@
-'
-' Written in VBScript
-'
-
-Function ReadLine()
-	Dim s
-	s = Split(WScript.StdIn.ReadLine, " ")
-	For i = LBound(s) To UBound(s)
-		s(i) = CLng(s(i))
-	Next
-	ReadLine = s
-End Function
-
-Function Gcd(a, b)
-	If a = 0 Then
-		Gcd = b
-	Else
-		Gcd = Gcd(b mod a, a)
-	End If
-End Function
-
-Function Solve(T)
-	Dim g, r
-	g = Abs(T(1) - T(2))
-	For i = 2 To UBound(T)
-		g = Gcd(g, Abs(T(1) - T(i)))
-	Next
-	r = T(1) mod g
-	If r = 0 Then
-		Solve = 0
-	Else
-		Solve = g - r
-	End If
-End Function
-
-C = ReadLine()(0)
-For CaseID = 1 To C
-	WScript.StdOut.WriteLine "Case #" & CaseID & ": " & Solve(ReadLine)
-Next

DELETED _lib/gcj/!gcj-tmpl.cpp
Index: _lib/gcj/!gcj-tmpl.cpp
==================================================================
--- _lib/gcj/!gcj-tmpl.cpp
+++ _lib/gcj/!gcj-tmpl.cpp
@@ -1,100 +0,0 @@
-//-----------------------------------------------------------------------------
-// >>Code Template<< (for Visual C++)
-
-#include <iostream>
-#include <sstream>
-#include <iomanip>
-#include <string>
-#include <vector>
-#include <set>
-#include <map>
-#include <algorithm>
-#include <numeric>
-#include <iterator>
-#include <complex>
-#include <functional>
-#include <queue>
-#include <stack>
-#include <cmath>
-#include <cassert>
-#include <cstring>
-#define  cout os
-using namespace std;
-typedef long long LL;
-typedef complex<double> CMP;
-void END_OF_INPUT_FOR_THIS_TEST_CASE(); // stub for multi-threading
-
-//-----------------------------------------------------------------------------
-// >>Main<<
-
-void case_main( ostream& os )
-{
-	END_OF_INPUT_FOR_THIS_TEST_CASE();
-}
-
-//-----------------------------------------------------------------------------
-// >>Code Template<< (Multi-Thread Solver)
-
-#if 1
-#undef cout
-#include <windows.h>
-#include <process.h>
-
-static const int THREAD_NUM = 2;
-volatile int     g_id;
-int              g_nCase;
-CRITICAL_SECTION g_cs;
-vector<string>   g_output;
-
-unsigned __stdcall thread_main( void* t_id ) {
-	for(;;) {
-		EnterCriticalSection(&g_cs);
-		const int id = ++g_id;
-		if(id>g_nCase) { LeaveCriticalSection(&g_cs); break; }
-		cerr << setw(4) << id << " @ " << (int)t_id << " start" << endl;
-
-		ostringstream ss;
-		ss << "Case #" << id << ": ";
-		case_main( ss );
-
-		EnterCriticalSection(&g_cs);
-		if(g_output.size()<id) g_output.resize(id);
-		g_output[id-1] = ss.str();
-		cerr << setw(4) << id << " @ " << (int)t_id << " end" << endl;
-		LeaveCriticalSection(&g_cs);
-	}
-	return 0;
-}
-
-void END_OF_INPUT_FOR_THIS_TEST_CASE() { LeaveCriticalSection(&g_cs); }
-
-int main() {
-	cin >> g_nCase;
-	string dummy; getline(cin, dummy);
-
-	InitializeCriticalSection(&g_cs);
-	vector<HANDLE> thread;
-	for(int i=0; i<THREAD_NUM; ++i)
-		thread.push_back( (HANDLE)_beginthreadex(0, 0, &thread_main, (void*)i, 0, 0) );
-	WaitForMultipleObjects( thread.size(), &thread[0], TRUE, INFINITE );
-	DeleteCriticalSection(&g_cs);
-
-	copy( g_output.begin(), g_output.end(), ostream_iterator<string>(cout) );
-}
-
-//-----------------------------------------------------------------------------
-// >>Code Template<< (Single-Thread Solver)
-
-#else
-#undef cout
-void END_OF_INPUT_FOR_THIS_TEST_CASE() {}
-int main() {
-	int nCase; cin >> nCase;
-	string dummy; getline(cin, dummy);
-
-	for(int id=1; id<=nCase; ++id) {
-		cout << "Case #" << id << ": ";
-		case_main( cout );
-	}
-}
-#endif

DELETED _lib/gcj/!gcj-tmpl.java
Index: _lib/gcj/!gcj-tmpl.java
==================================================================
--- _lib/gcj/!gcj-tmpl.java
+++ _lib/gcj/!gcj-tmpl.java
@@ -1,40 +0,0 @@
-import java.io.*;
-import java.util.*;
-
-public class <CLASSNAME>
-{
-	public static void main(String[] arg)
-	{
-		Scanner sc = new Scanner(System.in);
-		int T = sc.nextInt();
-		for(int C=1; C<=T; ++C)
-		{
-			System.out.printf("Case #%d: ", C);
-			(new <CLASSNAME>(sc)).caseMain();
-		}
-	}
-
-	Scanner sc;
-	<CLASSNAME>( Scanner sc ) { this.sc = sc; }
-
-	void caseMain()
-	{
-		int n = sc.nextInt();
-		long[] x = new long[n];
-		for(int i=0; i<n; ++i) x[i] = sc.nextInt();
-		long[] y = new long[n];
-		for(int i=0; i<n; ++i) y[i] = sc.nextInt();
-
-		System.out.println(solve(x, y, n));
-	}
-
-	long solve(long[] x, long[] y, int n)
-	{
-		Arrays.sort(x);
-		Arrays.sort(y);
-		long s = 0;
-		for(int i=0; i<n; ++i)
-			s += x[i] * y[n-1-i];
-		return s;
-	}
-}

DELETED _lib/gcj/!gcj-tmpl.lua
Index: _lib/gcj/!gcj-tmpl.lua
==================================================================
--- _lib/gcj/!gcj-tmpl.lua
+++ _lib/gcj/!gcj-tmpl.lua
@@ -1,19 +0,0 @@
-CN = io.stdin:read("*n")
-for C=1, CN, 1 do
-	io.stdout:write( string.format("Case #%d: ", C) )
-
-	N = io.stdin:read("*n")
-	M = io.stdin:read("*n")
-
-	fav = {}
-	for i=1, M, 1 do
-		T = io.stdin:read("*n")
-		fav[i] = {}
-		for j=1, T, 1 do
-			X, Y = io.stdin:read("*n", "*n")
-			fav[i][j] = {X, Y}
-		end
-	end
-
-	io.stdout:write( string.format("%d %d\n", N, M) )
-end

DELETED _lib/gcj/!gcj-tmpl.py
Index: _lib/gcj/!gcj-tmpl.py
==================================================================
--- _lib/gcj/!gcj-tmpl.py
+++ _lib/gcj/!gcj-tmpl.py
@@ -1,11 +0,0 @@
-import sys
-
-CN = int(sys.stdin.next())
-for C in range(1, CN+1):
-	print "Case #%d:"%(C),
-	N = int(sys.stdin.next())
-	M = int(sys.stdin.next())
-	for i in range(0, M):
-		a = map(int, sys.stdin.next().split(" "))
-		print a
-	print N, M

DELETED _lib/geo/area.cpp
Index: _lib/geo/area.cpp
==================================================================
--- _lib/geo/area.cpp
+++ _lib/geo/area.cpp
@@ -1,23 +0,0 @@
-//-------------------------------------------------------------
-// Area of a polygon
-//   O(N)
-//
-// Verified by
-//   - SRM 337 Div1 LV3
-//   - SRM 486 Div1 LV3
-//-------------------------------------------------------------
-
-double outer_prod( pt a, pt b )
-{
-	return (a.real()*b.imag() - b.real()*a.imag())/2;
-}
-
-double area( const vector<pt>& q )
-{
-	double a = 0.0;
-
-	pt o = q[0];
-	for(int i=1; i+1<q.size(); ++i)
-		a += outer_prod(q[i]-o, q[i+1]-o);
-	return abs(a) / 2;
-}

DELETED _lib/geo/ccw.cpp
Index: _lib/geo/ccw.cpp
==================================================================
--- _lib/geo/ccw.cpp
+++ _lib/geo/ccw.cpp
@@ -1,18 +0,0 @@
-//-------------------------------------------------------------
-// ccw
-//
-// Verified by
-//   - SRM 492 Div1 LV1
-//-------------------------------------------------------------
-
-double outer_prod(const CMP& a, const CMP& b) { return imag(conj(a)*b); }
-double inner_prod(const CMP& a, const CMP& b) { return real(conj(a)*b); }
-
-int ccw(const CMP& a, CMP b, CMP c) {
-	b -= a; c -= a;
-	if( outer_prod(b,c) > 0 ) return +1; // counter clockwise
-	if( outer_prod(b,c) < 0 ) return -1; // clockwise
-	if( inner_prod(b,c) < 0 ) return +2; // c--[a--b] on line
-	if( norm(b) < norm(c) )   return -2; // [a--b]--c on line
-	return 0; // [a--c--b] on line
-}

DELETED _lib/geo/circle_3pt.cpp
Index: _lib/geo/circle_3pt.cpp
==================================================================
--- _lib/geo/circle_3pt.cpp
+++ _lib/geo/circle_3pt.cpp
@@ -1,55 +0,0 @@
-
-//-------------------------------------------------------------
-// The circle passing three points
-//
-// Verified by
-//   - AOJ 0010
-//-------------------------------------------------------------
-
-vector<double> solve_linear_eq( int n, vector< vector<double> > M, const vector<double>& V )
-{
-	vector<double> A(V);
-	for(int i=0; i<n; ++i)
-	{
-		// pivot
-		if( M[i][i] == 0 )
-			for(int j=i+1; j<n; ++j)
-				if( M[j][i] != 0 )
-					{swap(M[i], M[j]); swap(A[i], A[j]); break;}
-		if( M[i][i] == 0 )
-			throw "no anser";
-
-		// M[i][i] <-- 1
-		double p = M[i][i];
-		for(int j=i; j<n; ++j)
-			M[i][j] /= p;
-		A[i] /= p;
-
-		// M[*][i] <-- 0
-		for(int j=0; j<n; ++j) if(j!=i)
-		{
-			double r = M[j][i];
-			for(int k=i; k<n; ++k)
-				M[j][k] -= M[i][k] * r;
-			A[j] -= A[i] * r;
-		}
-	}
-	return A;
-}
-
-void circle_3pt( CMP p1, CMP p2, CMP p3, CMP* c, double* r )
-{
-	p2 -= p1;
-	p3 -= p1;
-	// c == p2/2 + A p2 i == p3/2 + B p3 i
-
-	vector< vector<double> > M(2, vector<double>(2));
-	vector<double> V(2);
-	M[0][0] = -p2.imag(); M[0][1] = +p3.imag(); V[0] = (p3.real()-p2.real())/2.0;
-	M[1][0] = +p2.real(); M[1][1] = -p3.real(); V[1] = (p3.imag()-p2.imag())/2.0;
-	V = solve_linear_eq(2, M, V);
-
-	double A=V[0], B=V[1];
-	*c = p1 + p2/2.0 + A * p2 * CMP(0,1);
-	*r = abs(*c - p1);
-}

DELETED _lib/geo/convexHull.cpp
Index: _lib/geo/convexHull.cpp
==================================================================
--- _lib/geo/convexHull.cpp
+++ _lib/geo/convexHull.cpp
@@ -1,48 +0,0 @@
-
-//-------------------------------------------------------------
-// Convex Hull -- Andrew's Monotone Chain
-//   O(N log N)
-//
-// Verified by
-//   - SRM 336 Div1 LV3
-//   - TCO09 Round2 LV2
-//   - SRM 486 Div1 LV3
-//-------------------------------------------------------------
-
-double outer_prod(const CMP& a, const CMP& b) { return imag(conj(a)*b); }
-double inner_prod(const CMP& a, const CMP& b) { return real(conj(a)*b); }
-
-int ccw(const CMP& a, CMP b, CMP c) {
-	b -= a; c -= a;
-	if( outer_prod(b,c) > 0 ) return +1; // counter clockwise
-	if( outer_prod(b,c) < 0 ) return -1; // clockwise
-	if( inner_prod(b,c) < 0 ) return +2; // c--a--b on line
-	if( norm(b) < norm(c) )   return -2; // a--b--c on line
-	return 0;
-}
-
-bool byX( const CMP& a, const CMP& b ) {
-	if( a.real() != b.real() )
-		return a.real() < b.real();
-	return a.imag() < b.imag();
-}
-
-vector<CMP> convex_hull( vector<CMP> p )
-{
-	#define IS_RIGHT <0   // skip on-line verts
-	//#define IS_RIGHT ==-1 // take all
-
-	sort(p.begin(), p.end(), &byX);
-
-	vector<CMP> ans;
-	for(int i=0; i<p.size(); ans.push_back(p[i++])) // left-to-right
-		while( ans.size()>=2 && ccw(ans[ans.size()-2], ans[ans.size()-1], p[i]) IS_RIGHT )
-			ans.pop_back();
-	if( ans.size() == p.size() )
-		return ans;
-	for(int i=p.size()-2; i>=0; ans.push_back(p[i--])) // right-to-left
-		while( ans.size()>=2 && ccw(ans[ans.size()-2], ans[ans.size()-1], p[i]) IS_RIGHT )
-			ans.pop_back();
-	ans.pop_back();
-	return ans;
-}

DELETED _lib/geo/convexHull_old.cpp
Index: _lib/geo/convexHull_old.cpp
==================================================================
--- _lib/geo/convexHull_old.cpp
+++ _lib/geo/convexHull_old.cpp
@@ -1,49 +0,0 @@
-
-//-------------------------------------------------------------
-// Convex Hull -- Gift Wrapping
-//   O(N log N)
-//
-// ToDo : Rewrite by CCW...
-//
-// Verified by
-//   - SRM 336 Div1 LV3
-//-------------------------------------------------------------
-
-bool byY( CMP a, CMP b )
-{
-	if( a.imag() == b.imag() )
-		return a.real() < b.real();
-	return a.imag() < b.imag();
-}
-
-bool byArg( CMP a, CMP b )
-{
-	return arg(a) < arg(b);
-}
-
-bool isRight( CMP a, CMP b, CMP c )
-{
-	return arg((c-b) / (b-a)) < 0;
-}
-
-vector<CMP> convex_hull( vector<CMP> p )
-{
-	vector<CMP>::iterator it = min_element( p.begin(), p.end(), byY );
-	CMP o = *it;
-	p.erase(it);
-	for(int i=0; i<p.size(); ++i)
-		p[i] -= o;
-	sort( p.begin(), p.end(), byArg );
-
-	vector<CMP> q;
-	q.push_back( CMP(0,0) );
-	q.push_back( p[0] );
-	for(int i=1; i<p.size(); ++i) {
-		while( isRight(q[q.size()-2], q[q.size()-1], p[i]) )
-			q.pop_back();
-		q.push_back( p[i] );
-	}
-	for(int i=0; i<q.size(); ++i)
-		q[i] += o;
-	return q;
-}

DELETED _lib/geo/pt_in_poly.cpp
Index: _lib/geo/pt_in_poly.cpp
==================================================================
--- _lib/geo/pt_in_poly.cpp
+++ _lib/geo/pt_in_poly.cpp
@@ -1,28 +0,0 @@
-
-//-------------------------------------------------------------
-// The circle passing three points
-//
-// Verified by
-//   - AOJ 0012 (only triangles)
-//-------------------------------------------------------------
-
-double outer_prod( CMP a, CMP b )
-{
-	return (a.real()*b.imag() - b.real()*a.imag())/2;
-}
-
-bool point_in_polygon( vector<CMP>& ps, CMP p )
-{
-	bool in = false;
-	for(int i=0; i<ps.size(); ++i) {
-		CMP a = ps[i] - p;
-		CMP b = ps[(i+1)%ps.size()] - p;
-		if(a.imag() > b.imag()) swap(a,b);
-		if(a.imag()<=0 && 0<b.imag()) {
-			if( outer_prod(a,b) < 0 )
-				in = !in;
-		}
-		//if( outer_prod(a,b)==0 && inner_prod(a,b)<=0 ) return ON;
-	}
-	return in;
-}

DELETED _lib/graph/__goldberg.cpp
Index: _lib/graph/__goldberg.cpp
==================================================================
--- _lib/graph/__goldberg.cpp
+++ _lib/graph/__goldberg.cpp
@@ -1,126 +0,0 @@
-//----------------------------------------------------------------------------
-// Goldberg(-Tarjan?) �@ (�ő嗬)
-//
-// �e�ʕt���O���t G ��ł� Src ���� Dest �ւ̍ő嗬�ʂ����߂�
-//   - �O���t�͗אڍs��œn�����ƁBG �͔j�󂳂�܂��B
-//   - ���ʂ����łȂ�Flow���̂��~�����Ƃ��� G �̗e�ʂ������Ă镔��������
-//
-// �v�Z�ʂ�
-//   - O(V^3)
-//       # �ؖ��͂܂��������Ă��Ȃ�
-//
-// �A���S���Y���̊T��
-//   ��{�`�͈ȉ��̒ʂ�
-//      1. �Ƃ肠����Src����o�Ă�ӑS���ɗe�ʌ��E�܂ŗ���
-//      2. ���� > ���o �ȃm�[�h�����‚��ăo�����X�����
-//          - ��������Dest�܂œ��B����p�X������Ȃ�A���
-//            Dest�ɋ߂������֐i�ޕӂւ̗��o�ʂ𑝂₷
-//          - �Ȃ��Ȃ�A���Src�֋߂��������痈��ӂ����
-//            �����ʂ����炷
-//      3. 2.���A���o�����X�Ȓ��_���Ȃ��Ȃ�܂ŌJ��Ԃ�
-//   �u���� > ���o �ȃm�[�h�����‚��āv�̏����́A���o���o�����X��ς�������
-//   Queue�ɓ˂�����ł������Ƃ�FIFO�ł��B
-//   �uDest�ɋ߂������v�uSrc�ɋ߂������v����Ɍ����ɕ]������K�v�͂Ȃ��āA
-//   ������x�̋ߎ��ł��\���B���̃R�[�h�ł͂��̋ߎ��l�� d[] �Ŏ����Ă���B
-//   d[] �̕]����K�x�ɂ����Ƃ�����肷��heuristics�ő����Ȃ�炵���B
-//   �Ƃ��������Ȃ��Ǝ��ۏ�͎g�����ɂȂ�Ȃ��炵���B���Ƃł����Ƃ��B
-//----------------------------------------------------------------------------
-
-#include <vector>
-#include <queue>
-#include <limits>
-using namespace std;
-
-typedef int                        Vert;
-typedef int                        Capacity;
-typedef vector< vector<Capacity> > Graph;
-
-Capacity goldberg_tarjan( Graph& G, Vert Src, Vert Dest )
-{
-	vector<Capacity> e( G.size() ); // excess   : ������ - ���o��
-	vector<int>      d( G.size() ); // distance : Dest�ւ̋���||Src�ւ̋���-N�̉���
-	d[Src] = G.size();
-
-	queue<Vert> Q; // e[v]>0 �ȃm�[�h�����Ă����L���[
-
-	// �Ƃ肠����Src����S�͂ŗ���
-	for(int v=0; v!=G.size(); ++v)
-		if( G[Src][v] )
-		{
-			G[v][Src] += G[Src][v];
-			e[v]      += G[Src][v];
-			G[Src][v]  = 0;
-			Q.push(v);
-		}
-
-	// e[v]>0 �ȃm�[�h���Ȃ��Ȃ�܂ŌJ��Ԃ�
-	while( !Q.empty() )
-	{
-		Vert v = Q.front();
-
-		for(int u=0; u!=G.size() && e[v]; ++u)
-			if( G[v][u] && d[v]>d[u] ) // �K�؂ȕ����ɗ��� "PUSH"
-			{
-				Capacity f = min(e[v], G[v][u]);
-				G[v][u] -= f;
-				G[u][v] += f;
-				e[v] -= f;
-				e[u] += f;
-				if( e[u]==f && u!=Src && u!=Dest ) Q.push(u);
-			}
-
-		if( e[v] == 0 ) // �o�����X��ꂽ�̂ŃL���[���珜��
-			Q.pop();
-		else // �o�����X���ĂȂ��̂͂��������̂�d[v]�𒲐����Ă�蒼�� "RELABEL"
-		{
-			Capacity m = numeric_limits<Capacity>::max();
-			for(int u=0; u!=G.size(); ++u)
-				if( G[v][u] )
-					m = min(m, d[u]+1);
-			d[v] = m;
-		}
-	}
-
-	// Dest �ւ̗����ʂ�Ԃ�
-	return e[Dest];
-}
-
-//----------------------------------------------------------------------------
-// PKU 1459
-//----------------------------------------------------------------------------
-
-#include <iostream>
-#include <locale>
-
-int main()
-{
-	cin.imbue( std::locale("C") );
-	for(int n,np,nc,m; cin>>n>>np>>nc>>m;)
-	{
-		// 0: src, 1: dst, 2: ...
-		Graph g(n+2, vector<Capacity>(n+2));
-
-		while(m--)
-		{
-			int u, v, z; char _;
-			cin >> _ >> u >> _ >> v >> _ >> z;
-			g[u+2][v+2] = z;
-		}
-
-		while(np--)
-		{
-			int u, z; char _;
-			cin >> _ >> u >> _ >> z;
-			g[0][u+2] = z;
-		}
-
-		while(nc--)
-		{
-			int u, z; char _;
-			cin >> _ >> u >> _ >> z;
-			g[u+2][1] = z;
-		}
-
-		cout << goldberg_tarjan(g, 0, 1) << endl;
-	}
-}

DELETED _lib/graph/biMatch.cpp
Index: _lib/graph/biMatch.cpp
==================================================================
--- _lib/graph/biMatch.cpp
+++ _lib/graph/biMatch.cpp
@@ -1,40 +0,0 @@
-
-//-------------------------------------------------------------
-// Bipertite max-matching (by duality, min-vert-cover)
-//   O(V E)
-//
-// Verified by
-//   SRM 397 Div1 Lv3
-//-------------------------------------------------------------
-
-typedef int           vert;
-typedef vert          edge;
-typedef vector<edge>  edges;
-typedef vector<edges> graph;
-
-bool augment( graph& G, int v, vector<vert>& matchTo, bool visited[] )
-{
-	for(int i=0; i<G[v].size(); ++i) {
-		vert u = G[v][i];
-		if( visited[u] ) continue;
-		visited[u] = true;
-
-		if( matchTo[u]<0 || augment(G, matchTo[u], matchTo, visited) )
-			{ matchTo[v]=u, matchTo[u]=v; return true; }
-	}
-	return false;
-}
-
-template<int NV>
-int biMatch( graph& G, int L ) // [0,L):left, [L,?):right
-    // only left->right edges are used during computation
-{
-	vector<vert> matchTo(G.size(), -1);
-	int ans = 0;
-	for(vert v=0; v<L; ++v) {
-		bool visited[NV] = {};
-		if( augment(G, v, matchTo, visited) )
-			++ans;
-	}
-	return ans;
-}

DELETED _lib/graph/isGraphical.cpp
Index: _lib/graph/isGraphical.cpp
==================================================================
--- _lib/graph/isGraphical.cpp
+++ _lib/graph/isGraphical.cpp
@@ -1,34 +0,0 @@
-//-------------------------------------------------------------
-// Check the given list can be a degree list of some graph
-//   O(n^2)  (I suspect it can be made faster, though...)
-//
-// Verified by
-//   - SRM 398 Div1 LV3
-//
-//((
-// Havel-Hakimi
-//   If G[0], ..., G[n] (decreasing) is graphical,
-//   then G[1]-1, G[2]-1, ..., G[G[0]]-1, G[G[0]+1], .., G[n]
-//   is also graphical.
-//))
-//-------------------------------------------------------------
-
-bool isGraphical( vector<int> G )
-{
-	sort( G.begin(), G.end() );
-
-	vector<int>::iterator b = lower_bound( G.begin(), G.end(), 1 );
-	vector<int>::iterator e = G.end();
-
-	while( b < e )
-	{
-		int n = *(--e);
-		if( e-b < n )
-			return false;
-		for(vector<int>::iterator i=e-n; i!=e; ++i)
-			--*i;
-		inplace_merge( b, e-n, e );
-		b = lower_bound( G.begin(), G.end(), 1 );
-	}
-	return true;
-}

DELETED _lib/graph/maxFlow.cpp
Index: _lib/graph/maxFlow.cpp
==================================================================
--- _lib/graph/maxFlow.cpp
+++ _lib/graph/maxFlow.cpp
@@ -1,66 +0,0 @@
-
-//-------------------------------------------------------------
-// Dinic's Algorithm
-//   O(V E)
-//
-// G : bidirectional (G[i].has(j) <==> G[j].has(i))
-// F : flow-capacity F[i][j] = Capacity, F[j][i] = 0
-//
-// Verified by
-//   - SRM 399 Div1 LV3
-//   - PKU 1459
-//   - CodeCraft 09 CUTS
-//   - SRM 465 Div1 LV2
-//-------------------------------------------------------------
-
-static const int NV = 512;
-typedef int           flow;
-typedef int           vert;
-typedef vert          edge;
-typedef vector<edge>  edges;
-typedef vector<edges> graph;
-typedef flow          flow_graph[NV][NV];
-
-flow dinic_dfs( graph& G, flow_graph F, vert v, vert D,
-                int LV[], flow flow_in, int blocked[] )
-{
-	flow flow_out = 0;
-	for(int i=0; i!=G[v].size(); ++i) {
-		int u = G[v][i];
-		if( LV[v]+1==LV[u] && F[v][u] ) {
-			flow f = min(flow_in-flow_out, F[v][u]);
-			if( u==D || !blocked[u] && (f=dinic_dfs(G,F,u,D,LV,f,blocked)) ) {
-				F[v][u]  -= f;
-				F[u][v]  += f;
-				flow_out += f;
-				if( flow_in == flow_out ) return flow_out;
-			}
-		}
-	}
-	blocked[v] = (flow_out==0);
-	return flow_out;
-}
-
-flow maxFlow( graph& G, flow_graph F, vert S, vert D )
-{
-	for( flow total=0 ;; ) {
-		int LV[NV] = {0};
-		vector<int> Q(1, S);
-		for(int lv=1; !Q.empty(); ++lv) {
-			vector<int> Q2;
-			for(int i=0; i!=Q.size(); ++i) {
-				edges& ne = G[Q[i]];
-				for(int j=0; j!=ne.size(); ++j)
-					if( F[Q[i]][ne[j]] && !LV[ne[j]] && ne[j]!=S )
-						LV[ne[j]]=lv, Q2.push_back(ne[j]);
-			}
-			Q.swap(Q2);
-		}
-
-		if( !LV[D] )
-			return total;
-
-		int blocked[NV] = {};
-		total += dinic_dfs( G, F, S, D, LV, 0x7fffffff, blocked );
-	}
-}

DELETED _lib/graph/mincostFlow.cpp
Index: _lib/graph/mincostFlow.cpp
==================================================================
--- _lib/graph/mincostFlow.cpp
+++ _lib/graph/mincostFlow.cpp
@@ -1,108 +0,0 @@
-//-------------------------------------------------------------
-// MinCost-MaxFlow
-//   O(??)
-//
-// Verified by
-//   - SRM 487 Div2 LV3
-//   - SRM 491 Div1 LV3
-//-------------------------------------------------------------
-
-#include <iostream>
-#include <string>
-#include <vector>
-#include <map>
-#include <queue>
-using namespace std;
-
-template<typename T>
-class IdGen
-{
-	map<T, int> v2id_;
-	vector<T>   id2v_;
-public:
-	int v2id(const T& v) {
-		if( !v2id_.count(v) ) { v2id_[v] = size(); id2v_.push_back(v); }
-		return v2id_[v];
-	}
-	const T& id2v(int i) const { return id2v_[i]; }
-	int size() const { return id2v_.size(); }
-};
-
-template<typename Vert, typename Cost, typename Flow, int NV=256>
-class MinCostFlow
-{
-	IdGen<Vert> idgen;
-
-	vector<int> G[NV];
-	Flow F[NV][NV];
-	Cost C[NV][NV];
-
-public:
-	void addEdge( Vert s_, Vert t_, Cost c, Flow f )
-	{
-		int s = idgen.v2id(s_), t = idgen.v2id(t_);
-		G[s].push_back(t);
-		G[t].push_back(s);
-		C[s][t] = c;
-		C[t][s] = -c;
-		F[s][t] = f;
-		F[t][s] = 0;
-	}
-
-	pair<Cost, Flow> calc( Vert s_, Vert t_ )
-	{
-		const int N=idgen.size(), S=idgen.v2id(s_), T=idgen.v2id(t_);
-		static const Cost COST_INF = 1e+300; // !!EDIT HERE!!
-		static const Flow FLOW_INF = 0x7fffffff;
-
-		Cost total_cost = 0;
-		Flow total_flow = 0;
-		vector<Cost> h(N, 0); // potential
-		for(Flow RF=FLOW_INF; RF>0; ) // residual flow
-		{
-			// Dijkstra -- find the min-cost path
-			vector<Cost> d(N, COST_INF);  d[S] = 0;
-			vector<int>  prev(N, -1);
-
-			typedef pair< Cost, pair<int,int> > cedge;
-			priority_queue< cedge, vector<cedge>, greater<cedge> > Q;
-			Q.push( cedge(0, make_pair(S,S)) );
-			while( !Q.empty() ) {
-				cedge e = Q.top(); Q.pop();
-				if( prev[e.second.second] >= 0 )
-					continue;
-				prev[e.second.second] = e.second.first;
-
-				int u = e.second.second;
-				for(int i=0; i<G[u].size(); ++i) {
-					int  v = G[u][i];
-					Cost r_cost = C[u][v] + h[u] - h[v];
-					if( F[u][v] > 0 && d[v] > d[u]+r_cost )
-						Q.push( cedge(d[v]=d[u]+r_cost, make_pair(u,v)) );
-				}
-			}
-
-			if( prev[T] < 0 )
-				break; // Finished
-
-			// Run the flow as much as possible
-			Flow f = RF;
-			for(int u=T; u!=S; u=prev[u])
-				f = min(f, F[prev[u]][u]);
-			RF         -= f;
-			total_flow += f;
-
-			for(int u=T; u!=S; u=prev[u])
-			{
-				total_cost    += f * C[prev[u]][u];
-				F[prev[u]][u] -= f;
-				F[u][prev[u]] += f;
-			}
-
-			// Update the potential
-			for(int u=0; u<N; ++u)
-				h[u] += d[u];
-		}
-		return make_pair(total_cost, total_flow);
-	}
-};

DELETED _lib/graph/undAllPairMinCut.cpp
Index: _lib/graph/undAllPairMinCut.cpp
==================================================================
--- _lib/graph/undAllPairMinCut.cpp
+++ _lib/graph/undAllPairMinCut.cpp
@@ -1,99 +0,0 @@
-
-//-------------------------------------------------------------
-// Gomory-Hu
-// All Pair Minimum Cuts for an Undirected Graph
-//   O(V^4)
-//
-// G : undirected (G[i].has(j) <==> G[j].has(i))
-// F : flow-capacity F[i][j] = F[j][i] = Capacity
-//
-// Verified by
-//   - CodeCraft 09 CUTS
-//-------------------------------------------------------------
-
-static const int NV = 512;
-typedef int           flow;
-typedef int           vert;
-typedef vert          edge;
-typedef vector<edge>  edges;
-typedef vector<edges> graph;
-typedef flow          flow_graph[NV][NV];
-
-flow_graph GHF;
-void GomoryHu( graph& G, flow_graph F, flow_graph AllPairMinCut )
-{
-	int N = G.size();
-	vector<vert> parent(N, 0);
-	vector<flow> weight(N, 0);
-
-	// construct the Gomory-Hu Tree
-	for(vert s=1; s<N; ++s)
-	{
-		vert t = parent[s];
-
-		// mincut between s and t
-		memcpy(GHF, F, sizeof(GHF));
-		flow mincut = 0;
-		for(;;) {
-			// bfs
-			vector<vert> prev(N, -1); prev[s]=s;
-			queue<vert> Q; Q.push(s);
-			while( !Q.empty() ) {
-				vert v = Q.front(); Q.pop();
-				for(int i=0; i<G[v].size(); ++i) {
-					vert u = G[v][i];
-					if( prev[u]==-1 && GHF[v][u] )
-						prev[u]=v, Q.push(u);
-				}
-			}
-
-			if( prev[t] == -1 )
-			{
-				// done: separeted to {v|prev[v]!=-1} and {v|prev[v]==-1}
-				// split t's children
-				weight[s] = mincut;
-				for(vert v=0; v<N; ++v)
-					if( parent[v]==t && prev[v]!=-1 && v!=s )
-						parent[v] = s;
-				vert pt = parent[t];
-				if( prev[pt]!=-1 )
-					parent[s]=pt, parent[t]=s, weight[s]=weight[t], weight[t]=mincut;
-				break;
-			}
-
-			// flow...
-			flow cur = 0x7fffffff;
-			for(vert v=t; v!=s; v=prev[v])
-				cur = min(cur, GHF[prev[v]][v]);
-			for(vert v=t; v!=s; v=prev[v])
-				GHF[prev[v]][v] -= cur, GHF[v][prev[v]] += cur;
-			mincut += cur;
-		}
-	}
-
-	// AllPairMinCuts over the G-H tree
-	for(vert s=0; s<N; ++s) {
-		vector<vert> ps_s;
-		for(vert v=s ;; v=parent[v])
-			{ ps_s.push_back(v); if( v==parent[v] ) break; }
-		reverse(ps_s.begin(), ps_s.end());
-
-		for(vert t=s+1; t<N; ++t) {
-			vector<vert> ps_t;
-			for(vert v=t ;; v=parent[v])
-				{ ps_t.push_back(v); if( v==parent[v] ) break; }
-			reverse(ps_t.begin(), ps_t.end());
-
-			vert ca;
-			for(int i=0; i<ps_s.size() && i<ps_t.size() && ps_s[i]==ps_t[i]; ++i)
-				ca = ps_s[i];
-
-			flow s_to_root = 0x7fffffff, t_to_root = 0x7fffffff;
-			for(vert v=s; v!=ca; v=parent[v])
-				s_to_root = min(s_to_root, weight[v]);
-			for(vert v=t; v!=ca; v=parent[v])
-				t_to_root = min(t_to_root, weight[v]);
-			AllPairMinCut[s][t] = AllPairMinCut[t][s] = min(s_to_root, t_to_root);
-		}
-	}
-}

DELETED _lib/graph/undMinCut.cpp
Index: _lib/graph/undMinCut.cpp
==================================================================
--- _lib/graph/undMinCut.cpp
+++ _lib/graph/undMinCut.cpp
@@ -1,68 +0,0 @@
-
-//-------------------------------------------------------------
-// Stoer-Wagner
-// Minimum cost for making the graph unconnected
-//   O(V^3)
-//
-// G : bidirectional (G[v][u] == G[u][v])
-//
-// Verified by
-//   - SRM 340 Div1 LV3
-//-------------------------------------------------------------
-
-typedef int cost;
-typedef int vert;
-typedef vector< vector<cost> > graph;
-
-cost bidi_minCut( graph G, int N )
-{
-	vector<vert> V;
-	for(int v=0; v<N; ++v)
-		V.push_back(v);
-
-	cost result = 0x7fffffff;
-	for(int n=N; n>1; --n)
-	{
-		// invariant:
-		//   the current set of vertices = {V[0] .. V[n-1]}
-
-		// order the vertices
-		//   v0 = 0 (arbitrary),
-		//   v1 = (u that maximizes \Sigma G[v0][u])
-		//   v2 = (u that maximizes \Sigma G[v0][u]+G[v1][u])
-		//   ...
-		//   vn-2
-		//   vn-1 = (maximizes lastCut = \Sigma G[v0][u]+G[vn-2][u])
-		vector<int> vs;
-		cost lastCut = 0;
-		{
-			vector<cost> ws(n, 0);
-			ws[0] = 1;
-
-			for(int i=0; i<n; ++i) {
-				int j = max_element(ws.begin(), ws.end()) - ws.begin();
-				lastCut = ws[j];
-
-				vs.push_back(j);
-				ws[j] = -1;
-				for(int k=0; k<n; ++k)
-					if( ws[k] != -1 )
-						ws[k] += G[V[k]][V[j]];
-			}
-		}
-
-		// update mincut
-		//   lemma: "the min cost for separating vn-2 and vn-1"==lastCut
-		result = min(result, lastCut);
-
-		// reduce the graph (unify vn-2 and vn-1)
-		//   for testing the case vn-2 and vn-1 is connected
-		vert v2=vs[n-2], v1=vs[n-1];
-		for(int i=0; i<n; ++i) {
-			G[V[v2]][V[i]] += G[V[v1]][V[i]];
-			G[V[i]][V[v2]] += G[V[i]][V[v1]];
-		}
-		V.erase( V.begin() + v1 );
-	}
-	return result;
-}

DELETED _lib/numeric/bigintFake.cpp
Index: _lib/numeric/bigintFake.cpp
==================================================================
--- _lib/numeric/bigintFake.cpp
+++ _lib/numeric/bigintFake.cpp
@@ -1,14 +0,0 @@
-
-string add(const string& a, const string& b)
-{
-	int n = max(a.size(), b.size()), carry=0;
-	string c(n, '0');
-	for(int i=0; i<n; ++i) {
-		int x = (a.size()<=i ? 0 : a[a.size()-1-i]-'0')
-		      + (b.size()<=i ? 0 : b[b.size()-1-i]-'0') + carry;
-		c[n-1-i] = char('0'+x%10);
-		carry = x/10;
-	}
-	if( carry ) c = char('0'+carry)+c;
-	return c;
-}

DELETED _lib/numeric/erathos.cpp
Index: _lib/numeric/erathos.cpp
==================================================================
--- _lib/numeric/erathos.cpp
+++ _lib/numeric/erathos.cpp
@@ -1,9 +0,0 @@
-	static const int N = 999999;
-	vector<bool> isp(N+1, true);
-	vector<int> ps;
-	for(int p=2; p<=N; ++p)
-		if( isp[p] ) {
-			ps.push_back(p);
-			for(int q=p+p; q<=N; q+=p)
-				isp[q] = false;
-		}

DELETED _lib/numeric/fft.cpp
Index: _lib/numeric/fft.cpp
==================================================================
--- _lib/numeric/fft.cpp
+++ _lib/numeric/fft.cpp
@@ -1,39 +0,0 @@
-
-//-------------------------------------------------------------
-// Fast Discrete Fourier Transform
-//   O(n log n)
-//   n must be a power of 2.
-//
-// Verified by
-//   - SRM 436 LV3
-//-------------------------------------------------------------
-
-CMP tmp[65536*4];
-
-template<int F>
-void fft_impl( CMP a[], int n, int stride = 1 )
-{
-	if( n > 1 )
-	{
-		CMP *ev=a, *od=a+stride;
-		int s2=stride*2, n2=n/2;
-
-		fft_impl<F>( ev, n2, s2 );
-		fft_impl<F>( od, n2, s2 );
-
-		for(int i=0; i<n; ++i) tmp[i] = ev[s2*(i%n2)] + od[s2*(i%n2)]*polar(1.0, F*2*M_PI*i/n);
-		for(int i=0; i<n; ++i) a[stride*i] = tmp[i];
-	}
-}
-
-void fft( vector<CMP>& a )
-{
-	fft_impl<+1>(&a[0], a.size());
-}
-
-void ifft( vector<CMP>& a )
-{
-	fft_impl<-1>(&a[0], a.size());
-	for(int i=0; i<a.size(); ++i)
-		a[i] /= a.size();
-}

DELETED _lib/numeric/gcd_lcm.cpp
Index: _lib/numeric/gcd_lcm.cpp
==================================================================
--- _lib/numeric/gcd_lcm.cpp
+++ _lib/numeric/gcd_lcm.cpp
@@ -1,28 +0,0 @@
-LL gcd(LL a, LL b)
-{
-	while(a)
-		swap(a, b%=a);
-	return b;
-}
-
-LL lcm(LL a, LL b)
-{
-	return a/gcd(a,b)*b;
-}
-
-// Assumption: (a/b)*b+(a%b) == a
-//   For typical C/C++ compilers, % on negatives satisfies this.
-//   To be sure, replace a%b by a-a/b*b
-
-LL xgcd(LL a, LL b, LL* x, LL* y) // ax+by=g
-{
-	if(b) {
-		LL yy, g = xgcd(b,a%b,&yy,x);
-		*y = yy - a/b**x;
-		return g;
-	}
-	else {
-		*x=1, *y=0;
-		return a;
-	}
-}

DELETED _lib/numeric/linearEq.cpp
Index: _lib/numeric/linearEq.cpp
==================================================================
--- _lib/numeric/linearEq.cpp
+++ _lib/numeric/linearEq.cpp
@@ -1,40 +0,0 @@
-
-//-------------------------------------------------------------
-// Linear Equation Solver
-//   O(n^3)
-//
-// Verified by
-//   - SRM 398 Div1 LV3
-//   - AOJ 0004
-//-------------------------------------------------------------
-
-vector<double> solve_linear_eq( int n, vector< vector<double> > M, const vector<double>& V )
-{
-	vector<double> A(V);
-	for(int i=0; i<n; ++i)
-	{
-		// pivot
-		if( M[i][i] == 0 )
-			for(int j=i+1; j<n; ++j)
-				if( M[j][i] != 0 )
-					{swap(M[i], M[j]); swap(A[i], A[j]); break;}
-		if( M[i][i] == 0 )
-			throw "no anser";
-
-		// M[i][i] <-- 1
-		double p = M[i][i];
-		for(int j=i; j<n; ++j)
-			M[i][j] /= p;
-		A[i] /= p;
-
-		// M[*][i] <-- 0
-		for(int j=0; j<n; ++j) if(j!=i)
-		{
-			double r = M[j][i];
-			for(int k=i; k<n; ++k)
-				M[j][k] -= M[i][k] * r;
-			A[j] -= A[i] * r;
-		}
-	}
-	return A;
-}

DELETED _lib/numeric/linearEq2.cpp
Index: _lib/numeric/linearEq2.cpp
==================================================================
--- _lib/numeric/linearEq2.cpp
+++ _lib/numeric/linearEq2.cpp
@@ -1,71 +0,0 @@
-//-------------------------------------------------------------
-// Linear Equation Solver for Teplitz Matrix
-//   O(n^2)
-//
-// Faster solver for matrices satisfying
-//   M[i][j] == M[i+1][j+1]
-// for all i, j>0
-// Be careful that this implementation does not support the
-// case M[0][0]=M[1][1]=...=0. (zero-division)
-//
-// TODO: support M[0][0]=M[1][1]=...=0
-// QUICKHACK:
-//    find the least degenerated k and 
-//    calc initial vectors by gauss_jordan...
-//
-// Verified by
-//   - (Checked by random data to match with linearEq.cpp)
-//-------------------------------------------------------------
-
-vector<double> solve_teplitz( int n, vector< vector<double> > M, const vector<double>& V )
-{
-	vector<double> f, b, x;
-
-	// invariants
-	//   M|k f|k = (1 0 ... 0)
-	//   M|k b|k = (0 ... 0 1)
-	//   M|k x|k = V|k
-	f.push_back(    1 / M[0][0] );
-	b.push_back(    1 / M[0][0] );
-	x.push_back( V[0] / M[0][0] );
-
-	for(int k=2; k<=n; ++k)
-	{
-		// M|k (f|k-1 0) = (1 0 ... 0 ea)
-		double ea = 0;
-		for(int i=0; i<k-1; ++i)
-			ea += M[k-1][i] * f[i];
-
-		// M|k (0 b|k-1) = (eb 0 ... 0 1)
-		double eb = 0;
-		for(int i=1; i<k; ++i)
-			eb += M[0][i] * b[i-1];
-
-		// u(1 ea) + v(eb 1) = (1 0)
-		//==> u = 1 / (1-ea.eb),  v = -ea /(1-ea.eb)
-		double u = 1/(1-ea*eb), v = -ea/(1-ea*eb);
-		vector<double> ff;
-		for(int i=0; i<k; ++i)
-			ff.push_back( u*(i<k-1 ? f[i] : 0) + v*(i>0 ? b[i-1] : 0) );
-
-		// similarly...
-		u = -eb/(1-ea*eb), v = 1/(1-ea*eb);
-		vector<double> bb;
-		for(int i=0; i<k; ++i)
-			bb.push_back( u*(i<k-1 ? f[i] : 0) + v*(i>0 ? b[i-1] : 0) );
-
-		// update
-		f.swap(ff);
-		b.swap(bb);
-
-		// M|k (x|k 0) (V|k ec)
-		double ec = 0;
-		for(int i=0; i<k-1; ++i)
-			ec += M[k-1][i] * x[i];
-		x.push_back(0);
-		for(int i=0; i<k; ++i)
-			x[i] -= b[i] * (ec-V[k-1]);
-	}
-
-	return x;
-}

DELETED _lib/numeric/matrix.cpp
Index: _lib/numeric/matrix.cpp
==================================================================
--- _lib/numeric/matrix.cpp
+++ _lib/numeric/matrix.cpp
@@ -1,55 +0,0 @@
-
-//-------------------------------------------------------------
-// Matrix Operations
-//
-// Verified by
-//   - SRM342 Div1 LV3
-//   - SRM341 Div1 LV3
-//   - SRM338 Div1 LV2
-//-------------------------------------------------------------
-
-vector<LL> vMul( const vector< vector<LL> >& A, const vector<LL>& B )
-{
-	const int n = A.size();
-
-	vector<LL> C(n);
-	for(int i=0; i<n; ++i)
-		for(int j=0; j<n; ++j)
-			C[i] += A[i][j] * B[j];
-	return C;
-}
-
-vector< vector<LL> > mMul( const vector< vector<LL> >& A, const vector< vector<LL> >& B )
-{
-	const int n = A.size();
-
-	vector< vector<LL> > C(n, vector<LL>(n));
-	for(int i=0; i<n; ++i)
-		for(int j=0; j<n; ++j) {
-			LL Cij = 0;
-			for(int k=0; k<n; ++k)
-				Cij += A[i][k] * B[k][j];
-			C[i][j] = Cij;
-		}
-	return C;
-}
-
-vector< vector<LL> > mPow( vector< vector<LL> > M, LL t ) // t>0
-{
-	vector< vector<LL> > R;
-	for(; t; t>>=1, M=mMul(M,M))
-		if( t&1 )
-			R = (R.empty() ? M : mMul(R,M));
-	return R;
-}
-
-vector< vector<LL> > mAdd( const vector< vector<LL> >& A, const vector< vector<LL> >& B )
-{
-	const int n = A.size();
-
-	vector< vector<LL> > C(n, vector<LL>(n));
-	for(int i=0; i<n; ++i)
-		for(int j=0; j<n; ++j)
-			C[i][j] = A[i][j] + B[i][j];
-	return C;
-}

DELETED _lib/numeric/mebius.cpp
Index: _lib/numeric/mebius.cpp
==================================================================
--- _lib/numeric/mebius.cpp
+++ _lib/numeric/mebius.cpp
@@ -1,36 +0,0 @@
-
-//-------------------------------------------------------------
-// Mebius Function
-//   O( sqrt(N) )
-//
-// mebius(n) = 0      if n=0 mod p*p for some p
-// mebius(n) = (-1)^k otherwise, where k is the number of n's prime factors
-//
-// properties:
-//  1: mebius(nm) = mebius(n) mebius(m)   if gcd(n,m)=1
-//  2: f(n) = Sigma_{d|n} g(d)
-//       <=>
-//     g(n) = Sigma_{d|n} f(d)mebius(n/d)
-//
-// Verified by
-//   - CodeCraft '11 PSTR
-//
-// Further memoization or more might be needed.
-//-------------------------------------------------------------
-
-int mebius(int n)
-{
-	int f = 1;
-	for(int q=2; q*q<=n; ++q)
-		if( n%q == 0 )
-		{
-			if( n%(q*q) == 0 )
-				return 0;
-			n /= q;
-			f *= -1;
-		}
-	if(n>1)
-		f *= -1;
-	return f;
-}
-

DELETED _lib/numeric/modArith.cpp
Index: _lib/numeric/modArith.cpp
==================================================================
--- _lib/numeric/modArith.cpp
+++ _lib/numeric/modArith.cpp
@@ -1,93 +0,0 @@
-
-//-------------------------------------------------------------
-// Modulo Arithmetics
-//
-// Verified by
-//   - TCO10 R3 LV3
-//-------------------------------------------------------------
-
-static const int MODVAL = 1000000007; // op/ depends on primarity of the MODVAL
-struct mint
-{
-	int val;
-	mint():val(0){}
-	mint(int x):val(x%MODVAL) {}
-	mint(LL  x):val(x%MODVAL) {}
-};
-
-mint operator+(mint x, mint y) { return x.val+y.val; }
-mint operator-(mint x, mint y) { return x.val-y.val+MODVAL; }
-mint operator*(mint x, mint y) { return LL(x.val)*y.val; }
-mint POW(mint x, int e) {
-	mint v = 1;
-	for(;e;x=x*x,e>>=1)
-		if(e&1)
-			v=v*x;
-	return v;
-}
-mint operator/(mint x, mint y) { return x * POW(y, MODVAL-2); }
-
-vector<mint> FAC_(1,1);
-void FAC_INIT(int n) { for(int i=1; i<=n; ++i) FAC_.push_back( FAC_.back()*i ); }
-mint FAC(mint x)       { return FAC_[x.val]; }
-mint C(mint n, mint k) { return k.val<0 || n.val<k.val ? 0 : FAC(n) / (FAC(k) * FAC(n-k)); }
-
-/*
-// MODVAL must be a prime!!
-LL GSS(LL k, LL b, LL e) // k^b + k^b+1 + ... + k^e
-{
-	if( b >  e ) return 0;
-	if( k <= 1 ) return k*(e-b+1);
-	return DIV(SUB(POW(k, e+1), POW(k,b)), k-1);
-}
-
-LL Cpascal(LL n, LL k)
-{
-	vector< vector<LL> > c(n+1, vector<LL>(k+1));
-	for(LL nn=1; nn<=n; ++nn)
-		for(LL kk=0; kk<=min(nn,k); ++kk)
-			c[nn][kk] = kk==0 || kk==nn ? 1 : ADD(c[nn-1][kk-1], c[nn-1][kk]);
-	return c[n][k];
-}
-
-vector< vector<LL> > MATMUL(vector< vector<LL> >& a, vector< vector<LL> >& b)
-{
-   int N = a.size();
-   vector< vector<LL> > c(N, vector<LL>(N));
-   for(int i=0; i<N; ++i)
-   for(int j=0; j<N; ++j)
-   for(int k=0; k<N; ++k)
-      c[i][j] = ADD(c[i][j], MUL(a[i][k],b[k][j]));
-   return c;
-}
-
-// works for non-prime MODVAL
-LL GEO(LL x_, LL e) // x^0 + x^1 + ... + x^e-1
-{
-   vector< vector<LL> > v(2, vector<LL>(2));
-   vector< vector<LL> > x(2, vector<LL>(2));
-   v[0][0] = v[1][1] = 1;
-   x[0][0] = x_; x[0][1] = 0;
-   x[1][0] = 1 ; x[1][1] = 1;
-   for(;e;x=MATMUL(x,x),e>>=1)
-      if(e&1)
-         v = MATMUL(v, x);
-   return v[1][0];
-}
-
-// works for non-prime MODVAL
-LL HYP(LL x_, LL e) // e x^0 + (e-1) x^1 + ... + 1 x^e-1 = GEO(x,1)+GEO(x,2)+...+GEO(x,e)
-{
-   vector< vector<LL> > v(3, vector<LL>(3));
-   vector< vector<LL> > x(3, vector<LL>(3));
-   v[0][0] = v[1][1] = v[2][2] = 1;
-   x[0][0] = x_; x[0][1] = 0; x[0][2] = 0;
-   x[1][0] = 1 ; x[1][1] = 1; x[1][2] = 0;
-   x[2][0] = 0 ; x[2][1] = 1; x[2][2] = 1;
-   e++;
-   for(;e;x=MATMUL(x,x),e>>=1)
-      if(e&1)
-         v = MATMUL(v, x);
-   return v[2][0];
-}
-*/

DELETED _lib/numeric/modArithOld.cpp
Index: _lib/numeric/modArithOld.cpp
==================================================================
--- _lib/numeric/modArithOld.cpp
+++ _lib/numeric/modArithOld.cpp
@@ -1,94 +0,0 @@
-
-//-------------------------------------------------------------
-// Modulo Arithmetics
-//
-// Verified by
-//   - SRM 397 Div1 LV2
-//   - SRM 428 Div1 LV2
-//   - CodeCraft 2010 CNTINT DRAW
-//-------------------------------------------------------------
-
-static const LL MODVAL = 1000000007; // must fit in 32-bits
-
-LL ADD(LL x, LL y) { return (x+y)%MODVAL; }
-LL SUB(LL x, LL y) { return (x-y+MODVAL)%MODVAL; }
-LL MUL(LL x, LL y) { return (x*y)%MODVAL; }
-LL POW(LL x, LL e) {
-	LL v = 1;
-	for(;e;x=MUL(x,x),e>>=1)
-		if(e&1)
-			v = MUL(v, x);
-	return v;
-}
-
-// MODVAL must be a prime!!
-LL DIV(LL x, LL y) { return MUL(x, POW(y, MODVAL-2)); }
-
-// MODVAL must be a prime!!
-LL C(LL n, LL k) {
-	LL v = 1;
-	for(LL i=1; i<=k; ++i)
-		v = DIV(MUL(v, n-i+1), i);
-	return v;
-}
-
-// MODVAL must be a prime!!
-LL GSS(LL k, LL b, LL e) // k^b + k^b+1 + ... + k^e
-{
-	if( b >  e ) return 0;
-	if( k <= 1 ) return k*(e-b+1);
-	return DIV(SUB(POW(k, e+1), POW(k,b)), k-1);
-}
-
-
-
-
-LL Cpascal(LL n, LL k)
-{
-	vector< vector<LL> > c(n+1, vector<LL>(k+1));
-	for(LL nn=1; nn<=n; ++nn)
-		for(LL kk=0; kk<=min(nn,k); ++kk)
-			c[nn][kk] = kk==0 || kk==nn ? 1 : ADD(c[nn-1][kk-1], c[nn-1][kk]);
-	return c[n][k];
-}
-
-vector< vector<LL> > MATMUL(vector< vector<LL> >& a, vector< vector<LL> >& b)
-{
-   int N = a.size();
-   vector< vector<LL> > c(N, vector<LL>(N));
-   for(int i=0; i<N; ++i)
-   for(int j=0; j<N; ++j)
-   for(int k=0; k<N; ++k)
-      c[i][j] = ADD(c[i][j], MUL(a[i][k],b[k][j]));
-   return c;
-}
-
-// works for non-prime MODVAL
-LL GEO(LL x_, LL e) // x^0 + x^1 + ... + x^e-1
-{
-   vector< vector<LL> > v(2, vector<LL>(2));
-   vector< vector<LL> > x(2, vector<LL>(2));
-   v[0][0] = v[1][1] = 1;
-   x[0][0] = x_; x[0][1] = 0;
-   x[1][0] = 1 ; x[1][1] = 1;
-   for(;e;x=MATMUL(x,x),e>>=1)
-      if(e&1)
-         v = MATMUL(v, x);
-   return v[1][0];
-}
-
-// works for non-prime MODVAL
-LL HYP(LL x_, LL e) // e x^0 + (e-1) x^1 + ... + 1 x^e-1 = GEO(x,1)+GEO(x,2)+...+GEO(x,e)
-{
-   vector< vector<LL> > v(3, vector<LL>(3));
-   vector< vector<LL> > x(3, vector<LL>(3));
-   v[0][0] = v[1][1] = v[2][2] = 1;
-   x[0][0] = x_; x[0][1] = 0; x[0][2] = 0;
-   x[1][0] = 1 ; x[1][1] = 1; x[1][2] = 0;
-   x[2][0] = 0 ; x[2][1] = 1; x[2][2] = 1;
-   e++;
-   for(;e;x=MATMUL(x,x),e>>=1)
-      if(e&1)
-         v = MATMUL(v, x);
-   return v[2][0];
-}

DELETED _lib/numeric/nextComb.cpp
Index: _lib/numeric/nextComb.cpp
==================================================================
--- _lib/numeric/nextComb.cpp
+++ _lib/numeric/nextComb.cpp
@@ -1,15 +0,0 @@
-//-------------------------------------------------------------
-// Next Combination
-//
-// Verified by
-//   - SRM345 Div1 LV3
-//-------------------------------------------------------------
-
-LL next_combination(LL p)
-{
-	assert( p > 0 );
-	LL lsb = p & -p;
-	LL rem = p + lsb;
-	LL rit = rem & ~p;
-	return rem | (rit/lsb >> 1)-1;
-}

DELETED _lib/numeric/sternBrocot.cpp
Index: _lib/numeric/sternBrocot.cpp
==================================================================
--- _lib/numeric/sternBrocot.cpp
+++ _lib/numeric/sternBrocot.cpp
@@ -1,52 +0,0 @@
-
-//-------------------------------------------------------------
-// Enumerate all positive reduced fractions in a 
-// systematic manner.
-//   [pl/ql, pr/qr] --> [pl/ql, pm/qm] + [pm/qm, pr/qr]
-//                  --> ...
-//
-// Verified by
-//   - CodeCraft 09 FRACTION
-//-------------------------------------------------------------
-
-void sternBrocot(LL pl = 0, LL ql = 1, LL pr = 1, LL qr = 0)
-{
-	LL pm = pl + pr, qm = ql + qr;
-	sternBrocot(pl, ql, pm, qm);
-	sternBrocot(pm, qm, pr, qr);
-}
-
-/*
-bool fleq(LL p1, LL q1, LL p2, LL q2)
-{
-	return p1*q2 <= p2*q1;
-}
-
-void sternBrocot(LL pl = 0, LL ql = 1, LL pr = 1, LL qr = 0)
-{
-	for(;;) {
-		LL pm = pl + pr, qm = ql + qr;
-
-		if( fleq(c,d,pm,qm) ) {
-//			pr = pm, qr = qm; // [pl/ql, pm/qm]
-			LL X = c*ql - d*pl;
-			LL Y = d*pr - c*qr;
-			LL k = max(1LL, Y/X);
-			pr += k*pl;
-			qr += k*ql;
-		}
-		else if( fleq(pm,qm,a,b) ) {
-//			pl = pm, ql = qm; // [pm/qm, pr/qr]
-			LL X = b*pr - a*qr;
-			LL Y = a*ql - b*pl;
-			LL k = max(1LL, Y/X);
-			pl += k*pr;
-			ql += k*qr;
-		}
-		else {
-			cout << pm+qm*base << "/" << qm << endl;
-			break;
-		}
-	}
-}
-*/

DELETED _lib/typical/amap.cpp
Index: _lib/typical/amap.cpp
==================================================================
--- _lib/typical/amap.cpp
+++ _lib/typical/amap.cpp
@@ -1,72 +0,0 @@
-//
-// accumulative map
-//   basically a map<K,V>, but it eliminates duplication lazily
-//   first keep all inserted elements in a vector. only when the buffer
-//   overflows, it eliminates duplication
-//
-
-template<typename K, typename V, int LIMIT = (1<<26)/3/sizeof(pair<K,V>)>
-struct amap
-{
-	typedef typename vector< pair<K,V> >::iterator iterator;
-
-	vector< pair<K,V> > kv;
-	amap() { kv.reserve(LIMIT); }
-
-	void add(const K& k, const V& v)
-	{
-		kv.push_back( make_pair(k,v) );
-		if( kv.size() == LIMIT )
-			normalize();
-	}
-
-	iterator begin()       { return kv.begin(); }
-	iterator end()         { return kv.end(); }
-	void swap( amap& rhs ) { kv.swap(rhs.kv); }
-
-	// Tested: SRM 469 Lv3
-	void normalize()
-	{
-		sort(kv.begin(), kv.end());
-		int i=0;
-		for(int j=0; j<kv.size(); ++i)
-		{
-			int k = j;
-			kv[i] = kv[k];
-			while( ++j<kv.size() && kv[k].first==kv[j].first )
-				kv[i].second += kv[j].second;
-		}
-		kv.resize(i);
-	}
-/*
-	// Not Tested (Prefer First)
-	void normalize()
-	{
-		sort(kv.begin(), kv.end());
-		int i=0;
-		for(int j=0; j<kv.size(); ++i)
-		{
-			int k = j;
-			kv[i] = kv[k];
-			while( ++j<kv.size() && kv[k].first==kv[j].first )
-				{}
-		}
-		kv.resize(i);
-	}
-
-	// Not Tested (Prefer Last)
-	void normalize()
-	{
-		sort(kv.begin(), kv.end());
-		int i=0;
-		for(int j=0; j<kv.size(); ++i)
-		{
-			int k = j;
-			while( ++j<kv.size() && kv[k].first==kv[j].first )
-				{}
-			kv[i] = kv[j-1];
-		}
-		kv.resize(i);
-	}
-*/
-};

DELETED _lib/typical/bfs.cpp
Index: _lib/typical/bfs.cpp
==================================================================
--- _lib/typical/bfs.cpp
+++ _lib/typical/bfs.cpp
@@ -1,26 +0,0 @@
-// code template for BFS
-
-		State start = /*start*/;
-		vector<State> Q( 1, start );
-		set<State>    V; V.insert(start);
-
-		for(int step=0; !Q.empty(); ++step)
-		{
-			vector<State> Qold; Qold.swap(Q);
-			for(int qi=0; qi<Qold.size(); ++qi)
-			{
-				State& cur = Qold[qi];
-				if( /*isGoal(cur)*/ )
-				{
-				}
-
-				foreach(next)
-				{
-					if( !V.count(next) )
-					{
-						V.insert(next);
-						Q.push_back(next);
-					}
-				}
-			}
-		}

DELETED _lib/typical/bitop.cpp
Index: _lib/typical/bitop.cpp
==================================================================
--- _lib/typical/bitop.cpp
+++ _lib/typical/bitop.cpp
@@ -1,7 +0,0 @@
-int bitcnt(LL x)
-{
-	int c = 0;
-	for(; x; x>>=1)
-		c += x&1;
-	return c;
-}

DELETED _lib/typical/comb.cpp
Index: _lib/typical/comb.cpp
==================================================================
--- _lib/typical/comb.cpp
+++ _lib/typical/comb.cpp
@@ -1,18 +0,0 @@
-//-------------------------------------------------------------
-// number of combinations choosing k out of n
-// # you might better consider to use Pascal's triangle
-// # for comb modulo some number...
-//
-// Verified by
-//   - SRM 350 Div1 LV2
-//-------------------------------------------------------------
-
-LL comb(LL n, LL k)
-{
-	k = min(k, n-k);
-
-	LL c = 1;
-	for(LL i=0; i<k; ++i)
-		c *= n-i, c /= i+1;
-	return c;
-}

DELETED _lib/typical/dp.cpp
Index: _lib/typical/dp.cpp
==================================================================
--- _lib/typical/dp.cpp
+++ _lib/typical/dp.cpp
@@ -1,112 +0,0 @@
-// Tested: SRM 454 Lv2
-template<typename T>
-struct DP2
-{
-	const int N1, N2;
-	vector<T> data;
-	DP2(int N1, int N2, const T& t = T())
-		: N1(N1), N2(N2), data(N1*N2, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2)
-		{ return data[ (i1*N2)+i2 ]; }
-	void swap(DP2& rhs)
-		{ data.swap(rhs.data); }
-};
-
-// Tested: Codeforces #13 C
-template<typename T>
-struct DP2x
-{
-	const int N1, N2;
-	vector<T> data;
-	DP2x(int, int N2, const T& t = T())
-		: N1(2), N2(N2), data(N1*N2, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2)
-		{ i1&=1; return data[ (i1*N2)+i2 ]; }
-	void swap(DP2x& rhs)
-		{ data.swap(rhs.data); }
-};
-
-// Tested: SRM 452 Lv3
-template<typename T>
-struct DP3
-{
-	int N1, N2, N3;
-	vector<T> data;
-	DP3(int N1, int N2, int N3, const T& t = T())
-		: N1(N1), N2(N2), N3(N3), data(N1*N2*N3, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2, int i3)
-		{ return data[ ((i1*N2)+i2)*N3+i3 ]; }
-	void swap(DP3& rhs)
-		{ data.swap(rhs.data); }
-};
-
-// Tested: SRM 468 Lv2
-template<typename T>
-struct DP3x
-{
-	int N1, N2, N3;
-	vector<T> data;
-	DP3x(int, int N2, int N3, const T& t = T())
-		: N1(2), N2(N2), N3(N3), data(N1*N2*N3, t) { assert(data.size()*sizeof(T) < (1<<26)); }
-	T& operator()(int i1, int i2, int i3)
-		{ i1&=1; return data[ ((i1*N2)+i2)*N3+i3 ]; }
-	void swap(DP3x& rhs)
-		{ data.swap(rhs.data); }
-};
-
-// Not Tested
-template<typename T>
-struct DP4
-{
-	int N1, N2, N3, N4;
-	vector<T> data;
-	DP4(int N1, int N2, int N3, int N4, const T& t = T())
-		: N1(N1), N2(N2), N3(N3), N4(N4), data(N1*N2*N3*N4, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2, int i3, int i4)
-		{ return data[ (((i1*N2)+i2)*N3+i3)*N4+i4 ]; }
-	void swap(DP4& rhs)
-		{ data.swap(rhs.data); }
-};
-
-// Not Tested
-template<typename T>
-struct DP4x
-{
-	int N1, N2, N3, N4;
-	vector<T> data;
-	DP4x(int, int N2, int N3, int N4, const T& t = T())
-		: N1(2), N2(N2), N3(N3), N4(N4), data(N1*N2*N3*N4, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2, int i3, int i4)
-		{ i1&=1; return data[ (((i1*N2)+i2)*N3+i3)*N4+i4 ]; }
-	void swap(DP4x& rhs)
-		{ data.swap(rhs.data); }
-};
-
-
-// Tested: SRM 351 Lv2
-template<typename T>
-struct DP5
-{
-	int N1, N2, N3, N4, N5;
-	vector<T> data;
-	DP5(int N1, int N2, int N3, int N4, int N5, const T& t = T())
-		: N1(N1), N2(N2), N3(N3), N4(N4), N5(N5), data(N1*N2*N3*N4*N5, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2, int i3, int i4, int i5)
-		{ return data[ ((((i1*N2)+i2)*N3+i3)*N4+i4)*N5+i5 ]; }
-	void swap(DP5& rhs)
-		{ data.swap(rhs.data); }
-};
-
-// Not Tested
-template<typename T>
-struct DP5x
-{
-	int N1, N2, N3, N4, N5;
-	vector<T> data;
-	DP5x(int, int N2, int N3, int N4, int N5, const T& t = T())
-		: N1(2), N2(N2), N3(N3), N4(N4), N5(N5), data(N1*N2*N3*N4*N5, t) { assert(data.size()*sizeof(T)<(1<<26)); }
-	T& operator()(int i1, int i2, int i3, int i4, int i5)
-		{ i1&=1; return data[ ((((i1*N2)+i2)*N3+i3)*N4+i4)*N5+i5 ]; }
-	void swap(DP5x& rhs)
-		{ data.swap(rhs.data); }
-};

DELETED _lib/typical/idGen.cpp
Index: _lib/typical/idGen.cpp
==================================================================
--- _lib/typical/idGen.cpp
+++ _lib/typical/idGen.cpp
@@ -1,21 +0,0 @@
-//-------------------------------------------------------------
-// ID Assignment
-//
-// Verified by
-//   - ACM/ICPC Tokyo 2010 A
-//   - SRM 491 Div1 LV3
-//-------------------------------------------------------------
-
-template<typename T>
-class IdGen
-{
-	map<T, int> v2id_;
-	vector<T>   id2v_;
-public:
-	int v2id(const T& v) {
-		if( !v2id_.count(v) ) { v2id_[v] = size(); id2v_.push_back(v); }
-		return v2id_[v];
-	}
-	const T& id2v(int i) const { return id2v_[i]; }
-	int size() const { return id2v_.size(); }
-};

DELETED _lib/typical/inMap.cpp
Index: _lib/typical/inMap.cpp
==================================================================
--- _lib/typical/inMap.cpp
+++ _lib/typical/inMap.cpp
@@ -1,5 +0,0 @@
-
-bool inMap(int H, int W, int y, int x)
-{
-	return (0<=y && y<H && 0<=x && x<W);
-}

DELETED _lib/typical/maxRectArea.cpp
Index: _lib/typical/maxRectArea.cpp
==================================================================
--- _lib/typical/maxRectArea.cpp
+++ _lib/typical/maxRectArea.cpp
@@ -1,39 +0,0 @@
-//-------------------------------------------------------------
-// h[] : list of heights
-//   i-th rectangle is located at (i,0)--(i+1,h[i])
-//
-// calculate the area of maximum sub-rectangle
-//
-// Verified by
-//   - SRM 337 Div1 LV2
-//-------------------------------------------------------------
-
-		// solve
-		vector<int> left(n);
-		{
-			map<LL, int> h; h[-1] = -1;
-			for(int i=0; i<n; ++i) {
-				// position of the highest building < R[i]
-				map<LL,int>::iterator it = h.lower_bound(R[i]);
-				left[i] = (--it)->second+1;
-				h.erase( ++it, h.end() );
-				h[ R[i] ] = i;
-			}
-		}
-		vector<int> right(n);
-		{
-			map<LL, int> h; h[-1] = n;
-			for(int i=n-1; i>=0; --i) {
-				// position of the highest building < R[i]
-				map<LL,int>::iterator it = h.lower_bound(R[i]);
-				right[i] = (--it)->second-1;
-				h.erase( ++it, h.end() );
-				h[ R[i] ] = i;
-			}
-		}
-		LL ans = 0;
-		for(int i=0; i<n; ++i) {
-			LL area = R[i] * (right[i] - left[i] + 1);
-			ans = max(ans, area);
-		}
-		return ans;

DELETED _lib/typical/setCover.cpp
Index: _lib/typical/setCover.cpp
==================================================================
--- _lib/typical/setCover.cpp
+++ _lib/typical/setCover.cpp
@@ -1,69 +0,0 @@
-
-//-------------------------------------------------------------
-// (Minimum) Set Cover
-//   NP-complete
-//
-// Consider using SetCoverEasy, when mask is downward closed.
-//
-// Verified by
-//   - SRM345 Div1 LV3
-//-------------------------------------------------------------
-
-LL next_combination(LL p)
-{
-	LL lsb = p & -p;
-	LL rem = p + lsb;
-	LL rit = rem & ~p;
-	return rem|(((rit/lsb)>>1)-1);
-}
-
-// Is it possible to cover the goal by at most k elements from mask?
-// Assumption:  goal \subseteq \bigcup mask
-bool canCover( LL goal, set<LL>& mask, int k )
-{
-	// optimizer
-	for(bool update=true; update; ) {
-		update = false;
-
-		// if *it \subseteq *jt for some jt, then we NEVER use it
-		// if *it is the only mask covering some city, we ALWAYS use it
-		for(set<LL>::iterator it=mask.begin(); it!=mask.end(); ) {
-			bool isSubset = false;
-			LL   onlyByMe = *it & goal;
-			for(set<LL>::iterator jt=mask.begin(); jt!=mask.end(); ++jt)
-				if( it!=jt ) {
-					onlyByMe &= ~*jt;
-					isSubset |= (*it & *jt & goal) == (*it & goal);
-				}
-
-			update |= isSubset | !!onlyByMe;
-
-			if( isSubset )
-				mask.erase(it++);
-			else if( onlyByMe ) {
-				if( --k < 0 )
-					return false;
-				goal &= ~*it;
-				mask.erase(it++);
-			}
-			else
-				++it;
-		}
-
-		if( mask.size()<=k || goal==0 )
-			return true;
-	}
-
-	// exhaustive search
-	vector<LL> ms(mask.begin(), mask.end());
-	for(LL i=(1LL<<k)-1; i<(1LL<<ms.size()); i=next_combination(i))
-	{
-		LL gg = goal;
-		for(LL j=0; (1LL<<j)<=i; ++j)
-			if( i & (1<<j) )
-				gg &= ~ms[j];
-		if( gg == 0 )
-			return true;
-	}
-	return false;
-}

DELETED _lib/typical/setCoverEasy.cpp
Index: _lib/typical/setCoverEasy.cpp
==================================================================
--- _lib/typical/setCoverEasy.cpp
+++ _lib/typical/setCoverEasy.cpp
@@ -1,74 +0,0 @@
-
-//-------------------------------------------------------------
-// (Minimum) Set Cover
-//   NP-complete
-//
-// Assumption: "mask" must be downward closerd,
-//   i.e., X subset Y & Y in mask => X in mask
-//
-// Verified by
-//   - Google Code Jam 09 Round2 C
-//     (The "Only" part never tested, in the test cases of
-//      the problem, that optimization had never used.
-//      Be careful to use that; to be on the safer side,
-//      consider using the <true,false> mode.)
-//-------------------------------------------------------------
-
-int minCover_exhaustive( int goal, const set<int>& mask )
-{
-	typedef int STATE;
-
-	vector<STATE> Q( 1, 0 );
-	set<STATE> visited;
-	for(int step=0; step<=mask.size(); ++step)
-	{
-		vector<STATE> Q2;
-		for(int i=0; i<Q.size(); ++i)
-		{
-			STATE s = Q[i];
-			if( s == goal )
-				return step;
-			for(set<int>::const_iterator it=mask.begin(); it!=mask.end(); ++it)
-				if( !visited.count( s|*it ) )
-				{
-					visited.insert( s|*it );
-					Q2.push_back( s|*it );
-				}
-		}
-		Q.swap(Q2);
-	}
-	return -1;
-}
-
-template<bool Subs, bool Only>
-int minCover( int goal, set<int> mask )
-{
-	if( Subs )
-		for(set<int>::iterator it=mask.begin(); it!=mask.end(); )
-		{
-			bool needed = true;
-			for(int m=1; m<=goal; m<<=1)
-				if( (goal&m) && !(*it&m) && mask.count(*it|m))
-					{needed = false; break;}
-
-			if( needed )
-				++it;
-			else
-				mask.erase(it++);
-		}
-
-	if( Only )
-		for(set<int>::iterator it=mask.begin(); it!=mask.end(); )
-		{
-			int onlyByMe = goal & *it;
-			for(set<int>::iterator jt=mask.begin(); jt!=mask.end(); ++jt)
-				onlyByMe &= ~*jt;
-	
-			if( !onlyByMe )
-				++it;
-			else
-				mask.erase(it++), goal &= ~onlyByMe;
-		}
-
-	return minCover_exhaustive(goal, mask);
-}

DELETED _lib/typical/ternery.cpp
Index: _lib/typical/ternery.cpp
==================================================================
--- _lib/typical/ternery.cpp
+++ _lib/typical/ternery.cpp
@@ -1,18 +0,0 @@
-// ternery search (for shita-ni-totsu)
-{
-	double x = ~min~;
-	double w = ~max~;
-	for(int i=0; i<~iteration~; ++i) // or, while( w-x > ~eps~ ) be careful for precision!
-	{
-		double y = 2*x/3 + w/3;
-		double z = x/3 + 2*w/3;
-		double fx = f(x);
-		double fy = f(y);
-		double fz = f(z);
-		double fw = f(w);
-		if( fx < fy )      w = y;
-		else if( fz > fw ) x = z;
-		else if( fy < fz ) w = z;
-		else               x = y;
-	}
-}

DELETED _lib/typical/unionfind.cpp
Index: _lib/typical/unionfind.cpp
==================================================================
--- _lib/typical/unionfind.cpp
+++ _lib/typical/unionfind.cpp
@@ -1,33 +0,0 @@
-//-------------------------------------------------------------
-// Union-Find
-//   Balancing + Path-Compression : O(invAckermann)
-//
-// Verified by
-//   - SRM Member Pilot 2 Div1 LV2
-//-------------------------------------------------------------
-
-struct UnionFind
-{
-	vector<int> uf, sz;
-	int nc;
-
-	UnionFind(int N): uf(N), sz(N,1), nc(N)
-		{ for(int i=0; i<N; ++i) uf[i] = i; }
-	int size()
-		{ return nc; }
-	int Find(int a)
-		{ return uf[a]==a ? a : uf[a]=Find(uf[a]); }
-	bool Union(int a, int b)
-		{
-			a = Find(a);
-			b = Find(b);
-			if( a != b )
-			{
-				if( sz[a] >= sz[b] ) swap(a, b);
-				uf[a] = b;
-				sz[b] += sz[a];
-				--nc;
-			}
-			return (a!=b);
-		}
-};

ADDED   lib/dataStructure/fenwickTree.cpp
Index: lib/dataStructure/fenwickTree.cpp
==================================================================
--- lib/dataStructure/fenwickTree.cpp
+++ lib/dataStructure/fenwickTree.cpp
@@ -0,0 +1,30 @@
+//-------------------------------------------------------------
+// Fenwick Tree
+//   O(log N) per each operation
+//
+// Verified by
+//   - SRM424 Div1 LV3
+//-------------------------------------------------------------
+
+template<typename T = LL>
+struct FenwickTree
+{
+	vector<T> x;
+	FenwickTree(size_t n, const T& v = T()) : x(n, v) {}
+
+	void incr( int k, const T& a ) { // z[k] += a;
+		for(; k < x.size(); k|=k+1)
+			x[k] += a;
+	}
+
+	T sum(int i, int j) { // z[i]+...+z[j] : inclusive
+		if(i)
+			return sum(0, j) - sum(0, i-1);
+		else {
+			T v = T();
+			for(; j>=0; j=(j&(j+1))-1)
+				v += x[j];
+			return v;
+		}
+	}
+};

ADDED   lib/dataStructure/rmq.cpp
Index: lib/dataStructure/rmq.cpp
==================================================================
--- lib/dataStructure/rmq.cpp
+++ lib/dataStructure/rmq.cpp
@@ -0,0 +1,53 @@
+//-------------------------------------------------------------
+// Range Minimum Query
+//   O(N log N) construction
+//   O(1) per each query
+//   returns the LEFTMOST/ALL minimum index (I hope)
+//
+// Verified by
+//   - SRM337 Div1 LV2
+//   - SRM431 Div1 LV3
+//-------------------------------------------------------------
+
+template<typename T>
+struct RMQ
+{
+	vector< vector<int> > rm;
+	vector<T> d;
+
+	RMQ( const vector<T>& d ) : d(d) {
+		int n = d.size();
+
+		// rm[k][x] = i s.t. d[i] is the minimum in [x, x+2^k)
+		rm.push_back( vector<int>(n) );
+		for(int x=0; x<n; ++x)
+			rm[0][x] = x;
+		for(int k=1; (1<<k)<=n; ++k) {
+			rm.push_back( rm[k-1] );
+			for(int x=0; x+(1<<k-1)<n; ++x)
+				if( d[rm[k][x]] > d[rm[k-1][x + (1<<k-1)]] )
+					rm[k][x] = rm[k-1][x + (1<<k-1)];
+		}
+	}
+
+	int operator()(int L, int R) const { // inclusive
+		int k=0;
+		for(; L+(1<<k) < R-(1<<k)+1; ++k) {}
+		LL i = rm[k][L];
+		LL j = rm[k][R-(1<<k)+1];
+		return (d[i]<=d[j] ? i : j);
+	}
+
+	vector<int> all(int L, int R) const {
+		vector<int> ans;
+		int minValue = d[(*this)(L, R)];
+		while( L <= R ) {
+			int C = (*this)(L, R);
+			if( minValue < d[C] )
+				break;
+			ans.push_back(C);
+			L = C+1;
+		}
+		return ans;
+	}
+};

ADDED   lib/doc/hitori-sharp.pdf
Index: lib/doc/hitori-sharp.pdf
==================================================================
--- lib/doc/hitori-sharp.pdf
+++ lib/doc/hitori-sharp.pdf
cannot compute difference between binary files

ADDED   lib/doc/lb_ub.txt
Index: lib/doc/lb_ub.txt
==================================================================
--- lib/doc/lb_ub.txt
+++ lib/doc/lb_ub.txt
@@ -0,0 +1,21 @@
+lb = lower_bound(k)
+  ... k-1] lb [k ...
+
+
+ub = upper_bound(k)
+  .., k] ub [k+1 ...
+
+
+なので任意の開閉区間をイテレータの [) 区間に直すには
+
+  [A, B]
+    = [ lower_bound(A), upper_bound(B) )
+
+  [A, B)
+    = [ lower_bound(A), lower_bound(B) )
+
+  (A, B]
+    = [ upper_bound(A), upper_bound(B) )
+
+  (A, B)
+    = [ upper_bound(A), lower_bound(B) )

ADDED   lib/doc/nim.txt
Index: lib/doc/nim.txt
==================================================================
--- lib/doc/nim.txt
+++ lib/doc/nim.txt
@@ -0,0 +1,38 @@
+// SRM338 Div1 LV3
+
+定理のステートメント
+
+    対象となるゲーム
+      - 二人ゲーム
+      - 二人とも動かせる手の集合は同じ
+      - 動かせなくなった方が負け
+      - 無限ループしない
+
+    nim
+      - サイズ s1, s2, .., sk の山がある
+      - どれか一つ山を選んで 1 ～ 山サイズ の任意個取り除ける
+      - 最後の山を取った方が勝ち（山がなくなって打つ手が無くなった方が負け）
+
+    *n をサイズ n の山が1個だけの nim とする
+      - *0 は負け
+      - *n は勝ち n>=1
+
+    ゲーム
+      - 状態 G から打てる手で遷移する先を G1, ..., Gk とする
+        G = {G1, ..., Gk} と書く
+
+    等価
+      - 二つのゲームG, Fが等価なことを G ≡ F と書く
+        等価というのは勝ち負け一致よりもっと強い。二つのゲームの和を取ると必敗になること
+        *n と *n は同サイズの山2つのnimは必敗なので等価。同サイズでなければ必勝なので等価でない
+
+    定理
+      - G1≡*n1, ..., Gk≡*nk ならば、G={G1, ..., Gk} ≡ *mex(n1,...,nk)
+        where mex(S) = min(nat \setminus S)
+        帰納的に、全てのゲームはなんらかの *n と等価
+
+
+おまけ
+  @G を G のnim値とする。つまりG≡*@G
+  G = X + Y ならば @G = @X xor @Y
+  山が複数のnimの勝敗は全部xorとった1つ山のnimに等しいのと同じ原理。

ADDED   lib/gcj/!gcj-smpl.aml
Index: lib/gcj/!gcj-smpl.aml
==================================================================
--- lib/gcj/!gcj-smpl.aml
+++ lib/gcj/!gcj-smpl.aml
@@ -0,0 +1,40 @@
+(****************************************************************************
+ * Written in
+ *    Alice 1.4
+ *    http://www.ps.uni-saarland.de/alice/
+ ****************************************************************************)
+
+fun solve (n : int, k : int) : bool =
+	let
+		open Word
+		infix 5 << andb
+		val mask = (0wx1 << fromInt n) - 0wx1
+	in
+		(fromInt k) andb mask = mask
+	end
+
+fun main () =
+	let
+		fun readLine () =
+			case TextIO.inputLine TextIO.stdIn of
+			| NONE    => []
+			| SOME(s) => map Int.fromString (String.tokens (fn x => x = #" ") s)
+
+		fun parseOne c =
+			case readLine() of
+			| [SOME n, SOME k] => (c, (n, k))
+			| _                => assert false
+
+		fun spawn solveOne (c, inp) =
+			(c, solve inp)
+
+		fun printOne (c, ans) =
+			print ("Case #" ^ Int.toString (c+1) ^ ": " ^ (if ans then "ON" else "OFF") ^ "\n")
+	in
+		case readLine () of
+		| [SOME t] => List.app printOne (List.map solveOne (List.tabulate (t, parseOne)))
+		| _        => assert false
+	end
+
+do main ()
+do OS.Process.exit 0

ADDED   lib/gcj/!gcj-smpl.as
Index: lib/gcj/!gcj-smpl.as
==================================================================
--- lib/gcj/!gcj-smpl.as
+++ lib/gcj/!gcj-smpl.as
@@ -0,0 +1,33 @@
+/////////////////////////////////////////////////////////////////////////////
+// Written in
+//    ActionScript (using mtasc 1.14)
+//    http://www.mtasc.org/
+/////////////////////////////////////////////////////////////////////////////
+
+class A
+{
+	static function main(mc)
+	{
+		var input_string = "<paste the input data here>";
+		var input = input_string.split("\n");
+		var T = Number(input[0]);
+
+		var output_string = ""
+		for(var C=1; C<=T; ++C)
+		{
+			var theCase = input[C].split(" ");
+			var N = Number(theCase[0]);
+			var K = Number(theCase[1]);
+			output_string += "Case #" + C + ": " + (solve(N,K)?"ON":"OFF") + "\n"
+		}
+
+		_root.createTextField("tf",0,0,0,800,600);
+		_root.tf.text = output_string;
+	}
+
+	static function solve(N, K)
+	{
+		var mask = (1<<N) - 1;
+		return (K & mask) == mask;
+	}
+}

ADDED   lib/gcj/!gcj-smpl.bl
Index: lib/gcj/!gcj-smpl.bl
==================================================================
--- lib/gcj/!gcj-smpl.bl
+++ lib/gcj/!gcj-smpl.bl
@@ -0,0 +1,133 @@
+#############################################################################
+# Written in
+#    Blue 1.7.5
+#    http://www.lechak.info/blue/
+#############################################################################
+
+global FROM_STRING = func {
+	arg s;
+
+	a = [];
+	i = 0;
+	loop {
+		a = a.append( s.substr(i,i+1).num() );
+		(i=i+1) >= s.length() ? return;
+	};
+	return a;
+};
+
+global TO_STRING = func {
+	arg x;
+	lexical s = "";
+	x.map( func{ s = s + this; } );
+	return s;
+};
+
+global LESS = func {
+	arg x;
+	arg y;
+
+	xl = x.length();
+	yl = y.length();
+	(xl < yl) ? return 1;
+	(xl > yl) ? return 0;
+	i = 0;
+	return loop {
+		(i == xl)     ? return 0;
+		(x[i] < y[i]) ? return 1;
+		(x[i] > y[i]) ? return 0;
+		i = i+1;
+	};
+};
+
+global SUB = func {
+	arg x;
+	arg y;
+
+	x.length() > y.length() ? (y = [].resize(x.length()-y.length(),0).merge(y));
+	z = [].resize(x.length(), 0);
+	i = x.length()-1;
+	carry = 0;
+	loop {
+		z[i]  = x[i] - y[i] - carry;
+		carry = z[i] < 0 ? (z[i]=z[i]+10; 1) : 0;
+		(i=i-1) < 0 ? return;
+	};
+	head = 0;
+	loop {
+		head == z.length()-1 ? return;
+		z[head] != 0 ? return;
+		head = head + 1;
+	};
+	return z.slice(head, z.length());
+};
+
+global DBL = func {
+	arg x;
+
+	z = [].resize(x.length(), 0);
+	i = x.length()-1;
+	carry = 0;
+	loop {
+		z[i]  = x[i] + x[i] + carry;
+		carry = z[i] > 9 ? (z[i]=z[i]-10; 1) : 0;
+		((i=i-1) < 0) ? return;
+	};
+	return (carry==1 ? [1].merge(z) : z);
+};
+
+global DIF = func {
+	arg x;
+	arg y;
+	return LESS(x, y) ? SUB(y, x) : SUB(x, y);
+};
+
+global MOD = func {
+	arg x;
+	arg y;
+
+	LESS(x,y) ? return x;
+	z = MOD(x, DBL(y));
+	LESS(z,y) ? return z;
+	return SUB(z,y);
+};
+
+global ISZERO = func {
+	arg x;
+	return x[0] == 0;
+};
+
+global GCD = func {
+	arg x;
+	arg y;
+	return ISZERO(x) ? y : GCD(MOD(y,x),x);
+};
+
+#############################################################################
+
+solve = func {
+	arg N;
+	arg t;
+	t = t.map( func{return FROM_STRING(this);} );
+
+	lexical t0 = t[0];
+	lexical g  = DIF(t0, t[1]);
+	t.slice(2, t.length()).map(func{
+		g = GCD(g, DIF(t0, this));
+	});
+	r = MOD(t[0], g);
+	return TO_STRING( ISZERO(r) ? r : DIF(g,r) );
+};
+
+input = sys.library("streams").stdio().read(1000000000).split("\n").map(func{return this.split(" ");});
+C = input[0][0].num();
+
+CaseID = (args.length() >= 2 ? args[1].num() : 1);
+loop {
+	"Case #".print();
+	CaseID.print();
+	": ".print();
+	solve( input[CaseID][0].num(), input[CaseID].slice(1,input[CaseID].length()) ).print();
+	"\n".print();
+	((CaseID = CaseID+1) > C) ? return;
+};

ADDED   lib/gcj/!gcj-smpl.cl
Index: lib/gcj/!gcj-smpl.cl
==================================================================
--- lib/gcj/!gcj-smpl.cl
+++ lib/gcj/!gcj-smpl.cl
@@ -0,0 +1,150 @@
+/////////////////////////////////////////////////////////////////////////////
+// Writte in
+//    Claire v3.3.37
+//    http://www.claire-language.com/
+/////////////////////////////////////////////////////////////////////////////
+
+D :: 1000000000
+
+long <: ephemeral_object(hi:integer, lo:integer)
+
+[long!(x: integer) : long ->
+	long(hi = x / D, lo = x mod D)
+]
+
+[+(x:long, y:long) : long ->
+	let lolo := x.lo - D + y.lo in
+	(
+		if ( lolo < 0 )
+			long(hi = x.hi + y.hi,     lo = lolo + D)
+		else
+			long(hi = x.hi + y.hi + 1, lo = lolo)
+	)
+]
+
+[-(x:long, y:long) : long ->
+	let lolo := x.lo - y.lo in
+	(
+		if ( lolo < 0 )
+			long(hi = x.hi - y.hi - 1, lo = lolo + D)
+		else
+			long(hi = x.hi - y.hi,     lo = lolo)
+	)
+]
+
+[*(x:long, y:integer) : long ->
+	let r := long!(0) in
+	let c := x in
+		(while (y > 0) (
+			(if (y mod 2 = 1) r :+ c),
+			c := c + c,
+			y :/ 2
+		),
+		r)
+]
+
+[+(x:long, y:integer) : long ->
+	x + long!(y)
+]
+
+[<(x:long, y:long) : boolean ->
+	if (x.hi < y.hi)
+		true
+	else if (x.hi > y.hi)
+		false
+	else
+		x.lo < y.lo
+]
+[>=(x:long, y:long) : boolean -> not(x < y)]
+[<=(x:long, y:long) : boolean -> not(y < x)]
+[>(x:long, y:long) : boolean -> (y < x)]
+
+[>=(x:long, y:integer) : boolean -> not(x < long!(y))]
+[<=(x:long, y:integer) : boolean -> not(long!(y) < x)]
+[<(x:long, y:integer) : boolean -> (x < long!(y))]
+[>(x:long, y:integer) : boolean -> (long!(y) < x)]
+
+[>=(x:integer, y:long) : boolean -> not(long!(x) < y)]
+[<=(x:integer, y:long) : boolean -> not(y < long!(x))]
+[<(x:integer, y:long) : boolean -> (long!(x) < y)]
+[>(x:integer, y:long) : boolean -> (y < long!(x))]
+
+
+[string!(x: long) : string ->
+	if (x.hi > 0)
+		let los := string!(x.lo) in
+			string!(x.hi) /+ make_string(9 - length(los), '0') /+ los
+	else
+		string!(x.lo)
+]
+
+/////////////////////////////////////////////////////////////////////////////
+
+[readLine() : list<integer> ->
+	list<integer>{integer!(s) | s in explode(freadline(stdin), " ")}
+]
+
+[solve(R:integer, k:integer, N:integer, g:list<integer>) ->
+	let dest := make_list(N, 0) in
+	let earn := make_list(N, long!(0)) in
+	(
+		for s in (0 .. N - 1) (
+			let ride := long!(0) in
+			let i    := 0        in
+			(
+				while (i < N)
+				(
+					let ride2 := ride + g[((s + i) mod N) + 1] in
+					(if (ride2 > k)
+						break()
+					else
+						(ride := ride2, i :+ 1))
+				),
+				earn[s + 1] := ride,
+				dest[s + 1] := ((s + i) mod N) + 1
+			)
+		),
+		let firstVisitTime := make_list(N,       -1) in
+		let firstVisitEarn := make_list(N, long!(0)) in
+		let q         := 1        in
+		let totalEarn := long!(0) in
+		let i         := 0        in
+		(
+			(while (i < R) (
+				if (firstVisitTime[q] = -1)
+				(
+					firstVisitTime[q] := i,
+					firstVisitEarn[q] := totalEarn,
+					totalEarn         :+ earn[q],
+					q                 := dest[q],
+					i                 :+ 1
+				)
+				else
+				(
+					let loopSize := i - firstVisitTime[q] in
+					let loopEarn := totalEarn - firstVisitEarn[q] in
+					let rest     := R - i in
+					(
+						totalEarn :+ loopEarn * (rest / loopSize),
+						// clear
+						firstVisitTime := make_list(N, -1),
+						i              := R - (rest mod loopSize)
+					)
+				)
+			)),
+			princ(string!(totalEarn))
+		)
+	)
+]
+
+[main() ->
+	let T := readLine()[1] in
+		for C in (1 .. T)
+		(
+			printf("Case #~A: ", C),
+			let RkN := readLine() in solve(RkN[1], RkN[2], RkN[3], readLine()),
+			princ("\n")
+		)
+]
+
+(main())

ADDED   lib/gcj/!gcj-smpl.cs
Index: lib/gcj/!gcj-smpl.cs
==================================================================
--- lib/gcj/!gcj-smpl.cs
+++ lib/gcj/!gcj-smpl.cs
@@ -0,0 +1,60 @@
+using System;
+using System.Linq;
+
+class C
+{
+	static void Main()
+	{
+		long[] Ts = readLongArray();
+		long   T  = Ts[0];
+
+		for(long C=1; C<=T; ++C)
+		{
+			long[] RkN = readLongArray();
+			long   R = RkN[0];
+			long   k = RkN[1];
+			long   N = RkN[2];
+			long[] g = readLongArray();
+
+			Console.WriteLine("Case #{0}: {1}", C, solveSmall(R,k,N,g));
+		}
+	}
+
+	static long[] readLongArray()
+	{
+		return (from s in Console.ReadLine().Split(' ') select long.Parse(s)).ToArray();
+	}
+
+	static long solveSmall(long R, long k, long N, long[] g)
+	{
+		long totalEarn = 0;
+
+		for(int i=0; i<R; ++i)
+		{
+			long ride = 0;
+			for(int q=0; q<g.Length; ++q)
+				if( ride + g[q] > k )
+				{
+					g = rotate(g, q);
+					break;
+				}
+				else
+				{
+					ride += g[q];
+				}
+			totalEarn += ride;
+		}
+
+		return totalEarn;
+	}
+
+	static long[] rotate(long[] g, int s)
+	{
+		long[] g2 = new long[g.Length];
+		for(int i=s; i<g.Length; ++i)
+			g2[i-s] = g[i];
+		for(int i=0; i<s; ++i)
+			g2[i+g.Length-s] = g[i];
+		return g2;
+	}
+}

ADDED   lib/gcj/!gcj-smpl.vbs
Index: lib/gcj/!gcj-smpl.vbs
==================================================================
--- lib/gcj/!gcj-smpl.vbs
+++ lib/gcj/!gcj-smpl.vbs
@@ -0,0 +1,39 @@
+'
+' Written in VBScript
+'
+
+Function ReadLine()
+	Dim s
+	s = Split(WScript.StdIn.ReadLine, " ")
+	For i = LBound(s) To UBound(s)
+		s(i) = CLng(s(i))
+	Next
+	ReadLine = s
+End Function
+
+Function Gcd(a, b)
+	If a = 0 Then
+		Gcd = b
+	Else
+		Gcd = Gcd(b mod a, a)
+	End If
+End Function
+
+Function Solve(T)
+	Dim g, r
+	g = Abs(T(1) - T(2))
+	For i = 2 To UBound(T)
+		g = Gcd(g, Abs(T(1) - T(i)))
+	Next
+	r = T(1) mod g
+	If r = 0 Then
+		Solve = 0
+	Else
+		Solve = g - r
+	End If
+End Function
+
+C = ReadLine()(0)
+For CaseID = 1 To C
+	WScript.StdOut.WriteLine "Case #" & CaseID & ": " & Solve(ReadLine)
+Next

ADDED   lib/gcj/!gcj-tmpl.cpp
Index: lib/gcj/!gcj-tmpl.cpp
==================================================================
--- lib/gcj/!gcj-tmpl.cpp
+++ lib/gcj/!gcj-tmpl.cpp
@@ -0,0 +1,100 @@
+//-----------------------------------------------------------------------------
+// >>Code Template<< (for Visual C++)
+
+#include <iostream>
+#include <sstream>
+#include <iomanip>
+#include <string>
+#include <vector>
+#include <set>
+#include <map>
+#include <algorithm>
+#include <numeric>
+#include <iterator>
+#include <complex>
+#include <functional>
+#include <queue>
+#include <stack>
+#include <cmath>
+#include <cassert>
+#include <cstring>
+#define  cout os
+using namespace std;
+typedef long long LL;
+typedef complex<double> CMP;
+void END_OF_INPUT_FOR_THIS_TEST_CASE(); // stub for multi-threading
+
+//-----------------------------------------------------------------------------
+// >>Main<<
+
+void case_main( ostream& os )
+{
+	END_OF_INPUT_FOR_THIS_TEST_CASE();
+}
+
+//-----------------------------------------------------------------------------
+// >>Code Template<< (Multi-Thread Solver)
+
+#if 1
+#undef cout
+#include <windows.h>
+#include <process.h>
+
+static const int THREAD_NUM = 2;
+volatile int     g_id;
+int              g_nCase;
+CRITICAL_SECTION g_cs;
+vector<string>   g_output;
+
+unsigned __stdcall thread_main( void* t_id ) {
+	for(;;) {
+		EnterCriticalSection(&g_cs);
+		const int id = ++g_id;
+		if(id>g_nCase) { LeaveCriticalSection(&g_cs); break; }
+		cerr << setw(4) << id << " @ " << (int)t_id << " start" << endl;
+
+		ostringstream ss;
+		ss << "Case #" << id << ": ";
+		case_main( ss );
+
+		EnterCriticalSection(&g_cs);
+		if(g_output.size()<id) g_output.resize(id);
+		g_output[id-1] = ss.str();
+		cerr << setw(4) << id << " @ " << (int)t_id << " end" << endl;
+		LeaveCriticalSection(&g_cs);
+	}
+	return 0;
+}
+
+void END_OF_INPUT_FOR_THIS_TEST_CASE() { LeaveCriticalSection(&g_cs); }
+
+int main() {
+	cin >> g_nCase;
+	string dummy; getline(cin, dummy);
+
+	InitializeCriticalSection(&g_cs);
+	vector<HANDLE> thread;
+	for(int i=0; i<THREAD_NUM; ++i)
+		thread.push_back( (HANDLE)_beginthreadex(0, 0, &thread_main, (void*)i, 0, 0) );
+	WaitForMultipleObjects( thread.size(), &thread[0], TRUE, INFINITE );
+	DeleteCriticalSection(&g_cs);
+
+	copy( g_output.begin(), g_output.end(), ostream_iterator<string>(cout) );
+}
+
+//-----------------------------------------------------------------------------
+// >>Code Template<< (Single-Thread Solver)
+
+#else
+#undef cout
+void END_OF_INPUT_FOR_THIS_TEST_CASE() {}
+int main() {
+	int nCase; cin >> nCase;
+	string dummy; getline(cin, dummy);
+
+	for(int id=1; id<=nCase; ++id) {
+		cout << "Case #" << id << ": ";
+		case_main( cout );
+	}
+}
+#endif

ADDED   lib/gcj/!gcj-tmpl.java
Index: lib/gcj/!gcj-tmpl.java
==================================================================
--- lib/gcj/!gcj-tmpl.java
+++ lib/gcj/!gcj-tmpl.java
@@ -0,0 +1,40 @@
+import java.io.*;
+import java.util.*;
+
+public class <CLASSNAME>
+{
+	public static void main(String[] arg)
+	{
+		Scanner sc = new Scanner(System.in);
+		int T = sc.nextInt();
+		for(int C=1; C<=T; ++C)
+		{
+			System.out.printf("Case #%d: ", C);
+			(new <CLASSNAME>(sc)).caseMain();
+		}
+	}
+
+	Scanner sc;
+	<CLASSNAME>( Scanner sc ) { this.sc = sc; }
+
+	void caseMain()
+	{
+		int n = sc.nextInt();
+		long[] x = new long[n];
+		for(int i=0; i<n; ++i) x[i] = sc.nextInt();
+		long[] y = new long[n];
+		for(int i=0; i<n; ++i) y[i] = sc.nextInt();
+
+		System.out.println(solve(x, y, n));
+	}
+
+	long solve(long[] x, long[] y, int n)
+	{
+		Arrays.sort(x);
+		Arrays.sort(y);
+		long s = 0;
+		for(int i=0; i<n; ++i)
+			s += x[i] * y[n-1-i];
+		return s;
+	}
+}

ADDED   lib/gcj/!gcj-tmpl.lua
Index: lib/gcj/!gcj-tmpl.lua
==================================================================
--- lib/gcj/!gcj-tmpl.lua
+++ lib/gcj/!gcj-tmpl.lua
@@ -0,0 +1,19 @@
+CN = io.stdin:read("*n")
+for C=1, CN, 1 do
+	io.stdout:write( string.format("Case #%d: ", C) )
+
+	N = io.stdin:read("*n")
+	M = io.stdin:read("*n")
+
+	fav = {}
+	for i=1, M, 1 do
+		T = io.stdin:read("*n")
+		fav[i] = {}
+		for j=1, T, 1 do
+			X, Y = io.stdin:read("*n", "*n")
+			fav[i][j] = {X, Y}
+		end
+	end
+
+	io.stdout:write( string.format("%d %d\n", N, M) )
+end

ADDED   lib/gcj/!gcj-tmpl.py
Index: lib/gcj/!gcj-tmpl.py
==================================================================
--- lib/gcj/!gcj-tmpl.py
+++ lib/gcj/!gcj-tmpl.py
@@ -0,0 +1,11 @@
+import sys
+
+CN = int(sys.stdin.next())
+for C in range(1, CN+1):
+	print "Case #%d:"%(C),
+	N = int(sys.stdin.next())
+	M = int(sys.stdin.next())
+	for i in range(0, M):
+		a = map(int, sys.stdin.next().split(" "))
+		print a
+	print N, M

ADDED   lib/geo/area.cpp
Index: lib/geo/area.cpp
==================================================================
--- lib/geo/area.cpp
+++ lib/geo/area.cpp
@@ -0,0 +1,23 @@
+//-------------------------------------------------------------
+// Area of a polygon
+//   O(N)
+//
+// Verified by
+//   - SRM 337 Div1 LV3
+//   - SRM 486 Div1 LV3
+//-------------------------------------------------------------
+
+double outer_prod( pt a, pt b )
+{
+	return (a.real()*b.imag() - b.real()*a.imag())/2;
+}
+
+double area( const vector<pt>& q )
+{
+	double a = 0.0;
+
+	pt o = q[0];
+	for(int i=1; i+1<q.size(); ++i)
+		a += outer_prod(q[i]-o, q[i+1]-o);
+	return abs(a) / 2;
+}

ADDED   lib/geo/ccw.cpp
Index: lib/geo/ccw.cpp
==================================================================
--- lib/geo/ccw.cpp
+++ lib/geo/ccw.cpp
@@ -0,0 +1,18 @@
+//-------------------------------------------------------------
+// ccw
+//
+// Verified by
+//   - SRM 492 Div1 LV1
+//-------------------------------------------------------------
+
+double outer_prod(const CMP& a, const CMP& b) { return imag(conj(a)*b); }
+double inner_prod(const CMP& a, const CMP& b) { return real(conj(a)*b); }
+
+int ccw(const CMP& a, CMP b, CMP c) {
+	b -= a; c -= a;
+	if( outer_prod(b,c) > 0 ) return +1; // counter clockwise
+	if( outer_prod(b,c) < 0 ) return -1; // clockwise
+	if( inner_prod(b,c) < 0 ) return +2; // c--[a--b] on line
+	if( norm(b) < norm(c) )   return -2; // [a--b]--c on line
+	return 0; // [a--c--b] on line
+}

ADDED   lib/geo/circle_3pt.cpp
Index: lib/geo/circle_3pt.cpp
==================================================================
--- lib/geo/circle_3pt.cpp
+++ lib/geo/circle_3pt.cpp
@@ -0,0 +1,55 @@
+
+//-------------------------------------------------------------
+// The circle passing three points
+//
+// Verified by
+//   - AOJ 0010
+//-------------------------------------------------------------
+
+vector<double> solve_linear_eq( int n, vector< vector<double> > M, const vector<double>& V )
+{
+	vector<double> A(V);
+	for(int i=0; i<n; ++i)
+	{
+		// pivot
+		if( M[i][i] == 0 )
+			for(int j=i+1; j<n; ++j)
+				if( M[j][i] != 0 )
+					{swap(M[i], M[j]); swap(A[i], A[j]); break;}
+		if( M[i][i] == 0 )
+			throw "no anser";
+
+		// M[i][i] <-- 1
+		double p = M[i][i];
+		for(int j=i; j<n; ++j)
+			M[i][j] /= p;
+		A[i] /= p;
+
+		// M[*][i] <-- 0
+		for(int j=0; j<n; ++j) if(j!=i)
+		{
+			double r = M[j][i];
+			for(int k=i; k<n; ++k)
+				M[j][k] -= M[i][k] * r;
+			A[j] -= A[i] * r;
+		}
+	}
+	return A;
+}
+
+void circle_3pt( CMP p1, CMP p2, CMP p3, CMP* c, double* r )
+{
+	p2 -= p1;
+	p3 -= p1;
+	// c == p2/2 + A p2 i == p3/2 + B p3 i
+
+	vector< vector<double> > M(2, vector<double>(2));
+	vector<double> V(2);
+	M[0][0] = -p2.imag(); M[0][1] = +p3.imag(); V[0] = (p3.real()-p2.real())/2.0;
+	M[1][0] = +p2.real(); M[1][1] = -p3.real(); V[1] = (p3.imag()-p2.imag())/2.0;
+	V = solve_linear_eq(2, M, V);
+
+	double A=V[0], B=V[1];
+	*c = p1 + p2/2.0 + A * p2 * CMP(0,1);
+	*r = abs(*c - p1);
+}

ADDED   lib/geo/convexHull.cpp
Index: lib/geo/convexHull.cpp
==================================================================
--- lib/geo/convexHull.cpp
+++ lib/geo/convexHull.cpp
@@ -0,0 +1,48 @@
+
+//-------------------------------------------------------------
+// Convex Hull -- Andrew's Monotone Chain
+//   O(N log N)
+//
+// Verified by
+//   - SRM 336 Div1 LV3
+//   - TCO09 Round2 LV2
+//   - SRM 486 Div1 LV3
+//-------------------------------------------------------------
+
+double outer_prod(const CMP& a, const CMP& b) { return imag(conj(a)*b); }
+double inner_prod(const CMP& a, const CMP& b) { return real(conj(a)*b); }
+
+int ccw(const CMP& a, CMP b, CMP c) {
+	b -= a; c -= a;
+	if( outer_prod(b,c) > 0 ) return +1; // counter clockwise
+	if( outer_prod(b,c) < 0 ) return -1; // clockwise
+	if( inner_prod(b,c) < 0 ) return +2; // c--a--b on line
+	if( norm(b) < norm(c) )   return -2; // a--b--c on line
+	return 0;
+}
+
+bool byX( const CMP& a, const CMP& b ) {
+	if( a.real() != b.real() )
+		return a.real() < b.real();
+	return a.imag() < b.imag();
+}
+
+vector<CMP> convex_hull( vector<CMP> p )
+{
+	#define IS_RIGHT <0   // skip on-line verts
+	//#define IS_RIGHT ==-1 // take all
+
+	sort(p.begin(), p.end(), &byX);
+
+	vector<CMP> ans;
+	for(int i=0; i<p.size(); ans.push_back(p[i++])) // left-to-right
+		while( ans.size()>=2 && ccw(ans[ans.size()-2], ans[ans.size()-1], p[i]) IS_RIGHT )
+			ans.pop_back();
+	if( ans.size() == p.size() )
+		return ans;
+	for(int i=p.size()-2; i>=0; ans.push_back(p[i--])) // right-to-left
+		while( ans.size()>=2 && ccw(ans[ans.size()-2], ans[ans.size()-1], p[i]) IS_RIGHT )
+			ans.pop_back();
+	ans.pop_back();
+	return ans;
+}

ADDED   lib/geo/convexHull_old.cpp
Index: lib/geo/convexHull_old.cpp
==================================================================
--- lib/geo/convexHull_old.cpp
+++ lib/geo/convexHull_old.cpp
@@ -0,0 +1,49 @@
+
+//-------------------------------------------------------------
+// Convex Hull -- Gift Wrapping
+//   O(N log N)
+//
+// ToDo : Rewrite by CCW...
+//
+// Verified by
+//   - SRM 336 Div1 LV3
+//-------------------------------------------------------------
+
+bool byY( CMP a, CMP b )
+{
+	if( a.imag() == b.imag() )
+		return a.real() < b.real();
+	return a.imag() < b.imag();
+}
+
+bool byArg( CMP a, CMP b )
+{
+	return arg(a) < arg(b);
+}
+
+bool isRight( CMP a, CMP b, CMP c )
+{
+	return arg((c-b) / (b-a)) < 0;
+}
+
+vector<CMP> convex_hull( vector<CMP> p )
+{
+	vector<CMP>::iterator it = min_element( p.begin(), p.end(), byY );
+	CMP o = *it;
+	p.erase(it);
+	for(int i=0; i<p.size(); ++i)
+		p[i] -= o;
+	sort( p.begin(), p.end(), byArg );
+
+	vector<CMP> q;
+	q.push_back( CMP(0,0) );
+	q.push_back( p[0] );
+	for(int i=1; i<p.size(); ++i) {
+		while( isRight(q[q.size()-2], q[q.size()-1], p[i]) )
+			q.pop_back();
+		q.push_back( p[i] );
+	}
+	for(int i=0; i<q.size(); ++i)
+		q[i] += o;
+	return q;
+}

ADDED   lib/geo/pt_in_poly.cpp
Index: lib/geo/pt_in_poly.cpp
==================================================================
--- lib/geo/pt_in_poly.cpp
+++ lib/geo/pt_in_poly.cpp
@@ -0,0 +1,28 @@
+
+//-------------------------------------------------------------
+// The circle passing three points
+//
+// Verified by
+//   - AOJ 0012 (only triangles)
+//-------------------------------------------------------------
+
+double outer_prod( CMP a, CMP b )
+{
+	return (a.real()*b.imag() - b.real()*a.imag())/2;
+}
+
+bool point_in_polygon( vector<CMP>& ps, CMP p )
+{
+	bool in = false;
+	for(int i=0; i<ps.size(); ++i) {
+		CMP a = ps[i] - p;
+		CMP b = ps[(i+1)%ps.size()] - p;
+		if(a.imag() > b.imag()) swap(a,b);
+		if(a.imag()<=0 && 0<b.imag()) {
+			if( outer_prod(a,b) < 0 )
+				in = !in;
+		}
+		//if( outer_prod(a,b)==0 && inner_prod(a,b)<=0 ) return ON;
+	}
+	return in;
+}

ADDED   lib/graph/__goldberg.cpp
Index: lib/graph/__goldberg.cpp
==================================================================
--- lib/graph/__goldberg.cpp
+++ lib/graph/__goldberg.cpp
@@ -0,0 +1,126 @@
+//----------------------------------------------------------------------------
+// Goldberg(-Tarjan?) �@ (�ő嗬)
+//
+// �e�ʕt���O���t G ��ł� Src ���� Dest �ւ̍ő嗬�ʂ����߂�
+//   - �O���t�͗אڍs��œn�����ƁBG �͔j�󂳂�܂��B
+//   - ���ʂ����łȂ�Flow���̂��~�����Ƃ��� G �̗e�ʂ������Ă镔��������
+//
+// �v�Z�ʂ�
+//   - O(V^3)
+//       # �ؖ��͂܂��������Ă��Ȃ�
+//
+// �A���S���Y���̊T��
+//   ��{�`�͈ȉ��̒ʂ�
+//      1. �Ƃ肠����Src����o�Ă�ӑS���ɗe�ʌ��E�܂ŗ���
+//      2. ���� > ���o �ȃm�[�h�����‚��ăo�����X�����
+//          - ��������Dest�܂œ��B����p�X������Ȃ�A���
+//            Dest�ɋ߂������֐i�ޕӂւ̗��o�ʂ𑝂₷
+//          - �Ȃ��Ȃ�A���Src�֋߂��������痈��ӂ����
+//            �����ʂ����炷
+//      3. 2.���A���o�����X�Ȓ��_���Ȃ��Ȃ�܂ŌJ��Ԃ�
+//   �u���� > ���o �ȃm�[�h�����‚��āv�̏����́A���o���o�����X��ς�������
+//   Queue�ɓ˂�����ł������Ƃ�FIFO�ł��B
+//   �uDest�ɋ߂������v�uSrc�ɋ߂������v����Ɍ����ɕ]������K�v�͂Ȃ��āA
+//   ������x�̋ߎ��ł��\���B���̃R�[�h�ł͂��̋ߎ��l�� d[] �Ŏ����Ă���B
+//   d[] �̕]����K�x�ɂ����Ƃ�����肷��heuristics�ő����Ȃ�炵���B
+//   �Ƃ��������Ȃ��Ǝ��ۏ�͎g�����ɂȂ�Ȃ��炵���B���Ƃł����Ƃ��B
+//----------------------------------------------------------------------------
+
+#include <vector>
+#include <queue>
+#include <limits>
+using namespace std;
+
+typedef int                        Vert;
+typedef int                        Capacity;
+typedef vector< vector<Capacity> > Graph;
+
+Capacity goldberg_tarjan( Graph& G, Vert Src, Vert Dest )
+{
+	vector<Capacity> e( G.size() ); // excess   : ������ - ���o��
+	vector<int>      d( G.size() ); // distance : Dest�ւ̋���||Src�ւ̋���-N�̉���
+	d[Src] = G.size();
+
+	queue<Vert> Q; // e[v]>0 �ȃm�[�h�����Ă����L���[
+
+	// �Ƃ肠����Src����S�͂ŗ���
+	for(int v=0; v!=G.size(); ++v)
+		if( G[Src][v] )
+		{
+			G[v][Src] += G[Src][v];
+			e[v]      += G[Src][v];
+			G[Src][v]  = 0;
+			Q.push(v);
+		}
+
+	// e[v]>0 �ȃm�[�h���Ȃ��Ȃ�܂ŌJ��Ԃ�
+	while( !Q.empty() )
+	{
+		Vert v = Q.front();
+
+		for(int u=0; u!=G.size() && e[v]; ++u)
+			if( G[v][u] && d[v]>d[u] ) // �K�؂ȕ����ɗ��� "PUSH"
+			{
+				Capacity f = min(e[v], G[v][u]);
+				G[v][u] -= f;
+				G[u][v] += f;
+				e[v] -= f;
+				e[u] += f;
+				if( e[u]==f && u!=Src && u!=Dest ) Q.push(u);
+			}
+
+		if( e[v] == 0 ) // �o�����X��ꂽ�̂ŃL���[���珜��
+			Q.pop();
+		else // �o�����X���ĂȂ��̂͂��������̂�d[v]�𒲐����Ă�蒼�� "RELABEL"
+		{
+			Capacity m = numeric_limits<Capacity>::max();
+			for(int u=0; u!=G.size(); ++u)
+				if( G[v][u] )
+					m = min(m, d[u]+1);
+			d[v] = m;
+		}
+	}
+
+	// Dest �ւ̗����ʂ�Ԃ�
+	return e[Dest];
+}
+
+//----------------------------------------------------------------------------
+// PKU 1459
+//----------------------------------------------------------------------------
+
+#include <iostream>
+#include <locale>
+
+int main()
+{
+	cin.imbue( std::locale("C") );
+	for(int n,np,nc,m; cin>>n>>np>>nc>>m;)
+	{
+		// 0: src, 1: dst, 2: ...
+		Graph g(n+2, vector<Capacity>(n+2));
+
+		while(m--)
+		{
+			int u, v, z; char _;
+			cin >> _ >> u >> _ >> v >> _ >> z;
+			g[u+2][v+2] = z;
+		}
+
+		while(np--)
+		{
+			int u, z; char _;
+			cin >> _ >> u >> _ >> z;
+			g[0][u+2] = z;
+		}
+
+		while(nc--)
+		{
+			int u, z; char _;
+			cin >> _ >> u >> _ >> z;
+			g[u+2][1] = z;
+		}
+
+		cout << goldberg_tarjan(g, 0, 1) << endl;
+	}
+}

ADDED   lib/graph/biMatch.cpp
Index: lib/graph/biMatch.cpp
==================================================================
--- lib/graph/biMatch.cpp
+++ lib/graph/biMatch.cpp
@@ -0,0 +1,40 @@
+
+//-------------------------------------------------------------
+// Bipertite max-matching (by duality, min-vert-cover)
+//   O(V E)
+//
+// Verified by
+//   SRM 397 Div1 Lv3
+//-------------------------------------------------------------
+
+typedef int           vert;
+typedef vert          edge;
+typedef vector<edge>  edges;
+typedef vector<edges> graph;
+
+bool augment( graph& G, int v, vector<vert>& matchTo, bool visited[] )
+{
+	for(int i=0; i<G[v].size(); ++i) {
+		vert u = G[v][i];
+		if( visited[u] ) continue;
+		visited[u] = true;
+
+		if( matchTo[u]<0 || augment(G, matchTo[u], matchTo, visited) )
+			{ matchTo[v]=u, matchTo[u]=v; return true; }
+	}
+	return false;
+}
+
+template<int NV>
+int biMatch( graph& G, int L ) // [0,L):left, [L,?):right
+    // only left->right edges are used during computation
+{
+	vector<vert> matchTo(G.size(), -1);
+	int ans = 0;
+	for(vert v=0; v<L; ++v) {
+		bool visited[NV] = {};
+		if( augment(G, v, matchTo, visited) )
+			++ans;
+	}
+	return ans;
+}

ADDED   lib/graph/isGraphical.cpp
Index: lib/graph/isGraphical.cpp
==================================================================
--- lib/graph/isGraphical.cpp
+++ lib/graph/isGraphical.cpp
@@ -0,0 +1,34 @@
+//-------------------------------------------------------------
+// Check the given list can be a degree list of some graph
+//   O(n^2)  (I suspect it can be made faster, though...)
+//
+// Verified by
+//   - SRM 398 Div1 LV3
+//
+//((
+// Havel-Hakimi
+//   If G[0], ..., G[n] (decreasing) is graphical,
+//   then G[1]-1, G[2]-1, ..., G[G[0]]-1, G[G[0]+1], .., G[n]
+//   is also graphical.
+//))
+//-------------------------------------------------------------
+
+bool isGraphical( vector<int> G )
+{
+	sort( G.begin(), G.end() );
+
+	vector<int>::iterator b = lower_bound( G.begin(), G.end(), 1 );
+	vector<int>::iterator e = G.end();
+
+	while( b < e )
+	{
+		int n = *(--e);
+		if( e-b < n )
+			return false;
+		for(vector<int>::iterator i=e-n; i!=e; ++i)
+			--*i;
+		inplace_merge( b, e-n, e );
+		b = lower_bound( G.begin(), G.end(), 1 );
+	}
+	return true;
+}

ADDED   lib/graph/maxFlow.cpp
Index: lib/graph/maxFlow.cpp
==================================================================
--- lib/graph/maxFlow.cpp
+++ lib/graph/maxFlow.cpp
@@ -0,0 +1,66 @@
+
+//-------------------------------------------------------------
+// Dinic's Algorithm
+//   O(V E)
+//
+// G : bidirectional (G[i].has(j) <==> G[j].has(i))
+// F : flow-capacity F[i][j] = Capacity, F[j][i] = 0
+//
+// Verified by
+//   - SRM 399 Div1 LV3
+//   - PKU 1459
+//   - CodeCraft 09 CUTS
+//   - SRM 465 Div1 LV2
+//-------------------------------------------------------------
+
+static const int NV = 512;
+typedef int           flow;
+typedef int           vert;
+typedef vert          edge;
+typedef vector<edge>  edges;
+typedef vector<edges> graph;
+typedef flow          flow_graph[NV][NV];
+
+flow dinic_dfs( graph& G, flow_graph F, vert v, vert D,
+                int LV[], flow flow_in, int blocked[] )
+{
+	flow flow_out = 0;
+	for(int i=0; i!=G[v].size(); ++i) {
+		int u = G[v][i];
+		if( LV[v]+1==LV[u] && F[v][u] ) {
+			flow f = min(flow_in-flow_out, F[v][u]);
+			if( u==D || !blocked[u] && (f=dinic_dfs(G,F,u,D,LV,f,blocked)) ) {
+				F[v][u]  -= f;
+				F[u][v]  += f;
+				flow_out += f;
+				if( flow_in == flow_out ) return flow_out;
+			}
+		}
+	}
+	blocked[v] = (flow_out==0);
+	return flow_out;
+}
+
+flow maxFlow( graph& G, flow_graph F, vert S, vert D )
+{
+	for( flow total=0 ;; ) {
+		int LV[NV] = {0};
+		vector<int> Q(1, S);
+		for(int lv=1; !Q.empty(); ++lv) {
+			vector<int> Q2;
+			for(int i=0; i!=Q.size(); ++i) {
+				edges& ne = G[Q[i]];
+				for(int j=0; j!=ne.size(); ++j)
+					if( F[Q[i]][ne[j]] && !LV[ne[j]] && ne[j]!=S )
+						LV[ne[j]]=lv, Q2.push_back(ne[j]);
+			}
+			Q.swap(Q2);
+		}
+
+		if( !LV[D] )
+			return total;
+
+		int blocked[NV] = {};
+		total += dinic_dfs( G, F, S, D, LV, 0x7fffffff, blocked );
+	}
+}

ADDED   lib/graph/mincostFlow.cpp
Index: lib/graph/mincostFlow.cpp
==================================================================
--- lib/graph/mincostFlow.cpp
+++ lib/graph/mincostFlow.cpp
@@ -0,0 +1,108 @@
+//-------------------------------------------------------------
+// MinCost-MaxFlow
+//   O(??)
+//
+// Verified by
+//   - SRM 487 Div2 LV3
+//   - SRM 491 Div1 LV3
+//-------------------------------------------------------------
+
+#include <iostream>
+#include <string>
+#include <vector>
+#include <map>
+#include <queue>
+using namespace std;
+
+template<typename T>
+class IdGen
+{
+	map<T, int> v2id_;
+	vector<T>   id2v_;
+public:
+	int v2id(const T& v) {
+		if( !v2id_.count(v) ) { v2id_[v] = size(); id2v_.push_back(v); }
+		return v2id_[v];
+	}
+	const T& id2v(int i) const { return id2v_[i]; }
+	int size() const { return id2v_.size(); }
+};
+
+template<typename Vert, typename Cost, typename Flow, int NV=256>
+class MinCostFlow
+{
+	IdGen<Vert> idgen;
+
+	vector<int> G[NV];
+	Flow F[NV][NV];
+	Cost C[NV][NV];
+
+public:
+	void addEdge( Vert s_, Vert t_, Cost c, Flow f )
+	{
+		int s = idgen.v2id(s_), t = idgen.v2id(t_);
+		G[s].push_back(t);
+		G[t].push_back(s);
+		C[s][t] = c;
+		C[t][s] = -c;
+		F[s][t] = f;
+		F[t][s] = 0;
+	}
+
+	pair<Cost, Flow> calc( Vert s_, Vert t_ )
+	{
+		const int N=idgen.size(), S=idgen.v2id(s_), T=idgen.v2id(t_);
+		static const Cost COST_INF = 1e+300; // !!EDIT HERE!!
+		static const Flow FLOW_INF = 0x7fffffff;
+
+		Cost total_cost = 0;
+		Flow total_flow = 0;
+		vector<Cost> h(N, 0); // potential
+		for(Flow RF=FLOW_INF; RF>0; ) // residual flow
+		{
+			// Dijkstra -- find the min-cost path
+			vector<Cost> d(N, COST_INF);  d[S] = 0;
+			vector<int>  prev(N, -1);
+
+			typedef pair< Cost, pair<int,int> > cedge;
+			priority_queue< cedge, vector<cedge>, greater<cedge> > Q;
+			Q.push( cedge(0, make_pair(S,S)) );
+			while( !Q.empty() ) {
+				cedge e = Q.top(); Q.pop();
+				if( prev[e.second.second] >= 0 )
+					continue;
+				prev[e.second.second] = e.second.first;
+
+				int u = e.second.second;
+				for(int i=0; i<G[u].size(); ++i) {
+					int  v = G[u][i];
+					Cost r_cost = C[u][v] + h[u] - h[v];
+					if( F[u][v] > 0 && d[v] > d[u]+r_cost )
+						Q.push( cedge(d[v]=d[u]+r_cost, make_pair(u,v)) );
+				}
+			}
+
+			if( prev[T] < 0 )
+				break; // Finished
+
+			// Run the flow as much as possible
+			Flow f = RF;
+			for(int u=T; u!=S; u=prev[u])
+				f = min(f, F[prev[u]][u]);
+			RF         -= f;
+			total_flow += f;
+
+			for(int u=T; u!=S; u=prev[u])
+			{
+				total_cost    += f * C[prev[u]][u];
+				F[prev[u]][u] -= f;
+				F[u][prev[u]] += f;
+			}
+
+			// Update the potential
+			for(int u=0; u<N; ++u)
+				h[u] += d[u];
+		}
+		return make_pair(total_cost, total_flow);
+	}
+};

ADDED   lib/graph/undAllPairMinCut.cpp
Index: lib/graph/undAllPairMinCut.cpp
==================================================================
--- lib/graph/undAllPairMinCut.cpp
+++ lib/graph/undAllPairMinCut.cpp
@@ -0,0 +1,99 @@
+
+//-------------------------------------------------------------
+// Gomory-Hu
+// All Pair Minimum Cuts for an Undirected Graph
+//   O(V^4)
+//
+// G : undirected (G[i].has(j) <==> G[j].has(i))
+// F : flow-capacity F[i][j] = F[j][i] = Capacity
+//
+// Verified by
+//   - CodeCraft 09 CUTS
+//-------------------------------------------------------------
+
+static const int NV = 512;
+typedef int           flow;
+typedef int           vert;
+typedef vert          edge;
+typedef vector<edge>  edges;
+typedef vector<edges> graph;
+typedef flow          flow_graph[NV][NV];
+
+flow_graph GHF;
+void GomoryHu( graph& G, flow_graph F, flow_graph AllPairMinCut )
+{
+	int N = G.size();
+	vector<vert> parent(N, 0);
+	vector<flow> weight(N, 0);
+
+	// construct the Gomory-Hu Tree
+	for(vert s=1; s<N; ++s)
+	{
+		vert t = parent[s];
+
+		// mincut between s and t
+		memcpy(GHF, F, sizeof(GHF));
+		flow mincut = 0;
+		for(;;) {
+			// bfs
+			vector<vert> prev(N, -1); prev[s]=s;
+			queue<vert> Q; Q.push(s);
+			while( !Q.empty() ) {
+				vert v = Q.front(); Q.pop();
+				for(int i=0; i<G[v].size(); ++i) {
+					vert u = G[v][i];
+					if( prev[u]==-1 && GHF[v][u] )
+						prev[u]=v, Q.push(u);
+				}
+			}
+
+			if( prev[t] == -1 )
+			{
+				// done: separeted to {v|prev[v]!=-1} and {v|prev[v]==-1}
+				// split t's children
+				weight[s] = mincut;
+				for(vert v=0; v<N; ++v)
+					if( parent[v]==t && prev[v]!=-1 && v!=s )
+						parent[v] = s;
+				vert pt = parent[t];
+				if( prev[pt]!=-1 )
+					parent[s]=pt, parent[t]=s, weight[s]=weight[t], weight[t]=mincut;
+				break;
+			}
+
+			// flow...
+			flow cur = 0x7fffffff;
+			for(vert v=t; v!=s; v=prev[v])
+				cur = min(cur, GHF[prev[v]][v]);
+			for(vert v=t; v!=s; v=prev[v])
+				GHF[prev[v]][v] -= cur, GHF[v][prev[v]] += cur;
+			mincut += cur;
+		}
+	}
+
+	// AllPairMinCuts over the G-H tree
+	for(vert s=0; s<N; ++s) {
+		vector<vert> ps_s;
+		for(vert v=s ;; v=parent[v])
+			{ ps_s.push_back(v); if( v==parent[v] ) break; }
+		reverse(ps_s.begin(), ps_s.end());
+
+		for(vert t=s+1; t<N; ++t) {
+			vector<vert> ps_t;
+			for(vert v=t ;; v=parent[v])
+				{ ps_t.push_back(v); if( v==parent[v] ) break; }
+			reverse(ps_t.begin(), ps_t.end());
+
+			vert ca;
+			for(int i=0; i<ps_s.size() && i<ps_t.size() && ps_s[i]==ps_t[i]; ++i)
+				ca = ps_s[i];
+
+			flow s_to_root = 0x7fffffff, t_to_root = 0x7fffffff;
+			for(vert v=s; v!=ca; v=parent[v])
+				s_to_root = min(s_to_root, weight[v]);
+			for(vert v=t; v!=ca; v=parent[v])
+				t_to_root = min(t_to_root, weight[v]);
+			AllPairMinCut[s][t] = AllPairMinCut[t][s] = min(s_to_root, t_to_root);
+		}
+	}
+}

ADDED   lib/graph/undMinCut.cpp
Index: lib/graph/undMinCut.cpp
==================================================================
--- lib/graph/undMinCut.cpp
+++ lib/graph/undMinCut.cpp
@@ -0,0 +1,68 @@
+
+//-------------------------------------------------------------
+// Stoer-Wagner
+// Minimum cost for making the graph unconnected
+//   O(V^3)
+//
+// G : bidirectional (G[v][u] == G[u][v])
+//
+// Verified by
+//   - SRM 340 Div1 LV3
+//-------------------------------------------------------------
+
+typedef int cost;
+typedef int vert;
+typedef vector< vector<cost> > graph;
+
+cost bidi_minCut( graph G, int N )
+{
+	vector<vert> V;
+	for(int v=0; v<N; ++v)
+		V.push_back(v);
+
+	cost result = 0x7fffffff;
+	for(int n=N; n>1; --n)
+	{
+		// invariant:
+		//   the current set of vertices = {V[0] .. V[n-1]}
+
+		// order the vertices
+		//   v0 = 0 (arbitrary),
+		//   v1 = (u that maximizes \Sigma G[v0][u])
+		//   v2 = (u that maximizes \Sigma G[v0][u]+G[v1][u])
+		//   ...
+		//   vn-2
+		//   vn-1 = (maximizes lastCut = \Sigma G[v0][u]+G[vn-2][u])
+		vector<int> vs;
+		cost lastCut = 0;
+		{
+			vector<cost> ws(n, 0);
+			ws[0] = 1;
+
+			for(int i=0; i<n; ++i) {
+				int j = max_element(ws.begin(), ws.end()) - ws.begin();
+				lastCut = ws[j];
+
+				vs.push_back(j);
+				ws[j] = -1;
+				for(int k=0; k<n; ++k)
+					if( ws[k] != -1 )
+						ws[k] += G[V[k]][V[j]];
+			}
+		}
+
+		// update mincut
+		//   lemma: "the min cost for separating vn-2 and vn-1"==lastCut
+		result = min(result, lastCut);
+
+		// reduce the graph (unify vn-2 and vn-1)
+		//   for testing the case vn-2 and vn-1 is connected
+		vert v2=vs[n-2], v1=vs[n-1];
+		for(int i=0; i<n; ++i) {
+			G[V[v2]][V[i]] += G[V[v1]][V[i]];
+			G[V[i]][V[v2]] += G[V[i]][V[v1]];
+		}
+		V.erase( V.begin() + v1 );
+	}
+	return result;
+}

ADDED   lib/numeric/bigintFake.cpp
Index: lib/numeric/bigintFake.cpp
==================================================================
--- lib/numeric/bigintFake.cpp
+++ lib/numeric/bigintFake.cpp
@@ -0,0 +1,14 @@
+
+string add(const string& a, const string& b)
+{
+	int n = max(a.size(), b.size()), carry=0;
+	string c(n, '0');
+	for(int i=0; i<n; ++i) {
+		int x = (a.size()<=i ? 0 : a[a.size()-1-i]-'0')
+		      + (b.size()<=i ? 0 : b[b.size()-1-i]-'0') + carry;
+		c[n-1-i] = char('0'+x%10);
+		carry = x/10;
+	}
+	if( carry ) c = char('0'+carry)+c;
+	return c;
+}

ADDED   lib/numeric/erathos.cpp
Index: lib/numeric/erathos.cpp
==================================================================
--- lib/numeric/erathos.cpp
+++ lib/numeric/erathos.cpp
@@ -0,0 +1,9 @@
+	static const int N = 999999;
+	vector<bool> isp(N+1, true);
+	vector<int> ps;
+	for(int p=2; p<=N; ++p)
+		if( isp[p] ) {
+			ps.push_back(p);
+			for(int q=p+p; q<=N; q+=p)
+				isp[q] = false;
+		}

ADDED   lib/numeric/fft.cpp
Index: lib/numeric/fft.cpp
==================================================================
--- lib/numeric/fft.cpp
+++ lib/numeric/fft.cpp
@@ -0,0 +1,39 @@
+
+//-------------------------------------------------------------
+// Fast Discrete Fourier Transform
+//   O(n log n)
+//   n must be a power of 2.
+//
+// Verified by
+//   - SRM 436 LV3
+//-------------------------------------------------------------
+
+CMP tmp[65536*4];
+
+template<int F>
+void fft_impl( CMP a[], int n, int stride = 1 )
+{
+	if( n > 1 )
+	{
+		CMP *ev=a, *od=a+stride;
+		int s2=stride*2, n2=n/2;
+
+		fft_impl<F>( ev, n2, s2 );
+		fft_impl<F>( od, n2, s2 );
+
+		for(int i=0; i<n; ++i) tmp[i] = ev[s2*(i%n2)] + od[s2*(i%n2)]*polar(1.0, F*2*M_PI*i/n);
+		for(int i=0; i<n; ++i) a[stride*i] = tmp[i];
+	}
+}
+
+void fft( vector<CMP>& a )
+{
+	fft_impl<+1>(&a[0], a.size());
+}
+
+void ifft( vector<CMP>& a )
+{
+	fft_impl<-1>(&a[0], a.size());
+	for(int i=0; i<a.size(); ++i)
+		a[i] /= a.size();
+}

ADDED   lib/numeric/gcd_lcm.cpp
Index: lib/numeric/gcd_lcm.cpp
==================================================================
--- lib/numeric/gcd_lcm.cpp
+++ lib/numeric/gcd_lcm.cpp
@@ -0,0 +1,28 @@
+LL gcd(LL a, LL b)
+{
+	while(a)
+		swap(a, b%=a);
+	return b;
+}
+
+LL lcm(LL a, LL b)
+{
+	return a/gcd(a,b)*b;
+}
+
+// Assumption: (a/b)*b+(a%b) == a
+//   For typical C/C++ compilers, % on negatives satisfies this.
+//   To be sure, replace a%b by a-a/b*b
+
+LL xgcd(LL a, LL b, LL* x, LL* y) // ax+by=g
+{
+	if(b) {
+		LL yy, g = xgcd(b,a%b,&yy,x);
+		*y = yy - a/b**x;
+		return g;
+	}
+	else {
+		*x=1, *y=0;
+		return a;
+	}
+}

ADDED   lib/numeric/linearEq.cpp
Index: lib/numeric/linearEq.cpp
==================================================================
--- lib/numeric/linearEq.cpp
+++ lib/numeric/linearEq.cpp
@@ -0,0 +1,40 @@
+
+//-------------------------------------------------------------
+// Linear Equation Solver
+//   O(n^3)
+//
+// Verified by
+//   - SRM 398 Div1 LV3
+//   - AOJ 0004
+//-------------------------------------------------------------
+
+vector<double> solve_linear_eq( int n, vector< vector<double> > M, const vector<double>& V )
+{
+	vector<double> A(V);
+	for(int i=0; i<n; ++i)
+	{
+		// pivot
+		if( M[i][i] == 0 )
+			for(int j=i+1; j<n; ++j)
+				if( M[j][i] != 0 )
+					{swap(M[i], M[j]); swap(A[i], A[j]); break;}
+		if( M[i][i] == 0 )
+			throw "no anser";
+
+		// M[i][i] <-- 1
+		double p = M[i][i];
+		for(int j=i; j<n; ++j)
+			M[i][j] /= p;
+		A[i] /= p;
+
+		// M[*][i] <-- 0
+		for(int j=0; j<n; ++j) if(j!=i)
+		{
+			double r = M[j][i];
+			for(int k=i; k<n; ++k)
+				M[j][k] -= M[i][k] * r;
+			A[j] -= A[i] * r;
+		}
+	}
+	return A;
+}

ADDED   lib/numeric/linearEq2.cpp
Index: lib/numeric/linearEq2.cpp
==================================================================
--- lib/numeric/linearEq2.cpp
+++ lib/numeric/linearEq2.cpp
@@ -0,0 +1,71 @@
+//-------------------------------------------------------------
+// Linear Equation Solver for Teplitz Matrix
+//   O(n^2)
+//
+// Faster solver for matrices satisfying
+//   M[i][j] == M[i+1][j+1]
+// for all i, j>0
+// Be careful that this implementation does not support the
+// case M[0][0]=M[1][1]=...=0. (zero-division)
+//
+// TODO: support M[0][0]=M[1][1]=...=0
+// QUICKHACK:
+//    find the least degenerated k and 
+//    calc initial vectors by gauss_jordan...
+//
+// Verified by
+//   - (Checked by random data to match with linearEq.cpp)
+//-------------------------------------------------------------
+
+vector<double> solve_teplitz( int n, vector< vector<double> > M, const vector<double>& V )
+{
+	vector<double> f, b, x;
+
+	// invariants
+	//   M|k f|k = (1 0 ... 0)
+	//   M|k b|k = (0 ... 0 1)
+	//   M|k x|k = V|k
+	f.push_back(    1 / M[0][0] );
+	b.push_back(    1 / M[0][0] );
+	x.push_back( V[0] / M[0][0] );
+
+	for(int k=2; k<=n; ++k)
+	{
+		// M|k (f|k-1 0) = (1 0 ... 0 ea)
+		double ea = 0;
+		for(int i=0; i<k-1; ++i)
+			ea += M[k-1][i] * f[i];
+
+		// M|k (0 b|k-1) = (eb 0 ... 0 1)
+		double eb = 0;
+		for(int i=1; i<k; ++i)
+			eb += M[0][i] * b[i-1];
+
+		// u(1 ea) + v(eb 1) = (1 0)
+		//==> u = 1 / (1-ea.eb),  v = -ea /(1-ea.eb)
+		double u = 1/(1-ea*eb), v = -ea/(1-ea*eb);
+		vector<double> ff;
+		for(int i=0; i<k; ++i)
+			ff.push_back( u*(i<k-1 ? f[i] : 0) + v*(i>0 ? b[i-1] : 0) );
+
+		// similarly...
+		u = -eb/(1-ea*eb), v = 1/(1-ea*eb);
+		vector<double> bb;
+		for(int i=0; i<k; ++i)
+			bb.push_back( u*(i<k-1 ? f[i] : 0) + v*(i>0 ? b[i-1] : 0) );
+
+		// update
+		f.swap(ff);
+		b.swap(bb);
+
+		// M|k (x|k 0) (V|k ec)
+		double ec = 0;
+		for(int i=0; i<k-1; ++i)
+			ec += M[k-1][i] * x[i];
+		x.push_back(0);
+		for(int i=0; i<k; ++i)
+			x[i] -= b[i] * (ec-V[k-1]);
+	}
+
+	return x;
+}

ADDED   lib/numeric/matrix.cpp
Index: lib/numeric/matrix.cpp
==================================================================
--- lib/numeric/matrix.cpp
+++ lib/numeric/matrix.cpp
@@ -0,0 +1,55 @@
+
+//-------------------------------------------------------------
+// Matrix Operations
+//
+// Verified by
+//   - SRM342 Div1 LV3
+//   - SRM341 Div1 LV3
+//   - SRM338 Div1 LV2
+//-------------------------------------------------------------
+
+vector<LL> vMul( const vector< vector<LL> >& A, const vector<LL>& B )
+{
+	const int n = A.size();
+
+	vector<LL> C(n);
+	for(int i=0; i<n; ++i)
+		for(int j=0; j<n; ++j)
+			C[i] += A[i][j] * B[j];
+	return C;
+}
+
+vector< vector<LL> > mMul( const vector< vector<LL> >& A, const vector< vector<LL> >& B )
+{
+	const int n = A.size();
+
+	vector< vector<LL> > C(n, vector<LL>(n));
+	for(int i=0; i<n; ++i)
+		for(int j=0; j<n; ++j) {
+			LL Cij = 0;
+			for(int k=0; k<n; ++k)
+				Cij += A[i][k] * B[k][j];
+			C[i][j] = Cij;
+		}
+	return C;
+}
+
+vector< vector<LL> > mPow( vector< vector<LL> > M, LL t ) // t>0
+{
+	vector< vector<LL> > R;
+	for(; t; t>>=1, M=mMul(M,M))
+		if( t&1 )
+			R = (R.empty() ? M : mMul(R,M));
+	return R;
+}
+
+vector< vector<LL> > mAdd( const vector< vector<LL> >& A, const vector< vector<LL> >& B )
+{
+	const int n = A.size();
+
+	vector< vector<LL> > C(n, vector<LL>(n));
+	for(int i=0; i<n; ++i)
+		for(int j=0; j<n; ++j)
+			C[i][j] = A[i][j] + B[i][j];
+	return C;
+}

ADDED   lib/numeric/mebius.cpp
Index: lib/numeric/mebius.cpp
==================================================================
--- lib/numeric/mebius.cpp
+++ lib/numeric/mebius.cpp
@@ -0,0 +1,36 @@
+
+//-------------------------------------------------------------
+// Mebius Function
+//   O( sqrt(N) )
+//
+// mebius(n) = 0      if n=0 mod p*p for some p
+// mebius(n) = (-1)^k otherwise, where k is the number of n's prime factors
+//
+// properties:
+//  1: mebius(nm) = mebius(n) mebius(m)   if gcd(n,m)=1
+//  2: f(n) = Sigma_{d|n} g(d)
+//       <=>
+//     g(n) = Sigma_{d|n} f(d)mebius(n/d)
+//
+// Verified by
+//   - CodeCraft '11 PSTR
+//
+// Further memoization or more might be needed.
+//-------------------------------------------------------------
+
+int mebius(int n)
+{
+	int f = 1;
+	for(int q=2; q*q<=n; ++q)
+		if( n%q == 0 )
+		{
+			if( n%(q*q) == 0 )
+				return 0;
+			n /= q;
+			f *= -1;
+		}
+	if(n>1)
+		f *= -1;
+	return f;
+}
+

ADDED   lib/numeric/modArith.cpp
Index: lib/numeric/modArith.cpp
==================================================================
--- lib/numeric/modArith.cpp
+++ lib/numeric/modArith.cpp
@@ -0,0 +1,93 @@
+
+//-------------------------------------------------------------
+// Modulo Arithmetics
+//
+// Verified by
+//   - TCO10 R3 LV3
+//-------------------------------------------------------------
+
+static const int MODVAL = 1000000007; // op/ depends on primarity of the MODVAL
+struct mint
+{
+	int val;
+	mint():val(0){}
+	mint(int x):val(x%MODVAL) {}
+	mint(LL  x):val(x%MODVAL) {}
+};
+
+mint operator+(mint x, mint y) { return x.val+y.val; }
+mint operator-(mint x, mint y) { return x.val-y.val+MODVAL; }
+mint operator*(mint x, mint y) { return LL(x.val)*y.val; }
+mint POW(mint x, int e) {
+	mint v = 1;
+	for(;e;x=x*x,e>>=1)
+		if(e&1)
+			v=v*x;
+	return v;
+}
+mint operator/(mint x, mint y) { return x * POW(y, MODVAL-2); }
+
+vector<mint> FAC_(1,1);
+void FAC_INIT(int n) { for(int i=1; i<=n; ++i) FAC_.push_back( FAC_.back()*i ); }
+mint FAC(mint x)       { return FAC_[x.val]; }
+mint C(mint n, mint k) { return k.val<0 || n.val<k.val ? 0 : FAC(n) / (FAC(k) * FAC(n-k)); }
+
+/*
+// MODVAL must be a prime!!
+LL GSS(LL k, LL b, LL e) // k^b + k^b+1 + ... + k^e
+{
+	if( b >  e ) return 0;
+	if( k <= 1 ) return k*(e-b+1);
+	return DIV(SUB(POW(k, e+1), POW(k,b)), k-1);
+}
+
+LL Cpascal(LL n, LL k)
+{
+	vector< vector<LL> > c(n+1, vector<LL>(k+1));
+	for(LL nn=1; nn<=n; ++nn)
+		for(LL kk=0; kk<=min(nn,k); ++kk)
+			c[nn][kk] = kk==0 || kk==nn ? 1 : ADD(c[nn-1][kk-1], c[nn-1][kk]);
+	return c[n][k];
+}
+
+vector< vector<LL> > MATMUL(vector< vector<LL> >& a, vector< vector<LL> >& b)
+{
+   int N = a.size();
+   vector< vector<LL> > c(N, vector<LL>(N));
+   for(int i=0; i<N; ++i)
+   for(int j=0; j<N; ++j)
+   for(int k=0; k<N; ++k)
+      c[i][j] = ADD(c[i][j], MUL(a[i][k],b[k][j]));
+   return c;
+}
+
+// works for non-prime MODVAL
+LL GEO(LL x_, LL e) // x^0 + x^1 + ... + x^e-1
+{
+   vector< vector<LL> > v(2, vector<LL>(2));
+   vector< vector<LL> > x(2, vector<LL>(2));
+   v[0][0] = v[1][1] = 1;
+   x[0][0] = x_; x[0][1] = 0;
+   x[1][0] = 1 ; x[1][1] = 1;
+   for(;e;x=MATMUL(x,x),e>>=1)
+      if(e&1)
+         v = MATMUL(v, x);
+   return v[1][0];
+}
+
+// works for non-prime MODVAL
+LL HYP(LL x_, LL e) // e x^0 + (e-1) x^1 + ... + 1 x^e-1 = GEO(x,1)+GEO(x,2)+...+GEO(x,e)
+{
+   vector< vector<LL> > v(3, vector<LL>(3));
+   vector< vector<LL> > x(3, vector<LL>(3));
+   v[0][0] = v[1][1] = v[2][2] = 1;
+   x[0][0] = x_; x[0][1] = 0; x[0][2] = 0;
+   x[1][0] = 1 ; x[1][1] = 1; x[1][2] = 0;
+   x[2][0] = 0 ; x[2][1] = 1; x[2][2] = 1;
+   e++;
+   for(;e;x=MATMUL(x,x),e>>=1)
+      if(e&1)
+         v = MATMUL(v, x);
+   return v[2][0];
+}
+*/

ADDED   lib/numeric/modArithOld.cpp
Index: lib/numeric/modArithOld.cpp
==================================================================
--- lib/numeric/modArithOld.cpp
+++ lib/numeric/modArithOld.cpp
@@ -0,0 +1,94 @@
+
+//-------------------------------------------------------------
+// Modulo Arithmetics
+//
+// Verified by
+//   - SRM 397 Div1 LV2
+//   - SRM 428 Div1 LV2
+//   - CodeCraft 2010 CNTINT DRAW
+//-------------------------------------------------------------
+
+static const LL MODVAL = 1000000007; // must fit in 32-bits
+
+LL ADD(LL x, LL y) { return (x+y)%MODVAL; }
+LL SUB(LL x, LL y) { return (x-y+MODVAL)%MODVAL; }
+LL MUL(LL x, LL y) { return (x*y)%MODVAL; }
+LL POW(LL x, LL e) {
+	LL v = 1;
+	for(;e;x=MUL(x,x),e>>=1)
+		if(e&1)
+			v = MUL(v, x);
+	return v;
+}
+
+// MODVAL must be a prime!!
+LL DIV(LL x, LL y) { return MUL(x, POW(y, MODVAL-2)); }
+
+// MODVAL must be a prime!!
+LL C(LL n, LL k) {
+	LL v = 1;
+	for(LL i=1; i<=k; ++i)
+		v = DIV(MUL(v, n-i+1), i);
+	return v;
+}
+
+// MODVAL must be a prime!!
+LL GSS(LL k, LL b, LL e) // k^b + k^b+1 + ... + k^e
+{
+	if( b >  e ) return 0;
+	if( k <= 1 ) return k*(e-b+1);
+	return DIV(SUB(POW(k, e+1), POW(k,b)), k-1);
+}
+
+
+
+
+LL Cpascal(LL n, LL k)
+{
+	vector< vector<LL> > c(n+1, vector<LL>(k+1));
+	for(LL nn=1; nn<=n; ++nn)
+		for(LL kk=0; kk<=min(nn,k); ++kk)
+			c[nn][kk] = kk==0 || kk==nn ? 1 : ADD(c[nn-1][kk-1], c[nn-1][kk]);
+	return c[n][k];
+}
+
+vector< vector<LL> > MATMUL(vector< vector<LL> >& a, vector< vector<LL> >& b)
+{
+   int N = a.size();
+   vector< vector<LL> > c(N, vector<LL>(N));
+   for(int i=0; i<N; ++i)
+   for(int j=0; j<N; ++j)
+   for(int k=0; k<N; ++k)
+      c[i][j] = ADD(c[i][j], MUL(a[i][k],b[k][j]));
+   return c;
+}
+
+// works for non-prime MODVAL
+LL GEO(LL x_, LL e) // x^0 + x^1 + ... + x^e-1
+{
+   vector< vector<LL> > v(2, vector<LL>(2));
+   vector< vector<LL> > x(2, vector<LL>(2));
+   v[0][0] = v[1][1] = 1;
+   x[0][0] = x_; x[0][1] = 0;
+   x[1][0] = 1 ; x[1][1] = 1;
+   for(;e;x=MATMUL(x,x),e>>=1)
+      if(e&1)
+         v = MATMUL(v, x);
+   return v[1][0];
+}
+
+// works for non-prime MODVAL
+LL HYP(LL x_, LL e) // e x^0 + (e-1) x^1 + ... + 1 x^e-1 = GEO(x,1)+GEO(x,2)+...+GEO(x,e)
+{
+   vector< vector<LL> > v(3, vector<LL>(3));
+   vector< vector<LL> > x(3, vector<LL>(3));
+   v[0][0] = v[1][1] = v[2][2] = 1;
+   x[0][0] = x_; x[0][1] = 0; x[0][2] = 0;
+   x[1][0] = 1 ; x[1][1] = 1; x[1][2] = 0;
+   x[2][0] = 0 ; x[2][1] = 1; x[2][2] = 1;
+   e++;
+   for(;e;x=MATMUL(x,x),e>>=1)
+      if(e&1)
+         v = MATMUL(v, x);
+   return v[2][0];
+}

ADDED   lib/numeric/nextComb.cpp
Index: lib/numeric/nextComb.cpp
==================================================================
--- lib/numeric/nextComb.cpp
+++ lib/numeric/nextComb.cpp
@@ -0,0 +1,15 @@
+//-------------------------------------------------------------
+// Next Combination
+//
+// Verified by
+//   - SRM345 Div1 LV3
+//-------------------------------------------------------------
+
+LL next_combination(LL p)
+{
+	assert( p > 0 );
+	LL lsb = p & -p;
+	LL rem = p + lsb;
+	LL rit = rem & ~p;
+	return rem | (rit/lsb >> 1)-1;
+}

ADDED   lib/numeric/sternBrocot.cpp
Index: lib/numeric/sternBrocot.cpp
==================================================================
--- lib/numeric/sternBrocot.cpp
+++ lib/numeric/sternBrocot.cpp
@@ -0,0 +1,52 @@
+
+//-------------------------------------------------------------
+// Enumerate all positive reduced fractions in a 
+// systematic manner.
+//   [pl/ql, pr/qr] --> [pl/ql, pm/qm] + [pm/qm, pr/qr]
+//                  --> ...
+//
+// Verified by
+//   - CodeCraft 09 FRACTION
+//-------------------------------------------------------------
+
+void sternBrocot(LL pl = 0, LL ql = 1, LL pr = 1, LL qr = 0)
+{
+	LL pm = pl + pr, qm = ql + qr;
+	sternBrocot(pl, ql, pm, qm);
+	sternBrocot(pm, qm, pr, qr);
+}
+
+/*
+bool fleq(LL p1, LL q1, LL p2, LL q2)
+{
+	return p1*q2 <= p2*q1;
+}
+
+void sternBrocot(LL pl = 0, LL ql = 1, LL pr = 1, LL qr = 0)
+{
+	for(;;) {
+		LL pm = pl + pr, qm = ql + qr;
+
+		if( fleq(c,d,pm,qm) ) {
+//			pr = pm, qr = qm; // [pl/ql, pm/qm]
+			LL X = c*ql - d*pl;
+			LL Y = d*pr - c*qr;
+			LL k = max(1LL, Y/X);
+			pr += k*pl;
+			qr += k*ql;
+		}
+		else if( fleq(pm,qm,a,b) ) {
+//			pl = pm, ql = qm; // [pm/qm, pr/qr]
+			LL X = b*pr - a*qr;
+			LL Y = a*ql - b*pl;
+			LL k = max(1LL, Y/X);
+			pl += k*pr;
+			ql += k*qr;
+		}
+		else {
+			cout << pm+qm*base << "/" << qm << endl;
+			break;
+		}
+	}
+}
+*/

ADDED   lib/typical/amap.cpp
Index: lib/typical/amap.cpp
==================================================================
--- lib/typical/amap.cpp
+++ lib/typical/amap.cpp
@@ -0,0 +1,72 @@
+//
+// accumulative map
+//   basically a map<K,V>, but it eliminates duplication lazily
+//   first keep all inserted elements in a vector. only when the buffer
+//   overflows, it eliminates duplication
+//
+
+template<typename K, typename V, int LIMIT = (1<<26)/3/sizeof(pair<K,V>)>
+struct amap
+{
+	typedef typename vector< pair<K,V> >::iterator iterator;
+
+	vector< pair<K,V> > kv;
+	amap() { kv.reserve(LIMIT); }
+
+	void add(const K& k, const V& v)
+	{
+		kv.push_back( make_pair(k,v) );
+		if( kv.size() == LIMIT )
+			normalize();
+	}
+
+	iterator begin()       { return kv.begin(); }
+	iterator end()         { return kv.end(); }
+	void swap( amap& rhs ) { kv.swap(rhs.kv); }
+
+	// Tested: SRM 469 Lv3
+	void normalize()
+	{
+		sort(kv.begin(), kv.end());
+		int i=0;
+		for(int j=0; j<kv.size(); ++i)
+		{
+			int k = j;
+			kv[i] = kv[k];
+			while( ++j<kv.size() && kv[k].first==kv[j].first )
+				kv[i].second += kv[j].second;
+		}
+		kv.resize(i);
+	}
+/*
+	// Not Tested (Prefer First)
+	void normalize()
+	{
+		sort(kv.begin(), kv.end());
+		int i=0;
+		for(int j=0; j<kv.size(); ++i)
+		{
+			int k = j;
+			kv[i] = kv[k];
+			while( ++j<kv.size() && kv[k].first==kv[j].first )
+				{}
+		}
+		kv.resize(i);
+	}
+
+	// Not Tested (Prefer Last)
+	void normalize()
+	{
+		sort(kv.begin(), kv.end());
+		int i=0;
+		for(int j=0; j<kv.size(); ++i)
+		{
+			int k = j;
+			while( ++j<kv.size() && kv[k].first==kv[j].first )
+				{}
+			kv[i] = kv[j-1];
+		}
+		kv.resize(i);
+	}
+*/
+};

ADDED   lib/typical/bfs.cpp
Index: lib/typical/bfs.cpp
==================================================================
--- lib/typical/bfs.cpp
+++ lib/typical/bfs.cpp
@@ -0,0 +1,26 @@
+// code template for BFS
+
+		State start = /*start*/;
+		vector<State> Q( 1, start );
+		set<State>    V; V.insert(start);
+
+		for(int step=0; !Q.empty(); ++step)
+		{
+			vector<State> Qold; Qold.swap(Q);
+			for(int qi=0; qi<Qold.size(); ++qi)
+			{
+				State& cur = Qold[qi];
+				if( /*isGoal(cur)*/ )
+				{
+				}
+
+				foreach(next)
+				{
+					if( !V.count(next) )
+					{
+						V.insert(next);
+						Q.push_back(next);
+					}
+				}
+			}
+		}

ADDED   lib/typical/bitop.cpp
Index: lib/typical/bitop.cpp
==================================================================
--- lib/typical/bitop.cpp
+++ lib/typical/bitop.cpp
@@ -0,0 +1,7 @@
+int bitcnt(LL x)
+{
+	int c = 0;
+	for(; x; x>>=1)
+		c += x&1;
+	return c;
+}

ADDED   lib/typical/comb.cpp
Index: lib/typical/comb.cpp
==================================================================
--- lib/typical/comb.cpp
+++ lib/typical/comb.cpp
@@ -0,0 +1,18 @@
+//-------------------------------------------------------------
+// number of combinations choosing k out of n
+// # you might better consider to use Pascal's triangle
+// # for comb modulo some number...
+//
+// Verified by
+//   - SRM 350 Div1 LV2
+//-------------------------------------------------------------
+
+LL comb(LL n, LL k)
+{
+	k = min(k, n-k);
+
+	LL c = 1;
+	for(LL i=0; i<k; ++i)
+		c *= n-i, c /= i+1;
+	return c;
+}

ADDED   lib/typical/dp.cpp
Index: lib/typical/dp.cpp
==================================================================
--- lib/typical/dp.cpp
+++ lib/typical/dp.cpp
@@ -0,0 +1,112 @@
+// Tested: SRM 454 Lv2
+template<typename T>
+struct DP2
+{
+	const int N1, N2;
+	vector<T> data;
+	DP2(int N1, int N2, const T& t = T())
+		: N1(N1), N2(N2), data(N1*N2, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2)
+		{ return data[ (i1*N2)+i2 ]; }
+	void swap(DP2& rhs)
+		{ data.swap(rhs.data); }
+};
+
+// Tested: Codeforces #13 C
+template<typename T>
+struct DP2x
+{
+	const int N1, N2;
+	vector<T> data;
+	DP2x(int, int N2, const T& t = T())
+		: N1(2), N2(N2), data(N1*N2, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2)
+		{ i1&=1; return data[ (i1*N2)+i2 ]; }
+	void swap(DP2x& rhs)
+		{ data.swap(rhs.data); }
+};
+
+// Tested: SRM 452 Lv3
+template<typename T>
+struct DP3
+{
+	int N1, N2, N3;
+	vector<T> data;
+	DP3(int N1, int N2, int N3, const T& t = T())
+		: N1(N1), N2(N2), N3(N3), data(N1*N2*N3, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2, int i3)
+		{ return data[ ((i1*N2)+i2)*N3+i3 ]; }
+	void swap(DP3& rhs)
+		{ data.swap(rhs.data); }
+};
+
+// Tested: SRM 468 Lv2
+template<typename T>
+struct DP3x
+{
+	int N1, N2, N3;
+	vector<T> data;
+	DP3x(int, int N2, int N3, const T& t = T())
+		: N1(2), N2(N2), N3(N3), data(N1*N2*N3, t) { assert(data.size()*sizeof(T) < (1<<26)); }
+	T& operator()(int i1, int i2, int i3)
+		{ i1&=1; return data[ ((i1*N2)+i2)*N3+i3 ]; }
+	void swap(DP3x& rhs)
+		{ data.swap(rhs.data); }
+};
+
+// Not Tested
+template<typename T>
+struct DP4
+{
+	int N1, N2, N3, N4;
+	vector<T> data;
+	DP4(int N1, int N2, int N3, int N4, const T& t = T())
+		: N1(N1), N2(N2), N3(N3), N4(N4), data(N1*N2*N3*N4, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2, int i3, int i4)
+		{ return data[ (((i1*N2)+i2)*N3+i3)*N4+i4 ]; }
+	void swap(DP4& rhs)
+		{ data.swap(rhs.data); }
+};
+
+// Not Tested
+template<typename T>
+struct DP4x
+{
+	int N1, N2, N3, N4;
+	vector<T> data;
+	DP4x(int, int N2, int N3, int N4, const T& t = T())
+		: N1(2), N2(N2), N3(N3), N4(N4), data(N1*N2*N3*N4, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2, int i3, int i4)
+		{ i1&=1; return data[ (((i1*N2)+i2)*N3+i3)*N4+i4 ]; }
+	void swap(DP4x& rhs)
+		{ data.swap(rhs.data); }
+};
+
+
+// Tested: SRM 351 Lv2
+template<typename T>
+struct DP5
+{
+	int N1, N2, N3, N4, N5;
+	vector<T> data;
+	DP5(int N1, int N2, int N3, int N4, int N5, const T& t = T())
+		: N1(N1), N2(N2), N3(N3), N4(N4), N5(N5), data(N1*N2*N3*N4*N5, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2, int i3, int i4, int i5)
+		{ return data[ ((((i1*N2)+i2)*N3+i3)*N4+i4)*N5+i5 ]; }
+	void swap(DP5& rhs)
+		{ data.swap(rhs.data); }
+};
+
+// Not Tested
+template<typename T>
+struct DP5x
+{
+	int N1, N2, N3, N4, N5;
+	vector<T> data;
+	DP5x(int, int N2, int N3, int N4, int N5, const T& t = T())
+		: N1(2), N2(N2), N3(N3), N4(N4), N5(N5), data(N1*N2*N3*N4*N5, t) { assert(data.size()*sizeof(T)<(1<<26)); }
+	T& operator()(int i1, int i2, int i3, int i4, int i5)
+		{ i1&=1; return data[ ((((i1*N2)+i2)*N3+i3)*N4+i4)*N5+i5 ]; }
+	void swap(DP5x& rhs)
+		{ data.swap(rhs.data); }
+};

ADDED   lib/typical/idGen.cpp
Index: lib/typical/idGen.cpp
==================================================================
--- lib/typical/idGen.cpp
+++ lib/typical/idGen.cpp
@@ -0,0 +1,21 @@
+//-------------------------------------------------------------
+// ID Assignment
+//
+// Verified by
+//   - ACM/ICPC Tokyo 2010 A
+//   - SRM 491 Div1 LV3
+//-------------------------------------------------------------
+
+template<typename T>
+class IdGen
+{
+	map<T, int> v2id_;
+	vector<T>   id2v_;
+public:
+	int v2id(const T& v) {
+		if( !v2id_.count(v) ) { v2id_[v] = size(); id2v_.push_back(v); }
+		return v2id_[v];
+	}
+	const T& id2v(int i) const { return id2v_[i]; }
+	int size() const { return id2v_.size(); }
+};

ADDED   lib/typical/inMap.cpp
Index: lib/typical/inMap.cpp
==================================================================
--- lib/typical/inMap.cpp
+++ lib/typical/inMap.cpp
@@ -0,0 +1,5 @@
+
+bool inMap(int H, int W, int y, int x)
+{
+	return (0<=y && y<H && 0<=x && x<W);
+}

ADDED   lib/typical/maxRectArea.cpp
Index: lib/typical/maxRectArea.cpp
==================================================================
--- lib/typical/maxRectArea.cpp
+++ lib/typical/maxRectArea.cpp
@@ -0,0 +1,39 @@
+//-------------------------------------------------------------
+// h[] : list of heights
+//   i-th rectangle is located at (i,0)--(i+1,h[i])
+//
+// calculate the area of maximum sub-rectangle
+//
+// Verified by
+//   - SRM 337 Div1 LV2
+//-------------------------------------------------------------
+
+		// solve
+		vector<int> left(n);
+		{
+			map<LL, int> h; h[-1] = -1;
+			for(int i=0; i<n; ++i) {
+				// position of the highest building < R[i]
+				map<LL,int>::iterator it = h.lower_bound(R[i]);
+				left[i] = (--it)->second+1;
+				h.erase( ++it, h.end() );
+				h[ R[i] ] = i;
+			}
+		}
+		vector<int> right(n);
+		{
+			map<LL, int> h; h[-1] = n;
+			for(int i=n-1; i>=0; --i) {
+				// position of the highest building < R[i]
+				map<LL,int>::iterator it = h.lower_bound(R[i]);
+				right[i] = (--it)->second-1;
+				h.erase( ++it, h.end() );
+				h[ R[i] ] = i;
+			}
+		}
+		LL ans = 0;
+		for(int i=0; i<n; ++i) {
+			LL area = R[i] * (right[i] - left[i] + 1);
+			ans = max(ans, area);
+		}
+		return ans;

ADDED   lib/typical/setCover.cpp
Index: lib/typical/setCover.cpp
==================================================================
--- lib/typical/setCover.cpp
+++ lib/typical/setCover.cpp
@@ -0,0 +1,69 @@
+
+//-------------------------------------------------------------
+// (Minimum) Set Cover
+//   NP-complete
+//
+// Consider using SetCoverEasy, when mask is downward closed.
+//
+// Verified by
+//   - SRM345 Div1 LV3
+//-------------------------------------------------------------
+
+LL next_combination(LL p)
+{
+	LL lsb = p & -p;
+	LL rem = p + lsb;
+	LL rit = rem & ~p;
+	return rem|(((rit/lsb)>>1)-1);
+}
+
+// Is it possible to cover the goal by at most k elements from mask?
+// Assumption:  goal \subseteq \bigcup mask
+bool canCover( LL goal, set<LL>& mask, int k )
+{
+	// optimizer
+	for(bool update=true; update; ) {
+		update = false;
+
+		// if *it \subseteq *jt for some jt, then we NEVER use it
+		// if *it is the only mask covering some city, we ALWAYS use it
+		for(set<LL>::iterator it=mask.begin(); it!=mask.end(); ) {
+			bool isSubset = false;
+			LL   onlyByMe = *it & goal;
+			for(set<LL>::iterator jt=mask.begin(); jt!=mask.end(); ++jt)
+				if( it!=jt ) {
+					onlyByMe &= ~*jt;
+					isSubset |= (*it & *jt & goal) == (*it & goal);
+				}
+
+			update |= isSubset | !!onlyByMe;
+
+			if( isSubset )
+				mask.erase(it++);
+			else if( onlyByMe ) {
+				if( --k < 0 )
+					return false;
+				goal &= ~*it;
+				mask.erase(it++);
+			}
+			else
+				++it;
+		}
+
+		if( mask.size()<=k || goal==0 )
+			return true;
+	}
+
+	// exhaustive search
+	vector<LL> ms(mask.begin(), mask.end());
+	for(LL i=(1LL<<k)-1; i<(1LL<<ms.size()); i=next_combination(i))
+	{
+		LL gg = goal;
+		for(LL j=0; (1LL<<j)<=i; ++j)
+			if( i & (1<<j) )
+				gg &= ~ms[j];
+		if( gg == 0 )
+			return true;
+	}
+	return false;
+}

ADDED   lib/typical/setCoverEasy.cpp
Index: lib/typical/setCoverEasy.cpp
==================================================================
--- lib/typical/setCoverEasy.cpp
+++ lib/typical/setCoverEasy.cpp
@@ -0,0 +1,74 @@
+
+//-------------------------------------------------------------
+// (Minimum) Set Cover
+//   NP-complete
+//
+// Assumption: "mask" must be downward closerd,
+//   i.e., X subset Y & Y in mask => X in mask
+//
+// Verified by
+//   - Google Code Jam 09 Round2 C
+//     (The "Only" part never tested, in the test cases of
+//      the problem, that optimization had never used.
+//      Be careful to use that; to be on the safer side,
+//      consider using the <true,false> mode.)
+//-------------------------------------------------------------
+
+int minCover_exhaustive( int goal, const set<int>& mask )
+{
+	typedef int STATE;
+
+	vector<STATE> Q( 1, 0 );
+	set<STATE> visited;
+	for(int step=0; step<=mask.size(); ++step)
+	{
+		vector<STATE> Q2;
+		for(int i=0; i<Q.size(); ++i)
+		{
+			STATE s = Q[i];
+			if( s == goal )
+				return step;
+			for(set<int>::const_iterator it=mask.begin(); it!=mask.end(); ++it)
+				if( !visited.count( s|*it ) )
+				{
+					visited.insert( s|*it );
+					Q2.push_back( s|*it );
+				}
+		}
+		Q.swap(Q2);
+	}
+	return -1;
+}
+
+template<bool Subs, bool Only>
+int minCover( int goal, set<int> mask )
+{
+	if( Subs )
+		for(set<int>::iterator it=mask.begin(); it!=mask.end(); )
+		{
+			bool needed = true;
+			for(int m=1; m<=goal; m<<=1)
+				if( (goal&m) && !(*it&m) && mask.count(*it|m))
+					{needed = false; break;}
+
+			if( needed )
+				++it;
+			else
+				mask.erase(it++);
+		}
+
+	if( Only )
+		for(set<int>::iterator it=mask.begin(); it!=mask.end(); )
+		{
+			int onlyByMe = goal & *it;
+			for(set<int>::iterator jt=mask.begin(); jt!=mask.end(); ++jt)
+				onlyByMe &= ~*jt;
+	
+			if( !onlyByMe )
+				++it;
+			else
+				mask.erase(it++), goal &= ~onlyByMe;
+		}
+
+	return minCover_exhaustive(goal, mask);
+}

ADDED   lib/typical/ternery.cpp
Index: lib/typical/ternery.cpp
==================================================================
--- lib/typical/ternery.cpp
+++ lib/typical/ternery.cpp
@@ -0,0 +1,18 @@
+// ternery search (for shita-ni-totsu)
+{
+	double x = ~min~;
+	double w = ~max~;
+	for(int i=0; i<~iteration~; ++i) // or, while( w-x > ~eps~ ) be careful for precision!
+	{
+		double y = 2*x/3 + w/3;
+		double z = x/3 + 2*w/3;
+		double fx = f(x);
+		double fy = f(y);
+		double fz = f(z);
+		double fw = f(w);
+		if( fx < fy )      w = y;
+		else if( fz > fw ) x = z;
+		else if( fy < fz ) w = z;
+		else               x = y;
+	}
+}

ADDED   lib/typical/unionfind.cpp
Index: lib/typical/unionfind.cpp
==================================================================
--- lib/typical/unionfind.cpp
+++ lib/typical/unionfind.cpp
@@ -0,0 +1,33 @@
+//-------------------------------------------------------------
+// Union-Find
+//   Balancing + Path-Compression : O(invAckermann)
+//
+// Verified by
+//   - SRM Member Pilot 2 Div1 LV2
+//-------------------------------------------------------------
+
+struct UnionFind
+{
+	vector<int> uf, sz;
+	int nc;
+
+	UnionFind(int N): uf(N), sz(N,1), nc(N)
+		{ for(int i=0; i<N; ++i) uf[i] = i; }
+	int size()
+		{ return nc; }
+	int Find(int a)
+		{ return uf[a]==a ? a : uf[a]=Find(uf[a]); }
+	bool Union(int a, int b)
+		{
+			a = Find(a);
+			b = Find(b);
+			if( a != b )
+			{
+				if( sz[a] >= sz[b] ) swap(a, b);
+				uf[a] = b;
+				sz[b] += sz[a];
+				--nc;
+			}
+			return (a!=b);
+		}
+};