PSM

Pプログラミング S初心者の Mメモ書き

Python3で組み合わせ計算を高速化

初めに

Pythonで組み合わせ計算 {}_n\mathrm{C}_rを計算しようとすると,少々問題が発生します.
nやrが小さいうちはいいんですが,これが大きくなってくると計算が非常に長くなってくるのです.
実務や研究では少々計算が遅くても別に問題はないと思うのですが,競技プログラミングだとそれが命取りになってくるんですね……

組み合わせ計算のための元(と逆元)のテーブルの計算

階乗をメモしておく

というわけで,計算が重いなら予め計算した結果をテーブルにしといちゃいましょうというテクニックがあるそうです.

さて
 {}_n\mathrm{C}_r = \frac{n!}{r!(n-r)!}
からわかるように,単純に考えれば元 x! 1 \le x \le N(Nは考えられる最大)で計算してリストに書き込んでおけばよいことになります.
よって

#元
N  = 10
g =  [math.factorial(i) for i in range(N) ]

のように求めておけば,例えばあとで10!の値が欲しければg[10]を呼び出すことにより,O(1)で組み合わせ計算を行うことができます.

上のテーブルを使えば組み合わせ計算は

def combination(n,r):
    return int(g[n]/(g[r]*g[n-r]))

のように求めることができます.


出力できる数値に最大値の制限がある場合(競プロ特有の話)

上のようにテーブルに記憶させた場合でも,nやrが大きくなってくると,そもそもの階乗の計算量が膨大になります.

ここで,出力に制限がある状態にはさらに高速にテーブルを作成することができるようです.
具体的には「出力は 10^9+7で割った余りとせよ」のような場合のことです.


 {}_n\mathrm{C}_r = \frac{n!}{r!(n-r)!}
について,元 x!と逆元 (x!)^{-1} 1 \le x \le N(Nは考えられる最大)で計算しておくことを考えます.

 10^9+7素数なので,フェルマーの小定理より
 x ≡ x^{10^9 + 7} \ \ (\mathrm{mod} \ 10^9 + 7)
両辺を x^2で割って
 x^{-1} ≡ x^{10^9 + 5} \ \ (\mathrm{mod} \ 10^9 + 7)
となっています.

元と逆元のテーブルは,以下のように作成し,使うことができるようです.
(下記のコードの元は自分で考えたものではありません)

mod = 10**9+7  #出力の制限
N = 10**4 
g1 = [1, 1]  # 元テーブル
g2 = [1, 1]  #逆元テーブル
inverse = [0, 1]  #逆元テーブル計算用テーブル

for i in range( 2, N + 1 ):
    g1.append( ( g1[-1] * i ) % mod )
    inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
    g2.append( (g2[-1] * inverse[-1]) % mod )

#以下のcom関数を用いて組み合わせを計算
def com(n, r, mod):
    if ( r<0 or r>n ):
        return 0
    r = min(r, n-r)
    return g1[n] * g2[r] * g2[n-r] % mod

関連する定理など

二分累乗法

xをn乗したいとき,普通に計算するとn回累乗計算を行うことになるので,計算量は \mathrm{O}(n)になります.
これを高速にする二分累乗法というものがあり,Python3で以下のように実装できます.

# 二分累乗法
def power(x, n):
    y=1
    while( n>0 ):
        if ( n%2==0 ):
            x = x * x
            n = n/2
        else:
            y = y * x
            n = n-1
    return y

nを二進法展開するように2でどんどん割っています.
1.奇数のときは答えにxを掛けてnから1引く
2.偶数のときは答えはそのままにxを2乗してnを1/2にする
という操作を繰り返すのですが,最低でも2回に1回はnが1/2になるので計算量が  \mathrm{O}(\log_2 n)に減少するというわけです.

より数学的な(でも簡略な)説明は以下.

nを二進法で展開すると
 n  = 2^{s_1} + 2^{s_2} + \cdots + 2^{s_k}  ただし   0 \le s_1 < s_2 < \cdots < s_k
