행렬 곱셈 및 고속 거듭제곱

1. 행렬 곱셈

1.1. 정의

두 행렬 \\(A, B\\)가 주어졌을 때, 행렬 \\(A\\)의 크기가 \\(n \times m\\)이고 행렬 \\(B\\)의 크기가 \\(m \times p\\)이면 이 두 행렬은 곱셈이 가능합니다. 이 연산의 결과로 생성되는 행렬 \\(C\\)의 크기는 \\(n \times p\\)가 됩니다. 예를 들어 다음과 같은 행렬 \\(A\\)와 \\(B\\)가 있습니다: \\( A=\\begin{bmatrix} a\_{11} & a\_{12} & a\_{13}\\\\ a\_{21} & a\_{22} & a\_{23} \\end{bmatrix}\\) \\( B=\\begin{bmatrix} b\_{11} & b\_{12}\\\\ b\_{21} & b\_{22}\\\\ b\_{31} & b\_{32} \\end{bmatrix} \\) 이때 행렬 \\(C = AB\\)는 아래와 같이 계산됩니다:
\[C=AB=\begin{bmatrix} a_{11} & a_{12} & a_{13}\\ a_{21} & a_{22} & a_{23} \end{bmatrix}\begin{bmatrix} b_{11} & b_{12}\\ b_{21} & b_{22}\\ b_{31} & b_{32} \end{bmatrix}=\begin{bmatrix} a_{11}b_{11}+a_{12}b_{21}+a_{13}b_{31} & a_{11}b_{12}+a_{12}b_{22}+a_{13}b_{32}\\ a_{21}b_{11}+a_{22}b_{21}+a_{23}b_{31} & a_{21}b_{12}+a_{22}b_{22}+a_{23}b_{32} \end{bmatrix} \]

1.2. 주요 특징

  • 행렬 곱셈이 정의되기 위해서는 첫 번째 행렬의 열 개수가 두 번째 행렬의 행 개수와 일치해야 합니다.
  • 행렬 곱셈은 일반적으로 교환법칙이 성립하지 않습니다 (즉, 대부분의 경우 \\(AB \neq BA\\)입니다).

1.3. 코드 구현

아래는 C++에서 행렬 곱셈을 구현하기 위한 구조체 및 기본적인 입출력 코드입니다. 행렬의 원소는 `long long` 타입을 사용하여 큰 수를 처리할 수 있도록 하고, `std::vector`를 통해 동적으로 크기를 조절하는 행렬을 구현합니다.
#include <iostream>
#include <vector>
#include <stdexcept> // 예외 처리를 위해

// 행렬을 나타내는 구조체
struct Matrix {
    int rows, cols;
    std::vector<std::vector<long long>> elements;

    // 기본 생성자
    Matrix() : rows(0), cols(0) {}

    // 크기를 지정하여 초기화하는 생성자
    Matrix(int r, int c) : rows(r), cols(c), elements(r, std::vector<long long>(c, 0)) {}

    // 행렬 곱셈 연산자 오버로드
    Matrix operator*(const Matrix& other) const {
        // 곱셈 가능 여부 확인: 현재 행렬의 열 수와 다른 행렬의 행 수가 같아야 함
        if (cols != other.rows) {
            throw std::runtime_error("행렬의 차원이 곱셈에 적합하지 않습니다.");
        }

        Matrix result(rows, other.cols); // 결과 행렬 초기화 (현재 행의 수 x 다른 행렬의 열 수)

        for (int i = 0; i < rows; ++i) { // 결과 행렬의 행 인덱스
            for (int j = 0; j < other.cols; ++j) { // 결과 행렬의 열 인덱스
                for (int k = 0; k < cols; ++k) { // 내부 곱셈을 위한 인덱스 (현재 열 인덱스 == 다른 행렬의 행 인덱스)
                    result.elements[i][j] += elements[i][k] * other.elements[k][j];
                }
            }
        }
        return result;
    }

