format: text

text:
# Nuclear Mean Fields through Selfconsistent Semiclassical Calculations11footnote 1This work is supported in part by the French-Algerian CNRS/DEF scientific research agreement under Contract No. PNC 10090

## 1 Introduction

Mean-field calculations have been extremely successfull over the last 3 decades to describe the structure of stable as well as radioactive nuclei and this over a very wide range of nuclear masses. Especially effective nucleon-nucleon interactions of the type of Skyrme  and Gogny forces  have been particularly efficient in this context. Such phenomenological effective interactions can be understood as mathematically simple parametrisations of a density-dependent effective G-matrix (see  and  for a review on such effective forces).

Together with the exact treatment of the mean-field problem in the Hartree-Fock (HF) approach, semiclassical approximations thereof have proven very appropriate. Especially the approach known as the Extended Thomas-Fermi (ETF) method has been shown  to describe very accurately average nuclear properties in the sense of the Bethe-Weizsäcker mass formula . In their selfconsistent version the ETF calculations determine the structure of a given nucleus by minimizing the total energy with respect to a variation of the neutron and proton densities. Such calculations require, however, only integrated quantities as e.g. the total nuclear energy.

The aim of the present paper is, on the contrary, to investigate local quantities such as nuclear mean-field potentials, effective mass and spin-orbit form factors which are at the basis of the description of the nuclear structure and which can be obtained as function of the selfconsistent semiclassical densities. The convergence of these local (non integrated) semiclassical quantities and their comparison to the corresponding HF distributions has, to our knowledge, never been extensively investigated as this will be done here.

This paper is organized as follows. After specifying in section 2 the ETF expressions up to order $`\mathrm{}^4`$ for $`\tau [\rho ]`$ and $`\stackrel{}{J}[\rho ]`$ using the general but rather cumbersome form of these expressions given in , we show in section 3 that the semiclassical expansions which define these quantities converge rapidly for reasonable forms of the nuclear density $`\rho (\stackrel{}{r})`$. Once these expressions and their convergence established, we compare in section 4 the average mean fields obtained using these semiclassical densities to the central potentials resulting from a HF calculation. We finally conclude giving an outlook on how these calculations can be generalized to excited and rotating nuclei.

## 2 Form factors for Skyrme Interactions

For effective nucleon-nucleon interactions of the Skyrme type the total energy of a nucleus is a functional of the local densities $`\rho _q(\stackrel{}{r})`$, the kinetic energy densities $`\tau _q(\stackrel{}{r})`$ and the so called spin-orbit densities $`\stackrel{}{J}_q(\stackrel{}{r})`$

$$E=(\rho _q(\stackrel{}{r}),\tau _q(\stackrel{}{r}),\stackrel{}{J}_q(\stackrel{}{r}))d^3r,$$
(1)
where the subscript $`q=\{n,p\}`$ denotes the nucleon charge state. In the case of broken time-reversal symmetry the energy density depends, in addition, on other local quantities , such as the current density $`\stackrel{}{j}(\stackrel{}{r})`$ and the spin-vector density $`\stackrel{}{\rho }(\stackrel{}{r})`$. In what follows we will restrict ourselves to time-reversal symmetric systems leaving the case of broken time-reversal symmetry, particularly encountered in the case of rotating nuclei, to a subsequent publication.

The total energy determined in this way is exact within the HF formalism. A semiclassical approximation thereof is obtained when instead of using the exact quantum-mechanical densities $`\rho _q(\stackrel{}{r})`$, $`\tau _q(\stackrel{}{r})`$ $`\stackrel{}{J}_q(\stackrel{}{r})`$, etc. a semiclassical approximation for these quantities is used. The semiclassical densities $`\tau _q(\stackrel{}{r})`$ and $`\stackrel{}{J}_q(\stackrel{}{r})`$ are obtained in the so called Extended Thomas-Fermi model  as functions of the local density $`\rho (\stackrel{}{r})`$ and of its derivatives. The best known of these functionals is the Thomas-Fermi approximation for the kinetic energy density in the form

$$\tau _q^{(TF)}[\rho _q]=\frac{3}{5}(3\pi ^2)^{2/3}\rho _q^{5/3},q=\{n,p\}.$$
(2)
Once these functional expression are given, the total energy of the nuclear system is then uniquely determined by the knowledge of the local densities of protons and neutrons. That such a functional dependence of the total energy on the local densities $`\rho _q(\stackrel{}{r})`$ should exist is guaranteed by the famous theorem by Hohenberg and Kohn . In the most general quantum mechanical case this functional is, however, perfectly unknown. The great advantage of the semiclassical approach used here, consists in the fact that, in connection with effective interactions of the Skyrme type such an energy functional $``$ can be derived explicitly. In addition, it is to be noted that the semiclassical functionals obtained in the ETF formalism such as $`\tau [\rho ]`$ are completely general and valid for any local interaction and any nucleus, and can therefore be determined once and forever.

Once the functional of the total energy is known, one is able, due to the Hohenberg-Kohn theorem, to perform density variational calculations, where the local densities $`\rho _q(\stackrel{}{r})`$ are the variational quantities. One should, however, keep in mind that, as the ETF functionals are of semiclassical nature, the density functions $`\rho _q(\stackrel{}{r})`$ obtained as a result of the variational procedure can only be semiclassical in nature, i.e. of the liquid-drop type. Taking into account that in such a process the particle numbers $`N`$ and $`Z`$ should be conserved, one can formulate the variational principle in the form

$$\frac{\delta }{\delta \rho _q}\left\{[\rho _n(\stackrel{}{r}),\rho _p(\stackrel{}{r})]\lambda _n\rho _n(\stackrel{}{r})\lambda _p\rho _p(\stackrel{}{r})\right\}d^3r$$
(3)
with the Lagrange multipliers $`\lambda _n`$ and $`\lambda _p`$ to ensure the conservation of neutron and proton number.

This density variational problem has been solved in two different ways in the past: either by resolving the Euler-Lagrange equation  resulting from eq. (3) or by carrying out the variational calculation in an restricted subspace of functions adapted to the problem, i.e. being of semiclassical nature, free of shell oscillations in the nuclear interior and falling off exponentially in the nuclear surface. It has been shown  that modified Fermi functions which for spherical symmetry take the form

$$\rho _q(\stackrel{}{r})]=\rho _{0_q}[1+exp(\frac{rR_{0_q}}{a_q})]^{\gamma _q},q=\{n,p\}$$
(4)
are particularly well suited in this context and that the semiclassical energies obtained are, indeed, very close to those resulting from the exact variation .

* Fig. 1. Comparison of selfconsistent neutron and proton HF (solid line) and ETF densities (dashed line) for 208Pb calculated with the SkM Skyrme force.

As an example of the quality of the semiclassical density obtained by such a restricted variation in the subspace of modified Fermi functions, eq. (4), we show in Fig. 1 a comparison of the neutron and proton densities obtained in this way within the ETF approach with the corresponding Hartree-Fock (HF) densities for the nucleus 208Pb calculated with the Skyrme interaction SkM . It should be emphasized here that a similarly good agreement as the one shown in Fig. 1 is obtained for other nuclei or using other effective interactions such as the Skyrme forces SIII  and SLy4 .

The energy density $``$ appearing in eqs. (1) and (3) can be written for a Skyrme interaction as defined in ref.  in the compact form

$`(\stackrel{}{r})={\displaystyle \frac{\mathrm{}^2}{2m}}\tau +B_1\rho ^2+B_2(\rho _n^2+\rho _p^2)+B_3\rho \tau +B_4(\rho _n\tau _n+\rho _p\tau _p)`$
$`B_5(\stackrel{}{}\rho )^2B_6\left[(\stackrel{}{}\rho _n)^2+(\stackrel{}{}\rho _p)^2\right]+\rho ^\alpha [B_7\rho ^2+B_8(\rho _n^2+\rho _p^2)]`$
$`B_9\left[\stackrel{}{J}\stackrel{}{}\rho +\stackrel{}{J_n}\stackrel{}{}\rho _n+\stackrel{}{J_p}\stackrel{}{}\rho _p\right]+_{Coul}(\stackrel{}{r}),`$ (5)
given in terms of the coefficients $`B_1B_9`$ (see table 1) instead of the usual Skyrme force parameters $`t_0,t_1,t_2,t_3,x_0,x_1,x_2,x_3,W_0`$ :

|  | $`B_1`$ |  |  | $`\frac{1}{2}t_0(1+\frac{x_0}{2})`$ |  | $`B_6`$ |  |  | $`\frac{1}{16}\left[3t_1(x_1+\frac{1}{2})+t_2(x_2+\frac{1}{2})\right]`$ |
| --- | --- | --- | --- | --- | --- | --- | --- | --- | --- |
|  | $`B_2`$ |  |  | $`\frac{t_0}{2}(\frac{1}{2}+x_0)`$ |  | $`B_7`$ |  |  | $`\frac{1}{12}t_3(1+\frac{x_3}{2})`$ |
|  | $`B_3`$ |  |  | $`\frac{1}{4}\left[t_1(1+\frac{x_1}{2})+t_2(1+\frac{x_2}{2})\right]`$ |  | $`B_8`$ |  |  | $`\frac{1}{12}t_3(\frac{1}{2}+x_3)`$ |
|  | $`B_4`$ |  |  | $`\frac{1}{4}\left[t_1(x_1+\frac{1}{2})t_2(x_2+\frac{1}{2})\right]`$ |  | $`B_9`$ |  |  | $`\frac{1}{2}W_0`$ |
|  | $`B_5`$ |  |  | $`\frac{1}{16}\left[3t_1(1+\frac{x_1}{2})t_2(1+\frac{x_2}{2})\right]`$ |  |  |  |  |  |

* Tab. 1. Correspondence between the coefficients $`B_i`$ used in the text and the usual Skyrme force parameters.

In eq. (5) non-indexed quantities like $`\rho `$ correspond to the sum of neutron and proton densities as $`\rho =\rho _n+\rho _p`$ and $`_{Coul}`$ is the Coulomb energy density which can be written as the sum of the direct and an exchange contribution, the latter being taken into account in the Slater approximation

$$_{Coul}(\stackrel{}{r})=\frac{e^2}{2}\rho _p(\stackrel{}{r})d^3r^{}\frac{\rho _p(\stackrel{}{r^{}})}{|\stackrel{}{r}\stackrel{}{r^{}}|}\frac{3}{4}e^2(\frac{3}{\pi })^{1/3}\rho _p^{4/3}(\stackrel{}{r}).$$
(6)
The HF equation is obtained through the variational principle which states that the total energy of eq. (1) should be stationary with respect to any variation of the single-particle wave functions $`\phi _j^{(q)}`$ :

$$\widehat{}_q\phi _j^{(q)}=\left(\stackrel{}{}\frac{\mathrm{}^2}{2m_q^{}(\stackrel{}{r})}\stackrel{}{}+V_q(\stackrel{}{r})i\stackrel{}{W}_q(\stackrel{}{r})(\stackrel{}{}\times \stackrel{}{\sigma })\right)\phi _j^{(q)}=\epsilon _j^{(q)}\phi _j^{(q)}.$$
(7)
Here appear different form factors such as the central one-body potential $`V_q(\stackrel{}{r})`$, the effective mass $`m_q^{}(\stackrel{}{r})`$ and the spin-orbit potential $`\stackrel{}{W}_q(\stackrel{}{r})`$ which are all defined as functional derivatives of the total energy density. One obtains from eq. (5) :

$`V_q(\stackrel{}{r})={\displaystyle \frac{\delta (\stackrel{}{r})}{\delta \rho _q(\stackrel{}{r})}}=`$
$`=2(B_1\rho +B_2\rho _q)+B_3\tau +B_4\tau _q+2(B_5\mathrm{\Delta }\rho +B_6\mathrm{\Delta }\rho _q)+(2+\alpha )B_7\rho ^{\alpha +1}`$
$`+B_8\left[\alpha \rho ^{\alpha 1}{\displaystyle \underset{q}{}}\rho _q^2+2\rho ^\alpha \rho _q\right]+B_9(\mathrm{div}\stackrel{}{J}+\mathrm{div}\stackrel{}{J}_q)+V_{Coul}(\stackrel{}{r})\delta _{pq},`$ (8)
$`f_q(\stackrel{}{r})={\displaystyle \frac{m}{m_q^{}(\stackrel{}{r})}}={\displaystyle \frac{2m}{\mathrm{}^2}}{\displaystyle \frac{\delta (\stackrel{}{r})}{\delta \tau _q(\stackrel{}{r})}}=1+{\displaystyle \frac{2m}{\mathrm{}^2}}\left[B_3\rho (\stackrel{}{r})+B_4\rho _q(\stackrel{}{r})\right]`$ (9)
and

$$\stackrel{}{W}_q(\stackrel{}{r})=\frac{\delta (\stackrel{}{r})}{\delta \stackrel{}{J}_q(\stackrel{}{r})}=B_9\stackrel{}{}(\rho +\rho _q).$$
(10)
The Coulomb potential in eq. (8) is easily obtained as

$$V_{Coul}(\stackrel{}{r})=e\frac{\rho _p(\stackrel{}{r^{}})}{|\stackrel{}{r}\stackrel{}{r^{}}|}d^3r^{}(\frac{e}{\pi })^{1/3}\rho _p^{1/3}(\stackrel{}{r}).$$
(11)
It is noteworthy in this connection that for such an effective mass (9) and the spin-orbit potential (10) the energy density (5) takes the simple form

$`(\stackrel{}{r})={\displaystyle \frac{\mathrm{}^2}{2m}}{\displaystyle \underset{q}{}}f_q\tau _q+B_1\rho ^2+B_2(\rho _n^2+\rho _p^2)B_5(\stackrel{}{}\rho )^2B_6\left[(\stackrel{}{}\rho _n)^2+(\stackrel{}{}\rho _p)^2\right]`$
$`+\rho ^\alpha [B_7\rho ^2+B_8(\rho _n^2+\rho _p^2)]+{\displaystyle \underset{q}{}}\stackrel{}{J}_q\stackrel{}{W}_q+_{Coul}(\stackrel{}{r})`$ (12)
which simplifies somewhat the calculation.

