text:
# Stable quark stars beyond neutron stars: can they account for the missing matter ?

## 1 Introduction

Many authors have found solutions for strange quark stars with shells of hadronic composition . We investigate here the structure of a spherically symmetric stable dark ’star’ as outlined in the abstract. We use for the quark core the Tolman solution VI  with infinite central energy density and pressure. We also introduce a superfluid hydrogen shell which can be very large depending on its inner energy density. For normal hydrogen densities the hydrogen shell is very thin and negligible in mass. Objects with a mass of few solar masses can in principle be detected by gravitational lensing. However, in case of abnormally large initial energy density of the hydrogen shell, we arrive at very large masses of the star, which can not be easily detected by lensing. In this case, such stars are candidates for ’nonbaryonic’ dark matter. The three shells are discussed in the following section. For details of this analysis see reference .

## 2 The shells of the star

We consider 3 shells of radii r0 (quark shell), r1 (hadron shell), and r2 (hydrogen shell). The quark core of the star is described according to Tolman’s solution VI  and energy density and pressure are given by
$`\rho _e=\frac{3}{56\pi G_N}\frac{1}{r^2}`$ and $`p=\frac{1}{56\pi G_N}\frac{1}{r^2}\frac{(19(B/A)r)}{(1(B/A)r)}`$.
Where A, B are integration constants. We consider the core made up of u,d,s quarks and antileptons. allowing for the three chemical potentials associated with charge ($`\mu _Q`$), baryon number ($`\mu _B`$) and total lepton number ($`\mu _L`$). The hadronic phase is composed of baryons, leptons and antileptons. We use in the hadronic phase the vacuum pressure from . The TOV equations allow for a cosmological term, which in the hadronic phases must be chosen to cancel the positive vacuum pressure. It thus appears with negative sign in the deconfined core phase. The Gibbs equal pressure condition cannot be met, at zero temperature, neglecting a repulsive and very large second virial coefficient in the hadronic phase. Thus the equilibrium condition corresponds to equal energy density . The hydrodynamic and QCD equilibrium conditions between the two phases (hadronic and plasma) together with the condition of locally vanishing overall electric charge determine the core radius, mass and all chemical potentials for charge, baryon number and lepton number. We derive the mass $`M`$ and the radius $`r0`$ of the core: M(quark core) = 0.457 solar masses and r0(quark core) = 3.149 km. The baryon number of the core is $`N_{Bcore}0.287N_{Bsun}`$. Therefore the core itself is not bound, but stabilized by the outer hadron shell. We simplify the hadron shell composition structure by considering only neutrons, protons and electrons, though near the QCD transition zone also $`\mathrm{\Sigma }^{}`$, $`\mathrm{\Lambda }`$, $`\mu ^+`$, $`e^+`$ and all antineutrinos intervene. We obtain : $`M(\text{hadronic shell})=1.344\text{solar masses}`$ and $`r1(\text{hadronic shell})=9.047`$ km. Therefore the total mass of both phases inside the radius r1 is: $`M(\text{total inside r1})=1.801\text{solar masses}`$ and the baryon number inside r=r1 is 1.853 $`N_B(sun)`$. The nonrelativistic approximation of the TOV equation used here overestimates somewhat both mass and radius. We proceed to calculate the mass and radius of the hydrogen shell which starts at radius r1 and ends at radius r2. We use here the original approximate relations for superfluid hydrogen 111We do not include any further light elements, fermions or bosons, in the hydrogen shell., or alike boson-condensed non relativistic system, first described by N. N. Bogolubov . The results dependent on the parameter H, which is the initial energy density of hydrogen at r=r1, normalized to $`\alpha m_em_pm_H^2`$. The total radii, masses and baryon numbers are shown in table 1 together with H. The normal hydrogen shell corresponds to the first line.

## 3 The big annihilation

We follow the path of the early universe right after the QCD phase transition which is important for the genesis of quark stars. In particular fig. 1 shows the temperature as a function of the baryochemical potential. At about T $``$ 30 MeV the annihilation of matter and antimatter begins. We denote this ’the big annihilation’. Note that laboratory experiments measured a copious production of antibaryons at a temperature of $``$ 170 MeV see e.g. . The statistical fluctuations of net baryon number normalized to a volume of 1 fm3 following this path, are shown in fig. 2. They show maxima at the phase boundary and at temperatures below 10 MeV. The rest of the condensation of such a star needs large gravitational fluctuations too.

## 4 Conclusions

We investigated the structure of a spherically symmetric stable dark ’star’. For a low density hydrogen shell we find for the mass and the radius of the ’star’ $`\text{mass}=1.8\text{solar masses}`$ and $`\text{radius}=9.2\text{km}`$. For a high density hydrogen shell supporting densities of up to one atom / $`(10fm)^{3}`$ we find a large range of possible masses from 1.8 to 375 solar masses. The radii vary accordingly from 9 to 1200 km. Such objects are candidates for dark matter of any kind, despite their clearly baryonic nature, if they have been formed in such a way that they do not affect nucleosynthesis and the cosmic microwave background fluctuations at decoupling. We have calculated the path of the early universe after the QCD phase transition in the (T, $`\mu _B`$) plane, and the net baryon number fluctuations in the absence of gravitational fluctuations. We found that the statistical net baryon number fluctuations are large at the QCD phase transition at $`\mu _B0`$ and at T $``$ 10 MeV and less and at high $`\mu _B`$. These are regions where such stars could be formed. The mass dominating part of the hadron shell is very similar to a neutron star. This can be seen in our estimate of the mass : $`1.34`$ solar masses, slightly below the Chandrasekhar bound. The distinctive feature to a neutron star is that the hadron shell does not extend to the center of the ’star’. The quark core taken for itself and according to Tolmans universal solution, with $`\rho _er^2`$ towards the center is not stable. This means that it contains less nucleons counted by number than counted by mass. Yet this does not imply any instability of the ’star’ as a whole, as it would be the case for a one-phase neutron star. The phase boundary between quark core and hadron shell with its surface pressure stabilizes the ’star’. Despite genuine baryonic origin, besides the (anti)leptonic parts, such ’stars’ are indeed candidates for any kind of dark matter. If 1.8 solar masses indeed characterize the mass of these ’stars’ then they could be observed through gravitational microlensing of the galactic halo towards the Magellanic clouds . However masses much larger than few solar masses are hardly detectable by lensing. We also note, that the composition of the ’star’ considered here lacks heavy elements and hence the metallicity nearby is unusually low compared to burning stars of similar mass. They could possibly emit and absorb H lines, e.g. $`\text{Lyman}\alpha `$ lines taking into account gravitational redshift (which is z=0.56 for a mass of 1.8 solar masses) and high pressure compared to normal conditions. Another possibility is that they are part of a gravitating system including one or more shining stars.

Acknowledgments

We thank Dr. J. Schaffner-Bielich, Dr. R. von Steiger and Prof. K. Pretzl for discussions.


format: text

text:
# Electronic susceptibilities in systems with anisotropic Fermi surfaces

## I Introduction

The possible existence of quantum critical points in the phase diagrams of many materials has led to a detailed study of the low temperature behavior of the susceptibilities of electron systems. The critical properties of the system are determined by the energy and momentum dependence of the response function of the electron system associated to the order parameter in the ordered phaseH76 ; M85 ; M93 . It has been shown that the low temperature spin susceptibility of the isotropic electron liquid has an unexpected non analytic dependence on temperature, when high order perturbative corrections are consideredBKV97 . These corrections are irrelevant in the Renormalization Group senseS94 ; P92 ; MCC98 , and do not modify the basic properties of the electron liquid, as described by Landau’s theory. However, they can lead to unexpected power law dependences in many physical quantities at low temperatures, or change the order of the phase transitionsBKV99 . The origin of these non analiticities in homogeneous response functions has been traced back to the special properties of $`2k_F`$ scattering in the isotropic electron liquidCM01 .

It is well known that anisotropic Fermi surfaces can have regions where scattering becomes more singular than in the isotropic electron liquid, the so called “ hot spots ”. When two portions of the Fermi surface are flat and parallel, nesting occurs, and the susceptibilities diverge logarithmically, $`\mathrm{Re}\chi (\stackrel{}{𝐐},\omega )\mathrm{log}(\mathrm{\Lambda }/\omega )`$, where $`\stackrel{}{𝐐}`$ is the nesting vector. A saddle point in the density of states leads also to logarithmic divergences in two dimensions. The hot spots at the Fermi surface can be characterized by the frequency dependence of $`\mathrm{Im}\chi (\stackrel{}{𝐐},\omega )`$, where $`\stackrel{}{𝐐}`$ spans the hot spots. The usual behavior in a Fermi liquid is $`\mathrm{Im}\chi (\stackrel{}{𝐐},\omega )|\omega |`$, in any dimension $`D`$. For an isotropic Fermi surface, if $`|\stackrel{}{𝐐}|=2k_F`$, one has $`\mathrm{Im}\chi (\stackrel{}{𝐐},\omega )|\omega |^{(D1)/2}`$. For $`D=1`$ the imaginary part of the $`2k_F`$ susceptibility approaches a constant at low frequencies. By a Kramers Kronig transformation, it can be shown that the real part should diverge logarithmically, leading to the deviations from Landau’s theory which characterize a Luttinger liquid.

It is also possible to show that, when $`\stackrel{}{𝐐}`$ connects two saddle points in an anisotropic Fermi surface, $`\mathrm{Im}\chi (\stackrel{}{𝐐},\omega )|\omega |^{(D2)/2}`$. This result implies the existence of logarithmic divergences for $`D=2`$, which have been extensively studied in relation to high Tc superconductorsVH , and lead to deviations from Landau’s theoryGGV96 . In addition to saddle points, a generic anisotropic Fermi surface can show inflection points (see Fig.1). The existence of these points at the Fermi surface, which do not require any special fine tuning of the chemical potential, leads toGGV97  $`\mathrm{Im}\chi (\stackrel{}{𝐐},\omega )|\omega |^{(D2)/2+1/4}`$. For $`D=2`$, scattering between these points is more singular than the $`2k_F`$ scattering considered previously.

In the present work, we analyze scattering at inflection points in a two dimensional anisotropic Fermi surface. In the next section, we present the main features of the two loop calculation, extending the method used in referenceCM01 . The main results are obtained in section III, while the finer details of the calculation are deferred to the appendices. Applications to Fermi surfaces of different shapes are given in section IV, and section V discusses the main results of our work.

## II The method

We consider a system of two-dimensional (2D) fermions interacting through a generic short-ranged effective potential $`U(q)`$. For the sake of simplicity, we shall also assume that the interaction only affects electrons of opposite spins, which is a reasonable approximation when the momentum dependence of $`U(q)`$ is weak. It was shown in refs.BKV97 ; CM01  that while the lowest order ($`U`$) perturbative corrections are well behaved, higher order corrections can lead to an anomalous behaviour in the low-energy properties of the Fermi liquid. To be more precise, the uniform spin susceptibility of a 2D electron system shows a linear $`T`$ dependence, which contradicts the usual Sommerfeld expansion in powers of $`(T/E_F)^2`$. Such anomalous behaviour was traced back to the peculiarities of $`2k_F`$ scattering, i.e. the occurrence of particle-hole pairs lying on opposite sides of the Fermi surface (FS). This special wavevector plays a key role in the $`q`$-dependent susceptibility of electronic systems already in the non-interacting case, with the appearance of a square-root singularity around $`2k_F`$ which is directly related to the jump in the occupation number. If one considers the uniform susceptibility, though, the singularities associated with $`2k_F`$ scattering can only show up indirectly through the excitation of a virtual particle-hole pair, which explains the absence of anomalous corrections at lowest order in the interaction strength.

In the general (non circular) case, among all the wavevectors connecting opposite sides of the FS, the inflection points play a special role, due to the flatness of the Fermi surface (the extreme case being the one of a perfectly flat FS, or perfect nesting, which leads to strong instabilities). According to the previous discussion, the second-order diagrams which lead to non-standard behavior are the ones containing a particle-hole bubble whose transferred momentum can match the special value $`\stackrel{~}{𝐐}`$. Such diagrams are depicted in figure 2.

Diagrams a) and b) are vertex corrections. They have opposite signs and cancel in the case of a perfectly $`q`$-independent interaction potential: the fermion propagators involved are the same in both diagrams, the only difference being in the momentum carried by the interaction. To be specific, with the notations of figure 2, diagram a) is proportional to $`U(\stackrel{~}{𝐐})^2`$, while diagram b) involves some momentum average of the interaction, and there is no reason for a perfect cancellation in the general case. Diagram c) is a self-energy correction, and will be considered separately.

## III Uniform susceptibilities

After integration over Matsubara frequencies, the vertex correction a) of figure 2 in the zero-frequency, zero-momentum limit reads:

$$\chi (T)=\frac{d^2p}{(2\pi )^2}\frac{d^2Q}{(2\pi )^2}U(𝐐)^2\mathrm{\Delta }(\xi _𝐩)\mathrm{\Delta }(\xi _{𝐩+𝐐})L(\xi _{𝐩+𝐐}\xi _𝐩,𝐐)$$
(1)
where the Lindhard function in 2 space dimensions is defined as

$$L(i\mathrm{\Omega },Q)=\frac{d^2k}{(2\pi )^2}\frac{n(\xi _𝐤)n(\xi _{𝐤+𝐐})}{i\mathrm{\Omega }+\xi _𝐤\xi _{𝐤+𝐐}}$$
(2)
and 111This is the finite temperature generalization, within our perturbative treatment, of the well-known Fermi Liquid formula expressing the singularity of $`G(k)G(k+q)`$ lines in the vertex function, see e.g. A.A. Abrikosov, L.P. Gorkov and I.E. Dzyaloshinski, Methods of Quantum Field Theory in Statistical Physics (Dover, New York, 1975). $`\mathrm{\Delta }(\xi )=\beta /4\mathrm{cosh}^2(\beta \xi /2)`$. Such $`\mathrm{\Delta }`$-functions, which constrain momenta to lie within a shell of thickness $`T`$ from the Fermi surface, are strongly temperature dependent, and they are responsible for the leading temperature dependence of the susceptibility. Indeed, in (1) we have omitted terms proportional to $`\mathrm{\Delta }(\xi _{𝐩+𝐐})n(\xi _𝐩)`$ coming from the poles of the Lindhard function, which are less $`T`$-dependent since they receive contributions mainly from regions far from the Fermi surface. Taking advantage of time-reversal symmetry ($`\xi _𝐤=\xi _𝐤`$), we can write

$$L(\mathrm{\Delta }\xi _p,Q)=\frac{d^2k}{(2\pi )^2}n(\xi _𝐤)\left[\frac{1}{\mathrm{\Delta }\xi _𝐩\mathrm{\Delta }\xi _𝐤}\frac{1}{\mathrm{\Delta }\xi _𝐩+\mathrm{\Delta }\xi _𝐤}\right]$$
(3)
where we have defined $`\mathrm{\Delta }\xi _𝐤=\xi _{𝐤+𝐐}\xi _𝐤`$. As was pointed out in the previous section, the most singular contributions to (1) come from regions where the momentum $`𝐐`$ flowing through the Lindhard function $`L`$ connects parts of the FS which are almost parallel, since this makes the denominators in eq. (3) small on large regions of $`k`$-space. Otherwise stated, the scattering processes taking place within a particle-hole pair are enhanced around special wavevectors $`\stackrel{~}{𝐐}`$ due to the peculiar geometry of the FS. In the case of a spherically symmetric FS, any momentum $`\stackrel{~}{𝐐}`$ of modulus $`2k_F`$ is a source of enhanced scattering, but the deviation from parallelicity is quadratic as we move in the direction tangent to the surface (see fig. 1, left). More singular is the case of inflection points occurring when the curvature of the FS vanishes, leading to a cubic, or even quartic dispersion (see figure 1, right), which is a quite generic phenomenon when dealing with electrons on a lattice.

In the next subsections, we shall present the calculation of the leading $`T`$-dependence of the diagram a) in the simplest circular case as well as for more complex FS shapes. The result for diagram b) can obtained by replacing $`U(𝐐)U(𝐩𝐤)`$ in eq. (1). This can only lead to a change in the prefactors, but will not alter the leading temperature dependence of the susceptibility. The self-energy diagram c) has a different structure, and will be analysed at the end of the section.

### III.1 Isotropic Fermi surface

By choosing an appropriate coordinate system, the dispersion relation around any point on a spherical FS (and, generically, about non-special points of an anisotropic FS) can be expanded as:

$$\xi _𝐤/v=k_y+ak_x^2$$
(4)
$`v`$ being the Fermi velocity at that particular point (that we shall identify as $`\stackrel{~}{𝐐}/2`$), and $`a>0`$ being related to the FS curvature. The above expression is assumed to be valid up to a momentum cutoff $`\mathrm{\Lambda }`$ which is larger than the one imposed by the finite temperature $`\mathrm{\Delta }`$-functions 222The momentum cutoff should in principle be different in the $`x`$ and $`y`$ directions, but this is of no practical relevance here.. Introducing $`𝐪=𝐐\stackrel{~}{𝐐}`$, we can write

$$\xi _{𝐤+𝐐}/v=(k_y+q_y)+a(k_x+q_x)^2$$
(5)
were we have used the property $`\xi _{𝐤+\stackrel{~}{𝐐}}=\xi _𝐤`$ (reflection symmetry). We now change variables to

$`k_{}`$ $`=`$ $`k_y+ak_x^2`$ (6)
$`k_{}`$ $`=`$ $`k_x+q_x/2`$
such that $`\xi _𝐤=vk_{}`$ and

$$\xi _{𝐤+𝐐}/v=k_{}q_{}+2ak_{}^2+(3/2)aq_{}^2$$
(7)
Apart from the constant shift $`q_x/2`$ introduced for later convenience, this new coordinate system is locally equivalent to polar coordinates, which would be the natural choice when dealing with a perfectly symmetric FS. We shall focus on the first term in (3), which is independent of $`q_{}`$, and therefore turns out to be the most singular. Omitting unimportant multiplicative factors, we have for the real part of $`L`$:

