Algorithm And Tips For Competitive Programming

Tomoki Imai

2013

1 テンプレート

各種バッドノウハウを含む。

1.1 C++

// compile in C++11. use -std=c++11.
#include <iostream>
#include <iomanip>
#include <vector>
#include <valarray>
#include <map>
#include <set>
#include <list>
#include <queue>
#include <stack>
#include <bitset>
#include <utility>
#include <numeric>
#include <algorithm>
#include <functional>
#include <complex>
#include <string>
#include <sstream>

#include <cstdio>
#include <cstdlib>
#include <cctype>
#include <cstring>

// this require C++11
#include <unordered_set>
#include <unordered_map>
#include <random>

using namespace std;

#define all(c) c.begin(),c.end()
#define repeat(i,n) for(int i=0;i<static_cast<int>(n);i++)
#define debug(x) #x << "=" << (x)
#define dump(x) cerr << debug(x) << " (L:" << __LINE__ << ")"<< endl

typedef long long ll;
typedef vector<int> vi;

template<typename T>
ostream& operator<<(ostream& os,const vector<T>& vec){
    os << "[";
    for(const auto& v : vec){
        os << v << ",";
    }
    os << "]";
    return os;
}

template<typename T>
T input(){
    T t;cin >> t;
    return t;
}

template<typename T>
vector<T> input(const int N){
    vector<T> v(N);
    repeat(i,N) cin >> v[i];
    return v;
}



int main(){
    return 0;
}

1.2 Python

#!/usr/bin/env python
#coding:utf-8

from __future__ import division
from __future__ import print_function
from __future__ import unicode_literals

from math import *
from string import *
from fractions import *
from itertools import *

def main():
    pass

if __name__ == "__main__":
    main()

2 算術型

2.1 int

基本中の基本。\(10^{9}\)くらい。こわいときにはlong longを使うことを推奨。

2.2 long long

大きい整数。\(10^{18}\)?くらい。

2.2.1 ビット演算

ビットDPとかに使う。

// 下から2番目のビットを立てる。1llと書くことに注意
x |= (1ll << 1);

// 集合演算
// 和
ll z = x | y;
// 積
ll z = x & y;
// xor
ll z = x ^ y;
// ループ
for(int i=0;i<N;i++){
    // i番目が立っていれば
    if(x & (1ll << i)){
        //something.
    }
}

ll dp[1 << N][N];
dp[0][0] = 1;
// Bitのほうが先!
for(int i=0;i<(1<<N);i++){
    for(int j=0;j<N;j++){
        //something.
    }
}

2.3 double

floatは使ってはだめ。

2.4 char

-128 ~ 127くらい。ちいさい。基本的には文字を入れるのに使う。vector<char>を
vector<bool>の代わりに使ってもいい。

char c = 'a';
//cctypeが必要。大文字に変換する。すでに大文字のときは何も起こらない。
c = toupper(c);
//小文字に変換する。
c = tolower(c);
//vector<bool>の代わり。
vector<char> used(10,false);

2.5 bool

true(==1)とかfalse(==0)を入れるためだけに使う。ただしvector<bool>は使ってはいけない。

2.6 補助関数

上記の型に関する便利な関数。

床関数、天井関数、および四捨五入。返り値はdouble。cmathをincludeする。

double x = 0.3;
int f = floor(x); // -> 0
int c = ceil(x); // -> 1
int r = round(x) // -> 0

任意の場所で四捨五入したいときには、\(10^{n}\)をかけて、roundした後に、
\(10^{n}\)で割ればいい。

double x = 0.123456789;
//0.123
double r = round(x*1000) / 1000.0;

2.6.0.1


template<typename T>
int popcount(T x){
    int ret = 0;
    while(x){
        x &= x-1;
        ret++;
   }
    return ret;
}

3 入出力

基本はcin,coutを使おう。

3.1 cin,cout

iostream,iomanipをincludeしておくこと。

3.1.1 cin

基本的な使い方について。

int n;
cin >> n;
vector<int> V(n);
rep(i,n) cin >> V[i];

こうすると、短く書ける。

入力の最後まで読む。

int n;
while(cin >> n){
    //処理
}

n=0のとき終わりとかの場合は、条件に&&n!=0とかをつける。
数値をカンマ区切りで読み込む。

int x,y;char c;
//cにカンマが入る
cin >> x >> c >> y;

冗長かもだけど、一番楽。
空白とか含めて一行読む。

string s;
getline(cin,s);

改行文字は、sに入らず、かつ読み捨てされる。
cinでは、改行文字は読み捨てないことに注意しよう。
つまり、数値<改行>文字列<改行>を読みたいときには、

int n;string s;
// 数値
cin >> n;
// 改行よみとばし
cin.ignore();
// 文字列
getline(cin,s);

とする。cinは改行文字を残すので、ignoreでそれを読み捨てないといけない。
また、ignoreの引数は読み捨てる文字数。引数なしの場合は1を渡したのと同等の効果がある。

3.1.2 cout

有効数字等が設定されている問題は、必ず多めに出力すること。多めに出す分には大丈夫。

基本の使い方

int n;vector<int> V(n);
cout << n << endl;
rep(i,n) cout << V[i] << " ";
cout << endl;

以下主なiomanipの使い方

int n = 123;double d = 1.23;

//10進数 −> 123
cout << dec << n << endl ;

//8進数 −> 173
cout << oct << n << endl ;

//16進数 −> 7b
cout << hex << n << endl ;

//16進数かつ、大文字 −> 7B
cout << hex << uppercase << n << endl;

//10進数に戻す
cout << dec;

//幅が10になるようにする。デフォルトは右寄せ
// -> xxxxxxx123 (default)
cout << setfill('x') << setw(10) << right << n << endl;

// -> 123xxxxxxx
cout << setfill('x') << setw(10) << left << n << endl;
// -> 123yyyyyyy
cout << setfill('y') << setw(10) << n << endl;

//小数点以下10桁表示に。
cout << fixed << setprecision(10);

// -> 1.2300000000
cout << d << endl;
// -> 12.3000000000
cout << 10*d << endl;

//小数点の表示を元に戻す
cout.unsetf(ios::fixed);
// -> 1.23
cout << d << endl;

基本的には、引数のあるマニピュレータの効果は保存される。

3.2 scanf,printf

C++では、cstdioをinclude。複雑な書式とかが必要なときにはこっちを使うといいかもしれない。

3.2.1 scanf

基本的な使い方

int n;char tmp[256];
scanf("%d\n",&n);
gets(tmp);

stringに直接いれるのはだめ。scanfはcinと同様に改行を残す。getlineするなら
cin.ignore。getsするなら、直前のscanfで改行を読んでおく必要がある。
また、scanfで改行を読むのでなく、直後にgetc(stdin)してもいい。

3.2.2 printf

基本的な使い方

int n = 100;
printf("n is %d\n",n);

scanfとほとんど同様の使い方ができる。

書式指定子
指定子 出力書式
%c 文字
%s 文字列
%d 符号付き10進整数
%u 符号なし10進
%f 10進浮動小数点数
%o 符号なし8進
%x 符号なし16進(Xなら大文字)
%% %記号

3.3 高速化

以下のコードをmain関数の最初に書くことで、cin,coutの速度が2倍程度になる。

ios::sync_with_stdio(false);
cin.tie(0);

ただし、このコードはstdioとの同期を切るという意味なので、これを使うときにはprintfやscanfを使用してはだめ。

4 std::vector

ここでは、配列のSTL版である、vectorの使いかたについて書く。
ここに書かれている関数は、string等にも用いることができるものが多い。
ちなみに、vector<bool>は使ってはいけない。bitsetや、vector<char>をつかうこと。
また、all(vector)は、vector.begin(),vector.end()とdefineしている。

4.1 基本

//要素数10で、初期値は-1にする。
vector<int> vec(10,-1);
//vecの最初から3つの要素をコピーする。
vector<int> newvec(vec.begin(),vec.begin()+3);
//vecの最初から3つの要素を削除する。
vec.erase(vec.begin(),vec.begin()+3);

4.2 並び換え

4.2.1 sort

C++のsortは、\(O(n \log n)\)で、introsort。何も指定しない場合には昇順にソートされる。
注意すべきなのは、C++11では、最悪ケースで\(O(n \log n)\)となること。C++03では特に何も制限はないが、g++ならば\(O(n \log n)\)である。

//昇順 (sort(vec.begin(),vec.end())) (2,1,3) -> (1,2,3)
sort(all(vec));
//降順 (ただ単純にreverseしてもいい) (2,1,3) -> (3,2,1)
sort(all(vec),greater<int>());

第三引数には関数を渡すことができる。

bool comp(const int& a ,const int& b){
    return abs(a) < abs(b);
}
int main(){
    vector<int> vec(10);
    //絶対値が小さい順にソート
    sort(all(vec),comp);
    sort(all(vec),[](int l,int r) {return abs(l)<abs(r)});
}

4.2.2 stable_sort

sortとちがって、同じ優先順位の要素の順番は保存される。最悪計算量は\(O(n \log ^ 2 n)\)である。

stable_sort(all(vec),comp);

4.3 unique

隣あう同じ要素を一つにまとめる。eraseすることに注意。

int ints[] = {1,1,2,1,1};
vector<int> vec(ints,ints+5);
vec.erase(unique(all(vec)),vec.end());
// 1 2 1
repeat(i,vec.size()) cout << vec[i] << endl;

4.4 rotate

rotateは第二引数の場所を先頭にするように回転する。

int[] ints = {1,2,3,4,5};
vector<int> vec(ints,ints+5);
rotate(vec.begin(),vec.begin()+1),vec.end()); //2,3,4,5,1
rotate(vec.begin(),vec.end()-1,vec.end()); //5,1,2,3,4

4.5 next_permutation

順列をすべて列挙する。\(N!\)個なので、それなりの勢いで大きくなる。章末の付録参照。

vector<int> V = {3,2,1,4};
//ソートすること。
sort(all(V));

do{
   for(int i=0;i<V.size();i++){
       cout << V[i] << " ";
   }
   cout << endl;
}while(next_permutation(all(V)));

4.6 座標圧縮

#include <iostream>
#include <vector>
#include <iomanip>
#include <algorithm>
#define rep(x,n) for(int x=0;x<(n);x++)
#define all(vec) vec.begin(),vec.end()


using namespace std;
template<typename T>
ostream& operator<<(ostream& os,const vector<T>& vec){
    os << "[";
    for(const auto& v : vec){
        os << v << ",";
    }
    os << "]";
    return os;
}

// 重なる長方形がカバーしている面積
int main(){
    int k = 0;
    while(true){
        k++;
        int n; cin >> n;
        if(n == 0) break;
        // (X1,Y1) is upper left.
        // (X1,Y1) is lower right.
        vector<double> X1(n),X2(n),Y1(n),Y2(n);
        rep(i,n){
            double x,y,r;
            cin >> x >> y >> r;
            X1[i] = x-r;
            X2[i] = x+r;
            Y1[i] = y-r;
            Y2[i] = y+r;
        }
        // X,Yは全てのX座標,Y座標
        vector<double> X,Y;
        X.insert(X.end(),all(X1));
        X.insert(X.end(),all(X2));
        Y.insert(Y.end(),all(Y1));
        Y.insert(Y.end(),all(Y2));

        // sortして重複を削除
        sort(all(X));
        sort(all(Y));
        X.erase(unique(all(X)),X.end());
        Y.erase(unique(all(Y)),Y.end());

        // X1,X2,Y1,Y2の値をX,Y中のindexにする。この時点でintにすべきかも。
        rep(i,n){
            X1[i] = lower_bound(all(X),X1[i]) - X.begin();
            X2[i] = lower_bound(all(X),X2[i]) - X.begin();
            Y1[i] = lower_bound(all(Y),Y1[i]) - Y.begin();
            Y2[i] = lower_bound(all(Y),Y2[i]) - Y.begin();
        }
        vector<vector<char>> is_on(X.size(),vector<char>(Y.size()));
        for(int i=0;i<n;i++){
            for(int sx=X1[i];sx<X2[i];sx++){
                for(int sy=Y1[i];sy<Y2[i];sy++){
                    is_on[sx][sy] = true;
                }
            }
        }
        double ans = 0;
        for(int i=0;i<X.size()-1;i++){
            for(int j=0;j<Y.size()-1;j++){
                if(is_on[i][j]){
                    ans += (X[i+1]-X[i]) * (Y[j+1]-Y[j]);
                }
            }
        }
        cout << fixed << setprecision(2);
        cout << k << " " << ans << endl;
    }
}

5 探索

5.1 全探索

全部しらべる。

int linear_search(vector<int> V,int val){
    repeat(i,V.size()){
        if(V[i] == val) return i;
    }
    return -1;
}

5.2 二分探索

ある条件を満たす最小のものを探す。ただし単調増加な物にしかつかえない。
叙述関数をPとすると、

double lower = 0,upper = 1000000;
for(int i=0;i<200;i++){
    double m = (lower+upper) / 2;
    if(P(m)){
        upper = m;
    }else{
        lower = m;
    }
    cout << upper << endl;
}

とすると、upperに求めたい値がはいる。もしみつからなかった場合には、値は
変わらない。なので、lower,upperには極端な値を設定すること。200という回数
は、すこし多め。100で十分。対象がvectorの場合は以下のように書ける。

#include <iostream>
#include <vector>
#include <algorithm>

using namespace std;