All the expressions derived so far (eqs. (6) – (12)) are exact and when used with densities constructed from the single-particle wave functions $`\phi _j^{(q)}(\stackrel{}{r})`$, solutions of the HF equation (7) these quantities contain all the quantum effects of the system. If one is interested in the semiclassical approximation of these form factors one can immediately conclude, from the analytical form of eqs. (9) and (10) and the smooth behavior of the semiclassical densities, as demonstrated in Fig. 1, on the smooth behavior of the semiclassical effective mass and spin-orbit form factors. As the nuclear quantal density is well reproduced on the average by the ETF densities it appears evident that the same is going to be the case for the effective mass and the spin-orbit potential, when semiclassical, i.e. liquid-drop type densities are used in eqs. (9) and (10).

Things are, however, less evident for the central nuclear potentials. Not only is the functional expression, eq. (8), much more complicated then those for the effective mass and spin-orbit potential, but the central potential is also the only of the three functional derivatives that depends not only on the local densities $`\rho _q(\stackrel{}{r})`$ and their derivatives but also on the kinetic energy density $`\tau _q(\stackrel{}{r})`$ and the spin-orbit density $`\stackrel{}{J}_q(\stackrel{}{r})`$ which are the quantities for which the Extended Thomas-Fermi approach has written down functional expressions. We therefore choose to study the convergence of the semiclassical series corresponding to these functional expressions of $`\tau _q[\rho _q]`$ and $`\stackrel{}{J}_q[\rho _q]`$ before investigating the quality of the agreement between the HF central potential and the one obtained when using these semiclassical ETF functionals.

## 3 Convergence of ETF functionals

The semiclassical expansions of kinetic energy density $`\tau _q`$ and spin-orbit density $`\stackrel{}{J}_q`$ as functions of the local density $`\rho _q`$ are functional expressions with $`\mathrm{}`$ as order parameter. These expressions can be obtained for instance through the semiclassical $`\mathrm{}`$ expansions developed by Wigner  and Kirkwood  or through the semiclassical method of Baraff and Borowitz . In either of the two approaches one obtains functional expressions like

$$\tau _q^{(ETF)}[\rho _q]=\tau _q^{(TF)}[\rho _q]+\tau _q^{(2)}[\rho _q]+\tau _q^{(4)}[\rho _q]$$
(13)
written here for the kinetic energy density $`\tau _q(\stackrel{}{r})`$ where $`\tau _q^{(TF)}[\rho _q]`$ is the well known Thomas-Fermi expression already given in eq. (2), $`\tau _q^{(2)}[\rho _q]`$ the semiclassical correction of second order and $`\tau _q^{(4)}[\rho _q]`$ is of fourth order in $`\mathrm{}`$. The ETF expressions such as $`\tau _q^{(ETF)}[\rho _q]`$, eq. (13), up to order $`\mathrm{}^4`$ are to be understood as the converging part of an asymptotic series.

The second order term $`\tau _q^{(2)}[\rho _q]`$ has already been derived in ref.  for a Hamiltonian, eq. (7), with an effective mass $`m_q^{}=m/f_q`$ and a spin-orbit potential $`\stackrel{}{W}_q`$

$`\tau _q^{(2)}[\rho _q]={\displaystyle \frac{1}{36}}{\displaystyle \frac{(\stackrel{}{}\rho _q)^2}{\rho _q}}+{\displaystyle \frac{1}{3}}\mathrm{\Delta }\rho _q+{\displaystyle \frac{1}{6}}{\displaystyle \frac{\stackrel{}{}\rho _q\stackrel{}{}f_q}{f_q}}+{\displaystyle \frac{1}{6}}\rho _q{\displaystyle \frac{\mathrm{\Delta }f_q}{f_q}}{\displaystyle \frac{1}{12}}\rho _q({\displaystyle \frac{\stackrel{}{}f_q}{f_q}})^2`$
$`+{\displaystyle \frac{1}{2}}\left({\displaystyle \frac{2m}{\mathrm{}^2}}\right)^2\rho _q({\displaystyle \frac{\stackrel{}{W}_q}{f_q}})^2`$ (14)
the first term of which is known as the Weizsäcker correction . It was sometimes used in the past with an adjustable parameter (instead of $`\frac{1}{36}`$) in order to mock up the absence of the other second order and all of the fourth order terms. It has, however, been shown (see e.g. ref. ) that such a procedure is unable to correctly describe both the slope of the surface of the nuclear mass or charge density and, at the same time, the height of nuclear fission barriers in the actinide region. From this analysis we conclude that the inclusion of fourth order terms in the semiclassical functionals is, in fact, without credible alternatives.

In the following we are going to exploit the expressions for the 4th order functionals $`\tau _q`$ and $`\stackrel{}{J}_q`$ developed by Grammaticos and Voros . These authors have taken the convention “that any free-standing gradient operator acts only on the rightmost term” (see their remark after eq. (III.13) of ref. ). In our present work we prefer to rewrite these terms in a more conventional way and have any free standing gradient operator act, as usual, on all the terms that appear on its right hand side. We therefore write (subscripts GV refer to the Grammaticos/Voros convention)

$$\left[(\stackrel{}{}f_q\stackrel{}{})^2f_q\right]_{GV}=\frac{1}{2}(\stackrel{}{}f_q\stackrel{}{})(\stackrel{}{}f_q)^2,$$

$$\left[(\stackrel{}{}\rho _q\stackrel{}{})^2f_q\right]_{GV}=\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{}\rho _q)\frac{1}{2}\stackrel{}{}f_q\stackrel{}{}(\stackrel{}{}\rho _q)^2,$$

and

$$\left[(\stackrel{}{}f_q\stackrel{}{})(\stackrel{}{}\rho _q\stackrel{}{})f_q\right]_{GV}=\frac{1}{2}\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}f_q)^2,$$

plus terms that are obtained from these ones by interchanging the role of $`f_q`$ and $`\rho _q`$. One has also to keep in mind that Grammaticos and Voros use a slightly different definition of the effective-mass and spin-orbit form factors than the ones given in eqs. (9) and (10) above :

$$f_{_{GV}}=\frac{1}{m}f$$

and

$$\stackrel{}{S}_{_{GV}}=\frac{1}{\mathrm{}^2}\stackrel{}{W}.$$

Using these expressions we obtain the following form for the 4th order kinetic energy density, where contributions from the spin-orbit interaction have been, temporarily left out.

$`\tau _q^{(4)}[\rho _q]=(3\pi ^2)^{2/3}{\displaystyle \frac{\rho _q^{1/3}}{4320}}\{24{\displaystyle \frac{\mathrm{\Delta }^2\rho _q}{\rho _q}}60{\displaystyle \frac{\stackrel{}{}\rho _q\stackrel{}{}(\mathrm{\Delta }\rho _q)}{\rho _q^2}}28\left({\displaystyle \frac{\mathrm{\Delta }\rho _q}{\rho _q}}\right)^214{\displaystyle \frac{\mathrm{\Delta }(\stackrel{}{}\rho _q)^2}{\rho _q^2}}`$
$`+{\displaystyle \frac{280}{3}}{\displaystyle \frac{(\stackrel{}{}\rho _q)^2\mathrm{\Delta }\rho _q}{\rho _q^3}}+{\displaystyle \frac{184}{3}}{\displaystyle \frac{\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}\rho _q)^2}{\rho _q^3}}96({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}})^436{\displaystyle \frac{\mathrm{\Delta }^2f_q}{f_q}}+36{\displaystyle \frac{\mathrm{\Delta }(\stackrel{}{}f_q)^2}{f_q^2}}`$
$`18\left({\displaystyle \frac{\mathrm{\Delta }f_q}{f_q}}\right)^272{\displaystyle \frac{\stackrel{}{}f_q\stackrel{}{}(\stackrel{}{}f_q)^2}{f_q^3}}+54\left({\displaystyle \frac{\stackrel{}{}f_q}{f_q}}\right)^4+12{\displaystyle \frac{\mathrm{\Delta }(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{f_q\rho _q}}`$
$`+24{\displaystyle \frac{\stackrel{}{}f_q\stackrel{}{}(\mathrm{\Delta }\rho _q)}{f_q\rho _q}}36{\displaystyle \frac{\stackrel{}{}\rho _q\stackrel{}{}(\mathrm{\Delta }f_q)}{f_q\rho _q}}+24{\displaystyle \frac{\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}f_q)^2}{f_q^2\rho _q}}12{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)\mathrm{\Delta }f_q}{f_q^2\rho _q}}`$
$`24{\displaystyle \frac{\stackrel{}{}f_q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{f_q^2\rho _q}}44{\displaystyle \frac{\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{f_q\rho _q^2}}16{\displaystyle \frac{\stackrel{}{}f_q\stackrel{}{}(\stackrel{}{}\rho _q)^2}{f_q\rho _q^2}}`$
$`52{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)\mathrm{\Delta }\rho _q}{f_q\rho _q^2}}+30{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)^2}{f_q^2\rho _q^2}}+{\displaystyle \frac{260}{3}}{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)(\stackrel{}{}\rho _q)^2}{f_q\rho _q^3}}\}.`$ (15)
The interested reader, who might want to use the semiclassical functionals calculated here, will find in the Appendix the expression that the 4th order kinetic energy density takes in the case of spherical symmetry as well as all the other semiclassical functionals developed below.

Until now we have not taken into account the spin-orbit interaction. Its influence on the semiclassical ETF functionals is treated in ref.  and its contribution $`\tau _q^{(4)_{_{\mathrm{so}}}}`$ constitutes simply an additive term to the spin-orbit independent part of the kinetic energy density considered above. According to  :

$`\tau _q^{(4)_{_{\mathrm{so}}}}[\rho ]=(3\pi ^2)^{2/3}\left({\displaystyle \frac{2m}{\mathrm{}^2}}\right)^2{\displaystyle \frac{\rho _q^{1/3}}{4f_q^2}}\{[{\displaystyle \frac{1}{4}}\stackrel{}{W}_q\mathrm{\Delta }\stackrel{}{W}_q+{\displaystyle \frac{1}{2}}\stackrel{}{W}_q\stackrel{}{}(\mathrm{div}\stackrel{}{W}_q)+{\displaystyle \frac{1}{8}}\mathrm{\Delta }(\stackrel{}{W}_q^2)`$
$`+{\displaystyle \frac{1}{4}}(\mathrm{div}\stackrel{}{W}_q)^2]{\displaystyle \frac{1}{2f_q}}[2\stackrel{}{W}_q(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{W}_q+\mathrm{div}\stackrel{}{W}_q(\stackrel{}{}f_q\stackrel{}{W}_q)+\stackrel{}{}f_q(\stackrel{}{W}_q\stackrel{}{})\stackrel{}{W}_q`$
$`+\stackrel{}{W}_q^2\mathrm{\Delta }f_q+\stackrel{}{W}_q\stackrel{}{}(\stackrel{}{W}_q\stackrel{}{}f_q){\displaystyle \frac{1}{2}}\stackrel{}{}f_q\stackrel{}{}(\stackrel{}{W}_q^2)]+{\displaystyle \frac{3}{4f_q^2}}[(\stackrel{}{}f_q)^2\stackrel{}{W}_q^2+(\stackrel{}{W}_q\stackrel{}{}f_q)^2`$
$`({\displaystyle \frac{2m}{\mathrm{}^2}})^2\stackrel{}{W}_q^4]+{\displaystyle \frac{1}{6\rho _q}}[\stackrel{}{W}_q(\stackrel{}{}\rho _q\stackrel{}{})\stackrel{}{W}_q+(\stackrel{}{W}_q\stackrel{}{}\rho _q)\mathrm{div}\stackrel{}{W}_q]`$
$`{\displaystyle \frac{1}{6f_q\rho _q}}[(\stackrel{}{}f_q\stackrel{}{}\rho _q)\stackrel{}{W}_q^2+(\stackrel{}{}f_q\stackrel{}{W}_q)(\stackrel{}{}\rho _q\stackrel{}{W}_q)]\}.`$ (16)
It is now interesting to investigate the relative importance of the different contributions in eqs. (2), (14), (15) and (16) to the kinetic energy density obtained when using the selfconsistent semiclassical densities generated by a variational procedure restricted to functions of the type of eq. (4) as explained above. As can be seen on Fig. 2 (a) the Thomas-Fermi contribution to $`\tau [\rho ]`$ is largely dominant, at least in the nuclear bulk. Semiclassical corrections play, however, a significant role in the nuclear surface with a second-order correction which is much larger than the fourth-order term (multiplied for better visibility by a factor 10 in Fig. 2 (a)). The different contributions to the second and fourth order functional are given respectively in part (b) and (c) of the figure. We show the contributions coming form gradient terms of $`\rho `$ (term 1 and 2 in eq. (14)), of $`f`$ (terms 4 and 5) and the mixed term (term 3) as well as the spin-orbit contribution (last term) and similarly in part (c) of the figure for the fourth-order term. As can be seen, the gradient term of $`\rho `$ is dominant in 2nd order, whereas in 4nd order the spin-orbit contribution becomes also crucial. The selfconsistent HF neutron kinetic energy density is also shown on Fig. 2 (a). One notices that, except for quantum oscillations in the nuclear interior, the HF kinetic energy density is quite nicely reproduced if the semiclassical corrections $`\tau _2`$ and $`\tau _4`$ are taken into account.

* Fig. 2 Contributions from the different orders in the semiclassical expansion to the kinetic-energy density $`\tau \left[\rho \right]`$ for the selfconsistent neutron density distribution shown in Fig. 1 for 208Pb (TF (solid line), 2nd (dashed line) and 4th-order multiplied (dotted line)) are compared with the corresponding HF density (dash-dotted line) (part (a)). Different contributions to 2nd (part (b)) and 4th order (part (c)) coming from gradient terms of $`\rho `$ (solid line), of $`f`$ (dashed), of mixed terms containing gradient terms of $`\rho `$ and $`f`$ (dotted), and of the spin-orbit coupling (dash-dotted line).