    // 행렬 내용을 출력하는 함수
    void display() const {
        for (int i = 0; i < rows; ++i) {
            for (int j = 0; j < cols; ++j) {
                std::cout << elements[i][j] << (j == cols - 1 ? "" : " ");
            }
            std::cout << '\n';
        }
    }
};

int main() {
    std::ios_base::sync_with_stdio(false);
    std::cin.tie(NULL);

    int r1, c1; // 첫 번째 행렬 A의 행과 열 크기
    std::cin >> r1 >> c1;
    Matrix matrixA(r1, c1);
    for (int i = 0; i < r1; ++i) {
        for (int j = 0; j < c1; ++j) {
            std::cin >> matrixA.elements[i][j];
        }
    }

    int c2; // 두 번째 행렬 B의 열 크기 (행 크기는 c1과 같아야 함)
    std::cin >> c2;
    Matrix matrixB(c1, c2); // matrixB의 행은 matrixA의 열과 같아야 함
    for (int i = 0; i < c1; ++i) {
        for (int j = 0; j < c2; ++j) {
            std::cin >> matrixB.elements[i][j];
        }
    }
    
    try {
        Matrix matrixC = matrixA * matrixB;
        matrixC.display();
    } catch (const std::runtime_error& e) {
        std::cerr << "오류: " << e.what() << std::endl;
    }

    return 0;
}

2. 행렬 고속 거듭제곱

2.1. 정의

일반적인 숫자의 거듭제곱을 효율적으로 계산하는 고속 거듭제곱(Exponentiation by Squaring) 알고리즘처럼, 행렬의 거듭제곱 \\(A^k\\) 역시 동일한 원리로 빠르게 계산할 수 있습니다. 행렬 \\(A^k = A \cdot A \cdot \ldots \cdot A\\) (\\(k\\)번 곱함)는 이 알고리즘을 통해 \\(O(\log k)\\)의 시간 복잡도로 구해집니다. 행렬 거듭제곱을 수행하려면 해당 행렬 \\(A\\)가 반드시 정방 행렬(square matrix)이어야 합니다.

2.2. 활용 사례

행렬 고속 거듭제곱은 주로 복잡한 점화식(recurrence relations)의 계산을 가속화하는 데 사용됩니다. 몇 가지 대표적인 예시를 통해 그 활용법을 살펴보겠습니다.

2.2.1. 예시 1: 피보나치 수열의 N번째 항

문제: 피보나치 수열의 \\(n\\)번째 항을 \\(m\\)으로 나눈 나머지를 구하시오 (\\(n \le 10^9\\)). 해결: \\(n\\)의 값이 매우 크기 때문에 \\(O(n)\\) 시간 복잡도를 가지는 단순 반복문으로는 해결할 수 없습니다. 대신 행렬 고속 거듭제곱을 이용하여 \\(O(\log n)\\) 시간 내에 해결할 수 있습니다. 피보나치 수열의 점화식은 \\(F_k = F_{k-1} + F_{k-2}\\) 입니다. 이 점화식을 행렬 형태로 표현하면 다음과 같습니다: \\(①\\) 먼저 \\(F_n\\) 항을 포함하는 행렬을 만듭니다:
\[\begin{bmatrix} F_n \\ \ \ \end{bmatrix} = \begin{bmatrix} 1 & 1\\ \ & \ \end{bmatrix} \begin{bmatrix} F_{n-1} \\ F_{n-2} \end{bmatrix} \]
\\(②\\) 오른쪽 행렬의 원소들이 \\(n-1\\), \\(n-2\\)와 같이 인덱스 차이가 발생하므로, 왼쪽 행렬은 \\(\\begin{bmatrix} F\_n \\\\ F\_{n-1} \\end{bmatrix}\\)과 같이 구성하여 일관성을 유지합니다. \\(③\\) 이를 통해 완전한 점화 행렬식을 구성할 수 있습니다:
\[\begin{bmatrix} F_n \\ F_{n-1} \ \end{bmatrix}=\begin{bmatrix} 1 & 1\\ 1 & 0 \end{bmatrix}\begin{bmatrix} F_{n-1} \\ F_{n-2} \end{bmatrix} \]
이 행렬식을 바탕으로 변환 행렬 \\(T = \\begin{bmatrix} 1 & 1\\\\ 1 & 0 \\end{bmatrix}\\)에 대해 \\((n-1)\\)번 고속 거듭제곱을 수행하면 \\(T^{n-1}\\)을 얻습니다. 여기에 초기 상태 벡터 \\(\\begin{bmatrix} F\_{1} \\\\ F\_{0} \\end{bmatrix} = \\begin{bmatrix} 1 \\\\ 0 \\end{bmatrix}\\)을 곱하면 \\(\\begin{bmatrix} F\_n \\\\ F\_{n-1} \\end{bmatrix}\\)을 구할 수 있습니다 (일반적으로 \\(n \ge 1\\)에 대해 \\(F_0=0, F_1=1\\)을 가정). 결과 행렬의 첫 번째 원소가 \\(F_n\\)이 됩니다. 코드 구현 (피보나치 수열): 행렬 구조체를 활용하여 고속 거듭제곱 연산자를 추가합니다.
#include <iostream>
#include <vector>
#include <stdexcept> // 예외 처리를 위해