$`L`$ $`=`$ $`{\displaystyle \frac{1}{v}}{\displaystyle 𝑑k_{}n(vk_{})_\mathrm{\Lambda }^\mathrm{\Lambda }𝑑k_{}\frac{1}{A2k_{}+Bk_{}^2}}`$ (8)
$`=`$ $`{\displaystyle \frac{1}{v\sqrt{B}}}{\displaystyle 𝑑k_{}n(vk_{})\frac{\theta (A2k_{})}{\sqrt{A2k_{}}}}`$
with $`A=2p_{}2ap_{}^2`$ and $`B=2a`$ (the $`\theta `$-function ensures that the integral is real). In the second term of eq. (8) we have performed the $`k_{}`$ integration by pushing the momentum cutoff to infinity. The main point is that the former expression can now be integrated by parts to give a further $`\mathrm{\Delta }`$ constraint on $`k_{}`$ CM01 :

$$L=\frac{1}{\sqrt{B}}_\mathrm{\Lambda }^{A/2}𝑑k_{}\sqrt{A2k_{}}\mathrm{\Delta }(vk_{})+\mathrm{}$$
(9)
where the ellipsis stands for terms which are not confined to the region near the FS. By inspection of the results for $`AT/v`$, $`AT/v`$ and $`A0`$ respectively, we see that the $`\mathrm{\Delta }(vk_{})`$ function behaves qualitatively as a $`\delta (vk_{}+T)`$. Therefore, to study the temperature dependence of the susceptibility we can replace the previous expression by

$$L\frac{\sqrt{A/2+T/v}}{v\sqrt{a}}$$
(10)
where there is an implicit $`\theta `$-function of the argument of the square root. We are left with a tractable expression for the real part of the Lindhard function, that we shall use to evaluate the 2-loop diagram of figure 2.a.

We now perform the remaining integrals in (1) in the following order: $`dq_{},dp_{}`$ then $`dp_{}`$ and $`dq_{}`$. The first integral is trivial, since $`q_{}`$ only enters in $`\mathrm{\Delta }(\xi _{𝐩+𝐐})`$. Moreover, $`\xi _{𝐩+𝐐}`$ is linear in $`q_{}`$ (cf. eq. (7)) so that the integration just gives $`1/v`$. The $`p_{}`$ integral can also be performed straightforwardly, and we are left with an expression of the form

$`\chi `$ $``$ $`{\displaystyle \frac{\stackrel{~}{U}^2}{v^3\sqrt{a}}}{\displaystyle 𝑑q_{}𝑑p_{}\sqrt{T/vap_{}^2}}`$
$``$ $`{\displaystyle \frac{\mathrm{\Lambda }\stackrel{~}{U}^2}{v^3\sqrt{a}}}{\displaystyle 𝑑p_{}\sqrt{T/vap_{}^2}}{\displaystyle \frac{\mathrm{\Lambda }\stackrel{~}{U}^2}{v^4a}}T`$
where again we have neglected unimportant multiplicative factors and we have defined $`\stackrel{~}{U}=U(\stackrel{~}{𝐐})`$. Within our treatment, we have recovered the result that the spin susceptibility of an isotropic 2D Fermi liquid is intrinsically linear in temperatureCM01 . For a circular FS, this can be written as

$$\chi (T)=\chi _0+\chi _1T$$
(12)
Incidentally, our calculation suggests that the low-temperature correction to the susceptibility is positive, in agreement with refs. M99 ; HT98 ; CM01 .

### III.2 Anisotropic FS with inflection points

In the vicinity of an inflection point, the dispersion relation can be written as

$$\xi _k/v=k_ybk_x^3+gk_x^4$$
(13)
where $`b`$ and $`g`$ can be chosen to be positive. A change of variables similar to eq. (6) of the previous section leads to

$$L=\frac{1}{v}𝑑k_{}n(vk_{})𝑑k_{}\frac{1}{A2k_{}+Bk_{}^2+Ck_{}^4}$$
(14)
with $`A=2p_{}3bq_{}p_{}^22gp_{}^4`$, $`B=3bq_{}`$ and $`C=2g`$. The $`k_{}`$ integral can be rewritten in the form

$$I=\frac{1}{C^{1/4}(A2k_{})^{3/4}}\frac{dx}{1+\alpha x^2+x^4}$$
(15)
with $`\alpha =(B/\sqrt{C})/\sqrt{A2k_{}}`$. We shall be interested in the region close to the edge ($`k_{}A/2`$), where $`\alpha `$ is large and positive. We can then drop the quadratic term and perform the integration:

$$I=\frac{1}{[bq_{}(A2k_{})]^{1/2}}$$
(16)
This can be integrated by parts in $`dk_{}`$ to give

$$L=\frac{1}{\sqrt{bq_{}}}_{A/2b^2q_{}^2/g}^{A/2}𝑑k_{}\sqrt{A2k_{}}\mathrm{\Delta }(vk_{})+\mathrm{}$$
(17)
where the ellipsis stands for a term which is not confined close to the FS (the limits of integration account for the condition $`\alpha 1`$). Provided that the $`\mathrm{\Delta }`$ function lies entirely inside the domain of integration, i.e.

$$q_{}>q_{min}=\sqrt{\frac{gT}{b^2v}}$$
(18)
the result takes the form

$$L\frac{\sqrt{A/2+T/v}}{v\sqrt{bq_{}}}$$
(19)
The region of phase space we have just identified is the one which gives the leading temperature dependence in the susceptibility. Indeed, for $`k_{}`$ outside the range of integration considered above (implying $`\alpha 1`$), the result of the integral (15) is $`I(A2k_{})^{3/4}`$ instead of eq. (16), leading to a weaker (linear) temperature dependence in the final result. The same holds if we consider a negative $`q_{}`$ ($`\alpha <0`$).

The calculation now follows the same lines as in the previous case. The $`q_{}`$ integration yields a factor $`1/v`$, and the $`p_{}`$ integration can be performed by replacing $`\mathrm{\Delta }(vp_{})\delta (p_{}T/v)/v`$, which gives

$`\chi `$ $``$ $`\stackrel{~}{U}^2{\displaystyle _{q_{min}}^\mathrm{\Lambda }}𝑑q_{}{\displaystyle _\mathrm{\Lambda }^\mathrm{\Lambda }}𝑑p_{}{\displaystyle \frac{\sqrt{T/vbq_{}p_{}^2gp_{}^4}}{v^3\sqrt{bq_{}}}}`$ (20)
$``$ $`{\displaystyle \frac{\stackrel{~}{U}^2T}{bv^4}}{\displaystyle _{q_{min}}^\mathrm{\Lambda }}{\displaystyle \frac{dq_{}}{q_{}}}={\displaystyle \frac{\stackrel{~}{U}^2}{bv^4}}T\mathrm{log}\left({\displaystyle \frac{gT}{vb^2\mathrm{\Lambda }^2}}\right)`$
Taking into account the scattering from the regions of the FS far from the inflection points, whose behavior is given by eq. (12), the susceptibility reads

$$\chi (T)=\chi _0+\chi _1T\chi _1^{}T\mathrm{log}T$$
(21)
Once again, the sign of the correction is such that the susceptibility increases with temperature. However, the contribution coming from the other diagram b) has opposite sign. As a rule of thumb, one can argue that the overall vertex-correction is positive if the effective interaction is peaked around $`\stackrel{~}{𝐐}`$ and negative otherwise (it vanishes when the momentum dependence of $`U(q)`$ is flat, since in that case the two diagrams perfectly cancel).

### III.3 Special inflection points

The previous analysis assumes the existence of a generic inflection point along the Fermi surface. This is a situation which can be achieved, in an anisotropic system, for a finite range of values of the filling or the chemical potential. These points are characterized by the absence of a quadratic term in the expansion of the dispersion relation around the Fermi surface presented in eq. (13). For certain values of the parameters, however, which require a fine tuning of the filling or the chemical potential, the cubic coefficient, $`b`$, or the quartic one, $`g`$, in eq. (13) can be zero as well. Two such situations are schematically shown in Fig.3.

We first consider the case when the cubic term in the dispersion relation parallel to the Fermi surface vanishes ($`b=0`$), which is realized in the $`tt^{}`$ Hubbard model in a square lattice (see point $`G`$ in the top panel of figure 3), or in simple tight binding models on the triangular lattice, for instance. The susceptibility becomes more anomalous than in the generic case discussed previously, as can be seen by letting $`b0`$ in eq. (20). One has respectively $`A=2p_{}3gq_{}^2p_{}^22gp_{}^4`$ and $`B=3gq_{}^2`$. The condition $`\alpha 1`$ now corresponds to $`A2k_{}gq_{}^4`$, which modifies the limits of integration in eq. (17). Repeating the same arguments as before with $`q_{min}=(T/gv)^{1/4}`$, we obtain

$$L=\frac{\sqrt{A/2+T/v}}{v\sqrt{g}q_{}}$$
(22)
leading to

$$\chi =\frac{\stackrel{~}{U}^2}{v^3}\left(\frac{T}{vg}\right)^{3/4}\mathrm{log}\left(\frac{T}{gv\mathrm{\Lambda }^4}\right)$$
(23)
The other possibility is that the quartic term vanishes ($`g=0`$), which can occur in a tight-binding model with orthorombic symmetry, considering two different hopping parameters $`t_at_b`$ (see point $`D`$ in the bottom panel of figure 3). In that case, however, not only $`g`$ but all the even coefficients in the dispersion relation vanish. This leads to perfect nesting between opposite branches of the Fermi surface, giving rise to a much more singular behavior $`\chi (T)\mathrm{log}T`$.

