text:
# Coexistence of antiferromagnetism and superconductivity in the Anderson lattice

## 1 Introduction

It is a common feature of several strongly correlated electronic systems that the low temperature ordered phases compete with each other. In particular, competition between antiferromagnetic and superconducting orders is an important characteristic of heavy-fermion systems, which is also shared by high-$`T_c`$ materials  and low-dimensional systems . The closeness of the superconducting phase to the antiferromagnetic phase in heavy fermion compounds has lead to the conjecture that the attractive interaction leading to superconductivity is actually mediated by a magnetic excitation , instead of the traditional phonon mechanism. Also, the interplay between antiferromagnetism and superconductivity has been proposed to be described by a $`SO(5)`$ symmetry breaking . The complexity of these systems arises from the interplay between Kondo screening of the local moments, the antiferromagnetic (RKKY) interaction between the local moments and superconducting correlations between the heavy quasiparticles.

Heavy fermion systems that exhibit both superconductivity and antiferromagnetism exhibit ratios between the Néel temperature $`T_N`$ and the superconducting critical temperature $`T_c`$ that can vary substantially (of the order of $`T_N/T_c1100`$), with coexistence of both types of order below $`T_c`$. The coexistence of both types of order can be tuned by external parameters such as externally applied pressure or chemical pressure (involving changes in the stoichiometry). Examples of heavy-fermion materials which exhibit antiferromagnetic and superconducting order at low temperature are $`URu_2Si_2`$ and $`U_{0.97}Th_{0.03}Be_{13}`$. It has recently been found that $`UPd_2Al_3`$ ($`T_N=14.3`$ K and $`T_c=2`$ K) and $`UNi_2Al_3`$ ($`T_N=4.5`$ K and $`T_c=1.2`$ K) show coexistence of superconductivity and local moment antiferromagnetism.  However, in the $`Ce`$-based heavy-fermions magnetism tipically competes with superconductivity. In the prototype heavy-fermion system CexCu2Si2 both coexistence and competition between $`d`$wave superconductivity and magnetic order has been clearly observed in a small range of $`x`$ values around $`x0.99`$ pressure. This system exhibits a magnetic “A phase” at low temperature whose detailed nature is not yet known. Increasing pressure reduces the critical temperature $`T_A`$ of the A phase. Recent studies of $`\mathrm{CeCu}_2\mathrm{Si}_2`$ samples near stoichiometric composition have shown that a $`d`$-wave superconducting phase expels the magnetic “A phase” when $`T_A`$ approaches $`T_c`$ under increasing pressure.

In general terms, the local moments due to the $`f`$-electrons are progressively quenched as the temperature lowers. In dilute systems the picture is well understood as due to the Kondo screening by the conduction electrons. In dense systems however the picture is more involved. At low temperatures the local moments are not completely quenched. In $`U`$-based materials such as $`URu_2Si_2`$, $`UPt_3`$ or $`UPd_3`$, for instance, the remaining moments are quite small of the order of $`0.010.03\mu _B`$ but for other systems such as $`UPd_2Al_3`$ the local moment is quite large of the order of $`0.85\mu _B`$. The low temperature magnetic behaviour of $`Ce`$-based compounds such as $`CeCu_2Si_2`$ and $`CeCu_{6x}Au_x`$ has been interpreted as due to the vicinity to a quantum critical point  where the Néel temperature tends to zero. Two pictures arise however : in the first one the Kondo temperature is high (the moments are quenched at a finite temperature) and when the system approaches the quantum critical point there are no free moments (assuming that quenching is complete). Then the system has to order due to a Fermi surface instability of the spin density wave type. Another possible situation is one in which the moments are not completely quenched down to $`T=0`$ and are free to orient themselves leading to magnetism. In the case of $`CeCu_{6x}Au_x`$ evidence has been recently found that the second picture seems to hold  but a small but finite Kondo temperature has been quoted for this material (see the scond reference in ). On the other hand the high value of the Kondo temperature for the $`CeCu_2Si_2`$ compound  indicates possibly that the first scenario should hold. Furthermore, recent experiments  with $`CeRh_{1x}Ir_xIn_5`$ also reveal unusual coexistence of magnetism and superconductivity. It appears that in this system the $`f`$-electrons are more band-like than localized. On the other hand in $`UPd_2Al_3`$ it is the dual character of the $`5f`$ electrons that leads to the high value of the local moment in coexistence with the itinerant electrons and with the superconductivity .

On the theoretical side it is believed that the Anderson model and its extension to the lattice captures the basic physics of the heavy-fermions . In this model the conduction electrons, $`c`$, (usually regarded as free) hybridize with local states, $`f`$, where the electrons are strongly interacting in such a way that the Coulomb repulsion, $`U`$, between two $`f`$-electrons is the largest energy scale in the problem. Frequently the limit $`U\mathrm{}`$ is taken implying that double occupancy is forbidden. The limit $`U=\mathrm{}`$ has been studied using the slave boson technique. In particular, it has been shown that superconducting instabilities arise in the $`p`$ and $`d`$-wave channels because of the effective (RKKY) interaction between the $`f`$electrons. Recently, the magnetic and superconducting instabilities of the normal phase were studied in the random phase approximation (RPA)  by taking into account slave boson fluctuations above the condensate. However, the competition/coexistence between both types of ordering was not considered.

In this work we consider the $`U=\mathrm{}`$ Anderson lattice model in the slave-boson approach. Because our aim is to study the interplay between magnetism and superconductivity, we explicitly introduce antiferromagnetic and superconducting couplings phenomenologically. The coupling constants are taken as independent even though they are related if the superconducting mechanism is mediated via the RKKY interaction. Using a mean-field approach we study directly the ordered phases and determine regimes of coexistence or competition between the two ordered phases depending on the parameters of the model.

In this work we will take the conduction electrons to be non-interacting but we should also mention that attempts to include interactions between the conduction electrons have been carried out . The inclusion of this more realistic interaction has been found to be required in some systems to attain a better understanding of the experimental results. We focus our attention on the coexistence of superconducting correlations and magnetic ordering in heavy fermions which, to our knowledge, has not yet been studied theoretically despite the considerable recent experimental effort devoted to this subject.

## 2 Model Hamiltonian and quasiparticle spectrum

The microscopic description of superconductivity and magnetic order in the Anderson lattice model is a still unsolved problem. In the folllowing we shall consider an effective Hamiltonian which originates from the $`U=\mathrm{}`$ Anderson model with two additional phenomenological terms: one where superconducting correlations are explicitly assumed between the local $`f`$-electrons (since it is believed that pairing occurs between heavy quasi-particles which, therefore, have essentially $`f`$ character) and another term where local spins are coupled antiferromagnetically. The infinite Coulomb repulsion between the $`f`$ electrons is described within Coleman’s slave-boson technique with a condensation amplitude $`\sqrt{z}`$. The effective Hamiltonian is therefore:

$`H^{MF}`$ $`=`$ $`{\displaystyle \underset{k\sigma }{}}(ϵ_k\mu )c_{k\sigma }^{}c_{k\sigma }+{\displaystyle \underset{k\sigma }{}}(ϵ_f\mu )f_{k\sigma }^{}f_{k\sigma }`$ (1)
$`+`$ $`\sqrt{z}V{\displaystyle \underset{k\sigma }{}}\left(f_{k\sigma }^{}c_{k\sigma }+c_{k\sigma }^{}f_{k\sigma }\right)`$
$`+`$ $`z{\displaystyle \underset{k}{}}\left(\mathrm{\Delta }_f\eta _kf_k^{}f_{k,}^{}+\mathrm{\Delta }_f^{}\eta _kf_{k,}f_k\right)`$
$``$ $`2m{\displaystyle \underset{k\sigma }{}}\sigma \left(f_{k\sigma }^{}f_{kQ,\sigma }+f_{k\sigma }^{}f_{k+Q,\sigma }\right)`$
$`+`$ $`(ϵ_fϵ_0)(z1)N_s+{\displaystyle \frac{2N_sm^2}{J_m}}{\displaystyle \frac{N_s|\mathrm{\Delta }_f|^2}{J_f}}`$
The $`c`$ and $`f`$ operators refer to conduction and localized electrons and obey the usual anticommutation relations. For simplicity, the hybridization potential $`V`$ is assumed to be momentum independent, $`ϵ_0`$ and $`ϵ_f`$ denote the bare and renormalized $`f`$-level energies, and $`N_s`$ denotes the number of lattice sites.

Although it is known that conduction electrons provide an effective RKKY interaction between $`f`$ electrons, we remark that the superconducting and magnetic order parameters in (1) cannot be simultaneously derived from a Hubbard-Stratonovich decoupling of a single RKKY term of the form $`\stackrel{~}{J}_{ij}\stackrel{}{S}_i\stackrel{}{S}_j`$, as discussed, for instance, in ref.  This is why we have phenomenologically introduced those terms.

In writing down the pairing term in (1) we have in mind that in the slave boson formulation a slave-boson operator $`b_i`$ is associated with every $`f_i`$ operator to prevent double occupancy. Condensation of the slave bosons is described by the replacement $`b_i<b_i>=<b_i^{}>=\sqrt{z}`$, hence the factor $`z`$ in the superconducting term of (1). The superconducting order parameter is given by $`\mathrm{\Delta }_f=\frac{zJ_f}{N_s}_k\eta _k<f_{k,}f_k>`$, where $`\eta _k`$ denotes any of the possible pairing symmetries $`\eta _k^{(s)}=\mathrm{cos}k_x+\mathrm{cos}k_y`$, $`\eta _k^{(p,i)}=\sqrt{2}\mathrm{sin}k_i`$ and $`\eta _k^{(d)}=\mathrm{cos}k_x\mathrm{cos}k_y`$ for $`s`$, $`p`$ and $`d`$ waves, respectively. Here we consider two dimensions for simplicitly of the calculations. We take a square lattice even though several heavy-fermions have a complicated lattice structure since we want to capture the main features.

The magnetic order parameter is given by the mean-field equation $`m_i=J_m<f_i^{}f_if_i^{}f_i>`$ and $`m_i=2m\mathrm{cos}(\stackrel{}{r}_i\stackrel{}{Q})`$ where $`\stackrel{}{Q}=(\pi ,\pi )`$ is the antiferromagnetic ordering vector. We consider only commensurate antiferromagnetic correlations since it is particularly relevant to the systems referred.

In the calculations we assume a simple dispersion for the conduction electrons, of the form $`ϵ_k=2t_{i=x,y}\mathrm{cos}k_i`$. For a given particle density the chemical potential must be computed from the condition $`n_c+n_f=n`$. The mean-field equations obtained from minimization of the free energy of Hamiltonian (1) are solved numerically. The solution gives the interplay between the boson condensation, the magnetization, and superconducting pairing as function of band-filling and of the various model parameters. Throughout this work we will consider only the case of $`d`$-wave pairing. The other two pairings give qualitatively similar results in most regimes as stressed before for the case of no magnetism and only superconducting order . In different regimes different pairing symmetries become the most stable one but the overall behavior is similar. However we will get back to this point later. We will consider temperatures such that the slave bosons are condensed and therefore we quench the $`f`$-electron density fluctuations.

