Problem:


Score : 1600 points

Problem Statement

Joisino is planning on touring Takahashi Town.
The town is divided into square sections by north-south and east-west lines.
We will refer to the section that is the x-th from the west and the y-th from the north as (x,y).

Joisino thinks that a touring plan is good if it satisfies the following conditions:

Let (p,q) be the section where she starts the tour. Then, X_1 \leq p \leq X_2 and Y_1 \leq q \leq Y_2 hold.

Let (s,t) be the section where she has lunch. Then, X_3 \leq s \leq X_4 and Y_3 \leq t \leq Y_4 hold.

Let (u,v) be the section where she ends the tour. Then, X_5 \leq u \leq X_6 and Y_5 \leq v \leq Y_6 hold.

By repeatedly moving to the adjacent section (sharing a side), she travels from the starting section to the ending section in the shortest distance, passing the lunch section on the way.

Two touring plans are considered different if at least one of the following is different: the starting section, the lunch section, the ending section, and the sections that are visited on the way.
Joisino would like to know how many different good touring plans there are.
Find the number of the different good touring plans.
Since it may be extremely large, find the count modulo 10^9+7.

Constraints

1 \leq X_1 \leq X_2 < X_3 \leq X_4 < X_5 \leq X_6 \leq 10^6

1 \leq Y_1 \leq Y_2 < Y_3 \leq Y_4 < Y_5 \leq Y_6 \leq 10^6

Input

Input is given from Standard Input in the following format:

X_1 X_2 X_3 X_4 X_5 X_6
Y_1 Y_2 Y_3 Y_4 Y_5 Y_6

Output

Print the number of the different good touring plans, modulo 10^9+7.

Sample Input 1

1 1 2 2 3 4
1 1 2 2 3 3

Sample Output 1

10

The starting section will always be (1,1), and the lunch section will always be (2,2).
There are four good touring plans where (3,3) is the ending section, and six good touring plans where (4,3) is the ending section.
Therefore, the answer is 6+4=10.

Sample Input 2

1 2 3 4 5 6
1 2 3 4 5 6

Sample Output 2

2346

Sample Input 3

77523 89555 420588 604360 845669 973451
2743 188053 544330 647651 709337 988194

Sample Output 3

137477680



Sample input:
1 1 2 2 3 4
1 1 2 2 3 3


Sample output:
10

Submission status: Wrong Answer

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np
from functools import lru_cache

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

@lru_cache()
def make_comb(n):
    return fact[n] * fact_inv[:n+1] % MOD * fact_inv[:n+1][::-1] % MOD

answer = 0
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    make_comb(c)
    make_comb(d)

Submission status: Wrong Answer

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

N,M,*X = map(int,read().split())

MOD = 10 ** 9 + 7

def mult(a,b,c,d,e,f):
    # (a+bx+cx^2)(d+ex+fx^2) modulo 1-4x+2x^2-x^3
    a,b,c,d,e = a*d,a*e+b*d,a*f+b*e+c*d,b*f+c*e,c*f
    b += e; c -= 4*e; d += 2*e; e = 0
    a += d; b -= 4*d; c += 2*d; d = 0
    a %= MOD; b %= MOD; c %= MOD
    return a,b,c

# (1/x)^i modulo (1-4x+2x^2-x^3)
M = 10 ** 5
A1 = [0] * (M+1)
a,b,c = 1,0,0
for i in range(M+1):
    A1[i] = (a,b,c)
    a,b,c = b+4*a,c-2*a,a
    a %= MOD; b %= MOD; c %= MOD

# (1/x)^Mi modulo (1-4x+2x^2-x^3)
A2 = [0] * (M+1)
a,b,c = 1,0,0
d,e,f = A1[M]
for i in range(M+1):
    A2[i] = (a,b,c)
    a,b,c = mult(a,b,c,d,e,f)

def power(n):
    # (1/x)^n modulo (1-4x+2x^2-x^3)
    q,r = divmod(n,M)
    a,b,c = A1[r]
    d,e,f = A2[q]
    return mult(a,b,c,d,e,f)

X.append(N)

a,b,c = 0,1,1
prev_x = 0
for x in X:
    a,b,c = mult(a,b,c,*power(x - prev_x))
    b -= a; c -= a
    prev_x = x

answer = a
print(answer)

Submission status: Wrong Answer

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

