8. 発展的な事例¶
この章では，前章でまだ力尽きていない人のために，spring_mass.py
を利用するもう少しだけ発展的な例を 2 つほど駆け足で紹介します．
1 つめは，バネ・質点のメッシュ構造による布らしきもののシミュレーションです．2 つめは，描画ストラテジを入れ替えることによる 3 次元化です．
8.1. 布シミュレーション¶
- マウス移動: 速度に応じて周辺の質点へ外力を加える
- ESCキー: 終了
8.1.1. バネ・質点メッシュを配置¶
まず，前章で作成したインタラクティブ放物シミュレーションの ver 10.0 をコピーして cloth_main.py を作ります．これは ActorFactory
に create_spring を追加する直前のバージョンです．
1import pygame
2import spring_mass as spm
3
4
5class ActorFactory:
6 def __init__(self, world, actor_list):
7 self.world = world
8 self.actor_list = actor_list
9
10 def create_point_mass(self, pos, fixed=False):
11 vel = (0, 0)
12 mass = 1.0
13 radius = 5
14 viscous = 0.05
15 restitution = 0.0
16 if fixed:
17 PointMassClass = spm.FixedPointMass
18 color = "gray"
19 else:
20 PointMassClass = spm.PointMass
21 color = "green"
22 return PointMassClass(pos, vel, self.world, radius, mass, viscous,
23 restitution, spm.CircleDrawer(color, width=0))
24
25 def create_spring(self, p1, p2, natural_len=0):
26 spring_const = 0.5
27 return spm.Spring(p1, p2, self.world, spring_const, natural_len,
28 spm.LineDrawer("white", 1))
29
30 def create_collision_resolver(self):
31 return spm.CollisionResolver(self.world, self.actor_list)
32
33 def create_boundary(self, name):
34 width, height = self.world.size
35 geometry = {"top": ((0, -1), (0, 0)),
36 "bottom": ((0, 1), (0, height)),
37 "left": ((-1, 0), (0, 0)),
38 "right": ((1, 0), (width, 0))}
39 normal, point_included = geometry[name]
40 return spm.Boundary(normal, point_included, self.world, self.actor_list)
41
42
43class AppMain:
44 def __init__(self):
45 pygame.init()
46 width, height = 600, 400
47 self.screen = pygame.display.set_mode((width, height))
48
49 self.world = spm.World((width, height), dt=1.0, gravity_acc=(0, 0.05))
50 self.actor_list = []
51 self.factory = ActorFactory(self.world, self.actor_list)
52
53 self.actor_list.append(self.factory.create_boundary("top"))
54 self.actor_list.append(self.factory.create_boundary("bottom"))
55 self.actor_list.append(self.factory.create_boundary("left"))
56 self.actor_list.append(self.factory.create_boundary("right"))
57
58 self.add_cloth()
59
60 def add_cloth(self):
61 n_columns, n_rows = 16, 12
62 x_offset, y_offset = 100, 20
63 spacing = 25
64 is_fixed = lambda i, j: (j == 0 and i in [0, n_columns // 2, n_columns - 1])
65 is_leftmost = lambda i, j: (i == 0)
66 is_topmost = lambda i, j: (j == 0)
67
68 plist = []
69 for j in range(n_rows):
70 for i in range(n_columns):
71 pos = (x_offset + i * spacing, y_offset + j * spacing)
72 p = self.factory.create_point_mass(pos, is_fixed(i, j))
73 plist.append(p)
74 self.actor_list.append(p)
75
76 if not is_leftmost(i, j):
77 sp = self.factory.create_spring(plist[-1], plist[-2], spacing)
78 self.actor_list.append(sp)
79 if not is_topmost(i, j):
80 sp = self.factory.create_spring(plist[-1], plist[-1 - n_columns], spacing)
81 self.actor_list.append(sp)
82
83 def update(self):
84 for a in self.actor_list:
85 a.update()
86 self.actor_list[:] = [a for a in self.actor_list if a.is_alive]
87
88 def draw(self):
89 self.screen.fill(pygame.Color("black"))
90 for a in self.actor_list:
91 a.draw(self.screen)
92 pygame.display.update()
93
94 def run(self):
95 clock = pygame.time.Clock()
96
97 while True:
98 frames_per_second = 60
99 clock.tick(frames_per_second)
100
101 should_quit = False
102 for event in pygame.event.get():
103 if event.type == pygame.QUIT:
104 should_quit = True
105 elif event.type == pygame.KEYDOWN and event.key == pygame.K_ESCAPE:
106 should_quit = True
107 if should_quit:
108 break
109
110 self.update()
111 self.draw()
112
113 pygame.quit()
114
115
116if __name__ == "__main__":
117 AppMain().run()
その上で以下の変更を加えます．
- 11-15, 25-28, 49行目
ActorFactoryにcreate_springメソッドを作ります．重力加速度やバネマスのパラメータを適宜調整しておきます．- 52行目
CollisionResolverは使わないので外しておきます．- 60-81行目:
AppMainにadd_clothを追加します．このメソッドではファクトリのcreate_point_massとcreate_springで質点とバネを生成しながらメッシュ状に結合し，actor_listに追加していきます．基本的には，69-70行目に始まる 2 重の
forループで，最上行から最下行の順で，行内では左から右の順で質点を作り出していきます．さらに必要に応じてバネによる接続を行います．具体的には- 最左列でないなら，現在の質点とその左隣の質点をつなぐ (76-78行目)
- 最上行でないなら，現在の質点とその上隣の質点とつなぐ (79-81行目)
この「左隣」や「上隣」を簡単に表すため，
actor_listとは別に質点だけを追加しておくリストplistを用意しています．最後の要素plist[-1]が最新の質点で，plist[-2]はその左隣です．こういうとき，Python のリストの「負のインデックス」は便利です．上隣も同様にplist[-1 - n_columns]で表せます．最左列かどうか，最上行かどうかの判断に使う関数
is_leftmost，is_topmostはラムダ式で定義しています (65-66行)．「ラムダ式って何だっけ」という人は 6.6 節を確認してください．また，最上行の左端，中央，右端の場合は固定質点にします (64, 72行目)． 64 行目のようなときに
in演算子を使えるのは Python ならではかも知れません．
実行すると，2次元アレイ上にバネで結ばれた質点が，重力によって懸垂する様子が見えると思います．
わかりやすい例としては，物理パラメータ群はこのサンプルコードのままで，時間ステップ (49行目の dt) だけをちょっと増やしてみると，比較的すぐに計算が発散します．こういうことが起きにくいように，
5 章で説明した通り Symplectic 積分法を導入しているのですが，十分ではないことがわかります．
発散を避ける安直な方法は時間ステップ dt を小さくすることです．ただし，もちろんシミュレーションの表示はその分だけ「スローモーション」になります．表示を同じ速さにするためには，例えば dt を N
分の 1 にしたのなら，110行目を:
for _ in range(N):
self.update()
として，N 時刻分の計算を進めてから描画に進むという方法が考えられます．ただし N を大きくすると計算量が増えるので，60 フレーム/秒で動かすのが難しくなるでしょう．
根本的には，より安定な数値積分法を導入することが必要といえます．
8.1.2. マウスによるインタラクション¶
布がぶら下がっているだけでは寂しいので，マウスによるインタラクションを導入します．ここでは，マウスポインタが動くと，その周囲の一定の範囲にある質点に力を及ぼすことにします．
1import pygame
2import spring_mass as spm
3
4
5def blow_by_mouse(actor_list, mouse_pos, mouse_rel):
6 mpos = pygame.math.Vector2(mouse_pos)
7 mrel = pygame.math.Vector2(mouse_rel)
8 max_distance = 50
9 sensitivity = 0.002
10 plist = [a for a in actor_list if isinstance(a, spm.PointMass)]
11 for p in plist:
12 distance = (mpos - p.pos).magnitude()
13 if distance <= max_distance:
14 p.receive_force(sensitivity * (max_distance - distance) * mrel)
15
95 def update(self):
96 blow_by_mouse(self.actor_list, pygame.mouse.get_pos(), pygame.mouse.get_rel())
97 for a in self.actor_list:
98 a.update()
99 self.actor_list[:] = [a for a in self.actor_list if a.is_alive]
100
AppMain の update 内で，現在のマウスポインタの位置と速度を取得して関数 blow_by_mouse を呼び出します (96行目)．
blow_by_mouse では，for 文ですべての質点についてループしながら，マウスポインタ位置からの距離が max_distance 以下ならば，マウスポインタの速度に比例する力を及ぼします．比例係数は距離に応じて変わるようにしています．
8.1.3. バネの破断¶
Spring を FragileSpring に置き換えれば，布を「引きちぎる」かのような挙動を実現することができます．
37 def create_spring(self, p1, p2, natural_len=0):
38 spring_const = 0.5
39 break_threshold = 5.0
40 return spm.FragileSpring(p1, p2, self.world, spring_const, natural_len,
41 spm.LineDrawer("white", 1), break_threshold)
42
8.2. 3 次元シミュレーションとその可視化¶
2 つめの例では，バネ・質点系のシミュレーションを 2 次元 (2D) から 3 次元 (3D) に拡張します．
- 左クリック: 質点を生成
- 右クリック: 固定質点を生成
- スペースキーを押しながらクリック: バネでつながれた質点を生成
- Ctrlキーを押しながらマウス移動: 座標系回転
- ESCキー: 終了
8.2.1. spring_mass.py の 3D 化¶
3D 拡張と一口に言っても，その内容は大きく 2 つに分かれます．1 つは，力学シミュレーション自体を 2D の計算から 3D の計算に変更することです．2 つめは，得られた 3D の計算結果を 2D のスクリーン上に描画するように変更することです．
このうち前者は実は簡単です．spring_mass.py では，力や運動の計算には pygame.math.Vector2 を用いてすべてベクトル表記で行っています．pygame のリファレンスマニュアルをよく読んでいる人はお気づきと思いますが，pygame には 3 次元ベクトルを表すクラス pygame.math.Vector3
も用意されているので，Vector2 を Vector3 に置き換えて，あとは初期化さえ何とかすれば 3D 化できてしまいます．
実は，このことを見越して pygame.math.Vector2 という名前を直接使わずに PgVector という略名を導入していました．PgVector が
pygame.math.Vector3 を指すように変更すれば，修正は基本的に完了です．
1import pygame
2
3
4PGVECTOR_DIMENSION = 2
5
6
7def PgVector(v):
8 if PGVECTOR_DIMENSION == 3:
9 if len(v) == 3:
10 return pygame.math.Vector3(v)
11 else:
12 return pygame.math.Vector3((v[0], v[1], 0))
13 else:
14 return pygame.math.Vector2(v)
15
16
17class World:
18 def __init__(self, size, dt, gravity_acc):
19 self.size = size
20 self.dt = dt
21 self.gravity_acc = PgVector(gravity_acc)
22
23
24class CircleDrawer:
25 def __init__(self, color, width):
26 self.color = pygame.Color(color)
27 self.width = width
28
29 def __call__(self, screen, center, radius):
30 pygame.draw.circle(screen, self.color, center[0:2], radius, self.width)
31
32
33class LineDrawer:
34 def __init__(self, color, width):
35 self.color = pygame.Color(color)
36 self.width = width
37
38 def __call__(self, screen, pos1, pos2):
39 pygame.draw.line(screen, self.color, pos1[0:2], pos2[0:2], self.width)
40
7-14行目が PgVector の新しい定義です．これまでいろいろ作ってきた
2D のアプリケーションはそのまま動くようにしたいので，デフォルトでは 2D
で動くようにして，特に指定したときだけ 3D モードになるようにします．定数 PGVECTOR_DIMENSION が 2 か 3 かに応じて Vector2 または
Vector3 を生成するようにします．また，Vector3 なのに 2D の初期値が与えらえた場合は第 3 要素に 0 を入れることにします (12行目)．
描画クラス CircleDrawer と LineDrawer も対応が必要です．真面目に 3D 対応するのは後から考えるとして，ひとまず PgVector 型オブジェクトには [0:2] をつけてしまいます．これで PgVector が 2D の場合も 3D の場合も，最初の 2 要素だけ，つまり x 座標と y 座標のみを使取り出します (30, 39行)．2D の場合は元のまま変わらず，3D の場合は x-y 平面への正射影を表示することになります．
ここまで変更したら，これまで作った 2D アプリケーション interactive_main.py や slinky_drop_main.py
や cloth_main.py が問題なく動作することを確認しておいてください．
さて，いよいよ 3D 計算をするためのメインプログラムを用意します．
interactive_main.py をコピーして interactive3d_main.py
を作成し，以下の変更を加えます．
1import random
2import pygame
3import spring_mass as spm
4
5spm.PGVECTOR_DIMENSION = 3
6
7
8class ActorFactory:
9 def __init__(self, world, actor_list):
10 self.world = world
11 self.actor_list = actor_list
12
13 def create_point_mass(self, pos, fixed=False):
14 vel = (random.uniform(-10, 10), random.uniform(-10, 0))
15 mass = 10
16 radius = 10
17 viscous = 0.01
18 restitution = 0.95
19 if fixed:
20 PointMassClass = spm.FixedPointMass
21 color = "gray"
22 else:
23 PointMassClass = spm.PointMass
24 color = "green"
25 return PointMassClass(pos, vel, self.world, radius, mass, viscous,
26 restitution, spm.CircleDrawer(color, width=0))
27
28 def create_spring(self, p1, p2):
29 spring_const = 0.01
30 natural_len = 20
31 break_threshold = 5.0
32 return spm.FragileSpring(p1, p2, self.world, spring_const, natural_len,
33 spm.LineDrawer("white", width=2), break_threshold)
34
35 def create_collision_resolver(self):
36 return spm.CollisionResolver(self.world, self.actor_list)
37
38 def create_boundary(self, name):
39 width, height = self.world.size
40 geometry = {"top": ((0, -1), (0, 0)),
41 "bottom": ((0, 1), (0, height)),
42 "left": ((-1, 0), (0, 0)),
43 "right": ((1, 0), (width, 0))}
44 normal, point_included = geometry[name]
45 return spm.Boundary(normal, point_included, self.world, self.actor_list)
46
47
48class AppMain:
49 def __init__(self):
50 pygame.init()
51 width, height = 600, 400
52 self.screen = pygame.display.set_mode((width, height))
53
54 self.world = spm.World((width, height), dt=1.0, gravity_acc=(0, 0.1))
55 self.actor_list = []
56 self.factory = ActorFactory(self.world, self.actor_list)
57
58 self.actor_list.append(self.factory.create_collision_resolver())
59 self.actor_list.append(self.factory.create_boundary("top"))
60 self.actor_list.append(self.factory.create_boundary("bottom"))
61 self.actor_list.append(self.factory.create_boundary("left"))
62 self.actor_list.append(self.factory.create_boundary("right"))
63
64 self.point_mass_prev = None
65
66 def add_connected_point_mass(self, pos, button):
67 if button == 1:
68 fixed = False
69 elif button == 3:
70 fixed = True
71 else:
72 return
73 p = self.factory.create_point_mass(pos, fixed)
74 self.actor_list.append(p)
75
76 if self.point_mass_prev is not None:
77 sp = self.factory.create_spring(p, self.point_mass_prev)
78 self.actor_list.append(sp)
79 if pygame.key.get_pressed()[pygame.K_SPACE]:
80 self.point_mass_prev = p
81
82 def update(self):
83 for a in self.actor_list:
84 a.update()
85 self.actor_list[:] = [a for a in self.actor_list if a.is_alive]
86
87 def draw(self):
88 self.screen.fill(pygame.Color("black"))
89 for a in self.actor_list:
90 a.draw(self.screen)
91 pygame.display.update()
92
93 def run(self):
94 clock = pygame.time.Clock()
95
96 while True:
97 frames_per_second = 60
98 clock.tick(frames_per_second)
99
100 should_quit = False
101 for event in pygame.event.get():
102 if event.type == pygame.QUIT:
103 should_quit = True
104 elif event.type == pygame.KEYDOWN and event.key == pygame.K_ESCAPE:
105 should_quit = True
106 elif event.type == pygame.KEYUP and event.key == pygame.K_SPACE:
107 self.point_mass_prev = None
108 elif event.type == pygame.MOUSEBUTTONDOWN:
109 self.add_connected_point_mass(event.pos, event.button)
110 if should_quit:
111 break
112
113 self.update()
114 self.draw()
115
116 pygame.quit()
117
118
119if __name__ == "__main__":
120 AppMain().run()
まず，spring_mass モジュールをいつものように spm という名前でインポートしてから，spm.PGVECTOR_DIMENSION を 3 に書き換えます (5
行目)．このようにモジュール内のグローバル定数を書き換えるのはあまりお行儀のよいことではないですが，まあよしとしましょう．
14行目で質点の初速度を，54行目で重力加速度を設定していますが，
interactive_main.py でいろいろ試してもらったので人によって値が違うかも知れません．ここでは当初の値に戻しています．まあ好きにしてください．
走らせてみると，3D シミュレーションの結果が表示されるはずです!!
…え? 何も変わったように見えない? そうですね．PgVector は 2 要素だけで初期化された場合 z 成分は 0 になるようにしたので，すべての Actorの初期 z 位置は 0，初速度も 0 なので，結局ずっと z 位置は 0 のままです．結果として，x-y 平面内で運動するものを x-y 平面に正射影しているだけなので，見た目は何も変わりません．
試しに 14行目を変えて z 方向にもランダムな初速度を与えてみれば，質点間の衝突が起きなくなることがわかると思います．質点間で z 位置が一般に一致しなくなるからです．
しかしこのままではあまり面白くないので，描画処理の修正に移ることにします．
8.2.2. 3D 可視化モジュール vis3d.py¶
3D の可視化を本格的に行うなら，OpenGL や DirectX といった 3D 描画を行う低レベルのライブラリを使うか，それらを内部的に利用する何らかの高機能ライブラリを使うことになるのですが，ここでは新しいライブラリを学ぶことはせずに，pygame の機能だけを使って簡易的な 3D 化を実現します．
spring_mass モジュールでは，描画処理は Strategy パターンで入れ替えられるようになっています．デフォルトの
CircleDrawer と LineDrawer は x-y 平面への正射影を行うようにさっき変更しましたが，別の投影を行うような CircleDrawer と
LineDrawer を新たに作ってそれらを使うようにすれば，可視化結果が置き換わります．
interactive3d_main.py の中で定義してもよいのですが，他のアプリケーションからも使えるように新たに vis3d.py というモジュールのファイルを作成することにします．
1import pygame
2
3
4class View:
5 def __init__(self, world):
6 self.world = world
7
8 def project_point(self, pos):
9 width = self.world.size[0]
10 pos_screen = (pos[2] + width/2, pos[1])
11 return pos_screen
12
13
14class CircleDrawer:
15 def __init__(self, view, color, width):
16 self.view = view
17 self.color = pygame.Color(color)
18 self.width = width
19
20 def __call__(self, screen, center, radius):
21 center_screen = self.view.project_point(center)
22 pygame.draw.circle(screen, self.color, center_screen, radius, self.width)
23
24
25class LineDrawer:
26 def __init__(self, view, color, width):
27 self.view = view
28 self.color = pygame.Color(color)
29 self.width = width
30
31 def __call__(self, screen, pos1, pos2):
32 pos1_screen = self.view.project_point(pos1)
33 pos2_screen = self.view.project_point(pos2)
34 pygame.draw.line(screen, self.color, pos1_screen, pos2_screen, self.width)
このモジュールでは，新しい CircleDrawer と LineDrawer を定義しています．3D 空間から 2D スクリーンへの投影機能は両者から共通して使用しますので，View という別のクラスにしておきます．
3D から 2D への投影を project_point メソッドにより行います．最初から凝った投影をするとよくわからなくなるので，まずは簡単に，y-z 平面への正射影を実装しておきます．10行目で，与えられた3 次元座標 pos
のうち z 成分 (pos[2]) をスクリーンの x 座標として，y 成分
(pos[1]) をスクリーンの y 座標として使います．ただし z = 0 のときにスクリーンの中央に描画されるよう，pos[2] には width/2 を足しています．
CircleDrawer と LineDrawer はこの View のインスタンスを保持して，描画の要求が来るたびに 3D 空間座標を 2D スクリーン座標に変換し
(21, 32-33行目)，描画します．
あとは，この新しい描画クラスを使うように interactive3d_main.py
を修正する必要があります．
1import random
2import pygame
3import spring_mass as spm
4import vis3d
5
6spm.PGVECTOR_DIMENSION = 3
7
8
9class ActorFactory:
10 def __init__(self, world, view, actor_list):
11 self.world = world
12 self.view = view
13 self.actor_list = actor_list
14
15 def create_point_mass(self, pos, fixed=False):
16 vel = (random.uniform(-10, 10), random.uniform(-10, 0))
17 mass = 10
18 radius = 10
19 viscous = 0.01
20 restitution = 0.95
21 if fixed:
22 PointMassClass = spm.FixedPointMass
23 color = "gray"
24 else:
25 PointMassClass = spm.PointMass
26 color = "green"
27 return PointMassClass(pos, vel, self.world, radius, mass, viscous,
28 restitution, vis3d.CircleDrawer(self.view, color, width=0))
29
30 def create_spring(self, p1, p2):
31 spring_const = 0.01
32 natural_len = 20
33 break_threshold = 5.0
34 return spm.FragileSpring(p1, p2, self.world, spring_const, natural_len,
35 vis3d.LineDrawer(self.view, "white", width=2), break_threshold)
36
37 def create_collision_resolver(self):
38 return spm.CollisionResolver(self.world, self.actor_list)
39
40 def create_boundary(self, name):
41 width, height = self.world.size
42 geometry = {"top": ((0, -1), (0, 0)),
43 "bottom": ((0, 1), (0, height)),
44 "left": ((-1, 0), (0, 0)),
45 "right": ((1, 0), (width, 0))}
46 normal, point_included = geometry[name]
47 return spm.Boundary(normal, point_included, self.world, self.actor_list)
48
49
50class AppMain:
51 def __init__(self):
52 pygame.init()
53 width, height = 600, 400
54 self.screen = pygame.display.set_mode((width, height))
55
56 self.world = spm.World((width, height), dt=1.0, gravity_acc=(0, 0.1))
57 self.view = vis3d.View(self.world)
58 self.actor_list = []
59 self.factory = ActorFactory(self.world, self.view, self.actor_list)
60
61 self.actor_list.append(self.factory.create_collision_resolver())
62 self.actor_list.append(self.factory.create_boundary("top"))
63 self.actor_list.append(self.factory.create_boundary("bottom"))
64 self.actor_list.append(self.factory.create_boundary("left"))
65 self.actor_list.append(self.factory.create_boundary("right"))
66
67 self.point_mass_prev = None
68
69 def add_connected_point_mass(self, pos, button):
70 if button == 1:
71 fixed = False
72 elif button == 3:
73 fixed = True
74 else:
75 return
76 p = self.factory.create_point_mass(pos, fixed)
77 self.actor_list.append(p)
78
79 if self.point_mass_prev is not None:
80 sp = self.factory.create_spring(p, self.point_mass_prev)
81 self.actor_list.append(sp)
82 if pygame.key.get_pressed()[pygame.K_SPACE]:
83 self.point_mass_prev = p
84
85 def update(self):
86 for a in self.actor_list:
87 a.update()
88 self.actor_list[:] = [a for a in self.actor_list if a.is_alive]
89
90 def draw(self):
91 self.screen.fill(pygame.Color("black"))
92 for a in self.actor_list:
93 a.draw(self.screen)
94 pygame.display.update()
95
96 def run(self):
97 clock = pygame.time.Clock()
98
99 while True:
100 frames_per_second = 60
101 clock.tick(frames_per_second)
102
103 should_quit = False
104 for event in pygame.event.get():
105 if event.type == pygame.QUIT:
106 should_quit = True
107 elif event.type == pygame.KEYDOWN and event.key == pygame.K_ESCAPE:
108 should_quit = True
109 elif event.type == pygame.KEYUP and event.key == pygame.K_SPACE:
110 self.point_mass_prev = None
111 elif event.type == pygame.MOUSEBUTTONDOWN:
112 self.add_connected_point_mass(event.pos, event.button)
113 if should_quit:
114 break
115
116 self.update()
117 self.draw()
118
119 pygame.quit()
120
121
122if __name__ == "__main__":
123 AppMain().run()
まず ActorFactory は初期化時に View のインスタンスを受け取って保持するようにします (10, 12行目)．その上で，質点の描画ストラテジを，
spm.CircleDrawer ではなく vis3d.CircleDrawer のインスタンスに置き換えます (28行目)．バネも同様です (35行目)．self.view を渡す必要があることを忘れないようにしてください．
AppMain の __init__ メソッドでは，View のインスタンスを生成し (57行目)，ファクトリの生成時に渡しておきます (59行目)．
これを実行すると，x-y 平面内で運動する Actor たちを y-z 平面に正射影して眺めることができます．何やら画面中央の縦一直線上で動いているのがが見えると思います．x-y 平面を y-z 平面に正射影すると 1 本の直線になるので，これで正しい動作です．…わかりますかね?
せっかく z 成分を可視化できるようになったので，初速度に z 成分を与えてみます．この例では z 成分はランダムな正の値を持つようにしています．スクリーン上では右に向かって飛んでいくことになるはずです．こちらが z の正の方向です．
15 def create_point_mass(self, pos, fixed=False):
16 vel = (random.uniform(-10, 10), random.uniform(-10, 0), random.uniform(-0, 10))
17 mass = 10
18 radius = 10
19 viscous = 0.01
20 restitution = 0.95
このままだと z の正の方向に飛んで行ったまま帰ってこないので， z 方向の手前と奥にも壁を作ってみます．
40 def create_boundary(self, name):
41 width, height, depth = self.world.size
42 geometry = {"top": ((0, -1, 0), (0, 0, 0)),
43 "bottom": ((0, 1, 0), (0, height, 0)),
44 "left": ((-1, 0, 0), (0, 0, 0)),
45 "right": ((1, 0, 0), (width, 0, 0)),
46 "front": ((0, 0, -1), (0, 0, -depth/2)),
47 "back": ((0, 0, 1), (0, 0, depth/2))}
48 normal, point_included = geometry[name]
49 return spm.Boundary(normal, point_included, self.world, self.actor_list)
50
51
52class AppMain:
53 def __init__(self):
54 pygame.init()
55 width, height, depth = 600, 400, 500
56 self.screen = pygame.display.set_mode((width, height))
57
58 self.world = spm.World((width, height, depth), dt=1.0, gravity_acc=(0, 0.1))
59 self.view = vis3d.View(self.world)
60 self.actor_list = []
61 self.factory = ActorFactory(self.world, self.view, self.actor_list)
62
63 self.actor_list.append(self.factory.create_collision_resolver())
64 self.actor_list.append(self.factory.create_boundary("top"))
65 self.actor_list.append(self.factory.create_boundary("bottom"))
66 self.actor_list.append(self.factory.create_boundary("left"))
67 self.actor_list.append(self.factory.create_boundary("right"))
68 self.actor_list.append(self.factory.create_boundary("front"))
69 self.actor_list.append(self.factory.create_boundary("back"))
70
71 self.point_mass_prev = None
width と height はこれまで通りそれぞれ 600 と 400 として，新たに depth = 500 を定義し (55行目)，World 型オブジェクトに渡しておきます (58行目)．
7 章に入るときに World の定義を変えたのは，実はこうやって 2 要素でも 3 要素でも初期化できるようにするためでした．
ActorFactory の create_boundary では，この depth を利用して
"front" と "back" の壁を生成できるようにしておきます (41-47行目)．それぞれ z = -depth/2 と z = depth/2 の位置に置くことにします．2D のスクリーン座標系と，3D 空間座標系を図示すると以下のようになります．
AppMain の __init__ で，front と back の壁を生成するようにすれば (68-69行目)，それらの位置で跳ね返るようになるのが確認できるはずです．
Boundary も自然に 3D 拡張されていることに注目してください．
以降では，座標系を回転させたり，投影方法を変えたりしていきます．その様子がよくわかるように，壁が置かれる位置をワイヤフレームで表示するようにしておきます．
6spm.PGVECTOR_DIMENSION = 3
7
8
9def draw_wireframe(screen, size, view):
10 width, height, depth = size
11 line_drawer = vis3d.LineDrawer(view, "gray64", 1)
12 dot_drawer = vis3d.CircleDrawer(view, "gray64", 0)
13
14 x_list = [0, width]
15 y_list = [0, height]
16 z_list = [-depth/2, depth/2]
17 for y in y_list:
18 for z in z_list:
19 x1, x2 = x_list
20 line_drawer(screen, (x1, y, z), (x2, y, z))
21 for x in x_list:
22 for z in z_list:
23 y1, y2 = y_list
24 line_drawer(screen, (x, y1, z), (x, y2, z))
25 for x in x_list:
26 for y in y_list:
27 z1, z2 = z_list
28 line_drawer(screen, (x, y, z1), (x, y, z2))
29 dot_drawer(screen, (0, 0, 0), radius=5)
30
31
32class ActorFactory:
117 def draw(self):
118 self.screen.fill(pygame.Color("black"))
119 draw_wireframe(self.screen, self.world.size, self.view)
120 for a in self.actor_list:
121 a.draw(self.screen)
122 pygame.display.update()
draw_wireframe という関数を作り，AppMain の draw メソッドの中で呼び出します (119行目)．draw_wireframe の中では，6 面の壁によって構成される直方体の各辺を描画し (17-28行目)，また，原点に小さな円を描画しています (29行目)．
注目すべきは，これらの描画も vis3d.LineDrawer と
vis3d.CircleDrawer によって行っているという点です．これによって，
View の内容を変えると，Actor だけでなくワイヤフレームの表示も適切に視点変換されることになります．
今の時点では y-z 平面へ正射影なので，front と back の壁の位置に縦線が描画されるだけです．跳ね返りが正しい位置で起きていることが確認できれば OK です．
8.2.3. 座標変換¶
もう少し 3D らしい表示にするため，座標系の回転と拡大縮小を実装します． 6 面の壁が構成する直方体の真ん中を中心として回転・拡大縮小することにしましょう．
本来，座標変換の定式化には線形代数の知識が必要ですが，
pygame.math.Vector3 には回転を行うメソッドが用意されています．これを使えば自分で計算する必要はありません．有効活用しましょう．
ただし，Vector3 のメソッドとして用意されているのは原点まわりの回転だけです．今，原点は直方体の中心ではありませんので，あらかじめ平行移動してから回転する必要があります．
また，拡大縮小は質点の位置だけではなく，半径にも適用する必要があります．
1import pygame
2from spring_mass import PgVector
3
4
5class View:
6 def __init__(self, world):
7 self.world = world
8 self.mag = 0.5
9 self.rot_x = 15
10 self.rot_y = 15
11
12 def project_point(self, pos, radius):
13 w, h, _ = self.world.size
14 pos_translated = PgVector(pos) - PgVector((w/2, h/2, 0))
15 pos_rotated = pos_translated.rotate_y(self.rot_y).rotate_x(self.rot_x)
16 pos_scaled = pos_rotated * self.mag
17 pos_screen = (pos_scaled[0] + w/2, pos_scaled[1] + h/2)
18 radius_screen = radius * self.mag
19 return pos_screen, radius_screen
20
21
22class CircleDrawer:
23 def __init__(self, view, color, width):
24 self.view = view
25 self.color = pygame.Color(color)
26 self.width = width
27
28 def __call__(self, screen, center, radius):
29 center_screen, radius_screen = self.view.project_point(center, radius)
30 pygame.draw.circle(screen, self.color, center_screen, radius_screen, self.width)
31
32
33class LineDrawer:
34 def __init__(self, view, color, width):
35 self.view = view
36 self.color = pygame.Color(color)
37 self.width = width
38
39 def __call__(self, screen, pos1, pos2):
40 pos1_screen, _ = self.view.project_point(pos1, 0)
41 pos2_screen, _ = self.view.project_point(pos2, 0)
42 pygame.draw.line(screen, self.color, pos1_screen, pos2_screen, self.width)
まず View の属性として拡大率 mag，x 軸と y 軸周りの回転角度
rot_x と rot_y を用意しておきます (8-10行目)．
位置と半径の両方を変換するため，project_point メソッドには引数
radius を追加し，返り値も位置と半径のタプルに変更します．
座標変換の処理は一見ややこしいですが，順番に見ていけばそれほど難しいことはありません:
- 直方体の中心が原点になるよう平行移動 (14行目)
- 直方体の中心まわりで回転 (15行目)
- 拡大・縮小(16行目)
- x-y 平面に正射影し，原点をスクリーン左上に戻す (17行目)
- 半径の拡大・縮小 (18行目)
これで，少し斜めの視点から可視化した結果が描画されるはずです．回転角度や拡大率をいろいろ書き換えて試してみてください．
なお，マウスでクリックした位置と質点が生成される位置は，スクリーン上で一致しなくなることに注意してください．スクリーンの座標 (x, y) をクリックすると，質点は 3 次元位置 (x, y, 0) に生成されますが，それは平行移動，回転，拡大・縮小を経て表示されるからです．後で直します．
8.2.4. 透視投影¶
3D での運動計算と表示ができていることはわかりましたが，ちょっと物足りないと感じるかもしれません．近いものは大きく，遠いものは小さく見えるような透視投影を実装すると，より 3D らしく見えるようになります．
透視投影の考え方は難しくありません．視点位置に仮想的なカメラの光学中心があり，物体からの光線は，光学中心を通り，光学中心から焦点距離だけ離れた位置にあるスクリーン (フィルムあるいはイメージセンサ) に像を作ります．光学中心から物体までの距離 (奥行き) を \(d\)，焦点距離を \(f\) とすると，三角形の相似から，像と物体の大きさの比 (拡大率) は \(f/d\) になります．
正射影の際に定数として与えていた拡大率を，こうやって計算するように変えればよいことになります．
5class View:
6 def __init__(self, world):
7 self.world = world
8 self.camera_offset_z = 2000
9 self.focal_len = 1000
10 self.rot_x = 15
11 self.rot_y = 15
12
13 def project_point(self, pos, radius):
14 w, h, _ = self.world.size
15 pos_translated = PgVector(pos) - PgVector((w/2, h/2, 0))
16 pos_rotated = pos_translated.rotate_y(self.rot_y).rotate_x(self.rot_x)
17
18 depth_from_camera = self.camera_offset_z + pos_rotated[2]
19 mag = self.focal_len / depth_from_camera
20
21 pos_scaled = pos_rotated * mag
22 pos_screen = (pos_scaled[0] + w/2, pos_scaled[1] + h/2)
23 radius_screen = radius * mag
24 return pos_screen, radius_screen
25
26
27class CircleDrawer:
View の属性として，拡大率の代わりに，焦点距離 focal_len と，視点から回転中心までの距離 camera_offset_z を設けます (8-9行目)．ただし，視線は z の正の方向を向き，回転中心を通るものとします．
このとき奥行き \(d\) は，camera_offset_z に，回転中心から物体位置までの奥行き pos_rotated[2] を加えることで得られます．これにより拡大率を求めます (18-19行目)．
3D の球を透視投影すると，スクリーン上では一般に楕円になります．が，面倒なので円で近似してしまいましょう．位置と同じ拡大率で半径を拡大・縮小します (23行目)．
以上で透視投影の実装は完了です．ワイヤフレームもそうですが，質点の大きさも遠近法で変わって見えるはずです．質点が小さくてわかりにくいかも知れません．お好みで半径を大きくしたりしてみてください．
38 def create_point_mass(self, pos, fixed=False):
39 vel = (random.uniform(-10, 10), random.uniform(-10, 0), random.uniform(-0, 10))
40 mass = 10
41 radius = 30
42 viscous = 0.01
43 restitution = 0.95
8.2.5. マウスによるインタラクション¶
ここまで来ると，マウスで座標系を回転してみたりしたくなります．回転量を与えるメソッドを用意することにします．
ついでに，決め打ちにしていた focal_len や camera_offset_z も
View の初期化時に与えられるようにしておきます．意欲がある人は，マウスホイールでこれらを変更できるようにすると面白いかも知れません．
5class View:
6 def __init__(self, world, camera_offset_z, focal_len):
7 self.world = world
8 self.camera_offset_z = camera_offset_z
9 self.focal_len = focal_len
10 self.rot_x = 15
11 self.rot_y = 15
12
13 def rotate_x_by(self, deg):
14 self.rot_x += deg
15
16 def rotate_y_by(self, deg):
17 self.rot_y += deg
18
19 def project_point(self, pos, radius):
これを使って，AppMain 側から回転量を制御します．ここでは，Ctrl キーを押しながらマウスを移動させたとき，移動量に応じて座標系が回転するようにしています．
75class AppMain:
76 def __init__(self):
77 pygame.init()
78 width, height, depth = 600, 400, 500
79 self.screen = pygame.display.set_mode((width, height))
80
81 self.world = spm.World((width, height, depth), dt=1.0, gravity_acc=(0, 0.1))
82 self.view = vis3d.View(self.world, 2000, 1000)
83 self.actor_list = []
84 self.factory = ActorFactory(self.world, self.view, self.actor_list)
124 def run(self):
125 clock = pygame.time.Clock()
126
127 while True:
128 frames_per_second = 60
129 clock.tick(frames_per_second)
130
131 should_quit = False
132 for event in pygame.event.get():
133 if event.type == pygame.QUIT:
134 should_quit = True
135 elif event.type == pygame.KEYDOWN and event.key == pygame.K_ESCAPE:
136 should_quit = True
137 elif event.type == pygame.KEYUP and event.key == pygame.K_SPACE:
138 self.point_mass_prev = None
139 elif event.type == pygame.MOUSEBUTTONDOWN:
140 self.add_connected_point_mass(event.pos, event.button)
141 if should_quit:
142 break
143
144 mouse_dx, mouse_dy = pygame.mouse.get_rel()
145 if pygame.key.get_pressed()[pygame.K_LCTRL]:
146 rot_sensitivity = 0.1
147 self.view.rotate_x_by(mouse_dy * rot_sensitivity)
148 self.view.rotate_y_by(-mouse_dx * rot_sensitivity)
149
150 self.update()
151 self.draw()
一つやり残した仕事がありました．マウスクリック位置を正しく 3D 空間内の質点生成位置に正しく反映することです．
2D から 3D への変換なので，奥行きは何らかの事前情報として与えてやる必要があります．ここでは，視点位置からの奥行き camera_offset_z の面上に決め打ちにします．
この処理を View に unproject_point メソッドとして実装します．
32 def unproject_point(self, pos_screen):
33 w, h, _ = self.world.size
34 mag = self.focal_len / self.camera_offset_z
35
36 pos_scaled_x, pos_scaled_y = pos_screen[0] - w/2, pos_screen[1] - h/2
37 pos_rotated = PgVector((pos_scaled_x / mag, pos_scaled_y / mag, 0))
38 pos_translated = pos_rotated.rotate_x(-self.rot_x).rotate_y(-self.rot_y)
39 pos = pos_translated + PgVector((w/2, h/2, 0))
40 return pos
41
42
43class CircleDrawer:
基本的には
project_point の処理を逆順にたどるだけです．回転の順序も逆順にしなくてはならない点に注意してください (38行目)．
奥行きを決め打ちにする件は，34行目の分母が camera_offset_z であることと 37行目で pos_rotated の z 成分を 0 にしていることで実現されています．
131 should_quit = False
132 for event in pygame.event.get():
133 if event.type == pygame.QUIT:
134 should_quit = True
135 elif event.type == pygame.KEYDOWN and event.key == pygame.K_ESCAPE:
136 should_quit = True
137 elif event.type == pygame.KEYUP and event.key == pygame.K_SPACE:
138 self.point_mass_prev = None
139 elif event.type == pygame.MOUSEBUTTONDOWN:
140 self.add_connected_point_mass(self.view.unproject_point(event.pos),
141 event.button)
142 if should_quit:
143 break
AppMain 側で，質点を生成するときにこの変換を通すことで，スクリーン上のクリックした位置に質点が現れることになるはずです．
奥行きについては，例えば右クリックで複数の固定質点を出現させてから座標系を回転してみれば「奥行き camera_offset_z の面上に決め打ち」の意味がよくわかると思います．
8.2.6. 遮蔽関係の表現¶
すべての質点を同じ色で描画しているうちは，以上のような実装で特に違和感はないと思います．しかし，異なる色の質点を表示すると途端に何かおかしいことに気づくはずです．
ファクトリを修正して，生成する度に質点の色が変わるようにしてみます．
32class ActorFactory:
33 def __init__(self, world, view, actor_list):
34 self.world = world
35 self.view = view
36 self.actor_list = actor_list
37 self.color_list = ["red", "magenta", "blue", "cyan", "green", "yellow"]
38 self.color_index = 0
39
40 def create_point_mass(self, pos, fixed=False):
41 vel = (random.uniform(-10, 10), random.uniform(-10, 0), random.uniform(-0, 10))
42 mass = 10
43 radius = 30
44 viscous = 0.01
45 restitution = 0.95
46 if fixed:
47 PointMassClass = spm.FixedPointMass
48 color = "gray"
49 else:
50 PointMassClass = spm.PointMass
51 color = self.color_list[self.color_index]
52 self.color_index = (self.color_index + 1) % len(self.color_list)
53 return PointMassClass(pos, vel, self.world, radius, mass, viscous,
54 restitution, vis3d.CircleDrawer(self.view, color, width=0))
55
56 def create_spring(self, p1, p2):
わかるでしょうか．質点の前後関係による遮蔽が正しく表現できていません．
これに対する簡単な対策は，視点から遠い順に描画することです．ちょうど油絵を描くときに，奥の物体を手前の物体で上描きしていく様子に似ているため画家のアルゴリズムと呼ばれることがあります (アルゴリズムというほどのものかという気もしますが)．互いにめり込まない球の前後関係を表現するだけならこれで十分です．
42 def depth(self, pos):
43 w, h, _ = self.world.size
44 pos_translated = PgVector(pos) - PgVector((w/2, h/2, 0))
45 pos_rotated = pos_translated.rotate_y(self.rot_y).rotate_x(self.rot_x)
46 depth_from_camera = self.camera_offset_z + pos_rotated[2]
47 return depth_from_camera
48
49
50class CircleDrawer:
View クラスに，視点位置からの奥行きを返すメソッド depth を追加しておきます．内容は project_point の途中で計算したものと同じです．
120 def draw(self):
121 self.screen.fill(pygame.Color("black"))
122 draw_wireframe(self.screen, self.world.size, self.view)
123 plist = [a for a in self.actor_list if isinstance(a, spm.PointMass)]
124 rest = [a for a in self.actor_list if not isinstance(a, spm.PointMass)]
125 plist.sort(key=lambda p: self.view.depth(p.pos), reverse=True)
126 for a in rest + plist:
127 a.draw(self.screen)
128 pygame.display.update()
これを利用して，AppMain の draw では，奥行き順で描画を行います．まず actor_list のうち質点のみを plist に，その他を rest
に取り出します (123-124行目)．plist を奥行きの降順でソートします
(125行目)．その後改めて rest と plist を描画します．
これで違和感のない描画ができるようになったはずです．
画家のアルゴリズムは考え方も実装も簡単ですが，表現できることに限界があります．例えば，CollisionResolver を外して質点どうしが重なり得るようにしてみてください．球と球が互いにめり込んだような絵にならなくてはならないのですが，当然そんな描画結果は得られません．
より一般に，3D グラフィクスでは物体形状は小さな三角形や四角形などの多角形 (ポリゴン) をつなぎ合わせたものとして表現されるのが通常ですが，ポリゴン相互間で同じ問題が生じます．
そのため，本格的な 3D グラフィクスを実現する場合は奥行きの比較を画素ごとに行います．よく用いられるのはデプスバッファ法と呼ばれる方法で，図形を 2D のスクリーンに描画する際に，各画素の色だけでなく奥行き量も書き込みます．ただし，その位置に既に書き込まれている奥行き量が，今書き込もうとしている奥行きより近い場合は，色も奥行きも更新しません．
このような奥行き比較の他にも，3D グラフィクスでは座標変換，陰影計算などのため大量の計算が発生します．これらを並列演算により高速に行うためのハードウェアは GPU (Graphics Processing Unit) と呼ばれます．GPU は，現在では 3D グラフィクスの計算だけではなく，科学技術計算，機械学習などに幅広く応用されています．