long long FIB_MOD_VALUE; // 모듈러 연산을 위한 전역 변수

// 행렬을 표현하는 구조체 (모듈러 연산 포함)
struct AdvancedMatrix {
    int dim; // 정방 행렬의 차원 (rows == cols)
    std::vector<std::vector<long long>> elements;

    AdvancedMatrix() : dim(0) {}
    AdvancedMatrix(int d) : dim(d), elements(d, std::vector<long long>(d, 0)) {}

    // 행렬 곱셈 연산자 오버로드 (모듈러 연산 적용)
    AdvancedMatrix operator*(const AdvancedMatrix& other) const {
        if (dim != other.dim) {
            throw std::runtime_error("행렬의 차원이 곱셈에 적합하지 않습니다.");
        }

        AdvancedMatrix result(dim);
        for (int i = 0; i < dim; ++i) {
            for (int j = 0; j < dim; ++j) {
                for (int k = 0; k < dim; ++k) {
                    result.elements[i][j] = (result.elements[i][j] + elements[i][k] * other.elements[k][j]) % FIB_MOD_VALUE;
                }
            }
        }
        return result;
    }

    // 행렬 고속 거듭제곱 연산자 오버로드
    AdvancedMatrix operator^(long long power) const {
        if (dim == 0) return AdvancedMatrix(0);
        
        AdvancedMatrix res(dim);
        // 결과 행렬을 단위 행렬로 초기화
        for (int i = 0; i < dim; ++i) {
            res.elements[i][i] = 1;
        }

        AdvancedMatrix base = *this; // 현재 행렬을 밑으로 사용
        while (power > 0) {
            if (power % 2 == 1) { // power가 홀수이면 결과에 현재 base를 곱함
                res = res * base;
            }
            base = base * base; // base를 제곱하여 다음 단계로 이동
            power /= 2;
        }
        return res;
    }
};

int main() {
    std::ios_base::sync_with_stdio(false);
    std::cin.tie(NULL);

    long long n_term; // 구하고자 하는 피보나치 항의 번호
    std::cin >> n_term >> FIB_MOD_VALUE;

    if (n_term == 0) { // F_0 = 0
        std::cout << 0 << std::endl;
        return 0;
    }
    if (n_term == 1) { // F_1 = 1
        std::cout << 1 % FIB_MOD_VALUE << std::endl;
        return 0;
    }

    AdvancedMatrix fibTransform(2); // 2x2 변환 행렬
    fibTransform.elements = {{1, 1}, {1, 0}};

    // F_n을 구하려면 (n-1)번 변환 행렬을 곱해야 합니다.
    // T^(n-1)의 (0,0) 원소가 F_n이 됩니다.
    fibTransform = fibTransform ^ (n_term - 1);
    std::cout << fibTransform.elements[0][0] << std::endl;

    return 0;
}

