Solution (source code)

= Solution

Put $s=\delta Z>0$. For the <Fay solution>, $1/\sinh(ns)=2e^{-ns}[1+O(e^{-2ns})]$. At $s\gg1$ the first harmonic therefore gives
$$
\boxed{q=4\delta e^{-s}\sin\theta+O(\delta e^{-2s}),}
$$
which agrees with the late-time limit in (i).

For the opposite limit, expand the reciprocal <hyperbolic sine> as a <geometric series> and interchange absolutely convergent sums at every $s>0$:
$$
\begin{aligned}
q&=4\delta\sum_{j=0}^{\infty}\sum_{n=1}^{\infty}e^{-(2j+1)ns}\sin n\theta\\
&=4\delta\sum_{j=0}^{\infty}\frac{e^{-(2j+1)s}\sin\theta}{1-2e^{-(2j+1)s}\cos\theta+e^{-2(2j+1)s}}\\
&=2\delta\sin\theta\sum_{j=0}^{\infty}\frac1{\cosh((2j+1)s)-\cos\theta}.
\end{aligned}
$$
This exact positive-denominator form is useful near the shock. For small $s$ and $|\theta|$, expand the denominator of the terms with small $(2j+1)s$ and use $\sin\theta\sim\theta$. The leading sum is
$$
q\sim\frac{4\delta\theta}{s^2}\sum_{j=0}^{\infty}\frac1{(2j+1)^2+(\theta/s)^2}.
$$
For the shock-layer scaling $\theta=s\eta$, the passage to this sum follows from dominated convergence: $\cosh z-1\ge z^2/2$ and $1-\cos\theta\ge c\theta^2$ for small $|\theta|$ provide a summable bound proportional to $[(2j+1)^2+\eta^2]^{-1}$. The same leading expression matches the outer small-angle range $s\ll|\theta|\ll1$. More generally, split the exact sum at $(2j+1)s=b$ with $\max(s,|\theta|)\ll b\ll1$: the low terms have vanishing relative Taylor error, while the discarded exact and approximate tails are smaller than the leading expression by $O(\max(s,|\theta|)/b)$.

Now apply the partial-fraction identity for the <hyperbolic tangent> with $\eta=\theta/s$:
$$
\sum_{j=0}^{\infty}\frac1{(2j+1)^2+\eta^2}=\frac{\pi}{4\eta}\tanh(\pi\eta/2).
$$
The <Fay shock-layer asymptotics> are therefore
$$
\boxed{q(\theta,Z)\sim\frac\pi Z\tanh\!\left(\frac{\pi\theta}{2\delta Z}\right),\qquad |\theta|\ll1,\ \delta Z\ll1.}
$$
The angular shock thickness is $O(\delta Z)$, and its matched states are $\pm\pi/Z$. At $\theta=0$, both formulas vanish; the limiting slopes agree. The denominator throughout this derivation is $\sinh(n\delta Z)$ as printed in the PDF, not the erroneous $\sin(n\delta Z)$ in the converted TeX.