You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
I just wrote a function, based on a recurrence relation that can be found on Wikipedia.
Do you have any suggestion whether it is useful for this package, or else if this code could be a useful addition elsewhere?
Any (other) comment of course welcome too.
Cheers,
Evaluate the incomplete Bell polynomial B_{N,K}(x)
where x is an array containing the abscissa values [x_1, ... x{n-k+1}]
Definition and background information, see
https://en.wikipedia.org/wiki/Bell_polynomials#Definitions
https://en.wikipedia.org/wiki/Bell_polynomials#Properties
https://en.wikipedia.org/wiki/Bell_polynomials#Table_of_values
"""
evaluate_incomplete_Bell_poly = function(N,K,all_x)
# see formula C2
# let us construct the triangular table
if K > N
return 0.
end
all_B_old = zeros(N+1,1)
all_B_old[1] = 1.
all_B_new = zeros(N+1,1)
all_B_new[1] = 1.
for k=1:K
for n=1:N
all_B_new[n+1] = sum([binomial(n-1,i-1)*all_x[i]*all_B_old[n-i+1] for i=1:n-k+1])
end
all_B_old .= all_B_new
end
return all_B_new[end]
end
## test the function
x = randn(11,1)
@assert evaluate_incomplete_Bell_poly(0,0,x) ≈ 1.
@assert evaluate_incomplete_Bell_poly(1,0,x) ≈ 0.
@assert evaluate_incomplete_Bell_poly(3,0,x)≈ 0.
@assert evaluate_incomplete_Bell_poly(0,1,x) ≈ 0.
@assert evaluate_incomplete_Bell_poly(0,3,x) ≈ 0.
@assert evaluate_incomplete_Bell_poly(1,1,x) ≈ x[1]
@assert evaluate_incomplete_Bell_poly(2,2,x) ≈ x[1]^2
@assert evaluate_incomplete_Bell_poly(3,3,x) ≈ x[1]^3
@assert evaluate_incomplete_Bell_poly(6,4,x) ≈ 45*x[1]^2*x[2]^2 + 20*x[1]^3*x[3]
@assert evaluate_incomplete_Bell_poly(5,4,x) ≈ 10*x[1]^3*x[2]
@assert evaluate_incomplete_Bell_poly(6,3,x) ≈ 15*x[2]^3 + 60*x[1]*x[2]*x[3] + 15*x[1]^2*x[4]
# use some results in (Abbas, M.; Bouroubi, S. (2005). "On new identities for Bell's polynomial". Discrete Math. 293 (1–3)) as reference
# see https://vinar.vin.bg.ac.rs//bitstream/id/12791/4374.pdf
@assert evaluate_incomplete_Bell_poly(9,7,x) ≈ 378*x[1]^5*x[2]^2 + 84*x[1]^6*x[3]
@assert evaluate_incomplete_Bell_poly(10,7,x) ≈ 3150*x[1]^4*x[2]^3 + 2520*x[1]^5*x[2]*x[3] + 210*x[1]^6*x[4]
@assert evaluate_incomplete_Bell_poly(11,7,x) ≈ 17325*x[1]^3*x[2]^4 + 34650*x[1]^4*x[2]^2*x[3] + + 4620*x[1]^5*x[3]^2 + 6930*x[1]^5*x[2]*x[4] + 462*x[1]^6*x[5]
println("tests passed")
Dear contributors of the Combinatorics.jl package,
Would you please know where I could find an implementation of (the evaluation) of the incomplete Bell polynomials?
https://en.wikipedia.org/wiki/Bell_polynomials#Definitions#Exponential%20Bell%20polynomials
Thank you!
The text was updated successfully, but these errors were encountered: