Numerical integration: Difference between revisions

Content added Content deleted
m (→‎{{header|Sidef}}: Fix link: Perl 6 --> Raku)
(→‎{{header|Raku}}: Move some calculations out of hot loops, more efficient looping constructs, arguably less idiomatic, but less horribly slow)
Line 4,183: Line 4,183:
sub leftrect(&f, $a, $b, $n) {
sub leftrect(&f, $a, $b, $n) {
my $h = ($b - $a) / $n;
my $h = ($b - $a) / $n;
$h * [+] do f($_) for $a, $a+$h ... $b-$h;
my $end = $b-$h;
⚫
my $sum = 0;
loop (my $i = $a; $i <= $end; $i += $h) { $sum += f($i) }
$h * $sum;
}
}

sub rightrect(&f, $a, $b, $n) {
sub rightrect(&f, $a, $b, $n) {
my $h = ($b - $a) / $n;
my $h = ($b - $a) / $n;
⚫
my $sum = 0;
$h * [+] do f($_) for $a+$h, $a+$h+$h ... $b;
loop (my $i = $a+$h; $i <= $b; $i += $h) { $sum += f($i) }
$h * $sum;
}
}

sub midrect(&f, $a, $b, $n) {
sub midrect(&f, $a, $b, $n) {
my $h = ($b - $a) / $n;
my $h = ($b - $a) / $n;
my $sum = 0;
$h * [+] do f($_) for $a+$h/2, $a+$h+$h/2 ... $b-$h/2;
my ($start, $end) = $a+$h/2, $b-$h/2;
loop (my $i = $start; $i <= $end; $i += $h) { $sum += f($i) }
$h * $sum;
}
}

sub trapez(&f, $a, $b, $n) {
sub trapez(&f, $a, $b, $n) {
my $h = ($b - $a) / $n;
my $h = ($b - $a) / $n;
my $partial-sum += f($_) * 2 for $a+$h, $a+$h+$h ... $b-$h;
my $partial-sum = 0;
$h / 2 * [+] f($a), f($b), $partial-sum;
my ($start, $end) = $a+$h, $b-$h;
loop (my $i = $start; $i <= $end; $i += $h) { $partial-sum += f($i) * 2 }
$h / 2 * ( f($a) + f($b) + $partial-sum );
}
}

sub simpsons(&f, $a, $b, $n) {
sub simpsons(&f, $a, $b, $n) {
my $h = ($b - $a) / $n;
my $h = ($b - $a) / $n;
my $h2 = $h/2;
my $h2 = $h/2;
my ($start, $end) = $a+$h, $b-$h;
my $sum1 = f($a + $h2);
my $sum1 = f($a + $h2);
my $sum2 = 0;
my $sum2 = 0;
loop (my $i = $start; $i <= $end; $i += $h) {
for $a+$h, *+$h ... $b-$h {
$sum1 += f($i + $h2);
$sum1 += f($_ + $h2);
$sum2 += f($i);
⚫
$sum2 += f($_);
}
}
($h / 6) * (f($a) + f($b) + 4*$sum1 + 2*$sum2);
($h / 6) * (f($a) + f($b) + 4*$sum1 + 2*$sum2);
}
}

sub integrate($f, $a, $b, $n, $exact) {
sub integrate($f, $a, $b, $n, $exact) {
⚫
my @r0;
my $e = 0.000001;
my $e = 0.000001;
@r0.push: "$f\n in [$a..$b] / $n\n";
my $r0 = "$f\n in [$a..$b] / $n\n"
@r0.push: ' exact result: '~ $exact.round($e);
~ ' exact result: '~ $exact.round($e);


my (@r1,@r2,@r3,@r4,@r5);
my ($r1,$r2,$r3,$r4,$r5);
my &f;
my &f;
EVAL "&f = $f";
EVAL "&f = $f";
my $p1 = Promise.start( { @r1.push: ' rectangle method left: '~ leftrect(&f, $a, $b, $n).round($e) } );
my $p1 = Promise.start( { $r1 = ' rectangle method left: '~ leftrect(&f, $a, $b, $n).round($e) } );
my $p2 = Promise.start( { @r2.push: ' rectangle method right: '~ rightrect(&f, $a, $b, $n).round($e) } );
my $p2 = Promise.start( { $r2 = ' rectangle method right: '~ rightrect(&f, $a, $b, $n).round($e) } );
my $p3 = Promise.start( { @r3.push: ' rectangle method mid: '~ midrect(&f, $a, $b, $n).round($e) } );
my $p3 = Promise.start( { $r3 = ' rectangle method mid: '~ midrect(&f, $a, $b, $n).round($e) } );
my $p4 = Promise.start( { @r4.push: 'composite trapezoidal rule: '~ trapez(&f, $a, $b, $n).round($e) } );
my $p4 = Promise.start( { $r4 = 'composite trapezoidal rule: '~ trapez(&f, $a, $b, $n).round($e) } );
my $p5 = Promise.start( { @r5.push: ' quadratic simpsons rule: '~ simpsons(&f, $a, $b, $n).round($e) } );
my $p5 = Promise.start( { $r5 = ' quadratic simpsons rule: '~ simpsons(&f, $a, $b, $n).round($e) } );


await $p1, $p2, $p3, $p4, $p5;
await $p1, $p2, $p3, $p4, $p5;
@r0, @r1, @r2, @r3, @r4, @r5;
$r0, $r1, $r2, $r3, $r4, $r5;
}
}

.say for integrate '{ $_ ** 3 }', 0, 1, 100, 0.25; say '';
.say for integrate '{ $_ ** 3 }', 0, 1, 100, 0.25; say '';
.say for integrate '1 / *', 1, 100, 1000, log(100); say '';
.say for integrate '1 / *', 1, 100, 1000, log(100); say '';
Retrieved from "https://rosettacode.org/wiki/Numerical_integration"