# Exploding paraxial beams, vortex beams, and cylindrical beams of light with
finite power in linear media, and their enhanced longitudinal field

Miguel A. Porras  Grupo de Sistemas Complejos, ETSIME, Universidad Politécnica
de Madrid, Rios Rosas 21, 28003 Madrid, Spain

###### Abstract

We present a set of paraxial light beams with cylindrical symmetry, smooth and
localized transversal profile carrying finite power, that develop intensity
singularities when they are focused in a linear medium, such as vacuum. They
include beams with orbital angular momentum and with radial polarization, in
which case they develop punctual phase and polarization singularities
surrounded by infinitely bright rings, along with singular longitudinal
fields. In practice, these effects are manifested in focal intensities and
spot sizes, vortex bright ring intensities and radii, and strengths of the
longitudinal field, that strongly change with the lens aperture radius.
Continuous control of these focal properties is thus exercised without
changing the light incident on the lens, with substantially the same collected
power, and while maintaining paraxial focusing conditions. As solutions of the
Schrödinger equation, these exploding beams have analogues in other areas of
physics where this equation is the fundamental dynamical model.

##  I Introduction

Inspired by what happens to some wave functions in quantum mechanics PERES ,
Aiello has recently introduced a class of paraxial light beam with a localized
transversal profile and finite power that develops a singularity when its is
ideally focused AIELLO . Real, apertured versions of these beams have
subsequently been demonstrated in experiments AIELLO2 . These beams, even if
apertured, outperform standard Gaussian beams of similar intensity and power
in terms of focal intensity and resolution. In a sense, these beams are able
to reproduce with a finite amount of power the focusing behavior of plane
waves, though their finite-aperture realizations do not exceed the diffraction
limit AIELLO . The electric field of these beams is factorized in the
coordinates $x$ and $y$ in the transversal plane, i. e., it has a rectangular
geometry, which greatly facilitates its analytical treatment.

In many experimental settings, however, cylindrical symmetry is advisable, if
not mandatory; in particular, lenses, whose aperture plays a crucial role in
the behavior of these beams, are most often circular. In this paper we
describe cylindrical symmetric beams with finite power and similar exploding
behavior in their intensity when they are focused. Considering beams with
cylindrical geometry, we can also examine other types of beams that are of
great interest today. We describe exploding beams with orbital angular
momentum, and exploding radially polarized beams, both with finite power. They
form a punctual dark vortex surrounded by an infinitely bright ring,
accompanied by an infinitely strong longitudinal component.

Their apertured versions approach the above ideal behavior as the lens
aperture is increased, which offers a practical way to control these
properties without changing the illuminating beam, including the collected
power, e. g., to enhance the amplitude of the axial component up to $66$
percent of that of the transverse component under paraxial conditions. These
exploding beams, vortex beams and radially polarized beams with large, but
paraxial apertures, greatly outclass the performance of standard illuminating
fields such as Gaussian or Laguerre-Gaussian beams of similar power and
intensity. As their Cartesian counterparts, exploding beams of finite aperture
do not beat the diffraction limit for plane wave illumination, but approach it
closely. The practical advantages of exploding versus uniform illumination for
controlling the properties of the focused light is discussed at the end of
this paper.

Of course, the exploding behavior of these light beams has nothing to do with
blow-up light beams in self-focusing nonlinear media KIVSHAR ; FIBICH . Yet,
they are somewhat similar in that they both appear in paraxial fields, are
dissolved when more precise models are considered, and their mathematical
existence has physical repercussions. The results presented here can
straightforwardly extended to waves in other areas of physics such as
probability waves in quantum mechanics, or mechanical waves in acoustics, and
to matter waves and electron beams.

##  II Ideal exploding beams and vortex beams

In the Fresnel regime, the focused electric field $\psi(x,y,z)e^{i(kz-\omega
t)}$ of frequency $\omega$ and propagation constant $k$ (e. g., $k=\omega/c$
in vacuum), and the field $\psi(x,y)e^{-i\omega t}$ at the input plane of a
thin lens of focal length $f$, are related by

$\psi(x,y,z)=\frac{k}{2\pi
iz}\int\\!\\!dx^{\prime}dy^{\prime}\psi(x^{\prime},y^{\prime},0)e^{\frac{ik}{2z}[(x-x^{\prime})^{2}+(y-y^{\prime})^{2}]}\,,$
(1)

where $\psi(x,y,0)=\psi(x,y)e^{-ik(x^{2}+y^{2})/2f}$, and the integral extends
over the aperture of the lens. We are interested in an illuminating field of
the form $\psi(x,y)=\psi(r)e^{is\varphi}$, where $(r,\varphi)$ are polar
coordinates in the transversal plane, and $s=0,\pm 1,\pm 2,\dots$ i. e., in
light beams with revolution symmetry about the propagation axis without and
with orbital angular momentum associated with vortex of topological charge
$s$. For these fields, Fresnel integral in (1) can more conveniently be
written as

$\displaystyle\psi(r,z)e^{is\varphi}$ $\displaystyle=$
$\displaystyle\frac{ke^{is\varphi}}{i^{|s|+1}z}e^{\frac{ikr^{2}}{2z}}\int_{0}^{R}dr^{\prime}r^{\prime}\psi(r^{\prime})$
(2) $\displaystyle\times$ $\displaystyle e^{\frac{-ikr^{\prime
2}}{2}\left(\frac{1}{f}-\frac{1}{z}\right)}J_{|s|}\left(\frac{k}{z}rr^{\prime}\right)\,,$

where $J_{n}(\cdot)$ is the Bessel function of the first kind and order $n$
GRADS , and $R$ is the radius of the lens aperture.

We consider the illuminating, collimated field with transversal distribution
of amplitude

$\psi(r)=\sqrt{\frac{P}{A}}\frac{(r/\sigma)^{|s|}}{(1+r^{2}/\sigma^{2})^{\mu+1}}\,,$
(3)

where the length $\sigma$ scales the field transversally, and $P$ and $A$ are
constants to be conveniently fixed. With the real number $\mu$ satisfying
$\mu>(|s|-1)/2$, the power of the illuminating field in the entire transversal
plane, $2\pi\int_{0}^{\infty}drr|\psi(r)|^{2}$, is finite, and therefore this
field is in principle physically realizable. With the choice
$A=\pi\sigma^{2}\Gamma(|s|+1)\Gamma(2\mu+1-|s|)/\Gamma(2\mu+2)$, where
$\Gamma(\cdot)$ is the Gamma function, the beam power is just the constant $P$
appearing in (3). Examples of (3) with $s=0$ and $|s|=1$ are depicted as solid
curves in Figs. 1(a) and 2(a).

At the focal plane, $z=f$, the integral in (2) with (3) can be carried out
analytically (see 6.565.4 in Ref. GRADS ) when neglecting the finiteness of
the lens aperture, i. e., with $R=\infty$:

$\displaystyle\psi(r,f)e^{is\varphi}$ $\displaystyle=$
$\displaystyle\sqrt{\frac{P}{A}}\frac{k}{i^{|s|+1}f}\frac{\sigma^{2}}{2^{\mu}\Gamma(\mu+1)}$
(4) $\displaystyle\times$ $\displaystyle
e^{\frac{ikr^{2}}{2f}}\left(\frac{k}{f}\sigma
r\right)^{\mu}K_{|s|-\mu}\left(\frac{k}{f}\sigma r\right)e^{is\varphi}\,,$

for $r>0$, where $K_{\nu}(\cdot)$ is the modified Bessel function of the
second kind and order $\nu$ GRADS . Using that $K_{-\nu}(s)=K_{\nu}(s)$, and
the asymptotic behavior
$K_{\nu}(\alpha)\simeq(1/2)\Gamma(\nu)(\alpha/2)^{-\nu}$ for $\nu>0$ for small
argument GRADS , (4) is seen to approach infinity for $r\rightarrow 0$ when
$\mu<|s|/2$. For these values of $\mu$, the vortex-less beam ($s=0$) reaches
infinitely large amplitude at the beam center $r=0$, and the vortex beam
($s\neq 0$) vanishes at $r=0$, since $J_{|s|}(0)=0$ for $s\neq 0$, i.e., has
an punctual vortex surrounded by an infinitely bright ring. In short, the
unapertured focused field of illumination in (3) carries finite power and
produces infinitely intense field approaching the focal point if

$\frac{|s|-1}{2}<\mu<\frac{|s|}{2}\,,$ (5)

and this unbounded field surrounds an infinitely narrow, dark vortex when
$s\neq 0$. The singularity originates from the small but persistent tails of
(3) at large radius that are superposed with uniform phases at the focal
point, whose contribution to the beam power is nevertheless finite and small.

We point out that the family of light beams with amplitude
$\psi(r)e^{is\varphi}e^{-ir^{2}/2f}$ in (3) at $z=0$ and the singular field
$\psi(r,f)e^{is\varphi}$ in (4) at $z=f$ are solutions of the Schödinger
equation

$\frac{\partial\psi}{\partial z}=\frac{i}{2k}\Delta_{\perp}\psi\,,$ (6)

