数値計算についての助言

 〜 微分積分学演習II,複素関数論演習の数値計算レポートに 取り組んでいる人へ 〜

加藤 2006年11月25日
2007年 7月7日 例4の間違いを訂正.

数値計算は便利なものだ. しかし,慣れるまでの最初のうちは, 山のようなエラーに遭遇して 嫌になることがあるかもしれない. 数値計算ができなくて困っている人は, 以下のことを試してみて欲しい. なお,以下に記す内容は 「エラーの嵐から脱出し,正しい計算結果を得るための道筋」に関する 私の経験則である. 数値計算法やfortran文法については, 数値計算の講義資料などを参照されたい.

目次
  1. エラーメッセージをちゃんと読む
  2. 実行時エラーの源はwriteで見つける
  3. 計算結果が正しいことを確かめる
  4. おすすめ文献

1. コンパイルエラーが出たときは,エラーメッセージをちゃんと読む

苦労してプログラムを書き終わったあと,普通 まず最初に遭遇するのはコンパイルエラーだ. ここでめげてはいけない. まずは エラーメッセージが何を言っているのかちゃんと読む ことが必要. コンパイルエラーには様々な種類があるが, エラーメッセージが伝えようとしていることは,要するに

「あんたの書いたプログラムは文法的に間違っている (だから俺は何をどうしていいのかわからなくて, このままじゃ計算はできないよ〜).」

ということだ. エラーメッセージには,大概 プログラムのどこがどう間違っているのか書かれている. これを読み解けば,プログラムをどうなおせばよいかわかる訳だ.

メッセージは英語で書かれている.まずは わからない単語を全部辞書で調べる こと. 耳慣れない単語が多くて最初は戸惑うだろうし, 単語を全部調べても メッセージの内容がわからない場合だってあるだろう. でも,「どうやらこのあたりに文句を言われているようだ」 というあたりがつけばしめたもの. 怪しいと思う箇所をなおしながら試行錯誤することで, 間違いの箇所はたいていわかる.

例として, 2つの実変数 x, y に対して a=2x+y を出力するプログラムを考える.

例1: コンパイルエラーが出るプログラム
(一番左のカラムの数字は説明の都合のためにつけた 行番号である.プログラムには書かれていない.)
 1  program example1
 2  implicit none
 3  real::x,y,a
 4  
 5  x=1.0
 6  y=2.0
 7  a=2x+y
 8  write(*,*) 'x=',x,' y=',y
 9  write(*,*) 'answer=',answer
10  stop
11 
12  end program example1

このプログラムをコンパイルすると, 次のエラーメッセージが出る (コンパイラによって文言は多少違うかもしれない). 青字の部分は,エラーが出たところの行番号を示している.

0007) a=2x+y
 Variable X follows another operand
0009) write(*,*) ' answer=',answer
 ANSWER must appear in a type declaration

7行目のエラーメッセージは "変数 X が別の被演算子の次に来ている(続いている)" と言っている. もしも「被演算子」の意味がわからなくても variableとfollowの意味さえわかれば, 「変数 x が何かの後に来るのがいけないらしい (-> xの前が怪しそう)」 ということは想像がつく. 7行目の「a=2x+y」を「a=2*x+y」になおすとこのエラーはでなくなる.

9行目のエラーメッセージは "ANSWERは型宣言に現れねばいかん" と言っている. 宣言について文句を言われているようなので 9行目と宣言文(3行目)を見かえしてみると, 9行目の answer という変数が3行目で宣言されていないことがわかる. 3, 7行目の「a」を「answer」にかえるか, あるいは9行目を 「write(*,*) ' answer=',a」 にかえるかすれば,エラーは出なくなる.

エラーメッセージに指摘された行を調べると 文法ミスが見つかることが多い. ただし,上の例の9行目エラーのように, エラーの源が指摘された行よりはるか 前に存在することもある. エラーの源がわからないときは, 指摘された行より前にある文も疑ってかかる とよい.

例1ではエラーは2個しか出なかったが,長いプログラムの場合は 100個くらいのエラーがどばどば出ることもある. そういうとき,後ろの方のエラーは前のエラー (行番号が若いもの) の 副作用として出ていることが多い. したがって, エラーの嵐に出会ったときは上のエラーから順番につぶしていく と効率がよい.

目次に戻る

2. 実行時エラーの源はwriteで見つける

コンパイルがエラーなく通れば,実行ファイル (Windowsなら *.exe) ができる. これを実行すれば計算結果が得られるはずだが, 一発でうまく実行できず,エラーが出て計算が止まることもある. これはよくあることだ. こんなことで数値計算をあきらめてはいけない. 実行時エラーが出た場合は,まずは エラーメッセージをちゃんと読む こと. もしもエラーメッセージをちゃんと読んでも なぜエラーになるのかわからない場合は, write文を使ってどこまで計算が進んでいるのか確認するとよい.

再び例で説明する.次の例2を実行すると, 実行時エラーとして「Floating point divide by zero」という メッセージが出て, 計算結果は何も出力されない. このエラーメッセージは「ゼロ割り」を指摘しているのだが, 以下ではこのメッセージの内容がわからない場合の対応を説明する.

例2: 実行時エラーの出るプログラム
 1  program example2
 2  implicit none
 3  real::x,y,a,b,answer
 4  
 5  x=1.0
 6  y=2.0
 7  a=2*x+y
 9  b=1.0/(a**2-16.0)
10  answer=a*b
11  write(*,*) 'x=',x,' y=',y
12  write(*,*) 'Answer=',answer
13  stop
14
15  end program example2

以下の例では,エラーの原因を探すために 各行の間にwrite文を挿入した.

