As in the main editorial, consider the related problem with only operations of type $$$G$$$. Applying several operations of type $$$G$$$ to an array is equivalent to applying some permutation to the codomain $$$\{1 \ldots k\}$$$ of the array. Since these operations are invertible, and we need to count equivalence classes with respect to this operation, it makes sense to try Burnside's lemma. That gives an explicit formula:

For a given permutation $$$\sigma \in S_k$$$, we have that $$$a = \sigma \circ a$$$ if and only if $$$a$$$ maps $$$\{1 \ldots n\}$$$ into the set of fixed points of $$$\sigma$$$. So the formula can be simplified to

Each permutation on $$$\{1 \ldots k\}$$$ with exactly $$$i$$$ fixed points is uniquely determined by the set of its fixed points, and by the derangement (no-fixpoint-permutation) it induces on the other $$$k - i$$$ points in its domain. So, the formula further reduces to

Counting derangements is easy and well-known, so this formula allows direct calculation of $$$A_i$$$ in $$$O(k \log{i})$$$. But since it is impossible to use more than $$$i$$$ different values in an array of length $$$i$$$, this can be reduced to $$$O(i \log{i})$$$ by just replacing $$$k$$$ with $$$\min(i, k)$$$.

The solution to the original problem in terms of $$$A_i$$$ is the same as described in the main editorial.

Implementation: 147499995

The editorial describes how to solve the problem in terms of $$$\sum_{j=0}^k S(i, j)$$$ for various values of $$$i$$$. This is another way to compute only that quantity.

Let $$$g_j(t) = \sum_{i=0}^\infty S(i, j) \frac{t^i}{i!}$$$ as a formal power series. It is the EGF of the $$$j$$$-th column of the Stirling numbers. Since multiplication of EGFs corresponds to partitioning, it can be seen that $$$g_j(t) = \frac{(g_1(t))^j}{j!}$$$: The $$$j!$$$ in the denominator forgets the order in which the sets were used in the partition. It's also easy to see that $$$g_1(t) = e^t - 1$$$.

To solve the problem, we need to calculate $$$\sum_{j=0}^k S(i, j)$$$ for various values of $$$i$$$. So, we can consider the function $$$h(t) = \sum_{j=0}^k g_j(t) = \sum_{j=0}^k \frac{(g_1(t))^j}{j!}$$$ and calculate its coefficients. There is obviously not enough time to evaluate all terms of this sum; we will have to cancel most of them out and then undo this. One way to achieve this by noticing that $$$h(t)$$$ is almost the same as $$$\exp{(g_1(t))}$$$. So, $$$h$$$ satisfies $$$h'(t) \approx h(t) \cdot g_1'(t)$$$. We can calculate the exact value easily and then solve the resulting ODE as follows:

This formula can be directly implemented in $$$O(n \log n)$$$, at least if you have access to a library/template for the necessary formal power series/generating function operations. Implementation: 147380982