下記のうち、いずれか１つの課題を選択し、プログラムを作成してレポートにまとめ提出すること。レポートは、文章で記述した報告書としてまとめること。まとめ方および提出先については、来週指示します。たとえ、最終的にうまくいかない場合でも、出来たところまでで問題点を記述して、期限内に提出すること。期限内に提出されない場合には単位を認めません。
レポートには、最低限
１）プログラムリスト
２）V(40,i) i=1,50の値
３）V(j,40) j=1,50の値
を含めなさい。
(1,1) V=0 (99,1)
+---------------------+ 外側の電極は正方形で電位が零、 内側の電極は
| | 図に示す板状で電位が１である。
| | 形状の対称性より、1/4 の部分のみ計算すればよい。
| V=1 |
| +=========+ | 電位を二次元配列 V(50,50) で確保したとき、
| (25,50) (75,50) | 電極位置以外の電位 V(i,j) は
| | V(i,j) = ( V(i-1,j) + V(i+1,j)
| | +V(i,j-1) + V(i,j+1) )/4.0
+---------------------+ としてまわりの電位から計算できる。この関係を
(1,99) (99,99) 用いて適当な初期値 (たとえば V(i,j)=0.0)から
はじめて変化がなくなるまで逐次計算を繰り返す(緩和法)。
なお対称性より、 V(i,51) = V(i,49)、V(51,j) = V(49,j) の関係がある。
プログラミングのヒント：
電位分布V(x,y)は、空間に電荷が存在しなければ、
xおよびyに関する２次微分の値が0でなければなりません。
空間を、等間隔で区切って、その値をV(i,j)で表わすと、
注目する点とその上下左右の点での電位の間に、
V(i,j) = ( V(i-1,j) + V(i+1,j) +V(i,j-1) + V(i,j+1) )/4.0
の関係が成り立つことが数学的に示されます。
計算は、最初に適当な値（例えば0）を電位分布V(i,j)に与えておき、
上の式を用いて、修正していきます。具体的には、
（０）電圧を記憶しておく配列を確保する。例えば
real*8 V(50,50), Vt(50,50)
（１）全てのVt(i,j)に0.0を代入する
（２）電極での値を与える。例えば、
do i=1,50
Vt(1,i)=0.0
Vt(i,1)=0.0
end do
do i=25,50
Vt(i,50)=1.0
end do
（３）Vtの値をVにコピーする。
（４）電極以外の点での電圧を式にしたがってより正しい値Vtを計算する。
do j=2,49
do i=2,49
Vt(i,j) = ( V(i-1,j) + V(i+1,j) +V(i,j-1) + V(i,j+1) )/4.0
end do
end do
（５）Vt(50,i) i=2,49およびVt(i,50) i=2,24の計算を行う。
どのようにしたら良いかを考えよう。
（対称性より、 V(i,51) = V(i,49)、V(51,j) = V(49,j) の関係がある）
（６）３、４、５を繰り返す。
何回繰り返したら正しい答えが見つかるか考えよう。
（７）Vt(40,i) i=1,50の値、Vt(j,40) j=1,50の値を表示させ、
答えが正しいか確認する。
この方法は、分かりやすいが計算時間が非常に多くかかります。より速く計算結果が得られる方法を考えてみましょう。
課題２：
Monte Calro法を用いて円周率（π）の近似値を求めなさい。
考え方：
(a) ２桁の乱数XiとXi+1を作り、Xi*Xi+Xi+1*Xi+1≧10000を満たす乱数の数をp、Xi*Xi+Xi+1*Xi+1＜10000を満たす乱数の数をqとすれば、πの近似値は4q/(p+q)である。
(b) K桁の乱数の作り方は、Xi（K桁）に定数C（素数）を乗じて、下位K桁をとりだし、これを新しい乱数Xi+1とする。
(1) 乱数を作るプログラムを、Xiを与えてXi+1を返すサブルーチンとして作りなさい。
(2) p+q=500からはじめ、２倍つづ大きくして、１６０００までの近似値の変化を見なさい。
(3) 素数Cを変えて、結果の変化を調べなさい。
(4) 結果を次のような表にして示しなさい。
C p+q p pai
17 500 102 3.184