where $\Delta_{\perp}=\partial^{2}/\partial_{x}^{2}+\partial^{2}/\partial
y^{2}$, and as such these beams can be directly translated to other areas of
physics where the Schrödinger equation is the dynamical model, e. g., in the
dynamics of free electron wave packets in quantum mechanics BLIOKH .

Figure 1:  (a) Intensity $|\psi(r)|^{2}$ of the exploding illumination in (3)
with power $P=10$ W, $\sigma=1$ mm, $s=0$, and $\mu=-0.25$. (b) Intensity of
the unapertured exploded field at the focal plane as given by (4) when
$\omega=2.5$ rad/fs-1, $k=\omega/c=8.33$ mm-1, and $f=200$ mm (solid). (c)
Peak intensity at $r=0$, numerically evaluated from (2) with (3), versus $z$
about the focus. (d) and (e) Radial profiles at planes close to the focal
plane, evaluated from (2) with (3). For comparison, in all plots the dashed
curves represent the same quantities for Gaussian illumination of the same
power and peak intensity. The vertical lines in (c) delimit the focal region
of the focused Gaussian beam.

Figure 2:  (a) Intensity $|\psi(r)|^{2}$ of the exploding illumination in (3)
with power $P=10$ W, $\sigma=1$ mm, $s=1$, and $\mu=0.25$. (b) Intensity of
the unapertured exploded field at the focal plane as given by (4) when
$\omega=2.5$ rad/fs-1, $k=\omega/c=8.33$ mm-1 (vacuum), and $f=200$ mm. (c)
Gray curve: Peak intensity at the maximum of the bright ring surrounding the
vortex, numerically evaluated from (2) with (3), versus $z$ about the focal
point. Black curve: radius of the maximum of the bright ring versus $z$. (d)
and (e) Radial profiles at planes close to the focal plane evaluated from (2)
with (3). For comparison, in all plots the dashed curves represent the same
quantities for Laguerre-Gauss illumination with $s=1$ and zero radial order,
$\psi(r)\propto r\exp(-r^{2}/a^{2})e^{i\varphi}$, having the same power and
peak intensity.

Figures 1 (a) and (b) compare the exploding illuminating field in (3) at $z=0$
and exploded field in (4) at $z=f$ for $s=0$ and $\mu=-0.25$ (solid curves)
with illuminating Gaussian beam of the same power and peak intensity and its
focused field (dashed curves). Figures 1(c), (d) and (e) show peak intensities
at $r=0$ versus $z$ and radial profiles at values of $z$ close to $f$
illustrating how the beam profile approaches the singular profile in (4),
compared to the same quantities for the Gaussian illumination (dashed curves).
The singularity is only formed at the focal plane because the fast
oscillations $\exp\left[\frac{-ikr^{\prime
2}}{2}\left(\frac{1}{f}-\frac{1}{z}\right)\right]$ in (2) out of focus makes
the integral to converge to finite values. Compared to the Gaussian
illumination, the singularity develops explosively in a tiny axial region of
the standard depth of focus [vertical lines in (c)].

Thus, as already pointed out in AIELLO  and AIELLO2 , the exploding profile in
(3) reproduces, to a certain extent, what happens when a plane wave of
infinite lateral extent and power is ideally focused, namely, both focused
fields result in an infinitely bright point at the focus, with the substantial
difference that the exploding profile carries finite power, and is therefore
physically realizable, in the same sense that a standard Gaussian beam is.
This is the cylindrically symmetric counterpart of the concentrating beam
factorized in $x$ and $y$ described and realized in AIELLO  and AIELLO2 .

Also with cylindrical symmetric intensity, Figs. 2(a) and (b) show
illuminating exploding and exploded transversal profiles of vortex beams with
orbital angular mementum given by (3) and (4) with $s=1$ and $\mu=0.25$ (solid
curves), compared to illuminating Laguerre-Gauss beam with $s=1$ (and zero
radial order) having the same power and peak intensity, and its focused field
(dashed curves). The peak intensity at the radius of the bright ring
surrounding the vortex grows up explosively with $z$ to infinity at $f$, as
seen in Fig. 2(c) (solid gray curve), and the same time that the bright ring
shrinks acceleratedly down to zero to form a punctual vortex at $f$, as also
seen in Fig. 2(c) (solid black curve). For comparison, Fig. 2(c) also shows
the peak intensity and radius of the bright ring for the Laguerre-Gauss
illumination (dashed curves). The ripples in the explosively growing and
shrinking radial profiles at planes close to the focal plane, observed in Fig.
2(d) and (e), disappear in the smoother, singular profile at the focal plane.

##  III Observable effects of the exploding behavior

Figure 3:  Physical properties of finite-aperture exploding beams. (a) Focal
intensity when exploding profiles (3) with $s=0$, $\sigma=1$ mm, power $P=10$
W, frequency $\omega=2.5$ fs-1 and the indicated values of $\mu$ illuminate a
lens of focal length $f=200$ mm and increasing radius $R$. (b) The same in
vertical logarithmic scale to visualize the much smaller focal intensity
produced by Gaussian illumination of the same power and peak intensity on the
lens. (c) FWHM diameter of the focal spot produced by the exploding profile
(3), by the Gaussian profile, and by a plane wave, as functions of the radius
of the lens. (d) Radial profiles of intensity at the focal plane due to
exploding (solid), Gaussian of the same peak intensity and power (dashed), and
plane wave (dotted) illuminations for lens radii $R=51.65$ mm (0.25 NA) and
$R=115.47$ mm (0.5 NA). For the Gaussian illumination there is no change with
$R$.

The above description is an ideal behavior that originates from ignoring the
finite aperture of focusing systems, and that disappears when finite
transversal extents are considered. Still, the existence of these singular
beams in a linear medium has observable effects in experiments that are absent
with standard illuminating beams. If a Gaussian or any other standard light
beam carrying finite power is focused, the width of the focal spot and its
intensity do not change once the lens aperture is appreciably larger than the
beam spot size; by contrast, focusing the exploding field in (3) produces, in
principle, a brighter and narrower focal spot as the aperture radius $R$ of
the focusing system is opened. This behavior also mimics what happens when a
plane wave illuminates the lens of increasing radius $R$, with the important
difference that the focused power does not change significantly when opening
the aperture (once $R\gg\sigma$), since the power in the beam tails is
increasingly negligible. In addition, focusing the exploding vortex fields in
(3) produces increasingly brighter and narrower ring surrounding the vortex at
the focal plane when opening the aperture without significantly changing the
beam power. These effects are described below under conditions of paraxial
focusing, conditions under which they are supposed to be physically valid.

For beams without orbital angular momentum ($s=0$, hence $-1/2<\mu<0$), the
intensity at the focal point $(r,z)=(0,f)$ with a lens of aperture radius $R$
can easily be evaluated from (2) with (3) to be

$|\psi(0,f)|^{2}=\frac{P}{A}\left(\frac{k\sigma^{2}}{2f|\mu|}\right)^{2}\left[1-\left(1+\rho^{2}\right)^{|\mu|}\right]^{2}\,,$
(7)

where $\rho=R/\sigma$, which is seen in Fig. 3(a) to grow without bound with
$R$ for any allowed value of $\mu$ even when the aperture radius $R$ is much
larger than the beam spot size of (3) on the lens, measured by $\sigma$. The
curves end when the angle of a marginal ray as seen from the focus is
$30^{\circ}$, taken, according to Siegman SIEGMAN , as a limit of validity of
the paraxial approximation, i. e., for a numerical aperture of NA $=0.5$.
Evaluation of the derivative of (7) with respect to $\mu$ reveals more
pronounced intensity enhancements (at large values $R/\sigma$) when using
values of $\mu$ about $-0.42$, as can also be appreciated in Fig. 3(a). By
contrast, the intensity at the focal point of Gaussian illumination approaches
a constant value as soon as the aperture radius $R$ is slightly larger than
its spot size, as seen in Fig. 3(b). Note the logarithmic vertical scale to
visualize the tiny focal intensity of the focused Gaussian illumination of the
same power and peak intensity. Intensities two orders of magnitude higher than
those attainable with the Gaussian illumination are obtained while keeping
paraxial focusing conditions. Also, the diameter of the focal spot decreases
as the aperture radius is increased up to the limit imposed by the paraxial
approximation, in contrast with the constant and much larger focal diameter
with the Gaussian beam of the same power and peak intensity, as seen in Fig.
3(c). Bright spots of area two orders of magnitude smaller than with the
Gaussian illumination are obtained. Two examples of focused transversal
profiles with numerical apertures $0.5$ and $0.25$ are compared to the low-
intensity and much broader focused Gaussian beam in Fig. 3(d).

Figures 3(c) and (d) include also diffraction-limited focal diameters and
transversal profiles (dotted curves) corresponding to plane wave illumination
on the lens of radius $R$. It is clear that focusing the exploding profile
with finite power imitates focusing of a plane wave, in the sense that
increasing the numerical aperture results in tighter focal spots. The peak
intensities of the focused transversal profiles with illumination (3) and a
plane wave are equated in Fig. 3(d) to visualize that the width of the focal
spot is, for any aperture radius, only slightly above the diffraction limit,
in contrast to what happens with the Gaussian illumination.

