NumRoom01 の履歴(No.6)
図1.10 剛体球による微分散乱断面積の分析 (教科書11ページ)†
概要†
剛体球 (半径は指定可能) による散乱問題について、散乱角度 θ を指定したときの衝突径数 b, θ とその刻み幅 dθ を指定したときの b の変化量 |db|, およびそれらの積 b|db| を、θ の関数として出力するプログラムです。
プログラムファイル†
hsphere.f (ソースファイル)
hs.cnt (入力ファイル)
hs.dat (計算結果)
hs.outlist (標準出力ファイル)
hsphere.exe (実行形式ファイル: Windows7/10 64bit 版用)
実行方法†
実行ファイル名を hsphere.exe とします (違う名前のファイルを使用する場合は以下を適当に読み替えてください)。コマンドプロンプトあるいはターミナルで
hspere.exe < hs.cnt
とタイプして Enter キーを押せば、入力ファイル hs.cnt の内容に従って計算が行われます。環境によっては
./hspere.exe < ./hs.cnt
とします。作業ディレクトリは実行・入力ファイルを置いているところに変更 (指定) する必要があります。正常に終了したら、0 あるいは stop 0 とターミナルに表示されます。
入力ファイルの説明†
以下において、L* はファイルの *行目を指します。C* は桁数です。* は単一の数字あるいは数字の範囲を表します。たとえば L10 は「10行目」のことであり、C1-4 は「1桁目から4桁目まで」という意味です。
- L1
表題 (C1-60)。標準出力またはターミナル上にそのまま書き出されます。C61 以降は無視されます。 - L2-L3
プログラム中で使用する (OPEN の対象とする) ファイルの機番 (C2-4), ファイルの status (C6-13), ファイル名 (C16-65) を指定します。ファイル名はパス付きで指定することが可能です。サンプルでは 2つのファイルを指定していますが、機番およびファイルの重複がなければ、いくつ指定しても構いません。ただし機番として 5 を指定することは避けましょう。また、100 以上の機番を指定することはできません。C1 を空白以外にすると、その行の指定内容はスキップされます。 - L4
機番として 100 以上または 0 よりも小さい値を指定することで、ファイル指定処理の終了を宣言します。 - L5
人間用の見出し。プログラム中では何の処理もされません。 - L6
ユーザのコメント (C1-50)。C51 以降は無視されます。 - L7 剛体球の半径 (C1-10) を F 型の実数で指定します。
- L8
人間用の見出し。プログラム中では何の処理もされません。 - L9
計算結果を出力するファイル機番、角度の最小値、同最大値、同刻みを、それぞれ 10桁の F 型実数で指定します。機番は、プログラムの中で整数値 KIBO に変換されます。この機番のファイルが L2-3 (冒頭部) で指定されていれば、そのファイルに結果が出力され、指定されていない場合には、fort.[KIBO] というファイルが生成されます ([KIBO] には数字が入る)。角度に関する数値は、degree 単位で指定します。 - L10
桁数を見やすくするための補助。プログラム中では何の処理もされません。
出力の説明†
以下では、掲載している入力ファイル hs.cnt をそのまま用いた場合に生成されるファイル名に基づいて説明を行います。これらのファイル名は (入力ファイルの名前も)、ユーザが自由に書き換えることが可能です。計算内容を変更し、それ以前に作成したファイルを保持したい場合には、出力されるファイル名を変更しましょう。その際には、できるだけ計算内容が読みとれるようなファイル名にすることをお勧めします。
hs.outlist†
標準出力ファイル。計算内容が出力されます。意図した計算内容になっているか確認する癖をつけましょう。また、エラーが生じたときには、その内容が簡単に書き出されます。
たとえば剛体球の半径に負の値を指定してみると、終了コード 1 で計算がストップします。このとき、hs.outlist には
-- ERROR IN INPUT --:
HSRAD <= 0
と表示されるはずです。ただし、あらゆる入力エラーに対してこのような表示がなされるわけではありません。入力パラメータを正しく指定する (物理的に意味のある数字を与える) ことは、基本的にはユーザの責任でなすべきことです。
hs.dat†
左から順に、散乱角 θ (単位: degree), b, |db|, b|db| が出力されます。hs.cnt を用いて計算した結果 hs.dat が、図1.10 にプロットされているものです。
注意点と補足†
- このプログラムは、DO ループを利用して角度 θ を変化させ、θ ごとに b や |db| を解析的に計算し、ファイルに出力するという単純なものです。ただし、プログラム中で使用するファイル名の指定や、入力パラメータの与え方など、第2章以降で使用するプログラムでも使用する基本的な内容が含まれていますので、できればプログラムの中身を丁寧に見て、しっかりと習得してください。
- なお、ファイル名の指定は、井芹康統氏 (千葉経済短大) が作成されたサブルーチン FOPEN によって行っています。FOPEN を用いることで、計算の入出力に使用するファイル名を変更したい場合に、いちいちソースコードを書き換える必要がなくなります。FOPEN は、色々なケースに対して計算を繰り返すようなとき、絶大な威力を発揮します。
発展課題†
- hs.dat に、微分断面積を追加してみましょう。もちろん正解は教科書の式(1.17)にあるとおりですが、出力されている b|db| から微分断面積を計算してみてください。
- 角度の刻みを変えると |db| や b|db| は変化しますが、微分断面積は変化しないことを確認しましょう。
- 反射角度を入力ファイルで与えて、その角度以上に粒子が散乱する断面積 (反射断面積) を計算する機能を追加してください。反射断面積は解析的に与えられますが、微分断面積を積分することでも得られます (教科書 p.12 を参照)。2つの方法で算出した反射断面積が一致することを確認してください。
※数値積分については、第3章の数値計算で紹介するプログラムを必要に応じて参照してください。