In mathematics,   Faulhaber's formula,   named after Johann Faulhaber,   expresses the sum of the p-th powers of the first n positive integers as a (p + 1)th-degree polynomial function of n,   the coefficients involving Bernoulli numbers.



Task

Generate the first 10 closed-form expressions, starting with p = 0.



Related tasks

   Bernoulli numbers.
   evaluate binomial coefficients.


See also

   The Wikipedia entry:   Faulhaber's formula.
   The Wikipedia entry:   Bernoulli numbers.
   The Wikipedia entry:   binomial coefficients.


sub bernoulli_number($n) {
 
    return 1/2 if $n == 1;
    return 0/1 if $n % 2;
 
    my @A;
    for 0..$n -> $m {
        @A[$m] = 1 / ($m + 1);
        for $m, $m-1 ... 1 -> $j {
            @A[$j - 1] = $j * (@A[$j - 1] - @A[$j]);
        }
    }
 
    return @A[0];
}
 
sub binomial($n, $k) {
    $k == 0 || $n == $k ?? 1 !! binomial($n-1, $k-1) + binomial($n-1, $k);
}
 
sub faulhaber_s_formula($p) {
 
    my @formula = gather for 0..$p -> $j {
        take '('
            ~ join('/', (binomial($p+1, $j) * bernoulli_number($j)).Rat.nude)
            ~ ")*n^{$p+1 - $j}";
    }
 
    my $formula = join(' + ', @formula.grep({!m{'(0/1)*'}}));
 
    $formula .= subst(rx{ '(1/1)*' }, '', :g);
    $formula .= subst(rx{ '^1'» }, '', :g);
 
    "1/{$p+1} * ($formula)";
}
 
for 0..9 -> $p {
    say "f($p) = ", faulhaber_s_formula($p);
}