下記のプログラムは、perlという簡易言語で書いたものです。 中身を解読しつつ、 自分の知っている言語(FortranやCやJava)に翻訳して動かしてみてください。
#!/usr/bin/perl -w
# ----------------------------------------------------------------
# 振動子の方程式を解く (2000-11-16, 大信田 丈志)
# dx/dt = v, dv/dt = (1/m)F_spring(x) = -x, m=1
# ----------------------------------------------------------------
# 刻み幅の設定
$dt = 0.05;
# 初期条件を読み込む
if( @ARGV == 2 ) {
($x0,$v0) = @ARGV; $tmax = 6.5;
} elsif( @ARGV == 3 ) {
($x0,$v0,$tmax) = @ARGV;
} else {
die "usage: $0 x0 v0 [tmax]\n";
}
# 数値計算
($x,$v) = ($x0,$v0); $t = 0.0;
while( $t < $tmax ) {
printf("%6.3f %f %f\n", $t, $x, $v);
$xh = $x + $v *($dt/2);
$vh = $v + &F_spring($x) *($dt/2);
$xnew = $x + $vh * $dt;
$vnew = $v + &F_spring($xh) * $dt;
($x,$v) = ($xnew,$vnew);
$t += $dt;
}
exit 0;
# ----------------------------------------------------------------
# バネの力
sub F_spring {
my($x) = @_;
return (-$x); # [バネの力] = -kx, k=1
}
# ----------------------------------------------------------------