If $$$\gcd(x!, p) \ne 1$$$ then there is no modular inverse of $$$x! \pmod{p^e}$$$.

I have an idea that I am not sure of (As I never tested it before), but you can try it.

Let the modulo of this problem be $$$M$$$. Initially, we will first have to prime factorize $$$M$$$, this can be done in $$$O(\sqrt{M})$$$, so we know all the primes which appear in the prime factorization of $$$M$$$

When you have an integer $$$x$$$ you can represent it as the product of two parts, one part is coprime with the $$$M$$$ and the other part consists as a product of primes from the prime factorization of $$$M$$$

Notice that the part non-coprime with $$$M$$$ is equal to $$$gcd(x,M)$$$, so the remaining part is equal to $$$\dfrac{x}{gcd(x,M)}$$$, let these two parts be $$$q$$$ and $$$p$$$ respectively, we are going to store $$$p$$$ as a normal integer

As for $$$q$$$, notice $$$q$$$ will only consist of product of primes from set of primes appearing in the prime factorization of $$$M$$$, there are at most $$$8$$$ such primes when $$$M≤10^6$$$, so we are actually going to store the prime factorization of $$$q$$$ as an array of exponents, where the $$$i-th$$$ element donates the exponent of $$$i-th$$$ prime in the prime factorization of $$$q$$$, the array size will not exceed $$$8$$$

I don't know how can you apply addition or subtraction to numbers when they are in this form, but you can apply multiplication and division like this :

let's suppose we have two numbers $$$x,y$$$ in the above form and we want to find their product $$$z$$$ modulo $$$M$$$ (also in the above form), here is what we are going to do

First for the $$$p$$$ part, we simple multiply the two parts together taken under modulo $$$M$$$, and you have the $$$p$$$ part for their product $$$z$$$

As for the $$$q$$$ part, we are just going to add the corresponding elements of the exponent arrays of $$$x,y$$$, which will result in the exponent array of their product

Division is similar to multiplication, if we want to calculate $$$\dfrac{x}{y}$$$, instead of multiplying the $$$p$$$ parts, we multiply the $$$p$$$ part of the $$$x$$$ by the inverse of the $$$p$$$ part of $$$y$$$, finding the inverse is possible by using the Extended Euclidean Algorithm

And as for the $$$q$$$ part, we subtract the corresponding elements of the exponent array of $$$y$$$ from the exponent array of $$$x$$$

If now you want to find the actual value of the number taken modulo $$$M$$$, you will evaluate the prime factorization using Fast Power (AKA Binary Exponentiation) and finally multiply the results of each prime by the $$$p$$$ part (All operations are under modulo $$$M$$$ here), and you will have your desired value

This might be helpful — https://cp-algorithms.com/combinatorics/binomial-coefficients.html#mod-prime-pow

No idea how to solve it, but I know that Japanese can: Library-Checker

For nck mod m. If m is small, I feel like k cannot be that far away from n, otherwise binomial would become zero.

Suppose $$$n! = p^e \cdot x(x\perp p)$$$. We want to find $$$e$$$ and $$$x \bmod p^\alpha$$$.

$$$e$$$ is just $$$\sum_{i\geq 1}\lfloor{n\over p^i}\rfloor$$$. Let $$$M = \prod_{1\leq i < p^\alpha,i\perp p}i$$$, then $$$x(n) \equiv M^{\lfloor{n\over p^\alpha}\rfloor}x(\lfloor{n\over p}\rfloor)x(n \bmod p^\alpha)$$$, which can be calculated recursively.

Hey, I wrote a blog recently just about that! Essentially, you just compute a pair $$$(k, a)$$$ such that $$$n! = p^k a$$$, where $$$\gcd(a, p) = 1$$$. Then you can cancel out the powers of $$$p$$$, and for co-prime parts you just compute inverses as usual.