### III.4 Self-energy correction

After integration over Matsubara frequencies, which now requires some more attention due to the presence of two fermion lines of equal argument (a double pole in the complex-plane integrals), the anomalous part of the self-energy correction c) of figure 2 can be reduced to the form

$$\chi \stackrel{~}{U}^2Td^2Qd^2pd^2k\frac{\mathrm{\Delta }(\xi _𝐤)\mathrm{\Delta }(\xi _𝐩)}{(\mathrm{\Delta }\xi _𝐩\mathrm{\Delta }\xi _𝐤)^2}$$
(24)
with $`\mathrm{\Delta }\xi _𝐤=\xi _{𝐤+𝐐}\xi _𝐤`$. The $`Q`$-integration is now restricted to the region close to (within $`T/v`$ of) $`\stackrel{~}{𝐐}`$. We shall not go through all the calculations of the self-energy diagram, which can be performed following the same lines of the previous sections. The results for the temperature dependence are analogous to those given by equations (12), (21) and (23). This can be understood by noting that although the denominator in eq. (24) is more singular than the one of eq. (3), the additional anomalies that it carries with it are cancelled by the explicit $`T`$ factor in front of eq. (24), leading to the same temperature dependence as the vertex corerction. Its sign is also the same as the vertex diagram a).

The results of this section are summarized in table 1.

## IV Examples

### IV.1 Superconducting cuprates

It is often assumed that a tight-binding model on a square lattice with nearest neighbor ($`t`$) and next-nearest neighbor ($`t^{}`$) hopping reproduces well the band structure of the layered cuprates:

$$\epsilon (𝐤)=2t(\mathrm{cos}k_x+\mathrm{cos}k_y)4t^{}\mathrm{cos}k_x\mathrm{cos}k_y$$
(25)
where $`t^{}/t0.25`$. This case corresponds roughly to the top set of Fermi surfaces in fig.3. The dispersion relation above has a saddle point at a doping $`\delta _{VHS}`$ corresponding to a chemical potential is $`E_{VHS}=4|t^{}|`$. The curvature of the Fermi surface along the diagonals becomes negative at a higher doping $`\delta _c`$, where the chemical potential is $`E_c=8|t^{}|+16\frac{|t^{}|^3}{t^2}`$. For fillings such that $`E_cE_FE_{VHS}`$, the Fermi surface has 8 inflection points. From these values, and the previous analysis, one can obtain a qualitative picture of the temperature dependence of the susceptibilities, when the electron density is in this range:

i ) For $`T|t^{}|`$, the susceptibility is determined by $`t`$ only. As the doping is close to half filling, we expect $`\chi (T)|T|^0`$, the result for perfect nesting.

ii) For $`T|t^{}|`$ and $`T|E_FE_{VHS}|`$, the same behavior as in i) should be observed.

iii) For $`T|t^{}|`$ and $`T|E_FE_c|`$, the susceptibility is dominated by contribution from the area near the special point discussed in section IIIC. Hence, $`\chi (T)|T|^{3/4}\mathrm{log}T`$.

iv) For $`T|t^{}|`$ and $`T|E_FE_c|`$, and $`T|E_FE_{VHS}|`$, the contributions from the saddle point and from the special point in the previous paragraph are absent. Thus, $`\chi (T)|T|\mathrm{log}(T)`$, because of the presence of the inflection points.

We can make the estimates of the crossover region in the $`T`$doping plane more precise from the doping dependence of the coefficient of the cubic term $`b`$ in eq.(13). Expanding around the saddle point, we obtain $`b|E_FE_{VHS}|`$. Hence, the crossover between regions ii) and iv) takes place at a temperature $`T^{}|E_FE_{VHS}|`$. Performing a similar calculation around the situation $`E_F=E_c`$, we have $`b\sqrt{|E_FE_c|}`$, so that the crossover temperature is $`T^{}(E_cE_F)^2`$. For fillings $`E_FE_c`$ but with no inflection points in the Fermi surface, we obtain a crossover to the $`\chi (T)|T|`$ behavior due to $`2k_F`$ scattering, with a crossover temperature $`T^{}|E_FE_c|`$. The different regimes are schematically shown in Fig.4.

Taking realistic numbers for the dispersion relation, our analysis predicts anomalous low temperature behavior in all the region between $`\delta _{VHS}`$ and $`\delta _c`$, corresponding to the strongly overdoped region which is experimentally accessible. This shows that non-standard behaviour of the physical properties should be expected even in a regime which is usually believed to be well described by normal Fermi liquid theory.

### IV.2 Quasi-1D organic compounds

Organic conductors are often very anisotropic due to the planar structure of their molecules. For example, the salts of the family (TM)2X (TM=TMTTF,TMTSF and X=inorganic anion) are all isostructural and can be viewed as two-dimensional arrays of weakly coupled 1D chains, since the electronic overlaps in the transverse direction are 10 times smaller than in the chain direction (the transfer integrals in the third direction are 500 times smaller, and can be neglected, see for instance ref.BJ99 ). The band structure is well represented as:

$$\epsilon (𝐤)=2t_a\mathrm{cos}(k_aa)2t_b\mathrm{cos}(k_bb)$$
(26)
assuming an orthorombic structure with lattice parameters $`b2a`$. This case corresponds to the bottom set of Fermi surfaces in fig. 3. The parameter $`t_b1030meV`$ sets the scale below which the FS is modulated in the $`b`$ direction, so that the predicted enhancement of susceptibilities due to inflection points should be observable at and below room temperature. The value of the anisotropy ratio $`\tau =t_b/t_a`$ is large enough to ensure that the system is well described by a Fermi liquid picture down to very low temperatures. The filling factor $`\rho `$ is fixed by charge transfer and varies from compound to compound, ranging from $`1/2`$ to $`1`$ hole per TM site. The Fermi surface has two Van Hove singularities at $`E_F=\pm E_{VHS}=\pm 2t(1\tau )`$, and 4 inflection points in all the interval $`0<|E_F|<E_{VHS}`$ Taking $`\tau =0.1`$, this corresponds to the region of fillings $`0.3<\rho <1.7`$. In the absence of higher harmonics in eq.(26), $`E_F=0`$ corresponds to half-filling ($`\rho =1`$), and the Fermi surface has perfect nesting, as $`\epsilon (𝐤)=\epsilon (𝐤+𝐐)`$, where $`𝐐=(\pi ,\pi )`$ (hopping between more distant neighbors will suppress this effect). The points in the phase diagram where the topology of the Fermi surface changes, leading to different behaviors of the electronic susceptibility, are sketched in figure 5.

In the (TM)2X compounds, the spin susceptibility shows a sizeable increase in the metallic phase up to room temperature (see e.g. fig. 8 of reference BJ99 ), which cannot be explained by “standard” Fermi liquid theory (the latter predicts variations on the scale of the Fermi temperature). On the other hand, the presence of enhanced scattering close to inflection points could well be the underlying mechanism of this anomalous temperature dependence, and should be taken into account when studying the low temperature phase transitions of such compounds.

Following the same procedure used in the previous subsection, the electron susceptibility will undergo a succession of crossovers upon varying the filling, which can be achieved either by anion substitution or by applying pressure to the samples.

## V Conclusions

We have analyzed the corrections to Fermi liquid behavior in anisotropic interacting electronic systems in two dimensions, which arise from the existence of points in the Fermi surface where scattering is enhanced. Besides the extensively studied case of a saddle point, we have analyzed in detail the influence of inflection points, which do not require any special fine tuning of the chemical potential or the filling. The presence of these points enhance the anomalous dependence on temperature which arise from $`2k_F`$ scattering in isotropic Fermi surfacesCM01 . We find that the corrections which were linear in $`|T|`$ change into $`|T|\mathrm{log}|T|`$. The absence of symmetries also implies the lack of cancellation between different diagrams, so that these anomalies should be observed in both the spin and charge susceptibilities.

For special fillings, more singular behavior is expected. In the case of systems with tetragonal or hexagonal symmetry, when the Fermi surface is close to these fillings, the corrections to the susceptibilities go as $`|T|^{3/4}\mathrm{log}|T|`$, showing that the existence of non integer $`T`$ dependences does not need to violate Landau’s model for the low energy excitations of a Fermi liquid.

We have also discussed the possible crossovers between the different regimes analyzed, and the experimental consequences that they may lead to. In the case of the superconducting cuprates, anomalous susceptibilities should appear in the strongly overdoped region, above the doping $`\delta _{VHS}`$ characterized by Van Hove singularities in the density of states. On the other hand, all of the organic conductors of the family (TM)2X should fall in the region of fillings where anomalous corrections to the susceptibility are important. Of course, the analysis presented here should also apply to other classes of quasi two-dimensional systems (heavy fermion materials, Sr2RuO4, electrically doped 2D organic films, other organic conductors …).

Finally, let us point out that the breakdown of the Sommerfeld expansion for the spin susceptibility suggests that the free energy $``$ itself has a non analytic dependence on $`T`$, once that high order interactions ($`2k_F`$ scattering) are taken into account. If, as was proposeded inM88 ; M99 ; M01 , and numerically verified in HK99 , the role of temperature and magnetic field is interchangeable in the functional form of $`(T,H)`$, one can conclude that the anomalous $`T`$-dependences calculated here for the susceptibility are also expected in the specific-heat coefficient $`\gamma =C/T`$.

