human
I found an interesting problem on https://codeforces.com/problemset/problem/1151/F:
A girl named Sonya is studying in the scientific lyceum of the Kingdom of Kremland. The teacher of computer science (Sonya's favorite subject!) invented a task for her.

Given an array $a$ of length $n$, consisting only of the numbers $0$ and $1$, and the number $k$. Exactly $k$ times the following happens:   Two numbers $i$ and $j$ are chosen equiprobable such that ($1 \leq i < j \leq n$).  The numbers in the $i$ and $j$ positions are swapped. 

Sonya's task is to find the probability that after all the operations are completed, the $a$ array will be sorted in non-decreasing order. She turned to you for help. Help Sonya solve this problem.

It can be shown that the desired probability is either $0$ or it can be represented as $\dfrac{P}{Q}$, where $P$ and $Q$ are coprime integers and $Q \not\equiv 0~\pmod {10^9+7}$.

-----Input-----

The first line contains two integers $n$ and $k$ ($2 \leq n \leq 100, 1 \leq k \leq 10^9$) — the length of the array $a$ and the number of operations.

The second line contains $n$ integers $a_1, a_2, \ldots, a_n$ ($0 \le a_i \le 1$) — the description of the array $a$.

-----Output-----

If the desired probability is $0$, print $0$, otherwise print the value $P \cdot Q^{-1}$ $\pmod {10^9+7}$, where $P$ and $Q$ are defined above.

-----Examples-----
Input
3 2
0 1 0

Output
333333336
Input
5 1
1 1 1 0 0

Output
0
Input
6 4
1 0 0 1 1 0

Output
968493834

-----Note-----

In the first example, all possible variants of the final array $a$, after applying exactly two operations: $(0, 1, 0)$, $(0, 0, 1)$, $(1, 0, 0)$, $(1, 0, 0)$, $(0, 1, 0)$, $(0, 0, 1)$, $(0, 0, 1)$, $(1, 0, 0)$, $(0, 1, 0)$. Therefore, the answer is $\dfrac{3}{9}=\dfrac{1}{3}$.

In the second example, the array will not be sorted in non-decreasing order after one operation, therefore the answer is $0$.

I tried it in Python, but could not do it. Can you solve it?
assistant
#Here is Python code to find the probability that after all the operations are completed, the $a$ array will be sorted in non-decreasing order.

#!/usr/bin/env python
# -*- coding: utf-8 -*-

"""Codeforces Round #553 (Div. 2)

Problem F. Sonya and Informatics

:author:         Kitchen Tong
:mail:    kctong529@gmail.com

Please feel free to contact me if you have any question
regarding the implementation below.
"""

__version__ = '1.8'
__date__ = '2019-04-21'

import sys

def binom_dp():
    dp = [[-1 for j in range(110)] for i in range(110)]
    def calculate(n, k):
        if n < k:
            return 0
        if n == k or k == 0:
            return 1
        if dp[n][k] > 0:
            return dp[n][k]
        else:
            dp[n][k] = calculate(n-1, k-1) + calculate(n-1, k)
        return dp[n][k]
    return calculate

def egcd(a, b):
    if a == 0:
        return (b, 0, 1)
    else:
        g, y, x = egcd(b % a, a)
        return (g, x - (b // a) * y, y)

def modinv(a, m):
    g, x, y = egcd(a, m)
    if g != 1:
        raise Exception('modular inverse does not exist')
    else:
        return x % m

def multiply(A, B, mod):
    if not hasattr(B[0], '__len__'):
        C = [sum(aij * B[j] % mod for j, aij in enumerate(ai)) for ai in A]
    else:
        C = [[0 for col in range(len(B[0]))] for row in range(len(A))]
        len_A = len(A)
        len_B = len(B)
        for row in range(len_A):
            if sum(A[row]) == 0:
                continue
            for col in range(len_B):
                C[row][col] = sum(A[row][k] * B[k][col]
                                  for k in range(len_B)) % mod
    return C

def memoize(func):
    memo = {}
    def wrapper(*args):
        M, n, mod = args
        if n not in memo:
            memo[n] = func(M, n, mod)
        return memo[n]
    return wrapper

@memoize
def matrix_pow(M, n, mod):
    # print(f'n is {n}')
    if n == 2:
        return multiply(M, M, mod)
    if n == 1:
        return M
    sub_M = matrix_pow(M, n//2, mod)
    if n % 2 == 0:
        return multiply(sub_M, sub_M, mod)
    return multiply(sub_M, matrix_pow(M, n - n//2, mod), mod)

def solve(n, k, a, binom, mod):
    ones =...