2.2.2. 예시 2: 미로 찾기 (P207)

해결: 만약 모든 간선의 가중치가 1이었다면, 단순히 인접 행렬을 \\(t\\)번 거듭제곱하여 \\(t\\)번 이동 후의 경로 수를 계산할 수 있을 것입니다. 그러나 이 문제에서는 간선에 가중치가 부여되어 있습니다. 중요한 점은 간선 가중치가 최대 9로 매우 작다는 것입니다. 이러한 특성을 활용하여, 각 노드를 여러 "층(layer)"으로 분할하여 가중치 있는 간선을 가중치 1인 간선으로 변환하는 기법을 적용할 수 있습니다. 예를 들어, 각 노드 \\(u\\)를 \\(u\_0, u\_1, \ldots, u\_8\\)과 같이 9개의 서브 노드로 분할합니다 (0-based 인덱스).
  • 가중치 \\(w\\)인 간선 \\(u \to v\\)는 새로운 그래프에서 \\(u\_{w-1} \to v\_0\\)로 연결됩니다.
  • 또한, 각 노드 \\(u\\) 내에서 \\(u\_i \to u\_{i+1}\\) (\\(i=0, \ldots, 7\\)) 간선을 추가하여, \\(u\_0\\)에서 \\(u\_{w-1}\\)까지 \\(w-1\\) 단계를 이동할 수 있게 합니다. 이 간선들의 가중치는 모두 1입니다.
이 변환을 통해 모든 간선의 가중치를 1로 만들 수 있으며, 결과적으로 전체 그래프의 크기는 커지지만, 행렬 고속 거듭제곱을 적용할 수 있는 형태로 변경됩니다. 변환된 그래프의 인접 행렬을 \\(t\\)번 거듭제곱한 후, 시작 노드에 해당하는 서브 노드에서 목표 노드에 해당하는 서브 노드까지의 경로 수를 찾아 최종 답을 구합니다. 일반적으로 시작 노드 1 (0-based: 0)에서 목표 노드 \\(n\\) (0-based: \\(n-1\\))까지의 경로를 찾으려면, \\(1\_0\\)에서 \\(n\_0\\)까지의 경로를 계산합니다. 코드 구현 (미로 찾기): 최대 간선 가중치에 따라 행렬의 차원을 조절합니다.
#include <iostream>
#include <vector>
#include <string>
#include <stdexcept>

// 모듈러 값 정의
const long long MAZE_MOD_VALUE = 2009;
const int MAX_EDGE_WEIGHT = 9; // 최대 간선 가중치 (1부터 9까지)

// 계층형 그래프 행렬 구조체
struct LayeredMatrix {
    int dimension; // 행렬의 차원 (rows == cols)
    std::vector<std::vector<long long>> elements;

    LayeredMatrix() : dimension(0) {}
    LayeredMatrix(int d) : dimension(d), elements(d, std::vector<long long>(d, 0)) {}

    // 행렬 곱셈 (모듈러 연산 포함)
    LayeredMatrix operator*(const LayeredMatrix& other) const {
        if (dimension != other.dimension) {
            throw std::runtime_error("행렬의 차원이 곱셈에 적합하지 않습니다.");
        }
        LayeredMatrix result(dimension);
        for (int i = 0; i < dimension; ++i) {
            for (int j = 0; j < dimension; ++j) {
                for (int k = 0; k < dimension; ++k) {
                    result.elements[i][j] = (result.elements[i][j] + elements[i][k] * other.elements[k][j]) % MAZE_MOD_VALUE;
                }
            }
        }
        return result;
    }

