Problem:
 Score : 100 points    Problem Statement There is a simple directed graph G with N vertices, numbered 1, 2, \ldots, N .  For each i and j ( 1 \leq i, j \leq N ), you are given an integer a_{i, j} that represents whether there is a directed edge from Vertex i to j .
If a_{i, j} = 1 , there is a directed edge from Vertex i to j ; if a_{i, j} = 0 , there is not.  Find the number of different directed paths of length K in G , modulo 10^9 + 7 .
We will also count a path that traverses the same edge multiple times.      Constraints  All values in input are integers.  1 \leq N \leq 50  1 \leq K \leq 10^{18}  a_{i, j} is 0 or 1 .  a_{i, i} = 0         Input Input is given from Standard Input in the following format:  N  K  a_{1, 1}  \ldots  a_{1, N}  :  a_{N, 1}  \ldots  a_{N, N}       Output Print the number of different directed paths of length K in G , modulo 10^9 + 7 .        Sample Input 1 4 2
0 1 0 0
0 0 1 1
0 0 0 1
1 0 0 0      Sample Output 1 6  G is drawn in the figure below:   There are six directed paths of length 2 :   1 → 2 → 3  1 → 2 → 4  2 → 3 → 4  2 → 4 → 1  3 → 4 → 1  4 → 1 → 2        Sample Input 2 3 3
0 1 0
1 0 1
0 0 0      Sample Output 2 3  G is drawn in the figure below:   There are three directed paths of length 3 :   1 → 2 → 1 → 2  2 → 1 → 2 → 1  2 → 1 → 2 → 3        Sample Input 3 6 2
0 0 0 0 0 0
0 0 1 0 0 0
0 0 0 0 0 0
0 0 0 0 1 0
0 0 0 0 0 1
0 0 0 0 0 0      Sample Output 3 1  G is drawn in the figure below:   There is one directed path of length 2 :   4 → 5 → 6        Sample Input 4 1 1
0      Sample Output 4 0       Sample Input 5 10 1000000000000000000
0 0 1 1 0 0 0 1 1 0
0 0 0 0 0 1 1 1 0 0
0 1 0 0 0 1 0 1 0 1
1 1 1 0 1 1 0 1 1 0
0 1 1 1 0 1 0 1 1 1
0 0 0 1 0 0 1 0 1 0
0 0 0 1 1 0 0 1 0 1
1 0 0 0 1 0 1 0 0 0
0 0 0 0 0 1 0 0 0 0
1 0 1 1 1 0 1 1 1 0      Sample Output 5 957538352  Be sure to print the count modulo 10^9 + 7 .

Solution (Python):
N, K = map(int, input().split())
A = [list(map(int, input().split())) for _ in range(N)]

MOD = 10**9+7

def dot(A, B):
    N = len(A)
    K = len(A[0])
    assert len(B) == K
    M = len(B[0])
    C = [[0] * M for _ in range(N)]
    for i in range(N):
        for j in range(M):
            for k in range(K):
                C[i][j] += (A[i][k] * B[k][j]) % MOD
                C[i][j] %= MOD
    return C



def mp(A, n):
    assert len(A) == len(A[0])
    B = [[0] * len(A) for _ in range(len(A))]
    for i in range(len(A)):
        B[i][i] = 1

    while n:
        n, i = divmod(n, 2)
        if i:
            B = dot(A, B)
        A = dot(A, A)

    return B


ans = sum([sum(r) % MOD for r in mp(A, K)]) % MOD
print(ans)

Solution (Python):
N, K = map(int, input().split())
A = [list(map(int, input().split())) for _ in range(N)]

MOD = 10**9+7

def dot(A, B):
    N = len(A)
    K = len(A[0])
    assert len(B) == K
    M = len(B[0])
    C = [[0] * M for _ in range(N)]
    for i in range(N):
        for j in range(M):
            C[i][j] = sum(A[i][k] * B[k][j] for k in range(K)) % MOD
    return C



def mp(A, n):
    assert len(A) == len(A[0])
    B = [[0] * len(A) for _ in range(len(A))]
    for i in range(len(A)):
        B[i][i] = 1

    while n:
        n, i = divmod(n, 2)
        if i:
            B = dot(A, B)
        A = dot(A, A)

    return B


ans = sum([sum(r) % MOD for r in mp(A, K)]) % MOD
print(ans)

Solution (Python):
N, K = map(int, input().split())
A = [list(map(int, input().split())) for _ in range(N)]

MOD = 10**9+7

def dot(A, B):
    N = len(A)
    K = len(A[0])
    assert len(B) == K
    M = len(B[0])
    C = [[0] * M for _ in range(N)]
    for i in range(N):
        for j in range(M):
            for k in range(K):
                C[i][j] += A[i][k] * B[k][j]
            C[i][j] %= MOD
    return C



