text:
# A Model for the 𝜆-Transition of Helium†

## I Introduction

In 1938 London presented his famous suggestion about the connection between the $`\lambda `$-transition in 4He and the Bose-Einstein condensation in the ideal Bose gas (IBG). About ten years later this point of view was confirmed experimentally by showing the absence of such a transition in 3He. Theoretically it was supported by Feynman who argued that London’s view was essentially correct. Furthermore, the microscopic IBG provides a means of understanding some properties of He II. The actual properties of helium are, however, in many respects quite different from the IBG. In view of this we propose a model which is a modification of the IBG invented such that relevant equilibrium properties of the real system are reproduced.

The IBG can be defined by the microscopic states

$$\mathrm{\Psi }_{\mathrm{IBG}}(𝒓_j,n_k)=𝒮\underset{𝒌}{}\left[\phi _𝒌\right]^{n_k},$$
(1)
which are the symmetrized ($`𝒮`$) product states of single-particle functions (s.p.f.) $`\phi _𝒌`$ of momentum $`\mathrm{}𝒌`$. The states (1) depend on the coordinates $`𝒓_j`$ (where $`j=1,\mathrm{},N`$) of the atoms and on the parameters $`n_k=n_𝒌`$. The temperature dependent expectation values $`n_k`$ of the occupation numbers display the phase transition of the IBG.

Our model is based on two major assumptions:

1. It starts with an IBG-related guess for the form of the many-body states. As in Eq. (1), these states depend on occupation numbers $`n_k`$, but beyond Eq. (1) they contain correlations due to the use of localized and phase ordered single-particle functions.
2. The phase transition is introduced phenomenologically by using the IBG values for $`n_k`$. Due to this feature our model is closely related to the IBG and will therefore be called almost ideal Bose gas model, or AIBG.

The basic idea of the AIBG for introducing the relevant correlations is the following: Normally exponential functions $`\phi _𝒌\mathrm{exp}(\mathrm{i}𝒌𝒓_j)`$ are used for the s.p.f. in Eq. (1). Alternatively one may consider sinus functions which depend on the component $`x_j=\widehat{𝒙}𝒓_j`$ of $`𝒓_j`$ as follows:

$$\phi _𝒌\mathrm{sin}(qx_j+\varphi _j).$$
(2)
Here $`\widehat{𝒙}`$ is a unit vector in an arbitrary direction, and $`q=\widehat{𝒙}𝒌`$. For these s.p.f. we introduce the concept of phase ordering (p.o.):

$$\varphi _j=\varphi _0\text{(phase ordering)}.$$
(3)
Phase ordered (p.o.) s.p.f. with the same $`q`$ (but in general different $`𝒌`$) are correlated because their squares $`|\phi _𝒌(𝒓_j)|^2`$ exhibit common extrema. Thus p.o. leads to an extra spatial correlation between all atoms with the same $`q`$. Such correlations may be important if many atoms share the same momentum. Therefore p.o. is a potentially decisive correlation near the transition point in an IBG-like model.

A p.o. may be introduced by physical (or Dirichlet) boundary conditions at the walls of the macroscopic volume $`V`$. Then the resulting effects are surface effects and vanish like $`V^{1/3}`$ for $`V\mathrm{}`$. In contrast to this, we obtain finite effects by assuming s.p.f. $`\phi _𝒌(\mathrm{l}.\mathrm{p}.\mathrm{o}.)`$ which are not only phase ordered (p.o.) but also localized (l.). For this purpose the volume $`V`$ is thought to be divided into $`V/V_0`$ finite boxes of size $`V_0`$. The $`\phi _𝒌(\mathrm{l}.\mathrm{p}.\mathrm{o})`$ are then the s.p.f. which are confined to one of these boxes and subject to Dirichlet boundary conditions. For these s.p.f. we construct a state of the form (1):

$$\mathrm{\Psi }_\mathrm{M}(𝒓_j,n_k,N_0)=𝒮\underset{\mathrm{vol}}{}\underset{𝒌}{}\left[\phi _𝒌(\mathrm{l}.\mathrm{p}.\mathrm{o}.)\right]^{n_k}.$$
(4)
The product runs over all momenta and all finite volumes. The states $`\mathrm{\Psi }_\mathrm{M}`$ depend on the additional parameter $`N_0`$ which denotes the number of atoms in $`V_0`$. For $`N_0\mathrm{}`$ these states reduce to IBG states. As a surface effect of the finite boxes the l.p.o. introduces correlations of the order $`N_0^{1/3}`$.

As it stands, the states $`\mathrm{\Psi }_\mathrm{M}`$ cannot be used together with a realistic Hamiltonian $`H`$. Therefore, we multiply them by a suitable Jastrow factor,

$$\mathrm{\Psi }(𝒓_j,n_k,N_0)=F\mathrm{\Psi }_\mathrm{M},F=\underset{i<j}{}f(r_{ij}).$$
(5)
We use Jastrow functions $`f(|𝒓_i𝒓_j|)`$ which have been determined by minimizing the energy $`E_0=F|H|F`$.

We have now sketched the basic idea of the AIBG, namely the modification of the IBG due to a phase ordering of the s.p.f. The resulting model is rather close to the IBG; in spite of some formal similarity of the underlying states it is not related to the quasi-particle model (note that $`n_k=N`$ for Eq. (5)). In Sec. II the AIBG will be defined in detail. This includes the evaluation of the energy

$$E(n_k)=F\mathrm{\Psi }_\mathrm{M}\left|H\right|F\mathrm{\Psi }_\mathrm{M}$$
(6)
Evaluated with $`n_k_{\mathrm{IBG}}`$ this energy yields a logarithmic singularity (Sec. III).

Since our approach is based on the IBG (above point 2), the critical exponent $`\beta `$ of the model condensate density has the value $`1/2`$. This differs from the critical exponent $`\nu 1/3`$ of superfluid density $`\rho _\mathrm{s}`$. Therefore, $`\rho _\mathrm{s}`$ cannot be identified with the square of the condensate wave function. In the AIBG we argue that due to p.o. also noncondensed particles contribute to $`\rho _\mathrm{s}`$; in any IBG-like model ($`\beta =1/2`$) such a contribution will be required in order to account for the experimental value $`\nu 1/3`$. This contribution leads to the peculiar consequence that the superfluid component has a nonvanishing entropy $`S_\mathrm{s}0`$. A prediction for $`S_\mathrm{s}0`$ will either falsify the model right away, or provide a crucial test of the model. Therefore we extend the model in Sec. V to VI such that $`\rho _\mathrm{s}`$ and $`S_\mathrm{s}`$ can be calculated.

Section IV shows how the localization of the s.p.f. can be reconciled with a macroscopic p.o. and introduces generalized s.p.f. which are able to describe the coherent motion of noncondensed particles with the condensate. Section V presents the model expression for $`\rho _\mathrm{s}`$ and shows that it yields excellent fits to the experimental data. With the model parameters fixed by this fit we obtain in Sec. VI the crucial model prediction for $`S_\mathrm{s}0`$. This prediction is compared to the available experimental data; it is found that the prediction is at the border of present-day experimental detectability.

## II The almost ideal Bose gas model (AIBG)

For a complete definition of our model we define the l.p.o. s.p.f. (Sec. II A), calculate the energy $`E`$ (Sec. II B) and specify the assumptions about the expectation values of the parameters (Sec. II C).

### A Localized, phase ordered single-particle functions

The l.p.o. s.p.f. in Eq. (4) are formally defined as follows: The macroscopic volume $`V`$ is divided into $`V/V_0`$ cubic boxes of size $`V_0`$. Within each box we define orthonormalized s.p.f. which obey Dirichlet boundary conditions at the walls:

