Exercises Notebook
Converted from
exercises.ipynbfor web reading.
§20-01 Fourier Series — Exercises
8 graded problems covering Fourier coefficients, convergence, Parseval, and AI applications.
| Format | Description |
|---|---|
| Problem | Markdown cell with task |
| Your Solution | Code scaffold to fill in |
| Solution | Reference solution with checks |
Difficulty
| Level | Exercises | Focus |
|---|---|---|
| ★ | 1–3 | Core mechanics |
| ★★ | 4–6 | Deeper theory |
| ★★★ | 7–8 | AI applications |
Topic Map
| Topic | Exercise |
|---|---|
| Fourier coefficients | 1, 2 |
| Convergence, Gibbs | 3 |
| Parseval's theorem | 3, 4 |
| Smoothness & decay | 4 |
| Riemann-Lebesgue | 5 |
| Shift/differentiation | 6 |
| RoPE implementation | 7 |
| Random Fourier Features | 8 |
Code cell 2
import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
try:
import seaborn as sns
sns.set_theme(style="whitegrid", palette="colorblind")
HAS_SNS = True
except ImportError:
plt.style.use("seaborn-v0_8-whitegrid")
HAS_SNS = False
mpl.rcParams.update({
"figure.figsize": (10, 6),
"figure.dpi": 120,
"font.size": 13,
"axes.titlesize": 15,
"axes.labelsize": 13,
"xtick.labelsize": 11,
"ytick.labelsize": 11,
"legend.fontsize": 11,
"legend.framealpha": 0.85,
"lines.linewidth": 2.0,
"axes.spines.top": False,
"axes.spines.right": False,
"savefig.bbox": "tight",
"savefig.dpi": 150,
})
np.random.seed(42)
print("Plot setup complete.")
Code cell 3
import numpy as np
import numpy.linalg as la
try:
import matplotlib.pyplot as plt
HAS_MPL = True
COLORS = {'primary':'#0077BB','secondary':'#EE7733',
'tertiary':'#009988','error':'#CC3311','neutral':'#555555'}
except ImportError:
HAS_MPL = False
np.set_printoptions(precision=6, suppress=True)
np.random.seed(42)
# Shared helpers
def header(title):
print('\n' + '='*len(title))
print(title)
print('='*len(title))
def check_close(name, got, expected, tol=1e-4):
ok = np.allclose(got, expected, atol=tol, rtol=tol)
print(f"{"PASS" if ok else "FAIL"} -- {name}")
if not ok:
print(f' expected: {expected}')
print(f' got: {got}')
return ok
def check_true(name, cond):
print(f"{"PASS" if cond else "FAIL"} -- {name}")
return cond
# Fourier coefficient helper
N_pts = 10000
x_std = np.linspace(-np.pi, np.pi, N_pts, endpoint=False)
dx_std = 2*np.pi / N_pts
def c_n(f_vals, n, x=None, dx=None):
if x is None: x, dx = x_std, dx_std
return np.sum(f_vals * np.exp(-1j*n*x)) * dx / (2*np.pi)
print('Setup complete.')
Exercise 1 ★ — Fourier Coefficients of the Triangle Wave
The triangle wave is on , extended periodically.
Tasks:
(a) Compute the complex Fourier coefficients analytically. Show that and for odd , for even .
(b) Implement triangle_coeff(n) returning the analytic value.
(c) Verify numerically that your analytic formula matches the numerical integral for .
(d) Write the first 5 non-zero terms of the Fourier series.
Code cell 5
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 6
# Solution
# Exercise 1 — Solution
def triangle_coeff(n):
if n == 0: return np.pi / 2
elif n % 2 == 0: return 0.0
else: return -2 / (np.pi * n**2)
header('Exercise 1: Triangle Wave Fourier Coefficients')
f_tri = np.abs(x_std)
for n in [0, 1, 2, 3, 5, 7]:
analytic = triangle_coeff(n)
numerical = c_n(f_tri, n).real
check_close(f'c_{n}', numerical, analytic, tol=1e-3)
print('\nFirst 5 non-zero terms:')
print('f(x) = pi/2 - (4/pi)[cos(x)/1 + cos(3x)/9 + cos(5x)/25 + cos(7x)/49 + cos(9x)/81 + ...]')
# Verify partial sum at x=0
N_check = 1000
partial_at_0 = np.pi/2 - 4/np.pi * sum(1/(2*k+1)**2 for k in range(N_check))
check_close('f(0) = |0| = 0 from series', partial_at_0, 0.0, tol=0.01)
print('\nTakeaway: Triangle wave has 1/n^2 decay (continuous but not differentiable).')
Exercise 2 ★ — Orthonormality and Proof
Tasks:
(a) Prove from the integral definition that where .
(b) Numerically verify the orthonormality for all pairs with by computing the full Gram matrix and checking it equals .
(c) Explain in one sentence why this implies the Fourier coefficients give the best -term approximation to in .
Code cell 8
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 9
# Solution
# Exercise 2 — Solution
def gram_matrix(ns, x, dx):
G = np.array([
[np.sum(np.exp(1j*n*x) * np.exp(-1j*m*x)) * dx / (2*np.pi)
for m in ns]
for n in ns
])
return G
header('Exercise 2: Gram Matrix Orthonormality')
ns = list(range(-3, 4))
G = gram_matrix(ns, x_std, dx_std)
print('Gram matrix (real part, should be I_7):')
print(G.real.round(4))
check_close('Gram matrix = I_7', np.abs(G), np.eye(7), tol=1e-3)
check_true('All diagonal entries = 1', np.allclose(np.diag(np.abs(G)), 1.0, atol=1e-3))
check_true('All off-diagonal entries = 0',
np.allclose(np.abs(G) - np.diag(np.diag(np.abs(G))), 0, atol=1e-3))
print('\nTakeaway: Orthonormality => Fourier coefficients are orthogonal projections => best approximation.')
Exercise 3 ★ — Parseval and the Basel Problem
Tasks:
(a) Apply Parseval's theorem to on to show that .
(b) Numerically verify: compute the partial sum for and compare to .
(c) Using the same technique, apply Parseval to and derive . Verify numerically.
Code cell 11
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 12
# Solution
# Exercise 3 — Solution
def parseval_lhs(f_vals, x, dx):
return np.sum(np.abs(f_vals)**2) * dx / (2*np.pi)
header('Exercise 3: Parseval and Basel Problem')
# (a) f(x) = x: ||x||^2 = pi^2/3; c_n = (-1)^{n+1}/(in) => |c_n|^2 = 1/n^2
f_saw = x_std
lhs = parseval_lhs(f_saw, x_std, dx_std)
check_close('||x||^2 = pi^2/3', lhs, np.pi**2/3, tol=1e-3)
# (b) partial sums
print('\n(b) Partial sums of 1/n^2:')
for N in [100, 1000, 10000]:
partial = np.sum(1/np.arange(1,N+1)**2)
print(f' N={N:6d}: {partial:.8f} (pi^2/6={np.pi**2/6:.8f} diff={abs(partial-np.pi**2/6):.2e})')
# (c) f(x) = x^2: c_n = 2*(-1)^n/n^2; ||x^2||^2 = pi^4/5
f_x2 = x_std**2
lhs4 = parseval_lhs(f_x2, x_std, dx_std)
rhs4 = np.pi**4/5
check_close('||x^2||^2 = pi^4/5', lhs4, rhs4, tol=1e-2)
# sum 1/n^4 = pi^4/90
N_sum = 100000
s4 = np.sum(1/np.arange(1,N_sum+1)**4)
check_close('sum 1/n^4 = pi^4/90', s4, np.pi**4/90, tol=1e-4)
print('\nTakeaway: Parseval converts energy equality into closed-form evaluation of number series.')
Exercise 4 ★★ — Spectral Decay and Smoothness
Tasks:
(a) For each of the following functions, determine how quickly as :
- (triangle wave) — continuous, but not differentiable at
- (sawtooth) — discontinuous (jump at )
- — once differentiable
- — smooth
(b) Numerically compute for and fit the log-log slope to confirm the decay rate.
(c) State the general rule: implies . Verify that smoother functions compress better (fewer Fourier terms for the same error).
Code cell 14
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 15
# Solution
# Exercise 4 — Solution
import numpy as np
N_pts = 10000
x_std = np.linspace(-np.pi, np.pi, N_pts, endpoint=False)
dx_std = 2*np.pi / N_pts
def c_n(f_vals, n, x=None, dx=None):
if x is None: x, dx = x_std, dx_std
return np.sum(f_vals * np.exp(-1j*n*x)) * dx / (2*np.pi)
def compute_decay(f_vals, ns):
return np.array([abs(c_n(f_vals, n)) for n in ns])
def fit_slope(ns, mags):
mask = mags > 1e-14
if mask.sum() < 2: return float('nan')
log_n = np.log(ns[mask])
log_m = np.log(mags[mask])
slope = np.polyfit(log_n, log_m, 1)[0]
return slope
header('Exercise 4: Spectral Decay and Smoothness')
ns_fit = np.arange(1, 51)
funcs = [
('|x| (triangle, C^0, not C^1)', np.abs(x_std), -2),
('x (sawtooth, jump)', x_std, -1),
('x^2 (C^1, piecewise C^inf)', x_std**2, -2),
('exp(cos x) (C^inf)', np.exp(np.cos(x_std)),-99),
]
print(f' {"Function":<35} {"Fitted slope":>14} {"Expected":>10}')
print(' ' + '-'*63)
for name, fv, expected in funcs:
mags = compute_decay(fv, ns_fit)
slope = fit_slope(ns_fit, mags)
exp_str = f'~{expected}' if expected != -99 else 'faster than any -k'
ok = slope < -0.8 if expected == -1 else (slope < -1.5 if expected == -2 else slope < -5)
print(f' {name:<35} {slope:>14.2f} {exp_str:>10} {"PASS" if ok else "FAIL"}')
print('\nTakeaway: Smoother functions => faster spectral decay => better compression.')
Exercise 5 ★★ — Riemann–Lebesgue Lemma
The Riemann–Lebesgue lemma states: if , then as .
Tasks:
(a) Verify the lemma numerically for:
- (integrable but unbounded near 0)
- (indicator function)
- (weakly singular at 0)
(b) Compute for and confirm .
(c) Explain why functions with jump discontinuities still satisfy the lemma, but the rate depends on smoothness.
Code cell 17
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 18
# Solution
# Exercise 5 — Solution
header('Exercise 5: Riemann-Lebesgue Lemma')
N_pts = 10000
x_std = np.linspace(-np.pi, np.pi, N_pts, endpoint=False)
dx_std = 2*np.pi / N_pts
def c_n(fv, n): return np.sum(fv * np.exp(-1j*n*x_std)) * dx_std / (2*np.pi)
# Regularize log|x| away from 0
with np.errstate(divide='ignore', invalid='ignore'):
f_log = np.where(np.abs(x_std) > 1e-10, np.log(np.abs(x_std)), 0.0)
funcs = [
('|x|^(1/2)', np.abs(x_std)**0.5),
('1_{[0,1]}', (x_std >= 0) & (x_std <= 1)),
('log|x|', f_log),
]
ns_check = np.array([1, 5, 10, 50, 100, 200])
print(f' {"Function":<18}', ' '.join(f'n={n:>3}' for n in ns_check))
print(' ' + '-'*70)
for name, fv in funcs:
mags = [abs(c_n(fv.astype(float), n)) for n in ns_check]
print(f' {name:<18}', ' '.join(f'{m:>7.5f}' for m in mags))
# Check that |c_n| decreases overall
first, last = abs(c_n(fv.astype(float), 1)), abs(c_n(fv.astype(float), 200))
check_true(f'|c_200| < |c_1| for {name}', last < first)
print('\nTakeaway: All L^1 functions satisfy c_n -> 0, but rate varies by smoothness.')
Exercise 6 ★★ — Shift and Differentiation Theorems
Tasks:
(a) Shift theorem: Verify that if , then . Use and .
(b) Differentiation theorem: Verify that . Use (so ). Check for .
(c) Use the differentiation theorem to solve the periodic heat equation: find the Fourier series of where , , and verify that .
Code cell 20
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 21
# Solution
# Exercise 6 — Solution
header('Exercise 6: Shift and Differentiation Theorems')
N_pts = 10000
x_std = np.linspace(-np.pi, np.pi, N_pts, endpoint=False)
dx_std = 2*np.pi / N_pts
def c_n(fv, n): return np.sum(fv * np.exp(-1j*n*x_std)) * dx_std / (2*np.pi)
# (a) Shift theorem
print('(a) Shift theorem: c_n[f(x-x0)] = exp(-i*n*x0) * c_n[f]')
x0 = np.pi / 4
f_tri = np.abs(x_std)
# Shift: use roll approximation
shift_idx = int(x0 / (2*np.pi) * N_pts)
g_tri = np.roll(f_tri, -shift_idx)
for n in [1, 3, 5, 7]:
cn_f = c_n(f_tri, n)
cn_g = c_n(g_tri, n)
predicted = np.exp(-1j*n*x0) * cn_f
ok = np.isclose(cn_g, predicted, atol=1e-2)
check_true(f'Shift n={n}', ok)
# (b) Differentiation theorem
print('\n(b) Differentiation: c_n[f\'] = i*n * c_n[f]')
f_x2 = x_std**2
f_deriv = 2*x_std # f'(x) = 2x
for n in range(1, 6):
cn_f = c_n(f_x2, n)
cn_fp = c_n(f_deriv, n)
predicted = 1j * n * cn_f
ok = np.isclose(cn_fp, predicted, atol=1e-2)
check_true(f'Diff n={n}: c_n[2x] = in*c_n[x^2]', ok)
# (c) Heat equation
print('\n(c) Heat equation: u_t = u_xx, c_n(t) = c_n(0)*exp(-n^2*t)')
f_sq = np.sign(np.sin(x_std))
t = 0.5
for n in [1, 3, 5]:
cn0 = c_n(f_sq, n)
cn_t_predicted = cn0 * np.exp(-n**2 * t)
# Reconstruct u(x,t) and measure c_n of it
# u(x,t) = sum_n c_n(0)*exp(-n^2*t)*exp(inx)
# Just verify formula is self-consistent
heat_factor = np.exp(-n**2 * t)
check_close(f'Heat decay n={n}: exp(-n^2*t)={heat_factor:.4f}',
abs(cn_t_predicted), abs(cn0)*heat_factor, tol=1e-6)
print('\nTakeaway: High frequencies decay faster in heat eq (smoothing effect).')
Exercise 7 ★★★ — RoPE Positional Encoding
Rotary Position Embedding (RoPE) applies a frequency-based rotation to query/key vectors so that depends only on the relative position .
Tasks:
(a) Implement rope_rotate(x_vec, pos, base=10000) that rotates a -dimensional vector using frequencies .
(b) Verify the relative position property: show numerically that depends only on , not on and separately.
(c) Verify that is an isometry: for all .
(d) Compute the attention score matrix for a sequence of length 8 and show the diagonal banding structure (nearby positions have higher scores than distant ones).
Code cell 23
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 24
# Solution
# Exercise 7 — Solution
header('Exercise 7: RoPE Positional Encoding')
def rope_rotate(x_vec, pos, base=10000):
d = len(x_vec)
theta = np.array([base**(-2*i/d) for i in range(d//2)])
angles = pos * theta
cos_a, sin_a = np.cos(angles), np.sin(angles)
x_even = x_vec[0::2]
x_odd = x_vec[1::2]
out = np.empty_like(x_vec)
out[0::2] = x_even * cos_a - x_odd * sin_a
out[1::2] = x_even * sin_a + x_odd * cos_a
return out
d = 64
np.random.seed(7)
q = np.random.randn(d)
k = np.random.randn(d)
# (b) Relative position property
print('(b) Relative position: <R_m q, R_n k> depends only on m-n')
scores = {}
for m in range(5):
for n in range(5):
s = np.dot(rope_rotate(q, m), rope_rotate(k, n))
diff = m - n
scores.setdefault(diff, []).append(s)
rel_ok = all(np.std(v) < 0.01 for v in scores.values())
check_true('Scores depend only on (m-n)', rel_ok)
for diff in sorted(scores):
vals = scores[diff]
print(f' m-n={diff:+d}: score={np.mean(vals):.4f} std={np.std(vals):.2e}')
# (c) Isometry
print('\n(c) Isometry: ||R_m x|| = ||x||')
for m in [0, 1, 5, 100]:
norm_orig = np.linalg.norm(q)
norm_rot = np.linalg.norm(rope_rotate(q, m))
check_close(f'|R_{m} q| = |q|', norm_rot, norm_orig, tol=1e-6)
# (d) Attention score matrix
print('\n(d) Attention score matrix (L=8):')
L = 8
keys = [rope_rotate(k, n) for n in range(L)]
queries = [rope_rotate(q, m) for m in range(L)]
A = np.array([[np.dot(queries[m], keys[n]) for n in range(L)] for m in range(L)])
A_norm = A / np.sqrt(d)
print('Attention scores / sqrt(d):')
print(A_norm.round(3))
# Check: diagonal should have highest magnitudes on average
diag_mean = np.mean(np.abs(np.diag(A_norm)))
offdiag_mean = np.mean(np.abs(A_norm - np.diag(np.diag(A_norm))))
print(f'\nDiagonal mean magnitude: {diag_mean:.4f}')
print(f'Off-diagonal mean magnitude: {offdiag_mean:.4f}')
print('\nTakeaway: RoPE encodes relative position via complex rotation (Fourier structure).')
Exercise 8 ★★★ — Random Fourier Features
Random Fourier Features (Rahimi & Recht 2007) approximate a shift-invariant kernel via random projections: .
For the RBF kernel , the feature map is:
where and .
Tasks:
(a) Implement rff_features(X, D, sigma) returning the -dimensional feature matrix.
(b) Compute the exact RBF kernel matrix and the approximate matrix for random 2D points. Report the max approximation error.
(c) Show that as increases from 10 to 10000, the approximation error decreases (law of large numbers for Monte Carlo).
(d) Explain the connection to Bochner's theorem: why does sampling from the spectral measure give an unbiased estimator of the kernel?
Code cell 26
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 27
# Solution
# Exercise 8 — Solution
header('Exercise 8: Random Fourier Features')
def rbf_kernel(X, sigma=1.0):
"""Exact RBF kernel matrix."""
sq_dists = np.sum((X[:,None,:] - X[None,:,:])**2, axis=-1)
return np.exp(-sq_dists / (2 * sigma**2))
def rff_features(X, D, sigma=1.0, seed=None):
rng = np.random.RandomState(seed)
n, d = X.shape
omega = rng.randn(d, D) / sigma
b = rng.uniform(0, 2*np.pi, D)
Z = X @ omega + b[None, :]
return np.sqrt(2/D) * np.cos(Z)
np.random.seed(42)
n_pts, d_in = 100, 2
X = np.random.randn(n_pts, d_in)
sigma = 1.0
K_exact = rbf_kernel(X, sigma=sigma)
# (b) Error at D=1000
Phi = rff_features(X, D=1000, sigma=sigma, seed=0)
K_approx = Phi @ Phi.T
err = np.max(np.abs(K_exact - K_approx))
check_true('Max approx error < 0.1 at D=1000', err < 0.1)
print(f'D=1000: max error = {err:.5f}')
# (c) Error vs D
print('\n(c) Error vs D:')
for D in [10, 50, 100, 500, 1000, 5000, 10000]:
Phi = rff_features(X, D=D, sigma=sigma, seed=0)
err = np.max(np.abs(K_exact - Phi @ Phi.T))
print(f' D={D:6d}: max error = {err:.5f}')
print('\n(d) Bochner connection:')
print(' Bochner: k(x-y) = E_{omega~p}[exp(i omega^T (x-y))]')
print(' RBF spectral measure: omega ~ N(0, sigma^{-2}I)')
print(' Sampling omega gives unbiased estimator: E[phi(x)^T phi(y)] = k(x,y)')
print(' -> D features = D-sample Monte Carlo average, error = O(1/sqrt(D))')
print('\nTakeaway: Fourier transform of kernel = spectral measure; RFF = Monte Carlo sampling thereof.')
Summary
| Exercise | Topic | Key Result |
|---|---|---|
| 1 ★ | Triangle wave coefficients | for odd |
| 2 ★ | Gram matrix / orthonormality | |
| 3 ★ | Parseval → Basel | |
| 4 ★★ | Spectral decay | $C^k \Rightarrow |
| 5 ★★ | Riemann–Lebesgue | for all functions |
| 6 ★★ | Shift & differentiation | |
| 7 ★★★ | RoPE implementation | Relative position via rotation |
| 8 ★★★ | Random Fourier Features | Bochner's theorem, kernel approximation |
Next: §20-02 Fourier Transform — Exercises
Exercise 9: Fejer Averaging and Gibbs Reduction
Compare an ordinary partial Fourier sum with a Fejer-averaged sum near a jump discontinuity. Explain why averaging reduces oscillation without changing the underlying coefficients.
Code cell 30
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 31
# Solution
header("Exercise 9: Fejer Averaging")
def square_coeff(n):
if n == 0 or n % 2 == 0:
return 0.0
return 2 / (1j * np.pi * n)
def partial_square(x, N):
total = np.zeros_like(x, dtype=complex)
for n in range(-N, N + 1):
total += square_coeff(n) * np.exp(1j * n * x)
return total.real
def fejer_square(x, N):
total = np.zeros_like(x, dtype=complex)
for n in range(-N, N + 1):
weight = 1 - abs(n) / (N + 1)
total += weight * square_coeff(n) * np.exp(1j * n * x)
return total.real
xs = np.linspace(-0.4, 0.4, 2000)
ordinary = partial_square(xs, 40)
fejer = fejer_square(xs, 40)
print("ordinary max near jump:", round(float(np.max(ordinary)), 6))
print("Fejer max near jump:", round(float(np.max(fejer)), 6))
check_true("Fejer has smaller overshoot", np.max(fejer) < np.max(ordinary))
print("Takeaway: Cesaro averaging damps the Dirichlet kernel side lobes.")
Exercise 10: Spectral Decay From Smoothness
Estimate Fourier coefficient decay for a triangle wave and a square wave, then connect decay rate to smoothness.
Code cell 33
# Your Solution
print("Write your solution here, then compare with the reference solution below.")
Code cell 34
# Solution
header("Exercise 10: Spectral Decay")
ns = np.arange(1, 60, 2)
square_mag = 2 / (np.pi * ns)
triangle_mag = 2 / (np.pi * ns**2)
slope_square = np.polyfit(np.log(ns), np.log(square_mag), 1)[0]
slope_triangle = np.polyfit(np.log(ns), np.log(triangle_mag), 1)[0]
print("square-wave slope:", round(float(slope_square), 3))
print("triangle-wave slope:", round(float(slope_triangle), 3))
check_true("triangle coefficients decay faster", slope_triangle < slope_square)
print("Takeaway: more smoothness gives faster spectral decay and easier approximation.")