ゆるふわブログ

プログラミング,大学の勉強,日常生活で感じたことをゆるふわに書いていけたらなと思います.技術的に拙いところがあっても温かい目で見守っていただければ幸いです.

ポエム

※ ポエムなので真に受けないでください。

「そういう人は生きている価値がないと思います」とのたまう輩がいる。鈍感なのだ。鈍感だから人を見ただけで、人と話しただけで人を決めつけられるのである。侮れるのである。そういう輩は命の尊さを見ようとしない。人の魂の輝きを、人の可能性を、その誇り高さを理解していない。ある者は見せかけだけの美醜に拘泥し、ある者は自らの尊厳を貶め、皆輝く道を自ら閉ざしてしまう。断言しよう。「生きている価値がない」などということは断じてないのである。弱い君は、君の立ち位置から、君の才能に合わせて、他ならぬ君自身が着実に成長できる道を思い描いて、少しずつでも前に進む。生きる努力をする。それだけである。それを他人が貶めることなどできはしないのである。君の可能性を否定することなどできないのである。君はそんな妄言に惑わされてはいけない。たとえそれがいかに強烈なマインドコントロールをはらんでいてもである。生きてさえいれば、何にだってなれる。その道はいつだって見つけることができるのだ。生きてさえいれば。いかに絶望的な状況に思われても死んではいけない。自らの意志によらず、卑劣な連中に追い込まれて命を散らすことなどあってはならない。

なぜ、平然と人を傷つけられるのだろう?なぜ、人を絶望の底に叩き落とすことになんの躊躇もないのだろう?自らの卑劣さを恥もせず、あまつさえ自分たちに正義があると思っているのである。仕草が気にくわないとか、そんな些細な理由で殺すのである。私は悲しい。人の気持ちがわかるのなら、それを機敏に察知できるのなら、それを踏みにじるのではなく寄り添うこともできたはずではないか?それは彼らの選択なのである。彼らはそれを望んだのである。そんな悲しいことはないではないか!
きれいごとなのは百も承知であるが、「だから貴様らは永遠に争いを止められないのだ」と皮肉を言いたくもなる。
人はなぜ争うのか?お互いを憎み、血を流すことが嵯峨だとでもいうのか?絶望に突き落とされて、人の卑劣さ・狡猾さを骨の髄まで思い知って、その果てにそんなことばかり考える。そこはまさに地獄であった。そんな私は、「文学」に当てられた、正義に燃えた異常者なのである。

アニメや漫画、あの世界はいい。美しい。その一言に尽きる。そこに理想の世界がある。しかし、あれらは「宝石のような毒」なのだと思う。宝石のように輝かしい。涙を流しながら感動すらする。心が洗われ、清められるようである。しかし、同時に生きる力をもがれているような気がする。その高潔さに憧れるほど、その気高さを追い求めるほど、現実を、理不尽な世界を生き抜く冷徹な心を失ってゆく。毒なのだ。

思うに私はこの絶望を味わい尽くさなければならないのではないか?絶望を見つめ、理解しなければならない。たとえそれで自らが崩れ落ちれしまうのだとしても。それが、サイエンスにしがみつくものとしての意地である。使命である。

対等性の具現化、確率漸化式 その一

恩師、青木亮二氏がこのような本を刊行したようである。

入試のツボを押さえる重点学習/数学IAIIB

なかなかに良い本だ。もともと氏は大学への数学に記事を執筆していて、自らの記事をコピーして授業中に配るなどざらであった。(版権は大丈夫なのだろうか?) ゆえに本をパラパラとめくっていると見覚えのある問題が整然と並んでいる。(なお解けなくなっている模様)

この本に陳列されているのは、教科書や一般の各種問題集に掲載された有名かつ論理ギャップの少ない典型題とは一線を画す、珠玉の名問たちである。そのどれもが、一瞥しただけでは ("特殊な" 訓練を積まなければ) 解答が簡単には思い浮かばない、非自明さの芳香を放っている。これこそが数学であると、これこそが力であると、氏の叫ぶ声が聞こえる。ーそして何より重要なことであるがーその問題の主張自体はひどく単純明快なのである。

眠りにつく前にこの本の一節を鑑賞するのが最近の楽しみである。

ここまで、恩師を尊敬するいい生徒を演じてきたわけであるが、実は私はあの男が嫌いである。とても嫌いである。あの男は社会に適応することこそが唯一無二の正しい生き方であると信じて疑わない生粋の「善人」であり、私は随分と屈辱を味わった。私は元来あの大勢の人間が狭い教室へ押し込まれる予備校という空間が嫌いであるが、それに加えてこの恥辱である。それでも私が予備校へと毎度重い足を運んだのは、「この男から何としても数理的技術を盗みとってやる!」という執念、ただそれだけに依っていた。私はあの男を軽蔑していながらも、その数学的技量・センスは誰よりも認めていたのである。それに憧憬すら抱いていたのかもしれない。

さて、氏の悪口はこれぐらいにして本題に入ろう。氏の記事の一節、「何が対等かが問題だ」の軽い解説を書こうと思う。解説記事の解説という随分と不恰好な形にはなるが、"対等性" というものに対する直感の曖昧さ、また確率という分野、特に状態遷移に着目し漸化式を立てて解くという問題の難しさを鑑みれば、わかりやすい解説はいくらあっても多いということはなかろう。多くの受験生がつまづくポイントなのではないかと思う。かつての私もその一人である。

対等性という直感

かつての私は対等性に関するこの記事を読み、真っ先に「本当にそうか?」と思った。
ここでいう対等性とは以下のような "直感" である。

大、中、小の 3 つの玉が入った箱から、玉を順に 2 つ取り出すとき、2 つ目が中である確率は?

10 本中 3 本があたりのくじを次々と引いて行くとき、3 人目があたりを引く確率は?

前者は 1 つ目も 2 つ目も対等であるから  \frac{1}{3}、後者は 1 人目も 3 人目も対等であるから  \frac{3}{10} というのである。本当にそうだろうか???

思うに、直感というのは具体的な計算に裏打ちされていなければならない。ただ「なんとなくそう思うから」というのは究極的には直感でなくただの当てずっぽうである。我々は依拠するところのない思い込みや欺瞞の中で生きているのである。ゆえに往々にして人間の直感は実際とずれる。客観的議論には不向きである。きちんとその正当性が保証されているステップを踏んで、ギャップの小さい、「無理のない」推論こそが必要なのではないか。

ここでは、次のような計算を裏に添えると良いと思う。

「3 つ全てを引くものとして大中小の順列を考える。2 つ目を中で固定すると、順列は  2! 通り。これを全順列  3! 通りで割って、 \frac{2!}{3!} = \frac{1}{3} である。」

「10 本全てを引くものとする。10 本のくじを区別し、順列を考える。3 人目があたりである順列は、3 本のあたりから 1 つを選び ( \binom{3}{1} = 3 通り)、残りを並べる ( 9!)。よって、確率は  \frac{3\cdot 9!}{10!} = \frac{3}{10}。」

単元「場合の数」をきちんと履修した読者諸賢であれば、さして難しくはないと思う。しかし、このような裏付けを逐一考えてみるのが、直感の補正に役立つのではないであろうか。次節はそのような姿勢を持って見ていただきたい。

奇怪なトーナメント

次のような問題を考える。

次図のようなトーナメントを考える。
A〜H の 8 チームがくじ引きで対等にあ〜くに入る。
試合はどちらのチームも勝つ確率は  \frac{1}{2} である。

