Exercises

$\tilde{f} = \mathcal{F}\{f\}$, amplitude $= |\tilde{f}|$, phase $= \angle\tilde{f}$. Similarly for $g$.

$h_1 = \mathcal{F}^{-1}\{|\tilde{f}| \cdot e^{j\angle\tilde{g}}\}$ uses $f$'s amplitude and $g$'s phase. The result looks like $g$ (edges and structure from $g$'s phase) with $f$'s contrast profile.

SSIM$(h_1, g) >$ SSIM$(h_1, f)$ confirms phase dominance. Typical values: SSIM to phase source $\approx 0.6$--$0.8$, SSIM to amplitude source $\approx 0.1$--$0.3$.

$\mathbf{Y} = \frac{1}{M}\sum_i y_i \mathbf{a}_i\mathbf{a}_i^H$. Its expectation is $\mathbb{E}[\mathbf{Y}] = \mathbf{x}_0\mathbf{x}_0^H + \|\mathbf{x}_0\|^2\mathbf{I}$.

$\mathbf{z}^{(0)} = \sqrt{\lambda_1(\mathbf{Y})}\mathbf{v}_1$. After phase alignment: $\text{err} = \|\mathbf{z}^{(0)} - e^{j\hat{\phi}}\mathbf{x}_0\|/\|\mathbf{x}_0\|$.

At $M/N = 6$: error $\approx 0.25$--$0.35$. At $M/N = 12$: error $\approx 0.15$--$0.20$. The error decreases but does not vanish — Wirtinger flow iterations are needed to drive it below $10^{-2}$.

$\mathbf{h} = \mathbf{A}^T\mathbf{z}$ gives all inner products. Residuals: $\mathbf{r} = |\mathbf{h}|^2 - \mathbf{y}$. Gradient: $\nabla = \frac{1}{2M}\mathbf{A}^H(\mathbf{r} \odot \mathbf{h})$.

On a semilog plot, the error decreases linearly (straight line), confirming the linear convergence rate. Typical convergence: error $< 10^{-6}$ in $\sim$100 iterations for $\mu = 1.0$.

$\mu < 0.5$: slow convergence. $\mu \in [0.8, 2.0]$: optimal range. $\mu > 3.0$: oscillation or divergence.

$\|\hat{\mathbf{x}} - e^{j\phi}\mathbf{x}_0\|^2 = \|\hat{\mathbf{x}}\|^2 + \|\mathbf{x}_0\|^2 - 2\text{Re}(e^{-j\phi}\hat{\mathbf{x}}^H\mathbf{x}_0)$.

Maximizing $\text{Re}(e^{-j\phi}\hat{\mathbf{x}}^H\mathbf{x}_0)$ over $\phi$: set $e^{-j\phi} = \overline{\hat{\mathbf{x}}^H\mathbf{x}_0} /|\hat{\mathbf{x}}^H\mathbf{x}_0|$, giving $\hat{\phi} = \angle(\hat{\mathbf{x}}^H\mathbf{x}_0)$.

Since $e^{j\phi}\mathbf{x}_0$ is equally valid for any $\phi$, comparing $\hat{\mathbf{x}}$ to $\mathbf{x}_0$ directly would include a spurious phase error. Phase alignment removes this trivial ambiguity.

X = cp.Variable((N, N), hermitian=True)
constraints = [X >> 0]
constraints += [cp.trace(A_i @ X) == y_i for i in range(M)]
prob = cp.Problem(cp.Minimize(cp.trace(X)), constraints)
prob.solve(solver=cp.MOSEK)

Eigenvalues of $\hat{\mathbf{X}}$: $\lambda_1 \gg \lambda_2 \approx \ldots \approx \lambda_N \approx 0$. Ratio $\lambda_1/\lambda_2 > 10^8$ confirms rank-1 recovery.

Log-log slope $\approx 6$, confirming $O(N^6)$ scaling.

At each iteration: $\mathcal{T} = \{i : |\langle\mathbf{a}_i, \mathbf{z}\rangle| \leq 3\sqrt{y_i}\}$. Use only indices in $\mathcal{T}$ for the gradient sum.

Both methods converge, but truncated WF is $\sim$2$\times$ faster (fewer iterations to the same accuracy) due to removal of outlier gradient contributions.

Standard WF fails $\sim$60% of the time at $M = 2N$. Truncated WF succeeds $\sim$85% of the time. The phase transition (50% success) occurs at $M/N \approx 2.5$ for truncated WF vs. $M/N \approx 3.5$ for standard WF.

With $L = 1$ Fourier magnitude, alternating projections typically stagnate at relative error $\sim$0.3--0.5. The Fourier measurement structure has insufficient diversity for reliable recovery.

The error plateaus at $L = 3$--$4$. Theory predicts $L \geq 3$ coded patterns suffice for generic recovery. Additional masks reduce noise sensitivity but have diminishing returns for noiseless recovery.

Let $h_i = \langle\mathbf{a}_i, \mathbf{z}\rangle$. $\frac{\partial}{\partial \bar{\mathbf{z}}} |h_i| = \frac{h_i \mathbf{a}_i}{2|h_i|}$.

$\nabla_{\bar{\mathbf{z}}} g = \frac{1}{2M}\sum_i 2(|h_i| - \sqrt{y_i}) \cdot \frac{h_i\mathbf{a}_i}{2|h_i|} = \frac{1}{2M}\sum_i (1 - \sqrt{y_i}/|h_i|) h_i\mathbf{a}_i$.

When $|h_i| = 0$, the gradient term is $0/0$. Truncation removes indices where $|h_i|$ is very small, avoiding numerical instability.

$\langle\mathbf{a}, \mathbf{x}\rangle = (1 + j)/\sqrt{2}$, so $|\langle\mathbf{a}, \mathbf{x}\rangle|^2 = |1+j|^2/2 = 1$. $\text{tr}(\mathbf{a}\mathbf{a}^H \mathbf{x}\mathbf{x}^H) = \text{tr}\begin{pmatrix} 1 & -j \\ j & 1 \end{pmatrix}/2 = 1$.

$\mathbf{X} = \begin{pmatrix} 1 & -j \\ j & 1 \end{pmatrix}$. Eigenvalues: $\lambda_1 = 2$, $\lambda_2 = 0$. Rank 1, PSD (both eigenvalues $\geq 0$).

$\mathbf{X}_1 = \mathbf{e}_1\mathbf{e}_1^H$, $\mathbf{X}_2 = \mathbf{e}_2\mathbf{e}_2^H$. Average: $(\mathbf{X}_1 + \mathbf{X}_2)/2 = \mathbf{I}/2$, which has rank 2. Hence the rank-1 PSD set is non-convex.

All three methods show error $\propto 1/\text{SNR}$ at moderate-to-high SNR. Truncated WF and amplitude flow are slightly more robust at low SNR ($< 15$ dB).

To achieve the same relative error as linear recovery at SNR $= 20$ dB, phase retrieval needs SNR $\approx 24$--$26$ dB — a $4$--$6$ dB tax.

At low SNR ($\leq 10$ dB): amplitude flow is most robust (the amplitude loss is less sensitive to large noise). At high SNR ($\geq 30$ dB): all three methods converge to similar floor errors determined by measurement noise.

After the gradient step $\mathbf{z}' = \mathbf{z}^{(t)} - \mu\nabla f$, apply $\mathbf{z}^{(t+1)} = \mathcal{H}_s(\mathbf{z}')$.

Sparse WF converges $\sim$2$\times$ faster and achieves $\sim$3 dB lower NMSE than standard WF for $s = 10$. The improvement comes from the sparsity constraint preventing noise amplification in the zero components.

Sparse WF fails (NMSE $> -5$ dB) at $s \approx 30$--$40$ for the given measurement budget. This is consistent with the $M = O(s^2\log N)$ requirement: at $M = 6144$ and $N = 1024$, the maximum recoverable sparsity is $s \lesssim \sqrt{M/\log N} \approx 25$.

Alternate between Fourier domain (replace amplitude with $\sqrt{y_k}$, keep current phase) and spatial domain (apply support constraint).

Random init: converges in $\sim$30% of runs, median error $\sim$0.15. Spectral init: converges in $\sim$90% of runs, error $\sim$0.02 after 500 iterations.

Spectral init + 50 GS + 100 WF achieves error $\sim$0.005, matching pure WF (200 iter) at similar total cost. GS provides cheap early iterations; WF refines to high accuracy.

$L(\mathbf{X}, \boldsymbol{\lambda}, \mathbf{S}) = \text{tr}(\mathbf{X}) - \sum_i \lambda_i (\text{tr}(\mathbf{A}_i\mathbf{X}) - y_i) - \text{tr}(\mathbf{S}\mathbf{X})$. Dual: $\max_{\boldsymbol{\lambda}} \sum_i \lambda_i y_i$ s.t. $\mathbf{I} - \sum_i \lambda_i \mathbf{A}_i \succeq 0$.

For $\mathbf{X} = \alpha\mathbf{I}$ with $\alpha$ large enough, $\text{tr}(\mathbf{A}_i \mathbf{X}) = \alpha\|\mathbf{a}_i\|^2 > y_i$ is not necessarily feasible. But the constraint set is non-empty for the true $\mathbf{X}_0$, and we can perturb to strict feasibility.

If $\sum_i \lambda_i \mathbf{A}_i \preceq \mathbf{I}$ and $(\mathbf{I} - \sum_i \lambda_i \mathbf{A}_i)\mathbf{X}_0 = 0$, then $\boldsymbol{\lambda}$ is dual feasible and complementary slackness holds. By strong duality, $\mathbf{X}_0$ is primal optimal.

The success rate jumps from $\sim$20% at $L = 2$ to $\sim$95% at $L = 3$. This matches the theoretical prediction that $L \geq 3$ coded patterns suffice.

Random uniform phase: 95% success at $L = 3$. Random $\pm 1$: 80% success (less diverse). Limited range: 40% success (insufficient diversity).

The masks must make the effective measurement matrix $[\mathbf{D}_1\mathbf{F}; \ldots; \mathbf{D}_L\mathbf{F}]$ satisfy an RIP-like condition. Full $[0, 2\pi)$ phase randomness maximizes the incoherence.

$|F(\mathbf{k})|^2 = \mathcal{F}\{R(\boldsymbol{\tau})\}$ by the Wiener--Khinchin theorem. So from magnitude-only data, we can reconstruct $R$ by inverse Fourier transform.

$R(\boldsymbol{\tau})$ has bandwidth $2B$ (convolution doubles the support in frequency). But the autocorrelation of two point targets at distance $\Delta r$ has peaks at $0$, $+\Delta r$, and $-\Delta r$ — the peaks at $\pm\Delta r$ mix with sidelobes.

Successful phase retrieval recovers $\chi(\mathbf{r})$ at its native resolution $\delta_r = c/(2B)$. Two targets at $1.5\delta_r$ are cleanly resolved. The autocorrelation image shows merged peaks — the targets are not individually resolvable.

On a semilog plot, $\log e_t \approx \log e_0 + t\log\rho$. Linear regression gives $\rho$.

The estimated $c$ increases linearly with $M/N$, consistent with the theoretical bound $\rho = 1 - c_1(M/N)/N$ where $c_1$ is a universal constant.

$h(\mathbf{z}) = \frac{1}{4M_2}\sum_{j=1}^{M_2} (|\langle\mathbf{a}_j, \mathbf{z}\rangle|^2 - y_j)^2 + \frac{\alpha}{2M_1}\sum_{i=1}^{M_1} |\langle\mathbf{a}_i, \mathbf{z}\rangle - z_i|^2$.

The gradient adds a linear term from the complex measurements: $\nabla_{\bar{\mathbf{z}}} h = \nabla_{\bar{\mathbf{z}}} f + \frac{\alpha}{M_1}\sum_i (\langle\mathbf{a}_i, \mathbf{z}\rangle - z_i)\mathbf{a}_i$.

At $\sim$20% complex measurements, the hybrid method matches pure linear recovery quality. Even 5% complex measurements reduce the recovery error by $\sim$10 dB compared to pure magnitude-only.

For $y_i = |\langle\mathbf{a}_i, \mathbf{x}\rangle|^2 + \eta_i$ with $\eta_i \sim \mathcal{N}(0, \sigma^2)$: $\frac{\partial y_i}{\partial \bar{\mathbf{x}}} = \langle\mathbf{a}_i, \mathbf{x}\rangle\mathbf{a}_i$. $\mathbf{J} = \frac{1}{\sigma^2}\sum_i |\langle\mathbf{a}_i, \mathbf{x}\rangle|^2 \mathbf{a}_i\mathbf{a}_i^H$.

$\mathbb{E}[\mathbf{J}] \propto \frac{M}{\sigma^2} (\|\mathbf{x}\|^2\mathbf{I} + \mathbf{x}\mathbf{x}^H)$. The CRB is $\text{MSE} \geq \text{tr}(\mathbf{J}^{-1}) \propto \frac{N\sigma^2}{M\|\mathbf{x}\|^2}$.

Empirically, WF achieves MSE within a factor of 2--3 of the CRB for $M/N \geq 6$, confirming near-optimal statistical efficiency.

$\min_{\mathbf{V}} \sum_i (\text{tr}(\mathbf{A}_i \mathbf{V}\mathbf{V}^H) - y_i)^2$. Gradient: $\nabla_\mathbf{V} = 4\sum_i (\text{tr}(\mathbf{A}_i\mathbf{V}\mathbf{V}^H) - y_i) \mathbf{A}_i\mathbf{V}$.

$r = 1$: equivalent to amplitude flow; works but slow convergence near saddle points. $r = 2$: no spurious local minima; converges reliably. $r \geq 2$ is recommended.

For $N = 64$: BM with $r = 2$ is $\sim$30% slower than WF (gradient involves $N \times 2$ matrix operations). For $N > 256$: BM can be faster than SDP but slower than direct WF. BM is useful when you want the lifted solution $\mathbf{X}$ (e.g., for uncertainty quantification).

Build the $M \times N$ sensing matrix $\mathbf{A}$ with entries $A_{in} = e^{-j\kappa(\|\mathbf{s}_i - \mathbf{p}_n\| + \|\mathbf{r}_i - \mathbf{p}_n\|)}$ for each Tx-Rx-frequency triple. Apply $L = 4$ random phase masks.

The full pipeline ($\sim$11 s) achieves within 4.4 dB of the coherent baseline. The sparse WF stage provides the best marginal gain per second. ADMM-TV adds 1.8 dB at moderate cost. The "phase retrieval tax" is $\sim$4.4 dB.