### 2.1 Magnetic non-superconducting phase

We begin by discussing the magnetic but non-superconducting phase of the model. A detailed study of the superconducting non-magnetic phases has been presented elsewhere.

In order to diagonalize (1) for $`J_f=0`$ we introduce a quasi-particle operator which is a linear combination of the operators forming the basis $`(c_{k\sigma }^{};c_{k+Q,\sigma }^{};f_{k\sigma }^{};f_{k+Q,\sigma }^{})`$. We write the Hamiltonian in this basis. The mean-field equations are obtained varying the effective Hamiltonian in the usual way. Minimization of the free energy of (1) with respect to $`z`$ gives:

$$ϵ_fϵ_0=\frac{V}{2\sqrt{z}}\frac{1}{N_s}\underset{k\sigma }{}\left(<f_{k\sigma }^{}c_{k\sigma }>+<c_{k\sigma }^{}f_{k\sigma }>\right).$$
(2)
The condition $`n_f+z=1`$ reads:

$$z=1\frac{1}{N_s}\underset{k\sigma }{}<f_{k\sigma }^{}f_{k\sigma }>,$$
(3)
and the chemical potential is related to the total particle density as:

$$n=\frac{1}{N_s}\underset{k\sigma }{}<c_{k\sigma }^{}c_{k\sigma }>+\frac{1}{N_s}\underset{k\sigma }{}<f_{k\sigma }^{}f_{k\sigma }>.$$

The magnetization is given by:

$$m=\frac{J_m}{2N_s}\underset{k\sigma }{}\sigma <f_{k\sigma }^{}f_{kQ,\sigma }+f_{k\sigma }^{}f_{k+Q,\sigma }>.$$
(4)
For the sake of comparison, we show in Figure 1 the quasiparticle bands and the density of states for a situation where the magnetic order parameter is zero and nonzero, respectively. We see that in addition to the hybridization gap, there is an additional gap due to the magnetic order. The magnetic order with momentum $`\stackrel{}{Q}`$ reduces the size of the Brillouin zone to a half of that of a non-magnetic system.

### 2.2 Coexistence of antiferromagnetism and superconductivity

The Hamiltonian matrix $`H_k`$ can be written in the basis

$$(c_k,c_{k+Q,},f_{k,},f_{k+Q,},c_{k,}^{},c_{kQ,}^{},f_{k,}^{},f_{kQ,}^{})$$

as:

$$H_k=\left(\begin{array}{cc}A_+& D\\ D^{}& A_{}\end{array}\right)$$

where the matrices $`A_\pm `$ and $`D`$ are given by

$$A_\pm =\left(\begin{array}{cccc}\pm ϵ_k& 0& \pm \sqrt{z}V& 0\\ 0& \pm ϵ_{k+Q}& 0& \pm \sqrt{z}V\\ \pm \sqrt{z}V& 0& \pm ϵ_f& 2m\\ 0& \pm \sqrt{z}V& 2m& \pm ϵ_f\end{array}\right)$$

and

$$D=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& z\mathrm{\Delta }_f\eta _k& 0\\ 0& 0& 0& z\mathrm{\Delta }_f\eta _{k+Q}\end{array}\right)$$

The energies are measured with respect to the chemical potential. The eigenvectors of the matrix $`H_k`$ are the Bogolubov operators expressed in the same basis.

Variation of the free energy with respect to $`z`$ shows that (2) must now be replaced by

$`ϵ_fϵ_0`$ $`=`$ $`{\displaystyle \frac{V}{2\sqrt{z}N_s}}{\displaystyle \underset{\stackrel{}{k},\sigma }{}}\left(<f_{\stackrel{}{k},\sigma }^{}c_{\stackrel{}{k},\sigma }>+<c_{\stackrel{}{k},\sigma }^{}f_{\stackrel{}{k},\sigma }>\right)`$ (5)
$``$ $`{\displaystyle \frac{2N_s}{zJ_f}}\mathrm{\Delta }_f^2`$
while the other mean-field equations remain unaltered.

In Figure 2 we show a typical quasiparticle band structure for the case where there is coexistence of magnetism and superconductivity. The quasiparticle bands are symmetric around the chemical potential due to the particle-hole structure of the Bogoliubov operators. Due to the superconducting order the spectrum is gapless.

## 3 Phase diagrams

The phase diagram of the model is quite rich due to the various correlations and order parameters considered. In this work we focus our attention on regimes where magnetism and superconductivity coexist or compete.

In Figure 3 we plot the magnetization and the $`d`$-wave superconducting order parameter as functions of temperature for a typical case. The critical temperature for the antiferromagnetic order parameter $`m`$, $`T_m`$, is larger than the critical temperature, $`T_c`$, for $`d`$-wave superconductivity. At low temperatures the two phases coexist. While $`\mathrm{\Delta }_f0`$ the magnetization $`m`$ is increasing and then $`m`$ decreases to zero at $`T_m`$ in the usual mean-field like manner.

The left panel of Figure 4 shows the order parameters $`m`$ and $`\mathrm{\Delta }_f`$ as functions of band-filling at a fixed low temperature $`T=0.001`$. For low to intermediate band-fillings superconductivity exists. As the band-filling increases, the superconductivity is less favorable while the magnetic order parameter appears, as expected, since the magnetic order is more favored if the local electron density is higher. The value of $`\mathrm{\Delta }_f`$ vanishes at low densities because the $`f`$-level occupancy also becomes small in that limit ($`z1`$, $`n_f0`$) and Cooper pairing occurs only between the $`f`$-electrons in the model considered. In the high density limit, $`\mathrm{\Delta }_f`$ also tends to zero because the $`f`$-level occupancy is higher, approaching 1, and freezing of the charge fluctuations occurs due to the infinite on-site repulsion. Furthermore, a comparison with earlier results, where magnetic order was not considered, shows that magnetic order lowers the maximum value of the band-filling for which $`\mathrm{\Delta }_f0`$, indicating that the two effects compete with each other. At smaller band-fillings incommensurate antiferromagnetic order also stabilizes. Indeed it extends to lower band-fillings as compared to the commensurate case. However the superconducting order expels the incommensurate case as well and including both types of ordering the phase diagram is qualitatively similar for commensurate or incommensurate order if superconductivity is also allowed.

The behavior of the order parameters against magnetic coupling, $`J_m`$, at low temperature, is shown in the right panel of Figure 4. Increasing $`J_m`$ leads to a crossover from a region where $`m=0`$ to a regime where $`m`$ grows with $`J_m`$ while $`\mathrm{\Delta }_f`$ follows the opposite trend. Keeping $`J_m`$ fixed and decreasing $`J_f`$ leads to the opposite result where the superconductivity disappears in favor of magnetism.

In Figure 5 we plot the two critical temperatures as functions of the band-filling and of the hybridization. The behavior of the critical temperatures against density follows the same trend as in Figure 4: for low $`n`$ the superconducting temperature is higher while for higher densities $`T_m>T_c`$, indicating which phase is favored as one lowers the temperature from the disordered high temperature phase. Lowering $`ϵ_0`$ or $`V`$ has the tendency to increase $`f`$-occupancy favoring magnetic order over superconducting order. At the point where the two temperatures cross the magnetic temperature falls abruptly to zero if the Cooper pairing symmetry is $`d`$-wave. This does not happen if the pairing symmetry is either extended $`s`$-wave or $`p`$-wave. In these cases there is no abrupt expulsion of the magnetization when the magnetic critical temperature becomes lower than the superconducting critical temperature and the phase where AF and SC coexist extends for smaller band-filling values and for larger values of the hybridization.

In the light of the experimental results, it is clearly interesting to compare qualitatively these results with the experimental phase diagrams, where either external pressure or chemical pressure are varied. Increasing pressure is expected to increase both hopping and hybridization amplitudes, while probably keeping $`V/t`$ and other parameters approximately constant, to a first approximation. In Figure 6 we plot the mean field temperatures as functions of $`V`$ for fixed ratio $`V/t`$. As pressure increases the magnetic critical temperature decreases and $`T_c`$ increases. Once again at the point where the two temperatures cross the magnetic temperature falls abruptly to zero if the Cooper pairing symmetry is $`d`$-wave. This result, showing expulsion of a spin-density-wave by $`d`$-wave superconductivity but not in the case of the other symmetries, is very similar to that of a previous study where the problem of local moment formation in a superconducting phase was addressed: increasing pressure causes expulsion of local moments by $`d`$-wave superconductivity. Therefore, such a phenomenon occurs either in a scenario of (unscreened or partially screened) ordered local moments or in a spin-density-wave scenario, where magnetic order appears as a Fermi surface instability.

The abrupt decrease of $`m`$ as $`T_c`$ crosses $`T_m`$ signals a quantum phase transition that can be tuned using the external pressure as a parameter: as pressure is reduced (at zero $`T`$) the groundstate of the system changes abruptly from nonmagnetic but superconducting to magnetic and superconducting at a critical value of $`V_c`$. The transition appears to be first order. We note that the same behavior is found for the compound $`CeCu_2Si_2`$ .

Experimentally, most systems have quite large effective masses. Large mass enhancements are obtained when $`z`$ is small ($`n_f1`$). This is indeed observed in regimes where the AF order parameter is large but in these regimes superconductivity is absent. In Figure 7 we show results for the quasiparticle bands and densities of states where a strong mass enhancement is seen for a large magnetization in opposition to a situation where the magnetization is zero. In the calculations above, we have not been able to find regimes where there is a very large (say larger than $`200`$) mass enhancement. In the model considered, the appearance of superconductivity is restricted to the mixed valent regime. In this situation we do not expect very large densities of states at the chemical potential. In a more realistic approach $`U`$ is large but finite and larger band-fillings are allowed ($`0n_f2`$). We expect therefore that in this regime we might find large effective masses together with superconducting and magnetic order. The restriction to the mixed valent regime is a consequence of the slave boson approach . In this method the chemical potential is always pinned to the lowest quasiparticle band, the density is always smaller than one, and the superconductivity only appears in the intermediate valence regime with moderate effective masses . Nevertheless, our results reproduce qualitatively well the experimental phase diagrams and the competition between the $`d`$-wave superconducting and magnetic phases, such as that observed, for instance, in $`CeCu_2Si_2`$.

