Godel's Solution To Einstien's GR

godel.py

# Create a 3D animation (MP4) that visualizes Gödel's rotating universe features:
# - Light-cone tipping along rings (as null-direction "wedges")
# - The null ring at r = r_c
# - A sample time-traveler worldline (helical path) that goes out, loops in the CTC region, and returns
# - Camera rotation for spatial intuition
#
# If ffmpeg is not available, the code will fall back to GIF.
#
# Outputs:
#   /mnt/data/godel_3d.mp4 (preferred)
#   /mnt/data/godel_3d.gif (fallback if MP4 fails)
#   /mnt/data/godel_3d_still.png (key frame)
import numpy as np
import math
import matplotlib.pyplot as plt
from matplotlib import animation
from mpl_toolkits.mplot3d import Axes3D  # noqa: F401
import os

# -------------------- Gödel-type helpers --------------------
def D_of_r(r, m):
    return np.sinh(m*r)/m

def H_of_r(r, m, omega):
    return (4.0*omega/(m**2)) * (np.sinh(0.5*m*r)**2)

def rc_analytic(omega):
    return float(np.arccosh(3.0) / (math.sqrt(2.0) * omega))

# Null slopes along a ring at radius r (dr = 0): dt/dphi = -H(r) ± D(r)
def null_slopes(r, m, omega):
    H = H_of_r(r, m, omega)
    D = D_of_r(r, m)
    return -H + D, -H - D

# -------------------- Scene construction --------------------
omega = 1.0
m = math.sqrt(2.0)*omega
rc = rc_analytic(omega)

# Time range for drawing (arbitrary units)
t_min, t_max = -4.0, 4.0

# Radii to place "cone wedges" (inside, at, outside r_c)
r_list = [0.5*rc, rc, 1.5*rc]

# Build traveler worldline:
# Segment A: move radially out (phi=0), t increases
# Segment B: at r=r1>rc, loop phi from 0 -> 2π with timelike slope dt/dphi = -H(r1) (between null slopes)
# Segment C: move radially back (phi=2π), small positive dt
r0 = 0.3*rc
r1 = 1.35*rc
# Segment A
Na = 120
ra = np.linspace(r0, r1, Na)
phia = np.zeros(Na)
ta = np.linspace(-1.0, 0.0, Na)
# Segment B
Nb = 260
phib = np.linspace(0.0, 2.0*np.pi, Nb)
rb = np.full(Nb, r1)
slope_mid = -H_of_r(r1, m, omega)  # choose the mid between the two null slopes
tb = np.linspace(ta[-1], ta[-1] + slope_mid*(phib[-1]-phib[0]), Nb)
# Segment C
Nc = 120
rc_seg = np.linspace(r1, r0, Nc)
phic = np.full(Nc, 2.0*np.pi)
tc = np.linspace(tb[-1], tb[-1]+0.5, Nc)

r_path = np.concatenate([ra, rb, rc_seg])
phi_path = np.concatenate([phia, phib, phic])
t_path = np.concatenate([ta, tb, tc])
x_path = r_path * np.cos(phi_path)
y_path = r_path * np.sin(phi_path)
z_path = t_path

# Precompute cylinder rings at r = r_c for a few t-levels
phi_ring = np.linspace(0.0, 2.0*np.pi, 200)
x_rc = rc * np.cos(phi_ring)
y_rc = rc * np.sin(phi_ring)
t_levels = np.linspace(t_min, t_max, 6)

