本履歴

購入した古本の履歴と時々プログラミング

数学パズル本に敢てプログラムで挑戦(8)

マーチン・ガードナーサム・ロイド本の続き

69. The Canals on Mars

火星で発見された水路を通って各町を必ず一回回ってスタート地点に戻ってくる問題。火星の半球面全体に水路を張るのはすごい技術というか無駄な気がするけど、まあプログラム的にはハミルトン閉路を求めるものですね。各都市には文字が振られていて、正解路は何か意味を持つとのこと。これはチェックが大変だと思ったら、正解は実質一通りしかなかった。

通常ハミルトンがあるかは解けないですが、今回は辺の数が少ないので幅優先で探索したらさくっと解答が出ました。いつもこのくらいだと嬉しいんですが。実行時間がコーディング時間より長くなると悲しいよね。まあ、グラフの構造を入力するのが一番大変かな。 くだらないのが、本によると多くの読者(数千人)が”THERE IS NO POSSIBLE WAY”(そんな道存在しないよ)と言ってきたとのこと。ナンジャと思ったらそれが正解メッセージでしたとさ。本当にくだらない。

cities="THYAWNREIOIELPESSBOS".chars.to_a
edges=[
    [1,2],
    [0,5,9,10,11],
    [0,3,4,8,17],
    [2,4],
    [2,3,7,12,14],
    [1,9,13,16,18,19],
    [10,11,14],
    [4,10,12],
    [2,15,17],
    [1,5,13],
    [1,6,7,11,14,16],#10
    [1,6,10],
    [4,7,17],
    [5,9,18],
    [4,6,10],
    [8,17,19],
    [5,10,17],
    [2,8,12,15,16,19],
    [5,13,19],
    [5,15,17,18]
]

n=cities.size
ans=[]
s=s_bit=0
q=[]
init_state=[[s],s_bit]
q.push init_state
g_bit=(?1*n).to_i(2)

until q.empty?
    state=q.shift
    seq,bit=state
    u=seq.last
    if bit==g_bit
        ans << seq if u==s #hamiltonian cycle path
        next
    end
    edges[u].each{|v|
        if bit[v]==0
            next_state=[seq.dup << v, bit+(1<<v)]
            q.push next_state
        end
    }
end

puts ans.size
ans.each{|seq|puts seq[0..-2].map{|i|cities[i]}.join} #=> THEREISNOPOSSIBLEWAY

f:id:uru:20170430162424j:plain

数学パズル本に敢てプログラムで挑戦(7)

今日も今日とて、マーチン・ガードナーのMy Best Mathematical and Logic Puzzlesを眺めていたら見つけた問題。

39. The Game of Hip

6x6のチェッカーボードに先手後手が黒石、白石を置いていく。自分の石で正方形(斜めも可)を作ったら負け。この時、引き分けになる置き方は存在するか?

このゲーム何処かで見かけたことあるな・・・。いずれにせよ本によると、マーチン・ガードナーは当初存在しない、勝ちか負けになると予想していたらしいが、どこかの大学の数学科の生徒が引き分けになる最終図を見つけたとのこと。つまり、6x6が黒石18個、白石18個で埋められていて且つそれぞれの色で正方形が作られないパターンだ。36箇所から18箇所選ぶ組み合わせの数は90億なので、それぞれが正方形を含まない(6x6上の正方形のパターンは105通り)をチェックすると相当無理、ここでは小さい場合について解き、その結果からn=6を構築した。

まず、6x6の正解のパターンがあるとすると、その左上、4x4の部分チェッカーボードも正方形を含まない。同様に、右上、左下、右下も4x4の正方形ができないパターンになっている。そこでまず、4x4の正方形が出来ないパターンを全て求める。これは5,006通りあった。その中から一つ左上の候補を固定する。右上の部分の4x4だが重なっている部分は左上と同じパターンにならないといけないので5,006からかなり少なくなる。20通りくらい。同様に左下の候補も20通りほど。右下は右上と左下の4x4と重なっているので更に候補が絞れる、高々数通りだ。このようにして候補の6x6パターンを作る。当然黒石が18個のチェックも忘れずに。