answer = 0
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    # (1+A)^c(1+B)^d / (A-B)^2 における A^aB^b の係数
    # まずはa+b+2次式部分を抽出する：A側の次数で持つ
    D = a + b + 2
    # A^i B^j の寄与。
    L = max(0, D-d)
    R = min(c, D)
    if L > R:
        continue
    x = fact[c] * fact_inv[L:R+1] % MOD * fact_inv[c-R:c-L+1][::-1] % MOD
    L,R = D-R,D-L
    y = fact[d] * fact_inv[L:R+1] % MOD * fact_inv[d-R:d-L+1][::-1] % MOD
    x *= y[::-1]
    x %= MOD
    # B=1と見立てる。(1-A)^2 で割って、A^(a-L)の係数
    np.cumsum(x,out=x)
    x %= MOD
    np.cumsum(x,out=x)
    x %= MOD
    L,R = D-R,D-L
    if 0 <= a-L < len(x):
        answer += sgn * x[a-L]

print(answer)

Submission status: Time Limit Exceeded

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
 
mod = 10**9+7 #出力の制限
N = 2*(10**6)
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 )
def sihenkei(a,b,c,d):#a to c  b to d
  return (cmb(c+d+2,c+1,mod)-cmb(d+1+a,a,mod)-cmb(b+c+1,c+1,mod)+cmb(a+b,a,mod))%mod
x1,x2,x3,x4,x5,x6=map(int,input().split())
y1,y2,y3,y4,y5,y6=map(int,input().split())
ans=0
for i in range(x3,x4+1):
  for j in range(y3,y4+1):
    ans+=sihenkei(i-x2,j-y2,i-x1,j-y1)*sihenkei(x5-i,y5-j,x6-i,y6-j)
    ans%=mod
print(ans)

Submission status: Time Limit Exceeded

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

mod = 10**9+7 #出力の制限
N = 2*(10**6)
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 )
def sihenkei(a,b,c,d):#a to c  b to d
  return (cmb(c+d+2,c+1,mod)-cmb(d+1+a,a,mod)-cmb(b+c+1,c+1,mod)+cmb(a+b,a,mod))%mod
x1,x2,x3,x4,x5,x6=map(int,input().split())
y1,y2,y3,y4,y5,y6=map(int,input().split())
ans=0
for i in range(x3,x4+1):
  for j in range(y3,y4+1):
    ans+=sihenkei(i-x2,j-y2,i-x1,j-y1)*sihenkei(x5-i,y5-j,x6-i,y6-j)
    ans%=mod
print(ans)

Submission status: Time Limit Exceeded