となる.これより x^n
 \begin{split} x^n  &= x^{2^{s_1} + 2^{s_2} + \cdots + 2^{s_k}} \\ &= x^{2^{s_1}}x^{2^{s_2}} \cdots x^{2^{s_k}} \end{split}
となる.このとき,各々の x^{2^{s_i}} x^{2^{s_{i-1}}}から順に計算していけるので,計算量は  \mathrm{O}(\log_2 n)になる.

フェルマーの小定理

wikipedia:フェルマーの小定理

p を素数とし、a を p の倍数でない整数(a と p は互いに素)とするときに、
 a^{p-1} ≡ 1 \ \ (\mathrm{mod} \ p)
すなわち、a の p − 1 乗を p で割った余りは 1 であるというもの

Python3で最小木問題とプリム法

最小木問題を解くアルゴリズムのうち,プリム法をPython3で実装したのでメモ.
プリム法の実装の理解だけでなく,matplotlibを使ったので良い練習になりました.

wikipedia:プリム法

import matplotlib.pyplot as plt
import random
import numpy as np
import math
import sys

#最小木問題をプリム法で解く

#点群をランダム発生
p = [[random.uniform(0,100), random.uniform(0,100)] for i in range(100)]

#空のリストを用意
e = []  #部分最小木
av = []  #部分最小木に含まれる頂点
pv = list(p)  #その他の頂点

#始点は任意
av.append(p[0])
pv.remove(p[0])

while 1:
    #比較のための初期値として十分大きい数を用意
    shortestlen = sys.float_info.max
    #部分最小木に追加する辺を為す2点の初期化
    kav = None  #部分最小木に含まれる頂点からの候補
    kpv = None  #その他の頂点からの候補
    for v1 in av:
        for v2 in pv:
            #重みの計算(ここでは距離)
            dist = math.sqrt(  (v1[0]-v2[0])**2 + ( v1[1]-v2[1] )**2 )
            if shortestlen > dist:
                shortestlen = dist
                #候補の更新
                kav = v1
                kpv = v2
    if shortestlen > 0:
        #部分最小木に追加
        e.append([kav , kpv])
    #確定した頂点を移す
    av.append(kpv)
    pv.remove(kpv)
    #全頂点に対する判定を終えたら終了
    if len(pv) == 0:
        break


# 最小木の描画
for l in e:
    l = np.array(l).T
    plt.plot(l[0], l[1], 'k-', lw=2)
plt.show()


結果は
f:id:muromura:20170701180710p:plain
となりました.

Python3で任意長のリストを取ろうと思ったらいつのまにかリスト内包表記を勉強していた

Python3で任意長のリストを宣言しておくのってどうしたらいいんだと思って調べていたら,いつの間にかリスト内包表記を勉強していました.

これを使うと,
1.とりあえずfor文で回した結果を片っ端からappendなどでリストに突っ込む
      ↓↓↓
2.要素1つ1つ判定しながらリストを加工(削除,要素の変更など)する
というような操作が1行で書けます.
正確に言うと,先に判定をして加工も済んだ要素で,リストをつくるという感じのようです.

任意長の配列

[0, 0, 0, 0, 0, 0, 0]

例えば 上のような長さ7のリストが欲しい時は,

[0]*7

と書けばよいですが,リスト内包表記という”[]の中にfor文を書く”やり方でも作ることができます.

 [0 for i in range(7)]

リスト内包表記

基本形

 [i for i in range(7)]

出力は

[0, 1, 2, 3, 4, 5, 6]

一番左のiがリストの要素となっています.
この場合では,右側のfor文でiが0から1つずつ増えていくと定義されていることになります.
よって例えば以下のように一番左のiを(i+1)に変えると,

 [(i+1) for i in range(7)]

出力は

[1, 2, 3, 4, 5, 6, 7]

となります.

if(後置)