Figure 4:  Physical properties of finite-aperture exploding vortex beams. (a)
Focal intensity when exploding profiles (3) with $|s|=1$, $\sigma=1$ mm, power
$P=10$ W, frequency $\omega=2.5$ fs-1 and $\mu=0.08$ illuminate a lens of
focal length $f=200$ mm and increasing radius $R$ (solid curve) compared to
the same quantity for Laguerre-Gauss (LG) illumination of the same power and
peak intensity on the lens. (b) Radius of the maximum of the bright ring at
the focal plane produced by the exploding profile (3) (solid), by the
Laguerre-Gauss profile (dashed), and by an uniform plane wave with a punctual
vortex (dotted), as functions of the radius of the lens. The gray curves of
the same type represent the inner radius and outer radius (half width at half
maximum of the maximum intensity) for the corresponding illuminations. (c)
Radial profiles of intensity at the focal plane due to exploding (solid),
Laguerre-Gauss of the same peak intensity and power (dashed), and vortex plane
wave (dotted) illuminations for lens radii $R=51.65$ mm (0.25 NA) and
$R=115.47$ mm (0.5 NA). For the Gaussian illumination there is no change with
$R$. (d) Intensity of the bright ring and its radius as a function of
propagation distance about the focus for the exploding beam with radius
$R=115.47$ mm (solid curves) compared to the ideal exploding beam. The axial
region is so short that the same quantities for Laguerre-Gauss beams do not
appreciably change (dashed lines).

For exploding beams with orbital angular momentum ($|s|\neq 0$, hence
$(|s|-1)/2<\mu<|s|/2$) the intensity of the bright ring at the focal plane
grows similarly when opening the aperture, in contrast to what happens to
Laguerre-Gauss beams, as seen in Fig. 4(a). As with $s=0$, this enhancement is
further magnified by using the value $\mu\simeq(|s|-1)/2+0.08$, e. g., the
value $\mu=0.08$ for $|s|=1$ chosen in Fig. 4. Again focal intensities about
two orders of magnitude higher than with Laguerre-Gaussian illumination of the
same power, peak intensity and vorticity on the lens are easily attained.
Figure 4(b) shows the shrinking radius of the bright ring at the focal plane
when increasing the aperture radius $R$, in contrast to the constant radius
with the Laguerre-Gauss illumination, and mimicking the shrinking diffraction-
limited radius for plane wave illumination with a punctual vortex (solid,
dashed and dotted curves, respectively). The gray curves of the same type
locate the inner and outer radius of the corresponding bright rings, in order
to visualize the constant thickness of the bright ring at the focal plane for
the Laguerre-Gauss illumination, and the shrinking thicknesses for exploding
and plane wave illuminations. Two examples of transversal profiles at the
focal plane with numerical apertures $0.5$ and $0.25$ are seen in Fig. 4(c),
where they are compared to the much broader and less intense focused Laguerre-
Gauss, and are also seen to be slightly broader than the diffraction-limited
profiles of the same numerical apertures. Interestingly, the exploding beam
focused by a lens of finite radius behaves very approximately as the ideal
exploding beam from an infinite aperture out of the focal plane, as
illustrated in Fig. 4(d) plotting the intensity of the maximum of the bright
ring and its radius about the focus. Only in a tiny fraction (about one wave
length) of the focal region (of about $0.5$ mm) the bright ring intensity and
its radius depart from the ideal behavior, reaching at the focus finite and
non-zero values, respectively.

##  IV Enhanced longitudinal component

The ideally infinite gradients of these singular beams and vortex beams at the
focal plane suggests the presence of a strong longitudinal or axial component
of the electric field, even under paraxial conditions, since the axial
component is directly related to the gradient of the transversal part.
Enhancement of the axial component commonly relies on nonparaxial focusing and
the use of specific fields such as radially polarized beams QUABIS ; COLLIN .
While paraxial, standard, radially polarized fields present quite small
longitudinal components, exploding radially polarized fields have longitudinal
components comparable in amplitude to that of the transversal component.

For a paraxial beam of transversal components
$\bm{\psi}_{\perp}=\psi_{x}\bm{u}_{x}+\psi_{y}\bm{u}_{y}$, the axial component
can be evaluated from LAX

$\psi_{z}=\frac{i}{k}\nabla_{\perp}\cdot\bm{\psi}_{\perp}=\frac{i}{k}\left(\frac{\partial\psi_{x}}{\partial
x}+\frac{\partial\psi_{y}}{\partial y}\right)\,.$ (8)

Strictly speaking, and according to Gauss law for the electric field, the term
$(i/k)\partial\psi_{z}/\partial z$ should be included in the right hand side
of (8), but according to the perturbative method in LAX , the $x$ and $y$
derivatives are leading-order terms determining the axial component under
paraxial conditions, and the $z$ derivative is a small correction. In
addition, the axial component is evaluated below only at the focal plane,
where it takes maximum values and therefore $\partial\psi_{z}/\partial z=0$.

In order to maintain the cylindrical symmetry of the intensity of the
vectorial field, we consider left or right handed circular polarizations,
$\bm{\psi}_{\perp}=\psi(r,z)e^{is\varphi}\bm{u}_{l,r}$, where
$\bm{u}_{l,r}=(\bm{u}_{x}\pm i\bm{u}_{y})/\sqrt{2}$, and the upper and lower
signs stand for left and right handed polarizations, respectively, and where
$\psi(r,z)$ is given by Fresnel integral in (2) with (3). Some calculation of
derivatives leads to

$\psi_{z}(r,z)e^{i(s\pm
1)\varphi}=\frac{1}{\sqrt{2}}\frac{i}{k}\left[\frac{\partial\psi(r,z)}{\partial
r}\mp\frac{s}{r}\psi(r,z)\right]e^{i(s\pm 1)\varphi}\,.$ (9)

With $s=-1$ (hence $0<\mu<1/2$) and left handed circular polarization, and
with $s=+1$ (also $0<\mu<1/2$) and right handed polarization, the axial
components do not carry any vorticity, can then have maxima at the beam center
(where the transversal components vanish), which are equal in both cases and
given by $\psi_{z}(r,z)=(i/\sqrt{2}k)\left[\partial\psi(r,z)/\partial
r+\psi(r,z)/r\right]$.

We can now easily consider the radially or azimuthally polarized beams

$\bm{\psi}_{\perp}=\frac{\psi(r,z)}{\sqrt{2}}\left[e^{-i\varphi}\bm{u}_{l}\pm
e^{i\varphi}\bm{u}_{r}\right]\,,$ (10)

having the same amplitude associated with the transversal components,
$|\bm{\psi}_{\perp}|=|\psi(r,z)|$, as that of the individual vortex beams with
circular polarizations. The axial component vanishes for azimuthal
polarization [for the minus sign in 10], and is given by

$\psi_{z}=\frac{i}{k}\left[\frac{\partial\psi(r,z)}{\partial
r}+\frac{1}{r}\psi(r,z)\right]=\frac{i}{k}\frac{1}{r}\frac{\partial[r\psi(r,z)]}{\partial
r}\,,$ (11)

for radial polarization [for the plus sign in 10]. Writing the ideal focal
profile in (4) for $|s|=1$ in the compact form
$\psi(r,f)=C\alpha^{\mu}K_{1-\mu}(\alpha)$, where $\alpha=k\sigma r/f$ and $C$
is the constant in the first row in (4), and where we have neglected the
negligible curvature factor $e^{ikr^{2}/2f}$ in the Debye approximation
(symmetric focused field about the focal plane), one can readily evaluate the
axial component in (11) as

$\psi_{z}=\frac{i\sigma}{f}C\alpha^{\mu}\left[K_{2-\mu}(\alpha)-2\frac{K_{1-\mu}(\alpha)}{\alpha}\right]\,,$
(12)

which is singular at $\alpha=0$, i. e., at $r=0$. Indeed the ratio of
amplitudes of the axial and transversal components close to the origin behaves
as $|\psi_{z}|/|\psi|=2\sigma\mu/f\alpha=2\mu/kr$, implying that the
singularity in the axial field is stronger than that of the transversal field.
Of course, these are mathematical equations for ideal exploding radially
polarized beams, but they have physical manifestations in their apertured
versions. From (11), (2) with finite radius $R$ and illumination (3) with
$|s|=1$, the axial component of the radially polarized beam at the focal plane
can be evaluated to be

$\psi_{z}(0,f)=\frac{ik}{f^{2}}\sqrt{\frac{P}{A}}\frac{\sigma^{3}}{2\mu}\left[\frac{\rho^{2}}{(1+\rho^{2})^{\mu}}+\frac{1-(1+\rho^{2})^{1-\mu}}{1-\mu}\right]\,,$
(13)

where again $\rho=R/\sigma$. The amplitude of the axial component grows as
$(R/\sigma)^{2-2\mu}$ at large $R/\sigma$, while the amplitude of the
transversal component at the bright ring surrounding the polarization
singularity grows with the potential law $(R/\sigma)^{1-2\mu}$ (also at large
$R/\sigma$) of lower power, as can be seen in Fig. 5(a). This difference
implies a linear increase with $R/\sigma$ of the ratio of the axial and
transversal amplitudes, as can be appreciated in Fig. 5(b), that reaches a
maximum value of $0.665$ in the limit of validity of the paraxial
approximation (solid curve). As expected, for Laguerre-Gauss radially
polarized beams, this ratio reaches a constant maximum value that does not
exceed $0.049$ (dashed curve), and, again, the behavior of the axial component
of the exploding radially polarized beam mimics the behavior of radially
polarized uniform plane wave with a punctual polarization singularity in its
center, whose maximum ratio is about $0.757$ (dotted curve). For illustration
Figs. 5(c) and (d) depict transversal profiles of amplitude of the axial and
transversal components at the focal plane for NA $=0.25$ and $0.5$. The above
analysis indicates that the individual, circularly polarized exploding vortex
beams have also enhanced axial components, but we have ignored them because
they simply are $\sqrt{2}$ smaller than that of the radially polarized beam.

Figure 5:  Exploding radially polarized beams. (a) Amplitude at the focal
point of the axial component $|\psi_{z}(0,f)|$ (black curve), as given by
(13), and of the transversal component
$|\bm{\psi}_{\perp}(r_{m},f)|=|\psi(r_{m},f)|$ at its radial maximum $r_{m}$
(gray curve), when the radially polarized beam in (10), with the exploding
profile $\psi(r)$ in (3) with $\sigma=1$ mm, $P=10$ W, frequency $\omega=2.5$
fs-1, $\mu=0.08$ and $|s|=1$, is focused in vacuum with a lens of focal length
$f=200$ mm, both as functions of the radius $R$ of the lens aperture. (b)
Quotient of the axial and transversal amplitudes as a function of $R$ (solid
curve) for the above conditions, for radially polarized Laguerre-Gauss
illumination of the same power and peak intensity on the lens (dashed curve),
and for radially polarized plane wave illumination with a punctual
polarization singularity in its center (dotted curve). (c) and (d) Radial
profiles of the amplitude of the axial (black) and transversal (gray)
components for $R=51.65$ mm (NA= $0.25$) and $R=115.47$ mm (NA $=0.5$).

##  V Exploding versus uniform illumination

Figure 6:  Solid curves: For the focusing conditions and exploding
illumination as in Fig. 3, (a) peak intensity at focus and (b) FWHM of the
focal spot as functions of the lens aperture $R$. Dotted curves: The same but
for top-hat illumination of the same power of radius $R$, i. e., always
filling the aperture. Dashed curves: The same for top-hat illumination of the
same power and fixed radius $R_{0}=30$ mm.

From what is exposed so far, it arises the practical question of the possible
advantages of using exploding illumination over uniform illumination (of
finite transversal extent). Of course the response depends on the particular
application. Smooth profiles are more convenient in nonlinear optics
applications to avoid instabilities and collapse. Reshaping Gaussian-like beam
profiles into uniform flat-top profiles involves large power losses if based
on truncation; otherwise it requires using sophisticated diffractive and
refractive optics DICKEY , or more recently developed techniques FRIES , but
uniform illumination provides better resolution MAZNEV  in applications such
as microscopy.

Suppose here we dispose of the fixed amount of power $P$ supplied by a laser
system which can be shaped as the exploding profile or as a flat-top profile,
and wish to increase the intensity and/or diminish the size of the focal spot
by opening the aperture, or simply to control these properties. As seen in
Fig. 6 (a) and (b) the intensity is considerably higher and the spot size
slightly smaller (the diffraction limit) with a flat-top profile that always
fills the aperture (dotted curves) than with exploding illumination (solid
curves), but this arrangement requires reshaping the flat-top profile each
time the aperture radius is increased. The situation is different in a setting
where reshaping of the illumination is difficult or not possible (dashed
curves). The peak intensity of the flat-top illumination of the power $P$ and
fixed radius $R_{0}$ is at first lower than with the exploding illumination
because of truncation losses, then higher when more power is collected, and
finally constant once $R_{0}$ is smaller than the aperture radius $R$,
becoming eventually smaller than the focal intensity with the exploding
illumination. Similarly, the focal spot is tighter but then constant and wider
than with the exploding illumination of increasing aperture radius. Thus, the
exploding profile offers the possibility of a smoother control of the
intensity and width of the focused spot over wider ranges without the
necessity of reshaping the input beam. Similar conclusion holds for the
intensity and thickness of the bright ring of exploding vortex beams, and for
the strength of the axial component of exploding radially polarized beams.

##  VI Conclusions

In conclusion, we have reported the existence of cylindrically symmetric,
paraxial beams, vortex beams, and radially polarized beams of light that
produce singular (infinite) intensities when they are ideally focused. They
mimic the behavior of ideally focused plane waves but with a localized
transversal profile that carries a finite amount of power.

With real focusing systems of finite lateral extent, this peculiar behavior
manifests as focused beams whose intensity and spot size, vortex intensity and
radius, and strength of the longitudinal component strongly change with the
aperture size, even if the focused power is not significantly altered. With
the same exploding illumination, focal intensities exceeding by two orders of
magnitude, focal spots and vortex radii smaller by two orders of magnitude,
and longitudinal component one order of magnitude higher than the same
properties for standard illuminating beams of similar power and intensity can
be achieved by increasing the aperture radius while keeping paraxial focusing.

These effects are absent with standard beams, and can find application in
linear or nonlinear optics experiments where a precise control of the width
and intensity of the focal spot, of the radius and intensity of the vortex
bright ring, or the strength of the axial component of the electric field, are
crucial, e. g., in second harmonic generation, also with radially polarized
light YEW  and in surfaces BISS , acceleration of electrons GUPTA ; VARIN ,
particle trapping ZHAN , laser material processing DREVINSKAS , etc. These
exploding beams can also be used as alternative beams to standard Laguerre-
Gauss beams for quantum entanglement of states of orbital angular momentum
FICKLER . Outside the field of optics, given the generality of the Schrödinger
equation, exploding free electron wave packets BLIOKH  and vortex electron
beams MCMORRAN  do exist and could find application in electron microscopy, as
well as exploding acoustic vortices for Mie particle trapping MARZO .

The author acknowledges support from Projects of the Spanish Ministerio de
Economía y Competitividad No. MTM2015-63914-P, and No. FIS2017-87360-P.

## References

  * (1) A. Peres, Quantum theory: Concepts and methods (Kluver Academic Publishers, San Diego, CA, 1994), 5th ed., pp. 81–82.
  * (2) A. Aiello, “Spontaneous generation of singularities in paraxial optical fields,” Opt. Lett. 41, 1668–1671 (2016).
  * (3) A. Aiello, M. Paúr, B. Stoklasa, Z. Hradil, J. Rehacek, and L. L. Sánchez-Soto, “Observation of concentrating paraxial beams,” OSA Continuum 3, 2387–2394 (2020).
  * (4) Y. S. Kivshar and D. E. Pelinovsky, “Self-Focusing and Transverse Instabilities of Solitary Waves,” Phys. Rep. 331, 117–195 (2000).
  * (5) G. Fibich, The nonlinear Schrödinger Equation: Singular solutions and optical collapse (Springer: Heidelberg, 2015).
  * (6) I. S. Gradshteyn and I. M. Ryzhik,Table of Integrals, Series, and Products, Fifth Edition (Academic, 1994).
  * (7) K.Y. Bliokh, Y.P. Bliokh, S. Savel’ev, and F. Nori, “Semiclassical Dynamics of Electron Wave Packet States with Phase Vortices,” Phys. Rev. Lett. 99, 190404 (2007).
  * (8) A. E. Siegman, Lasers, (Mill Valley, 1986). See pp. 628-630.
  * (9) S. Quabis, R. Dorn, M. Eberler, O. Glöckl, and G. Leuchs, “Focusing light to a tighter spot,” Opt. Commun. 179 1–7 (2000).
  * (10) H Wang, L. Shi, B. Lukyanchuk, and C. Sheppard, “Creation of a needle of longitudinally polarized light in vacuum using binary optics,” Nature Photonics 2, 501–505 (2008).
  * (11) M. Lax, W. H. Louisell, and W. B. McKnight, “From Maxwell to paraxial wave optics,” Phys. Rev. A 11, 1365–1370 (1975).
  * (12) F. M. Dickey, L. S. Weichman, and R. N. Shagam, “Laser beam shaping techniques,” Proc. SPIE 4065, 338–348 (2000).
  * (13) F. Fries, M Fröbel, P. Y. Ang, S. Lenk, and S. Reineke, “Real-time beam shaping without additional optical elements,” Light Sci. Appl. 7, 18 (2018).
  * (14) A. A. Maznev and O. B. Wright, “Upholding the diffraction limit in the focusing of light and sound,” Wave Motion 68, 182–189 (2017).
  * (15) E.Y.S. Yew and C. J. R. Sheppard, “Second harmonic generation polarization microscopy with tightly focused linearly and radially polarized beams,” Opt. Commun. 275, 453–457 (2007).
  * (16) D. P. Biss and T. G. Brown, “Polarization-vortex-driven second-harmonic generation,” Opt. Lett. 28, 923–925 (2003).
  * (17) D. N. Gupta, N. Kant, D. E. Kim, and H. Suk, “Electron acceleration to GeV energy by a radially polarized laser,” Phys. Lett. A 368, 402–407 (2007).
  * (18) C. Varin, M. Piché, M. A. Porras, “Acceleration of electrons from rest to GeV energies by ultrashort transverse magnetic laser pulses in free space,” Phys. Rev. E 71, 026603 (2005).
  * (19) Q. Zhan, “Trapping metallic Rayleigh particles with radial polarization,” Opt. Express 12, 3377–3382 (2004).
  * (20) R. Drevinskas, J. Zhang, M. Beresna, M. Gecevicius, A. G. Kazanskii, Y. P. Svirko, and P. G. Kazansky, “Laser material processing with tightly focused cylindrical vector beams,” Appl. Phys. Lett. 108, 221107 (2016).
  * (21) R. Fickler, R. Lapkiewicz, W. N. Plick, M. Krenn, C. Schaeff, S. Ramelov, and A. Zeilinger, Science 338, 640 (2012).
  * (22) B. J. McMorran, A. Agrawal, I. M. Anderson, A. A. Herzing, H. J. Lezec, J. J. McClelland, and J. Unguris, “Electron vortex beams with high quanta of orbital angular momentum, ” Science 331, 192–195 (2011).
  * (23) A. Marzo, M. Caleap, and B. W. Drinkwater, “Acoustic Virtual Vortices with Tunable Orbital Angular Momentum for Trapping of Mie Particles,” Phys. Rev. Lett. 120, 044301 (2018).


