from math import gcd


TARGET_N = 123
TARGET_K = 4_567_891
TARGET_E = TARGET_K // 2
MOD = 1_234_567_891


def binomial_table(n, r, mod):
    c = [[0] * (r + 1) for _ in range(n + 1)]
    c[0][0] = 1
    for i in range(1, n + 1):
        c[i][0] = 1
        hi = min(i, r)
        for j in range(1, hi + 1):
            value = c[i - 1][j - 1] + c[i - 1][j]
            c[i][j] = value if mod == 0 else value % mod
    return c


def choose_from_table(c, n, r):
    if n < 0 or r < 0 or n < r:
        return 0
    if n >= len(c) or r >= len(c[0]):
        return 0
    return c[n][r]


def run_count_table(n, emax, mod):
    comb = binomial_table(2 * emax + 1, n, mod)
    runs = [[0] * (emax + 1) for _ in range(n + 1)]

    for length in range(2, n + 1):
        for e in range(1, emax + 1):
            all_positive = choose_from_table(comb, 2 * e - 1, length - 1)
            too_large = choose_from_table(comb, e - 1, length - 1)
            if mod == 0:
                runs[length][e] = all_positive - length * too_large
            else:
                runs[length][e] = (all_positive - length * too_large) % mod
    return runs


def cumulative_counts(n, emax, mod):
    runs = run_count_table(n, emax, mod)
    total = [[0] * (emax + 1) for _ in range(n + 1)]
    zero = [[0] * (emax + 1) for _ in range(n + 1)]
    total[0][0] = 1
    zero[0][0] = 1

    for i in range(1, n + 1):
        zero[i] = total[i - 1][:]
        row = zero[i][:]

        for length in range(2, i + 1):
            source = zero[i - length]
            run = runs[length]
            for before, source_value in enumerate(source):
                if source_value == 0:
                    continue
                for used in range(1, emax - before + 1):
                    run_value = run[used]
                    if run_value == 0:
                        continue
                    if mod == 0:
                        row[before + used] += source_value * run_value
                    else:
                        row[before + used] = (row[before + used] + source_value * run_value) % mod

        total[i] = row

    cumulative = []
    running = 0
    for e in range(emax + 1):
        running += total[n][e]
        if mod:
            running %= mod
        cumulative.append(running)
    return cumulative


def binomial_large_mod(n, r):
    if r < 0 or n < r:
        return 0
    if r > n - r:
        r = n - r

    numerator = [n - r + 1 + i for i in range(r)]
    for d in range(2, r + 1):
        remaining = d
        for i, value in enumerate(numerator):
            g = gcd(value, remaining)
            if g == 1:
                continue
            numerator[i] = value // g
            remaining //= g
            if remaining == 1:
                break
        assert remaining == 1

    result = 1
    for value in numerator:
        result = result * (value % MOD) % MOD
    return result


def polynomial_tail_value(n, e):
    start = n - 1
    degree = n
    emax = start + degree
    cumulative = cumulative_counts(n, emax, MOD)

    current = [cumulative[start + i] for i in range(degree + 1)]
    differences = []
    for _ in range(degree + 1):
        differences.append(current[0])
        current = [(current[i + 1] - current[i]) % MOD for i in range(len(current) - 1)]

    x = e - start
    result = 0
    for r, diff in enumerate(differences):
        result = (result + diff * binomial_large_mod(x, r)) % MOD
    return result


def brute_count(n, emax):
    identity = tuple(range(n))
    states = {(identity, (0,) * n)}
    results = {(0,) * n}

    for _ in range(1, 2 * emax + 1):
        nxt = set()
        for perm, count in states:
            for p in range(n - 1):
                next_perm = list(perm)
                next_count = list(count)
                next_count[next_perm[p + 1]] += 1
                next_perm[p], next_perm[p + 1] = next_perm[p + 1], next_perm[p]
                next_perm_t = tuple(next_perm)
                next_count_t = tuple(next_count)
                if next_perm_t == identity:
                    results.add(next_count_t)
                nxt.add((next_perm_t, next_count_t))
        states = nxt

    return len(results)


def run_checkpoints():
    for n in range(2, 6):
        values = cumulative_counts(n, 4, 0)
        for e in range(5):
            assert values[e] == brute_count(n, e)

    assert cumulative_counts(3, 2, 0)[2] == 8
    assert cumulative_counts(12, 17, 0)[17] == 2_457_178_250

    for n in range(2, 9):
        start = n - 1
        degree = n
        emax = start + degree + 5
        values = cumulative_counts(n, emax, 0)

        current = [values[start + i] for i in range(degree + 1)]
        differences = []
        for _ in range(degree + 1):
            differences.append(current[0])
            current = [current[i + 1] - current[i] for i in range(len(current) - 1)]

        for e in range(start, emax + 1):
            predicted = 0
            choose = 1
            for r, diff in enumerate(differences):
                if r > 0:
                    choose = choose * (e - start - r + 1) // r
                predicted += diff * choose
            assert predicted == values[e]


def main():
    run_checkpoints()
    print(polynomial_tail_value(TARGET_N, TARGET_E))


if __name__ == "__main__":
    main()