要素をifで判定してリストに加えるかどうか決めることもできます.
この場合,for文のあとにif文を書きます.

 [i for i in range(7) if i<5]

5より小さいiのみが要素として追加され,

[0, 1, 2, 3, 4]

となります.

if else

if elseを使いたい時は最初にもってくる必要があります.

 [i if i<5 else 4 for i in range(7)]

iが5より小さい時はそのままiが要素となり,iが5以上のときは4が要素となります.

[0, 1, 2, 3, 4, 4, 4]

集計

リストなどで集計関数がそのまま使えます.例えば,

sum( [(i+1) for i in range(7)])

について,全要素の和となるので

28

が得られます.28は完全数ですね.

2重for

forを2回重ねることもできます.

[[i,j] for i in range(2) for j in range(5)]

この場合はiが2通り,jが5通り動くので,長さ10のリストができます.

[[0, 0],
 [0, 1],
 [0, 2],
 [0, 3],
 [0, 4],
 [1, 0],
 [1, 1],
 [1, 2],
 [1, 3],
 [1, 4]]

Python3で二分探索

wikipediaに載っているc++スクリプトpythonで書き直しただけですが,”忘れた時にとりあえずコピペしたら動く”がモットーなので自分用にメモ.

[wiki:二分探索]

実装は以下.

#二分探索
def binarySearch(ary, key, imin, imax):
    if ( imax < imin ):
        return "KEY_NOT_FOUND"
    else:
        imid = imin + (imax - imin)//2
        if (ary[imid] > key):
            return binarySearch(ary, key, imin, imid - 1)
        elif (ary[imid] < key):
            return binarySearch(ary, key, imid + 1, imax)
        else:
            return imid

試してみたらちゃんと動きました.

a = [1,2,2,4,5,6,7,10,20,40]
binarySearch(a, 40, 0, len(a)-1)

再帰ってプログラミング初心者が苦しむ難関の1つですよね……

Pythonで2進法方のリストを作る

Pythonで二進法の各桁を各要素とするようなリストをつくって,降順に出力することを考えました.
その際,文字列の各要素のリスト化や文字列のゼロパディングなど,いくつか覚えておきたい操作があったのでメモしておきます.

操作

数値型のリストのソート

a = [0,3,4,2,1]

# 昇順
sorted(a)
#降順
sorted(a, reverse = True)

文字列の各文字を要素に持つリストをつくる

s  = "abcde"

#リスト化
list(a)

これを応用して,何桁かの数字を各桁を要素に持つリストにバラしたい場合は

a = 12345

# リスト化
map(int, list(str(a)))

二進数に変換

a = 31
format(a, 'b')

出力はstr型になるのに注意が必要です.

文字列の0埋め

"20"
"123"
"02"

#4桁にそろえる
"0020"
"0123"
"0002"

にしたいという話です.

a = "20"
a.zfill(4)

サンプルスクリプト

l  = 4
f = 2**l-1
for i in range(2**l):
    ff = list(map(int,list(format(f,'b').zfill(l))))
    print(ff)
    f-=1
[1, 1, 1, 1]
[1, 1, 1, 0]
[1, 1, 0, 1]
[1, 1, 0, 0]
[1, 0, 1, 1]
[1, 0, 1, 0]
[1, 0, 0, 1]
[1, 0, 0, 0]
[0, 1, 1, 1]
[0, 1, 1, 0]
[0, 1, 0, 1]
[0, 1, 0, 0]
[0, 0, 1, 1]
[0, 0, 1, 0]
[0, 0, 0, 1]
[0, 0, 0, 0]

こんな感じにでてきます.
なんか無駄が多い気がしますが……


あとがき

競技プログラミング初心者(2日)です.
あくまで練習のために始めたので,競プロに取り組む上ではどうも邪道らしいんですが,Pythonでやっていきたいと思います.

今回は上で作ったリストをフラグとして使ったのですが,結局,TLEが出てしまって方向転換を余儀なくされました.

python3での標準入出力