f:id:Ysmr_Ry:20181004221103j:plain

(1) B が決勝戦に残る確率  p を求めよ

まず、見ての通りトーナメントは二分木である。枝分かれごとに確率  \frac{1}{2} が乗ぜられていく形である。ゆえに全ての葉の確率を足すと 1 である。("決勝戦" は根の 1 つ後であるから和は 2 となる)
愚直に計算する。
あ〜くのうち、B が固定したある 1 箇所に確率は、等しく  \frac{7!}{8!} = \frac{1}{8} である。
各々について、それに葉に対応する確率が乗ぜられるので、 p = \frac{1}{8}\left(\frac{1}{2^4}+\frac{1}{2^4}+\frac{1}{2^4}+\frac{1}{2^4}+\frac{1}{2^3}+\frac{1}{2^3}+\frac{1}{2}+1\right) = \frac{1}{4} である。(長々しいカッコの中は計算するまでもなく 2 なのであるが)

これを、"対等性の直感" でいうと、

「各チームは対等であり、決勝に残る確率は  \frac{1}{8}。決勝戦は 2 チーム参加するため、 p=\frac{2}{8}=\frac{1}{4}

ありうる誤答

このような誤答が考えられる。

「各チームは対等なので、決勝戦に残る確率は  \frac{1}{8} である。」

要するに、「2 をかけるのを忘れませんか? なぜ 2 をかけるのかわかりますか?」ということである。
少なくとも私はー未熟さゆえかもしれないがー初めてこれを見たときにそう感じた。対等性という直感の空漠たる様を感じたのである。(直感の感覚なのでメタ感覚である)

先に見たように、具体的な計算による推論では 2 をかけることはむしろ明らかである。(ー"決勝戦"とは、全てを束ねる根から 1 つ降りた 2 つの広がりのことを指すのだからー)

この推論による裏付けを元に、模糊たる "対等性の直感" の特性をあぶり出してみよう。

<特性 1>
対象が異なるものは "対等" ではない

どういうことであろうか。「決勝戦にたどりつく」という目的において、実は「く」とそれ以外では対等ではないのである。
「く」に B が当たったとして、「決勝戦にたどりつく」とはとりもなおさずトーナメント上右側のノードにたどりつくということである。それに対し、「く」以外ではトーナメント上左側に限定される。つまり、B が当たった場所によって (それが「く」であるか否かによって)、目的地たるノードは都合よく変更されるのである。これは対等ではない。対象は固定されていなければならない。つまり、正しくはこうである。

「決勝戦のうち、トーナメント上左側のノードに到達する確率は、各チームの対等性から  \frac{1}{8}。同じくトーナメント上右側のノードに到達するのも  \frac{1}{8}。ゆえに 2 をかける。」

すなわち、左ノードに着目した時は、「く」を選んだ時の確率 0 を含めての  \frac{1}{8} であり、右ノードに着目した時は、「く」以外を選んだ時の確率 0 を含めての  \frac{1}{8} なのである。ここまでの行程を先は「決勝戦は 2 チーム参加するため」の一言でまとめてしまったのである。(決して同じことを繰り返して言っただけではない。この違い、それこそ「直感」が通ずるであろうか?)

(2) A が あ に入ることだけが決まった。この時、B が決勝戦に残る確率  q を求めよ

これもまず愚直な計算をする。

い〜くのうち、B が特定の 1 箇所に入る確率は、 \frac{6!}{7!} = \frac{1}{7} である。
・い〜きの場合
葉の確率の和はあ〜きで 1 であるから、 1-\frac{1}{2^4} = \frac{15}{16}。よって、 \frac{1}{7}\times\frac{15}{16}
・くの場合
とりもなおさず決勝戦なので、 \frac{1}{7}

よってその和は、 \frac{31}{112}

これを "対等性の直感" では、

「A が決勝に残る場合は、決勝戦右側ノードへ到達する確率であり、対等性から  \frac{1}{7}。また、A が残らない場合は、左右のノード分で  \frac{2}{7}。ゆえに、 \frac{1}{2^4}\times\frac{1}{7}+\left(1-\frac{1}{2^4}\right)\times\frac{2}{7}=\frac{31}{112}。」

これにより、

<特性 2>
固定化されたものの動きを考える

が見える。

すると、

(3) A が あ に入ることだけが決まった。この時、A が B に敗れる確率  r を求めよ

は見えやすくなり、

「A が敗れるのは、優勝しない場合であるから、確率は  1-\frac{1}{2^5} = \frac{31}{32}。そのうち、残りのチームのどれが A に勝つかは対等であるから、 \frac{31}{32}\times\frac{1}{7} = \frac{31}{224}

愚直に計算して裏付けをするのは読者への練習問題とする。

終わりに

いかがであったであろうか。長くなりそうなので、今回は対等性の直感の具現化に的を絞った。
次回は状態遷移と確率漸化式について書く (予定である)

AtCoder Beginner Contest 106 敗戦記 (Haskell/C++)

はい…

A Garden

やるだけ

import Control.Applicative ((<$>))

main :: IO ()
main = do
  [a,b] <- map read . words <$> getLine :: IO [Int]
  print $ a*b-a-b+1
B 105

やるだけ。約数列挙 O(N) でも間に合いますね。

main :: IO ()
main = do
  n <- readLn :: IO Int
  print $ length $ filter odd8 [1..n]

factors :: Int -> [Int]
factors n = filter ((==0).(n`mod`)) [1..n]

odd8 :: Int -> Bool
odd8 n
  | even n    = False
  | otherwise = length (factors n) == 8
C To Infinity

1 でないなら 5000 兆乗するので指数関数的に爆発していきます。
よって、先頭から見ていって最初に現れる 1 でないものが基本的には答えです。
k 文字 1 が続くなら 1 です。(ここら辺のコーナーケースで何回か失敗した記憶) (水色失格)

import Data.List
import Data.Maybe

main :: IO ()
main = do
  s <- getLine
  k <- readLn :: IO Integer
  let none = filter (/='1') s
      noneId = findIndex (/='1') s

  putChar $ if null none
    then '1'
    else if fromJust (noneId) < fromIntegral k 
      then head none
      else '1'

  putStrLn ""
D AtCoder Express 2

できませんでした()
まあ結論をいうと区間を 2 次元座標と捉えて 2 次元累積和をします (この発想がなかった)
[p, q] のクエリがきた時は座標が (p,p), (p,p+1), ..., (p, q), (p+1, p), (p+1, p+1), ... (q,q) のところ ([(p,p),(q,q)] の矩形) の和がわかれば良いです。

#include <cstdio>
#include <algorithm>
#include <vector>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

int N, M, Q;
int I[510][510], S[510][510];

int main()
{
  scanf( "%d%d%d", &N, &M, &Q );
  rep( i, M )
  {
    int L, R;
    scanf( "%d%d", &L, &R );

    ++S[L][R];
  }

  rep( i, N+1 ) rep( j, N+1 )
    S[i][j+1] += S[i][j];

  rep( j, N+1 ) rep( i, N+1 )
    S[i+1][j] += S[i][j];

  rep( t, Q )
  {
    int p, q;
    scanf( "%d%d", &p, &q );

    printf( "%d\n", S[q][q]-S[q][p-1]-S[p-1][q]+S[p-1][p-1] );
  }
  
  return 0;
}
おわりに

