$q$-multinomial series

Question

I have the following code, which produces a $q$-multinomial coefficient, but selected randomly according to a Poisson distribution.

Consider a 3D lattice path selected uniformly at random from a cuboidal domain with random dimensions $m1,m2,m3$, then average over all dimensions, where each value has mean $l1,l2,l3$, respectively:

qmultinomial[q_, m1_, m2_, m3_] := 
 Normal@Series[QFactorial[m1 + m2 + m3, q]/(
   QFactorial[m1, q] QFactorial[m2, q] QFactorial[m3, q]), {q, 0, 
    m1 m2 + m2 m3}]
l1 = 3;
l2 = 3;
l3 = 3;
maxi = 10;
maxj = 10;
maxk = 10;
poissonqmult = 
   N@Series[
     Total@Flatten@
       Table[1/Multinomial[i, j, k] (Exp[-l1] l1^i)/i! (
         Exp[-l2] l1^j)/j! (Exp[-l3] l1^k)/
         k! qmultinomial[q, i, j, k], {i, 0, maxi}, {j, 0, maxj}, {k, 
         0, maxk}], {q, 0, 64}]; // AbsoluteTiming
ListPlot@Transpose[{Table[
    i, {i, 1, Length@CoefficientList[poissonqmult, q]}], 
   CoefficientList[poissonqmult, q]}]

which gives

{32.9553, Null}

The problem is the code is slow, and I need to increase $maxi,maxj,maxk$ to finite values which are quite large (technically infinite). Is there some way I can speed it up? Perhaps using ParallelTable?

Answer 1

Drop the series expansion altogether and define the qmultinomial directly:

ClearAll[qmultinomial]
qmultinomial[q_, m1_, m2_, m3_] := 
 QFactorial[m1 + m2 + m3, q]/(QFactorial[m1, q] QFactorial[m2, q] QFactorial[m3, q])

Then retain the rest of the code as written and you will obtain ostensibly the same result, but in 1.5 s (your original code took ca. vs. 20 s on my machine):

l1 = 3;
l2 = 3;
l3 = 3;
maxi = 10;
maxj = 10;
maxk = 10;
poissonqmult = 
   N@Series[
     Total@Flatten@
       Table[1/Multinomial[i, j, k] (Exp[-l1] l1^i)/
          i! (Exp[-l2] l1^j)/j! (Exp[-l3] l1^k)/k! qmultinomial[q, i, 
          j, k], {i, 0, maxi}, {j, 0, maxj}, {k, 0, maxk}], {q, 0, 
      64}]; // AbsoluteTiming

ListPlot@CoefficientList[poissonqmult, q]

(* Out: {1.58515, Null} *)