AOJ 2858 Prime-Factor Prime

解法

区間篩の要領で解ける。
[1, 1e5] ぐらいまでの素数を普通に篩で検出し、素数が見つかれば [L, R] にカウントを割り振っていく。
今回の場合は、素数 p が検出できたら、[L, R] には p, p^2, ... の倍数のところに1をカウントしていくことになる。
しかしこれだけだと不十分で、素因数分解でよくある (小さい素数)*(でかい素数)のでかい素数のほうがカウントされない。
なので、カウントする配列の他に、[L, R] の値を直に持っている配列を用意し、割り振るたびにその値を p で割っていくとよい。
最後に、カウンタの値と、割られて残った数が1であるか(1でないならそれも素因数なので、追加で1カウント)によって、その値が prime factor prime であるか判定する。

計算量は、篩の中で p, p^2, ... を見ていくから、結局 O(PlogRloglogP + (R - L))
(ただし、P は篩のテーブルのサイズ)で、体感はもっと早い。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using ll = long long;

int main() {
    const int max_p = 1 << 18;

    int l, r; cin >> l >> r;

    vector<int> ns(r - l + 1);
    for(int i = 0; i <= r - l; ++i) {
        ns[i] = i + l;
    }
    vector<int> cnt(r - l + 1);
    vector<bool> is_prime(max_p, true);
    is_prime[0] = is_prime[1] = false;
    for(int i = 2; i < max_p; ++i) {
        if(is_prime[i]) {
            for(int j = i + i; j < max_p; j += i) is_prime[j] = false;
            ll x = i;
            while(x <= r) {
                for(int j = (l + x - 1) / x * x; j <= r; j += x) {
                    cnt[j - l] += 1;
                    ns[j - l] /= i;
                }
                x *= i;
            }
        }
    }
    int ans = 0;
    for(int i = 0; i < r - l + 1; ++i) {
        if(ns[i] > 1) cnt[i] += 1;
        ans += is_prime[cnt[i]];
    }
    cout << ans << endl;
}

AOJ 2593 Square in Circles

解法

辺の長さをまず決め打つと、それぞれの円と交わる部分が区間として求まる。
得られた区間の連結な部分の総長さで、決めうった値よりも長い所があれば、そのような正方形を配置することができる。
これは O(n(logn)^2) で解ける。

しゃくとりとかスタックとか使ってうまくやれば O(nlogn) でも解けそうだけど、まあ定数倍もそこまでじゃないし楽な方を選択すれば良さそう。

ありえそうな間違いといえば、隣接する2つの円だけみて判定とか?
(3つ以上の円がかなり近い位置に並んでいると、横幅で3つ以上の円を利用できることがあるのでNG)

ソースコード

#include <bits/stdc++.h>
using namespace std;

using pdd = pair<double, double>;

double solve(vector<double> x, vector<double> r) {
    const int n = x.size();
    auto check = [&] (double h) {
        vector<pdd> ps(n);
        for(int i = 0; i < n; ++i) {
            double x1 = x[i], x2 = x[i];
            if(h <= r[i]) {
                x1 -= sqrt(r[i] * r[i] - h * h);
                x2 += sqrt(r[i] * r[i] - h * h);
            }
            ps[i] = make_pair(x1, x2);
        }
        sort(begin(ps), end(ps));
        double lx = -1e9, rx = -2e9;
        for(int i = 0; i < n; ++i) {
            if(rx < ps[i].first) {
                if(rx - lx >= 2 * h) return true;
                lx = ps[i].first, rx = ps[i].second;
            } else {
                rx = max(rx, ps[i].second);
            }
        }
        return rx - lx >= 2 * h;
    };
    double lb = 0, ub = 2e5;
    for(int lp = 0; lp < 60; ++lp) {
        const auto mid = (lb + ub) / 2;
        (check(mid) ? lb : ub) = mid;
    }
    return lb * 2;
}