It is interesting to note in this connection that despite the fact that the 2nd -order contribution $`\tau _q^{(2)}(\stackrel{}{r})`$ is one order of magnitude larger than the 4th-order contribution, after integration the 2nd -order contribution of $`\tau [\rho ]`$ to the total energy, i.e. the integral $`_qf_q(\stackrel{}{r})\tau _q^{(2)}(\stackrel{}{r})d^3r`$ is of the same order of magnitude than the corresponding 4th-order contribution (see e.g. ), which seems to indicate that there is a stronger cancellation taking place in the 2nd-order than in 4th-order contribution.

We have done the same study for the proton distribution, for other nuclei and used other effective interactions of the Skyrme type, namely the Skyrme SIII force  and the SLy4 Skyrme force . The conclusions made above remain valid in all these cases, only the relative importance of the semiclassical corrections $`\tau ^{(2)}[\rho ]`$ and $`\tau ^{(4)}[\rho ]`$ increases slightly when one goes from heavy to light nuclei.

As already discussed above, ETF functionals like $`\tau ^{(ETF)}[\rho ]`$, eq. (13), up to order $`\mathrm{}^4`$ constitute the converging part of an asymptotic expansion which needs to be truncated. Comparing, indeed, the $`\rho `$ dependence of the different orders of the semiclassical functional ($`\rho ^{5/3}`$ for $`\tau ^{(TF)}[\rho ]`$, $`\rho `$ for $`\tau ^{(2)}[\rho ]`$ and $`\rho ^{1/3}`$ for $`\tau ^{(4)}[\rho ]`$) one concludes that a term $`\tau ^{(6)}[\rho ]`$ in the ETF functional would show a $`\rho `$ dependence of the form $`\rho ^{1/3}`$ and would therefore diverge in the limit $`r\mathrm{}`$ for densities that fall off exponentially at large distances.

Let us now turn to the spin-orbit density $`\stackrel{}{J}`$. It is given in ref.  in the form of second rank tensor which is related to the components of the vector $`\stackrel{}{J}`$ by the relation

$$J_\lambda =\underset{\mu \nu }{}ϵ_{\lambda \mu \nu }J_{\mu \nu }$$

where $`ϵ_{\lambda \mu \nu }`$ is the Levi-Civita symbol. The spin being a purely quantal property with no classical analogon, there is no contribution to the semiclassical functional of $`\stackrel{}{J}`$ in lowest order, i.e. at the level of the Thomas-Fermi approach whereas one obtains for the 2nd and 4th-order contributions to the semiclassical expansion of the spin-orbit density

$$\stackrel{}{J}_q^{(2)}=\frac{2m}{\mathrm{}^2}\frac{\rho _q}{f_q}\stackrel{}{W}_q$$
(17)
and

$`\stackrel{}{J}_q^{(4)}=(3\pi ^2)^{2/3}{\displaystyle \frac{2m}{\mathrm{}^2}}{\displaystyle \frac{\rho _q^{1/3}}{8f_q}}\{[\mathrm{\Delta }\stackrel{}{W}_q+\stackrel{}{}(\mathrm{div}\stackrel{}{W}_q)]+{\displaystyle \frac{1}{f_q}}[\stackrel{}{W}_q\mathrm{\Delta }f_q+(\stackrel{}{W}_q\stackrel{}{})\stackrel{}{}f_q`$
$`+(\stackrel{}{}f_q\times \mathrm{rot}\stackrel{}{W}_q)+2(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{W}_q]{\displaystyle \frac{1}{f_q^2}}[(\stackrel{}{}f_q)^2\stackrel{}{W}_q+(\stackrel{}{}f_q\stackrel{}{W}_q)\stackrel{}{}f_q2({\displaystyle \frac{2m}{\mathrm{}^2}})^2\stackrel{}{W}_q^3]`$
$`{\displaystyle \frac{1}{3\rho _q}}[(\stackrel{}{}\rho _q\stackrel{}{})\stackrel{}{W}_q+\mathrm{div}\stackrel{}{W}_q\stackrel{}{}\rho _q{\displaystyle \frac{1}{f_q}}((\stackrel{}{}f_q\stackrel{}{}\rho _q)\stackrel{}{W}_q+(\stackrel{}{}f_q\stackrel{}{W}_q)\stackrel{}{}\rho _q)]\}`$ (18)
Let us again investigate the convergence of the semiclassical expansion associated this time with the vector field $`\stackrel{}{J}_q[\rho ]`$ and compare it with the corresponding HF density. We show in Fig. 3 the contributions to $`\stackrel{}{J}_q[\rho ]`$ from 2nd and from 4th order as well as $`\stackrel{}{J}_{HF}`$.

* Fig. 3. Contributions from the different orders in the semiclassical expansion (2nd (solid line) and 4th-order (dashed line)) to the radial part of the vector field $`\stackrel{}{J}_q\left(\stackrel{}{r}\right)`$ shown here for the selfconsistent spherical neutron distribution of 208Pb as compared to the HF spin-orbit density (dash-dotted line).

We would like to check now that the semiclassical functionals which we have written down up to order $`\mathrm{}^4`$ are indeed correct. We perform this test numerically in the following way :

One notices that when calculating the total energy through eq. (12) the kinetic energy density $`\tau _q`$ does not appear by itself, but only in connection with the form factor $`f_q`$. The $`\tau `$-dependent part of the total energy is simply obtained through the integral $`d^3r_qf_q\tau _q`$. In the contribution at order $`\mathrm{}^4`$ to this integral one can then perform integrations by parts to obtain an expression which contains only second order derivatives of the density $`\rho _q`$ and the effective-mass form factor $`f_q`$ .

$`{\displaystyle }f_q\tau _q^{(4)}d^3r=(3\pi ^2)^{2/3}{\displaystyle }d^3r\rho _q^{1/3}\{{\displaystyle \frac{1}{270}}f_q\left({\displaystyle \frac{\mathrm{\Delta }\rho _q}{\rho _q}}\right)^2{\displaystyle \frac{1}{240}}f_q{\displaystyle \frac{\mathrm{\Delta }\rho _q}{\rho _q}}\left({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}}\right)^2+{\displaystyle \frac{1}{810}}f_q\left({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}}\right)^2`$
$`{\displaystyle \frac{1}{240}}{\displaystyle \frac{(\mathrm{\Delta }f_q)^2}{f_q}}+{\displaystyle \frac{1}{120}}{\displaystyle \frac{\mathrm{\Delta }f_q(\stackrel{}{}f_q)^2}{f_q^2}}{\displaystyle \frac{1}{240}}{\displaystyle \frac{(\stackrel{}{}f_q)^4}{f_q^3}}+{\displaystyle \frac{1}{360}}\mathrm{\Delta }f_q{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{f_q\rho _q}}`$
$`{\displaystyle \frac{1}{360}}\mathrm{\Delta }\rho _q{\displaystyle \frac{(\stackrel{}{}f_q)^2}{f_q\rho _q}}{\displaystyle \frac{7}{2160}}\mathrm{\Delta }f_q({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}})^2+{\displaystyle \frac{1}{540}}({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}})^2{\displaystyle \frac{(\stackrel{}{}f_q)^2}{f_q}}+{\displaystyle \frac{7}{2160}}{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{f_q\rho _q^2}}`$
$`{\displaystyle \frac{11}{3240}}\left({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}}\right)^2{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{\rho _q}}+{\displaystyle \frac{7}{1080}}\mathrm{\Delta }\rho _q{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{}\rho _q)}{\rho _q^2}}+{\displaystyle \frac{1}{180}}\mathrm{\Delta }f_q{\displaystyle \frac{\mathrm{\Delta }\rho _q}{\rho _q}}\},`$
(19)
where $`\tau _q^{(4)}`$ is the spin-independent part of the kinetic energy density.

The dependence of the total energy density on the spin degrees of freedom enters in two different ways : through the spin-orbit part of the kinetic energy density $`\tau _q^{so}`$ and through a term $`\stackrel{}{J}_q\stackrel{}{W}_q`$. This total spin-dependence is then given by $`_q\left(\frac{\mathrm{}^2}{2m}f_q\tau _q^{so}[\rho ]+\stackrel{}{W}_q\stackrel{}{J}_q[\rho ]\right)d^3r`$. It can be shown after some integration by parts with eqs. (14), (16), (17) and (18) that in the different orders of the semiclassical expansion this integral takes on the simple form

$$\left[\frac{\mathrm{}^2}{2m}f_q\tau _q^{(2)_{_{so}}}[\rho ]+\stackrel{}{W}_q\stackrel{}{J}_q^{(2)}[\rho ]\right]d^3r=\frac{m}{\mathrm{}^2}\underset{q}{}\frac{\rho _q}{f_q}\stackrel{}{W}_q^2d^3r$$
(20)
and that

$`{\displaystyle \left[\frac{\mathrm{}^2}{2m}f_q\tau _q^{(4)_{_{so}}}[\rho ]+\stackrel{}{W}_q\stackrel{}{J}_q^{(4)}[\rho ]\right]d^3r}=`$
$`=(3\pi ^2)^{2/3}{\displaystyle \frac{m}{\mathrm{}^2}}{\displaystyle }{\displaystyle \frac{\rho _q^{1/3}}{f_q}}\{{\displaystyle \frac{1}{4}}(\mathrm{div}\stackrel{}{W}_q)^2{\displaystyle \frac{3}{8}}\mathrm{div}\stackrel{}{W}_q{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{W}_q)}{f_q}}`$
$`+{\displaystyle \frac{1}{16}}\stackrel{}{W}_q^2{\displaystyle \frac{\mathrm{\Delta }f_q}{f_q}}+{\displaystyle \frac{1}{8}}{\displaystyle \frac{(\stackrel{}{}f_q\stackrel{}{W}_q)^2}{f_q^2}}+{\displaystyle \frac{1}{24}}\mathrm{div}\stackrel{}{W}_q{\displaystyle \frac{(\stackrel{}{}\rho _q\stackrel{}{W}_q)}{\rho _q}}`$
$`+{\displaystyle \frac{1}{48}}\stackrel{}{W}_q^2{\displaystyle \frac{\mathrm{\Delta }\rho _q}{\rho _q}}{\displaystyle \frac{1}{72}}\stackrel{}{W}_q^2\left({\displaystyle \frac{\stackrel{}{}\rho _q}{\rho _q}}\right)^2+{\displaystyle \frac{1}{2}}\left({\displaystyle \frac{m}{\mathrm{}^2}}\right)^2{\displaystyle \frac{\stackrel{}{W}_q^4}{f_q^2}}\}d^3r.`$ (21)
We have tested the semiclassical functionals given above by verifying numerically that these integral relations (19), (20) and (21) hold true.

We have also evaluated different integrals involving these functionals and which have been calculated in ref. . We obtain agreement with their results of the order of 1 to 2 %, which is of the same order as their agreement between the results of a full variational calculation and one in the restricted subspace of modified Fermi functions.

As can be seen on eq. (8) only the divergence of the vector field $`\stackrel{}{J}`$ is present in the expression of the central one-body potential. One obtains from eqs. (17) and (18) respectively the contributions to 2nd order

$$\mathrm{div}\stackrel{}{J}_q^{(2)}=\frac{2m}{\mathrm{}^2}\frac{1}{f_q}\left[\rho _q\mathrm{div}\stackrel{}{W}_q+(\stackrel{}{}\rho _q\stackrel{}{W}_q)\frac{\rho _q}{f_q}(\stackrel{}{}f_q\stackrel{}{W}_q)\right]$$
(22)
and, after some lengthy but straightforward calculation, to $`4^{\mathrm{th}}`$ order

$`\mathrm{div}\stackrel{}{J}_q^{(4)}=(3\pi ^2)^{2/3}{\displaystyle \frac{2m}{\mathrm{}^2}}{\displaystyle \frac{\rho _q^{1/3}}{8f_q}}\{2\mathrm{d}\mathrm{i}\mathrm{v}\left[\stackrel{}{}(\mathrm{div}\stackrel{}{W}_q)\right]+{\displaystyle \frac{1}{2f_q}}[2\mathrm{\Delta }f_q\mathrm{div}\stackrel{}{W}_q+(\stackrel{}{}^3f_q\stackrel{}{W}_q)`$
$`+3\mathrm{\Delta }(\stackrel{}{}f_q\stackrel{}{W}_q)+(\stackrel{}{}f_q\mathrm{\Delta }\stackrel{}{W}_q)+4\stackrel{}{}f_q\stackrel{}{}(\mathrm{div}\stackrel{}{W}_q)]{\displaystyle \frac{1}{f_q^2}}[(\stackrel{}{}f_q)^2\mathrm{div}\stackrel{}{W}_q`$
$`+3\mathrm{\Delta }f_q(\stackrel{}{}f_q\stackrel{}{W}_q)+5\stackrel{}{}f_q(\stackrel{}{W}_q\stackrel{}{})\stackrel{}{}f_q+5\stackrel{}{}f_q(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{W}_q)`$
$`+2\stackrel{}{}f_q(\stackrel{}{}f_q\times \mathrm{rot}\stackrel{}{W}_q)2({\displaystyle \frac{2m}{\mathrm{}^2}})^2(\stackrel{}{W}_q^2\mathrm{div}\stackrel{}{W}_q+\stackrel{}{W}_q\stackrel{}{}(\stackrel{}{W}_q^2))]`$
$`+{\displaystyle \frac{6}{f_q^3}}[(\stackrel{}{}f_q)^2(\stackrel{}{}f_q\stackrel{}{W}_q)({\displaystyle \frac{2m}{\mathrm{}^2}})^2\stackrel{}{W}_q^2(\stackrel{}{}f_q\stackrel{}{W}_q)]{\displaystyle \frac{1}{6\rho _q}}[\mathrm{\Delta }(\stackrel{}{}\rho _q\stackrel{}{W}_q)+(\stackrel{}{}\rho _q\mathrm{\Delta }\stackrel{}{W}_q)`$
$`(\stackrel{}{}^3\rho _q\stackrel{}{W}_q)+2\mathrm{\Delta }\rho _q\mathrm{div}\stackrel{}{W}_q+6\stackrel{}{}\rho _q\stackrel{}{}(\mathrm{div}\stackrel{}{W}_q)]+{\displaystyle \frac{1}{3f_q\rho _q}}[2(\stackrel{}{}f_q\stackrel{}{}\rho _q)\mathrm{div}\stackrel{}{W}_q`$
$`+\stackrel{}{W}_q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{}\rho _q)+\mathrm{\Delta }\rho _q(\stackrel{}{}f_q\stackrel{}{W}_q)+2\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{W}_q)+\mathrm{\Delta }f_q(\stackrel{}{}\rho _q\stackrel{}{W}_q)`$
$`+\stackrel{}{}\rho _q(\stackrel{}{}f_q\times \mathrm{rot}\stackrel{}{W}_q)+2\stackrel{}{}\rho _q(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{W}_q)]{\displaystyle \frac{1}{3f_q^2\rho _q}}[5(\stackrel{}{}f_q\stackrel{}{}\rho _q)(\stackrel{}{}f_q\stackrel{}{W}_q)`$
$`+(\stackrel{}{}f_q)^2(\stackrel{}{}\rho _q\stackrel{}{W}_q)2({\displaystyle \frac{2m}{\mathrm{}^2}})^2\stackrel{}{W}_q^2(\stackrel{}{}\rho _q\stackrel{}{W}_q)]+{\displaystyle \frac{2}{9\rho _q^2}}[\stackrel{}{}\rho _q(\stackrel{}{}\rho _q\stackrel{}{})\stackrel{}{W}_q)`$
$`+(\stackrel{}{}\rho _q)^2\mathrm{div}\stackrel{}{W}_q]{\displaystyle \frac{2}{9f_q\rho _q^2}}[(\stackrel{}{}f_q\stackrel{}{}\rho _q)(\stackrel{}{}\rho _q\stackrel{}{W}_q)+(\stackrel{}{}\rho _q)^2(\stackrel{}{}f_q\stackrel{}{W}_q)]\}.`$ (23)
## 4 Average nuclear potentials