pythonでの標準入出力を求められる機会って稀によくありますよね.
pythonって2だったり3だったり紛らわしいし,調べてみてもコピペで動いてくれるスクリプトが載っている記事を見つけられなかったので,自分の検索能力のなさを嘆きつつ苦しんだ結果をメモ書き程度に残しておきます.
(環境はpython3.6.0を想定しています)

入力

1行に1つの入力を受け取る

28

とか

test

とかを受け取りたい時です.
たぶんどの言語でもそうなんでしょうけど,pythonでも入力関数を一回呼び出すごとに1行分の入力を読み込んでくれるようです.
というわけで,

a = input()

とすればaに入力を代入できます.

ただしこれだとintだろうとfloatだろうとstrだろうと何でもstr型で読み込んでしまうので,型変換までついでにやっちゃうなら

a = int(input())

のように書けます.

1行に複数の入力を受け取る

a b c d e

1 2 3 4 5

とかを受け取りたい時です.

文字列で受け取るなら

a = input().split()

とすると,半角スペースのところで区切ってリストにして受け取ってくれます.

型を指定して読み込むときには,変数が1つのときの場合のようには書けず,map関数を使う必要があります.
intで受け取りたいなら

a = list(map(int, input().split()))

となります.
注意点としては,python3ではmap()を使って読み込むとmap objectというものになってしまうので,list()で変換する必要があります.

また,受け取りたい入力が数個程度で,その数がわかっているならこんな風にもできました
たとえば

3 4

に対して

a,b = list(map(int, input().split()))

とすればaに3が,bに4がint型で代入されます.


複数行に複数の入力を受け取る

3
1 2 3
4 5 6
7 8 9

競プロでよくある,みたいなのを受け取りたい時です.

input()を一回呼び出すごとに1行分の入力を受け取ってくれるので,for文などで行数分呼び出すか,whileなどを使えば良いようです.
というわけで,

n = int(input())
a = []
for i in range(n):
    a.append(list(map(int, input().split())))

最初に空のリストを用意して,appendでどんどん突っ込んでいくという方針です.
numpyを使うならここでついでにnp.arrayにしておいても良いかもしれません.


出力

基本

python3では標準出力printが関数となっていました(python2系との違い).

print("output")

のように書きます.
勘違いしていたのですが,ここではint型も勝手に型変換して出力してくれるようなので,例えば

print(1234)

とint型を直接書いても,勝手に変換して"1234"と出力してくれます.
複数の引数をとることもできて,

n = 2
print("春が",n, "階から降ってきた")

とすれば

春が 2 階から降ってきた

リストを空白などの文字区切りで出力する

join()を用いると,リストの要素を任意の文字を糊にしてくっつけた文字列に変換することができます.
["君","羊","青"]を"と"でくっつくて"君と羊と青"にするイメージです.

s = ["君","羊","青"]
joined_s = "と".join(map(str, s))
print(joined_s)

とすれば

君と羊と青

"と"の部分を半角スペース" "などにすれば半角スペース区切りで出力できます.

プログラミング初心者がUnityチュートリアルで苦しんだポイント

1 はじめに

プログラミング初心者がUnityのチュートリアルに挑戦してつまづいた部分をメモしていきます.
(逐次更新する予定です)

「初心者ゆうてもピンキリやがな」という話ですが,私がどのくらい初心者かというと,for文やif文くらいは自信を持ってわかると言えるレベルです.
クラス?宣言?オブジェクト指向
なんとなく説明できないこともないこともないこともないよ!ってなものです.

さて.
現在のUnityのバージョン(5.6.1)に対してチュートリアルで使われているバージョンが若干古く,UIなどの仕様変更が多々あっていちいちつまづきました.
プログラミングを勉強し始めてから,
あっこれ死にゲーってやつだ
と割と早く悟ることができたのですが,盲点なのは自分の記憶力が恐ろしく悪かったことです.
いつまでも序盤の罠でしに続けているわけにはいかないので,主に自分の為にそれらをメモしていこうと思います.

