...

高知大学教育学部の情報数学のテキスト 1 チューリング・マシン

by user

on
Category: Documents
10

views

Report

Comments

Transcript

高知大学教育学部の情報数学のテキスト 1 チューリング・マシン
高知大学教育学部の情報数学のテキスト
1
チューリング・マシン
チューリング・マシンは計算模型のひとつで計算機を数学的に議論するための、単純化・理想化
された仮想機械である。1936 年にイギリスの数学者アラン・チューリングの論文「計算可能数に
ついて──決定問題への応用」で発表された。
チューリング・マシンは 1 本のテープと、このテープ上の記号を読み書きするためのヘッドを一
つ持つだけの仮想機械である。チューリング・マシンには有限個の「状態」があって、各時刻でそ
のうち何れか一つの状態をとることができる。
チューリング・マシンのプログラムは次の五つの記号で表される命令のリストで表現されます。
一組の命令は
<現在の状態、読み込んだカラー、書き込むカラー、新しい状態、動く方向>
の形で表されます。 このチューリング・マシンではスペース(または空白)のカラーはBで表
現します。また、動く方向はLとRで示します。 チューリング・マシンの処理は、次の六つのス
テップの繰り返しです。
■ ヘッドの現在いるセルの記号を読む
■ プログラムの全体を調べて、現在の状態である<現在の状態、読み込んだカラー>で始ま
る命令セットを見つけだす。該当する命令セットがなければ、停止する。見つかったら、その命令
セットの<書き込むカラー、新しい状態、動く方向>を確認する。
■ ヘッドの状態を<新しい状態>にする。
■ 現在のセルの内容を<書き込むカラー>の記号に変える。
■ <動く方向>によって指定された位置にヘッドを移動する。
■ 最初のステップに戻る。
TURING.EXE でプログラミングしてみましょう。
TURING.EXE を立ち上げます。
セルとヘッドの初期状態はメニューの set の table で指定します。メニューの set の table をク
リックします。
1
このダイアログボックスで、初期状態と初期データを入力します。通常、初期状態は 0 ですの
で、デフォルトのままで良いです。初期データは、1行に一文字入力します。例えば、次のように,、
001101011 と入力します。
OK のボタンを押すと次のようになります。
チューリング・マシンのプログラムは、小さな子供の遊んでいる様子をイメージすればいいです。
2
上図のような初期データが与えられたとき 0 をすべて 1 に変えるプログラムを考えて見ます。
先頭から右に一つずつ調べて、0 があれば 1 に変え、1 であればそのままにします。一つ処理し
た後も同じことをすればいいので、状態を変える必要はありません。従って、
0
0
1
0
R
0
1
1
0
R
とプログラミングすればいいです。大して知的なことをしていないので、プログラムも簡単です。
このプログラムを左側のメモにプログラムを打ち込みます。プログラムは直接編集することも保存
することも読み込むこともできます。命令リストの最後は何も書き込まないで下さい。改行だけの
行を作らないという意味です。
Start のボタンを押せば、チューリング・マシンが動き出します。
すべての0が1になりました。
同様に、0 と 1 をすべて入れ替えるプログラムは、
0
0
0
1
1
0
0
0
R
R
3
とプログラミングすればいいです。初期データを入力し直して、実行してみて下さい。
<例題> SAMPLE DATA : XXXXX が与えられた時、XXXXXX のように、最後に一つ X を
追加するプログラムを作れ。
解答例:
0 X X 0 R
0 B X 1 R
これを
0 X X 0 R
0 B X 0 R
とすると何が起こるか確かめなさい。
問題:例えば SAMPLE DATA : XXXXX が与えられた時、XXXXXXX のように、最後に2つ
X を追加するプログラムを作れ。
<例題> SAMPLE DATA : XXX が与えられた時、BBBBXXX のように、X 達を移動するプ
ログラムを作れ。
解答例:
0 X B 1 R
1 X X 1 R
1 B B 2 R
2 B X 3 L
3 B B 4 L
4 X X 4 L
4 B B 0 R
2 X X 2 R
3 X X 3 L
これは次のように考えればよい。初期状態 0 で X を読めば、空白 (B) を書き込み、「X を抱えて
置く場所を見つけるために区切りの空白を見つけようと右に移動している状態 1 」として、右に
移動する。状態 1 で X を読めば、X を読み飛ばして、同じ状態のまま右に移動する。状態 1 で
B (空白) を読めば、B (空白) を読み飛ばして、
「X を抱えて置く場所を見つけるために最初の空白
を見つけようと右に移動している状態 2 」として、右に移動する。状態 2 で B (空白) を読めば、
X を書き込み、「次の X を見つけるために区切りの空白を見つけようと左に移動している状態 3 」
として、左に移動する。状態 3 で B (空白) を読めば、B (空白) を読み飛ばして、「次の X を見
つけるために最初の空白を見つけようと左に移動している状態 4 」として、左に移動する。状態
4 で X を読めば、X を読み飛ばして、同じ状態のまま左移動する。状態 4 で B (空白) を読めば、
B (空白) を読み飛ばして、状態を 0 にして右に移動する。
0 X B 1 R
1 X X 1 R
1 B B 2 R
2 B X 3 L
3 B B 4 L
4
4 X X 4 L
4 B B 0 R
これでプログラムが出来たと思って実行してみると
と、「X を抱えて置く場所を見つけるために最初の空白を見つけようと右に移動している状態 2 」
で X を読んだところで停止します。すでに移動した X がある場合の考察を忘れていました。状態
2 で X を読めば、X を読み飛ばして、同じ状態のまま右移動するという命令
2 X X 2 R
を追加して、再び実行します。
今度は、
「次の X を見つけるために区切りの空白を見つけようと左に移動している状態 3 」で、X
を読んだところで停止します。ここでも、すでに移動した X がある場合の考察を忘れていました。
3 X X 3 L
を追加して、再び実行します。
5
今度は旨くいきました。このようにプログラムを作るときは、最初からすべてが完全に出来ること
はまずないです。何度も何度も失敗し、あらゆる失敗を経験した後、プログラムが作れるようにな
ります。以下で学ぶどのプログラミング言語を使う場合でも、プログラミングが出来るようになる
かどうかは、何度も何度も失敗し、あらゆる失敗を経験するまで我慢出来るかどうかです。
問題: SAMPLE DATA : XXX が与えられた時、BBBBXXXXXX のように、X の個数の2倍
の X を描くプログラムを作れ。
例えば SAMPLE DATA : BXXXX が与えられた時、
0 B B 1 R
1 X X 1 R
1 B B 2 L
2 X X 2 L
2 B B 1 R
のプログラムを実行すると何が起こるか確かめなさい。ここで、B は空白を示しています。
例えば SAMPLE DATA : BXXXX が与えられた時、
0 B B 1 R
1 X X 1 R
1 B X 2 L
2 X X 2 L
2 B B 1 R
のプログラムを実行すると何が起こるか確かめなさい。ここで、B は空白を示しています。
<例題> SAMPLE DATA : ()(()()) や (()(() が与えられた時、括弧がうまく対応しているかど
うか判定できるプログラムを作れ。
解答例:
0 ( Y 1 R
0 X X 0 R
0 ) Z 3 R
0 B B 5 R
1 ( ( 1 R
6
1 ) X 2 L
1 X X 1 R
1 B B 4 R
2 ( X 1 R
2 X X 2 L
2 Y Z 0 R
2 Z Z 3 R
5 B O 6 R
6 B K 7 R
3 ( ( 3 R
3 ) ) 3 R
3 B B 8 R
8 B N 9 R
9 B O 10 R
4 B B 8 R
実行すると
となります。
少し長いですが、2 進数の足し算をするプログラムを作ってみます。
7
0 0 Z 1 R
0 1 W 1 R
1 0 0 1 R
1 1 1 1 R
1 + + 1 R
1 = = 2 L
2 0 Y 4 L
2 1 Y 5 L
3 0 Y 5 L
3 1 Y 6 L
4 + + 7 L
5 + + 8 L
6 + + 9 L
7 X X 7 L
8 X X 8 L
9 X X 9 L
7 0 X 10 R
7 1 X 11 R
8 0 X 11 R
8 1 X 12 R
9 0 X 12 R
9 1 X 13 R
10 0 0 10 R
10 1 1 10 R
10 + + 10 R
10 = = 10 R
10 X X 10 R
10 Y Y 10 R
11 0 0 11 R
11 1 1 11 R
11 + + 11 R
11 = = 11 R
11 X X 11 R
11 Y Y 11 R
12 0 0 12 R
12 1 1 12 R
12 + + 12 R
12 = = 12 R
12 X X 12 R
12 Y Y 12 R
13 0 0 13 R
13 1 1 13 R
13 + + 13 R
8
13 = = 13 R
13 X X 13 R
13 Y Y 13 R
10 B B 14 L
11 B B 15 L
12 B B 16 L
13 B B 17 L
14 0 B 18 R
14 1 B 19 R
15 0 B 20 R
15 1 B 21 R
16 0 B 22 R
16 1 B 23 R
17 0 B 24 R
17 1 B 25 R
18 B 0 14 L
19 B 1 14 L
20 B 0 15 L
21 B 1 15 L
22 B 0 16 L
23 B 1 16 L
24 B 0 17 L
25 B 1 17 L
14 B B 14 L
15 B B 15 L
16 B B 16 L
17 B B 17 L
14 = = 26 R
15 = = 27 R
16 = = 28 R
17 = = 29 R
26 B 0 2 L
27 B 1
28 B 0
29 B 1
2 = =
3 = =
2 L
3 L
3 L
2 L
3 L
2 Y Y 2 L
3 Y Y 3 L
2 + + 30 L
3 + + 31 L
30 X X 30 L
31 X X 31 L
9
30 Z X 32 R
30 W X 33 R
31 Z X 33 R
31 W X 34 R
7 Z X 32 R
7 W X 33 R
8 Z X 33 R
8 W X 34 R
9 Z X 34 R
9 W X 35 R
32 X X 32 R
32 + + 32 R
32 Y Y 32 R
32 = = 32 R
32 0 0 32 R
32 1 1 32 R
33 X X 33 R
33 + + 33 R
33 Y Y 33 R
33 = = 33 R
33 0 0 33 R
33 1 1 33 R
34 X X 34 R
34 + + 34 R
34 Y Y 34 R
34 = = 34 R
34 0 0 34 R
34 1 1 34 R
35 X X 35 R
35 + + 35 R
35 Y Y 35 R
35 = = 35 R
35 0 0 35 R
35 1 1 35 R
32 B B 36 L
33 B B 37 L
34 B B 38 L
35 B B 39 L
36 0 B 40 R
36 1 B 41 R
37 0 B 42 R
37 1 B 43 R
38 0 B 44 R
10
38 1 B 45 R
39 0 B 46 R
39 1 B 47 R
40 B 0 36 L
41 B 1 36 L
42 B 0 37 L
43 B 1 37 L
44 B 0 38 L
45 B 1 38 L
46 B 0 39 L
47 B 1 39 L
36 B B 36 L
37 B B 37 L
38 B B 38 L
39 B B 39 L
36 = = 48 R
37 = = 49 R
38 = = 50 R
39 = = 51 R
48 B 0 -1 R
49 B 1 -1 R
50 B 0 52 R
51 B 1 52 R
52 0 0 52 R
52 1 1 52 R
52 B B 53 L
53 0 B 54 R
53 1 B 55 R
54 B 0 53 L
55 B 1 53 L
53 B B 53 L
53 = = 56 R
56 B 1 -1 R
30 0 X 57 R
30 1 X 58 R
31 0 X 58 R
31 1 X 59 R
57 X X 57 R
57 + + 57 R
57 Y Y 57 R
57 = = 57 R
57 0 0 57 R
57 1 1 57 R
11
58 X X 58 R
58 + + 58 R
58 Y Y 58 R
58 = = 58 R
58 0 0 58 R
58 1 1 58 R
59 X X 59 R
59 + + 59 R
59 Y Y 59 R
59 = = 59 R
59 0 0 59 R
59 1 1 59 R
57 B B 60 L
58 B B 61 L
59 B B 62 L
60 0 B 63 R
60 1 B 64 R
61 0 B 65 R
61 1 B 66 R
62 0 B 67 R
62 1 B 68 R
63 B 0 60 L
64 B 1 60 L
65 B 0 61 L
66 B 1 61 L
67 B 0 62 L
68 B 1 62 L
60 B B 60 L
61 B B 61 L
62 B B 62 L
60 = = 69 R
61 = = 70 R
62 = = 71 R
69 B 0 2 L
70 B 1 2 L
71 B 0 3 L
4 0 0 72 L
4 1 1 72 L
5 0 0 73 L
5 1 1 73 L
6 0 0 74 L
6 1 1 74 L
72 0 0 72 L
12
72 1 1 72 L
73 0 0 73 L
73 1 1 73 L
74 0 0 74 L
74 1 1 74 L
72 + + 75 L
73 + + 76 L
74 + + 77 L
75 X X 75 L
76 X X 76 L
77 X X 77 L
75 0 X 10 R
75 1 X 11 R
76 0 X 11 R
76 1 X 12 R
77 0 X 12 R
77 1 X 13 R
75 Z X 78 R
75 W X 79 R
76 Z X 79 R
76 W X 80 R
77 Z X 80 R
77 W X 81 R
78 X X 78 R
78 + + 78 R
78 0 0 78 R
78 1 1 78 R
78 Y Y 78 R
78 = = 78 R
79 X X 79 R
79 + + 79 R
79 0 0 79 R
79 1 1 79 R
79 Y Y 79 R
79 = = 79 R
80 X X 80 R
80 + + 80 R
80 0 0 80 R
80 1 1 80 R
80 Y Y 80 R
80 = = 80 R
81 X X 81 R
81 + + 81 R
13
81 0 0 81 R
81 1 1 81 R
81 Y Y 81 R
81 = = 81 R
78 B B 82 L
79 B B 83 L
80 B B 84 L
81 B B 85 L
82 0 B 86 R
82 1 B 87 R
83 0 B 88 R
83 1 B 89 R
84 0 B 90 R
84 1 B 91 R
85 0 B 92 R
85 1 B 93 R
86 B 0 82 L
87 B 1 82 L
88 B 0 83 L
89 B 1 83 L
90 B 0 84 L
91 B 1 84 L
92 B 0 85 L
93 B 1 85 L
82 B B 82 L
83 B B 83 L
84 B B 84 L
85 B B 85 L
82 = = 94 R
83 = = 95 R
84 = = 96 R
85 = = 97 R
94 B 0 98 L
95 B 1 98 L
96 B 0 99 L
97 B 1 99 L
98 = = 98 L
99 = = 99 L
98 Y Y 98 L
99 Y Y 99 L
98 0 Y 78 R
98 1 Y 79 R
99 0 Y 79 R
14
99 1 Y 80 R
98 + + -1 R
99 + + 100 R
100 Y Y 100 R
100 = = 100 R
100 0 0 100 R
100 1 1 100 R
100 B B 101 L
101 0 B 102 R
101 1 B 103 R
102 B 0 101 L
103 B 1 101 L
101 B B 101 L
101 = = 104 R
104 B 1 -1 R
このプログラムをコピーします。
初期データをセットします。
15
実行すると
となります。
分かり易い別の例をやってみます。初期データをセットします。
16
実行すると
となります。
このように数値計算も努力すればプログラミングできます。しかし複雑なプログラムを作るのは
大変です。私は大きな模造紙にプログラムの流れを描いて、チューリング・マシンの動きを図示し
ながら、丸一日かかってこのプログラムを作り上げました。チューリング・マシンは実際にプログ
ラミングする機械ではなく、計算機の原理を示す機械です。
問題: 2 進数の 2 の補数を計算するプログラムを作れ。
次は現実の計算機を考察します。
17
2
機械語
和田秀男著「計算数学」朝倉書店 第二章「機械語」に載っている単純なコンピュータの仕組み
をとうして、機械語とアセンブリ言語を見てみよう。
このコンピュータにはメモリ(記憶装置)として8192語(1語=16ビット)があり、演算
装置があり、演算装置で計算した結果を記憶する16ビットからなるアキュムレータと次にメモリ
の何処を実行するかを示す16ビットからなるプログラムカウンタがある。上のソフトには演算装
置は表示されていない。アキュムレータとプログラムカウンタが左側に、メモリが中央にある。
メモリにはコンピュータのプログラムとデータを入れる。コンピュータのプログラムをメモリに
入れる方式がプログラム内蔵方式といわれ、フォン・ノイマンが考案したと言われている。
このコンピュータの命令はオペコードと呼ばれる命令の種類を表す数とオペランドと呼ばれるメ
モリの番地を示す数の組で表される。
例えば15と−7の和をこのコンピュータで計算したいときは次のようにする。
TinyCom を起ち上げる。
「メモリーの修正」ボタンを使って、メモリにデータをセットしていく。「メモリーの修正」ボ
タンの右のエディットに「0000:0010000000000100」と入力する。ここで、: の左の4桁が10進
数で番地を示し、: の右の16桁が2進数でその番地の内容を示している。
18
「メモリーの修正」ボタンをクリックする。
メモリの一番上番地が「0000:0010000000000100」と変化した。ここで、: の左の4桁が10進数
で番地を示し、: の右の16桁が2進数でその番地の内容を示している。
「メモリーの修正」ボタンの右のエディットに「0001:0110000000000101」と入力する。
「メモリーの修正」ボタンをクリックする。
19
メモリの2番目番地が「0001:0110000000000101」と変化した。
この操作を繰り返し、メモリが
0000:0010000000000100
0001:0110000000000101
0002:0100000000000110
0003:1010000000000000
0004:0000000000001111
0005:1111111111111001
0006:0000000000000000
となるようにする。
最初の4個が命令の列で、残りがデータです。データは16ビットで整数を表し、正の整数は最
初のビットを 0 とし残り15ビットを使って2進数で表します。負の整数は2の補数表示という
方式を使っていて、例えば -7 を表すには次のように変形して表現を求めます。
0000000000000111
1111111111111000
7 を2進数表示する
0 と 1 を入れ替える
1111111111111001
1 を足す
20
負の整数はこの操作で表現を求めます。これを2の補数表示といいます。2の補数表示の利点は
数の表示が一意になることです。つまり、最初のビットが0なら正、最初のビットが1なら負と
し、残りの15ビットで絶対値を表現することにすると0の表現がの表現が 0000000000000000 と
1000000000000000 と二通り出来てしまい不便です。
このプログラムの意味は
0000:0010000000000100
4 番地のデータをアキュムレータに入れる
0001:0110000000000101
0002:0100000000000110
0003:1010000000000000
5 番地のデータをアキュムレータに加える
6 番地にアキュムレータの内容を保存
停止する
0004:0000000000001111
0005:1111111111111001
15 を表す
-7 を表す
0006:0000000000000000
0 を表す
です。
このプログラムを実行するには、プログラムカウンタが 0000000000000000 になっていることを
確かめて、スタートのボタンをクリックします。
6 番地が 0000000000001000 となりました。すなわち、15 + −7 = 8 です。プログラムカウンタ
が 0000000000000011 になりました。停止しろという命令で止まりました。
6 番地を 0000000000000000 にもどし、プログラムカウンタも 0000000000000000 に戻し、
21
ステップのボタンをクリックします。
アキュムレータが 0000000000001111 になり、プログラムカウンタが 0000000000000001 になりま
した。0番地の命令が実行されました。ステップのボタンをクリックするとプログラムカウンタが
指している番地の命令が1つ実行されます。
このコンピュータの命令は7つしかないから1から7までの通し番号が付けられていて、2進法
で表せば3ビットで表せ、メモリの番地は0から8191(= 213 − 1)までなので13ビットで
表せる。それで、1語16ビットのうち左3ビットに命令の番号を入れ、右の残り13ビットにメ
モリの番地を入れることにより、1つの命令を1語で表現している訳です。
コンピュータにとってはこれでよいが、人間には2進数の数字の羅列を読むのは辛いものがある
ので、人間に分かり易い言語としてアセンブリ言語が作られた。オペコードやオペランドを英単語
や文字で表します。オペコードをアセンブリ言語で表したのが次の表の load, store, add, subtract,
stop, jump および jump minus である。番地は文字や文字列で表す。
このコンピュータのオペコードは7つで、次の表の様になっている。
10進数
2進数
アセンブリコード
説明
1
2
001
010
load
store
load M で M 番地の内容をアキュムレータにコピー
store M でアキュムレータの内容と M 番地にコピー
3
011
add
add M で M 番地の内容とアキュムレータの内容を
加え、結果をアキュムレータにコピー
4
100
subtract
subtract M でアキュムレータの内容から M 番地
の内容を引き、結果をアキュムレータにコピー
5
6
101
110
stop
jump
stop でコンピュータの動作は停止する
jump M で無条件で M の番地へ移動
7
111
jump minus
jump minus M でアキュムレータが負の時
M の番地へ移動
アセンブリ言語では上のプログラムを次のように表現する。
22
ラベル
オペコード
オペランド
対応する機械語
load
A
0000:0010000000000100
add
store
stop
B
C
0001:0110000000000101
0002:0100000000000110
0003:1010000000000000
A
B
15
-7
0004:0000000000001111
0005:1111111111111001
C
0
0006:0000000000000000
アセンブリ言語の設計はある程度の自由度があるが、和田秀男著「計算数学」朝倉書店 第二章
「機械語」に載っている単純なコンピュータのアセンブリ言語、および TinyCom のアセンブリ言
語の命令は
ラベル オペコード オペランド
からなり、ラベルやオペランドは空白のこともある。ラベルが空白でもラベル欄には半角スペース
を入れる。ラベルを記入するときはその前に空白を入れてはいけません!
!
!
TinyCom を起ち上げ、
右側のメモを見るとアセンブリ言語のプログラムらしきものが書いてある。これはこのコンピュー
タの7つの命令をすべて使ったプログラムを書いているだけで無意味なプログラムです。これを消
して、上のアセンブリ言語のプログラムを打ち込む。文字も空白もすべて半角文字や半角スペース
を使います。ラベルを記入するときはその前に空白を入れてはいけません!
!
!
23
右下にある「アセンブル」ボタンをクリックします。
メモリに右側のアセンブリ言語のプログラムを機械語に翻訳されたプログラムがセットされまし
た。アセンブリ言語のプログラムを機械語に翻訳するプログラムをアセンブラと言います。コン
ピュータを使うのが随分楽になりました。今後は TinyCom のアセンブリ言語のプログラムを調べ
ます。
足し算をしましたが、引き算をするには subtract 命令を使います。subutract M と命令すれば、
アキュムレータから M とラベル付けされた番地の内容を引き、その結果をアキュムレータに入れ
ます。例えば、A + B − C を D に入れる (D ← A + B − C) にはアセンブリ言語で次のようにプ
ログラミングします。
24
ラベル
オペコード
オペランド
対応する機械語
load
A
0000:0010000000000101
add
subtract
store
B
C
D
0001:0110000000000110
0002:1000000000000111
0003:0100000000001000
A
stop
6
0004:1010000000000000
0005:0000000000000110
B
C
D
9
-8
0
0006:0000000000001001
0007:1111111111111000
0008:0000000000000000
ラベル A, B, C, D を付けられた番地には好きな値を入れて下さい。今の場合は、6 + 9 − (−8)
を計算しています。アセンブルし
スタートボタンをクリックすると
ラベル D で指示された8番地が 0000000000010111 になりました。6 + 9 − (−8) = 23 を計算し
ました。
問題:A ← 2A − B − C のプログラムを作れ。
A番地の内容を16倍して、その結果を A 番地に保存するには次のようにプログラミングします。
25
ラベル
A
オペコード
オペランド
load
A
add
store
add
A
A
A
store
add
A
A
store
add
store
A
A
A
stop
-3
A には好きな数字を入れます。add を15回使うよりはプログラムが短くなっています。
A 番地の内容と B 番地の内容を入れ替えるには補助の番地 C を使って、次のようにします。
ラベル
オペコード
オペランド
load
A
store
load
C
B
store
load
store
A
C
B
A
stop
27
B
C
-18
0
A, B, C には好きな数字を入れます。これを
ラベル
A
B
オペコード
オペランド
load
store
A
B
load
store
B
A
stop
27
-18
としては駄目です。
問題:A ← B ← C ← A のプログラムを作れ。つまり、3つの値を入れ替えるプログラムを
作れ。
jump minus M という命令は、アキュムレータが負の時、プログラムの行の先頭に M と書か
れているところ(ラベルが M である命令の行)を次に実行しなさいと言う、条件分岐命令です。
26
jump M という命令は、無条件にプログラムの行の先頭に M と書かれているところ(ラベルが M
である命令の行)を次に実行しなさいと言う、無条件分岐命令です。これらの命令を使うと A 番
地の内容の絶対値を B 番地に保存する(B ← |A|)プログラムは次のように書くことが出来ます。
ラベル
オペコード
オペランド
L1
load
jump minus
store
A
L0
B
stop
load
subtract
Zero
A
L0
Zero
jump
0
A
B
-95
0
L1
ラベル Zero 番地は 0 でなければなりませんが、A, B 番地は好きな値で良いです。
A番地の内容が -95 の様に負の数であれば、load A を実行すれば、アキュムレータに負の数が
セットされるので次の jump minus L0 の命令で、L0 とラベル付けされた番地の命令 load Zero を
次に実行します。一方、A番地の内容が 95 の様に正の数であれば、load A を実行すれば、アキュ
ムレータに正の数がセットされるので次の jump minus L0 の命令は無視され、次の store B の命
令が次に実行されます。また jump L1 の命令を実行すると次は L1 とラベル付けされた番地の命
令 store B が次に実行されます。このように jump minus と jump だけがプログラムの流れを変え
ます。
A 番地の値と B 番地の値の小さい方を C 番地に保存するプログラムは次のようになります。
ラベル
オペコード
オペランド
load
A
subtract
jump minus
load
B
L0
B
L1
store
stop
C
L0
A
L1
A
load
jump
-9
B
-12
問題:A 番地の値と B 番地の値の大きい方を C 番地に保存するプログラムを作れ。
このように色々なことが出来ます。もっと複雑なプログラミングの例が知りたければ和田秀男著
「計算数学」朝倉書店 第二章「機械語」を読んで下さい。かけ算や割り算やユークリッドの互除
法による最大公約数を求めるプログラムやサブルーチンのプログラムまで載っています。
現実のコンピュータは構造が非常に複雑です。Knuth の THE ART OF COMPUTER PRO-
GRAMMING で説明されている仮想のコンピュータは次の図の様になっています。
27
このコンピュータで 500 個の素数を求めるプログラムは
* 500 個の素数の表の印刷 PRIME.PRO
*
L
EQU 500
PRINTER EQU 18
これから探す素数の個数
PRIME
BUF0
BUF1
EQU -1
EQU 800
EQU BUF0+25
素数表の記憶領域
BUFFER[0] の記憶領域
BUFFER[1] の記憶領域
START
ORIG 600
IOC 0(PRINTER)
改頁
2H
LD1 =1-L=
LD2 =3=
INC1 1
P1. 表の作成開始. J<--1.
N<--3.
P2. N は素数, J<--J+1.
ST2 PRIME+L,1
J1Z 2F
PRIME[J]<--N.
P3. すでに 500 個見つけたか.
INC2 2
ENT3 2
ENTA 0
P4. N を進める.
P5. K<--2.
P6. PRIME[K] | N か.
ENTX 0,2
DIV PRIME,3
JXZ 4B
R=0 か.
CMPA PRIME,3
INC3 1
P7. PRIME[K] の方が大きいか.
P8. K を進める.
JG 6B
JMP 2B
OUT TITLE(PRINTER)
さもなければ N は素数.
P9. 見出しの印刷
ENT4 BUF1+10
ENT5 -50
Set B<--1.
Set M<--0.
INC5 L+1
LDA PRIME,5
M を進める.
P10. 一行分の印刷内容を入れる.
4H
6H
2H
2H
4H
印字機の機番
Q>PRIME[K] なら飛び越し,
28
CHAR
STX 0,4(1:4)
DEC4 1
DEC5 50
(正の範囲で rI5 の内容が 50 個ずつ減る)
J5P 4B
OUT 0,4(PRINTER)
LD4 24,4
P11. 一行印刷.
緩衝領域の切替.
J5N 2B
HLT
rI5=0 なら完了.
* INITIAL CONTENTS OF TABLES AND BUFFERS
ORIG PRIME+1
CON 2
最初の 2.
TITLE
ORIG BUF0-5
ALF FIRST
表題の文字情報.
ALF FIVE
ALF HUND
ALF RED P
ALF RIMES
ORIG BUF0+24
CON BUF1+10
ORIG BUF1+24
CON BUF0+10
各緩衝領域が他の緩衝領域を参照.
END START
終わり.
となります。これを MIX にロードし
アセンブルし
29
実行し
Printer を表示すると
となります。別のプログラムの例は
ORIG 0
30
X
Y
CON 24
CON 16
ORIG 20
BUFFER CON 0
START
L1
ORIG 100
IOC 0(18)
ENTA 0
LDX X
DIV Y
JXZ L2
LDA Y
STA X
STX Y
JMP L1
L2
LDA Y
CHAR
STX BUFFER
OUT BUFFER(18)
HLT
END START
です。これは24と16の最大公約数を計算しています。
とプログラムを打ち込み、
31
アセンブルし、
実行し、Printer を表示すると
となります。
これも今となっては古いコンピュータです。MIX のプログラミングは、Donard E. Knuth 「The
Art of Computer Programming Volumes 1-4A」 にあります。この本はコンピュータの世界では
32
最も有名な本です。翻訳もあります。アルゴリズムを集大成した本で、計算機科学者のバイブルで
す。Donard E. Knuth は、世界中の人が使っている複雑な数式や化学式を含む文章や楽譜や色々
な言語の文章を作るためのソフト TeX の作者としても有名です。
現実に作られたパソコンの CPU の解説は、「はじめて読む8086」アスキー出版局。198
7や「はじめて読む486」アスキー出版局。1994がありますが、これも古いです。昔はこん
な本を読んで、私もコンピュータの勉強をしましたが、今では機械語やアセンブリ言語の勉強をす
る必要はほとんどなくなりました。唯、機械語を知っていれば、C 言語のポインタやコアやスレッ
ドの仕組みがそれなりに理解できます。
次は、昔私が作った Borland の Turbo Assembler TASM 5.0J による 8 人の女王の問題を解
くプログラムです。左側が 8086 の機械語で右側がアセンブリ言語によるプログラムです (勿論
PENTIUM でも同じアセンブリ言語のプログラムで動きます) 。対称性は考慮していません。 8
人の女王の問題を拡張した N 人の女王の問題を 8 人以下で解けます。
1 0000
2 0000
.MODEL
.STACK
SMALL
200H
3 0000
4 0000
5 0009
09*(??)
09*(00)
MM
JJ
.DATA
DB
DB
9 DUP(?)
9 DUP(0)
6 0012
7 0022
10*(00)
10*(00)
KK
LL
DB
DB
16 DUP(0)
16 DUP(0)
????
????
????
I
N
M
DW
DW
DW
?
?
?
8 0032
9 0034
10 0036
11 0038
12 0000
.CODE
QUEENSSTART:
13 0000
14 0003
15 0005
B8 0000s
8E D8
B4 02
MOV
MOV
MOV
AX,@DATA
DS,AX
AH,2
16 0007
17 0009
18 000B
B2 3F
CD 21
B4 01
MOV
INT
MOV
DL,’?’
21H
AH,1
19 000D
20 000F
CD 21
2C 30
INT
SUB
21H
AL,’0’
21 0011
22 0013
23 0016
B4 00
A3 0036r
B4 02
MOV
MOV
MOV
AH,0
[M],AX
AH,2
24 0018
25 001A
B2 0D
CD 21
MOV
INT
DL,0DH
21H
26 001C
27 001E
28 0020
B2 0A
CD 21
C7 06 0032r 0001
MOV
INT
MOV
DL,0AH
21H
[I],1
29 0026
30 002C
C7 06 0034r 0000
FF 06 0034r
MOV
INC
[N],0
[N]
L3:
L4:
33
31 0030
32 0033
A1 0034r
3B 06 0036r
MOV
CMP
AX,[N]
AX,[M]
33 0037
34 0039
77 7F
BB 0009r
JA
MOV
L6
BX,OFFSET JJ
35 003C
36 0040
37 0043
03 1E 0034r
80 3F 01
74 E7
ADD
CMP
JE
BX,[N]
BYTE PTR [BX],1
L4
38 0045
39 0048
BB 0012r
03 1E 0032r
MOV
ADD
BX,OFFSET KK
BX,[I]
40 004C
41 0050
42 0054
2B 1E 0034r
03 1E 0036r
80 3F 01
SUB
ADD
CMP
BX,[N]
BX,[M]
BYTE PTR [BX],1
43 0057
44 0059
74 D3
BB 0022r
JE
MOV
L4
BX,OFFSET LL
45 005C
46 0060
47 0064
03 1E 0032r
03 1E 0034r
4B
ADD
ADD
DEC
BX,[I]
BX,[N]
BX
48 0065
49 0068
80 3F 01
74 C2
CMP
JE
BYTE PTR [BX],1
L4
50 006A
51 006D
52 0070
A1 0034r
BB 0000r
03 1E 0032r
MOV
MOV
ADD
AX,[N]
BX,OFFSET MM
BX,[I]
53 0074
54 0076
55 0079
88 07
BB 0009r
03 1E 0034r
MOV
MOV
[BX],AL
BX,OFFSET JJ
56 007D
57 0080
C6 07 01
BB 0012r
ADD
MOV
MOV
BX,[N]
BYTE PTR [BX],1
BX,OFFSET KK
58 0083
59 0087
60 008B
03 1E 0032r
2B 1E 0034r
03 1E 0036r
ADD
SUB
ADD
BX,[I]
BX,[N]
BX,[M]
61 008F
62 0092
C6 07 01
BB 0022r
MOV
MOV
BYTE PTR [BX],1
BX,OFFSET LL
63 0095
64 0099
65 009D
03 1E 0032r
03 1E 0034r
4B
ADD
ADD
DEC
BX,[I]
BX,[N]
BX
66 009E
67 00A1
C6 07 01
FF 06 0032r
MOV
INC
BYTE PTR [BX],1
[I]
68 00A5
69 00A8
70 00AC
A1 0032r
3B 06 0036r
77 03
MOV
CMP
JA
AX,[I]
AX,[M]
L5
71 00AE
72 00B1
E9 FF75
8B 0E 0036r
JMP
MOV
L3
CX,[M]
73 00B5
E8 004B
CALL
PRINTRESULT
L5:
34
74 00B8
75 00BC
FF 0E 0032r
83 3E 0032r 00
76 00C1
77 00C3
DEC
CMP
[I]
[I],0
74 3C
BB 0000r
JE
MOV
ENDQUEENS
BX,OFFSET MM
78 00C6
79 00CA
80 00CC
03 1E 0032r
8A 07
32 E4
ADD
MOV
XOR
BX,[I]
AL,[BX]
AH,AH
81 00CE
82 00D1
A3 0034r
BB 0009r
MOV
MOV
[N],AX
BX,OFFSET JJ
83 00D4
84 00D8
85 00DB
03 1E 0034r
C6 07 00
BB 0012r
ADD
MOV
MOV
BX,[N]
BYTE PTR [BX],0
BX,OFFSET KK
86 00DE
87 00E2
03 1E 0032r
2B 1E 0034r
ADD
SUB
BX,[I]
BX,[N]
88 00E6
89 00EA
90 00ED
03 1E 0036r
C6 07 00
BB 0022r
ADD
MOV
MOV
BX,[M]
BYTE PTR [BX],0
BX,OFFSET LL
91 00F0
92 00F4
03 1E 0032r
03 1E 0034r
ADD
ADD
BX,[I]
BX,[N]
93 00F8
94 00F9
95 00FC
4B
C6 07 00
E9 FF2D
DEC
MOV
JMP
BX
BYTE PTR [BX],0
L4
96 00FF
97 00FF
98 0101
L6:
ENDQUEENS:
MOV
AH,4CH
BB 0000r
INT
PRINTRESULT
MOV
21H
PROC
BX,OFFSET MM
101 0106
102 0106
103 0107
43
B4 02
PRINTLOOP:
INC
MOV
BX
AH,2
104 0109
105 010B
8A 17
80 C2 30
MOV
ADD
DL,[BX]
DL,’0’
106 010E
107 0110
108 0112
CD 21
B2 20
CD 21
INT
MOV
INT
21H
DL,’ ’
21H
109 0114
110 0116
E2 F0
B4 02
LOOP
MOV
PRINTLOOP
AH,2
111 0118
112 011A
113 011C
B2 0D
CD 21
B2 0A
MOV
INT
MOV
DL,0DH
21H
DL,0AH
114 011E
115 0120
CD 21
C3
INT
RET
21H
99 0103
100 0103
116 0121
B4 4C
CD 21
PRINTRESULT
35
ENDP
117
END
QUEENSSTART
プログラミングの歴史を「闘うプログラマー」[上] G ・パスカル・ザカりー著山岡洋一訳 日経
BP 出版センターから引用してみよう。
「初期のコンピュータ・プログラムは、特定のコンピュータを補助するものにすぎなかった。
第二次世界大戦以前には、機械式の計算機がほとんどで、スイッチを切り替えたり、配線を変えた
り、歯車をうごかすのがプログラミングの仕事にあたっていた。1930 年代には、当時の最強の機
械式計算機、微分解析機で新しい問題を解く準備をするのに、何日もかかった。それから十年たっ
て、デジタル・コンピュータができたが、むずかしい問題を解くには、セットアップに何日もか
かった。 パンチ・カードや紙テープに記録された要求を読みとれるコンピュータができたが、カー
ドやテープは、人間の手で機械に送り込まなければならなかった。プログラミングはまだ初期の段
階にあり、進歩が必要だった。 画期的な進歩がおとずれたのは、1944 年であった。ハンガリー生
まれのアメリカの数学者、ジョン・フォン・ノイマンが、プログラム記憶方式という考えを発展さ
せた。プログラムを記憶させるという考え方自体は、この分野ではだれでも考えていたことだが、
ノイマンはとくにその重要性をはっきりととらえていた。プログラム記憶方式では、コンピュータ
への命令は、コンピュータ自体のメモリーに保存され、データとおなじようにあつかわれる。こ
うすれば、ごく短時間でプログラムを立ち上げられ、プログラムを修正したり、いくつかのプログ
ラムを切り替えるのも簡単になる。 プログラム記憶方式は、初期のコンピュータ文化の中に広が
り、プログラミングが爆発的に増え、この方式を支持する人が急増した。しかし、道のりはけわ
しかった。デジタル・コンピュータには、オンかオフの二つの状態しかないため、 1 (オン)と
0 (オフ)でできた二進法のメッセージにしか応答できない。 プログラムはすべて、この二つの
数字だけで表さなければならないので、 ごく普通の数式が、すぐに、目がくらむほど複雑になる。
1940 年代後半には、 コンピュータのプログラミングは「気が狂うほどむずかしい」と言われた。
まもなく、二進法の文字列を、もっと簡単につくる方法が考え出された。最初に発明されたのは、
二進コードを自動的に打ち出す特殊なタイプライターだ。つぎに、もっと応用範囲の広い「アセ
ンブリー」言語ができた。これは、文字や記号で 0 と 1 を表すものだ。アセンブリーでプログラ
ムを書くようになったのは進歩だが、それでも、融通のきかないコンピュータの命令セットに忠実
に従って書く必要がある。アセンブリー・コードを書きこなすには、命令セットを完璧に知ってい
なければならない。それに、命令セットは、プロセッサーの設計に依存するため、コンピュータの
機種によってちがっている。つまり、苦労して覚えたアセンブリー言語の知識も、そのコンピュー
タが使われなくなれば、無駄になる可能性がある。 1950 年代初期には、コンピュータの利用が多
い組織、とくにアメリカの陸海空軍は、ソフトウエアが頭痛のタネであり、コストがかかることに
気づくようになった。一流のプログラマーは、プログラムを書きやすくする方法を模索していた。
1951 年、アメリカ海軍予備役補給局の数学者、グレース・マレー・ホッパーが、コンパイラを開
発した。コンパイラは、プログラマーの命令を、コンピュータを制御する 1 と 0 の文字列、つま
りマシン語に翻訳するプログラムである。コンパイラで、ハードウエアの束縛や頭の痛い二進コー
ドから、プログラマは解放されるようになった。 ホッパーが突破口を開いた後、プログラミング
を簡単にするために、さまざまな努力が行われるようになった。その中でおそらくもっとも重要な
ものは、 IBM が開発したコンパイラ、FORTRAN だろう。FORTRAN には、PUNCH, READ
DRUM, IF DIVIDE CHECK など三十二種類の命令があり、それぞれに、コンピュータに必要な
正確なマシン語が対応している。1950 年代後半、FORTRAN の影響力は絶大であった。あるコン
ピュータ史研究家はこう述べている。「論理的に考えられる人なら、だれでも、やる気があればコ
36
ンピュータのプログラミングを学べるようになった。コンピュータの内部構造やむずかしいアセン
ブリー言語を熟知したスペシャリストである必要はなくなった。FORTRAN の簡単なコマンドだ
けを使えば、コンピュータを命令どおりにうごかせる。コンパイラーが自動的に、命令を正しい
マシン・コードに翻訳してくれる」 FORTRAN によって、おなじ命令セットを使って、何種類も
のコンピュータ用にプログラミングできるようになったが、種類のちがうマシンで使うには、プロ
グラムを変更しなければならないことが多かった。それに、FORTRAN は科学技術の問題を解決
する事を目的としている。他の用途のために、COBOL などの言語が開発された。やがて、コン
ピュータ言語のジャングルが、プログラマにとって、新たな頭痛のタネになった。どの言語を得意
とするかによって、キャリアが決まるようになった。」
私がコンピュータを勉強しだした40年ぐらい前は、コンピュータを勉強することは FORTRAN
(FORmulaTRANslater)を勉強することでした。最初に実用化されたプログラミング言語は FOR-
TRAN ですが、これは数値計算のための言語で、私が最初に勉強し、20 年間使っていた言語です
が、今では古い言語です。宇宙開発のソフトなど開発している「真のプログラマ」は FORTRAN
とアセンブリ言語であらゆるプログラムを作るといいます。
FORTRAN は科学技術計算に使われ、事務計算には COBOL という言語が使われていました。
LISP (LISt Processor) は 1960 年ごろ J.McCarthy によって開発されたリストを基本データ構造と
する関数型の言語です。LISP は人工知能のためのアセンブラといわれリスト処理に優れています。
論理型のプログラミング言語 PROLOG (PROgramming in LOGic) は 1972 年に A. Colmerauer
によってフランスで開発されます。他の言語のように手続きを記述するのではなく、事実(物事の
論理的関係)を記号論理に基づいて記述し、パターン・マチィングとバックトラッキングでプログ
ラムを実行します。LISP と同じくリスト処理に優れ人工知能のための言語で、日本の第五世代コ
ンピュータの開発用言語としても使われました。しばらく忘れ去られていましたが、Bruce A. Tate
著「7 つの言語 7 つの世界」 Ruby, Io, Prolog, Scala, Erlang, Clojure, and Haskell 平成 23 年
3200 円の影響で再び注目されています。
1975 年 1 月、米「ポピュラー・エレクトロニクス」誌で紹介された8ビット・マイコンキットで
ある「アルテア」が経営難に陥っていた MITS 社から通信販売されます。397 ドルのキットに購入
申し込みが殺到し、1年分と見込んだ個数を初日で突破します。非常に高価で、専門教育を受けた
一部のエリートしか触ることのできなかったコンピュータを個人が我が手にすることができるとい
うことで、医者や弁護士といった人たちが競って購入したそうです。若き日のビル・ゲイツと友人
のポール・アレンがこのアルテア用の BASIC インタープリタを開発します。これが圧倒的人気を
得ます。これより一般の人がマイコン(当時はパソコンのことをマイコンと言っていました)で使
うプログラミング言語は BASIC になりました。高等学校の数学 A,B,C の教科書にも BASIC
のプログラムが大量に載っていましたが、BASIC の処理系もなくなり、センター試験のコンピュー
タの問題も実際の BASIC の処理系でなく、 JIS 規格を基準として問題を作るようになり、情報と
いう教科の導入も相まって、指導要領が改訂され、数学の教科書から BASIC は無くなりました。
BASIC の次にマイコンのプログラミング言語として出てきたのは、GOTO 文有害論が出て、そ
れの対策として開発された構造的プログラミング言語 PASCAL のマイコン版 PASCAL(TURBO
PASCAL) です。大学では、PASCAL でプログラミング教育をしていた時代もありますが、純粋
な PASCAL の処理系も現在では無く(FREE PASCAL と言うのがありますが、オブジェクト指
向言語の DELPHI の代用品です)、ワークステーションの言語として生まれ、すべてのコンピュー
タに普及した C を最初の言語として教えているところが日本の大学では多いです。アメリカの大
学では、例えば、 MIT では LISP の方言である Scheme で教育しています。構造的プログラミン
グ言語の次に、オブジェクト指向言語と言われる SMALLTALK や C++ や Java が出てきます。
37
計算機でできることはすべて C または C++ でプログラミングできますが、コンピュータの性能
が良くなったので、インターネットがらみで、Perl や Ruby や Python や Erlang や Clojure や
Scala や Haskell といた言語ができます。
この授業では、中学校・高等学校の数学の先生になるためのプログラミング教育なので、本格的
なプログラミング教育ではなく、代表的なプログラミング言語達の紹介と数学 A, B, C の教科書
に載っていた BASIC のプログラムのアルゴリズムが理解でき、どれかのプログラミング言語でプ
ログラミングできることを目指します。面白いプログラミング言語が見つかれば自分で勉強して下
さい。
現在の主流のプログラミング言語には手続き型プログラミング言語、関数型プログラミング言
語、論理型プログラミング言語、オブジェクト指向プログラミング言語があります。オブジェクト
指向プログラミング言語は難しいので、手続き型プログラミング言語の代表として C 言語、関数
型プログラミング言語の代表として Scheme、論理型プログラミング言語の代表として Prolog の
基本を学びます。
昔、ネットワークについて話を聞いたり、Java によるネットワークプログラミングなどを勉強
したときは、いくらやっても全然頭に勉強したことが残りませんでした。そのときには、私の環
境では物理的に複数のコンピュータをネットワークに繋げるわけでもなく、大学のコンピュータで
ホームページを作っても極端に制限された使い方しかできなかったので、ネットワークの勉強は私
には本当には必要なかったからだと思います。最近、私のパソコンがマルチスレッドに対応しだし
(CPU が i3, i5 や i7 になった)、家のネットワーク環境にルータが導入され、複数のパソコンを有
線や無線でインターネットに繋げるように知らない間になっていて(pikara がつながるようにし
てくれていた)、やっとマルチスレッドやネットワークのプログラミングを実際に実験しながら勉
強する物理的な勉強の環境が揃ってきました。金高堂でたまたま手に取った「C 言語逆引きハンド
ブック」という本にマルチスレッドやネットワークのプログラミングの基本が書いてありました。
何の勉強でもそうですが、特にコンピュータの勉強は、唯本を読むだけの抽象的な勉強ではなく、
「実際にプログラムを打ち込み、動かし、失敗し、何処が間違いか調べ、修正し、実際にプログラ
ムを打ち込み、動かし、失敗し、何処が間違いか調べ、修正し、
・
・
・」の試行錯誤を繰り返しなが
ら、集中的に時間を掛けて勉強する必要があります。そうすれば、次に何を勉強すればよいか自然
に分かります。現在は色々なことを勉強するのに本当に良い時代になりました。ただ、猛烈な勢い
で世の中変わっています。勉強をし続けないと取り残されてしまいます。
38
3
Windows フォーム・アプリケーション
通常の C++ の解説と異なり、結果が正しいかどうか見れば分かる、関数のグラフを描くことか
ら始めます。
まず、Microsoft Visual Studio 2010 を立ち上げる。
スタートページがなければ、「表示」のメニューでスタートページをクリックすればよい。
上図は Professional 版であることを示していますが、それは私が幾何学の研究の為のプログ
ラムや数独など各種パズルを解いたり問題を作ったりするプログラムやオセロや囲碁などの各種
ゲームのプログラムを作るためにメモリーを沢山使う64ビットのプログラミングがしたいから、
Professional 版を使っているからです。学生さんはアカデミック版を購入できますから、1 万円位で
生協で Professional 版を購入出来ますが、以下のプログラムは Professional 版固有の機能は使って
いないので、学生さんが勉強のために、以下のプログラムを実行するだけなら、無料の Microsoft
Visual Studio 2010 Express 版で十分です。
新しいプロジェクト・
・
・をクリックする。
Windows フォームアプリケーションをクリックし、図では場所 (L):が D:\jyugyou\となっていま
す(自分で適当なディレクトリ、例えば D:\jyugyou\を作り、ここを指定すれば、プログラムが
D:\jyugyou\に作られます)が、何処にプログラムが保存されるか気にしなければ、そのままで良
いですので、名前 (N): のエディットボックスに例えば en2 と入力し、
39
OK ボタンをクリックする。
マウスで Form1 を適当な大きさにする。
ツールボックスのメニューが無ければ、表示のメニューから
40
ツールボックスを選択して、ツールボックスを表示する。
ツールボックスにある PictureBox をクリックし、マウスを Form1 に移動し、ドラッグして、Form1
の中央に Picture Box を配置する。
プロパティが表示されていなければ
41
または、ツールボックスの欄の右端のアイコンを押して
で、プロパティを探し、プロパティを表示する。
プロパティで、雷のアイコンをクリックして、イベントを表示する。
42
一番下にある paint をダブルクリックする。
必要なら Form1.h をドラッグして、左側に移す。
Form1.h の一番下にある
43
private: System::Void pictureBox1_Paint(System::Object^
System::Windows::Forms::PaintEventArgs^ e) {
sender,
}
の { と } の間に、次のようになるようキーボードから打ち込む。
これは picturebox1 のイベント paint が発生したとき、どのように反応すべきかを指示するプロ
グラムを作ることになる。
private: System::Void pictureBox1_Paint(System::Object^
System::Windows::Forms::PaintEventArgs^ e) {
sender,
Graphics^ g = e->Graphics;
Pen^ pen = gcnew Pen(Color::Blue);
double pi = Math::PI;
g->DrawEllipse(pen, 50, 50, 260, 260);
for (double t=0; t<2*pi; t += pi/80) {
g->DrawLine(pen, (int)(180+130*cos(t)), (int)(180-130*sin(t)),
(int)(180+130*cos(2*t)), (int)(180-130*sin(2*t)));
}
}
44
Form1.h ファイルの先頭
#include <math.h>
を追加する。
コンパイルし、実行する(Debug の左の赤い三角形のアイコンをクリックする)と
45
となる。図が欠けていれば、Form1.h[デザイン] のタブをクリックし、Form や PictureBox を
適当な大きさに修正し、再度 Debug の左の赤い三角形のアイコンをクリックする.。
別のプログラムを作ってみます。
Microsoft Visual Studio 2010 を立ち上げ、上と同じように適当な名前でプロジェクトを作り、
Form に PictureBox と Button 2個を配置する。
Form1.h ファイルの先頭に
#include <math.h>
を追加する。
button1 をダブルクリックし、
private: System::Void button1_Click(System::Object^
System::EventArgs^ e) {
sender,
Graphics^ g = pictureBox1->CreateGraphics();
Brush^ brush = gcnew SolidBrush(Color::White);
g->FillRectangle(brush, 0, 0, 360, 360);
Pen^ pen = gcnew Pen(Color::Blue, 1);
g->DrawEllipse(pen, 50, 50, 260, 260);
double pi = Math::PI;
for (double t=0; t<2*pi; t += pi/80) {
g->DrawLine(pen, (int)(180+130*cos(t)), (int)(180-130*sin(t)),
(int)(180+130*cos(2*t)), (int)(180-130*sin(2*t)));
}
}
と打ち込む。さらに、button2 をダブルクリックし、
private: System::Void button2_Click(System::Object^
System::EventArgs^ e) {
Graphics^ g = pictureBox1->CreateGraphics();
Brush^ brush = gcnew SolidBrush(Color::Yellow);
46
sender,
g->FillRectangle(brush, 0, 0, 360, 360);
Pen^ pen = gcnew Pen(Color::Red, 1);
g->DrawEllipse(pen, 50, 50, 260, 260);
double pi = Math::PI;
for (double t=0; t<2*pi; t += pi/80) {
g->DrawLine(pen, (int)(180+130*cos(t)), (int)(180-130*sin(t)),
(int)(180+130*cos(3*t)), (int)(180-130*sin(3*t)));
}
}
と打ち込む。実行すると
となります。button1 をクリックすると
となり、button2 をクリックすると
47
となります。
(−5π ≤ x ≤ 5π) のグラフを描け。
Form の中央に PictureBox を配置する。
例題:y = x sin x
Form1.h ファイルの先頭に
#include <math.h>
を追加する。
picturebox1 のイベント paint を次のように定義する。
private: System::Void pictureBox1_Paint(System::Object^
System::Windows::Forms::PaintEventArgs^ e) {
Graphics^ g = e->Graphics;
Pen^ pen = gcnew Pen(Color::Blue);
double pi = Math::PI;
g->DrawLine(pen, 0, 180, 360, 180);
g->DrawLine(pen, 180, 0, 180, 360);
int ox = (int)(180-180.0*5*pi/15.7079);
int oy = (int)(180-10*(-5*pi)*sin(-5*pi));
48
sender,
for (double t=-5*pi; t<=5*pi; t += pi/80) {
g->DrawLine(pen, ox, oy,
(int)(180+180.0*t/15.7079), (int)(180-10*t*sin(t)));
ox = (int)(180+180.0*t/15.7079);
oy = (int)(180-10*t*sin(t));
}
}
コンパイルし、実行すると
となる。
例題:(1/2, 0) を中心とする直径1の円上の点を動点として、直径 OP の円群を描け。輪郭線に
カージオイドが浮かび上がる。
Form の中央に PictureBox を配置する。
Form1.h ファイルの先頭に
#include <math.h>
を追加する。
picturebox1 のイベント paint を次のように定義する。
49
private: System::Void pictureBox1_Paint(System::Object^
System::Windows::Forms::PaintEventArgs^ e) {
sender,
Graphics^ g = e->Graphics;
Pen^ pen1 = gcnew Pen(Color::Black);
g->DrawLine(pen1, 0, 180, 360, 180);
g->DrawLine(pen1, 180, 0, 180, 360);
double pi = Math::PI;
int K = 150;
for (double t = 0; t<=2*pi; t += pi/20) {
double x = 0.25+cos(t)/4.0;
double y = sin(t)/4.0;
int x1 = (int)(K*(x-sqrt(0.5+cos(t)/2.0)/2.0)+180);
int y1 = (int)(K*(-y-sqrt(0.5+cos(t)/2.0)/2.0)+180);
int x2 = (int)(K*sqrt(0.5+cos(t)/2.0));
int y2 = (int)(K*sqrt(0.5+cos(t)/2.0));
g->DrawEllipse(pen1, x1, y1, x2, y2);
}
}
コンパイルし、実行すると
となる。
例題:原点を中心とし半径が a の円 O の外側を半径が b の円 C が円 O に外接しながら滑るこ
となく転がるとき、円 C 上の点 P の軌跡を考える。ただし、点 P のはじめの位置は、円 O と x
軸の正の部分との交点 A とする。
円 C が転がるとき、動径 OC が表す角を θ とし、そのときの外接する点を Q、点 P の座標を
(x,y) とすると、P の軌跡は θ を媒介変数として
x = (a + b) cos θ − b cos
a+b
a+b
θ, y = (a + b) sin θ − b sin
θ
b
b
で表される。この曲線を エピサイクロイド という。エピサイクロイドを描くプログラムを作れ。
Form に PictureBox と Button と Label 2個と textBox 2個を配置する。
50
Form1.h ファイルの先頭に
#include <math.h>
を追加する。
button1 をダブルクリックし、
private: System::Void button1_Click(System::Object^
sender,
System::EventArgs^ e) {
Graphics^ g = pictureBox1->CreateGraphics();
double a = System::Convert::ToDouble(textBox1->Text);
double b;
Double::TryParse(textBox2->Text, b);
double K = 170/(a+2*b);
Brush^ brush = gcnew SolidBrush(Color::White);
g->FillRectangle(brush, 0, 0, 360, 360);
Pen^ pen2 = gcnew Pen(Color::Red, 2);
g->DrawEllipse(pen2, 180-(int)(K*a), 180-(int)(K*a),
(int)(2*K*a), (int)(2*K*a));
Pen^ pen = gcnew Pen(Color::Blue, 2);
double pi = Math::PI;
int ox = (int)(180+K*((a+b)*cos(0.0)-b*cos(0.0)));
int oy = (int)(180-K*((a+b)*sin(0.0)-b*sin(0.0)));
for (double t=0; t<=2*pi*b; t += pi/80) {
int nx = (int)(180+K*((a+b)*cos(t)-b*cos((a+b)*t/b)));
int ny = (int)(180-K*((a+b)*sin(t)-b*sin((a+b)*t/b)));
g->DrawLine(pen, ox, oy, nx, ny);
ox = nx;
oy = ny;
}
}
51
コンパイルし、実行すると
となる。button1 をクリックすると
となる。
レポート問題:原点を中心とし半径が a の円 O の内側を半径が b a > b > 0) の円 C が円 O
に内接しながら滑ることなく転がるとき、円 C 上の点 P の軌跡を考える。ただし、点 P のはじ
めの位置は、円 O と x 軸の正の部分との交点 A とする。
円 C が転がるとき、動径 OC が表す角を θ とし、そのときの外接する点を Q、点 P の座標
を (x,y) とすると、P の軌跡は θ を媒介変数として
x = (a − b) cos θ + b cos
a−b
a−b
θ, y = (a − b) sin θ − b sin
θ
b
b
で表される。この曲線を ハイポサイクロイド という。ハイポサイクロイドを描くプログラムを
作れ。
52
レポート問題:極方程式
r = a + b cos θ
で表される曲線を リマソン という。特に、a = b のとき、極方程式
r = a(1 + cos θ)
で表される曲線を カージオイド という。リマソンを表示するように、上のプログラムを変更せよ。
レポート問題:
x = (a + b) cos θ − c cos
a+b
a+b
θ, y = (a + b) sin θ − c sin
θ
b
b
で表される曲線を描くプログラムを作れ。
レポート問題:
x = (a − b) cos θ + c cos
a−b
a−b
θ, y = (a − b) sin θ − c sin
θ
b
b
で表される曲線を描くプログラムを作れ。
53
4
コンソール・プログラミング
自分のためにちょこちょこと計算するには、コンソールアプリケーションを使います。
プログラムの作成と実行
手順は以下の通り
1. 新しいプロジェクトを作成する。
2. C++ ソースファイルをプロジェクトに追加する。
3. ソースコードを入力する。
4. 実行可能ファイルをビルドする。
5. プログラムを実行する。
6. プログラムを保存する。
新しいプロジェクトの作成
1. Visual Studio アイコンをクリックして Visual C++ IDE を開く。
2. 「ファイル」メニューを開き、「新規作成」→ 「プロジェクト」をクリックする。
3. 「プロジェクトの種類」で「Visual C++」を選択する。
4. 「テンプレート」または「インストールされたテンプレート」セクションで「Win32 コン
ソールアプリケーション」を選択する。
5. 「プロジェクト名」または「名前」テキストボックスにプロジェクトの名前(例えば hello)
を入力する。
6. プロジェクトのディレクトリを選択する。通常はデフォルトでよい。
7. 「OK」をクリックする。
8. 「Win32 アプリケーションウィザード」が起動する。
9. ページの左側で「アプリケーションの設定」をクリックするか、または「次へ」をクリック
する。
10. 「追加のオプション」から「空のプロジェクト」を選択する。
11. 「完了」をクリックする。これで、すべてのコンパイラ設定がコンソールプロジェクトのた
めに初期化されるはずです。
プロジェクトへの C++ ソースファイルの追加
プログラムには、ソースファイルが少なくとも1つ(たいていは複数)必要です。
1. メニューバーの「新しい項目の追加」アイコン(通常は左から二つ目のアイコン)をクリッ
クする。「新しい項目の追加」ダイアログボックスが表示される。
2. 「Visual C++」カテゴリから「コード」を選択する。
54
3. テンプレートセクションで「C++ ファイル (.cpp)」アイコンを選択する。
4. 「ファイル名」テキストボックスにプログラムファイルの名前(hello)を入力し、
「追加」を
クリックする。
これで、空のソースコードファイルが作成され、ソースコードプログラムを入力する準備が整った。
ソースコードの入力
この時点で、ソースコードを IDE に直接入力するか、別のソースからコピーして貼り付けるか
を選択できる。
実行可能プログラムのビルド
プログラムのソースコードを正しく入力できたという自信がある場合は、「ビルド」メニューか
ら「ソリューションのビルド」を選択するか、IDE ウィンドウの上部にあるアイコンリストで右
向き三角形アイコンをクリックする。そうすると、IDE がプログラムのコンパイルとリンクを試
みる。それらが正常終了した場合は、「出力」ウィンドウに次のメッセージが表示される。
ビルド:1 正常終了、0 失敗、0 更新不要、0 スキップ
正常終了しなかった場合は、エラーメッセージが表示される。プログラムをデバッグしてエラー
を修正し、再び「ソリューションのビルド」を実行する。
右向き三角形アイコンを使用した場合、エラーがなければプログラムが実行を自動的に開始す
る。「ソリューションのビルド」メニュー項目を使用した場合は、プログラムを自分で開始する必
要がある。
プログラムの実行
すべてのエラーを取り除いたら、「デバッグ」メニューで「デバッグなしで開始」を選択し、プ
ログラムを実行する。
プログラムの保存
「ファイル」メニューの「すべてを保存」をクリックする。保存するのを忘れて IDE を閉じよ
うとすると、IDE がそれを知らせてくれる。
ごちゃごちゃ書いていますが無視してともかくやってみます。
まず、Microsoft Visual Studio 2010 を立ち上げる。
このようにならなければ、「表示」のメニューで「スタートページ」をクリックします。
新しいプロジェクト・
・
・をクリックする。
55
Win32 コンソール アプリケーションをクリックし、図では場所 (L):が D:\jyugyou\となって
います(自分で適当なディレクトリ、例えば D:\jyugyou\を作り、ここを指定すれば、プログラム
が D:\jyugyou\に作られます)が、何処にプログラムが保存されるか気にしなければ、そのまま
で良いですので、名前 (N): のエディットボックスに例えば hello と入力し、
OK ボタンをクリックする。
次へのボタンをクリックし、空のプロジェクトをチェックする。
56
完了ボタンをクリックする。
クラスビューになってなければ、クラスビューをクリックする。
一段目の左から二つ目のアイコン(新しい項目の追加)をクリックする。
57
C++ ファイルをクリックし、名前の欄に hello と打ち込む。
追加のボタンをクリックする。
hello.cpp のファイルが作られ、開いている。ここにプログラムを打ち込む。
58
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char *argv[])
{
cout << "Hello, world" << endl;
getchar();
return 0;
}
コンパイルし、実行する(Debug の横の赤い三角をクリックする)と
と Hello,world と表示する。
以下では、このプログラムを雛形として
cout << "Hello, world" << endl;
の部分を作りかえて、プログラミングしていきます。
59
4.1
式と計算
いろいろな計算を行うには、変数と呼ばれるデータの記憶場所に数値を保存させる。変数は名前
(変数名)を付けて区別する。変数名はアルファベット、数字、_(アンダーバー)からなる文字
列で、先頭の文字には数字は使えない。大文字と小文字は区別するので abc ABC Abc AbC はす
べて異なる変数を表す。次の識別子(名前)はキーワード(予約語)として使われるので、別の用
途に使えない。予約語は文字がゴチックになるので、すぐわかります。
auto
double int
struct
break
else
long
case
const
default
enum
float
goto
register
short
sizeof
typedef
unsigned
volatile
char
continue
do
extern
for
if
return
signed
static
switch
union
void
while
asm
_cs
_ds
_es
_ss
cdecl
far
huge
interrupt near
pascal
C++ では更に次の識別子もキーワードです。
class inline overload public friend new virtual operator
this
delete
次は三つの数を入力し、その和と平均を出力する C のプログラムです。
Microsoft Visual Studio 2010 で、 C のプログラムを作るには次のようにします。基本的には
C++ のプログラミングと同じですが、Microsoft Visual Studio 2010 を起ち上げます。
のような画面になることもあります。上でやったように、
「表示」のメニューで「スタートページ」
をクリックし、新しいプロジェクト・
・
・をクリックしても良いですが、次のようにすることも出来
ます。「ファイル(F)」メニューの「新規作成(N)」の「プロジェクト(P)」をクリックします。
60
Win32 コンソール アプリケーションをクリックし、名前 (N): のエディットボックスに例えば
heikin3 と入力し、
OK ボタンをクリックする。
次へのボタンをクリックし、空のプロジェクトをチェックする。
61
完了ボタンをクリックする。
クラスビューになってなければ、クラスビューをクリックし、一段目の左から二つ目のアイコ
ン(新しい項目の追加)をクリックし、C++ ファイルをクリックし、名前の欄に heikin.c と打ち
込む。
追加のボタンをクリックする。
62
heikin.c のファイルが作られ、開いている。
#include <stdio.h>
int main()
{
double a, b, c, wa, heikin;
printf("a b c = ");
scanf("%lf %lf %lf", &a, &b, &c);
wa = a + b + c;
heikin = wa / 3;
printf("和 = %lf\n", wa);
printf("平均 = %lf\n", heikin);
getchar();
return 0;
}
と入力する。
C のプログラムも C++ のプログラムも main() 関数から実行する。
int main(int argc, char *argv[])
{
}
になっていたり、
int main()
{
}
となっているが、あまり気にしなくて良い。上の場合は、実行プログラムに引数を与えて、引数に
応じて、プログラムの実行を変化したいときに使うが、上のようにして、引数を与えなくても良い。
#include <stdio.h>
63
は、
printf("a b c = ");
scanf("%lf %lf %lf", &a, &b, &c);
printf("和 = %lf\n", wa);
printf("平均 = %lf\n", heikin);
getchar();
で使われているライブラリ関数 printf(), scanf(), getchar() を使うために必要なインクルードファイ
ルです。ライブラリ関数とは、メーカーが作って、ライブラリに納められている関数です。printf(),
scanf(), getchar() のプロトタイプ(どのような引数を取り、どのような値を戻すか)が stdio.h に
書かれています。これを include 命令で取り込んでいます。関数は使う前にプロトタイプ(どのよ
うな引数を取り、どのような値を戻すか)を書いておかなくてはなりません。
double a, b, c, wa, heikin;
は、使用する変数の宣言です。プログラムで使用する変数は記憶するための場所を確保するために、
前もって宣言する必要がある。C では変数の宣言は関数定義の先頭で行う。double は基本データ
型の一つ倍精度浮動小数点数型を示す。
double a, b, c, wa, heikin;
は変数 a, b, c, wa, heikin は倍精度浮動小数点数のための変数であると宣言している。double は
4Byte で +- 2.225073858507201e-308 から+-1.797693134862316e+308 の範囲の実数を記憶する
ことができます。よく使う他の基本データ型には char と int と float がある。
printf("a b c = ");
の printf(”a b c = ”) は、コンソール(画面)に a b c = と表示します。printf() 関数は ” と ” で
囲われた文字列を一部の例外を除いてそのまま表示します。
scanf("%lf %lf %lf", &a, &b, &c);
の scanf() は標準ライブラリ関数で書式付きでキーボードからデータを入力します。scanf() 関数
の最初の引数が書式用の文字列です。%lf がキーボードから入力された数値(文字列)を倍精度浮
動小数点数に変換することを示します。したがって、この場合はスペースで区切られた数値を倍精
度浮動小数点数として三個入力することを指示しています。scanf() 関数はキーボードから文字を
読み込み、書式で指定された形式に従ってそれを解釈し、残りの引数へ結果を格納する。二番目以
降の引数は、それに対応する変換入力を格納できるようにポインタでなければならない。ポインタ
にするためには変数名の前に & を付ければよい。今の場合、変数 a, b, c にこの順で倍精度浮動小
数点数が入力される。ポインタは説明するのがやっかいなので、気になって仕方がなければ、自分
でインターネットで調べるか、C または C++ の入門書を読んで下さい。ここでは、scanf() では
変数の頭に & を付ければよいと覚えておけばいいです。C++ では、以下で見るようにもっと簡単
です。
wa = a + b + c;
64
は代入文です。代入文は 変数=式; の形式で式の値を計算し、その値を変数に代入します。もと
在った変数の値はなくなり新しい値で置き換わります。今の場合、変数 a, b, c の値の和を取りそ
れを変数 wa の値とします。
heikin = wa / 3;
は変数 wa の値を3で割り変数 heikin に代入します。+ と / は算術演算子です。よく使う算術演
算子には +, -, *, /, % があります。+ は加法を、 - は単項のマイナス演算子と減法を、 * は
乗法を、/ は除法(a / b は a を b で割った値になります。ただし、整数の除法は小数点以下切り
捨てて整数値にします)を、また % は整数の余り (a % b は a を b で割った余りですが、a が負の
数のときには、負の余りが計算されます) を求めるのに使われる。+ - を加法演算子* / % を乗法
演算子といい、乗法演算子の方が加法演算子より先に計算される。乗法演算子、加法演算子同士の
場合は左から右に順に計算される。
最後に
printf("和 = %lf\n", wa);
printf("平均 = %lf\n", heikin);
の printf() 関数で wa と heikin の値が出力される。printf() も最初の引数が書式で %lf の所に wa
および heikin の値が代入されて出力される。 wa および heikin が倍精度浮動小数点数型なので
%lf が使われます。 %\n は改行を指示します。
getchar();
は、これがなければ、Microsoft Visual Studio 2010 の統合開発環境では、直ちにコンソール画面
が閉じますから、これを実行して、何か文字が入力されるまで待ちます。例えば、エンターキーを
押すとプログラムが終了します。
return 0;
は
int main()
と main() 関数が int (整数値)を戻す関数と定義されているので、正常に終了した事を示す値 0
を返しています。
C++ では、次のようにプログラミングする事もできます。
#include <iostream.h>
using namespace std;
int main()
{
double a, b, c, wa, heikin;
cout << "a b c = ";
cin >> a >> b >>c;
wa = a + b + c;
heikin = wa / 3;
65
cout << "和 = " << wa << "\n";
cout << "平均 = " << heikin << endl;
getchar(); getchar();
return 0;
}
このプラグムの
#include <iostream.h>
は
cout << "a b c = ";
cin >> a >> b >>c;
cout << "和 = " << wa << "\n";
cout << "平均 = " << heikin << endl;
getchar(); getchar();
の cout, cin, getchar() のために必要なインクルードファイルです。
using namespace std;
は、 cout, cin などが std という namespace で宣言、定義されているので、std という namespace
を使いますという宣言です。即ち、cin と cout は名前が重複しないように、 namespace の std と
いうところで宣言されているので、
using namespace std;
とするか、std::cin と std::cout の形で使う必要があります。std::cin と std::cout とするのはめん
どくさいので、
using namespace std;
としています。
C++ ではいつでも
#include <iostream.h>
using namespace std;
または
#include <iostream>
using namespace std;
を最初に書くものだと思って下さい。
cout << "a b c = ";
で、コンソールに ”a b c = ” と表示します。
cin >> a >> b >>c;
66
で、キーボードから入力されたスペースで区切られた 3 つの数字を a, b, c にセットします。a, b,
c 1つづつに >> を使います。入力の時、アドレスを渡す(ポインタを使う)必要がありません。
書式の指定も必要ありません。
cout << "和 = " << wa << "\n";
cout << "平均 = " << heikin << endl;
で、 文字列”和 = ” と wa の値と "\n"(改行)と 文字列”平均 = ” と heikin の値と endl (改
行)を表示します。改行はこのように endl でも可能です。このように数値の桁数の指定が必要に
ないときは、C++ の方が楽です。
次は三桁の整数を入力して百、十、一の位の数を表示する C のプログラムです。
#include <stdio.h>
void main()
{
int n, hyaku, juu, iti;
printf("三桁の数 ");
scanf("%d", &n);
hyaku = n / 100;
n %= 100;
juu = n / 10;
iti = n % 10;
printf("百の位 = %d\n", hyaku);
printf("十の位 = %d\n", juu);
printf("一の位 = %d\n", iti);
}
int n, hyaku, juu, iti;
は変数 n, hyaku, juu, iti を整数のための変数であると宣言している。int は 2Byte で -32768 か
ら 32767 までの整数値を記憶する事ができます。これの値はコンパイラに依ります。4Byte で
-2147483648 から 2147483647 までの整数値であることもあります。もっと大きい整数値を記憶し
たければ long int を使います。long int は 4Byte で -2147483648 から 2147483647 までの整数値
を記憶できます。これより大きな整数値を使いたければ、後で出てくる配列を使います。
scanf() および printf() の書式中の %d は整数値に変換することを指示している。
hyaku = n / 100;
は n および 100 が整数値であるから、n の値を 100 で割った商(小数点以下切り捨て)を変数
hyaku に代入することを指示している。
n %= 100;
は n = n % 100; の省略形で n の値を 100 で割った余りを n の新しい値とするよう指示してい
る。このような省略形は任意の算術演算子に適用できる。たとえば、
67
a += 10; b -= c;
はそれぞれ
n *= a+b;
a = a + 10; b = b - c;
と同じである。
n /= s;
n = n * (a+b);
n = n / s;
C++ のプログラムは次のようになります。
#include <iostream>
using namespace std;
void main()
{
int n, hyaku, juu, iti;
cout << "三桁の数 ";
cin >> n;
hyaku = n / 100;
n %= 100;
juu = n / 10;
iti = n % 10;
cout << "百の位 = " << hyaku << endl;
cout << "十の位 = " << juu << endl;
cout << "一の位 = " << iti << endl;
}
次のプログラムは atan() 関数を用いて円周率をもとめる C のプログラムである。
#include <stdio.h>
#include <math.h>
void main()
{
printf("円周率=%20.15f\n", 4 * atan(1));
getchar();
}
ここでは
void main()
{
}
と main() 関数の戻り値を void (値を戻さない)と宣言しているが、このようにしても良いです。
この場合は return 文は必要ありません。
atan() は倍精度浮動小数点数の引数を取り、その arctangent を倍精度浮動小数点数として返す
関数である。数学関数を用いるためには
#include <math.h>
68
のインクルード文が必要である。以下にあげるのはよく使われる数学関数である。いずれも一つま
たは二つの double の引数をとり、double の値を返す。
sin(x)
x の sine, 単位はラジアン
cos(x)
tan(x)
x の cosine, 単位はラジアン
x の tangent, 単位はラジアン
exp(x)
log(x)
log10(x)
指数関数 ex
x の自然対数(基数は e)(x > 0)
x の常用対数(基数は10)(x > 0)
pow(x, y)
sqrt(x)
xy
x の平方根(x ≥ 0)
fabs(x)
x の絶対値
floor(x)
x を越えない最大の整数
ceil(x)
x 以上の最小の整数
printf() 関数の書式の %20.15f は double の値を20桁取り小数点以下15桁表示せよとの指示
である。
C++ のプログラムは次のようになります。
#include <iostream>
#include <math.h>
using namespace std;
void main()
{
cout << "円周率=" << 4 * atan(1) << endl;
getchar();
}
C++ で、値を20桁取り小数点以下15桁表示せよと言った指示はめんどくさいのでここでは示
さない。気になるなら、自分でインターネットで調べなさい。私のホームページにも載せてあり
ます。
例題:球の半径 r の値を入力し、表面積と体積を求める C++ のプログラムを作れ。
解答例
次のようなプログラムを作ればいいです。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main(int argc, char* argv[])
{
double PI = 4 * atan(1.0);
double r, V, S;
69
cout << "半径:";
cin >> r;
S = 4 * PI * r * r;
V = 4 * PI * r * r* r / 3;
cout << "S=" << S << " V=" << V << endl;
getchar();
getchar();
return 0;
}
例題:標準体重は 身長 (m)× 身長 (m) ×22 の ±10% 以内です。身長を入力し、その標準体重
を出力するプログラムを作れ。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
double A, W, L, H;
cout << "身長 (m)=";
cin >> A;
getchar();
W = A * A * 22;
L = W * 0.9;
H = W * 1.1;
cout << "標準体重(kg)は" << L << "から" << H << "の範囲です" << endl;
getchar();
return 0;
}
例題: 24 を、 2 ∗ 2 ∗ 2 ∗ 2 ではなく、底が e の対数関数 log() と底が e の指数関数 exp()
を利用して、計算するプログラムを作れ。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main(int argc, char* argv[])
{
double a, b;
70
cout << "底=";
cin >> a;
cout << "指数=";
cin >> b;
getchar();
cout << a << "の" << b << "乗は" << exp(b*log(a)) << "です" << endl;
getchar();
return 0;
}
例題: 余弦定理 c2 = a2 + b2 − 2ab cos C を利用して、a,b,c の値から三角形の角 C を求める
プログラムを作れ。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main(int argc, char* argv[])
{
double a, b, c, C;
double PI = 4*atan(1.0);
cout << "a="; cin >> a;
cout << "b="; cin >> b;
cout << "c="; cin >> c;
getchar();
C = acos((a*a+b*b-c*c)/(2*a*b))/PI*180;
cout << "C=" << C << endl;
getchar();
return 0;
}
演習問題
1. 自然数 n を入力し n を 3 で割ったときの商と余りを表示するプログラムを作れ。
2. 4桁の数を入力して千、百、十、一の位の数を出力するプログラムを書け。
3. 対数関数 LOG() または log10(x) を利用して入力した数の桁数を出力するプログラムを書け。
4. 余弦定理 c2 = a2 + b2 − 2ab cos C を利用して、a,b,C の値から三角形の辺 AB の長さ c を
求めるプログラムを作れ。
5. 三角形 ABC の辺 BC の長さ a と角 A から、外接円の半径 R を求めるプログラムを作れ
71
旺文社の昔の数学 A の教科書には、次のような最大公約数を求めるプログラムが載っています。
10 REM 最大公約数
20 INPUT "N, M"; N, M
30 G=N
40 IF N>INT(N/G)*G THEN 80
50 IF M>INT(M/G)*G THEN 80
60 PRINT "最大公約数=";G
70 GOTO 100
80 G=G-1
90 GOTO 40
100 END
これを C++ に翻訳すると次のようになります。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int n, m, g;
cout << "n m:";
cin >> n >> m;
getchar();
g = n;
L0: if (n % g != 0) goto L1;
if (m % g != 0) goto L1;
cout << "最大公約数=" << g << endl;
goto L2;
L1: g = g-1;
goto L0;
L2: getchar();
return 0;
}
これは昔のプログラミングスタイルで、読みにくいプログラムになります。
goto L2;
は goto 文と言われ、この場合は L2 というラベルの付いた文に制御を移します。つまり、
L2: getchar();
72
を実行します。ラベルは文字列にコロン : を付けて表します。BASIC の場合は行番号で分岐先を
指示しています。
40 IF N>INT(N/G)*G THEN 80
は、N が G で割り切れないなら80行にジャンプしなさいという IF 文です。
このプログラムのように goto 文を多用すると読みにくいプログラムになるので、構造的プログ
ラミングと言って、 PASCAL や C と言ったプログラミング言語が作られました。n と m をス
ペースで区切って、続けて入力します。
これらプログラムのアルゴリズムは
1. N を入力する
2. M を入力する
3. G を N から N-1, N-2, N-3, · · ·, 1 と、G が N と M も共に割り切るまで、一つずつ減らし
ていく。G が N と M も共に割り切るとき、G は N と M の最大公約数である。
4. G を N と M の最大公約数として表示する。
です。
次のようにプログラミングする方が良いです。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int n, m, g;
cout << "n m:";
cin >> n >> m;
getchar();
for (g=n; g>0; g--) {
if (n % g == 0 && m % g == 0)
break;
}
cout << "最大公約数=" << g << endl;
getchar();
return 0;
}
for 文を使っています。g-- は g = g - 1 と同じです。
for 文は
73
for (初期設定; 条件式; 後処理) {
命令の並び
}
の形式で、まず、初期設定を実行し、条件式が真であれば、命令の並びを実行し、後処理をし、条
件式が真であれば、命令の並びを実行し、後処理をし、を条件式が偽になるまで繰り返します。条
件式の一回目の判定が偽であれば、命令の並びは一回も実行されません。
for (初期設定; 条件式; 後処理) {
命令の並び
}
で、命令の並びが1つの命令の場合は
for (int i=1; i<= n; i++)
r *= i;
や
for (int i=0; i<=n; i++)
if (i % 2 == 0)
s += pow(x, 2*i) / fact(2*i);
else
s -= pow(x, 2*i) / fact(2*i);
のように、文が1つの場合は、命令の並びを囲む { と } は省略できます。
for (g=n; g>0; g--) {
if (n % g == 0 && m % g == 0)
break;
}
は、g を n とし 0 より大きい間は g を1つづつ減らしながら
if (n % g == 0 && m % g == 0)
break;
もし n が g で割り切れ(n % g == 0)且つ(&&)m が g で割り切れる(m % g == 0)なら、for
文を抜ける( break 文)ことを指示しています。&& は「かつ」を表す論理演算子です。
if 文は、 Scheme の場合と異なり、2つの形式があります。
if (条件式) {
命令の並び
}
の形で、条件式が真 (true) の時、命令の並びを実行するか
if (条件式) {
命令の並び
} else {
命令の並び
}
74
の形で、条件式が真 (true) の時、上の命令の並びを実行し、そうでないとき下の命令の並びを実
行します。
if (条件式) {
命令の並び
}
や
if (条件式) {
命令の並び
} else {
命令の並び
}
で、命令の並びが一個の命令の場合は、
if (i % 2 == 0)
s += pow(x, 2*i) / fact(2*i);
else
s -= pow(x, 2*i) / fact(2*i);
のように、命令の並びを囲む { と } は省略できます。
これは次のようにプログラミングすることもできます。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int n, m;
cout << "n m:";
cin >> n >> m;
getchar();
while (m != 0) {
int temp = n;
n = m;
m = temp % m;
}
cout << "最大公約数=" << n << endl;
getchar();
return 0;
}
75
while 文を使っています。このアルゴリズムを「ユークリッドの互除法」といいます。
while 文は
while (条件式) {
命令の並び
}
の形式で、条件式が真 (true) の間は{ と } の間の命令達を繰り返します。今の場合、m が 0 でな
い間は { と } の間の命令達を繰り返します。
次の例は素因数分解の BASIC のプログラムです。これも旺文社の昔の数学 A の教科書で見つ
けたプログラムです。
10 REM 素因数分解
20 INPUT N
30 M=1
40 M=M+1
50 IF INT(N/M)*M<N THEN 40
60 PRINT M;
70 N=N/M
80 IF N>1 THEN 50
90 END
これも次のようにプログラミングするのが良いです。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int n, k;
cout << "n=";
cin >> n;
getchar();
k = 2;
while (n > 1) {
while (n % k == 0) {
cout << k << " ";
n = n / k;
}
k++;
}
cout << endl;
76
getchar();
return 0;
}
while 文が入れ子になっています。
例題:ヘロンの公式を使って、三角形の面積を計算するプログラムを作れ。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main(int argc, char* argv[])
{
double a, b, c;
cout << "a="; cin >> a;
cout << "b="; cin >> b;
cout << "c="; cin >> c;
getchar();
if (a+b<=c || b+c<=a || c+a<=b) {
cout << "三角形ではありません" << endl;
getchar();
return 1;
}
double S = (a+b+c)/2;
S = sqrt(S*(S-a)*(S-b)*(S-c));
cout << "面積は" << S << "です" << endl;
getchar();
return 0;
}
if (a+b<=c || b+c<=a || c+a<=b) {
cout << "三角形ではありません" << endl;
getchar();
return 1;
}
の if 文は、a + b が c 以下(a+b<=c) か b + c が a 以下(b+c<=a) か c + a が b 以下
(c+a<=b) の時という意味です。
例題:Zeller の公式を使って、西暦の年月日を入力して曜日を計算するプログラムを作れ。
#include <iostream>
#include <stdio.h>
77
using namespace std;
int main(int argc, char* argv[])
{
int y, m, d, w;
cout << "西暦年:"; cin >> y;
cout << "月:"; cin >> m;
cout << "日:"; cin >> d;
getchar();
m = m - 2;
if (m < 1) {
y--; m += 12;
}
int b = y % 100;
int c = y / 100;
w = ((26*m-2)/10+d+b+b/4+5*c+c/4) % 7;
switch(w) {
case 0: cout << "日曜日です" << endl; break;
case 1: cout << "月曜日です" << endl; break;
case 2: cout << "火曜日です" << endl; break;
case 3: cout << "水曜日です" << endl; break;
case 4: cout << "木曜日です" << endl; break;
case 5: cout << "金曜日です" << endl; break;
case 6: cout << "土曜日です" << endl; break;
}
getchar();
return 0;
}
switch(w) {
case 0: cout << "日曜日です" << endl; break;
case 1: cout << "月曜日です" << endl; break;
case 2: cout << "火曜日です" << endl; break;
case 3: cout << "水曜日です" << endl; break;
case 4: cout << "木曜日です" << endl; break;
case 5: cout << "金曜日です" << endl; break;
case 6: cout << "土曜日です" << endl; break;
}
は、w が 0 なら日曜日ですと表示し、w が 1 なら月曜日ですという様に w の値で場合分けします。
演習問題
78
1. ある都市の中型タクシーの料金は、2 km までは 500 円、あとは 300 m までごとに 80 円増
すことになっている。乗車距離 x (km) を入力し、料金 y (円) を表示するプログラムを作れ。
2. 西暦年を入力して閏年かどうかを判定するプログラムを作れ。1583年以降グレゴリオ暦
が使われていて、閏年は、4で割り切れるが100の倍数で無いときである。ただし400
で割り切れれば閏年である。
3. a,b,c が定数、a ̸= 0 のとき a,b,c を入力して、二次方程式 ax2 + bx + c = 0 の解の個数を求
めるプログラムを作れ。
4. 三角形 ABC の3辺の長さ a, b, c を与えこれから三角形 ABC が鋭角三角形、直角三角形、
鈍角三角形のいずれであるかを判定するプログラムを作れ。
例題:整数を入力し、素数かどうかを判定するプログラムを作れ。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main(int argc, char* argv[])
{
int n;
bool flag;
cout << "n="; cin >> n; getchar();
if (n % 2 == 0) {
cout << "2で割り切れます" << endl;
getchar();
return 0;
}
flag = true;
int i;
for (i=3; i<=sqrt((double)n); i += 2) {
if (n % i == 0) {
flag = false;
break;
}
}
if (flag == true)
cout << "素数です" << endl;
else
cout << i << "で割り切れます" << endl;
getchar();
79
return 0;
}
これは次のように関数を使ってもいいです。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
bool primeP(int n)
{
for (int i=2; i<=sqrt((double)n); i++) {
if (n % i == 0)
return false;
}
return true;
}
int main(int argc, char* argv[])
{
int n;
cout << "n="; cin >> n; getchar();
if (primeP(n))
cout << "素数です" << endl;
else
cout << "合成数です" << endl;
getchar();
return 0;
}
bool は true か false の値をとる真理値の型です。
関数は
bool primeP(int n)
{
for (int i=2; i<=sqrt((double)n); i++) {
if (n % i == 0)
return false;
}
return true;
}
のように、
80
戻り値の型 関数名 (引数の並び)
{
関数の本体
}
からなります。上の関数は整数値である1つの引数 n を取り、真理値 true か false を戻す関数
です。関数は使う前に宣言しておく必要があります。今の場合、primeP() が使われるのは main()
関数で、 main() 関数の前に primeP() が定義されているので、これで良いですが、もし
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
bool primeP(int n);
int main(int argc, char* argv[])
{
int n;
cout << "n="; cin >> n; getchar();
if (primeP(n))
cout << "素数です" << endl;
else
cout << "合成数です" << endl;
getchar();
return 0;
}
bool primeP(int n)
{
for (int i=2; i<=sqrt((double)n); i++) {
if (n % i == 0)
return false;
}
return true;
}
のように、 main() 関数の後ろに primeP() が定義されている場合は、main() 関数の前に
bool primeP(int n);
のように primeP() がどのような引数を取り、どのような戻り値を返すか宣言しておく必要があり
√
ます。bool primeP(int n) 関数は、整数値 n が i=2 から n まで i を1つづつ増やしながら、ど
81
れかの i が n を割り切るなら false を返し、すべて割り切らなければ true を返す関数として定義
しています。
例題: 解析学で学ぶように、cos(x) や sin(x) は
cos(x) =
∞
∑
(−1)n 2n
x ,
(2n)!
n=0
sin(x) =
∞
∑
(−1)n 2n+1
x
(2n + 1)!
n=0
と級数を使って表すことができます。x はラジアンです。この級数はすべての x に対して定義さ
れていて非常に計算効率が良いです。これを使って、cos(x) を計算するプログラムを作れ。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main(int argc, char* argv[])
{
double x, y, new_c, old_c;
double t;
int n;
cout << "X="; cin >> x; getchar();
double PI = 4 * atan(1.0);
y = x * PI / 180;
n = 0;
t = 1;
new_c = 1;
do {
cout << n << "回目:" << new_c << endl;
old_c = new_c;
n++;
t = t / ((2*n - 1)*2*n)*y*y;
if (n % 2 == 0)
new_c += t;
else
new_c -= t;
} while (fabs(new_c - old_c)>0.0000001);
cout << n << "回目: cos(" << x << ")=" << new_c << endl;
getchar();
return 0;
}
do while 文を使っています。do while 文は
82
do {
命令の並び
} while (条件式) ;
の形式で、まず、命令の並びを実行した後、条件式が真 (true) の間は{ と } の間の命令達を繰り
返します。
このプログラムはもっと分かり易いプログラムを私のホームページに載せています。
演習問題
1. 整数を入力しその約数をすべて表示するプログラムを作れ。
2. cos(x) を計算するプログラムを真似て、 sin(x) を計算するプログラムを作れ。
3. ライプニッツの級数を利用して、円周率を計算するプログラムを作れ。
π
1 1 1 1
= 1 − + − + − ···
4
3 5 7 9
がライプニッツの級数です。本当はグレゴリーの級数という方が正しいそうで、歴史的には
余りにも有名ですが、円周率の計算には向かない。
4. 自然数 n について、入力された n の各桁の数字の和を求めるプログラムを作れ。
5. a2 + b2 = c2 となる自然数 a, b, c の組みとその個数を求めるプログラムを作れ。但し c < 100
とする。
6. 入力された自然数 a, b の公約数の総和を求めるプログラムを作れ。
4.2
配列
いままで扱った変数には一つの数値しか記憶できない。多数の数値をまとめて効率よく扱うには
配列を使う。配列の宣言は変数名の後ろに [ と ] で配列の大きさを囲む。
int a[10];
は整数型のデータ記憶場所を 10 個準備し、a[0],a[1],a[2], ... ,a[9] でアクセスできる。a[0]
から a[9] であることに注意。次はフィボナッチ数列の値を求める配列を用いたプログラムである。
#include "stdafx.h"
#include <iostream>
#include <stdio.h>
using namespace std;
int main()
{
int f[21], i;
f[0] = 1;
83
f[1] = 1;
for (i=2; i <= 20; i++)
f[i] = f[i-1] + f[i-2];
for (i=0; i <= 20; i++)
cout << i << " " <<
getchar();
return 0;
f[i] << endl;
}
次はベクトルの外積を計算するプログラムである。
#include <iostream>
#include <stdio.h>
using namespace std;
void vp(double a[], double b[], double c[])
{
c[0] = a[1] * b[2] - a[2] * b[1];
c[1] = a[2] * b[0] - a[0] * b[2];
c[2] = a[0] * b[1] - a[1] * b[0];
}
int main()
{
double a[3], b[3], c[3];
cout << "ベクトル1の x, y, z?";
cin >> a[0] >> a[1] >>a[2];
cout << "ベクトル2の x, y, z?";
cin >>b[0] >> b[1] >>b[2];
vp(a, b, c);
cout << "外積=" << c[0] << " " << c[1] << " " << c[2] << endl;
getchar();
getchar();
return 0;
}
配列は
void vp(double a[], double b[], double c[])
{
c[0] = a[1] * b[2] - a[2] * b[1];
c[1] = a[2] * b[0] - a[0] * b[2];
c[2] = a[0] * b[1] - a[1] * b[0];
}
84
のように関数に渡すこともできる。このときアドレス渡しとなる。すなわち、関数内で値を変更す
ると実引数も値がその値になる。もう少し説明すると
#include <iostream>
#include <stdio.h>
using namespace std;
void func(int a)
{
a = 6;
}
int main()
{
int n = 10;
func(n);
cout << n << endl;
getchar();
}
を実行すると 6 ではなく 10 と表示する。これは
void func(int a)
{
a = 6;
}
が値渡しになっていて、func() の中で a の値を変えても、元の n の値は変わりません。これを
#include <iostream>
#include <stdio.h>
using namespace std;
void func(int *a)
{
*a = 6;
}
int main()
{
int n = 10;
func(&n);
cout << n << endl;
85
getchar();
}
とポインタを使って、アドレス渡しにすると 6 と表示します。
例題:素数を小さい方から100個見つけるプログラムを作れ。
旺文社数学Aには次のプログラムがある。
10 REM 素数の発見
20 INPUT "求める素数の個数";X
30 DIM P(X)
40 J=1:P(J)=2
50 PRINT J,P(J)
60 N=3
70 FOR M=1 TO J
80
IF INT(N/P(M))*P(M)=N THEN 130
90 NEXT M
100 J=J+1
110 P(J)=N
120 PRINT J,P(J)
130 N=N+1
140 IF J<X THEN 70
150 END
30 行目の DIM P(X) が配列の宣言で P(0) から P(X) までの X+1 個の実数型の配列を宣言して
いる。C に翻訳すると次のようになる。
#include <stdio.h>
#include <stdlib.h>
int main()
{
int x;
int *p;
int j, n, m;
printf("求める素数の個数"); scanf("%d", &x);
if ((p = (int *)malloc((x+1)*sizeof(int))) == NULL) {
printf("malloc error!\n");
exit(1);
}
j = 1;
p[j] = 2;
printf("(%d %d) ", j, p[j]);
n = 3;
L1: for (m=1; m<=j; m++)
if (n % p[m] == 0)
86
goto L2;
p[++j] = n;
printf("(%d %d) ", j, p[j]);
L2: n++;
if (j < x)
goto L1;
free(p);
getchar();
return 0;
}
BASIC と異なり C では実行時に配列の大きさを決めることができないのでポインタを宣言し
(今の場合 int *p)、malloc() 関数を用いて必要な領域を割り当てる。
if ((p = (int *)malloc((x+1)*sizeof(int))) == NULL) {
printf("malloc error!\n");
exit(1);
}
その領域が不要になれば
free(p);
で、領域を開放する。そうすれば配列と同様な使い方ができる。
p[++j] = n は ++j; p[j] = n; と同じでまずjを一つ増やしそれをjの値として使いnの値
を代入する。。もしこれが p[j++] = n; であれば、元のjの値を使い代入をすました後、jを一
つ増やす。このプログラムは次のようにすると分かりやすい。
#include <stdio.h>
#include <stdlib.h>
int main()
{
int x;
int *p;
int j, n, m;
int flag;
printf("求める素数の個数"); scanf("%d", &x);
if ((p = (int *)malloc((x+1)*sizeof(int))) == NULL) {
printf("malloc error!\n");
exit(1);
}
j = 1;
p[j] = 2;
printf("(%d %d) ", j, p[j]);
n = 3;
87
while (j < x) {
flag = 1;
for (m=1; m<=j; m++)
if (n % p[m] == 0) {
flag = 0;
break;
}
if (flag == 1) {
p[++j] = n;
printf("(%d %d) ", j, p[j]);
}
n++;
}
free(p);
getchar();
return 0;
}
しかし、素数を小さい方から100個見つけるのなら、C++ で次のようにするのが良い。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int p[101];
int ind, n;
p[1] = 2;
ind = 1;
n = 3;
while (ind < 100) {
bool flag = true;
for (int i = 1; i<= ind; i++) {
if (n % p[i] == 0) {
flag = false; break;
}
}
if (flag) {
ind++; p[ind] = n;
}
n += 2;
88
}
for (int i=1; i<=ind; i++) {
cout << p[i] << ’\t’;
if (i % 10 == 0) cout << endl;
}
getchar();
return 0;
}
私のホームページに vector を使った別のプログラムを載せてあります。
例題:配列の要素を並べ替える(ソートと言います)プログラムを作れ。次は入力した数値を大
きい順に並べ替える BASIC のプログラムである。
10 REM データの並べ替え
20 INPUT "データの個数";N
30 DIM A(N)
40 FOR J=1 TO N
50
INPUT A(J)
60 NEXT J
70 FOR J=1 TO N-1
80
FOR K=J+1 TO N
90
100
IF A(J) > A(K) THEN 110
X=A(J):A(J)=A(K):A(K)=X
110
NEXT K
120
PRINT A(J)
130 NEXT J
140 PRINT A(N)
150 END
これは最小値選択法と呼ばれるアルゴリズムです。データの個数が多くなると色々なアルゴリズム
を使う必要があります。この分野は最も研究されている分野の1つです。
これを C に翻訳すると次のようになる。
#include <stdio.h>
int main()
{
double *a, x;
int n, j, k;
printf("データの個数"); scanf("%d", &n);
if ((a = (double *)malloc((n+1)*sizeof(double))) == NULL) {
printf("malloc error!\n");
exit(1);
}
for (j=1; j<=n; i++)
89
scanf("%lf", &a[j]);
for (j=1; j<n; j++) {
for (k=j+1; k<=n; k++) {
if (a[j]>a[k]) continue;
x=a[j]; a[j]=a[k]; a[k]=x;
}
printf("%f\n", a[j]);
}
printf("%f\n", a[n]);
free(a);
getchar();
return 0;
}
continue 文が使われていることに注意。
for (k=j+1; k<=n; k++) {
if (a[j]>a[k]) continue;
x=a[j]; a[j]=a[k]; a[k]=x;
}
で、if 文で a[j] > a[k] であれば、以下の命令達を無視し、すぐ for 文の k++ を実行します。
C++ だと次のようになります。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int n;
double *a;
cout << "データの個数 N="; cin >> n;
a = new double[n];
for (int i=0; i<n; i++) {
cout << "a[" << i+1 << "]="; cin >> a[i];
}
getchar();
for (int i=0; i<n; i++)
for (int k=i+1; k<n; k++)
if (a[i] > a[k]) {
double temp = a[i];
a[i] = a[k];
90
a[k] = temp;
}
for (int i=0; i<n; i++) {
cout << a[i] << ’\t’;
if ((i+1) % 5 == 0) cout << endl;
}
delete[] a:
getchar();
return 0;
}
C では実行時に配列の大きさを決めることができないのでポインタを宣言し(今の場合 double
*a)、malloc() 関数を用いて必要な領域を割り当てましたが、C++ では
double *a;
とポインタを宣言し、new 命令を用いて
a = new double[n];
で、必要な領域を割り当てます。領域が不要になれば
delete[] a;
で、領域を開放します。そうすれば配列と同様な使い方ができる。現在は C++ では、配列の代わ
りに、もっと便利なコンテナクラス vector が使えます。興味があれば、私のホームページを見て
下さい。
例題:5000未満の素数をすべて求めるプログラムを作れ。ここではエラトステネスの篩とい
う有名なギリシャ時代から知られているアルゴリズムを使った例を示します。例えば、C で作ると
次のようになります。
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
int main()
{
int i, k;
int *p;
if ((p = (int *)malloc(5001*sizeof(int))) == NULL) {
printf("malloc error!\n");
exit(1);
}
for (i=1; i<=5000; i++)
p[i] = 1;
p[1] = 0;
91
for (i=4; i<=5000; i += 2)
p[i] = 0;
for (i=3; i <= sqrt((double)5000); i += 2) {
if (!p[i])
continue;
for (k=2*i; k<=5000; k += i)
p[k] = 0;
}
for (i=1; i<=5000; i++)
if (p[i])
printf("%4d ", i);
free(p);
getchar();
return 0;
}
ここでは、malloc() による動的な配列を使っています。C++ で次のようにすることも出来ます。
#include <iostream>
#include <stdio.h>
using namespace std;
int main(int argc, char* argv[])
{
int limit = 5000;
bool p[5000];
for (int i=2; i<limit; i++)
p[i] = true;
p[0] = false;
p[1] = false;
int ind = 2;
while (ind < limit) {
while (ind < limit && p[ind] == false)
ind++;
for (int k=2*ind; k<limit; k += ind)
p[k] = false;
ind++;
}
int n = 0;
for (int i=2; i<limit; i++)
if (p[i] == true) {
cout << i << ’\t’;
92
n++;
if (n % 10 == 0) cout << endl;
}
getchar();
return 0;
}
次は二次正方行列を入力し、その行列式を計算するプログラムである。
#include <stdio.h>
double det(double a[][3])
{
return a[1][1]*a[2][2]-a[2][1]*a[1][2];
}
int main()
{
double a[3][3];
int i;
for (i=1; i<=2; i++)
for (j=1; j<=2; j++) {
printf("a[%d][%d]=", i, j); scanf("%lf", &a[i][j]);
}
printf("行列式 = %f\n", det(a));
getchar();
getchar();
return 0;
}
二次元の配列は double a[3][3] のように宣言する。この場合、
a[0][0], a[0][1], a[0][2], a[1][0], a[1][1], a[1][2], a[2][0], a[2][1], a[2][2]
の 3 × 3 = 9 個にアクセスできるようになる。メモリが無駄になるので上のプログラムは次のよ
うに書く方がよい。引数として関数に渡すときは最初の大きさだけが省略できる。すなわち、関
数 det の定義において、double det(double a[3][3]) と宣言してもよいが、最初の 3 だけは省略で
きる。
#include <stdio.h>
double det(double a[][2])
{
return a[0][0]*a[1][1]-a[1][0]*a[0][1];
}
int main()
{
double a[2][2];
93
int i;
for (i=0; i<2; i++)
for (j=0; j<2; j++) {
printf("a[%d][%d]=", i+1, j+1); scanf("%lf", &a[i][j]);
}
printf("行列式 = %f\n", det(a));
getchar();
getchar();
return 0;
}
配列を使えば、大きな桁数の計算が出来るようになります。ここでは「C -言語とプログラミン
グ-」米田信夫編 産業図書にある π の計算のプログラムを研究する事にします。
#include <iostream>
#include <stdio.h>
#define RADIX 10
#define W 1100
using namespace std;
void init(int a[], int n)
{
int i;
a[0] = n;
for (i=1; i<=W; i++)
a[i] = 0;
}
void add(int a[], int b[], int c[])
{
int t, x, i;
t = 0;
for (i=W; i>=0; i--) {
x = b[i] + c[i] + t;
a[i] = x % RADIX;
t = x / RADIX;
}
if (t != 0) {
cout << "overflow\n";
94
}
}
void sub(int a[], int b[], int c[])
{
int t, x, i;
t = 1;
for (i=W; i>=0; i--) {
x = b[i] + (RADIX - 1 - c[i]) + t;
a[i] = x % RADIX;
t = x / RADIX;
}
if (t != 1)
cout << "overflow\n";
}
int top(int a[], int p)
{
while (p <= W && a[p] == 0)
p++;
return p;
}
void div(int a[], int b[], int n)
{
int t, x, i;
t = 0;
for (i=0; i<=W; i++) {
x = t * RADIX + b[i];
a[i] = x / n;
t = x % n;
}
}
void marctan(int a[], int n, int d)
{
int e[W+1], f[W+1];
int p, i;
init(a, 0);
init(e, n);
95
div(e, e, d);
p = top(e, 0);
add(a, a, e);
i = 3;
while (p <=W) {
div(e, e, d);
div(e, e, d);
div(f, e, i);
if (i % 4 == 1)
add(a, a, f);
else
sub(a, a, f);
p = top(e, p);
i += 2;
}
}
int main()
{
int a[W+1], b[W+1], pi[W+1];
int i;
marctan(a, 16, 5);
marctan(b, 4, 239);
sub(pi, a, b);
printf("%d.", pi[0]);
for (i=1; i<=1000; i++) {
if (i%50==1 && i != 1)
printf(" ");
printf("%d", pi[i]);
if (i % 5 == 0)
putchar(’ ’);
if (i % 50 == 0)
putchar(’\n’);
}
getchar();
return 0;
}
π の小数点以下1000桁は
3.
14159 26535 89793 23846 26433 83279 50288 41971 69399 37510
58209 74944 59230 78164 06286 20899 86280 34825 34211 70679
96
82148 08651 32823 06647 09384 46095 50582 23172 53594 08128
48111 74502 84102 70193 85211 05559 64462 29489 54930 38196
44288 10975 66593 34461 28475 64823 37867 83165 27120 19091
45648 56692 34603 48610 45432 66482 13393 60726 02491 41273
72458 70066 06315 58817 48815 20920 96282 92540 91715 36436
78925 90360 01133 05305 48820 46652 13841 46951 94151 16094
33057 27036 57595 91953 09218 61173 81932 61179 31051 18548
07446 23799 62749 56735 18857 52724 89122 79381 83011 94912
98336 73362 44065 66430 86021 39494 63952 24737 19070 21798
60943 70277 05392 17176 29317 67523 84674 81846 76694 05132
00056 81271 45263 56082 77857 71342 75778 96091 73637 17872
14684 40901 22495 34301 46549 58537 10507 92279 68925 89235
42019 95611 21290 21960 86403 44181 59813 62977 47713 09960
51870 72113 49999 99837 29780 49951 05973 17328 16096 31859
50244 59455 34690 83026 42522 30825 33446 85035 26193 11881
71010 00313 78387 52886 58753 32083 81420 61717 76691 47303
59825 34904 28755 46873 11595 62863 88235 37875 93751 95778
18577 80532 17122 68066 13001 92787 66111 95909 21642 01989
となります。
π
の計算の歴史 アルキメデス(紀元前 283-212 年)は n = 6, 12, 24, 48, 96 の正 n 角形の
長さを計算し、半角の公式
x
sin( ) = ±
2
√
1 − cos x
x
, cos( ) = ±
2
2
√
1 + cos x
2
を繰り返し使うことにより
10
1
<π<3
71
7
という評価を得ました。この値を改良しようとする中世のあらゆる試みは実を結びませんでした
3
が、ついに、アドリアン・ファン・ルーメンが(1580 年に)、アルキメデスの方法を使い、何年も
計算した末に小数点以下 20 桁までの値を得ることに成功しました。ルドルフ・ファン・ケーレン
は(1596, 1616 年に)35 桁まで計算しました。この精度に達するために、ルドルフは n = 62̇60 ま
で計算を続けなければなりませんでした。ドイツでは π のことをルドルフの数と言います。
ライプニッツが公式
arctan x = x −
x5
x7
x9
x11
x3
+
−
+
−
+ ···
3
5
7
9
11
を発見します。この式で x = 1 と置けば、有名なライプニッツの公式
π
1 1 1 1
1
1
=1− + − + −
+
− ···
4
3 5 7 9 11 13
(1682 年)が得られます。この式は美しさとは裏腹に、あまり役に立ちません。なにしろ、50 桁
欲しければ 1050 項の計算を「永遠に続ける努力で」
(オイラー 1737)、行わねばなりません。ずっ
√
と有効なのは tan(π/6) = 1/ 3 を使うことです。すると
√
1
1
1
+
−
+ · · ·)
π = 2 3(1 −
2
3·3 5·3
7 · 33
97
という公式が得られます。この式を用いて、Th.F. ド・ラニーは、「信じられないほどの仕事量を
示している」210 項を加えることによって、127 桁(113 桁目に間違いがあった)を 1719 年に計
算しました。
tan x + tan y
1 − tan x tan y
tan(x + y) =
の式に u = tan x, v = tan y を代入すれば、| arctan u + arctan v| < π/2 のとき、
arctan u + arctan v = arctan(
u+v
)
1 − uv
が得られます。ここで、u = 1/2, v = 1/3 を代入すれば、オイラーの公式 (1737 年)
π
1
1
= arctan + arctan
4
2
3
が得られます。有名なジョン・マチンの公式は次のようにして求めます。u = v = 1/5 と置けば
2 · arctan
1
2/5
5
= arctan(
) = arctan
5
1 − 1/25
12
が得られます。u = v = 5/12 と置けば
2 · arctan
5
5/6
120
= arctan(
) = arctan
12
1 − 25/144
119
となります。そこで、u = 120/119 と置いて、v を
1−u
1
u+v
= 1 となるように、したがって v =
=−
1 − uv
1+u
239
とおきます。これらの式をまとめると、
π
1
1
= 4 · arctan − arctan
4
5
239
となります。この公式を使ってマチンは 100 桁の値を求めたことを W. ジョーンズ(1706 年)が
細かい点は述べずに紹介しています。
演習問題
1. 三次元のベクトルを配列に入力し、その長さ(ノルム)を計算するプログラムを書け。
2. 三次元のベクトルを2個、配列に入力し、その内積を計算するプログラムを書け。
3. 三次正方行列を配列に入力し、その行列式を計算するプログラムを書け。
4. 二次正方行列を二個入力し、その積を計算するプログラムを書け。
4.3
方程式の根の近似解法
昔の数学 B の教科書には正の数 a の平方根
√
a を求める次のようなアルゴリズムが載っている。
まず、縦の長さが 1、横の長さが a の長方形からはじめ面積 a を変えないで、ある規則にした
がって縦と横の長さだけ変え、正方形に近づけていく。正方形が得られたとき、その辺の長さを x
√
とすれば、 x2 = a から平方根 a を求めることができる。長方形を正方形に近づける手順におい
て、k 回目の縦、横の長さをそれぞれ yk , xk とする。縦と横の長さの平均
98
xk+1 =
xk + yk
2
によって k+1 回目の横の長さを定める。xk+1 · yk+1 = a から縦の長さは
yk+1 =
a
xk+1
である。
x0 = a, y0 = 1 からはじめ、上の計算を続け、数列 {xk }, {yk } を求め k を大きくしていくと、
√
それぞれの数列が平方根 a に限りなく近づく。
この算法によるプログラムは、次のようになる。
PRINT "長方形から正方形による平方根の値"
INPUT "縦=1, 横の長さ A を入力 A=";A
PRINT "
X=A : Y = 1
X
Y"
FOR I=1 TO 5
X=(X+Y)/2
Y=A/X
PRINT USING "0.0000000000000000" X, USING "0.0000000000000000" Y
NEXT I
END
実行すると
長方形から正方形による平方根の値
縦=1, 横の長さ A を入力 A=3
X
Y
2.0000000000000000
1.7500000000000000
1.7321428571428572
1.5000000000000000
1.7142857142857142
1.7319587628865978
1.7320508100147274
1.7320508075688772
1.7320508051230272
1.7320508075688774
となる。わずか5回の繰り返しで、小数大15位まで一致し、良い近似値が得られている。
なお、この算法では、平方根を求めるわけであるから、xk , yk の両方を求める必要はない。yk =
を xk+1 =
xk + yk
に代入して
2
xk+1 =
1
a
(xk +
)
2
xk
が得られる。これを使えばよい。これは数学 C にあるニュートン法の特別な場合である。
この算法による C++ のプログラムは、次のようになる。
#include <iostream>
#include <stdio.h>
99
a
xk
using namespace std;
int main(int argc, char* argv[])
{
double a, x, y;
cout << "a="; cin >> a;
getchar();
x = a; y = 1.0;
for (int i=0; i<5; i++) {
x = (x + y) / 2;
y = a / x;
cout << x << "\t" << y << endl;
}
getchar();
return 0;
}
一般的に方程式の近似解を求めるには二分法とニュートン法がよく使われる。これは昔の数学 C
で解説されている。一般的に方程式の近似解を求めるには二分法とニュートン法がよく使われる。
これは数学 C で解説されている。次は二分法により平方根を求める旺文社数学 C の BASIC のプ
ログラムである。
100 REM 二分法
110 EPS = 0.00001
120 INPUT "S=";S
130 DEF FNY(X)=X^2-S
140 IF S<1 THEN P=0:Q=1 ELSE P=1:Q=S
150 FOR N=1 TO 20
160
170
R=(P+Q)/2:NI=N
IF FNY(R)*FNY(Q)<0 THEN P=R ELSE Q=R
180
IF ABS(Q-P)<=EPS THEN 200
190 NEXT N
200 PRINT NI;"回反復 近似値=";R
210 END
実行すると
S=3
18.0 回反復 近似値= 1.7320480346679688
となる。
s の平方根を求めるには f (x) = x2 −s = 0 の解を求めればよい。s の値により f (p) < 0, f (q) > 0
なる p, q を初期値として与える。
100
C++ で二分法により平方根を求めてみよう。s の平方根を求めるには f (x) = x2 − s = 0 の解
を求めればよい。s の値により f (p) < 0, f (q) > 0 なる p, q を初期値として与える。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
double fn(double x, double s)
{
return x * x - s;
}
int main(int argc, char *argv[])
{
double s, eps = 0.0000001;
double p, q, r;
cout << "s="; cin >> s; getchar();
if (s < 1) {
p = 0; q = 1;
} else {
p = 1; q = s;
}
int n;
for (n=0; n<200; n++) {
r = (p + q) / 2;
if (fn(r, s) * fn(q, s) < 0)
p = r;
else
q = r;
if (fabs(q - p) < eps)
break;
}
cout << n+1 << "回反復 近似値=" << p << endl;
cout << "p*p=" << p * p << endl;
getchar();
return 0;
}
次は f (x) = x2 − cos x = 0 の近似解を求める二分法のプログラムである。二分法ではグラフを
描くとか何らかの方法で f (a) × f (b) < 0 となる初期値 a, b を見つけなければならない。
#include <iostream>
101
#include <stdio.h>
#include <math.h>
using namespace std;
double fn(double x)
{
return x * x - cos(x);
}
int main(int argc, char* argv[])
{
double a, b, m, eps = 0.00000000001;
a = 0.5; b = 2.0; m = (a+b)/2;
cout << "fn(" << a << ")=" << fn(a) << endl;
cout << "fn(" << b << ")=" << fn(b) << endl;
while (fabs(fn(a)-fn(b)) > eps) {
m = (a+b)/2;
if (fn(m) > 0) b = m; else a = m;
}
cout << "SOLUTION=" << m << endl;
getchar();
return 0;
}
次に Newton-Raphson 法で f (x) = x2 − cos x の解を求める。
f (x) = 0 の解 x の近似解を xk とし、Taylor 級数展開の x に関して一次までの項で f (x) を近
似すると、
f (x) = f (xk ) + f ′ (xk )(x − xk ) = 0
を得る。これより xk+1 = xk − f (xk )/f ′ (xk ) の逐次計算を行うと、ある条件のもとで真の解に収
束する。
x1
102
x0
方程式によっては、ニュートン法では、求めたい解の近くに初期値を取っても、別の解の近似値
が求まったり、近似解がまったく求められない場合もある。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
double func(double x)
{
return x * x - cos(x);
}
double diff(double x)
{
return 2 * x + sin(x);
}
int main(int argc, char* argv[])
{
double x, dx;
int i, max = 100;
cout << "初期値?"; cin >> x; getchar();
i = 0;
do {
dx = - func(x) / diff(x);
x = x + dx;
} while (i <= max && fabs(dx) > 0.0000001);
cout << "近似解=" << x << endl;
getchar();
return 0;
}
同じことを Java でプログラミングすると次のようになる。
import java.io.BufferedReader;
import java.io.IOException;
import java.io.InputStreamReader;
public class Newton {
static double f(double x) {
return x * x - Math.cos(x);
103
}
static double df(double x) {
return 2 * x + Math.sin(x);
}
public static void main(String[] args) {
int i, max = 100;
double x, dx;
BufferedReader rd =
new BufferedReader(new InputStreamReader(System.in));
try {
String line;
System.out.print("初期値?");
line = rd.readLine();
x = Double.parseDouble(line);
}
catch(IOException e) {
System.out.println("入力エラーが発生しました。");
return; // 入力エラーなら終了する
}
catch(NumberFormatException e) {
System.out.println("実数を入力して下さい。");
return; // 入力ミスなら終了する
}
System.out.println("
i
x
f(x)
i=0;
do {
dx = - f(x) / df(x);
x += dx;
System.out.println(i+","+x+","+f(x)+","+df(x));
i++;
} while((i <= max) && (Math.abs(dx) > 0.00000001));
System.out.println("解="+x);
}
}
実行すると
初期値?2
i
x
f(x)
df(x)
0,1.1004523758499203,0.75780251717421,3.0923172164951502
1,0.8553926151203572,0.07577432278402518,2.465613729520581
2,0.8246601764269862,0.0012578634935394017,2.383637503039626
3,0.8241324688825565,3.730086470810079E-7,2.382223774388821
104
df(x)");
4,0.8241323123025361,3.26405569239796E-14,2.382223354880568
5,0.8241323123025225,2.220446049250313E-16,2.3822233548805314
解=0.8241323123025225
となる。
注
上の図は
prologues := 1;
defaultfont := "rptmr";
defaultscale := 1.5;
beginfig(1);
u=1cm;
vardef f(expr x)=x*x-cosd(x) enddef;
vardef g(expr x)=2*x+sind(x) enddef;
drawarrow (-2.5u,0)--(2.5u,0);
drawarrow (0,-1.1u)--(0,3.5u);
numeric a, b, c;
a := -2.2;
draw (a*u,f(a)*u)
for i=-2.2 step 0.2 until 2.3: ..(i*u, f(i)*u) endfor;
b := 2-f(2)/g(2);
c := b-f(b)/g(b);
draw (2u,0)--(2u,f(2)*u)--(b*u,0)--(b*u,f(b)*u)--(c*u,0);
label.bot(btex $x_0$ etex, (2u,0));
label.bot(btex $x_1$ etex, (b*u,0));
endfig;
end.
というプログラムで、 MetaPost を使って作画した。MetaPost は John Hobby が作った META-
FONT 類似のソフトで、 PostScript(EPS)形式のファイルを出力します。 1995 年から AT&T ライ
センス不要のフリーソフトになりました。METAPOST は METAFONT をお手本にして作られた
プログラミング言語で、 METAFONT は白黒の図形しか作れないのに対して、METAPOST なら
フルカラーの図形が作れるなどの細かい違いはありますが、共通している部分が多いです。Meta-
Post は LATEX を導入すれば、多分自動的にインストールされている。METAFONT は Donald E.
Knuth が作ったフォント作成ツールで、 METAFONT を使って色々な gf 形式のフォントが作れ
ます。 METAFONT を使い、色々の人の無料奉仕で、アラビア語やロシア語やギリシャ語や漢文
など色々な言語が LATEX で使えるようになっています。更に、LATEX は数式や楽譜や化学式など
が綺麗に書ける無料のワープロソフトで、数学の先生は全員使っています。このテキストも LATEX
で作っています。
発展問題:ルート2の計算の歴史
√
√
√
2 のもっと上手な計算法 a + b において、b が a に比べて小さいものと考えると、 a + b ≈
√
√
√
a となります。 a + b と a との差である未知量 δ を考えて近似を良くしようとすると、
√
√
a+b= a+δ
105
より
√
√
2
a + b = ( a + δ) = a + 2 aδ + δ 2
√
となります。δ が小さいから、δ 2 を無視すれば、δ = b/(2 a) が得られ、したがって、
√
a+b≈
√
b
a + √ (|b| ≪ |a| のとき)
2 a
√
が得られます。さて、 2 の近似値として v = 1.4 から始めることにして、a = v 2 , b = 2−a = 2−v 2
とおきます。すると、上の式より新しい近似値として
v+
2 − v2
1 2
= ( + v)
2v
2 v
が得られ、この公式を次々と使うことが出来て、
1.4
1.414285
1.4142135642
1.4142135623730950499
1.4142135623730950488016887242096980790
が得られます。まったく同じ計算を 60 進法で行うことができ、1, 25 (即ち 1 +
24
60
51
602
25
60 )
10
603 )
で始めると
1, 24, 51, 10 (即ち 1 + +
+
が得られます。この値が紀元前 1900 年のバビロニアの粘土
板で見つかるのです。
1450 年頃アルカルサーディーが上の式を改良します。
√
√
b
a+b= a+ √ +δ
2 a
を考えて、平方を計算すると、
a+b=a+b+
√
b2
bδ
+ 2 aδ + √ + δ 2
4a
a
が得られます。最後の2項を無視すれば、
√
√
b
b2
a+b≈ a+ √ − √
2 a 8 a3
√
が得られます。 2 に対して今度は
v+
2 − v2
4 − 4v 2 + v 4
3v
3
1
−
=
+
−
2v
8v 3
8
2v 2v 3
が得られ、v = 1.4 に次々と使っていけば、速く収束する近似列
1.4
1.4142128
1.41421356237309504870
1.41421356237309504880168872420969807856967187537694807317643
√
が得られます。上の近似式は a で割って、 ab を x で置き換えると極めて簡潔な式
√
x √
x x2
1+x≈1+ , 1+x≈1+ −
2
2
8
106
になります。もっと精密な近似式
√
1
1
1
5 4
1 + x = 1 + x − x2 + x3 −
x + ···
2
8
16
128
をニュートンが 1665 年に見つけます。
発展問題:複素関数のニュートン法
Newton-Raphson 法は複素関数に対してもまったく同様に使えます。
f (z) = 0 の解 z の近似解を zk とし、Taylor 級数展開の z に関して一次までの項で f (z) を近
似すると、
f (z) = f (zk ) + f ′ (zk )(z − zk ) = 0
を得る。これより zk+1 = zk − f (zk )/f ′ (zk ) の逐次計算を行うと、ある条件のもとで真の解に収束
する。
方程式によっては、ニュートン法では、求めたい解の近くに初期値を取っても、別の解の近似値が
√
1
3
求まったり、近似解がまったく求められない場合もある。f (z) = z 3 −1 の根は z = 1, z = − ±
i
2
2
です。初期値がどの根に収束するか、色分けすると面白い図ができます。
この図を描く C++ のプログラムは Form に PictureBox を貼り付け、サイズを 400 * 400 にし
ます。
picturebox1 のイベント paint を次のように定義する。
private: System::Void pictureBox1_Paint(System::Object^ sender,
System::Windows::Forms::PaintEventArgs^ e) {
Graphics^ g = e->Graphics;
Pen^ pen = gcnew Pen(Color::Blue);
int x0 = 200, y0 = 200;
int gx, gy;
double limit = 0.000001;
for (double x=-1; x<=1; x+=0.005) {
for (double y=-1; y<=1; y+=0.005) {
double a=x, b=y;
while (fabs((a-1)*(a-1)+b*b)>limit
107
&& fabs((a+0.5)*(a+0.5)+(b-sqrt(3.0)/2)*(b-sqrt(3.0)/2))>limit
&& fabs((a+0.5)*(a+0.5)+(b+sqrt(3.0)/2)*(b+sqrt(3.0)/2))>limit) {
double c1 = a*a*a-3*a*b*b-1;
double c2 = 3*a*a*b-b*b*b;
double d1 = 3*(a*a-b*b);
double d2 = 6*a*b;
a = a - (c1*d1+c2*d2)/(d1*d1+d2*d2);
b = b - (c2*d1-c1*d2)/(d1*d1+d2*d2);
}
if (fabs((a-1)*(a-1)+b*b)<=limit)
pen = gcnew Pen(Color::Blue);
else if (fabs((a+0.5)*(a+0.5)+(b-sqrt(3.0)/2)*(b-sqrt(3.0)/2))<=limit)
pen = gcnew Pen(Color::Red);
else
pen = gcnew Pen(Color::Yellow);
gx = (int)(x0 + 200*x);
gy = (int)(y0 - 200*y);
g->DrawLine(pen, gx, gy, gx+1, gy+1);
}
}
}
そして、Form1.h の先頭に
#include <math.h>
と打ち込みます。
実行すると
108
となります。
演習問題
1. 立方根の近似値を計算するプログラムを書け。
2. ニュートン法によって 1, 2, ... , 20 の平方根の近似値を計算し、結果を数表の形に見やすく
出力するプログラムを書け。
3. 関数 f (x) = x2 − cos x のグラフを描くプログラムを書け。
4.4
数値積分
定積分を数値的に求める方法が数学 C で説明されている。区分求積法、数値積分法として台形
公式とシンプソンの公式である。区分求積法は定積分の定義であるが、普通は収束が速い台形公式
∫b
やシンプソンの公式が使われる。 a f (x)dx の値を近似的に計算することを考える。
いま、区間 [a, b] を n 等分して、幅が h = (b − a)/n の小区間に分割する。そして、分点を小さ
い方から順に P0 , P1 , P2 , · · · , Pn とし、その x 座標の値を a = x0 , x1 , x2 , · · · , xn = b とする。さら
に、各々の分点 Pk における関数の値を yk = f (xk ) (k = 0, 1, 2, · · · , n) とし、座標 (xk , yk ) の点
を Qk とする。曲線の弧 Qk Qk+1 と三線分 Qk Pk , Pk Pk+1 , Pk+1 Qk+1 で囲まれる図形の面積を弧
Qk Qk+1 を線分 Qk Qk+1 で近似し、台形の面積
h
(yk + yk+1 )
2
で近似する。
これらを加えると、つぎの公式が得られる。
∫ b
h
f (x)dx ≈ {y0 + 2(y1 + y2 + · · · + yn−1 ) + yn }
2
a
Qn
Q0
Q1
Q2
P0
P1
P2
109
Pn
次は台形公式公式による
∫2
1
dx
1 x
の近似値を求める BASIC のプログラムである。
100 REM 台形公式
110 A=1:B=2
120 DEF FNY(X)=1/X
130 N=4
140 H=(B-A)/N
150 S=0
160 FOR K=1 TO N-1
170
S=S+FNY(A+K*H)
180 NEXT K
190 S=(FNY(A)+2*S+FNY(B))*H/2
200 PRINT "S=";S
210 PRINT "log(2)=";LOG(2)
220 END
実行すると
S= 0.6970238095238095
log(2)= 0.6931471805599453
です。
210 PRINT "log(2)=";LOG(2)
は、真の値と比較するために入れています。分割数 N が小さすぎるので、余り良い近似値ではな
いです。
上の BASIC のプログラムを C++ に翻訳すると次のようになります。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
double func(double x) {
return 1.0/x;
}
int main()
{
double a = 1, b = 2;
int n = 4;
double h = (b-a)/n;
double s = 0;
110
for (int k=1; k<n; k++)
s = s + func(a+k*h);
s = (func(a)+2*s+func(b))*h/2;
cout << "s=" << s << endl;
cout << "log(2)=" << log(2.0) << endl;
getchar();
return 0;
}
シンプソンの公式の公式は次のようにして求める。台形公式で n 等分していたのを今度は 2n 等
分する。曲線 y = f (x) の弧 Q2k Q2k+1 Q2k+2 の代わりに点 Q2k , Q2k+1 , Q2k+2 を通る放物線で近
似し面積を求める。
Qn
Q0
Q1
Q2
P0
c=
a+b
2
P1
P2
Pn
とし、点 (a, f (a)), (c, f (c)), (b, f (b)) を通る放物線を y = ux2 + vx + w とする。

2

 f (a) = ua + va + w
2
f (b) = ub + vb + w


f (c) = uc2 + vc + ww
が成り立つ。したがって、
w = f (a) − ua2 − va
v(b − a) = f (b) − f (a) + u(a2 − b2 )
2
(f (a) + f (b) − 2f (c))
u=
(a − b)2
が成り立つ。それらを次の式に代入する。
111
∫b
(ux2 + vx + w)dx
= (b − a){ 13 u(b2 + ab + a2 ) + 12 v(b + a) + w}
a
= (b − a){ 13 u(b2 + ab + b2 ) + 12 v(b + a) + f (a) − ua2 − va}
= (b − a){ 13 u(b2 + ab − 2a2 ) + 12 v(b − a) + f (a)}
= (b − a){ 13 u(b2 + ab − 2a2 ) + 12 (f (b) − f (a) + u(a2 − b2 )) + f (a)}
2
−a2 +2ab−b2
= (b − a){ 31 (a−b)
+ 21 (f (b) + f (a))}
2 (f (a) + f (b) − 2f (c))
2
=
b−a
6 (f (b)
+ f (a) + 4f (c))
したがって、 h = (b − a)/(2n) とすれば、近似値は
h
(y2k+2 + 4y2k+1 + y2k )
3
である。これらを加えると、次の近似公式がえられる。これをシンプソンの公式という。
∫ b
h
f (x)dx ≈ {y0 + 4(y1 + y3 + · · · + y2n−1 ) + 2(y2 + y4 + · · · + y2n−2 ) + y2n }
3
a
つぎはシンプソンの公式による
∫2
1
dx
1 x
の近似値を求める BASIC のプログラムである。
100 REM シンプソンの公式
110 A=1:B=2
120 DEF FNY(X)=1/X
130 N=4
140 H=(B-A)/N
150 S1=0:S2=0
160 FOR K=1 TO N-1 STEP 2
170
S1=S1+FNY(A+K*H):S2=S2+FNY(A+(K+1)*H)
180 NEXT K
190 S=(FNY(A)+4*S1+2*S2-FNY(B))*H/3
200 PRINT "S=";S
210 END
実行すると
S= 0.6931545306545307
となる。
上の BASIC のプログラムを C++ に翻訳すると次のようになります。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
double f(double x) {
112
return 1.0/x;
}
int main()
{
double a = 1, b = 2;
double n = 4;
double h = (b-a)/n;
double s1 = 0, s2 = 0;
for (int k=1; k<n; k += 2) {
s1 = s1 + f(a+k*h);
s2 = s2 + f(a+(k+1)*h);
}
double s = (f(a) + 4*s1 + 2*s2 - f(b))*h/3;
cout << "s=" << s << endl;
cout << "log(2)=" << log(2.0) << endl;
getchar();
return 0;
}
n の値を大きくすると良い近似が得られます。
演習問題
∫1√
1. 定積分 0 1 − x2 dx を、区間の分割数が10、20、
・
・
・、60の場合についてシンプソン
の公式で計算し、誤差を調べるプログラムを書け。ただし、この定積分の値は π4 であり、有
効桁数6桁で表すと0.785396である。
2. 放物線 y = x2 の x = 0 から x = 1 までの弧長を求めるプログラムを書け。
4.5
連立一次方程式の解法
数学 C では連立一次方程式を解くガウスの消去法が説明され BASIC のプログラムが与えられ
ている。つぎの BASIC のプログラムはガウスの消去法で連立一次方程式を解くものである。
100 REM ガウスの消去法
110 INPUT "N=";N
120 DIM A(N,N),B(N)
130 REM 方程式の入力
140 FOR I=1 TO N
150
FOR J=1 TO N
160
170
180
PRINT "A(";USING "#" I;USING "#" J;")=";:INPUT A(I,J)
NEXT J
PRINT "B(";USING "#" I;")=";:INPUT B(I)
190 NEXT I
113
200 REM 前進消去
210 FOR K=1 TO N-1
220
230
REM 行の入れ換え
AMAX=ABS(A(K,K)):KMAX=K
240
250
260
FOR I=K+1 TO N
IF ABS(A(I,K))>AMAX THEN KMAX=I:AMAX=ABS(A(I,K))
NEXT I
270
280
IF AMAX<0.00001 THEN PRINT "解けない":END
IF KMAX=K THEN 340
290
300
310
FOR J=K TO N
T=A(K,J):A(K,J)=A(KMAX,J):A(KMAX,J)=T
NEXT J
320
330
T=B(K):B(K)=B(KMAX):B(KMAX)=T
REM 前進消去の中心部
340
350
360
FOR I=K+1 TO N
A(I,K)=A(I,K)/A(K,K)
FOR J=K+1 TO N
370
380
A(I,J)=a(I,J)-A(I,K)*A(K,J)
NEXT J
390
B(I)=B(I)-A(I,K)*B(K)
400
NEXT I
410 NEXT K
420 IF ABS(A(N,N))<0.00001 THEN PRINT "解けない" : END
430 REM 後退代入
440 FOR K=N TO 1 STEP -1
450
460
FOR I=K+1 TO N
B(K)=B(K)-A(K,I)*B(I)
470
NEXT
480
B(K)=B(K)/A(K,K)
490 NEXT K
500 REM 解の出力
510 PRINT "X ";
520 FOR K=1 TO N
530
PRINT B(K);
540 NEXT K
550 END
実行すると
N=3
A( 1 1 )= ? 1
A( 1 2 )= ? 2
A( 1 3 )= ? -1
B( 1 )= ? -3
114
A( 2 1 )= ? 2
A( 2 2 )= ? 4
A( 2 3 )= ? 1
B( 2 )= ? 0
A( 3 1 )= ? 3
A( 3 2 )= ? 8
A( 3 3 )= ? 1
B( 3 )= ? -3
X 1.0000000000000007 -1.0000000000000002 2.0
となる。上の BASIC のプログラムを C++ に翻訳すると次のようになります。
#include <iostream>
#include <stdio.h>
#include <math.h>
using namespace std;
int main() {
int n;
cout << "n="; cin >> n;
double **A, *B;
A = new double*[n];
for (int i=0; i<n; i++)
A[i] = new double[n];
B = new double[n];
for (int i=0; i<n; i++) {
for (int j=0; j<n; j++) {
cout << "A[" << i << "][" << j << "]="; cin >> A[i][j];
}
cout << "B[" << i << "]="; cin >> B[i];
}
for (int k=0; k<n; k++) {
double AMAX = fabs(A[k][k]); int KMAX = k;
for (int i=k+1; i<n; i++) {
if (fabs(A[i][k]) > AMAX) {
KMAX = i;
AMAX = fabs(A[i][k]);
}
}
if (AMAX < 0.00001) {
cout << "解けない" << endl;
for (int i=0; i<n; i++)
115
delete A[i];
delete A;
getchar();
getchar();
return 1;
}
if (KMAX != k) {
for (int j=k; j<n; j++) {
double T = A[k][j];
A[k][j]= A[KMAX][j];
A[KMAX][j] = T;
}
double T = B[k];
B[k] = B[KMAX];
B[KMAX] = T;
}
for (int i=k+1; i<n; i++) {
A[i][k] = A[i][k] / A[k][k];
for (int j=k+1; j<n; j++) {
A[i][j] = A[i][j] - A[i][k] * A[k][j];
}
B[i] = B[i] - A[i][k] * B[k];
}
}
if (fabs(A[n-1][n-1]) < 0.00001) {
cout << "解けない" << endl;
for (int i=0; i<n; i++)
delete A[i];
delete A;
delete B;
getchar();
getchar();
return 1;
}
for (int k=n-1; k>=0; k--) {
for (int i=k+1; i<n; i++) {
B[k] = B[k] - A[k][i] * B[i];
}
B[k] = B[k] / A[k][k];
}
cout << " X ";
for (int k=0; k<n; k++) {
cout << B[k] << " ";
116
}
cout << endl;
for (int i=0; i<n; i++)
delete A[i];
delete A;
delete B;
getchar();
getchar();
return 0;
}
n次正方行列の逆行列を求めるには連立一次方程式をn個同時に解けばよい。次は逆行列を求め
る BASIC プログラムである。
100 REM ガウスの消去法
110 INPUT "N=";N
120 DIM A(N,N),B(N,N)
130 REM 行列の入力
140 FOR I=1 TO N
150
FOR J=1 TO N
160
170
PRINT "A(";USING "#" I;USING "#" J;")=";:INPUT A(I,J)
IF I=J THEN B(I,J)=1 ELSE B(I,J)=0
180
NEXT J
190 NEXT I
200 REM 前進消去
210 FOR K=1 TO N-1
220
REM 行の入れ換え
230
AMAX=ABS(A(K,K)):KMAX=K
240
250
FOR I=K+1 TO N
IF ABS(A(I,K))>AMAX THEN KMAX=I:AMAX=ABS(A(I,K))
260
270
280
NEXT I
IF AMAX<0.00001 THEN PRINT "解けない":END
IF KMAX=K THEN 340
290
300
FOR J=K TO N
T=A(K,J):A(K,J)=A(KMAX,J):A(KMAX,J)=T
310
317
320
NEXT J
FOR J=1 TO N
T=B(K,J):B(K,J)=B(KMAX,J):B(KMAX,J)=T
323
330
NEXT J
REM 前進消去の中心部
340
350
360
FOR I=K+1 TO N
A(I,K)=A(I,K)/A(K,K)
FOR J=K+1 TO N
370
A(I,J)=a(I,J)-A(I,K)*A(K,J)
117
380
387
NEXT J
FOR J=1 TO N
390
393
B(I,J)=B(I,J)-A(I,K)*B(K,J)
NEXT J
400
NEXT I
410 NEXT K
420 IF ABS(A(N,N))<0.00001 THEN PRINT "解けない" : END
430 REM 後退代入
440 FOR K=N TO 1 STEP -1
450
457
460
FOR I=K+1 TO N
FOR J=1 TO N
B(K,J)=B(K,J)-A(K,I)*B(I,J)
463
470
NEXT J
NEXT
477
480
483
FOR J=1 TO N
B(K,J)=B(K,J)/A(K,K)
NEXT J
490 NEXT K
500 REM 解の出力
510 PRINT "逆行列 "
520 FOR K=1 TO N
527
FOR J=1 TO N
530
533
536
PRINT B(K,J);
NEXT J
PRINT
540 NEXT K
550 END
実行すると
N=2
A( 1 1 )= ? 2
A( 1 2 )= ? 1
A( 2 1 )= ? 3
A( 2 2 )= ? 4
逆行列
0.8000000000000002 -0.20000000000000004
-0.6000000000000001 0.4
となる。
上の BASIC のプログラムを C++ に翻訳すると次のようになります。
#include <iostream>
#include <stdio.h>
#include <math.h>
118
using namespace std;
int main() {
int n;
cout << "n="; cin >> n;
double **A, **B;
A = new double*[n];
for (int i=0; i<n; i++)
A[i] = new double[n];
B = new double*[n];
for (int i=0; i<n; i++)
B[i] = new double[n];
for (int i=0; i<n; i++) {
for (int j=0; j<n; j++) {
cout << "A[" << i << "][" << j << "]="; cin >> A[i][j];
if (i == j)
B[i][j] = 1;
else
B[i][j] = 0;
}
}
for (int k=0; k<n; k++) {
double AMAX = fabs(A[k][k]); int KMAX = k;
for (int i=k+1; i<n; i++) {
if (fabs(A[i][k]) > AMAX) {
KMAX = i;
AMAX = fabs(A[i][k]);
}
}
if (AMAX < 0.00001) {
cout << "解けない" << endl;
for (int i=0; i<n; i++)
delete A[i];
delete A;
for (int i=0; i<n; i++)
delete B[i];
delete B;
getchar();
getchar();
return 1;
}
119
if (KMAX != k) {
for (int j=k; j<n; j++) {
double T = A[k][j];
A[k][j]= A[KMAX][j];
A[KMAX][j] = T;
}
for (int j=k; j<n; j++) {
double T = B[k][j];
B[k][j]= B[KMAX][j];
B[KMAX][j] = T;
}
}
for (int i=k+1; i<n; i++) {
A[i][k] = A[i][k] / A[k][k];
for (int j=k+1; j<n; j++) {
A[i][j] = A[i][j] - A[i][k] * A[k][j];
}
for (int j=k+1; j<n; j++) {
B[i][j] = B[i][j] - A[i][k] * B[k][j];
}
}
}
if (fabs(A[n-1][n-1]) < 0.00001) {
cout << "解けない" << endl;
for (int i=0; i<n; i++)
delete A[i];
delete A;
for (int i=0; i<n; i++)
delete B[i];
delete B;
getchar();
getchar();
return 1;
}
for (int k=n-1; k>=0; k--) {
for (int i=k+1; i<n; i++) {
for (int j=0; j<n; j++) {
B[k][j] = B[k][j] - A[k][i] * B[i][j];
}
}
for (int j=0; j<n; j++) {
B[k][j] = B[k][j] / A[k][k];
}
120
}
cout << " 逆行列: " << endl;
for (int k=0; k<n; k++) {
for (int j=0; j<n; j++) {
cout << B[k][j] << " ";
}
cout << endl;
}
for (int i=0; i<n; i++)
delete A[i];
delete A;
for (int i=0; i<n; i++)
delete B[i];
delete B;
getchar();
getchar();
return 0;
}
演習問題
1. 基本的な考え方はガウスの消去法と同じであるが、前進消去の段階で、対角線より下にある
成分だけでなく、対角線より上にある成分も同時に消去してしまう方法もある。これをガウ
ス・ジョルダンの消去法あるいは掃出し法という。ガウス・ジョルダンの消去法により連立
一次方程式を解くプログラムを書け。
挑戦問題
1. バイオリズムは、身体の周期が23日、感情の周期が28日、知性の周期が33日で、次の
ように計算されます。
physical = sin(T /23 ∗ 2 ∗ P I)
sensitivity = sin(T /28 ∗ 2 ∗ P I)
intellectual = sin(T /33 ∗ 2 ∗ P I)
ここで、T は誕生日から計算する日までの日数で、PI は円周率です。各々の値が正ならば好
調、負ならば不調、零ならば好不調の変わり目の要注意日となります。ある日の前後のバイ
オリズム曲線を描く場合には、上の式で T をその日の前後に変化させて、値を線で結べばよ
いです。バイオリズムを表示するプログラムを作れ。
121
参考: sage によるプログラミング
sage という数式処理ソフトを使うとこのテキストで述べた計算が簡単にできます。


−x2 − x3 + x4 = 0



 x +x +x +x =6
1
2
3
4

2x1 + 4x2 + x3 − 2x4 = −1



 3x + x − 2x + 2x = 3
1
2
3
4
という連立方程式を解くには、sage の notebook で
A
=
Matrix(QQ, [[0, -1, -1, 1],[1, 1, 1, 1],[2, 4, 1, -2][3, 1, -2, 2]])
B = vector([0, 6, -1, 3])
solution = A.solve_right(B)
show(solution)
と入力し、
evaluate のボタンをクリックすれば良いです。

2 5

A= 3 1
5 4
122

4

2 
6
(1)
の逆行列は
A = matrix(QQ, [[2, 5, 4],[3, 1, 2],[5, 4, 6]])
show(A.inverse())
と入力し、evaluate のボタンをクリックすれば良いです。
C や C++ で不定積分を計算するのは大変難しいですが、sage では簡単で
var(’x’)
f(x) = e^x * cos(x)
f_int(x) = integrate(f, x)
show(f_int)
√
と入力すればよく、 (1 − x2 ) の 0 から 1 の定積分は
f(x) = sqrt(1 - x^2)
f_integral = integrate(f, (x, 0, 1))
show(f_integral)
と入力すればいいです。
123
関数のグラフを描くのも簡単で
f(x) = e^x * cos(x)
plot(f, (x, -2, 8))
と入力すれば
を描きます。三次元空間内の曲面のグラフも簡単に描けます。
sage は Windows でも使えるとインターネットには書いていますが、その方法をいくつかやって
みましたが、私はいずれも成功しませんでした。Windows にインストールするのは難しいです。
4万円の怪しげなノートパソコンをパソコン工房で購入し、それに無料の ubuntu という linux
の OS をインストールし、それに sage をインストールして使っています。sage は100種類ぐら
いのフリーのソフトを Python でまとめ上げたソフトだそうです。
linux の世界には無数のフリーソフトがあります。例えば、めちゃくちゃ強い無料の将棋ソフト
もあります。使えるプリンタが限られているとか、ソフトをインストールするのに分からなくて、
インターネットで情報を探し回ったりしなければなりませんが、十年前と比べると ubuntu 12.04
は隔世の感があるくらい使いやすくなっています。昔は linux では、ソフトをインストールし、そ
れが動くことで満足していましたが、現在は各種フリーソフトを駆使し、実際に仕事をする時代に
なりました。光ファイバにインターネットががっていれば、すごく快適です。Windows と随分勝
手が違い、教えてくれる友人がいなければ、時間が無駄に過ぎていき、私のように苦労しますが、
得られるものは多いです。
コンピュータの世界は、いかに楽をするかの戦いの歴史です。テレビはどうして映るか知らなく
ても見えますが、知っている人、作れる人がいなくなれば、テレビは消滅します。コンピュータの
ソフトも作り方を知らなくても使えますが、プログラミングの方法やアルゴリズムを知っていて、
ソフトを自由自在に作れる人がいなくなれば、コンピュータが作り上げた文明も消滅します。子供
達に聞かれれば、この授業で述べたくらいの簡単な仕組みを答えられるようになっていることが教
師になる人達には必要なことだと思います。
124
Fly UP