  We present an exact analytical solution to the problem of shear dispersion
given a general initial condition. The solution is expressed as an infinite
series expansion involving Mathieu functions and their eigenvalues. The
eigenvalue system depends on the imaginary parameter $q=2ik$Pe, with $k$ the
wavenumber that determines the tracer scale in the initial condition and Pe the
P\'{e}clet number. Solutions are valid for all Pe, $t>0$, and $k>0$ except at
specific values of $q=q_{\ell}^{EP}$ called Exceptional Points (EPs), the first
occurring at $q_{0}^{EP}\approx1.468i$. For values of $q \lessapprox 1.468i$,
all the eigenvalues are real, different and eigenfunctions decay with time,
thus shear dispersion can be represented as a diffusive process. For values of
$q \gtrapprox 1.468i$, pairs of eigenvalues coalesce at EPs becoming complex
conjugates, the eigenfunctions propagate and decay with time, and so shear
dispersion is no longer a purely diffusive process.
  The limit $q\rightarrow0$ is approached by the small P\'{e}clet number limit
for all finite $k>0$, or equally by the large P\'{e}clet number limit as long
as $2k \ll 1/$Pe. The latter implies $k\rightarrow0$, strong separation of
scales between the tracer and flow. The limit $q\rightarrow\infty$ results from
large P\'{e}clet number for any $k>0$, or from large $k$ and non-vanishing Pe.
  We derive an exact closure that is continuous in wavenumber space. At small
$q$, the closure approaches a diffusion operator with an effective diffusivity
proportional to $U_0^2/\kappa$, for flow speed $U_0$ and diffusivity $\kappa$.
At large $q$, the closure approaches the sum of an advection operator plus a
half-derivative operator (differential operator of fractional order), the
latter with coefficient proportional to $\sqrt{\kappa U_0}$.
