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.


with javascript_semantics
include builtins\pfrac.e -- (0.8.0+)
 
function bernoulli(integer n)
    sequence a = {}
    for m=0 to n do
        a = append(a,{1,m+1})
        for j=m to 1 by -1 do
            a[j] = frac_mul({j,1},frac_sub(a[j+1],a[j]))
        end for
    end for
    if n!=1 then return a[1] end if
    return frac_uminus(a[1])
end function
 
function binomial(integer n, k)
    if n<0 or k<0 or n<k then ?9/0 end if
    if n=0 or k=0 then return 1 end if
    integer num = 1,
            denom = 1
    for i=k+1 to n do
        num *= i
    end for
    for i=2 to n-k do
        denom *= i
    end for
    return num / denom
end function
 
procedure faulhaber(integer p)
    string res = sprintf("%d : ", p)
    frac q = {1, p+1}
    for j=0 to p do
        frac bj = bernoulli(j)
        if frac_ne(bj,frac_zero) then
            frac coeff = frac_mul({binomial(p+1,j),p+1},bj)
            string s = frac_sprint(coeff)
            if j=0 then
                if s="1" then
                    s = ""
                end if
            else
                if s[1]='-' then
                    s[1..1] = " - "
                else
                    s[1..0] = " + "
                end if
            end if
            res &= s&"n"
            integer pwr = p+1-j
            if pwr>1 then
                res &= sprintf("^%d", pwr)
            end if
        end if
    end for
    printf(1,"%s\n",{res})
end procedure
 
for i=0 to 9 do
    faulhaber(i)
end for