Model (1) is an effective Hamiltonian for the interacting $`c`$ and $`f`$ particles. Since the quasiparticles are heavy close to the top of the lowest band, as evidenced by the high values of the specific heat jump, we have modelled the superconducting interaction as taking place between the local $`f`$ electrons, as usual. One may wonder however the effect of adding a pairing term, with coupling constant $`J_c`$, between the $`c`$-electrons, since these are present at low energies close to the chemical potential. Considering a pairing term in the Hamiltonian in the standard way, with amplitude $`\mathrm{\Delta }_c`$, we can solve the mean-field equations as before. The effect of this added term is shown in Figure 8. In general, the hybridization exchanges electrons between the $`c`$ and the $`f`$ levels. Due to the restriction on the level occupancy due to the infinite $`U`$ repulsion, for a fixed $`V`$ as the density increases the number of f-electrons decreases and the superconductivity is destroyed, as discussed above. However, if the $`c`$ electrons pair then the order parameter $`\mathrm{\Delta }_c`$ increases. As a consequence in Figure 8 we see that as the density increases the ordering in the c-electrons increases. Actually, if we consider both types of pairing then the pairing in the f electrons extends to higher densities. So in that sense adding the pairing between the c electrons favors superconductivity at higher densities. However, we have found that the adddition of this type of pairing inhibits the magnetism in the f-electrons probably due to the greater stability in the pairing channel.

We remark however, that had we written a BCS pairing term only among the $`c`$-electrons in the Hamiltonian (1), the resulting phase diagram would be quite different: the superconducting temperature would be maximum at zero hybridization and rapidly decrease with increasing $`V`$. This is easily understood because increasing $`V`$ makes $`c`$-electrons heavier hence reducing their pairing amplitude. Clearly this is the opposite trend to that found experimentally in the phase diagrams as a function of pressure where, at low hybridizations, superconductivity is absent giving place to the magnetic order. Even if we do not include the possibility of magnetic order then it was shown before  that for small values of the hybridization the superconducting critical temperature is an increasing function of pressure. On the other hand, by having chosen a pairing interaction between $`f`$-electrons, we obtained increasing $`T_c`$ since the $`f`$-electrons become more mobile upon increasing $`V`$. Therefore it is justified to consider only a pairing term between the heavy f-electrons.

## 4 Optical and dynamic conductivity

The study of the particle-hole excitations of a system can be probed by studying the finite frequency conductivity. Studies of the Anderson lattice have been carried out previously  in the disordered phases and experimental results for the heavy fermions have been reviewed in Ref . In traditional superconductors the optical conductivity clearly shows a threshold at twice the gap energy. Since the d.c. conductivity of a perfect lattice is infinite, most studies are carried out taking into account scattering off impurities  considering moderate scattering up to the dirty limit . In this work we will study the optical conductivity in the clean limit and we shall study its behavior in the ordered phases as well.

The real part of the dynamic conductivity is given by

$$\sigma _{\alpha \beta }(\stackrel{}{q},\omega )=\frac{1e^{\beta \omega }}{2\omega V}_{\mathrm{}}^{\mathrm{}}𝑑te^{i\omega t}<j_\alpha ^{}(\stackrel{}{q},t)j_\beta (\stackrel{}{q},0)>$$
(6)
Writing the trace, inserting a decomposition of the identity in terms of the exact many-body energy states and integrating over time we obtain

$`\sigma _{\alpha \beta }(\stackrel{}{q},\omega )`$ $`=`$ $`{\displaystyle \frac{1e^{\beta \omega }}{2\omega V}}2\pi {\displaystyle \underset{n,m}{}}{\displaystyle \frac{e^{\beta E_n}}{Z}}\delta (\omega +E_nE_m)`$ (7)
$`<n|j_\alpha ^{}(\stackrel{}{q})|m><m|j_\beta (\stackrel{}{q})|n>`$
Consider now $`\alpha =\beta `$ and $`\stackrel{}{q}=0`$. The $`\alpha =x,y`$ component of the current operator can be written as

$$j_\alpha =it\underset{\stackrel{}{k},\sigma }{}\left(\stackrel{}{\chi }_\stackrel{}{k}\right)_\alpha c_{\stackrel{}{k}\sigma }^{}c_{\stackrel{}{k}\sigma }$$
(8)
where $`\stackrel{}{\chi }_\stackrel{}{k}=_\stackrel{}{\delta }\stackrel{}{\delta }e^{i\stackrel{}{k}\stackrel{}{\delta }}`$ where $`\sigma `$ is the spin component. Writing the electronic operators in terms of the Bogolubov operators it is straightforward to calculate the matrix elements of the current operator.

We can also calculate the finite-momentum and finite-frequency conductivity $`\sigma _{xx}(\stackrel{}{q},\omega )`$. Starting from eq. (8) and taking similar steps to those followed for the calculation of the optical conductivity it is easy to calculate the dynamical conductivity.

In Figures 9 and 10 we show the optical conductivity in the clean limit for two typical cases. Figure 9 refers to a regime where there is superconductivity but no magnetization at low temperature. Figure 10 refers to a regime where there is coexistence of superconductivity and antiferromagnetism at low temperature.

The experimental study of the optical conductivity, over a large range of frequencies, of some heavy fermion systems has only recently become available.  In reference  it was found that effective mass $`m^{}`$ of the quasiparticles increased about 50 times when the temperature decreases below the Néel temperature. Furthermore, the authors found for $`T<T_N`$ that a two peak structure developed at finite frequencies. At temperatures higher than $`T_N`$ the $`\sigma (\omega )`$ presents a single peak, separated from the Drude weight by a finite gap (if disorder is included the Drude peak broadens to finite energies). When the temperature drops below $`T_N`$ a second peak, respecting to the magnetic gap, shows up in $`\sigma (\omega )`$ at lower, but finite, energies. The features reproduced by our calculation are in qualitative agreement with the data presented in Ref. . Also they are in qualitative agreement with the results of . Note however that our model is not appropriate to describe this material since one needs to take into account the dual nature of the f-electrons, which is not included in our model. This shows that the results obtained are qualitative general trends that are captured by our simplified model.

In order to preserve momentum the optical conductivity probes the transitions between different bands coupled by the current operator. If we allow that momentum is interchanged then the dynamic conductivity also probes excitations along the same band in addition to across bands. In Figs. 11, 12, 13 we present the dynamic conductivity for different points in the Brillouin zone in three typical situations. In Fig. 11 the set of parameters and temperature are such that the system is superconducting, in Fig. 12 the system is both superconducting and magnetic and in Fig. 13 the system is only magnetic. The gap structure evident in the optical conductivity ($`\stackrel{}{q}=0`$) is now absent because low energy excitations are in general allowed since there is coupling between states in the same band. In the superconducting case (but not magnetic) the dynamic conductivity extends to zero frequency but vanishes in this limit due to the energy dispersion of the $`d`$-wave symmetry. If there is coexistence of superconductivity and magnetism then the same structure appears with the double-peak feature, characteristic of the magnetic phase, enhanced, as observed in the optical conductivity. If the superconducting order parameter is zero then the dynamic conductivity is finite at zero energy.

## 5 Summary

The interplay between magnetic correlations, the Kondo effect and superconducting correlations in heavy fermion systems is a difficult problem to solve. While previous studies on the $`U=\mathrm{}`$ Anderson lattice have identified the instabilities towards magnetic or superconducting order by taking into account the slave-boson fluctuations, the description of the ordered phases has not previously been carried out. In this work we have phenomenologically studied the interplay between the superconducting and antiferromagnetic ordered phases by studying their dependence on the model parameters. We have found that Cooper pairing symmetry plays an important role in determining the regimes of either coexistence or competition between superconductivity and antiferromagnetism: indeed, we have found that $`d`$-wave superconductivity coexists with antiferromagnetic order but the former expels the latter abruptly as the hybridization between the $`f`$-level and the conduction electrons increases. We also found that the $`f`$-level occupancy affects the two types of order in different ways: higher occupancy favors antiferromagnetism while a lower occupancy favors superconductivity. Hence, the quasiparticle mass enhancement relative to the bare electron was found not to exceed a few tens in the superconducting phase. The optical conductivity has also been computed, reflecting the band structure of the ordered phases.

The author wants to thank Miguel Araújo and Nuno Peres for many discussions on this subject.

## References


format: text

text:
# Supercooled Liquids under Shear: Computational Approach

## Abstract

We have performed molecular dynamics simulations for a model two-dimensional soft-core mixture in a supercooled state. The mixture exhibits a slow structural relaxation in a quiescent state, however, the relaxation is much enhanced in sheared states. There observed surprisingly small anisotropy both in the coherent and incoherent density correlation functions even under extremely strong shear which is $`10^3`$ times faster than the structural relaxation rate. The present simulation results agree well with predictions of the recently developed mode-coupling theory in shear.

As liquids are cooled toward the glass transition, the dynamics is drastically slowed down, while only small changes can be detected in the static properties. One of the main targets in theoretical investigations on the glass transition is to identify the mechanism of the drastic slowing-down. Beside this fundamental problem, a striking example occurs when one brings glassy materials away from equilibrium, for instance, by changing temperature rapidly or applying shear flow to them. There appears a variety of unique phenomena such as aging or shear thinning. Although these phenomena are not only conceptually but also practically important, physical properties of glassy materials in nonequilibrium conditions has not yet been understood well. A couple of years ago, we performed extensive molecular dynamics (MD) simulations in two dimensions (2D) and three dimensions (3D) on binary soft core mixtures with and without shear flow. It was found that the dynamical properties of the mixtures under shear can be fairly mapped onto those at quiescent states at higher temperatures yamamoto1998 ; yamamoto2002 . We found also the surprising isotropy in the tagged particle motions even under extremely high shear, which may justify the simple mapping idea. In the present study, we calculate intermediate scattering functions by using the method proposed by Onuki onuki1979 ; onuki1997  to investigate microscopic dynamics of glassy materials in shear flow. Simulations have been done in 2D to compare the present computational results directry with the theory developed recently in 2D miyazaki2002 .

Our model system is composed of two different particle spices $`1`$ and $`2`$, which interact via the soft-core potential

$$v_{\alpha \beta }(r)=ϵ(\sigma _{\alpha \beta }/r)^{12},\sigma _{\alpha \beta }=(\sigma _\alpha +\sigma _\beta )/2,$$
(1)
where $`r`$ is the distance between two particles, and $`\alpha ,\beta 1,2`$. We take the mass ratio $`m_2/m_1=2`$, the size ratio $`\sigma _2/\sigma _1=1.4`$, and the numbers of particles $`N_1=N_2=5000`$. Simulations were performed in the absence and presence of shear flow with being particle density and temperature fixed at $`n=(N_1+N_2)/V=0.8/\sigma _1^2`$ and $`k_BT=0.526ϵ`$, respectively. In the sheared case, we kept the temperature at a constant using the Gaussian constraint thermostat to eliminate the viscous heating effect. Our method of applying shear is as follows: The system was at rest for $`t<0`$ for a very long equilibration time and was then sheared for $`t>0`$. Here we added the average velocity $`\dot{\gamma }y`$ to the velocities of all the particles in the $`x`$ direction at $`t=0`$ and afterwards maintained the shear flow by using the Lee-Edwards boundary condition. Simulation data have been taken and accumulated in steady states which can be realized after transient states.

