• 変数と型について
• for文による繰り返し処理
• if文による条件分岐

Pythonに限らず、プログラミング言語には「変数(Variable)」という概念が出てくる。変数は簡単なように見えて、意外に理解が難しい概念である。例え話を多用するととっつきやすが、後で混乱するもとなので、ここでは変数の実装に近いところから説明する。ただし、簡略化してあり、実装そのものではないので注意して欲しい。

計算機とは、メモリにあるデータをCPUで処理して、またメモリに書き戻す機械である。メモリには「番地」という通し番号がついており、例えば足し算なら「0番地と1番地のデータを読み込んで足してから2番地に書き込め」といったことを指示する。しかし、いちいちどの値がどこにあるか番地で覚えるのは面倒だ。そこで変数というラベルを使うことにする。

以下のプログラムを見てみよう。

a = 10

これは、aという変数を用意し、そこに10という値を代入せよ、という意味だ。実際には、メモリ上にaというラベルを貼り、そこに10を書き込む、という処理をする。同様にb = 20とすると、bというラベルが作られ、そこに20が書き込まれる。

a = 10
b = 20

printという命令を使うと、その値を確認することができる。

print(a) # => 10

今後、# =>という記号は、「左を実行すると、右の結果が得られるよ」という意味だと約束する。

また、Jupyter Notebookでは、変数名だけを含む、もしくはセルの最後に変数名だけを含むセルを実行すると、その値を確認することができる。

a #=> 10

さて、この状態でc = a + bを実行してみよう。まず、右のa+bが評価される。その結果得られた値30をcの指す場所に書きこもうとするが、まだcという変数は作られていないので、まずcの場所が作られてからそこに値が書き込まれる。

同様に、a = a + bを実行した場合、aの場所はすでに存在するので、その値が更新される。

このようにPythonでは、代入文があった時、その左辺にある変数が未定義なら作成され、定義済みなら値を更新するので覚えておこう。

また、その時に定義されていない変数を使おうとするとエラーになる。

a = 10
b = 10
a = c #=> NameError: name 'c' is not defined

上記の例では、aとbだけ定義された状態で、aにcというラベルの指す内容を代入せよと指示しているが、この時点ではcというラベルがないためにエラーになっている。

変数とは、メモリにつけられた「ラベル」である。メモリには整数しか保存できないが、プログラムでは小数点を含む数や文字列なども表現したい。そこで、「このメモリの値をどう解釈するか」を指定する必要がでてくる。これが「型(Type)」である。全ての値や変数には型がある。

例えば、1が代入された変数は「整数」の型を持っている。変数の型はtypeという命令で知ることができる。

a = 1
type(a) # => int

小数点を含む数字は浮動小数点数(float型)になる。

b = 1.0
type(b) # => float

ダブルクォーテーション" "やシングルクォーテーションマーク' 'で囲まれた文字は「文字列」として扱われる。

c = "test"
type(c) # => str

例えば文字列の場合、メモリに0x74という数値があった場合、それをアルファベット小文字のtと解釈しましょう、という約束を決めておく。変数cのstrという型は、cが指す先のメモリを「文字列として解釈する」という意味を持つ。すると、今cが指している「0x74 65 73 74」という数値列がtestという文字列として解釈される。

同様に、例えば「1.0」という浮動小数点数は、メモリ上では「0x3FF0000000000000」という8バイトの数字で表現されている。これをある約束(IEEE754という規格)に従って解釈すると「1.0」という数値となる。

整数や浮動小数点数は四則演算が可能である。同じ型同士の演算は原則としてその型になる。

type(1+2) # => int
type(1.0+2.0) # => float

ただし、整数の割り算だけ注意を要する。Pythonでは、整数同士の割り算は、たとえ割り切れる場合でも、値は浮動小数点数になる。

4 / 2 # => 2.0

割り算のあまりを切り捨てた整数値が欲しければ//と、/記号を二つ続けた演算子を用いる。

