数学特別講義i：くらしの中のコンピュータと数学（2009年4月17日）

平成２１年度 筑波大学 理工学群数学類 数学特別講義Ⅰ

くらしの中のコンピュータと 数学 照井 章 （筑波大学 数理物質科学研究科 数学専攻） 2009年4月17日

誤入学おめでとう ございます！

自己紹介 • 1991年：筑波大学自然学類に入学 • 1995年：筑波大学大学院数学研究科に入 学

• 1999年：筑波大学助手 （2007年から助教に職名が変わる）

• 専門：計算機数学（計算代数、数式処理）

きょうのお題

•

身近なくらしの中で接している情報 をもたらすコンピュータの技術、そ の基盤となる数学 気象予報

• 身近なくらしの中で使っているコン ピュータの技術、その基盤となる 数学

暗 号

気象予報のはなし

気象予報

•

私達のくらしに欠かせない情報の一 つ 「お弁当屋さんからロケットまで」

• 昔：予報官の経験と勘 • 現在：コンピュータが大活躍

数値予報 大気の状態変化を数値的に計算して、将来の状態を 予測する天気用法の手法

（気象庁ホーム ページより）

数値予報の歴史 • 1922年: リチャードソン（英） 初めての数値予報（手計算！残念ながら失敗）

• 「64,000人の計算員が巨大なホールに集ま り、指揮者のもとで整然と計算を行えば、実 際の天候の変化と同程度の速さで予報が行え るだろう」 （「リチャードソンの夢」）

数値予報の歴史 • 1930∼40年代: • 気象観測技術の発達 （ラジオゾンデ、レーダー等の電子機器）

• 電子計算機の開発 • 1950年: 最初期の電子計算機ENIACを用いて 最初の数値予報に成功（米国）

数値予報の歴史

•

1955年: 米国で数値予報業 務を開始

•

1959年: 日本（気象庁）で 数値予報業務を開始

数値予報をどう行うか？

実際の数値予報の手法はかなり複雑なので より易しい例題を用いて解説します

板に熱を加えた時の 温度分布を求めよう • 四角い板の境界と下部 から熱を加える

• 熱平衡状態に達した時 の、板の温度分布 u(x, y) を求めよう

• 境界の温度分布は一定 • 下部から加える熱量は 一定

熱力学: フーリエの基本則 • 厚さdの板 • 面Bは面Aより温度が

A ∆x

B

∆uだけ高い

•

熱はBからAへと伝わ る

• AからBへ伝わる熱量 qは？

∆y

フーリエの基本則 A

• AからBへ伝わる熱 量qは

q = −(∆y) dk

• k: 熱伝導率 • ∆x −→ 0

!

∆u ∆x

"

∆x

B

.

とすると q = −(∆y) dk

∆y

!

∂u ∂x

"

.

温度勾配

フーリエの基本則 •時間Δtの間に

A

伝わる熱量を Δqとすると ! " ∂u ∆q = −(∆y) dk ∆t. ∂x

∆x

B

∆y

今度はこれを板の上で 考える（時間: Δt） y

∆q4

y + ∆y

熱の入り: ΔQ

∆q1

y

∆q3

x

x + ∆x

∆q2

x

それぞれの熱の出入りは ∆q1 = −(∆y) dk ∆q2 = −(∆y) dk ∆q3 = −(∆y) dk ∆q3 = −(∆y) dk

! ∂u !! ∆t ∂x !x ! ∂u !! ∆t ! ∂x x+∆x ! ∂u !! ∆t ! ∂y y ! ∂u !! ∆t ! ∂y y+∆y

熱の入りの合計ΔQは ∆Q = ∆q1 − ∆q2 + ∆q3 − ∆q4 = (∆x)(∆y)d  $ ∂u $  k ∂x x+∆x − k × ∆x

.

$ ∂u $

∂x x

+

$ ∂u $ k ∂y $

$  ∂u $ − k ∂y $ y+∆y y . ∆y

