« 2010年1月 | トップページ | 2010年3月 »

2010年2月

2010年2月21日 (日)

Project Euler : Problem 16 ~ 整数のリスト化

 問題はこちらをご覧ください。

 

 方針としては「2 の 1000乗を計算して、それをリストに変換して合計を出す」ということになると思います。
 それでは、どうやって整数をリストに変換するか?

 十進法限定ならば、かなり簡単にできます。

import Data.Char (digitToInt) -- 整数をリストに変換 dexToList :: Integral a => a -> [Int] dexToList = (map digitToInt) . show problem016 :: Integral a => a -> Int problem016 = sum . dexToList main = print $ problem016 (2 ^ 1000)

 これだけだと面白くないので、基数を指定して整数をリストに変換する関数を考えてみました。

-- 整数をリストに変換 -- ex : intToList 10 123 => [1,2,3] -- ex : intToList 2 123 => [1,1,1,1,0,1,1] intToList :: Integral a => a -> a -> [a] intToList radix n = iter n [] where iter n prod | n < radix = n : prod | otherwise = iter q (r : prod) where (q, r) = quotRem n radix problem016 :: Integral a => a -> a problem016 = sum . (intToList 10) main = print $ problem016 (2 ^ 1000)

Project Euler : Problem 15 ~ 組み合わせ

 問題はこちらをご覧ください。

 

 組み合わせの問題。

-- 下降階乗冪 fallingFactorial :: Integral a => a -> a -> a fallingFactorial m n = product [m - n + 1 .. m] -- 組み合わせの数 combSize :: Integral a => a -> a -> a combSize n r = div (fallingFactorial n r) (fallingFactorial r r) problem015 :: Integral a => a -> a -> a problem015 m n = combSize (m + n) n main = print $ problem015 20 20

Project Euler : Problem 14 ~ Collatz 問題とメモ化

 問題はこちらをご覧ください。

 

 お馴染みの「Collatz 問題」ですね。
 まずは、Brute force 版です。

import Data.List (maximumBy) import Data.Ord (comparing) collatz :: Integral a => a -> Int -> Int collatz 1 c = c collatz n c | even n = collatz (div n 2) (c + 1) | otherwise = collatz (n * 3 + 1) (c + 1) problem014 :: Integral a => a -> (a, Int) problem014 n = maximumBy (comparing snd) xs where xs = [(x, collatz x 1) | x <- [2 .. n]] main = print $ problem014 1000000
 単純に定義どおりにステップ数を数えただけです。最後に "maximumBy" でステップ数の一番多いものを探しました。
 何の工夫もしていないので約 14 秒もかかってしまいました。

 そこで Ruby 1.9 で解いたときと同じように、配列によるメモ化の手法を使ってみました。
 Haskell の配列は書き換えにコストがかかるという話を聞いたことがあったのですが、今回コードでは書き換えは行われないので大丈夫だと思います。

import Data.Array import Data.List (maximumBy) import Data.Ord (comparing) problem014 :: Integer -> (Integer, Int) problem014 n = maximumBy (comparing snd) (assocs arr) where arr = listArray (1, n) (1 : [collatz i i 0 | i <- [2 .. n]]) collatz n m c | m < n = (arr ! m) + c | even m = collatz n (div m 2) (c + 1) | otherwise = collatz n (m * 3 + 1) (c + 1) main = print $ problem014 1000000
 今回は約 2.7 秒で答えが出ました。やはり、かなりのスピードアップになりました。
 Haskell の配列はリストと違っていろいろと制約があったりして、使いこなすにはまだまだ勉強が必要そうです。

2010年2月20日 (土)