    // 행렬 고속 거듭제곱
    LayeredMatrix operator^(long long power) const {
        if (dimension == 0) return LayeredMatrix(0);
        
        LayeredMatrix res(dimension);
        // 단위 행렬로 초기화
        for (int i = 0; i < dimension; ++i) {
            res.elements[i][i] = 1;
        }

        LayeredMatrix base = *this;
        while (power > 0) {
            if (power % 2 == 1) {
                res = res * base;
            }
            base = base * base;
            power /= 2;
        }
        return res;
    }
};

int main() {
    std::ios_base::sync_with_stdio(false);
    std::cin.tie(NULL);

    int num_original_nodes; // 원래 그래프의 노드 수
    long long num_steps;    // 이동 횟수 (t)
    std::cin >> num_original_nodes >> num_steps;

    // 변환된 그래프의 총 차원: (원래 노드 수) * (최대 간선 가중치)
    int total_dimension = num_original_nodes * MAX_EDGE_WEIGHT;
    LayeredMatrix adjacencyMatrix(total_dimension);

    // 각 노드 u 내에서 u_i -> u_{i+1} 연결 (가중치 1)
    // 0-based original node index: i, 0-based layer index: j
    for (int i = 0; i < num_original_nodes; ++i) {
        for (int j = 0; j < MAX_EDGE_WEIGHT - 1; ++j) {
            // (i * MAX_EDGE_WEIGHT + j)는 현재 노드 u의 j번째 층 서브 노드 인덱스
            // (i * MAX_EDGE_WEIGHT + j + 1)는 현재 노드 u의 (j+1)번째 층 서브 노드 인덱스
            adjacencyMatrix.elements[i * MAX_EDGE_WEIGHT + j][i * MAX_EDGE_WEIGHT + j + 1] = 1;
        }
    }

    // 원래 그래프의 간선 추가: u -> v (가중치 'w')는 u_{w-1} -> v_0 로 변환
    for (int i = 0; i < num_original_nodes; ++i) {
        for (int j = 0; j < num_original_nodes; ++j) {
            char weight_char;
            std::cin >> weight_char;
            int weight_val = weight_char - '0';
            if (weight_val > 0) {
                // 출발 노드 i의 (weight_val - 1)번째 층 서브 노드에서 (0-based)
                // 도착 노드 j의 0번째 층 서브 노드로 연결
                adjacencyMatrix.elements[i * MAX_EDGE_WEIGHT + (weight_val - 1)][j * MAX_EDGE_WEIGHT] = 1;
            }
        }
    }

    // num_steps번 이동 후의 경로 수를 계산
    LayeredMatrix finalStateMatrix = adjacencyMatrix ^ num_steps;

    // 시작 노드 1 (0-based: 0)의 첫 번째 층 (0)에서
    // 끝 노드 n (0-based: num_original_nodes-1)의 첫 번째 층 (0)까지의 경로 수
    std::cout << finalStateMatrix.elements[0][(num_original_nodes - 1) * MAX_EDGE_WEIGHT] << std::endl;

    return 0;
}

2.2.3. 예시 3: 佳佳의 피보나치 (\\(F_1+2F_2+\cdots+nF_n\\))

문제: \\((F\_1+2F\_2+3F\_3+\\cdots+nF\_n)\\bmod{m}\\)의 값을 구하시오. 해결: 이 합계 문제를 행렬 고속 거듭제곱으로 해결하기 위해 새로운 점화식을 유도해야 합니다. 먼저, 피보나치 수열의 부분합을 \\(S_k = F_1 + F_2 + \cdots + F_k\\)로 정의합니다. 주어진 합계 \\(T_n = F\_1+2F\_2+3F\_3+\cdots+nF\_n\\)는 다음과 같이 변형할 수 있습니다:
\[\begin{aligned}&T_n = n\times F_n + (n-1)F_{n-1} + \dots + 1F_1 \\&= n S_n - (S_{n-1} + S_{n-2} + \dots + S_1) \\&= n S_n - \sum_{i=1}^{n-1}S_i \end{aligned} \]
여기서 \\(P_n = \sum_{i=1}^{n-1}S_i\\)로 두면, 원래 합계는 \\(n \times S_n - P_n\\)으로 간소화됩니다. 이제 \\(P_n, S_n\\) 및 피보나치 항들을 포함하는 상태 벡터에 대한 점화 행렬을 구성합니다. 관련된 점화식은 다음과 같습니다:
  • \\(P_k = P_{k-1} + S_{k-1}\\) (\\(P_k\\)는 \\(S_1 \dots S_{k-1}\\)의 합)
  • \\(S_k = S_{k-1} + F_k\\) (\\(S_k\\)는 \\(F_1 \dots F_k\\)의 합)
  • \\(F_k = F_{k-1} + F_{k-2}\\)