一方、入る熱量と温度の 上昇率の関係は

∆Q = c ρ (∆x)(∆y)d 比熱

密度

!

∂u ∂t

"

温度の 上昇率

∆t

以上の関係をまとめると ∂u cρ = ∂t

! ! k ∂u ∂x x+∆x − k ∆x

! ∂u !

∂x x

+

! ∂u ! k ∂y !

y+∆y

! ∂u ! − k ∂y !

∆y

∆x −→ 0, ∆y −→ 0 として極限をとると ! " ! " ∂u ∂ ∂u ∂ ∂u cρ = k + k . ∂t ∂x ∂x ∂y ∂y

熱伝導率 k: 一定とすると ∂u cρ =k ∂t

!

2

2

∂ u ∂ u + ∂x2 ∂y 2

"

.

y

板の温度分布 u(x, y) は 以下の微分方程式をみたす

• 熱平衡状態の場合 2

2

∂ u ∂ u + = 0 ∂x2 ∂y 2

（ラプラス方程式）

f(x, y) が与えられる場合 • 外部から熱量 ! "  

2

−

2

∂ u ∂ u + ∂x2 ∂y 2

=f

（ポアソン方程式）

楕円型偏微分方程式

境界条件 • 板の外周（境界）でどのように熱を与 えるかで、板の温度分布が決まる

• 境界条件: g(x, y)  楕円型偏微分方程式の 境界値問題

偏微分方程式の数値解法: 差分法 • 板の温度分布 u(x, y) を求める ! 2 " 2 板の内部で ∂ u ∂ u • −k + =f ∂x2

∂y 2

• 熱伝導率 k = 1 • f(x, y): 外部から与える熱量 

• 板の境界で u(x, y) = g(x, y) をみたす • g(x, y): 境界条件 

差分法

温度を求める点 境界条件

y

• 領域を x 方 向、y 方向に それぞれ m+1 分割する (m = 3の例)

g14

g24

g34

g03

u13

u23

u33

g43

g02

u12

u22

u32

g42

g01

u11

u21

u31

g41

g10

g20

g30

1

• h = 1/(m+1) とおく

0

1

x

各点で、微分を差分で 近似する（中心差分） ! 2 ! ∂ u! − " − ! 2 ∂x i,j

! ∂u !

∂x i+1/2,j

! ∂u ! − ∂x i−1/2,j h

に、1階微分（差分）を代入する ! ∂u !! ui+1,j − ui,j − " − ∂x !i+1/2,j h

! ∂u !! ui,j − ui−1,j − "− ! ∂x i−1/2,j h

y

g14

g24

g34

g03

u13

u23

u33

g43

g02

u12

u22

u32

g42

g01

u11

u21

u31

g41

g10

g20

g30

1

0

1

x

各点で、微分を差分で 近似する（中心差分） 代入すると

! 2 ! ∂ u! 1 − [u − 2u + u ] " − i+1,j i,j i−1,j ! 2 ∂x i,j h2

y についても同様に差分を計算すると

! 2 ! ∂ u! 1 − [u − 2u + u ] " − i,j+1 i,j i,j−1 ! 2 2 ∂y i,j h

格子点 (i, j) における 微分方程式の差分による近似は −ui,j−1 − ui−1,j + 4ui,j − ui+1,j − ui,j+1 = h fi,j 2

ui,j •   に関する1次方程式 • すべての内部格子点についてまとめる と、連立1次方程式が得られる

連立1次方程式 Au = b u = (u11 , u21 , u31 , u12 , u22 , u32 , u13 , u23 , u33 ) t





4 −1 −1 −1 4 −1  −1     −1 4 −1   −1  4 −1 −1   . −1 −1 4 −1 −1 A=     −1 −1 4 −1     −1 4 −1    −1 −1 4 −1 −1 −1 4





