修正欧拉法计算
返回
###############################
# 修正欧拉法计算
# (dx/dt)2=-x
# 分解为 dv/dt = -x
# dx/dt = v
###############################
use strict;
my(%t,$t1,$x,$v,$ref);
$x = 1.0;
$v = 0.0;
$t{dt} = 0.01;
$t{tmax} = 100;
open(OUT,">euler1.csv");
print OUT 'NO,x,v,',"\n";
for ($t1=0; $t1<$t{tmax}; $t1+=$t{dt}) {
$t{k00} = $v;
$t{k01} = -$x;
$t{k10} = $v+$t{dt}*$t{k01};
$t{k11} = -($x+$t{dt}*$t{k00});
$x = $x + $t{dt}*($t{k00}+$t{k10})/2.0;
$v = $v + $t{dt}*($t{k01}+$t{k11})/2.0;
print OUT $t1,',';
print OUT $x,',';
print OUT $v,',';
print OUT "\n";
}
close(OUT);
__END__;
返回