transform_charts.py — src/26-pnp

1#!/usr/bin/env python32# -*- coding: utf-8 -*-3"""4transform_charts.py  (26-pnp, item 2 anchor: V1 + V3)5Representation-relativity made concrete on a structured (convolution) problem:6the same target is FAR in the position chart and NEAR in the Fourier chart,7because the Fourier transform diagonalises convolution and is forward-computable (FFT).8  V1  cyclic convolution: position chart O(n^2) vs Fourier chart O(n log n).9  V3  the diagonalising transform (FFT) is itself in P: T, T^{-1} poly and exact.10"""11import time, numpy as np12rng = np.random.default_rng(0)1314def direct_conv(a, b):                       # position chart: O(n^2)15    return np.convolve(a, b)1617def fourier_conv(a, b):                      # Fourier chart: O(n log n)18    n = len(a) + len(b) - 119    m = 1 << (n-1).bit_length()20    return np.fft.irfft(np.fft.rfft(a, m) * np.fft.rfft(b, m), m)[:n]2122print(f"{'n':>7} {'position (s)':>13} {'Fourier (s)':>12} {'pos/Fourier':>12} {'max|diff|':>11}")23ns, tpos, tfou = [], [], []24for n in (2048, 4096, 8192, 16384):25    a = rng.standard_normal(n); b = rng.standard_normal(n)26    t0=time.perf_counter(); c1=direct_conv(a,b); tp=time.perf_counter()-t027    t0=time.perf_counter(); c2=fourier_conv(a,b); tf=time.perf_counter()-t028    ns.append(n); tpos.append(tp); tfou.append(tf)29    print(f"{n:>7} {tp:>13.4f} {tf:>12.5f} {tp/tf:>12.1f} {np.max(np.abs(c1-c2)):>11.2e}")3031bp=np.polyfit(np.log(ns), np.log(tpos),1)[0]32bf=np.polyfit(np.log(ns), np.log(tfou),1)[0]33print(f"\nposition chart ~ n^({bp:.2f}) (the O(n^2) direct sum)")34print(f"Fourier  chart ~ n^({bf:.2f}) (the O(n log n) FFT) -- same target, near in one chart, far in the other")3536# V3: the transform itself is in P and exact (T and T^{-1} polynomial)37x = rng.standard_normal(1<<16)38err = np.max(np.abs(np.fft.ifft(np.fft.fft(x)).real - x))39print(f"\nFFT round-trip on n=2^16: max|T^-1 T x - x| = {err:.2e}  -- the diagonalising transform is forward-computable (P).")40