def mp(A, n):
    assert len(A) == len(A[0])
    B = [[0] * len(A) for _ in range(len(A))]
    for i in range(len(A)):
        B[i][i] = 1

    while n:
        n, i = divmod(n, 2)
        if i:
            B = dot(A, B)
        A = dot(A, A)

    return B


ans = sum([sum(r) % MOD for r in mp(A, K)]) % MOD
print(ans)

Solution (Python):
n, k = map(int, input().split())
a = [list(map(int, input().split())) for _ in range(n)];

import numpy as np





graph = [[] for _ in range(n)]
for i, row in enumerate(a):
    for j, value in enumerate(row):
        if value:
            graph[i].append(j)






















import numpy as np
a = np.array(a, dtype=int)

M = 10**9 + 7
def powmod(a, n):
    if n == 0:
        return np.eye(np.alen(a), dtype=int)
    if n % 2:
        return (powmod(a, n - 1).dot(a)) % M
    r = powmod(a, n // 2)

    calc = []
    for row in r.T:
        calc.append(((row * r) % M).sum(axis=-1))
    return np.array(calc).T % M

print(powmod(a, k).sum() % M)

Solution (Python):
import sys
input = sys.stdin.readline

def matmul(a, b, mod):
    m = len(a)
    l = len(a[0])
    n = len(b[0])
    ret = [[0] * n for _ in range(m)]
    for i in range(m):
        for j in range(n):
            e = 0
            for k in range(l):
                e = (e + a[i][k] * b[k][j]) % mod
            ret[i][j] = e
    return ret

mod = 10**9+7
N, K = map(int, input().split())
a = [list(map(int, input().split())) for _ in range(N)]
ans = [[i == j for j in range(N)] for i in range(N)]
while K > 0:
    if K & 1 == 1:
        ans = matmul(ans, a, mod)
    a = matmul(a, a, mod)
    K >>= 1
print(sum(sum(row) % mod for row in ans) % mod)

Solution (Python):
import numpy as np
n, k = map(int, input().split())
a = [list(map(int , input().split())) for _ in range(n)]
a = np.array(a, dtype=object)
p = 10**9 + 7

def mat_pow(a, n, mod):
    res = np.eye(len(a), dtype=object)
    while(n > 0):
        if(n & 1 == 1):
            res = res.dot(a) % mod
        a = a.dot(a) % mod
        n = n >> 1
    return res
ak = mat_pow(a, k, p)

print(ak.sum() % p)
"""
ans = 0
for i in range(n):
    for j in range(n):
        ans += ak[i, j]
        ans %= p

print(ans)
"""

Solution (Python):
import numpy as np
n, k = map(int, input().split())
a = [list(map(int , input().split())) for _ in range(n)]
a = np.array(a, dtype=object)
p = 10**9 + 7

def mat_pow(a, n, mod):
    res = np.eye(len(a), dtype=object)
    while(n > 0):
        if(n & 1 == 1):
            res = np.dot(res, a) % mod
        a = np.dot(a, a) % mod
        n = n >> 1
    return res
ak = mat_pow(a, k, p)


ans = 0
for i in range(n):
    for j in range(n):
        ans += ak[i, j]
        ans %= p

print(ans)

Solution (Python):
n, k = map(int, input().split())
mod = 10**9+7
A = [list(map(int, input().split())) for _ in range(n)]


def mul(a, b):
    c = [[0]*len(b[0]) for i in range(len(a))]
    for i in range(len(a)):
        for k in range(len(b)):
            for j in range(len(b[0])):
                c[i][j] = (c[i][j] + a[i][k]*b[k][j])%mod

    return c


def pow(a, n):
    b = [[0]*len(a) for i in range(len(a))]
    for i in range(len(a)):
        b[i][i] = 1
    while n > 0:
        if n & 1 == 1:
            b = mul(a, b)
        a = mul(a, a)
        n = n>>1
    return b

B = pow(A, k)

ans = 0
for i in range(n):
    ans += sum(B[i])
    ans %= mod
print(ans)

Solution (Python):
def Mmul(a, b):  
    
    m = len(a)
    n = len(a[0])
    if len(b) != n:
        raise ValueError('列と行の数が一致しません')
    p = len(b[0])
    ans = [[0]*p for _ in range(m)]
    for i in range(m):
        for j in range(p):
            for k in range(n):
                ans[i][j] += a[i][k]*b[k][j] % MOD
    return ans


def Mfactorial(n, p):  
    if (p == 1):
        return n
    if (p % 2 == 0):
        t = Mfactorial(n, p // 2)
        return Mmul(t, t)
    return Mmul(Mfactorial(n, p - 1), n)


n, k = map(int, input().split())
a = [list(map(int, input().split())) for _ in range(n)]
MOD = 10**9+7
print(sum([sum(i) % MOD for i in Mfactorial(a, k)]) % MOD)

Solution (Python):

import sys, math, itertools, collections, bisect
input = lambda: sys.stdin.buffer.readline().rstrip().decode('utf-8')
inf = float('inf') ;mod = 10**9+7
mans = inf ;ans = 0 ;count = 0 ;pro = 1

n,k = map(int,input().split())
A = [list(map(int,input().split())) for i in range(n)]
def calc(A,B):
  res = [[0]*n for i in range(n)]
  for i in range(n):
    for j in range(n):
      for k in range(n):
        res[i][j] += A[i][k] * B[k][j]
      res[i][j] %= mod
  return res


def power(A,k):
  if k == 1:
    return A
  elif k % 2 == 0:
    tmp = power(A,k//2)
    return calc(tmp,tmp)
  elif k % 2 == 1:
    return calc(A,power(A,k-1))

result = power(A,k)
print(sum(map(sum,result))%mod)

Solution (Python):
N, K = map(int, input().split())
A = [list(map(int, input().split())) for i in range(N)]
MOD = 10**9 + 7


def dot(M1, M2, mod=10**9+7):
  n = len(M1)
  res = [[0]*n for _ in range(n)]
  for i in range(n):
    for j in range(n):
      for k in range(n):
        res[i][j] = (res[i][j] + M1[i][k] * M2[k][j]) % MOD
  
  return res


def power(M, k, mod=10**9+7):
  if k == 1:
    return M
  elif k%2 == 0:
    return power(dot(M, M, mod), k//2, mod)
  else:
    return dot(M, power(dot(M, M, mod), k//2), mod)


trans = power(A, K) 

ans = 0
for t in trans:
  ans += sum(t)%MOD
  ans %= MOD

print(ans)

Solution (Python):
mod = 10**9+7


def mat_mul(A, B):
    
    C = [[0]*N for _ in range(N)]
    for i in range(N):
        for j in range(N):
            C[i][j] = sum(A[i][k] * B[k][j] for k in range(N)) % mod
    return C


def matpow(A, n):
    B = [[0]*N for _ in range(N)]
    for i in range(N):
        B[i][i] = 1
    while n > 0:
        if n & 1:
            B = mat_mul(A, B)
        A = mat_mul(A, A)
        n //= 2
    return B


if __name__ == "__main__":
    N, K = map(int, input().split())
    A = [list(map(int, input().split())) for _ in range(N)]
    A = matpow(A, K)
    ans = 0
    for a in A:
        ans = (ans + sum(a)) % mod
    print(ans)

Solution (Python):
import sys
reader = (s.rstrip() for s in sys.stdin)
input = reader.__next__

mod = 10**9+7

def mul(a, b):
    n = len(a)
    c = [[0]*n for i in range(n)]
    for i in range(n):
        for j in range(n):
            for k in range(n):
                c[i][j] += (a[i][k]*b[k][j])%mod
    return c

def pow(a, n):
    b = [[0]*len(a) for i in range(len(a))]
    for i in range(len(a)):
        b[i][i] = 1
    while n:
        if n&1:
            b = mul(a, b)
        a = mul(a, a)
        n >>=1
    return b

n,k = map(int, input().split())
a = [list(map(int, input().split())) for i in range(n)]

b = pow(a, k)
ans = sum(sum(i) for i in b) % mod
print(ans)

Solution (Python):
N, K = map(int, input().split())
A = [list(map(int, input().split())) for i in range(N)]
MOD = 10**9 + 7


def dot(M1, M2, mod=10**9+7):
  n = len(M1)
  res = [[0]*n for _ in range(n)]
  for i in range(n):
    for j in range(n):
      for k in range(n):
        res[i][j] += M1[i][k] * M2[k][j] % MOD
  
  return res


def power(M, k, mod=10**9+7):
  if k == 1:
    return M
  elif k%2 == 0:
    return power(dot(M, M, mod), k//2)
  else:
    return dot(M, power(dot(M, M, mod), k//2), mod)


trans = power(A, K) 

ans = 0
for t in trans:
  ans += sum(t)%MOD
  ans %= MOD

print(ans)

Solution (Python):
import numpy as np
def main():
    mod=10**9+7
    n,k=map(int,input().split())
    A=np.matrix([list(map(int,input().split()))for _ in range(n)],dtype='O')
    D=A.copy()
    for i in range(k.bit_length()-1,0,-1):
        D=(D*D)%mod
        if k&1<<(i-1)>0:
            D=(D*A)%mod
    print(int(((D*np.ones(n,dtype='O')[:,np.newaxis])%mod).sum()%mod))
main()

Solution (Python):
def examA():
    N = I()
    H = LI()
    dp = [inf]*(N)
    dp[0] = 0
    for i in range(N-1):
        dp[i+1] = min(dp[i+1],dp[i]+abs(H[i+1]-H[i]))
        if i<N-2:
            dp[i + 2] = min(dp[i + 2], dp[i]+abs(H[i + 2] - H[i]))
    ans = dp[-1]
    print(ans)
    return

def examB():
    N, K = LI()
    H = LI()
    dp = [inf]*N
    dp[0] = 0
    for i in range(N):
        for k in range(1,K+1):
            if i+k>=N:
                break
            if dp[i+k]>dp[i]+abs(H[i+k]-H[i]):
                dp[i+k] = dp[i]+abs(H[i+k]-H[i])
    ans = dp[-1]
    print(ans)
    return

def examC():
    N = I()
    A = [0]*N; B = [0]*N; C = [0]*N
    for i in range(N):
        A[i],B[i],C[i] = LI()
    dp = [[0]*3 for _ in range(N+1)]
    for i in range(N):
        dp[i+1][0] = max(dp[i][1],dp[i][2]) + A[i]
        dp[i + 1][1] = max(dp[i][0], dp[i][2]) + B[i]
        dp[i + 1][2] = max(dp[i][1], dp[i][0]) + C[i]
    ans = max(dp[N])
    print(ans)
    return

def examD():
    ans = 0
    print(ans)
    return

def examE():
    ans = 0
    print(ans)
    return

def examF():
    
    S = SI(); N = len(S)
    T = SI(); M = len(T)
    dp = [[0] * (M + 1) for _ in range(N + 1)]
    for i in range(N):
        for j in range(M):
            if S[i] == T[j]:
                dp[i + 1][j + 1] = max(dp[i + 1][j + 1], dp[i][j] + 1)
            dp[i+1][j+1] = max(dp[i+1][j+1], dp[i+1][j], dp[i][j+1])

    i = N; j = M
    ans = []
    while i > 0 and j > 0:
        if dp[i][j] == dp[i-1][j]:
            i -= 1
        elif dp[i][j] == dp[i][j-1]:
            j -= 1
        else:
            ans.append(S[i-1])
            i -= 1
            j -= 1
    print("".join(ans[::-1]))
    return

def examG():
    
    
    def topological_sort(n,outs,ins):
        q = deque(v1 for v1 in range(n) if ins[v1] == 0)
        res = []
        while q:
            v1 = q.popleft()
            res.append(v1)
            for v2 in outs[v1]:
                ins[v2] -= 1
                if ins[v2] == 0:
                    q.append(v2)
        return res
    N, M = LI()
    V = [[]for _ in range(N)]
    start = [0]*N
    for _ in range(M):
        x,y = LI()
        x -= 1; y -= 1
        V[x].append(y)
        start[y] += 1
    Order = topological_sort(N, V, start)

    dp = [0]*N
    for i in Order:
        for v in V[i]:
            dp[v] = max(dp[v],dp[i]+1)
    ans = max(dp)
    print(ans)
    return

def examH():
    H, W = LI()
    A = [SI() for _ in range(H)]
    dp = [[0]*W for _ in range(H)]
    dp[0][0] = 1
    for i in range(H):
        for j in range(W):
            if A[i][j]=="
                continue
            dp[i][j] %=mod
            if i+1<H and dp[i+1][j]!="
                dp[i+1][j] +=dp[i][j]
            if j+1<W and dp[i][j+1]!="
                dp[i][j+1] +=dp[i][j]
    ans = dp[-1][-1]
    print(ans)
    return

def examI():
    N = I()
    P = LFI()
    dp = [[0]*(N+1) for _ in range(N+1)]
    dp[0][0] = 1
    for i in range(N):
        for j in range(i+1):
            dp[i+1][j+1] += dp[i][j]*P[i]
            dp[i+1][j] += dp[i][j]*(1-P[i])
    ans = 0
    for v in dp[N][(1+N//2):]:
        ans += v
    print(ans)
    return

def examJ():
    N = I()
    A = LI()
    dp = [[[0]*(N+1) for _ in range(N+1)]for _ in range(N+1)]
    S = Counter(A)
    print(S)
    dp[S[1]][S[2]][S[3]] = 0
    que = deque()
    que.append((S[1],S[2],S[3]))
    used = defaultdict(set)
    while(que):
        n1,n2,n3 = que.popleft()
        if used[(n1,n2,n3)]:
            continue
        used[(n1,n2,n3)] = True
        now = n1+n2+n3
        if n1>0:
            dp[n1-1][n2][n3] += dp[n1][n2][n3] + N/now * n1/now
            que.append((n1-1, n2, n3))
        if n2>0:
            dp[n1+1][n2-1][n3] += dp[n1][n2][n3] + N/now * n2/now
            que.append((n1+1,n2-1,n3))
        if n3>0:
            dp[n1][n2+1][n3-1] += dp[n1][n2][n3] + N/now * n3/now
            que.append((n1,n2+1,n3-1))

    ans = dp[0][0][0]
    print(ans)
    return

def examK():
    N, K = LI()
    A = LI()
    dp = [-1]*(K+1)
    for i in range(1,K+1):
        flag = False
        for a in A:
            if a>i:
                break
            if dp[i-a]==-1:
                flag = True
                break
        if flag:
            dp[i] = 1

    ans = dp[K]
    if ans==1:
        print("First")
    else:
        print("Second")
    return

def examM():
    N, K = LI()
    A = LI()
    dp = [[0]*(K+1) for _ in range(N+1)]
    dp[0][0] = 1
    for i in range(N):
        cum = [0]*(K+1)
        cum[0] = dp[i][0]
        for c in range(K):
            cum[c+1] = cum[c] + dp[i][c+1]
        for j in range(K+1):
            if j<=A[i]:
                dp[i+1][j] = cum[j]
            else:
                dp[i+1][j] = cum[j] - cum[j-A[i]-1]
            dp[i+1][j] %=mod

    ans = dp[N][K]%mod
    print(ans)
    return

def examO():
    N = I()
    A = [LI() for _ in range(N)]
    loop = 2**N
    dp = [[0]*loop for _ in range(N+1)]
    dp[0][0] = 1
    for mask in range(loop):
        cur = 0
        for j in range(N):
            if mask&(1<<j)==(1<<j):
                cur +=1
        for i in range(N):
            if mask&(1<<i)==0 and A[i][cur]==1:
                dp[cur+1][mask|(1<<i)] += dp[cur][mask]
                dp[cur + 1][mask | (1 << i)] %=mod

    ans = dp[N][loop-1]
    print(ans)
    return

def examP():
    def dfs_dp(n,s,fords,visited):
        W = 1; B = 1
        visited[s] = True
        for v in fords[s]:
            if visited[v]:
                continue
            cur = dfs_dp(n, v, fords,visited)
            W *= sum(cur)
            B *= cur[0]
            W %= mod
            B %= mod

        return (W,B)
    N = I()
    V = [[]for _ in range(N)]
    for _ in range(N-1):
        x,y = LI()
        x -=1; y -=1
        V[x].append(y)
        V[y].append(x)
    visited = [False]*N
    ans = dfs_dp(N,0,V,visited)
    print(sum(ans)%mod)
    return

def examQ():
    N = I()
    H = LI()
    maxH = max(H)
    A = LI()
    dp = [0](maxH+1)
    cur = 0
    for i in range(N):
        dp[B[i]] = cur+dp[A[i]]
    return

def examR():
    def list_dot(A1,A2):
        res = [[0]*len(A2) for _ in range(len(A1))]
        for i in range(len(A1)):
            for j in range(len(A2)):
                cur = 0
                for k in range(len(A1[0])):
                    cur += A1[i][k]*A2[k][j]
                res[i][j] = cur%mod
        return res
    N, K = LI()
    A = [LI() for _ in range(N)]
    now = [[1]*N]
    while(K>0):
        if K%2==1:
            now = list_dot(now,A)
        K //=2
        A = list_dot(A,A)
    ans = 0
    for i in now:
        for k in i:
            ans += k
            ans %= mod
    print(ans)
    return

import sys,copy,bisect,itertools,heapq,math
from heapq import heappop,heappush,heapify
from collections import Counter,defaultdict,deque
def I(): return int(sys.stdin.readline())
def LI(): return list(map(int,sys.stdin.readline().split()))
def LFI(): return list(map(float,sys.stdin.readline().split()))
def LSI(): return list(map(str,sys.stdin.readline().split()))
def LS(): return sys.stdin.readline().split()
def SI(): return sys.stdin.readline().strip()
global mod,mod2,inf,alphabet
mod = 10**9 + 7
mod2 = 998244353
inf = 10**18
alphabet = [chr(ord('a') + i) for i in range(26)]

sys.setrecursionlimit(10**6)

if __name__ == '__main__':
    examR()

"""

"""

Solution (Python):
mod = 10**9 +7

n, k = map(int, input().split())
a = [list(map(int, input().split()))  for _ in range(n)]

def matmul_nn(a,b):
    
    n = len(a)
    ret = [[0]*n for _ in range(n)]
    for i in range(n):
        for j in range(n):
            temp = 0
            for k in range(n):
                temp += a[i][k]*b[k][j]
                temp = temp % mod
            ret[i][j] = temp
    return ret

k_bin = str(bin(k))[2:]
k_bin = k_bin[::-1]
k_len = len(k_bin)

a_2 = [0]* k_len
a_2[0] = a

for i in range(1,k_len):
    a_2[i] = matmul_nn(a_2[i-1],a_2[i-1])

a_k = a_2[k_len-1]
for i in range(k_len-1):
    if( k_bin[i] == '1'):
        a_k = matmul_nn(a_k, a_2[i])

ans = 0
for i in range(n):
    ans += sum(a_k[i])
    ans = ans % mod

print(ans)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
X,M,S=matrix(A,"O").reshape(N,N),10**9+7,1
while K:S*=~K%2or X;K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
X,M,S=matrix(A,"O").reshape(N,N),10**9+7,1
while K:S*=~K%2or X;K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
import numpy as np
N, K = map(int, input().split())
M = 10**9+7
xs = []
for _ in range(N):
    xs.append(list(map(int, input().split())))
X = np.array(xs, dtype="object")
S = np.ones(N, dtype="object")

for _ in range(60):

    if K%2:
        S = np.dot(X, S)%M
    K = K//2
    X = np.dot(X, X)%M
S = np.array(S, dtype="int64")
r = 0
for s in S:
    r = (r + s)%M
print(r)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
X=matrix(A,"O").reshape(N,N)
M=10**9+7
S=1
while K:S*=~K%2or X;K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
M,S,X=10**9+7,1,matrix(A,"O").reshape(N,N)
while K:S*=[1,X][K%2];K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
X=matrix(A,"O").reshape(N,N)
M=10**9+7
S=1
while K:S*=~K%2or X;K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
X=matrix(A,'O').reshape(N,N)
M=10**9+7
S=1
while K:S*=~K%2or X;K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
X=matrix(A,"O").reshape(N,N)
M=10**9+7
S=1
while K:S*=1-K%2or X;K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
from numpy import*
N,K,*A=map(int,open(0).read().split())
M,S,X=10**9+7,1,matrix(A,"O").reshape(N,N)
while K:S*=[1,X][K%2];K//=2;X=X*X%M
print(sum(S)%M)

Solution (Python):
import numpy
N,K,*A=map(int,open(0).read().split())
M,S,X=10**9+7,1,numpy.matrix(A,"O").reshape(N,N)
while K:S*=[1,X][K&1];K//=2;X=X*X%M
print(S.sum()%M)

Solution (Python):
import numpy as u
N,K,*A=map(int,open(0).read().split())
M,S,X=10**9+7,1,u.matrix(A,"O").reshape(N,N)
exec("S*=[1,X][K&1];K//=2;X=X*X%M;"*60)
print(S.sum()%M)

Solution (Python):
import numpy as np
N, K = map(int, input().split())
M = 10**9+7
xs = []
for _ in range(N):
    xs.append(list(map(int, input().split())))
X = np.array(xs, dtype="object")
S = np.ones(N, dtype="object")

for _ in range(60):

    if K%2:
        S = np.dot(X, S)%M
    K = K//2
    X = np.dot(X, X)%M
S = np.array(S, dtype="int64")
r = 0
for s in S:
    r = (r + s)%M
print(r)

Solution (Python):
N,K = map(int, input().split())
A = [[int(a) for a in input().split()] for _ in range(N)]

mod = 10**9+7

import numpy as np

a = np.array(A, dtype=object)
Mat = np.eye(N, dtype=object)

while K:
    if K%2:
        Mat = np.dot(Mat, a)%mod
    a = np.dot(a, a)%mod
    K //= 2
    
ans = 0
for i in range(N):
    for j in range(N):
        ans += Mat[i, j]
        ans %= mod
        
print(ans)

Solution (Python):
mod = 10**9+7
n, K = map(int, input().split())
a = [[int(x) for x in input().split()] for _ in range(n)]

A = [a]

i = 1
while 2**i <= K:
    mtrix = []
    for j in range(n):
        raw = []
        for l in range(n):
            tmp = 0
            for k in range(n):
                tmp += A[-1][j][k] * A[-1][k][l]
            raw.append(tmp % mod)
        mtrix.append(raw)
    A.append(mtrix)
    i += 1

bit_max = i

ans = [[[0]*n for _ in range(n)]]
for i in range(n):
    for j in range(n):
        if i == j:
            ans[-1][i][j] = 1

for i in range(bit_max):
    if (K >> i) & 1:
        mtrix = []
        for j in range(n):
            raw = []
            for l in range(n):
                tmp = 0
                for k in range(n):
                    tmp += ans[-1][j][k] * A[i][k][l]
                raw.append(tmp % mod)
            mtrix.append(raw)
        ans.append(mtrix)

ans = sum(map(lambda x: sum(x) % mod, ans[-1]))

print(ans % mod)

Solution (Python):

import sys
INF = float("inf")
import numpy as np

MOD = 1000000007  


def solve(N: int, K: int, a: "List[List[int]]"):
    a = np.asarray(a, dtype=object)
    ans = np.identity(N, dtype=object)
    while K > 0:
        if K & 1:
            ans = np.mod(np.dot(ans, a), MOD)
        a = np.mod(np.dot(a, a), MOD)
        K >>= 1
    print(ans.sum() % MOD)

    return


def main():

    def iterate_tokens():
        for line in sys.stdin:
            for word in line.split():
                yield word
    tokens = iterate_tokens()
    N = int(next(tokens))  
    K = int(next(tokens))  
    a = [[int(next(tokens)) for _ in range(N)]
         for _ in range(N)]  
    solve(N, K, a)


if __name__ == '__main__':
    main()

Solution (Python):
import numpy as np


def main():
    N, K = map(int, input().split())
    mod = 10 ** 9 + 7
    A = [[int(x) for x in input().split()] for _ in range(N)]
    A = np.asarray(A, dtype=object)
    ans = np.identity(N, dtype=object)
    while K > 0:
        if K & 1:
            ans = np.mod(np.dot(ans, A), mod)
        A = np.mod(np.dot(A, A), mod)
        K >>= 1
    print(ans.sum() % mod)


if __name__ == '__main__':
    main()

Solution (Python):
n, k = map(int, input().split())
matrix = [list(map(int, input().split())) for i in range(n)]
MOD = 10**9 + 7

def _mul(A, B, MOD):
    C = [[0] * len(B[0]) for i in range(len(A))]
    for i in range(len(A)):
        for k in range(len(B)):
            for j in range(len(B[0])):
                C[i][j] = (C[i][j] + A[i][k]*B[k][j]) % MOD
    return C

def pow_matrix(A, n, MOD):
    """A**nをダブリングによって求める"""
    B = [[0] * len(A) for i in range(len(A))]
    for i in range(len(A)):
        B[i][i] = 1
    while n > 0:
        if n & 1:
            B = _mul(A, B, MOD)
        A = _mul(A, A, MOD)
        n  = n // 2
    return B 
          
ans_matrix = pow_matrix(matrix, k, MOD)
ans = 0
for i in ans_matrix:
    ans += sum(i)
    ans %= MOD
print(ans)

Solution (Python):
def mmult(A, B):
    global mod
    n, m, l = len(A), len(B), len(B[0])
    ret = [[0]*l for _ in range(n)]
    for i in range(n):
        for j in range(m):
            for k in range(l):
                ret[i][k] = (ret[i][k]+A[i][j]*B[j][k])%mod
    return ret
def mpow(A, n):
    if n == 0: return [[1 if i==j else 0 for j in range(len(A))] for i in range(len(A))]
    if n % 2: return mmult(mpow(A, n-1), A)
    return mpow(mmult(A, A), n//2)

mod = 10**9+7
N, K = map(int, input().split())
X = [[int(a) for a in input().split()] for _ in range(N)]
print(sum([sum(x) for x in mpow(X, K)]) % mod)

Solution (Python):
import numpy as np


def main():
    MOD = int(1e9 + 7)

    n, k = map(int, input().split())
    can = [[int(c) for c in input().split()] for _ in range(n)]
    can = np.array(can, dtype=object)

    answer = np.eye(n, dtype=object)
    while k:
        if k % 2:
            answer = np.dot(answer, can) % MOD
        can = np.dot(can, can) % MOD
        k //= 2

    total = 0
    for i in range(n):
        for j in range(n):
            total += answer[i, j]
            if total >= MOD:
                total -= MOD

    print(total)


if __name__ == "__main__":
    main()

Solution (Python):
import sys
import numpy as np


def main():
    MOD = int(1e9 + 7)

    n, k = map(int, input().split())
    can = [[int(c) for c in input().split()] for _ in range(n)]
    

    can = np.array(can, 'O')

    
    
    
    
    

    

    answer = np.eye(n, dtype=object)
    while k:
        if k % 2:
            answer = np.dot(answer, can) % MOD
        can = np.dot(can, can) % MOD
        k //= 2

    total = 0
    for i in range(n):
        for j in range(n):
            total += answer[i,j]
            if total >= MOD:
                total -= MOD

    print(total)


if __name__ == "__main__":
    main()

Solution (Python):
import sys
import numpy as np


def main():
    MOD = int(1e9 + 7)

    n, k = map(int, input().split())
    can = [[int(c) for c in input().split()] for _ in range(n)]
    can = np.matrix(can, 'O')

    
    
    

    
    
    
    
    

    

    answer = np.eye(n, dtype=int)
    while k:
        if k % 2:
            answer = answer * can % MOD
        can = can * can % MOD
        k //= 2

    total = 0
    for i in range(n):
        for j in range(n):
            total += answer[i,j]
            if total >= MOD:
                total -= MOD

    print(total)


if __name__ == "__main__":
    main()

Solution (Python):
import sys
import numpy as np


def main():
    MOD = int(1e9 + 7)

    n, k = map(int, input().split())
    can = [[int(c) for c in input().split()] for _ in range(n)]
    can = np.matrix(can, 'O')
    
    
    
    

    b = np.eye(n, n, 0, 'O')
    while k:
        b = [b, b*can % MOD][k & 1]
        can = can*can % MOD
        k >>= 1

    print(b.sum() % MOD)

    
    
    
    
    
    

    
    
    
    
    
    

    


if __name__ == "__main__":
    main()

Solution (Python):
import sys

n, k = [int(i) for i in sys.stdin.readline().split()]
a_ls = []
MOD = 10**9+7
for i in range(n):
    a_ls.append([int(i) for i in sys.stdin.readline().split()])


def matmul(a, b):
    x = len(a)
    res = [[0 for i in range(x)]for j in range(x)]
    for i in range(x):
        for j in range(x):
            for k in range(x):
                res[i][j] += (a[i][k] * b[k][j])
            res[i][j] %= MOD
    return res

def power(a, n):
    x = len(a)
    res = a
    for i in range(n - 1):
        res = matmul(res, a)
    return res

def plus(a, b):
    x = len(a)
    res = [[0 for i in range(x)]for j in range(x)]
    for i in range(x):
        for j in range(x):
            res[i][j] = (a[i][j] + b[i][j]) % MOD
    return res

k_2 = bin(k)[2:][::-1]
base = a_ls
res = [[1 if i == j else 0 for i in range(n)]for j in range(n)]
for ind, _k in enumerate(k_2):
    if _k == "1":
        res = matmul(res, base)
    base = power(base, 2)
print(sum([sum(i) % MOD for i in res]) % MOD)

Solution (Python):
mod = 10**9+7


def matmul(A, B):
    C = [[0] * len(B[0]) for _ in range(len(A))]
    for i in range(len(A)):
        for j in range(len(B[0])):
            for k in range(len(B)):
                C[i][j] = (C[i][j] + A[i][k] * B[k][j]) % mod
    return C


def matpow(A, p):
    n = len(A)
    B = [[0] * n for _ in range(n)]
    for i in range(n):
        B[i][i] = 1
    while p > 0:
        if p & 1:
            B = matmul(B, A)
        A = matmul(A, A)
        p >>= 1
    return B


N, K = map(int, input().split())
adj = []
for _ in range(N):
    adj.append(list(map(int, input().split())))
mat = []
for col in range(N):
    tmp = []
    for row in range(N):
        tmp.append(adj[row][col])
    mat.append(tmp)

vec = [[1] * N]
ans = matmul(vec, matpow(mat, K))
print(sum(ans[0]) % mod)

Solution (Python):
import numpy as np
mod = 10**9 + 7

N, K = map(int,input().split())
A = [list(map(int,input().split())) for _ in range(N)]
A = np.array(A, dtype=object)
B = np.eye(N, dtype=object)

while K > 0:
    if K & 1:
        B = np.dot(B, A) % mod
    A = np.dot(A, A) % mod
    K >>= 1

ans = B.sum() % mod
print(ans)

Solution (Python):
p = 10**9+7
def pow(X,k):
    if k==0:
        return [[1 if i==j else 0 for j in range(N)] for i in range(N)]
    if k==1:
        return X
    if k%2==0:
        Y = pow(X,k//2)
        Z = [[0 for _ in range(N)] for _ in range(N)]
        for i in range(N):
            for j in range(N):
                for m in range(N):
                    Z[i][j] = (Z[i][j]+Y[i][m]*Y[m][j])%p
        return Z
    else:
        Y = pow(X,(k-1)//2)
        Z = [[0 for _ in range(N)] for _ in range(N)]
        for i in range(N):
            for j in range(N):
                for m in range(N):
                    Z[i][j] = (Z[i][j]+Y[i][m]*Y[m][j])%p
        Z1 = [[0 for _ in range(N)] for _ in range(N)]
        for i in range(N):
            for j in range(N):
                for m in range(N):
                    Z1[i][j] = (Z1[i][j]+Z[i][m]*X[m][j])%p
        return Z1
N,K = map(int,input().split())
A = [list(map(int,input().split())) for _ in range(N)]
X = pow(A,K)
cnt = 0
for i in range(N):
    for j in range(N):
        cnt = (cnt+X[i][j])%p
print(int(cnt))