5/2 # => 2

整数型と浮動小数点数型の演算結果は浮動小数点数型になる。

print(1+1.0) # => 2.0
type(1+1.0) # => float

文字列同士は足し算ができる。

"Hello " + "World!" # => "Hello World!"

文字列と数値の演算はできない。

"1" + 2  # => Type Error

intやfloatで囲むと文字列から整数や浮動小数点数に変換できるので、演算も可能になる。

int("1") # => 1
int("1") + 2 # => 3
float("1") + 2 # => 3.0

真偽値(bool)とは「真(True)」であるか「偽(False)」であるかの二値だけを取る型で、条件分岐などで使われる。

値を比較すると真偽値になる。条件分岐や、ループの終了条件等に用いる。

1 == 1 #=> True
1 == 2 #=> False
1 != 1 #=> False
1 != 2 #=> True
1 < 2  #=> True
1 > 2  #=> False
1.0 < 2.0 #=> True

文字列の比較もできる。

"test" == "test" #=> True
"hoge" < "piyo" #=> True

真偽値は、notをつけると真偽値が逆転する。

not True # => False
not False # => True

二つの真偽値を使って論理演算もできる。andは「かつ」、orは「または」を意味する。例えば「真かつ偽」は「偽」に、「真または偽」は「真」となる。

True and False # => False
True or False # => True
not True # => False
not False # => True

「0.1」「123.45」といった、整数ではない値を表現する。整数と同様に四則演算や比較ができるし、負の数も扱える。

0.5 + 0.5 #=> 1.0
0.5 * 0.5 #=> 0.25
0.5 < 1.0 #=> True

ただし、浮動小数点数は内部的には その数値に最も近い近似値 を扱っているため、誤差が存在する。例えば0.1を3回足しても0.3にならないことに注意。

0.1 + 0.1 + 0.1 #=> 0.30000000000000004
0.1 + 0.1 == 0.2 # => True
0.1 + 0.1 + 0.1 == 0.3 #=> False

したがって、浮動小数点数同士の等号比較は信頼できない。浮動小数点数同士の等号比較は、意図通りに動作したりしなかったりする。初心者がよく入れるバグなので注意すること。

Pythonは複素数も扱うことができる。虚数単位がjなので注意。実部、虚部はreal、imagで取り出すことができる。 また、整数で記述しても浮動小数点数として扱われることに注意。複素数の宣言は1+2jのように書くか、complex(1,2)のように書く。

1 + 2j # => (1+2j)
complex(1,2) # => (1+2j)
(1+2j) + (2+4j) #=> (3+6j)
(1+2j)*(1-2j) => (5+0j)
(1+2j).real #=> 1.0
(1+2j).imag #=> 2.0

繰り返し処理は以下のように書ける。

for i in range(10):
  print(i)

これは、iの値を0から9まで変化させながら、print(i)を実行しなさい、という意味である。なお、inの前の変数はなんでも良い。例えば以下のようにしても同じ結果となる。

for j in range(10):
  print(j)

また、単に10回繰り返したいだけで、現在何回目かは不要である場合には_を使う。

for _ in range(10):
  print("Hello") # Helloが10回表示される

「もし〜なら」という処理を条件分岐と呼び、if文で書く。ifの後ろには真偽値を与えるような式を書く。

if 5 > 3:
  print("5>3")

「もし〜ならAをせよ、そうでなければBをせよ」という場合にはelse:を使う。

if 5 < 3:
  print("A") # 実行されない
else:
  print("B") # => B

「もし〜ならAをせよ、そうでない場合で〜ならBをせよ」という場合にはelifを使う。

a, b = 1, 2
if a == b:
  print("a == b")
elif a > b:
  print("a > b")
else:
  print("a < b")

if文は入れ子構造にできる。