One now controls all the ingredients which enter into the calculation of the nuclear central potentials, eq. (8), for effective interactions of the Skyrme type. It is now interesting to look at the convergence of the expressions which define these average fields calculated with the semiclassical functionals $`\tau _q[\rho ]`$ and $`\stackrel{}{J}_q[\rho ]`$ and to check how these potentials compare with the ones obtained in the HF approach.

We therefore show in Fig. 4 the neutron and proton nuclear central potentials for the nucleus 208Pb obtained for the Skyrme interaction SkM. The Coulomb potential $`V_{Coul}`$ for the proton field has been left out in this investigation, because it is directly given through the proton density (see eq. (11)). As the latter is very well reproduced, except for quantum oscillations in the nuclear interior, we already know that the Coulomb potential calculated through this semiclassical density will, indeed, reproduce on the average the exact one calculated from the quantum-mechanical densities.

We show a comparison between the HF neutron and proton central potentials with the ones obtained using the selfconsistent semiclassical densities $`\rho _n`$ and $`\rho _p`$ but restricting ourselves to the TF approximations for the functional $`\tau [\rho ]`$, eq. (2), and $`\stackrel{}{J}[\rho ]`$ (which is zero as explained above). We do not want to call this the Thomas-Fermi approximation to the nuclear central fields since even if we have used the above mentioned functionals in their Thomas-Fermi approximation, the nuclear structure has been determined through a full variational calculation including the functionals up to order $`\mathrm{}^4`$.

As seen on the figure the reproduction of the HF selfconsistent fields is already quite remarkable at the lowest (TF) order in the semiclassical expansion. Apart from shell oscillations in the nuclear interior and a small wiggle in the TF potential in the surface region the agreement seems very satisfactory.

* Fig. 4. Comparison of the Hartree-Fock central nuclear potentials (solid line) for protons and neutrons with the corresponding semiclassical potentials (dashed line) obtained using the Thomas-Fermi approximation for the functionals $`\tau \left[\rho \right]`$ and $`\stackrel{}{J}\left[\rho \right]`$.

It is now interesting to study the contributions to the nuclear central fields coming from higher orders in the semiclassical expansion. For this reason we show in Fig. 5 the corrections $`\delta V_n^{(2)}`$ and $`\delta V_n^{(4)}`$ defined as (see eq. (8))

$$\delta V_n^{(2)}=(B_3+B_4)\tau _n^{(2)}+B_3\tau _p^{(2)}+B_9(2\mathrm{div}\stackrel{}{J}_n^{(2)}+\mathrm{div}\stackrel{}{J}_p^{(2)}),$$

(24)
$$\delta V_n^{(4)}=(B_3+B_4)\tau _n^{(4)}+B_3\tau _p^{(4)}+B_9(2\mathrm{div}\stackrel{}{J}_n^{(4)}+\mathrm{div}\stackrel{}{J}_p^{(4)})$$

together with the semiclassical TF potential already shown in Fig. 4. It turns out that these corrections are rather small and we have to multiply $`\delta V_n^{(2)}`$ by a factor of 10 and $`\delta V_n^{(4)}`$ by a factor of 100 to make their relative importance better visible on figure 5.

* Fig. 5. Semiclassical neutron potential using the TF approximation (solid line) of the semiclassical functionals $`\tau \left[\rho \right]`$, and $`\stackrel{}{J}\left[\rho \right]`$ and corrections coming from second (dashed) and fourth order (dotted line) in the semiclassical expansion. For better visibility the second order correction has been multiplied by a factor of 10 and the fourth order by a factor of 100.

It can be seen that both these terms give a contribution in the nuclear surface where the lowest order term showed some deviation from the HF potentials. It is therefore to be expected that potential using the second order functionals $`\tau ^{(2)}[\rho ]`$ and $`\stackrel{}{J}^{(2)}[\rho ]`$ will partially correct for this deficiency and be quite close to the selfconsistent HF potentials. Due to the smallness of the 4th order term in Fig. 5 we can expect the semiclassical potentials obtained using the full semiclassical functionals up to 4th order to be practically indistinguishable from the ones using the 2nd order corrections only. This conclusion is, indeed, confirmed on Fig. 6.

* Fig. 6. Comparison of the Hartree-Fock central nuclear potentials (solid line) for protons and neutrons with the semiclassical potentials (dashed line) obtained by using the semiclassical functionals $`\tau \left[\rho \right]`$ and $`\stackrel{}{J}\left[\rho \right]`$ up to order $`\mathrm{}^2`$ (dashed line). The semiclassical potentials including the functionals up to order $`\mathrm{}^4`$ are indistinguishable from the latter ones.

The same kind of calculations have been performed also for lighter nuclei down to 40Ca. The quality of the agreement between HF and semiclassical potentials is the same as the one obtained for the nucleus 208Pb studied above. As already mentioned, the effective-mass form factor $`f_q(\stackrel{}{r})`$, eq. (9) and the spin-orbit potential $`\stackrel{}{W}_q(\stackrel{}{r})`$, eq. (10) are, except for shell oscillations in the nuclear interior, very well reproduced, since the nuclear densities which directly determine these quantities are well reproduced. For this reason we do not explicitly show these quantities. These conclusions remain valid when other effective interactions of the Skyrme type are used.

## 5 Summary and conclusions

We have demonstrated that using the Extended Thomas-Fermi approach one is not only able to give a very precise description of average nuclear properties but that this method is also able to reproduce quite nicely local quantities, not only neutron and proton density distributions but also the corresponding nuclear central potentials, effective-mass form factors and spin-orbit potentials. These are precisely the ingredients of the Schrödinger like Hartree-Fock equation, eq. (7), which arises from the variational principle. Within this semiclassical approach which relies on a density variational calculation one should therefore be able to solve in an approximate way the quantum mechanical problem without having to go through the full selfconsistency problem of the HF approach. This is the essential idea of an approximate solution of the HF problem known as the “expectation value method” (EVM) . It consists in constructing the ground state Slater determinant from the eigenfunctions of eq. (7) using the ETF fields $`V_q(\stackrel{}{r})`$, $`m_q^{}(\stackrel{}{r})`$ and $`\stackrel{}{W}_q(\stackrel{}{r})`$ (to second order in the ETF expansions) and calculating with this Slater determinant the expectation value of the total Skyrme Hamiltonian. The rapid convergence of the ETF functionals demonstrated above explains, a posteriori, the success of this approach.

The interested reader might object that nowadays, where computational power has been increased tremendously, there is no real need for semiclassical approximations, but all calculations of nuclear structure should be directly performed at the level of the Hartree-Fock model (or beyond). The point, however, is that as soon as one is interested in nuclear systematics where one is looking at the behavior of nuclei over a wide range of the nuclear chart, semiclassical approximations are without any credible substitute.

A point that might e.g. be interesting to study is the variation of the diffuseness of the nuclear densities and central potentials when increasing the nuclear excitation and/or when going to rotating nuclei. The approach we have developed here is, indeed, easily generalized to the description of excited or rotating nuclei. If one is interested in hot nuclear systems one simply need to replace the semiclassical functionals for $`\tau [\rho ]`$ and $`\stackrel{}{J}[\rho ]`$ which we have developed by those derived in ref.  for the ETF approach at finite temperature. In this case the coefficients in the semiclassical expansions, eqs. (2), (14), (15), (16) etc. are to be replaced by a combination of Fermi integrals .

If on the other hand one is interested in the description of systems breaking time reversal symmetry, the energy density, eq. (5), will be changed  and some additional densities will appear, such as the current density $`\stackrel{}{j}(\stackrel{}{r})`$ or the spin-vector density $`\stackrel{}{\rho }(\stackrel{}{r})`$ which, in the case of rotations, are a manifestation of the time-odd part of the density matrix generated by the cranking piece of the Hamiltonian. This causes not only a change in the analytical form of quantities such as the average potentials $`V_q(\stackrel{}{r})`$ which are now going to depend on these additional densities but also leads to the appearance of additional terms in the functionals $`\tau [\rho ]`$ and $`\stackrel{}{J}[\rho ]`$. These functionals have already been determined in  which makes it quite straightforward to calculate selfconsistent semiclassical fields such as $`V_q(\stackrel{}{r})`$ for rotating nuclei very similarly to what we have done here in the static case and investigate the dependence of these quantities with increasing angular momentum. Investigations along these lines are currently in progress.

The present method can also be profitably exploited to establish some easy to use parametrization (as function of mass number A and isospin parameter I=(N-Z)/A) of central and spin-orbit potential and effective mass and this over a wide region of the nuclear chart . Such an investigation could be the ideal starting point for Strutinski shell correction calculations as e.g. formulated in the so-called Extended Thomas-Fermi plus Strutinski integral (ETFSI) method (see ref.  and references given therein) where average mean fields like the ones investigated here are used to determine shell corrections.

Acknowledgments

The authors are grateful to K. Pomorski for having suggested the investigation undertaken here and for his constant interest in this work. The authors have also benefitted from fruitful discussions with P. Quentin. One of us (K.B.) gratefully acknowledges the hospitality extended to him during his stay at the Institut de Recherches Subatomiques of Strasbourg.

Appendix

As mentioned in the text we gather here the expressions the ETF functionals take in the case of spherical symmetry.

For the 4th order spin independent part of the kinetic energy density, eq. (15), one obtains (primes denoting derivatives with respect to the radial variable $`r`$) :

$`\tau _q^{(4)}[\rho ]=(3\pi ^2)^{2/3}{\displaystyle \frac{\rho _q^{1/3}}{4320}}\{{\displaystyle \frac{24}{\rho _q}}[\rho _q^{(4)}+{\displaystyle \frac{4}{r}}\rho _q^{^{\prime \prime \prime }}]{\displaystyle \frac{8}{\rho _q^2}}[11\rho _q^{^{\prime \prime \prime }}\rho _q^{^{}}+7(\rho _q^{^{\prime \prime }})^2+{\displaystyle \frac{36}{r}}\rho _q^{^{\prime \prime }}\rho _q^{^{}}{\displaystyle \frac{1}{r^2}}(\rho _q^{^{}})^2]`$
$`+{\displaystyle \frac{8}{3\rho _q^3}}[81\rho _q^{^{\prime \prime }}(\rho _q^{^{}})^2+{\displaystyle \frac{70}{r}}(\rho _q^{^{}})^3]{\displaystyle \frac{96}{\rho _q^4}}(\rho _q^{^{}})^4{\displaystyle \frac{36}{f_q}}[f_q^{(4)}+{\displaystyle \frac{4}{r}}f_q^{^{\prime \prime \prime }}]+{\displaystyle \frac{18}{f_q^2}}[4f_q^{^{\prime \prime \prime }}f_q^{^{}}+`$
$`+3(f_q^{^{\prime \prime }})^2+{\displaystyle \frac{4}{r}}f_q^{^{\prime \prime }}f_q^{^{}}{\displaystyle \frac{4}{r^2}}(f_q^{^{}})^2]{\displaystyle \frac{144}{f_q^3}}f_q^{^{\prime \prime }}(f_q^{^{}})^2+{\displaystyle \frac{54}{f_q^4}}(f_q^{^{}})^4+{\displaystyle \frac{12}{f_q\rho _q}}[3f_q^{^{}}\rho _q^{^{\prime \prime \prime }}+2f_q^{^{\prime \prime }}\rho _q^{^{\prime \prime }}`$
$`2f_q^{^{\prime \prime \prime }}\rho _q^{^{}}+{\displaystyle \frac{6}{r}}f_q^{^{}}\rho _q^{^{\prime \prime }}{\displaystyle \frac{4}{r}}f_q^{^{\prime \prime }}\rho _q^{^{}}+{\displaystyle \frac{2}{r^2}}f_q^{^{}}\rho _q^{^{}}]+{\displaystyle \frac{12}{f_q^2\rho _q}}[f_q^{^{}}f_q^{^{\prime \prime }}\rho _q^{^{}}2(f_q^{^{}})^2\rho _q^{^{\prime \prime }}{\displaystyle \frac{2}{r}}(f_q^{^{}})^2\rho _q^{^{}}]`$
$`{\displaystyle \frac{4}{f_q\rho _q^2}}[11f_q^{^{\prime \prime }}(\rho _q^{^{}})^2+32f_q^{^{}}\rho _q^{^{}}\rho _q^{^{\prime \prime }}+{\displaystyle \frac{26}{r}}f_q^{^{}}(\rho _q^{^{}})^2]+{\displaystyle \frac{30}{f_q^2\rho _q^2}}(f_q^{^{}})^2(\rho _q^{^{}})^2+{\displaystyle \frac{260}{3f_q\rho _q^3}}f_q^{^{}}(\rho _q^{^{}})^3\}`$ (A.1)
When giving the spin dependent part of $`\tau _q^{(4)}[\rho ]`$ we take advantage of the fact that the spin-orbit potential $`\stackrel{}{W}_q`$ has for Skyrme forces the simple form of eq. (10) which allows us to introduce the quantity

