この文書の URL は http://www.cc.kyoto-su.ac.jp/~mtkg/lecture/comp_B/2012/report1.html です。

前半のまとめの問題

問1. Zeller の公式

Zeller の公式を利用して西暦2112年9月3日の曜日を求めよ。 西暦年月日 (year, month, day) を引数とする関数 zeller(year, month, day) を定義し、それを使うこと。

問2. 覆面算

次の覆面算の解を求めるプログラムを作成し、実際に解を求めよ。

　　　ＫＹＯＴＯ
＋）　ＯＳＡＫＡ
────────
　　　ＴＯＫＹＯ

ただし、解は次の条件を満たすものであること。

ヒント： 文字ごとに do ループを回せばよい。 異なる文字に同じ数字が入っている場合は cycle 文でループをスキップする。 cycle 文が実行されると処理がその do ループの最後に移動し、結果的にループの次の繰り返しに移る。

! cycle 文の例 -- 奇数だけ画面に表示する do ループ
do i = 1, 10
   if (mod(i, 2) == 0) cycle    ! i が偶数の場合は以下をスキップし、次の i に移る。
   write(*,*) i
end do

問3. ルジャンドル多項式

ルジャンドル多項式 Pn(x) は |x| <= 1 を満たす x に対して次の漸化式で定義される。

P0(x) = 1
P1(x) = x
Pn(x) = [(2*n-1) * x * Pn-1(x) - (n-1) * Pn-2(x)] / n   (n = 2, 3, ...)

P5(x) の値を -1 から 1 までの x について 0.1 刻みで求め、 Gnuplot でグラフをプロットせよ。 任意の n と x についてルジャンドル多項式 Pn(x) の値を返す関数 lp(n, x) を定義し、それを用いること。

■いくつかのルジャンドル関数をプロットした図

問4. 多数桁計算

「整数÷整数」の計算を小数の任意の桁まで行うプログラムを作成せよ。 2つの整数 a, b と桁数 n はキーボードから入力する。 ただし、桁数の上限はあらかじめ決めておいてよい（以下のヒントでは 1000 としている）。

ヒント： 筆算をプログラムに直せばよい。 商の整数部を表す整数 p と小数部を表す整数配列 q(:) を用意し、筆算をループで表現する。

! nmax    桁数の最大
! p       商の整数部
! q(i)    商の小数部, q(i) は商の小数第 i 桁を表す
program main
  implicit none
  integer, parameter :: nmax = 1000
  integer :: p, q(nmax)
  integer :: a, b, n
    :

  ! 結果の表示
  write(*,*) p
  write(*,'(5(X,10I1))') q(1:n)

  stop
end program main

問5. 測地線

地球上の2地点の緯度・経度から2地点間の（大円に沿った）最短距離を求めるプログラムを作成し、京都−ローマ間の距離を求めよ。 大円とは球と球の中心を含む平面との交線のことである。 地球は回転楕円体ではなく半径 6370 km の球であるとしてよい。 京都とローマの緯度・経度は次の値を用いること。

京都　　　北緯35度00分30秒　　東経135度46分16秒
ローマ　　北緯41度53分35秒　　東経 12度28分60秒

ヒント： 球上の異なる2点を A, B、球の中心を O とする。 角 AOB のなす角がわかれば AB 間の最短距離（弧 AB の長さ）が計算できる。 角 AOB のなす角は OA ベクトルと OB ベクトルの内積から求められる（cos の逆関数 acos を利用せよ）。 OA ベクトルと OB ベクトルの内積は A と B の座標がわかれば計算可能である。 半径 r の球面中の緯度θ・経度φにある点の座標 (x, y, z) は r, θ, φ を使って次のように表される。

x = r cosθ cosφ
y = r cosθ sinφ
z = r sinθ

問6. 多数桁計算（Napier 数）

自然対数の底 e （Napier 数）を任意の桁数（小数第 K 位）まで計算するプログラムを作成せよ。 Napier 数の求め方は以下の通り（京大冨田先生の解説より）。

e = 1 + 1/1! + 1/2! + 1/3! + 1/4! + ...     [exp(x) の x=1 でのテイラー展開]

を次のように変形する。

        1      1      1           1       1
e = 1 + - (1 + - (1 + - (1 + ... --- (1 + -) ...)))
        1      2      3          M-1      M

小数第 K 位までの多数桁小数を考え、「多数桁小数÷整数」の計算を次のように繰り返す。

このようにして任意の精度で e が求められる。 小数第 K 位まで求める場合、M は M > K*log(10.0) であれば十分である。 例えば K=100 桁の計算をするには K*log(10.0) = 230.25… より M=231 とすればよい。