次に、そのパターンが正方形を含むかどうかだが、正方形を含まない4x4から作ってあるので6x6上でも幾つかの正方形構成条件は満たさないことが分かっている。結局チェックする条件は105から減って36通り。ここまでやってあとはプログラムをぶん回した。結果、数十分かかった・・・。最終結果は本に載っていた解答も出てきたので問題なさそうだ。それにしても、回転、反転、bit反転含めてで24通りしかないんだから、マーチン・ガードナーが見つけられなくてもしようがない。目視で確認してみたけど人間は斜めの正方形見つけるのは不可能だね(笑)

プログラムを改良して速くしたいんだけど3x3に分割して合成するくらいしか思いつかない。それも速いかどうか定かでないし。結局これって36変数の210節の4-SATを解いているわけだし、構造になんか良い性質ないと難しそう。

#generate square position on nxn board
def gene_square_pos(n)
    ret=[]
    r=0...n
    n.times{|x|1.upto(n-1){|y|
        ul=[x,0]
        ur=[0,y]
        next unless r===x+y
        dl=[x+y,x]
        dr=[y,x+y]
        v=[ul,ur,dl,dr]
        d=n-x-y
        d.times{|dx|d.times{|dy|ret << v.map{|_x,_y|[_x+dx,_y+dy]}}}
    }}
    ret
end

def anti_pos(pos)
    pos.map{|p4|p4.map{|p|!p}}
end

def move_square_pos(pos,dx,dy)
    pos.map{|p4|p4.map{|x,y|[x+dx,y+dy]}}
end

#turn some portion of position into code
def code(pos,xrange,yrange)
    ret=0
    xrange.each{|x|yrange.each{|y|ret=2*ret+(pos[x][y] ? 1 : 0)}}
    ret
end

def lcode(pos)
    n=pos.size
    code(pos,0...n,0...n/2)
end

def rcode(pos)
    n=pos.size
    code(pos,0...n,n/2...n)
end

def ucode(pos)
    n=pos.size
    code(pos,0...n/2,0...n)
end

def dcode(pos)
    n=pos.size
    code(pos,n/2...n,0...n)
end

#check whether the given position make a square?
def pos_make_square?(pos,sq_pos)
    sq_pos.any?{|pts|
        pts.all?{|p|pos.dig(*p)}||pts.all?{|p|!pos.dig(*p)}
    }
end

