#1169609 (C++17) No.541 3 x N グリッド上のサイクルの個数

• 問題: No.541 3 x N グリッド上のサイクルの個数
• ユーザー: hly1204
• 提出日時: 2026-06-24 00:45:13
• 言語: C++17 (gcc 15.2.0 + boost 1.89.0)
• 結果: AC
• 実行時間: 2 ms / 2,000 ms
• コード長: 4,598 bytes
• コンパイル時間: 839 ms
• コンパイル使用メモリ: 101,148 KB
• 実行使用メモリ: 6,400 KB
• 最終ジャッジ日時: 2026-06-24 00:45:17
• 合計ジャッジ時間: 2,834 ms
• ジャッジサーバーID （参考情報）: judge2_0 / judge1_0

ファイルパターン	結果
other	AC * 62

ソースコード

// competitive-verifier: PROBLEM https://yukicoder.me/problems/no/2883

#include <algorithm>
#include <array>
#include <cassert>
#include <iostream>
#include <utility>
#include <vector>

using uint         = unsigned;
using ull          = unsigned long long;
constexpr uint MOD = 1000000007;

constexpr uint PowMod(uint a, ull e) {
    for (uint res = 1;; a = (ull)a * a % MOD) {
        if (e & 1) res = (ull)res * a % MOD;
        if ((e /= 2) == 0) return res;
    }
}

constexpr uint InvMod(uint a) { return PowMod(a, MOD - 2); }

// returns */Q in F((x^-1)) such that */Q = A[0]*x^-1 + A[1]*x^-2 + ...
// with deg(Q) minimized.
std::vector<uint> BerlekampMassey(std::vector<uint> A) {
    reverse(begin(A), end(A));
    int k = size(A), degA = k - 1, degB = k;
    std::vector<uint> B(k + 1), Q0, Q1  = {1U};
    for (B[k] = 1;; swap(Q0, Q1), swap(A, B), std::swap(degA, degB)) {
        while (degA >= 0 && A[degA] == 0) --degA;
        if (degA < 0 || degA - degB < -k) return Q1;
        Q0.resize(size(Q1) + (degB - degA));
        k -= (degB - degA) * 2;
        const uint a = InvMod(A[degA]);
        for (int i = degB - degA; i >= std::max(-k, 0); --i) {
            const uint d = (ull)B[degB--] * a % MOD;
            for (int j = 0; j <= degA; ++j) B[i + j] = (B[i + j] + (ull)A[j] * (MOD - d)) % MOD;
            for (int j = 0; j < (int)size(Q1); ++j) Q0[i + j] = (Q0[i + j] + (ull)Q1[j] * d) % MOD;
        }
    }
}

// returns [x^[-deg(Q), 0)] x^k/Q, x^k/Q in F((x^-1))
std::vector<uint> BostanMoriT(const std::vector<uint> &Q, long long k) {
    assert(k >= 0);
    const int degQ = size(Q) - 1;
    std::vector<uint> U(degQ);
    if (k == 0) return U[0] = InvMod(Q.back()), U;
    std::vector<uint> V(size(Q));
    for (int i = 0; i < (int)size(Q); ++i)
        for (int j = i % 2; j < (int)size(Q); j += 2)
            V[(i + j) / 2] = (V[(i + j) / 2] + (ull)Q[i] * (i % 2 == 0 ? Q[j] : MOD - Q[j])) % MOD;
    const auto T = BostanMoriT(V, k / 2);
    for (int i = 0; i < (int)size(T); ++i)
        for (int j = 0; j < (int)size(Q); ++j) {
            const int l = i * 2 + (int)(k % 2) + j;
            if (l >= degQ && l < degQ * 2)
                U[l - degQ] = (U[l - degQ] + (ull)T[i] * (j % 2 == 0 ? Q[j] : MOD - Q[j])) % MOD;
        }
    return U;
}

// returns x^k mod Q
std::vector<uint> MonomialModPoly(long long k, const std::vector<uint> &Q) {
    const auto invQ = BostanMoriT(Q, k);
    std::vector<uint> R(size(Q) - 1);
    for (int i = 0; i < (int)size(invQ); ++i)
        for (int j = 0; j < (int)size(Q); ++j)
            if (i + j >= (int)size(invQ))
                R[i + j - (int)size(invQ)] =
                    (R[i + j - (int)size(invQ)] + (ull)invQ[i] * Q[j]) % MOD;
    return R;
}

uint BMBM(const std::vector<uint> &A, long long k) {
    const auto c = MonomialModPoly(k, BerlekampMassey(A));
    uint res     = 0;
    for (int i = 0; i < (int)size(c); ++i) res = (res + (ull)c[i] * A[i]) % MOD;
    return res;
}

int main() {
    std::ios::sync_with_stdio(false);
    std::cin.tie(nullptr);
    long long N;
    std::cin >> N;

    uint dp[50][0b1000] = {};

    dp[1][0b1]   = 1;
    dp[1][0b10]  = 1;
    dp[1][0b100] = 1;
    dp[1][0b11]  = 1;
    dp[1][0b110] = 1;
    dp[1][0b111] = 1;

    for (int i = 2; i < 50; ++i) {
        for (int j = 0b1; j <= 0b111; ++j) {
            if (j == 0b101) {
                dp[i][j] = ((ull)dp[i][j] + dp[i - 1][0b111] + dp[i - 1][0b101]) % MOD;
                continue;
            }
            for (int k = 0b1; k <= 0b111; ++k) {
                if (k == 0b101) {
                    if (j == 0b11 || j == 0b110) continue;
                    if (j == 0b111) {
                        dp[i][j] = (dp[i][j] + 1) % MOD;
                        if (i > 2) {
                            for (int k = 1; k <= i - 2; ++k) {
                                dp[i][j] = ((ull)dp[i][j] + dp[k][0b1] + dp[k][0b100]) % MOD;
                            }
                        }
                        continue;
                    }
                }
                if ((j & k) > 0) dp[i][j] = (dp[i][j] + dp[i - 1][k]) % MOD;
            }
        }
    }

    auto sum = [&](int r) {
        uint res = 0;
        for (int i = 0b1; i <= 0b111; ++i) res = (res + dp[r][i]) % MOD;
        return res;
    };

    std::vector<uint> A(49);

    for (int i = 1; i < 50; ++i) {
        for (int j = 1; j <= i; ++j) {
            A[i - 1] = (A[i - 1] + (ull)sum(j) * (i - j + 1)) % MOD;
        }
        // std::cout << A[i - 1] << '\n';
    }

    std::cout << BMBM(A, N - 1) << '\n';

    return 0;
}