Haskell で「順列」と「組み合わせ」

 "Project Euler" などの数学的な問題を解くときに「辞書順に生成された順列」を使いたい時があります。( ここでいう「辞書順」とは、例えば "[1,2,3]" を関数に渡したときに "[[1,2,3],[1,3,2],[2,1,3],[2,3,1],[3,1,2],[3,2,1]]" が返ってくることを言います )
 Hasekell には Data.List モジュールに "permutations" という順列を生成する関数があるのですが、これが辞書順の答えを返してくれないんです。"[1,2,3]" を渡すと "[[1,2,3],[2,1,3],[3,2,1],[2,3,1],[3,1,2],[1,3,2]]" が返ってきます。

 そこで、Haskell の練習も兼ねて自分で作ってみることにしました。
 紙と鉛筆で順列を考える場合のように、樹形図を作っていくイメージです。

permutation :: [a] -> Int -> [[a]] permutation [] _ = [[]] permutation ns size = perm size (length ns - 1) [(ns, [])] where perm 0 _ xs = [a | (_, a) <- xs] perm c n xs = perm (c - 1) (n - 1) $ concatMap (f n) xs f n (xs, ys) = [(as ++ bs, ys ++ [b]) | x <- [0 .. n], let (as, b : bs) = splitAt x xs]
*Main> permutation [1 .. 4] 2 [[1,2],[1,3],[1,4],[2,1],[2,3],[2,4],[3,1],[3,2],[3,4],[4,1],[4,2],[4,3]]

 ちょっと行数が多い気はしますが、ネットで見つけたいくつかの「辞書順に順列を生成する関数」と比べても 4 〜 5 倍の速さで答えを出せるようなので自分としては満足かなと……。

 

 さらにこの関数をもとに「組み合わせを生成する関数」もつくってみました。

combination :: [a] -> Int -> [[a]] combination [] _ = [[]] combination ns size = comb size [(ns, [])] where comb 0 xs = [a | (_, a) <- xs] comb c xs = comb (c - 1) $ concatMap comb' xs comb' (x : xs, ys) = (xs, ys ++ [x]) : comb' (xs, ys) comb' _ = []
*Main> combination [1 .. 4] 2 [[1,2],[1,3],[1,4],[2,3],[2,4],[3,4]]

2010年2月 9日 (火)

Project Euler : Problem 13

 問題はこちらをご覧ください。

 

 今回は、かなり手抜きです。解説も省略。