例3: 例2のプログラムに write文を追加したもの.
 1  program example3
 2  implicit none
 3  real::x,y,a,b,answer
 4  
 5  x=1.0
 6  y=2.0
 7 write(*,*) 'ここまではOKかな? (1) x=', x, ' y=',y 
 8  a=2*x+y
 9  write(*,*) 'ここまではOKかな? (2) a=',a 
10  b=1.0/(a**2-16.0)
11 write(*,*) 'ここまではOKかな? (3) b=',b 
12  answer=a*b
13  write(*,*) 'x=',x,' y=',y
14  write(*,*) 'Answer=',answer
15  stop
16
17  end program example3

例3を実行すると,例2と同様に実行時エラーが出て計算は止まる. ただし,出力として次のものが得られる.

  ここまではOKかな? (1) x=1.00000     y=2.00000    
  ここまではOKかな? (2) a=4.00000    

7,9行目のwriteは出力されたが, 11行目のwriteの結果「ここまではOKかな? (3) b=...」が出力されない. このことから,9行目までは計算が進んでいて, 10行目か11行目のどちらかで計算が止まったのだということがわかる. このプログラムでは10行目で 1 を a**2-16.0=0 で割っている. これがエラーの原因である.

例3のように write文を使って次の2点を確認 すれば,大概の実行時エラーの原因はわかる.

ただし,いつもこの方法でうまくいく訳ではない. やっかいなのは,実行時エラーが出たり出なかったりする場合だ. よくあることとして 「write文がないときは実行時エラーが出るのに, write文を追加すると実行時エラーが出なくなってしまう」 現象がある. こういうときには subroutine間でのデータの受け渡しのところで 変なミスをしているとか, 配列をどこかで壊しているとか, そんな類の間違いをしていることが多い. このような間違いはwrite文の周辺を探しても発見できない. プログラム全体を疑ってかかる必要がある.

エラーの源を探す際には, コンパイラオプションを活用するのも有効だ. コンパイルするときに使う命令(f77, g77, ftn95, fortなど)には, いろんなオプションが用意されている. 「警告を沢山だしてくれ (文法ミスじゃなくても怪しそうな所は教えてくれ)」とか, コンパイラによっては「実行時に配列を壊したときは教えてくれ」とか, そういったことをお願いできる機能がある. オプションの機能と使い方はコンパイラによって異なる. 使ってみようという人は各自コンパイラの help 等で調べよ.

目次に戻る

3. 計算結果が正しいことを確かめる

コンパイルエラーもなくなった,実行時エラーもなくなった... ここでよかったよかったと安心して そのままレポートを提出してはいけない. なぜなら,

エラーが出なくても計算結果が正しいとは限らない

から. 数値計算の間違いをなくす作業は, エラーがなくなった段階から始まると言ってもよい.

次の例4は,入力 x に対して x の3乗根 (1/3乗)を出力するプログラムである. この例ではコンパイルエラーも実行時エラーも出ない. しかし,計算結果は間違っている.

例4: エラーは出ないが間違った答えを出すプログラム
 1  program example4
 2  implicit none
 3  real::x,answer
 4  
 5  x=8.0
 6  answer=x**(1/3)
 7  write(*,*) 'x=',x,' answer=',answer
 8  stop
 9
10  end program example4

実行結果
 x=8.00000 answer=1.00000

計算結果が間違っているという表現は正しくないかもしれない. 計算機はプログラムに書かれた命令どおり正しく計算をしている. しかしその計算は人間の意図したものとは異なっている. 人間の意図通り計算すれば,answerは 「x**(1/3)=8**(1/3)=2」となるはずだが, 計算機の出す答えは「x**(1/3)=8**(0)=1」である. この例のように, エラーはでないが プログラムには間違いがある というのはよくあることだ.

上の例は数値計算しなくても正解がわかる問題だったので, 計算結果が間違っていることはすぐにわかった. しかし,普通は数値計算で扱うのは答えがわからない問題である. 答えがわからない問題について 計算結果が正しいことを確認するのはそう簡単ではない. 難しい作業ではあるが, 計算結果が正しいかどうか確認することは プログラムを書くこと以上に大切な作業 である.

目次に戻る

3-1 計算結果が正しいかどうか調べるときのヒント

計算結果が正しいことを確認するための万能な方法はない. 問題に応じてもし可能なら,次のことを試すとよい.

目次に戻る

3-2 どこかにあるはずの間違いを探す方法

計算結果はどうも正しくないようだ... プログラムにはどこかに間違いがあるみたいだ... でもどこが間違っているのかわからない. そんなとき,私なら次のことを試す. 間違いを発見するための万能な方法は無い. 基本となる 出発点は全てを疑ってかかること だ.  「今目の前にあるプログラムを書いた奴 (=過去の自分) は, 自分よりも100倍まぬけでおっちょこちょいで信用できない奴なのだ, 考えられんようなお馬鹿な間違いがどこに潜んでいても不思議はない」  そう思ってプログラムを見返すと,間違いが見つかりやすいような気がする.

目次に戻る

4. おすすめ文献

  1. 数値計算I, II, IIIの講義資料. わからないことがあれば,まずは講義資料をちゃんと読むこと.
  2. 浦昭二, Fortran77入門, 培風館, 1990. 情報量が多く辞書的に使えるので便利.
  3. For the Better Fortran Programming (Programming Style Guide for Fortran User)
    北大電子研,一宮さんによる解説 (日本語). どのようにプログラムを書けば, バグ (プログラミングでの間違い) が入りにくく,とりやすくなるか? 教科書には載っていないような様々な実践的なポイントについて解説されている. 学部生には難しい箇所があるかもしれないが, 将来のため一読しておいて損はない.
目次に戻る