int main(){
    vector<int> v = {1,2,3,1,2,3,1,2,3};
    // ソートする必要あり。
    //  i  0 1 2 3 4 5 6 7 8
    // -> [1,1,1,2,2,2,3,3,3]
    sort(v.begin(),v.end());
    // 2以上の数値が初めて現れる場所 -> 3
    int lower = lower_bound(v.begin(),v.end(),2) - v.begin();
    // 2より大きい数値が初めて表われる場所 -> 6
    int upper = upper_bound(v.begin(),v.end(),2) - v.begin();
    // 2の個数
    int number_of_2 = upper - lower;
}

5.3 三分探索

凸関数の極値を求める

// 凸関数の極大な点をもとめる
template<typename F,typename T>
T ternary_search(F f,T left,T right){
    for(int i=0;i<1000;i++){
        T l = (2*left + right) / 3;
        T r = (left + 2*right) / 3;
        if(f(l) < f(r)){
            left = l;
        }else{
            right = r;
        }
    }
    return (left+right)/2;
}

// 凹関数の極小な・を求める
template<typename F,typename T>
T ternary_search_concave(F f,T left,T right){
    return ternary_search([f](T x){return -f(x);},left,right);
}

5.3.1 合計の重さXを達成するような組み合わせの数え上げ

#include <vector>
using namespace std;

// if we have objects that each has w_{i} weight.
//  we calc how many combination that has total weight X.
//  O(2^(N/2) logN)?
template<typename T>
int howmany_combination_such_that_total_weight_is_X(vector<T> W,int X){
    // sort(all(W));
    int N = W.size();
    vector<vector<T>> parts(2);
    copy(W.begin(),W.begin()+N/2,back_inserter(parts[0]));
    copy(W.begin()+N/2,W.end(),back_inserter(parts[1]));

    int ret = 0;
    vector<vector<T>> sums(2,{0});
    for(size_t p=0;p<parts.size();p++){
        for(size_t i=0;i<parts[p].size();i++){
            vector<T> nex;
            for(size_t j=0;j<sums[p].size();j++){
                nex.push_back(sums[p][j]);
                nex.push_back(sums[p][j]+parts[p][i]);
            }
            sums[p].swap(nex);
        }
    }
    for(vector<T>& sum : sums){
        sort(sum.begin(),sum.end());
    }

    for(T s : sums[0]){
        ret += upper_bound(sums[1].begin(),sums[1].end(),X-s)
             - lower_bound(sums[1].begin(),sums[1].end(),X-s);
    }
    return ret;
}

5.4 ぐるぐる

いつかつかうかも。 http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=1189

// WH/2,WH/2 を中心としたぐるぐる。
// まんなかは1
vector<vector<int>> guruguru(int N){
    const int WH = sqrt(N) + 2;
    vector<vector<int>> field(WH,vector<int>(WH,N + 100));
    {
        int cx=WH/2,cy=WH/2;
        field[cy][cx] = 1;
        int c = 2;
        for(int i=0;;i++){
            int sx = cx + i+1;
            int sy = cy + i;
            for(int j=0;j<2*i+2;j++){
                field[sy-j][sx] = c;
                c++;
                if(c >= MAXN) goto end;
            }
            vector<int> dx = {-1,0,1};
            vector<int> dy = {0,1,0};
            // left
            int ny = sy-2*i-1;
            int nx = sx;
            for(int k=0;k<3;k++){
                for(int j=0;j<2*(i+1);j++){
                    ny = ny + dy[k];
                    nx = nx + dx[k];
                    field[ny][nx] = c;
                    c++;
                    if(c >= MAXN) goto end;
                }
            }
        }
end:;
    }
    return field;
}

6 文字列操作

stringをincludeする。cctypeもいるかも。

6.1 std::string

charをラップしたテンプレートクラス。 基本的な使い方について

6.1.1 部分列

          //0123456789
string str("abcdefghij");
// 5番目以降
str.substr(5);    // "fghij"
// 5番目から3つ
str.substr(5,3); // "fgh"
//全部小文字にする
transform(s.begin(),s.end(),s.begin(),::tolower);

substrは一つの引数でそこから先全部、二つの場合は第一引数の位置から、第二
引数の数だけ持ってくる。

6.1.2 検索

stringには、いくつかのfindが定義されている。線形検索なので、早い検索が必
要なときには後述するKMP法やBM法を用いること。

第二引数として探すための最初の位置を指定できる。

6.1.2.1 Boyer Moore法

6.1.2.2 KMP法

6.2 stringstream

cinやcoutのようなstreamをstringを元に作成したりする。基本的には、string
をフォーマットしたり、intやlongに、intやlongから変換するために使用する。

stringstream ss;
ss << 102;
string s;
ss >> s;

6.3 再帰下降構文解析

BNFを書いて、それにしたがっていく。左再帰除去を忘れずに。

#include <iostream>
#include <string>
#include <cctype>
#include <sstream>

using namespace std;

typedef string::const_iterator Cursor;
class ParseError{};

// <四則演算の式> ::= <乗算除算の式> (+ or -) <乗算除算の式> (+ or -) ...
// <乗算除算の式> ::= <括弧か数> (* or /) <括弧か数> (* or /) ...
// <括弧か数>     ::= '(' <四則演算の式> ')' or <数>
// <数>           ::= (0|1|2|3|4|5|6|7|8|9)+

int expression(Cursor&);
int term(Cursor&);
int factor(Cursor&);
int number(Cursor&);

// <四則演算の式> ::= <乗算除算の式> (+ or -) <乗算除算の式> (+ or -) ...
int expression(Cursor &c){
    int ret = term(c);
    while(*c == '+' or *c == '-'){
        if(*c == '+'){
            c++;
            ret += term(c);
        }else{
            c++;
            ret -= term(c);
        }
    }
    return ret;
}

// <乗算除算の式> ::= <括弧か数> (* or /) <括弧か数> (* or /) ...
int term(Cursor &c){
    int ret = factor(c);
    while(*c == '*' or *c == '/'){
        if(*c == '*'){
            c++;
            ret *= factor(c);
        }else{
            c++;
            ret /= factor(c);
        }
    }
    return ret;
}

// <括弧か数>     ::= '(' <四則演算の式> ')' or <数>
int factor(Cursor &c){
    if(*c == '('){
        c++;
        int ret = expression(c);
        c++; // ')'
        return ret;
    }else{
        return number(c);
    }
}

// <数>           ::= (0|1|2|3|4|5|6|7|8|9)+
int number(Cursor &c){
    stringstream ss;
    while(isdigit(*c)){
        ss << *c;
        c++;
    }
    int ret;
    ss >> ret;
    return ret;
}

int main(){
    int test_case;
    cin >> test_case;
    cin.ignore();
    for(int i=0;i<test_case;i++){
        string s;
        getline(cin,s);
        Cursor c = s.begin();
        cout << expression(c) << endl;
    }
    return 0;
}

6.4 LL1文法(<-> 再帰下降で書けるかどうか)

予測型構文解析表を埋めて重複をチュックする。テストしたい。
一般に、"X -> Xb,X -> a"なる生成規則は、"X->aX',X'->bX',X'->"とできる。
これは前者がab*なる形の文字列を導出するからである。

#include <iostream>
#include <vector>
#include <set>
#include <map>

using namespace std;

struct Rule{
    char from;
    vector<char> to;
    Rule(char from_,vector<char> to_)
        : from(from_),to(to_){
    }
};

bool operator==(const Rule& lhs, const Rule& rhs){
    if(lhs.from != rhs.from or lhs.to.size() != rhs.to.size()){
        return false;
    }
    for(size_t i=0;i<lhs.to.size();i++){
        if(lhs.to[i] != rhs.to[i]){
            return false;
        }
    }
    return true;
}

bool operator<(const Rule& lhs, const Rule& rhs){
    if(lhs.from != rhs.from){
        return lhs.from < rhs.from;
    }
    if(lhs.to.size() != rhs.to.size()){
        return lhs.to.size() < rhs.to.size();
    }
    for(size_t i=0;i<lhs.to.size();i++){
        if(lhs.to[i] != rhs.to[i]){
            return lhs.to[i] < rhs.to[i];
        }
    }
    return false;
}

bool operator>(const Rule& lhs, const Rule& rhs){
    return not(lhs == rhs) and not(lhs < rhs);
}

ostream& operator<<(ostream& os,Rule r){
    os << r.from << "->";
    for(char c : r.to){
        os << " " << c;
    }
    return os;
}

template<typename T>
ostream& operator<<(ostream& os,set<T> s){
    os << "{";
    for(T v : s){
        os << v << ",";
    }
    os << "}";
    return os;
}

// can be written in recursive descent parsing?
//  see tiger book(p.46).
bool can_be_written_in_recursive_descent_parsing(const set<char>& terminal,
                                                 const set<char>& non_terminal,
                                                 const vector<Rule>& rules){
    for(Rule r : rules){
        cerr << r << endl;
    }
    map<char,bool> nullable;
    map<char,set<char>> first,follow;
    for(char c : non_terminal){
        nullable[c] = false;
        first[c] = {};
        follow[c] = {};
    }
    for(char c : terminal){
        nullable[c] = false;
        first[c] = {c};
        follow = {};
    }

    {
        bool changed;
        do{
            changed = false;
            for(Rule r : rules){
                bool all_nullable = true;
                for(char t : r.to){
                    all_nullable = all_nullable and nullable[t];
                }
                if(all_nullable and not nullable[r.from]){
                    nullable[r.from] = true;
                    changed = true;
                }
                for(size_t i=0;i<r.to.size();i++){
                    bool nul = true;
                    for(size_t k=0;k<i;k++){
                        nul = nul and nullable[r.to[k]];
                    }
                    if(nul){
                        size_t bef = first[r.from].size();
                        first[r.from].insert(first[r.to[i]].begin(),
                                             first[r.to[i]].end());
                        changed = changed or (bef != first[r.from].size());
                    }
                }
                for(size_t i=0;i<r.to.size();i++){
                    // nullable[y[k]] for k <- [i+1,]
                    bool nul = true;
                    for(int k=i+1;k<(int)r.to.size();k++){
                        nul = nul and nullable[r.to[k]];
                    }
                    if(nul){
                        size_t bef = follow[r.to[i]].size();
                        follow[r.to[i]].insert(follow[r.from].begin(),
                                               follow[r.from].end());
                        changed = changed or (bef != follow[r.to[i]].size());
                    }
                }
                for(size_t i=0;i<r.to.size();i++){
                    for(size_t j=i+1;j<r.to.size();j++){
                        // nullable[y[k]] for k <- [i+1,j-1]
                        bool nul = true;
                        for(size_t k=i+1;k<j;k++){
                            nul = nul and nullable[r.to[k]];
                        }
                        if(nul){
                            size_t bef = follow[r.to[i]].size();
                            follow[r.to[i]].insert(first[r.to[j]].begin(),
                                                   first[r.to[j]].end());
                            changed = changed or (bef != follow[r.to[i]].size());
                        }
                    }
                }
            }
        }while(changed);
    }
    for(char c : non_terminal){
        cerr << c << endl
             << "nullable? " << nullable[c] << endl
             << "first " << first[c] << endl
             << "follow " << follow[c] << endl;
    }
    map<char,map<char,set<Rule>>> maps;
    for(Rule r : rules){
        // check r.to is nullable or not.
        //  and collect first.
        set<char> to_first;
        bool to_nullable = true;
        for(char t : r.to){
            to_first.insert(first[t].begin(),first[t].end());
            to_nullable = to_nullable and nullable[t];
            if(not to_nullable) break;
        }

        for(char tf : to_first){
            maps[r.from][tf].insert(r);
        }
        if(to_nullable){
            for(char t : follow[r.from]){
                maps[r.from][t].insert(r);
            }
        }
    }

    cerr << "Map" << endl;
    bool dup = false;
    for(const pair<char,map<char,set<Rule>>> m : maps){
        int from = m.first;
        for(const pair<char,set<Rule>> p : m.second){
            int to = p.first;
            cerr << char(from) << " " << char(to) << " " << p.second << endl;
            dup = dup or p.second.size() > 1;
        }
    }
    return not dup;
}


int main(){
    set<char> terminal;
    set<char> non_terminal;

    {
        // terminal;
        int n;cin >> n;
        for(int i=0;i<n;i++){
            char c;cin >> c;
            terminal.insert(c);
        }
    }

    {
        // non_terminal;
        int n;cin >> n;
        for(int i=0;i<n;i++){
            char c;cin >> c;
            non_terminal.insert(c);
        }
    }

    // number of rule.
    int n;
    cin >> n;
    vector<Rule> rules;
    for(int i=0;i<n;i++){
        int k;cin >> k;
        char from;cin >> from;
        vector<char> to(k);
        for(int j=0;j<k;j++){
            cin >> to[j];
        }
        rules.push_back(Rule(from,to));
    }
    cout << can_be_written_in_recursive_descent_parsing(terminal,non_terminal,rules) << endl;
    return 0;
}

              // 9                  // 9
// 3          // $ + - * / i n ( )  // $ + - * / i n ( )
// a c d      // 3                  // 6
// 3          // S E T              // S E T D U F
// X Y Z      // 10                 // 12
// 6          // 2 S E $            // 2 S E $
// 1 X Y      // 3 E E p T          // 2 E T D
// 1 X a      // 3 E E m T          // 3 D + T D
// 0 Y        // 1 E T              // 3 D - T D
// 1 Y c      // 3 T T * F          // 0 D
// 1 Z d      // 3 T T / F          // 2 T F U
// 3 Z X Y Z  // 1 T F              // 3 U * F U
//  -> 0      // 1 F i              // 3 U / F U
              // 1 F n              // 0 U
              // 3 F ( E )          // 1 F i
              //  -> 0              // 1 F n
                                    // 3 F ( E )
                                    //  -> 1

