Past exam of the mathematics course of the University of Cambridge / 2019 / iii / Paper 335 (Direct and Inverse Scattering of Waves)

Use the outgoing Green function for the three-dimensional Helmholtz equationLet ${\mathbf{n}}^{\prime }$ point out of the obstacle $D$ into the exterior. Apply Green second identity to the field and $G_k$ in a large exterior region. Their Sommerfeld radiation conditions make the outer-sphere contribution vanish. The normal on the exterior region's inner boundary is $-{\mathbf{n}}^{\prime }$, so the Kirchhoff–Helmholtz representation isThe first integral is the incident source field ${\psi }_i$. Therefore the scattered field is the surface term. Since ${\psi }_i$ solves the homogeneous Helmholtz equation within the obstacle, its surface term is zero at exterior observation points; the same representation can be written using only scattered-field traces:The unknown boundary condition determines these traces, but the representation itself does not require choosing that condition.

For large $r=|\mathbf{r}|$, the outgoing kernel hasThus the surface representation of an obstacle far-field pattern, with ${\psi }_s=e^{ikr}{f}_{\infty }/r+O(r^{-2})$, givesConjugating the Herglotz wave function gives $\overline{v_g}({\mathbf{r}}^{\prime })={\int }_{S^2}{e}^{-ik\hat{\mathbf{r}}\cdot {\mathbf{r}}^{\prime }}\overline{g(\hat{\mathbf{r}})}dS$. Interchanging the angular and obstacle-surface integrals yields the Herglotz pairing with an obstacle far field:This contains only the scattered-field boundary traces and the Herglotz function, with the outward-obstacle normal convention fixing the sign.

The boundary identity gives a numerical route to shape recovery from a Herglotz boundary identity. First recover an admissible density $g$ from measured far-field patterns, enforcing the supplied Herglotz pairing with an obstacle far field for available incident directions:Then evaluate its Herglotz wave function throughout a search region and find a smooth positive radial profile satisfyingThe same $h(\theta )$ must work for every azimuth $\varphi$, expressing the assumed axial symmetry. In practice expand $g$ and $h$ in basis functions, enforce these equations at collocation points, and minimize the complex residual together with the far-field-data residual. The exclusion of an interior Dirichlet eigenvalue ensures uniqueness of the interior Dirichlet problem for the trial boundary data. Use regularization of an inverse problem to control measurement error and unstable density components; neither uniqueness nor numerical stability follows merely from writing the residual equations.

For rich incident-direction data, the conjugated pairing constraints form an adjoint far-field integral equation for $g$. For only one fixed incidence, the displayed pairing is a single complex scalar constraint and cannot by itself determine an arbitrary density. Additional data or a joint constrained shape fit is then needed; the trace condition supplies the actual geometric information.

Let $V=k_0^2(n^2-1)$ be the scattering potential supported in $D$. The background Helmholtz equation and its outgoing kernel give the scalar Lippmann-Schwinger equationReplacing the unknown field in the integral by the incident field gives the Born approximation for scalar wave scattering:For the Rytov approximation, work where ${\psi }_i\ne 0$ and set $\psi ={\psi }_i{e}^{\varphi }$. Dividing the wave equation by the total field and subtracting the incident equation givesNeglect the quadratic gradient term. If ${\varphi }_1={\psi }_{s,B}/{\psi }_i$, direct differentiation shows $\Delta {\varphi }_1+2\nabla \log {\psi }_i\cdot \nabla {\varphi }_1=-V$, soExpanding the exponential in the potential gives ${\psi }_R={\psi }_i+{\psi }_{s,B}+O(V^2)={\psi }_B+O(V^2)$. Thus both approximations have the same first-order field, although exponentiation retains a particular family of higher powers.

For unit plane wave incidence ${\psi }_i=e^{i{k}_0{\hat{\mathbf{r}}}_0\cdot \mathbf{r}}$, define the scattering wavevector difference $\mathbf{q}=k_0({\hat{\mathbf{r}}}_0-\hat{\mathbf{r}})$ andThe outgoing Green function for the three-dimensional Helmholtz equation givesThe Born far-field pattern is $A$, a Fourier transform sample of the potential. Expanding the far-field Rytov exponential in $1/r$ gives the same leading outgoing term, while the logarithmic representation preserves its accumulated phase and amplitude correction.

Use the leading far-field kernel from the previous part. DefineThen ${\varphi }_1={\int }_{D}KV$ and, because the incident intensity is one,To first order in weak contrast, $V=2\mu k_0^{2}W+O({\mu }^2)$. Its linearized part is a centered stationary Gaussian random field, with autocorrelation function of a random field $C_V(\mathbf{s})=4{\mu }^2{k}_0^4{C}_W(\mathbf{s})$. The exponential moment of a Gaussian linear functional therefore gives the usual Gaussian intensity in the first Rytov approximation:Equivalently the exponent is $\mathbb{E}|{\varphi }_1{|}^2+\mathrm{Re}\mathbb{E}{\varphi }_1^2$. The second term generally matters because a finite real random medium need not generate a circular complex Gaussian amplitude.

If the full stated potential $V=k_0^2(2\mu W+{\mu }^2{W}^2)$ is retained, it is not Gaussian and has mean $m_V={\mu }^2{k}_0^2$. To retain every term through order ${\mu }^2$ in the intensity of this first Rytov field, the cumulant expansion instead givesHere $C_V^{\mathrm{cov}}=\mathbb{E}[(V-m_V)(V^{\prime }-m_V)]$; the raw autocorrelation differs by $m_V^2=O({\mu }^4)$. The stated remainder uses Gaussian finite-dimensional statistics of $W$, for which odd total moments vanish. This mean term is absent if one first truncates the potential to $2\mu k_0^{2}W$. Both formulas concern the first Rytov approximation, and exponentiating it does not supply the omitted higher-order scattering corrections.

For admissible data, the Moore–Penrose inverse of an operator selects the minimum-norm least-squares solution:Its domain is $\mathrm{ran}A\oplus (\mathrm{ran}A{)}^{\perp }$. When the equation is consistent it is the solution of $Ax=y$ orthogonal to the kernel; for inconsistent admissible data it solves the normal equation for a linear inverse problem and projects away the component orthogonal to the range.

Take a singular system of a compact operator satisfying $A{v}_j={\sigma }_j{u}_j$, $A^{\ast }{u}_j={\sigma }_j{v}_j$. The generalized inverse isThe Picard criterion requires the squared coefficients to be summable. For an infinite-rank compact operator, ${\sigma }_j\to 0$, so $\Vert A^{†}{u}_j\Vert =1/{\sigma }_j\to \infty$ even though $\Vert u_j\Vert =1$. Hence $A^{†}$ is unbounded. Finite-rank compact operators do not have this obstruction. Data failing the Picard condition need not have a least-squares minimizer at all, even when the infimum residual is zero.

A regularization of an inverse problem replaces the unbounded generalized inverse by bounded operators $R_{\alpha }:Y\to X$ satisfying $R_{\alpha }y\to A^{†}y$ as $\alpha \downarrow 0$ on admissible exact data. A spectral regularization method can be writtenwhere $f_{\alpha }(\sigma )\to 1$ for fixed $\sigma >0$ while $f_{\alpha }(\sigma )/\sigma$ remains bounded for each fixed $\alpha$. A uniform bound on $f_{\alpha }$ supplies the usual dominated-convergence argument on exact data.

This defines $f_{\alpha }$ as the dimensionless damping factor. In the alternative convention calling the entire inverse coefficient the spectral filter, that filter is $f_{\alpha }(\sigma )/\sigma$. For data error bounded by $\delta$, a sufficient regularization parameter choice also requires $\delta \Vert R_{\alpha (\delta )}\Vert \to 0$; consistency alone does not control noise.

Minimize the data-misfit functional $J(x)=\frac{1}{2}\Vert Ax-y{\Vert }^2$. For a perturbation $h$, $DJ(x)[h]=\mathrm{Re}⟨A^{\ast }(Ax-y),h⟩$, so its gradient is $A^{\ast }(Ax-y)$. Gradient descent gives Landweber iteration:On a singular component, the iteration error is multiplied by $1-\tau {\sigma }_j^2$, whose modulus is less than one. Starting from zero preserves orthogonality to $\ker A$ and makes the exact-data iterates converge to $x^{†}$.

For the unit-step formulas and noise bound used below, rescale so that $\tau =1$ and $\Vert A\Vert \le 1$. Finite iteration count bounds amplification of small singular-value components; the stopping count provides the regularization parameter.

Use early stopping of Landweber iteration. For each finite $n$, the Landweber spectral filter gives a bounded reconstruction operator; as $n\to \infty$, it approaches the generalized inverse on exact admissible data. Taking $\alpha =1/n$, a sufficient rule for a convergent regularization of an inverse problem isor equivalently $\alpha (\delta )\to 0$ and $\delta /\sqrt{\alpha (\delta )}\to 0$. For example, $n(\delta )=\lfloor {\delta }^{-1}\rfloor$ satisfies both in the normalized problem. The first condition removes exact-data bias; the second prevents arbitrarily small measurement errors from being amplified by excessive iteration.

Add and subtract the exact-data iterate. The noise-bias decomposition for linear regularization givesso the supplied Landweber noise amplification bound yieldsWith unit step and zero initial iterate, the exact-data error is $x^{†}-x_n=(I-A^{\ast }A{)}^n{x}^{†}$. Its squared norm isby the dominated convergence theorem. Stopping with $n\to \infty$ but $n{\delta }^2\to 0$ therefore makes both errors vanish.

On each left singular vector $u_j$, use $A^{\ast }{u}_j={\sigma }_j{v}_j$ and $(I-A^{\ast }A)v_j=(1-{\sigma }_j^2)v_j$. The finite geometric series givesThus with $\alpha =1/n$, the dimensionless Landweber spectral filter isThe full inverse-coefficient convention instead uses $[1-(1-{\sigma }^2{)}^{1/\alpha }]/\sigma$. With general step $\tau$, replace ${\sigma }^2$ inside the power by $\tau {\sigma }^2$.

In the normalized problem, $0\le f_n(\sigma )\le \min (n{\sigma }^2,1)$, so the inverse coefficient is at most $\min (n\sigma ,1/\sigma )\le \sqrt{n}$. This verifies directly why finite iteration is stable and why admitting progressively smaller singular values eventually amplifies noise.