a, b, c = 1, 2, 3
if a < b:
  if b < c:
    print("a < b < c")
  elif c < a:
    print("c < a < b")
  else:
    print( "a <= c <= b")

以上の知識で、何かコードを書いてみよう。ある方程式を解きたいが、その解が厳密にはわからないとする。その場合でも数値計算で必要な精度で解を求めることができる。そのような場合によく用いられるのがニュートン法である。

いま、

$$ f(x) = 0 $$

という方程式を解きたいとする。もし、真の解を$x$として、それに近い値$\tilde{x} = x+\epsilon$があったとする。$f(x)$を$\tilde{x}$の周りでテイラー展開すると、

$$ f(\tilde{x} - \varepsilon) = f(\tilde{x}) - \varepsilon f'(\tilde{x}) + O(\varepsilon^2) $$

$\varepsilon$の二次の項を無視した状態で、$f(\tilde{x} - \varepsilon)$が0となるように$\varepsilon$の値を選ぶと、

$$ \varepsilon = \frac{f(\tilde{x})}{f'(\tilde{x})} $$

$\tilde{x} = x + \varepsilon$であったから、$\varepsilon$を引けば、より真の値に近づくはずである。以上から、

$$ x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)} $$

という数値列を得る。この数列が収束するということは$x_{n+1} = x_n$なので、$f(x_n)=0$が満たされなければならず、それはすなわち$x_n$が解に収束したことを示す。これを確認してみよう。

いま、$x^3 = 1$の解を知りたいとする。この時、$f(x) = 0$の形に書きたいので、$f(x) = x^3 - 1$である。$f'(x) = 3 x^2$であるから、対応するニュートン法のアルゴリズムは、

$$ x_{n+1} = x_n - \frac{x^3 -1}{3 x^2} $$

である。

新しいノートブックを開き、最初のセルに以下のように入力し、実行せよ。

def newton(x):
  for _ in range(10):
    x = x - (x**3-1)/(3*x**2)
    print(x)

ここでdef newton(x):とあるのは、「newtonという関数を定義し、xという名前で値を受け入れるよ」という宣言である。関数については次回紹介するが、ここでは

• def 関数名(入力):という形で関数が宣言できる
• 定義した関数名(入力)という形で呼び出すことができる
• 関数内でreturn 値とすると、値を返すことができる

ということを覚えておけば良い。

初期値として2.0を入れてみよう。2つ目のセルに、以下のように入力、実行せよ。

newton(2.0)

以下のような実行結果が出れば成功である。

1.4166666666666665
1.1105344098423684
1.0106367684045563
1.0001115573039492
1.0000000124431812
1.0000000000000002
1.0
1.0
1.0
1.0

急速に真の値である1.0に近づいていくことがわかるだろう。実際、ニュートン法の収束は非常に早く、一度繰り返すごとに精度が倍になっていく。

$x^3 -1 = 0$の解は、実数の範囲では$x=1$しかないが、複素数まで考えれば三つ存在し、それぞれ$x=1, -1/2 \pm \sqrt{3}/2 i$である。複素数の場合にもニュートン法が機能することを確認しよう。 3つ目のセルに、以下の入力し、実行せよ。

newton(2+0j)

先程と同じ値だが、複素数として入力している。1.0に収束するが、表示が複素数になることを確認せよ。

次は-1 + 1jを入力してみよう。4つ目のセルに以下のように入力し、実行せよ。

newton(-1+1j)

$x=-1/2 + \sqrt{3}/2 i$に収束するはずである。

同様に、5つ目のセルに初期値として-1-1jを入力してみよう。

newton(-1-1j)

今度は$x=-1/2 - \sqrt{3}/2 i$に収束するはずである。

先程、ニュートン法が複素数の場合にも機能し、$x^3 - 1 = 0$の解を三つとも見つけられることを確認した。 ナイーブに考えると、複数の解がある場合、初期値に近い解に収束すると考えられる。では、複素平面のどこからスタートしたらどこに収束するだろうか？単純に三分割になるだろうか？