nums :: [Integer] nums = [37107287533902102798797998220837590246510135740250, 46376937677490009712648124896970078050417018260538, 74324986199524741059474233309513058123726617309629, 91942213363574161572522430563301811072406154908250, 23067588207539346171171980310421047513778063246676, 89261670696623633820136378418383684178734361726757, 28112879812849979408065481931592621691275889832738, 44274228917432520321923589422876796487670272189318, 47451445736001306439091167216856844588711603153276, 70386486105843025439939619828917593665686757934951, 62176457141856560629502157223196586755079324193331, 64906352462741904929101432445813822663347944758178, 92575867718337217661963751590579239728245598838407, 58203565325359399008402633568948830189458628227828, 80181199384826282014278194139940567587151170094390, 35398664372827112653829987240784473053190104293586, 86515506006295864861532075273371959191420517255829, 71693888707715466499115593487603532921714970056938, 54370070576826684624621495650076471787294438377604, 53282654108756828443191190634694037855217779295145, 36123272525000296071075082563815656710885258350721, 45876576172410976447339110607218265236877223636045, 17423706905851860660448207621209813287860733969412, 81142660418086830619328460811191061556940512689692, 51934325451728388641918047049293215058642563049483, 62467221648435076201727918039944693004732956340691, 15732444386908125794514089057706229429197107928209, 55037687525678773091862540744969844508330393682126, 18336384825330154686196124348767681297534375946515, 80386287592878490201521685554828717201219257766954, 78182833757993103614740356856449095527097864797581, 16726320100436897842553539920931837441497806860984, 48403098129077791799088218795327364475675590848030, 87086987551392711854517078544161852424320693150332, 59959406895756536782107074926966537676326235447210, 69793950679652694742597709739166693763042633987085, 41052684708299085211399427365734116182760315001271, 65378607361501080857009149939512557028198746004375, 35829035317434717326932123578154982629742552737307, 94953759765105305946966067683156574377167401875275, 88902802571733229619176668713819931811048770190271, 25267680276078003013678680992525463401061632866526, 36270218540497705585629946580636237993140746255962, 24074486908231174977792365466257246923322810917141, 91430288197103288597806669760892938638285025333403, 34413065578016127815921815005561868836468420090470, 23053081172816430487623791969842487255036638784583, 11487696932154902810424020138335124462181441773470, 63783299490636259666498587618221225225512486764533, 67720186971698544312419572409913959008952310058822, 95548255300263520781532296796249481641953868218774, 76085327132285723110424803456124867697064507995236, 37774242535411291684276865538926205024910326572967, 23701913275725675285653248258265463092207058596522, 29798860272258331913126375147341994889534765745501, 18495701454879288984856827726077713721403798879715, 38298203783031473527721580348144513491373226651381, 34829543829199918180278916522431027392251122869539, 40957953066405232632538044100059654939159879593635, 29746152185502371307642255121183693803580388584903, 41698116222072977186158236678424689157993532961922, 62467957194401269043877107275048102390895523597457, 23189706772547915061505504953922979530901129967519, 86188088225875314529584099251203829009407770775672, 11306739708304724483816533873502340845647058077308, 82959174767140363198008187129011875491310547126581, 97623331044818386269515456334926366572897563400500, 42846280183517070527831839425882145521227251250327, 55121603546981200581762165212827652751691296897789, 32238195734329339946437501907836945765883352399886, 75506164965184775180738168837861091527357929701337, 62177842752192623401942399639168044983993173312731, 32924185707147349566916674687634660915035914677504, 99518671430235219628894890102423325116913619626622, 73267460800591547471830798392868535206946944540724, 76841822524674417161514036427982273348055556214818, 97142617910342598647204516893989422179826088076852, 87783646182799346313767754307809363333018982642090, 10848802521674670883215120185883543223812876952786, 71329612474782464538636993009049310363619763878039, 62184073572399794223406235393808339651327408011116, 66627891981488087797941876876144230030984490851411, 60661826293682836764744779239180335110989069790714, 85786944089552990653640447425576083659976645795096, 66024396409905389607120198219976047599490197230297, 64913982680032973156037120041377903785566085089252, 16730939319872750275468906903707539413042652315011, 94809377245048795150954100921645863754710598436791, 78639167021187492431995700641917969777599028300699, 15368713711936614952811305876380278410754449733078, 40789923115535562561142322423255033685442488917353, 44889911501440648020369068063960672322193204149535, 41503128880339536053299340368006977710650566631954, 81234880673210146739058568557934581403627822703280, 82616570773948327592232845941706525094512325230608, 22918802058777319719839450180888072429661980811197, 77158542502016545090413245809786882778948721859617, 72107838435069186155435662884062257473692284509516, 20849603980134001723930671666823555245252804609722, 53503534226472524250874054075591789781264330331690] problem013 :: Integer problem013 = iter $ sum nums where iter n | n < 10 ^ 10 = n | otherwise = iter (div n 10) main = print problem013

2010年2月 1日 (月)

既約分数クイズ、ファレイ数列、フォードの円

 こちらのブログ「既約分数クイズ」の存在を知りました。

 ちょっと面白そうだったので、Haskell でコードを書いてみました。答えは (分子, 分母) の組がリストの形で出てきます。