11institutetext:  1: Department of Physics, University of California, Davis,
CA 9516, USA; 11email<EMAIL_ADDRESS>
2: Department of Electrophysics, National Chiao Tung University, Hsinchu
30010, Taiwan; 11email<EMAIL_ADDRESS>
3: National High Magnetic Field Laboratory, Florida State University,
Tallahassee, Florida 32306-4005, USA; 11email<EMAIL_ADDRESS>
4: Theoretical Division, Los Alamos National Laboratory, Los Alamos, New
Mexico 87545, USA; 11email<EMAIL_ADDRESS>
# Hyperfine fields and magnetic structure in the B phase of CeCoIn5

Nicholas J. Curro1 Ben-Li Young2 Ricardo R. Urbano3 Matthias J. Graf4

###### Abstract

We re-analyze Nuclear Magnetic Resonance (NMR) spectra observed at low
temperatures and high magnetic fields in the field-induced B-phase of CeCoIn5.
The NMR spectra are consistent with incommensurate antiferromagnetic order of
the Ce magnetic moments. However, we find that the spectra of the In(2) sites
depend critically on the direction of the ordered moments, the ordering
wavevector and the symmetry of the hyperfine coupling to the Ce spins.
Assuming isotropic hyperfine coupling, the NMR spectra observed for
$\mathbf{H}~{}||~{}[100]$ are consistent with magnetic order with wavevector
$\mathbf{Q}=\pi(\frac{1+\delta}{a},\frac{1}{a},\frac{1}{c})$ and Ce moments
ordered antiferromagnetically along the [100] direction in real space. If the
hyperfine coupling has dipolar symmetry, then the NMR spectra require Ce
moments along the [001] direction. The dipolar scenario is also consistent
with recent neutron scattering measurements that find an ordered moment of
0.15$\mu_{B}$ along [001] and
$\mathbf{Q_{n}}=\pi(\frac{1+\delta}{a},\frac{1+\delta}{a},\frac{1}{c})$ with
incommensuration $\delta=0.12$ for field $\mathbf{H}~{}||~{}[1\bar{1}0]$.
Using these parameters, we find that a hyperfine field with dipolar
contribution is consistent with findings from both experiments. We speculate
that the B phase of CeCoIn5 represents an intrinsic phase of modulated
superconductivity and antiferromagnetism that can only emerge in a highly
clean system.

###### Keywords:

NMR superconductivity heavy fermion magnetism

###### pacs:

P76.60.-k 75.30.Fv 74.10.+v

††journal: Journal of Low Temperature Physics

##  1 Introduction

The heavy-fermion superconductor CeCoIn5 has attracted considerable attention
since its discovery in 2001.1 Not only does this unconventional d-wave
superconductor exhibit non-Fermi liquid behavior associated with proximity to
a proposed quantum critical point, but it is unique among the heavy-fermion
superconductors in that it also exhibits a new thermodynamic phase (B phase)
that exists only within the superconducting phase near $H_{c2}$.2, 3 Initially
this B phase was identified as the elusive Fulde-Ferrell-Larkin-Ovchinnikov
(FFLO) superconducting phase first predicted to exist in Pauli-limited
superconductors over 40 years ago.4, 5, 6, 7, 8 In fact, recent NMR work by
Young and coworkers9 identified the presence of incommensurate
antiferromagnetic order in the B phase in contrast to the standard predictions
for the FFLO phase.10, 11 Signatures of magnetism were also seen in other NMR
experiments.12, 13 The NMR spectra of the In sites in the B phase do not
reveal additional paramagnetic resonances that might be associated with
macroscopic phase separation, but rather are consistent with superconducting
order coexisting with a homogeneous modulation of antiferromagnetic order.9
Despite initial arguments to the contrary,14 recent neutron scattering results
by Kenzelmann and coworkers now provide conclusive proof for long-range static
incommensurate antiferromagnetic order.15

Figure 1: (Color online) The phase diagram of CeCoIn5 in high field as
determined by specific heat.5 Solid points represent second order phase
transitions and open points are first order transitions. The solid blue
squares are the points at which the spectra in Fig. 3 were obtained.