複素平面の様々な場所を初期として、

• 1+0jに収束したら赤
• -0.5+87jに収束したら緑
• -0.5-87jに収束したら青

に塗ることで、「どの場所からスタートしたらどこに収束するか」という「収束地図」を作ってみよう。

新しいノートブックを開き、以下のプログラムを4つのセルに分けて入力せよ。

最初のセルで、必要なライブラリをインポートしよう。入力したら実行するのを忘れないこと。

from PIL import Image, ImageDraw
import IPython

2つ目のセルに、関数newtonを実装しよう。先ほどと異なり、ニュートン法による反復を10回繰り返したのちに、収束した値を返す。

def newton(x):
  for _ in range(10):
    x = x - (x**3-1)/(3*x**2)
  return x

ここで、return xのインデントに注意。x = x - (x**3-1)/(3*x**2)ではなく、forと同じ高さにしなければならない。

3つ目のセルに、以下を入力せよ。

def plot(draw, s):
    hs = s/2
    red = (255, 0, 0)
    green = (0, 255, 0)
    blue = (0, 0, 255)
    for x in range(s):
        for y in range(s):
            z = complex(x-hs, y-hs)/s*4 + 0.01
            z = newton(z)
            # ここを埋めよ
            draw.rectangle([x, y, x+1, y+1], fill=c)

ただし、上記の「ここを埋めよ」の箇所に

• zの実部が正ならc = red
• zの実部が負かつ虚部が正ならc = green
• zの実部が負かつ虚部も負ならc = blue

を実行するようにプログラムを書くこと。

ここまで入力したプログラムを用いて、収束地図を作ってみよう。

size = 512
im = Image.new("RGB", (size, size))
draw = ImageDraw.Draw(im)
plot(draw, size)
im.save("test.png")
IPython.display.Image("test.png")

正しく入力できていれば、上記を実行した際に「収束地図」が描けたはずだ。どのような地図になっただろうか？

エラーが表示されたら、エラー内容をよく見てどこが間違っているかを確認せよ。

先程は$x^3 - 1 = 0$の解を考えた。次は$x^4 - 1 = 0$の解を考えてみよう。この方程式には$x = \pm 1$、$x \pm i$の4つの解が存在する。この解をニュートン法で探し、「収束地図」を描こう。

まず、先程の課題で作成したプログラムに名前をつけて保存し、複製して、そちらを修正することにしよう。ファイルメニューから「名前を変更」をクリック、もしくは上の「Untitled.ipynb」をクリックし、「newton3.ipynb」と名前を変えよう。その後、「ファイル」メニューの「保存」を選ぶと保存される。

この状態で「ファイル」メニューから「ドライブで探す」をクリックせよ。「マイドライブ」の「Colab Notebooks」に「newton3.ipynb」があるはずである。それを右クリックし「コピーを作成」を選ぶと「newton3.ipynbのコピー」というファイルが作成されるので、それを右クリックして「名前を変更」を選び、「newton4.ipynb」という名前にしよう。「newton4.ipynb」が作成されたら、右クリック→「アプリで開く」→「Colaboratory」を選ぶことで開くことができる。

「newton4.ipynb」が開かれたら、まず2つ目のセルにある関数newtonを次のように修正しよう。

def newton(x):
  for _ in range(10):
    x = x - (x**4-1)/(4*x**3)
  return x

$f(x) = x^4 -1$とし、$f(x)$と$f'(x)$を公式に代入しただけである。

3つ目のセルにある関数plotを次のように修正しよう。

def plot(draw, s):
    hs = s/2
    red = (255, 0, 0)
    green = (0, 255, 0)
    blue = (0, 0, 255)
    purple = (255, 0 ,255)  # ここを追加
    for x in range(s):
        for y in range(s):
            z = complex(x-hs, y-hs)/s*4 + 0.01
            z = newton(z)
            # ここを埋めよ
            draw.rectangle([x, y, x+1, y+1], fill=c)

