⚖️ Physics with Python · ฟิสิกส์ด้วย Python · 用 Python 学物理

1 🎯 A formula answers. A program argues back. · สูตรให้คำตอบ โปรแกรมเถียงกลับได้ · 公式给你答案，程序会跟你争

2 🪨 Drop a stone, two ways · ปล่อยก้อนหิน สองวิธี · 扔一块石头，用两种方法

terminal

python falling.py

falling.py

"""falling.py - drop a stone, two ways: with the formula, and step by step.

    python falling.py

The formula says where the stone is at any time. The simulation does not know the
formula: it only knows that gravity changes the speed, and speed changes the position,
and it repeats that in small steps.

They should agree. Where they disagree tells you something useful about simulations.
"""
G = 9.81          # metres per second, per second - how fast "how fast" changes
HEIGHT = 45.0     # a stone dropped from a tenth-floor balcony


def exact(t):
    """Where the stone is after t seconds, from the formula in the textbook."""
    return HEIGHT - 0.5 * G * t * t


def simulate(dt, stop):
    """Step through time in jumps of dt seconds. No formula, just repetition."""
    y, v, t = HEIGHT, 0.0, 0.0
    while t < stop - 1e-12:
        y = y - v * dt          # v is how fast it is falling, so y goes down
        v = v + G * dt          # then gravity makes it a little faster
        t = t + dt
    return y


if __name__ == "__main__":
    t = 2.0
    print(f"a stone dropped from {HEIGHT} m, where is it after {t} s?")
    print(f"  the formula says      {exact(t):8.4f} m above the ground")
    print()
    print("  step size     simulation      difference")
    for dt in (0.5, 0.1, 0.01, 0.001):
        y = simulate(dt, t)
        print(f"  {dt:>7} s   {y:8.4f} m     {y - exact(t):+8.4f} m")
    print()
    print("halve the step, halve the error. The simulation is not wrong,")
    print("it is approximate, and you can see exactly how approximate.")

what it printed

a stone dropped from 45.0 m, where is it after 2.0 s?
  the formula says       25.3800 m above the ground

  step size     simulation      difference
      0.5 s    30.2850 m      +4.9050 m
      0.1 s    26.3610 m      +0.9810 m
     0.01 s    25.4781 m      +0.0981 m
    0.001 s    25.3898 m      +0.0098 m

halve the step, halve the error. The simulation is not wrong,
it is approximate, and you can see exactly how approximate.

3 🎾 Now let the air push back · ทีนี้ให้อากาศดันกลับบ้าง · 现在让空气推回来

air.py

"""air.py - the same drop, but now the air pushes back.

    python air.py

Air resistance grows with speed. The faster the ball falls, the harder the air pushes,
until the push equals the weight and the ball stops speeding up. That speed has a name:
terminal velocity.

There is a formula for the final speed. There is no simple formula for the whole fall,
which is exactly why you simulate it.
"""
import math

G = 9.81
MASS = 0.057          # a tennis ball, in kilograms
RADIUS = 0.033        # metres
RHO = 1.2             # density of air, kg per cubic metre
DRAG = 0.5            # drag coefficient of a sphere, roughly
K = 0.5 * RHO * DRAG * math.pi * RADIUS ** 2      # the air-resistance constant
DT = 0.001


def terminal_speed():
    """The speed where the air pushes up exactly as hard as gravity pulls down."""
    return math.sqrt(MASS * G / K)


def fall(height, with_air=True):
    """Drop the ball and report (time to land, speed when it lands)."""
    y, v, t = height, 0.0, 0.0
    while y > 0:
        push = K * v * v if with_air else 0.0         # air pushes back, harder when fast
        a = G - push / MASS
        v = v + a * DT
        y = y - v * DT
        t = t + DT
    return t, v


if __name__ == "__main__":
    print(f"tennis ball: {MASS} kg, air constant K = {K:.6f}")
    print(f"terminal speed from the formula: {terminal_speed():.2f} m/s")
    print()
    print("  height      no air: time, speed        with air: time, speed")
    for h in (10, 45, 100, 500):
        t0, v0 = fall(h, with_air=False)
        t1, v1 = fall(h, with_air=True)
        print(f"  {h:>4} m     {t0:5.2f} s  {v0:6.2f} m/s        {t1:5.2f} s  {v1:6.2f} m/s")
    print()
    print("from 10 m the air hardly matters. From 500 m it decides everything:")
    print(f"the ball never goes faster than about {terminal_speed():.0f} m/s, however far it falls.")

what it printed

tennis ball: 0.057 kg, air constant K = 0.001026
terminal speed from the formula: 23.34 m/s

  height      no air: time, speed        with air: time, speed
    10 m      1.43 s   14.01 m/s         1.47 s   12.84 m/s
    45 m      3.03 s   29.71 m/s         3.45 s   20.91 m/s
   100 m      4.51 s   44.29 m/s         5.92 s   23.02 m/s
   500 m     10.10 s   99.04 m/s        23.07 s   23.34 m/s

