LeetCode daily problem P3539: Sum of array products of magic sequences

---

Timeline

Timeline

2025-10-12

init

A combinatorial counting problem that uses DP to simulate the carry rule of binary addition; advanced dynamic programming.__builtin_popcount

---

Problem:

P3539 Sum of array products of magic sequences

https://leetcode.cn/problems/find-sum-of-array-product-of-magical-sequences/?envType=daily-question&envId=2025-10-12

Try exhaustive enumeration

At first I only thought of brute force, but I didn’t even enumerate completely. It was frustrating.

#include <algorithm>
#include <stack>
#include <vector>

using std::stack;
using std::vector;

const int MOD = 1000000007; // 10^9 + 7

class Solution {
    public:
        bool magicalSeq(vector<int> &seq, int m, int k, int length)
        {
                // Magic sequence: the length of the sequence seq is m.
                // 0 <= seq[i] < nums.length
                // The binary form of 2^seq[0] + 2^seq[1] + ... + 2^seq[m - 1] has m set bits.
                // That is, if two numbers in {seq[0], seq[1] .. seq[M-1]} are the same
                // Recursively merge identical numbers into one with the original value + 1, until there are no identical numbers left.
                int val;
                int set = 0;
                if (seq.size() != m) {
                        return false;
                }

                std::sort(seq.begin(), seq.end());

                if (seq[0] < 0 || seq[m - 1] >= length) {
                        return false;
                }

                stack<int, std::vector<int> > st(seq);

                while (!st.empty()) {
                        val = st.top();
                        st.pop();
                        set++;
                        if (!st.empty() && val == st.top()) {
                                st.top()++;
                                set--;
                        }
                }

                if (set == k) {
                        return true;
                } else {
                        return false;
                }
        }

        int magicalSum(int m, int k, vector<int> &nums)
        {
                int i, j;
                int n = nums.size();
                unsigned long magical_sum = 0;
                unsigned long array_mul;

                vector<vector<int> > result;
                int res_size;

                if (m <= n) {
                        // Generate all subsequences of nums of length m.
                        // predict array, similar to SVE instructions.
                        vector<int> predict(n, 0);
                        // Start from the lexicographically smallest one.
                        std::fill(predict.end() - m, predict.end(), 1);

                        do {
                                vector<int> tmp;
                                for (i = 0; i < n; i++) {
                                        if (predict[i] != 0) {
                                                tmp.push_back(i);
                                        }
                                }

                                result.push_back(tmp);
                        } while (std::next_permutation(predict.begin(), predict.end()));
                }

                res_size = result.size();
                // Determine whether it is a magic sequence.
                for (i = 0; i < res_size; i++) {
                        if (magicalSeq(result[i], m, k, n)) {
                                // If it is a magic sequence, compute the array product of all its permutations.
                                // But if the numbers are the same, different permutations have the same array product.
                                // Therefore, directly multiply by the factorial of m.
                                array_mul = 1;
                                for (j = 0; j < m; j++) {
                                        array_mul = (array_mul * nums[result[i][j]]) % MOD;
                                }
                                for (j = 1; j <= m; j++) {
                                        array_mul = (array_mul * j) % MOD;
                                }
                                magical_sum = (magical_sum + array_mul) % MOD;
                        }
                }
                return magical_sum;
        }
};

Dynamic Programming

The solution to this problem still uses dynamic programming, but it is advanced dynamic programming.

Problem

You are given two integersMandK, and an integer arraynums。

An integer sequenceseqIf the following conditions are satisfied, it is called magic sequence:

• seqthe sequence length ofM。
• 0 <= seq[i] < nums.length
• 2^seq[0] + 2^seq[1] + ... + 2^seq[M - 1]of binary form haveKpiece set。

of this sequence array product defined asprod(seq) = (nums[seq[0]] * nums[seq[1]] * ... * nums[seq[M - 1]])。

Return all valid magic sequential array product of sum 。

Since the answer may be very large, return the result modulo10^9 + 7modulo。

set Refers to a bit whose value is 1 in the binary representation of a number.

What we need to do is:

• Given arraynums, lengthn
• take the length asMthe sequenceseq(You can repeatedly choose elements from nums.)
• Condition:2^seq[0] + ... + 2^seq[M-1]In binary, there isKItem 1
• Array product:prod(seq) = nums[seq[0]] * ... * nums[seq[M-1]]
• Goal: the sum of the array products of all magic sequences (mod 1e9+7)

Number of sequence permutations

Suppose wenumsindices0..n-1take numbers:

• each numberitakenr_ipiece
• The total satisfies:r_0 + r_1 + ... + r_{n-1} = M

Number of permutations of length M:

$$\frac{M!}{r_0! \, r_1! \, ... \, r_{n-1}!}$$

Formula for the number of permutations of length M

Number of permutations without considering repeated numbers

• If all numbers are distinct, then the number of permutations of a sequence of length M is M!
  • Because M distinct elements can have M! permutations

Considering repeated numbers

• If some numbers appear repeatedly, for example:

$$seq = [a, a, b, b, b, c]$$

• Swapping repeated numbers does not produce a new sequence
• In combinatorics, we need to divide by the factorial of each repeated number, eliminating the effect of repeated permutations
• Suppose number i appearsr_itimes, the formula for the total number of permutations is:

$$\text{排列数} = \frac{M!}{r_0! \, r_1! \, \dots \, r_{n-1}!}$$

This is the classicmultiset permutation formula: number of permutations of length M = factorial of total elements / factorial of each repeated element

array product

$$\prod_{i=0}^{n-1} nums[i]^{r_i} \cdot \frac{M!}{r_0! \, r_1! \, ... \, r_{n-1}!}$$

Array product formula

The array product of a sequence is defined as:

$$\text{prod(seq)} = nums[seq[0]] \cdot nums[seq[1]] \cdot \dots \cdot nums[seq[M-1]]$$

• If number i appears in the sequencer_itimes, then the product can be written as:

$$\prod_{i=0}^{n-1} nums[i]^{r_i}$$

• Combining with the number of permutations, a certain sequence { r0, r_1, r_2, … ,r } the total contribution of the sequence is:

$$\text{总贡献} = \frac{M!}{r_0! \, r_1! \, ... \, r_{n-1}!} \cdot \prod_{i=0}^{n-1} nums[i]^{r_i}$$

• Explanation:
  • The array product corresponding to each sequence of length M is∏ nums[i]^r_i
  • There areM! / ∏ r_i!kinds of permutations (because permutations of repeated numbers are not counted as new sequences)
  • Therefore multiplying them gives the total contribution of that combination

Dynamic programming approach

Use DP to simulate the carry rule of binary addition

• Enumerating all possibler_i(the occurrence count of each number) is very complex
• The idea of dynamic programming:Consider the numbers in nums one by one

State definition

• f[i][j][p][q]represents:
  • Already considered the firsticount (nums[0…i])
  • total takenj
  • the sum of the low i bits of the current mask (carry state) isp
  • the number of set bits in the low i bits isq
  • value = corresponding to all sequences∏ r_t! * nums[t]^r_tthe sum, i.e., the total contribution of all these cases

transition formula

When considering thei+1number:

• Suppose we selectrpiecenums[i+1]
• New state:

$$f[i+1][j+r][\lfloor p/2 \rfloor + r][q + (p \bmod 2)] += f[i][j][p][q] * nums[i+1]^r * r!$$

• p/2andp%2the role of:
  • Simulate splitting the binary of mask from low to high bits
  • p % 2= the number of set bits in the low bits
  • p / 2= the mask after shifting left by one

Initialize

• fori=0, i.e., the first number:

$$f[0][j][j][0] = nums[0]^j * j!$$

• the number of set bits in the low 0 bits is 0

---

Result

• total set bits of mask = set bits of high mask + low q

if (__builtin_popcount(p) + q == K)
    res += f[n-1][M][p][q] * M!

Code implementation

The code uses many techniques; here we first introduce the following

Fast exponentiation

Core idea

Using the exponent’sbinary representationandThe splitting property of powers：

$$x^y = x^{(b_0·2^0 + b_1·2^1 + ... + b_k·2^k)} = \prod_{i | b_i=1} x^{2^i}$$

In other words:

• holdyView it as binary;
• Each time you square the base, it is equivalent to doubling the exponent;
• If the current binary bit is 1, multiply this term into the result.

Example

Suppose we want to compute:

$$3^{13}$$

13 in binary is:1101₂ = 8 + 4 + 0 + 1

So:

$$3^{13} = 3^{8} × 3^{4} × 3^{1}$$

We can obtain these powers by repeatedly squaring:

Step	Power	Result
3¹	3
3²	9
3⁴	81
3⁸	6561

So:

$$3^{13} = 3⁸ × 3⁴ × 3¹ = 6561 × 81 × 3 = 1594323$$

Corresponding algorithm logic

long long quickPow(long long x, long long y) {
    long long res = 1;
    long long cur = x;  // Current base

    while (y > 0) {
        if (y & 1) {
            res *= cur;  // Current bit is 1 → multiply the current power in
        }
        y >>= 1;         // Right shift by one bit is equivalent to dividing by 2.
        cur *= cur;      // Square the base → double the exponent
    }
    return res;
}

Fermat’s Little Theorem

According toFermat’s Little Theorem：

$$a^{p-1} ≡ 1 \pmod{p}$$

So:

$$a^{p-2} ≡ a^{-1} \pmod{p}$$

Fermat’s Little Theorem requires p to be a prime number.

In other words:

The modular inverse of a number = this number raised to the (mod - 2) power, modulo mod.

Final code implementation

#include <vector>
using std::vector;

class Solution {
    public:
        long long quickmul(long long x, long long y, long long mod)
        {
                long long res = 1, cur = x % mod;
                while (y) {
                        if (y & 1) {
                                res = res * cur % mod;
                        }
                        y >>= 1;
                        cur = cur * cur % mod;
                }
                return res;
        };

        int magicalSum(int m, int k, vector<int> &nums)
        {
                int n = nums.size();
                const long long mod = 1e9 + 7;

                // Compute the factorials from 1 to m
                vector<long long> fac(m + 1, 1);
                for (int i = 1; i <= m; i++) {
                        fac[i] = fac[i - 1] * i % mod;
                }
                // ifac[i]=(i!)^(−1) mod mod
                vector<long long> ifac(m + 1, 1);
                for (int i = 2; i <= m; i++) {
                        // Use Fermat's Little Theorem to compute the modular inverse of each i
                        ifac[i] = quickmul(i, mod - 2, mod);
                }
                // Use the recurrence relation to convert the inverse of a single integer into factorial inverses.
                for (int i = 2; i <= m; i++) {
                        ifac[i] = ifac[i - 1] * ifac[i] % mod;
                }

                // numsPower[i][j] represents nums[i] to the j-th power modulo mod.
                // If i appears j times in the sequence, then the product of the sequence is, by iterating over all i in the sequence,
                // Multiply numsPower[i][j] together.
                vector numsPower(n, vector<long long>(m + 1, 1));
                for (int i = 0; i < n; i++) {
                        for (int j = 1; j <= m; j++) {
                                numsPower[i][j] = numsPower[i][j - 1] * nums[i] % mod;
                        }
                }

                vector f(n, vector(m + 1, vector(m * 2 + 1, vector<long long>(k + 1, 0))));

                for (int j = 0; j <= m; j++) {
                        f[0][j][j][0] = numsPower[0][j] * ifac[j] % mod;
                }

                for (int i = 0; i + 1 < n; i++) {
                        for (int j = 0; j <= m; j++) {
                                for (int p = 0; p <= m * 2; p++) {
                                        for (int q = 0; q <= k; q++) {
                                                int q2 = p % 2 + q;
                                                if (q2 > k) {
                                                        break;
                                                }
                                                for (int r = 0; r + j <= m; r++) {
                                                        int p2 = p / 2 + r;
                                                        f[i + 1][j + r][p2][q2] +=
                                                                f[i][j][p][q] *
                                                                numsPower[i + 1][r] % mod *
                                                                ifac[r] % mod;
                                                        f[i + 1][j + r][p2][q2] %= mod;
                                                }
                                        }
                                }
                        }
                }
                long long res = 0;
                for (int p = 0; p <= m * 2; p++) {
                        for (int q = 0; q <= k; q++) {
                                if (__builtin_popcount(p) + q == k) {
                                        res = (res + f[n - 1][m][p][q] * fac[m] % mod) % mod;
                                }
                        }
                }
                return res;
        }
};