Overview

This problem asks us to simulate the diffusion of particles among \(N\) connected vessels and find the amount of positively charged particles in a specific vessel \(G\) after \(T\) steps. Since the number of steps \(T\) is very large (up to \(10^{18}\)), a naive simulation is too slow. Instead, we speed up the process using matrix exponentiation (binary exponentiation).

Analysis

1. Limitations of Naive Simulation

The transition of particles in each step is represented by a linear combination (addition and scalar multiplication) of the current amount of particles (positive and negative) in each vessel. However, since the number of steps \(T\) is up to \(10^{18}\), a naive simulation that updates the state step-by-step will not run within the time limit (resulting in TLE).

2. Linearity and State Expansion

For each vessel \(i\), we need to manage two states: the amount of positive particles \(x_i^{+}\) and the amount of negative particles \(x_i^{-}\). Since there are \(N\) vessels, there are \(2N\) states in total.

We can represent these states as a single column vector \(\mathbf{x}\): $\(\mathbf{x} = \begin{pmatrix} x_1^{+} \\ \vdots \\ x_N^{+} \\ x_1^{-} \\ \vdots \\ x_N^{-} \end{pmatrix}\)$

A one-step transition can be expressed as a linear mapping on this vector \(\mathbf{x}\), which is a matrix-vector product \(\mathbf{y} = A \mathbf{x}\) using a \((2N) \times (2N)\) transition matrix \(A\).

3. Speedup using Matrix Exponentiation

If a one-step transition is represented as \(\mathbf{x}^{(1)} = A \mathbf{x}^{(0)}\), then the state after \(T\) steps is as follows: $\(\mathbf{x}^{(T)} = A^T \mathbf{x}^{(0)}\)$

The \(T\)-th power of the matrix (\(A^T\)) can be computed in \(O((2N)^3 \log T)\) time complexity using binary exponentiation. This allows for fast computation even when \(T = 10^{18}\).

Algorithm

1. Calculation of Transition Probabilities and Modular Inverse

The transition probability \(P_{i,j}\) from vessel \(i\) to vessel \(j\) is defined as follows: $\(P_{i,j} = \frac{R_{i,j}}{\sum_{j'=1}^{N} R_{i,j'}}\)$

Since we perform calculations on rational numbers modulo \(10^9+7\), we compute the modular inverse of the denominator \(D_i = \sum_{j'=1}^{N} R_{i,j'}\) as \(D_i^{-1} \pmod{10^9+7}\) using Fermat’s Little Theorem. Thus, the transition probability can be expressed as a multiplication: $\(P_{i,j} \equiv R_{i,j} \times D_i^{-1} \pmod{10^9+7}\)$

2. Construction of Transition Matrix \(A\)

We construct a transition matrix \(A\) of size \((2N) \times (2N)\). We map indices \(0\) to \(N-1\) to “positive charges” and indices \(N\) to \(2N-1\) to “negative charges”.

For a transition from vessel \(i\) to vessel \(j\): * If there is no inversion rule (\(F_{i,j} = 0\)): * Positive particles remain positive: \(A[j][i] \leftarrow A[j][i] + P_{i,j}\) * Negative particles remain negative: \(A[j+N][i+N] \leftarrow A[j+N][i+N] + P_{i,j}\) * If there is an inversion rule (\(F_{i,j} = 1\)): * Positive particles invert to negative: \(A[j+N][i] \leftarrow A[j+N][i] + P_{i,j}\) * Negative particles invert to positive: \(A[j][i+N] \leftarrow A[j][i+N] + P_{i,j}\)

3. Matrix Exponentiation and Outputting the Answer

In the initial state \(\mathbf{x}^{(0)}\), vessel \(S\) has a positive particle amount of \(1\), and all other states are \(0\). Therefore, the amount of positive particles in vessel \(G\) after \(T\) steps corresponds to the \((G, S)\) entry of the cumulative transition matrix \(A^T\), which is \(A^T[G][S]\) (or \(A^T[G-1][S-1]\) in 0-indexed notation).

Complexity

• Time Complexity: \(O(N^3 \log T)\)

  • The size of the transition matrix is \(2N \times 2N\).
  • Each matrix multiplication takes \(O((2N)^3) = O(N^3)\) operations.
  • Since we perform matrix multiplication \(O(\log T)\) times using binary exponentiation, the overall time complexity is \(O(N^3 \log T)\).
  • When \(N \le 30\) and \(T \le 10^{18}\), the number of operations is approximately \(60^3 \times 60 \approx 1.3 \times 10^7\), which easily runs within the time limit.