6.5 Z algorithm

文字列のすべての部分文字列に対し、それと全体の共通接頭文字列の長さのテーブルを作る。 O(N)

// z-algorithm (O(|S|)
// Return V where V[i] = length of the largest common prefix of S and S[i]
// Example:
//  S: aaabaaaab
//  V: 921034210
// Ref: https://ivanyu.me/blog/2013/10/15/z-algorithm/
vector<int> zAlgorithm(const std::string& s)
{
    std::vector<int> z(s.size());
    z[0] = s.size();

    // left-bound, right-bound of current z-box.
    int lt = 0, rt = 0;

    for (int k = 1; k < s.size(); k++) {
        if (k > rt) {
            // If k is outside the current z-box, then we do naive computation.
            int n = 0;
            while (n + k < s.size() && s[n] == s[n + k])
                n++;
            z[k] = n;
            if (n > 0) {
                lt = k;
                rt = k + n - 1;
            }
        } else {
            // If k is inside the current z-box, then consider two cases
            int p = k - lt; // pair index
            // rightPartLen is a length of k..rt 
            int rightPartLen = rt - k + 1;

            if (z[p] < rightPartLen) {
                z[k] = z[p];
            } else {
                // s[i - k] is a A[i] where A starts from k of S.
                int i = rt + 1;
                while (i < s.size() && s[i] == s[i - k])
                    i++;
                z[k] = i - k;
                lt = k;
                rt = i - 1;
            }
        }
    }
    return z;
}

7 整数論

7.1 最大公約数,最小公倍数

ユークリッドの互除法を使う。intをlong longに置換してもいい。 \(O(\log n)\)

7.1.1 最大公約数

int gcd(int a,int b){
    return b==0 ? a : gcd(b,a%b);
}

7.1.2 最小公倍数

int lcm(int a,int b){
    return a*b / gcd(a,b);
}

7.2 mod

long longに入らないような答えのときにmodが登場する。

7.2.1 modの計算式について

\[\begin{aligned} a \equiv c & \pmod m \\ b \equiv d & \pmod m\end{aligned}\]

の条件下では以下の式が成り立つ。

\[\begin{aligned} a+b \equiv c+d & \pmod m \\ a-b \equiv c-d & \pmod m \\ a \times b \equiv c \times d & \pmod m\end{aligned}\]

さらに、mが素数の場合、以下の関係が成り立つ。

\[\begin{aligned} a ^ m \equiv a \pmod m \\ a ^ {m-1} \equiv 1 \pmod m \\ a ^ {m-2} \equiv \frac{1}{a} \pmod m\end{aligned}\]

つまり、\(a\)で割ることと、\(a^{m-2}\)を掛けることは同じである。
これは、\(C(10000,5000) \pmod p\)といった式を計算する際、次の冪乗の演算と組みあわせて用いる。

7.2.2 冪乗のmod

いわゆるmod_pow。計算量は\(O(\log n)\)

template<typename T,typename F>
T apply_doubling(const T& x,
                 const long long& n,
                 const T& id_elem,
                 const F& binary_op){
    if(n == 0) return id_elem;
    T ret = apply_doubling(binary_op(x,x),n/2,id_elem,binary_op);
    if(n % 2 == 1) ret = binary_op(ret,x);
    return ret;
}

long long mod_pow(long long x,long long n,long long mod){
    return apply_doubling(x,n,1ll,[mod](ll a,ll b){return a*b%mod;});
}

ちなみにC++のpowを使うときに、引数が整数で、返り値も整数であることを期待
するときには、上記のpowを使うべき。なぜならC++のpowは
double,double->doubleな関数であるから。

7.3 素数

7.3.1 エラトステネスの篩

#include <iostream>
#include <vector>
using namespace std;

//上限より余裕を取ること。
vector<bool> sieve(const int M){
    vector<bool> isPrime(M);
    for(int i=2;i<M;i++) isPrime[i] = true;
    for(int i=2;i*i<M;i++){
        if(not isPrime[i]) continue;
        for(int j=i*i;j<M;j+=i){
            isPrime[j] = false;
        }
    }
    return isPrime;
}

// リストがほしいときはこっち
vector<int> sieve(const int M){
    vector<bool> isPrime(M);
    vector<int> primes;
    for(int i=2;i<M;i++) isPrime[i] = true;
    for(int i=2;i*i<M;i++){
        if(not isPrime[i]) continue;
        for(int j=i*i;j<M;j+=i){
            isPrime[j] = false;
        }
    }
    for(int i=2;i<M;i++){
        if(isPrime[i]) primes.push_back(i);
    }
    return primes;
}
from math import *

def sieve(N):
    primes = set()
    for i in range(2,N):
        primes.add(i)

    for i in range(2,ceil(sqrt(N))):
        if i in primes:
            for j in range(i*i,N,i):
                primes.discard(j)

    return primes

素数のリストが欲しかったら、適当に突っ込むこと。
実際には\(O(n \log \log n)\)だけれど、大体\(O(n)\)だと思っていい。

7.3.2 素因数分解

// primes is primes under sqrt(N).
vector<int> prime_decomposition(int N,const vector<int>& primes){
    vector<int> ret;
    for(int p : primes){
        while(N % p == 0){
            N /= p;
            ret.push_back(p);
        }
    }
    if(N != 1) ret.push_back(N);
    return ret;
}

7.4 コンビネーション

くみあわせ。

#include <iostream>
#include <vector>
#include <algorithm>

using namespace std;
typedef long long ll;

//いい感じのやつ (n=61まで大丈夫)
//  modにしても良い。O(N)
ll nCr(int n,int r){
    if(n < r) return 0;
    r = min(r,n-r);
    ll ret = 1;
    for(int i=0;i<r;i++){
        ret *= n-i;
        ret /= i+1;
    }
    return ret;
}


ll mod_pow(ll x,ll n,ll mod){
    if(n==0) return 1;
    ll ret = mod_pow(x * x % mod,n/2,mod);
    if(n%2 == 1)  ret = ret * x % mod;
    return ret;
}
ll nCr_mod(ll n,ll r,ll mod){
    if(n < r) return 0;
    r = min(r,n-r);
    ll ret = 1;
    for(int i=0;i<r;i++){
        ret = (ret * (n-i)) % mod;
        // ret /= i+1;
        ret = (ret * mod_pow(i+1,mod-2,mod)) % mod;
    }
    return ret;
}

// http://www.ravco.jp/cat/view.php?cat_id=6179
// 区別するn 個のものから,繰り返し用いることを許して,r 個取り出して作った組
// <-> 区別しないr 個のボールを,区別するn 個の箱に配る場合の総数
ll nHr(ll n,ll r){
    return nCr(r+n-1,r);
}

ll nHr_mod(ll n,ll r,ll mod){
    return nCr_mod(r+n-1,r,mod);
}

//パスカルの三角形 (n=66まで大丈夫)
//たくさん必要になるときはこっちのほうがいい。
ll combi_tri(int n,int r){
    int N = n+1;
    vector<vector<ll> > memo(N,vector<ll>(N,0));
    memo[0][0] = 1;
    for(int i=1;i<N;i++){
        memo[i][0] = memo[i-1][0];
        for(int j=1;j<N;j++){
            memo[i][j] = memo[i-1][j-1] + memo[i-1][j];
        }
    }
    return memo[n][r];
}

// testing...
template <typename BidirectionalIterator>
bool next_combination(BidirectionalIterator first1,
                      BidirectionalIterator last1,
                      BidirectionalIterator first2,
                      BidirectionalIterator last2){
    if((first1 == last1) or (first2 == last2)){
        return false;
    }
    BidirectionalIterator m1 = last1;
    BidirectionalIterator m2 = last2;
    --m2;
    while(--m1 != first1 && not (*m1 < *m2)){
    }
    bool result = (m1 == first1) and not (*first1 < *m2);
    if(not result){
        while(first2 != m2 and not(*m1 < *first2)){
            ++first2;
        }
        first1 = m1;
        iter_swap(first1,first2);
        ++first1;++first2;
    }
    if((first1 != last1) and (first2 != last2)){
        m1 = last1; m2 = first2;
        while((m1 != first1) and (m2 != last2)){
            iter_swap(--m1,m2);
            ++m2;
        }
        reverse(first1,m1);
        reverse(first1,last1);
        reverse(m2,last2);
        reverse(first2,last2);
    }
    return not result;
}

template<typename BidirectionalIterator>
bool next_combination(BidirectionalIterator first,
                      BidirectionalIterator middle,
                      BidirectionalIterator last){
    return next_combination(first,middle,middle,last);
}

template<class BidirectionalIterator>
inline bool prev_combination(BidirectionalIterator first,
                             BidirectionalIterator middle,
                             BidirectionalIterator last){
  return next_combination(middle,last,first,middle);
}

template<typename T>
ostream& operator<<(ostream& os,const vector<T>& vec){
    os << "[";
    for(const auto& v : vec){
        os << v << ",";
    }
    os << "]";
    return os;
}

typedef unsigned long long ull;

void print_bit(ull S, int n=64){
    for(int i=n-1; i>=0; i--){
        if(S>>i & 1) std::cout << 1;
        else std::cout << 0;
    }
    std::cout << std::endl;
}

void subset_combination(int n, int k){
    ull S = (1ULL << k) - 1ULL;
    ull E = ~((1ULL << n) - 1ULL);
    while(!(S & E)){
        print_bit(S, n);
        ull smallest = S & -S;
        ull ripple = S + smallest;
        ull nsmallest = ripple & -ripple;
        S = ripple | (((nsmallest / smallest) >> 1) - 1);
    }
}


int main(){
    // poyoから2つえらぶ
    vector<int> poyo = {1,2,3,4,5};
    int N = poyo.size(),R = 3;
    do{
        vector<int> c(poyo.begin(),poyo.begin()+R);
        cerr << c  << endl;
    }while(next_combination(poyo.begin(),poyo.begin()+R,poyo.end()));

    return 0;
}

7.5 確率的なアレ

テスト中。

// n< 341550071728321 , ok.
#include <iostream>
#include <vector>
#include <cmath>
using namespace std;

typedef long long ll;

// will not overflow even if mod mod is too big..
ll mod_mul(ll a,ll b,ll mod){
    if(b == 0) return 0;
    ll res = mod_mul((a+a)%mod,b/2,mod);
    if(b % 2 == 1) res = (res + a)%mod;
    return res;
}

// use mod_mul if mod is too big.
ll mod_pow(ll x,ll n,ll mod){
    if(n==0) return 1;
    ll ret = mod_pow(x * x % mod,n/2,mod);
    if(n%2== 1)  ret = ret * x % mod;
    return ret;
}

// return probably prime.
//  if n < 341550071728321 return correct answer.
bool miller_rabin(ll n){
    if(n == 2) return true;
    if(n % 2 == 0 or n <= 1) return false;
    ll s=0,d=n-1;
    static const ll a[] = {2,3,5,7,11,13,17,n+1};

    while(d % 2 == 0){
        s++;d/=2;
    }
    for(int i=0;a[i]<n;i++){
        ll x = mod_pow(a[i],d,n);
        if(x != 1){
            ll r;
            for(r=0;r<s;r++){
                if(x == n-1) break;
                x = mod_mul(x,x,n);
            }
            if(r == s) return false;
        }
    }
    return true;
}

ll gcd(ll a,ll b){
    return b==0 ? a : gcd(b,a%b);
}

inline ll random(ll x,ll c,ll m){
    return (mod_mul(x,x,m)+c)%m;
}

vector<bool> sieve(const int M){
    vector<bool> isPrime(M);
    for(int i=2;i<M;i++) isPrime[i] = true;
    for(int i=2;i*i < M;i++){
        if(not isPrime[i]) continue;
        for(int j=i*i;j<M;j+=i){
            isPrime[j] = false;
        }
    }
    return isPrime;
}

// primes must contatin 2.
vector<ll> pollard_rho(ll n,const vector<ll>& primes,bool precheck=true){
    if(n == 0 or n == 1) return vector<ll>(1,n);
    vector<ll> ret;
    static const int constants[] = {1,51,73,0};
    if(precheck){
        for(size_t i=0;i<primes.size();i++){
            while(n % primes[i] == 0){
                n /= primes[i];
                ret.push_back(primes[i]);
            }
            if(n == 1) return ret;
        }
    }
    if(miller_rabin(n)){
        ret.push_back(n);
        return ret;
    }
    for(int i=0;constants[i];i++){
        ll x = 2,y = 2,d = 1;
        while(d == 1){
            x = random(x,constants[i],n);
            y = random(random(y,constants[i],n),constants[i],n);
            d = gcd(abs(x-y),n);
        }
        if(d == n) continue;

        vector<ll> dp = pollard_rho(d,primes,false);
        vector<ll> ndp = pollard_rho(n/d,primes,false);
        for(size_t j=0;j<dp.size();j++){
            ret.push_back(dp[j]);
        }
        for(size_t j=0;j<ndp.size();j++){
            ret.push_back(ndp[j]);
        }
        return ret;
    }
    return ret;
}


int main(){
    int N = 100;
    vector<bool> isprime = sieve(5);
    vector<ll> primes;
    for(int i=0;i<isprime.size();i++){
        if(isprime[i]){
            primes.push_back(i);
        }
    }
    for(int i=0;i<N;i++){
        //cout << i << " " << miller_rabin(i) << " " << endl;
        vector<ll> ps = pollard_rho(i,primes,true);
        for(int j=0;j<ps.size();j++){
            cout << ps[j] << " ";
        }
        cout << endl;
    }
    return 0;
}

7.6 乱数

XORShiftをつかったらうれしいかもしれない。

