Overview

Particles diffuse across \(N\) containers, and the amount of “positive charge” and “negative charge” particles in each container is tracked. This problem asks for the amount of positive charge particles in container \(G\) after \(T\) steps. We can model this as repeated linear transformations and solve it efficiently using matrix exponentiation.

Discussion

States Can Be Represented by Linear Transformations

The key point of this problem is that one step of the update can be represented as a linear transformation (matrix-vector multiplication).

Each container \(i\) has two types of values: “positive charge amount \(x_i^+\)” and “negative charge amount \(x_i^-\)”. We can represent the state as a vector of length \(2N\) by listing all of these values.

Looking at the update equations:

\[y_j^{+} = \sum_{i} \Bigl( (1 - F_{i,j}) x_i^{+} P_{i,j} + F_{i,j} x_i^{-} P_{i,j} \Bigr)\]

\[y_j^{-} = \sum_{i} \Bigl( (1 - F_{i,j}) x_i^{-} P_{i,j} + F_{i,j} x_i^{+} P_{i,j} \Bigr)\]

we can see that each component of \(y\) is a linear combination (sum of scaled terms) of the components of \(x\). In other words, using some \(2N \times 2N\) matrix \(A\), we can write:

\[\mathbf{y} = A \mathbf{x}\]

Repeating this \(T\) times yields:

\[\mathbf{x}_T = A^T \mathbf{x}_0\]

Issues with the Naive Approach

Updating step-by-step naively would take \(O(T \cdot N^2)\) time, which will definitely TLE since \(T \leq 10^{18}\). Therefore, we use matrix exponentiation (binary exponentiation) to find \(A^T\) in \(O(\log T)\) matrix multiplications.

Handling Rational Numbers

\(P_{i,j} = R_{i,j} / (\sum_{j'} R_{i,j'})\) involves division. The answer is a rational number, but the output must be modulo \(10^9+7\). Since \(10^9+7\) is a prime, we replace division with the modular multiplicative inverse (using Fermat’s Little Theorem: \(a^{-1} \equiv a^{p-2}\)). By precomputing the inverse of each row sum, we can perform all calculations using integer arithmetic modulo \(10^9+7\).

Algorithm

State Indexing

In the vector of length \(2N\), we define the indices as follows: - Index \(i\) (\(0 \leq i < N\)): positive charge amount in container \(i\) - Index \(N+i\): negative charge amount in container \(i\)

Constructing the Transition Matrix \(A\)

For each \((i, j)\), calculate \(P_{i,j} = R_{i,j} \cdot (\text{row sum})^{-1}\): - No inversion (\(F_{i,j}=0\)): - \(A[j][i] \mathrel{+}= P_{i,j}\) (positive \(\to\) positive) - \(A[N+j][N+i] \mathrel{+}= P_{i,j}\) (negative \(\to\) negative) - With inversion (\(F_{i,j}=1\)): - \(A[j][N+i] \mathrel{+}= P_{i,j}\) (negative \(\to\) positive) - \(A[N+j][i] \mathrel{+}= P_{i,j}\) (positive \(\to\) negative)

Here, \(A[\text{row}][\text{col}]\) represents the matrix in \(\mathbf{y} = A\mathbf{x}\).

Extracting the Answer

The initial vector \(\mathbf{x}_0\) has \(1\) for the “positive charge of container \(S\)” and \(0\) for all other components. Thus, \(\mathbf{x}_0\) is a unit vector where only index \(S\) is \(1\).

The positive charge amount in container \(G\) in \(\mathbf{x}_T = A^T \mathbf{x}_0\) is the component at index \(G\), which is:

\[\mathbf{x}_T[G] = (A^T)[G][S]\]

This is exactly the \((G, S)\) element of \(A^T\). We can obtain the answer simply by extracting this single element from the result of the matrix exponentiation.

Complexity

The size of the matrix is \(D = 2N \leq 60\).

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

Since \(N \leq 30\) and \(T \leq 10^{18}\), the number of operations is around \(60^3 \times 60 \approx 1.3 \times 10^7\), which is fast enough.