• Space Complexity: \(O(N^2)\)

  • The memory required to store a \(2N \times 2N\) matrix is \(O(N^2)\). Since \(N \le 30\), the memory usage is extremely small.

Implementation Points

1. 0-indexed Conversion: The input vessel indices \(S, G\) and inversion rules \(u_k, v_k\) are 1-based, so we subtract \(1\) to convert them to 0-based in our implementation.

2. Optimizing Matrix Multiplication (Pruning): When computing the matrix product \(C_{i,j} = \sum_k A_{i,k} B_{k,j}\), skipping the calculation (using continue) when \(A_{i,k} == 0\) in the inner loop provides a constant-factor speedup.

3. Modular Inverse: Since \(10^9+7\) is prime, modInverse can be computed in \(O(\log(\text{MOD}))\) time using binary exponentiation.

   Source Code

#include <iostream>
#include <vector>

using namespace std;

const int MOD = 1000000007;

// Function to compute (base^exp) % MOD
long long power(long long base, long long exp) {
    long long res = 1;
    base %= MOD;
    while (exp > 0) {
        if (exp % 2 == 1) res = (res * base) % MOD;
        base = (base * base) % MOD;
        exp /= 2;
    }
    return res;
}

// Function to compute modular inverse
long long modInverse(long long n) {
    return power(n, MOD - 2);
}

typedef vector<vector<long long>> Matrix;

// Matrix multiplication under modulo MOD
Matrix multiply(const Matrix& A, const Matrix& B) {
    int n = A.size();
    Matrix C(n, vector<long long>(n, 0));
    for (int i = 0; i < n; ++i) {
        for (int k = 0; k < n; ++k) {
            if (A[i][k] == 0) continue;
            for (int j = 0; j < n; ++j) {
                C[i][j] = (C[i][j] + A[i][k] * B[k][j]) % MOD;
            }
        }
    }
    return C;
}

// Binary exponentiation for Matrix
Matrix matrixPower(Matrix A, long long p) {
    int n = A.size();
    Matrix res(n, vector<long long>(n, 0));
    for (int i = 0; i < n; ++i) res[i][i] = 1;
    while (p > 0) {
        if (p & 1) res = multiply(res, A);
        A = multiply(A, A);
        p >>= 1;
    }
    return res;
}

int main() {
    // Optimize standard I/O operations for performance
    ios_base::sync_with_stdio(false);
    cin.tie(NULL);

    int N, M;
    long long T;
    if (!(cin >> N >> M >> T)) return 0;

    int S, G;
    cin >> S >> G;
    --S; --G; // Convert to 0-based indexing

    vector<vector<long long>> R(N, vector<long long>(N));
    vector<long long> D(N, 0);
    for (int i = 0; i < N; ++i) {
        for (int j = 0; j < N; ++j) {
            cin >> R[i][j];
            D[i] = (D[i] + R[i][j]) % MOD;
        }
    }

    vector<vector<int>> F(N, vector<int>(N, 0));
    for (int k = 0; k < M; ++k) {
        int u, v;
        cin >> u >> v;
        --u; --v; // Convert to 0-based indexing
        F[u][v] = 1;
    }

    vector<long long> D_inv(N);
    for (int i = 0; i < N; ++i) {
        D_inv[i] = modInverse(D[i]);
    }

    int sz = 2 * N;
    Matrix A(sz, vector<long long>(sz, 0));

    // Construct the transition matrix A
    for (int i = 0; i < N; ++i) {
        for (int j = 0; j < N; ++j) {
            long long P = (R[i][j] * D_inv[i]) % MOD;
            if (P == 0) continue;
            int f = F[i][j];
            if (f == 0) {
                // No inversion: positive to positive, negative to negative
                A[j][i] = (A[j][i] + P) % MOD;
                A[j + N][i + N] = (A[j + N][i + N] + P) % MOD;
            } else {
                // Inversion: positive to negative, negative to positive
                A[j][i + N] = (A[j][i + N] + P) % MOD;
                A[j + N][i] = (A[j + N][i] + P) % MOD;
            }
        }
    }

    // Compute A^T
    Matrix AT = matrixPower(A, T);

    // The answer is the positive charge in container G, starting from container S
    long long ans = AT[G][S];
    cout << ans << "\n";

    return 0;
}

---

This editorial was generated by gemini-3.5-flash-high.