Figure 2: (Color online) The unit cell of CeCoIn5. The Ce atoms (yellow) sit
at the eight corners. The In(1) atoms sit in the center of the top and bottom
faces (orange). The Co atoms are grey and the In(2) atoms are green. For the
field oriented in the $ab$ plane, there are two inequivalent In(2) atoms,
depending on whether the field is parallel (In(2a) or perpendicular (In(2b))
to the unit cell face.

The antiferromagnetism in CeCoIn5 was first identified by Young et al. due to
the presence of a broad spectrum observed at the In(2) sites in this material
(see Figs. 2 and 3). The In(1) and Co sites, in contrast, showed no splitting.
Young et al. pointed out that these observations place constraints on the
possible magnetic structure, but do not uniquely identify the structure. They
proposed a minimal model where the magnetic structure consists of ordered
local Ce spins with moments $\mathbf{S}_{0}$ along the applied magnetic field
direction (along [100]), with an ordering wavevector of the form
$\mathbf{Q}=\pi(\frac{1+\delta}{a},\frac{1}{a},\frac{1}{c})$. The structure of
the NMR spectra revealed the incommensurate nature, but the value of the
modulation $\delta$ remained undetermined since the hyperfine field at the
In(2) site depends on the product of the size of the ordered moment and the
incommensuration. In the neutron diffraction experiment, Kenzelmann et al.
oriented the field along the [1$\bar{1}$0] direction, and observed
$\mathbf{Q_{n}}=\pi(\frac{1+\delta}{a},\frac{1+\delta}{a},\frac{1}{c})$ with
$\delta=0.12$ and moments along [001]. A crucial observation was that $\delta$
is independent of the applied field in the B phase, in contrast to the
predictions for the FFLO phase. By proposing a Ginzburg-Landau model for the
coupling of antiferromagnetism and superconductivity, they showed that the
superconducting order parameter in the B phase acquires a component with
finite momentum because of the strong coupling between the incommensurate
magnetism and superconductivity. The neutron data confirm the NMR observation
that this exotic state disappears immediately above $H_{c2}$, where the system
returns to a fully homogeneous normal phase, yet with strong deviations from
conventional Fermi liquid theory.3

A priori, the NMR and neutron scattering results suggest different magnetic
structures. In order to address this discrepancy, we investigate several
possible magnetic structures allowed by the NMR results. The neutron
diffraction results suggest that the applied field, $\mathbf{H}$, the moments,
$\mathbf{S}_{0}$, and the incommensuration wavevector
$\mathbf{Q}_{i}=\frac{\pi}{a}(\delta,\delta,0)$ are all mutually orthogonal.
In contrast, the proposed NMR scenario suggested
$\mathbf{S}_{0}~{}||~{}\mathbf{H}~{}||~{}\mathbf{Q}_{0}$. As we show below,
this scenario is the most likely for isotropic transferred hyperfine couplings
between the Ce spins and the In(2) nuclei. On the other hand, if the coupling
tensors are anisotropic, then other magnetic structures are possible as argued
by Koutroulakis et al.16, 17 If the hyperfine tensor has purely dipolar
symmetry, then we find that for $\mathbf{H}~{}||~{}[100]$ the most likely
magnetic structures satisfy $\mathbf{S}_{0}~{}||~{}[001]$ and
$\mathbf{Q}_{i}\perp\mathbf{S}_{0}$. The spectra of the In(1) and In(2) differ
slightly depending on the orientation of $\mathbf{Q}_{i}$ in the plane, but
the data are most consistent with $\mathbf{Q}_{i}||[010]$ or
$\mathbf{Q}_{i}||[110]$. The anisotropic coupling scenario offers a picture
that is both physically more reasonable and consistent with neutron
diffraction observations of $\mathbf{Q}_{i}\perp\mathbf{S}_{0}$ and
$\mathbf{H}_{0}\perp\mathbf{S}_{0}$, though for fields
$\mathbf{H}_{0}~{}||~{}[1\bar{1}0]$.

##  2 Analysis

Figure 3: (Color online) Fixed field NMR spectra at 11.1 T in CeCoIn5 showing
how the In(1), In(2) and Co sites evolve from the normal state (top) to the B
phase (bottom), adapted from Young et al.9

###  2.1 Spectra

In order to explain the broad double-peak structure of the In(2) spectrum in
Fig. 3, there must be a distribution of local fields both parallel and
antiparallel to the applied field $\mathbf{H}_{0}||[100]$ with values ranging
up to 1.3 kOe.9 The resonance frequency is given by
$f=\gamma|\mathbf{H_{0}}+\mathbf{H}_{hf}|+f_{Q}$, where $\gamma$ is the
gyromagnetic ratio and $f_{Q}$ is the contribution from the quadrupolar
interaction at the nucleus. Since the electric field gradient (EFG) at all
three nuclear sites is unaffected by the onset of superconductivity or
magnetism, the dramatic line broadening effects observed in the B phase can be
attributed entirely to the onset of the static hyperfine field,
$\mathbf{H}_{hf}$. $f_{Q}$ is a temperature independent constant that depends
on the particular site and we will not address it further. Experimentally, we
find no significant broadening in the B phase for the In(1) or the Co, but a
broad, double-peak spectrum for the In(2a) site (previously referred to as
In(2)||, the In(2) site on the unit cell face that lies parallel to the
field). Independent measurements show no significant broadening at the In(2b)
(previously referred to as In(2)⟂). These results put stringent constraints on
any candidate magnetic structure.

The double-peak structure observed for the In(2a) arises because there is a
distribution of local hyperfine fields that lie either parallel or
antiparallel to the applied field. If $\mathbf{H}_{hf}$ is parallel to
$\mathbf{H}_{0}$ and is modulated along a direction $\hat{r}$, then
$f(r)=f_{0}+\gamma h_{hf}^{0}\cos(qr)$, where $h_{hf}^{0}$ is the magnitude of
the modulation, $r$ is the distance along the modulation, and $q$ is the
wavevector. In this case, the spectrum will then be given by

$\mathcal{P}_{||}(f)\propto|df/dr|^{-1}=\frac{q^{-1}}{\sqrt{\gamma^{2}(h_{hf}^{0})^{2}-(f-f_{0})^{2}}}.$
(1)

On the other hand, if $\mathbf{H}_{hf}(r)~{}\perp~{}\mathbf{H}_{0}$ then
$f(r)=\sqrt{f_{0}^{2}+\gamma^{2}(h_{hf}^{0})^{2}\cos^{2}(qr)}$ and the
spectrum is given by

$\mathcal{P}_{\perp}(f)\propto\frac{q^{-1}f}{\sqrt{f^{2}-f_{0}^{2}}\sqrt{\gamma^{2}(h_{hf}^{0})^{2}-f^{2}-f_{0}^{2}}}.$
(2)

These spectra are shown in Fig. 4, and clearly show that for the parallel
case, there are two double peaks at frequencies both below and above $f_{0}$,
whereas for the perpendicular case, there is only a single peak at higher
frequency. A key result of the re-analysis of the NMR spectra is, in agreement
with our earlier analysis, that the In(2a) sites require a hyperfine field
parallel to the applied magnetic field,
$\mathbf{H}_{hf}~{}||~{}\mathbf{H}_{0}$, to account for the broadening and
double-peak spectrum. On the other side, the In(2b) sites show little or no
broadening, consistent with a vanishing or perpendicular hyperfine field,
i.e., at the In(2b) sites either $\mathbf{H}_{hf}=0$ or
$\mathbf{H}_{hf}\perp\mathbf{H}_{0}$.

Figure 4: (Color online) The calculated spectra of $\mathcal{P}_{||}(f)$
(blue) and $\mathcal{P}_{\perp}(f)$ (red) assuming $f_{0}$= 100 MHz and
$\gamma h_{hf}^{0}$ = 1 MHz. For the perpendicular case, the spectrum is only
weakly affected by the hyperfine fields, whereas for the parallel case it
broadens dramatically resulting in a double-peak structure.

###  2.2 Hyperfine couplings

The hyperfine interaction is given by the Hamiltonian

$\mathcal{H}_{hf}=\mathbf{\hat{I}}\cdot\mathbb{A}\cdot\mathbf{S}(\mathbf{r}=\mathbf{0})+\sum_{i}\mathbf{\hat{I}}\cdot\mathbb{B}_{i}\cdot\mathbf{S}(\mathbf{r}_{i}),$
(3)

where the hyperfine coupling tensor $\mathbb{A}$ represents the on-site
coupling to an electron spin, $\mathbf{S}$, at the nuclear site
$\mathbf{r}=0$, $\mathbb{B}$ represents a transferred hyperfine coupling to an
electron spin on a distant (ligand) site at $\mathbf{r}_{i}$.18 In CeCoIn5,
these sites are the nearest neighbor Ce 4$f$ electrons, and the sum is over
the nearest neighbors. For static ordering of the Ce spins, Eq. (3) can be re-
written as:
$\mathcal{H}_{hf}=\gamma\hbar\mathbf{\hat{I}}\cdot\mathbf{H}_{hf}$, where the
magnitude and direction of the hyperfine field $\mathbf{H}_{hf}$ depend
critically on the hyperfine tensors for the particular site and the magnetic
structure. The tensorial $\mathbb{A}$ term represents hyperfine coupling to
the itinerant conduction electrons, which we will ignore, since we are only
concerned with static contributions to $\mathbf{H}_{hf}$ from the static local
Ce ordering.

The transferred hyperfine tensor $\mathbb{B}$ is generally not diagonal in the
crystal axis basis,19 and may be written as the sum of isotropic and dipolar
contributions.20 To lowest order, the tensor can be approximated by a scalar
(isotropic) interaction, since the transferred hyperfine interaction is
typically at least one order of magnitude greater than the direct dipolar
interaction. However, in the CeCoIn5 compound, there is evidence that the
hyperfine interaction is not purely isotropic, and therefore we must consider
dipolar symmetries as well. Indeed, the magnitude of the dipolar portion is
found to be enhanced by the delocalized nature of the electrons in the
solid.21 Therefore, we write the hyperfine fields at the ligand sites as:

$\displaystyle\mathbf{H}_{hf}(\mathbf{r})$ $\displaystyle=$
$\displaystyle\sum_{i=1}^{4}\mathbb{B}_{i}\cdot\mathbf{S}(\mathbf{r}+\mathbf{r_{i}})/\gamma\hbar\mbox{~{}~{}~{}~{}~{}~{}~{}~{}~{}
at In(1) }$ (4) $\displaystyle\mathbf{H}_{hf}(\mathbf{r})$ $\displaystyle=$
$\displaystyle\sum_{j=1}^{2}\mathbb{B}_{i}\cdot\mathbf{S}(\mathbf{r}+\mathbf{r_{j}})/\gamma\hbar\mbox{~{}~{}~{}~{}~{}~{}~{}~{}~{}
at Co }$ (5) $\displaystyle\mathbf{H}_{hf}(\mathbf{r})$ $\displaystyle=$
$\displaystyle\sum_{k=1}^{2}\mathbb{B}_{i}\cdot\mathbf{S}(\mathbf{r}+\mathbf{r_{k}})/\gamma\hbar\mbox{~{}~{}~{}~{}~{}~{}~{}~{}~{}
at In(2a) }$ (6) $\displaystyle\mathbf{H}_{hf}(\mathbf{r})$ $\displaystyle=$
$\displaystyle\sum_{l=1}^{2}\mathbb{B}_{i}\cdot\mathbf{S}(\mathbf{r}+\mathbf{r_{l}})/\gamma\hbar,\mbox{~{}~{}~{}~{}~{}~{}~{}~{}~{}
at In(2b) }$ (7)

where $\mathbf{r}_{i}=(\pm\frac{a}{2},\pm\frac{a}{2},0)$ for the In(1) nearest
neighbor Ce sites, $\mathbf{r}_{j}=(0,0,\pm\frac{c}{2})$ for the Co nearest
neighbor Ce sites, and $\mathbf{r}_{l}=(\pm\frac{a}{2},0,z_{0})$ for the
In(2a) nearest neighbor Ce sites, and
$\mathbf{r}_{k}=(0,\pm\frac{a}{2},z_{0})$ for the In(2b) nearest neighbor Ce
sites. The couplings $\mathbb{B}_{i}=\mathbb{B}_{\rm iso}+\mathbb{B}_{\rm
dip}$, where:

$\mathbb{B}_{\rm iso}=B_{\rm iso}\left(\begin{array}[]{ccc}1&0&0\\\ 0&1&0\\\
0&0&1\end{array}\right)$ (8)

and

$\mathbb{B}_{\rm dip}=\frac{B_{\rm
dip}}{r^{2}}\left(\begin{array}[]{ccc}2x^{2}-y^{2}-z^{2}&3xy&3xz\\\
3xy&-x^{2}+2y^{2}-z^{2}&3yz\\\ 3xz&3yz&-x^{2}-y^{2}+2z^{2}\end{array}\right)$
(9)

are evaluated for each site. Here $\mathbf{r}=(x,y,z)$ is the vector joining
the particular site to the Ce atom in question.

Table 1: The hyperfine fields at the In(1) and In(2) sites for various magnetic structures. Here we assumed that the applied magnetic field is $\mathbf{H}_{0}~{}||~{}[100]$, the moment is 0.15 $\mu_{B}$ and its modulation is $\delta=0.12$. The column at the far right indicates magnetic structures in agreement with the NMR spectra. Case | $\mathbf{Q}_{i}$ | $\mathbf{S}_{0}$ | $\mathbb{B}$ | In(1) | In(2a) | In(2b) | Agreement?
---|---|---|---|---|---|---|---
(1.1) | [100] | [100] | iso | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}=0$ | yes (a)
(1.2) | [100] | [100] | dip | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}=0$ | no
(1.3) | [100] | [001] | iso | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | no
(1.4) | [100] | [001] | dip | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | yes (b)
(2.1) | [010] | [100] | iso | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | no
(2.2) | [010] | [100] | dip | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}\approx 0$ | no
(2.3) | [010] | [001] | iso | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}=0$ | no
(2.4) | [010] | [001] | dip | $\mathbf{H}_{\rm hf}=0$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | yes (c)
(3.1) | [110] | [100] | iso | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | no
(3.2) | [110] | [100] | dip | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | no
(3.3) | [110] | [001] | iso | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | no
(3.4) | [110] | [001] | dip | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | yes (d)
(4.1) | [1$\bar{1}$0] | [100] | iso | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | no
(4.2) | [1$\bar{1}$0] | [100] | dip | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | no
(4.3) | [1$\bar{1}$0] | [001] | iso | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | $\mathbf{H}_{\rm hf}~{}||~{}[001]$ | no
(4.4) | [1$\bar{1}$0] | [001] | dip | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | $\mathbf{H}_{\rm hf}~{}||~{}[100]$ | $\mathbf{H}_{\rm hf}~{}||~{}[010]$ | yes (e)