言い訳を述べますが、競技プログラミングにおいてこういう瞬間はあるんだと思います。
自分の中にそのアイデアというか概念がなければ、それが即わからない・解けないにつながってしまうのです。
〜色といっても(知らんですが)〜点問題を〜%の正答率で解けるみたいな感じだと思うので、それを 100% に近づけるのは点数低めであっても難しいものなのだろうと。
"自分がいかにわかっていないか" を思い知る瞬間というのは何を勉強していてもあると思うし、それが非常に重要であると思うので、こういう瞬間に巡り会う機会を増やす (コンテストになるべく出る)、巡り合った時に確実に潰す、この 2 つを大事にしていきたいと思いました。
水色といってもまだまだ ABC を確実に全完できるレベルには達していないことが今回わかったし、初心者を早く脱っさなければならないとお尻に火がつきました。
これからも ABC といっても油断せず参加したいと思います。

AtCoder Beginner Contest 105 (by Haskell)

はい。

A AtCoder Crackers

問題よく読んでなかったですが 0 か 1 しかない感。

import Control.Applicative ((<$>))

fi = fromIntegral

main :: IO ()
main = do
  [n,k] <- map read . words <$> getLine :: IO [Int]
  print $ if min (ceiling (fi n/fi k)*k-n) (n-(n`div`k)*k) > 0
    then 1
    else 0

B Cakes and Donuts

ループ回して探索すれば良いと思います。

main :: IO ()
main = do
  n <- readLn :: IO Int
  putStrLn $ if null $ filter (\i -> i*4 <= n && (n-i*4)`mod`7 == 0) [0..n]
    then "No"
    else "Yes"

C Base -2 Number

なんかよく考えたら 2 進数もどきなのでどうとでもなった気がするが、頭が死んでいたので愚直にやった。divMod と quotRem 的な罠にはまった() anamorphism です。

import Data.List
import Data.Char

div' :: Int -> Int -> Int
div' a b
  | a`mod`b == 0 = a`div`b
  | otherwise    = a`div`b+1
  
mod' :: Int -> Int -> Int
mod' a b
  | a`mod`b==0 = 0
  | otherwise  = a`mod`b-b

f :: Int -> Maybe (Char, Int)
f i
  | i == 0    = Nothing
  | i == 1    = Just ('1',0)
  | otherwise = Just (chr (ord '0'+(i`mod'`(-2))), i`div'`(-2))

main :: IO ()
main = do
  n <- readLn :: IO Int
  let cs = unfoldr f n

  putStrLn $ if null cs
    then "0"
    else foldl' (\ans c -> c : ans)  "" cs

D Candy Distribution

一瞬アレっとなりました。C++ で AC してから Haskell で書き直した…
まあ累積和を取るところ、mod M で考えるところまではパッと出てくると思います。
制約をみると O(N^2) で探索させる気はないので、線形オーダーでどうにかできればうれしい。
区間 [i,j) の和が S[j]-S[i] で表せることを思えば、S[i] を mod M の値で分類し、頻度を数えれば良いことがわかります。(S[i] と S[j] の mod M が同じ時その差は M の倍数。mod M が同じものから異なる 2 つを選ぶ場合の数)
mod M が 0 の場合は自身だけを選ぶ場合も含まれることに注意。

import Control.Applicative ((<$>))
import Data.List

main :: IO ()
main = do
  [n,m] <- map read . words <$> getLine :: IO [Integer]
  as <- map read . words <$> getLine :: IO [Integer]
  let
    ss = map (`mod` m) $ scanl1 (+) as
    cs = map (\i -> (head i,genericLength i)) $ group $ sort ss

  print $ sum $ map (\(i,c) -> c*(c-1)`div`2 + (if i == 0 then c else 0)) cs