int main() {
    int n;
    while(cin >> n, n) {
        vector<double> x(n), r(n);
        for(int i = 0; i < n; ++i) {
            cin >> x[i] >> r[i];
        }
        cout << fixed << setprecision(10) << solve(move(x), move(r)) << endl;
    }
}

AOJ 1620 Boolean Expression Compressor

解法

前もってすべての式を生成してしまえばOK。
真理値表は 16bit 分の情報があればよく、bitDP 的な要領で生成していけば良い。
つまり、キーが真理値表で、値がその最小長さである。
生成するのにはダイクストラが楽だと思う。
ダイクストラで遷移するときは、今見ている式(これは真理値表と長さのペアのことを指している)を、これまでに作ったすべての式と組み合わせて新しい式を作る。
これまでに作ったすべての式を見る部分の計算量がヤバそう(式としてありえるのが 2^16 個あるから)だが、式の長さが高々16文字までしか使えないことを考えると、実は生成できる真理値表はそもそもかなり少ない。
なので、なんかいい感じに作ろうと頑張らなくてもサボって通せる。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using pii = pair<int, int>;

constexpr int inf = 1e9;

int eval_impl(string const& s, int& p, const int S) {
    if(s[p] == '0' || s[p] == '1') {
        return s[p++] == '1';
    }
    if(isalpha(s[p])) {
        return (S >> (s[p++] - 'a')) & 1;
    }
    if(s[p] == '-') {
        return !eval_impl(s, ++p, S);
    }
    if(s[p] == '(') {
        const bool t1 = eval_impl(s, ++p, S);
        const char op = s[p++];
        const bool t2 = eval_impl(s, p, S);
        p++; // )
        if(op == '^') return t1 ^ t2;
        else          return t1 & t2;
    }
    assert(false);
    return -1;
}

int eval(string const& s) {
    int truth = 0;
    for(int S = 0; S < 1 << 4; ++S) {
        int p = 0;
        truth |= (eval_impl(s, p, S) << S);
    }
    return truth;
}

int main() {
    constexpr int all_mask = (1 << 16) - 1;
    vector<int> ans(1 << 16, inf);
    ans[0] = ans[(1 << 16) - 1] = 1;
    ans[eval("a")] = ans[eval("b")] = ans[eval("c")] = ans[eval("d")] = 1;
    vector<int> ts = {0, (1 << 16) - 1, eval("a"), eval("b"), eval("c"), eval("d")};
    vector<int> len = {1, 1, 1, 1, 1, 1};
    priority_queue<pii, vector<pii>, greater<>> que;
    for(int i = 0; i < (int)ts.size(); ++i) {
        que.emplace(len[i], ts[i]);
    }
    while(!que.empty()) {
        const auto l = que.top().first;
        const auto t = que.top().second;
        que.pop();
        if(ans[t ^ all_mask] > l + 1) {
            ans[t ^ all_mask] = l + 1;
            que.emplace(l + 1, t ^ all_mask);
            ts.push_back(t ^ all_mask), len.push_back(l + 1);
        }
        const int m = ts.size();
        for(int i = 0; i < m; ++i) {
            if(len[i] + l + 3 > 16) continue;
            if(ans[t & ts[i]] > len[i] + l + 3) {
                ans[t & ts[i]] = len[i] + l + 3;
                ts.push_back(t & ts[i]), len.push_back(len[i] + l + 3);
                que.emplace(len[i] + l + 3, t & ts[i]);
            }
            if(ans[t ^ ts[i]] > len[i] + l + 3) {
                ans[t ^ ts[i]] = len[i] + l + 3;
                ts.push_back(t ^ ts[i]), len.push_back(len[i] + l + 3);
                que.emplace(len[i] + l + 3, t ^ ts[i]);
            }
        }
    }

    string exp;
    while(cin >> exp, exp != ".") {
        const auto t = eval(exp);
        cout << ans[t] << endl;
    }
}

感想

一昨年は結局本番で通せなかったなあ。

他の人のコードを読むと頑張って賢く生成しようとしていて大変そうだった。
とりあえず愚直生成で遅かったら考える、ぐらいの気持ち(運良く早くてよかった)。

