Problem:


Score : 500 points

Problem Statement

We have a grid of squares with N rows and M columns. Let (i, j) denote the square at the i-th row from the top and j-th column from the left. We will choose K of the squares and put a piece on each of them.

If we place the K pieces on squares (x_1, y_1), (x_2, y_2), ..., and (x_K, y_K), the cost of this arrangement is computed as:

\sum_{i=1}^{K-1} \sum_{j=i+1}^K (|x_i - x_j| + |y_i - y_j|)

Find the sum of the costs of all possible arrangements of the pieces. Since this value can be tremendous, print it modulo 10^9+7.

We consider two arrangements of the pieces different if and only if there is a square that contains a piece in one of the arrangements but not in the other.

Constraints

2 \leq N \times M \leq 2 \times 10^5

2 \leq K \leq N \times M

All values in input are integers.

Input

Input is given from Standard Input in the following format:

N M K

Output

Print the sum of the costs of all possible arrangements of the pieces, modulo 10^9+7.

Sample Input 1

2 2 2

Sample Output 1

8

There are six possible arrangements of the pieces, as follows:

((1,1),(1,2)), with the cost |1-1|+|1-2| = 1

((1,1),(2,1)), with the cost |1-2|+|1-1| = 1

((1,1),(2,2)), with the cost |1-2|+|1-2| = 2

((1,2),(2,1)), with the cost |1-2|+|2-1| = 2

((1,2),(2,2)), with the cost |1-2|+|2-2| = 1

((2,1),(2,2)), with the cost |2-2|+|1-2| = 1

The sum of these costs is 8.

Sample Input 2

4 5 4

Sample Output 2

87210

Sample Input 3

100 100 5000

Sample Output 3

817260251

Be sure to print the sum modulo 10^9+7.



Sample input:
2 2 2


Sample output:
8

Submission status: Accepted

Solution code:
#!usr/bin/env python3
from collections import defaultdict
from collections import deque
from heapq import heappush, heappop
import sys
import math
import bisect
import random
import itertools
sys.setrecursionlimit(10**5)
stdin = sys.stdin
bisect_left = bisect.bisect_left
bisect_right = bisect.bisect_right
def LI(): return list(map(int, stdin.readline().split()))
def LF(): return list(map(float, stdin.readline().split()))
def LI_(): return list(map(lambda x: int(x)-1, stdin.readline().split()))
def II(): return int(stdin.readline())
def IF(): return float(stdin.readline())
def LS(): return list(map(list, stdin.readline().split()))
def S(): return list(stdin.readline().rstrip())
def IR(n): return [II() for _ in range(n)]
def LIR(n): return [LI() for _ in range(n)]
def FR(n): return [IF() for _ in range(n)]
def LFR(n): return [LI() for _ in range(n)]
def LIR_(n): return [LI_() for _ in range(n)]
def SR(n): return [S() for _ in range(n)]
def LSR(n): return [LS() for _ in range(n)]
mod = 1000000007
inf = float('INF')