$$A_q=\rho +\rho _q\text{ }\text{ }\stackrel{}{W}_q=B_9\stackrel{}{}A_q$$
(A.2)
Using the vector identity

$$\mathrm{rot}(\mathrm{rot}\stackrel{}{a})=\stackrel{}{}(\mathrm{div}\stackrel{}{a})\mathrm{\Delta }\stackrel{}{a}$$

one shows that because of the form of the spin-orbit potential, eq. (10), one simply has

$$\stackrel{}{W}\mathrm{\Delta }\stackrel{}{W}=\stackrel{}{W}\stackrel{}{}(\mathrm{div}\stackrel{}{W})$$

which then allows us to write the spin dependent part of $`\tau _q^{(4)_{_{\mathrm{so}}}}`$ in the form of the local densities $`\rho _n`$ and $`\rho _p`$ and of their derivatives :

$`\tau _q^{(4)_{_{\mathrm{so}}}}[\rho ]=(3\pi ^2)^{2/3}({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2{\displaystyle \frac{\rho _q^{1/3}}{4f_q^2}}\{[{\displaystyle \frac{3}{4}}\stackrel{}{}A_q\stackrel{}{}^3A_q+{\displaystyle \frac{1}{8}}\mathrm{\Delta }(\stackrel{}{}A_q)^2+{\displaystyle \frac{1}{4}}(\mathrm{\Delta }A_q)^2]`$
$`{\displaystyle \frac{1}{2f_q}}[\mathrm{\Delta }A_q(\stackrel{}{}f_q\stackrel{}{}A_q)+2\stackrel{}{}A_q(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{}A_q+\stackrel{}{}f_q(\stackrel{}{}A_q\stackrel{}{})\stackrel{}{}A_q+\mathrm{\Delta }f_q(\stackrel{}{}A_q)^2`$
$`+\stackrel{}{}A_q\stackrel{}{}(\stackrel{}{}A_q\stackrel{}{}f_q){\displaystyle \frac{1}{2}}\stackrel{}{}f_q\stackrel{}{}(\stackrel{}{}A_q)^2]+{\displaystyle \frac{3}{4f_q^2}}[(\stackrel{}{}f_q)^2(\stackrel{}{}A_q)^2+(\stackrel{}{}f_q\stackrel{}{}A_q)^2`$
$`({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(\stackrel{}{}A_q)^4]+{\displaystyle \frac{1}{6\rho _q}}[\stackrel{}{}A_q(\stackrel{}{}\rho _q\stackrel{}{})\stackrel{}{}A_q+\mathrm{\Delta }A_q(\stackrel{}{}\rho _q\stackrel{}{}A_q)]`$
$`{\displaystyle \frac{1}{6f_q\rho _q}}[(\stackrel{}{}f_q\stackrel{}{}\rho _q)(\stackrel{}{}A_q)^2+(\stackrel{}{}f_q\stackrel{}{}A_q)(\stackrel{}{}\rho _q\stackrel{}{}A_q)]\}`$ (A.3)
which for spherical symmetry reads

$`\tau _q^{(4)_{_{\mathrm{so}}}}[\rho ]=(3\pi ^2)^{2/3}\left({\displaystyle \frac{mW_0}{\mathrm{}^2}}\right)^2{\displaystyle \frac{\rho _q^{1/3}}{4f_q^2}}\{{\displaystyle \frac{1}{2}}[2A_q^{}A_q^{\prime \prime \prime }+(A_q^{\prime \prime })^2+{\displaystyle \frac{6}{r}}A_q^{}A_q^{\prime \prime }{\displaystyle \frac{1}{r^2}}(A_q^{})^2]`$
$`{\displaystyle \frac{A_q^{}}{f_q}}\left[f_q^{\prime \prime }A_q^{}+2f_q^{}A_q^{\prime \prime }+{\displaystyle \frac{2}{r}}f_q^{}A_q^{}\right]+{\displaystyle \frac{3(A_q^{})^2}{4f_q^2}}\left[2(f_q^{})^2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(A_q^{})^2\right]`$
$`+{\displaystyle \frac{\rho _q^{}A_q^{}}{3f_q\rho _q}}[f_q(A_q^{\prime \prime }+{\displaystyle \frac{1}{r}}A_q^{})f_q^{}A_q^{}]\}.`$ (A.4)
The 4th order spin-orbit density $`\stackrel{}{J}_q^{(4)}[\rho ]`$, eq. (18) written in terms of the function $`A_q`$ defined above is given by

$`\stackrel{}{J}_q^{(4)}[\rho ]=(3\pi ^2)^{2/3}{\displaystyle \frac{mW_0}{\mathrm{}^2}}{\displaystyle \frac{\rho _q^{1/3}}{8f_q}}\{2\stackrel{}{}^3A_q+{\displaystyle \frac{1}{f_q}}[\mathrm{\Delta }f_q\stackrel{}{}A_q+(\stackrel{}{}A_q\stackrel{}{})\stackrel{}{}f_q`$
$`+2(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{}A_q]{\displaystyle \frac{1}{f_q^2}}[(\stackrel{}{}f_q)^2\stackrel{}{}A_q+(\stackrel{}{}f_q\stackrel{}{}A_q)\stackrel{}{}f_q`$
$`2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(\stackrel{}{}A_q)^3]{\displaystyle \frac{1}{3\rho _q}}[(\stackrel{}{}\rho _q\stackrel{}{})\stackrel{}{}A_q+\mathrm{\Delta }A_q\stackrel{}{}\rho _q]`$
$`+{\displaystyle \frac{1}{3f_q\rho _q}}[(\stackrel{}{}f_q\stackrel{}{}\rho _q)\stackrel{}{}A_q+(\stackrel{}{}f_q\stackrel{}{}A_q)\stackrel{}{}\rho _q]\}`$ (A.5)
which in the case of a spherically symmetric system takes the form

$`\stackrel{}{J}_q^{(4)}[\rho ]=(3\pi ^2)^{2/3}{\displaystyle \frac{mW_0}{\mathrm{}^2}}{\displaystyle \frac{\rho _q^{1/3}}{4f_q}}\{[A_q^{\prime \prime \prime }+{\displaystyle \frac{2}{r}}A_q^{\prime \prime }{\displaystyle \frac{2}{r^2}}A_q^{}]+{\displaystyle \frac{1}{f_q}}[f_q^{\prime \prime }A_q^{}+f_q^{}A_q^{\prime \prime }+{\displaystyle \frac{1}{r}}f_q^{}A_q^{}]`$
$`{\displaystyle \frac{A_q^{}}{f_q^2}}[(f_q^{})^2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(A_q^{})^2]{\displaystyle \frac{\rho _q^{}}{3f_q\rho _q}}[f_q(A_q^{\prime \prime }+{\displaystyle \frac{1}{r}}A_q^{})f_q^{}A_q^{}]\}\stackrel{}{u}_r`$ (A.6)
with the unit vector in radial direction $`\stackrel{}{u}_r`$.

The corresponding expressions for $`\mathrm{div}\stackrel{}{J}_q^{(4)}[\rho ]`$ are the following :

$`\mathrm{div}\stackrel{}{J}_q^{(4)}[\rho ]=(3\pi ^2)^{2/3}{\displaystyle \frac{mW_0}{\mathrm{}^2}}{\displaystyle \frac{\rho _q^{1/3}}{8f_q}}\{2\mathrm{\Delta }^2A_q+{\displaystyle \frac{1}{2f_q}}[2\mathrm{\Delta }f_q\mathrm{\Delta }A_q+\stackrel{}{}^3f_q\stackrel{}{}A_q`$
$`+3\mathrm{\Delta }(\stackrel{}{}f_q\stackrel{}{}A_q)+5\stackrel{}{}f_q\stackrel{}{}^3A_q]{\displaystyle \frac{1}{f_q^2}}[(\stackrel{}{}f_q)^2\mathrm{\Delta }A_q+3\mathrm{\Delta }f_q(\stackrel{}{}f_q\stackrel{}{}A_q)`$
$`+5\stackrel{}{}f_q(\stackrel{}{}A_q\stackrel{}{})\stackrel{}{}f_q+5\stackrel{}{}f_q(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{}A_q)2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(\mathrm{\Delta }A_q(\stackrel{}{}A_q)^2`$
$`+\stackrel{}{}A_q\stackrel{}{}(\stackrel{}{}A_q)^2)]+{\displaystyle \frac{6}{f_q^3}}[(\stackrel{}{}f_q)^2(\stackrel{}{}f_q\stackrel{}{}A_q)({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(\stackrel{}{}A_q)^2(\stackrel{}{}f_q\stackrel{}{}A_q)]`$
$`{\displaystyle \frac{1}{6\rho _q}}\left[\mathrm{\Delta }(\stackrel{}{}\rho _q\stackrel{}{}A_q)+7\stackrel{}{}\rho _q\stackrel{}{}^3A_q\stackrel{}{}^3\rho _q\stackrel{}{}A_q+2\mathrm{\Delta }\rho _q\mathrm{\Delta }A_q\right]`$
$`+{\displaystyle \frac{1}{3f_q\rho _q}}[2(\stackrel{}{}f_q\stackrel{}{}\rho _q)\mathrm{\Delta }A_q+\stackrel{}{}A_q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{}\rho _q)+\mathrm{\Delta }\rho _q(\stackrel{}{}f_q\stackrel{}{}A_q)+`$
$`+2\stackrel{}{}\rho _q\stackrel{}{}(\stackrel{}{}f_q\stackrel{}{}A_q)+\mathrm{\Delta }f_q(\stackrel{}{}\rho _q\stackrel{}{}A_q)+2\stackrel{}{}\rho _q(\stackrel{}{}f_q\stackrel{}{})\stackrel{}{}A_q]`$
$`{\displaystyle \frac{1}{3f_q^2\rho _q}}[5(\stackrel{}{}f_q\stackrel{}{}\rho _q)(\stackrel{}{}f_q\stackrel{}{}A_q)+(\stackrel{}{}f_q)^2(\stackrel{}{}\rho _q\stackrel{}{}A_q)`$
$`2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(\stackrel{}{}A_q)^2(\stackrel{}{}\rho _q\stackrel{}{}A_q)]+{\displaystyle \frac{2}{9\rho _q^2}}[\stackrel{}{}\rho _q(\stackrel{}{}\rho _q\stackrel{}{})\stackrel{}{}A_q+(\stackrel{}{}\rho _q)^2\mathrm{\Delta }A_q]`$
$`{\displaystyle \frac{2}{9f_q\rho _q^2}}[(\stackrel{}{}f_q\stackrel{}{}\rho _q)(\stackrel{}{}\rho _q\stackrel{}{}A_q)+(\stackrel{}{}\rho _q)^2(\stackrel{}{}f_q\stackrel{}{}A_q)]\}.`$
(A.7)
and in the case of spherical symmetry

$`\mathrm{div}\stackrel{}{J}_q^{(4)}[\rho ]=(3\pi ^2)^{2/3}{\displaystyle \frac{mW_0}{\mathrm{}^2}}{\displaystyle \frac{\rho _q^{1/3}}{4f_q}}\{(A_q^{(4)}+{\displaystyle \frac{4}{r}}A_q^{\prime \prime \prime })+{\displaystyle \frac{1}{f_q}}[f_q^{\prime \prime \prime }A_q^{}+2f_q^{\prime \prime }A_q^{\prime \prime }+2f_q^{}A_q^{\prime \prime \prime }`$
$`+{\displaystyle \frac{5}{r}}f_q^{}A_q^{\prime \prime }+{\displaystyle \frac{3}{r}}f_q^{\prime \prime }A_q^{}{\displaystyle \frac{1}{r}}f_q^{}A_q^{}]{\displaystyle \frac{1}{f_q^2}}[3(f_q^{})^2A_q^{\prime \prime }+4f_q^{}f_q^{\prime \prime }A_q^{}+{\displaystyle \frac{4}{r}}(f_q^{})^2A_q^{}`$
$`({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(A_q^{})^2(3A_q^{\prime \prime }+{\displaystyle \frac{2}{r}}A_q^{})]+{\displaystyle \frac{3}{f_q^3}}f_q^{}A_q^{}[(f_q^{})^2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(A_q^{})^2]`$
$`{\displaystyle \frac{1}{3\rho _q}}[2\rho _q^{}A_q^{\prime \prime \prime }+\rho _q^{\prime \prime }A_q^{\prime \prime }+{\displaystyle \frac{1}{r}}\rho _q^{\prime \prime }A_q^{}+{\displaystyle \frac{5}{r}}\rho _q^{}A_q^{\prime \prime }+{\displaystyle \frac{1}{r^2}}\rho _q^{}A_q^{}]+{\displaystyle \frac{1}{3f_q\rho _q}}[3f_q^{}\rho _q^{}A_q^{\prime \prime }`$
$`+f_q^{}\rho _q^{\prime \prime }A_q^{}+2f_q^{\prime \prime }\rho _q^{}A_q^{}+{\displaystyle \frac{4}{r}}f_q^{}\rho _q^{}A_q^{}]{\displaystyle \frac{1}{3f_q^2\rho _q}}\rho _q^{}A_q^{}[3(f_q^{})^2({\displaystyle \frac{mW_0}{\mathrm{}^2}})^2(A_q^{})^2]`$
$`+{\displaystyle \frac{2}{9\rho _q^2}}(\rho _q^{})^2(A_q^{\prime \prime }+{\displaystyle \frac{1}{r}}A_q^{}){\displaystyle \frac{2}{9f_q\rho _q^2}}f_q^{}(\rho _q^{})^2A_q^{}\}`$
(A.8)
which can also be obtained directly from eq. (A.6), remembering that the divergence of a vector field $`\stackrel{}{a}`$ that has only a radial component $`a_r`$ is given by