AOJ 2373 HullMarathon

解法

自分は数式で捉えて、ラグランジュの未定乗数法で解いた。
ラグランジュの未定乗数法については以下のサイトがわかりやすい(ただし厳密な証明はない)。
mathtrain.jp

与えられた問題を定式化する。
今、r[0], r[1], ..., r[n - 1] をこの順に反時計回りに使うとして、それぞれの間の成す角を θ[0], ..., θ[n - 1] とする。
すると、以下の問題になる。

maximize sum(r[i] * r[(i + 1) % n] * sin(θ[i]) / 2) (for i = 0, ..., n - 1)
s.t. sum(θ[i]) = 2π (for i = 0, ..., n - 1)

ラグランジュの未定乗数法の通りに方程式を立てると、
r[0] * r[1] / (r[i] * r[(i + 1) % n]) * cos(θ[0]) = cos(θ[i])  (forall i = 0, ..., n - 1)
を満たすようなθを求める問題になる。これは二分探索でできる。

凸包パートが必要そうに見えるが、全探索するとサボれる。
r のうちどれを使うか最初に決めて、それらだけで next_permutation を行うと、単に三角形の面積の和を求めるだけでよくなる。
ただし、next_permutation をするときに先頭を固定してしわないと定数倍で落ちるので注意(最初を固定して定数倍落とすのが嬉しい場面はたまにあるので覚えて損はない)。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using ld = long double;

const ld pi = acos(-1.0);
constexpr ld eps = 1e-9;

int main() {
    int n; cin >> n;
    vector<ld> rr(n);
    for(int i = 0; i < n; ++i) {
        cin >> rr[i];
    }
    sort(begin(rr), end(rr));

    ld ans = 0;
    for(int S = 1; S < (1 << n); ++S) {
        vector<ld> r;
        for(int i = 0; i < n; ++i) {
            if(S & (1 << i)) {
                r.push_back(rr[i]);
            }
        }
        const int m = r.size();
        if(m < 3) continue;
        do {
            ld tans = 0, ang_sum = 0;
            auto check = [&] (ld theta) {
                ang_sum = theta;
                tans = 0.5 * r[0] * r[1] * sin(theta);
                for(int i = 1; i < m; ++i) {
                    const auto ang = acos(max(-1.0L, min(1.0L, r[0] * r[1] / r[i] / r[(i + 1) % m] * cos(theta))));
                    tans += 0.5 * r[i] * r[(i + 1) % m] * sin(ang);
                    ang_sum += ang;
                }
                return ang_sum < 2 * pi;
            };
            ld lb = 0, ub = pi;
            for(int lp = 0; lp < 100; ++lp) {
                const auto mid = (lb + ub) / 2;
                (check(mid) ? lb : ub) = mid;
            }
            if(abs(ang_sum - 2 * pi) < eps) {
                ans = max(ans, tans);
            }
        } while(next_permutation(begin(r) + 1, end(r)));
    }
    cout << fixed << setprecision(10) << ans << endl;
}

VK Cup 2015 - Finals A. Matching Names

問題概要

n 個の名前と、n 個のペンネームが与えられる。これらを1対1対応させることを考える。
名前 a とペンネーム b をペアにした時の満足度を lcp(a, b) と定義する。
lcp(a, b) の総和が最大となるような、1対1対応を出力せよ。

1 <= n <= 10^5
文字列の長さの総和は 8 * 10^5

解法