Figure 1 (a) shows the geometry of shear flow in the present simulation. As shown in Figure 1 (b), shear flow with the rate $`\dot{\gamma }`$ advect a positional vector r as

$$𝐫(t)=𝐫+\dot{\gamma }tr_y𝐞_x$$
(2)
in the time duration $`t`$, where $`𝐞_\alpha `$ is a unit vector in $`\alpha x,y`$ axis. A similar advection can be defined in Fourier space for a wave vector k as

$$𝐤(t)=𝐤\dot{\gamma }tk_x𝐞_y.$$
(3)
The above definition enable us to calculate Fourier component $`𝐤`$ of the time correlation function

$$C(𝐤,t)A_𝐤(0)A_{𝐤(t)}(t)$$
(4)
in shear flow. We calculated the incoherent and coherent scattering functions by using the definitions

$$Fs(𝐤,t)=\frac{1}{N}\underset{i=1}{\overset{N}{}}e^{[i\{𝐤(t)𝐫_i(t)𝐤𝐫_i(0)\}]}$$
(5)
and

$$F(𝐤,t)=\frac{1}{N}\underset{i=1}{\overset{N}{}}e^{[i𝐤(t)𝐫_i(t)]}\underset{j=1}{\overset{N}{}}e^{[i𝐤𝐫_j(0)]},$$
(6)
respectively. To examine anisotropy in the dynamics, the wave vector is taken in four different directions $`𝐤_{10}`$, $`𝐤_{01}`$, $`𝐤_{11}`$, and $`𝐤_{11}`$, where

$$𝐤_{\alpha \beta }=\frac{k_0}{\sqrt{\alpha ^2+\beta ^2}}(\alpha 𝐞_x+\beta 𝐞_y),$$
(7)
$`k_0=2\pi /\sigma _1`$, and $`\alpha ,\beta 0,1`$ as shown in Figure 1 (d).

Figure 2 (a) and (b) show $`Fs_1(𝐤,t)`$ and $`F_{11}(𝐤,t)`$, respectively, at $`T=0.526`$ with and without shear flow. Here the subscript $`1`$ denotes the smaller particle component. The so-called $`\alpha `$ relaxation time $`\tau _\alpha `$ of the present mixture is defined by

$$F_{11}(k_0,\tau _\alpha )Fs_1(k_0,\tau _\alpha )=e^1$$
(8)
in the quiescent state. The followings have been found; i) $`Fs_1(𝐤,t)`$ and $`F_{11}(𝐤,t)`$ behave quite similarly both in quiescent and sheared states. ii) Shear accelerates drastically the microscopic structural relaxation in the supercooled state. The structural relaxation time $`\tau _\alpha `$ decreases strongly with increasing shear rate as $`\tau _\alpha \dot{\gamma }^\nu `$ with $`\nu 1`$. iii) The acceleration in the dynamics due to shear occurs almost isotropically. There observed surprisingly small anisotropy in the correlation functions even in extremely strong shear $`\dot{\gamma }\tau _\alpha 10^3`$. This simplicity in the dynamics is quite different from behaviors of other complex fluids such as critical fluids or polymers in shear. Finally we note that the recent mode-coupling theory in shear flow miyazaki2002  predicts almost the same behaviors.


text:
# Nuclear pasta structure in hot neutron stars

## References

Hashimoto, M., Seki, H., & Yamada, M. 1984, Prog. Theor. Phys., 71, 320

Lassaut, M., Flocard, H., Bonche, P., Heenen, P. H., & Suraud, E. 1987, A&A, 183, L3

Maruyama, T., Niita, K., Oyamatsu, K., Maruyama, T., Chiba, S., & Iwamoto, A. 1998, Phys.Rev.C, 57, 655

Ravenhall, D. G., Pethick, C. J., & Wilson, J. R. 1983, Phys.Rev.Lett, 50, 2066

Watanabe, G., Iida, K., Sato, K. 2000, Nucl. Phys. A676, 455; 2001, Nucl. Phys. A687, 512; 2003, Nucl. Phys. A726, 357

Watanabe, G., Sato, K., Yasuoka, K., & Ebisuzaki, T. 2002, Phys.Rev.C, 66, 012801(R); 2003, Phys.Rev.C, 68, 035806 (nucl-th/0308007)


text:
# References

Physical origin of ”laser induced” molecule alignment effect.

V.A.Kuz’menko

Troitsk Institute for Innovation and Fusion Research,

Troitsk, Moscow region, 142190, Russian Federation.

## Abstract

The alternative to dynamic alignment explanation of experimental results on spatial-asymmetric dissociation of molecules in a laser field is proposed. The concept of geometrical alignment is sufficient for explanation of these results. In this case the spatial anisotropy of interaction of molecules with laser radiation is transferred from one field to another through the ordinary mechanism of nonlinear optical interactions. Thus laser radiation does not create alignment, but only registers it. A physical basis of such nonlinear processes is inequality of forward and reversed optical transitions that corresponds to a concept of time invariance violation in electromagnetic interactions. Directions of the further researches in the field of alignment spectroscopy are discussed.

PACS number: 42.50.Gy

Last years the impressing successes are achieved in the field of control of atoms motion with the help of laser radiation . Laser radiation allows to keep the atoms long time and to cool its up to ultra low temperature.

For molecules the situation is quit different. In the field of keeping and cooling of molecules by laser radiation only minimal successes are achieved . Here the effects, connected to orientation of molecules in space are usually studied. These effects are divided into two types: passive (geometrical) and active (dynamic) alignment . The concept of geometrical alignment is connected with the generally recognized opinion, that the probability of optical transition depends on orientation of the molecule in space relative to the direction and polarization of the laser beam . How strong is this dependence? What is it shape? Such information is absent for present day. However, it has a fundamental meaning. At the same time the opportunity to obtain such information exists and will be discussed below.

The dynamic alignment effect usually is connected with the assumption that intense laser field is capable to create the certain orientation of molecules in space. It is supposed, that the dipole moment induced in a polarizable molecule by the laser field can hinder free rotation of molecules. Molecules appear in the so-called pendular states, in which the molecular axis librates about the electric field vector. When the laser radiation is switched off, the molecules continue free rotation . Such assumption seems fantastic. The situation is interesting, because of any reliable experimental proofs of existence of the dynamic alignment effect are absent. At the same time this assumption is widely used for an explanation of experimental results. In turn, these experimental results are used as the proofs of the existence of the effect of dynamic alignment.

Some doubts in accuracy of such proof arose earlier . The purpose of the present article is to add the big piece of such doubts and to show, that for explanation of experimental results there are enough conception about the geometrical alignment which anisotropy can be transferred from one laser field to another by ordinary nonlinear interactions.

We shall consider in more detail typical experiment on dynamic alignment of iodine molecules . Cooled to low temperature in jet the iodine molecules interact with the powerful nonresonant laser pulse with linear polarization. Simultaneously with it or later the molecules are dissociated by the second powerful laser pulse. Then the fragments of molecules are ionized and their angular distributions are studied. Thus without the first laser field the dissociation of molecules is isotropic. However, when the first laser pulse is present, the strong spatial anisotropy of molecule dissociation is observed. The effect is explained as a result of the dynamic alignment of molecules and, in turn, it serves as the proof of existence of dynamic alignment.

The alternative explanation assumes the lack of any dynamic alignment. There exist only the dependence of optical transition probability from orientation of molecules in space (geometric alignment effect) and nonlinear optical phenomena. Laser radiation does not create orientation of molecules, but only registers it. The situation can be discussed in a popular now terms of optical information storage \[9-11\]. Radiation of the first laser writes down the information and radiation of the second laser reads out this information. Reading of the information occurs only at a certain orientation of the molecules in space, which is set by a direction of the first laser beam. Thus the effect of geometrical alignment is transferred from one laser field to another.

The physical mechanism of process can be described as follows. Under action of the first laser field the Raman transitions occur between various rotational sublevels of ground state of iodine molecules (Fig.1). Thus the laser radiation interacts most effectively only with the molecules, which have appropriate orientation in space. So, the effect of geometric alignment takes place. The essence of the explanation takes into consideration, that the spectroscopic properties of molecules, which participated in Raman transition, differ from those for molecules, which did not interact with laser radiation (usually similar problems are discussed within the framework of the ”coherency” concept). The molecules have property of memory about the initial state and aspiration to return in this condition.

To get precisely in the initial state the molecule should have certain orientation in space concerning the direction of the second laser beam. Under action of the second laser field the reversed in the initial state Raman transition occurs. The experiments show, that the cross-section of the reversed in the initial state optical transition can considerably exceed cross-section of the forward optical transition . Therefore radiation of the second laser field interacts mainly with the molecules, which participated in forward optical transitions. If the quantum of the second laser field is great enough, or a multiphoton process takes place, the transition into the initial state is accompanied by dissociation of grater or smaller part of molecules. Thus the effect of geometric alignment is transferred from one laser field to another as a result of ordinary process of four-photon mixing.

Representation about inequality of forward and reversed optical transitions lays in a basis of the given alternative explanation of ”induced” anisotropy of dissociation of molecules. The same representation underlies a physical explanation of the majority nonlinear optical phenomena . This representation corresponds to a principle of time invariance violation in electromagnetic interactions. Such conclusion is very strong, but it has a physical basis  and the direct experimental proof . Recently the new experimental data about T-invariance violation in optics are received . Authors believe that they find unambiguous experimental evidence of a time non-reversal interaction between light and metallic planar chiral nanostructures.

One of the most important problems of experimental researches in optics is independent measurement parameters of forward and reversed optical transitions. It is an uneasy task . The certain information on the given problem can be received in the field of discussed experiments on alignment spectroscopy. We keep in mind a more detailed studying the transition of process from nonadiabatic to adiabatic mode. Under adiabatic regime the case is understood, when the duration of the front and back edges of the first laser pulse exceeds the duration of the rotational period of molecules. In this case all molecules have an opportunity to interact effectively with laser radiation. The alignment effect is observed during action of the first laser pulse, but is absent after its ending. The physical essence of the phenomenon, probably, can be explained as follow: the front edge of a laser pulse carries out the forward Raman transitions and the low intensive back edge of a laser pulse carries out mainly the reversed Raman transitions in the initial state. As a result the memory of molecules is deleted and the ”coherency” disappears (in nonlinear optics the concept of inequality of forward and reversed transitions is, obviously, a physical equivalent, base of the concept of ”coherency”).

