Let's try to figure out what's the value of $(x_{1},x_{2},x_{3},...)$.
We can discover that the value is the number of different paths from $(0,0,0,...)$ to $(x_{1},x_{2},x_{3},...)$, one step in the path is defined as increase one of the co-ordinates by one.
So we can know the value of $(x_{1},x_{2},x_{3},...)$ is $\frac{(\sum_{k=1}^{N} x_{k})!}{\prod_{k=1}^{N} (x_{k}!)}$
If you don't know why, image that there are $x_{1}$ $vectors$ of $(1,0,0,...)$, $x_{2}$ $vectors$ of $(0,1,0,...)$, and so on. Therefore, no matter how you permutate them, they always form a path from $(0,0,0,...)$ to $(x_{1},x_{2},x_{3},...)$, the number of such permutations is $\frac{(\sum_{k=1}^{N} x_{k})!}{\prod_{k=1}^{N} (x_{k}!)}$ (Let's call it $M$ later)
The problem is, in which condition $M$ is not divisible by P ?
We can just take only $P$'s multiples into consideration.
The number of factor $P$ in $n!$ is $\left[\frac{n}{P}\right]+\left[\frac{n}{P^{2}}\right]+\left[\frac{n}{P^{3}}\right]+...$
For simplicity, we call the value $f(n)$ (regard $P$ as constant)
Then M is not divisible by P if and only if $f(\sum_{k=1}^{N} x_{k})=\sum_{k=1}^{N} f(x_{k})$
i.e., for every $r\in\mathbb{N}$, $\sum_{k=1}^{N} (x_{k}\mod P^{r})<P^{r}$ must hold.
If we transform every $n$ into $P-dicimal$ numbers, call it $T(n)$, then the above equation holds if and only if $\sum_{k=1}^{N} (k_{th}$ digit of $T(x_{k}))<P$.
So be question becomes, how many set, $(x_{1},x_{2},x_{3},...)$, are there, such that for each $i\in\mathbb{N}$, $\sum_{k=1}^{N} (i_{th}$ digit of $T(x_{k}))<P$ ?
Let's try to design a Dynamic Programming method to get the answer.
Suppose that we can enumerate from high digit to low digit.
A state includes the following information:
1. $i$, represent we are considering the $i_{th}$ digit now.
2. For every $T(x_{k})$, the range of its $(i-1)_{th}$ digit should be.(whether the $(i-1)_{th}$ digit has low bound or up bound)
For example, if the $i_{th}$ digit has a low bound of $a$, a up bound of $b$, the $(i-1)_{th}$ digit has no bounds if the $i_{th}$ digit$\in[a+1,b-1]$, has a low bound if the $i_{th}$ digit is $a$, has a up bound if the $i_{th}$ digit is $b$
Thus, we can update the values by enumerate the subsets of its set of bounds.
It's convenient to use bitmask to represent the sets.
Partial time complexity: $O(3^{2N}\log_{P}\max(x_{k}))$
However, we can optimize it by skipping if the value is 0.
Use another DP to calculate the number of cases, such that $\sum_{k=1}^{N} \{ i_{th}$ digit of $x_{k}$, within the given bound$\}<P$.
Partial time complexity: $O(NP)$.
Total time complexity: $O(3^{2N}NP\log_{P}\max(x_{k}))$.
code: