Part III › Chapter 9

MIMO & Water-Filling

Multiple-input multiple-output (MIMO) channels are decomposed via SVD into independent sub-channels, then power is allocated optimally using the water-filling algorithm.

9.1 The MIMO Channel Model

A MIMO system with \(n_T\) transmit and \(n_R\) receive antennas has the input-output relation:

\[ \mathbf{Y} = H\mathbf{X} + \mathbf{Z}, \qquad H \in \mathbb{C}^{n_R \times n_T} \]

where \(\mathbf{Z} \sim \mathcal{CN}(0, \sigma^2 I)\) is complex Gaussian noise and\(H\) is the channel matrix whose entries represent the complex gain from each transmit to each receive antenna.

The capacity with known \(H\) at both ends and power constraint\(\operatorname{tr}(Q) \le P\) where \(Q = E[\mathbf{X}\mathbf{X}^\dagger]\):

\[ C = \max_{Q \succeq 0,\;\mathrm{tr}(Q)\le P} \log_2\det\!\left(I + \frac{1}{\sigma^2}HQH^\dagger\right) \]

9.2 SVD Decomposition into Parallel Channels

The singular value decomposition of \(H\):

\[ H = U\Lambda V^\dagger, \qquad \Lambda = \operatorname{diag}(\lambda_1, \ldots, \lambda_r) \]

By transmitting along right singular vectors (columns of \(V\)) and receiving along left singular vectors (columns of \(U\)), the MIMO channel decomposes into\(r = \operatorname{rank}(H)\) independent scalar sub-channels:

\[ \tilde{Y}_i = \lambda_i \tilde{X}_i + \tilde{Z}_i, \quad i = 1,\ldots,r \]

Each sub-channel has gain \(\lambda_i\) (a singular value) and independent noise\(\tilde{Z}_i \sim \mathcal{CN}(0, \sigma^2)\). The capacity becomes:

\[ C = \sum_{i=1}^{r} \log_2\!\left(1 + \frac{P_i\lambda_i^2}{\sigma^2}\right) \]

subject to \(\sum_i P_i \le P\), \(P_i \ge 0\)

9.3 Water-Filling Algorithm

We must maximize \(\sum_i \log_2(1+P_i\gamma_i)\) where \(\gamma_i = \lambda_i^2/\sigma^2\)subject to \(\sum P_i \le P\). Using Lagrange multipliers with KKT conditions:

\[ P_i^* = \left(\mu - \frac{1}{\gamma_i}\right)^+, \qquad \text{where } (x)^+ = \max(0, x) \]

The water level \(\mu\) is chosen so that\(\sum_i P_i^* = P\). The “water-filling” metaphor: imagine pouring water into vessels of different heights (noise floors \(1/\gamma_i\)) — water fills the lower vessels first, and the surface is the water level \(\mu\).

Python: Water-Filling Simulation

Implements the water-filling algorithm for 4 parallel channels with gain-to-noise ratios 4, 1, 0.5, 0.25. Four plots: the water-filling bar diagram, capacity vs total power budget (water-filling vs equal allocation), per-channel power vs budget, and a MIMO 4×4 capacity distribution compared to SISO.

Python

script.py193 lines

import numpy as np
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches

plt.style.use('dark_background')

def water_filling(channel_gains, total_power, n_iter=1000):
    """
    Water-filling power allocation.
    channel_gains: array of lambda_i^2 / sigma_i^2  (gain-to-noise ratio per channel)
    Returns allocated powers P_i >= 0 summing to total_power.
    """
    n = len(channel_gains)
    noise_floors = 1.0 / np.array(channel_gains)  # sigma_i^2 / lambda_i^2

# Sort by noise floor
    idx = np.argsort(noise_floors)
    nf_sorted = noise_floors[idx]

# Binary search for water level mu
    lo, hi = np.max(nf_sorted), np.max(nf_sorted) + total_power
    for _ in range(200):
        mu = (lo + hi) / 2
        P = np.maximum(0, mu - nf_sorted)
        if P.sum() < total_power:
            hi = mu * 2
            lo, hi = lo, mu * 2
            # simple linear search instead
            break

# Linear search for mu
    mu_vals = np.linspace(np.min(nf_sorted), np.min(nf_sorted) + total_power * 2, 10000)
    best_mu = mu_vals[0]
    for mu in mu_vals:
        P_try = np.maximum(0, mu - nf_sorted)
        if P_try.sum() >= total_power:
            best_mu = mu
            break

P_sorted = np.maximum(0, best_mu - nf_sorted)
    # Scale to exactly match total power (numerical precision)
    if P_sorted.sum() > 0:
        P_sorted *= total_power / P_sorted.sum()

# Unsort
    P = np.zeros(n)
    P[idx] = P_sorted
    return P, best_mu

fig, axes = plt.subplots(2, 2, figsize=(13, 10))
fig.patch.set_facecolor('#0a0a0f')
for ax in axes.flat:
    ax.set_facecolor('#0d0d1a')

# ── Channel parameters (4 parallel sub-channels) ─────────────────────────
n_channels = 4
# Gain-to-noise ratios for each sub-channel (lambda_i^2 / sigma_i^2)
gnr = np.array([4.0, 1.0, 0.5, 0.25])  # descending: ch1 best, ch4 worst
noise_floors = 1.0 / gnr
total_power = 2.0

P_opt, mu_opt = water_filling(gnr, total_power)

print(f"Water-filling results (P_total = {total_power}):")
for i in range(n_channels):
    cap_i = 0.5 * np.log2(1 + P_opt[i] * gnr[i]) if P_opt[i] > 0 else 0
    print(f"  Ch {i+1}: noise_floor={noise_floors[i]:.3f}  P*={P_opt[i]:.4f}  C_i={cap_i:.4f} bits")
C_total = sum(0.5 * np.log2(1 + P_opt[i] * gnr[i]) for i in range(n_channels) if P_opt[i] > 0)
print(f"  Total capacity: {C_total:.4f} bits")

# Equal power allocation for comparison
P_equal = np.ones(n_channels) * total_power / n_channels
C_equal = sum(0.5 * np.log2(1 + P_equal[i] * gnr[i]) for i in range(n_channels))
print(f"  Equal allocation capacity: {C_equal:.4f} bits")
print(f"  Water-filling gain: {(C_total - C_equal)/C_equal*100:.1f}%")

# ── 1. Water-filling diagram ───────────────────────────────────────────────
ax = axes[0, 0]
x_positions = np.arange(n_channels)
bar_width = 0.5
colors_ch = ['#818cf8', '#34d399', '#fbbf24', '#f87171']

# Bottom: noise floors
bars_noise = ax.bar(x_positions, noise_floors, bar_width,
                    color=[c + '55' for c in ['#818cf8','#34d399','#fbbf24','#f87171']],
                    edgecolor='#374151', label='Noise floor 1/γᵢ')

# Top: allocated power
for i in range(n_channels):
    if P_opt[i] > 1e-6:
        ax.bar(x_positions[i], P_opt[i], bar_width,
               bottom=noise_floors[i], color=colors_ch[i],
               alpha=0.85, edgecolor='white', lw=0.7,
               label=f'Power P*_{i+1}={P_opt[i]:.3f}' if P_opt[i] > 1e-4 else f'P*_{i+1}≈0')

# Water level line
ax.axhline(mu_opt, color='#60a5fa', lw=2, ls='--', label=f'Water level μ = {mu_opt:.3f}')
ax.set_xticks(x_positions)
ax.set_xticklabels([f'Ch {i+1}\nγ={gnr[i]}' for i in range(n_channels)], color='#9ca3af')
ax.set_ylabel('Power level', color='#9ca3af')
ax.set_title('Water-Filling Power Allocation', color='#a5b4fc', fontweight='bold')
ax.legend(fontsize=7, facecolor='#0d0d1a', edgecolor='#374151', labelcolor='white', loc='upper right')
ax.tick_params(colors='#6b7280')
for spine in ax.spines.values(): spine.set_edgecolor('#1e293b')

# ── 2. Capacity vs total power ─────────────────────────────────────────────
ax = axes[0, 1]
P_totals = np.linspace(0.01, 10, 200)
C_wf_list = []
C_eq_list = []

for Pt in P_totals:
    P_wf, _ = water_filling(gnr, Pt)
    C_wf_list.append(sum(0.5*np.log2(1+P_wf[i]*gnr[i]) for i in range(n_channels) if P_wf[i]>1e-9))
    Peq = Pt / n_channels
    C_eq_list.append(sum(0.5*np.log2(1+Peq*gnr[i]) for i in range(n_channels)))

ax.plot(P_totals, C_wf_list, color='#818cf8', lw=2.5, label='Water-filling (optimal)')
ax.plot(P_totals, C_eq_list, color='#f87171', lw=2, ls='--', label='Equal power allocation')
ax.fill_between(P_totals, C_eq_list, C_wf_list, alpha=0.15, color='#818cf8', label='WF gain')
ax.scatter([total_power], [C_total], color='#34d399', s=80, zorder=5,
           label=f'Current P={total_power}: C={C_total:.2f}')
ax.set_xlabel('Total Power P', color='#9ca3af')
ax.set_ylabel('Total Capacity C  [bits]', color='#9ca3af')
ax.set_title('Water-Filling vs Equal Allocation', color='#a5b4fc', fontweight='bold')
ax.legend(fontsize=8, facecolor='#0d0d1a', edgecolor='#374151', labelcolor='white')
ax.tick_params(colors='#6b7280')
for spine in ax.spines.values(): spine.set_edgecolor('#1e293b')

# ── 3. Per-channel allocation across budgets ───────────────────────────────
ax = axes[1, 0]
for i in range(n_channels):
    pi_list = []
    for Pt in P_totals:
        P_wf, _ = water_filling(gnr, Pt)
        pi_list.append(P_wf[i])
    ax.plot(P_totals, pi_list, color=colors_ch[i], lw=2,
            label=f'Ch {i+1}  (γ={gnr[i]})')

ax.set_xlabel('Total Power Budget P', color='#9ca3af')
ax.set_ylabel('Allocated Power Pᵢ*', color='#9ca3af')
ax.set_title('Per-Channel Power vs Budget', color='#a5b4fc', fontweight='bold')
ax.legend(fontsize=9, facecolor='#0d0d1a', edgecolor='#374151', labelcolor='white')
ax.tick_params(colors='#6b7280')
ax.text(6, 0.2, 'Weak channels (low γ) activated later\nas water level rises',
        color='#9ca3af', fontsize=8,
        bbox=dict(boxstyle='round,pad=0.3', facecolor='#0d0d1a', edgecolor='#374151'))
for spine in ax.spines.values(): spine.set_edgecolor('#1e293b')

# ── 4. SVD decomposition of a random MIMO channel ─────────────────────────
ax = axes[1, 1]
rng = np.random.default_rng(7)
n_t, n_r = 4, 4  # 4x4 MIMO
n_trials = 400
capacities_mimo = []
capacities_single = []

for _ in range(n_trials):
    H = (rng.standard_normal((n_r, n_t)) + 1j * rng.standard_normal((n_r, n_t))) / np.sqrt(2)
    _, sv, _ = np.linalg.svd(H)
    gains = sv**2  # eigenvalues of H H†
    # Water-fill over parallel channels
    Pt_mimo = 4.0
    P_wf_mimo, _ = water_filling(gains, Pt_mimo)
    sigma2 = 1.0
    C_mimo = sum(np.log2(1 + P_wf_mimo[i]*gains[i]/sigma2) for i in range(len(gains)) if P_wf_mimo[i]>1e-9)
    capacities_mimo.append(C_mimo)
    # Single antenna baseline
    snr_single = Pt_mimo / sigma2
    C_single = np.log2(1 + snr_single)
    capacities_single.append(C_single)

bins = np.linspace(0, max(max(capacities_mimo), max(capacities_single)) + 1, 40)
ax.hist(capacities_mimo, bins=bins, color='#818cf8', alpha=0.75, label=f'MIMO 4×4  (mean={np.mean(capacities_mimo):.1f})')
ax.hist(capacities_single, bins=bins, color='#f87171', alpha=0.75, label=f'SISO  (C={capacities_single[0]:.1f})')
ax.axvline(np.mean(capacities_mimo), color='#818cf8', lw=1.5, ls='--')
ax.axvline(capacities_single[0], color='#f87171', lw=1.5, ls='--')
ax.set_xlabel('Capacity  [bits/channel use]', color='#9ca3af')
ax.set_ylabel('Count', color='#9ca3af')
ax.set_title('MIMO 4×4 vs SISO Capacity Distribution', color='#a5b4fc', fontweight='bold')
ax.legend(fontsize=9, facecolor='#0d0d1a', edgecolor='#374151', labelcolor='white')
ax.tick_params(colors='#6b7280')
for spine in ax.spines.values(): spine.set_edgecolor('#1e293b')

fig.suptitle('Chapter 9: MIMO & Water-Filling', color='#a5b4fc',
             fontsize=14, fontweight='bold', y=1.01)
plt.tight_layout()
plt.savefig('output.png', dpi=130, bbox_inches='tight',
            facecolor='#0a0a0f', edgecolor='none')
print("\nSaved output.png")

Click Run to execute the Python code

Code will be executed with Python 3 on the server

← Gaussian Channels Next: Kolmogorov Complexity →