h f11 + g10 + g01  h2 f21 + g20   2  h f31 + g30 + g41    2  h f12 + g02    2 . h f22 b=    h2 f32 + g42   2  h f13 + g03 + g14    2  h f23 + g24  h2 f33 + g43 + g34 2

連立1次方程式 Au = b を解くことにより、 温度分布 u = t (u11 , u21 , u31 , u12 , u22 , u32 , u13 , u23 , u33 )

が得られる！

y g14

g24

g34

g03

u13

u23

u33

g43

g02

u12

u22

u32

g42

g01

u11

u21

u31

g41

g10

g20

g30

1

0

1

x

実際の数値予報 数理モデル の構成

問題を 差分化

初期値の設定 （観測など）

y 2

2

∂ u ∂ u + = 0 ∂x2 ∂y 2

g14

g24

g34

g03

u13

u23

u33

g43

g02

u12

u22

u32

g42

g01

u11

u21

u31

g41

g10

g20

g30

1

0

1

x

計算

計算の難しさ • しばしば、大規模な連立1次方程式の 解法が必要

• それを予報時間に合わせて反復計算 • 超高性能計算（能力・機器）が必要 ハイパフォーマンスコンピューティング (High Performance Computing, HPC)

筑波大学における ハイパフォーマンスコン ピューティング

CP-PACS (1996∼2005) （計算物理学研究センター（当時））

• 理論物理学の計算 を主目的に構築

• 1996年11月、世 界のスーパーコン ピュータTop 500 リスト第1位！ （出典: 筑波大学計算科学研究センター）

T2K筑波 (2008∼) （計算科学研究センター） • T2Kオープンスー パーコンピュータ

!

• 筑波・東大・京大 で仕様の基本部分 を共通化

• 2008年11月、Top 500リスト第32位 （国内第3位）

! ! ! ! （出典: !

筑波大学プレスリリース）

""#$ #$%&'()*+,-./,0123 %&' 456,78456 123 '(( (0>?@,AB?C&'()D! !

暗号のはなし

暗号の用途

•

秘匿 他人に知られずに情報を送る

•

認証 情報が真正であることを保証 する

（ほぼ）毎日使う暗号: SSL (Secure Socket Layer)

• 秘匿 すべての入出力を暗 号化して情報を送る

• 認証 TWINSのサーバが 本物であることを 保証する https://twins.tsukuba.ac.jp/

暗号の歴史

• 情報あるところに暗号あり • 人類の歴史＝暗号の歴史？

シーザー暗号

（出典: Wikipedia）

エニグマ暗号機 • 第二次世界大戦にかけ て、ナチス・ドイツが 使用

• イギリスの数学者アラ ン・チューリング率い るグループ（「ウルト ラ」Ultra）が解読 英国王立戦争博物館（ロンドン）蔵 （出典: Wikipedia）

暗号化方式と「かぎ」 平文

暗号文

平文 復号

暗号化 暗号化鍵

復号鍵 解 読

太郎 (Bob) さん

盗聴者

花子 (Alice) さん

共有鍵暗号: 暗号化と復号の鍵が同一 平文

暗号文

平文 復号

暗号化 暗号化鍵

復号鍵 解 読

太郎 (Bob) さん

花子 (Alice) さん

盗聴者

鍵をどうやって安全に渡すか？

公開鍵暗号

•1976年、DiffieとHellmanに よって発表された概念

•2つの鍵のペアを用意 •1つは秘密（秘密鍵） •もう1つは公開（公開鍵）

公開鍵暗号: 暗号化の場合 平文

暗号文 暗号化 公開鍵

太郎 (Bob) さん

平文 復号 秘密鍵

花子 (Alice) さん

安全に暗号化して情報を伝達可能

公開鍵暗号: 認証の場合 平文

暗号文 復号 公開鍵

太郎 (Bob) さん

平文 暗号化 秘密鍵

花子 (Alice) さん

花子さんが送った情報であり、 かつ改ざんされていないことが保証される

RSA公開鍵暗号 • 1977年、Rivest, Shamir, Adelmanに よって公表された、初の実用的な公開 鍵暗号