$$\mathrm{div}\stackrel{}{a}=\frac{1}{r^2}\frac{}{r}(r^2a_r).$$


text:
# Antiferromagnetism and hole pair checkerboard in the vortex state of high Tc superconductors

## Abstract

We propose a microscopic state for the vortex phase of BSCO superconductors. Around the vortex core or above $`H_{c2}`$, the $`d`$ wave hole pairs form a checkerboard localized in the antiferromagnetic background. We discuss this theory in connection with recent STM experiments.

Since the theoretical prediction of the antiferromagnetic (AF) order in the vortex state of the high $`T_c`$ superconductors (SC)zhang1997 ; arovas1997 , tremendous experimentallake2001 ; lake2002 ; miller2001 ; mitrovic2001 ; hoffman2002 ; khaykovich2002 ; mitrovic2002  and theoreticalbruus1998 ; ogata1999 ; hu1999 ; demler2001 ; hu2001 ; lee2001 ; chen2001 ; sachdev2002 ; zhu2002 ; polkovnikov2002  progress has been made in recent years. Neutron scattering experimentskatano2000 ; lake2001 ; lake2002 ; khaykovich2002  in the $`La_{2x}Sr_xCuO_4`$ compound first reported enhanced magnetic scattering induced by the magnetic field in the vortex state. The field induced moment scales linearly (with possible logarithmic corrections) with the external field, or the number of vortices in the system, in accordance with theoretical predictionsarovas1997 ; demler2001 . The $`\mu `$SR experimentmiller2001  in the underdoped $`YBa_2Cu_3O_{7+\delta }`$ compound reveals static AF order in the vortex core. The NMR experiment.mitrovic2001 ; mitrovic2002  in the $`YBCO`$ compound under high magnetic field shows enhanced AF ordering, predominantly from the region inside the vortex core. These experiments have important implications on our theoretical understanding of high $`T_c`$ superconductivity. In the $`SO(5)`$ theory, the transition from the AF Mott insulator phase to the SC phase is understood in terms of the rotation of the $`SO(5)`$ superspin vector as the chemical potential or doping level is variedzhang1997 . Recent experiments in the AF vortex state show that the $`SO(5)`$ rotation picture can be validated in real space, as a function of the distance from the center of the vortex core, where the superspin rotates continuously from AF to SC directions.

While the salient features of the recent AF vortex experiments are in overall agreement with the theoretical predictions, there are additional experimental features worth exploring. A striking feature is revealed in the STM experiments by Hoffman et al.hoffman2002 , where the local density of states (DOS) near the vortex core show a two dimensional checkerboard-like modulation with a $`4a\times 4a`$ charge unit cell. Here $`a`$ is the dimension of the $`CuO_2`$ unit cell. This modulation decays exponentially away from the center of the vortex core, with a decay length of about $`10a`$. More recently, a similar pattern has also been seen in the absence of the applied magnetic fieldhowald2002 , possibly induced by the impurities at the surface. Understanding the four unit cell checkerboard modulation structure around the vortex core is the main focus of our present work.

One possible explanation of the STM experiments is the formation of the stripe AF orderzaanen2000 ; emery1999  around the vortex core or impurities. If the stripes are static, the question is how the 1D stripes on two different layers of the bi-layer $`Bi_2Sr_2CaCu_2O_{8+\delta }`$ material are oriented relative to each other. If we assume a strong coupling between the bi-layers, one would naturally expect the stripes on two different layers to be oriented along the same direction. This state breaks the rotational symmetry and is inconsistent with the 2D checkerboard pattern observed in the STM experiments. On the other hand, it is possible that the stripes on two different layers are oriented along the orthogonal directions. In this case, the STM tunneling should still predominantly probe the upper layer, and observe a strongly anisotropic quasi-1D pattern. Therefore, the static 1D stripe picture can be experimentally falsified as the instrumental resolution improves. Current observations could be consistent with a picture of dynamically fluctuation stripes. The question is the time scale of the fluctuation. As the magnetic field increases, these fluctuations are expected to slow down, and eventually, a static 1D ordering pattern should emerge. This scenario could be tested in the underdoped materials or in the $`Nd_{2x}Ce_xCuO_{4\delta }`$ materials, where $`H_{c2}`$ is within experimental reach. Besides the stripe state, one could also try to construct a Wigner crystal state of the doped holes in a AF background. Such an insulating state can only be realized for rational doping fractions. The optimally doped system is close to $`x=1/8`$ doping, with one hole per eight sites. In this case, the holes can form a $`\sqrt{8}a\times \sqrt{8}a`$ superlattice uniformly distributed over the AF spin background. Such a state has the full four-fold rotational symmetry as observed in the STM experiment. Unfortunately, the charge unit cell in this case is smaller than the experimentally observed $`4a\times 4a`$ structure. Finally, the simplest explanation of the checkerboard charge structure around the impurities could be the Friedel oscillation associated with the fermi surface nestinglee2002 . The observed $`4a\times 4a`$ charge pattern could be produced by the quasi-particles at the anti-node. However, these quasi-particles are gapped, and it is not clear if they can explain the low energy feature at $`7meV`$ observed in the vortex core experiment.

In this work, we propose a specific microscopic state around the vortex core and above $`H_{c2}`$ for the $`YBCO`$ and $`BSCO`$ class of superconductors. The insulating version of this state at doping $`x=1/8`$ is depicted in Fig. 1. In this state, the $`d`$ wave hole pairs form a Wigner crystal, dubbed a $`d`$ wave hole pair checkerboard. The center-of-mass of the $`d`$ wave hole pair is located at the center of every four plaquette on the original lattice. This state has the $`4a\times 4a`$ checkerboard symmetry as observed in the experiment, and the doping level for the insulating $`x=1/8`$ state is reasonably close to the optimal doping level of the cuprates. Therefore, by forming the Wigner crystal state of the hole pairs rather than the holes themselves, the doping level can be compatible with the observed size of the charge unit cell. Here the $`d`$ wave pair is drawn on the plaquette for illustrative purposes. In the actual system, the pair size can be more extended. In this Letter, we discuss the case where the AF regions are connected, therefore, antiphase domain walls would be energetically unfavorable. However, the antiphase domain walls would be possible if the hole pairs are more extended along the four-fold symmetry direction, and if there is an effective coupling among the triplet magnons on next-nearest-neighbor plaquettes. We propose that the hole pair checkerboard state is energetically competitive in the microscopic $`tJ`$ model with short ranged Coulomb interaction. This state takes advantage of the underlying AF exchange energy in the background, while the holes form singlet $`d`$ wave bonds to minimize their exchange energy as well. Finally, the kinetic energy of the hole pair competes with the Coulomb interaction. At $`x=1/8`$, if the Coulomb interaction is stronger than the kinetic energy of the hole pair, the hole pair checkerboard state could win. We believe that the $`YBCO`$ and $`BSCO`$ materials are in the regime where the kinetic energy of the hole pairs favor the uniform SC state. However, when the applied magnetic field impedes the coherent hole pair motion either around the vortex core or in the state above $`H_{c2}`$, we argue that the hole pair checkerboard state can be realized, either in the static or dynamically fluctuating form. Therefore, the central idea behind this work is that the magnetic field destroys the phase coherence of the hole pair by localizing them into a crystal, without breaking the pair.

The hole pair checkerboard state arises naturally from the projected $`SO(5)`$ model on a latticezhang1999 ; altman2001 ; dorneich2002 . In this approach, the original lattice is divided into non-overlapping plaquettes, depicted by solid squares in Fig. 1. In the spirit of real-space renormalization, we retain only five low energy states on the plaquette, the half-filled singlet state, denoted by $`|\mathrm{\Omega }`$, the triplet magnon state $`t_\alpha ^{}|\mathrm{\Omega }`$ and $`d`$ wave hole pair state $`t_h^{}|\mathrm{\Omega }`$ with 2 holes on the plaquette. The $`d`$ wave hole pair on a plaquette can be justified from the microscopic modelsscalapino1996 , but a more extended hole pair can also be incorporated into our model. We define the center of the non-overlapping plaquettes to be the lattice sites of the projected $`SO(5)`$ model, which has the spacing of $`2a`$. Here we follow the notations and the conventions of reference zhang1999 . The projected $`SO(5)`$ model is defined in terms of the following Hamiltonian.

$`H`$ $`=`$ $`\mathrm{\Delta }_s{\displaystyle \underset{x}{}}t_\alpha ^{}t_\alpha (x)+\stackrel{~}{\mathrm{\Delta }}_c{\displaystyle \underset{x}{}}t_h^{}t_h(x)`$ (1)
$``$ $`{\displaystyle \frac{J_s}{2}}{\displaystyle \underset{xx^{}}{}}(t_\alpha ^{}+t_\alpha )(x)(t_\alpha ^{}+t_\alpha )(x^{})`$
$``$ $`{\displaystyle \frac{J_c}{2}}{\displaystyle \underset{xx^{}}{}}e^{i\varphi _{ij}}t_h^{}(x)t_h(x^{})+h.c.`$
$`+`$ $`V_1{\displaystyle \underset{xx^{}}{}}\rho (x)\rho (x^{})+V_2{\displaystyle \underset{xx^{}}{}}\rho (x)\rho (x^{}),`$
where $`\rho (x)=t_h^{}t_h(x)`$ is defined as the hole-pair density.

As described in reference zhang1999 , the $`\mathrm{\Delta }_s`$ term describes the effective spin gap energy for the magnons, the $`\stackrel{~}{\mathrm{\Delta }}_c`$ term describe the effective chemical potential for the hole pairs, while the $`J_s`$ and $`J_c`$ terms describe the effective AF and SC Josephson couplings respectively. Here $`e^{i\varphi _{ij}}`$ corresponds to the phase introduced by the magnetic field in the vortex state, and we have also included the nearest neighbor and next nearest neighbor Coulomb interactions $`V_1`$ and $`V_2`$. This model can be approximately derived from the microscopic Hubbard and $`tJ`$ model by a renormalization group transformationaltman2001 . The model without the Coulomb terms have been studied extensively by numerical calculations, and the phase diagram has been determineddorneich2002 .

If one neglects the spin degree of freedom, the $`SO(5)`$ model reduces to the hard core boson model on a two dimensional square lattice. Even this simple model has a rich phase diagram as the density, the kinetic and the potential energies are variedbruder1993 ; hebert2002 . Besides the superfluid ground state, the system can also develop charge-density-wave order and supersolid order in the ground state. In the strong coupling limit at half-filling, the ground state develops checker-board charge order where every other site is occupied by a boson. At quarter-filling, a similar checker-board pattern develops where every four sites is occupied by a boson. This state is stabilized by the nearest neighbor and the next-nearest neighbor interactions $`V_1`$ and $`V_2`$. Including the AF spin order and translating the hole density and lattice size from the $`SO(5)`$ lattice model to the original lattice, we find that the $`1/4`$ filled checkerboard charge ordering of the hard-core bosons corresponds exactly to the $`d`$ wave pair checkerboard state illustrated in Fig. 1. However, we shall see that the AF exchange energies also play an important role in determining the competition between the various states.

We consider the following mean field state for the Hamiltonian (1),

$$|\mathrm{\Psi }=\underset{x}{}(\mathrm{cos}\theta +\mathrm{sin}\theta (\mathrm{cos}\alpha t_\alpha ^{}+\mathrm{sin}\alpha e^{i\varphi }t_h^{}))(x)|\mathrm{\Omega }$$
(2)
and minimize the energy $`E=\mathrm{\Psi }|H|\mathrm{\Psi }`$ as a function of the variational parameters $`\theta (x),\alpha (x)`$ and $`\varphi (x)`$. At quarter filling of the hole pair bosons, (which corresponds to $`x=1/8`$ in the actual system), the mean field theory predicts a quantum phase transition at $`V_1+V_2=V_c6J`$ for $`J_c=2J_s2J`$. For $`V_1+V_2<V_c`$, the system is in an uniform SC state. On the other hand, for $`V_1+V_2>V_c`$, a hole pair checkerboard state depicted in Fig. 1 is favored. We shall choose our parameters such that the system is in an uniform SC state, but it is close to the hole pair checkerboard state.

In order to reduce the number of variables, we have assumed four-fold rotational symmetry as well as symmetry by reflection along the diagonals. We can further simplify the problem by fixing the SC phase $`\varphi `$ assuming that the magnetic field is uniform in a region of radius equal to the London penetration depth (here, we present data with $`30a`$ but we have checked that this effect is very small).

The boundary sites are taken to be in a bulk SC phase with a fixed hole-pair density, which is taken to be the same as the averaged one. We fixed the center of the vortex core to be fully anti-ferromagnetic, i.e., $`\alpha =0`$ and $`\theta =\pi /4`$.

In Fig.2, we present the result of charge distribution which shows a clear checkerboard pattern. The contrast decreases very quickly away from the vortex core over a distance $`10a`$ which is roughly consistent with the experimental value hoffman2002 . We would like to emphasize that this checkerboard pattern persists for smaller interactions $`V_1`$ and $`V_2`$, but with a smaller contrast. For sufficiently small contrast, the modulation of the charge density may be hard to detect in the integrated LDOShoffman2002 .