実際に作ったものや感想についてはまた別の記事で.


2 つまづいた部分

2.1 Assetについて

Asset Tabを開き,導入したいassetを検索してインポートします.
今回はspace shooterで検索するとチュートリアル用の一式を手に入れることができます.
Asset storeのタブは,ショートカット⌘9か,Window > Asset Storeで開くことができます.


2.2 Unity Web PlayerとWebGL

チュートリアルの動画ではWeb Playerが紹介されていますが,現在のUnityでは非推奨としており,どうもそのうち使えなくなるようです.
そこで代わりにWebGLを使用することにしました.
作業環境や対象とする環境は,⌘+Shift+Bでビルトの設定画面を開き,そこで設定することができます.


2.3 Gameウィンドウの設定

ビルトの設定画面で左下のplayer settingsをクリックすると,Inspectorタブで設定ができるようになります.
ここではResolution and Presantationのところでサイズを変更し,そのあとGameタブの左上のFree Aspectとなっているところをクリックして+を選択し,新しく600×900に設定し直しました.


2.4 スクリプトエディターについて

Unityをインストールすると勝手についてくるMonoDevelopを使えばとりあえずスクリプトは書けます.
サジェスト機能など,使い心地などがちょっとアレな場合はVisual Studioなどを使うこともできます.
www.visualstudio.com


2.5 RIgidbodyの取得

Unity公式のチュートリアル動画では

rigidbody.velocity = movement

のように,rigidbodyのvelocityメソッドに直接movementの情報を与えられると書いていますが,これはUnity4のための書き方であって,現在のUnity5ではエラーとなってしまいます.
Unity5では,まずvoid FixedUpdate()の前にvoid Start()を宣言し,その中でGetComponent()を使ってrigidbodyを取得してやる必要があるようです.

つまりスクリプト

void Start()
{
	rb = GetComponent<Rigidbody> ();
}

void FixedUpdate()
{
	float moveHorizontal = Input.GetAxis ("Horizontal");
	float moveVertical = Input.GetAxis ("Vertical");

	Vector3 movement = new Vector3 (moveHorizontal, 0.0f, moveVertical);
	rb.velocity = movement;
}

と書き換えれば大丈夫です.


2.6 Mathf

チュートリアルの中で,公式にリファレンスがあるので逐一チェックすると良いよと言われました.
試しにMathfで検索してみたら公式の(しかも日本語の!)リファレンスがヒットしました.
docs.unity3d.com
これによるとMathfはUnityEngineの中に入っている,標準的な数学関数を扱うためのモジュールのようです.


2.7 public float speedなど

public class PlayerController : MonoBehaviour
{
	public float speed;
	public float xMin, xMax, zMin, zMax;
}

のようにクラスの先頭で[public+型+変数名]で宣言した変数は,UnityのInspectorで後から与えることができます.
こんな感じに新しい入力窓が勝手にできます.
f:id:muromura:20170612194751p:plain


また詳しいことはわかりませんが,class宣言時に[System.Serializable]を使うと,以下のように折り畳んだ形で変数を与えられるようです.
f:id:muromura:20170612195538p:plain





3 まとめ

そもそもCとC++C#の違いすらわかっていない人間なので,コードの意味もなんとなくしか追えず,ひたすら辛い時間が続いています.
しかしその辛さを乗り越えて余りあるのが,ゲームを作る楽しさです.

とりあえずチュートリアルの範疇では面倒なキャラクターの3Dモデリングやテクスチャー作成をする必要がないのも,素晴らしい点です.
システムとしてはそんな複雑なことをやっていなくても,とりあえずモデリングとテクスチャーがしっかりしていたら見れるものになりますもんね.

それでは.




4 この記事を書くにあたって参考にさせていただいた記事

tech-camp.in
idol-ch.net
docs.unity3d.com
brian.hatenablog.jp
curious-boys.hatenablog.com