結果

問題 No.3672 Volume 3D
コンテスト
ユーザー harurun
提出日時 2026-08-05 03:33:27
言語 Python3
(3.14.7 + numpy 2.5.2 + scipy 1.18.0 + ACL)
コンパイル:
python3 -mpy_compile _filename_
実行:
python3 _filename_
結果
AC  
実行時間 889 ms / 2,000 ms
+ 249µs
コード長 2,581 bytes
記録
記録タグの例:
初AC ショートコード 純ショートコード 純主流ショートコード 最速実行時間
コンパイル時間 572 ms
コンパイル使用メモリ 21,472 KB
実行使用メモリ 117,944 KB
最終ジャッジ日時 2026-09-04 22:02:25
合計ジャッジ時間 7,597 ms
ジャッジサーバーID
(参考情報)
judge2_0 / judge1_0
このコードへのチャレンジ
(要ログイン)
ファイルパターン 結果
sample AC * 1
other AC * 20
権限があれば一括ダウンロードができます

ソースコード

diff #
raw source code

#!/usr/bin/env python3
from __future__ import annotations

import sys
from decimal import Decimal, getcontext, localcontext
from typing import Iterable, Sequence

PRECISION = 90
getcontext().prec = PRECISION


def compute_pi() -> Decimal:
    """Gauss-Legendre法でDecimalの円周率を計算する。"""
    with localcontext() as ctx:
        ctx.prec = PRECISION + 20
        one = Decimal(1)
        two = Decimal(2)
        four = Decimal(4)

        a = one
        b = one / two.sqrt()
        t = one / four
        p = one

        for _ in range(9):
            next_a = (a + b) / two
            b = (a * b).sqrt()
            t -= p * (a - next_a) * (a - next_a)
            a = next_a
            p *= two

        return +((a + b) * (a + b) / (four * t))


PI = compute_pi()
D4 = Decimal(4)
D12 = Decimal(12)
D3 = Decimal(3)
D2 = Decimal(2)


def intersection_volume(values: Sequence[int]) -> Decimal:
    xa, ya, za, ra_i, xb, yb, zb, rb_i = values

    dx = xa - xb
    dy = ya - yb
    dz = za - zb
    d2_i = dx * dx + dy * dy + dz * dz
    sum_i = ra_i + rb_i
    diff_i = abs(ra_i - rb_i)

    if d2_i >= sum_i * sum_i:
        return Decimal(0)

    if d2_i <= diff_i * diff_i:
        r = Decimal(min(ra_i, rb_i))
        return +(D4 * PI * r * r * r / D3)

    ra = Decimal(ra_i)
    rb = Decimal(rb_i)
    d2 = Decimal(d2_i)
    d = d2.sqrt()
    sum_r = ra + rb

    # sum_r - dを有理化し、外接直前の桁落ちを避ける。
    overlap_depth = (sum_r * sum_r - d2) / (sum_r + d)
    radius_difference = ra - rb
    bracket = (
        d2
        + D2 * d * sum_r
        - D3 * radius_difference * radius_difference
    )
    return +(PI * overlap_depth * overlap_depth * bracket / (D12 * d))


def solve_tokens(tokens: Iterable[bytes]) -> str:
    values = [int(token) for token in tokens]
    if not values:
        raise ValueError("empty input")

    t = values[0]
    if not 1 <= t <= 100000:
        raise ValueError(f"T is out of range: {t}")
    expected = 1 + 8 * t
    if len(values) != expected:
        raise ValueError(f"expected {expected} integers, got {len(values)}")

    output: list[str] = []
    offset = 1
    for _ in range(t):
        answer = intersection_volume(values[offset:offset + 8])
        output.append(f"{answer:.60E}")
        offset += 8
    return "\n".join(output) + "\n"


def solve_text(text: str) -> str:
    return solve_tokens(token.encode() for token in text.split())


def main() -> None:
    sys.stdout.write(solve_tokens(sys.stdin.buffer.read().split()))


if __name__ == "__main__":
    main()
0