2.3 三角形與畢氏定理：用 SymPy 驗證（公元前 500 年 / 公元 2020 年）#

import sympy
from sympy import Integer

triples = [(3, 4, 5), (5, 12, 13), (8, 15, 17), (119, 120, 169), (3367, 3456, 4825)]
for a, b, c in triples:
    a, b, c = Integer(a), Integer(b), Integer(c)
    print(a**2 + b**2 == c**2)
# 預期輸出：True True True True True（Plimpton 322 的數據精确成立）

import sympy
from sympy import Point, Matrix, Polygon
from sympy.geometry import Triangle

def perp_offset(P, Q, d, side=1):
    """從線段 PQ 出發，沿其法向（side=1 一側，side=-1 另一側）擴張 d。"""
    v = Matrix([Q.x - P.x, Q.y - P.y])
    n = Matrix([-v[1], v[0]]) * side      # 法向量 = 旋轉 90°
    n = n / n.norm() * d
    return Point(Q.x + n[0], Q.y + n[1]), Point(P.x + n[0], P.y + n[1])

def square_on(P, Q, side=1):
    """以 PQ 為一邊的正方形，其餘兩頂點向外（side=±1）。"""
    R, S = perp_offset(P, Q, P.distance(Q), side)
    return Polygon(P, Q, R, S)

A = Point(0, 0); B = Point(3, 0)
C = Point(0, 4)   # 直角在 A

# 三個正方形
sq_c = square_on(B, C, side=-1)   # 邊 c（斜邊）
sq_a = square_on(A, B, side=-1)   # 邊 a = 3
sq_b = square_on(A, C, side=1)    # 邊 b = 4（注意方向）

print(sq_a.area)      # 預期輸出：9
print(sq_b.area)      # 預期輸出：16
print(sympy.simplify(sq_c.area))    # 預期輸出：25
print(sq_a.area + sq_b.area == sq_c.area)
# 預期輸出：True —— I.47 對 (3,4,5) 成立

import sympy
from sympy import Point, symbols, simplify

u, v = symbols('u v', positive=True)   # 兩直角邊長
A = Point(0, 0)
B = Point(u, 0)
C = Point(0, v)
a2 = A.distance(B) ** 2
b2 = A.distance(C) ** 2
c2 = B.distance(C) ** 2
print(sympy.simplify(c2), sympy.simplify(a2 + b2), sympy.simplify(c2 - a2 - b2) == 0)
# 預期輸出：u**2 + v**2  u**2 + v**2  True

import sympy
from sympy import Point, Rational, simplify
from sympy.geometry import Triangle

# 直角三角形，直角在 A=(0,0)：B=(4,0)、C=(0,3)
A = Point(0, 0)
B = Point(4, 0)
C = Point(0, 3)

# 頂點 A 到斜邊 BC 的高，腳點 D
BC = sympy.Line(B, C)
D = sympy.Line(A, A + __import__('sympy').Point(0,1)).intersection(BC)

import sympy
from sympy import Point, Rational, simplify
from sympy.geometry import Triangle, Line

A = Point(0, 0)
B = Point(4, 0)
C = Point(0, 3)

# 斜邊 BC 的直線；過 A 作 BC 的垂線，取交點 D
BC = Line(B, C)
alt = BC.perpendicular_line(A)
D = alt.intersection(BC)[0]                 # 交點（高腳）
print(D)                                    # 預期輸出：Point2D(48/25, 36/25)

# 三個三角形：原三角形 + 兩個小三角形
t0 = Triangle(A, B, C)
t1 = Triangle(A, D, B)
t2 = Triangle(A, C, D)

# 相似性驗證：三對互相似
print(t0.is_similar(t1), t0.is_similar(t2), t1.is_similar(t2))
# 預期輸出：True True True

# 比例式驗證：AB^2 = BD * BC，AC^2 = CD * BC，兩式相加得 BC^2
BD = B.distance(D)
DC = C.distance(D)
BC_ = B.distance(C)
ratio1 = B.distance(A) ** 2
ratio2 = A.distance(C) ** 2
print(simplify(BD * BC_ - ratio1))    # 預期輸出：0
print(simplify(DC * BC_ - ratio2))    # 預期輸出：0
print(simplify(BD * BC_ + DC * BC_ - BC_ ** 2))
# 預期輸出：0 —— 因 BD + DC = BC，故 a^2 + b^2 = c^2。證畢。

import sympy
from sympy import symbols, simplify, Rational

a, b = symbols('a b', positive=True)
little = (b - a) ** 2          # 中間小正方形
four_tri = 4 * Rational(1, 2) * a * b   # 四個三角形
total = simplify(little + four_tri)     # 大正方形 c^2
print(total)                            # 預期輸出：a**2 + b**2

import sympy
from sympy import symbols, simplify, Eq, solve

a, b, c = symbols('a b c', positive=True)
trap_area = (a + b) * (a + b) / 2          # 梯形面積
tri_sum  = a * b / 2 + a * b / 2 + c ** 2 / 2   # 兩繞射三角形+ 直角三角形
rel = simplify(Eq(trap_area, tri_sum))