irreducibleFraction :: Int -> [(Int, Int)] irreducibleFraction n = [(p, q) | q <- [1 .. n], p <- [0 .. q], gcd p q == 1]
*Main> irreducibleFraction 4 [(0,1),(1,1),(1,2),(1,3),(2,3),(1,4),(3,4)]
 Haskell では条件を書くだけで、何の工夫もせずに答えが出てしまいました。

 「他の解法はないものか?」と調べていたら、既約分数に関連して「ファレイ数列」というものがあることが分かりました。これは既約分数が大きさの順に並んでいるものです。
 そこで、先ほどの関数を使って「ファレイ数列」を作ってみることにしました。ついでにより分数らしく見えるようにしてみました。

import Ratio import Data.List irreducibleFraction :: Int -> [(Int, Int)] irreducibleFraction n = [(p, q) | q <- [1 .. n], p <- [0 .. q], gcd p q == 1] farey :: Int -> [(Int, Int)] farey n = sortBy comp (irreducibleFraction n) where comp (p1, q1) (p2, q2) = compare (p1 % q1) (p2 % q2) fareyToStr :: [(Int, Int)] -> String fareyToStr ns = init $ concat $ map toStr ns where toStr (a, b) = show a ++ "/" ++ show b ++ " "
*Main> fareyToStr $ farey 4 "0/1 1/4 1/3 1/2 2/3 3/4 1/1"
 さらに調べると、「ファレイ数列」には「中関数」を使って求める方法があるようなので、自分でも考えて見ました。
farey2 :: Int -> [(Int, Int)] farey2 n = farey' [(0, 1), (1, 1)] where farey' ((x1, y1) : xs@((x2, y2) : _)) | y1 + y2 > n = (x1, y1) : farey' xs | otherwise = farey' ((x1, y1) : (x1 + x2, y1 + y2) : xs) farey' xs = xs
*Main> fareyToStr $ farey2 4 "0/1 1/4 1/3 1/2 2/3 3/4 1/1"

 

 ところで Wikipedia の「ファレイ数列」の項には「フォードの円」という図が載っています。
 「これって、自作の "Turtle Graphics on Ruby/SDL" で描けるかな?」ということでコードを書いてみました。

# -*- coding: utf-8 -*- require 'turtle_graphics' # # == ファレイ数列 # def farey(n) q_arr = (1 .. n).to_a p_arr = Array.new f_arr = Array.new farey = q_arr.each do |q| p_arr = (0 .. q).to_a f_arr.push(p_arr.map {|p| [p, q]}) end f_arr = f_arr.flatten(1).select{|p, q| p.gcd(q) == 1} return f_arr.sort_by{|p, q| p / q.to_f} return end # # == フォードの円を描く # class Screen def ford_circle(p, q, c) q = q.to_f r = 1 / (2 * q * q) draw_circle(p / q, r, r, c).update # draw_circle(p / q, 1 - r, r, c).update end end sc = Screen.create(600, 600) sc.set_range(0, 1, 0, 1) # 描画範囲 (x_min, x_max, y_min, y_max) farey(10).each do |p, q| sc.ford_circle(p, q, White) end sc.main_loop

Screenshotturtle_graphics_on_rubysd

 Wikipedia には F5 の図が載っていますが、今回のコードでは F10 の図を描いてみました。
 当然、コードを書き換えれば F50 だろうが F100 だろうが表示してくれます(相当時間がかかると思われますが……)
 また、描画範囲を書き換えて拡大表示させれば、小さい円もそれぞれしっかり他の円に接しているのがわかります。
 ただし、あまり大きく拡大表示すると、0/1 の円と 1/1 の円の表示が乱れるようです。これは "turtle_graphics.rb" の "Screen#draw_circle" 内で使用している、"draw_aa_circle" というメソッドの問題のようです。この部分を "draw_circle" に書き換えれば、どんなに拡大しても表示は乱れないと思われます。

« 2010年1月 | トップページ | 2010年3月 »

2016年7月
          1 2
3 4 5 6 7 8 9
10 11 12 13 14 15 16
17 18 19 20 21 22 23
24 25 26 27 28 29 30
31            
フォト

最近のトラックバック

無料ブログはココログ