When the short laser pulse is used (much less than the duration of the molecular rotational period), the nonadiabatic mode is realized. In this case only small part of molecules, which have proper orientation in space, interacts effectively with laser radiation. These molecules can not return in the initial state and the alignment effect is observed for a long time as short periodic revivals . So, the detailed studying of efficiency of the back edge of the pulse on suppression of alignment effect can give certain information about ratio of cross-sections of forward and reversed optical transitions.

The major direction of alignment spectroscopy should become experimental studying of dependence of probability of optical transition from orientation of molecules in space relative to the direction of a laser beam. The basic way for getting of such information will consist in studying the consecutive optical transitions in crossed laser beams. The simplest example of such transitions is a Raman transition under action of two laser fields with different frequency. In work  the directions of two laser beams coincided and the authors observed, that the character of spectrum was changed when the polarization direction of second laser field was rotated. The effect was explained as the result of dynamic alignment and creation of pendular states. Obviously, it is wrong interpretation. We have two consecutive optical transitions. On the first transition the radiation of the first laser selects molecules with the certain orientation. That is the geometrical alignment takes place. Radiation of the second laser induces transition in these already selected and oriented molecules. So, the character of this transition should depend on the direction and polarization of radiation of the second laser beam.

The next step should be experimental studying of optical transitions in crossed laser beams. Dependence of efficiency of Raman transition on the angle between the laser beams will give the information about the dependence of probability of optical transition from orientation of molecules in space relative to a direction of the laser beam. It is possible to use various optical schemes on the basis of a four-photon mixing process. It is possible to try to receive the information separately about forward and separately about reversed optical transitions. It is an extensive work for many years.

In conclusion, the alternative explanation of the molecule dynamic alignment effect is proposed. It takes into account only the geometric alignment effect and nonlinear character of a molecule-field interaction. A basis of nonlinear character of such interactions is inequality of forward and reversed optical transitions that corresponds to a principle of time invariance violation in electromagnetic interactions. Directions of the further researches in the field of alignment spectroscopy are discussed.

Fig.1 Four-photon mixing scheme for alternative explanation of molecule alignment experiments.


text:
# Polarization control of optical transmission of a periodic array of elliptical holes in a metal film

## Abstract

Spectral dependencies of polarized optical transmission of a metal film with a periodic array of elliptical nanoholes have been studied. Such nanostructured metal films exhibit the enhanced broadband optical transmission which can be controlled by selecting polarization of incident and/or transmitted light.

Optical components for applications in extreme ambient conditions and high level of integration require novel design approaches. Recent progress in understanding of the optical properties of nanostructured metal films has allowed the design of spectral filters and devices which exploit surface plasmon polariton excitations (SPPs) in metal nanostructures.ebbessen-98 ; jopa  Such nanostructures could act as band-pass spectral filtersebbessen-98 ; krishnan-01  or, if combined with nonlinear materials, could perform all-optical stitching and limiting functions at low control light intensities.prl-02 ; apl-02 ; prb-02 ; dykhne-03  Polarisation properties of light can also be controlled with metallic nanostructures.zheludev-03 ; sanders

Functionality of metallic nanostructures are determined by the interaction of photons with collective excitations of conduction electrons at a metal surface.agran  The coupling and mutual transformations between photons and surface plasmons depend upon the frequencies and wave vectors of the excitations, as well as on the polarization state of incident electromagnetic field. The polarization state determines the excitation of SPPs in one or another generally non-equivalent directions of the Brillouin zone related to the periodic structure. By choosing the light polarization one could not only provide the optimal intensity transmittance through the structure but could also introduce various polarization-sensitive effects related to the chirality of the structure and non-locality of its electromagnetic response. zheludev-03 ; schwanecke-03

In this letter we report the studies of polarization properties of the enhanced light transmission through a metal film with a periodic arrangement of oriented elliptical nanoholes. Such metallic structure exhibits a broadband transmission stemming from the multiple SPP resonances permitted by the low-symmetry arrangement of elliptical holes in a square array. The transmission spectrum strongly depends on the polarization state of the incident light even at normal incidence, the property of striking contrast with a square array of circular holes where transmission spectrum has well-defined resonances independent on the polarization of the incident light.prl-01  In the case of elliptical holes, the incident linearly polarized light is converted into elliptically polarized in transmission. Using this feature we can continuously tune the transmission spectrum of the structure from blue to red wavelengths by changing the polarization state of the incident light and selecting appropriate polarization component of the transmitted light.

The metallic nanostructures of overall size of $`20\times 20`$ $`\mu ^2`$ studied here were fabricated in a 40 nm gold film on a glass substrate using ion beam milling. The structure consists of elliptical holes with the main axes of the ellipse of 250 and 500 nm (Fig. 1). The holes are arranged in a square array with period $`D=`$500 nm. The transmission spectra of the nanostructures were recorded using a spectrograph equipped with a CCD-camera coupled to a long-working-distance optical microscope. Collimated white light from a stabilized tungsten-halogen source passed through a polarizer and illuminated the sample at normal incidence. The light transmitted through the nanostructure was collected by an objective lens. The beam position at the nanostructure was monitored by an imaging CCD camera. The light exiting the structure passed through the analyzer and was launched into a fiber bundle connected to the spectrograph. Polarization mode scrambling in the fiber bundle ensured that polarization-related effects do not influence spectral measurements.

The transmission spectra of the structure obtained at normal incidence without polarization selectivity of the transmitted light are shown in Fig. 2 for different polarization states of the incident light with respect to the array orientation. For the incident light polarization parallel to the principal lattice axes ($`x`$ or $`y`$) the transmission spectrum is observed with broad features in the 450–550 nm and 750–800 nm spectral ranges. With the increase of the polarization azimuth $`\varphi `$ up to about 20–25o the transmission increases in the 550–750 nm wavelength range, where the spectrum becomes rather flat. With further increase of $`\varphi `$ when the light is polarized along the short principal axis of the elliptical holes, the transmission in the long-wavelength spectral range (750–800 nm) is significantly suppressed. For the polarisation azimuth along the long principal axis of ellipse, the transmission at the far red part of the spectrum is enhanced, but becomes smaller at around 600–700 nm. Thus, with rotation of the polarization of the incident light, the intensity of the spectral components of the transmitted light exhibit complex, oscillatory behavior different in different spectral ranges. The relatively well-defined long-wavelength band at around 780–800 nm has transmission minimum for the polarization of light parallel to the short principal axis of an ellipse and maximum for the polarization along the long principal axis. For all wavelengths, the dependencies show a 2-fold symmetry corresponding to the symmetry of the array: a square lattice with basis elements of a 2-fold rotational symmetry.

To understand polarization sensitivity of transmission, the spectrum of the SPP Bloch waves on a periodically structured film should be considered. The strict consideration of the symmetry properties of the eigenstates of a periodic structure with an asymmetrical lattice basis would require the approach analogous to the linear combination of atomic orbitals (LCAO) method or the Wannier fanctions conventionally used in solid state physics.kittel  Nevertheless, it will be instructive to consider a simplified model using the SPP Bloch modes.

The field of SPP modes on a periodically structured surface can be described by

$$E(x,y)=U_{𝐤_{SP}}(\xi ,\eta )exp[i(k_{SP}^{(x)}x+k_{SP}^{(y)}y)],$$
(1)
where $`U_{𝐤_{SP}}(\xi ,\eta )`$ possesses the periodicity of the array and is the SPP Bloch function and $`𝐤_{SP}`$ is the wave vector of the Bloch wave. Here the coordinate systems ($`\xi `$,$`\eta `$) of the basis and ($`x`$,$`y`$) of the lattice are not generally independent and introduced to simplify symmetry considerations. At normal incidence of light, the allowed SPP wavevectors on the periodically structured surface are determined by

$$𝐤_{SP}=\pm p\frac{2\pi }{D}𝐮_x\pm q\frac{2\pi }{D}𝐮_y,$$
(2)
where $`𝐮_x`$ and $`𝐮_y`$ are the unit reciprocal lattice vectors of the periodic structure, $`D`$ is its periodicity (assumed to be the same in both $`x`$\- and $`y`$-directions), and $`p`$ and $`q`$ are integer numbers corresponding to the different directions of the SPP Brillouin zone. At normal incidence, SPP can be excited if the electric field of the incident light has a component in the direction of SPP propagation: $`(𝐄𝐤_{SP})0`$. Thus, the polarization dependence of the coupling efficiency in different directions of the Brillouin zone is proportional to $`|pcos\varphi +qsin\varphi |`$, where $`\varphi `$ is the polarization azimuth angle with respect to the $`x`$ principal axis of the lattice. It is clear from this analysis that for a square lattice with basis elements of circular symmetry no polarization dependencies of the SPP excitation and, thus, of the enhanced transmission can be expected at normal incidence.

However, in the case of an elliptical basis of the lattice, polarization dependencies are significant. Taking into account different scattering properties of the ellipse in different directions, one can understand this by considering the mixing of the SPP Bloch states of the square lattice due to the lower symmetry of the basis elements. Elliptical holes with different size in the different directions modify the Brillouin zone structure by introducing size-dependent anizotropic scattering. This leads to the reduction of the Brillouin zone symmetry showing a 2-fold rotation axis instead of 4-fold, appropriate for a square lattice, effectively creating a medium with pronounced linear birefringence in dichroism. Linearly polarized waves with polarization azimuth along the main axes of the elliptical basis will be the polarization eigenstates of the structure. As a result of the anizotropic SPP Brillouin zone, the transmission spectrum of the structure is dependent on the polarization of the incident light.

Mixing and re-excitation of surface plasmon resonances in different directions of the Brillouin zone are responsible for anizotropic retardation and eventually for the elliptization of the polarization state of the transmitted light. Polarized transmission spectra show that the changes of the polarization state of the transmission light are strongly wavelength dependent (Fig. 3). This allows identification of transmission resonances contributing to the broad transmission spectrum. Several resonances can be identified by the different degree of elliptization at around 750, 650, and 550 nm. Thus, by controlling the polarization state of the incident and exercising polarization selection of the transmitted light, transmission spectrum of the elliptical hole array can be tuned at specific resonant band corresponding to one or another SPP mode of the system, quasi-continuously from blue to red wavelengths of the visible spectrum (Fig. 4).

In conclusion, we have studied polarization properties of the transmission spectra of a periodic array of elliptical holes in a metal film. Due to the reduced symmetry of the lattice, elliptical holes allow to achieve the enhanced broadband transmission which can be controlled by choosing the polarization of the incident and/or transmitted light, offering numerous applications in passive and active optical components. Combined with non-linear effects which are significantly enhanced in metallic structures due to surface plasmons, the observed polarization properties will allow further applications in integrated photonic circuits.

###### Acknowledgements.

We are indebt to C. Cassidi for sputtering thin gold films used in these experiments. This work was supported in part by EPSRC.