Solution code:
X = list(map(int, input().split()))
Y = list(map(int, input().split()))
mod = int(1e+9) + 7
def extgcd(a, b):
  if b == 0:
    return 1, 0
  else:
    x, y, u, v, k, l = 1, 0, 0, 1, a, b
    while l != 0:
      x, y, u, v = u, v, x - u * (k // l), y - v * (k // l)
      k, l = l, k % l
    return x
def inved(x):
  a = extgcd(x, mod)
  return a % mod

F = [
  [
    [1, 0, 0, Y[5]-Y[4]+1, 0], 
    [-1, 0, 0, 0, 0]
  ],
  [
    [1, X[5]-X[4]+1, 0, 0, 0], 
    [-1, 0, 0, X[5]-X[4]+1, 0]
  ],
  [
    [1, 0, 0, Y[3]-Y[2]+1, 0], 
    [-1, 0, 0, 0, Y[3]-Y[2]+1]
  ],
  [
    [1, 0, X[3]-X[2]+1, X[3]-X[2]+1, 0], 
    [-1, X[3]-X[2]+1, 0, 0, X[3]-X[2]+1]
  ],
  [
    [1, 0, 0, 0, Y[1]-Y[0]+1], 
    [-1, 0, 0, 0, 0]
  ],
  [
    [1, 0, 0, 0, X[1]-X[0]+1],
    [-1, 0, X[1]-X[0]+1, 0, 0]
  ],
  [
    [1, X[2]-X[5]-1, X[0]-X[3]-1, X[4]+Y[4]-X[3]-Y[3], X[2]+Y[2]-X[1]-Y[1]]
  ]
]
fact = [1 for i in range(2000001)]
invf = [0 for i in range(2000001)]
for i in range(2000000):
  fact[i+1] = ((i+1) * fact[i]) % mod
invf[-1] = inved(fact[2000000])
for i in range(2000000, 0, -1):
  invf[i-1] = (invf[i] * i) % mod
def encryptic(left, mat):
  cont = []
  for i in range(len(left)):
    for j in range(len(mat)):
      tmp = [left[i][0] * mat[j][0]]
      for k in range(1, 5):
        tmp.append(left[i][k] + mat[j][k])
      cont.append(tmp)
  return cont
Poly = [[1, 0, 0, 0, 0]]
for p in F:
  Poly = encryptic(p, Poly)
S = 0
D = {}
for p in Poly:
  a, b, c, d = -p[1], -p[2], p[3], p[4]
  s = p[0]
  M = min(b, c-a-2)
  tmp = 0
  for i in range(M+1):
    pix = (i+1) * invf[2+i+a] * invf[c-i-a-2]
    pix *= invf[b-i] * invf[d+i-b]
    pix %= mod
    tmp += pix
    tmp %= mod
  tmp *= fact[c] * fact[d]
  tmp %= mod
  S += s * tmp
  S %= mod
print(S)

Submission status: Memory Limit Exceeded

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np
from functools import lru_cache
from operator import itemgetter

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

@lru_cache(16)
def make_comb(n):
    return fact[n] * fact_inv[:n+1] % MOD * fact_inv[:n+1][::-1] % MOD

answer = 0
tasks = []
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    tasks.append((a,b,c,d,sgn))

tasks.sort(key = itemgetter(2))

for a,b,c,d,sgn in tasks:
    # (1+A)^c(1+B)^d / (A-B)^2 における A^aB^b の係数
    # まずはa+b+2次式部分を抽出する：A側の次数で持つ
    D = a + b + 2
    # A^i B^j の寄与。
    L = max(0, D-d)
    R = min(c, D)
    if L > R:
        continue
    x = make_comb(c)[L:R+1]
    L,R = D-R,D-L
    y = make_comb(d)[L:R+1]
    x = x * y[::-1] % MOD
    # B=1と見立てる。(1-A)^2 で割って、A^(a-L)の係数
    np.cumsum(x,out=x)
    x %= MOD
    np.cumsum(x,out=x)
    x %= MOD
    L,R = D-R,D-L
    if 0 <= a-L < len(x):
        answer += sgn * x[a-L]

answer %= MOD
print(answer)

Submission status: Memory Limit Exceeded

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np
from functools import lru_cache

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

@lru_cache()
def make_comb(n):
    return fact[n] * fact_inv[:n+1] % MOD * fact_inv[:n+1][::-1] % MOD

answer = 0
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    # (1+A)^c(1+B)^d / (A-B)^2 における A^aB^b の係数
    # まずはa+b+2次式部分を抽出する：A側の次数で持つ
    D = a + b + 2
    # A^i B^j の寄与。
    L = max(0, D-d)
    R = min(c, D)
    if L > R:
        continue
    x = make_comb(c)[L:R+1]
    L,R = D-R,D-L
    y = make_comb(d)[L:R+1]
    x = x * y[::-1] % MOD
    # B=1と見立てる。(1-A)^2 で割って、A^(a-L)の係数
    np.cumsum(x,out=x)
    x %= MOD
    np.cumsum(x,out=x)
    x %= MOD
    L,R = D-R,D-L
    if 0 <= a-L < len(x):
        answer += sgn * x[a-L]

answer %= MOD
print(answer)

Submission status: Accepted

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np
from functools import lru_cache
from operator import itemgetter

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

@lru_cache(8)
def make_comb(n):
    return fact[n] * fact_inv[:n+1] % MOD * fact_inv[:n+1][::-1] % MOD

answer = 0
tasks = []
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    tasks.append((a,b,c,d,sgn))

tasks.sort(key = itemgetter(2))

for a,b,c,d,sgn in tasks:
    # (1+A)^c(1+B)^d / (A-B)^2 における A^aB^b の係数
    # まずはa+b+2次式部分を抽出する：A側の次数で持つ
    D = a + b + 2
    # A^i B^j の寄与。
    L = max(0, D-d)
    R = min(c, D)
    if L > R:
        continue
    x = make_comb(c)[L:R+1]
    L,R = D-R,D-L
    y = make_comb(d)[L:R+1]
    x = x * y[::-1] % MOD
    # B=1と見立てる。(1-A)^2 で割って、A^(a-L)の係数
    np.cumsum(x,out=x)
    x %= MOD
    np.cumsum(x,out=x)
    x %= MOD
    L,R = D-R,D-L
    if 0 <= a-L < len(x):
        answer += sgn * x[a-L]

answer %= MOD
print(answer)

Submission status: Accepted

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np
from functools import lru_cache

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

@lru_cache(5)
def make_comb(n):
    return fact[n] * fact_inv[:n+1] % MOD * fact_inv[:n+1][::-1] % MOD

answer = 0
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    # (1+A)^c(1+B)^d / (A-B)^2 における A^aB^b の係数
    # まずはa+b+2次式部分を抽出する：A側の次数で持つ
    D = a + b + 2
    # A^i B^j の寄与。
    L = max(0, D-d)
    R = min(c, D)
    if L > R:
        continue
    x = make_comb(c)[L:R+1]
    L,R = D-R,D-L
    y = make_comb(d)[L:R+1]
    x = x * y[::-1] % MOD
    # B=1と見立てる。(1-A)^2 で割って、A^(a-L)の係数
    np.cumsum(x,out=x)
    x %= MOD
    np.cumsum(x,out=x)
    x %= MOD
    L,R = D-R,D-L
    if 0 <= a-L < len(x):
        answer += sgn * x[a-L]

answer %= MOD
print(answer)

Submission status: Accepted

Solution code:
import sys
read = sys.stdin.buffer.read
readline = sys.stdin.buffer.readline
readlines = sys.stdin.buffer.readlines

import itertools
import numpy as np

X = list(map(int,readline().split()))
Y = list(map(int,readline().split()))

for i in [1,3,5]:
    X[i] += 1
    Y[i] += 1

X1 = X[:2]; X2 = X[2:4]; X3 = X[4:]
Y1 = Y[:2]; Y2 = Y[2:4]; Y3 = Y[4:]

def cumprod(arr,MOD):
    L = len(arr); Lsq = int(L**.5+1)
    arr = np.resize(arr,Lsq**2).reshape(Lsq,Lsq)
    for n in range(1,Lsq):
        arr[:,n] *= arr[:,n-1]; arr[:,n] %= MOD
    for n in range(1,Lsq):
        arr[n] *= arr[n-1,-1]; arr[n] %= MOD
    return arr.ravel()[:L]

def make_fact(U,MOD):
    x = np.arange(U,dtype=np.int64); x[0] = 1
    fact = cumprod(x,MOD)
    x = np.arange(U,0,-1,dtype=np.int64); x[0] = pow(int(fact[-1]),MOD-2,MOD)
    fact_inv = cumprod(x,MOD)[::-1]
    return fact,fact_inv

U = 2 * 10 ** 6 + 10
MOD = 10**9 + 7
fact, fact_inv = make_fact(U,MOD)

answer = 0
for p in itertools.product([0,1],repeat=6):
    x1,x2,x3,y1,y2,y3 = [A[i] for A,i in zip([X1,X2,X3,Y1,Y2,Y3],p)]
    sgn = (-1) ** sum(p)
    a,b,c,d = x2-x1, x3-x2, x2-x1+y2-y1, x3-x2+y3-y2
    c += 2; d += 2; sgn = -sgn
    # (1+A)^c(1+B)^d / (A-B)^2 における A^aB^b の係数
    # まずはa+b+2次式部分を抽出する：A側の次数で持つ
    D = a + b + 2
    # A^i B^j の寄与。
    L = max(0, D-d)
    R = min(c, D)
    if L > R:
        continue
    x = fact[c] * fact_inv[L:R+1] % MOD * fact_inv[c-R:c-L+1][::-1] % MOD
    L,R = D-R,D-L
    y = fact[d] * fact_inv[L:R+1] % MOD * fact_inv[d-R:d-L+1][::-1] % MOD
    x *= y[::-1]
    x %= MOD
    # B=1と見立てる。(1-A)^2 で割って、A^(a-L)の係数
    np.cumsum(x,out=x)
    x %= MOD
    np.cumsum(x,out=x)
    x %= MOD
    L,R = D-R,D-L
    if 0 <= a-L < len(x):
        answer += sgn * x[a-L]

answer %= MOD
print(answer)