이 점화식들을 기반으로 상태 벡터를 \\(\\begin{bmatrix} P\_k \\\\ S\_k \\\\ F\_{k+1} \\\\ F\_k \\end{bmatrix}\\)로 설정하고, 이 벡터를 \\(\\begin{bmatrix} P_{k+1} \\\\ S_{k+1} \\\\ F_{k+2} \\\\ F_{k+1} \\end{bmatrix}\\)로 변환하는 행렬을 찾아봅시다: \\(①\\) \\(P_{k+1} = P_k + S_k\\): \\(\\begin{bmatrix} P_{k+1} \\\\ S_{k+1} \\\\ F_{k+2} \\\\ F_{k+1} \end{bmatrix}=\begin{bmatrix} 1 & 1 & 0 & 0 \\\\ \cdot & \cdot & \cdot & \cdot \\\\ \cdot & \cdot & \cdot & \cdot \\\\ \cdot & \cdot & \cdot & \cdot \end{bmatrix}\begin{bmatrix} P_k \\\\ S_k \\\\ F_{k+1} \\\\ F_k \end{bmatrix}\\) \\(②\\) \\(S_{k+1} = S_k + F_{k+1}\\): \\(\\begin{bmatrix} P_{k+1} \\\\ S_{k+1} \\\\ F_{k+2} \\\\ F_{k+1} \end{bmatrix}=\begin{bmatrix} 1 & 1 & 0 & 0 \\\\ 0 & 1 & 1 & 0 \\\\ \cdot & \cdot & \cdot & \cdot \\\\ \cdot & \cdot & \cdot & \cdot \end{bmatrix}\begin{bmatrix} P_k \\\\ S_k \\\\ F_{k+1} \\\\ F_k \end{bmatrix}\\) \\(③\\) \\(F_{k+2} = F_{k+1} + F_k\\) 및 \\(F_{k+1}\\)은 그대로 전달:
\[\begin{bmatrix} P_{k+1} \\ S_{k+1} \\ F_{k+2} \\ F_{k+1} \end{bmatrix}=\begin{bmatrix} 1 & 1 & 0 & 0 \\ 0 & 1 & 1 & 0 \\ 0 & 0 & 1 & 1 \\ 0 & 0 & 1 & 0 \end{bmatrix}\begin{bmatrix} P_k \\ S_k \\ F_{k+1} \\ F_k \end{bmatrix} \]
따라서 최종 변환 행렬은 위와 같습니다. 초기 상태 벡터는 \\(k=0\\)일 때를 기준으로 설정합니다. \\(P_0 = 0\\) (\\(S_i\\)의 합은 없음) \\(S_0 = 0\\) (\\(F_i\\)의 합은 없음) \\(F_1 = 1\\) \\(F_0 = 0\\) 그러므로 초기 상태 벡터는 \\(\\begin{bmatrix} 0 \\\\ 0 \\\\ 1 \\\\ 0 \end{bmatrix}\\)입니다. 이 초기 벡터에 변환 행렬을 \\(n\\)번 거듭제곱하여 곱한 후, 결과 벡터에서 \\(P_n\\)과 \\(S_n\\) 값을 추출하여 \\(n \times S_n - P_n\\)을 계산합니다. 코드 구현 (합계 계산):
#include <iostream>
#include <vector>
#include <stdexcept>