text:
# Measurements of Event and Jet Shapes at HERA

## 1 Introduction

The HERA collider is a precision tool for the study of QCD. The H1 and ZEUS Collaborations have performed a wealth of analyses which investigate the hadronic final state in $`ep`$ collisions. In addition to the measurement of plain jet cross sections, it is rewarding to study the substructure of jets, or in the case of event shapes, the structure of the whole hadronic final state. In this case, uncertainties connected to the determination of cross sections for jet production do not compromise the results.
In the regime of high transverse momenta, experimental uncertainties and non-perturbative effects are small, which makes possible precision tests of perturbative QCD and the extraction of the strong coupling constant.
On the other hand, understanding the nature of hadronisation is a challenge, which is important for all high energy experiments involving hadrons. In this area, the concept of power corrections offers prospects beyond phenomenological models implemented in Monte Carlo generators.

## 2 Subjet Multiplicities

The internal structure of jets is determined by QCD radiation processes. Non-perturbative fragmentation effects become small at high transverse jet energies $`E_T`$ and values of the resolution scale $`y_{\mathrm{cut}}`$ that are not too low. Perturbative QCD calculations can be applied. In this region of phasespace they depend only weakly on parton density functions, but are directly proportional to the strong coupling constant $`\alpha _s`$.
The mean subjet multiplicity is defined as

$$n_{\mathrm{sbj}}(y_{\mathrm{cut}})=\frac{1}{N_{\mathrm{jets}}}\underset{i=1}{\overset{N_{\mathrm{jets}}}{}}n_{\mathrm{sbj}}^i(y_{\mathrm{cut}}),$$

where $`n_{\mathrm{sbj}}^i(y_{\mathrm{cut}})`$ is the number of subjets in jet $`i`$ and $`N_{\mathrm{jets}}`$ is the total number of jets in the sample.
The ZEUS collaboration has measured the mean subjet multiplicity in neutral current deep inelastic scattering (NC DIS) as a function of $`E_T^{\mathrm{jet}}`$Chekanov:2002ux . Jets were defined by using the longitudinally invariant $`k_T`$ cluster algorithm in the laboratory frame. In Fig. 1 the measurements of $`n_{\mathrm{sbj}}`$ are compared with the NLO prediction from DISENTCatani:1997vz , corrected for hadronisation efects. The overall description of the data is good, so that the measurements can be used to make a determination of $`\alpha _s`$. The region used for the fit was restricted to $`E_T^{\mathrm{jet}}>25\text{GeV}`$, where hadronisation corrections differ by less than $`17\%`$ from unity.
The QCD fit of the mean subjet multiplicity for $`25<E_T^{\mathrm{jet}}<71\text{GeV}`$ at $`y_{\mathrm{cut}}=10^2`$ yields

$$\alpha _s(m_Z)=0.1187\pm 0.0017(\text{stat.})_{0.0009}^{+0.0024}(\text{syst.})_{0.0076}^{+0.0093}(\text{th.}).$$

This result is consistent with recent determinations by the H1 Adloff:2000tq ; Adloff:2000qk  and ZEUS Chekanov:2002be ; Chekanov:2002pv  Collaborations and with the PDG value Hagiwara:2002fs . The total error is dominated by the theoretical uncertainty due to missing terms beyond NLO.

## 3 Substructure Dependence of Jet Cross Sections

Quarks and gluons carry different colour charges, therefore the internal structure of a jet is expected to depend on the type of the primary parton from which it originated. Obervables which probe the internal structure are the subjet multiplicity (section 2) and the integrated jet shape $`\mathrm{\Psi }(r)`$, which is defined as

$$\mathrm{\Psi }(r)=\frac{E_T(r)}{E_T^{\mathrm{jet}}},$$

where $`E_T(r)`$ is the transverse energy within a cone of radius $`r=\sqrt{\mathrm{\Delta }\varphi ^2+\mathrm{\Delta }\eta ^2}`$ about the jet axis.
In photon-proton reactions, two types of processes contribute to jet production. In the direct process the photon interacts with a parton in the proton as a whole, whereas in the resolved process it behaves as a source of partons which subsequently interact with partons from the proton. Jets from direct and resolved processes exhibit different angular distributions and $`E_T`$ spectra. Consequently, dijets from both processes differ in their distributions of invariant mass $`M^{jj}`$ and the angle between the jet-jet axis and the beam direction in the dijet centre-of-mass system $`\mathrm{cos}\theta ^{}`$. The cross section for subprocesses with gluons in the final state is larger for the resolved than for the direct process, if an appropriate region of phase space is chosen. Therefore, the experimental investigation of jet substructure as a function of jet kinematics allows one to study the characteristics of quark- and gluon-initiated jets.
The ZEUS Collaboration has measured differential inclusive jet and dijet cross sections in photoproduction ZEUS:2003eps518 . The internal structure of jets, in terms of the integrated jet shape and subjet multiplicity, has been used to select samples enriched in quark- and gluon-initiated jets, classified as “thin” and “thick” samples, respectively. Fig. 2 shows the measured inclusive jet cross section as a function of $`E_T`$ and $`\eta _{\mathrm{jet}}`$, and the dijet cross section as a function of $`|\mathrm{cos}\theta ^{}|`$ and $`M^{jj}`$. The distributions are in general well described by calculations using PYTHIA, which have been normalised to the total measured cross section. The “thick” data sample behaves according to expectations for gluon initiated jets in the ivestigated phase space: these jets are more forward, have a softer $`E_t`$ spectrum (inclusive jets), lower invariant mass and the jet axis is closer to the beam direction (dijets). This behaviour is consistent with the expected dominance of the resolved subprocess for gluon initiated jets.

## 4 Event Shapes

Event shape variables are designed to study QCD by measuring properties of the flow of hadronic energy-momentum. A suitable frame of reference in DIS is the Breit frame, which divides an event into hemispheres corresponding to the proton remnant and the hadronic final state evolving from the struck parton. Event shape observables are then usually defined for the particles in the current hemisphere (CH) alone, thereby avoiding non-perturbative effects connected with the proton remnant. A commonly used event shape variable is thrust $`\tau 1T`$:

$$T=\frac{_{hϵ\mathrm{CH}}|p_{zh}|}{_{hϵ\mathrm{CH}}|p_h|},$$

where $`p_h`$ is the momentum of particle $`h`$ in the current hemisphere and $`p_{zh}`$ this momentum projected on the axis of the exchanged virtual photon. The description of such infrared-collinear safe observables by fixed order calculations faces two difficulties. Firstly, non-perturbative hadronisation effects can be large, even at formally perturbative scales. In the approach by Dokshitzer, Webber et al. Dokshitzer:1997ew  these effects can be addressed using power corrections proportional to $`1/Q`$, governed by one universal free parameter $`\alpha _0`$. Secondly, the convergence of the perturbative series at low values of the event shape variables is very poor, making a resummation to all orders in $`\alpha _s`$ necessary. Resummed NLL calculations, matched to NLO, have recently become available for two kinds of thrust, $`\tau _c`$ and $`\tau `$, the jet broadening $`B`$, the jet mass $`\rho _0`$ and the C-parameterDasgupta:2002dc . While these event shapes are mainly sensitive to the transition from 1-jet like to 2-jet like topologies, there are now also measurements and calculations available for the out-of-plane momentum $`K_{\mathrm{out}}`$ and the azimuthal correlation $`\chi `$, which probe the transition from 2-jet to 3-jet topologies Banfi:2001ci .
Mean values of event shapes as a function of $`Q`$ have been studied by the H1 Adloff:1999gn  and ZEUS Chekanov:2002xk  Collaborations. The H1 Collaboration has also measured spectra of event shape variables H1:2003eps111 . Fig. 3 shows the distributions of thrust and out-of-plane momentum in bins of momentum transfer $`Q`$. The thrust spectra are compared with the results of a fit based on resummed calculations matched to NLO and including power corrections. For the out-of-plane momentum the matching is not yet available, therefore a comparison with the RAPGAP Monte Carlo event generator Jung:1995gf  is presented. This includes LO matrix elements and parton showers. Both observables are well described over the whole region of phase space.
The ZEUS Collaboration has made a fit to the mean values of event shape variables using DISASTER++ Graudenz:1997gv  and the Dokshitzer-Webber power corrections, as shown in Fig. 4 (left). The contours correspond to 1-$`\sigma `$ statistical errors. Consistent values of $`\alpha _s`$ are obtained for four of the shape variables with $`\alpha _0`$ values that agree within $`10\%`$. However, the lack of overall consistency in $`\alpha _s`$ and $`\alpha _0`$ suggests that higher order processes which are not included may be important. On the right hand side of Fig. 4 results of a fit to event shape spectra by the H1 Collaboration are shown, including resummed calculations. Here the contours indicate 1-$`\sigma `$ statistical and experimental systematic errors, added in quadrature. The results for the five studied event shapes are consistent with each other and with the world average of $`\alpha _s`$ Bethke:2000ai . The universal non-perturbative parameter $`\alpha _0`$ is confirmed to be 0.5 at the $`10\%`$ level.

## 5 Conclusion

Differential distributions of jet and event shapes have been measured by the H1 and ZEUS Collaborations with high precision. Calculations based on perturbative and non-perturbative QCD describe the measurements over a large region of phasespace.


text:
# Einstein, Podolsky, Rosen, and Shannon*footnote **footnote *Dedicated to the memory of James T. Cushing.

## Abstract

The Einstein-Podolsky-Rosen paradox (1935) is reexamined in the light of Shannon’s information theory (1948). The EPR argument did not take into account that the observers’ information was localized, like any other physical object.

I had the privilege of meeting Jim Cushing in 1986, during a conference that Daniel Greenberger had organized in the World Trade Center in New York City . Some time later Cushing sent me a copy of his book , where our encounter is related in a footnote:

> “I mentioned to Peres that his position appeared to be an instrumentalist one. He replied with no apparent discomfort that others had told him that before. For a physicist’s statement on an instrumentalist interpretation of quantum mechanics, see Peres (1988) .”

Jim had asked me that question with the same tone as if he were asking whether I was a cannibal. For a more recent discussion on this subject, see .

On the other hand, I never had the privilege of meeting Albert Einstein. He died when I was an undergraduate. I had always been fascinated by Einstein, as any normal Jewish boy would be, and later in my life I even got the impression that I came to know him personally. This is because my PhD thesis advisor was Nathan Rosen, who had been a close collaborator of Einstein. Together they built the Einstein-Rosen bridge in General Relativity , and together with Podolsky they formulated the famous EPR paradox . Rosen told me many anecdotes about Einstein and his reactions to various events. (Rosen’s wife Hanna, who was an accomplished pianist, gave piano accompaniment to Einstein who played the violin.)