#A
def A():
    a, b = LI()
    if a <= 5:
        print(0)
    elif 6 <= a <= 12:
        print(b // 2)
    else:
        print(b)
    return

#B
def B():
    r, d, x = LI()
    for i in range(10):
        print(x * r - d)
        x = x * r - d
    return

#C
def C():
    n, m = LI()
    imo = [0 for i in range(n+1)]
    for _ in range(m):
        l, r = LI_()
        imo[l] += 1
        imo[r + 1] -= 1
    for i in range(1,n):
        imo[i] += imo[i - 1]
    ans = 0
    #print(imo)
    for i in imo:
        if i == m:
            ans += 1
    print(ans)
    return

#D
def D():
    n, m = LI()
    a = LI()
    bc = LIR(m)
    bc.sort(key=lambda x: x[1], reverse=True)
    a.sort()
    index = 0
    suma = 0
    for b, c in bc:
        if index + b >= n:
            if a[n - 1] <= c:
                suma += (n - index) * c
            else: 
                br = bisect_left(a, c)
                if br > index:
                    suma += (br - index) * c
                    index = br
                    break
                else:
                    break
            index = n
            break
        if a[index + b - 1] < c:
            suma += b * c
            index += b
        else:
            if a[index] < c:
                bl = bisect_left(a, c)
                suma += (bl - index) * c
                index = bl
            else:
                break
    suma += sum(a[index:])
    print(suma)

    return

#E
def E():
    n, m, k = LI()
    def cmb(n, r, mod):
        if ( r<0 or r>n ):
            return 0
        r = min(r, n-r)
        return g1[n] * g2[r] * g2[n-r] % mod

    g1 = [1, 1] 
    g2 = [1, 1]
    inverse = [0, 1] 

    for i in range( 2, n*m + 1 ):
        g1.append( ( g1[-1] * i ) % mod )
        inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
        g2.append((g2[-1] * inverse[-1]) % mod)
        
    ans = 0
    for d in range(1,n+1):
        ans += d*cmb(n * m - 2, k - 2, mod) * (n - d) * m ** 2 
    for d in range(1,m+1):
        ans += d*cmb(n * m - 2, k - 2, mod) * (m - d) * n ** 2
    print(ans%mod)
    

    
    return

#F
def F():
    return

#G
def G():
    return

#H
def H():
    return

#Solve
if __name__ == '__main__':
    E()

Submission status: Accepted

Solution code:
def factorial_mod(n, mod):
    a = 1
    for i in range(1, n + 1):
        a *= i
        a %= mod
    return a


def comb_mod(n, k, mod):
    fact_n = factorial_mod(n, mod)
    fact_k = factorial_mod(k, mod)
    fact_n_k = factorial_mod(n - k, mod)
    return (fact_n * pow(fact_k, mod - 2, mod) * pow(fact_n_k, mod - 2, mod)) % mod


n, m, k = map(int, input().split())
s = 0
mod = 10 ** 9 + 7
for x in range(m):
    s += (m - x - 1) * (m - x) // 2 * n **2
for y in range(n):
    s += (n - y - 1) * (n - y) // 2 * m ** 2


print(s * comb_mod(n * m - 2, k - 2, mod) % mod)

Submission status: Accepted

Solution code:
import math
import sys
sys.setrecursionlimit(10 ** 6 + 10)
from collections import defaultdict

from math import factorial

def combinations_count(n, r):
    return factorial(n) // (factorial(n - r) * factorial(r))

def permutations_count(n, r):
    return math.factorial(n) // math.factorial(n - r)



def main():
    mod = 10 ** 9 + 7
    high, width, num = map(int, input().split())
    toori = combinations_count(high * width - 2, num - 2)

    kyori = high + width - 2

    ans = 0
    for now_kyori in range(1, kyori + 1):
        kari = 0
        if now_kyori < high:
            kari += width * width * (high - now_kyori)
        if now_kyori < width:
            kari += high * high * (width - now_kyori)
        ans += now_kyori * kari

        if now_kyori >= high and now_kyori >= width:
            break

    ans = ans * toori % mod

    print(ans)



if __name__ == '__main__':
    main()

Submission status: Accepted

Solution code:
import sys
import bisect
from collections import defaultdict

ans=0
#n=int(input())
n,m,k=map(int,input().split())
#a=[(1,int(i)) for i in input().split()]

# http://drken1215.hatenablog.com/entry/2018/06/08/210000
max_value = n*m+2
mod = 10 ** 9 + 7

fac = [0 for _ in range(max_value)]
finv = [0 for _ in range(max_value)]
inv = [0 for _ in range(max_value)]

def com_init():
    fac[0] = 1
    fac[1] = 1
    finv[0] = 1
    finv[1] = 1
    inv[1] = 1

    for i in range(2, max_value):
        fac[i] = fac[i - 1] * i % mod
        inv[i] = mod - inv[mod % i] * (mod // i) % mod
        finv[i] = finv[i - 1] * inv[i] % mod

def com(n, k):
    if n < k:
        return 0
    if n < 0 or k < 0:
        return 0
    return fac[n] * (finv[k] * finv[n - k] % mod) % mod

def dist(a,b):
    if a<=0 or b<=0:
        return 0
    
    a,b=min(a,b),max(a,b)
    a,b=a%mod,b%mod
    '''
    rt=((a*(a+1))%mod)*(a+b)//2
    rt+=((a*(a+b))%mod)*(b-a-1) if a!=b else 0
    rt+= - a*a if a==b else 0
    '''
    rt=((b*a*(a+b))//2)%mod
    #rt//=2
    return rt%mod

com_init()
inv2=500000004
#print(ans,dist(2,2),dist(1,0))
d=0
for i in range(n):
    for j in range(m):
        d+=dist(j,i+1)+dist(i,m-j)+dist(j+1,n-i-1)+dist(n-i,m-j-1)
        d%=mod
ans=d*com(n*m-2,k-2)
ans*=inv2

print(ans%mod)

Submission status: Accepted

Solution code:
N, M, K = map(int, input().split())

MOD = 10**9 + 7
R = 1
for i in range(K - 2):
    R = (R * (N * M - 2 - i) * pow(i + 1, MOD - 2, MOD)) % MOD

ans = 0
for i in range(1, N + 1):
    for j in range(1, M + 1):
        u = (N - i + 1)
        d = i
        l = j
        r = (M - j + 1)
        D = 0
        D += u * l * (u + l - 2) // 2
        D += u * r * (u + r - 2) // 2
        D += d * l * (d + l - 2) // 2
        D += d * r * (d + r - 2) // 2
        D -= u * (u - 1) // 2
        D -= d * (d - 1) // 2
        D -= l * (l - 1) // 2
        D -= r * (r - 1) // 2
        ans += D

ans = (ans * R * pow(2, MOD - 2, MOD)) % MOD
print(ans)

Submission status: Accepted

Solution code:
N, M, K = map(int, input().split())

MOD = 10**9 + 7
L = N*M

fact = [1]*(L+1)
rfact = [1]*(L+1)
r = 1
for i in range(1, L+1):
  fact[i] = r = r * i % MOD
rfact[L] = r = pow(fact[L], MOD-2, MOD)
for i in range(L, 0, -1):
  rfact[i-1] = r = r * i % MOD

L = max(N, M)
D = [0]*(L+1)
for i in range(L):
    D[i+1] = D[i] + (i+1)

ans = 0
for i in range(N):
    for j in range(M):
        ans += D[j] * (i+1) + D[M-1-j] * i + D[i] * M
    ans %= MOD
print(ans * fact[N*M-2] * rfact[N*M-K] * rfact[K-2] % MOD)

Submission status: Accepted

Solution code:
def cmb(n, r, mod):
    if ( r<0 or r>n ):
        return 0
    r = min(r, n-r)
    return g1[n] * g2[r] * g2[n-r] % mod

A=2*10**5
mod = 10**9+7
g1 = [1, 1]
g2 = [1, 1]
inverse = [0, 1]

for i in range( 2, A + 1 ):
    g1.append( ( g1[-1] * i ) % mod )
    inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
    g2.append( (g2[-1] * inverse[-1]) % mod )

N,M,K=map(int,input().split())
keisuu=cmb(N*M-2,K-2,mod)

sum_=0
for i in range(N):
    a=min(abs(i),abs(N-i-1))
    b=max(abs(i),abs(N-i-1))
    sum_+=(M**2)*((a*(a+1)//2)+(b*(b+1)//2))
for i in range(M):
    a=min(abs(i),abs(M-i-1))
    b=max(abs(i),abs(M-i-1))
    sum_+=(N**2)*((a*(a+1)//2)+(b*(b+1)//2))
print((keisuu * (sum_//2))%mod)

Submission status: Accepted

Solution code:
import math
import sys

sys.setrecursionlimit(10000)
INF = float('inf')


def inv(x):
    return pow(x, MOD - 2, MOD)


def comb(n, r, mod):
    # https://fushime2.hatenablog.com/entry/2016/11/12/155543
    return (math.factorial(n) * inv(math.factorial(n - r)) * inv(math.factorial(r))) % mod


N, M, K = list(map(int, input().split()))


def pyramid(n):
    return n * (n + 1) * (n + 2) // 6


MOD = 10 ** 9 + 7

# K ==2 の場合
# print(pyramid(N - 1) * M ** 2 + pyramid(M - 1) * N ** 2)
k2 = (pyramid(N - 1) * M ** 2 + pyramid(M - 1) * N ** 2) % MOD
# k2 を何回使うか
c = comb(N * M - 2, K - 2, mod=MOD)
print(k2 * c % MOD)

Submission status: Accepted

Solution code:
def getInv(N):#Qはmod
    inv = [0] * (N + 1)
    inv[0] = 1
    inv[1] = 1
    for i in range(2, N + 1):
        inv[i] = (-(Q // i) * inv[Q%i]) % Q
    return inv

Q = 10**9+7
modinv = getInv(4*10**5)
N, M, K = map( int, input().split())
T = N*M
comb = 1
for i in range(K-2):
    comb *= T-2-i
    comb %= Q
    comb *= modinv[K-2-i]
    comb %= Q
ans = 0
for i in range(1, M):
    ans += N**2*(M-i)*i
    ans %= Q
for i in range(1, N):
    ans += M**2*(N-i)*i
    ans %= Q
ans %= Q
print( ans*comb%Q)

Submission status: Accepted

Solution code:
def mod_inv(n: int, mod: int)->int:
    b, u, v = mod, 1, 0
    while b > 0:
        t = n // b

        n -= t * b
        u -= t * v

        n, b = b, n
        u, v = v, u

    return (u+mod) % mod


def cell_distance(N: int, M: int, K: int)->int:
    MOD = 10**9 + 7

    C = 1
    for k in range(K-2):
        C *= (N*M - K + k + 1)
        C %= MOD
        C *= mod_inv(k+1, MOD)
        C %= MOD

    S = 0
    for d in range(N):
        S += d * (N-d) * M * M
        S %= MOD

    for d in range(M):
        S += d * (M-d) * N * N
        S %= MOD

    S *= C
    S %= MOD

    return S


if __name__ == "__main__":
    N, M, K = map(int, input().split())

    ans = cell_distance(N, M, K)
    print(ans)

Submission status: Accepted

Solution code:
n,m,k=map(int,input().split())
ans=0
INF=10**9+7
if n<=m:
    M=(m-1)*m*(2*m-1)//3
    for i in range(1,n+1):
        for j in range(i,n+1):
            if i==j:
                ans=ans+m*m*(m-1)//2-(m-1)*m*(2*m-1)//6
            else:
                d=j-i
                ans=ans+d*m+(m+d)*m*(m-1)-M
            ans=ans%INF
else:
    (n,m)=(m,n)
    M=(m-1)*m*(2*m-1)//3
    for i in range(1,n+1):
        for j in range(i,n+1):
            if i==j:
                ans=ans+m*m*(m-1)//2-(m-1)*m*(2*m-1)//6
            else:
                d=j-i
                ans=ans+d*m+(m+d)*m*(m-1)-M
            ans=ans%INF
import math
ans=(math.factorial(n*m-2)//(math.factorial(n*m-k)*math.factorial(k-2)))*ans
ans=ans%INF
print(ans)

Submission status: Accepted

Solution code:
def dist(arr, n):
	res = 0; su = 0
	for i in range(n):
		res = (res + (arr[i] * i) - su) % mod
		su = (su + arr[i]) % mod

	return res % mod

n,m,k = map(int,input().split())
mod = int(1e9+7)
fact = [1]
for i in range(1,n*m+1):
	fact.append((fact[-1] * i)%mod)

x = []
y = []

for i in range(1,n+1):
	for j in range(1,m+1):
		x.append(i)
		y.append(j)

y.sort()

tot = (dist(x,n*m) + dist(y,n*m))%mod
N = n*m - 2
Need = k-2
num = fact[N]
den = (fact[Need] * fact[N-Need])%mod
ncr = (num * pow(den,mod-2,mod)) % mod
ans = (tot * ncr) % mod
print(ans)

Submission status: Accepted

Solution code:
n,m,k = map(int,input().split())
x = m*n
def cmb(n, r):
    if n - r < r: r = n - r
    if r == 0: return 1
    if r == 1: return n

    numerator = [n - r + k + 1 for k in range(r)]
    denominator = [k + 1 for k in range(r)]

    for p in range(2,r+1):
        pivot = denominator[p - 1]
        if pivot > 1:
            offset = (n - r) % p
            for k in range(p-1,r,p):
                numerator[k - offset] /= pivot
                denominator[k] /= pivot

    result = 1
    for k in range(r):
        if numerator[k] > 1:
            result *= int(numerator[k])

    return result
ss = cmb(x,k)
s = (ss*k*(k-1))//2
r = 10**9 +7
av = (((m**3-m) *n**2)//6 + ((n**3-n)*m**2)//6)%r
d = (x*(x-1))//2
print(int(((av*s)//d)%r))

Submission status: Accepted

Solution code:
import math
n, m, k = map(int, input().split())

mod = 10**9 + 7

def mpow(x, n):
  ans = 1
  while n != 0:
    if n&1:
      ans = ans * x % mod
    x = x * x % mod
    n = n // 2
  return ans

fac = [0] * (n*m+2)
finv = [0] * (n*m+2)
inv = [0] * (n*m+2)
fac[0] = 1
fac[1] = 1
finv[0] = 1
finv[1] = 1
inv[1] = 1

for i in range(2, n*m+1):
  fac[i] = fac[i-1] * i % mod
  inv[i] = mod - inv[mod%i] * (mod // i) % mod
  finv[i] = finv[i-1] * inv[i] % mod

def choose(n,k):
  if n == 0 and k == 0:
    return 1
  if n < 0:
    return 0
  if n < k:
    return 0
  return fac[n] * (finv[k] * finv[n-k] % mod) % mod


sum = 0
for i in range(1, n+1):
  for j in range(1, m+1):
    score = i-1 + j-1
    cnt = (n-i+1) * (m-j+1)
    sum += (score * cnt) % mod
    #print((i,j), score, cnt)

for i in range(2, n+1):
  for j in range(1, m):
    score = i-1 + m-j
    cnt = (n-i+1) * j
    sum += (score * cnt) % mod
    #print((i,j), score, cnt)

sum = sum * choose(m*n-2,k-2)
sum = sum % mod
print(sum)

Submission status: Accepted

Solution code:
import math
N, M, K = map(int, input().split())
sosu = 10 ** 9 + 7

def P(n, r):
    return math.factorial(n)//math.factorial(n-r)

def C(n, r):
    return P(n, r)//math.factorial(r)


dis = M * M * N * N * (N+1) // 2 % sosu
dis = dis - M * M * N * (N+1) * ( 2* N + 1) // 6 % sosu
dis = dis + N * N * M * M * (M+1) //2 % sosu
dis = dis - N * N * M * (M+1) * (2*M + 1) // 6 % sosu

def cmb(n, r, mod):
    if ( r<0 or r>n ):
        return 0
    r = min(r, n-r)
    return g1[n] * g2[r] * g2[n-r] % mod

mod = sosu #出力の制限
N = N * M
g1 = [1, 1] # 元テーブル
g2 = [1, 1] #逆元テーブル
inverse = [0, 1] #逆元テーブル計算用テーブル

for i in range( 2, N + 1 ):
    g1.append( ( g1[-1] * i ) % mod )
    inverse.append( ( -inverse[mod % i] * (mod//i) ) % mod )
    g2.append( (g2[-1] * inverse[-1]) % mod )


if K >= 3:
    dis = dis * cmb(N-2, K-2, mod) % sosu
if dis >= 0:
    print(dis)
else:
    print(dis + sosu)

Submission status: Accepted

Solution code:
N, M, K = map(int,input().split())
MOD = 10**9+7

def framod(n, mod, a=1):
    for i in range(1,n+1):
        a=a * i % mod
    return a

def power(n, r, mod):
    if r == 0: return 1
    if r%2 == 0:
        return power(n*n % mod, r//2, mod) % mod
    if r%2 == 1:
        return n * power(n, r-1, mod) % mod

def comb(n, k, mod):
    a=framod(n, mod)
    b=framod(k, mod)
    c=framod(n-k, mod)
    return (a * power(b, mod-2, mod) * power(c, mod-2, mod)) % mod

# 0~N*M-1までの1次元上と考える
ans = 0

def f(x):
    return x*(x+1)//2

c = comb(N*M-2,K-2,MOD)

for i in range(N*M):
    tate, yoko = i//M, i%M
    migi = f(M-1-yoko)
    hidari = f(yoko)
    sita = (hidari+migi)*(N-1-tate) + M*f(N-1-tate)
    coef = migi + sita
    ans += coef*c
    ans %= MOD
print(ans)

Submission status: Accepted

Solution code:
def init():
    global factorial
    MOD = 10 ** 9 + 7
    factorial = [0] * 200001
    factorial[0] = 1
    for i in range(1, 200001):
        factorial[i] = (factorial[i - 1] * i) % MOD

def main():
    init()
    n, m, k = map(int, input().split())
    s = 0
    MOD = 10 ** 9 + 7
    for i in range(n):
        s += (((((n - 1 - i) * (n - i)) // 2) % MOD) * ((m * m) % MOD)) % MOD
    for i in range(m):
        s += (((((m - 1 - i) * (m - i)) // 2) % MOD) * ((n * n) % MOD)) % MOD
    print((((s * ((factorial[n * m - 2] * pow(factorial[k - 2], MOD - 2, MOD)) % MOD)) % MOD) * pow(factorial[n * m - k], MOD - 2, MOD)) % MOD)

main()

Submission status: Accepted

Solution code:
# -*- coding: utf-8 -*-
import bisect
import heapq
import math
import random
import sys
from collections import Counter, defaultdict, deque
from decimal import ROUND_CEILING, ROUND_HALF_UP, Decimal
from functools import lru_cache, reduce
from itertools import combinations, combinations_with_replacement, product, permutations
from operator import add, mul, sub

sys.setrecursionlimit(100000)


def read_int():
    return int(input())


def read_int_n():
    return list(map(int, input().split()))


def read_float():
    return float(input())


def read_float_n():
    return list(map(float, input().split()))


def read_str():
    return input().strip()


def read_str_n():
    return list(map(str, input().split()))


def error_print(*args):
    print(*args, file=sys.stderr)


def mt(f):
    import time

    def wrap(*args, **kwargs):
        s = time.time()
        ret = f(*args, **kwargs)
        e = time.time()

        error_print(e - s, 'sec')
        return ret

    return wrap


class Mod:
    def __init__(self, m):
        self.m = m

    def add(self, a, b):
        return (a + b) % self.m

    def sub(self, a, b):
        return (a - b) % self.m

    def mul(self, a, b):
        return ((a % self.m) * (b % self.m)) % self.m

    def div(self, a, b):
        return self.mul(a, self.pow(b, self.m-2))

    def pow(self, x, y):
        if y == 0:
            return 1
        elif y == 1:
            return x % self.m
        elif y % 2 == 0:
            return self.pow(x, y//2)**2 % self.m
        else:
            return self.pow(x, y//2)**2 * x % self.m


class Combination:
    def __init__(self, n, mod):

        g1 = [1, 1]
        g2 = [1, 1]
        inverse = [0, 1]
        for i in range(2, n + 1):
            g1.append((g1[-1] * i) % mod)
            inverse.append((-inverse[mod % i] * (mod//i)) % mod)
            g2.append((g2[-1] * inverse[-1]) % mod)
        self.MOD = mod
        self.N = n
        self.g1 = g1
        self.g2 = g2
        self.inverse = inverse

    def __call__(self, n, r):
        if (r < 0 or r > n):
            return 0
        r = min(r, n-r)
        return self.g1[n] * self.g2[r] * self.g2[n-r] % self.MOD


@mt
def slv(N, M, K):
    m = Mod(10**9+7)
    d = Combination(N*M, 10**9+7)(N*M-2, K-2)

    ans = 0
    for i in range(1, N+1):
        for j in range(1, M+1):
            if i == N and j == M:
                continue
            if i == N or j == M:
                a = 1
            else:
                a = 2
            ans = m.add(ans, m.mul(i*j*a, m.mul(i + j, d)))
    ans = m.div(ans, 2)

    return ans


def main():
    N, M, K = read_int_n()
    print(slv(N, M, K))


if __name__ == '__main__':
    main()

Submission status: Accepted

Solution code:
from functools import reduce

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)

# mを法とするaの乗法的逆元
def modinv(a, m):
    g, x, y = egcd(a, m)
    if g != 1:
        raise Exception('modular inverse does not exist')
    else:
        return x % m

MOD = 10**9+7

N,M,K = map(int,input().split())

patterns = reduce(lambda a,b: (a*b)%MOD, range(N*M-K+1, N*M+1))
patterns *= modinv(reduce(lambda a,b: (a*b)%MOD, range(1, K+1)), MOD)

patterns *= (N+M)*K*(K-1)
patterns *= modinv(6, MOD)
patterns %= MOD

print(patterns)

Submission status: Accepted

Solution code:
p =10**9 +7
def comb(n,k,p):
  """power_funcを用いて(nCk) mod p を求める"""
  from math import factorial
  if n<0 or k<0 or n<k: return 0
  if n==0 or k==0: return 1
  a=factorial(n) %p
  b=factorial(k) %p
  c=factorial(n-k) %p
  return (a*power_func(b,p-2,p)*power_func(c,p-2,p))%p
def power_func(a,b,p):
  """a^b mod p を求める"""
  if b==0: return 1
  if b%2==0:
    d=power_func(a,b//2,p)
    return d*d %p
  if b%2==1:
    return (a*power_func(a,b-1,p ))%p
N,M,K =map(int,input().split())
def keisan(a,b):
    foo = a*(a**2 +3*a +2)//6
    foo *= b**2
    foo = foo%p
    return foo
ans = (keisan(N-1,M) + keisan(M-1,N))%p
foo = comb(N*M-2,K-2,p)
bar =9*17
bar =bar%p
ans = (ans*foo)%p
print(ans)

Submission status: Accepted

Solution code:
# Enter your code here. Read input from STDIN. Print output to STDOUT

MOD=10**9+7
fact=[1]
infact=[1]
temp=1
intemp=1
for i in xrange(1,10**6):
    temp*=i
    temp%=MOD
    fact+=[temp]
    intemp*=pow(i,MOD-2,MOD)
    intemp%=MOD
    infact+=[intemp]
    

def binom(a,b):
    up=fact[a]
    down=(infact[b]*infact[a-b])%MOD
    ans=(up*down)%MOD
    return ans


N,M,k=list(map(int,raw_input().strip().split(' ')))
if N==1 or M==1:
    if M<N:
        N,M=M,N
    ans=(binom(M*N,k)*binom(k,2))%MOD
    ans*=binom(M+1,3)
    ans%=MOD
    ans*=pow(binom(M*N,2),MOD-2,MOD)
    ans%=MOD
    print(ans)

else:       


    ans=(binom(M*N,k)*binom(k,2))%MOD
    ans*=((M*M*binom(N+1,3))%MOD+(N*N*binom(M+1,3))%MOD)%MOD
    ans*=pow(binom(M*N,2),MOD-2,MOD)
    ans%=MOD
    print(ans)

Submission status: Accepted

Solution code:
def mod_exp(b,e,mod):
    r = 1
    while e > 0:
        if (e&1) == 1:
            r = (r*b)%mod
        b = (b*b)%mod
        e >>= 1

    return r

# get degree of p in n! (exponent of p in the factorization of n!)
def fact_exp(n,p):
    e = 0
    u = p
    t = n
    while u <= t:
        e += t//u
        u *= p
    return e

# convert given number n into array of its base b representation
# most significant digit is at rightmost position in array
def get_base_digits(n,b):
    d = []
    while n > 0:
        d.append(n % b)
        n  = n // b
    return d

# Using Fermat's little theorem to compute nCk mod p
# considering cancelation of p in numerator and denominator
# Note: p must be prime
def fermat_binom_advanced(n,k,p):
    # check if degrees work out
    num_degree = fact_exp(n,p) - fact_exp(n-k,p)
    den_degree = fact_exp(k,p)
    if num_degree > den_degree:
        return 0
    if k > n:
        return 0
    # calculate numerator and cancel out occurrences of p
    num = 1
    for i in range(n,n-k,-1):
        cur = i
        while cur%p == 0:
            cur //= p
        num = (num*cur)%p
    # calculate denominator and cancel out occurrences of p
    denom = 1
    for i in range(1,k+1):
        cur = i
        while cur%p == 0:
            cur //= p
        denom = (denom*cur)%p
    # numerator * denominator^(p-2) (mod p)
    return (num * mod_exp(denom,p-2,p))%p

# Using Lucas' theorem to split the problem into smaller sub-problems
# In this version assuming p is prime
def lucas_binom(n,k,p):
    # get n and k in base p representation
    np = get_base_digits(n,p)
    kp = get_base_digits(k,p)
    # calculate (nCk) = (n0 choose k0)*(n1 choose k1) ... (ni choose ki) (mod p)
    binom = 1
    for i in range(len(np)-1,-1,-1):
        ni = np[i]
        ki = 0
        if i < len(kp):
            ki = kp[i]
        binom = (binom * fermat_binom_advanced(ni,ki,p)) % p
    return binom

mod = 10**9+7

n,m,k = map(int,input().split())
temp = 0
for d in range(1,m):
	temp += d*n*n*(m-d)
for d in range(1,n):
	temp += d*m*m*(n-d)
l = lucas_binom(n*m-2,k-2,mod)
print(temp*l%mod)

Submission status: Accepted

Solution code:
def cmb(n, r):
    if n - r < r: r = n - r
    if r == 0: return 1
    if r == 1: return n
 
    numerator = [n - r + k + 1 for k in range(r)]
    denominator = [k + 1 for k in range(r)]
 
    for p in range(2,r+1):
        pivot = denominator[p - 1]
        if pivot > 1:
            offset = (n - r) % p
            for k in range(p-1,r,p):
                numerator[k - offset] /= pivot
                denominator[k] /= pivot
 
    result = 1
    for k in range(r):
        if numerator[k] > 1:
            result *= int(numerator[k])
 
    return result
 
list=input().split(" ")
n=int(list[0])
m=int(list[1])
k=int(list[2])
 
niko1=int((n*(n+1)*(2*n+1)/6 - (n+1)*n*(n+1)/2 + n*n*(n+1)/2)/2)
niko2=int((m*(m+1)*(2*m+1)/6 - (m+1)*m*(m+1)/2 + m*m*(m+1)/2)/2)
niko1=niko1%(10**9+7)
niko2=niko2%(10**9+7)
niko1=niko1*m
niko2=niko2*n
niko1=niko1%(10**9+7)
niko2=niko2%(10**9+7)
niko1=niko1*m
niko2=niko2*n
niko1=niko1%(10**9+7)
niko2=niko2%(10**9+7)
niko=niko1+niko2
kosu=cmb(m*n-2, k-2)
kosu=kosu%(10**9+7)
kotae=int(niko*kosu)
kotae=kotae%(10**9+7)

print(kotae)

Submission status: Accepted

Solution code:
import sys
range = xrange
input = raw_input

MOD = 10**9+7

def modinv(x):
    return pow(x,MOD-2,MOD)

def nck(n,k):
    a = 1
    b = 1
    c = 1
    for i in range(1,n+1):
        a = a*i%MOD
    for i in range(1,k+1):
        b = b*i%MOD
    for i in range(1,n-k+1):
        c = c*i%MOD
    return a*modinv(b*c%MOD)%MOD
        
n,m,k = [int(x) for x in input().split()]

inv = modinv(n*m)
m1inv = modinv(n*m-1)

ans = 0

def summer(x):
    return x*(x-1)//2

p = (k*inv%MOD)*(m1inv*(k-1)%MOD)%MOD

xfac = n*n*p%MOD
for xi in range(m):
    ans += xfac * summer(m-xi)%MOD

yfac = m*m*p%MOD
for yi in range(n):
    ans += yfac * summer(n-yi)%MOD
#for xi in range(m):
#    for yi in range(n):
#        ans += n*p%MOD * summer(m-xi)%MOD
#        ans += m*p%MOD * summer(n-yi)%MOD
        #for xj in range(m-xi):
        #    s = xj
        #    ans += n*s*p%MOD
        
        #for yj in range(n-yi):
        #    s = yj
        #    ans += m*s*p%MOD
        
        #for xj in range(m):
        #    for yj in range(n):
        #        s = abs(xi-xj) + abs(yi - yj)     
        #        ans += p*s%MOD
print nck(n*m,k)*ans%MOD

Submission status: Accepted

Solution code:
N,M,K = list(map(int, input().split()))
mod = 10**9 + 7

tn = 0
tm = 0
for n in range(N):
  tn = (tn + n*(n+1)//2)%mod
tn = tn*M*M %mod
for m in range(M):
  tm = (tm + m*(m+1)//2)%mod
tm = tm*N*N %mod

#N*M C K-2
class combinationmod():
  def __init__(self,MAX, mod):
    self.mod =mod
    self.fac=[0]*(MAX+1)
    self.finv=[0]*(MAX+1)
    self.inv=[0]*(MAX+1)
    self.fac[0]=self.fac[1]=1
    self.finv[0]=self.finv[1]=1
    self.inv[1]=1
    for i in range(2,MAX+1):
      self.fac[i] = self.fac[i-1]*i%mod
      self.inv[i] = (-self.inv[mod%i]*(mod//i))%mod
      self.finv[i] = self.finv[i-1]* self.inv[i] %mod

  def nCk(self,n,k):
    if n<k or n<0 or k<0:
      return 0
    return self.fac[n]*(self.finv[k]*self.finv[n-k]%self.mod)%self.mod

c = combinationmod(200000,int(1e9+7))


ans = (tn + tm)*c.nCk(N*M-2, K-2) %mod
print(ans)

Submission status: Accepted

Solution code:
# coding: utf-8
# Your code here!
# coding: utf-8
# Your code here!

import sys

n,m,k = [int(i) for i in input().split()]
MOD=10**9+7
SIZE = 2*(10**5)
inv = [0]*(SIZE+1)# inv[j] = j^{-1} mod MOD
fac = [0]*(SIZE+1)# fac[j] = j! mod MOD
finv = [0]*(SIZE+1)# finv[j] = (j!)^{-1} mod MOD
inv[1] = 1
fac[0] = fac[1] = 1
finv[0] = finv[1] = 1
for i in range(2,SIZE+1):
    inv[i] = MOD -(MOD//i)*inv[MOD%i]%MOD
    fac[i] = fac[i-1]*i%MOD
    finv[i]= finv[i-1]*inv[i]%MOD

def choose(n,r): # nCk mod MOD の計算
    if 0 <= r <= n:
        return (fac[n]*finv[r]%MOD)*finv[n-r]%MOD
    else:
        return 0


ans0 = 0
ans1 = 0
for i in range(n):
    ans0 += i*(n-i)
    ans0 %= MOD
ans0 *= m
ans0 %= MOD
ans0 *= m
ans0 %= MOD

for i in range(m):
    ans1 += i*(m-i)
    ans1 %= MOD
ans1 *= n
ans1 %= MOD
ans1 *= n
ans1 %= MOD

ans = (ans0+ans1)%MOD


print(ans*choose(n*m-2,k-2)%MOD)

Submission status: Accepted

Solution code:
N, M, K = list(map(int, input().split()))

MOD = int(1e9 + 7)
NM = N * M

def su(n):
  return n * (n + 1) // 2
class comb():
  F = [1, 1]
  Fi = [1, 1]
  I = [0, 1]
  def __init__(self, num, mod):
    self.MOD = mod
    for i in range(2, num + 1):
      self.F.append((self.F[-1] * i) % mod)
      self.I.append(mod - self.I[mod % i] * (mod // i) % mod)
      self.Fi.append(self.Fi[-1] * self.I[i] % mod)
  def com(self, n, k):
    if n < k: return 0
    if n < 0 or k < 0: return 0
    return self.F[n] * (self.Fi[k] * self.Fi[n - k] % self.MOD) % self.MOD

co = comb(NM, MOD)
times = co.com(NM - 2, K - 2)
Ans = 0
for i in range(N):
  for j in range(M):
    k = (su(i) + su(N - i - 1)) * M
    k = k + (su(j) + su(M - j - 1)) * N
    k = k * times
    Ans += k

print(Ans // 2 % MOD)

Submission status: Accepted

Solution code:
import math,string,itertools,fractions,heapq,collections,re,array,bisect,sys,copy,functools
import random

sys.setrecursionlimit(10**7)
inf = 10**20
eps = 1.0 / 10**10
mod = 10**9+7
dd = [(-1,0),(0,1),(1,0),(0,-1)]
ddn = [(-1,0),(-1,1),(0,1),(1,1),(1,0),(1,-1),(0,-1),(-1,-1)]

def LI(): return list(map(int, sys.stdin.readline().split()))
def LLI(): return [list(map(int, l.split())) for l in sys.stdin.readlines()]
def LI_(): return [int(x)-1 for x in sys.stdin.readline().split()]
def LF(): return [float(x) for x in sys.stdin.readline().split()]
def LS(): return sys.stdin.readline().split()
def I(): return int(sys.stdin.readline())
def F(): return float(sys.stdin.readline())
def S(): return input()
def pf(s): return print(s, flush=True)
def pe(s): return print(str(s), file=sys.stderr)

def inv(x):
    return pow(x, mod - 2, mod)

# 1回だけの場合
def comb(n, b):
    if b > n - b:
        b = n - b
    r = 1
    for k in range(n, n-b, -1):
        r = r * k % mod
    d = 1
    for k in range(1, b+1):
        d = d * k % mod
    return r * inv(d) % mod

def main():
    n,m,k = LI()

    r = comb(n*m,k)
    cc = comb(n*m,2)
    tt = 0
    m2 = m * m
    n2 = n * n
    for i in range(n):
        v = (n-i) * (n-i-1) // 2
        tt += v * m2 % mod
        tt %= mod
    for i in range(m):
        v = (m-i) * (m-i-1) // 2
        tt += v * n2 % mod
        tt %= mod

    rr = r * k * (k-1) // 2 * tt
    rr = rr * inv(cc) % mod

    return rr%mod


print(main())

Submission status: Accepted

Solution code:
import sys
from collections import deque
from collections import defaultdict
import math
sys.setrecursionlimit(20000000)
input = sys.stdin.readline
n,m,k = map(int,input().split())
mod = 10**9 + 7
fact = [1]*(n*m+10000+1)
rfact = [1]*(n*m+10000+1)
for i in range(n*m+10000):
    fact[i+1] = rr = ((i+1) * fact[i]) % mod
def comb(N, K,mod):
    return fact[N] * pow(fact[K], mod-2, mod) * pow(fact[N-K], mod-2, mod) % mod
ans = 0
ni = comb(n*m-2,k-2,mod)
for i in range(n):
    for j in range(m):
        ko = (1+j)*j*n//2
        ko += (m-j)*(m-1-j)*n//2
        ko += (1+i)*i*m//2
        ko += (n-i)*(n-1-i)*m//2
        ans += ni*ko
        ans %= mod
print(ans*pow(2,mod-2,mod)%mod)