pmultinom is a library for calculating multinomial probabilities. The probabilities that can be calculated include the multinomial cumulative distribution function:
P(N₁ ≤ u₁, N₂ ≤ u₂, ⋯, N_k ≤ u_k)
In this case the usage would be

pmultinom(upper=us, size=n, probs=ps, method="exact")

where us is the vector containing u₁, u₂, ⋯, u_k, and n and ps are the parameters of the multinomial distribution. This usage is analogous to the use of pbinom. Another important case is the probability of seeing more than some minimum number of observations in each category:
P(N₁ > l₁, N₂ > l₂, ⋯, N_k > l_k)
In this case the usage would be

pmultinom(lower=ls, size=n, probs=ps, method="exact")

where this time ls is the vector containing l₁, l₂, ⋯, l_k. Notice that in this case these are greater than signs, not greater than or equal signs. This is analogous to the usage of pbinom with lower.tail=FALSE. With some creativity, these can be adapted to calculate the probability that the maximum or minimum of a multinomial random vector is a given number, or that a given category will be the most or least observed. pmultinom also supports a more general usage, in which both lower and upper bounds are specified:
P(l₁ < N₁ ≤ u₁, l₂ < N₂ ≤ u₂, ⋯, l_k < N_k ≤ u_k)
In this case the usage would be

pmultinom(lower=ls, upper=us, size=n, probs=ps, method="exact")

See vignette("pmultinom") for the above text in Latex plus an example application. Many thanks to Aislyn Schalck for advice and encouragement.