def display(pos)
    pos.each{|ary|puts ary.map{|v|v ? ?# : ?.}*""}
end

#doable up to 4
def gene_nxn_pos(n)
    ret=[]
    sq_pos=gene_square_pos(n)
    1.upto(n*n/2){|k|
        (n*n).times.map{|i|i.divmod n}.combination(k){|ary|
            pos=Array.new(n){Array.new(n,false)}
            ary.each{|x,y|pos[x][y]=true}
            next if pos_make_square?(pos,sq_pos)
            ret << pos
            ret << anti_pos(pos) unless 2*k==n*n
        }
    }
    ret
end

pos4x4=gene_nxn_pos(4)
sq_pos4x4=gene_square_pos(4)
sq_pos6x6=gene_square_pos(6)
pos6x6=[]
u={}
l={}
pos4x4.each{|pos|
    uc,lc=ucode(pos),lcode(pos)
    [[u,uc],[l,lc]].each{|h,c|
        if h[c]
            h[c] << pos
        else
            h[c]=[pos]
        end
    }
}

#eliminate the condition already examined
sq_pos6x6-=sq_pos4x4
sq_pos6x6-=move_square_pos(sq_pos4x4,0,2)
sq_pos6x6-=move_square_pos(sq_pos4x4,2,0)
sq_pos6x6-=move_square_pos(sq_pos4x4,2,2)

pos4x4.each_with_index{|ulpos,idx|
    p idx
    rc=rcode(ulpos)
    dc=dcode(ulpos)
    l[rc].each{|urpos|u[dc].each{|dlpos|
        dl_cand1=l[rcode(dlpos)]
        next unless dl_cand1
        dl_cand2=u[dcode(urpos)]
        next unless dl_cand2
        dl_cand=dl_cand1&dl_cand2
        next if dl_cand.empty?
        dl_cand.each{|drpos|
            pos=Array.new(6){Array.new(6,false)}
            #copy
            num=0
            4.times{|x|4.times{|y|
                if ulpos[x][y]
                    pos[x][y]=true
                    num+=1
                end
            }}
            4.times{|x|4.upto(5){|y|
                if urpos[x][y-2]
                    pos[x][y]=true
                    num+=1
                end
            }}
            4.times{|y|4.upto(5){|x|
                if dlpos[x-2][y]
                    pos[x][y]=true
                    num+=1
                end
            }}
            4.upto(5){|x|4.upto(5){|y|
                if drpos[x-2][y-2]
                    pos[x][y]=true
                    num+=1
                end
            }}
            next unless num==18
            #square check
            unless pos_make_square?(pos,sq_pos6x6)
                pos6x6 << pos
                display(pos)
                puts
            end
        }
    }}
}

puts pos6x6.size # => 24

My Best Mathematical and Logic Puzzles (Dover Recreational Math)

My Best Mathematical and Logic Puzzles (Dover Recreational Math)

続きを読む

数学パズル本に敢てプログラムで挑戦(6)

マーチン・ガードナーのサム・ロイド本を眺めていたら見つけた問題。解けそうだったので頭の中で組み合わせ考えたけど結局解けなかった。悔しいからRuby書いてしまいました(笑)

56. How can you score exactly 50 points?

ようは、与えられた数の集合{25, 27, 3, 12, 6, 15, 9, 30, 21, 19}から選んで合計が50にできるかというもの。O(n*k)で解けますね。っていうか、解けないから解無しだと思ったんですけどね。こういうテキ屋の景品って当たり無しが当然じゃないですか。 ちなみに、このパズル作家サム・ロイドですが(音がメトロイドっぽいよね(笑))相当古いみたいでpublic domain入りしてましたので、実際の紙面(MathPuzzle.comより)も載せておきます。マーチン・ガードナーの本と結構説明が違う感じですね。差別的な表現が抑えられているような。まあ、合法的にオリジナル本を参照できるなら、マーチン・ガードナー本を買う必要はあまりないかな。

points=[25, 27, 3, 12, 6, 15, 9, 30, 21, 19].sort

N = 50

dp = Array.new(N + 1, false)
dp[0] = []

points.each{|p|
    (N + 1).downto(0){|i|
        next if i + p > N
        if dp[i]
            dp[i + p] = dp[i].dup
            dp[i + p] << p
        end
    }
}

p dp[N] # => [6, 19, 25]

f:id:uru:20170429002831j:plain

Mathematical Puzzles of Sam Loyd (Dover Recreational Math)

Mathematical Puzzles of Sam Loyd (Dover Recreational Math)

マジ数学本に敢てプログラムで挑戦(1):包除原理(Inclusion-Exclusion Principle)

SpringerのGTMと言えば黄色のカバーで有名ですが、238 A Course in Enumerationに包除原理を使った面白い問題が

Ex. 5.7.

How many integer solutions \( x_1+x_2+x_3+x_4=30 \) exist with the restriction -10 \le x_{i} \le 20 ?

上限が無ければよくある問題でcombinationで計算できます。

O(n4)

数字が小さいので四重ループで解けます、がO(n4)なのでNが大きいとお手上げです。

N=30
L=-10
R=20

def f_bf
    ret=0
    L.upto(R){|x|
        L.upto(R){|y|
            L.upto(R){|z|
                L.upto(R){|w|
                    ret+=1 if x+y+z+w==N
        }}}}
    ret
end

O(n2)

次に考えるのは4変数一括ではなく2変数ずつで考える。具体的には\( x_1+x_2 \) の取りうる組み合わせ計算しておき、 x_1+x_2=a となる個数と x_3+x_4=30-a を満たす個数をかけ合わせて行く方法。これはO(n2)で計算できる。

def f_dc
    memo=Array.new((R-L)*2+N,0)
    L.upto(R){|x|
        (x+1).upto(R){|y|
            memo[x+y]+=2
    }}
    L.upto(R){|x|memo[2*x]+=1}
    ret=0
    (2*L).upto(2*R){|i|
        ret+=memo[i]*memo[N-i]
    }
    ret
end

O(1)

最後の方法が包除原理を使ったもので、これが一番早い。まず初めに、次の事実に注意:

\( x_1+x_2+x_3+x_4=30 \)を満たす正整数の組み合わせ個数は  _{29} C_3 である

これは、ボールを30個横に並べて3つの衝立で分割すること(◯◯|◯・・◯◯|◯・・◯|◯◯◯)に方程式の解が一対一に対応し、衝立を入れる場合の数は29箇所から3つ選ぶことになるから。次に\( y_i=x_i+10 \) と変形することで同値な問題

How many integer solutions \( (\ast) \ y_1+y_2+y_3+y_4=70 \) exist with the restriction (\star) \ 0 \le y_{i} \le 30 ?

にできます。すると求めるのは包除原理より、

(\ast) を満たす正整数の組み合わせ数-(\ast) を満たすけど1変数だけ(\star)を満たさない+(\ast) を満たすけど2変数だけ(\star)満たさない・・・

となりますね。j個の変数だけ(\star)を満たさないということはj個の変数の値は31以上。つまりj個の変数を\( z_i=y_i-31 \) とすることで元と同じ問題、ただし総和が70-31*jを解くことに帰着されます。4変数のうちj個を選ぶのは _4 C_j ですから結局次のようなアルゴリズムになります。

def f_inc_exc
    n=4
    k=N-L*n
    s=R-L+1
    #new problem: find the # of sol x+y+z+w=70 with the restriction 0<= each<=30?
    ret=0
    0.upto(n){|j|
        sgn=j.odd? ? -1 : 1
        temp=sgn*combination(n,j)*combination(n+k-j*s-1,n-1)
        ret+=temp
    }
    ret
end

小さいNでは差が分かりませんが大きいと違いは明らか。ブルートフォース版は10分待っても帰ってきませんでした。

N=10000
L=-5000
R=7000

       user     system      total        real
inclusion & exclusion : 828252025001
  0.000000   0.000000   0.000000 (  0.000093)
divide & conquer      : 828252025001
  7.710000   0.020000   7.730000 (  7.735661)
bruteforce            : ^Cex05_07.rb:12:in `block (4 levels) in f_bf': Interrupt

A Course in Enumeration (Graduate Texts in Mathematics)

A Course in Enumeration (Graduate Texts in Mathematics)

数学パズル本に敢てプログラムで挑戦(5)

チェス盤のタイリング問題に挑戦。この問題はALGORITHMIC PUZZLESではチュートリアルの例題として載っていた。 まあ有名ですからね。

Domino Tiling of Deficient Chessboards

8x8のチェス盤をドミノ(1*2)で隙間なく埋めることを考える。a)左上が欠けているチェス盤で可能か? b)左上、右下が欠けている場合は?

a)は偶奇性が異なる(空きは偶数である必要がある)からNGは明らか。b)はちょっとトリッキーだが、チェス盤が市松模様なのとドミノをどう置いても2色をカバーする必要があることからチェス盤の2色は同じセル数ある必要があるからNG。 念のため、深さ優先でプログラムを組んでドミノを置いてみたが、理屈通りできませんでした(笑) まあ、ロジックで出来ないを判定できるのはすごいことだ。

#domino tiling on N*N board

N = 8
board = Array.new(N + 1){Array.new(N + 1, false)}
(N + 1).times{|j|board[N][j] = true} #sentinel
(N + 1).times{|i|board[i][N] = true}
board[0][0] = true #deficient cell
board[N - 1][N - 1] = true

def find_next(b, x, y, n)
    x.upto(N){|i|return [i, y] unless b[i][y]}
    (y + 1).upto(n - 1){|j|N.times{|i|return [i, j] unless b[i][j]}}
    nil
end

def f(b)
    def dfs(b, x, y, rest, n)
        ret = 0

        #horizontal
        if b[x][y + 1] == false
            return 1 if rest == 2
            b[x][y] = b[x][y + 1] = true
            nx, ny = find_next(b, x, y, n)
            if nx
                ret += dfs(b, nx, ny,rest - 2, n)
            end
            b[x][y] = b[x][y + 1] = false
        end

        #vertical
        if b[x + 1][y] == false
            return 1 if rest == 2
            b[x][y] = b[x + 1][y] = true
            nx, ny = find_next(b, x, y, n)
            if nx
                ret += dfs(b, nx, ny, rest - 2, n)
            end
            b[x][y] = b[x + 1][y] = false
        end
        ret
    end

    #body
    n = b.size - 1
    rest = 0
    n.times{|i|n.times{|j|rest += 1 unless b[i][j]}}
    x,y = find_next(b, 0, 0, n)
    dfs(b, x, y, rest, n)
end

p f(board) # => 0

Cracking the Code Interview 16.19 Pond Sizes

もう1題は深さ優先探索の練習問題。長方形の区画の高さ情報(0は水面)が与えられるので全ての池を求めるというもの。池は上下左右斜めで繋がっているという定義。深さ優先探索再帰で書くと簡単だけど、たまにスタックがいっぱいになったりするので注意が必要。いつも思うのは、チェック先セルが範囲内かをチェックするのにRangeを使っているけど遅い気がするんだよなあ。不等号の方が速いはず。Ruby配列には範囲外を指定するとnilが返る仕様なんだけど、多次元だとnil[4]とかでmethod missingエラーになって処理が美しくないんだよね。Nilクラスにメソッド追加してnil返すようにすればいいんだけど。まあ、周りを壁で囲えば気にしなくていいんだけど、配列操作がRubyは非常に面倒くさいのでやりたくないしindex -1で最後にアクセスできるので、与えられる配列の最後に2つ番兵入れるだけだけどね。不等号と言えば、Enumerable#minも2値の場合は遅いと思う。ループの最深では不等号にしないとTLE。せっかく短く書けるからRuby使っているのに、冗長な書き方しないといけないのは本末転倒なんだけど。時間制限がある場合はしようがない。

def pond_sizes(pond)
    def dfs(x,y,pond)
        ret = 1
        pond[x][y]=-1
        [[1,0],[0,1],[-1,0],[0,-1],[1,1],[1,-1],[-1,1],[-1,-1]].each{|dx,dy|
            nx,ny = x+dx, y+dy
            next unless (0...pond.size)===nx
            next unless (0...pond.first.size)===ny
            ret += dfs(nx,ny,pond) if pond[nx][ny]==0
        }
        ret
    end

    h=pond.size
    w=pond.first.size
    ret=[]
    h.times{|i|w.times{|j|ret << dfs(i,j,pond) if pond[i][j]==0}}
    ret
end

if __FILE__==$0
    #text example
    pond=[[0,2,1,0],[0,1,0,1],[1,1,0,1],[0,1,0,1]]
    p pond_sizes(pond) #=> [2,4,1]
    
    pond =[[0]]
    p pond_sizes(pond) #=> [1]

    pond =[[1,0]*50,[0,1]*50]
    p pond_sizes(pond) #=> [100]
end

Cracking the Code Interview 17.21 Volume Of Histogram

「世界で闘うプログラミング力を鍛える本」の最新版の翻訳が出たとのことで暇つぶしに眺めていたら面白い問題を見つけた。ヒストグラム地形が与えられるので雨が降った時に貯まる水の量を求めるというもの。例えば、3,1,5 という入力だと、3と5の間に高さ2の水たまりが出来るから2を出力すると言った具合。2,0,1,0,3だと5、基本0以上の数値が与えられる。 山頂を境として右と左に水たまりが出来ることは明らかだ。漢字の「山」という字の右左のくぼみに水が貯まるイメージ。つまり、左から山頂までの水たまりの体積が分かったら、今度は右から同じことをすればOK。従って、まず山頂を求め、そこで左右に分割し、右は反転し、vol_of_slopeメソッドにあとは計算を任せる。vol_of_slopeは一番最後に最高峰の壁が出てくることになっているので、配列を読みながら逐次現在までの最高点の更新と、最高点以下の地形が出てきたらその差分が水で埋まるので加えていく。

def volume_of_histgram(a)
    max_hight = max_idx = 0
    a.each_with_index{|h, idx|
        if h > max_hight
            max_hight = h
            max_idx = idx
        end
    }
    vol_of_slope(a[0..max_idx]) + vol_of_slope(a[max_idx..-1].reverse)
end

#the rightmost must be the highest
def vol_of_slope(a)
    ret = 0
    prev_highest = 0
    a.each{|h|
        if h > prev_highest
            prev_highest = h
        elsif h < prev_highest
            ret += prev_highest - h
        end
    }
    ret
end