long long SUM_MOD_VALUE; // 모듈러 연산을 위한 전역 변수

// 합계 계산을 위한 행렬 구조체
struct SumMatrix {
    int rows, cols;
    std::vector<std::vector<long long>> data;

    SumMatrix(int r, int c) : rows(r), cols(c), data(r, std::vector<long long>(c, 0)) {}
    // 정방 행렬 생성을 위한 생성자
    SumMatrix(int d) : rows(d), cols(d), data(d, std::vector<long long>(d, 0)) {}


    // 행렬 곱셈 연산자
    SumMatrix operator*(const SumMatrix& other) const {
        if (cols != other.rows) { // this.cols != other.rows
            throw std::runtime_error("행렬의 차원이 곱셈에 적합하지 않습니다.");
        }
        SumMatrix result(rows, other.cols); // rows x other.cols
        for (int i = 0; i < rows; ++i) {
            for (int j = 0; j < other.cols; ++j) {
                for (int k = 0; k < cols; ++k) {
                    result.data[i][j] = (result.data[i][j] + data[i][k] * other.data[k][j]) % SUM_MOD_VALUE;
                }
            }
        }
        return result;
    }

    // 행렬 고속 거듭제곱 연산자 (정방 행렬에만 적용)
    SumMatrix operator^(long long power) const {
        if (rows != cols) { // 반드시 정방 행렬이어야 함
            throw std::runtime_error("행렬은 거듭제곱을 위해 정방 행렬이어야 합니다.");
        }
        SumMatrix res(rows); // 단위 행렬로 초기화
        for (int i = 0; i < rows; ++i) {
            res.data[i][i] = 1;
        }

        SumMatrix base = *this;
        while (power > 0) {
            if (power % 2 == 1) {
                res = res * base;
            }
            base = base * base;
            power /= 2;
        }
        return res;
    }
};

int main() {
    std::ios_base::sync_with_stdio(false);
    std::cin.tie(NULL);

    long long n_val; // 합계의 상한값
    std::cin >> n_val >> SUM_MOD_VALUE;

    // 변환 행렬 (4x4)
    SumMatrix transformMatrix(4);
    transformMatrix.data = {
        {1, 1, 0, 0},
        {0, 1, 1, 0},
        {0, 0, 1, 1},
        {0, 0, 1, 0}
    };

    // 초기 상태 벡터 (4x1)
    SumMatrix initialStateVector(4, 1);
    // [P_0, S_0, F_1, F_0]^T = [0, 0, 1, 0]^T
    initialStateVector.data = {
        {0}, // P_0
        {0}, // S_0
        {1}, // F_1
        {0}  // F_0
    };
    
    // n이 0인 경우, 합계는 0
    if (n_val == 0) {
        std::cout << 0 << std::endl;
        return 0;
    }

    // 변환 행렬을 n번 거듭제곱
    SumMatrix poweredMatrix = transformMatrix ^ n_val;

    // 거듭제곱된 행렬과 초기 상태 벡터 곱셈
    SumMatrix finalState = poweredMatrix * initialStateVector;

    // 결과 추출: finalState 벡터의 첫 번째 원소가 P_n, 두 번째 원소가 S_n
    long long P_n = finalState.data[0][0];
    long long S_n = finalState.data[1][0];

    // 최종 결과 계산: (n * S_n - P_n) % m
    long long result = (n_val % SUM_MOD_VALUE * S_n % SUM_MOD_VALUE - P_n + SUM_MOD_VALUE) % SUM_MOD_VALUE;
    std::cout << result << std::endl;

    return 0;
}

태그: MatrixMultiplication MatrixExponentiation RecurrenceRelation DynamicProgramming ComputationalMathematics

8월 1일 02:25에 게시됨