unsigned long xor128(){
    static unsigned long x=123456789,y=362436069,z=521288629,w=88675123;
    unsigned long t;
    t=(x^(x<<11));x=y;y=z;z=w; return( w=(w^(w>>19))^(t^(t>>8)) );
}

7.7 基数変換

vector<int> int_to_digits(int base,int N){
    vector<int> rev_ret;
    while(N != 0){
        rev_ret.push_back(N % base);
        N /= base;
    }
    reverse(all(rev_ret));
    return rev_ret;
}

int digits_to_int(int base,const vector<int>& digits){
    int ret = 0;
    for(int i=0;i<(int)digits.size();i++){
        ret *= base;
        ret += digits[i];
    }
    return ret;
}

7.8 二次方程式

// solve ax^2 + bx + c = 0.
vector<double> solve_quadratic_equation(double a,double b,double c){
    vector<double> ret;
    if(abs(a-EPS) < 0){
        ret.push_back(-c/b);
    }else{
        double d = b*b-4*a*c;
        if(d >= 0){
            if(b >= 0){
                ret.push_back(-2*c/(b+sqrt(d)));
                ret.push_back((-b-sqrt(d))/(2*a));
            }else{
                ret.push_back(-2*c/(b-sqrt(d)));
                ret.push_back((-b+sqrt(d))/(2*a));
            }
            if(abs(ret[0]-ret[1]) < EPS){
                ret.pop_back();
            }
        }
    }
    return ret;
}

7.9 連立方程式

テスト中

const int none = 0; // if no answer
const int one = 1;  // if there is exactly one answer
const int many = 2; // many answer.
// answer will be inserted in b.
//  test中
int normal_gauss_jordan(vector<vector<double>> A,vector<double>& b){
    const double EPS = 1e-8;
    int n = A.size();
    int m = A[0].size();
    int pi = 0,pj = 0;
    while(pi < n and pj < m){
        for(int i=pi+1;i<n;i++){ // pivot
            if(abs(A[i][pj]) > abs(A[pi][pj])){
                swap(A[i],A[pi]);
                swap(b[i],b[pi]);
            }
        }
        if(abs(A[pi][pj]) > EPS){
            double d = A[pi][pj];
            for(int j=0;j<m;j++){
                A[pi][j] /= d;
            }
            b[pi] /= d;
            for(int i=pi+1;i<n;i++){
                double k = A[i][pj];
                for(int j=0;j<m;j++){
                    A[i][j] -= k * A[pi][j];
                }
                b[i] -= k*b[pi];
            }
            pi++;
        }
        pj++;
    }
    for(int i=pi;i<n;i++){
        if(abs(b[i]) > EPS){
            return none;
        }
    }
    if(pi < m or pj < m){
        return many;
    }
    for(int j=m-1;j>=0;j--){
        for(int i=0;i<j;i++){
            b[i] -= b[j] * A[i][j];
        }
    }
    return one;
}

7.10 nと互いに素な数の個数

EulerのTotient関数による。

// nと互いに素な数の個数
// n (1-1/p1) ... (1-1/pn)
template<typename T>
T totient(T n){
    if(n == 0) return 0;
    T ans = n;
    for(T x=2;x*x <= n;x++){
        if(n % x == 0){
            ans -= ans/x;
            while(n % x == 0) n /= x;
        }
    }
    if(n > 1){
        ans -= ans/n;
    }
    return ans;
}

8 行列

8.1 基本要素

正方行列用 //いつかなおす。

// 適宜intにしたりすること。
typedef vector<vector<ll> > ll_mat;

8.2 基本演算

かけ算とmod。たしざんはcoming soon.

namespace std{
    ll_mat operator*(const ll_mat& lhs,const ll_mat& rhs){
        int N = lhs.size();
        ll_mat ret(N,vector<ll>(N));
        for(int row=0;row<N;row++){
            for(int col=0;col<N;col++){
                for(int k=0;k<N;k++){
                    ret[row][col] += rhs[row][k] * lhs[k][col];
                }
            }
        }
        return ret;
    }

    ll_mat operator%(ll_mat lhs,ll rhs){
        int N = lhs.size();
        ll_mat ret = lhs;
        for(int i=0;i<N;i++){
            for(int j=0;j<N;j++){
                ret[i][j] = ret[i][j] % rhs;
            }
        }
        return ret;
    }
};

8.3 基本操作

8.3.1 累乗

vector<vector<ll> >  mod_pow(vector<vector<ll> > x,ll n,ll mod){
    if(n==0){
        vector<vector<ll> > E(x.size(),vector<ll>(x.size()));
        for(int i=0;i<x.size();i++){
            E[i][i] = 1;
        }
        return E;
    }
    vector<vector<ll> > ret = mod_pow(x * x % mod,n/2,mod);
    if(n%2 == 1)  ret = ret * x % mod;
    return ret;
}

8.3.2 表示

void display_matrix(vector<vector<ll> > mat){
    for(int i=0;i<mat.size();i++){
        for(int j=0;j<mat[0].size();j++){
            cerr << mat[i][j] << " ";
        }
        cerr << endl;
    }
}

8.3.3 ベクトルとのかけ算

一次元列ベクトルとのかけ算

vector<ll> mat_multi(vector<vector<ll> > lhs,vector<ll> rhs,int mod){
    vector<ll> ret(rhs.size());
    int N = lhs.size();
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            ret[i] = (ret[i] + lhs[i][j] * rhs[j]) % mod;;
        }
    }
    return ret;
}

8.4 Gauss-Jordan

GF上で連立方程式を解く

// Gauss-Jordan. Solve equation on GF.
int invert(int x){
    int ret[2] = {0,1};
    return ret[x];
}

int modulo(int x){
    x %= 2;
    while(x < 0){
        x += 2;
    }
    return x;
}

const int none = 0; // if no answer
const int one = 1;  // if there is exactly one answer
const int many = 2; // many answer.
// answer will be inserted in b.
int gauss(matrix A,vector<int>& b){
    int n = A.size();
    int m = A[0].size();
    int pi = 0,pj = 0;
    while(pi < n and pj < m){
        for(int i=pi+1;i<n;i++){ // pivot
            if(abs(A[i][pj]) > abs(A[pi][pj])){
                swap(A[i],A[pi]);
                swap(b[i],b[pi]);
            }
        }
        if(abs(A[pi][pj]) > 0){
            int d = invert(A[pi][pj]);
            for(int j=0;j<m;j++){
                A[pi][j] = modulo(A[pi][j] * d);
            }
            b[pi] = modulo(b[pi]*d);
            for(int i=pi+1;i<n;i++){
                int k = A[i][pj];
                for(int j=0;j<m;j++){
                    A[i][j] = modulo(A[i][j] - k * A[pi][j]);
                }
                b[i] = modulo(b[i] - k*b[pi]);
            }
            pi++;
        }
        pj++;
    }
    for(int i=pi;i<n;i++){
        if(abs(b[i]) > 0){
            return none;
        }
    }
    if(pi < m or pj < m){
        return many;
    }
    for(int j=m-1;j>=0;j--){
        for(int i=0;i<j;i++){
            b[i] = modulo(b[i] - b[j] * A[i][j]);
        }
    }
    return one;
}

9 動的計画法およびそれに似たやつら。(TODO)

9.1 LCS

Longest common sequence.O(NM)

string lcs(const string& a,const string& b){
    const int n = a.size(),m = b.size();
    vector<vector<int> > x(n+1,vector<int>(m+1));
    vector<vector<int> > y(n+1,vector<int>(m+1));

    for(int i=0;i<n;i++){
        for(int j=0;j<m;j++){
            if(a[i] == b[j]){
                x[i+1][j+1] = x[i][j]+1;
                y[i+1][j+1] = 0;
            }else if(x[i+1][j] < x[i][j+1]){
                x[i+1][j+1] = x[i][j+1];
                y[i+1][j+1] = 1;
            }else{
                x[i+1][j+1] = x[i+1][j];
                y[i+1][j+1] = -1;
            }
        }
    }
    string ret;
    for(int i=n,j=m;i>0 and j>0;){
        if(y[i][j] > 0){
            i--;
        }else if(y[i][j] < 0){
            j--;
        }else{
            ret.push_back(a[i-1]);
            i--;j--;
        }
    }
    reverse(all(ret));
    return ret;
}

9.2 LIS

9.2.1 O(NlogN)

// O(Nlog(N))
//  Sphagetthi source より.
vector<int> longest_increasing_subsequence(const vector<int>& a){
    const int n = a.size();
    vector<int> A(n,1<<30);
    vector<int> id(n);
    for(int i=0;i<n;i++){
        id[i] = distance(A.begin(),lower_bound(A.begin(),A.end(),a[i]));
        A[id[i]] = a[i];
    }
    int m = *max_element(id.begin(),id.end());
    vector<int> b(m+1);
    for(int i=n-1;i>=0;i--){
        if(id[i] == m) b[m--] = a[i];
    }
    return b;
}

9.3 巡回セールスマン問題

bit演算をする。bitのループを先に回すこと。\(O(N^2\times2^{N})\)

ll dp[1 << N][N];
dp[0][0] = 1;
// Bitのほうが先!
for(int i=0;i<(1<<N);i++){
    for(int j=0;j<N;j++){
        //something.
    }
}

9.4 ナップサック問題

個数制限なしのとき

#include <iostream>
#include <vector>

using namespace std;

#define rep(i,n) for(int i=0;i<(int)n;i++)

typedef vector<int> vi;

template<class T>
inline int len(const T& t){
    return t.size();
}

struct Treasure{
    int value,weight;
    Treasure():
        value(0),weight(0) {};
    Treasure(int _value,int _weight)
        : value(_value),weight(_weight) {}
};

int main(){
    int W,N;
    cin >> W;
    cin >> N;
    vector<Treasure> v(N);
    repeat(i,N){
        cin >> v[i].value >> v[i].weight;
    }
    // もしも複数入れて良いなら、直接valueを更新する。
    vi value(W+1);
    repeat(n,N){
        vi tmp(W+1);
        repeat(i,W+1){
            tmp[i] = max(tmp[i],value[i]);
            int in = i+v[n].weight;
            if(in <= W){
                tmp[in] = max(tmp[in],value[i]+v[n].value);
            }
        }
        value.swap(tmp);
    }

    int retw=0;
    int retv=0;
    repeat(i,len(value)){
        if(value[i] > retv){
            retv = value[i];
            retw = i;
        }
    }
    cout << retv << endl << retw << endl;
    return 0;
}

9.5 ヒストグラム内の最大長方形のサイズ

#define repeat(i,n) for(int i=0;i<static_cast<int>(n);i++)
typedef vector<int> vi;
typedef vector<vi> vvi;

// see http://algorithms.blog55.fc2.com/blog-entry-132.html
//  max rectangle area in histogram.
template<typename T>
T max_rectangle_size_in_histogram(const vector<T>& v){
    // where,height. in this function,use index as where.
    typedef pair<int,T> piT;
    stack<piT> sta;
    T ret = 0;
    for(int i=0;i<v.size();i++){
        if(sta.empty() or sta.top().second < v[i]){
            sta.push(make_pair(i,v[i]));
        }else if(sta.top().second == v[i]){
            // pass
        }else if(sta.top().second > v[i]){
            piT t;
            while(not sta.empty() and sta.top().second >= v[i]){
                t = sta.top();
                ret = max(ret,t.second*(i-t.first));
                sta.pop();
            }
            sta.push(make_pair(t.first,v[i]));
        }
    }
    int total_width = v.size();
    while(not sta.empty()){
        piT t = sta.top();
        sta.pop();
        ret = max(ret,t.second*(total_width-t.first));
    }
    return ret;
}

// square area 
int solve(const vvi &matrix){
    int h = matrix.size();
    int w = matrix[0].size();
    vvi up(h,vi(w));

    {
        for(int i=0;i<h;i++){
            for(int j=0;j<w;j++){
                if(i == 0){
                    up[i][j] = matrix[i][j] == 1;
                }else {
                    if(matrix[i][j] == 0){
                        up[i][j] = 0;
                    }else{
                        up[i][j] = up[i-1][j]+1;
                    }
                }
            }
        }
    }

    int ret = 0;
    for(int y=0;y<h;y++){
        ret = max(ret,max_rectangle_size_in_histogram(up[y]));
    }
    return ret;
}

9.6 なんか長方形の和をもとめるやつ

#include <vector>
#include <iostream>

using namespace std;

// ret[i][j] = sum of all field in [0,i) x [0,j).
template<typename T>
vector<vector<T>> calc_sum(const vector<vector<T>>& field){
    int h = field.size();
    int w = field[0].size();
    vector<vector<T>> ret(h+1,vector<T>(w+1));
    for(int y=1;y<=h;y++){
        for(int x=1;x<=w;x++){
            ret[y][x] = ret[y][x-1] + field[y-1][x-1];
        }
    }
    for(int x=0;x<=w;x++){
        for(int y=1;y<=h;y++){
            ret[y][x] = ret[y][x] + ret[y-1][x];
        }
    }
    return ret;
}

// calc [sy,gy] x [sx,gx]
template<typename T>
T calc(const vector<vector<T>>& v,int sy,int sx,int gy,int gx){
    // for close section
    gy++;gx++;
    return v[gy][gx] - v[sy][gx] - v[gy][sx] + v[sy][sx];
}

int main(){
    vector<vector<int>> field = {{0,1,2},
                                 {3,4,5},
                                 {6,7,8}};
    auto s = calc_sum(field);
    // calc 3+5+6 + 6+7+8
    cout << calc(s,1,0,2,2) << endl;

    return 0;
}