###  2.3 Magnetic structure

The magnetic structure is given by
$\mathbf{S}=\mathbf{S_{0}}\cos[(\mathbf{Q}_{0}+\mathbf{Q}_{i})\cdot\mathbf{r}]$,
where the antiferromagnetic wavevector
$\mathbf{Q}_{0}=(\frac{\pi}{a},\frac{\pi}{a},\frac{\pi}{c})$ is commensurate
with the lattice. Neutron diffraction reports an incommensurate wavevector
$\mathbf{Q}_{i}=\frac{\pi}{a}(\delta,\delta,0)$ with spatial modulation
$\sqrt{2}a/\delta\approx 12a\approx 5.4$ nm in the $ab$ plane, and ordered
moment $\mathbf{S_{0}}$ at the Ce site. This modulation is significantly
shorter than the average inter-vortex distance of $\sim 14$ nm in a field of
10 T.13 So it does not support a picture of overlapping extended states from
the vortex cores leading to magnetic ordering. Similarly this specific
$\mathbf{Q}_{i}$ would be inconsistent with a field along [100], if it was
related to the alignment of vortices along the [100] direction. In our re-
analysis of the NMR spectra, we consider four cases: (1)
$\mathbf{Q}_{i}~{}||~{}[100]$, (2) $\mathbf{Q}_{i}~{}||~{}[010]$, and (3)
$\mathbf{Q}_{i}~{}||~{}[110]$, (4) $\mathbf{Q}_{i}~{}||~{}[1\bar{1}0]$. Case
(1) was proposed by Young et al.9 for NMR measurements under the condition
$\mathbf{H}_{0}~{}||~{}[100]$. Case (2) should be equally likely as case (1)
because of the tetragonal symmetry of the crystal structure. Case (3) was
proposed by Kenzelmann et al.15 for neutron diffraction measurements under the
condition $\mathbf{H}_{0}~{}||~{}[1\bar{1}0]$, and case (4) should be equally
likely as case (3). We have calculated the hyperfine fields for each of these
cases for both purely isotropic and purely dipolar couplings, for
$\mathbf{H}_{0}~{}||~{}[100]$ with moments along both [100] and [001], and the
results are summarized in Table (1).

The cases that are most consistent with the NMR observations are (1.1)
$\mathbf{Q}_{i}~{}||~{}[100]$, $\mathbf{S}_{0}~{}||~{}[100]$ and isotropic
coupling, (1.4) $\mathbf{Q}_{i}~{}||~{}[100]$, $\mathbf{S}_{0}~{}||~{}[001]$
and dipolar coupling, and (2.4) $\mathbf{Q}_{i}~{}||~{}[010]$,
$\mathbf{S}_{0}||[001]$ and dipolar coupling. Cases (3.4) and (4.4) are
consistent with the In(2) spectra, but give rise to an internal field at the
In(1) site that at first sight is inconsistent with experiment. We will
discuss this in more detail below in the discussion section. Figures 6 and 6
show the magnetic structure and hyperfine fields for cases (1.1) and (2.4).
Case (1.1) is identical to the one originally proposed by Young and
coworkers,9 which most likely will not minimize the magnetic contribution to
the free energy, as the moments are either parallel or antiparallel to the
applied field. Cases (1.4), (2.4), (3.4) and (4.4), in which the moments are
perpendicular to the applied field, are physically more reasonable for
antiferromagnetic ordering and agree with the neutron diffraction results of
$\mathbf{Q}_{i}\perp\mathbf{S}_{0}$ and $\mathbf{H}_{0}\perp\mathbf{S}_{0}$.15

##  3 Discussion

Figure 5: (Color online) Magnetic structure and hyperfine fields for isotropic
hyperfine couplings to the In(2) (case (1.1)). The in-plane tetragonal
structure is outlined in gray. The Ce atoms are yellow, and their moments are
indicated by red arrows pointing along [100]. The In(1) atoms are orange, and
the Co are not shown. The In(2a) are green and the In(2b) are blue. The
hyperfine fields at the In(2a) sites are indicated by blue arrows. Here the
hyperfine fields vanish at the In(1), Co and the In(2b) sites; the direction
of $\mathbf{H}_{0}$ is shown by the black arrow.

Figure 6: (Color online) Magnetic structure and hyperfine fields for dipolar
hyperfine couplings to the In(2) (case (2.4)). Same notation as in Fig. 6. The
Ce moments point along [001]. Here the hyperfine fields vanish at the In(1)
and Co sites, but not at the In(2) sites; the direction of $\mathbf{H}_{0}$ is
shown by the black arrow.

###  3.1 Hyperfine constants from the Knight shift

In CeCoIn5 measurements of the Knight shift in the normal state show that the
hyperfine couplings for the In(2a) site
($\mathbf{r}_{k}=(\pm\frac{a}{2},0,z_{0})$) for the applied field along (100),
(010) and (001) are 10.3 kOe/$\mu_{B}$, 0.0 kOe/$\mu_{B}$, and 32.4
kOe/$\mu_{B}$, respectively.22 For this site, Equations (8) and (9) yield:

$\displaystyle K_{a}$ $\displaystyle=$ $\displaystyle(2B_{iso}+B_{dip}(1-3\cos
2\theta_{z}))\chi_{a}$ (10) $\displaystyle K_{b}$ $\displaystyle=$
$\displaystyle(2B_{iso}-2B_{dip})\chi_{b}$ (11) $\displaystyle K_{c}$
$\displaystyle=$ $\displaystyle(2B_{iso}+B_{dip}(1+3\cos
2\theta_{z}))\chi_{c},$ (12)

where $\chi_{\alpha}$ is the susceptibility in the $\alpha$ direction. Using
the experimental numbers,22 we find $B_{iso}=B_{dip}=7.1$ kOe/$\mu_{B}$ and
$\theta_{z}=29^{\circ}$. The difference between this angle and that of the
crystal structure ($\theta_{z}=45^{\circ}$) probably is related to details of
the bonding of the In $4p$ orbitals, and will need further investigations.

###  3.2 Hyperfine fields in the B phase

Figure 7: (Color online) Real-space map (15 x 15 unit cells) of the hyperfine
field at the In(2a) (Eq. 6) (upper row) and the In(2b) (Eq. 7) (middle row) in
the $ab$ plane. Shown along the horizontal are the components of the hyperfine
field along the ($\hat{\mathbf{x}}$, $\hat{\mathbf{y}}$, $\hat{\mathbf{z}}$)
directions in (red, green blue) shading, for case (2.4) with
$\mathbf{Q}_{0}~{}||~{}[010]$ and $\mathbf{S}_{0}~{}||~{}[001]$. Black
corresponds to zero hyperfine field. The lower row shows the spin density
along the [001] direction (cyan), and the histogram of resonant frequencies
for the In(2a) (orange) and the In(2b) (yellow).

Figure 8: (Color online) The same notation as in Fig. 8, but for case (3.4)
with $\mathbf{Q}_{0}~{}||~{}[110]$ and $\mathbf{S}_{0}~{}||~{}[001]$. The
corresponding histogram of resonant frequencies is in better agreement with
experiment.

For cases (1.4), (2.4), (3.4) and (4.4) the hyperfine field at the In(2a) and
In(2b) sites varies spatially
$\sim\cos[(\mathbf{Q}_{0}+\mathbf{Q}_{i})\cdot\mathbf{r}]$ with modulus:

$\displaystyle\mathbf{H}_{hf}(2a)$ $\displaystyle=$ $\displaystyle
S_{0}B_{dip}[3\sin(2\theta_{z})\cos(\frac{\pi\delta_{x}}{2})\hat{\mathbf{a}}+(1-\cos(2\theta_{z})\sin(\frac{\pi\delta_{x}}{2})\hat{\mathbf{c}}]$
(13) $\displaystyle\mathbf{H}_{hf}(2b)$ $\displaystyle=$ $\displaystyle
S_{0}B_{dip}[3\sin(2\theta_{z})\cos(\frac{\pi\delta_{y}}{2})\hat{\mathbf{b}}+(1-\cos(2\theta_{z})\sin(\frac{\pi\delta_{y}}{2})\hat{\mathbf{c}}]$
(14)

In each case, the hyperfine field at the In(2a) site oscillates along the
modulation direction with a component along [100] and the resulting spectrum
is described by Figs. 3, 8 and 8. Using the values
$\delta_{x}=\delta_{y}\approx 0.12$ and $S_{0}\approx 0.15\mu_{B}$ as reported
by Kenzelmann et al.,15 we find $h_{hf}^{0}\approx 2.6$ kOe, which is about
twice the experimental value of 1.3 kOe. The difference may be related to
uncertainties in the hyperfine coupling itself 18 or changes in the magnetic
structure for the field along [100]. Note that the modulation and ordered
moment may differ for field oriented along [100], in which case Eq. (13) will
give a different value. In fact the NMR result is consistent with the locus of
points given by $\delta_{x}=\delta_{x0}\cos^{-1}(S_{0}^{0}/S_{0}-1)$, where
$\delta_{x0}=0.3$ and $S_{0}^{0}=0.14\mu_{B}$.

For cases (3.4) and (4.4), where the modulation is along [110] or
[1$\bar{1}$0], the hyperfine field at the In(1) site does not cancel but has a
component along the [001] direction. This field can give rise to a minor shift
and/or broadening of the In(1) line. The spectra (Fig. 3) clearly show that
the In(1) line shifts and is only slightly broadened. However, the shift and
broadening may come from the onset of spin shift suppression in the
superconducting state and the presence of superconducting vortices. Therefore,
we cannot distinguish the presence of a hyperfine field from the
antiferromagnetic structure at the In(1) site within experimental error. As
seen in Figs. 8 and 8, the calculated spectra for case (3.4) is closer to the
experimental one. We speculate that the true magnetic structure for the field
$\mathbf{H}_{0}$ along [100] is best described either by case (2.4), (3.4) or
(4.4) with $\mathbf{H}_{0}\perp\mathbf{S}_{0}$, which will minimize the free
energy of the antiferromagnet. In each case the hyperfine field at the In(2a)
and In(2b) have components perpendicular to $\mathbf{H}_{0}$, but these
components only give rise to small shifts of the resonant frequency that are
difficult to distinguish from the Knight shift. The crucial element is that
the hyperfine field at the In(2a) is along $\mathbf{H}_{0}$.

###  3.3 Nature of the B phase

The fact that the antiferromagnetism exists only in field and only within the
superconducting phase indicates a strong coupling between these order
parameters. Kenzelmann and coauthors15 analyzed the symmetry of the
superconducting state for such a coupling and concluded that the
superconducting order parameter, $\Delta_{\mathbf{Q}}$ acquires the finite
momentum $\mathbf{Q}$ of the antiferromagnetic order parameter,
${\mathbf{M_{Q}}}$. This corresponds to a modulation of the order parameter in
real space that presumably is out-of-phase with the antiferromagnetism. In
other words, the antiferromagnetic order is maximum at the nodes of
$\Delta_{\mathbf{Q}}$. Since then various microscopic models have been
proposed to explain the field induced antiferromagnetic order.23, 24
Curiously, this scenario is similar to that observed in the ferropnictide
SrFe2As2 under pressure.25 In this compound, a novel hybrid state of
coexisting superconductivity and antiferromagnetism emerges above 5 GPa. We
speculate that these two novel states may in fact be the same. Although
superconductivity and antiferromagnetism are known to coexist inhomogeneously
in a number of doped high Tc, heavy fermion, and ferropnictide systems, the
highly clean undoped CeCoIn5 and SrFe2As2 materials support the emergence of
this fragile but intrinsic thermodynamic phase of modulated antiferromagnetism
and superconductivity. In CeCoIn5, this phase is quickly destroyed by doping
and is replaced by a commensurate order at zero field for sufficiently high Cd
doping.26, 27, 28, 29 Clearly many questions about this new state of matter
remain unexplained, such as the driving mechanism(s), the origin of the
incommensurate wavevector, and the nature of the excitations.

##  4 Conclusions

In summary, we have shown that both NMR and neutron diffraction measurements
in the field-induced B phase of CeCoIn5 are consistent with magnetic
structures where $\mathbf{Q}_{i}\perp\mathbf{S}_{0}$ and
$\mathbf{H}_{0}\perp\mathbf{S}_{0}$. The incommensurate modulation
$\textbf{Q}_{i}$ lies possibly along either [010] or [110] direction for
magnetic field pointing along [100] in real space. Tetragonal equivalent
directions for $\textbf{Q}_{i}$, [100] and [1$\bar{1}$0], are also possible.
Based on our analysis of the NMR spectra, we speculate that the B phase of
CeCoIn5 represents an intrinsic phase of modulated superconductivity and
antiferromagnetism that can only emerge in a highly clean system. Further NMR
and neutron diffraction measurements are necessary for the same field
orientations to unravel the origin of the field-induced antiferromagnetic
structure.

###### Acknowledgements.

We would like to thank R. Movshovich, V. Mitrović, and M. Kenzelmann for
valuable discussions and sharing their results. Work at Los Alamos National
Laboratory was performed under the auspices of the US Department of Energy
under grant no. DE-AC52-06NA25396.

## References

  * 1 C. Petrovic et al.,  Europhys. Lett. 53, 354 (2001).
  * 2 V.A. Sidorov et al.,  Phys. Rev. Lett. 89, 157004 (2002).
  * 3 A. Bianchi et al.,  Phys. Rev. Lett. 91, 257001 (2003).
  * 4 A. Bianchi et al.,  Phys. Rev. Lett. 89, 137002 (2002).
  * 5 A. Bianchi et al.,  Phys. Rev. Lett. 91, 187004 (2003).
  * 6 H.A. Radovan et al.,  Nature 425, 51 (2003).
  * 7 K. Kakuyanagi et al.,  Phys. Rev. Lett. 94, 047602 (2005).
  * 8 K. Kumagai et al.,  Phys. Rev. Lett. 97, 227002 (2006).
  * 9 B.-L. Young et al.,  Phys. Rev. Lett. 98, 036402 (2007).
  * 10 A.B. Vorontsov, J.A. Sauls, and M.J. Graf,  Phys. Rev. B 72, 184501 (2005).
  * 11 A.B. Vorontsov and M.J. Graf,  Phys. Rev. B 74, 172504 (2006).
  * 12 V. F. Mitrović et al.,  Phys. Rev. Lett. 97, 117002 (2006).
  * 13 G. Koutroulakis et al.,  Phys. Rev. Lett. 101, 047004 (2008).
  * 14 Y. Matsuda and H. Shimahara,  J. Phys. Soc. Jpn. 76, 051005 (2007).
  * 15 M. Kenzelmann et al.,  Science 321, 1652 (2008).
  * 16 G. Koutroulakis et al., 2009 APS March Meeting, BAPS.2009.MAR.B41.6.
  * 17 V. Mitrović, (private communications).
  * 18 N.J. Curro, B.-L. Young, J. Schmalian, and D. Pines,  Phys. Rev. B 70, 235117 (2004).
  * 19 M.-A. Vachon et al.,  J. Phys.: Condens. Matter 20, 29225 (2008).
  * 20 N.J. Curro,  New. J. Phys. 8, 173 (2006).
  * 21 C.P. Slichter, Principles of Nuclear Magnetic Resonance, 3rd edn. Springer-Verlag (1992).
  * 22 N.J. Curro et al.,  Phys. Rev. B 64, 180514 (2001).
  * 23 Y. Yanase and M. Sigrist,  JPCS 150, 052287 (2009).
  * 24 H.-Y. Kee and D. Podolsky,  Europhysics Letters, 86, 57005 (2009); K. Miyake, J. Phys. Soc. Jpn. 77, 123703 (2008); Y. Yanase J. Phys. Soc. Jpn. 77, 063705 (2008); A Aperis et al., J. Superconductivity and Novel Magnetism 22 115 (2009)
  * 25 K. Kitagawa et al.,  arXiv:0906.4740 (2009).
  * 26 Y. Tokiwa et al.,  Phys. Rev. Lett. 101, 037001 (2008).
  * 27 R.R. Urbano et al.,  Phys. Rev. Lett. 99, 146402 (2007).
  * 28 R.R. Urbano, B.-L. Young, N.J. Curro,  Physica B 403, 1056 (2008).
  * 29 C.F. Miclea et al.,  Phys. Rev. Lett. 96, 117001 (2006).