なんとなく行列を思い浮かべて、対角成分と上三角成分のようなイメージが湧いた()
直交行列の群 O(n) が Cayley 変換により n(n-1)/2 次元の Lie 群とみなせることを思い出して終了しました((((

Mujin Programming Contest

直した

A コンテスト名

やるだけ

#include <iostream>
#include <string>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

std::string S;

int main()
{
  std::cin >> S;
  
  std::cout << (S.substr(0,5) == "MUJIN" ? "Yes" : "No") << std::endl;
  
  return 0;
}
B セキュリティ

やるだけ

#include <iostream>
#include <string>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

int A;
std::string S;

int main()
{
  std::cin >> A >> S;

  if( A == 0 )
  {
    puts("Yes");

    return 0;
  }

  rep( i, S.size() )
  {
    if( S[i] == '+' )
      ++A;
    else
      --A;

    if( A == 0 )
    {
      puts("Yes");

      return 0;
    }
  }

  puts("No");
  
  return 0;
}
C 右折

やるだけ (これに結構時間かかった) (色々と考えて通せたのはうれしかった) (でも 400 点)

#include <cstdio>
#include <vector>
#include <map>
#include <utility>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

using P = std::pair<int, int>;
using ll = long long;

int N, M;
char fld[2010][2010];
int horB[2010][2010], horE[2010][2010];
bool fl[2010][2010];
ll ans;

bool safe( int i, int j )
{ return i >= 0 && i < N && j >= 0 && j < M; }

int main()
{
  scanf( "%d%d", &N, &M );
  rep( i, N )
    scanf( "%s", fld[i] );

  rep( i, N )
  {
    int b = -1;

    rep( j, M )
    {
      if( fld[i][j] == '#' )
        b = -1;
      if( b == -1 && fld[i][j] != '#' )
        b = j;

      horB[i][j] = b;
    }

    b = -1;

    for( int j = M-1; j >= 0; --j )
    {
      if( fld[i][j] == '#' )
        b = -1;

      if( b == -1 && fld[i][j] != '#' )
        b = j;

      horE[i][j] = b;
    }
  }  

  /*rep( i, M ) rep( j, N )
    printf( "%d%c", horB[i][j], j==N-1?'\n':' ' );
  rep( i, M ) rep( j, N )
    printf( "%d%c", horE[i][j], j==N-1?'\n':' ' );
*/
  rep( j, M )
  {
    int b = -1;
    ll sum = 0;

    rep( i, N ) 
    {
      //printf( "(%d,%d)\n", i, j );

      if( b == -1 && fld[i][j] != '#' )
        b = i;

      if( fld[i][j] != '#' && safe( i, j+1 ) && fld[i][j+1] != '#' )
        sum += horE[i][j]-j;//, printf( "+sum: %d\n", horE[i][j]-j );

      if( fld[i][j] != '#' && safe( i, j-1 ) && fld[i][j-1] != '#' )
        sum += j-horB[i][j];//, printf( "+sum: %d\n", j-horB[i][j] );

      if( b != -1 && ( i == N-1 || fld[i+1][j] == '#' ) )
      {
        //printf( "i-b: %d, i: %d, b: %d\n", i-b, i, b );
        //printf( "sum: %lld\n", sum );
        ans += (i-b)*sum;
        sum = 0;
        b = -1;
      }
    }
  }

  printf( "%lld\n", ans );
  
  return 0;
}
D うほょじご

互除法の亜型。規則はない。制約が 1 <= M,N < 1000 なので全て試せる。手順を逆順に探索。(4 つの候補が出る) 手順の途中でも < 1000 は保たれることに注意。rev の逆写像に注意。(‘0’ を前から足せるだけ複数解のある多価関数。逆に rev^{-1}(90) など後ろに 0 があると解はない。(後ろに 0 があるということは rev する前に前に 0 があるということだが、10 進法ではそのようには見ない。)

#include <algorithm>
#include <iostream>
#include <string>
#include <sstream>
#include <queue>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

using P = std::pair<int, int>;

std::vector<int> rev( int x )
{
  std::vector<int> ret;
  std::ostringstream oss;
  oss << x;
  std::string t = oss.str();
    
  //printf( "rev(%d)\n", x );

  if( x % 10 == 0 )
    return ret;

  do
  {
    int y = 0;

    for( int i = t.size()-1; i >= 0; --i )
      y *= 10, y += t[i]-'0';

    ret.push_back( y );

    //printf( "y: %d\n", y );

    if( !x )
      break;

    t = "0"+t;
  } while( t.size() <= 3 );

  return ret;
}

int N, M;
std::queue<P> que;
bool used[1010][1010];
int ans;

int main()
{
  scanf( "%d%d", &N, &M );

  rep( i, 999 )
    que.push( P( i+1, 0 ) );
  rep( j, 999 )
    que.push( P( 0, j+1 ) );

  while( !que.empty() )
  {
    P p = que.front(); que.pop();

    int x = p.first, y = p.second;
    std::vector<int> rys = rev(y), rxs = rev(x);

    if( used[x][y] ) 
      continue;

    //printf( "(%d,%d) %d\n", x, y, que.size() );

    used[x][y] = true;

    if( x+y < 1000 )
    {
      auto rxys = rev(x+y);
      
      for( auto rxy : rxys ) if( rxy < 1000 && y < 1000 && rxy <= y && !used[rxy][y] )
        que.push( P( rxy, y ) );//, printf( "push (%d,%d)\n", rxy, y );
      
      for( auto rxy : rxys ) if( x < 1000 && rxy < 1000 && rxy <= x && !used[x][rxy] )
        que.push( P( x, rxy ) );//, printf( "push (%d,%d)\n", x, rxy );
    }

    for( auto ry : rys ) if( ry < 1000 && x+y < 1000 && ry <= x+y && !used[x+y][ry] )
      que.push( P( x+y, ry ) );//, printf( "push (%d,%d)\n", x+y, ry );
    
    for( auto rx : rxs ) if( rx < 1000 && x+y < 1000 && rx <= x+y && !used[rx][x+y] )
      que.push( P( rx, x+y ) );//, printf( "push (%d,%d)\n", rx, x+y );
  }

  repi( i, 1, N+1 ) repi( j, 1, M+1 ) if( !used[i][j] )
    ++ans;//, printf( "(%d,%d)\n", i, j );

  printf( "%d\n", ans );

  return 0;
}
E 迷路

迷路の問題。bfs でなく dijkstra でやる方が高速。”待ち” の選択がポイント。そこまでの最短路を dist として、dist%K から方向 d (=‘U’ とか) で移動したいとして、どれだけ待てば良いかを nxt[d][dist%K] でもっておけば良い。(nxt[s[i]][i]=0 として、遡って 1 を足す。末尾から最初の 0 にかけては最初まで遡った後 cyclic に通る) INF は 1<<30 だと 1 つテストケース通らない。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <string>
#include <queue>
#include <limits>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()
#define clr(a,v) memset((a),(v),sizeof(a))

using ll = long long;
using P = std::pair<int, int>;
using State = std::pair<ll, P>;

const std::string URDL = "URDL";
constexpr ll dx[4] = { 0, 1, 0, -1 }, dy[4] = { -1, 0, 1, 0 };
constexpr ll INF = std::numeric_limits<ll>::max()>>2;

ll N, M, K;
char d[100010], s[1010][1010];
ll si, sj, gi, gj;
bool used[1010][1010];
std::priority_queue<State, std::vector<State>, std::greater<State> > pque;
ll nxt[4][100010];
ll di[1010][1010];

bool safe( int i, int j )
{ return i >= 0 && i < N && j >= 0 && j < M; }

int main()
{
  scanf( "%lld%lld%lld", &N, &M, &K );
  scanf( "%s", d );
  rep( i, N )
  {
    scanf( "%s", s[i] );

    rep( j, M )
    {
      if( s[i][j] == 'S' )
        si = i, sj = j;
      if( s[i][j] == 'G' )
        gi = i, gj = j;
    }
  }

  std::fill( (ll*)di, (ll*)(di+1010), INF );
  
  clr( nxt, -1 );

  rep( i, K )
    nxt[URDL.find(d[i])][i] = 0;
  
  rep( i, 4 ) 
  {
    bool fl = false;
    ll p = 1;

    for( int j = K-1; j >= 0; --j )
    {
      if( !nxt[i][j] )
        fl = true, p = 0;

      if( fl )
        nxt[i][j] = p++;
    }

    for( int j = K-1; j >= 0; --j )
    {
      if( !nxt[i][j] )
        break;

      if( fl )
        nxt[i][j] = p++;
    }
  }

  pque.push( State( 0, P( si, sj ) ) );
  di[si][sj] = 0;

  while( !pque.empty() )
  {
    State st = pque.top(); pque.pop();

    ll dist = st.first, ci = st.second.first, cj = st.second.second;

    if( di[ci][cj] < dist )
      continue;

    if( ci == gi && cj == gj )
    {
      printf( "%lld\n", dist );

      return 0;
    }

    rep( d, 4 ) if( nxt[d][dist%K] != -1 )
    {
      ll cost = 1+nxt[d][dist%K];
      ll ni = ci+dy[d], nj = cj+dx[d];

      if( safe( ni, nj ) && s[ni][nj] != '#' && di[ni][nj] > dist+cost )
      {
        di[ni][nj] = dist+cost;

        pque.push( State( di[ni][nj], P( ni, nj ) ) );
      }
    }
  }

  puts("-1");

  return 0;
}
F チーム分け

こういう dp しゅき (初見で解けると言ってない)
dp[i][j] : i 人以上のグループを作り、j 人余っている。a[i]==k なる個数 c[k] を求めておくと、dp[i][c[i]+j-k*i] += comb(c[i]+j,k*i)*(k*i)!/(i!)^k/k!*dp[i+1][j] (c[i]+j から k*i を選び、k*i を i グループ k 個に分ける。k*i を並び替え、前から i 個ずつグループを作り、グループ内外の並び替えを考える)

#include <cstdio>
#define repi(i,a,b) for(ll i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

using ll = long long;

constexpr int MAX_N = 1010;
constexpr ll mod = 998244353;

ll N;
ll a[MAX_N], c[MAX_N+1];
ll dp[MAX_N+1][MAX_N+1];
ll fact[MAX_N+1], inv[MAX_N+1];

ll powMod( ll x, ll y )
{
  ll ret = 1;

  while( y > 0 )
  {
    if( y & 1 )
      ret = (ret*x) % mod;

    y >>= 1;
    x = (x*x) % mod;
  }

  return ret;
}

int main()
{
  scanf( "%lld", &N );
  rep( i, N )
    scanf( "%lld", a+i ), ++c[a[i]];
  
  fact[0] = 1;
  rep( i, N )
    fact[i+1] = ((i+1)*fact[i]) % mod;

  rep( i, N+1 )
    inv[i] = powMod( fact[i], mod-2 );

  dp[N+1][0] = 1;

  for( ll i = N; i >= 1; --i ) rep( j, N+1 )
  {
    ll f = 1;

    // dp[1][0] = k=3i=1 dp[2][3]
    for( ll k = 0; k*i <= j+c[i]; ++k )
    {
      //printf( "dp[%d][%d](%lld) += %lld*dp[%d][%d](%lld)\n", i, j+c[i]-k*i, dp[i][j+c[i]-k*i], fact[j+c[i]]*f%mod*inv[j+c[i]-k*i]%mod*inv[k]%mod, i+1, j, dp[i+1][j] );
      dp[i][j+c[i]-k*i] = (dp[i][j+c[i]-k*i]+fact[j+c[i]]*f%mod*inv[j+c[i]-k*i]%mod*inv[k]%mod*dp[i+1][j]%mod) % mod;
      f = (f*inv[i])%mod;
    }
  }

  printf( "%lld\n", dp[1][0] );

  return 0;
}
G 移動

荷重重心。av_1+bv_2+cv_3 (0 <= a,b,c, a+b+c <= K) の作る点は Σ_{k=0}^K 3_H_k = (K+1)(K+2)(K+3)/6 である。それから重複を引く。重複は、xv_1+yv_2+zv_3 = 0 なる (x,y,z) に対して、0 <= a+x, b+y, c+z, (a+x)+(b+y)+(c+z) <= K であるような (a,b,c) の個数。max(0,-x) <= a, max(0,-y) <= b, max(0,-z) <= c であり、max(0,-x)+max(0,-y)+max(0,-z) <= a+b+c <= min(k,k-x-y-z) なので、K = min(k,k-x-y-z)-max(0,-x)-max(0,-y)-max(0,-z) として計算すれば良い。K < 0 の時は、当然和は 0。gcd 取る時は abs 取る。

#include <cstdio>
#include <cstdlib>
#include <algorithm>
#define repi(i,a,b) for(int i=(a);i<(b);++i)
#define rep(i,a) repi(i,0,a)
#define all(a) (a).begin(), (a).end()

using ll = long long;

constexpr ll mod = 998244353;

ll powMod( ll x, ll y )
{
  ll ret = 1;

  while( y > 0 )
  {
    if( y & 1 )
      ret = (ret*x) % mod;

    y >>= 1;
    x = (x*x)%mod;
  }

  return ret;
}

ll func( ll K )
{ return K < 0 ? 0 : (K+1)*(K+2)%mod*(K+3)%mod*powMod(6,mod-2)%mod; }

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

int Q;

int main()
{
  scanf( "%d", &Q );
  rep( i, Q )
  {
    ll a, b, c, d, e, f, k;

    scanf( "%lld%lld%lld%lld%lld%lld%lld", &a, &b, &c, &d, &e, &f, &k );

    ll x = d*e-c*f, y = a*f-b*e, z = -a*d+b*c;
    ll g = gcd(llabs(x),gcd(llabs(y),llabs(z)));

    x /= g, y /= g, z /= g;

    printf( "%lld\n", (func(k)-func(std::min(k,k-x-y-z)-std::max(0ll,-x)-std::max(0ll,-y)-std::max(0ll,-z))+mod)%mod );
  }
  
  return 0;
}
H タイル張り

わからん。とりあえずわかりやすそうなこの辺りのコードを眺めていたがよくわからん (情報量が不足)

mujin-pc-2018.contest.atcoder.jp

mujin-pc-2018.contest.atcoder.jp

最初のやつをメインに見てたけど、集合の冪集合が最大 91 通りしかないらしく、それの隣接行列 (?) (列を 1 列経る毎の状態遷移) を W 乗している??? (繰り返し 2 乗法で高速化)
なんで集合の集合を考えなきゃいけないのかしっくりこない。dfs では そのまま or 2 個連続でフラグが立っているところを潰している っぽいが、そもそもこのフラグが editorial に沿っているのかもわからない (つじつまが合わない?)

教えて強い人()

pV 図上の吸発熱の判定

こんにちは,Ysmr-Ry です.
今回は僕が高校のとき曖昧だったことをまとめてみたいと思います.
意外とコレがなかなか教えてくれなかったりする.

話題は高校の熱力学です.(高校の熱力学の基礎的な知識を仮定します)
pV 図というのを書くじゃないですか.囲む面積が仕事になるので…
大学に行くとエントロピーという状態量を習って,TS 図とかいうのもかけたりするんですね.
それは囲む面積が熱を表すんです.pV と TS どちらも次元がエネルギーになってるからなんですけど.

つまり高校ではエネルギーを供給する熱と仕事のうち仕事にフォーカスを当てて,そちらをわかりやすいようにしてるのです.物事というのは往々にしてトレードオフであり,そうすると熱がわかりにくくなってしまうんですね.

今回はそのわかりにくい熱のお話です.

熱力学とは

熱力学というのは系の熱平衡状態を仮定して,状態量というマクロな少数パラメータで熱力学的状態を記述し,その上での法則を考えようというものです.モル数を固定すると,状態量としては,圧力  p,体積  V絶対温度  T,内部エネルギー  U などがあります.(高校で出てくるのはたぶんこれだけ)
大学に行くともっとたくさん出てくるのですが… そのうちで最も重要なのがエントロピー  S です.まぁざっくり言うと熱を司り,エントロピーが変化しないのなら熱がなく断熱過程です.(ここではその程度の認識で十分)
で,状態量の自由度は 2 なのです.状態量としてはいっぱいあるのですが,そのうちの 2 つの量がわかれば,他の状態量は既知の状態量との関係式をもとに原理的に決定できるのです.(僕はこういう理解なのですが間違ってたら教えてください)
冒頭で pV 図や TS 図が出てきましたが,これらは 2 次元平面です.平面上の 1 点を与えれば,実数 2 つ分の情報量を与えることができ,それが熱力学的な状態と 1 対 1 に対応しているというわけです.
で,高校での主役は pV 図です.その上での熱の話をしたいのです.熱を司るのはエントロピーなのです.
あらわに目にする状態量はこの場合  p, V なわけですが,その裏で  S が働いているという構図です.
 p, V S の関係を熱力学的な法則から考えたいのです.
まずは,"あらわでない状態量が裏で働いている" 様を体感してもらいましょう.

等温過程

あらわでない状態量として,絶対温度  T を考えます.
例えば, T=298K の状態を集めてくると pV 図上ではどうなるのでしょうか.熱力学的状態は pV 図上の 1 点に対応しているので,pV 図上の領域になることは容易にわかります.
それに加えて,自由度が 2 であることも考えると,そのうちの 1 つである  T を固定しているため,実質的な自由度は 1 になります.この場合自由度を次元と置き換えることができ,1 次元,すなわち曲線になることが予想されます*1
高校物理の熱力学で考える作業物質は理想気体がほとんどですから,具体的な曲線は Boyle の法則が与えます.
Boyle の法則により,等温であるとき, pV が一定となります.
 xy 座標平面で言えば  xy=const. であり,これは反比例のグラフです.(直角双曲線) (曲線にちゃんとなった)
また,理想気体の状態方程式から, pV の値は  T に比例します.
よって, T をとびとびに増やしてやると,このように平面を覆っているのがわかります.(一つ一つが  T が一定な状態を集めた双曲線)

右上に行くほど一定である  T の値が大きくなります.
これが "あらわでない状態量が裏で働いている" 例です.

断熱過程

では,裏で働いている状態量がエントロピー  S である場合を考えましょう.
熱の議論をするためには基本的に熱力学第一法則から出発しなければなりません.
熱力学第一法則は内部エネルギーの変化を  dU,気体がする仕事を  d'W,熱を  d'Q として, d'Q=dU+d'W となります.高校では  \Delta を使って書いていましたが,厳密には微小量です*2.ここで  d'W = pdV です.作業物質は理想気体なので, dU = nc_VdT です.( c_V は定積モル比熱)
すると,
\begin{align*}
\frac{d'Q}{T} &= nc_V\frac{dT}{T}+\frac{p}{T}dV \\
&= nc_V\frac{dT}{T}+nR\frac{dV}{V}
\end{align*}
1 行目から 2 行目へは状態方程式  pV=nRT を使いました.
両辺を積分すると,
 \begin{align*}
\int \frac{d'Q}{T} &= nc_V\log T+nR\log V \\
&= nc_V\log TV^\frac{R}{c_V} \\
&= nc_V\log \frac{pV^\frac{c_V+R}{c_V}}{nR}
\end{align*}
これが断熱で 0 ならば,Mayer の関係  c_P = c_V+R と比熱比  \gamma=\frac{c_P}{c_V} を使って, pV^\gamma = const. です.
Poisson の式として教科書にぽつりと載っていたりします.この式の表す pV 図上の曲線を断熱曲線といいますが,そのある接点での接線の傾きは, \frac{dp}{dV} = -\gamma\frac{p}{V} になります.接点の表す状態から接線に沿って進むことは断熱過程を表します.その点からの過程においてこのラインが吸発熱の境界線になります.
ちなみにエントロピーの定義が  \int \frac{d'Q}{T} だったりするので,いわば "裏でエントロピーが働いている" 曲線なわけです.

変形の仕方を少し変えると,
\begin{align*}
d'Q &= nc_VdT+pdV \\
&= nc_V\frac{Vdp+pdV}{nR}+pdV \\
&= \frac{c_V}{R}Vdp+\left(\frac{c_V}{R}+1\right)pdV
\end{align*}
となります.1 行目から 2 行目へは,状態方程式から, dT = \frac{d(pV)}{nR} = \frac{Vdp+pdV}{nR} を使いました.
 d'Q > 0 で吸熱のときは,
 \begin{align*}
\frac{c_V}{R}Vdp+\left(\frac{c_V}{R}+1\right)pdV &> 0 \\
\frac{dp}{dV} &> -\gamma\frac{p}{V}
\end{align*}
となります.
 d'Q < 0 で発熱のときは  \frac{dp}{dV} < -\gamma\frac{p}{V} です.
つまり,ある 1 点の状態からの過程は,そこを接点とする断熱曲線の接線よりも右上ならば吸熱,左下ならば発熱なのです.

図の紫の曲線が断熱曲線です.オレンジの点の状態から,紫の接線に沿って進めば断熱過程,青い線に沿って進めば吸熱,赤い線に沿って進めば発熱です.

つまり,

このように,青い吸熱の領域と赤い発熱の領域に分けられているイメージです.(この場合はオレンジの 1 点に関しての話です.)

おまけ

以下のような熱サイクルを考えてみましょう.

f:id:Ysmr_Ry:20171231030407p:plain

先ほどの話から,左,上に進む過程はそれぞれ発熱,吸熱と決まりますが,斜めの過程はそうとは限りません.
特にこの場合,途中で切り替わるようにできています*3.切り替わるのはこの斜めの線と断熱曲線が接している点です.

右下の方のオレンジの点が切り替わる点です.試しにそれよりも左上のもう 1 つのオレンジの点を見てみましょう.
この点を通る断熱曲線の接線がオレンジの線です.この線よりも熱サイクルの青の斜めの線の方が右上にあるので,この点では吸熱です.接点をすぎた先の点では発熱に切り替わります.
これは, \frac{d^2p}{dV^2} = \gamma\frac{p}{V^2} > 0 より断熱曲線が下に凸なことによります.

理想気体の場合  \gamma=\frac{5}{3} なので*4,この斜め線の方程式と Poisson の式を連立すれば切り替わる接点の座標がわかります.計算すると  (\frac{15}{8}V_0,\frac{9}{8}p_0) であり,この熱サイクルの熱効率は  \frac{48}{97} となります.(間違ってたら教えてください)

*1:ここで直線といいたくなるかもしれませんが,途中折れ曲がっていてもよく,折れ線を限りなく細かくしていくと曲線になるので一般には曲線です.

*2:僕の先生は熱力学第一法則を  \Delta で書くやつを信用するなとまで言っていました.

*3:過去に東工大がこのサイクルを出題したそうです.途中で切り替わるのは出題ミスではと僕の先生が言ってました.

*4:ここらへんは統計力学の結論です.

高速フーリエ変換 (FFT)

お久しぶりです,Ysmr-Ry です.
試験範囲なので高速フーリエ変換 (FFT) について僕の理解を記しておこうと思います.
高速フーリエ変換とは離散フーリエ変換を高速に計算するためのアルゴリズムで,コンピューターで (離散) フーリエ変換をするために重要な手法です.
フーリエ変換・離散フーリエ変換については,下記のようなサイトで調べればわかると思います.たぶん.僕読んでないけど.

www.ic.is.tohoku.ac.jp

成果物

見て楽しくないとあれかなぁと思ったので一応 Visualizer を作ってみました.
矩形波 (もどき) を Fourier 変換して表示してます.Enter を押すとマイクから音を拾ってきます.
もう一度 Enter を押すとその波形で固定して Fourier 変換してくれます.

CodePen - Fourier Visualizer

また,高速フーリエ変換の理解に欠かせないバタフライダイアグラムを自動で作ってくれるのも作ったので記事と合わせてみると良いかもしれません.↑↓ で数を増減できます.

CodePen - Butterfly Diagram Maker

4 つのポイント

以下の 4 点を押さえていればおそらく FFT を理解できるであろうと思います.

  • 分配法則
  • 約分
  • 二進数
  • 再帰的 (漸化式的) 定義

上の 3 つはすぐわかると思います.4 つ目が曲者なのでそのつもりで.
この 4 点を順に確認した後,離散フーリエ変換の高速化を考えます.

分配法則

ご存知だと思いますが一応.分配法則とは,
 2 \times (1+3) = 2\times1+2\times3 = 2+6 = 8 という風に,かっこの中の各項に掛け算を分配することを言います.

約分

これも大丈夫でしょう.
 \frac{2}{8} = \frac{1}{4} のようにすることです.

二進数

コンピューターと相性のいい数の表し方です.
我々が普段使っているのは十進数で,桁が上がるごとに桁の担っている数が 10 倍になっていきます.一の位,十の位,百の位,千の位…という風に.
二進数はこれが 2 倍になっていきます.一の位,ニの位,四の位,八の位…という風です.
例えば,11 を表したいなら,1+2+8 より,1011 となります.

再帰的 (漸化式的) 定義

ここまでは「何を当たり前のことを」という感じだったかもしれませんが,これはちょっと難しいです.
これは高校で習ったであろう漸化式を思い出してもらえるとわかりやすいです.
ご存知フィボナッチ数列が良さそうですね.
フィボナッチ数列とは,1, 1 から始めて直前の 2 つの数を足し合わせて次を作っていく数列です.
 1, 1, 2, 3, 5, 8, 13, \cdots
左から数えて n 番目のフィボナッチ数列の数を  f_n と書くと,「直前の 2 つの数を足し合わせて次を作っていく」は,
 f_n = f_{n-1} + f_{n-2}\quad(n\geq3)
と表せます.これが漸化式です.
フィボナッチ数列の作り方のように,最初 ( f_1, f_2) から順に作っていけば  f_n を得られますが,今は何もない状態からぽっと  f_n を聞かれた状況を考えてみましょう.最初から数列を作っていないので何もわかりません.
その時この漸化式を見ると,「n 番目を求めるには n-1 番目と n-2 番目を求めて足せばいいのさ!」と語りかけてきます.
そう,左から n 番目を求める問題を,n-1 番目と n-2 番目を求めると言い換えられるのです.
n-1 番目は n-2 番目と n-3 番目がわかればよく,n-2 番目は…と繰り返していけば,いつか 1 番目と 2 番目がわかればわかるというところまできます.1, 2 番目は 1 と決まっていますから,問題が解けたことになります.
このように,より小さな (この場合は左側) な問題の組み合わせに言い換えていくことができるのです.
よって,問題を解くのに必要なのは,より 小さな問題に言い換えるルール (漸化式)もっとも単純な問題の答え (1, 2 番目が 1) だけです.
これを覚えておいてください.

離散フーリエ変換 (DFT)

 W_N = e^{-i\frac{2\pi}{N}} とし,これを回転子と言います.(いきなり何)
難しそうですがなんのことはなく,名前の通り 360° を N 等分した分時計回りに回ることを表します.
例えば  W_4 によって円周を一周すると図のようになります.

f:id:Ysmr_Ry:20171108232907p:plain

はい.ちなみに円の右端を 0° として時計回りに回ってます.番号の順に回ってます.
 W_N^k W_N の回転を  k 回連続で行うことを表します.
ここでは簡単のため, W_N^k を円の中心から回転で行き着いた点を結ぶ矢印で表します.
例えば, W_4^2 ならば,円の右端から半周するので左端にいます.図の 3 の位置です.
よって,矢印 ← で表します.
すると,大きさ 8 の離散フーリエ変換は与えられた  x_i\,(i=1,2,\cdots,8) に対し,次の行列計算を行うことを指します.

 \begin{align*}
\begin{pmatrix}
X_0 \\ X_1 \\ X_2 \\ X_3 \\ X_4 \\ X_5 \\ X_6 \\ X_7
\end{pmatrix} = 
\begin{pmatrix}
\mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} \\
\mbox{→} & \mbox{↘} & \mbox{↓} & \mbox{↙} & \mbox{←} & \mbox{↖} & \mbox{↑} & \mbox{↗} \\
\mbox{→} & \mbox{↓} & \mbox{←} & \mbox{↑} & \mbox{→} & \mbox{↓} & \mbox{←} & \mbox{↑} \\
\mbox{→} & \mbox{↙} & \mbox{↑} & \mbox{↘} & \mbox{←} & \mbox{↗} & \mbox{↓} & \mbox{↖} \\
\mbox{→} & \mbox{←} & \mbox{→} & \mbox{←} & \mbox{→} & \mbox{←} & \mbox{→} & \mbox{←} \\
\mbox{→} & \mbox{↖} & \mbox{↓} & \mbox{↗} & \mbox{←} & \mbox{↘} & \mbox{↑} & \mbox{↙} \\
\mbox{→} & \mbox{↑} & \mbox{←} & \mbox{↓} & \mbox{→} & \mbox{↑} & \mbox{←} & \mbox{↓} \\
\mbox{→} & \mbox{↗} & \mbox{↑} & \mbox{↖} & \mbox{←} & \mbox{↙} & \mbox{↓} & \mbox{↘} \\
\end{pmatrix}\begin{pmatrix}
x_0 \\ x_1 \\ x_2 \\ x_3 \\ x_4 \\ x_5 \\ x_6 \\ x_7
\end{pmatrix}
\end{align*}

行列の 2 行目は 8 等分を 1 個飛ばしで,k 行目は 8 等分を k-1 個飛ばしで回っていることがわかります.
これは通常  O(n^2) の計算量になるのですが,隠れた対称性を見出すことでこの計算量を  O(n\log n) まで落とすのが高速フーリエ変換です.
バブルソートクイックソートぐらい違いますね.

高速フーリエ変換 (FFT)

高速フーリエ変換を一言でいうと,「 \frac{2}{8} って約分して  \frac{1}{4} にできんじゃん!」です.たぶん.
つまり, W_8^2 で 8 等分を 2 個飛ばしで回るのと, W_4^1 で 4 等分を順に回るのとでは区別がつきません.

f:id:Ysmr_Ry:20171109002725p:plain

f:id:Ysmr_Ry:20171108232907p:plain

これです.
これによって,行を半分に分けた時,回転していく様子が "対等" となり分配法則によってくくれるのです.
サイズの小さい場合から具体的に見てみましょう.

n = 2 の場合

離散フーリエ変換は,

 \begin{align*}
\begin{pmatrix} X_0 \\ X_1 \end{pmatrix} = \begin{pmatrix}
\mbox{→} & \mbox{→} \\ \mbox{→} & \mbox{←}
\end{pmatrix}\begin{pmatrix}
  x_0 \\ x_1
\end{pmatrix}
\end{align*}

となります.
ここで,行を半分 (1 個ずつ) に分けてみましょう.
1 行目はどちらも同じです.
2 行目は右側は左側を半周 (←) 回したものです.
よって,
 \begin{align*}
\begin{pmatrix} X_0 \\ X_1 \end{pmatrix} &= \begin{pmatrix}
\mbox{→} x_0 + \mbox{→} x_1 \\
\mbox{→} x_0 + \mbox{←} x_1 \\
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→}(x_0 +x_1) \\
\mbox{→}(x_0 + \mbox{←} x_1) \\
\end{pmatrix} \\
&= \mbox{→}
\begin{pmatrix}
x_0 +x_1 \\ x_0 + \mbox{←} x_1
\end{pmatrix} \\
&=
\begin{pmatrix}
x_0 +x_1 \\ x_0 + \mbox{←} x_1
\end{pmatrix}
\end{align*}

 \mbox{→} = W_2^0 = 1 ですからこうなります.
これを図にしてみるとこうなります.

f:id:Ysmr_Ry:20171109004714p:plain

左から線をたどっていきます.合流したら足します.
三角があったら流れてきた数に傍に書いてある回転子を書きます.
この図は, X_0 = x_0 + W_2^0 x_1 = x_0 + x_1 X_1 = x_0 + W_2^1 x_1 = x_0 + \mbox{←} x_1 を表しています.

この図を,蝶のようなのでバタフライダイアグラムといいます.

n = 4 の場合

離散フーリエ変換は,

 \begin{align*}
\begin{pmatrix}
X_0 \\ X_1 \\ X_2 \\ X_3
\end{pmatrix} = \begin{pmatrix}
\mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} \\
\mbox{→} & \mbox{↓} & \mbox{←} & \mbox{↑} \\
\mbox{→} & \mbox{←} & \mbox{→} & \mbox{←} \\
\mbox{→} & \mbox{↑} & \mbox{←} & \mbox{↓} \\
\end{pmatrix} \begin{pmatrix}
x_0 \\ x_1 \\ x_2 \\ x_3
\end{pmatrix}
\end{align*}

今度は偶数行目と奇数行目で考えてみましょう.偶数行目は行を半分に分けた時に同じものの繰り返しですが,奇数行目は半周分ずれています.
よって.

 \begin{align*}
\begin{pmatrix}
X_0 \\ X_2
\end{pmatrix} &= \begin{pmatrix}
\mbox{→} & \mbox{→} & \mbox{→} & \mbox{→} \\
\mbox{→} & \mbox{←} & \mbox{→} & \mbox{←} \\
\end{pmatrix} \begin{pmatrix}
x_0 \\ x_1 \\ x_2 \\ x_3
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→} x_0 + \mbox{→} x_1 + \mbox{→} x_2 + \mbox{→} x_3 \\
\mbox{→} x_0 + \mbox{←} x_1 + \mbox{→} x_2 + \mbox{←} x_3
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→} (x_0+x_2) + \mbox{→} (x_1+x_3) \\
\mbox{→} (x_0+x_2) + \mbox{←} (x_1+x_3) \\
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→} & \mbox{→} \\
\mbox{→} & \mbox{←} \\
\end{pmatrix} \begin{pmatrix}
x_0+x_2 \\ x_1+x_3
\end{pmatrix} \\
\begin{pmatrix}
X_1 \\ X_3
\end{pmatrix} &= \begin{pmatrix}
\mbox{→} & \mbox{↓} & \mbox{←} & \mbox{↑} \\
\mbox{→} & \mbox{↑} & \mbox{←} & \mbox{↓} \\
\end{pmatrix} \begin{pmatrix}
x_0 \\ x_1 \\ x_2 \\ x_3
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→} x_0 + \mbox{↓} x_1 + \mbox{←} x_2 + \mbox{↑} x_3 \\
\mbox{→} x_0 + \mbox{↑} x_1 + \mbox{←} x_2 + \mbox{↓} x_3
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→} (x_0+\mbox{←}x_2) + \mbox{↓} (x_1+\mbox{←}x_3) \\
\mbox{→} (x_0+\mbox{←}x_2) + \mbox{↑} (x_1+\mbox{←}x_3) \\
\end{pmatrix} \\
&= \begin{pmatrix}
\mbox{→} & \mbox{↓} \\
\mbox{→} & \mbox{↑} \\
\end{pmatrix} \begin{pmatrix}
x_0+\mbox{←}x_2 \\ x_1+\mbox{←}x_3
\end{pmatrix} \\
\end{align*}