// 座標圧縮
int main(){
    int n,m;cin >> n >> m;

    vector<int> index_of_x(n),index_of_y(n);
    for(int i=0;i<n;i++){
        cin >> index_of_x[i] >> index_of_y[i];
    }
    vector<int> value_of_x = index_of_x,
                value_of_y = index_of_y;

    sort(all(value_of_x));
    sort(all(value_of_y));
    value_of_x.erase(unique(all(value_of_x)),value_of_x.end());
    value_of_y.erase(unique(all(value_of_y)),value_of_y.end());
    for(int i=0;i<n;i++){
        index_of_x[i] = lower_bound(all(value_of_x),index_of_x[i]) - value_of_x.begin();
        index_of_y[i] = lower_bound(all(value_of_y),index_of_y[i]) - value_of_y.begin();
    }
    vector<vector<int>> p(value_of_y.size(),vector<int>(value_of_x.size()));
    for(int i=0;i<n;i++){
        p[index_of_y[i]][index_of_x[i]] += 1;
    }
    auto v = calc_sum(p);

    for(int i=0;i<m;i++){
        int x1,y1,x2,y2;
        cin >> x1 >> y1 >> x2 >> y2;

        int sx = lower_bound(all(value_of_x),x1) - value_of_x.begin();
        int gx = upper_bound(all(value_of_x),x2) - value_of_x.begin();
        int sy = lower_bound(all(value_of_y),y1) - value_of_y.begin();
        int gy = upper_bound(all(value_of_y),y2) - value_of_y.begin();

        int ans = v[gy][gx] - v[sy][gx] - v[gy][sx] + v[sy][sx];
        cout << ans << endl;
    }

    return 0;
}

10 データ構造

10.1 Union-Find

struct UnionFind{
    const int is_root = -1;
    vector<int> data;
    UnionFind(int n){
        data = vector<int>(n,-1);
    }
    // 親を探す
    int root(int x){
        if(data[x] < 0){
            return x;
        }else{
            return data[x] = root(data[x]);
        }
    }
    // x,yの含まれる集合を併合
    bool unite(int x,int y){
        x = root(x);
        y = root(y);
        if(x != y){
            // 大きいほうに追加する。
            if(data[y] < data[x]) swap(x,y);
            data[x] += data[y];
            data[y] = x;
            return true;
        }else{
            return false;
        }
    }
    // 同じ集合にいるかどうか
    bool same(int x,int y){
        return root(x) == root(y);
    }
    int size(int x){
        return -data[root(x)];
    }
};

10.2 ヒープ

#include <queue>
#include <iostream>

using namespace std;

typedef pair<int,int> pii;

struct Comp{
    bool operator()(pii left,pii right){
        if(left.second < right.second) return true;
        else if(left.second == right.second and left.first > right.first) return true;
        else return false;
    };
};

struct Robot{
    int y,x,dir,step;
    Robot(int y,int x,int dir,int step) : y(y),x(x),dir(dir),step(step) {};
};

// <,>を定義すればless<Robot>みたいに扱える。
bool operator<(const Robot& lhs,const Robot& rhs){
    return lhs.step < rhs.step;
}
bool operator>(const Robot& lhs,const Robot& rhs){
    return lhs.step > rhs.step;
}

int main(){
    // 何も書かないと降順。(おっきい方からでてくる。)
    // これは昇順(ちいさいほうから出てくる)にしたもの。
    priority_queue<int,vector<int>,greater<int> > Qi;
    //関数オブジェクトを使っていい感じにもできる。
    priority_queue<pii,vector<pii>,Comp> Q;
    // 自作クラスの場合はこんな感じ
    priority_queue<Robot,vector<Robot>,greater<Robot> > que;

    Q.push(make_pair(1,2));
    Q.push(make_pair(2,2));
    Q.push(make_pair(3,2));
    while(not Q.empty()){
        cout << Q.top().first << " " << Q.top().second << endl;
        Q.pop();
    }
}

10.3 bitset

限られた大きさのvector<bool>を使いたいときに、bitsetを使うことができる。

#include <iostream>
#include <bitset>

using namespace std;
const int N = 10;

struct bit_cmp{
    bool operator() (const bitset<N> &left,const bitset<N> &right) {
        for(int i=N-1;i>=0;i--){
            if(left[i] < right[i]) return true;
            if(left[i] > right[i]) return false;
        }
        return false;
    }
};


int main(){
    //定数じゃないとダメ。最初は全部false
    bitset<N> bits;
    // すべての要素をtrue -> 1111111111
    bits.set();
    if(bits.all()) cout << "all" << endl;
    // 立っているbitの数 -> 10
    cout << bits.count() << endl;
    // すべての要素をfalse -> 0000000000
    bits.reset();
    if(bits.none()) cout << "none" << endl;

    //1番目の要素をtrue -> 0100000000
    bits.set(1);
    if(bits.any()) cout << "any" << endl;

    // 0110000000
    bits.set(2);
    //1番目の要素をfalse -> 0010000000
    bits.reset(1);

    if(bits[2]) cout << 2 << endl;
    cout << bits << endl;

    bitset<N> newbits;
    // 和を取る
    bits |= newbits;
    // 積を取る
    bits &= newbits;

    // 関数オブジェクトを作る必要アリ
    map<bitset<N>,int,bit_cmp> M;
}

10.4 分数

テストちゅう.

typedef long long ll;
ll gcd(ll a,ll b){
    return b==0 ? a : gcd(b,a%b);
}

ll lcm(ll a,ll b){
    if(a < 0) a *= -1;
    if(b < 0) b *= -1;
    return a*b / gcd(a,b);
}

struct Fraction{
    ll n,d;
    Fraction(ll _n,ll _d){
        ll c = lcm(_n,_d);
        n = c / _d;
        d = c / _n;
        if(d < 0){
            n *= -1;
            d *= -1;
        }
    }

    Fraction(ll _n){
        Fraction(_n,1);
    }

    bool operator<(const Fraction& r) const{
        ll c_d = lcm(d,r.d);
        return n*(c_d/d)< r.n*(c_d/r.d);
    }
    bool operator>(const Fraction& r) const{
        return not ((*this) < r or (*this) == r);
    }
    bool operator<=(const Fraction& r) const{
        return (*this) < r or (*this) == r;
    }
    bool operator>=(const Fraction& r) const{
        return (*this) > r or (*this) == r;
    }
    bool operator==(const Fraction& r) const{
        return n == r.n and d == r.d;
    }
    Fraction operator+(const Fraction& r) const{
        ll c_d = lcm(d,r.d);
        return Fraction(n*(c_d/d)+r.n*(c_d/r.d),c_d);
    }
    Fraction operator-(const Fraction& r) const{
        return (*this) + (-r);
    }
    Fraction operator*(const Fraction& r) const{
        return Fraction(n*r.n,d*r.d);
    }
    Fraction operator/(const Fraction& r) const{
        return (*this) * Fraction(r.d,r.n);
    }
    Fraction operator+() const{
        return Fraction(n,d);
    }
    Fraction operator-() const{
        return (*this) * -1;
    }
    Fraction operator*(const ll& a) const{
        return Fraction(a*n,d) ;
    }
};

ostream& operator<<(ostream &os,const Fraction& f){
    os << f.n << "/" << f.d;
    return os;
}

10.5 セグメント木

10.5.1 普通の

RMQ (POJ 3264)

#include <vector>
using namespace std;

struct SegmentTree{
    // edit this things! Range Sum Query
    static const int INIT_VALUE = 0;
    inline int calc(int a,int b){
        return a+b;
    }

    // Range Minimum Query
    // static const int INIT_VALUE = 1 << 30;
    // inline int calc(int a,int b){
    //     return min(a,b);
    // }


    vector<int> data;
    int n;
    explicit SegmentTree(int _n){
        n = 1;
        while(n < _n) n *= 2;
        data = vector<int>(2*n-1,INIT_VALUE);
    }

    void update(int i,int v){
        int k = i + n-1;
        data[k] = v;
        while(k > 0){
            k = parent(k);
            data[k] = calc(data[get_left(k)],data[get_right(k)]);
        }
    }

    // return data of [a,b).
    inline int get(int a,int b){
        return get(a,b,0,0,n);
    }
private:
    // l,r is segment of node k.
    int get(int a,int b,int k,int l,int r){
        if(r <= a or b <= l) return INIT_VALUE;
        if(a <= l and r <= b) return data[k];
        int vl = get(a,b,get_left(k),l,(l+r)/2);
        int vr = get(a,b,get_right(k),(l+r)/2,r);
        return calc(vl,vr);
    }

    inline int get_left(int x){
        return 2*x+1;
    }
    inline int get_right(int x){
        return 2*x+2;
    }
    inline int parent(int x){
        return (x-1)/2;
    }
};

10.5.2 Lazy

テスト中 参考: http://d.hatena.ne.jp/kyuridenamida/20121114/1352835261

typedef pair<int,int> pii;
// update [l,r) v -> +v to all element in [l,r).
// get [l,r) -> return sum of elements in [l,r).
struct LazySegmentTree{
    // also used for out of range value.
    static const int DATA_INIT_VALUE = 0;
    static const int LAZY_INIT_VALUE = 0;
    vector<int> data;
    vector<int> lazy;
    int n;
    LazySegmentTree(int _n){
        n = 1;
        while(n < _n) n *= 2;
        data = vector<int>(2*n-1,DATA_INIT_VALUE);
        lazy = vector<int>(2*n-1,LAZY_INIT_VALUE);
    }

    int get(int a,int b){
        return get(a,b,0,0,n);
    }
    void update(int a,int b,int v){
        return update(a,b,v,0,0,n);
    }


private:
    // edit update_lazy and update_data.
    inline void update_lazy(int k,int l,int r){
        data[k] += lazy[k]*(r-l);
        if(k < n-1){ // node k has children
            lazy[get_left(k)] += lazy[k];
            lazy[get_right(k)] += lazy[k];
        }
        lazy[k] = 0;
        return;
    }
    inline void update_data(int k){
        data[k] = data[get_left(k)] + data[get_right(k)];
    }

    // l,r is segment of node k.
    int get(int a,int b,int k,int l,int r){
        update_lazy(k,l,r);
        // [a,b) and [l,r) are not crossed
        if(r <= a or b <= l) return DATA_INIT_VALUE;
        // [a,b) contains [l,r)
        if(a <= l and r <= b) return data[k];
        int vl = get(a,b,get_left(k),l,(l+r)/2);
        int vr = get(a,b,get_right(k),(l+r)/2,r);
        update_data(k);
        return vl+vr;
    }
    void update(int a,int b,int v,int k,int l,int r){
        update_lazy(k,l,r);
        // [a,b) and [l,r) are not crossed
        if(r <= a or b <= l) return;
        // [a,b) contains [l,r)
        if(a <= l and r <= b) {
            lazy[k] += v;
            update_lazy(k,l,r);
            return;
        }
        update(a,b,v,get_left(k),l,(l+r)/2);
        update(a,b,v,get_right(k),(l+r)/2,r);
        update_data(k);
        return;
    }

    inline int get_left(int x){
        return 2*x+1;
    }
    inline int get_right(int x){
        return 2*x+2;
    }
    // parent node is (n-1)/2.
    inline int parent(int x){
        return (x-1)/2;
    }
};

ostream& operator<<(ostream& os,const LazySegmentTree& seg){
    os << seg.data << endl << seg.lazy;
    return os;
}

10.6 定数個のみを保持するpriority_queue

ビームサーチとかに使える?

#include <vector>
#include <queue>
#include <iostream>


template<class T,
         class Container=std::vector<T>,
         class Compare=std::less<typename Container::value_type>>
struct lens_queue : public std::priority_queue<T,Container,Compare>{
    typedef std::priority_queue<T,Container,Compare> super;
    lens_queue(int _cap)
        : cap(_cap){
    }
    inline std::size_t capacity() const {return cap;}
    void push(const T& x){
        if(should_push(x)){
            super::push(x);erase_overflow();
        }
    }
    void push(T&& x){
        if(should_push(x)){
            super::push(std::move(x));erase_overflow();
        }
    }
    template <class... Args>
    void emplace(Args&&... args){
        this->c.emplace_back(std::forward<Args>(args)...);
        if(should_push(this->c.back())){
            std::push_heap(this->c.begin(),this->c.end(),this->comp);
            erase_overflow();
        }else{
            this->c.pop_back();
        }
    }
    // this crash lens_queue.
    Container&& take(){
        return std::move(this->c);
    }
    template<class T2,class Container2,class Compare2>
    friend std::ostream& operator<<(std::ostream&,
                                    const lens_queue<T2,Container2,Compare2>&);
private:
    std::size_t cap;
    void erase_overflow(){
        while(this->size() > capacity()){
            this->pop();
        }
    }
    inline bool should_push(const T& x){
        return this->size()<capacity() || !this->comp(this->top(),x);
    }
};

template<class T,class Container,class Compare>
std::ostream& operator<<(std::ostream& os,const lens_queue<T,Container,Compare>& que){
    os << "{";
    for(const auto& v : que.c){
        os << v << ",";
    }
    os << "}";
    return os;
}



signed main() {
    int n,k;char c;
    std::cin >> n >> k >> c;
    if(c == 'g'){
        // k = 3
        // 4 2 3 1 5 -> 3 4 5
        lens_queue<int,std::vector<int>,std::greater<int>> que(k);
        for(int i=0;i<n;i++){
            int v;std::cin >> v;
            que.push(v);
            // que.emplace(v);
        }
        vector<int> taken = que.take();
        // while(not que.empty()){
        //     std::cout << que.top() << std::endl;
        //     que.pop();
        // }
    }else if(c == 'l'){
        // k = 3
        // 4 2 3 1 5 -> 3 2 1
        lens_queue<int> que(k);
        for(int i=0;i<n;i++){
            int v;std::cin >> v;
            que.push(v);
            // que.emplace(v);
        }
        while(not que.empty()){
            std::cout << que.top() << std::endl;
            que.pop();
        }
    }
    return 0;
}

