Faulhaber's formula: Difference between revisions
→{{header|Maxima}}
m (→{{header|zkl}}: comment) |
|||
Line 1,410:
9 : 1/10n^10 + 1/2n^9 + 3/4n^8 - 7/10n^6 + 1/2n^4 - 3/20n^2
</pre>
=={{header|Lua}}==
{{trans|C}}
<lang lua>function binomial(n,k)
if n<0 or k<0 or n<k then return -1 end
if n==0 or k==0 then return 1 end
local num = 1
for i=k+1,n do
num = num * i
end
local denom = 1
for i=2,n-k do
denom = denom * i
end
return num / denom
end
function gcd(a,b)
while b ~= 0 do
local temp = a % b
a = b
b = temp
end
return a
end
function makeFrac(n,d)
local result = {}
if d==0 then
result.num = 0
result.denom = 0
return result
end
if n==0 then
d = 1
elseif d < 0 then
n = -n
d = -d
end
local g = math.abs(gcd(n, d))
if g>1 then
n = n / g
d = d / g
end
result.num = n
result.denom = d
return result
end
function negateFrac(f)
return makeFrac(-f.num, f.denom)
end
function subFrac(lhs, rhs)
return makeFrac(lhs.num * rhs.denom - lhs.denom * rhs.num, rhs.denom * lhs.denom)
end
function multFrac(lhs, rhs)
return makeFrac(lhs.num * rhs.num, lhs.denom * rhs.denom)
end
function equalFrac(lhs, rhs)
return (lhs.num == rhs.num) and (lhs.denom == rhs.denom)
end
function lessFrac(lhs, rhs)
return (lhs.num * rhs.denom) < (rhs.num * lhs.denom)
end
function printFrac(f)
io.write(f.num)
if f.denom ~= 1 then
io.write("/"..f.denom)
end
return nil
end
function bernoulli(n)
if n<0 then
return {num=0, denom=0}
end
local a = {}
for m=0,n do
a[m] = makeFrac(1, m+1)
for j=m,1,-1 do
a[j-1] = multFrac(subFrac(a[j-1], a[j]), makeFrac(j, 1))
end
end
if n~=1 then
return a[0]
end
return negateFrac(a[0])
end
function faulhaber(p)
io.write(p.." : ")
local q = makeFrac(1, p+1)
local sign = -1
for j=0,p do
sign = -1 * sign
local coeff = multFrac(multFrac(multFrac(q, makeFrac(sign, 1)), makeFrac(binomial(p + 1, j), 1)), bernoulli(j))
if not equalFrac(coeff, makeFrac(0, 1)) then
if j==0 then
if not equalFrac(coeff, makeFrac(1, 1)) then
if equalFrac(coeff, makeFrac(-1, 1)) then
io.write("-")
else
printFrac(coeff)
end
end
else
if equalFrac(coeff, makeFrac(1, 1)) then
io.write(" + ")
elseif equalFrac(coeff, makeFrac(-1, 1)) then
io.write(" - ")
elseif lessFrac(makeFrac(0, 1), coeff) then
io.write(" + ")
printFrac(coeff)
else
io.write(" - ")
printFrac(negateFrac(coeff))
end
end
local pwr = p + 1 - j
if pwr>1 then
io.write("n^"..pwr)
else
io.write("n")
end
end
end
print()
return nil
end
-- main
for i=0,9 do
faulhaber(i)
end</lang>
{{out}}
<pre>0 : n
1 : 1/2n^2 + 1/2n
2 : 1/3n^3 + 1/2n^2 + 1/6n
3 : 1/4n^4 + 1/2n^3 + 1/4n^2
4 : 1/5n^5 + 1/2n^4 + 1/3n^3 - 1/30n
5 : 1/6n^6 + 1/2n^5 + 5/12n^4 - 1/12n^2
6 : 1/7n^7 + 1/2n^6 + 1/2n^5 - 1/6n^3 + 1/42n
7 : 1/8n^8 + 1/2n^7 + 7/12n^6 - 7/24n^4 + 1/12n^2
8 : 1/9n^9 + 1/2n^8 + 2/3n^7 - 7/15n^5 + 2/9n^3 - 1/30n
9 : 1/10n^10 + 1/2n^9 + 3/4n^8 - 7/10n^6 + 1/2n^4 - 3/20n^2</pre>
=={{header|Maxima}}==
| |||