$$\phi _𝒌(\mathrm{l}.\mathrm{p}.\mathrm{o}.)=\{\begin{array}{ccc}\sqrt{\frac{8}{V_0}}\underset{i=1}{\overset{3}{}}\mathrm{sin}(k_ix_i)& & \text{for }𝒓V_0,\hfill \\ 0& & \text{elsewhere}.\hfill \end{array}$$
(7)
The Cartesian components $`k_i`$ of $`𝒌`$ are restricted to the discrete values $`q_n`$,

$$q_n=n\mathrm{\Delta }k,n=1,2,\mathrm{},\mathrm{\Delta }k=\pi /V_0^{\mathrm{\hspace{0.17em}1}/3}.$$
(8)
The many-body states of the AIBG are now defined by Eq. (5) with Eqs. (4) and (7).

This construction is chosen for simplicity. It contains artificial aspects like the cubic shape and the identical size of the boxes. The underlying physical picture is that there are finite regions (of various shapes and sizes) in which the atoms lower their free energy (see Sec. IV A) by assuming p.o. in just one direction. For macroscopic physical quantities we have to average (explicitly or implicitly) over such regions, or over the boxes used for Eq. (7). Formally we perform this averaging by replacing the discrete momentum sums by integrals:

$$\underset{q_n}{}\mathrm{}\underset{q}{}\mathrm{}=\frac{1}{\mathrm{\Delta }k}_0^{\mathrm{}}𝑑q\mathrm{}.$$
(9)
The justification of the lower integral bound will be discussed in Sec. IV B. The averaging over boxes of finite size and shape restores also the translational and rotational invariance of the system as a whole. The occupation number for $`\phi _𝒌`$ will therefore depend only on $`k`$ but not on the direction of $`𝒌`$.

### B Energy

We evaluate the energy Eq. (6) for the Hamiltonian

$$H=\underset{i}{}\frac{\mathrm{}^2}{2m}\mathrm{\Delta }_i+\underset{i<j}{}u(r_{ij})$$
(10)
with a realistic atom-atom potential $`u(r)`$, for example of Lennard-Jones type.

The ground state (g.s.) may be approximated by $`\mathrm{\Psi }_0=F`$. The corresponding g.s. energy $`E_0`$ can be written as

$$E_0(V,N)=F|H|F=\frac{N^2}{2V}d^3r\stackrel{~}{u}(r)g_0(r),$$
(11)
where $`g_0(r)`$ is the pair correlation function for $`\mathrm{\Psi }_0`$ and $`\stackrel{~}{u}(r)`$ is given by

$$\stackrel{~}{u}(r)=u(r)\frac{\mathrm{}^2}{2m}\mathrm{ln}f(r).$$
(12)
Sensible results for $`E_0`$ and $`g_0(r)`$ are obtained by a simple ansatz for $`f(r)`$.

The following evaluation of the energy $`E`$ for $`\mathrm{\Psi }=F\mathrm{\Psi }_\mathrm{M}`$ is based on the fact that the main correlations are determined by $`F`$. The additional correlations due to $`\mathrm{\Psi }_\mathrm{M}1`$ are small, they are of the order $`N_0^{1/3}1`$ for the relevant p.o. effects. The pair correlation function $`g(r)`$ for $`\mathrm{\Psi }=F\mathrm{\Psi }_\mathrm{M}`$ can therefore be written as

$$g(r,n_k)=g_0(r)+\mathrm{\Delta }g(r,n_k),\mathrm{\Delta }g1.$$
(13)
We relate $`\mathrm{\Delta }g`$ to the pair correlation function $`g_\mathrm{M}`$ of $`\mathrm{\Psi }_\mathrm{M}`$. (An explicit expression for $`g_\mathrm{M}`$ and details of the following derivation are given in Ref. ). Writing $`g_\mathrm{M}=1+\mathrm{\Delta }g_\mathrm{M}`$ it follows that $`\mathrm{\Delta }g0`$ for $`\mathrm{\Delta }g_\mathrm{M}0`$, or $`\mathrm{\Delta }g\mathrm{\Delta }g_\mathrm{M}`$. To lowest order in $`f(r)`$ Eq. (5) yields

$$\mathrm{\Delta }g(r,n_k)f(r)^2\mathrm{\Delta }g_\mathrm{M}(r,n_k).$$
(14)
An alternative approximation is $`\mathrm{\Delta }gg_0\mathrm{\Delta }g_\mathrm{M}`$; the difference to Eq. (14) marks the uncertainty in this approximation. For our purpose both approximations will turn out to be nearly equivalent. The relation between $`\mathrm{\Delta }g`$ and $`\mathrm{\Delta }g_\mathrm{M}`$ cannot be evaluated exactly; for a discussion of this problem see Ref. .

The potential energy for $`F\mathrm{\Psi }_\mathrm{M}`$ is $`E_{\mathrm{pot}}=(N^2/V)`$ $`d^3rgu`$. The kinetic energy can be divided into two parts: The Jastrow factors are taken into account by replacing $`u`$ in $`E_{\mathrm{pot}}`$ by $`\stackrel{~}{u}`$, as in Eq. (11). The kinetic contributions due to $`\mathrm{\Psi }_\mathrm{M}1`$ lead to $`\epsilon _kn_k`$ where $`\epsilon _k=\mathrm{}^2k^2/2m`$. Here $`m`$ is an effective mass (Ref. ) but for simplicity we do not introduce a new symbol. We emphasize that $`m0`$ for $`k0`$. This can be related to the fact that the static structure function $`S(k,T)`$ is finite for $`k0`$ and $`TT_\lambda `$. The complete result reads

$`E(V,N,n_k)`$ $`=`$ $`F\mathrm{\Psi }_\mathrm{M}\left|H\right|F\mathrm{\Psi }_\mathrm{M}`$ (15)
$`=`$ $`{\displaystyle \frac{N^2}{2V}}{\displaystyle d^3rg(r,n_k)\stackrel{~}{u}(r)}+{\displaystyle \underset{𝒌}{}}\epsilon _kn_k`$ (16)
$`=`$ $`E_0(V,N)+E_\mathrm{M}(V,N,n_k).`$ (17)
The energy $`E_\mathrm{M}=EE_0`$ is obtained by inserting Eqs. (13), (14) and (11):

$`E_\mathrm{M}(V,N,n_k)`$ $``$ $`{\displaystyle \underset{𝒌}{}}\epsilon _kn_k`$ (19)
$`+{\displaystyle \frac{N^2}{2V}}{\displaystyle d^3rf(r)^2\mathrm{\Delta }g_\mathrm{M}(r,n_k)\stackrel{~}{u}(r)}.`$
The decisive contributions in $`\mathrm{\Delta }g_\mathrm{M}`$ are long-ranged because the low momenta dominate near the transition. For the intended evaluation of the critical behaviour we may therefore use the approximation

$$\mathrm{\Delta }g_\mathrm{M}(r,n_k)\mathrm{\Delta }g_\mathrm{M}(0,n_k).$$
(20)
This yields

$$E_\mathrm{M}(V,N,n_k)=\underset{𝒌}{}\epsilon _kn_k+\frac{N}{2}\frac{w_0}{v}\mathrm{\Delta }g_\mathrm{M}(0,n_k),$$
(21)
where $`v=V/N`$ and

$`{\displaystyle \frac{w_0}{v}}`$ $``$ $`{\displaystyle \frac{1}{v}}{\displaystyle d^3rf(r)^2\stackrel{~}{u}(r)}{\displaystyle \frac{1}{v}}{\displaystyle d^3rg_0(r)\stackrel{~}{u}(r)}`$ (22)
$``$ $`{\displaystyle \frac{2E_{\mathrm{g}.\mathrm{s}.}}{N}}=14.3k_\text{B}\mathrm{K}.`$ (23)
Using McMillan’s parameters, the first integral yields $`14.8k_\text{B}\mathrm{K}`$. Within the uncertainty of Eq. (14) it may be replaced by the second integral. This integral is the theoretical g.s. energy (11). By equating this energy to the experimental ground state energy $`E_{\mathrm{g}.\mathrm{s}.}`$ we fix the value of $`w_0/v`$.

Up to a constant, the energy $`E_\mathrm{M}`$ (21) can be expressed by

$`E_\mathrm{M}`$ $`=`$ $`\mathrm{\Psi }_\mathrm{M}\left|H\right|\mathrm{\Psi }_\mathrm{M}\text{with}`$ (24)
$`H_\mathrm{M}`$ $`=`$ $`{\displaystyle \underset{i}{}}{\displaystyle \frac{\mathrm{}^2}{2m}}\mathrm{\Delta }_i+w_0{\displaystyle \underset{i<j}{}}\delta (𝒓_i𝒓_j).`$ (25)
This means that within sensible approximations the energy $`E_\mathrm{M}`$ in $`E=E_0+E_\mathrm{M}`$ can be calculated as the expectation value of a simple product state $`\mathrm{\Psi }_\mathrm{M}`$ with a simple Hamiltonian $`H_\mathrm{M}`$. It should be noted that $`H_\mathrm{M}`$ is solely introduced for evaluating $`E_\mathrm{M}`$, it cannot serve as a model Hamiltonian in other respects. In particular the stability of the system is guaranteed by the part $`E_0`$ of the total energy $`E=E_0+E_\mathrm{M}`$; the corresponding part of the Hamiltonian is not contained in $`H_\mathrm{M}`$.

The attractive zero-range force in $`H_\mathrm{M}`$ has the following meaning: The l.p.o. implies an enhanced probability of finding two atoms near together. Due to the Jastrow factors in Eq. (5) this enhancement is effective only in the range where the realistic interaction is attractive; the strength of this attraction is measured by $`w_0/v`$. A zero-range force can be used because for $`g_\mathrm{M}`$ the enhancement is approximately $`r`$-independent up to a few Å.

The evaluation of Eq. (25) is straightforward. The two-body matrix elements $`w_{𝒌𝒒}`$ of the $`\delta `$-force with the s.p.f. (7) are given by

$`w_{𝒌𝒒}`$ $`=`$ $`\phi _𝒌\phi _𝒒\left|w_0\delta (𝒓_1𝒓_2)\right|\phi _𝒌\phi _𝒒`$ (26)
$`=`$ $`{\displaystyle \frac{w_0}{V_0}}{\displaystyle \underset{i=1}{\overset{3}{}}}\left[\mathrm{\hspace{0.17em}1}+{\displaystyle \frac{1}{2}}\delta _{k_iq_i}\right]`$ (27)
$`=`$ $`{\displaystyle \frac{w_0}{V_0}}\left[1+{\displaystyle \frac{1}{2}}{\displaystyle \underset{i=1}{\overset{3}{}}}\delta _{k_iq_i}+\mathrm{}\right].`$ (28)
The somewhat artificial construction Eq. (7) implies p.o. in 3 orthogonal directions. In contrast to this we assume that the p.o. is actually realized in just one direction. Therefore we omit products of $`\delta `$-functions in the last expression of Eq. (26) which originate from simultaneous p.o. in 2 or 3 directions. The remaining $`\delta `$-terms describe the mutual correlation of all particle with the same momentum $`q`$ in one direction. The number $`\nu _q`$ of these particles is

$$\nu _q=\underset{k_2,k_3}{}n_k\text{where }k=|𝒌|=\sqrt{q^2+k_2^{\mathrm{\hspace{0.17em}2}}+k_3^{\mathrm{\hspace{0.17em}2}}}.$$
(29)
The obvious form of the correlation energy is $`E_{\mathrm{corr}}\nu _q^{\mathrm{\hspace{0.17em}2}}`$. The actual result follows from Eq. (25) with Eqs. (26) and (29):

$$E_{\mathrm{corr}}(V,N,n_k)=\frac{3w_0}{2v}\frac{1}{N_0^{\mathrm{\hspace{0.17em}1}/3}}\left[\frac{1}{N^{2/3}}\left(\underset{q}{}\nu _q^{\mathrm{\hspace{0.17em}2}}+2n_0\nu _0\right)\right].$$
(30)
For the replacement (9) of the sum by an integral we took into account a possible finite condensate fraction $`n_0/N`$. The index zero of $`n_0`$ and $`\nu _0`$ stands for the lowest possible momentum value.

The complete expression for the energy $`E`$ of Eq. (15) becomes

$$E(V,N,n_k)=E_0+E_\mathrm{M}=E_0+E_{\mathrm{IBG}}+E_{\mathrm{corr}}.$$
(31)
Here $`E_0=E_0(V,N)`$ is given by Eq. (11) and $`E_{\mathrm{IBG}}`$ by

$$E_{\mathrm{IBG}}(V,N,n_k)=\underset{𝒌}{}\epsilon _kn_k.$$
(32)
The critical contribution in Eq. (31) is the correlation energy $`E_{\mathrm{corr}}`$ (30). In the evaluation of Eq. (25) leading to Eq. (31) we dismissed a constant term stemming from the constant in Eq. (26), the higher order terms ($`N_0^{2/3}`$ or $`N_0^1`$) because they come from p.o. in more than one direction, and a term $`n_0^{\mathrm{\hspace{0.17em}2}}`$ because it is noncritical (that means eventually negligible near the transition point).

For the finite boxes the momentum sums have to be evaluated with $`\mathrm{\Delta }k=\pi /V_0^{\mathrm{\hspace{0.17em}1}/3}`$; the subsequent summation over all boxes yields, however, a factor $`V/V_0`$. For plain momentum sums (without $`\delta _{k_iq_i}`$-terms) this procedure can be abbreviated by performing the common summation with $`\mathrm{\Delta }k=\pi /V^{1/3}`$. The final momentum sums ins Eqs. (30) and (32) are performed in this common way. This implies that the occupation numbers $`n_k`$ do no longer refer to a finite box (as in Eq. (4)) but to the macroscopic system (with $`V`$ and $`N`$). In Eq. (30) this implies $`\nu _q^{\mathrm{\hspace{0.17em}2}}N^{5/3}`$ and an $`N`$-independent $`E_{\mathrm{corr}}/N`$.

### C Statistical assumptions

In order to calculate the thermodynamic energy $`E(T,V,N)`$ from Eq. (31) we need the temperature dependent expectation values of the parameters of $`\mathrm{\Psi }`$. Our model uses the IBG expression for $`n_k`$ and treats $`N_0`$ as an adjustable constant. We discuss and specify these assumptions.

Any model of the $`\lambda `$-transition except the IBG introduces the phase transition phenomenologically. In the AIBG this phenomenological assumption is the use of the expectation values $`n_k_{\mathrm{IBG}}`$. This assumption can be made plausible (to some extent) by observing that the additional contribution $`E_{\mathrm{corr}}`$ in $`E_\mathrm{M}=E_{\mathrm{IBG}}+E_{\mathrm{corr}}`$ is small in the sense

$$\frac{E_{\mathrm{corr}}}{E_{\mathrm{IBG}}}=𝒪(y)\text{where}y=\frac{|w_0/v|N_0^{1/3}}{k_\text{B}T_\lambda }1.$$
(33)
The actual parameter values will yield $`y=0.13`$, Eq. (58). For $`N_0\mathrm{}`$ or $`y0`$ one obtains $`E_\mathrm{M}E_{\mathrm{IBG}}`$ and thus $`n_k_{\mathrm{IBG}}`$. (The $`n_k`$-independent contribution $`E_0(V,N)`$ in $`E`$ is without influence on the expectation values). Therefore, we expect $`n_k=n_k_{\mathrm{IBG}}[\mathrm{\hspace{0.17em}1}+𝒪(y)]`$. The evaluation of $`E_{\mathrm{corr}}`$ with $`n_k_{\mathrm{IBG}}`$ is then valid up to first order in $`y`$; it is this contribution which yields the logarithmic singularity. Another argument for using the $`n_k_{\mathrm{IBG}}`$ is that the additional term $`\nu _q^{\mathrm{\hspace{0.17em}2}}`$ correlates many s.p. states with many others, and that it will be without much influence on the occupation of a particular s.p. state.

The AIBG assumes expectation values $`n_k`$ of the IBG-form,

$`n_k`$ $`=`$ $`{\displaystyle \frac{1}{\mathrm{exp}[(\epsilon _k\mu )/k_\text{B}T]1}}`$ (34)
$`=`$ $`{\displaystyle \frac{1}{\mathrm{exp}(x^2+\tau ^2)1}}.`$ (35)
Here $`\mu `$ is the chemical potential and $`k_\text{B}`$ is Boltzmann’s constant. We have introduced the dimensionless quantities $`\tau ^2=\mu /k_\text{B}T`$ and

$$x=\frac{\lambda |𝒌|}{\sqrt{4\pi }}\text{with}\lambda =\frac{2\pi \mathrm{}}{\sqrt{2\pi mk_\text{B}T}}.$$
(36)
The transition temperature of the IBG is given by the following condition for the thermal wave length $`\lambda =\lambda (T)`$:

$$\lambda (T_\lambda )=\left[v\zeta (3/2)\right]^{1/3},$$
(37)
where $`\zeta (3/2)=2.6124`$ denotes Riemann’s zeta function. In applying the AIBG to the real system we will identify $`T_\lambda `$ with the actual transition temperature; formally this can be achieved by inserting a suitable effective mass $`m=m(v)`$ in $`\lambda `$ of Eq. (36), Ref. . In the following we use the relative temperature

$$t=\frac{TT_\lambda }{T_\lambda }.$$
(38)
For $`t0`$ the chemical potential is determined by the particle number condition $`\mathrm{\Sigma }n_k=N`$. For $`t0`$ this condition yields $`\mu \epsilon _00`$ and thus $`n_0\mathrm{}`$. For simplicity we rename $`\mu \epsilon _0`$ by $`\mu `$; it is this new $`\mu `$ which vanishes at $`t=0`$ also for $`\epsilon _00`$.

We are interested in the critical properties. Therefore we expand $`\mu `$ or, equivalently, $`\tau `$ for $`|t|1`$:

$$\tau (t)=\sqrt{\frac{\mu }{k_\text{B}T}}=\{\begin{array}{ccc}at+bt^2+\mathrm{}& & (t>0),\\ a^{}|t|+b^{}t^2+\mathrm{}& & (t<0).\end{array}$$
(39)
This expansion is directly related to the model condensate fraction

$`{\displaystyle \frac{n_0}{N}}`$ $`=`$ $`1{\displaystyle \frac{1}{N}}{\displaystyle n_k}`$ (40)
$`=`$ $`\{\begin{array}{ccc}0& & (t>0),\\ f|t|+gt^2+\mathrm{}& & (t<0).\end{array}`$ (43)
The coefficients $`f,g,\mathrm{}`$ are determined by $`a^{},b{}_{}{}^{},\mathrm{}`$; in particular $`f=3/2+2\sqrt{\pi }a^{}/\zeta (3/2)`$. The IBG yields finite values for $`a,b,\mathrm{}`$, and $`a^{}=b^{}=\mathrm{}=0`$.

Deviating from the IBG, the AIBG admits a coefficient $`a^{}0`$ in Eq. (39) or $`f>3/2`$ in Eq. (40). This deviation is introduced phenomenologically for the following reason: For $`t>0`$ the evaluation of $`\nu _q^{\mathrm{\hspace{0.17em}2}}`$ with $`n_k_{\mathrm{IBG}}`$ yields straightforwardly a logarithmic singularity. This singularity can be directly continued to $`t<0`$ by admitting $`a^{}0`$; the singularity itself does not depend on the value for $`a^{}`$. The form of the expansion (39) is to some extent plausible because of its simplicity and symmetry. Physically $`a^{}0`$ means an energy gap of $`a^{}{}_{}{}^{2}t_{}^{2}k_\text{B}T`$ for the condensed s.p. state for $`t<0`$.

For evaluating the thermodynamic energy we need also the expectation values of products of $`n_k`$’s. (We do not need, however, the full statistical information contained in the density matrix). For such products the IBG yields

$$n_𝒌n_𝒒=n_𝒌n_𝒒(𝒌𝒒).$$
(44)
Compared to the IBG, the AIBG states depend on one additional parameter, namely $`N_0`$. In principle, this parameter (as well as the occupation numbers) should be determined from the condition of minimal free energy $`F`$ yielding a temperature dependent expectation value $`N_0`$. Instead of this, we treat $`N_0`$ as an adjustable constant. Simplified estimates (Sec. IV A and Ref. ) yield a sensible value for $`N_0`$ at $`t=0`$.

## III Logarithmic singularity

We evaluate the thermodynamic energy $`E(T,V,N)=E(V,N,n_k)`$. According to Eq. (31) it is of the form

$$E(T,V,N)=E_0(V,N)+E_{\mathrm{IBG}}(T,V,N)+E_{\mathrm{corr}}(T,V,N).$$
(45)
The term $`E_0(V,N)`$ is essential for a realistic compressibility but it does not contribute to the specific heat. The energy $`E_{\mathrm{IBG}}=\epsilon _kn_k`$ is calculated as in the IBG. Equations (30) and (44) yield for the decisive contribution $`E_{\mathrm{corr}}`$,

$$\frac{E_{\mathrm{corr}}(T,V,N)}{N}=\frac{3w_0}{2v}\frac{1}{N_0^{\mathrm{\hspace{0.17em}1}/3}}\frac{1}{N^{5/3}}[\nu _q^2+2n_0\nu _0].$$
(46)
Leaving away unessential constants we sketch how $`E_{\mathrm{corr}}`$ yields a logarithmic singularity: For $`|t|0`$ and $`x0`$ the occupation numbers (34) behave like $`n_x1/[x^2+t^2]`$ and $`\nu _q𝑑kk(k^2+q^2+t^2)^1\mathrm{ln}(q^2+t^2)`$. It follows that $`n_0\nu _0|t|\mathrm{ln}|t|`$ and $`\mathrm{\Sigma }\nu _q^2𝑑q[\mathrm{ln}(q^2+t^2)]^2\text{const.}+|t|\mathrm{ln}|t|`$.

We present now the results of the detailed calculation. The $`\nu _q`$ of Eq. (29) with Eq. (34) can be evaluated analytically,

$`{\displaystyle \frac{\nu _q}{N^{2/3}}}`$ $`=`$ $`{\displaystyle \frac{1}{N^{2/3}}}{\displaystyle \underset{k_2,k_3}{}}{\displaystyle n_k}`$ (47)
$`=`$ $`{\displaystyle \frac{v^{2/3}}{\lambda ^2}}\mathrm{ln}\left[1\mathrm{exp}\left(x^2\tau ^2\right)\right].`$ (48)
We used the dimensionless quantities $`\tau `$ of Eq. (39) and $`x=q\lambda /\sqrt{4\pi }`$. The integral in Eq. (46) becomes

$$\frac{1}{N^{5/3}}\nu _q^2=\frac{2}{\sqrt{\pi }}\zeta (3/2)^{5/3}\left(\frac{T}{T_\lambda }\right)^{5/2}L(\tau ),$$
(49)
where $`L(0)=8.30`$ and

$`L(\tau )`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}𝑑x\left(\mathrm{ln}\left[1\mathrm{exp}\left(x^2\tau ^2\right)\right]\right)^2`$ (50)
$`=`$ $`L(0)+4\pi \tau \mathrm{ln}\tau +4\pi \tau \left[\mathrm{ln}(2)1\right]+𝒪(\tau ^2).`$ (51)
The expectation value $`\nu _0`$ is given by Eq. (47) with $`x=0`$ because $`\epsilon _0`$ has been absorbed in $`\tau ^2=(\mu \epsilon _0)/k_\text{B}T`$ (remark after Eq. (38)). Using $`\nu _0\mathrm{ln}|t|`$, $`n_0|t|`$ and $`\tau |t|`$ we see that both terms in Eq. (46) yield a contribution $`t\mathrm{ln}|t|`$. The critical ($`|t|1`$) behaviour of the specific heat reads therefore

$`c_V(T,V)`$ $`=`$ $`\left({\displaystyle \frac{(E/N)}{T}}\right)_V`$ (52)
$`=`$ $`\{\begin{array}{ccc}A\mathrm{ln}|t|+B+\mathrm{}\hfill & & (t>0),\\ A^{}\mathrm{ln}|t|+B^{}+\mathrm{}\hfill & & (t<0).\end{array}`$ (55)
The detailed evaluation of Eq. (46) yields

$$A=A^{}=9\zeta (3/2)^{2/3}\frac{|w_0|}{v}\frac{1}{N_0^{\mathrm{\hspace{0.17em}1}/3}}\frac{1}{T_\lambda }.$$
(56)
The terms depending on the parameter $`a^{}`$ cancel. We insert $`|w_0|/v`$ of Eq. (21) and $`T_\lambda =2.17\mathrm{K}`$ in the result (56).

Experimentally a logarithmic singularity has been measured for $`c_P`$ over several decades by Ahlers. For saturated vapour pressure we may use $`c_Vc_P`$. (The consistency of a logarithmic form for both, $`c_V`$ and $`c_P`$, has been demonstrated by Lee and Puff ). The comparison of our result (56) with the experimental values $`AA^{}0.63k_\text{B}`$ (Eq. (48) of Ref. ) determines the model parameter $`N_0`$,

$$N_0^{\mathrm{\hspace{0.17em}1}/3}50\text{(adjusted)}.$$
(57)
This fixes also the small parameter $`y`$ of Eq. (33) which is a quantitative measure for the deviations of the AIBG from the IBG,

$$y=\frac{|w_0/v|N_0^{1/3}}{k_\text{B}T_\lambda }0.13.$$
(58)
Experimentally one finds a small difference between $`A`$ and $`A^{}`$ (Ahlers reports $`5\%`$) and an indication of a cut of the logarithmic singularity (slightly negative critical exponents for $`c_P`$). In the AIBG a difference between $`A`$ and $`A^{}`$ could be obtained if corrections to approximation (20) are taken into account which differ for $`t>0`$ and $`t<0`$. A possible cut of the singularity is connected with the lowest possible values of the momenta.

The coefficients $`B`$ and $`B^{}`$ in Eq. (52) can be calculated from Eq. (46), too. They depend on the value of the parameter $`a^{}`$ which will be determined only later (Sec. V C). Moreover a generalization of the s.p.f. to be introduced in Sec. IV C will yield another contribution to the jump in $`B^{}B`$.

## IV Phase ordering

We review various aspects of the p.o. assumed in the AIBG. First we show that this p.o. is favoured by the free energy (Sec. IV A). Then we review the necessity of the localization of the s.p.f. for obtaining finite correlation effects due to p.o. (Sec. IV B); this leads to the notion that the coherence range of p.o. approaches infinity at the transition point. Finally (Sec. IV C) we connect the p.o. with a potential superfluid flow.

### A Free energy due to phase ordering

The l.p.o. leads to a lower energy, $`E_{\mathrm{corr}}/N(w_0/v)`$ $`N_0^{1/3}<0`$. At the same time such an ordering reduces the entropy, $`\mathrm{\Delta }S(\mathrm{p}.\mathrm{o}.)<0`$. We show that $`\mathrm{\Delta }S(\mathrm{p}.\mathrm{o}.)`$ is relatively small and that p.o. is therefore favoured by the free energy.

Arbitrary phases in Eq. (2) are statistically equivalent to an equal weight of the two independent choices for the phase, $`\varphi _0`$ as in Eq. (3) or $`\varphi _0+\pi /2`$. We consider p.o. in just one direction. For each $`𝒌`$ we have then two s.p. states. The number of possible distributions of $`n_k`$ bosons on 2 states is $`n_k+1`$. The p.o. requires that all $`n_k`$ bosons go into the same s.p. state; compared to a statistical distribution this implies the entropy change $`\mathrm{\Delta }s(\mathrm{p}.\mathrm{o}.)=k_\text{B}\mathrm{ln}(n_k+1)`$. On the other hand, the $`n_k`$ atoms gain by p.o. the energy $`\mathrm{\Delta }e(\mathrm{p}.\mathrm{o}.)(E_{\mathrm{corr}}/N)n_k`$. The free energy change $`\mathrm{\Delta }f`$ is then

$$\mathrm{\Delta }f(\mathrm{p}.\mathrm{o}.)\frac{|w_0|}{vN_0^{1/3}}n_k+k_\text{B}T\mathrm{ln}(n_k+1).$$
(59)
In spite of the smallness (58) of the considered correlations the p.o. is favoured ($`\mathrm{\Delta }f(\mathrm{p}.\mathrm{o}.)<0`$) as soon as $`n_k1`$. The boson property of the atoms is decisive for this conclusion: It implies that $`|\mathrm{\Delta }s|`$ is not proportional to $`n_k`$ but only to $`\mathrm{ln}n_k`$; this behaviour may be paraphrased by ‘bosons like to go into the same state’. For the relevant lowest s.p. states ($`k0`$) the condition $`n_k1`$ is fulfilled for $`t1`$. Thus p.o. will indeed be adopted when the transition point is approached ($`t0^+`$).

As a simplification our ansatz (4) assumes p.o. for all s.p.f. Summing over all momenta and including numerical factors one finds then $`E_{\mathrm{corr}}T\mathrm{\Delta }S(\mathrm{p}.\mathrm{o}.)`$ for $`TT_\lambda `$ and for the $`N_0`$ of Eq. (57). Reversely, the condition $`E_{\mathrm{corr}}T_\lambda \mathrm{\Delta }S(\mathrm{p}.\mathrm{o}.)`$ constitutes a rough theoretical estimate for the model parameter $`N_0`$, leading to an $`N_0`$ of the order (57).

### B Phase coherence volume

In the Introduction we argued that for finite correlation effects due to p.o. the s.p.f. have to be localized. We review this point and find that the localization is not required for the lowest s.p. state. We present then an estimate for the extension of this lowest s.p.f.; within this extension the p.o. will be coherent.

The finiteness of the correlation energy ($`N_0^{1/3}`$) is due to the finite spacing $`\mathrm{\Delta }k=\pi /V_0^{\mathrm{\hspace{0.17em}1}/3}`$ of the $`q`$-values in the sum $`E_{\mathrm{corr}}\nu _q^{\mathrm{\hspace{0.17em}2}}`$. Nonlocalized s.p.f. would lead to the same $`E_{\mathrm{corr}}`$ provided that only s.p. states with

$$q_n=q_0+n\mathrm{\Delta }k,(n=0,1,2,\mathrm{}).$$
(60)
are occupied. In the macroscopic system the possible $`q`$-values are, however, dense, and the entropy drives the particles into the occupation of all available states. This is the reason why a finite $`\mathrm{\Delta }k`$ can be realized only for s.p.f. localized within a volume $`V_0=(\pi /\mathrm{\Delta }k)^3`$.

This argument shows that we may either start from finite volumes or from finite $`\mathrm{\Delta }k`$, the second feature follows from the first one. There is, however, one difference: Starting from the finite spacing $`\mathrm{\Delta }k`$ in Eq. (60), the natural condition for the lowest $`q`$-value is $`q_0<\mathrm{\Delta }k`$. The corresponding s.p.f. can then not be localized within $`V_0`$; as just stated there is also no such need for a localization for the required correlation effect.

In view of this we introduce the following modification of the states (4): Only the s.p.f. with $`q_{n>0}`$ are localized within $`V_0`$, the s.p.f. with $`q_0`$ have a larger extension. We present a crude estimate for the volume $`V_\mathrm{c}`$ of the lowest s.p.f. $`\phi _0`$ which becomes the condensate state for $`t<0`$: Let $`V_\mathrm{c}`$ be some multiple of $`V_0`$, that means $`V_\mathrm{c}=WV_0`$. In a volume $`V_\mathrm{c}`$ there are $`W`$ s.p. states below $`\mathrm{\Delta }k`$ out of which only one is occupied. A redistribution of the $`n_0`$ atoms over these $`W`$ states would increase the entropy by $`\mathrm{\Delta }sk_\text{B}\mathrm{ln}(n_0)^W`$. At the same time these atoms would loose their correlation energy, $`\mathrm{\Delta }en_0|w_0/v|N_0^{1/3}`$. The stability condition $`T_\lambda \mathrm{\Delta }s<\mathrm{\Delta }e`$ yields an upper bound for $`W`$. Using $`n_0(t>0)t^2`$ and $`n_0(t>0)=𝒪(N)`$ we obtain

$$V_\mathrm{c}(t)=WV_0V_0\{\begin{array}{ccc}t^2& & (t>0),\\ \mathrm{}& & (t<0).\end{array}$$
(61)
This result means the following: Without changing our previous calculations we can modify the states (4) in such a way that the lowest s.p.f. $`\phi _0`$ has the large and eventually infinite volume $`V_\mathrm{c}`$, whereas the s.p.f. with $`q_n>0`$ are localized. This modification of Eq. (4) has important consequences:

1. Within the volume $`V_\mathrm{c}`$ the direction of p.o. is defined by $`\phi _0`$. Therefore, $`V_\mathrm{c}(t)`$ of Eq. (61) is the phase coherence volume of the considered p.o.
2. The p.o. constitutes a symmetry breaking because in an infinite system the atoms are free to adopt arbitrary phases $`\varphi _j`$ in Eq. (2). (The average over the size and shape of finite boxes (assumed in Sec. II A) does not restore this symmetry). The result (61) means that this symmetry breaking changes its character from local to global at the transition point. Approaching the transition point ($`t0^+`$) the directions of p.o. of neighbouring $`V_0`$’s get aligned, and for $`t<0`$ the coherence is potentially infinite.
3. Strictly finite and separated boxes imply a momentum cut at $`\mathrm{\Delta }k`$ and consequently a cut of the logarithmic singularity. This cut would be described by a lower bound $`\mathrm{\Delta }k`$ in the integral in Eq. (9). For Eq. (60) we revise the replacement Eq. (9) of the momentum sums by integrals:

$$\underset{q_n=q_0+n\mathrm{\Delta }k}{}\mathrm{}_{q_0}\frac{1}{\mathrm{\Delta }k}_{q_{0,\mathrm{min}}}^{\mathrm{}}𝑑q\mathrm{}.$$
(62)
The brackets indicate an average over possible $`q_0`$’s with $`q_0\mathrm{\Delta }k`$. The lowest value for $`q_0`$ is determined by Eq. (61). A lower bound $`q_{0,\mathrm{min}}t^2`$ (for $`t>0`$) of the integral does not lead to a cut of the logarithmic singularity.

### C Complex phase ordering

The phenomenon of superfluidity can be explained by assuming a macroscopic wave function $`\psi `$. The phase $`\mathrm{\Phi }(𝒓)`$ of the complex $`\psi `$ yields the velocity $`v_\mathrm{s}=(\mathrm{}/m)\mathbf{}\mathrm{\Phi }`$ of a potential superfluid flow. In the AIBG the condensed particles for a macroscopic wave function for $`t<0`$, too, Eq. (61). In order to join the common description we replace the real s.p.f. $`\phi _0`$ by the complex s.p.f.:

$$\overline{\phi }_0=\phi _0\mathrm{exp}[\mathrm{i}\mathrm{\Phi }(𝒓)].$$
(63)
This replacement does not affect the previous calculations; moreover, the Jastrow factors in Eq. (5) are without influence on a flow due to $`\mathbf{}\mathrm{\Phi }0`$. With Eq. (63) all results obtained from the assumption of a macroscopic wave function apply to our model, too. Besides this basic conformity the AIBG leads, however, also to a peculiar modification which will be discussed in the following.

In an IBG-like model (with $`n_0|t|`$ or $`\beta =1/2`$) we have to assume that noncondensed particles contribute to the superfluid density $`\rho _\mathrm{s}`$ in order to get agreement with the experiment ($`\rho _\mathrm{s}|t|^{2/3}`$ or $`\nu =1/3`$). The considered p.o. should in some way imply that noncondensed particles move coherently with the condensate. We introduce the possibility of a net current of the noncondensed particles by replacing the original real s.p.f. $`\phi _𝒌`$ by

$$\overline{\phi }_𝒌=\phi _𝒌\mathrm{exp}\left[\mathrm{i}\mathrm{\Phi }_𝒌(𝒓)\right].$$
(64)
A coherent motion with the condensate can now be described by the condition $`\mathrm{\Phi }_𝒌=\mathrm{\Phi }`$. We call this condition complex phase ordering (c.p.o.) named after the complex phase factors in Eqs. (63) and (64). In contrast to this the p.o. (3) will be called real phase ordering (r.p.o.). We note two points: The spatial correlations of the s.p.f. are unchanged by the replacement (63), (64) because $`|\phi |^2=|\overline{\phi }|^2`$. Secondly, for a coherent macroscopic flow it is sufficient that the field $`\mathrm{\Phi }(𝒓)`$ is macroscopic; it is then not necessary that all contributing s.p.f. are of macroscopic range.

The replacement (63) is suggested by the familiar picture of a superfluid, and Eq. (64) with $`\mathrm{\Phi }_𝒌=\mathrm{\Phi }`$ is necessary to reconcile $`\beta =1/2`$ (model) with $`\nu =1/3`$ (experiment). The fairly straightforward introduction of the phase factors (63) and (63) determines the crucial model prediction for $`S_\mathrm{s}0`$. The fit to the experimental $`\rho _\mathrm{s}`$ determines the extent of the c.p.o. and consequently $`S_\mathrm{s}`$.

Presenting a peculiar prediction like $`S_\mathrm{s}0`$ we feel obliged to make it quantitative. For this purpose we have to write down explicit expressions for $`\rho _\mathrm{s}`$ and $`S_\mathrm{s}`$. This requires additional assumptions which are introduced by plausibility arguments. As explained in the last section, the quantitative result for $`S_\mathrm{s}`$ will not depend sensitively on these assumptions because the model expression for $`\rho _\mathrm{s}`$ is fitted to the experiment.

We require that the new modes $`\mathrm{\Phi }_k`$ do not destroy the correlation energy due to the r.p.o.; this determines the extent of the c.p.o. The field $`\mathrm{\Phi }(𝒓)`$ contains degrees of freedom which will be thermally excited. Let $`k_\mathrm{c}(t)`$ be the average amount of the momenta of this field

$$k_\mathrm{c}(t)=\left|\mathbf{}\mathrm{\Phi }(𝒓)\right|.$$
(65)
Because of $`|\phi |^2=|\overline{\phi }|^2`$ the fields $`\mathrm{\Phi }`$ and $`\mathrm{\Phi }_𝒌`$ do not directly influence the spatial correlations. However, for $`\mathrm{\Phi }_𝒌0`$ the noncondensed states $`𝒌`$ with $`\epsilon _k<\mathrm{}^2k_\mathrm{c}^{\mathrm{\hspace{0.17em}2}}/2m`$ would have a lower energy than the condensate state. This would invalidate the assumption about the occupation pattern and destroy the correlation energy. On the other hand, for $`kk_\mathrm{c}`$ the $`\mathrm{\Phi }_𝒌`$ will have only minor influence on the occupation pattern. Simplifying and quantifying this qualitative argument leads to the following condition for c.p.o.:

$$\mathrm{\Phi }_k(𝒓)=\{\begin{array}{ccc}\mathrm{\Phi }(𝒓)& & (k<k_\mathrm{c}),\\ \text{arbitrary}& & (k>k_\mathrm{c}).\end{array}(\mathrm{c}.\mathrm{p}.\mathrm{o}.)$$
(66)
This condition means that the localized s.p.f. with $`k<k_\mathrm{c}`$ adopt within their range the (macroscopically defined) phase $`\mathrm{\Phi }`$ of the condensate. The physical reason for this c.p.o. is that the correlations due to the r.p.o. (favoured by the free energy) must not be destroyed; in this way the real and complex p.o. are connected to each other.

## V Superfluid density

### A AIBG expression

Using the new s.p.f. (63) and (64) the states (5) become for $`t<0`$:

$$\mathrm{\Psi }=F𝒮\left(\phi _0\mathrm{exp}[\mathrm{i}\mathrm{\Phi }(𝒓)]\right)^{n_0}\underset{𝒌,\mathrm{vol}}{}\left(\phi _𝒌\mathrm{exp}[\mathrm{i}\mathrm{\Phi }_𝒌(𝒓)]\right)^{n_k}.$$
(67)
For these states we evaluate the quantum mechanical expectation value of the current operator. Taking into account the c.p.o. (66) we obtain

$`𝒋(𝒓,n_k)`$ $`=`$ $`\mathrm{\Psi }|{\displaystyle \frac{\mathrm{}}{2\mathrm{i}}}{\displaystyle \underset{j=1}{\overset{N}{}}}\mathbf{}_j\delta (𝒓𝒓_j)+\mathrm{c}.\mathrm{c}.|\mathrm{\Psi }`$ (68)
$`=`$ $`{\displaystyle \frac{\rho }{N}}{\displaystyle \frac{\mathrm{}}{m}}\left[\left(n_0+{\displaystyle \underset{k<k_\mathrm{c}}{}}{\displaystyle n_k}\right)\mathbf{}\mathrm{\Phi }(𝒓)+{\displaystyle \underset{k>k_\mathrm{c}}{}}{\displaystyle n_k\mathbf{}\mathrm{\Phi }_𝒌(𝒓)}\right].`$ (69)
All derivatives of the real functions ($`F`$, $`\phi _0`$ and $`\phi _𝒌`$) are cancelled by the complex conjugate (c.c.) term; the only surviving contributions are the phase derivatives. Each such derivative is accompanied by a factor $`\mathrm{\Psi }|\delta (𝒓𝒓_j)|\mathrm{\Psi }=\rho (𝒓)/(Nm)`$. The mass density $`\rho `$ is assumed to be a constant in the following.

The contributions of the nonordered phases are averaged out statistically. This is also the case for the thermal excitations of the $`\mathrm{\Phi }`$-field, $`\mathbf{}\mathrm{\Phi }_{\mathrm{th}}=0`$. A coherent nonvanishing flow can only be obtained by an additional nonequilibrium contribution $`\mathrm{\Phi }_\mathrm{s}`$ to $`\mathrm{\Phi }`$,

$$\mathrm{\Phi }=\mathrm{\Phi }_{\mathrm{th}}+\mathrm{\Phi }_\mathrm{s}.$$
(70)
The flux due to $`\mathrm{\Phi }_\mathrm{s}`$ might be relatively stable if it is small enough. This implies $`|\mathbf{}\mathrm{\Phi }_\mathrm{s}|k_\mathrm{c}`$ where $`k_\mathrm{c}=|\mathbf{}\mathrm{\Phi }_{\mathrm{th}}|`$.

Using Eq. (70) and $`\mathbf{}\mathrm{\Phi }_{\mathrm{th}}=0`$ we evaluate the statistical expectation value of Eq. (68). Writing $`𝒋_\mathrm{s}=\rho _\mathrm{s}𝒗_\mathrm{s}`$ with $`𝒗_\mathrm{s}=(\mathrm{}/m)\mathbf{}\mathrm{\Phi }_\mathrm{s}`$ we obtain for $`\rho _\mathrm{s}`$:

$$\frac{\rho _\mathrm{s}}{\rho }=\frac{1}{N}\left(n_0+\underset{k<k_\mathrm{c}}{}n_k\right)=\frac{\rho _0+\rho _\mathrm{c}}{\rho }.$$
(71)
This is the AIBG expression for the superfluid density. In a number of points we summarize its implications:

1. The density $`\rho _\mathrm{s}`$ is composed by all atoms adopting the phase $`\mathrm{\Phi }(𝒓)`$. It consists of the condensate density $`\rho _0`$ and the comoving density $`\rho _\mathrm{c}`$. The temperature dependence of this composition (Fig. 1) is obtained by fitting Eq. (71) to the data. The comoving density $`\rho _\mathrm{c}`$ has an internal structure which implies a nonvanishing entropy $`S_\mathrm{s}0`$ of the superfluid component.

1. Expression (71) is the simplest possibility to account in an IBG-like model with $`n_0|t|`$ for $`\rho _\mathrm{s}|t|^{2/3}`$. Fitting Eq. (71) to the experimental $`\rho _\mathrm{s}`$ fixes $`k_\mathrm{c}(t)`$, and consequently the difference between the densities $`\rho _0`$ and $`\rho _\mathrm{s}`$ and the model prediction for $`S_\mathrm{s}0`$.
2. The reproduction of the data with Eq. (71) requires asymptotically $`k_\mathrm{c}(t)|t|^{2/3}`$. Theoretically this asymptotic form is made plausible in Sec. IV B. It leads to a certain functional form of $`\rho _\mathrm{s}`$ which yields excellent fits to the experimental data (Sec. V C).
3. The result (71) contains a specific picture for the relation between the condensate fraction $`n_\mathrm{c}`$ and the superfluid fraction $`\rho _\mathrm{s}/\rho `$. The explanation of superfluidity is based on the assumption of a macroscopic wave function. This connection implies, however, an open question: What is the quantitative relation between $`n_\mathrm{c}(T)`$ and $`\rho _\mathrm{s}(T)/\rho `$, and in particular the relation between the $`T=0`$ values $`n_\mathrm{c}0.1`$ (see Ref. ) and $`\rho _\mathrm{s}/\rho =1`$ ? The AIBG proposes the following picture: The model condensate $`n_{\mathrm{c},\mathrm{M}}=n_0/N=\rho _0/\rho `$ is the fundamental quantity. It is depleted by the Jastrow factor $`F`$ in Eq. (67) to the real condensate $`n_\mathrm{c}`$, for $`T=0`$ the value $`n_{\mathrm{c},\mathrm{M}}=1`$ is reduced to $`n_\mathrm{c}0.1`$. The factor $`F`$ in Eq. (67) does, however, not deplete the current due to phase factors. Therefore $`n_{\mathrm{c},\mathrm{M}}1`$ implies $`\rho _\mathrm{s}/\rho 1`$. On the other hand, for $`|t|1`$ the contribution $`\rho _\mathrm{c}\rho _0`$ dominates $`\rho _\mathrm{s}`$.

### B Effective Ginzburg-Landau model

Ginzburg and Sobyanin have proposed an effective Ginzburg-Landau functional for the free energy $`F_{\mathrm{GL}}`$ as a function of the order parameter field $`\psi `$. In this ansatz singular coefficients (like $`|t|^{4/3}`$ for the $`|\psi |^2`$ term) are introduced in order to yield the right critical exponents. We follow this kind of approach for investigating the relation between the critical exponents of the order parameter and the superfluid density. This detour is, however, not required for the final result for $`S_\mathrm{s}`$ (see point 2 of Sec. V A).

For discussing the fluctuations of the $`\mathrm{\Phi }`$-field we define the order parameter field by

$$\psi =\sqrt{\frac{n_0}{V}}\mathrm{exp}[\mathrm{i}\mathrm{\Phi }(𝒓)].$$
(72)
For a definite phase the state (67) must be replaced by the appropriate coherent state. This means that $`n_0`$ in Eq. (72) has to be understood as a quantum mechanical expectation value of the occupation number $`n_0`$ in such a coherent state.

The statistical expectation value $`n_0|t|`$ can be obtained by minimizing the common Landau energy $`F_L/V=Rt|\psi |^2+U|\psi |^4`$ (with regular coefficients $`R`$ and $`U`$). Adding a naive kinetic energy term $`(\mathrm{}^2/2m)|\mathbf{}\psi |^2`$ leads to a violation of scaling invariance and to wrong critical exponents. The decisive feature of the AIBG expression (71) is that the mass density comoving with $`\mathbf{}\mathrm{\Phi }`$ is $`\rho _0+\rho _\mathrm{c}`$ instead of $`\rho _0=m|\psi |^2`$ alone. This suggests the following effective Ginzburg-Landau ansatz

$$\frac{F_{\mathrm{GL}}}{V}=\frac{\mathrm{}^2}{2m}\frac{\rho _\mathrm{s}}{\rho _0}\left|\mathbf{}\psi \right|^2+Rt\left|\psi \right|^2+U\left|\psi \right|^4.$$
(73)
Leaving away unessential constants, the asymptotic behaviour of the mass density $`\rho _\mathrm{c}`$ in Eq. (71) is

$$\rho _\mathrm{c}_{k<k_\mathrm{c}}\frac{d^3k}{k^2+t^2}k_\mathrm{c}.$$
(74)
This, $`\rho \rho _\mathrm{c}`$ for $`|t|0`$ and Eq. (65) determine the asymptotic kinetic energy in Eq. (73),

$$\frac{\rho _\mathrm{s}}{\rho _0}|\mathbf{}\psi |^2=\frac{\rho _\mathrm{c}}{m}|\mathbf{}\mathrm{\Phi }|^2k_\mathrm{c}^{\mathrm{\hspace{0.17em}3}}.$$
(75)
Scaling this kinetic part of $`F_{\mathrm{GL}}`$ with $`F_Lt^2`$ yields

$$k_\mathrm{c}|t|^{2/3}.$$
(76)
This implies $`\rho _\mathrm{s}\rho _\mathrm{c}|t|^{2/3}`$ and a singular mass coefficient $`\rho _\mathrm{s}/\rho _0|t|^{1/3}`$ in Eq. (73). Therefore, Eq. (73) is an effective Ginzburg-Landau ansatz in the same sense as the one by proposed by Ginzburg and Sobyanin. Equation (73) demonstrates how the critical exponent $`\beta =1/2`$ (from $`n_0|\psi |^2|t|^{2\beta }`$) can be connected to $`\nu =1/3`$ (from $`\rho _\mathrm{s}|t|^{2\nu }`$) on account of the comoving density $`\rho _\mathrm{c}`$. The divergent mass coefficient damps the critical fluctuations such that Eq. (73) becomes scaling invariant (this is discussed in more detail in Ref. ). The scaling invariance implies that Eq. (73) might be used down to $`|t|=0`$ and that the critical exponent of $`k_\mathrm{c}`$ and $`\rho _\mathrm{s}`$ might be indeed exactly $`2/3`$. This possibility is supported by the excellent fits obtained from Eq. (71) with Eq. (76).

### C Fit to the experimental data

We insert the leading $`\tau =a^{}|t|`$ and $`k_\mathrm{c}=k_1|t|^{2/3}`$ in Eq. (71). This yields

$$\frac{\rho _\mathrm{s}}{\rho }=a_1|t|^{2/3}+a_2|t|+a_3|t|^{4/3}+\mathrm{}\text{(MAS)}.$$
(77)
With 3 parameters we take this as our asymptotic model fit, called MAS. Figure 2 shows that MAS yields an excellent reproduction of the data by Greywall and Ahlers (GA), for saturated vapour pressure.

The 3-parameter standard fit (SF) used by GA is

$$\frac{\rho _\mathrm{s}}{\rho }=a_1|t|^{a_2}\left(1+a_3|t|^{1/2}+\mathrm{}\right)\text{(SF)}.$$
(78)
We compare both fits, MAS and SF, by calculating the sums $`\chi ^2`$ of the quadratic deviations:

$$\frac{\chi _{\mathrm{SF}}^2}{\chi _{\mathrm{MAS}}^2}10(|t|0.03).$$
(79)
This shows that the reproduction shown in Fig. 2 is not a trivial result; in the considered range ($`|t|0.03`$) the data cannot be reproduced by SF. We remark that SF (but probably also other 3-parameter fits) fit the data in the considerably smaller range $`|t|0.004`$; this range is used by GA for the fit. The range $`|t|0.03`$ seems to be appropriate for a 3-parameter fit because already a 1-parameter fit ($`a_1|t|^{2/3}`$) roughly reproduces the data in this range (see Fig. 1).

The applicability range of Eq. (71) can be considerably extended by using the following 4-parameter expansion for $`\tau `$ and $`k_\mathrm{c}`$:

$$\tau (t)=a^{}|t|,x_\mathrm{c}(t)=x_1|t|^{2/3}+x_2|t|+x_3|t|^{4/3}.$$
(80)
The form of $`x_\mathrm{c}=\lambda k_\mathrm{c}(t)/\sqrt{4\pi }`$ suggests itself because it leaves the expansion (77) unchanged. The 4-parameter model fit, Eq. (71) with Eq. (80), reproduces the data down to $`1\mathrm{K}`$ within the experimental errors. The fit yields the parameter values $`a^{}=3.019`$, $`x_1=2.7028`$, $`x_2=0.837`$ and $`x_3=3.842`$; these values have been used for Fig. 1. More details of this fit and a discussion with respect to the quasi-particle picture are given in Ref. .

The fit to the experimental data fixes the decomposition of $`\rho _\mathrm{s}`$ into $`\rho _0`$ and $`\rho _\mathrm{c}`$. This decomposition is the starting point for the evaluation of the superfluid entropy $`S_\mathrm{s}`$.

## VI Superfluid entropy

### A Model prediction

In $`\rho _\mathrm{s}=\rho _0+\rho _\mathrm{c}`$ only the condensate part $`\rho _0`$ corresponds to a macroscopic wave function and has thus zero entropy content. The comoving part $`\rho _\mathrm{c}`$ is made up by different s.p.f. (however with the same phase factor in Eq. (64)) and has therefore a nonvanishing entropy content. We determine this entropy.

We start with the well-known entropy expression $`S(n_k)`$ for a Bose gas with occupation numbers $`n_k`$. The IBG equilibrium entropy $`S_{\mathrm{IBG}}`$ is obtained from this expression by $`S_{\mathrm{IBG}}=S(n_k)=S(n_k_{\mathrm{IBG}})`$. In Sec. IV A we have discussed that for r.p.o. in one direction only every second s.p. state is occupied. Since the total particle number $`N=n_k`$ is fixed we have to put twice as many atoms in every second state (as compared to the plain IBG expression). This leads to the following model entropy $`S_\mathrm{M}`$,

$`S_\mathrm{M}(T,V,N)`$ $`=`$ $`{\displaystyle \frac{k_\text{B}}{2}}{\displaystyle \underset{𝒌}{}}[(1+2n_k)\mathrm{ln}(1+2n_k)`$ (82)
$`2n_k\mathrm{ln}(2n_k)].`$
For $`t=0`$ this yields $`S_\mathrm{M}/N=0.96k_\text{B}`$ instead of $`S_{\mathrm{IBG}}/N=1.28k_\text{B}`$. The difference $`\mathrm{\Delta }S(\mathrm{p}.\mathrm{o}.)=S_{\mathrm{IBG}}S_\mathrm{M}`$ has been compared to $`\mathrm{\Delta }E(\mathrm{p}.\mathrm{o})`$ in Sec. IV A. For $`t<0`$ the $`n_k`$ with $`\tau (t)`$ of Eq. (80) will be used in Eq. (82). The resulting overall behaviour of $`S_\mathrm{M}(T)`$ is then —in contrast to $`S_{\mathrm{IBG}}(T)`$— similar to that of the experimental entropy $`S(T)`$.

The entropy of the superfluid component is $`S_\mathrm{c}=S(\rho _\mathrm{c})`$. For calculating $`S_\mathrm{c}`$ we have to restrict the sum in Eq. (82) by $`k<k_\mathrm{c}`$, and to account for the c.p.o. (66). Prescribing the phase field $`\mathrm{\Phi }_𝒌`$ for each s.p.f. in Eq. (64) implies a two to one restriction for each of the three directions of $`\mathbf{}\mathrm{\Phi }_𝒌`$. For each s.p. state $`𝒌`$ this reduces the entropy by $`k_\text{B}\mathrm{ln}8`$ (as long as $`n_k1`$). With these specifications we obtain

$`S_\mathrm{c}(T,V,N)`$ $`=`$ $`{\displaystyle \frac{k_\text{B}}{2}}{\displaystyle \underset{k<k_\mathrm{c}}{}}[(1+2n_k)\mathrm{ln}(1+2n_k)`$ (84)
$`2n_k\mathrm{ln}(2n_k)\mathrm{ln}8].`$
The ratio $`S_\mathrm{c}/S_\mathrm{M}`$ is the AIBG prediction for the superfluid entropy fraction $`S_\mathrm{s}/S`$,

$$\frac{S_\mathrm{s}}{S}=\frac{S_\mathrm{c}}{S_\mathrm{M}}.$$
(85)
The calculated result is shown in Fig. 3.

### B Experimental detectability

We show that the AIBG prediction for $`S_\mathrm{s}0`$ is at the border of present-day experimental detectability. The most prominent experiment showing that $`S_\mathrm{s}=0`$ (or at least $`S_\mathrm{s}0`$) is the fountain pressure (FP) measurement. For two containers of liquid helium connected by a superleak a temperature difference $`dT`$ produces a pressure difference $`dP`$. Admitting $`S_\mathrm{s}0`$ this FP is given by

$$\left(\frac{dP}{dT}\right)_{\mathrm{FP}}=\frac{S}{V}\left(1\frac{S_\mathrm{s}/N_\mathrm{s}}{S/N}\right).$$
(86)
For $`S_\mathrm{s}=0`$ this reduces to the well-known London relation. The correction term is $`S_\mathrm{s}/S`$ multiplied by $`N/N_\mathrm{s}=\rho /\rho _\mathrm{s}`$; the model prediction for it is shown in Fig. 3.

The most accurate FP measurements are that by Singsaas and Ahlers . These authors assume the validity of the London relation and interpret their experiment as an entropy measurement. For this entropy $`S_{\mathrm{FP}}`$,

$$\frac{S_{\mathrm{FP}}}{V}=\left(\frac{dP}{dT}\right)_{\mathrm{FP}}$$
(87)
they find within the errors no deviation from the true (caloric) entropy $`S`$. The absolute values of $`S`$ near $`T_\lambda `$ are, however, uncertain by about 2% (Refs. ).

TABLE I. Comparison between the caloric ($`S`$) and the fountain pressure ($`S_{\mathrm{FP}}`$) entropy. The input is the experimental $`S_{\mathrm{FP}}(t)`$ and $`C_p(t)`$; the quantity $`\mathrm{\Delta }S=S_\lambda S`$ is calculated from $`C_p`$. The last column shows the resulting marginal evidence for $`SS_{\mathrm{FP}}`$ or, equivalently, for $`S_\mathrm{s}0`$.

| 1 | 2 | 3 | 4 | 5 |
| --- | --- | --- | --- | --- |
| $`|t|`$ | $`{\displaystyle \frac{\rho S_{\mathrm{FP}}}{J\mathrm{cm}^3\mathrm{K}^1}}`$ | $`{\displaystyle \frac{S_{\mathrm{FP}}}{Nk_\text{B}}}`$ | $`{\displaystyle \frac{S_{\mathrm{FP}}+\mathrm{\Delta }S}{Nk_\text{B}}}`$ | $`{\displaystyle \frac{100(SS_{\mathrm{FP}})}{S}}`$ |
| 0.0006906 | 0.2299 | 0.7571 | $`\mathrm{\hspace{0.17em}0.7619}\pm 0.0008`$ | $`\mathrm{\hspace{0.17em}0.07}\pm 0.11`$ |
| 0.000790 | 0.2296 | 0.7561 | $`\mathrm{\hspace{0.17em}0.7615}\pm 0.0008`$ | $`\mathrm{\hspace{0.17em}0.11}\pm 0.11`$ |
| 0.001013 | 0.2294 | 0.7555 | $`\mathrm{\hspace{0.17em}0.7623}\pm 0.0008`$ | $`\mathrm{\hspace{0.17em}0.02}\pm 0.11`$ |
| 0.001794 | 0.2275 | 0.7493 | $`\mathrm{\hspace{0.17em}0.7607}\pm 0.0009`$ | $`\mathrm{\hspace{0.17em}0.23}\pm 0.11`$ |
| 0.003177 | 0.2250 | 0.7411 | $`\mathrm{\hspace{0.17em}0.7601}\pm 0.0010`$ | $`\mathrm{\hspace{0.17em}0.30}\pm 0.13`$ |
| 0.005662 | 0.2208 | 0.7274 | $`\mathrm{\hspace{0.17em}0.7592}\pm 0.0011`$ | $`\mathrm{\hspace{0.17em}0.43}\pm 0.14`$ |
| 0.01007 | 0.2145 | 0.7069 | $`\mathrm{\hspace{0.17em}0.7596}\pm 0.0013`$ | $`\mathrm{\hspace{0.17em}0.39}\pm 0.17`$ |
| 0.03338 | 0.1854 | 0.6116 | $`\mathrm{\hspace{0.17em}0.7578}\pm 0.0021`$ | $`\mathrm{\hspace{0.17em}0.74}\pm 0.34`$ |
| 0.07973 | 0.1407 | 0.4648 | $`\mathrm{\hspace{0.17em}0.7547}\pm 0.0034`$ | $`\mathrm{\hspace{0.17em}1.63}\pm 0.71`$ |

The theoretical prediction of Fig. 3 suggests that one should compare the temperature dependences of $`S`$ and $`S_{\mathrm{FP}}`$ near $`T_\lambda `$ rather than the absolute values. This is done in Table I in a number of steps:

1. We start with the experimental values for $`S_{\mathrm{FP}}`$ in column 1 and 2. By using $`\rho (t)`$ from equation (A1) of Ref.  we relate the entropy to the particle number rather than to the density (column 3).
2. For column 4 we calculated $`\mathrm{\Delta }S=S_\lambda S(T)=𝑑TC_P/T`$ from the experimental specific heat $`C_P`$ of Ref. . If $`S_{\mathrm{FP}}`$ were equal to $`S`$ then column 4 should show the temperature independent $`S_\lambda =S(T_\lambda )`$.
3. From the temperature dependence of column 4 we deduce $`S_\lambda /N=0.7624k_\text{B}`$ as the limit of $`S_{\mathrm{FP}}+\mathrm{\Delta }S`$ for $`t0^{}`$. Then

$$\frac{SS_{\mathrm{FP}}}{S}=\frac{S_\lambda [S_{\mathrm{FP}}(t)+\mathrm{\Delta }S(t)]}{S_\lambda \mathrm{\Delta }S(t)}=\frac{S_\mathrm{s}/N_\mathrm{s}}{S/N},$$
(88)
may be calculated from column 4. The last equality in Eq. (88) follows from Eqs. (86) and (87). Column 5 displays the resulting experimental evidence (or nonevidence) for $`S_{\mathrm{FP}}S`$.
4. The errors included in the Table are the statistical errors of 0.1% for $`S_{\mathrm{FP}}(t)`$, and a 1% error for $`C_P`$. The systematic errors of $`S_{\mathrm{FP}}(t)`$ are not included; they should be less important because we consider only the $`t`$-dependence (and not the absolute values) in a relatively small interval.

Figure 4 compares column 5 of the Table with our theoretical prediction. There seems to be some indication of a deviation $`SS_{\mathrm{FP}}(t)`$. The error bars (from the Table) show, however, that the experimental evidence for such a deviation is at most marginally significant. In any case, the compilation of Table I and Fig. 4 shows that, and in which way, our theoretical prediction is within the reach of experimental detectability.

## VII Concluding remarks

We have proposed a microscopic model for the $`\lambda `$-transition of liquid helium. The model uses phenomenological assumptions, in particular the transition itself is introduced by the analogy to the IBG. The basic idea of the AIBG is that phase ordering leads to an extra energy $`\nu _q^{\mathrm{\hspace{0.17em}2}}`$ which —evaluated with the IBG occupation numbers— yields a logarithmic singularity. In this model the superfluid density is not identical to the square of the order parameter field ($`\rho _\mathrm{s}m|\psi |^2`$). A fit to the experimental $`\rho _\mathrm{s}`$ leads to the model prediction for a nonvanishing superfluid entropy $`S_\mathrm{s}`$.

The model is hardly related to standard models of liquid helium which are either microscopic variational approaches (like Refs. ) or the Landau-Wilson renormalization group theory. Because of its novelty and originality our model cannot compete in quality of foundation, exactness and completeness with those other theories. It leads, however, to specific predictions about the structure of $`\rho _\mathrm{s}`$. These predictions are without competition by other approaches and they can and should be tested experimentally. Therefore the presentation of the model seems appropriate in spite of a number of unsolved questions concerning the model assumptions.

In contrast to other approaches the AIBG offers also a possible solution of the so-called microscopic problem of liquid helium. This problem was formulated by Uhlenbeck: If $`\rho _\mathrm{s}`$ is identified with a single quantum state ($`\rho _\mathrm{s}=m|\psi |^2`$) then the approach to equilibrium ($`\rho _\mathrm{s}=\rho _\mathrm{s}(T)`$) cannot be understood. In the AIBG this problem is solved by the contribution of noncondensed particles to the superfluid density.


format: text

text:
# The origin of [Fe II] emission in NGC 4151

## 1 Introduction

Seyfert galaxies exhibit strong emission lines from forbidden electronic transitions—notably those of oxygen, nitrogen and sulphur at visible wavelengths \[Seyfert 1943\] and iron in the near-infrared \[Rieke & Lebofsky 1981, Simpson et al. 1996\]. These transitions can occur where ions are excited by collisions with electrons, yet can decay radiatively before subsequent collisions lead to de-excitation \[Peterson 1997, Osterbrock 1989\]. Hence the lines are associated with low-density regions of ionized hydrogen, where free electrons are abundant but the timescale for collisions is sufficiently long.

Some species, including neutral oxygen and Fe+, have similar ionization potentials to hydrogen, so do not survive where it is fully ionized. The most favourable conditions for \[O I\] and \[Fe II\] emission arise in partially ionized zones, the extent of which influences the observed intensity. Since H II regions around stars normally give way sharply to neutral gas beyond the Strömgren radius, their \[O I\] and \[Fe II\] lines are relatively weak. The central spectra of Seyferts, in contrast, indicate the presence of extensive partially-ionized hydrogen clouds as well as very highly-ionized species. These differences are explained somewhat naturally if the emission line clouds in Seyferts are photoionized by the ‘power-law’ continuum of the central active galactic nucleus (AGN), which is spectrally flatter than a black body \[Osterbrock 1989, Mouri et al. 1990\]. Numerical simulations \[Ferland & Netzer 1983, Mouri et al. 2000\] lend credence to this scenario.

The second mechanism which can generate enhanced zones of \[O I\] and \[Fe II\] emission is shock excitation of gas \[Osterbrock 1989\]. This appears to be the dominant cause in starburst galaxies, where shocks due to supernova remnants are important and supernovae have been associated directly with enhanced \[Fe II\] \[Lester et al. 1990, Greenhouse et al. 1991, Forbes et al. 1993, Greenhouse et al. 1997\]. Comparison of observed spectral line ratios with recent numerical models helps to distinguish between shocking and photoionization \[Mouri et al. 2000\], and indeed proves consistent with the former dominating in starbursts and the latter in Seyfert nuclei. Doppler broadening can also be used to check for kinematic disturbances.

Long-slit spectroscopy of Seyfert galaxies has revealed a spatially extended narrow line region (ENLR) in the visible (at radius greater than 4 arcsec in the case of NGC 4151) beyond the central narrow-line region (NLR). Unger et al. (1987) traced the bright \[O III\] lines at 4960Å and 5007Å, generated in fully-ionized regions, along with H$`\beta `$ (4861Å) in a sample of galaxies. The features were shown to have minimal line-of-sight velocity structure consistent with simple rotation, and line ratios characteristic of high excitation levels, pointing once again to photoionization. Subsequent studies \[Penston et al. 1990, Robinson et al. 1994\] have tended to support this conclusion. The elongated, sometimes visibly conical morphology of ENLR (and NLR resolved with the Hubble Space Telescope: HST) points to a link with collimated radiation from the nucleus \[Unger et al. 1987, Mulchaey et al. 1996, Evans et al. 1993\].

In the near-infrared, the origin of strong, extended \[Fe II\] emission in Seyfert galaxies has been somewhat controversial, with evidence for both shock excitation and photoionization occurring. It has been argued \[Greenhouse et al. 1991\] that shocks must be important because iron in our interstellar medium is mainly condensed onto grains which shocks can destroy via sputtering \[Savage & Sembach 1996\]. However, recent studies \[Simpson et al. 1996, Mouri et al. 2000\] indicate that the ratio \[Fe II\]/\[O I\] in both Seyferts and starbursts is consistent with the usual depletion levels of iron and oxygen from the gaseous phase, and that certain other signs of grain destruction are absent. Even if metal enhancement is not the dominant factor, it is likely that shocks due to outflows along the radio jets of AGN do make some contribution to emission from partially-ionized zones \[Simpson et al. 1996, Morse et al. 1996\]. In particular, Knop et al. 1996 find that whilst \[Fe II\] in NGC 4151 seems kinematically compatible with the ENLR, the 1.257 $`\mu m`$ feature is broader than the nearby Pa$`\beta `$ line

NGC 4151 is both a natural candidate in which to study \[Fe II\] excitation and an important test case for understanding the structure of active galaxies. It is one of the brightest nearby Seyfert galaxies ($`\mathrm{V}11.5`$, $`cz=995`$km s-1), yet eludes firm classification in the framework of unified AGN models \[Antonucci 1993\]. The highly extended narrow line emission and radio jet \[Perez et al. 1989, Pedlar et al. 1993\] are characteristic of a Seyfert 2 nucleus, whose collimation axis is almost perpendicular to the line of sight. However, the object also exhibits the strong broad lines and variability of a Seyfert 1 \[Maoz et al. 1991\], suggesting quite a different orientation. NGC 4151 has therefore been the subject of many past projects involving slit spectroscopy and narrow-band imaging.

An unusual feature which facilitates the present investigation is the pronounced misalignment between the ENLR at position angle $``$50 \[Perez et al. 1989, Evans et al. 1993, Kaiser et al. 2000\], and radio jet at $`77^{}`$ \[Pedlar et al. 1993, Mundell et al. 1995\]. Whatever its origin, this provides an opportunity to associate \[Fe II\] directly with the optical ENLR and/or radio jet, by tracing the emission along both axes. If a clear separation can be made, it will be possible to determine the relative influence of ionization by collimated X-rays and shocking by outflowing radio plasma. Such a task is inherently suited to integral field spectroscopy, which can form a complete, homogeneous picture without prior assumptions, slit alignment problems or the waveband limitations of a filter.

The intrumentation is discussed in the next section, followed, in Section 3, by details of the observations. In Section 4, we present the method of data reduction and analysis which presents particular challenges since the technological limitations required the fibre outputs to overlap at the slit. However this does not degrade the data provided that the field is critically sampled at the IFU input \[Allington-Smith & Content 1998\]. Finally, the results are discussed in Section 5 and compared with other work.

## 2 Instrumentation for integral field spectroscopy

The SMIRFS-IFU is an experimental 72-fibre integral field unit \[Haynes et al. 1999, Lee 1998\] feeding the CGS4 spectrograph\[Wright et al. 1993\] at the UK Infrared Telescope (UKIRT). It provides integral field spectroscopy in the J and H bands via lensed fused-silica fibres. It was designed and built by the Astronomical Instrumentation Group at Durham University as a prototype for much larger integral field units: TEIFU on the William Herschel Telescope\[Murray et al. 2000\] and the GMOS-IFU on Gemini \[Allington-Smith et al. 2000\]. Futher details may be found in Haynes et al. (1999).

SMIRFS is designed to intercept the telescope beam without modification to existing instrumentation, but consequently requires removal of the CGS4 calibration unit. In place of this module, the SMIRFS slit projection unit feeds light from the fibre bundle into the long-slit, circumventing the usual beam path; the IFU input is connected to another port of the UKIRT instrument support structure. As a result of this arrangement, any reference observations must be taken using the main telescope aperture.

At the image plane, a hexagonal microlens array provides optimal coupling to the slow telescope beam and avoids dead space between fibres (see Haynes et al. for more information). The field of view is $`6\times 4`$ arcsec, with 0.62 arcsec spacing between adjacent fibre centres. A second set of microlenses at the output restore the correct focal ratio for CGS4. When using the long camera, each fibre projects a spot whose FWHM is $``$1.8 pixels at the detector, giving a comparable spectral resolution to the standard 2 pixel wide slit. These spots are spaced two pixels apart, placing the instrument on the borderline between the ‘maximally packed’ and ‘resolved peaks’ categories defined by Turner (2001), although throughput variations obscure the regular pattern in practice. Since there is only one field, sky subtraction is performed by nodding the telescope off source.

The capability of the SMIRFS-IFU was previously demonstrated in a June 1997 commissioning run. Its throughput relative to CGS4 alone was measured to be $``$50%; since then, minor adjustments have been made to improve performance.

## 3 Observations

### 3.1 Science programme

NGC 4151 was observed with the SMIRFS-IFU on 15–17 February 1998, as part of a sample of Seyfert galaxies. Details of the measurement parameters are presented in table 1. The selected region of the J-band includes both \[Fe II\] and Pa$`\beta `$, at rest wavelengths of 1.2567$`\mu `$m and 1.2818$`\mu `$m respectively. Previous studies have sometimes used \[Fe II\] at 1.6435$`\mu `$m and Br$`\gamma `$ at 2.1655$`\mu `$m; such measurements are equivalent to ours, because the same upper energy level of Fe+ is involved. Since, however, the longer-wavelength lines are much further apart, they cannot be observed together at high dispersion and their ratio is more sensitive to differential extinction along the line of sight. The J-band \[Fe II\] feature is also brighter, and the uncooled SMIRFS better suited to working in that regime.

The spatial resolution quoted in table 1 is for the final, combined data set; it was estimated by producing images of the unresolved \[S IX\] line and Pa$`\beta `$ broad line. These FWHM values correspond to the Nyquist scale, suggesting that the seeing disc was undersampled—the implications are discussed briefly in §4.

The long axis of the IFU was orientated at $`50^{}`$ east of north, along the optical ENLR of the galaxy. The field of view, however, only covers the inner NLR at $`<4\text{′′}`$ \[Winge et al. 1999\], which is orientated less favourably at $``$60 \[Kaiser et al. 2000\]. With arcsecond resolution, this gives $``$1′′ FWHM separation from the radio axis. Hence a mosaic was made of several telescope pointings around the nucleus.

Target centering was performed with the aid of a dedicated program, which derives a rough map of the illuminated input elements from the data. For the peripheral mosaic positions, it was necessary to rely on specified telescope offsets from the centre. Each position was observed in turn, following the repeated sequence object–sky–object; this maximizes the time on source without giving unequal object and sky exposures at individual pointings. The 180s integration time was just long enough for the background to dominate the read noise and dark current. The whole mosaic was repeated until the full exposure was reached.

In the absence of the calibration unit, sky emission lines were used as a wavelength reference. Some arc lamp observations were also taken through the telescope. The IFU was flat-fielded using the illuminated dome and the twilight sky, whilst a detector flat was made with the calibration unit, prior to the installation of SMIRFS. The bright F and A stars BS4572 and HD105601 were observed in order to calibrate the spectral response of the instrument and measure sky absorption features.

### 3.2 IFU set-up and characteristics

The IFU was aligned at the focal plane by projecting a laser beam back through the fibres onto the UKIRT secondary mirror. The CGS4 slit was set parallel to the row of fibre outputs, before adjusting the optics of the slit projection unit so that fibres were imaged two pixels apart in the spectrograph. The two pixel magnification was verified by illuminating the IFU with a mask in place, generating an image of the fibre slit at the detector with a gap for every input row. Since the mapping of the input to the slit reverses every row, to ensure that elements are adjacent both at the input and output (Fig. 1), the distance between illuminated band centres at the detector alternates about the row separation. By fitting a straight line to the centroids, the magnification was measured as $`2.002\pm 0.005`$ pixels/fibre.

The J-band throughput of the IFU was measured by comparing observations of the standard stars HD 84800 and BS 4069, taken with the CGS4 slit and the SMIRFS-IFU respectively. Examination of the point-spread function (PSF) indicates that truncation of the seeing disc by the 4′′ slit did not cause significant light loss. The difference in airmass was small, giving an expected error due to atmospheric extinction of $`<1\%`$. It was found that the throughput is 50%, consistent with the figure of 49% from the 1997 commissioning run \[Lee 1998\].

Variations in the fibre transmissions were quantified by extracting a flat-field spectrum for each element (see §4.2) and integrating in wavelength (Fig. 2). The values therefore include crosstalk from neighbouring fibres. For Gaussian profiles with FWHM 1.8 pixels, the contribution from either neighbour is $``$10% of the total. Ignoring further smoothing due to interpolation and any correlation between throughputs, the measured and true RMS values are thus related as follows:

$$\sigma _{\mathrm{data}}=\sqrt{(0.8\sigma _{\mathrm{IFU}})^2+2\times (0.1\sigma _{\mathrm{IFU}})^2}=0.8\sigma _{\mathrm{IFU}}.$$

The measured RMS, $`\sigma _{\mathrm{data}}`$, indicates the signal-to-noise fluctuation; it is 8% overall or 6% excluding the dead fibre and some vignetted elements. The corresponding values of $`\sigma _{\mathrm{IFU}}`$, which characterizes the fibre bundle, are $``$10% and $``$7% respectively. It is concluded that whilst careful flat-fielding is important, inherent differences in the signal-to-noise ratio are reasonably modest.

The point spread function does not vary substantially between fibres. Whilst a previous estimate \[Haynes et al. 1999\] indicates that some profiles are up to 20% broader than others, the PSF is only known to an accuracy of $``$10%. Any variations at this level will contribute slightly to the Nyquist-scale noise in line width maps. The effect on line flux noise is negligible, since the broadening is a fraction of a pixel. Any inhomogeneities introduced by the fibres can be smoothed out to some extent by combining observations with slight spatial offsets.

## 4 Data reduction and analysis

### 4.1 Overview

Data were reduced in IRAF, using several purpose-written scripts and programs alongside standard tasks. As discussed in Turner (2001), the best method of dealing with fibre spectra depends largely on their separation and FWHM at the detector. For the SMIRFS-IFU, spectra are much more densely packed than in a multi-object design, precluding the use of existing extraction software. New tasks were also needed for reconstructing and working with the observed $`x`$$`y`$$`\lambda `$ volumes in the form of datacubes, and for bad pixel correction. Because the slit maps to the field of view in a complex way, flexure in the spectrograph cannot be compensated by re-centring the target (this would be true anyway for off-centre mosaic positions). Observations must therefore be co-added after conversion to datacubes, rather than in their raw format. This requires effective bad pixel correction for integration pairs.

Turner (2001) presents an alternative reduction method for fibre-IFU data, inspired by the very dense packing of spectra, one per row, in TEIFU \[Murray et al. 2000\]. With the long CGS4 camera, the raw format of SMIRFS-IFU data is intermediate between that of TEIFU and multi-object designs. Since fibres are imaged with a FWHM of 1.8 pixels, the output pattern is almost critically sampled and is acceptable to interpolate in pixel coordinates. Moreover, the two pixel fibre separation is large enough that moderate smoothing by the interpolant will be small in relation to the image resolution. At the same time, flexure-related errors and undersampling at the IFU input rule out much gain in precision through TEIFU-style reduction. Hence, as for most fibre-based instruments, the reduction is based on extracting 1D spectra.

Finding absolute fibre positions on the detector after flexure is complicated by vignetting at the slit ends, which prevents direct location of the first and last elements. With an already limited field of view, it was decided to rely on the single broken fibre for measuring offsets between frames, rather than masking off row ends at the input to create regular gaps. This is usually adequate, if somewhat awkward and instrument-specific.

With 0.62′′ elements at the IFU input, undersampling is potentially an important issue. Inspection of various observations (with and without SMIRFS) indicates that the seeing during the February 1998 run was $``$0.9′′, including integration over the input microlenses; this figure is relatively favourable for such coarse sampling. A very basic simulation of spatial frequency propagation through the instrument shows that, depending on the alignment of fibres, pixels and image features, aliased frequencies in the data can have amplitudes up to 10% of the mean \[Turner 2001\]. However, the Fourier spectrum is only affected significantly above $``$80% of the Nyquist frequency, so moderate smoothing during extraction and reconstruction can easily eliminate much of the contamination. In principle, offset observations can be combined with phase shifts such that aliasing is removed \[Lauer 1999\], but in this case the main effect would be to amplify noise because the offsets are not known or controlled precisely (and the seeing varies). Nevertheless, simply co-adding datacubes after alignment on a common grid does suppress residual alias components, which add out of phase like noise. At the end of the reduction process, aliasing artifacts are not a major source of error.

### 4.2 Calibration and extraction

Each observation consists of two integrations, offset by one detector pixel in the spectral direction to account for bad pixels. These raw frames are largely obscured by non-uniform dark current. The first step was to subtract a sky pair from an object pair, removing dark current in the process. However, this double use means that sky frames cannot be smoothed spatially to reduce their noise level before subtraction. The result was then divided by a detector flat-field image, obtained with the calibration unit before the IFU was installed.

A less accurate dark frame, taken earlier in the night, was also subtracted from the object separately. This reveals the fibre throughput pattern so that the dead fibre can be located. Ideally, this would be done using a sky minus dark image, but flexure of a fraction of a pixel can occur even between consecutive pointings. This is important for image reconstruction, but less so for sky subtraction—assuming Gaussian fibre profiles and an offset of 0.3 pixels, the sky residual due to throughput variations is likely to be around 3%, falling below the noise level.

A program was written to combine the pairs of offset, sky subtracted frames with bad pixel removal. This makes a crude noise estimate based on the RMS residuals between the two images; any input pixel which differs from the median of its neighbours by more than a specified number of standard deviations is excluded. The procedure can sometimes truncate narrow arc lamp emission lines, whose gradients can be much greater than the noise level, but is safe for the galaxy spectra (which have had sky lines subtracted already). Any obvious remaining blemishes were removed manually.

Spectra were straightened with respect to detector rows by applying a small pre-determined shift to each column. A curvature model was constructed by measuring the centroid of a suitable image at each wavelength and fitting a low-order polynomial to the values. Interpolation was performed using a cubic spline, which is good for dealing with slightly undersampled data. Inspection of corrected frames verifies that the curvature is practically constant along the slit and hence independent of flexure. The maximum difference between columns is 0.7 rows.

Absolute fibre positions were found by summing the straightened object minus dark image in wavelength and fitting the dead fibre position using standard IRAF routines. The corresponding sky subtracted image was resampled so that fibres fell at the mid-points between odd pairs of rows. These pairs were combined with the *blkavg* task, producing a stack of extracted fibre spectra (boundary rows were discarded). Sometimes spatial gradients make the dead fibre difficult to identify, but the procedure can be repeated if a reconstruction appears distorted.

The stack of extracted spectra was divided by a fibre flat, created in the same way but averaged over its wavelength range. Sky emission lines in the object minus dark image were used for wavelength calibration. These allow correction of grating offsets between observations and also provide an absolute calibration. The spectra from every observation were interpolated onto a suitable standard grid, truncating the end few pixels which fell outside the common range. The error in determining the relative offsets is around 0.2Å. The error in the linear wavelength fit is 0.7Å, which is more than adequate for the present work. Separate solutions were not determined for the spectra within each observation; variation along the slit ($`\pm 0.3`$Å) is smaller than the final random velocity errors in the mosaic and much smaller than the true gradient in the galaxy.

Finally, the stack was divided by a continuum-normalized stellar spectrum, to remove telluric absorption features. A Pa$`\beta `$ line in one of the stars was interpolated over. Unfortunately, the mosaicing process makes it difficult to re-observe the standard stars frequently, so there was sometimes a non-negligible difference in airmass from the galaxy. Although the main absorption band at 1.269$`\mu `$m has been removed effectively in the final data, there may be some absorption residuals at around 2% of the continuum level.

### 4.3 Construction of datacubes and mosaicing

Sets of extracted spectra were interpolated spatially at each wavelength step, onto a square grid of 0.15′′ pixels.111We used the method of Renka and Cline\[Renka 1983\] for convenience. The constituent datacubes of each mosaic cycle were combined with appropriate offsets using the program *mosaic*, which has been included as a prototype task in the *imspec* package \[Turner 2001\]. This calculates the smallest common grid from a list of pointings and resamples the images at each wavelength using a tapered sinc function. The precise interpolation method is not important, since the datacube grids are much finer than the resolution of the data. A mask image is also produced, recording the number of input pixels contributing at each point.

Unfortunately, telescope pointing errors of a fraction of a fibre are evident. This is partly because a dichroic mirror is used to guide on the target at visible wavelengths whilst observing in the infrared; the two images drift apart as atmospheric dispersion changes with the airmass along the line of sight. Modelling this effect may systematically reduce errors, but was outside the scope of the current project. There may consequently be some image distortion at the points of overlap; mosaicing nevertheless reveals any wider-scale extension of spectral features within the noise limit of the data. For the NGC 4151 observations, emission turns out to be confined mainly to the central field; three mosaic positions with appreciable signal-to-noise, offset diagonally east and west, were finally used for mapping.

After constructing each separate mosaic, the resulting datacubes were combined using *comcubes*, also included in *imspec*. This is very similar to *mosaic*, but calculates the offsets between datacubes by locating the centroids of the intensity peaks after integrating in wavelength. The individual overlap masks are also combined to show the total number of data points contributing at each output pixel. The final datacube is divided by the mask, to scale the pixels correctly.

The spectral response of the instrument was determined by fitting piecewise polynomials to the continua of the standard star observations, dividing the results by black-body curves appropriate for their spectral types. This was achieved using the *standard*, *sensfunc* and *calibrate* tasks in the IRAF *onedspec* package. The response varies at the level of 5%, and with atmospheric absorption features already removed it is not strongly dependent on airmass. The accuracy to which the curve can be determined is limited mainly by the wings of absorption lines in the stars and atmosphere, which are difficult to distinguish from modulation in the continuum. The final datacube was divided by the response function using *calibrate* directly and the values normalized to a unit mean.

### 4.4 Mapping the flux and velocity field

A velocity map was obtained from the combined mosaic. Spectra with mean fluxes below 2.5% of the central value (about $`3\sigma `$) were not included, since line emission is barely detected at the corresponding distance. The NOAO task *rvidlines* was used to measure a combined \[Fe II\] and narrow Pa$`\beta `$ redshift for each spectrum. This is more reliable for low signal-to-noise peripheral points than using \[Fe II\] alone. The results show that differences in the velocity fields of the two lines are comparable with the measurement errors and small compared to the line widths. Hence any deviations introduced by Pa$`\beta `$ will not affect line strength measurements significantly. The plain text results from *rvreidlines* were redirected into a file and converted to a 2D image using *listoim* in *imspec*.

Line fluxes were measured using the *specmap* task in *imspec*. This reads the datacube, velocity map and a list of wavebands for line and continuum estimates, producing a flux map as output. The measurement bands and linear continuum estimates can be overplotted on the spectra for verification. Over the wavelength range observed, there is little pure continuum to provide a suitable reference for \[Fe II\], but experimenting with different bands shows that the derived spatial profiles are robust.

An excitation map, showing the strength of \[Fe II\] relative to Pa$`\beta `$, was also produced. Since this involves dividing one image by the other, the result is sensitive to noise where emission is weak. The map was therefore created using line images from a spatially smoothed version of the datacube. Rough line width estimates were made by examining individual spectra with *implot* and *splot*.

## 5 Results and comparison with other data

A wavelength-integrated image of NGC 4151 is shown in Fig. 3(a). Whilst some horizontal broadening is caused by the IFU, because spectra overlap, there is also a real NW–SE extension of the continuum, following the stellar distribution seen in optical images. Fig. 3(b) is an image of the Pa$`\beta `$ broad line (including the narrow line at $``$6% of the total), which gives an indication of the spatial resolution.

A central spectrum is presented in Fig. 4; this includes the broad and narrow Pa$`\beta `$ components, \[Fe II\] and \[S IX\]. There may also be contributions from a weaker 1.2703$`\mu `$m \[Fe II\] feature, superimposed on the blue wing of the Pa$`\beta `$ broad line, and He I at 1.2791$`\mu `$m, to the left of the Pa$`\beta `$ narrow line. The small bump to the right of \[S IX\] is probably due to imperfect removal of atmospheric absorption features.

The combined \[Fe II\] and Pa$`\beta `$ velocity map (Fig. 5) traces the rotation of gas about an axis projected at a position angle of $``$110. It is broadly consistent both with previous slit measurements of the same features \[Knop et al. 1996\] and with optical spectroscopy of the NLR \[Mediavilla & Arribas 1995, Winge et al. 1999\]. Estimates of the line widths are complicated by the Pa$`\beta `$ broad line shape, lack of clear continuum and modest signal-to-noise ratio off centre. The line profiles have not been modelled with multiple components, but their overall deconvolved widths are close to the values found by Knop et al, who note that the width of \[Fe II\] relative to Pa$`\beta `$ suggests a contribution from shock excitation. At the centre, the FWHM values are approximately 400 km s-1 for \[Fe II\] and 250km s-1 for Pa$`\beta `$. Both lines are broader in the brightest parts of the surrounding extended emission (see below), reaching $``$600 km s-1 and $``$400 km s-1 respectively. Additional broadening in the direction of the radio jet is not evident.

As an alternative to Knop et al.’s suggestion that the greater velocity width of the \[Fe II\] emission is due to a contribution from shock excitation, we note that our result is consistent with that for a sample of Seyfert 2 galaxies \[Veilleux et al. 1997\]. Veilleux et al. cite this as evidence that the two lines are emitted in different volumes of gas. We also note that Kaiser et al. (2000) found different components in the velocity field of NGC4151 revealed by \[OIII\] emission. They suggested that the low velocity dispersion component arises from smooth flow in the gravitational field of the galaxy while the high dispersion component results from an outflowing wind. Our result could arise if the Pa$`\beta `$/\[Fe II\] ratio was higher in the low dispersion component than in the high dispersion component.

The spatial distribution of \[Fe II\] is mapped in Fig. 6, alongside examples of individual spectra. The image profile is extended along both axes of the IFU with an intrinsic RMS of 0.3–0.4′′ horizontally and 0.9–1.0′′ vertically (FWHM 0.8′′ and 2.2′′ respectively assuming a gaussian distribution). The narrowest point is at the intensity peak, as might be expected for a biconical shape. Although these scales are comparable to the image resolution, it is clear that the \[Fe II\] is aligned closer to the optical line emission than to the radio axis. The best fit to the slope of the horizontal centroid has a PA of $``$54, part way between the NLR and ENLR axes. For the data to fit the radio axis equally well, either end of the outside ($``$25%) contour would have to be displaced by $``$0.8′′.

Fig. 7 compares the \[Fe II\] map directly with a narrow-band HST image of \[O III\] line emission \[Kaiser et al. 2000\] and a contour map of the radio jet at 8 MHz \[Mundell et al. 1995\]. It can be seen that the region over which \[Fe II\] is detected in the SMIRFS data corresponds to the extent of the narrow line region. If there were a strong correspondence with the radio jet, one might expect a more elongated distribution, resolving some knot structure at $`77^{}`$, where the outflow interacts with gas. Although detection of the larger-scale ENLR would have provided a better separation, the spatial association of \[Fe II\] with photoionized gas is more convincing than with shocking by radio plasma. The F502N filter used for imaging \[O III\] covers almost all of the NLR velocity range \[Hutchings et al. 1999, Kaiser et al. 2000\], so represents the total emission well. It is perhaps interesting to note that \[Fe II\] appears better correlated with the positions of low velocity dispersion clouds, plotted in \[Kaiser et al. 2000\], than the NLR as a whole. However, this is at variance with our previous suggestion that the greater velocity width of the \[Fe II\] line compared with Pa$`\beta `$ suggests an association with the high dispersion component.

The Pa$`\beta `$ narrow line is slightly extended in both directions, and has a faint wing out to $``$1.5′′ South-West of the centre, where the NLR is brightest. Fig. 8 shows both the Pa$`\beta `$ flux and \[Fe II\]/Pa$`\beta `$ ratio. Since Pa$`\beta `$ falls off fairly rapidly from the centre, the ratio is highly sensitive to noise artifacts. The datacube was therefore smoothed with a 0.9′′ wide Gaussian kernel before producing the excitation map, preserving minimal resolution of spatial structure. The line ratio is $``$0.7 at the centre and $`2`$ along the radio jet at the edges. A dip near the centre is to be expected for photoionization, since the strong incident flux creates large fully ionized zones, favouring Pa$`\beta `$ over \[Fe II\] \[Mouri et al. 2000\]. The high areas could reflect an enhancement of \[Fe II\] emission by shocks, but the line images are more indicative of a drop in Pa$`\beta `$ at those points, with respect to \[O III\], than of strong \[Fe II\]. To the South-West, a few highly-ionized clouds may just happen to give stronger Pa$`\beta `$ away from the radio axis, whilst to the North-East, the NLR is angled slightly more towards the radio jet anyway.

## 6 Conclusions

The distribution of \[Fe II\] emission at the centre of NGC 4151 has been mapped in two dimensions, using the innovative SMIRFS-IFU. Although the line is only detected over a few spatial resolution elements, it is evident that its image profile is better aligned with the visible narrow line region than with the radio jet. The study therefore suggests that \[Fe II\] arises primarily through photoionization of gas, by collimated X-rays from the Seyfert nucleus. The mean velocity field is consistent with this interpretation, indicating that \[Fe II\], Pa$`\beta `$ and the NLR all have similar bulk kinematics. Likewise, the variation in \[Fe II\]/Pa$`\beta `$ is compatible with photoionization dominating. As suggested by previous authors, shock excitation by outflowing radio plasma may make a secondary contribution to the emission. However, the evidence for this comes mainly from line widths, whereas tentative measurements suggest, if anything, that \[Fe II\] is broader along the NLR than the radio jet.

Despite the small field and prototypical nature of the IFU, these results demonstrate the utility of integral field spectroscopy in investigating the properties of active galaxies. Without a two dimensional field, it would have been very difficult to determine the spatial orientation of the line emission. Narrow-band imaging can also fulfil this requirement, given appropriate filters, but follow-up spectroscopy is then required to form a complete picture. The experience gained with the SMIRFS-IFU and associated data analysis has paved the way for subsequent projects, both technical and observational.

Reduction of the IFU data has highlighted a few important observational issues. In particular, location of spectra on the detector is difficult, and it is recommended that more than one reference gap between fibres should be available in future work. With a number of sources of error and inhomogeneity, the availability of multiple observations at each position has been very valuable.

Whilst the results presented here support the argument for \[Fe II\] originating in photoionized gas, there is no doubt that observations at high spatial resolution, as well as of the wider ENLR, will provide further insight into the nature of NGC 4151 and other Seyfert galaxies. With a bright, point-like nucleus, NGC 4151 is an ideal target for J-band observations with adaptive optics. This would provide a useful counterpart to the Hubble Space Telescope observations of \[O III\] emission, which have revealed the cloud structure and detailed kinematics of the NLR. In conjunction with recent models of photoionized and shock excited gas, comparison of optical and infrared lines in individual clouds would provide important information about the physical conditions and kinematic structure of the NLR. The misalignment with the radio jet will again help determine its influence on the ambient gas.

## Acknowledgements

The authors thank George Dodsworth and Ian Lewis for their work on the design and construction of the SMIRFS IFU and Ray Sharples, Ian Parry and Reynier Peletier for their valuable contributions. The United Kingdom Infrared Telescope is operated by the Joint Astronomy Centre on behalf of the U.K. Particle Physics and Astronomy Research Council.