My first encounter with the EPR paper occurred when I was a graduate student, circa 1958. The subject of my research was gravitational radiation. At that time, I was rather ignorant of quantum theory, having graduated in mechanical and nuclear engineering, not in physics. One day, I came into Rosen’s office, and I found him sorting out his papers. On the ground, there were cartons full of old reprints, with the characteristic green covers of The Physical Review. One of them read: “Can Quantum-Mechanical Description of Physical Reality Be Considered Complete?” The authors were A. Einstein, B. Podolsky, and N. Rosen. I thought that this would be a nice item to have in my collection of reprints, and politely asked:

“Professor Rosen, may I take one of these reprints?”

He looked concerned.

“Ha, how many are left?”

We counted them, which was not difficult: there were two. Then he said, with some hesitation,

“Well, if there are two, you may have one.”

This is how I acquired the last available reprint of the famous article of Einstein, Podolsky, and Rosen.

Let us have a look at that wonderful paper. You will immediately notice that Eqs. (7) and (8) involve entangled wave-functions, and indeed the whole issue is about the physical consequences of such an entanglement. Entangled wave-functions were not new at that time: you can find one, for example, in Eq. (10) of Rosen’s 1931 seminal paper on the ground state of the hydrogen molecule , which is probably more famous among chemists than the EPR paper is among physicists.

Some time after that work, Rosen became a post-doc of Einstein at the Institute of Advanced Studies in Princeton. One day, at the traditional 3 o’clock tea, Rosen mentioned to Einstein a fundamental issue of interpretation related to entangled wave-functions. Einstein immediately saw the implications for his long standing disagreement with Bohr. As they discussed the problem, Boris Podolsky joined the conversation, and later proposed to write an article. Einstein acquiesced. When he later saw the text, he disliked the formal approach, but agreed to its publication. Then, as soon as the EPR article appeared, Podolsky relased its contents to the New York Times (4 May 1935, page 11) in a way implying that the authors had found that quantum mechanics was faulty. This infuriated Einstein, who after that no longer spoke with Podolsky.

The EPR “paradox” drew immediate attention. Niels Bohr  found the reasoning faulty, because it contradicted his complementarity principle. Bell, in his first article on hidden variables and contextuality , wrote “the Einstein-Podolsky-Rosen paradox is resolved in the way which Einstein would have liked least.” Actually, the example given by Bell in the proof of his celebrated theorem  is based on a much simpler entangled system: two spin-$`\frac{1}{2}`$ particles in a singlet state . In 1991, David Mermin came to Technion to give the annual Wunsch Lecture. He said that the EPR paper was wrong. After the talk, there was as usual a discussion period. Nathan Rosen politely commented: the paper is not wrong, it makes some assumptions, and then draws the logical conclusions; the assumptions were wrong.

The EPR article was not wrong, but it had been written too early. Only some years later, in 1948, Claude Shannon published his theory of information  (and it took many more years before the latter was included in the physicist’s toolbox). Shannon was employed by the Bell Telephone Company and his problem was to make communication more efficient. Shannon showed that information could be given a quantitative measure, that he called entropy. It was later proved that Shannon’s entropy is fully equivalent to ordinary thermodynamical entropy . Information can be converted to heat and can perform work. Information is not just an abstract notion . It requires a physical carrier, and the latter is (approximately) localized. After all, it was the business of the Bell Telephone Company to transport information from one telephone to another telephone, in a different location.

In the EPR article, the authors complain that “it is possible to assign two different wave functions to … the second system,” and then, in the penultimate paragraph, they use the word simultaneous no less than four times, a surprising expression for people who knew very well that this term was undefined in the theory of relativity. Let us examine this issue with Bohm’s singlet model. One observer, conventionally called Alice, measures the $`z`$-component of the spin of her particle and find $`+\mathrm{}/2`$. Then she immediately knows that if another distant observer, Bob, measures (or has measured, or will measure) the $`z`$-component of the spin of his particle, the result is certainly $`\mathrm{}/2`$. One can then ask: when does Bob’s particle acquire the state with $`s_z=\mathrm{}/2`$?

This question has two answers. The first answer is that the question is meaningless — this is undoubtedly true. The second answer is that, although the question is meaningless, it has a definite answer: Bob’s particle acquires this state instantaneously. This then raises a new question: in which Lorentz frame is it instantaneous? Here, there is also a definite answer: it is instantaneous in the Lorentz frame that we arbitrarily choose to perform our calculations . Lorentz frames are not material objects: they exist only in our imagination.

When Alice measures her spin, the information she gets is localized at her position, and will remain so until she decides to broadcast it. Absolutely nothing happens at Bob’s location. From Bob’s point of view, all spin directions are equally probable, as can be verified experimentally by repeating the experiment many times with a large number of singlets without taking in consideration Alice’s results. Thus, after each one of her measurements, Alice assigns a definite pure state to Bob’s particle, while from Bob’s point of view the state is completely random ($`\rho `$ is proportional to the unit matrix). It is only if and when Alice informs Bob of the result she got (by mail, telephone, radio, or by means of any other material carrier, which is naturally restricted to the speed of light) that Bob realizes that his particle has a definite pure state. Until then, the two observers can legitimately ascribe different quantum states to the same system. For Bob, the state of his particle suddenly changes, not because anything happens to that particle, but because Bob receives information about a distant event. Quantum states are not physical objects: they exist only in our imagination.

In summary, the question raised by EPR “Can quantum-mechanical description of physical reality be considered complete?” has a positive answer. However, reality may be different for different observers.

This work was supported by the Gerard Swope Fund.


text:
# Dissipationless shock waves in repulsive Bose-Einstein condensates

## Abstract

We consider formation of dissipationless shock waves in Bose-Einstein condensates with repulsive interaction between atoms. It is shown that big enough initial inhomogeneity of density leads to wave breaking phenomenon followed by generation of a train of dark solitons. Analytical theory is confirmed by numerical simulations.

Experiments on free expansion of Bose-Einstein condensate (BEC) have shown  that evolution of large and smooth distributions of BEC is described very well by hydrodynamic approximation  where dispersion and dissipation effects are neglected. At the same time, it is well known from classical compressible gas dynamics (see, e.g. ) that typical initial distributions of density and velocity can lead to wave breaking phenomenon when formal solution of hydrodynamical equations becomes multiple-valued. It means that near the wave breaking point one cannot neglect dispersion and/or dissipation effects which prevent formation of a multiple-valued region of a solution. If the dissipation effects dominate the dispersion ones, then the multiple-valued region is replaced by the classical shock, i.e. narrow layer with strong dissipation within, which separates smooth regions with different values of density, fluid velocity and other physical parameters. This situation was studied in classical gas dynamics and found many practical applications. If, however, the dispersion effects dominate dissipation ones, then the region of strong oscillations is generated in vicinity of wave breaking point . Observation of dark solitons in BEC  shows that the main role in dynamics of BEC is played by dispersion and nonlinear effects taken into account by standard Gross-Pitaevskii (GP) equation , and dissipation effects are relatively small and can be considered as perturbation. Hence, there are initial distributions of BEC which can lead to formation of dissipationless shock waves. Here we shall consider such possibility.

The starting point of our consideration is the fact that the sound velocity in BEC is proportional to the square root from its density (see, e.g.  and references therein). Thus, if we create inhomogeneous BEC with high density hump (with density $`\rho _1`$) in the center of lower density distribution (with density $`\rho _0`$), and after that release this central part of BEC, then the high density hump will tend to expand with velocity $`\sqrt{\rho _1}`$ greater than the sound velocity $`\sqrt{\rho _0}`$ of propagation of disturbance in lower density BEC. As a result, wave breaking and formation of dissipationless shock wave can occur in this case. Note that initial distribution of this type was realized in the recent experiment  where generation of oscillations was also observed. In this experiment the lower density disk-shaped BEC was confined by magnetic trap and density distribution had standard Thomas-Fermi (TF) parabolic form. Additional potential was applied to BEC in the central part of above TF profile which lead to narrow parabolic hump in the density distribution. After the central potential was switched off, the hump started to expand against wide lower density TF profile leading to generation of oscillations in the transition region between high and low density condensates. The theory of dissipationless shock waves in media described by one-dimensional (1D) nonlinear Schrödinger (NLS) equation was developed in . Since the GP equation in some cases can be reduced to the 1D NLS equation, this theory can be applied to description of dissipationless shock waves in BEC. Here we suggest such description and confirm it by numerical simulations.

We consider BEC confined in a “pancake” trap with the axial frequency $`\omega _z`$ much greater than the transverse one $`\omega _x=\omega _y=\omega _{}`$. Let the density of atoms in the central part of BEC be of order of magnitude $`n_0`$ and satisfy the condition $`n_0a_sa_z^21`$, where $`a_s>0`$ is the $`s`$-wave scattering length and $`a_z=(\mathrm{}/ma_z)^{1/2}`$ is the amplitude of quantum oscillations in the axial trap. Then the condensate wavefunction $`\psi `$ can be factorized as $`\psi =\varphi (z)\mathrm{\Psi }(x,y)`$, where $`\varphi (z)=\pi ^{1/4}a_z^{1/2}\mathrm{exp}(z^2/2a_z^2)`$ is the ground state wavefunction of axial motion, and $`\mathrm{\Psi }(x,y,t)`$ satisfies the reduced 2D GP equation

$$i\mathrm{}\mathrm{\Psi }_t=\frac{\mathrm{}^2}{2m}(\mathrm{\Psi }_{xx}+\mathrm{\Psi }_{yy})+V(x,y)\mathrm{\Psi }+g_{2D}|\mathrm{\Psi }|^2\mathrm{\Psi },$$
(1)
where $`V(x,y)`$ is the potential of a transverse trap, $`g_{2D}=2\sqrt{2\pi }\mathrm{}^2a_s/(ma_z)`$ is the effective nonlinear interaction constant, and $`\mathrm{\Psi }`$ is normalized to the number of atoms, $`|\mathrm{\Psi }|^2𝑑x𝑑y=N`$. The initial distribution of density is determined by the potential $`V(x,y)`$ and consists of wide background with a hump in its center. We assume that the background width is much greater than the hump’s width, so that at the initial stage of evolution we can consider an expansion of the central part against the constant background. In a similar way, at the initial stage of evolution, when the radius of the central part does not change considerably, we can neglect the curvature of axially symmetrical distribution and consider its 1D cross section. It means that we can neglect the dependence of $`\mathrm{\Psi }`$ on $`y`$ coordinate and consider $`\mathrm{\Psi }`$ as a function of $`x`$ and $`t`$ only. As a result, we arrive at 1D NLS equation with inhomogeneous initial distribution of density. To simplify the notation, we introduce dimensionless variables $`t^{}=\omega _zt/2`$, $`x^{}=x/a_z,`$ $`u=(2\sqrt{2\pi }a_sa_z/n_0)^{1/2}\mathrm{\Psi }`$. Then the initial stage of evolution of the wavefunction profile in the $`x`$ axis cross section is governed by the NLS equation,