Key Points of Implementation

• Division using modular inverse: Precompute the inverse of the row sum \(\sum_{j} R_{i,j}\) as inv, and let \(P_{i,j} = R_{i,j} \cdot \text{invrow}[i]\). Maintain all values modulo \(10^9+7\).

• Skip \(P_{i,j}=0\): Skip transitions where \(R_{i,j}=0\) as they do not contribute (even if an inversion rule is set, it has no effect).

• Be careful with matrix orientation: When treating the state as “matrix \(\times\) column vector”, store the coefficients of the update equations directly into \(A[\text{destination}][\text{source}]\). The final answer is the \((G, S)\) element.

• Simultaneous updates: One step of the update is calculated using “only the values before the update”. Since the matrix-vector multiplication format automatically represents this simultaneous update, no special consideration is needed.

• Large \(T\): Since \(T \leq 10^{18}\), read it as long long and use binary exponentiation.

  Source Code

#include <bits/stdc++.h>
using namespace std;

const long long MOD = 1000000007;

long long pw(long long a, long long b){
    a%=MOD; if(a<0)a+=MOD;
    long long r=1;
    while(b>0){
        if(b&1) r=r*a%MOD;
        a=a*a%MOD;
        b>>=1;
    }
    return r;
}
long long inv(long long a){ return pw(a, MOD-2); }

int main(){
    int N, M; long long T;
    scanf("%d %d %lld", &N, &M, &T);
    int S, G;
    scanf("%d %d", &S, &G);
    S--; G--;
    vector<vector<long long>> R(N, vector<long long>(N));
    for(int i=0;i<N;i++)for(int j=0;j<N;j++) scanf("%lld",&R[i][j]);
    vector<vector<int>> F(N, vector<int>(N,0));
    for(int k=0;k<M;k++){
        int u,v; scanf("%d %d",&u,&v); u--;v--;
        F[u][v]=1;
    }
    // P[i][j] = R[i][j] * inv(rowsum_i)
    vector<long long> invrow(N);
    for(int i=0;i<N;i++){
        long long s=0;
        for(int j=0;j<N;j++) s+=R[i][j];
        invrow[i]=inv(s%MOD);
    }
    int D=2*N;
    // state: positive i -> i, negative i -> N+i
    vector<vector<long long>> A(D, vector<long long>(D,0));
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            long long P = (R[i][j]%MOD)*invrow[i]%MOD;
            if(P==0) continue;
            if(F[i][j]==0){
                // no flip
                A[j][i] = (A[j][i] + P)%MOD;      // y_j+ from x_i+
                A[N+j][N+i] = (A[N+j][N+i] + P)%MOD; // y_j- from x_i-
            } else {
                A[j][N+i] = (A[j][N+i] + P)%MOD;  // y_j+ from x_i-
                A[N+j][i] = (A[N+j][i] + P)%MOD;  // y_j- from x_i+
            }
        }
    }
    // matrix power A^T
    auto mul=[&](vector<vector<long long>>&X, vector<vector<long long>>&Y){
        vector<vector<long long>> Z(D, vector<long long>(D,0));
        for(int i=0;i<D;i++){
            for(int k=0;k<D;k++){
                if(X[i][k]==0) continue;
                long long v=X[i][k];
                for(int j=0;j<D;j++){
                    Z[i][j]=(Z[i][j]+v*Y[k][j])%MOD;
                }
            }
        }
        return Z;
    };
    vector<vector<long long>> result(D, vector<long long>(D,0));
    for(int i=0;i<D;i++) result[i][i]=1;
    long long e=T;
    vector<vector<long long>> base=A;
    while(e>0){
        if(e&1) result=mul(result,base);
        base=mul(base,base);
        e>>=1;
    }
    // initial vector: e[S]=1 (positive at S)
    // after T steps: y = result * e, y[G] = result[G][S]
    printf("%lld\n", result[G][S]%MOD);
    return 0;
}

---

This editorial was generated by claude4.8opus-high.