Buggy python submission (verdict: Runtime Error):
import sys
input = sys.stdin.readline

from scipy.sparse import csr_matrix
from scipy.sparse.csgraph import connected_components

import time
start = time.time()

MOD = 10**9 + 7

N = int(input())
ball = (tuple(int(x) for x in row.split()) for row in sys.stdin.readlines())

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = [False] * (N+N+1)
components = []
for x in range(1,N+N+1):
    if visited[x]:
        continue
    V = set([x])
    E = []
    q = [x]
    visited[x] = True
    while q:
        y = q.pop()
        for z in graph[y]:
            if y < z:
                E.append((y,z))
            if visited[z]:
                continue
            V.add(z)
            visited[z] = True
            q.append(z)
    components.append((V,E))

def make_get_pattern(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[(x,y)] = 0
        else:
            get[(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[(x,y)] = 2
    while graph[x]:
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[(x,y)] = 3
        else:
            get[(y,x)] = 2
        x = y
    return get

def F(V,E):
    # V is connected
    if len(E) != len(V):
        return 0
    ret = 0
    E.sort()
    get = make_get_pattern(V)
    den1,den2 = 1,1
    dp1 = {x:0 for x in V}
    dp2 = {x:0 for x in V}
    for x,y in E:
        if get[(x,y)] == 0:
            k1 = dp1[x] + 1; k2 = dp2[x] + 1
        elif get[(x,y)] == 1:
            k1 = dp1[y] + 1; k2 = dp2[y] + 1
        elif get[(x,y)] == 2:
            k1 = dp1[x] + 1; k2 = dp2[y] + 1
        else:
            k1 = dp1[y] + 1; k2 = dp2[x] + 1
        dp1[x] += k1; dp1[y] += k1
        dp2[x] += k2; dp2[y] += k2
        den1 *= k1; den2 *= k2
        den1 %= MOD; den2 %= MOD
    return sum(pow(x,MOD-2,MOD) for x in (den1,den2))

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

time.time()-start

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

from scipy.sparse import csr_matrix
from scipy.sparse.csgraph import connected_components

MOD = 10**9 + 7

N = int(input())
ball = (tuple(int(x) for x in row.split()) for row in sys.stdin.readlines())

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = [False] * (N+N+1)
components = []
for x in range(1,N+N+1):
    if visited[x]:
        continue
    V = set([x])
    E = []
    q = [x]
    visited[x] = True
    while q:
        y = q.pop()
        for z in graph[y]:
            if y < z:
                E.append((y,z))
            if visited[z]:
                continue
            V.add(z)
            visited[z] = True
            q.append(z)
    components.append((V,E))

def make_get_pattern(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[(x,y)] = 0
        else:
            pass
            get[(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[(x,y)] = 2
    while graph[x]:
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[(x,y)] = 3
        else:
            get[(y,x)] = 2
        x = y
    return get

def F(V,E):
    # V is connected
    if len(E) != len(V):
        return 0
    ret = 0
    E.sort()
    get = make_get_pattern(V)
    den1,den2 = 1,1
    dp1 = {x:0 for x in V}
    dp2 = {x:0 for x in V}
    for x,y in E:
        if get[(x,y)] == 0:
            k1 = dp1[x] + 1; k2 = dp2[x] + 1
        elif get[(x,y)] == 1:
            k1 = dp1[y] + 1; k2 = dp2[y] + 1
        elif get[(x,y)] == 2:
            k1 = dp1[x] + 1; k2 = dp2[y] + 1
        else:
            k1 = dp1[y] + 1; k2 = dp2[x] + 1
        dp1[x] += k1; dp1[y] += k1
        dp2[x] += k2; dp2[y] += k2
        den1 *= k1; den2 *= k2
        den1 %= MOD; den2 %= MOD
    return sum(pow(x,MOD-2,MOD) for x in (den1,den2))

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

from scipy.sparse import csr_matrix
from scipy.sparse.csgraph import connected_components

MOD = 10**9 + 7

N = int(input())
ball = (tuple(int(x) for x in row.split()) for row in sys.stdin.readlines())

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = [False] * (N+N+1)
components = []
for x in range(1,N+N+1):
    if visited[x]:
        continue
    V = set([x])
    E = []
    q = [x]
    visited[x] = True
    while q:
        y = q.pop()
        for z in graph[y]:
            if y < z:
                E.append((y,z))
            if visited[z]:
                continue
            V.add(z)
            visited[z] = True
            q.append(z)
    components.append((V,E))

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while graph[x]:
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V,E):
    # V is connected
    if len(E) != len(V):
        return 0
    ret = 0
    E.sort()
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

from scipy.sparse import csr_matrix
from scipy.sparse.csgraph import connected_components

MOD = 10**9 + 7

N = int(input())
ball = (tuple(int(x) for x in row.split()) for row in sys.stdin.readlines())

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = [False] * (N+N+1)
components = []
for x in range(1,N+N+1):
    if visited[x]:
        continue
    V = set([x])
    q = [x]
    visited[x] = True
    while q:
        y = q.pop()
        for z in graph[y]:
            if visited[z]:
                continue
            V.add(z)
            visited[z] = True
            q.append(z)
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while graph[x]:
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

from scipy.sparse import csr_matrix
from scipy.sparse.csgraph import connected_components

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
row = []
col = []
for x,y in ball:
    row.append(x)
    col.append(y+N)
    graph[x].add(y+N)
    graph[y+N].add(x)

Ncomp,comp = connected_components(csr_matrix(([1]*(N+N),(row,col)),(N+N+1,N+N+1)))

components = [[] for _ in range(Ncomp-1)]
for i,x in enumerate(comp[1:],1):
    components[x-1].append(i)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while graph[x]:
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = [False] * (N+N+1)
components = []
for x in range(1,N+N+1):
    if visited[x]:
        continue
    V = [x]
    q = [x]
    visited[x] = True
    while q:
        y = q.pop()
        for z in graph[y]:
            if visited[z]:
                continue
            V.append(z)
            visited[z] = True
            q.append(z)
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = [False] * (N+N+1)
components = []
for x in range(1,N+N+1):
    if visited[x]:
        continue
    V = set([x])
    q = [x]
    visited[x] = True
    while q:
        y = q.pop()
        for z in graph[y]:
            if visited[z]:
                continue
            V.add(z)
            visited[z] = True
            q.append(z)
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = set(x for x in V if len(graph[x]) == 1)
    get = {},{}
    while deg1:
        x = deg1.pop()
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.add(y)
        elif len(se) == 0: deg1.remove(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

U = N+N+10
fact = [1] * (U+1); fact_inv = [1] * (U+1)
for n in range(1,U+1):
    fact[n] = fact[n-1] * n % MOD
fact_inv[U] = pow(fact[U],MOD-2,MOD)
for n in range(U,0,-1):
    fact_inv[n-1] = fact_inv[n] * n % MOD

X_to_Y = [[] for _ in range(N+1)]
Y_to_X = [[] for _ in range(N+1)]
for x,y in ball:
    X_to_Y[x].append(y)
    Y_to_X[y].append(x)

# とりあえず連結成分に分けよう
visited_x = set()
visited_y = set()
components = []
for x,y in ball:
    if x in visited_x:
        continue
    V = set()
    V_x = set()
    V_y = set()
    q = [(x,y)]
    while q:
        x,y = q.pop()
        V.add((x,y))
        if x not in V_x:
            V_x.add(x)
            for t in X_to_Y[x]:
                q.append((x,t))
        if y not in V_y:
            V_y.add(y)
            for t in Y_to_X[y]:
                q.append((t,y))
        V_x.add(x)
        V_y.add(y)
    visited_x |= V_x
    visited_y |= V_y
    components.append((V,V_x,V_y))

def make_get_patterns(V_x,V_y):
    rest_xy = {x:set(X_to_Y[x]) for x in V_x}
    rest_yx = {y:set(Y_to_X[y]) for y in V_y}
    deg1_x = set(x for x in V_x if len(rest_xy[x]) == 1)
    deg1_y = set(y for y in V_y if len(rest_yx[y]) == 1)
    get = {},{}
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.remove(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 1
    for x in V_x:
        if rest_xy[x]:
            y = rest_xy[x].pop()
            break
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    rest_yx[y].remove(x); deg1_y.add(y)
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 1
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.discard(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
    return get

def F(V,V_x,V_y):
    # 連結成分上の順序がうまくいっている確率を返す
    if len(V) != len(V_x) + len(V_y):
        return 0
    V = sorted(V)
    ret = 0
    for get in make_get_patterns(V_x,V_y):
        den = 1
        dp_x = {x:0 for x in V_x}
        dp_y = {y:0 for y in V_y}
        for x,y in V:
            if get[(x,y)] == 0: # x軸でとる
                k = dp_x[x] + 1
            else:
                k = dp_y[y] + 1
            dp_x[x] += k
            dp_y[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD
    if prob == 0:
        break
answer = prob * fact[N+N] % MOD
print(answer)



Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Time Limit Exceeded):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

U = N+N+10
fact = [1] * (U+1); fact_inv = [1] * (U+1)
for n in range(1,U+1):
    fact[n] = fact[n-1] * n % MOD
fact_inv[U] = pow(fact[U],MOD-2,MOD)
for n in range(U,0,-1):
    fact_inv[n-1] = fact_inv[n] * n % MOD

X_to_Y = [[] for _ in range(N+1)]
Y_to_X = [[] for _ in range(N+1)]
for x,y in ball:
    X_to_Y[x].append(y)
    Y_to_X[y].append(x)

# とりあえず連結成分に分けよう
visited_x = set()
visited_y = set()
components = []
for x,y in ball:
    if x in visited_x:
        continue
    V = set()
    V_x = set()
    V_y = set()
    q = [(x,y)]
    while q:
        x,y = q.pop()
        V.add((x,y))
        if x not in V_x:
            V_x.add(x)
            for t in X_to_Y[x]:
                q.append((x,t))
        if y not in V_y:
            V_y.add(y)
            for t in Y_to_X[y]:
                q.append((t,y))
        V_x.add(x)
        V_y.add(y)
    visited_x |= V_x
    visited_y |= V_y
    components.append((V,V_x,V_y))

def make_get_patterns(V_x,V_y):
    rest_xy = {x:set(X_to_Y[x]) for x in V_x}
    rest_yx = {y:set(Y_to_X[y]) for y in V_y}
    deg1_x = set(x for x in V_x if len(rest_xy[x]) == 1)
    deg1_y = set(y for y in V_y if len(rest_yx[y]) == 1)
    get = {},{}
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.remove(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 1
    for x in V_x:
        if rest_xy[x]:
            y = rest_xy[x].pop()
            break
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    rest_yx[y].remove(x); deg1_y.add(y)
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 1
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.discard(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
    return get

def F(V,V_x,V_y):
    # 連結成分上の順序がうまくいっている確率を返す
    if len(V) != len(V_x) + len(V_y):
        return 0
    V = sorted(V)
    ret = 0
    for get in make_get_patterns(V_x,V_y):
        den = 1
        dp_x = {x:0 for x in V_x}
        dp_y = {y:0 for y in V_y}
        for x,y in V:
            if get[(x,y)] == 0: # x軸でとる
                k = dp_x[x] + 1
            else:
                k = dp_y[y] + 1
            dp_x[x] += k
            dp_y[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD
answer = prob * fact[N+N] % MOD
print(answer)



Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Runtime Error):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

U = N+N+10
fact = [1] * (U+1); fact_inv = [1] * (U+1)
for n in range(1,U+1):
    fact[n] = fact[n-1] * n % MOD
fact_inv[U] = pow(fact[U],MOD-2,MOD)
for n in range(U,0,-1):
    fact_inv[n-1] = fact_inv[n] * n % MOD

X_to_Y = [[] for _ in range(N+1)]
Y_to_X = [[] for _ in range(N+1)]
for x,y in ball:
    X_to_Y[x].append(y)
    Y_to_X[y].append(x)

# とりあえず連結成分に分けよう
visited_x = set()
visited_y = set()
components = []
for x,y in ball:
    if x in visited_x:
        continue
    V = set()
    V_x = set()
    V_y = set()
    q = [(x,y)]
    while q:
        x,y = q.pop()
        V.add((x,y))
        if x not in V_x:
            V_x.add(x)
            for t in X_to_Y[x]:
                q.append((x,t))
        if y not in V_y:
            V_y.add(y)
            for t in Y_to_X[y]:
                q.append((t,y))
        V_x.add(x)
        V_y.add(y)
    visited_x |= V_x
    visited_y |= V_y
    components.append((V,V_x,V_y))

def make_get_patterns(V_x,V_y):
    assert len(V) == len(V_x) + len(V_y)
    rest_xy = {x:set(X_to_Y[x]) for x in V_x}
    rest_yx = {y:set(Y_to_X[y]) for y in V_y}
    deg1_x = set(x for x in V_x if len(rest_xy[x]) == 1)
    deg1_y = set(y for y in V_y if len(rest_yx[y]) == 1)
    get = {},{}
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.remove(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 1
    for x in V_x:
        if rest_xy[x]:
            y = rest_xy[x].pop()
            break
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    rest_yx[y].remove(x); deg1_y.add(y)
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 1
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.discard(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
    return get

def F(V,V_x,V_y):
    # 連結成分上の順序がうまくいっている確率を返す
    if len(V) != len(V_x) + len(V_y):
        return 0
    V = sorted(V)
    ret = 0
    for get in make_get_patterns(V_x,V_y):
        den = 1
        dp_x = {x:0 for x in V_x}
        dp_y = {y:0 for y in V_y}
        for x,y in V:
            if get[(x,y)] == 0: # x軸でとる
                k = dp_x[x] + 1
            else:
                k = dp_y[y] + 1
            dp_x[x] += k
            dp_y[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD
answer = prob * fact[N+N] % MOD
print(answer)



Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Runtime Error):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

U = N+N+10
fact = [1] * (U+1); fact_inv = [1] * (U+1)
for n in range(1,U+1):
    fact[n] = fact[n-1] * n % MOD
fact_inv[U] = pow(fact[U],MOD-2,MOD)
for n in range(U,0,-1):
    fact_inv[n-1] = fact_inv[n] * n % MOD

X_to_Y = [[] for _ in range(N+1)]
Y_to_X = [[] for _ in range(N+1)]
for x,y in ball:
    X_to_Y[x].append(y)
    Y_to_X[y].append(x)

# とりあえず連結成分に分けよう
visited_x = set()
visited_y = set()
components = []
for x,y in ball:
    if x in visited_x:
        continue
    V = set()
    V_x = set()
    V_y = set()
    q = [(x,y)]
    while q:
        x,y = q.pop()
        V.add((x,y))
        if x not in V_x:
            V_x.add(x)
            for t in X_to_Y[x]:
                q.append((x,t))
        if y not in V_y:
            V_y.add(y)
            for t in Y_to_X[y]:
                q.append((t,y))
        V_x.add(x)
        V_y.add(y)
    visited_x |= V_x
    visited_y |= V_y
    components.append((V,V_x,V_y))

def make_get_patterns(V_x,V_y):
    assert len(V) == len(V_x) + len(V_y)
    rest_xy = {x:set(X_to_Y[x]) for x in V_x}
    rest_yx = {y:set(Y_to_X[y]) for y in V_y}
    deg1_x = set(x for x in V_x if len(rest_xy[x]) == 1)
    deg1_y = set(y for y in V_y if len(rest_yx[y]) == 1)
    get = {},{}
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.remove(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 1
    for x in V_x:
        if rest_xy[x]:
            y = rest_xy[x].pop()
            break
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    rest_yx[y].remove(x); deg1_y.add(y)
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 1
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.discard(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
    return get

def F(V,V_x,V_y):
    # 連結成分上の順序がうまくいっている確率を返す
    if len(V) != len(V_x) + len(V_y):
        return 0
    V = sorted(V)
    ret = 0
    for get in make_get_patterns(V_x,V_y):
        den = 1
        dp_x = {x:0 for x in V_x}
        dp_y = {y:0 for y in V_y}
        for x,y in V:
            if get[(x,y)] == 0: # x軸でとる
                k = dp_x[x] + 1
            else:
                k = dp_y[y] + 1
            dp_x[x] += k
            dp_y[y] += k
            den *= k
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD
answer = prob * fact[N+N] % MOD
print(answer)



Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)

Buggy python submission (verdict: Runtime Error):
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

U = N+N+10
fact = [1] * (U+1); fact_inv = [1] * (U+1)
for n in range(1,U+1):
    fact[n] = fact[n-1] * n % MOD
fact_inv[U] = pow(fact[U],MOD-2,MOD)
for n in range(U,0,-1):
    fact_inv[n-1] = fact_inv[n] * n % MOD

X_to_Y = [[] for _ in range(N+1)]
Y_to_X = [[] for _ in range(N+1)]
for x,y in ball:
    X_to_Y[x].append(y)
    Y_to_X[y].append(x)

# とりあえず連結成分に分けよう
visited_x = set()
visited_y = set()
components = []
for x,y in ball:
    if x in visited_x:
        continue
    V = set()
    V_x = set()
    V_y = set()
    q = [(x,y)]
    while q:
        x,y = q.pop()
        V.add((x,y))
        if x not in V_x:
            V_x.add(x)
            for t in X_to_Y[x]:
                q.append((x,t))
        if y not in V_y:
            V_y.add(y)
            for t in Y_to_X[y]:
                q.append((t,y))
        V_x.add(x)
        V_y.add(y)
    visited_x |= V_x
    visited_y |= V_y
    components.append((V,V_x,V_y))

def make_get_patterns(V_x,V_y):
    assert len(V) == len(V_x) + len(V_y)
    rest_xy = {x:set(X_to_Y[x]) for x in V_x}
    rest_yx = {y:set(Y_to_X[y]) for y in V_y}
    deg1_x = set(x for x in V_x if len(rest_xy[x]) == 1)
    deg1_y = set(y for y in V_y if len(rest_yx[y]) == 1)
    get = {},{}
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.remove(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 1
    for x in V_x:
        if rest_xy[x]:
            y = rest_xy[x].pop()
            break
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    rest_yx[y].remove(x); deg1_y.add(y)
    while deg1_x or deg1_y:
        while deg1_x:
            x = deg1_x.pop()
            y = rest_xy[x].pop()
            se = rest_yx[y]; se.remove(x)
            if len(se) == 1: deg1_y.add(y)
            elif len(se) == 0: deg1_y.remove(y)
            get[0][(x,y)] = 0; get[1][(x,y)] = 1
        while deg1_y:
            y = deg1_y.pop()
            x = rest_yx[y].pop()
            se = rest_xy[x]; se.remove(y)
            if len(se) == 1: deg1_x.add(x)
            elif len(se) == 0: deg1_x.discard(x)
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
    return get

def F(V,V_x,V_y):
    # 連結成分上の順序がうまくいっている確率を返す
    if len(V) != len(V_x) + len(V_y):
        return 0
    V = sorted(V)
    ret = 0
    for get in make_get_patterns(V_x,V_y):
        den = 1
        dp_x = {x:0 for x in V_x}
        dp_y = {y:0 for y in V_y}
        for x,y in V:
            if get[(x,y)] == 0: # x軸でとる
                k = dp_x[x] + 1
            else:
                k = dp_y[y] + 1
            dp_x[x] += k
            dp_y[y] += k
            den *= k
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(*c)
    prob %= MOD
answer = prob * fact[N+N] % MOD
print(answer)



Fixed submission:
import sys
input = sys.stdin.readline

MOD = 10**9 + 7

N = int(input())
ball = [tuple(int(x) for x in input().split()) for _ in range(N+N)]

# x座標を1,2,...,N
# y座標をN+1,N+2,...,N+N

graph = [set() for _ in range(N+N+1)]
for x,y in ball:
    graph[x].add(y+N)
    graph[y+N].add(x)

visited = set()
components = []
for x in range(1,N+N+1):
    if x in visited:
        continue
    V = set([x])
    q = [x]
    while q:
        y = q.pop()
        for z in graph[y]:
            if z in V:
                continue
            V.add(z)
            q.append(z)
    visited |= V
    components.append(V)

def make_get_patterns(V):
    deg1 = [x for x in V if len(graph[x]) == 1]
    get = {},{}
    while deg1:
        x = deg1.pop()
        if not graph[x]:
            continue
        y = graph[x].pop()
        se = graph[y]; se.remove(x)
        if len(se) == 1: deg1.append(y)
        if x < y:
            get[0][(x,y)] = 0; get[1][(x,y)] = 0
        else:
            pass
            get[0][(y,x)] = 1; get[1][(y,x)] = 1
    for x in V:
        if graph[x]:
            y = graph[x].pop()
            break
    # 残りはサイクル
    graph[y].remove(x)
    if x > y: x,y = y,x
    get[0][(x,y)] = 0; get[1][(x,y)] = 1
    while True:
        if not graph[x]:
            break
        y = graph[x].pop()
        graph[y].remove(x)
        if x < y:
            get[0][(x,y)] = 1; get[1][(x,y)] = 0
        else:
            get[0][(y,x)] = 0; get[1][(y,x)] = 1
        x = y
    return get

def F(V):
    # V is connected
    E = sorted((x,y) for x in V if x <= N for y in graph[x])
    if len(E) != len(V):
        return 0
    ret = 0
    for get in make_get_patterns(V):
        den = 1
        dp = {x:0 for x in V}
        for x,y in E:
            if get[(x,y)] == 0:
                k = dp[x] + 1
            else:
                k = dp[y] + 1
            dp[x] += k
            dp[y] += k
            den *= k
            den %= MOD
        ret += pow(den,MOD-2,MOD)
    return ret % MOD

prob = 1
for c in components:
    prob *= F(c)
    prob %= MOD

answer = prob
for n in range(1,N+N+1):
    answer *= n
    answer %= MOD
print(answer)