解がそれぞれ1, -1, j, -jなので、

• zの実部が正かつ虚部が正ならc = red
• zの実部が正かつ虚部が負ならc = green
• zの実部が負かつ虚部が正ならc = blue
• zの実部が負かつ虚部も負ならc = purple

となるように「# ここを埋めよ」の部分を書け。

4つ目のセルはそのままで良い。すべてのセルを上から順番に実行、もしくは「ランタイム」メニューから「すべてのセルを実行」を実行せよ。もしおかしなことになったら、「ランタイム」から「再起動してすべてのセルを実行」せよ。

ニュートン法の繰り返し数が10だと原点付近の収束が甘い。20くらいにして再実行してみよ。逆に5に減らすとどうなるだろうか？

プログラムが何か意図しない動作をする場合、その原因となる箇所を「バグ」と呼ぶ。バグの語源については諸説あるようだが、詳しいことはわからない。「i」と「l」と「1」など、似ている文字を誤入力してしまったり、考慮すべきケースを忘れていたり、バグの原因は様々である。単純なバグについては、コンパイラや検査ツールの充実、テスト手法の向上などにより事前に検出できるようになってきた。そんななか、未だによく見かけるバグに「オーバーフローバグ」がある。コンピュータが扱える数字には上限がある。たとえば整数は32ビットで表現されることが多い。符号無し整数の場合、表現できる最大の数は4294967295、つまり43億ちょっとである。符号付きの場合は、符号に1ビット使うので最大の数はその半分になる。この数字を超える、すなわち最大値を取っている変数に1を足すと、またゼロに戻ってしまう。オーバーフローバグは、よくタイマー周りに潜む。最近だと、ボーイング787という飛行機の電源制御システムが、連続して248日動作させると不具合を起こすことが報告された。慣れたプログラマなら、「248日」と聞いた瞬間に「あ、オーバーフローやったな」と気がつく。248日とは21427200秒である。31ビットで表現できる最大の数は2147483647であるから10ミリ秒を単位に動作する31ビットのクロックが、248日でオーバーフローしたと考えられる。同根のバグに「497日問題」とか「49.7日問題」があるので、興味があれば調べられたい。

値があふれるのとは逆に、値が0になってるのに引き算をすることで値が大きくなるバグもあり、こちらもオーバーフローバグと呼ばれている。有名なのは「突然キレるガンジー」であろう。Cibilizationという、文明を発展させるゲームがある。歴史上の有名人をプレイヤーとして選び、世界制覇などを目指すゲームである。この中のプレイヤーにガンジーがいた。ガンジーは「非暴力、不服従」の提唱者であり、平和主義者なのであるが、文明がある程度発展すると突如として核攻撃をしかけてくるようになる。これはオーバーフローバグであった。Cibilizationでは各プレイヤーには攻撃性が設定されており、文明が民主主義を採用すると攻撃性が2下がる、という仕組みがあった。さて、ガンジーの攻撃性は1なのだが、インド文明が民主主義を採用し2を引くと-1になる。しかし、攻撃性は「符号なし8ビット整数」で表現されていたため、1から2を引くと攻撃性最大の255になってしまった。こうして「突然キレるガンジー」が誕生したのである。

整数の表現できる数値に最大値があることに起因するバグは根深く、発見が難しい。かくいう私も、4294967295回に一度意図しない動作をするというエグいバグを入れたことがある。43億回に一度発生するため、研究室のPCでは再現せず、スパコンを使ったときにたまに発生する、という感じだったので原因究明に時間がかかった。バグという概念が生まれてかなりの時間がたったが、まだ人類はバグを根絶できていない。それなら僕がバグを入れるのもやむを得まい(言い訳である)。

参考：オーバーフローが引き起こした面白いバグの話