結果

問題 No.8031 （物理学）長距離相互作用
ユーザー Heartbreak
提出日時 2025-04-10 23:47:05
言語 PyPy3
(7.3.15)
結果
AC  
実行時間 92 ms / 10,000 ms
コード長 2,989 bytes
コンパイル時間 381 ms
コンパイル使用メモリ 82,828 KB
実行使用メモリ 77,248 KB
最終ジャッジ日時 2025-04-10 23:47:08
合計ジャッジ時間 2,535 ms
ジャッジサーバーID
（参考情報） 		judge2 / judge1
このコードへのチャレンジ
（要ログイン）
ファイルパターン 結果
other AC * 11
権限があれば一括ダウンロードができます

ソースコード

diff # 

import sys
import math

def solve():
    try:
        a_str = sys.stdin.readline().split()
        b_str = sys.stdin.readline().split()
        if not a_str or not b_str or len(a_str) != 8 or len(b_str) != 8:
            sys.exit(1)
    except:
        sys.exit(1)

    try:
        a = [float(s) for s in a_str]
        b = [float(s) for s in b_str]
    except ValueError:
        sys.exit(1)

    pts = []
    chg = []

    idx = 0
    for i in [0, 1]:
        for j in [0, 1]:
            for k in [0, 1]:
                pts.append((i/2.0, j/2.0, k/2.0))
                chg.append(a[idx])
                idx += 1

    idx = 0
    for i in [1, 3]:
        for j in [1, 3]:
            for k in [1, 3]:
                pts.append((i/4.0, j/4.0, k/4.0))
                chg.append(b[idx])
                idx += 1

    N = 16
    q0 = 0.0
    o_idx = 0
    tol = 1e-9
    if abs(pts[0][0]) < tol and abs(pts[0][1]) < tol and abs(pts[0][2]) < tol:
        q0 = chg[0]
    else:
        sys.exit(1)

    eta = 2.0
    r_cut = 7.0
    g_cut = 15.0

    v_r = 0.0
    r_cut_sq = r_cut**2
    n_max = int(math.ceil(r_cut)) + 1

    for i in range(-n_max, n_max + 1):
        for j in range(-n_max, n_max + 1):
            for k in range(-n_max, n_max + 1):
                nvx, nvy, nvz = float(i), float(j), float(k)
                for m in range(N):
                    p = pts[m]
                    q = chg[m]
                    if abs(q) < tol:
                        continue
                    if m == o_idx and i == 0 and j == 0 and k == 0:
                        continue
                    px = p[0] + nvx
                    py = p[1] + nvy
                    pz = p[2] + nvz
                    d_sq = px**2 + py**2 + pz**2
                    if d_sq < tol**2:
                        continue
                    if d_sq <= r_cut_sq:
                        d = math.sqrt(d_sq)
                        v_r += q * math.erfc(eta * d) / d

    v_k = 0.0
    pi = math.pi
    v_cell = 1.0
    g_max = int(math.sqrt(g_cut)) + 1

    for h in range(-g_max, g_max + 1):
        for k_ in range(-g_max, g_max + 1):
            for l in range(-g_max, g_max + 1):
                if h == 0 and k_ == 0 and l == 0:
                    continue
                hkl_sq = float(h*h + k_*k_ + l*l)
                if hkl_sq > g_cut:
                    continue
                g_sq = (2*pi)**2 * hkl_sq
                if g_sq < tol**2:
                    continue
                s_cos = 0.0
                for j in range(N):
                    p = pts[j]
                    q = chg[j]
                    if abs(q) < tol:
                        continue
                    g_dot = 2*pi * (h*p[0] + k_*p[1] + l*p[2])
                    s_cos += q * math.cos(g_dot)
                v_k += (4*pi / g_sq) * math.exp(-g_sq / (4*eta**2)) * s_cos

    v_s = 0.0
    if abs(q0) > tol:
        v_s = 2 * eta * q0 / math.sqrt(pi)

    v = v_r + v_k - v_s
    print(f"{v:.8f}")

solve()
0