各質問に $O(1)$ とかその辺で答える必要がある。結構場合分けが面倒くさい。
まず、$N \lt A+B$ の場合、どうあがいても重ならずに置けないので0。以下それ以外とする。
プロペラみたいに分割すると重複無く数えられる。
●●○○○○
●●○○○○
●●○○○○
●●赤赤■■
●●赤赤■■
□□□□■■
全ての赤の位置に対して、「■の中に青が置かれる」「○と■にまたがって青が置かれる」パターン数を数えると、その合計の4倍が答えとなる。
(●、○、□に置かれるパターン数は、これを回転させたものと合致し、それぞれ同数だけ存在する)
$A \lt B$ だと $A$ より前に $B$ が下端の制約を受けるので、場合分けが増えてしまう。$A \ge B$ とする。
赤を上から $p$、左から $q$ の位置に置くとすると、$(p=0~N-A,q=0~N-A-B)$
青は、縦 $N$、横 $N-q-A$ の範囲からはみ出さないような位置に置ける。
$$\displaystyle \sum_{p=0}^{B-1}\sum_{q=0}^{N-A-B}(N-B+1)(N-q-A-B+1)$$
Wolfram Alpha先生で式変形すると、以下のようになる。
$$\dfrac{1}{2}B(N-B+1)(N-A-B+2)(N-A-B+1)$$
青は、縦 $N-p+B-1$、横 $N-q-A$ の範囲からはみ出さないような位置に置ける。
\begin{eqnarray}
& \sum_{p=B}^{N-A}\sum_{q=0}^{N-A-B}(N-p)(N-q-A-B+1) \\
=& \frac{1}{4}(N-A-B+2)(N-A-B+1)^2(N+A-B)
\end{eqnarray}
これを計算して合計して4倍すればよい。
Python3
t = int(input())
buf = []
MOD = 10 ** 9 + 7
for _ in range(t):
n, a, b = map(int, input().split())
if n < a + b:
buf.append(0)
continue
if a < b:
a, b = b, a
ans1 = (n - a - b + 2) * (n - a - b + 1) ** 2 * (n + a - b) % MOD
ans2 = 2 * b * (n - b + 1) * (n - a - b + 2) * (n - a - b + 1) % MOD
buf.append((ans1 + ans2) % MOD)
print('\n'.join(map(str, buf)))
典型的な、“縦に見るものを横に見る” 問題。
「このマスが照らされるような置き方の個数」を、全ての空きマスについて合計すると、それが答えとなる。
あるマスから上下左右に繋がる空きマスの個数を $L$ とすると、この中のどれか一つさえ照明が置かれていれば、そのマスは照らされることになる。
│
│
──○─
│
少なくともどれか1つに置かれる場合は $2^L-1$ 通り。また、これに関わらない他の空きマスはどうなっていてもよいので、$2^{K-L}$ 通り。
これを掛け合わせた $(2^{L}-1)2^{K-L}$ が、あるマスが照らされる置き方の個数となる。
$L$ の求め方は、空きマスが横にいくつ連続するか、縦にいくつ連続するかをそれぞれ個別に求めておくと、その合計-1が $L$ となる。
$2000 \times 2000$ のグリッドを何回か走査することになるので、素のPythonだと厳しい。NumbaやPyPyなどを使った方がよい。
Python3
import os
import sys
import numpy as np
def solve(h, w, f):
hor = np.zeros_like(f, dtype=np.int64)
ver = np.zeros_like(f, dtype=np.int64)
MOD = 10 ** 9 + 7
def mod_pow(x, a, MOD):
ret = 1
cur = x
while a:
if a & 1:
ret = ret * cur % MOD
cur = cur * cur % MOD
a >>= 1
return ret
for i in range(h):
suc = 0
for j in range(w):
if f[i][j] == 1:
suc = 0
else:
suc += 1
hor[i][j] = suc
for j in range(w):
suc = 0
for i in range(h):
if f[i][j] == 1:
suc = 0
else:
suc += 1
ver[i][j] = suc
for i in range(h):
for j in range(w - 2, -1, -1):
if hor[i][j] == 0 or hor[i][j + 1] == 0:
continue
else:
hor[i][j] = hor[i][j + 1]
for j in range(w):
for i in range(h - 2, -1, -1):
if ver[i][j] == 0 or ver[i + 1][j] == 0:
continue
else:
ver[i][j] = ver[i + 1][j]
blank = h * w - f.sum()
tot = hor + ver - 1
ans = 0
for i in range(h):
for j in range(w):
if tot[i][j] == -1:
continue
ans = (ans + (mod_pow(2, tot[i][j], MOD) - 1) * mod_pow(2, blank - tot[i][j], MOD)) % MOD
return ans
if sys.argv[-1] == 'ONLINE_JUDGE':
from numba.pycc import CC
cc = CC('my_module')
cc.export('solve', '(i8,i8,i1[:,:],)')(solve)
cc.compile()
exit()
if os.name == 'posix':
# noinspection PyUnresolvedReferences
from my_module import solve
else:
from numba import njit
solve = njit('(i8,i8,i1[:,:],)', cache=True)(solve)
print('compiled', file=sys.stderr)
h, w = map(int, input().split())
field = np.array([['.#'.index(c) for c in input()] for _ in range(h)], dtype=np.int8)
ans = solve(h, w, field)
print(ans)
微分積分
変曲点毎に区間を分けて区間毎に処理
プログラム上でのMODを取りながらの微分積分
このあたりがポイントとなる。
以下の2つの確率変数があったとき、$L,R$ に出現する値で区切る。
i 0 1 2 3 4 5 6 7 8 9 10
1 |--------------|
2 |--------------------|
↓
|-----|--------|-----------|
区間毎に、「全ての確率変数が、その区間内の実数 $x$ 以下の値を取る場合の確率」を、多項式で表現できる。($f(x)$ とする)
1つの確率変数についてみれば、
これを掛け合わせたのが、全体の確率となる。
i 0 1 2 3 4 5 6 7 8 9 10
1 |--------------|
2 |--------------------|
|-----|--------|-----------|
x-1 x-1
変数1 ----- ----- 1
6-1 6-1
x-3 x-3
変数2 0 ----- -----
10-3 10-3
かけあわせる↓
(x-1)(x-3) x-3
各区間の 0 ------------ -----
f(x) (6-1)(10-3) 10-3
最初に $L,R$ で区切ったのは、各確率変数がそこで滑らかでなくなり、それを掛け合わせた $f(x)$ も滑らかでなくなるから(滑らかでないと微分できない)。
さて、これはあくまで $x$ 以下の確率なので、$x$ ちょうどになる確率1)を求めたい。
これは微分した式 $f'(x)$ がそれに当たる。
求めるのは期待値なので、$xf'(x)$ が、「最大値が $x$ だった場合に答えに寄与する値」となる。
これを、今度は積分を取ると、「区間全体が答えに寄与する値」を求められる。
よって、区間の両端を $[P,Q)$ とすると、以下の式を区間毎に計算し、合計したもの(に、かけろと言われている数をかけたもの)が答えとなる。
実装上は、多項式を表現できるクラスを作っておいて、
$(a_0+a_1x+a_2x^2+...)$ と $b_0+b_1x$ をかける
$(a_0+a_1x+a_2x^2+...)$ を $b_0+b_1x$ で割る(割り切れることは保証される)
$(a_0+a_1x+a_2x^2+...)$ を微分する
$(a_0+a_1x+a_2x^2+...)$ を積分する
$(a_0+a_1x+a_2x^2+...)$ の $x$ に具体的な値を代入する
等ができるようにしておくと、$f(x)=1$ で初期化しておいて、区間 $[P, Q)$ を昇順に
もし $L_{max} \lt Q$ なら、その区間は答えに寄与するので、$\displaystyle \int_P^Q xf'(x)dx$ を答えに加算する
$Q$ が $L_i$ 由来なら、$f(x)$ にその確率変数の $\dfrac{x-L_i}{R_i-L_i}$ をかける
$Q$ が $R_i$ 由来なら、$f(x)$ をその確率変数の $\dfrac{x-L_i}{R_i-L_i}$ で割る
とすればよい。
Python3
import os
import sys
import numpy as np
def solve(inp):
MOD = 10 ** 9 + 7
def mod_pow(x, a, MOD):
ret = 1
cur = x
while a:
if a & 1:
ret = ret * cur % MOD
cur = cur * cur % MOD
a >>= 1
return ret
INVERSE_CACHE = {-1: -1}
def inv(x):
if x not in INVERSE_CACHE:
INVERSE_CACHE[x] = mod_pow(x, MOD - 2, MOD)
return INVERSE_CACHE[x]
# [x^0, x^1, ...] の係数で多項式を表現
def poly_mul(p1, p2):
# 今のところ、p2は ax+b という形のみ
ret = np.zeros(p1.size + 1, np.int64)
ret[:p1.size] = p1 * p2[0] % MOD
ret[1:] += p1 * p2[1]
ret %= MOD
return ret
def poly_div(p1, p2):
# 今のところ、p2は ax+b という形のみ
ret = np.zeros(p1.size - 1, np.int64)
a = 0
iv = inv(p2[0])
for i in range(p1.size - 1):
ret[i] = a = (p1[i] - a * p2[1]) % MOD * iv % MOD
return ret
def poly_derivative(p):
ret = np.zeros(p.size - 1, np.int64)
for i in range(p.size - 1):
ret[i] = p[i + 1] * (i + 1) % MOD
return ret
def poly_integral(p):
ret = np.zeros(p.size + 1, np.int64)
for i in range(p.size):
ret[i + 1] = p[i] * inv(i + 1) % MOD
return ret
def poly_substitute(p, x):
ret = 0
y = 1
for i in range(p.size):
ret = (ret + p[i] * y) % MOD
y = y * x % MOD
return ret
n = inp[0]
lr = inp[1:]
slr = [(-1, 0, -1)]
for i in range(2 * n):
slr.append((lr[i], i & 1, i >> 1))
slr.sort()
lll = lr[0::2]
rrr = lr[1::2]
max_l = lll.max()
poly_x = np.array([0, 1], np.int64)
poly = np.ones(1, np.int64)
ans = 0
for i in range(2 * n):
s, f, j = slr[i]
t, g, k = slr[i + 1]
if t > max_l and s < t:
f = poly_integral(poly_mul(poly_derivative(poly), poly_x))
tmp = poly_substitute(f, t) - poly_substitute(f, s)
ans = (ans + tmp) % MOD
l = lll[k]
r = rrr[k]
if g == 0:
mul = np.array([-l, 1], np.int64)
poly = poly_mul(poly, mul) * inv(r - l) % MOD
else:
div = np.array([-l, 1], np.int64)
poly = poly_div(poly, div) * (r - l) % MOD
for i in range(2, n + 2):
ans = ans * i % MOD
for i in range(n):
l = lll[i]
r = rrr[i]
ans = ans * (r - l) % MOD
return ans
if sys.argv[-1] == 'ONLINE_JUDGE':
from numba.pycc import CC
cc = CC('my_module')
cc.export('solve', '(i8[:],)')(solve)
cc.compile()
exit()
if os.name == 'posix':
# noinspection PyUnresolvedReferences
from my_module import solve
else:
from numba import njit
solve = njit('(i8[:],)', cache=True)(solve)
print('compiled', file=sys.stderr)
inp = np.fromstring(sys.stdin.read(), dtype=np.int64, sep=' ')
ans = solve(inp)
print(ans)