## VI Acknowledgements

We are thankful to R. Markiewicz and M. A. H. Vozmediano for helpful discussions. This work was financially supported by MEC (Spain) through grant PB96-0875, and the European Union through grant FMRXCT980183.

## Appendix A Inflection points in the $`tt^{}`$ model

We shall determine here the parameters of the dispersion relation (13) in the case of a tight-binding model on a square lattice with nearest ($`t`$) and next-nearest ($`t^{}`$) neighbor hopping. Let us focus to the doping levels close to $`\delta _c`$, the point where the inflection points of the Fermi surface merge in pairs on the diagonals of the BZ, leading to the most singular corrections to the susceptibility. It is then natural to rewrite the dispersion relation in a basis rotated by $`45^o`$:

$$\xi =4t(\mathrm{cos}p\mathrm{cos}q)+4t^{}(\mathrm{cos}^2p\mathrm{sin}^2q)E_F$$
(27)
where $`p=(k_x+k_y)/2`$ and $`q=(k_xk_y)/2`$. The Fermi surface crosses the diagonal ($`q=0`$) at a momentum $`p_F`$ given by $`E_F=4t\mathrm{cos}p_F+4t^{}\mathrm{cos}^2p_F`$. The dispersion relation can then be expanded as

$$\xi =A(pp_F)+Bq^2+Cq^4$$
(28)
with $`A=4[t\mathrm{sin}p_Ft^{}\mathrm{sin}2p_F]`$, $`B=2t\mathrm{cos}p_F4t^{}`$ and $`C=(4/3)[t^{}(t/8)\mathrm{cos}p_F]`$. The topology of the Fermi surface changes at two well-defined doping levels:

* the curvature changes sign at a doping $`\delta =\delta _c`$ given by the condition $`B=0`$. The corresponding Fermi energy is $`E_c=8t^{}+16(t^{})^3/t^2`$ and the coordinates of the inflection point are $`(p_c,q_c)=(\mathrm{arccos}2t^{}/t,0)`$, corresponding to the point $`G`$ of figure 3;
* Van Hove singularities arise at a doping $`\delta _{VHS}`$ given by $`E_{VHS}=4t^{}`$ (M-points in figure 3).

Inflection points appear in all the region of dopings $`\delta _{VHS}<\delta <\delta _c`$, following the dashed curve of figure 3 (top panel). The consequences on the physical properties of the system are summarized in figure 4.

#### Dispersion around inflection points.

The equation of the Fermi surface is $`\xi =0`$, which implicitly defines a function $`p=p(q)`$. Putting the second derivative $`p^{\prime \prime }(q)=0`$ yields the coordinates $`(p_0,q_0)`$ of the inflection points. For $`\delta \delta _c`$, setting $`u=14(t^{}/t)^2`$, we can write

$$p_0=p_c+\frac{E_FE_c}{4tu^{3/2}};q_0=\left(\frac{E_FE_c}{12ut^{}}\right)^{1/2}$$
(29)
so that the trajectory of the inflection points is parabolic around $`G`$. By expanding around $`(p_0,q_0)`$, we obtain an equation of the form (13) with

$$v=4tu^{3/2};b=\frac{t^{}}{t}\left(\frac{E_FE_c}{12t^{}}\right)^{1/2};g=\frac{t^{}}{4tu^{3/2}}$$
(30)
## Appendix B Inflection points in the $`t_at_b`$ model

We shall now derive the parameters of eq. (13) for a tight binding model on an orthorombic lattice, with anisotropic hopping ($`\tau =t_b/t_a1`$). Let us rewrite for simplicity the dispersion relation (26) as

$$\xi =2t_a[\mathrm{cos}k+\tau \mathrm{cos}p+\nu ]$$
(31)
with $`k=k_aa`$, $`p=k_bb`$ and $`\nu =E_F/2t`$. The equation of the Fermi surface is $`k=\mathrm{arccos}(\nu \tau \mathrm{cos}p)`$. The number of electrons per site is approximately given by $`\rho =2\pi ^1\mathrm{arccos}(\nu )`$. The Fermi surface has 2 Van Hove singularities at $`E_{VHS}=\pm 2t(1\tau )`$, and 4 inflection points for any $`0<|E_F|<E_{VHS}`$. At half filling ($`E_F=0`$), the two branches of the open Fermi surface are perfectly nested. The physical consequences of the changes in the Fermi surface topology occurring at those special fillings are sketched in figure 5.

#### Dispersion around inflection points.

By setting $`k^{\prime \prime }(p)=0`$ we find that the inflection points are located at

$$p_0=\mathrm{arccos}\left(\frac{\tau \nu }{1\nu ^2}\right);k_0=\mathrm{arccos}\left(\nu \frac{\tau ^2\nu }{1\nu ^2}\right)$$
(32)
For filling levels close to $`\nu =0`$ (half filling), the location of the inflection points describes a straight line of slope $`\tau /\sqrt{1\nu ^2}`$ in the $`(k,p)`$ plane (see bottom panel of figure 3). After a rotation of the coordinate axes, we obtain an equation of the form (13) with

$$v=2t\sqrt{1\nu ^2};b=\frac{\tau }{6\sqrt{1\nu ^2}};g=\frac{\nu \tau ^2}{24(1\nu ^2)^{3/2}}$$
(33)
The approach to the perfect nesting situation at half filling is signalled by a vanishing $`g`$, the coefficient of the quartic term in the dispersion relation.


text:
# On the theory of magnetic impurities in integrable correlated electron chains

(April 15, 2002)

In a recent preprint  Ge et al. reply to our recent response  to some concerns raised in . Several obvious remarks are necessary to clarify the situation, because the authors of  certainly misunderstand (or misinterpret) some of our statements.

1. In  and in  the authors imply that integrable magnetic impurities cannot exist in closed $`tJ`$ and Hubbard chains (cf. P. 795 of  and P. 8544 of ). It turns out that in  they admit that integrable magnetic impurities can exist in exactly solvable closed correlated electron chains.

2. The approach (i) in our answer  has been named the “quantum inverse scattering method” (also known as the algebraic Bethe ansatz cf. ) long before our work on impurities in correlated electron systems, see, e.g., . This standard definition has been used in our reply . The approach (ii) in our reply  is known as the “graded quantum inverse scattering method” . Does the strong emphasis in  on the difference between the co-ordinate Bethe ansatz and the quantum inverse scattering method imply that the authors of  believe that those two methods could yield different answers? We are not aware of any such contradictions between the two methods.

3. In  we pointed out that the impurity matrix changes the commutation relations in the spin sector (see below). This is absolutely correct, keeping in mind that two parameters, $`\theta `$ and $`S`$, which distinguish the impurity site from other sites of the chain are nonzero (note that the impurity scattering matrix used in our papers  mixes the states with $`S`$ and $`S+\frac{1}{2}`$; this hybridization is sometimes misunderstood).

4. The magnetic impurities we studied in our papers have an essentially different structure than those of ; hence, it is no wonder that the solutions do not coincide with ours. Ref.  considers the supersymmetric $`tJ`$ model with a different grading than the one considered by us and does not consider magnetic impurities. Since this represents a very different situation, it does not contradict our results. Actually, the special case of the impurity of  coincides with our results .

5. In the approach (ii) for the supersymmetric $`tJ`$ model the operators $`\widehat{A}_{12}`$, $`\widehat{A}_{13}`$, $`\widehat{A}_{21}`$ and $`\widehat{A}_{23}`$ acting on the vacuum state do indeed yield zero (cf. Eq. (3.27) of ) in the FFB grading, contrary to the statements in . One can see that the results of  for the special case of $`\theta =0`$ and spin, equal to the ones of the host, coincide with those of  (cf. Eq. (3.50) of  and (A1) of ). The operators $`\widehat{A}_{12}`$ and $`\widehat{A}_{21}`$ in  do not contain any characteristics of the impurity (i.e., $`\theta `$ and $`S`$), and are then equivalent to those studied in . This way, the argumentation of  can be applied to paper , and the criticism presented in  actually concerns  rather than  (in which we essentially used the method developed in ). However, the criticism presented in  is incorrect, because the authors do not take into account the fact that the eigenvalue of the transfer matrix is determined up to some multiplier . Moreover, the important commutation relations for the spin sector are those between the $`\widehat{A}_{ij}`$ ($`i,j=1,2`$) and $`\widehat{A}_{31}`$ and $`\widehat{A}_{32}`$ (or $`C_{1,2}`$, cf. ), and those between the latter two operators, which indeed are used as “creation operators” in . Namely, the changes in these commutation relations, but not in those between $`\widehat{A}_{12}`$, $`\widehat{A}_{13}`$, $`\widehat{A}_{21}`$ and $`\widehat{A}_{23}`$ (which are mentioned in ), determine the changes in the spin sector of Bethe ansatz equations due to the magnetic impurity.

6. The change of the “class of the representation” ($`l`$) implies the change of the symmetry in the considered model (we did not discuss the symmetry of the Lax operator in , however, it turns out that the symmetries of our impurity $`L`$-operators and those of the host are the same, unlike the case of Refs. ). Hence, our statement in response to  is correct.