11 グラフ

11.1 構成要素

隣接行列を使うか、vector<Edge>みたいのを使うかの二択。場合によってはNode
みたいなのも使うかも。隣接行列を使うとメモリとか探すのとか重い。

struct Edge{
    int to,cost;
    Edge(int to,int cost) : to(to),cost(cost) {};
};

11.2 ベルマンフォード

O(N|E|)

const int INF = 1 << 30;
// s:始点,dist:距離,prev:最短経路木
bool bellman(const vector<vector<Edge> >& graph,int s,vector<int> &dist,vector<int> &prev){
    int n = graph.size();
    for(int i=0;i<n;i++) dist[i] = INF;
    dist[s] = 0;
    for(int i=0;i<n;i++) prev[i] = -1;

    bool neg_cycle = false;
    for(int i=0;i<n;i++){
        for(int j=0;j<n;j++){
            for(size_t k=0;k<graph[j].size();k++){
                const Edge &e = graph[j][k];
                if(dist[e.to] > dist[j] + e.cost){
                    dist[e.to] = dist[j] + e.cost;
                    prev[e.to] = j;
                    if(i == n-1){
                        dist[e.to] = -INF;
                        neg_cycle = true;
                    }
                }
            }
        }
    }
    return !neg_cycle;
}

11.3 ダイクストラ

負の経路があったらダメ

#include <iostream>
#include <vector>
#include <queue>
#include <utility>
#include <algorithm>
#include <functional>

using namespace std;

typedef pair<int,int> pii;

struct Edge{
    int from,to,cost;
    Edge(int from,int to,int cost)
        : from(from),to(to),cost(cost) {};
};

int dijkstra(const int start,const int goal,
             const vector<vector<Edge> > &graph){
    typedef pair<int,int> pii;
    int N = graph.size();
    vector<char> visited(N,false);
    // cost where
    priority_queue<pii,vector<pii>,greater<pii> > que;
    que.push(make_pair(0,start));
    while(not que.empty()){
        int cost,where;
        cost = que.top().first;
        where = que.top().second;
        que.pop();
        if(visited[where]) continue;
        if(where == goal) return cost;
        visited[where] = true;
        for(int j=0;j<(int)graph[where].size();j++){
            que.push(make_pair(graph[where][j].cost+cost,graph[where][j].to));
        }
    }
    return -1;
}


int main(){
    int n,m;
    cin >> n >> m;
    vector<vector<Edge> > V(m);

    for(int i=0;i<n;i++){
        int a,b,cost;
        cin >> a >> b >> cost;
        V[a].push_back(Edge(a,b,cost));
        V[b].push_back(Edge(b,a,cost));
    }

    int ret = -1;
    int p,q;
    cin >> p >> q;
    vector<char> visited(m,false);
    //                 cost,where
    priority_queue<pii,vector<pii>, greater<pii> > Q;
    Q.push(make_pair(0,p));
    while(!Q.empty()){
        int cost,where;
        cost = Q.top().first;
        where = Q.top().second;
        Q.pop();
        if(visited[where]) continue;
        if(where == q){
             ret = cost;
             break;
        }
        visited[where] = true;
        for(int j=0;j<(int)V[where].size();j++){
            Q.push(make_pair(V[where][j].cost+cost,V[where][j].to));
        }
    }
    // 到達不能なときは-1
    cout << ret << endl;
    return 0;
}

11.4 ワーシャルフロイド

負の経路があってもOK。すべてのノードに対してすべてのノードへの距離を求め
る。もし負の閉路があったらiからiはマイナスになる。\(O(|V|^{3})\)

// m はノードの個数。NOはでかい数
vector<vector<int> > V(m,vector<int>(m,NO));
// i->iは0にする。
rep(i,m){
    V[i][i] = 0;
}

for(int k=0;k<m;k++){
    for(int i=0;i<m;i++){
        for(int j=0;j<m;j++){
            V[i][j] = min(V[i][j],V[i][k]+V[k][j]);
        }
    }
}

11.5 最小全域木

プリム法による。\(O(N^2log(N))\)だと思う。最小コストを求めるコードが以下。

struct Edge{
    int to, cost;
    Edge(int to_,int cost_)
        : to(to_),cost(cost_) {}
};


vector<vector<Edge>> minimum_spanning_tree(const vector<vector<Edge>> &graph){
    typedef tuple<int,int,int> cft;  // cost,from,to
    int N = graph.size();
    vector<vector<Edge>> ret(N);
    vector<char> used(N,false);
    priority_queue<cft,vector<cft>,greater<cft>> que;
    que.push(make_tuple(0,-1,0));
    while(not que.empty()){
        int cost = get<0>(que.top());
        int from = get<1>(que.top());
        int to = get<2>(que.top());
        que.pop();
        if(used[to]) continue;
        used[to] = true;
        // ignore first.
        if(from != -1){
            ret[from].push_back(Edge(to,cost));
        }
        for(const Edge& e : graph[to]){
            que.push(make_tuple(e.cost,to,e.to));
        }
    }
    return ret;
}

11.6 最大流

Dinic法による。また、最大流最小カット定理より、最大流と最小カットは一致する。

struct Edge{
    int to,cap,rev;
    Edge(int _to,int _cap,int _rev) : to(_to),cap(_cap),rev(_rev) {};
};

void add_edge(vector<vector<Edge> >& E,int from,int to,int cap){
    E[from].push_back(Edge(to,cap,E[to].size()));
    E[to].push_back(Edge(from,0,E[from].size()-1));
}

vector<int> levels(vector<vector<Edge> > &E,int s){
    vector<int> level(E.size(),-1);
    level[s] = 0;
    queue<int> Q;
    Q.push(s);
    while(!Q.empty()){
        int v = Q.front();
        Q.pop();
        for(size_t i=0;i<E[v].size();i++){
            Edge &e = E[v][i];
            if(e.cap > 0 and level[e.to] == -1){
                level[e.to] = level[v]+1;
                Q.push(e.to);
            }
        }
    }
    return level;
}

int good_path(vector<vector<Edge> > &E,
        vector<int> &iter,
        vector<int> &level,
        int v,int t,int f){
    if(v == t) return f;
    for(int &i=iter[v];i<(int)E[v].size();i++){
        Edge &e = E[v][i];
        if(e.cap > 0 and level[v] < level[e.to]){
            int d = good_path(E,iter,level,e.to,t,min(f,e.cap));
            if(d > 0){
                e.cap -= d;
                E[e.to][e.rev].cap += d;
                return d;
            }
        }
    }
    return 0;
}

int max_flow(vector<vector<Edge> > E,int s,int t){
    int flow = 0;
    const int INF = 1 << 30;
    while(true){
        vector<int> level = levels(E,s);
        if(level[t] < 0) return flow;
        vector<int> iter(E.size());
        int f;
        while((f=good_path(E,iter,level,s,t,INF)) > 0){
            flow += f;
        }
    }
}

int main(){
    int N,M;
    while(cin >> N >> M){
        // [0,N) is cow,[N,N+M) is barn.
        vector<vector<Edge> > E(N+M+2);
        int s = N+M;
        int t = N+M+1;
        for(int i=0;i<N;i++){
            add_edge(E,s,i,1);
        }
        for(int i=0;i<M;i++){
            add_edge(E,N+i,t,1);
        }

        for(int i=0;i<N;i++){
            int S;
            cin >> S;
            for(int j=0;j<S;j++){
                int k;
                cin >> k;
                k--;
                add_edge(E,i,N+k,1);
            }
        }

        cout << max_flow(E,s,t) << endl;
    }
    return 0;
}

11.7 最小費用流

Primal-Dual法による。

typedef pair<int,int> pii;

struct Edge{
    int to,cap,cost,rev;
    Edge(int to,int cap,int cost,int rev)
        : to(to),cap(cap),cost(cost),rev(rev) {};
};

void add_edge(vector<vector<Edge> > &E,int from,int to,int cap,int cost){
    E[from].push_back(Edge(to,cap,cost,E[to].size()));
    E[to].push_back(Edge(from,0,-cost,E[from].size()-1));
}

// s -> t (flow f)
//  if cant, return -1.
int min_cost_flow(vector<vector<Edge> > E,int s,int t,int f){
    const int INF = 1 << 30;
    int ret = 0;
    // potential
    vector<int> h(E.size());
    vector<int> prevv(E.size());
    vector<int> preve(E.size());

    while(f > 0){
        vector<int> dist(E.size(),INF);
        dist[s] = 0;
        priority_queue<pii,vector<pii>,greater<pii> > que;
        que.push(make_pair(0,s));
        while(!que.empty()){
            pii p = que.top();
            que.pop();
            int pf = p.first,ps = p.second;
            if(dist[ps] < pf) continue;
            for(int i=0;i<E[ps].size();i++){
                Edge &e = E[ps][i];
                if(e.cap > 0 and dist[e.to] > dist[ps] + e.cost + h[ps] - h[e.to]){
                    dist[e.to] = dist[ps] + e.cost + h[ps] - h[e.to];
                    prevv[e.to] = ps;
                    preve[e.to] = i;
                    que.push(make_pair(dist[e.to],e.to));
                }
            }
        }
        if(dist[t] == INF){
            return -1;
        }
        for(int v=0;v<E.size();v++){
            h[v] += dist[v];
        }
        int d = f;
        for(int v=t;v!=s;v=prevv[v]){
            d = min(d,E[prevv[v]][preve[v]].cap);
        }
        f -= d;
        ret += d * h[t];
        for(int v=t;v!=s;v=prevv[v]){
            Edge &e = E[prevv[v]][preve[v]];
            e.cap -= d;
            E[v][e.rev].cap += d;
        }
    }
    return ret;
}

int main(){
    while(true){
        int N,M;
        cin >> N >> M;
        if(N == 0 and M == 0) break;
        vector<pair<int,int> > men;
        vector<pair<int,int> > houses;
        for(int h=0;h<N;h++){
            for(int w=0;w<M;w++){
                char x;
                cin >> x;
                if(x == 'H') houses.push_back(make_pair(h,w));
                else if(x == 'm') men.push_back(make_pair(h,w));
            }
        }

        vector<vector<Edge> > E(men.size()+houses.size()+2);
        int s = men.size()+houses.size();
        int t = s+1;

        for(int i=0;i<men.size();i++){
            add_edge(E,s,i,1,0);
            for(int j=0;j<houses.size();j++){
                int dist = abs(men[i].first - houses[j].first)
                    + abs(men[i].second - houses[j].second);
                add_edge(E,i,men.size()+j,1,dist);
            }
        }

        for(int i=0;i<houses.size();i++){
            add_edge(E,men.size()+i,t,1,0);
        }

        cout << min_cost_flow(E,s,t,men.size()) << endl;
    }
    return 0;
}

11.8 無向グラフにおける最小カット

Nagamochi-Ibaraki/Storer-Wagnerの方法によってO(V^3)で計算できる。

struct Edge{
    int to,cap,rev;
    Edge(int _to,int _cap,int _rev) : to(_to),cap(_cap),rev(_rev) {};
};

void add_edge(vector<vector<Edge> >& E,int from,int to,int cap){
    E[from].push_back(Edge(to,cap,E[to].size()));
    E[to].push_back(Edge(from,0,E[from].size()-1));
}

// Nagamochi-Ibaraki/Stoer-Wagner
//  http://www.prefield.com/algorithm/graph/minimum_cut.html
int minimum_cut_of_undirected_graph(const vector<vector<Edge> > &graph){
    const int INF = 1 << 30;
    int n = graph.size();
    vector<vector<int> > adj(n,vector<int>(n));
    for(int i=0;i<n;i++){
        for(size_t j=0;j<graph[i].size();j++){
            Edge e = graph[i][j];
            adj[i][e.to] += e.cap;
        }
    }
    vector<int> h(n);
    for(int i=0;i<n;i++){
        h[i] = i;
    }
    int cut = INF;
    for(int m = n;m > 1;m--){
        vector<int> ws(m,0);
        int u,v,w;
        for(int k=0;k<m;k++){
            u = v;
            v = max_element(ws.begin(),ws.end()) - ws.begin();
            w = ws[v];
            ws[v] = -1;
            for(int i=0;i<m;i++){
                if(ws[i] >= 0) {
                    ws[i] += adj[h[v]][h[i]];
                }
            }
        }
        for(int i=0;i<m;i++){
            adj[h[i]][h[u]] += adj[h[i]][h[v]];
            adj[h[u]][h[i]] += adj[h[v]][h[i]];
        }
        h.erase(h.begin()+v);
        cut = min(cut,w);
    }
    return cut;
}