We also present the profiles of the AF order parameter $`n_2`$, the SC order parameter $`n_1`$ and the $`\pi `$ order parameter $`iL_{25}/2`$zhang1997  in Fig.3. Notice that the $`\pi `$ order parameter is only non-zero in the region where the AF and SC order parameters overlap. In our projected SO(5) model, the $`SO(5)`$ orthogonality conditionzhang1997 ; hu2001

$$L_{ab}n_c+L_{bc}n_a+L_{ca}n_b=0$$
(3)
is satisfied at every site. This relation predicts a quantitative relation among these three order parameters and the hole density, $`\rho n_2+L_{25}n_1=0`$,which is satisfied by our numerical results of the order parameters. The relative size of the AF and SC regions is determined by the ratio of $`J_s/J_c`$, and can be varied to account for the experimental observations at low temperature.

We shall now discuss the experimental consequences of our theoretical proposal. Our proposal can be tested most directly in the vortex state, or above $`H_{c2}`$, by measuring the static charge structure factor corresponding to a $`4a\times 4a`$ unit cell with $`XY`$ rotational symmetry. The magnetic structure factor could either be a commensurate $`(\pi ,\pi )`$ peak with an amplitude modulation of four unit cells or a shifted peak corresponding to an anti-phase modulation. Since the energy difference between these two types of magnetic order is very small, we conjecture that both types could occur in difference compounds of the high $`T_c`$ cuprates. Our proposal differs in detail from the stripe picture from the fact that the ordering pattern has square symmetry and from the fact that the ordering charge is the hole pair rather than single holes. Our proposed state is also consistent with the $`\mu `$sR and the NMR experiments in the vortex state of the $`YBCO`$ systemmiller2001 ; mitrovic2002 . The $`\mu `$sR experiment reports commensurate $`(\pi ,\pi )`$ AF order inside the vortex coremiller2001 . Because of the modulation of the AF order parameter amplitude, oxygen sites, for examples those located in between the singlet hole pair and the AF plaquettes, can still detect a spin signalmiller2001 ; mitrovic2002 . In this work we primarily focused on the charge ordering in the vortex state. Similar arguments may be applied to the impurity states as wellhowald2002 . However, we caution here that if the impurity is pair breaking, the idea of a order state of hole pairs may not be directly applicable. We would also like to call attention to ref. kim2001 , where experimental evidence for Wigner lattice was reported.

We also encourage the detailed study  of this hole pair checkerboard state in the numerical calculations of the $`tJ`$ model, including effects of the Coulomb interactionDagotto1995 . While we do not propose this state to be realized in the ground state for physical range of parameters in the absence of the magnetic field, it is interesting to test if this state is realized in the vortex state or above $`H_{c2}`$. It is also instructive to determine parameters for which this state is realized in the ground state, and see how far these parameters are from the physically reasonable values. This way, we can estimate the perturbations required for this state to reveal itself.

We would like to acknowledge useful discussions with Drs. J. Berlinsky, S. Davis, E. Dagotto, E. Demler, W. Hanke, C. Kallin, A. Kapitulnik, S. Kivelson and D.H. Lee. This work is supported by the NSF under grant numbers DMR-9814289. HDC and JPH are also supported by a SGF, and E.A. by a Heisenberg fund of the DFG (AR 324/3-1) and by KONWIHR CUHE.


text:
# Superconductivity Induced by the Proximity Effect of Singlet Resonating Valence Bond Order in High-𝑇_𝑐 Superconductors

## Abstract

Effects of inhomogeneous doping on the high-$`T_c`$ cuprate superconductors are studied within the framework of the $`t`$-$`J`$ model. Especially, the boundary between two non-superconducting regions with doping rates much higher and lower than the optimal one is examined. It is found that, although there is no superconductivity in the bulk, a superconducting region appears at the boundary because singlet resonating valence bond order and holon condensation occur simultaneously in this region. The critical temperature of the induced superconductivity can be higher than that of the optimally doped sample and the critical current density can be the same order as that in the superconducting state below the critical temperature. We also point out an experimental possibility to observe this phenomenon.

preprint: APS/123-QED

It has been recognized that there are two important ingredients for the superconductivity in high-$`T_c`$ cuprates, which are the spin correlation, leading to the spin gap, and the doped charge carriers with a sufficiently high density. In the slave-boson mean-field theory of the $`t`$-$`J`$ model, the former is realized as the singlet resonating valence bond (s-RVB) order of spinons and the latter the bose condensation of holons Anderson ; Baskaran-Zou-Anderson ; Anderson-Zou ; Isawa-Maekawa-Ebisawa ; Anderson-Baskaran-Zou-Hsu ; Suzumura-Hasegawa-Fukuyama ; Kotlier-Liu . In this framework, it is assumed that the spinons and the holons are well-defined low lying excitations, and the spin-charge separation is a good picture. Although there still remains some controversy on this point, it has been clarified that some physical properties of the high-$`T_c`$ cuprates can be well described by the $`t`$-$`J`$ model, for example, the transport properties Nagaosa-Lee  and the magnetic properties Tanamoto-Kuboki-Fukuyama ; Tanamoto-Kohno-Fukuyama . In this Letter, we take the $`t`$-$`J`$ model to be our starting point.

We consider the situation where the doping rate is spatially varying. We especially concentrate on the case where two non-superconducting regions, one is too over-doped and the other is too under-doped, are in contact with each other within the same CuO2 plane. We expect that the s-RVB order and the holon condensate can exist simultaneously at the boundary, thus giving rise to superconductivity. Since the s-RVB order and the holon condensation can occur at temperatures higher than the critical temperature of the optimally doped sample, $`T_{\mathrm{opt}}`$, this “boundary superconductivity”may also occur above $`T_{\mathrm{opt}}`$.

Our idea is similar with the “spin-gap proximity effect”previously introduced by Emery, Kivelson and Zacher Emery-Kivelson-Zacher  based on a somewhat different model for high $`T_c`$ superconductors, discussing the stripe order and its significance for the superconductivity. Recently, Kivelson has proposed a possibility to raise critical temperature using spin-gap proximity effect Kivelson . He proposes making a stack of spin-gap material and hole-rich material in a layered structure. In contrast to this, we consider the situation where the spin-gap and hole-rich regions exist within a CuO2 plane. This kind of situation has not been studied within the framework of the $`t`$-$`J`$ model and we consider this is of interest from both theoretical and experimental points of view.

First we introduce the $`t`$-$`J`$ model which describes the correlated electrons in each CuO2 plane. The Hamiltonian is given by

$``$ $`=`$ $`t{\displaystyle \underset{i,j\sigma }{}}(f_{i\sigma }^{}f_{j\sigma }b_j^{}b_i+\mathrm{h}.\mathrm{c}.)+J{\displaystyle \underset{i,j}{}}\stackrel{}{S}_i\stackrel{}{S}_j`$ (1)
$`{\displaystyle \underset{i}{}}\lambda _i({\displaystyle \underset{\sigma }{}}f_{i\sigma }^{}f_{i\sigma }+b_i^{}b_i1),`$
where $`i`$ and $`j`$ denote lattice points in the CuO2 plane, $`i,j`$ means the nearest neighbors, $`f_{i\sigma }`$ ($`\sigma =,`$) and $`b_i`$ are the annihilation operators of the spinon and the holon, respectively, $`\stackrel{}{S}_i`$ denotes $`\frac{1}{2}_{\alpha \beta }f_{i\alpha }^{}\stackrel{}{\sigma }_{\alpha \beta }f_{i\beta }`$ with $`\stackrel{}{\sigma }_{\alpha \beta }`$ being Pauli matrices, and $`\lambda _i`$ is the Lagrange multiplier for the constraint, $`_\sigma f_{i\sigma }^{}f_{i\sigma }+b_i^{}b_i=1`$. In this Letter the lattice spacing is taken to be unity. In addition to Eq. (1), we later introduce the coupling to the electromagnetic field.

The Eq. (1) is treated by the so-called slave-boson mean-field theory. The local constraint represented by $`\lambda _i`$ is replaced by the global one, which is described by the constant chemical potential for spinons $`\mu _\mathrm{F}`$ and that for holons $`\mu _\mathrm{B}`$. We define the statistical averages, $`\xi _\sigma f_{i\sigma }^{}f_{j\sigma }`$ and $`\eta b_i^{}b_j`$, which are assumed to be independent of $`i`$ and $`j`$ average . The order parameter for the s-RVB order is denoted as $`\mathrm{\Delta }_{ij}=\frac{3J}{8}f_if_jf_if_j`$, which is equal to $`\mathrm{\Delta }_x(\stackrel{}{r}_i)`$ if $`j=i+\widehat{x}`$ and $`\mathrm{\Delta }_y(\stackrel{}{r}_i)`$ if $`j=i+\widehat{y}`$ where $`\widehat{x}`$ and $`\widehat{y}`$ indicate the shift by one lattice spacing in positive $`x`$\- and $`y`$-direction, respectively. Although $`\xi `$, $`\eta `$ and $`\mathrm{\Delta }_{ij}`$ should be determined self-consistently in principle, we skip this process for simplicity by just assuming $`\xi `$ and $`\eta `$ to be constant. Including also the scalar and vector potential, $`\phi `$ and $`\stackrel{}{A}`$, and the fictitious gauge fields, $`a_0`$ and $`\stackrel{}{a}`$, which represent the fluctuation around the mean-field solution Nagaosa-Lee , we obtain the following Lagrangian,

$``$ $`=`$ $`{\displaystyle \underset{i\sigma }{}}f_{i\sigma }^{}\left\{\mathrm{}_\tau \mu _\mathrm{F}+ia_0(\stackrel{}{r}_i)\right\}f_{i\sigma }`$ (2)
$`+{\displaystyle \underset{i}{}}b_i^{}\left\{\mathrm{}_\tau \mu _\mathrm{B}+i(e\phi (\stackrel{}{r}_i)+a_0(\stackrel{}{r}_i))\right\}b_i`$
$`+{\displaystyle \underset{i}{}}i\phi (\stackrel{}{r}_i)\rho (\stackrel{}{r}_i)`$
$`(t\eta +{\displaystyle \frac{3}{8}}J\xi ){\displaystyle \underset{ij\sigma }{}}\{\mathrm{e}^{i\theta _{ji}^\mathrm{F}}f_{i\sigma }^{}f_{j\sigma }+\mathrm{h}.\mathrm{c}.\}`$
$`t\xi {\displaystyle \underset{ij}{}}\{\mathrm{e}^{i\theta _{ji}^\mathrm{B}}b_i^{}b_j+\mathrm{h}.\mathrm{c}.\}`$
$`{\displaystyle \underset{i}{}}[\mathrm{\Delta }_x^{}(\stackrel{}{r}_i)\chi _{ii+\widehat{x}}+\mathrm{\Delta }_y^{}(\stackrel{}{r}_i)\chi _{ii+\widehat{y}}+\mathrm{h}.\mathrm{c}.`$
$`{\displaystyle \frac{8}{3J}}(|\mathrm{\Delta }_x(\stackrel{}{r}_i)|^2+|\mathrm{\Delta }_y(\stackrel{}{r}_i)|^2)],`$
where $`\chi _{ij}=f_if_jf_if_j`$ and

$$\theta _{ji}^\mathrm{B}=\frac{1}{\mathrm{}}_{\stackrel{}{r}_j}^{\stackrel{}{r}_i}\stackrel{}{a}d\stackrel{}{l},\theta _{ji}^\mathrm{F}=\frac{1}{\mathrm{}}_{\stackrel{}{r}_j}^{\stackrel{}{r}_i}\left(\stackrel{}{a}+\frac{e}{c}\stackrel{}{A}\right)d\stackrel{}{l}.$$
(3)
The chemical potentials $`\mu _\mathrm{F}`$ and $`\mu _\mathrm{B}`$ are determined from

$$f_{i\sigma }^{}f_{i\sigma }=1\overline{\delta },b_i^{}b_i=\overline{\delta },$$
(4)
where $`\overline{\delta }`$ is the average holon density in the whole system. The background charge density, originating from doping, is denoted by $`\rho (\stackrel{}{r}_i)`$, whose spatial average equals $`e\overline{\delta }`$. In this Letter $`a_0`$ and $`\stackrel{}{a}`$ are treated perturbatively on the same line with Refs. Nagaosa-Lee  and Nagaosa-Lee2 .

First we study how the distribution of the holons, namely $`\delta (\stackrel{}{r}_i)`$, is determined under an inhomogeneous doping, which is described by $`\rho (\stackrel{}{r}_i)`$. For this purpose, we neglect the terms including $`\mathrm{\Delta }_{ij}`$ and $`\stackrel{}{a}`$ and $`\stackrel{}{A}`$. Then the spinon and the holon degrees of freedom can be integrated out, giving rise to terms, $`\pi _\mathrm{F}a_0^2+\pi _\mathrm{B}(a_0+e\phi )^2`$, where $`\pi _\mathrm{F}`$ and $`\pi _\mathrm{B}`$, respectively, are the polarization functions of spinons and holons. Now the magnitudes of $`\pi _\mathrm{F}`$ and $`\pi _\mathrm{B}`$ determine which of the spinon and the holon couples to $`\phi `$. In this Letter we assume $`\pi _\mathrm{F}\pi _\mathrm{B}`$, which may be valid in lower doping region. Then $`a_0`$ is almost fixed to zero and $`\phi `$ couples to the holons. Now $`\delta `$ and $`\phi `$ are determined from $`\rho `$ in a self-consistent way.

