Windows10でFORTRANを動かさねばならない(下)

出力結果は非常に単純なものである。
座標値（位置）に対する変位電流など電磁気学的な量を計算するだけである。
要するに、横軸が位置ｘで、縦軸が f(x) というだけの簡単なグラフを描きたいのだが、1980年初期にはグラフを描くのはなかなか大変なので、単に数値だけがだらだらと出力ファイルになって出てくるものである。
もちろん、今は、計算結果をExcelに読み込んで、グラフを描いてシミュレーション結果を評価しているようだった。

このプログラム、1980年代のコンピュータを想定しているので、配列などのサイズが非常に小さく制限されている。
当時は、大型コンピュータといっても、せいぜい1MB程度しかメモリが使えないのが普通であった。
そういう制限の下で、ごちゃごちゃと計算していたのだ。
せっかくWindowsになるのだから、グラフをもっと細かく精密に出したいという当然の要望があった。

実際のプログラムでは、サイズの上限の数値がプログラム中に散在していたのだが、それを１つの定数変数として宣言することで、1箇所直せばサイズを変更できるようにした。
グラフの細かさは100倍くらい楽勝と思われたが、何とか現状の2倍にして欲しいというので、10倍に設定した。

さて、この機能強化により、より細かい精密な評価ができるようになる……….はずだった。
x軸上の測定間隔を半分にしたら、出発点が同じなら、１つおきに元の値と完全一致する筈だ。
でも、一致しなかったのだ。
まあ、計算誤差というのもあるから、ある程度誤差が出るのは仕方がないのであるが、明らかにかなり違う。

想定外の問題が出てしまった。どうしよう。
ということで、デバッガを２つ立ち上げて、計算間隔を丁度倍違うようにして、同じｘ値に対する計算が始まる箇所で止めを、変数の値を全部調べ比較したのだ。
どこで何が違うのかを調べた。
違っても問題ないのかもしれないし、プログラムの意味を知らずに調べているので苦しい。
それでも、デバッガでステップ動作で調べていると、徐々に分かってくるものである。

やっと分かったのは、やはり初期化の問題だった。
内部ループの最初で行うべき初期化が、なぜか省略されていた。
電磁波の反射などの影響を延々と足しこんでいるらしい配列だったのだが、1回前の測定点の値にそのまま足しこんでいるようだった。
恐る恐るループの最初で疑わしい複素配列を初期化すると、結果がまるで変わった。

数字の羅列では良く分からないし、でもExecelに結果を読み込むのも面倒なので、Pythonで比較グラフを表示するツールを作った。
青線が元の計算結果、赤が内部ループにゼロクリアを入れた場合の計算結果である。
計算は、右から左に向かって進んでいく。
振動が消えてしまったが、これで良いのだろうか、悪化したのか、情報がないので分からないので、お客さんに聞いてみた。
昔は、赤い線のようなグラフになっていたそうだ。
このグラフが出るようになってから、急に資料の発掘が進み、30年以上前の資料でグラフが載っているものがあり、振動がないことがわかり、昔の状態に復元できたようで、魔法使いのように扱われ、とても感謝された。

ということで、一件落着になったのだが、こんなに上手くいくことは珍しい。
非常に古いプログラムで、資料が無い、作った人が居ない、その他問題が多いプログラムを扱うのは拒絶するのが正しい。
それでも諸般の事情、大人の事情などで作業しないといけないこともあるが、その場合は細心の注意で望まれたい。

このような作業をおこなうと、非常に力が身に付くことも確かである。
ハッカー養成ギブスとしては良いかも知れない。

（完）