diff = simplify(trap_area - tri_sum)     # = (a**2 + b**2 - c**2)/2
print(diff)
from sympy import Eq
print(Eq(diff, 0))                        # 預期輸出：a**2 + b**2 == c**2

import sympy
from sympy import Point, Rational, simplify, sqrt
from sympy.geometry import Triangle

# 直角三角形 3-4-5，直角在 A
t = Triangle(Point(0, 0), Point(4, 0), Point(0, 3))

print(t.is_right())                 # 預期輸出：True
print(t.is_scalene())               # 預期輸出：True（不等腰）
print(t.is_isosceles())             # 預期輸出：False
print(t.area)                       # 預期輸出：6
print(t.perimeter)                  # 預期輸出：12
print(t.inradius)                   # 預期輸出：1（內切圓半徑）
print(t.circumradius)               # 預期輸出：5/2（外接圓半徑）

# 中線（medians）— 對應 Euclid I.以後的結構輔助
m = t.medians[Point(0, 0)]
print(m)                            # 預期輸出：Segment2D(Point2D(0, 0), Point2D(2, 3/2))

# 角：t.angles 是「頂點 → 角」的字典
angles = t.angles
for k, v in angles.items():
    print(k, v)
# 預期輸出：
# Point2D(0, 3) pi/2
# Point2D(0, 0) pi/2  ← 對！！直角在 A = (0,0)？

print(angles)
# {Point2D(0, 3): ..., Point2D(0, 0): pi/2, Point2D(4, 0): ...}

import sympy
from sympy import pi
from sympy.geometry import Triangle

def vertex_of_right_angle(t):
    """回傳哪個頂點的角等於 pi/2"""
    if not t.is_right():
        return None
    for v, ang in t.angles.items():
        if sympy.simplify(ang - pi) == 0 or sympy.simplify(ang - pi/2) == 0:
            if sympy.simplify(ang - pi/2) == 0:
                return v
    return None

t = Triangle(Point(0, 0), Point(4, 0), Point(0, 3))
print(vertex_of_right_angle(t))     # 預期輸出：Point2D(0, 0)

import sympy
from sympy import sqrt, Rational, pi
from sympy.geometry import Triangle

e = Triangle(Point(0, 0), Point(1, 0), Point(Rational(1, 2), sqrt(3)/2))
print(e.is_equilateral())                  # 預期輸出：True
print(all(sympy.simplify(v - pi/3) == 0 for v in e.angles.values()))
# 預期輸出：True（每個角 = pi/3）
print(e.area)                              # 預期輸出：sqrt(3)/4

import sympy
from sympy import Point, simplify
from sympy.geometry import Triangle

t = Triangle(Point(0, 0), Point(4, 0), Point(0, 3))

# 中線：頂點到對邊中點
print(t.medians[Point(4, 0)])
# 預期輸出：Segment2D(Point2D(4, 0), Point2D(0, 3/2))

# 高：BC.perpendicular_line(A) 也建基於 perp，Triangle 有 altitude 屬性
alt = t.altitude(0)     # 第 0 個頂點 A 的高
print(alt)
# 這回傳 ABC 頂點 A 的「高」對象其實依版本不同。更可靠的方式：

from sympy.geometry import Line
mBC = Line(Point(4, 0), Point(0, 3))
hA = mBC.perpendicular_line(Point(0, 0))
print(hA.intersection(mBC)[0])
# 預期輸出：Point2D(48/25, 36/25) —— 與 3.2 節 D 相同

import math
import sympy

a, b = 3e200, 4e200
# 樸素畢氏：溢位為 inf
try:
    print(math.sqrt(a*a + b*b))
except OverflowError:
    print("Overflow!")     # 頂部的浮點工業表現：3e200**2 = inf
# 穩健版本
print(math.hypot(a, b))    # 預期輸出：5e+200
# 精確算術不會有這問題
print(sympy.sqrt(sympy.Integer(3*10**200)**2 + sympy.Integer(4*10**200)**2))
# 預期輸出：5000000000000000000000000000000...（200 位精確）

import sympy
from sympy import Point, simplify, pi
from sympy.geometry import Triangle

def pythagorean_report(A, B, C):
    t = Triangle(A, B, C)
    s1, s2, s3 = A.distance(B), B.distance(C), A.distance(C)
    subs = {s1: 'a', s2: 'b', s3: 'c'}
    print("邊長:", simplify(s1), simplify(s2), simplify(s3))
    # 逐一檢查 c^2 = a^2 + b^2 的三種排列
    ok = []
    for hyp, others in [
        (s3, (s1, s2)), (s2, (s1, s3)), (s1, (s2, s3)),
    ]:
        lhs = hyp ** 2
        rhs = others[0] ** 2 + others[1] ** 2
        ok.append(sympy.simplify(lhs - rhs) == 0)
    print("畢氏成立?", True in ok)
    if t.is_right():
        for v, ang in t.angles.items():
            if sympy.simplify(ang - pi/2) == 0:
                print("直角頂點:", v)
    return True in ok

print(pythagorean_report(Point(0,0), Point(4,0), Point(0,3)))
# 預期輸出：
# 邊長: 4 5 3
# 畢氏成立? True
# 直角頂點: Point2D(0, 0)
# True