We especially consider the following situation. The CuO2 planes, described by the above Lagrangian, are laid parallel to $`x`$-$`y`$ plane and stacked in $`z`$-direction. The region $`x<0`$ ($`x>0`$) is filled with a s-RVB ordered (hole-rich) material with a low (high) enough doping rate. Both of the regions are non-superconducting. The system is uniform in $`y`$\- and $`z`$-direction. We assume that the doping rate, which is described by $`\rho (\stackrel{}{r}_i)`$, changes at $`x=0`$ abruptly. In order to preserve the charge neutrality in the bulk, namely at $`x\pm \mathrm{}`$, the statistical average $`\phi (\stackrel{}{r})`$ must take constant imaginary values, which differ between $`x\mathrm{}`$ and $`x\mathrm{}`$. (Note that the chemical potential $`\mu _\mathrm{B}`$ is a constant in the whole system.) Near the boundary $`\phi (\stackrel{}{r})`$ changes smoothly connecting these two bulk limits. As in the case of semiconductor junctions Sze , the smooth connection is enabled by the polarization charge appearing at the boundary. The depth of the charged region is given by the screening length of the holons. This length is extremely long (short) for $`x<0`$ ($`x>0`$) and there still remains a sharp drop of the holon density at the boundary. Therefore, in this Letter, we approximate $`\delta `$ by a step function (see Fig. 1).

Next we study the s-RVB order by introducing the Ginzburg-Landau (GL) free energy of the order parameter $`\mathrm{\Delta }_d=\frac{1}{2}(\mathrm{\Delta }_x\mathrm{\Delta }_y)`$. Here we assume that the s-RVB order has $`d`$-wave symmetry. The GL free energy is obtained by the perturbative expansion of the free nergy with respect to $`\mathrm{\Delta }_d`$ Feder-Kallin . We introduce the half band width $`D=4(t\eta +3J\xi /8)`$ and the strength of attractive interaction $`V=3J/8`$. The dispersion of spinons, $`(D/2)(\mathrm{cos}k_x+\mathrm{cos}k_y)`$, is approximated by a parabolic one, $`D(1+|\stackrel{}{k}|^2/4)`$, and then the density of states at Fermi energy becomes $`N(0)=1/(\pi D)`$. At a temperature $`TD`$, the solution for $`_\sigma f_{i\sigma }^{}f_{i\sigma }=1\delta `$ is approximately given by $`\mu _\mathrm{F}=\delta /N(0)`$ and the critical temperature $`T_d`$ becomes $`T_dT^{}\mathrm{exp}[1/\{2VN(0)(1\delta )^2\}]`$ with $`T^{}=2D\mathrm{e}^\gamma /\pi `$ where $`\gamma `$ is the Euler’s constant. This expression is valid in the “weak coupling regime”, $`VN(0)1`$. The GL free energy density of $`\mathrm{\Delta }_d`$ is given by

$`F_{\mathrm{s}\mathrm{RVB}}=F_0+\alpha _d|\mathrm{\Delta }_d|^2+\beta |\mathrm{\Delta }_d|^4+\gamma _d|\stackrel{}{\mathrm{\Pi }}\mathrm{\Delta }_d|^2,`$ (5)
where $`\stackrel{}{\mathrm{\Pi }}=i+(2/\mathrm{})\stackrel{}{a}+(2\pi /\varphi _0)\stackrel{}{A}`$ with $`\varphi _0=hc/(2e)`$. The GL coefficients are given by

$`\alpha _d`$ $`=`$ $`2N(0)(1\delta )^2\mathrm{ln}{\displaystyle \frac{T}{T_d}},`$ (6)
$`\beta `$ $`=`$ $`{\displaystyle \frac{21\zeta (3)N(0)}{2\pi ^2T^2}}(1\delta )^4{\displaystyle \frac{c_1}{T^2D}},`$ (7)
$`\gamma _d`$ $`=`$ $`{\displaystyle \frac{7\zeta (3)N(0)D^2}{16\pi ^2T^2}}{\displaystyle \frac{c_2D}{T^2}},`$ (8)
where $`\zeta (x)`$ is the zeta function and $`c_1`$ and $`c_2`$ are constants. The doping dependence of the GL coefficients comes from $`T_d`$, $`\alpha _d`$ and $`\beta `$. Since the largest dependence appears from $`T_d`$, we neglect the $`(1\delta )^n`$-factors in $`\alpha _d`$ and $`\beta `$ in the following. Then the doping rate affects $`T_d`$ only.

Here we consider the situation where the transition temperatures in the regions $`x<0`$ and $`x>0`$, denoted by $`T_d`$ and $`T_d^{}`$, respectively, satisfies the condition, $`T_d>TT_d^{}`$. It is also assumed that holons are condensed in the region $`x>0`$. Note that the quantity with (without) dash () is defined in $`x>0`$ ($`x<0`$). We first determine the spatial variation of $`\mathrm{\Delta }_d`$, disregarding $`\stackrel{}{a}`$ and $`\stackrel{}{A}`$. Here $`\mathrm{\Delta }_d`$ is assumed to be real. The expectation value of $`\mathrm{\Delta }_d`$ is $`\sqrt{|\alpha _d|/(2\beta )}\overline{\mathrm{\Delta }}_d`$ at $`x\mathrm{}`$ and $`0`$ at $`x\mathrm{}`$. In $`x<0`$ we introduce a new variable $`\delta \mathrm{\Delta }_d=\mathrm{\Delta }_d\overline{\mathrm{\Delta }}_d`$. The spatial variation of $`\mathrm{\Delta }_d`$ is then governed by the following equations,

$`2|\alpha _d|\delta \mathrm{\Delta }_d(x)\gamma _d_x^2\delta \mathrm{\Delta }_d(x)=0(x<0),`$
$`\alpha _d^{}\mathrm{\Delta }_d(x)\gamma _d_x^2\mathrm{\Delta }_d(x)=0(x>0).`$ (9)
We can easily see that $`\mathrm{\Delta }_d`$ behaves as $`\overline{\mathrm{\Delta }}_d+g_1\mathrm{exp}(x/\xi _\mathrm{F})`$ for $`x<0`$ and $`g_2\mathrm{exp}(x/\xi _\mathrm{F}^{})`$ for $`x>0`$, where $`\xi _\mathrm{F}=\sqrt{\gamma _d/2|\alpha _d|}`$, $`\xi _\mathrm{F}^{}=\sqrt{\gamma _d/\alpha _d^{}}`$, and $`g_1`$ and $`g_2`$ are constants. The solution in the whole region is obtained by connecting these two solutions at $`x=0`$ with boundary condition,

$$_x\mathrm{\Delta }_d(0)=_x\mathrm{\Delta }_d(+0),\mathrm{\Delta }_d(0)=\mathrm{\Delta }_d(+0),$$
(10)
(see Fig. 1). This condition is in contrast with that for the proximity effect in ordinary superconductor-normal metal junction, where the change of the BCS interaction causes a discontinuity of the gap function deGennes . (Note that the quantity which is continuous at the boundary is not the gap function but the anomalous amplitude.) The constants $`g_1`$ and $`g_2`$ are obtained as

$$g_1=\frac{\xi _\mathrm{F}}{\xi _\mathrm{F}+\xi _\mathrm{F}^{}}\overline{\mathrm{\Delta }}_d,g_2=\frac{\xi _\mathrm{F}^{}}{\xi _\mathrm{F}+\xi _\mathrm{F}^{}}\overline{\mathrm{\Delta }}_d.$$
(11)
Note that non-zero $`\mathrm{\Delta }_d`$ in $`x>0`$ is due to the proximity effect of s-RVB order and $`\xi _\mathrm{F}^{}`$ gives the proximity length. The temperature dependence of the coherence lengths at $`T_d^{}TT_d`$ are as follows:

$`\xi _\mathrm{F}\sqrt{{\displaystyle \frac{c_2}{2}}}{\displaystyle \frac{D}{T_d}}{\displaystyle \frac{1}{\sqrt{1T/T_d}}},\xi _\mathrm{F}^{}={\displaystyle \frac{\sqrt{c_2}}{\sqrt{\mathrm{ln}T/T_d^{}}}}{\displaystyle \frac{D}{T}}.`$
It is interesting to see that $`\xi _\mathrm{F}^{}`$ is, except for the logarithmic factor, analogous to the ordinary proximity length in the clean limit $`\mathrm{}v_\mathrm{F}/T`$ where $`v_\mathrm{F}`$ is the Fermi velocity. If $`T`$ is close to $`T_d`$, $`\xi _\mathrm{F}`$ becomes much larger than $`\xi _\mathrm{F}^{}`$ and the amplitude of $`\mathrm{\Delta }_d`$ in the proximity region, i.e., $`g_2`$, becomes much smaller than $`\overline{\mathrm{\Delta }}_d`$. However if $`T`$ is not too close to $`T_d`$, $`\xi _\mathrm{F}`$ and $`\xi _\mathrm{F}^{}`$ can be the same order and $`g_2`$ becomes comparable to $`\overline{\mathrm{\Delta }}_d`$. Noting that the holons are already bose condensed in the proximity region, this result means the induced superconductivity in the boundary area.

The strength of the induced superconductivity becomes more apparent by estimating the supercurrent which can flow parallel to the boundary. In this region the gauge field $`\stackrel{}{a}`$ is massive because of the holon condensation, and $`\stackrel{}{a}`$ in Eq. (5) can be set to zero. Then the response to the actual vector potential $`\stackrel{}{A}(\stackrel{}{r})`$ is dominated by the spinons. Here we take the vector potential to be constant and parallel to $`y`$-direction. The total current in $`y`$-direction is given by

$`J`$ $`=`$ $`2c\gamma _dA_y{\displaystyle _0^{\mathrm{}}}\mathrm{\Delta }_d^2(x)dx`$ (12)
$`=`$ $`c\gamma _dA_y{\displaystyle \frac{\{\stackrel{~}{\xi }_\mathrm{F}^{}(A_y)\}^3}{\{\stackrel{~}{\xi }_\mathrm{F}^{}(A_y)+\xi _d\}^2}}\overline{\mathrm{\Delta }}_d^2,`$
where the $`A_y`$-dependent proximity length $`\stackrel{~}{\xi }_\mathrm{F}^{}(A_y)=\xi _\mathrm{F}^{}/\sqrt{1+p^{}A_y^2}`$ decreases with increasing $`A_y`$, with $`p^{}`$ being $`(2\pi \xi _\mathrm{F}^{}/\varphi _0)^2`$. The maximum of $`J`$ as a function of $`A_y`$ is given approximately by

$$J_{\mathrm{max}}=2c\left(\frac{2\pi }{\varphi _0}\right)^2\gamma _d\frac{1}{3\sqrt{3p^{}}}\frac{\xi _{\mathrm{F}}^{}{}_{}{}^{3}}{\xi _{\mathrm{F}}^{}{}_{}{}^{2}}\overline{\mathrm{\Delta }}_d^2.$$
(13)
For clarity, here we introduce the critical current density $`j_c`$ in the bulk superconducting phase. This quantity (in overdoped region) is calculated from Eq. (5) by assuming that holons are condensed and $`\stackrel{}{a}`$ is massive in the whole region. We obtain

$$j_c=2c\left(\frac{2\pi }{\varphi _0}\right)^2\gamma _d\frac{1}{3\sqrt{3p}}\overline{\mathrm{\Delta }}_d^2,$$
(14)
where $`p=2(2\pi \xi _\mathrm{F}/\varphi _0)^2`$. It is clear that $`J_{\mathrm{max}}=\sqrt{2}j_c\xi _{\mathrm{F}}^{}{}_{}{}^{2}/\xi _\mathrm{F}`$. If $`\xi _\mathrm{F}^{}\xi _\mathrm{F}`$, the critical current density in the proximity region, defined by $`J_{\mathrm{max}}/\xi _\mathrm{F}^{}`$, can be the same order as $`j_c`$.

Finally we point out some experimental set-ups to observe the effect predicted in this Letter. The simplest way may be to utilize the field effect transistor (FET) recently developed by Shön et al. Shoen1 ; Shoen3 . We propose two types of experiments, which we call “half-gate”experiment and “side-gate”experiment. The experimental set-ups are depicted in Fig. 2. In the former case (Fig. 2(a)) the gate electrode is attached parallel to the CuO2 plane, but it covers only half of it. As a result, a boundary between the s-RVB ordered region (uncovered by the gate) and the hole-rich region (covered by the gate) is formed under the edge of the gate. (Here we are assuming that the gate induces doping to the undoped compound.) However, in order to observe the effect studied in this Letter, the insulating layer separating the gate and the sample must be smaller than the s-RVB coherence length (approximately the superconducting coherence length below $`T_c`$ $``$ several nm.), which may be hardly realized.

The other configuration is the “side-gate”as depicted in Fig. 2(b). In this case the gate is applied to the CuO2 plane from the side and by applying a voltage a hole-rich layer is formed under the gate. This situation is a little different from one considered in this Letter where the hole-rich region is infinitely wide. However the same mechanism studied by us may also work. In this case, the thickness of the superconducting region is limited not by the s-RVB proximity length but by the thickness of the hole-rich region. This configuration can be realized more easily.

There may be other ways. One might introduce the spatial variation of doping chemically, for example, by using the epitaxial growth Terashima . If this can be realized, it may fit the situation studied in this Letter best.

In either case, the width of the superconducting region is very narrow and the fluctuation is essentialKumagai . Therefore it is not easy to observe good superconducting behavior. However we expect that whether the superconducting correlation exists or not is an experimentally detectable fact.

In summary, we have studied the effect of inhomogeneous doping in the high-$`T_c`$ superconductors. We have shown, based on the $`t`$-$`J`$ model, a possibility of superconductivity at the boundary between the s-RVB ordered and the hole-rich regions, caused by the proximity effect of the s-RVB order. We have also suggested experimental set-ups to observe this effect. Since our proposal is based on the spin-charge separation picture, these experiments may also shed a renewed light on the mechanism of high-$`T_c`$ superconductivity. On the other hand, from a theoretical point of view, the extension of the present study to more general inhomogeneous cases is also important in understanding recently found inhomogeneous phenomena in cuprates Iguchi ; Xu ; Lang .

The authors are grateful to K. Kuboki, H. Kohno and M. Mori for many useful discussions. This research is supported by a Grant-In-Aid for Scientific Research (1246207) from the Ministry of Education, Science, Sports and Culture, Japan.