from 10 m the air hardly matters. From 500 m it decides everything:
the ball never goes faster than about 23 m/s, however far it falls.

4 🎯 The 45 degree rule is wrong · กฎ 45 องศานั้นผิด · 45 度那条规则是错的

throw.py

"""throw.py - everyone is taught to throw at 45 degrees. Test it.

    python throw.py

With no air, 45 degrees really is the best angle, and you can prove it with algebra.
With air, the algebra stops being easy, so ask the computer instead: try every angle
and see which one actually goes furthest.
"""
import math

G = 9.81
MASS = 0.057
RADIUS = 0.033
K = 0.5 * 1.2 * 0.5 * math.pi * RADIUS ** 2
SPEED = 30.0          # metres per second, a hard throw
DT = 0.0005


def distance(angle_degrees, with_air=True):
    """How far the ball lands, thrown from the ground at this angle."""
    angle = math.radians(angle_degrees)
    x, y = 0.0, 0.0
    vx, vy = SPEED * math.cos(angle), SPEED * math.sin(angle)
    while y >= 0:
        speed = math.hypot(vx, vy)
        if with_air:
            ax = -K * speed * vx / MASS               # drag always opposes the motion
            ay = -G - K * speed * vy / MASS
        else:
            ax, ay = 0.0, -G
        vx, vy = vx + ax * DT, vy + ay * DT
        x, y = x + vx * DT, y + vy * DT
    return x


def best_angle(with_air):
    """Try every angle and keep the winner. No algebra required."""
    return max(range(20, 66), key=lambda a: distance(a, with_air))


if __name__ == "__main__":
    for air in (False, True):
        label = "with air" if air else "no air  "
        best = best_angle(air)
        print(f"{label}  best angle: {best} degrees, {distance(best, air):.2f} m")
    print()
    print("  angle    no air      with air")
    for a in (30, 35, 40, 45, 50, 55):
        print(f"  {a:>3}     {distance(a, False):6.2f} m    {distance(a, True):6.2f} m")
    print()
    print("45 degrees is the right answer to a question about a world with no air.")

what it printed

no air    best angle: 45 degrees, 91.74 m
with air  best angle: 39 degrees, 44.29 m

  angle    no air      with air
   30      79.45 m     42.43 m
   35      86.21 m     43.89 m
   40      90.34 m     44.29 m
   45      91.74 m     43.68 m
   50      90.34 m     42.12 m
   55      86.21 m     39.65 m

45 degrees is the right answer to a question about a world with no air.

5 ⏱ Where the pendulum formula gives up · จุดที่สูตรลูกตุ้มยอมแพ้ · 钟摆公式在哪里认输

pendulum.py

"""pendulum.py - the pendulum formula is a lie, and you can measure how big a lie.

    python pendulum.py

Every textbook gives the period of a pendulum as T = 2*pi*sqrt(L/g). Almost every
textbook then adds, quietly, "for small swings". This program measures how wrong
that formula gets as the swing grows.
"""
import math

G = 9.81
LENGTH = 1.0          # a one-metre pendulum
DT = 0.00002


def textbook_period():
    """The formula everyone memorises. It does not mention the angle at all."""
    return 2 * math.pi * math.sqrt(LENGTH / G)


def measured_period(start_degrees):
    """Let it swing and time it: release from rest, wait until it comes back."""
    theta = math.radians(start_degrees)
    omega, t = 0.0, 0.0
    # swing away from the start until it stops, that is half a period
    while True:
        alpha = -(G / LENGTH) * math.sin(theta)      # the real equation, no shortcut
        omega = omega + alpha * DT
        theta = theta + omega * DT
        t = t + DT
        if omega >= 0 and theta <= 0:                # back through the bottom, going up
            break
    return 2 * t


if __name__ == "__main__":
    t0 = textbook_period()
    print(f"a {LENGTH} m pendulum")
    print(f"  the textbook formula says every swing takes {t0:.4f} s,")
    print("  no matter how far you pull it back. Test that.")
    print()
    print("  pulled back     real period     the formula is wrong by")
    for angle in (1, 5, 10, 30, 60, 90, 120):
        t = measured_period(angle)
        print(f"  {angle:>4} deg       {t:.4f} s        {100 * (t - t0) / t0:+5.2f}%")
    print()
    print("at small angles the formula is excellent. At 90 degrees it is out by")
    print("almost a fifth, which a clock would notice within the hour.")

what it printed

a 1.0 m pendulum
  the textbook formula says every swing takes 2.0061 s,
  no matter how far you pull it back. Test that.

  pulled back     real period     the formula is wrong by
     1 deg       2.0061 s        +0.00%
     5 deg       2.0070 s        +0.05%
    10 deg       2.0099 s        +0.19%
    30 deg       2.0410 s        +1.74%
    60 deg       2.1529 s        +7.32%
    90 deg       2.3679 s        +18.04%
   120 deg       2.7541 s        +37.29%

at small angles the formula is excellent. At 90 degrees it is out by
almost a fifth, which a clock would notice within the hour.

6 👉 What to try next · ลองอะไรต่อดี · 接下来试什么