How can I compute the $n$-th complete Bell polynomial?

Question

I'm interested in computing the n-th complete Bell polynomial $ B_n(x_1,..., x_n) $ using the formula given as the last equation in the "Exponential Bell polynomials" section here (I tried posting the equation here but the math code would keep getting marked as an "error" in formatting, apologies), which is written as a sum in terms of the incomplete Bell polynomial $B_{n,k}(x_1, ..., x_{n-k+1})$.

Looking at the Mathematica documentation, I found that the incomplete Bell polynomial $B_{n,k}$ is implemented as

BellY[n,k,{x₁,…,xₙ₋ₖ₊₁}]

So I was thinking that I could compute the complete Bell polynomial using the equation from the article. The problem is that each summand has a different number of inputs, and I don't know how to implement a sum with this condition. Moreover, I would ideally like the polynomials to be outputted with variables (not evaluated at any specific values), and I also don't know how to do this.

Any and all ideas or suggestions are very much appreciated. Thanks for reading!

Answer 1

Clear["Global`*"]

Format[x[n_]] := Subscript[x, n]

bell[0] := 1

bell[n_Integer?Positive, var_Symbol : x] :=
 Sum[BellY[n, k, Array[var, n - k + 1]], {k, 1, n}]

Column[bell /@ Range[0, 7], Frame -> All] // TraditionalForm

Example:

bell[4] /. Thread[Array[x, 4] -> {1, 7, 5, 3}]

(* 213 *)

or

Evaluate[Array[y, 4]] = {1, 7, 5, 3};

bell[4, y]

(* 213 *)

EDIT: To also accept a list as an argument, extend the definition of bell with

bell[{}] = 1;

bell[lst_List] := Module[{n = Length@lst},
  Sum[BellY[n, k, lst[[;; n - k + 1]]], {k, 1, n}]]

Comparing with the original form

bell[Array[x, #]] === bell[#] & /@ Range[7]

(* {True, True, True, True, True, True, True} *)

bell[{1, 7, 5, 3}]

(* 213 *)