Buggy python submission (verdict: Time Limit Exceeded):
import numpy as np
N = int(input())
A = list(map(int, input().split()))
R = [0] * (N + 1)
P = [1] * (N + 1)
#R = np.zeros(N + 1, dtype=int)
#P = np.ones(N + 1, dtype=int)
mod = 10 ** 9 + 7

def rev(a):
    return pow(a, mod - 2, mod)
    
for n in range(N + 1):
    R[n] = (R[n - 1] + rev(n)) % mod
    if n != 0:
        P[n] = P[n - 1] * n
    
C = [(R[j] + R[N - j + 1] - 1) % mod for j in range(1, N + 1)]

print(sum([a * c for a, c in zip(C, A)]) * P[N] % mod)


Fixed submission:
import numpy as np
N = int(input())
A = list(map(int, input().split()))
R = np.zeros(N + 1, dtype="int64")
P = np.ones(N + 1, dtype="int64")
mod = 10 ** 9 + 7

def rev(a):
    return pow(a, mod - 2, mod)
    
for n in range(N + 1):
    r = rev(n)
    R[n] = (R[n - 1] + r) % mod
    if n != 0:
        P[n] = (P[n - 1] * n) % mod

C = [(R[j] + R[N - j + 1] - 1) % mod for j in range(1, N + 1)]

res = 0
for i in range(N):
    res = (res + (A[i] * C[i]) % mod) % mod

print(res * P[N] % mod)
#print((sum([a * c % mod for a, c in zip(C, A)]) * P[N]) % mod)

Buggy python submission (verdict: Time Limit Exceeded):
N = int(input())
A = list(map(int, input().split()))
R = [0] * (N + 1)
P = [1] * (N + 1)
mod = 10 ** 9 + 7

def power(x, y):
    if y == 0:
        return 1
    elif y == 1: 
        return x % mod
    elif y % 2 == 0: 
        return power(x, y//2)**2 % mod
    else: 
        return power(x, y//2)**2 * x % mod

def rev(a):
    return power(a, mod - 2)
    
for n in range(N + 1):
    R[n] = R[n - 1] + rev(n)
    if n != 0:
        P[n] = P[n - 1] * n

res = 0
for i in range(N):
    res += ((R[i + 1] + R[N - i] - 1) % mod * A[i] * P[N]) % mod

print(res % mod)


Fixed submission:
import numpy as np
N = int(input())
A = list(map(int, input().split()))
R = np.zeros(N + 1, dtype="int64")
P = np.ones(N + 1, dtype="int64")
mod = 10 ** 9 + 7

def rev(a):
    return pow(a, mod - 2, mod)
    
for n in range(N + 1):
    r = rev(n)
    R[n] = (R[n - 1] + r) % mod
    if n != 0:
        P[n] = (P[n - 1] * n) % mod

C = [(R[j] + R[N - j + 1] - 1) % mod for j in range(1, N + 1)]

res = 0
for i in range(N):
    res = (res + (A[i] * C[i]) % mod) % mod

print(res * P[N] % mod)
#print((sum([a * c % mod for a, c in zip(C, A)]) * P[N]) % mod)

Buggy python submission (verdict: Time Limit Exceeded):
N = int(input())
A = list(map(int, input().split()))
R = [0] * (N + 1)
P = [1] * (N + 1)
mod = 10 ** 9 + 7

def power(x, y):
    if y == 0:
        return 1
    elif y == 1: 
        return x % mod
    elif y % 2 == 0: 
        return power(x, y//2)**2 % mod
    else: 
        return power(x, y//2)**2 * x % mod

def rev(a):
    return power(a, mod - 2)
    
for n in range(N + 1):
    R[n] = R[n - 1] + rev(n)
    if n != 0:
        P[n] = P[n - 1] * n

res = 0
for i in range(N):
    res += ((R[i + 1] + R[N - i] - 1) * A[i] * P[N]) % mod

print(res % mod)


Fixed submission:
import numpy as np
N = int(input())
A = list(map(int, input().split()))
R = np.zeros(N + 1, dtype="int64")
P = np.ones(N + 1, dtype="int64")
mod = 10 ** 9 + 7

def rev(a):
    return pow(a, mod - 2, mod)
    
for n in range(N + 1):
    r = rev(n)
    R[n] = (R[n - 1] + r) % mod
    if n != 0:
        P[n] = (P[n - 1] * n) % mod

C = [(R[j] + R[N - j + 1] - 1) % mod for j in range(1, N + 1)]

res = 0
for i in range(N):
    res = (res + (A[i] * C[i]) % mod) % mod

print(res * P[N] % mod)
#print((sum([a * c % mod for a, c in zip(C, A)]) * P[N]) % mod)

Buggy python submission (verdict: Time Limit Exceeded):
N = int(input())
A = list(map(int, input().split()))
R = [0] * (N + 1)
P = [1] * (N + 1)
mod = 10 ** 9 + 7

def power(x, y):
    if y == 0:
        return 1
    elif y == 1: 
        return x % mod
    elif y % 2 == 0: 
        return power(x, y//2)**2 % mod
    else: 
        return power(x, y//2)**2 * x % mod

def rev(a):
    return power(a, mod - 2)
    
for n in range(N + 1):
    R[n] = R[n - 1] + rev(n)
    if n != 0:
        P[n] = P[n - 1] * n
    
C = [(R[j] + R[N - j + 1] - 1) % mod for j in range(1, N + 1)]

print(sum([a * c for a, c in zip(C, A)]) * P[N] % mod)

Fixed submission:
import numpy as np
N = int(input())
A = list(map(int, input().split()))
R = np.zeros(N + 1, dtype="int64")
P = np.ones(N + 1, dtype="int64")
mod = 10 ** 9 + 7

def rev(a):
    return pow(a, mod - 2, mod)
    
for n in range(N + 1):
    r = rev(n)
    R[n] = (R[n - 1] + r) % mod
    if n != 0:
        P[n] = (P[n - 1] * n) % mod

C = [(R[j] + R[N - j + 1] - 1) % mod for j in range(1, N + 1)]

res = 0
for i in range(N):
    res = (res + (A[i] * C[i]) % mod) % mod

print(res * P[N] % mod)
#print((sum([a * c % mod for a, c in zip(C, A)]) * P[N]) % mod)