ここまでが第 1 段階です.これは n = 2 の場合とやっていることが全く同じです.
n = 2 の場合の  x_0 +x_1 が偶数行目に  x_0 + \mbox{←} x_1 が奇数行目に対応しています.(n = 2 の場合はそれぞれ 1 個ずつだったのが,今は 2 個ずつ)
しかし,添字の番号が違っていますね. x_0 + x_1 に対応するのが  x_0+x_2, x_1+x_3 です.ここにどのような対応があるのでしょう?
答えは二進数にあります.添字を二進数で表してみると,(0,2) が (00,10), (1,3) が (01,11) です.二進数の一番上の桁に注目すると,どちらの組も (0,1) となっていることがわかります.これがルールです.
これが第 k 段階ならば上から k 桁だけに注目します.あとで詳しく考えます.
次に,
 \begin{align*}
\begin{pmatrix}
X_0 \\ X_2
\end{pmatrix} &= \begin{pmatrix}
x_0+x_2 + \mbox{→}(x_1+x_3) \\ x_0+x_2 + \mbox{←}(x_1+x_3)
\end{pmatrix} \\
\begin{pmatrix}
X_1 \\ X_3
\end{pmatrix} &= \begin{pmatrix}
x_0+\mbox{←}x_2 + \mbox{↓}(x_1+\mbox{←}x_3) \\ x_0+\mbox{←}x_2 + \mbox{↑}(x_1+\mbox{←}x_3)
\end{pmatrix}
\end{align*}