7. Point (5) of our answer to  pertains to approach (i), but not to approach (ii). However, in  the approach (ii) was used. It is, naturally, correct  that in the FFB grading of the approach (ii) one cannot use $`\widehat{A}_{21}`$ as a “raising operator”. But the authors of Ref.  misunderstand our statements  and incorrectly mix the two approaches.

8. Obviously, the statements of our answer  (and the results of our previous papers) do not contradict .

9. The claim in  that Ge et al. studied a spin impurity, without additional charge degrees of freedom, contradicts the fact that to according their Bethe ansatz equations, e.g., derived in , the valence (the occupation number at the impurity site) varies with external parameters (such as the chemical potential, a global (non-local) magnetic field), cf. . This is impossible if one studies a pure magnetic impurity, which has only spin degrees of freedom (in this case the valence should be one and not vary with the external parameters, even for $`q=1`$).

Discussions with P. Schlottmann are acknowledged.


text:
# 1 Introduction

## 1 Introduction

Analysis of various aspects of spontaneous emission by a system of two-level atoms, initiated in the classical paper of Dicke  was further developed by several authors (see e.g. ). In particular, in the case of spontaneous emission by two trapped atoms separated by a distance small compared to the radiation wavelength, where is a substantial probability that a photon emitted by one atom will be absorbed by the other, there are states of the system in which photon exchange can enhance or diminish spontaneous decay rates. The states with enhanced decay rate are called superradiant and analogously states with diminished decay rate are called subradiant . It was also shown by Dicke, that the system of two coupled two-level atoms can be treated as a single four-level system with modified decay rates. Note also that such model can be realized in a laboratory by two laser-cooled trapped ions, where the observation of superradiance and subradiance is possible .

Another aspects of the model of the spontaneous emission are studied in the present paper. When the compound system of two atoms is in an entangled state, the irreversible process of radiative decay usually destroys correlations and the state becomes unentangled. In the model studied here, the photon exchange produces correlations between atoms which can partially overcome decoherence caused by spontaneous radiation. As a result, some amount of entanglement can survive, and moreover there is a possibility that this process can entangle separable states of two atoms. The idea that dissipation can create entanglement in physical systems, was recently developed in several papers . In the present paper we show that the dissipation due to spontaneous emission can entangle two atoms that are initially prepared in a separable state. We study the dynamics of this process. In the Markovian approximation it is given by the semi-group $`\{T_t\}`$ of completely positive linear mappings acting on density matrices . We consider time evolution of initial state of the system as well as the evolution of its entanglement, measured by so called concurrence , in the case when the photon exchange rate $`\gamma `$ is close to spontaneous emission rate of the single atom $`\gamma _0`$ and we can use the approximation $`\gamma _0=\gamma `$ (similar model was also considered in ). We calculate asymptotic stationary states $`\rho _{\mathrm{as}}`$ for the semi-group $`\{T_t\}`$ and show that they depend on initial conditions (i.e. $`\{T_t\}`$ is relaxing but not uniquely relaxing). The concurrence of $`\rho _{\mathrm{as}}`$ also depends on initial state and can be non zero for some of them. We discuss in details some classes of initial states. In particular, we show that there are pure separable states evolving to entangled mixed states and such which remain separable during evolution. The first class contains physically interesting initial state when one atom in in excited state and the other is in ground state. The relaxation process given by the semi-group $`\{T_t\}`$ produces in this case the states with entanglement monotonically increasing in time to the maximal value. The class of pure maximally entangled initial states is also discussed. Similar ”production” of entanglement is shown to be present for some classes of mixed states. On the other hand, when the photon exchange rate $`\gamma `$ is smaller then $`\gamma _0`$, the relaxation process brings all initial states to the unique asymptotic state when both atoms are in its ground states. As we show, even in that case the dynamics can entangle two separable states, but the amount of entanglement is decreasing to zero.

## 2 Pair of two-level atoms

Consider two-level atom $`A`$ with ground state $`|0`$ and excited state $`|1`$. This quantum system can be described in terms of the Hilbert space $`_A=^2`$ and the algebra $`𝔄_A`$ of $`2\times 2`$ complex matrices. If we identify $`|1`$ and $`|0`$ with vectors $`\left(\begin{array}{c}1\\ 0\end{array}\right)`$ and $`\left(\begin{array}{c}0\\ 1\end{array}\right)`$ respectively, then the raising and lowering operators $`\sigma _+,\sigma _{}`$ defined by

$$\sigma _+=|10|,\sigma _{}=|01|$$
(1)
can be expressed in terms of Pauli matrices $`\sigma _1,\sigma _2`$

$$\sigma _+=\frac{1}{2}(\sigma _1+i\sigma _2),\sigma _{}=\frac{1}{2}(\sigma _1i\sigma _2)$$
(2)
For a joint system $`AB`$ of two two-level atoms $`A`$ and $`B`$, the algebra $`𝔄_{AB}`$ is equal to $`4\times 4`$ complex matrices and the Hilbert space $`_{AB}=_A_B=^4`$. Let $`_{AB}`$ be the set of all states of the compound system i.e.

$$_{AB}=\{\rho 𝔄_{AB}:\rho 0\text{and}\mathrm{tr}\rho =1\}$$
(3)
The state $`\rho _{AB}`$ is separable , if it has the form

$$\rho =\underset{k}{}\lambda _k\rho _k^A\rho _k^B,\rho _k^A_A,\rho _k^B_B,\lambda _k0\text{and}\underset{k}{}\lambda _k=1$$
(4)
The set $`_{AB}^{\mathrm{sep}}`$ of all separable states forms a convex subset of $`_{AB}`$. When $`\rho `$ is not separable, it is called inseparable or entangled. Thus

$$_{AB}^{\mathrm{ent}}=_{AB}_{AB}^{\mathrm{sep}}$$
(5)
If $`P_{AB}`$ is a pure state i.e. $`P`$ is one-dimensional projector, then $`P`$ is separable iff partial traces $`\mathrm{tr}_AP`$ and $`\mathrm{tr}_BP`$ are also projectors. For mixed states, the separability problem is much more involved. Fortunately, in the case of $`4`$ – level compound system there is a simple necessary and sufficient condition for separability: $`\rho `$ is separable iff its partial transposition $`\rho ^{T_A}`$ is also a state . Another interesting question is how to measure the amount of entanglement a given quantum state contains. For a pure state $`P`$, the entropy of entanglement

$$E(P)=\mathrm{tr}[(\mathrm{tr}_AP)\mathrm{log}_2(\mathrm{tr}_AP)]$$
(6)
is essentailly a unique measure of entanglement . For mixed state $`\rho `$ it seems that the basic measure of entanglement is the entanglement of formation

$$E(\rho )=\mathrm{min}\underset{k}{}\lambda _kE(P_k)$$
(7)
where the minimum is taken over all possible decompositions

$$\rho =\underset{k}{}\lambda _kP_k$$
(8)
Again, in the case of $`4`$ – level system, $`E(\rho )`$ can be explicitely computed and it turns out that $`E(\rho )`$ is the function of another useful quantity $`C(\rho )`$ called concurrence, which also can be taken as a measure of entanglement . Since in the paper we use concurrence to quantify entanglement, now we discuss its definition. Let

$$\rho ^{}=(\sigma _2\sigma _2)\overline{\rho }(\sigma _2\sigma _2)$$
(9)
where $`\overline{\rho }`$ is the complex conjugation of the matrix $`\rho `$. Define also

$$\widehat{\rho }=(\rho ^{1/2}\rho ^{}\rho ^{1/2})^{1/2}$$
(10)
Then the concurrence $`C(\rho )`$ is given by

$$C(\rho )=\mathrm{max}(\mathrm{\hspace{0.17em}0},2p_{\mathrm{max}}(\widehat{\rho })\mathrm{tr}\widehat{\rho })$$
(11)
where $`p_{\mathrm{max}}(\widehat{\rho })`$ denotes the maximal eigenvalue of $`\widehat{\rho }`$. The value of the number $`C(\rho )`$ varies from $`0`$ for separable states, to $`1`$ for maximally entangled pure states.

## 3 Decay in a system of closely separated atoms

We study the spontaneous emission of two atoms separated by a distance $`R`$ small compared to the radiation wavelength . At such distances there is a substantial probability that the photon emitted by one atom will be absorbed by the other. Thus the dynamics of the system is given by the master equation

$$\frac{d\rho }{dt}=L\rho ,\rho _{AB}$$
(12)
with the following generator $`L`$

$$\begin{array}{cc}\hfill L\rho =& \frac{\gamma _0}{2}[2\sigma _{}^A\rho \sigma _+^A+2\sigma _{}^B\rho \sigma _+^B(\sigma _+^A\sigma _{}^A+\sigma _+^B\sigma _{}^B)\rho \rho (\sigma _+^A\sigma _{}^A+\sigma _+^B\sigma _{}^B)]+\hfill \\ & \frac{\gamma }{2}[2\sigma _{}^A\rho \sigma _+^B+2\sigma _{}^B\rho \sigma _+^A(\sigma _+^A\sigma _{}^B+\sigma _+^B\sigma _{}^A)\rho \rho (\sigma _+^A\sigma _{}^B+\sigma _+^B\sigma _{}^A)]\hfill \end{array}$$
(13)
where