# Precompute null wedge lines at selected radii
def wedge_lines_for_radius(r, t0=0.0, span=np.pi/8.0, samples=50):
    # two short lines representing the two null directions at phi0=0 over small |Δφ|≤span
    s_plus, s_minus = null_slopes(r, m, omega)
    dphi = np.linspace(-span, span, samples)
    # place wedges centered at phi0=0 and height t0
    # First wedge (s_plus)
    phi1 = dphi
    t1 = t0 + s_plus * dphi
    r1 = np.full_like(dphi, r)
    x1 = r1 * np.cos(phi1)
    y1 = r1 * np.sin(phi1)
    # Second wedge (s_minus)
    phi2 = dphi
    t2 = t0 + s_minus * dphi
    r2 = np.full_like(dphi, r)
    x2 = r2 * np.cos(phi2)
    y2 = r2 * np.sin(phi2)
    return (x1, y1, t1), (x2, y2, t2)

wedges = []
for idx, rr in enumerate(r_list):
    # stagger t0 so wedges don't overlap visually
    t0 = -2.0 + idx*2.0
    wedges.append(wedge_lines_for_radius(rr, t0=t0))

# -------------------- Animation setup --------------------
fig = plt.figure(figsize=(7, 6))
ax = fig.add_subplot(111, projection='3d')

# Set scene bounds
Rmax = 1.8*rc
ax.set_xlim(-Rmax, Rmax)
ax.set_ylim(-Rmax, Rmax)
ax.set_zlim(t_min, t_max)

ax.set_xlabel("x")
ax.set_ylabel("y")
ax.set_zlabel("t")

# Draw static elements
# Critical cylinder rings
for tz in t_levels:
    ax.plot(x_rc, y_rc, np.full_like(x_rc, tz), linewidth=1.0)

# Null wedges
wedge_artists = []
for (line1, line2) in wedges:
    x1, y1, t1 = line1
    x2, y2, t2 = line2
    h1, = ax.plot(x1, y1, t1, linewidth=2.0)
    h2, = ax.plot(x2, y2, t2, linewidth=2.0)
    wedge_artists.extend([h1, h2])

# Traveler worldline (animated)
path_line, = ax.plot([], [], [], linewidth=2.5)

# A null ring (closed photon curve) drawn at t=0 on r=rc
ax.plot(x_rc, y_rc, np.zeros_like(x_rc), linewidth=2.0)

# A horizontal slice (plane) hint via a thin ring grid (optional)
for rgrid in np.linspace(0.3*rc, 1.6*rc, 5):
    ax.plot(rgrid*np.cos(phi_ring), rgrid*np.sin(phi_ring), np.zeros_like(phi_ring), linewidth=0.75)

# Save a still
still_path = "/mnt/data/godel_3d_still.png"
plt.savefig(still_path, dpi=160, bbox_inches="tight")

# Animation function
total_frames = 360
def init():
    path_line.set_data([], [])
    path_line.set_3d_properties([])
    return [path_line] + wedge_artists

def animate(i):
    # Reveal the path progressively
    k = max(2, int((i+1)/total_frames * len(x_path)))
    path_line.set_data(x_path[:k], y_path[:k])
    path_line.set_3d_properties(z_path[:k])
    # Rotate camera
    az = 30 + 360.0 * (i / total_frames)
    el = 22 + 5.0*np.sin(2*np.pi * i/total_frames)
    ax.view_init(elev=el, azim=az)
    return [path_line] + wedge_artists

anim = animation.FuncAnimation(fig, animate, init_func=init,
                               frames=total_frames, interval=30, blit=True)

mp4_path = "/mnt/data/godel_3d.mp4"
gif_path = "/mnt/data/godel_3d.gif"

saved = {}
try:
    Writer = animation.FFMpegWriter
    writer = Writer(fps=30, metadata=dict(artist="UDK"), bitrate=2400)
    anim.save(mp4_path, writer=writer, dpi=160)
    saved["mp4"] = mp4_path
except Exception as e:
    # Fallback to GIF if ffmpeg is unavailable
    try:
        from matplotlib.animation import PillowWriter
        anim.save(gif_path, writer=PillowWriter(fps=20))
        saved["gif"] = gif_path
    except Exception as e2:
        saved["error"] = f"Failed to save MP4 and GIF: {e}; {e2}"

saved, still_path