• 一方の鍵で暗号化 もう一方の鍵で復号

• 3人は、この功績により、ACM（米国 計算機学会）チューリング賞を受賞

RSA公開鍵暗号を 使ってみよう

下準備

通信に使う文字と符号を 決める

• •符号：自然数（の集まり）

文字：通信文（平文）に使う情報 文字と1対1対応

アルファベット 大文字で！

空 白

文字

A

B

C

・・・

Y

符号

1

2

3

・・・

25 26 27

Z

p, q, n を決める • p, q: 素数 •n = p q • n が符号の最大値を上回るようにする • p = 3, q = 11, n=3

11 = 33

Lを決める

•

L: p - 1 と q - 1 の最小公 倍数 (LCM)

•

L = LCM(2, 10) = 10

eを決める

•

e: Lと互いに素、かつLよ りも小さい正数

•

L = 10 より e = 3

dを決める

•

d: e d 1 (mod 10) を満たす数

•

e(=3) d を10で割って、 余りが1になるような数

dを決める

•

3と10は互いに素なので、 3 d = 10 b + 1 をみたすd, bが存在する

•

d = 7, b = 2 とおけばよい

秘密鍵, 公開鍵を定める

• 公開鍵: (n, e) = (33, 3) •

秘密鍵: (n, d) = (33, 7)

実際の通信

•

太郎さんから花子さんに TSUKUBA と送ろう

公開鍵を相手に送る

了解★

公開鍵は (33, 3) よ！

太郎さん: 暗号化のしかた • 文字を符号に直す（符号化） • 例: T → 20 • 花子さんの公開鍵 (n, e) に対し、符号をe乗 してnで割った余りを求める

• 例: 20^3 = 8000

14 (mod 33)

• 計算した数を文字に直す（暗号化） • 例: 14 → N

太郎さん: 暗号化のしかた 平文

符号

33で割った 剰余

3乗

暗号文

T

20

8000

14

N

S

19

6859

28

(28)

U

21

9261

21

U

K

11

1331

11

K

U

21

9261

21

U

B

2

8

8

H

A

1

1

1

A

花子さん: 復号のしかた • 文字を符号に直す（符号化） • 例: N → 14 • 自分の秘密鍵 (n, d) に対し、符号をd乗してn で割った余りを求める

• 例: 14^7 = 105413504

20 (mod 33)

• 計算した数を文字に直す（復号） • 例: 20 → T

花子さん: 復号のしかた 暗号文

符号

33で割った 剰余

7乗

復号文

N

14

105413504

20

T

(28)

28

13492928512

19

S

U

21

1801088541

21

U

K

11

19487171

11

K

U

21

1801088541

21

U

H

8

2097152

2

B

A

1

1

1

A

実際の通信

平文 TSUKUBA

暗号文 N(28)UKUHA

復号文 TSUKUBA

なぜ正しく動くか？

•M: 平文の符号 E: Mに対応する暗号文の符号 • E M^e (mod n) • （EはMをe乗してnで割った 剰余）

なぜ正しく動くか？

• F: Eを復号した符号 • F E^d (mod n) (M^e)^d (mod n) M^(ed) (mod n)

• Fは、Mをed乗したものをnで割っ た剰余に等しい

なぜ正しく動くか？

•

n, e, dの選び方から M

M^(ed) (mod n)

が成り立つ

•すなわち、Mをed乗したものをn で割った剰余は、M自身に等し い！

RSA公開鍵暗号は SSLにも使われている

ここまでで・・・

たぶん湧き上がる疑問

•

これってホントに計算機の 話？

•

解析や代数の話じゃない？

計算機数学の諸問題は、他の数学 を使ってアタックするものが多い 代数学 ・暗号 ・計算代数 ・数値線形代数 数理統計学 ・機械学習

解析学 ・シミュレーション ・数理計画法

計算機数学

幾何学 ・グラフ理論