で完成です.(第 2 段階)
これをバタフライダイアグラムで表すと,

f:id:Ysmr_Ry:20171109093939p:plain

となります.左側の列の蝶の集まりが第 1 段階,右側が第 2 段階に対応しています.
第 1 段階は n = 2 の場合に対応すると言ったように,n = 2 の時のバタフライが 2 つ入っています.
サイズ n の FFT を実行することを, fft_n(x_i) と表すとすると, fft_4(x_0,x_2,x_1,x_3) = fft_2(x_0,x_2) + fft_2(x_1,x_3) + (\mbox{一番右のバタフライ}) です.これが漸化式です.

n = 8 の場合

これも同様です.バタフライダイアグラムだけ示すと,

f:id:Ysmr_Ry:20171109094547p:plain

となります.1 列目 (第 1 段階) には n = 2 の場合のバタフライが 4 つ,1, 2 列目を合わせると n = 4 の場合のバタフライが 2 つ隠れています.
この場合も, fft_8(x_0,x_4,x_2,x_6,x_1,x_5,x_3,x_7) = fft_4(x_0,x_4,x_2,x_6) + fft_4(x_1,x_5,x_3,x_7) + (\mbox{一番右のバタフライ}) が漸化式になります.

先ほど確認したように,問題を解くには漸化式の他にもっとも単純な問題の答えが必要です.
もっとも単純なのはこの場合 n = 1 の場合を考えると簡単です.n = 1 の場合は,フーリエ変換する前もなくそのままになります.
 fft_1(x_i) = x_i です.これで問題が解けます.