int main(){
    int N,M;
    while(scanf("%d %d\n",&N,&M) != EOF){
        vector<vector<Edge> > graph(N);
        for(int i=0;i<M;i++){
            int a,b,c;
            scanf("%d %d %d\n",&a,&b,&c);
            add_edge(graph,a,b,c);
            add_edge(graph,b,a,c);
        }
        printf("%d\n",minimum_cut_of_undirected_graph(graph));
    }
    return 0;

11.9 強連結成分分解

struct Edge{
    int to,cost;
    Edge(int t,int c)
        :to(t),cost(c) {
    }
};

void check_back_number(const vector<vector<Edge>>& graph,
                       int v,vector<char>& used,vector<int>& back_number){
    used[v] = true;
    for(const auto& e : graph[v]){
        if(not used[e.to]){
            check_back_number(graph,e.to,used,back_number);
        }
    }
    back_number.push_back(v);
}

void collect_nodes(const vector<vector<Edge>> &graph,
                   int v,vector<char>& used,vector<int>& s){
    used[v] = true;
    s.push_back(v);
    for(const auto& e : graph[v]){
        if(not used[e.to]) {
            collect_nodes(graph,e.to,used,s);
        }
    }
}

vector<vector<int>> strongly_connected_components(const vector<vector<Edge>>& graph){
    const int N = graph.size();
    vector<vector<int>> scc;
    vector<int> back_number;
    {
        vector<char> used(N);
        for(int i=0;i<N;i++){
            if(not used[i]){
                check_back_number(graph,i,used,back_number);
            }
        }
    }

    {
        vector<char> used(N);
        vector<vector<Edge>> reversed_graph(N);
        for(int from=0;from<N;from++){
            for(const auto& e : graph[from]){
                reversed_graph[e.to].push_back(Edge(from,e.cost));
            }
        }
        reverse(all(back_number));
        for(int k : back_number){
            if(not used[k]){
                scc.push_back(vector<int>());
                collect_nodes(reversed_graph,k,used,scc.back());
            }
        }
    }
    return scc;
}
namespace SAT{
    // variable is integer .
    struct Literal{
        bool is_not;
        int var;
        Literal(bool _is_not,int _var)
            : is_not(_is_not),var(_var){}
        Literal not_() const{
            return Literal(not is_not,var);
        }
    };
    Literal make_yes(int var){
        return Literal(false,var);
    }
    Literal make_no(int var){
        return Literal(true,var);
    }
    using Clause = std::vector<Literal>;
    bool solve_2SAT(const int number_of_variable,
                    const vector<Clause>& problem,
                    vector<char>& output){
        // variable is [0..N).
        // "not n" will be assigned to (n+N).
        auto index = [number_of_variable](const bool& is_not,const int& var){
            return number_of_variable*is_not+var;
        };
        auto to_int = [number_of_variable,&index](const Literal& lit){
            return index(lit.is_not,lit.var);
        };

        vector<vector<Edge>> graph(number_of_variable*2);
        for(const Clause& clause : problem){
            // clause = {Literal,Literal}
            // a \/ b -> ~a -> b /\ ~b -> a
            for(int i=0;i<2;i++){
                Literal c = clause[i].not_(),
                        d = clause[not i];
                graph[to_int(c)].push_back(Edge(to_int(d),1));
            }
        }
        vector<vector<int>> scc = strongly_connected_components(graph);
        // where[i] = which component holds i
        vector<int> where(number_of_variable*2);
        for(int i=0;i<scc.size();i++){
            for(int node : scc[i]){
                where[node] = i;
            }
        }
        for(int var=0;var<number_of_variable;var++){
            if(where[index(false,var)] == where[index(true,var)]){
                return false;
            }
        }
        for(int var=0;var<number_of_variable;var++){
            output[var] = where[index(false,var)] < where[index(true,var)];
        }
        return true;
    }
}


// Codeforces 228E
signed main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    int n,m;
    cin >> n >> m;
    vector<SAT::Clause> problem;
    for(int i=0;i<m;i++){
        int a,b,c;
        cin >> a >> b >> c;
        a--;b--;
        if(c){
            // already
            // ~(A/\~B) /\ ~(~A/\B)
            //  = (~A \/ B) /\ (A \/ ~B)
            problem.push_back({SAT::make_no(a),SAT::make_yes(b)});
            problem.push_back({SAT::make_yes(a),SAT::make_no(b)});
        }else{
            // ~(A/\B) /\ ~(~A/\~B)
            //  = (~A \/ ~B) /\ (A \/ B)
            problem.push_back({SAT::make_no(a),SAT::make_no(b)});
            problem.push_back({SAT::make_yes(a),SAT::make_yes(b)});
        }
    }
    vector<char> output(n);
    bool is_solvable = SAT::solve_2SAT(n,problem,output);
    if(is_solvable){
        cout << count(all(output),true) << endl;
        for(int i=0;i<n;i++){
            if(output[i]){
                cout << i+1 << " ";
            }
        }
        cout << endl;
    }else{
        cout << "Impossible" << endl;
    }

    return 0;
}

11.10 橋の列挙

#include <iostream>
#include <vector>
#include <algorithm>

using namespace std;

#define all(c) c.begin(),c.end()
#define repeat(i,n) for(int i=0;i<static_cast<int>(n);i++)

struct Edge{
    int to,cost;
    Edge(int to_,int cost_)
        :to(to_),cost(cost_){
    }
};

int find_bridge_sub(int cur,int from,int number,
                    const vector<vector<Edge>>& graph,
                    vector<int> &pre_order,
                    vector<int> &low,
                    vector<vector<Edge> > &bridges){
    pre_order[cur] = low[cur] = number;
    number++;
    for(const Edge& e : graph[cur]){
        // not visited
        if(pre_order[e.to] == -1){
            number = find_bridge_sub(e.to,cur,number,graph,pre_order,low,bridges);
            if(pre_order[e.to] == low[e.to]){
                bridges[cur].push_back(e);
            }
        }
        if(from != e.to){ //ignore parent.
            low[cur] = min(low[cur],low[e.to]);
        }
    }
    return number;
}

// assume it is connected graph.
vector<vector<Edge>> find_bridge(const vector<vector<Edge>>& graph){
    int N = graph.size();
    vector<int> pre_order(N,-1),lowest_number_of_connected(N,N+1);
    vector<vector<Edge>> ret(N);
    find_bridge_sub(0,-1,0,graph,pre_order,lowest_number_of_connected,ret);

    // optional: make sure from < e.to
    //  You can use following if graph is bidirectional.
    // vector<vector<Edge>> swapped(N);
    // repeat(f,graph.size()){
    //     for(const Edge& e: ret[f]){
    //         if(f < e.to){
    //             swapped[f].push_back(e);
    //         }else{
    //             swapped[e.to].push_back(Edge(f,e.cost));
    //         }
    //     }
    // }
    // ret = swapped;
    return ret;
}


int main(){
    int V_SIZE,E_SIZE;
    cin >> V_SIZE >> E_SIZE;
    vector<vector<Edge>> graph(V_SIZE);
    repeat(i,E_SIZE){
        int s,t;
        cin >> s >> t;
        graph[s].push_back(Edge(t,1));
        graph[t].push_back(Edge(s,1));
    }

    vector<vector<Edge>> bridges = find_bridge(graph);
    repeat(f,graph.size()){
        vector<Edge> &es = bridges[f];
        sort(all(es),[](Edge lhs,Edge rhs){return lhs.to < rhs.to;});
        for(const Edge& e : es){
            cout << f << " " << e.to << endl;
        }
    }
    return 0;
}

11.11 Lowest Common Ancestor

木において、根から最も遠い、u,vの共通の祖先をLCAと呼ぶ。

#include <vector>
#include <cmath>
using namespace std;

// lowest common ancestor.
//  queue version.
struct LCASolver{
    vector<vector<int>> children;
    int root;
    // if root,parent is -1.
    vector<int> parent;
    vector<int> depth;

    vector<vector<int>> parent_pow2;

    LCASolver(vector<vector<int>> children_,int root_)
        : children(children_),root(root_){
        calc_parent_and_depth();
        calc_parent_pow2();
    };
    
    int lca(int u,int v){
        // make sure depth(u) > depth(v).
        if(depth[u] < depth[v]) swap(u,v);
        for(size_t k=0;k<parent_pow2.size();k++){
            if(((depth[u] - depth[v]) >> k) & 1){
                u = parent_pow2[k][u];
            }
        }

        if(u == v) return u;
        for(int k=(int)parent_pow2.size()-1;k>=0;k--){
            if(parent_pow2[k][u] != parent_pow2[k][v]){
                u = parent_pow2[k][u];
                v = parent_pow2[k][v];
            }
        }
        return parent_pow2[0][u];
    }

private:
    void calc_parent_and_depth(){
        parent = vector<int>(children.size(),-1);
        depth = vector<int>(children.size(),-1);
        sub_calc_parent_and_depth(root,-1,0);
    }
    void sub_calc_parent_and_depth(int cur,int par,int dep){
        parent[cur] = par;
        depth[cur] = dep;
        for(int child : children[cur]){
            if(child != par){
                sub_calc_parent_and_depth(child,cur,dep+1);
            }
        }
    }
    void calc_parent_pow2(){
        // parent_pow2[k][i] = 2^k parent of node i.
        parent_pow2 = vector<vector<int>>(ceil(log(children.size())/log(2)),
                                          vector<int>(children.size(),-1));
        parent_pow2[0] = parent;
        for(size_t k=0;k+1<parent_pow2.size();k++){
            for(size_t v=0;v<children.size();v++){
                if(parent_pow2[k][v] >= 0){
                    parent_pow2[k+1][v] = parent_pow2[k][parent_pow2[k][v]];
                }
            }
        }
    }
};

12 幾何

#include <complex>
using namespace std;

const double EPS = 1e-9;
typedef complex<double> point;
typedef vector<point> vertex;

namespace std{
    bool operator < (const point &lhs,const point &rhs){
        if(real(lhs) == real(rhs)){
            return imag(lhs) < imag(rhs);
        }else{
            return real(lhs) < real(rhs);
        }
    }
}

// ベクタの長さ
double vector_length(point a){
    return abs(a);
}

// 二点間距離
double point_distance(point a,point b){
    return abs(a-b);
}

// 単位ベクトル
point unit_vector(point a){
    return a / abs(a);
}

// 法線ベクトル
pair<point,point> normal_vector(point a){
    point n1 = a * point(0,1);
    point n2 = a * point(0,-1);
    return make_pair(n1,n2);
}

// 点が一緒かどうか
bool point_eq(point a,point b){
    return abs(a-b) < EPS;
}
// 内積 (dot product) : a・b = |a||b|cosΘ
double dot(point a,point b){
    return real(conj(a)*b);
}

// 外積 (cross product) : |a×b| = |a||b|sinΘ
double cross(point a,point b){
    return imag(conj(a)*b);
}

// rotate by theta (-pi,+pi)
// arg will return angle of point
point rotate(point p,double theta){
    double cos_theta = cos(theta),
           sin_theta = sin(theta);
    return point(p.real()*cos_theta - p.imag()*sin_theta,
                 p.real()*sin_theta + p.imag()*cos_theta);
}
// 二つの円の共通接線 (たぶんおかしい)
//  http://www.ftext.org/text/subsubsection/1427
//  - 離れている → 4本
//  - 外接している → 3本
//  - 交わっている → 2本
//  - 内接している → 1本
//  - 含む → 0本
vector<pair<point,point> > common_tangent_of_two_circles(point c1,double r1,point c2,double r2){
    // shift c2 by c1.
    point nc = c2 - c1;
    // c.img == 0
    double ar = arg(nc);
    point c = rotate(nc,-ar);

    vector<pair<point,point> > ret;
    vector<double> xs;
    for(int sign=-1;sign<=1;sign+=2){
        xs.push_back((r1*r1 + sign*r1*r2)/c.real());
    }
    for(double x1 : xs){
        double y2 = r1*r1 - x1*x1;
        if(y2 < EPS) continue; // maybe two circle is crossed
        for(int sign=-1;sign<=1;sign+=2){
            double y1 = sign*sqrt(y2);
            // x1*x + y1*y = r1*r1 is such line
            point h = point(x1,y1);
            point t = point(x1+1,(-x1*(x1+1)+r1*r1)/y1); // tekito
            if(abs(y1) < EPS){
                t = point(x1,y1);
            }

            ret.push_back(make_pair(h,t));
        }
    }
    transform(ret.begin(),ret.end(),ret.begin(),[ar,c1](pair<point,point> p){
            return make_pair(rotate(p.first,ar)+c1,
                             rotate(p.second,ar)+c1);
        });
    return ret;
}


// a1,a2を端点とする線分(la)とb1,b2を端点(lb)とする線分の交差判定
bool is_intersected_linesegment(point a1,point a2,point b1,point b2){
    if(max(a1.real(),a2.real()) + EPS < min(b1.real(),b2.real())) return false;
    if(max(b1.real(),b2.real()) + EPS < min(a1.real(),a2.real())) return false;
    if(max(a1.imag(),a2.imag()) + EPS < min(b1.imag(),b2.imag())) return false;
    if(max(b1.imag(),b2.imag()) + EPS < min(a1.imag(),a2.imag())) return false;
    return (cross(a2-a1,b1-a1)*cross(a2-a1,b2-a1) < EPS) and
           (cross(b2-b1,a1-b1)*cross(b2-b1,a2-b1) < EPS);

}

// 与えられた4点が正方形を成すかどうかを判定する
// ただし,ほぼ同じ点が含まれたときにはfalseとする
//  4辺および対角線が1:1:1:1:sqrt(2):sqrt(2)ならばyes
bool is_square(vector<point> ps){
    // assert ps.size() == 4
    sort(ps.begin(),ps.end());
    vector<double> dist;
    for(int i=0;i<4;i++){
        for(int j=i+1;j<4;j++){
            // there is almost same point
            if(abs(ps[i] - ps[j]) < EPS) return false;
            dist.push_back(abs(ps[i]-ps[j]));
        }
    }
    sort(dist.begin(),dist.end());
    for(int i=0;i<4;i++){
        if(abs(dist[0] - dist[i]) > EPS) return false;
    }
    for(int i=4;i<6;i++){
        if(abs(sqrt(2)*dist[1]-dist[i]) > EPS) return false;
    }
    return true;
}


