Numerical integration: Difference between revisions
Content added Content deleted
SqrtNegInf (talk | contribs) m (→{{header|Perl 6}}: moved commentary to top, fixed language version tag) |
SqrtNegInf (talk | contribs) (Added Perl example) |
||
Line 3,338: | Line 3,338: | ||
integrate := integral |
integrate := integral |
||
end;</lang> |
end;</lang> |
||
=={{header|Perl}}== |
|||
{{trans|Perl 6}} |
|||
<lang perl>use feature 'say'; |
|||
sub leftrect { |
|||
my($func, $a, $b, $n) = @_; |
|||
my $h = ($b - $a) / $n; |
|||
my $sum = 0; |
|||
for ($_ = $a; $_ < $b; $_ += $h) { $sum += $func->($_) } |
|||
$h * $sum |
|||
} |
|||
sub rightrect { |
|||
my($func, $a, $b, $n) = @_; |
|||
my $h = ($b - $a) / $n; |
|||
my $sum = 0; |
|||
for ($_ = $a+$h; $_ < $b+$h; $_ += $h) { $sum += $func->($_) } |
|||
$h * $sum |
|||
} |
|||
sub midrect { |
|||
my($func, $a, $b, $n) = @_; |
|||
my $h = ($b - $a) / $n; |
|||
my $sum = 0; |
|||
for ($_ = $a + $h/2; $_ < $b; $_ += $h) { $sum += $func->($_) } |
|||
$h * $sum |
|||
} |
|||
sub trapez { |
|||
my($func, $a, $b, $n) = @_; |
|||
my $h = ($b - $a) / $n; |
|||
my $sum = $func->($a) + $func->($b); |
|||
for ($_ = $a+$h; $_ < $b; $_ += $h) { $sum += 2 * $func->($_) } |
|||
$h/2 * $sum |
|||
} |
|||
sub simpsons { |
|||
my($func, $a, $b, $n) = @_; |
|||
my $h = ($b - $a) / $n; |
|||
my $h2 = $h/2; |
|||
my $sum1 = $func->($a + $h2); |
|||
my $sum2 = 0; |
|||
for ($_ = $a+$h; $_ < $b; $_ += $h) { |
|||
$sum1 += $func->($_ + $h2); |
|||
$sum2 += $func->($_); |
|||
} |
|||
$h/6 * ($func->($a) + $func->($b) + 4*$sum1 + 2*$sum2) |
|||
} |
|||
# round where needed, display in a reasonable format |
|||
sub sig { |
|||
my($value) = @_; |
|||
my $rounded; |
|||
if ($value < 10) { |
|||
$rounded = sprintf '%.6f', $value; |
|||
$rounded =~ s/(\.\d*[1-9])0+$/$1/; |
|||
$rounded =~ s/\.0+$//; |
|||
} else { |
|||
$rounded = sprintf "%.1f", $value; |
|||
$rounded =~ s/\.0+$//; |
|||
} |
|||
return $rounded; |
|||
} |
|||
sub integrate { |
|||
my($func, $a, $b, $n, $exact) = @_; |
|||
my $f = sub { local $_ = shift; eval $func }; |
|||
my @res; |
|||
push @res, "$func\n in [$a..$b] / $n"; |
|||
push @res, ' exact result: ' . rnd($exact); |
|||
push @res, ' rectangle method left: ' . rnd( leftrect($f, $a, $b, $n)); |
|||
push @res, ' rectangle method right: ' . rnd(rightrect($f, $a, $b, $n)); |
|||
push @res, ' rectangle method mid: ' . rnd( midrect($f, $a, $b, $n)); |
|||
push @res, 'composite trapezoidal rule: ' . rnd( trapez($f, $a, $b, $n)); |
|||
push @res, ' quadratic simpsons rule: ' . rnd( simpsons($f, $a, $b, $n)); |
|||
@res; |
|||
} |
|||
say for integrate('$_ ** 3', 0, 1, 100, 0.25); say ''; |
|||
say for integrate('1 / $_', 1, 100, 1000, log(100)); say ''; |
|||
say for integrate('$_', 0, 5_000, 5_000_000, 12_500_000); say ''; |
|||
say for integrate('$_', 0, 6_000, 6_000_000, 18_000_000);</lang> |
|||
{{out}} |
|||
<pre>$_ ** 3 |
|||
in [0..1] / 100 |
|||
exact result: 0.25 |
|||
rectangle method left: 0.245025 |
|||
rectangle method right: 0.255025 |
|||
rectangle method mid: 0.249988 |
|||
composite trapezoidal rule: 0.250025 |
|||
quadratic simpsons rule: 0.25 |
|||
1 / $_ |
|||
in [1..100] / 1000 |
|||
exact result: 4.60517 |
|||
rectangle method left: 4.654991 |
|||
rectangle method right: 4.556981 |
|||
rectangle method mid: 4.604763 |
|||
composite trapezoidal rule: 4.605986 |
|||
quadratic simpsons rule: 4.60517 |
|||
$_ |
|||
in [0..5000] / 5000000 |
|||
exact result: 12500000 |
|||
rectangle method left: 12499997.5 |
|||
rectangle method right: 12500002.5 |
|||
rectangle method mid: 12500000 |
|||
composite trapezoidal rule: 12500000 |
|||
quadratic simpsons rule: 12500000 |
|||
$_ |
|||
in [0..6000] / 6000000 |
|||
exact result: 18000000 |
|||
rectangle method left: 17999997 |
|||
rectangle method right: 18000003 |
|||
rectangle method mid: 18000000 |
|||
composite trapezoidal rule: 18000000 |
|||
quadratic simpsons rule: 18000000</pre> |
|||
=={{header|Perl 6}}== |
=={{header|Perl 6}}== |