$$iu_t+u_{xx}2|u|^2u=0,$$
(2)
where the primes in $`t^{}`$ and $`x^{}`$ are omitted for convenience of the notation.

Evolution of smooth pulses before the wave breaking point can be described in the hydrodynamic approximation which can be achieved by substitution

$$u(x,t)=\sqrt{\rho (x,t)}\mathrm{exp}\left(i^xv(x^{},t)𝑑x^{}\right)$$
(3)
and separation of the real and imaginary parts. As a result we obtain the system

$$\frac{1}{2}\rho _t+(\rho v)_x=0,\frac{1}{2}v_t+vv_x+\rho _x=0,$$
(4)
where we have neglected the so-called “quantum pressure” term with higher space derivatives what is correct until the density distribution has smooth enough profile. To get simple qualitative picture of the wave breaking of BEC density distribution, let us consider idealized case with a box-like hump in the initial distribution,

$$\rho (x,0)=\{\begin{array}{c}\rho _0,|x|>a,\hfill \\ \rho _1,|x|a.\hfill \end{array}$$
(5)
Although this distribution has a parameter with dimension of length—width of the hump—it does not play any role for $`ta/(2\sqrt{\rho _1})`$, i.e. until two waves propagating inward the hump with sound velocity $`c_1=2\sqrt{\rho _1}`$ meet at $`x=0`$. Hence, for this initial period of evolution the solution for $`x>0`$ can only depend on the self-similar variable $`\xi =(xa)/t`$, centered at $`x=a`$, and for $`x<0`$ on $`\xi =(x+a)/t`$ centered at $`x=a`$. Since the picture is symmetrical, it is enough to consider only a half of the solution corresponding to $`x>0`$. It is easy to find that the density is given by the formulae

$$\rho (x,t)=\{\begin{array}{cc}\rho _1\hfill & \mathrm{for}0<x<a2\sqrt{\rho _1}t,\hfill \\ \left[2\sqrt{\rho _1}(xa)/2t\right]^2/9\hfill & \mathrm{for}a2\sqrt{\rho _1}t<x<a+2(\sqrt{\rho _1}\sqrt{\rho _0})t,\hfill \\ \left[2\sqrt{\rho _0}+(xa)/2t\right]^2/9\hfill & \mathrm{for}a+2\sqrt{\rho _0}t<x<2+2(\sqrt{\rho _1}\sqrt{\rho _0})t,\hfill \\ \rho _0\hfill & \mathrm{for}x>a+2\sqrt{\rho _0}t,\hfill \end{array}$$
(6)
and similar formulae can be written for the velocity field $`v(x,t)`$. This solution describes the wave breaking phenomenon which takes place if $`\rho _1>4\rho _0`$. The density profile shown in Fig. 1 clearly illustrates the origin of the multiple-valued region which should be replaced by the oscillatory shock wave when the dispersion effects are taken into account.

For more realistic initial pulses the density profile is a smooth function without cusp points. In vicinity of the wave breaking point the solution can be approximated by a cubic function for one Riemann invariant $`\lambda _+=v/2+\sqrt{\rho }`$ of the system (4) and by constant value for another one $`\lambda _{}=v/2\sqrt{\rho }`$ (see ). After Galileo and scaling transformations the hydrodynamic solution can be written in the form

$$x(3\lambda _++\lambda _{})t=\lambda _+^3,\lambda _{}=\mathrm{const},$$
(7)
and again for $`t>0`$ it has a multiple-valued region of $`\lambda _+`$. It means that dispersion effects have to be taken into account which lead to formation of dissipationless shock wave after wave breaking point.

In framework of Whitham theory of modulations  one can obtain an approximate solution of the NLS equation (2) in analytic form where the dissipationless shock wave is presented as a modulated periodic nonlinear wave solution of the NLS equation. The density is expressed in terms of Jacobi elliptic function

$$\rho (x,t)=|u(x,t)|^2=\frac{1}{4}(\lambda _1\lambda _2\lambda _3+\lambda _4)^2+(\lambda _1\lambda _2)(\lambda _3\lambda _4)\mathrm{sn}^2[\sqrt{(\lambda _1\lambda _3)(\lambda _2\lambda _4)}\xi ,m],$$
(8)
where

$`\xi =x(\lambda _1+\lambda _2+\lambda _3+\lambda _4)t,m={\displaystyle \frac{(\lambda _1\lambda _2)(\lambda _3\lambda _4)}{(\lambda _1\lambda _3)(\lambda _2\lambda _4)}},`$ (9)
and parameters $`\lambda _i`$, $`i=1,2,3,4,`$ change slowly along the dissipationless shock. Their dependence on $`x`$ and $`t`$ is determined implicitly by the solution

$$xv_i(\lambda )t=w_i(\lambda ),i=1,2,3;\lambda _4=\overline{\lambda }=\mathrm{const}$$
(10)
of Whitham equations, where Whitham velocities $`v_i`$ and $`w_i`$ are given by quite complicated expressions in terms of elliptic integrals (see ):

$$w_i=\frac{8}{35}w_i^{(3)}+\frac{4}{5}\overline{\lambda }w_i^{(2)}\frac{1}{35}\overline{\lambda }^2v_i(\lambda )+\frac{1}{35}\overline{\lambda }^3,i=1,2,3,$$
(11)
$$w_i^{(k)}=W^{(k)}+(v_is_1)_iW^{(k)},$$
(12)
$$W^{(1)}=V=s_1,W^{(2)}=\frac{3}{8}s_1^2\frac{1}{2}s_2,W^{(3)}=\frac{5}{16}s_1^3\frac{3}{4}s_1s_2+\frac{1}{2}s_3,$$
(13)
$$v_i(\lambda )=\left(1\frac{L}{_iL}_i\right)V,_i/\lambda _i,i=1,2,3,4,$$
(14)
where

$$L=\frac{\mathrm{K}(m)}{\sqrt{(\lambda _1\lambda _3)(\lambda _2\lambda _4)}}$$
(15)
is a wavelength, $`\mathrm{K}(m)`$ is the complete elliptic integral of the first kind, and $`s_1,`$ $`s_2,`$ $`s_3`$ are determined by the expressions

$$s_1=\underset{i}{}\lambda _i,s_2=\underset{i<j}{}\lambda _i\lambda _j,s_3=\underset{i<j<k}{}\lambda _i\lambda _j\lambda _k.$$
(16)
Equations (11) can be solved with respect to $`\lambda _i,`$ $`i=1,2,3,`$ giving them as functions of $`x`$ and $`t`$. Subsequent substitution of these functions $`\lambda _i(x,t)`$, $`i=1,2,3,`$ into (8) yields the modulated periodic wave which represents the dissipationless shock wave. The resulting profile of density in dissipationless shock wave is shown in Fig. 2. At one its edge it consists of the train of dark solitons, and at another edge describes small amplitude oscillations propagating with local sound velocity into unperturbed region described by smooth solution of hydrodynamical equations. The modulated periodic wave replaces the multiple-valued region shown by dashed line which was obtained in hydrodynamic approximation after the wave breaking point. This multiple-valued region is analogous to that in Fig. 1 with account of change of variables due to Galileo and scaling transformations.

To check the described above picture of formation of dissipationless shock wave and to extend it to real 2D situation, we have solved numerically the 2D generalization of Eq. (2),

$$iu_t+u_{xx}+u_{yy}2|u|^2u=0,$$
(17)
with the initial condition

$$\rho (r,0)=\{\begin{array}{cc}\rho _0+(\rho _1\rho _0)(1r^2/a^2),\hfill & |r|a,\hfill \\ \rho _0,\hfill & |r|>a,\hfill \end{array}$$
(18)
(where $`r=(x^2+y^2)^{1/2}`$) similar to one studied experimentally . Plot of two-dimensional initial density distribution is shown in Fig. 3 and density distribution after time evolution $`t=2`$ is shown in Fig. 4. As we see, the parabolic hump expands with formation of dissipationless shock wave in the transition layer between high density region to low density one. To see more clearly the evolution of the hump, its cross section profiles are shown in Fig. 5 at different values of time $`t`$. Slowly propagating dark solitons are clearly seen as well as small amplitude sound waves propagating into undisturbed low density region. Dissipationless shock wave generated at the right side of the profile coincides qualitatively with results of analytic theory shown in Fig. 2. One may assume that oscillations in BEC density profile observed in experiment  have the same origin.

Further evolution of the hump will ultimately lead to its spreading over large area with small amplitude oscillations, that is the shock wave does not persist permanently. Rather, it is a transient phenomenon caused by different values of characteristic velocities in high density hump and low density background. Slowly propagating dark solitons are generated in the transient layer between these two regions in order to reconcile two different values of velocities of disturbance propagation. This mechanism of dark soliton generation is quite general and can manifest itself in various geometries and initial BEC distributions 111Formation of dark soliton trains in 1D geometry was considered in recent preprint .

In conclusion, we have studied theoretically and numerically the process of formation of dissipationless shock waves in the density distribution of BEC. We believe that oscillations in BEC density profile observed in recent experiment  can be explained by this theory.

This work was supported by FAPESP (Brazil) and CNPq (Brazil). A.M.K. thanks also RFBR (grant 01–01–00696) for partial support.

Figures captions

Fig. 1. Wave breaking of a box-like initial density distribution (shown by dashed line) in hydrodynamic approximation. The plot corresponds to the solution (6) with $`a=10,`$ $`\rho _0=1,`$ $`\rho _1=10,`$ $`t=1`$. Point $`A`$ propagates inward the box with local sound velocity $`2\sqrt{\rho _1}`$, point $`C`$ propagates outward along the background density with local sound velocity $`2\sqrt{\rho _0}`$, and point $`B`$ corresponds to intersection of two simple wave solutions with profiles $`AB`$ and $`CB`$.

Fig. 2. Formation of dissipationless shock wave after wave breaking point according to Whitham modulation theory applied to 1D NLS equation. Dashed line corresponds to multiple-valued region arising in hydrodynamic approximation given by Eq. (7) (it is analogous to the region between the points B and C in Fig. 1), and solid line represents modulated periodic wave given by Eqs. (8) and (10). Both profiles are calculated for $`\overline{\lambda }=10`$ at time $`t=1`$.

Fig. 3. Two-dimensional initial distribution of BEC density with paraboloid hump on constant background given by Eq. (18) with $`a=10,`$ $`\rho _0=1,`$ $`\rho _1=10`$.

Fig. 4. Two-dimensional density distribution of BEC after time $`t=2`$ of evolution from initial paraboloid density on constant background according to numerical solution of 2D NLS equation.

Fig. 5. Cross sections of density profile at different evolution time according to numerical solution of 2D NLS equation (17) with initial condition (18).