// a1,a2を端点とする線分(la)とb1,b2を端点とする線分(lb)の交点計算
point intersection_point_linesegment(point a1,point a2,point b1,point b2) {
    if(a1 == b1 or a1 == b2) return a1;
    if(a2 == b1 or a2 == b2) return a2;
    point b = b2-b1;
    double d1 = abs(cross(b, a1-b1));
    double d2 = abs(cross(b, a2-b1));
    double t = d1 / (d1 + d2);
    return a1 + (a2-a1) * t;
}

// 線分同士の最短距離
double dist_linesegment_and_linesegment(point a1,point a2,point b1,point b2){
    if(is_intersected_linesegment(a1,a2,b1,b2)){
        return 0;
    }
    return min(min(dist_linesegment_and_point(a1,a2,b1),
               dist_linesegment_and_point(a1,a2,b2)),
           min(dist_linesegment_and_point(b1,b2,a1),
               dist_linesegment_and_point(b1,b2,a2)));
}


// 2直線の直交判定 : a⊥b <=> dot(a, b) = 0
// みけんしょう
bool is_orthogonal(point a1,point a2,point b1,point b2) {
    return dot(a1-a2,b1-b2) < EPS;
}

// 2直線の平行判定 : a//b <=> cross(a, b) = 0
bool is_parallel(point a1,point a2,point b1,point b2) {
    return abs(cross(a2-a1,b2-b1)) < EPS;
}

// a1,a2を通る直線とb1,b2を通る直線の交差判定
bool is_intersected_line(point a1,point a2,point b1,point b2) {
    return not is_parallel(a1,a2,b1,b2);
}

// a1,a2を通る直線とb1,b2を通る直線の交点計算
point intersection_line(point a1,point a2,point b1,point b2) {
    point a = a2 - a1,b = b2 - b1;
    return a1 + a * cross(b, b1-a1) / cross(b, a);
}

// 直線と点との距離
double dist_line_and_point(point a1,point a2,point b){
    return abs(cross(a2-a1,b-a1)) / abs(a2-a1);
}

// 線分と点との距離
double dist_linesegment_and_point(point a1,point a2,point b){
    if(dot(a2-a1,b-a1) < EPS) return abs(b-a1);
    if(dot(a1-a2,b-a2) < EPS) return abs(b-a2);
    return dist_line_and_point(a1,a2,b);
}

// 直線と点の最短距離を実現する直線の点(すいせんの足)(みけんしょう)
point nearest_point_line_and_point(point a1,point a2,point b){
    return a1 + (a2-a1) * dot((a2-a1),(b-a1)) / norm(a2-a1);
}

// 線分と点の最短距離を実現する線分嬢の点(みけんしょう)
point nearest_point_linesegment_and_point(point a1,point a2,point b){
    if(dot(a2-a1,b-a1) < EPS) return a1;
    if(dot(a1-a2,b-a2) < EPS) return a2;
    return nearest_point_line_and_point(a1,a2,b);
}
// 円と線分の交差判定
bool is_cross_linesegment_and_circle(point c,double r,point a1,point a2){
    return (dist_linesegment_and_point(a1,a2,c) < r+EPS and
             (r < abs(c-a1) + EPS or r < abs(c-a2) + EPS));
}


// 点の進行方向
int ccw(point a,point b,point c){
    b -= a;c -= a;
    if(cross(b,c) > 0) return +1;    // counter clockwise
    if(cross(b,c) < 0) return -1;    // clockwise
    if(dot(b,c) < 0) return +2;      // c -- a -- b
    if(norm(b) < norm(c)) return -2; // a -- b -- c
    return 0;
}

// 点が真に多角形(凸?)の中にはいっているか
bool is_inner_point_vertex(const vector<point> &ps,point a){
    int cc = ccw(ps[0],ps[1],a);
    if(not(cc == 1 or cc == -1)) return false;
    for(size_t i=0;i<ps.size();i++){
        if(cc != ccw(ps[i],ps[(i+1)%ps.size()],a)) return false;
    }
    return true;
}

// 点が辺上、もしくは内部にある。(未検証)
bool is_inner_point_vertex_or_line(const vector<point> &ps,point a){
    for(size_t i=0;i<ps.size();i++){
        if(dist_linesegment_and_point(ps[i],ps[(i+1)%ps.size()],a) < EPS){
            return true;
        }
    }
    return is_inner_point_vertex(ps,a);
}


// 凸包 (UVA 109)
vector<point> convex_hull(vector<point> ps){
    int n = ps.size();
    int k = 0;
    sort(ps.begin(),ps.end());
    vector<point> ch(2*n);
    for(int i=0;i<n;ch[k++] = ps[i++]){
        while(k >= 2 and ccw(ch[k-2],ch[k-1],ps[i]) <= 0) --k;
    }
    for(int i=n-2,t=k+1;i>=0;ch[k++]=ps[i--]){
        while(k >= t and ccw(ch[k-2],ch[k-1],ps[i]) <= 0) --k;
    }
    ch.resize(k-1);
    return ch;
}

// remember,pts make convex.
// (http://judge.u-aizu.ac.jp/onlinejudge/cdescription.jsp?cid=ACAC002&pid=C)
double convex_diameter(const vector<point>& pts){
    const int n = pts.size();
    int is=0,js=0; // initial antipode.
    for(int i=1;i<n;i++){
        if(pts[i].imag() > pts[is].imag()) is = i;
        if(pts[i].imag() < pts[js].imag()) js = i;
    }
    double maxd = abs(pts[is]-pts[js]);
    int i,j,maxi,maxj;
    i = maxi = is;
    j = maxj = js;
    do{
        if(cross(pts[(i+1)%n]-pts[i],
                 pts[(j+1)%n]-pts[j]) >= 0){
            j = (j+1)%n;
        }else{
            i = (i+1)%n;
        }
        if(abs(pts[i]-pts[j]) > maxd){
            maxd = abs(pts[i]-pts[j]);
            maxi = i;maxj = j;
        }
    } while(not(i == is and j == js));
    // pts[maxi],pts[maxj] is pair of max diff.
    return maxd;
}

// 円と円の交点(2点ある前提)
vector<point> circles_point(point c1,double r1,point c2,double r2){
    double d = abs(c1-c2);
    double s = (r1+r2+d) / 2;
    double S = sqrt(s*(s-r1)*(s-r2)*(s-d));
    double h = 2 * S / d;
    point v = (c2-c1) / (abs(c2-c1));

    double m = sqrt(r1*r1 - h*h);

    vector<point> ret;
    ret.push_back(c1 + m*v+h*v*point(0,1));
    ret.push_back(c1 + m*v-h*v*point(0,1));
    return ret;
}

// clockwiseだと負
double triangle_area(point a,point b,point c){
    return cross(b-a,c-a)/2;
}

// clockwiseだと負
double vertex_area(vector<point> v){
    double ret = 0;
    for(int i=1;i<v.size()-1;i++){
        ret += triangle_area(v[0],v[i],v[i+1]);
    }
    return ret;
}

12.1 最近点対

#include <complex>
#include <algorithm>
#include <iostream>
#include <vector>


using namespace std;

const double INF = 1e9;
typedef complex<double> point;
double sub_closest_pair_in_field(const vector<point> &ps){
    int N = ps.size();
    if(N <= 3){
        double ret = INF;
        for(size_t i=0;i<ps.size();i++){
            for(size_t j=i+1;j<ps.size();j++){
                ret = min(ret,abs(ps[i]-ps[j]));
            }
        }
        return ret;
    }

    point middle = ps[N/2];
    vector<point> left(ps.begin(),ps.begin()+N/2);
    vector<point> right(ps.begin()+N/2+1,ps.end());

    double d = min(sub_closest_pair_in_field(left),
                   sub_closest_pair_in_field(right));

    vector<point> strip;
    for(const point& p : ps){
        if(abs(middle.real()-p.real()) < d){
            strip.push_back(p);
        }
    }
    sort(strip.begin(),strip.end(),[](const complex<double>& lhs,
                                      const complex<double>& rhs){
                                          return lhs.imag() < rhs.imag();
                                      });

    for(size_t i=0;i<strip.size();i++){
        for(size_t j=i+1;j<strip.size() and
                         (strip[j].imag() - strip[i].imag()) < d;j++){
            d = min(d,abs(strip[i]-strip[j]));
        }
    }
    return d;
}

// find closest pair in field.return distance.
//  O(nlogn) http://www.geeksforgeeks.org/closest-pair-of-points/
double closest_pair_in_field(vector<point> ps){
    sort(ps.begin(),ps.end(),[](const complex<double>& lhs,
                                const complex<double>& rhs){
                                    return lhs.real() < rhs.real();
                                });
    return sub_closest_pair_in_field(ps);
}

int main(){
    int n;cin >> n;
    vector<point> ps(n);
    repeat(i,n){
        double x,y;
        cin >> x >> y;
        ps[i].real(x);
        ps[i].imag(y);
    }

    cout << fixed << setprecision(10);
    cout << closest_pair_in_field(ps) << endl;
    return 0;
}

13 ゲーム

13.1 Nim

いくつかのコインの山がある。この中からプレイヤーは山を一つ選び、1個以上の任意の
数のコインを取る。最後のコインを取ったプレイヤーが勝ちである。この問題に
対しては以下のことが知られている。すべての山のxorをとったとき、それが0で
あるとき、後攻の勝ち、それ以外のときは先攻の勝ち。

#include <iostream>
#include <vector>

using namespace std;

int main(){
    int xor_sum = 0;
    vector<int> coins = {1,2,4,5,6};
    for(size_t i=0;i<coins.size();i++){
        xor_sum = xor_sum ^ coins[i];
    }
    if(xor_sum != 0) cout << "First Player Win" << endl;
    else cout << "Second Player Win" << endl;
}

14 なんか最適化とか

14.1 焼き鈍し

http://shindannin.hatenadiary.com/entry/20121224/1356364040

#include <chrono>
#include <iostream>
#include <algorithm>
using namespace std;
#define debug(x) #x << "=" << (x)
#define dump(x) std::cerr << debug(x) << " (L:" << __LINE__ << ")" << std::endl

unsigned long xor128(){
    static unsigned long x=123456789,y=362436069,z=521288629,w=88675123;
    unsigned long t;
    t=(x^(x<<11));x=y;y=z;z=w; return( w=(w^(w>>19))^(t^(t>>8)) );
}

signed main() {
    ios::sync_with_stdio(false);
    cin.tie(0);

    const chrono::milliseconds annealing_time(2000);
    const auto start_time = chrono::system_clock::now();
    const auto end_time = start_time + annealing_time;
    const int precision_for_force_next = 1000000;
    auto current_time = start_time;

    const auto score = [](double x){
        return max(-pow(-x+100,2)+100,-pow(-x+333,2)+333);
    };

    auto current_state = 0.0;     // initial state
    auto current_score = score(current_state);
    auto best_state = current_state;
    auto best_score = current_score;

    // chrono::system_clock should NOT be called many times...
    while((current_time=chrono::system_clock::now()) < end_time){
        const auto progress = 1.0*(current_time-start_time)/(end_time-start_time);
        const bool force_next = progress*precision_for_force_next
                              < xor128()%precision_for_force_next;

        // [0..100] -> [-50..50] balance is IMPORTANT.
        const auto next_state = current_state + xor128()%101-50;
        const auto next_score = score(next_state);
        if(force_next or (next_score > current_score)){
            current_state = next_state;
            current_score = next_score;
            dump(current_state);
            dump(current_score);
        }else{
            // undo change.
        }

        if(current_score > best_score){
            best_score = current_score;
            best_state = current_state;
            dump(best_score);
        }
    }
    dump(best_score);
    dump(best_state);

    return 0;
}

15 いろんなデータ

15.1 階乗

\(N\) \(N!\)
0 1
1 1
2 2
3 6
4 24
5 120
6 720
7 5040
8 40320
9 362880
10 3628800

15.2 数単位変換

N 日本語 英語
\(10^{0}\) one
\(10^{1}\) ten
\(10^{2}\) hundred
\(10^{3}\) thousand
\(10^{4}\) ten thousand
\(10^{5}\) 十万 hundred thousand
\(10^{6}\) 百万 million
\(10^{7}\) 千万 ten million
\(10^{8}\) hundred million
\(10^{9}\) 十億 billion
\(10^{10}\) 百億 ten billion
\(10^{11}\) 千億 hundred billion

15.3 bit

N \(2^N\) 備考
0 1
1 2 boolの大きさ
2 4
3 8
4 16
7 128 charの最大値+1
8 256 unsigned charの最大値+1
16 65,534
24 16,777,216
31 2,147,483,648 intの最大値+1(about 2*10^9)
32 4,294,967,296 unsigned intの最大値+1
52 4,503,599,627,370,496 doubleのprecision
63 9,223,372,036,854,775,808 long longの最大値+1
64 18,446,744,073,709,551,616 unsigned long longの最大値+1(about 10^19)

15.4 最低限の設定ファイル

15.4.1 vim用

最低限のもの。ホームにおく。

set fileencoding=utf-8
set nocompatible
set t_Co=256
set ambiwidth=double

syntax on
filetype plugin indent on

set nobackup
set noswapfile

set completeopt=menuone
set wildmode=list:longest

set smartindent
set autoindent
set tabstop=4
set softtabstop=4
set shiftwidth=4

set smarttab
set expandtab

set incsearch ignorecase hlsearch
set showmatch
set wildmenu

set listchars=tab:>-,trail:-
set list

set backspace=indent,eol,start