数学基礎論 ・機械（援用）証明 ・計算機の正当性

いろいろな数学をベースに、計算機 を使って新しい問題にチャレンジ 代数学 ・暗号 ・計算代数 ・数値線形代数 数理統計学 ・機械学習

解析学 ・シミュレーション ・数理計画法

計算機数学

幾何学 ・グラフ理論

数学基礎論 ・機械（援用）証明 ・計算機の正当性

数学の素養が重要！ 代数学 ・暗号 ・計算代数 ・数値線形代数 数理統計学 ・機械学習

解析学 ・シミュレーション ・数理計画法

計算機数学

幾何学 ・グラフ理論

数学基礎論 ・機械（援用）証明 ・計算機の正当性

誤 入学にならないために （旅にたとえて）

パルテノン神殿（ギリシャ、2005年）

高校までの数学（学問）: ガイドツアー • 授業＝ 先生のガイ ド付き

• 目的地まで 連れて行っ てくれる

•

解説完備

パルテノン神殿の裏側では 復元工事が続いている

大学からの数学（学問）: 個人旅行 • 授業＝ 旅のガイドブックの 紹介

• 自分で読み、 自分で計画を立て、 自分で行ってみよう

ギリシャの首都アテネにて、 自分の誕生日に乾杯！ （2005年9月16日）

御 入学おめでとうございます

皆さんのこれからの 大学生活が よい「旅」に なりますよう・・・

参考文献：数値予報 • 岩崎俊樹, 数値予報̶スーパーコンピュータ を利用した新しい天気予報, 情報フロンティ アシリーズ（情報処理学会 編）, 共立出版, 1993.

• 二宮洸三, 数値予報の基礎知識, オーム社, 2004.

• 時岡達志, 山岬正紀, 佐藤信夫, 気象の数値シ ミュレーション, 気象の教室 5, 東京大学出 版会, 1993.

参考文献：微分方程式の 数値解法

• 名取亮, 数値解析とその応用, コン ピュータ数学シリーズ 15, コロナ社, 1990.

参考文献：RSA公開鍵暗号 •

伊藤正史, サルにもわかるRSA暗号, まいとう情報通信 研究会 (webサイト). http://www.maitou.gr.jp/rsa/

•

山本和彦, はやわかりRSA, IIJ技術研究所 (webサイト). http://pgp.iijlab.net/crypt/rsa.html

•

T.H. Cormen, C.E. Leiserson, R.L. Rivest, C. Stein, Introduction to Algorithms, 2nd ed., MIT Press, 2001.

•

V. Shoup, A Computational Introduction to Number Theory and Algebra, Cambridge University Press, 2005.

イラストの出典 •

可知 豊, people biz - male smile (p. 41, 42, 44, 45, 56, 57, 62), Open Clip Art Library. http://openclipart.org/media/files/yyycatch/11305

•

可知 豊, people biz - female smile (p. 41, 42, 44, 45, 56, 57, 62), Open Clip Art Library. http://openclipart.org/media/files/yyycatch/11308

•

Anonymous, Man Silhouette (p. 41, 42), Clker.com - public domain royalty free clip art. http://www.clker.com/clipart-24011.html

•

Carsten Schmitz, Closed envelope (p. 41, 42, 44, 45), Open Clip Art Library. http://openclipart.org/media/files/c_schmitz/3117

イラストの出典 •

Juliane Krug, Key (p. 44, 45), Clker.com - public domain royalty free clip art. http://www.clker.com/clipart-3770.html

•

Barret Ruttan, Key (p. 41, 42, 44, 45), Open Clip Art Library. http://openclipart.org/media/files/barretr/2948

•

Nicu Buculei, Message in a Bottle (p. 41, 42, 44, 45, 62), Open Clip Art Library. http://openclipart.org/media/files/nicubunu/3342

•

Lars Sundström, Steel Floor Texture (p. 12, 13, 14, 15), Stock.XCHNG. http://www.sxc.hu/photo/871271