lcp なので脳死で SA から攻めたくなりそうだが、lcp を Trie の上で考える視点も必要。
Trie における2文字列間の lcp は、LCA になる。
また、lcp(a, b) = (|a| + |b| - (|a| + |b| - 2 * lcp(a, b)) / 2 であることを考えると、lcp(a, b) の総和の最大化は(|a|, |b| が定数なので)以下のように言い換えられる。

木の上に n 頂点ずつ位置が与えられるので、それらをいい感じ組分けして最短距離の和の最小化をせよ

これは Trie を構築した後に dfs で下の方から貪欲的にペアを作っていけば達成できる。形式的に証明しなくても直感的にわかると思う。

これで計算量は O(total_length) となる。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using pii = pair<int, int>;

struct alphabet_t { // default: 'a' - 'z'
    static constexpr int size = 26;
    static constexpr int char_to_index(char c) {
        assert('a' <= c && c <= 'z');
        return c - 'a';
    }
    static constexpr char index_to_char(int idx) {
        assert(0 <= idx && idx < size);
        return 'a' + idx;
    }
};

template <typename Alpha>
class trie {
    static constexpr int alpha_sz = Alpha::size;

public:
    void insert(std::string const& s, int tag, int idx) {
        auto cur_node = this;
        for(int p = 0; p < static_cast<int>(s.size()); ++p) {
            const auto c = Alpha::char_to_index(s[p]);
            if(!cur_node->next[c]) cur_node->next[c] = std::make_unique<trie<Alpha>>();
            cur_node = cur_node->next[c].get();
        }
        cur_node->idxs[tag].push_back(idx);
    }

    // 本質
    pair<int, array<vector<int>, 2>> solve(std::vector<pii>& ans) {
        auto res = make_pair(0, idxs);
        for(int i = 0; i < alpha_sz; ++i) {
            if(!next[i]) continue;
            auto ch = next[i]->solve(ans);
            res.second[0].insert(end(res.second[0]), begin(ch.second[0]), end(ch.second[0]));
            res.second[1].insert(end(res.second[1]), begin(ch.second[1]), end(ch.second[1]));
            res.first += ch.first;
        }
        while(!res.second[0].empty() && !res.second[1].empty()) {
            ans.emplace_back(res.second[0].back(), res.second[1].back());
            res.second[0].pop_back();
            res.second[1].pop_back();
        }
        assert(res.second[0].empty() || res.second[1].empty());
        res.first += res.second[0].size() + res.second[1].size();
        return res;
    }

private:
    std::array<std::unique_ptr<trie<Alpha>>, alpha_sz> next;
    array<vector<int>, 2> idxs;
};


int main() {
    int n; cin >> n;
    trie<alphabet_t> t;
    int total_len = 0;
    for(int i = 0; i < 2 * n; ++i) {
        string s; cin >> s;
        t.insert(s, i >= n, (i >= n ? i - n + 1 : i + 1));
        total_len += s.size();
    }
    vector<pii> ans;
    cout << (total_len - t.solve(ans).first) / 2 << endl;
    for(auto& p : ans) {
        cout << p.first << " " << p.second << endl;
    }
}

AOJ 0353 Sort

解法

とりあえず観察しないとよくわからない。
僕は観察してもよくわからなかったので、確実に言えることを探すことにした。
まず、一番小さい数は、左端に到達したら以降はもう動かない。その後に、二番目に小さい数が左から二番目の位置に到達したら以降は動かない。以下同様。
この事実を使うなら、まだ位置が確定していない値のうちで、最も小さい値が目的の位置に移動するまでにどのような変化が起こるか観察すべきだ。
そこで、次に位置を確定させる値に注目したとして、その値が適切な位置に移るまでを考える。
1 から i - 1 番目までは確定しているとする。次に注目するのは i 番目に小さい数である。その数が位置 j にあるとしよう。
この数が位置 i に到達するまでに発生するスワップ回数は、j - i 回である。これは考えてみると当然である。
これは、i から j 番目までの数字の並びは、スワップ回数に全く関係ないことを意味する。したがって、次に注目したい数が今どの位置にあるのかさえわかればよい、ということにほかならない。
この位置を効率よく探索するにはどうすればよいだろうか?要素がごっそり移り変わるので、単純なセグツリを使うということはできない(平衡二分木なら可能かも)。
これを考えるためには、後ろに追いやった数の並びについて、確実に言えることはなにかに注目しなければならない。
そして、いろいろ試すとわかるが、後ろに追いやられた値の全体としての並びは複雑でつらい気持ちになる。
しかしよく観察すると、後ろに追いやられた数の中で最小の値は、常に末尾に存在することに気がつく。
我々が必要なのは次に確定させたい値、それはつまり(選んでない中で)一番小さな値であり、その候補となる値以外の位置はどうでも良い。
そしてその候補となるのが、後ろに追いやった数の集合で最も小さい値である。
それらの位置は、前から位置を確定させるごとに後ろに追いやった数の集合の要素数さえわかれば、累積和の要領で求めることができる。

これまでの考察から、以下の解法を得る。

小さい値から順に注目し、前から位置を確定させていく。
このとき、注目した値を目的の位置に移るまでに後ろに追いやる数の集合をまとめ上げてしまう。この集合は、含まれる中で最も小さいものを得られるようなものにする(priority_queue など)。
注目した値が適切な位置に至るまでのスワップ回数は、まとめ上げた集合の要素数(の合計)に等しい。ただし注目している値自身は除く。
以降、この集合は分解することなく管理する。次に確定させる値のときも同様に、自身の位置から目的の位置までにある集合をまとめ上げ、自分以外を後ろに追いやるようにする。
これを最後の要素になるまで繰り返すと、答えが得られる。
集合のまとめ上げにマージテクを行えば、O(N(logN)^2) で解ける。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using ll = long long;

int main() {
    int n; cin >> n;
    vector<int> a(n);
    for(int i = 0; i < n; ++i) {
        cin >> a[i];
    }

    { // compress
        auto b = a;
        sort(begin(b), end(b));
        for(int i = 0; i < n; ++i) {
            a[i] = lower_bound(begin(b), end(b), a[i]) - begin(b);
        }
    }

    using S = priority_queue<int, vector<int>, greater<>>;
    function<void(S&, S&)> merge = [] (S& s1, S& s2) {
        if(s1.size() < s2.size()) swap(s1, s2);
        while(!s2.empty()) {
            s1.push(s2.top()); s2.pop();
        }
    };

    queue<shared_ptr<S>> que;
    for(int i = 0; i < n; ++i) {
        auto q = make_shared<S>(); q->push(a[i]);
        que.push(move(q));
    }
    ll ans = 0;
    for(int i = 0; i < n; ++i) {
        auto s = que.front(); que.pop();
        while(s->top() != i) {
            auto v = que.front(); que.pop();
            merge(*s, *v);
        }
        s->pop();
        ans += s->size();
        que.push(s);
    }

    cout << ans << endl;
}

感想

難しい。解説するのも(言葉だけだと)難しい。
かといって図を載せて解説する元気もなく…(公式解説を読んでください)。
よくこんな問題思いつくなあ。

AOJ 0323 Ruins

解法

二分探索するとよい。
海岸からの距離を y とすると、半円 i との交点の x 座標は、x[i] ± sqrt(r[i] * r[i] - y * y) となる。
これを区間と見て、すべての半円に対応する区間に共通区間があることが、その高さで共通領域がある条件である。
このギリギリを二分探索で調べて終わり。

ソースコード

#include <bits/stdc++.h>
using namespace std;

using ld = long double;

int main() {
    int n; cin >> n;
    vector<ld> x(n), r(n);
    for(int i = 0; i < n; ++i) {
        cin >> x[i] >> r[i];
    }

    auto check = [&] (ld y) {
        ld mini = -1e18, maxi = 1e18;
        for(int i = 0; i < n; ++i) {
            const auto d = sqrt(r[i] * r[i] - y * y);
            mini = max(mini, x[i] - d);
            maxi = min(maxi, x[i] + d);
        }
        return mini <= maxi;
    };
    ld lb = 0, ub = *min_element(begin(r), end(r));
    for(int lp = 0; lp <= 100; ++lp) {
        const auto mid = (lb + ub) / 2;
        (check(mid) ? lb : ub) = mid;
    }
    cout << setprecision(10) << fixed << lb << endl;
}