$$\sigma _\pm ^A=\sigma _\pm 𝕀,\sigma _\pm ^B=𝕀\sigma _\pm ,\sigma _\pm =\frac{1}{2}(\sigma _1\pm i\sigma _2)$$
(14)
Here $`\gamma _0`$ is the single atom spontaneous emission rate, and $`\gamma =g\gamma _0`$ is a relaxation constant of photon exchange. In the model, $`g`$ is the function of the distance $`R`$ between atoms and $`g1`$ when $`R0`$. In this section we investigate the time evolution of the initial density matrix $`\rho `$ of the compound system, governed by the semi - group $`\{T_t\}_{t0}`$ generated by $`L`$. In particular, we will study the time development of entanglement of $`\rho `$, measured by concurrence.

Assume that the distance between atoms is so small that the exchange rate $`\gamma `$ is close to $`\gamma _0`$ and we can use the approximation $`g=1`$. Under this condition we study evolution of the system and in particular we consider asymptotic states. Direct calculations show that the semi - group $`\{T_t\}`$ generated by $`L`$ with $`g=1`$ is relaxing but not uniquely relaxing i.e. there are as many stationary states as there are initial conditions. More precisely, for a given initial state $`\rho =(\rho _{jk})`$, the state $`\rho (t)`$ at time $`t`$ has the following matrix elements

$$\begin{array}{cc}& \rho _{11}(t)=e^{2\gamma _0t}\rho _{11}\hfill \\ & \rho _{12}(t)=\frac{1}{2}[e^{2\gamma _0t}(\rho _{12}+\rho _{13})+e^{\gamma _0t}(\rho _{12}\rho _{13})]\hfill \\ & \rho _{13}(t)=\frac{1}{2}[e^{2\gamma _0t}(\rho _{12}+\rho _{13})+e^{\gamma _0t}(\rho _{13}\rho _{12})]\hfill \\ & \rho _{14}(t)=e^{\gamma _0t}\rho _{14}\hfill \\ & \rho _{22}(t)=\frac{1}{4}e^{2\gamma _0t}(\rho _{22}+\rho _{33}+2\mathrm{R}\mathrm{e}\rho _{23})+\frac{1}{2}e^{\gamma _0t}(\rho _{22}\rho _{33})+\gamma _0te^{2\gamma _0t}\rho _{11}+\frac{1}{4}(\rho _{22}+\rho _{33}2\mathrm{R}\mathrm{e}\rho _{23})\hfill \\ & \rho _{23}(t)=\frac{1}{4}e^{2\gamma _0t}(\rho _{22}+\rho _{33}+2\mathrm{R}\mathrm{e}\rho _{23})+\frac{1}{2}e^{\gamma _0t}(\rho _{23}\rho _{32})+\gamma _0te^{2\gamma _0t}\rho _{11}\frac{1}{4}(\rho _{22}+\rho _{33}2\mathrm{R}\mathrm{e}\rho _{23})\hfill \\ & \rho _{24}(t)=e^{2\gamma _0t}(\rho _{12}+\rho _{13})+\frac{1}{2}e^{\gamma _0t}(2\rho _{12}+2\rho _{13}+\rho _{24}+\rho _{34})+\frac{1}{2}(\rho _{24}\rho _{34})\hfill \\ & \rho _{33}(t)=\frac{1}{4}e^{2\gamma _0t}(\rho _{22}+\rho _{33}+2\mathrm{R}\mathrm{e}\rho _{23})\frac{1}{2}e^{\gamma _0t}(\rho _{22}\rho _{33})+\gamma _0te^{2\gamma _0t}\rho _{11}+\frac{1}{4}(\rho _{22}+\rho _{33}2\mathrm{R}\mathrm{e}\rho _{23})\hfill \\ & \rho _{34}(t)=e^{2\gamma _0t}(\rho _{12}+\rho _{13})+\frac{1}{2}e^{\gamma _0t}(2\rho _{12}+2\rho _{13}+\rho _{24}+\rho _{34})\frac{1}{2}(\rho _{24}\rho _{34})\hfill \\ & \rho _{44}(t)=\frac{1}{2}e^{2\gamma _0t}(1+\rho _{11}\rho _{44}+2\mathrm{R}\mathrm{e}\rho _{23})2\gamma _0te^{2\gamma _0t}\rho _{11}+\frac{1}{2}(1+\rho _{11}+\rho _{44}+2\mathrm{R}\mathrm{e}\rho _{23})\hfill \end{array}$$

and remaining matrix elements can be obtained by hermiticity condition $`\rho _{kj}=\overline{\rho }_{jk}`$. In the limit $`t\mathrm{}`$ we obtain asymptotic (stationary) states parametrized as follows

$$\rho _{\mathrm{as}}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \alpha & \alpha & \beta \\ 0& \alpha & \alpha & \beta \\ 0& \overline{\beta }& \overline{\beta }& 12\alpha \end{array}\right)$$
(15)
where

$$\alpha =\frac{1}{4}(\rho _{22}+\rho _{33}2\mathrm{R}\mathrm{e}\rho _{23}),\beta =\frac{1}{2}(\rho _{24}\rho _{34})$$
(16)
We can also compute concurrence of the asymptotic state and the result is:
Concurrence of asymptotic state of the semi \- group $`\{T_t\}`$ generated by $`L`$ with $`g=1`$ equals to

$$C(\rho _{\mathrm{as}})=2|\alpha |=\frac{1}{2}|\rho _{22}+\rho _{33}2\mathrm{R}\mathrm{e}\rho _{23}|$$
(17)
where $`\rho _{jk}`$ are the matrix elements of the initial state.

## 4 Some examples

In this section we consider examples of initial states and its evolution.
I. Pure separable states
Let

$$\rho =P_{\mathrm{\Psi }\mathrm{\Phi }}=P_\mathrm{\Psi }P_\mathrm{\Phi }$$
(18)
where

$$\mathrm{\Psi }=\left(\begin{array}{c}\mathrm{\Psi }_1\\ \mathrm{\Psi }_2\end{array}\right)_A,\mathrm{\Phi }=\left(\begin{array}{c}\mathrm{\Phi }_1\\ \mathrm{\Phi }_2\end{array}\right)_B$$

are normalized. Then one can check that

$$\begin{array}{cc}\hfill \alpha & =\frac{1}{4}(1|\mathrm{\Psi },\mathrm{\Phi }|^2)\hfill \\ \hfill \beta & =\frac{1}{2}(|\mathrm{\Phi }_2|^2\mathrm{\Psi }_1\overline{\mathrm{\Psi }}_2|\mathrm{\Psi }_2|^2\mathrm{\Phi }_1\overline{\mathrm{\Phi }}_2)\hfill \end{array}$$
(19)
where $`,`$ is the inner product in $`^2`$. So

$$C(\rho _{\mathrm{as}})=\frac{1}{2}(1|\mathrm{\Psi },\mathrm{\Phi }|^2)$$
(20)
From the formula (21) we see that there are separable initial states for which asymptotic states are entangled. In particular, the asymptotic state has a maximal concurrence if vectors $`\mathrm{\Psi }`$ and $`\mathrm{\Phi }`$ are orthogonal and its concurrence is zero (the state remains separable) if $`|\mathrm{\Psi },\mathrm{\Phi }|=1`$.

Now we discuss some special cases.
a.When one atom is in excited state and the other is in ground state

$$\mathrm{\Psi }=|1,\mathrm{\Phi }=|0$$

the asymptotic (mixed) state is given by

$$\rho _{\mathrm{as}}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \frac{1}{4}& \frac{1}{4}& 0\\ 0& \frac{1}{4}& \frac{1}{4}& 0\\ 0& 0& 0& 0\end{array}\right)$$

It can also be shown that in this case the relaxation process produces the state $`\rho _t`$ with concurrence

$$C(\rho _t)=\frac{1e^{\gamma _0t}}{2}$$

increasing to the maximal value equal to $`1/2`$. Thus two atoms initially in separable state become entangled for all $`t`$ and the asymptotic (steady) state attains the maximal amount of entanglement.
b. When two atoms are in excited states

$$\mathrm{\Psi }=\mathrm{\Phi }=|1$$

the asymptotic state equals to

$$|0|0$$

Thus the relaxation process brings two atoms into ground states.
c. The state $`|0|0`$ is stationary state for semi - group $`\{T_t\}`$.
II. Pure maximally entangled states
Let

$$\rho =Q(a,\theta _1,\theta _2)=\left(\begin{array}{cccc}\frac{a^2}{2}& \frac{a\sqrt{1a^2}}{2}e^{i\theta _1}& \frac{a\sqrt{1a^2}}{2}e^{i\theta _2}& \frac{a^2}{2}e^{i(\theta _1+\theta _2)}\\ \frac{a\sqrt{1a^2}}{2}e^{i\theta _1}& \frac{1a^2}{2}& \frac{1a^2}{2}e^{i(\theta _1\theta _2)}& \frac{a\sqrt{1a^2}}{2}e^{i\theta _2}\\ \frac{a\sqrt{1a^2}}{2}e^{i\theta _2}& \frac{1a^2}{2}e^{i(\theta _1\theta _2)}& \frac{1a^2}{2}& \frac{a\sqrt{1a^2}}{2}e^{i\theta _1}\\ \frac{a^2}{2}e^{i(\theta _1+\theta _2)}& \frac{\sqrt{1a^2}}{2}e^{i\theta _2}& \frac{a\sqrt{1a^2}}{2}e^{i\theta _1}& \frac{a^2}{2}\end{array}\right)$$