解けるんですが,あと問題なのは,バタフライダイアグラムの左側の  x_i の並びです.これにはある規則があります.

なぜビットリバースの順なのか

バタフライダイアグラムの左側の  x_i の右に書いてあるかっこの中の数は添字を二進数で表したものです.
これをそれぞれ逆から読んでみると,000, 001, 010, 011, 100, … で,0, 1, 2, 3, 4, … と順に増えていることがわかります.これをビットリバースの順と言います.
なぜこうなるのか考えてみましょう.漸化式に沿って考えてみると,右から 2 番目のバタフライが 2 つある段階 (第  \log_2 n-1 段階) において,添字を二進数で表した時の上から  \log_2 n-1 桁だけ見ています.言い換えると,一番右の桁だけ無視しています.なので,( \log_2 n-1 桁)0 と ( \log_2 n-1 桁)1 の 2 つのセットがあります.これが 2 つのバタフライに対応します.
漸化式から一番右のバタフライと 2 つの半分のサイズの  fft に言い換えられるため,半分のサイズのバタフライダイアグラムに今と同じことを順に適用していくと,n = 2 の場合,一番上の桁が 0, 1 なので,ビットリバースの順に増えていくことに対応します.
例えば,n = 4 の場合,第 1 段階は 0, 1 ですが,第 2 段階まで考えると,その後ろに 0 をつけるか,1 をつけるかの 2 セットあります.これを順に書くと,00, 10, 01, 11 でビットリバースの順になっています.
n = 8 ならば,00, 10, 01, 11 の後ろに 0 をつけるか 1 をつけるかで,000, 100, 010, 110, 001, 101, 011, 111 です.
ビットリバースの順に増えていくのは,そこまでの列の後ろに 0 を加えたものと 1 を加えたものをくっつけて拡張するということを繰り返していくと作れます.
ここにも再帰的定義が隠れているのです.
(伝わる?)

まとめ

個人的にネットに転がってる FFT の説明がわかりにくかったので自分なりにまとめてみました.
徐々に難しくしていったつもりなんですが,ビットリバースのところがうまく伝わったか心配です…
不明な点があったらコメントください.間違っていた場合もお願いします.