where $`a[0,1],\theta _1,\theta _2[0,2\pi ]`$. Pure states $`Q(a,\theta _1,\theta _2)`$ are maximally entangled and form a family of all maximally entangled states of the $`4`$ \- level system . It turns out that $`\rho _{\mathrm{as}}`$ is defined by

$$\begin{array}{cc}\hfill \alpha & =\frac{1}{4}(1a^2)(1\mathrm{cos}(\theta _1\theta _2))\hfill \\ \hfill \beta & =\frac{1}{4}a\sqrt{1a^2}(e^{i\theta _1}e^{i\theta _2})\hfill \end{array}$$
(21)
and

$$C(\rho _{\mathrm{as}})=\frac{1}{2}(1a^2)(1\mathrm{cos}(\theta _1\theta _2))$$
(22)
From the formula (23) we see that there are initial maximally entangled states which asymptotically become separable ($`a=1`$ or $`\theta _1\theta _2=2k\pi `$) and such that the asymptotic concurrence is greater then $`0`$. States with $`a=0`$ and $`\theta _1\theta _2=(2k+1)\pi `$ remain maximally entangled. For example the state

$$\frac{1}{\sqrt{2}}(|0|1|1|0)$$
(23)
is stable. On the other hand, the concurrence of

$$\frac{1}{\sqrt{2}}(|0|1+|1|0)$$
(24)
goes to zero faster then the concurrence of

$$\frac{1}{\sqrt{2}}(|0|0+|1|1)$$
(25)
as shown on Fig. 1. below. In Dicke’s theory of spontaneous radiation processes the state (24) is called subradiant whereas the state (25) has half the lifetime of a single atom and therefore is called superradiant . We see that the time-dependence of concurrence reflects the relaxation properties of those states.

Fig. 1. Concurrence as the function of time for initial states: $`\frac{1}{\sqrt{2}}(|0|0+|1|1)`$

(dotted line) and $`\frac{1}{\sqrt{2}}(|0|1+|1|0)`$ (solid line).

III. Some classes of mixed states
a. Bell - diagonal states. Let

$$\rho _\mathrm{B}=p_1|\mathrm{\Phi }^+\mathrm{\Phi }^+|+p_2|\mathrm{\Phi }^{}\mathrm{\Phi }^{}|+p_3|\mathrm{\Psi }^+\mathrm{\Psi }^+|+p_4|\mathrm{\Psi }^{}\mathrm{\Psi }^{}|$$
(26)
where Bell states $`\mathrm{\Phi }^\pm `$ and $`\mathrm{\Psi }^\pm `$ are given by

$$\mathrm{\Phi }^\pm =\frac{1}{\sqrt{2}}(|0|0\pm |1|1),\mathrm{\Psi }^\pm =\frac{1}{\sqrt{2}}(|1|0\pm |0|1)$$
(27)
It is known that all $`p_i[0,1/2],\rho _\mathrm{B}`$ is separable, while for $`p_1>1/2,\rho _\mathrm{B}`$ is entangled with concurrence equal to $`2p_11`$ (similarly for $`p_2,p_3p_4`$) . Now the asymptotic state has the form

$$\rho _{\mathrm{as}}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \frac{p_4}{2}& \frac{p_4}{2}& 0\\ 0& \frac{p_4}{2}& \frac{p_4}{2}& 0\\ 0& 0& 0& 1p_4\end{array}\right)$$
(28)
with concurrence $`C(\rho _{\mathrm{as}})=p_4`$. So even when the initial state is separable, the asymptotic state becomes entangled.
b. Werner states . Let

$$\rho _\mathrm{W}=(1p)\frac{𝕀_4}{4}+p|\mathrm{\Phi }^+\mathrm{\Phi }^+|$$
(29)
If $`p>1/3,\rho _\mathrm{W}`$ is entangled with concurrence equal to $`(3p1)/2`$. On the other hand

$$\rho _{\mathrm{as}}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \frac{1p}{8}& \frac{p1}{8}& 0\\ 0& \frac{p1}{8}& \frac{1p}{8}& 0\\ 0& 0& 0& \frac{3+p}{4}\end{array}\right)$$
(30)
has the concurrence $`C(\rho _{\mathrm{as}})=\frac{1p}{4}`$, so the asymptotic states are entangled for all $`p1`$. Notice that even completely mixed state $`\frac{𝕀_4}{4}`$ evolves to entangled asymptotic state.
c. Maximally entangled mixed states . The states

$$\rho _\mathrm{M}=\left(\begin{array}{cccc}h(\delta )& 0& 0& \delta /2\\ 0& 12h(\delta )& 0& 0\\ 0& 0& 0& 0\\ \delta /2& 0& 0& h(\delta )\end{array}\right),h(\delta )=\{\begin{array}{cc}1/3\hfill & \delta [0,2/3]\hfill \\ \delta /2\hfill & \delta [2/3,1]\hfill \end{array}$$
(31)
are conjectured to be maximally entangled for a given degree of inpurity measured by $`\mathrm{tr}\rho ^2`$ . According to (18) the concurrence of the asymptotic state is given by

$$C(\rho _{\mathrm{as}})=\frac{1}{2}(12h(\delta ))$$
(32)
Fig. 2. Concurrence of $`\rho _\mathrm{M}`$ (solid line) and $`\rho _{\mathrm{as}}`$ (dotted line) as the function of $`\mathrm{tr}\rho _\mathrm{M}^2`$

Even in that case, there are initial states (for sufficiently small $`\mathrm{tr}\rho _\mathrm{M}^2`$) such that the asymptotic state is more entangled (see Fig. 2.).

## 5 Remarks on general case

In the case of arbitrary distance between the atoms i.e. when $`g[0,1)`$, semi-group generated by $`L`$ is uniquely relaxing, with the asymptotic state $`|0|0`$. Thus, for any initial state $`\rho `$, the concurrence $`C(\rho _t)`$ approaches $`0`$ when $`t\mathrm{}`$. But it does not mean that the function $`tC(\rho _t)`$ is always monotonic. The general form of $`C(\rho _t)`$ is rather involved, so we consider only some special cases.
1. Let the initial state of the compound system be equal to $`|0|1`$. This states evolves to

$$\rho _t=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \frac{1}{2}e^{\gamma _0t}(\mathrm{cosh}\gamma t+1)& \frac{1}{2}e^{\gamma _0t}\mathrm{sinh}\gamma t& 0\\ 0& \frac{1}{2}e^{\gamma _0t}\mathrm{sinh}\gamma t& \frac{1}{2}e^{\gamma _0t}(\mathrm{cosh}\gamma t1)& 0\\ 0& 0& 0& 1e^{\gamma _0t}\mathrm{cosh}\gamma t\end{array}\right)$$
(33)
with concurrence

$$C(\rho _t)=e^{\gamma _0t}\mathrm{sinh}\gamma t$$
(34)
In the interval $`[0,t_\gamma ]`$, where

$$t_\gamma =\frac{1}{2\gamma }\mathrm{ln}\frac{\gamma _0+\gamma }{\gamma _0\gamma }$$

the function (34) is increasing to its maximal value

$$C_{\mathrm{max}}=\frac{\gamma }{\gamma _0\gamma }\left(\frac{\gamma _0+\gamma }{\gamma _0\gamma }\right)^{\frac{\gamma _0+\gamma }{2\gamma }}$$

whereas for $`t>t_\gamma `$, $`C(\rho _t)`$ decreases to $`0`$. Thus for any nonzero photon exchange rate $`\gamma `$, dynamics given by the semi - group $`\{T_t\}`$ produces some amount of entanglement between two atoms which are initially in the ground state and excited state. Note that the maximal value of $`C(\rho _t)`$ depends only on emission rates $`\gamma _0`$ and $`\gamma `$.
2. For the initial states

$$\mathrm{\Psi }^\pm =\frac{1}{\sqrt{2}}(|0|1\pm |1|0)$$

the relaxation to the asymptotic state $`|0|0`$ is given by density matrices

$$\rho _t^\pm =\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& \frac{1}{2}e^{(\gamma _0\pm \gamma )t}& \frac{1}{2}e^{(\gamma _0\pm \gamma )t}& 0\\ 0& \frac{1}{2}e^{(\gamma _0\pm \gamma )t)}& \frac{1}{2}e^{(\gamma _0\pm \gamma )t}& 0\\ 0& 0& 0& 1e^{(\gamma _0\pm \gamma )t}\end{array}\right)$$
(35)
with the corresponding concurrence

$$C(\rho _t^\pm )=e^{(\gamma _0\pm \gamma )t}$$

The state $`\mathrm{\Psi }^{}`$ is no longer stable (as in the case of $`\gamma =\gamma _0`$), but during the evolution its concurrence goes to zero slower than $`C(\rho _t^+)`$ (Fig. 3.). For $`\gamma `$ close to $`\gamma _0`$, $`\mathrm{\Psi }^{}`$ is almost stable.

Fig. 3. $`C(\rho _t^+)`$ (dotted line) and $`C(\rho _t^{})`$ (solid line) for $`\gamma /\gamma _0=0.99`$